MSOcc_mod: Fit multi-scale occupancy model.

In StrattonCh/msocc: Package for Fitting and Analyzing Computationally Efficient Multi-scale Occupancy Models

Description

This function fits the Bayesian multi-scale occupancy model described by Dorazio and Erickson (2017) using the polya-gamma data augmentation strategy described by Polson et al. (2012). Note that this documentation assumes there are M sites, J_i samples within each site, and K_{ij} replicates from each sample.

Usage

msocc_mod(
  wide_data,
  site = list(model = ~1, cov_tbl),
  sample = list(model = ~1, cov_tbl),
  rep = list(model = ~1, cov_tbl),
  priors = list(site = list(mu0 = 0, Sigma0 = 9), sample = list(mu0 = 0, Sigma0 = 9),
    rep = list(mu0 = 0, Sigma0 = 9), a0 = 1, b0 = 1),
  num.mcmc = 1000,
  progress = T,
  print = NULL,
  seed = NULL,
  beta_bin = T
)

Arguments

wide_data	object of class data.frame containing site, sample, and PCR replicates in wide format. Column names should be site, sample, PCR1, PCR2, ... and contain no other columns.
site	object of class list containing the following elements: model formula describing the within site model. cov_tbl object of class data.frame containing site specific covariates; cov_tbl should have exactly one row for each site and have a column named 'site'.
sample	object of class list containing the following elements: model formula describing the within sample model. cov_tbl object of class data.frame containing sample specific covariates; cov_tbl should have exactly one row for each sample and have columns named 'site' and 'sample'.
rep	object of class list containing the following elements: model formula describing the within replicate model. cov_tbl object of class data.frame containing replicate specific covariates; cov_tbl should have exactly one row for each sample if replicate specific covariates are aggregated at the sample level and contain columns named 'site' and 'sample'. Otherwise, cov_tbl should have exactly one row for each replicate and contain columns name 'site', 'sample', and 'rep'.
priors	object of class list containing the following elements: site object of class list containing the following elements: mu0 prior mean for site-level regression coefficients Sigma0 prior covariance matrix for site-level regression coefficients sample object of class list containing the following elements: mu0 prior mean for sample-level regression coefficients Sigma0 prior covariance matrix for sample-level regression coefficients rep object of class list containing the following elements: mu0 prior mean for replicate-level regression coefficients Sigma0 prior covariance matrix for replicate-level regression coefficients a0,b0 numeric shape parameters for beta prior on probability parameters if sampled directly (see beta_bin argument)
num.mcmc	number of MCMC samples
progress	should sampling progress be printed?
print	interval for printing; defaults to 5 percent of num.mcmc of left NULL
seed	optional seed for reproducible samples
beta_bin	optional; should a beta-binomial sampler be used when possible? This option is considerably faster.

Details

This function fits the multi-scale occupancy model described by Dorazio and Erickson (2017). However, this function implements a fully Bayesian sampler based on the data augmentation strategy described by Polson et al. (2012)

Value

object of class list containing the following elements:

• beta an object of class matrix of samples from the joint posterior distribution of the regression coefficients at the site level

• psi an object of class numeric of samples from the joint posterior distribution of the site-level presence probability, psi. Note that this is only returned if beta_bin is TRUE and site$model = ~ 1

• alpha an object of class matrix of samples from the joint posterior distribution of the regression coefficients at the sample level

• theta an object of class matrix of samples from the joint posterior distribution of the sample-level presence probability, theta. Note that this is only returned if beta_bin is TRUE and sample$model = ~ 1 or sample$model = ~ site

• delta an object of class matrix of samples from the joint posterior distribution of the regression coefficients at the replicate level

• p an object of class numeric of sample from the posterior distribution of the replicate level detection probability. NOte that this is only returned if beta_bin is TRUE and rep$model = ~ 1.

• model.info an object of class list containing the following elements:
  

  • X design matrix for site level predictors

  • W design matrix for sample level predictors

  • V design matrix for replicate level predictors

  • M number of sites

  • J vector of number of samples per site

  • K vector of number of replicates per site-sample combination

  • z matrix of posterior samples of latent site-presence z

  • z.vec matrix of posterior samples of latent site-presence stretched across samples

  • A matrix of posterior samples of latent sample-presence A

  • Y vector of binomial responses (aggregated at sample level)

  • y vector of Bernoulli responses (stretched across site-sample combination)

  • site_mod model statement for site-level predictors

  • samp_mod model statement for sample-level predictors

  • rep_mod model statement for replicate-level predictors

  • num.mcmc number of MCMC samples run

  • beta_bin was beta-binomial sampler implemented if possible?

  • df empty data.frame of design in long format

Examples

data(fung)

# prep data
fung.detect <- fung %>%
  dplyr::select(1:4)

site.df <- fung %>%
  dplyr::select(-sample, -pcr1, -pcr2) %>%
  dplyr::distinct(site, .keep_all = TRUE) %>%
  dplyr::arrange(site)

sample.df <- fung %>%
  dplyr::select(-pcr1, -pcr2) %>%
  dplyr::arrange(site, sample)

# fit intercept model at all three levels use beta-binomial sampler
fung_mod1 <- msocc_mod(wide_data = fung.detect, progress = T,
                       site = list(model = ~ 1, cov_tbl = site.df),
                       sample = list(model = ~ 1, cov_tbl = sample.df),
                       rep = list(model = ~ 1, cov_tbl = sample.df), # covariates aggregated at sample level
                       num.mcmc = 1000, beta_bin = T)

# model sample level occurence by frog density
fung_mod2 <- msocc_mod(wide_data = fung.detect, progress = T,
                       site = list(model = ~ 1, cov_tbl = site.df),
                       sample = list(model = ~ frogs, cov_tbl = sample.df),
                       rep = list(model = ~ 1, cov_tbl = sample.df),
                       num.mcmc = 1000, beta_bin = T)

StrattonCh/msocc documentation built on Dec. 22, 2020, 2:51 a.m.