Skip to contents

Function for performing posterior predictive checks on spOccupancy model objects.

Usage

ppcOcc(object, fit.stat, group, ...)

Arguments

object

an object of class PGOcc, spPGOcc, msPGOcc, spMsPGOcc, intPGOcc, spIntPGOcc, lfMsPGOcc, sfMsPGOcc, tPGOcc, or stPGOcc.

fit.stat

a quoted keyword that specifies the fit statistic to use in the posterior predictive check. Supported fit statistics are "freeman-tukey" and "chi-squared".

group

a positive integer indicating the way to group the detection-nondetection data for the posterior predictive check. Value 1 will group values by row (site) and value 2 will group values by column (replicate).

...

currently no additional arguments

Details

Standard GoF assessments are not valid for binary data, and posterior predictive checks must be performed on some sort of binned data.

Author

Jeffrey W. Doser doserjef@msu.edu,
Andrew O. Finley finleya@msu.edu

Value

An object of class ppcOcc that is a list comprised of:

fit.y

a numeric vector of posterior samples for the fit statistic calculated on the observed data when object is of class PGOcc or spPGOcc. When object is of class msPGOcc, spMsPGOcc, lfMsPGOcc, or sfMsPGOcc, this is a numeric matrix with rows corresponding to posterior samples and columns corresponding to species. When object is of class intPGOcc or spIntPGOcc, this is a list, with each element of the list being a vector of posterior samples for each data set. When object is of class tPGOcc or stPGOcc, this is a numeric matrix with rows corresponding to posterior samples and columns corresponding to primary sampling periods.

fit.y.rep

a numeric vector of posterior samples for the fit statistic calculated on a replicate data set generated from the model when object is of class PGOcc or spPGOcc. When object is of class msPGOcc, spMsPGOcc, lfMsPGOcc, or sfMsPGOcc, this is a numeric matrix with rows corresponding to posterior samples and columns corresponding to species. When object is of class intPGOcc or spIntPGOcc, this is a list, with each element of the list being a vector of posterior samples for each data set. When object is of class tPGOcc or stPGOcc, this is a numeric matrix with rows corresponding to posterior samples and columns corresponding to primary sampling periods.

fit.y.group.quants

a matrix consisting of posterior quantiles for the fit statistic using the observed data for each unique element the fit statistic is calculated for (i.e., sites when group = 1, replicates when group = 2) when object is of class PGOcc or spPGOcc. When object is of class msPGOcc, spMsPGOcc, lfMsPGOcc, or sfMsPGOcc, this is a three-dimensional array with the additional dimension corresponding to species. When object is of class intPGOcc or spIntPGOcc, this is a list, with each element consisting of the posterior quantile matrix for each data set. When object is of class tPGOcc or stPGOcc, this is a three-dimensional array with the additional dimension corresponding to primary sampling periods.

fit.y.rep.group.quants

a matrix consisting of posterior quantiles for the fit statistic using the model replicated data for each unique element the fit statistic is calculated for (i.e., sites when group = 1, replicates when group = 2) when object is of class PGOcc or spPGOcc. When object is of class msPGOcc, spMsPGOcc, lfMsPGOcc, or sfMsPGOcc, this is a three-dimensional array with the additional dimension corresponding to species. When object is of class intPGOcc or spIntPGOcc, this is a list, with each element consisting of the posterior quantile matrix for each data set. When object is of class tPGOcc or stPGOcc, this is a three-dimensional array with the additional dimension corresponding to primary sampling periods.

The return object will include additional objects used for standard extractor functions.

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 = 0, var = 2.72),
                   alpha.normal = list(mean = 0, var = 2.72))
# Initial values
inits.list <- list(alpha = 0, beta = 0,
                   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%

# Posterior predictive check
ppc.out <- ppcOcc(out, fit.stat = 'chi-squared', group = 1)
summary(ppc.out)
#> 
#> Call:
#> ppcOcc(object = out, fit.stat = "chi-squared", group = 1)
#> 
#> Samples per Chain: 5000
#> Burn-in: 4000
#> Thinning Rate: 1
#> Number of Chains: 1
#> Total Posterior Samples: 1000
#> 
#> Bayesian p-value:  0.378 
#> Fit statistic:  chi-squared