Model-averaged predictions, optionally with standard errors.
# S3 method for averaging
predict(object, newdata = NULL, se.fit = FALSE,
interval = NULL, type = NA, backtransform = FALSE, full = TRUE, ...)If se.fit = FALSE, a vector of predictions, otherwise a list
with components: fit containing the predictions, and se.fit with
the estimated standard errors.
an object returned by model.avg.
optional data.frame in which to look for variables
with which to predict. If omitted, the fitted values are used.
logical, indicates if standard errors should be returned.
This has any effect only if the predict methods for each of the
component models support it.
currently not used.
the type of predictions to return (see documentation for
predict appropriate for the class of used component models). If
omitted, the default type is used. See ‘Details’.
if TRUE, the averaged predictions are
back-transformed from link scale to response scale. This makes sense
provided that all component models use the same family, and the
prediction from each of the component models is calculated on the link
scale (as specified by type. For glm, use type =
"link").
See ‘Details’.
if TRUE, the full model-averaged coefficients are used
(only if se.fit = FALSE and the component objects are a result of
lm).
arguments to be passed to respective predict
method (e.g. level for lme model).
Kamil Bartoń
predicting is possible only with averaging objects with
"modelList" attribute, i.e. those created via model.avg
from a model list, or from model.selection object with argument fit
= TRUE (which will recreate the model objects, see model.avg).
If all the component models are ordinary linear models, the prediction can be
made either with the full averaged coefficients (the argument full =
TRUE this is the default) or subset-averaged coefficients. Otherwise the
prediction is obtained by calling predict on each component model and
weighted averaging the results, which corresponds to the assumption that all
predictors are present in all models, but those not estimated are equal zero
(see ‘Note’ in model.avg). Predictions from component models
with standard errors are passed to par.avg and averaged in the same way
as the coefficients are.
Predictions on the response scale from generalized models can be calculated by
averaging predictions of each model on the link scale, followed by inverse
transformation (this is achieved with type = "link" and
backtransform = TRUE). This is only possible if all component models use
the same family and link function. Alternatively, predictions from each model on
response scale may be averaged (with type = "response" and
backtransform = FALSE). Note that this leads to results differing from
those calculated with the former method. See also
predict.glm.
model.avg, and par.avg for details of model-averaged
parameter calculation.
# \dontshow{
oop <- options(na.action="na.fail")
if(require(graphics)) { # }
# Example from Burnham and Anderson (2002), page 100:
fm1 <- lm(y ~ X1 + X2 + X3 + X4, data = Cement)
ms1 <- dredge(fm1)
confset.95p <- get.models(ms1, subset = cumsum(weight) <= .95)
avgm <- model.avg(confset.95p)
nseq <- function(x, len = length(x)) seq(min(x, na.rm = TRUE),
max(x, na.rm=TRUE), length = len)
# New predictors: X1 along the range of original data, other
# variables held constant at their means
newdata <- as.data.frame(lapply(lapply(Cement[, -1], mean), rep, 25))
newdata$X1 <- nseq(Cement$X1, nrow(newdata))
n <- length(confset.95p)
# Predictions from each of the models in a set, and with averaged coefficients
pred <- data.frame(
model = sapply(confset.95p, predict, newdata = newdata),
averaged.subset = predict(avgm, newdata, full = FALSE),
averaged.full = predict(avgm, newdata, full = TRUE)
)
opal <- palette(c(topo.colors(n), "black", "red", "orange"))
matplot(newdata$X1, pred, type = "l",
lwd = c(rep(2,n),3,3), lty = 1,
xlab = "X1", ylab = "y", col=1:7)
# For comparison, prediction obtained by averaging predictions of the component
# models
pred.se <- predict(avgm, newdata, se.fit = TRUE)
y <- pred.se$fit
ci <- pred.se$se.fit * 2
matplot(newdata$X1, cbind(y, y - ci, y + ci), add = TRUE, type="l",
lty = 2, col = n + 3, lwd = 3)
legend("topleft",
legend=c(lapply(confset.95p, formula),
paste(c("subset", "full"), "averaged"), "averaged predictions + CI"),
lty = 1, lwd = c(rep(2,n),3,3,3), cex = .75, col=1:8)
palette(opal)
}
options(oop)
Run the code above in your browser using DataLab