www.phidot.org

Hi, I'm trying to run some CJS models in marked to look at apparent survival in a bird population. The population has ~190 nest boxes where we band only juveniles at fledgling stage at every breeding season.

I ran some models in Program MARK and it seems like we have a 5 age classes structure. I am now trying to run the models in marked for reproducibility:

Here is the code I used to create 5 age classes in the design matrix:

Code: Select all

ddl <- make.design.data(cjs_proc, parameters = list(Phi = list(age.bins = c(0, 1, 2, 3, 4, Inf))), right = FALSE)
levels(ddl$Phi$age) <- c("1","2","3","4","5+")


And then I created the 6 age classes for detection p:

Code: Select all

new_age_categories <- cut(age_levels_numeric, breaks = c(0, 1, 2, 3, 4, 5, Inf), labels = c("1", "2", "3", "4", "5", "6+"), right = FALSE)
levels(ddl$p$age) <- new_age_categories

Is this the correct way of defining age-specific classes for Phi and p?

And then I run the CJS models:

Code: Select all

cjs_model <- mark(cjs_proc, ddl, model.parameters = list(
  Phi = list(formula=~time:age),
  p = list(formula=~time:age)
), output = FALSE)

1. Now I want to create a more "complex" age structure by regrouping age classes. Juveniles = 0 to 1 / adults = 2-4 / older adults = 5+. I want to test the hypothesis that only juvenile survival fluctuates in time and adults have different survival but stable (constant in time). How should I proceed?

2. For the detection (p), I know there is a monitoring issue from 2017. How can I make detection age-dependent but different between before and after 2017?

Thanks for your advice

Adrien

Clearly you didn't try any of that code because it all failed. You are confusing marked and RMark. Did you meant to post this under RMark sub-forum? The syntax close but there are differences like the model fitting is crm and not mark.

The binning is included in some functions but it is more flexible to use R functions to do the same thing. I do it that way and create new variables so you can see if is working correctly. I show this below using the tail function.
Note that survival is an interval parameter and age and time are labelled based on age/time at beginning of interval. Survival from age 0 to 1 is labelled 0. Whereas, p is an occasion parameter and for cjs model if they are released at age 0, that initial capture is not modelled and they can only be caught as 1 year olds. Now 0 is still in the ddl$p in marked because it uses the equivalent of square (rectangular) pims that all start at the beginning time and until they are released say at time 2, the initial values are not used.

Here is code for marked doing what you requested using the dipper data.

library(marked)
data(dipper)

dp=process.data(dipper,model="cjs",begin.time=2014)
ddl <- make.design.data(dp)
ddl$Phi$agebin=cut(ddl$Phi$Age,c(0, 1, 2, 3, 4, 5, Inf), right =FALSE)
ddl$p$agebin=cut(ddl$p$Age,c(-Inf,2,3,4,5,Inf),right=FALSE)
levels(ddl$Phi$agebin) <- c("0","1","2","3","4","5+")
levels(ddl$p$agebin) <- c("1","2","3","4","5+")

tail(ddl$p)
tail(ddl$Phi)

cjs_model <- crm(dp, ddl, model.parameters = list(
Phi = list(formula=~time:agebin),
p = list(formula=~time:agebin)))
cjs_model

ddl$Phi$agebin2=cut(ddl$Phi$Age,c(0, 2, 5, Inf), right =FALSE)
levels(ddl$Phi$agebin2)=c("Juveniles","Adults","Older Adults")

tail(ddl$Phi)

ddl$p$after2017=ifelse(ddl$p$Time>2,1,0)
tail(ddl$p)

cjs_model2 <- crm(dp, ddl, model.parameters = list(
Phi = list(formula=~time:agebin2),
p = list(formula=~after2017*agebin)))
cjs_model2

Hi jlaake,

Thank you very much for the quick reply.

Sorry for the confusion I was working with both a script in Rmark and marked simultaneously.

The code now works and computes apparent survival for 6 age classes (0-5+).

for the part:

Code: Select all

ddl$Phi$agebin2=cut(ddl$Phi$Age,c(0, 2, 5, Inf), right =FALSE)
levels(ddl$Phi$agebin2)=c("Juveniles","Adults","Older Adults")

If I want my "juveniles" to be only survival at age 0, would it be:

Code: Select all

ddl$Phi$agebin2=cut(ddl$Phi$Age,c(0, 1, 5, Inf), right =FALSE)

Thank you!

Yes but you shouldn't have to ask that question. Just use it and look at the results in the ddl dataframe as I suggested. --jeff