Compute Widely Applicable Information Criterion for spOccupancy Model Objects
waicOcc.RdFunction for computing the Widely Applicable Information Criterion
(WAIC; Watanabe 2010) for spOccupancy model objects.
Arguments
- object
an object of class
PGOcc,spPGOcc,msPGOcc,spMsPGOcc,intPGOcc,spIntPGOcc,lfJSDM,sfJSDM,lfMsPGOcc,sfMsPGOcc,tPGOcc,stPGOcc,svcPGBinom,svcPGOcc,svcTPGBinom,svcTPGOcc, orintMsPGOcc,svcMsPGOcc,tMsPGOcc,stMsPGOcc,svcTMsPGOcc.- by.sp
a logical value indicating whether to return a separate WAIC value for each species in a multi-species occupancy model or a single value for all species.
- ...
currently no additional arguments
References
Watanabe, S. (2010). Asymptotic equivalence of Bayes cross validation and widely applicable information criterion in singular learning theory. Journal of Machine Learning Research, 11:3571-3594.
Gelman, A., J. B. Carlin, H. S. Stern, D. B. Dunson, A. Vehtari, and D. B. Rubin. (2013). Bayesian Data Analysis. 3rd edition. CRC Press, Taylor and Francis Group
Gelman, A., J. Hwang, and A. Vehtari (2014). Understanding predictive information criteria for Bayesian models. Statistics and Computing, 24:997-1016.
Author
Jeffrey W. Doser doserjef@msu.edu,
Andrew O. Finley finleya@msu.edu
Value
When object is of class PGOcc, spPGOcc, msPGOcc, spMsPGOcc, lfJSDM, sfJSDM, lfMsPGOcc, sfMsPGOcc, tPGOcc, stPGOcc, svcPGBinom, svcPGOcc, svcTPGOcc, svcTPGBinom, svcMsPGOcc, tMsPGOcc, stMsPGOcc, svcTMsPGOcc
returns a vector with three elements corresponding to
estimates of the expected log pointwise predictive density (elpd), the
effective number of parameters (pD), and the WAIC. When by.sp = TRUE for multi-species models, object is a data frame with each row corresponding to a different species. When object is
of class intPGOcc or spIntPGOcc, returns a data frame with
columns elpd, pD, and WAIC, with each row corresponding to the estimated
values for each data source in the integrated model.
Details
The effective number of parameters is calculated following the recommendations
of Gelman et al. (2014). Note that when fitting multi-species occupancy models with the
range.ind tag, it is not valid to use WAIC to compare a model that
uses range.ind (i.e., restricts certain species to a subset of the locations)
with a model that does not use range.ind (i.e., assumes all species can occur at
all locations in the data set) or that uses different range.ind values.
Examples
set.seed(400)
# Simulate Data -----------------------------------------------------------
J.x <- 8
J.y <- 8
J <- J.x * J.y
n.rep <- sample(2:4, J, replace = TRUE)
beta <- c(0.5, -0.15)
p.occ <- length(beta)
alpha <- c(0.7, 0.4)
p.det <- length(alpha)
dat <- simOcc(J.x = J.x, J.y = J.y, n.rep = n.rep, beta = beta, alpha = alpha,
sp = FALSE)
occ.covs <- dat$X[, 2, drop = FALSE]
colnames(occ.covs) <- c('occ.cov')
det.covs <- list(det.cov = dat$X.p[, , 2])
# Data bundle
data.list <- list(y = dat$y,
occ.covs = occ.covs,
det.covs = det.covs)
# Priors
prior.list <- list(beta.normal = list(mean = rep(0, p.occ),
var = rep(2.72, p.occ)),
alpha.normal = list(mean = rep(0, p.det),
var = rep(2.72, p.det)))
# Initial values
inits.list <- list(alpha = rep(0, p.det),
beta = rep(0, p.occ),
z = apply(data.list$y, 1, max, na.rm = TRUE))
n.samples <- 5000
n.report <- 1000
out <- PGOcc(occ.formula = ~ occ.cov,
det.formula = ~ det.cov,
data = data.list,
inits = inits.list,
n.samples = n.samples,
priors = prior.list,
n.omp.threads = 1,
verbose = TRUE,
n.report = n.report,
n.burn = 4000,
n.thin = 1)
#> ----------------------------------------
#> Preparing to run the model
#> ----------------------------------------
#> ----------------------------------------
#> Model description
#> ----------------------------------------
#> Occupancy model with Polya-Gamma latent
#> variable fit with 64 sites.
#>
#> Samples per Chain: 5000
#> Burn-in: 4000
#> Thinning Rate: 1
#> Number of Chains: 1
#> Total Posterior Samples: 1000
#>
#> Source compiled with OpenMP support and model fit using 1 thread(s).
#>
#> ----------------------------------------
#> Chain 1
#> ----------------------------------------
#> Sampling ...
#> Sampled: 1000 of 5000, 20.00%
#> -------------------------------------------------
#> Sampled: 2000 of 5000, 40.00%
#> -------------------------------------------------
#> Sampled: 3000 of 5000, 60.00%
#> -------------------------------------------------
#> Sampled: 4000 of 5000, 80.00%
#> -------------------------------------------------
#> Sampled: 5000 of 5000, 100.00%
# Calculate WAIC
waicOcc(out)
#> elpd pD WAIC
#> -110.839502 3.482209 228.643422