R: Create Model Selection Tables from Bayesian Analyses
dictab
R Documentation
Create Model Selection Tables from Bayesian Analyses
Description
This function creates a model selection table based on the deviance
information criterion (DIC). The table ranks the models based on the
DIC and also provides delta DIC and DIC weights. dictab selects
the appropriate function to create the model selection table based on
the object class. The current version works with objects of bugs
and rjags classes.
Usage
dictab(cand.set, modnames = NULL, sort = TRUE, ...)
## S3 method for class 'AICbugs'
dictab(cand.set, modnames = NULL, sort = TRUE, ...)
## S3 method for class 'AICrjags'
dictab(cand.set, modnames = NULL, sort = TRUE, ...)
Arguments
cand.set
a list storing each of the models in the candidate model set.
modnames
a character vector of model names to facilitate the identification of
each model in the model selection table. If 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., Mod1, Mod2) are
supplied in the table in the same order as in the list of candidate
models.
sort
logical. If TRUE, the model selection table is ranked according
to the DIC values.
...
additional arguments passed to the function.
Details
dictab internally creates a new class for the cand.set
list of candidate models, according to the contents of the list. The
current function is implemented for bugs and jags
classes. The function constructs a model selection table based on the
DIC (Spiegelhalter et al. 2002). Note that DIC might not be appropriate
to select among a set of hierarchical models and that modifications to
the information criterion have been proposed (Millar 2009).
Value
dictab creates an object of class dictab with the
following components:
Modname
the names of each model of the candidate model set.
pD
the effective number of estimated parameters for each
model.
DIC
the deviance information criterion for each model.
Delta_DIC
the delta DIC of each model, measuring the difference
in DIC between each model and the top-ranked model.
ModelLik
the relative likelihood of the model given the
data (exp(-0.5*delta[i])). This is not to be confused with the
likelihood of the parameters given the data. The relative likelihood
can then be normalized across all models to get the model probabilities.
DICWt
the DIC weights, sensu Burnham and Anderson (2002) and
Anderson (2008). These measures indicate the level of support (i.e.,
weight of evidence) in favor of any given model being the most
parsimonious among the candidate model set.
Cum.Wt
the cumulative DIC weights. These are only meaningful
if results in table are sorted in decreasing order of DIC weights
(i.e., sort = TRUE).
Deviance
the deviance of each model.
Author(s)
Marc J. Mazerolle
References
Anderson, D. R. (2008) Model-based Inference in the Life Sciences:
a primer on evidence. Springer: New York.
Burnham, K. P., Anderson, D. R. (2002) Model Selection and
Multimodel Inference: a practical information-theoretic
approach. Second edition. Springer: New York.
Spiegelhalter, D. J., Best, N. G., Carlin, B. P., van der Linde,
A. (2002). Bayesian measures of complexity and fit. Journal of the
Royal Statistical Society, Series B64, 583–639.
See Also
aictabCustom, aictab, confset,
DIC, evidence
Examples
##from ?jags example in R2jags package
## Not run:
require(R2jags)
model.file <- system.file(package="R2jags", "model", "schools.txt")
file.show(model.file)
##data
J <- 8.0
y <- c(28.4,7.9,-2.8,6.8,-0.6,0.6,18.0,12.2)
sd <- c(14.9,10.2,16.3,11.0,9.4,11.4,10.4,17.6)
jags.data <- list (J = J, y = y, sd = sd)
jags.inits <- function(){
list(theta=rnorm(J, 0, 100), mu=rnorm(1, 0, 100),
sigma=runif(1, 0, 100))
}
jags.parameters <- c("theta", "mu", "sigma")
##run model
schools.sim <- jags(data = jags.data, inits = jags.inits,
parameters = jags.parameters,
model.file = model.file,
n.chains = 3, n.iter = 10)
#note that n.iter should be higher
##set up in list
Cand.mods <- list(schools.sim)
Model.names <- "hierarchical model"
##other models can be added to Cand.mods
##to compare them to the top model
##model selection table
dictab(cand.set = Cand.mods, modnames = Model.names)
detach(package:R2jags)
## End(Not run)