I want to investigate three fixed group effects as covariates in a Cormack-Jolly-Seber style survival analysis using Jags in R. I have successfully developed a separate model for each covariate; however, I don't know whether it is appropriate to compare the apparent survival estimates between models, so I want to create one model with all three covariates.

My dataset:

-1100 wild birds over 4 years (I know it is a relatively short timespan for such a study).

-For each bird, I have 4 data points (one per year). Either the bird was seen (1) or unseen (0).

-I currently only have sex data for a subset of my dataset (411 birds right now, this will be a higher number in a few months)

I want to know how survival rates vary by sex (male or female), location (site A or site B), and age (juvenile or adult). My main prediction is that survival will vary between sexes at one site but not at the other. So if I cannot incorporate age into a model with the other covariates, it is not the end of the world.

Specific questions:

  1. How can I compare specific parameter estimates between models? Or even within models (I created a few models where I allowed the apparent survival (phi) to vary by year, and so I have 3 values of apparent survival and do not know how to compare them even there).

  2. How can I incorporate multiple categorical covariates into CJS in jags?

Any advice? I am also open to other methods of analysis, I know there are many methods out there. Jags just seemed the most appropriate/easiest to learn when I began this process.

Below is the model I developed to look at sex as a fixed group effect and the code I used to prepare priors. Not shown below: I created two capture-recapture histories (one for female and one for male) and then bound them together and made a separate matrix that coded whether each entry was male or female. That is what I call "dat" and use as my data matrix.

model {

# Priors and constraints 
for (i in 1:nind){ 
    for (t in f[i]:(n.occasions-1)){ 
        logit(phi[i,t]) <- beta[group[i]] + gamma[t] 
        p[i,t] <- p.g[group[i]] 
        } #t 
    } #i 
# for survival parameters 
for (t in 1:(n.occasions-1)){
    gamma[t] ~ dnorm(0, 0.01)I(-10, 10)         # Priors for time effects 
    phi.g1[t] <- 1 / (1+exp(-gamma[t]))         # Back-transformed survival of males 
    phi.g2[t] <- 1 / (1+exp(-gamma[t]-beta[2])) # Back-transformed survival of females 
beta[1] <- 0                                    # Corner constraint 
beta[2] ~ dnorm(0, 0.01)I(-10, 10)              # Prior for difference in male and female survival 

# for recapture parameters 
for (u in 1:g){ 
    p.g[u] ~ dunif(0, 1)                        # Priors for group-spec. recapture 

# Likelihood 
for (i in 1:nind){ 
    # Define latent state at first capture 
    z[i,f[i]] <- 1 
    for (t in (f[i]+1):n.occasions){
        # State process 
        z[i,t] ~ dbern(mu1[i,t]) 
        mu1[i,t] <- phi[i,t-1] * z[i,t-1] 

        # Observation process 
        y[i,t] ~ dbern(mu2[i,t]) 
        mu2[i,t] <- p[i,t-1] * z[i,t] 
        } #t 
    } #i 


jags.data <- list(y = dat, f = f, nind = dim(dat)[1], 
                  n.occasions = dim(dat)[2], 
                  g = length(unique(group)), group = group) 
# Initial values
known.state.cjs <- function(ch){
   state <- ch
   for (i in 1:dim(ch)[1]){
      n1 <- min(which(ch[i,]==1))
      n2 <- max(which(ch[i,]==1))
      state[i,n1:n2] <- 1
      state[i,n1] <- NA
   state[state==0] <- NA
n.occasions = dim(dat)[2]
inits <- function(){list(z = known.state.cjs(dat), beta = c(NA, rnorm(1)), 
                         gamma=rnorm(n.occasions-1), p.g = runif(length(unique(group)), 0, 1) )}
# Parameters monitored
parameters <- c("phi.g1", "phi.g2", "p.g", "beta") 

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.