model.avg: Model averaging

Description

Model averaging based on an information criterion.

Usage

model.avg(object, ..., revised.var = TRUE)
## S3 method for class 'default':
model.avg(object, ..., beta = FALSE, rank = NULL, rank.args = NULL,
    revised.var = TRUE, dispersion = NULL, ct.args = NULL)
## S3 method for class 'model.selection':
model.avg(object, subset, fit = FALSE, ..., revised.var = TRUE)

Arguments

object

a fitted model object or a list of such objects, or a "model.selection" object. See Details.

...

for default method, more fitted model objects. Otherwise, arguments that are passed to the default method.

beta

logical, should standardized coefficients be returned?

rank

optionally, a rank function (returning an information criterion) to use instead of AICc, e.g. BIC or QAIC, may be omitted if object is a model list returned by get.models or a "m

rank.args

optional list of arguments for the rank function. If one is an expression, an x within it is substituted with a current model.

revised.var

logical, indicating whether to use revised formula for standard errors. See par.avg.

dispersion

the dispersion parameter for the family used. See summary.glm. This is used currently only with glm, is silently ignored otherwise.

ct.args

optional list of arguments to be passed to coefTable (besides dispersion).

subset

see subset method for "model.selection" object.

fit

if TRUE, the component models are fitted using get.models. See Details.

Value

An object of class "averaging" is a list with components:
msTablea data.frame with log-likelihood, IC, latex{$\Delta_{IC}$}{Δ{Delta}_IC} and Akaike weights for the component models. "df", "logLik", "AICc", "delta", "weight". Its attribute "term.codes" is a named vector with numerical representation of the terms in the row names of msTable.
avg.modelthe model averaged parameters. A matrix with rows for model terms, columns for averaged coefficient, unconditional standard error, adjusted SE (NA if dfs are not available), z-test statistic and associated p-value.
coef.shrinkagea vector of full model-averaged coefficients, see Note.
coefArraya 3-dimensional array of component models' coefficients, their standard errors and degrees of freedom.
importanceobject of class importance containing relative importance values of each term (including interactions), calculated as a sum of the Akaike weights over all of the models in which the term appears.
formulaa formula corresponding to the one that would be used in a single model. The formula contains only the averaged (fixed) coefficients.
callthe matched call.
The object has following attributes:
rankthe rank function used.
modelLista list of all component model objects.
betaCorresponds to the function argument.
nobsnumber of observations.
revised.varCorresponds to the function argument.

encoding

utf-8

Details

model.avg may be used either with a list of models, or directly with a model.selection object (e.g. returned by dredge). In the latter case, the models from the model selection table are not evaluated unless the argument fit is set to TRUE or some additional arguments are present (such as rank or dispersion). This results in much faster calculation, but has certain drawbacks, because the fitted component model objects are not stored, and some methods (e.g. predict, fitted, model.matrix or vcov) would not be available with the returned object. Otherwise, get.models is called prior to averaging, and ...are passed to it.

For a list of model types that are accepted see list of supported models.

rank is found by a call to match.fun and typically is specified as a function or a symbol or a character string specifying a function to be searched for from the environment of the call to lapply. rank must be a function able to accept model as a first argument and must always return a numeric scalar.

Several standard methods for fitted model objects exist for class averaging, including summary, predict, coef, confint, formula, and vcov. coef, vcov, confint and coefTable accept argument full that if set to TRUE, the full model-averaged coefficients are returned, rather than subset-averaged ones (when full = FALSE, being the default). logLik returns a list of logLik objects for the component models.

References

Burnham, K. P. and Anderson, D. R. (2002) Model selection and multimodel inference: a practical information-theoretic approach. 2nd ed. New York, Springer-Verlag.

Lukacs, P. M., Burnham K. P. and Anderson, D. R. (2009) Model selection bias and Freedman’s paradox. Annals of the Institute of Statistical Mathematics 62(1): 117–125.

Examples

Run this code

# Example from Burnham and Anderson (2002), page 100:
fm1 <- lm(y ~ ., data = Cement, na.action = na.fail)
(ms1 <- dredge(fm1))

#models with delta.aicc < 4
summary(model.avg(ms1, subset = delta < 4))

#or as a 95\% confidence set:
avgmod.95p <- model.avg(ms1, cumsum(weight) <= .95)
confint(avgmod.95p)

# The same result, but re-fitting the models via 'get.models'
confset.95p <- get.models(ms1, cumsum(weight) <= .95)
model.avg(confset.95p)

# Force re-fitting the component models
model.avg(ms1, cumsum(weight) <= .95, fit = TRUE)
# Models are also fitted if additional arguments are given
model.avg(ms1, cumsum(weight) <= .95, rank = "AIC")

# using BIC (Schwarz's Bayesian criterion) to rank the models
BIC <- function(x) AIC(x, k = log(length(residuals(x))))
model.avg(confset.95p, rank = BIC)
# the same result, using AIC directly, with argument k
# 'x' in a quoted 'rank' argument is substituted with a model object
# (in this case it does not make much sense as the number of observations is
# common to all models)
model.avg(confset.95p, rank = AIC, rank.args = alist(k = log(length(residuals(x)))))

Run the code above in your browser using DataLab