model.avg
and dredge
QAIC(object, ..., chat)
QAICc(object, ..., chat)
AICc
quasi
family used for models with over-dispersion.
AIC
and BIC
may also be used as a custom rank
function in dredge
and model.avg
.}
}
[object Object]
budworm <- data.frame( ldose = rep(0:5, 2), numdead = c(1, 4, 9, 13, 18, 20, 0, 2, 6, 10, 12, 16), sex = factor(rep(c("M", "F"), c(6, 6)))) budworm$SF = cbind( numdead = budworm$numdead, numalive = 20 - budworm$numdead)
budworm.lg <- glm(SF ~ sex*ldose, data = budworm, family = quasibinomial)
dd <- dredge(budworm.lg, rank = "QAIC", chat = summary(budworm.lg)$dispersion) # Average all models budworm.avg <- model.avg(get.models(dd, seq(nrow(dd))), method="NA") #model.avg(mod[[1]], mod[[2]], rank = "QAIC", rank.args = list(chat = 1))
plot(c(1,32), c(0,1), type = "n", xlab = "dose", ylab = "prob", log = "x") text(2^budworm$ldose, budworm$numdead/20, as.character(budworm$sex)) ld <- seq(0, 5, 0.1)
newdata <- data.frame(ldose=ld, sex=factor(rep("M", length(ld)), levels=levels(budworm$sex)))
# Predictions from global model / Males pred.lg <- predict(budworm.lg, newdata, se.fit=TRUE, type="response") matplot(2^ld, cbind(pred.lg$fit, pred.lg$fit - (2 * pred.lg$se.fit), pred.lg$fit + (2 * pred.lg$se.fit)), add=TRUE, type="l", col=1)
# Predictions from averaged model / Males pred.avg <- predict(budworm.avg, newdata, se.fit=TRUE, type="response") matplot(2^ld, cbind(pred.avg$fit, pred.avg$fit - (2 * pred.avg$se.fit), pred.avg$fit + (2 * pred.avg$se.fit)), add=TRUE, type="l", col=2)
newdata$sex[] <- "F"
# Predictions from global model / Females pred.lg <- predict(budworm.lg, newdata, se.fit=TRUE, type="response") matplot(2^ld, cbind(pred.lg$fit, pred.lg$fit - (2 * pred.lg$se.fit), pred.lg$fit + (2 * pred.lg$se.fit)), add=TRUE, type="l", col=1)
# Predictions from averaged model / Females pred.avg <- predict(budworm.avg, newdata, se.fit=TRUE, type="response") matplot(2^ld, cbind(pred.avg$fit, pred.avg$fit - (2 * pred.avg$se.fit), pred.avg$fit + (2 * pred.avg$se.fit)), add=TRUE, type="l", col=2)
legend("bottomright", legend=c("full", "averaged"), title="Model", col=1:2, lty=1)