modavgCustom(logL, K, modnames = NULL, estimate, se, second.ord = TRUE,
nobs = NULL, uncond.se = "revised", conf.level = 0.95,
c.hat = 1)
NULL
, the function
uses the names in the cand.set list of candidate models. If no names
appear in the list, generic names (e.g.estimate
vector.TRUE
, the function returns the second-order
Akaike information criterion (i.e., AICc)."old"
, or "revised"
, specifying the equation
used to compute the unconditional standard error of a model-averaged
estimate. With uncond.se = "old"
, computations are based on
equation 4.9 of Burnham andc_hat
. Note that values of c.hat
different from 1 are only appropriate for binomial GLM's with trials > 1
(i.e., success/trial or cbimodavgCustom
creates an object of class modavgCustom
with
the following components:modavgCustom
computes a model-averaged estimate from the vector
of parameter estimates specified in estimate
. Estimates and
their associated standard errors must be specified in the same order as
the log-likelihood, number of estimated parameters, and model names.
Estimates provided may be for a parameter of interest (i.e., beta
estimates) or predictions from each model. This function is most useful
when model input is imported into R from other software (e.g., Program
MARK, PRESENCE) or for model classes that are not yet supported by the
other model averaging functions such as modavg
or
modavgPred
.Buckland, S. T., Burnham, K. P., Augustin, N. H. (1997) Model selection: an integral part of inference. Biometrics 53, 603--618.
Burnham, K. P., Anderson, D. R. (2002) Model Selection and Multimodel Inference: a practical information-theoretic approach. Second edition. Springer: New York.
Dail, D., Madsen, L. (2011) Models for estimating abundance from repeated counts of an open population. Biometrics 67, 577--587.
Lebreton, J.-D., Burnham, K. P., Clobert, J., Anderson, D. R. (1992) Modeling survival and testing biological hypotheses using marked animals: a unified approach with case-studies. Ecological Monographs 62, 67--118.
MacKenzie, D. I., Nichols, J. D., Lachman, G. B., Droege, S., Royle, J. A., Langtimm, C. A. (2002) Estimating site occupancy rates when detection probabilities are less than one. Ecology 83, 2248--2255.
Royle, J. A. (2004) N-mixture models for estimating population size from spatially replicated counts. Biometrics 60, 108--115.
AICcCustom
, aictabCustom
modavg
, modavgShrink
,
modavgPred
##model averaging parameter estimate (natural average)
##vector with model LL's
LL <- c(-38.8876, -35.1783, -64.8970)
##vector with number of parameters
Ks <- c(7, 9, 4)
##create a vector of names to trace back models in set
Modnames <- c("Cm1", "Cm2", "Cm3")
##vector of beta estimates for a parameter of interest
model.ests <- c(0.0478, 0.0480, 0.0478)
##vector of SE's of beta estimates for a parameter of interest
model.se.ests <- c(0.0028, 0.0028, 0.0034)
##compute model-averaged estimate and unconditional SE
modavgCustom(logL = LL, K = Ks, modnames = Modnames,
estimate = mode.ests, se = model.se.ests, nobs = 121)
##model-averaging with shrinkage
##set up candidate models
data(min.trap)
Cand.mod <- list( )
##global model
Cand.mod[[1]] <- glm(Num_anura ~ Type + log.Perimeter,
family = poisson, offset = log(Effort),
data = min.trap)
Cand.mod[[2]] <- glm(Num_anura ~ Type + Num_ranatra, family = poisson,
offset = log(Effort), data = min.trap)
Cand.mod[[3]] <- glm(Num_anura ~ log.Perimeter + Num_ranatra,
family = poisson, offset = log(Effort), data = min.trap)
Model.names <- c("Type + log.Perimeter", "Type + Num_ranatra",
"log.Perimeter + Num_ranatra")
##model-averaged estimate with shrinkage (glm model type is already supported)
modavgShrink(cand.set = Cand.mod, modnames = Model.names,
parm = "log.Perimeter")
##equivalent manual version of model-averaging with shrinkage
##this is especially useful when model classes are not supported
##extract vector of LL
LLs <- sapply(Cand.mod, FUN = function(i) logLik(i)[1])
##extract vector of K
Ks <- sapply(Cand.mod, FUN = function(i) attr(logLik(i), "df"))
##extract betas
betas <- sapply(Cand.mod, FUN = function(i) coef(i)["log.Perimeter"])
##second model does not include log.Perimeter
betas[2] <- 0
##extract SE's
ses <- sapply(Cand.mod, FUN = function(i) sqrt(diag(vcov(i))["log.Perimeter"]))
ses[2] <- 0
##extract
modavgCustom(logL = LLs, K = Ks, modnames = Model.names,
nobs = nrow(min.trap), estimate = betas, se = ses)
Run the code above in your browser using DataLab