probitCJS: Perform MCMC analysis of a CJS model

In marked: Mark-Recapture Analysis for Survival and Abundance Estimation

Perform MCMC analysis of a CJS model

Description

Takes design data list created with the function make.design.data for model "probitCJS" and draws a sample from the posterior distribution using a Gibbs sampler.

Usage

probitCJS(
  ddl,
  dml,
  parameters,
  design.parameters,
  burnin,
  iter,
  initial = NULL,
  imat = NULL
)

Arguments

ddl	list of dataframes for design data; created by call to make.design.data
dml	list of design matrices created by create.dm from formula and design data
parameters	A model specification list with a list for Phi and p containing a formula and optionally a prior specification which is a named list. See 'Priors' section below.
design.parameters	Specification of any grouping variables for design data for each parameter
burnin	number of iteration to initially discard for MCMC burnin
iter	number of iteration to run the Gibbs sampler for following burnin
initial	A named list (Phi,p). If null and imat is not null, uses cjs.initial to create initial values; otherwise assigns 0
imat	A list of vectors and matrices constructed by process.ch from the capture history data

Value

A list with MCMC iterations and summarized output:

beta.mcmc	list with elements Phi and p which contain MCMC iterations for each beta parameter
real.mcmc	list with elements Phi and p which contain MCMC iterations for each real parameter
beta	list with elements Phi and p which contain summary of MCMC iterations for each beta parameter
reals	list with elements Phi and p which contain summary of MCMC iterations for each real parameter

Prior distribution specification

The prior distributions used in probitCJS are multivariate normal with mean mu a and variance tau*(X'X)^(-1) on the probit scale for fixed effects. The matrix X is the design matrix based on the model specification (located in parameters$Phi$formula and parameters$p$formula respectively). Priors for random effect variance components are inverse gamma with shape parameter 'a' and rate parameter 'b'. Currently, the default values are mu=0 and tau=0.01 for phi and p parameters. For all randome effects deault values are a=2 and b=1.0E-4. In addition to the variance component each random effect can be specified with a known unscaled covariance matrix, Q, if random effects are correlated. For example, to obtain a random walk model for a serially correlated effect see Examples section below. To specify different hyper-parameters for the prior distributions, it must be done in the parameters list. See the Examples section for changing the prior distributions. Note that a and b can be vectors and the Qs are specified via a list in order of the random effects specified in the model statements.

Author(s)

Devin Johnson

Examples

# This example is excluded from testing to reduce package check time
# Analysis of the female dipper data
data(dipper)
dipper=dipper[dipper$sex=="Female",]
# following example uses unrealistically low values for burnin and 
# iteration to reduce package testing time
fit1 = crm(dipper,model="probitCJS",model.parameters=list(Phi=list(formula=~time),
 p=list(formula=~time)), burnin=100, iter=1000)
fit1
# Real parameter summary
fit1$results$reals

# Changing prior distributions:
fit2 = crm(dipper,model="probitCJS",
   model.parameters=list(
     Phi=list(formula=~time, prior=list(mu=rep(0,6), tau=1/2.85^2)),
     p=list(formula=~time, prior=list(mu=rep(0,6), tau=1/2.85^2))
   ), burnin=100, iter=1000)
fit2
# Real parameter summary
fit2$results$reals

# Creating a Q matrix for random walk effect for 6 recapture occasions
A=1.0*(as.matrix(dist(1:6))==1)
Q = diag(apply(A,1,sum))-A

# Fit a RW survival process
fit3 = crm(dipper,model="probitCJS",
   model.parameters=list(
     Phi=list(
       formula=~(1|time), 
       prior=list(mu=0, tau=1/2.85^2, re=list(a=2, b=1.0E-4, Q=list(Q)))
       ),
     p=list(formula=~time, prior=list(mu=rep(0,6), tau=1/2.85^2))
   ), burnin=100, iter=1000)
fit3
# Real parameter summary (Not calculated for random effects yet)
fit3$results$reals

marked documentation built on Oct. 19, 2023, 5:06 p.m.