boral (version 0.4)

get.more.measures: Additional Information Criteria for boral models

Description

Calculates some information criteria beyond those from get.measures for an boral model, although this set of criteria takes much longer to compute!

Usage

get.more.measures(y, X = NULL, family, trial.size = NULL, num.lv, fit.mcmc, 
	site.eff, verbose = TRUE)

Arguments

y

The response matrix that the boral model was fitted to.

X

The model matrix used in the boral model. Defaults to NULL, in which case it is assumed no model matrix was used.

family

Either a single element, or a vector of length equal to the number of columns in $y$. The former assumes all columns of $y$ come from this distribution. The latter option allows for different distributions for each column of $y$. Elements can be one of "b

trial.size

Either equal to NULL, a single element, or a vector of length equal to the number of columns in $y$. If a single element, then all columns assumed to be binomially distributed will have trial size set to this. If a vector, different trial sizes are allowe

num.lv

The number of latent variables used in the fitted boral model.

fit.mcmc

All MCMC samples for the fitted boral model, as obtained from JAGS. These can be extracted by fitting an boral model using boral with save.model = TRUE, and then applying as.mcmc on

site.eff

A logical value indicating whether to include row effects were included in the model.

verbose

If TRUE, a notice is printed every 100 samples indicating progress in calculation of the marginal log-likelihood. Defaults to TRUE.

Value

• A list with the following components:
• marg.aic
  AIC (using on the marginal log-likelihood) evaluated at posterior mode.
• marg.bic
  BIC (using on the marginal log-likelihood) evaluated at posterior mode.
• marg.dic
  DIC based on the marginal log-likelihood.
• marg.waic
  WAIC based on the marginal log-likelihood.
• all.marg.logLik
  The marginal log-likelihood evaluated at all MCMC samples. This is done via repeated application of calc.marglogLik.
• num.params
  Number of estimated parameters used in the fitted model.

Warning

Using information criteria for variable selection should be done with extreme caution, for two reasons: 1) The implementation of these criteria are both heuristic and experimental. 2) Deciding what model to fit for ordination purposes should be driven by the science. For example, it may be the case that criteria suggests a model with 3 or 4 latent variables. However, if we interested in visualizing the data for ordination purposes, then models with 1 or 2 latent variables are far more appropriate. As an another example, whether or not we include row (site) effects when ordinating multivariate abundance data depends on if we are interested in differences between sites in terms of abundance (site.eff = FALSE) or in terms of species composition (site.eff = TRUE).

Details

Currently, four information criteria has been implemented in this function: 1) AIC (using the marginal likelihood) evaluated at the posterior mode; 2) BIC (using the marginal likelihood) evalulated at the posterior mode; 3) Deviance information criterion (DIC) based on the marginal log-likelihood; 4) Widely Applicable Information Criterion (WAIC, Watanabe, 2010) based on the marginal log-likelihood. Since flat priors are used in fitting boral models, then the posterior mode should be approximately equal to the maximum likelihood estimates.

All four criteria require computing the marginal log-likelihood across all MCMC samples. This takes a very long time to run, since Monte Carlo integration needs to be performed for all MCMC samples. Consequently, this function is currently not implemented as an argument in main boral fitting function, unlike get.measures which is available via the calc.ics = TRUE argument.

The two main differences between the criteria and those returned from get.measures are:

• The AIC and BIC computed here are based on the log-likelihood evalulated at the posterior mode, whereas the AIC and BIC fromget.measuresare evalulated at the posterior median. The posterior mode and median will be quite close to one another if the component-wise posterior distributions are unimodal and symmetric. Furthermore, given uninformative priors are used, then both will be approximate maximum likelihood estimators.
• The DIC and WAIC computed here are based on the marginal log-likelihood, whereas the DIC and WAIC fromget.measuresare based on the conditional log-likelihood. Criteria based on the two types of log-likelihood are equally valid, and to a certain extent, which one to use depends on the question being answered i.e., whether to condition on the latent variables or treat them as "random effects" (see discussions in Spiegelhalter et al. 2002, and Vaida and Blanchard, 2005). Having said that, there is evidence to suggests, for models with overdispersed count data, conditional DIC/WAIC may not perform as well as than marginal DIC/WAIC for overdispered abundance data (Millar, 2009).

In our very limited experience, we found BIC evalulated at the posterior mode tends to be quite stable, whereas marginal DIC and WAIC tend to overfit the number of latent variables.

References

• Millar, R. B. (2009). Comparison of hierarchical Bayesian models for overdispersed count data using DIC and Bayes' factors. Biometrics, 65, 962-969.
• Spiegelhalter, D. J., Best, N. G., Carlin, B. P., & Van Der Linde, A. (2002). Bayesian measures of model complexity and fit. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 64, 583-639.
• Vaida, F., & Blanchard, S. (2005). Conditional Akaike information for mixed-effects models. Biometrika, 92, 351-370.
• Watanabe, S. (2010). Asymptotic equivalence of Bayes cross validation and widely applicable information criterion in singular learning theory. The Journal of Machine Learning Research, 11, 3571-3594.

See Also

get.measures for several information criteria which take less time to compute, and are automatically implemented in boral with calc.ics = TRUE.

Examples

library(mvabund) ## Load a dataset from the mvabund package
data(spider)
y <- spider$abun
n <- nrow(y); p <- ncol(y); 
    
spider.fit.nb <- boral(y, family = "negative.binomial", num.lv = 2, 
     site.eff = TRUE, save.model = TRUE, calc.ics = TRUE)

all.mcmc <- as.mcmc(spider.fit.nb$jags.model)

## WATCH OUT! The following takes a very long time to run!
get.more.measures(y, family = "negative.binomial", 
     num.lv = 2, fit.mcmc = all.mcmc, site.eff = TRUE)