RMark Multistate Age Issues

[EmmaRigby]

Hello,
I've recently moved from Mark to RMark but have a couple of questions I'd be very grateful for help with!

I'm doing a multistate analysis of 651 individuals over 10 capture ocassions (Spring and Autumn for the years 2004-2008) and 6 locations. All individuals are grouped in the following categories (adult male, juvenile male, adult female, juvenile female), for example:

00A0BG0000 0 0 1 0;

All juveniles are born at roughly the same time and survival is thought to vary in their first year. This means if juveniles were caught in spring their first two capture occassions (spring and autumn) would be different but if they were not caught until autumn only their first capture occassion would be different:

1 1 2 2 2 2 2 2 2
...1 2 2 2 2 2 2 2
......1 1 2 2 2 2 2
.........1 2 2 2 2 2
............1 1 2 2 2
...............1 2 2 2
..................1 1 2
.....................1 2
........................1

I want to code juveniles as the above PIM and adults as just number 2, so I can then combine this with season and sex in my models to create more complicated analyses. The coding I've used so far is:

-------------------------------------------------------------------------------------
#CONVERT .INP FILE TO R FORMAT AND SPECIFY GROUPS:
wharfedale=convert.inp("Multi-state (seasons).inp", group.df=data.frame(sex=c(rep("Male",2),rep("Female",2)),age=rep(c("Adult","Juvenile"),2)))

#CREATE THE PROCESSED DATAFRAME:
wharfedale.process=process.data(wharfedale,model="Multistrata", begin.time=0, time.intervals=c(3.27,9.27,2.3,8.73,2.37,9.7,2.5,9.9,2.03), groups=c("sex","age"), age.var=2, initial.age=c(1,0))

#CREATE THE DESIGN DATA:
wharfedale.ddl=make.design.data(wharfedale.process)
-------------------------------------------------------------------------------------

but I'm unsure of how to proceed (age bins or something more complicated etc). I'm also a bit unsure of my begin.time and time.intervals. I've coded my time intervals in months as I want seasonal parameter estimates but this seems to affect the way my ages are displayed (i.e. in fractions of months!). Hopefully this won't be an issue though if I can get my PIM coding right?


I hope that made sense and thank you in advance for any help!

[jlaake]

You have made a great start and there is only a few things you haven't understood and mostly because I need to add some material to the appendix to clarify some of this.

Here is how it works. Animals are given an initial.age which is their age on the occasion of their first capture. Let's say that first occasion is k. Their age at occasion k+1 is initial.age+ time.interval[k]*age.unit and on occasion k+2 is initial.age+ time.interval[k]*age.unit + time.interval[k+1]*age.unit, etc. The default for age.unit is 1 and in your case you would want it to be 1/12 since your time intervals are in months and your aging is in years. Alternatively, you could leave age.unit=1 and express age in months such that juveniles are 0 and adults are 12+. Your choice and no difference other than labelling.

Now it is hopefully obvious that you have an issue with juveniles because they do not all have the same initial.age as those initially caught in spring are not the same age as those initially caught in fall. Therefore you need to add some more structure to your data to split these groups of juveniles to correctly assign an initial.age. I can work up some code that will help you do this and will send it off-list.

Then the only other issue you have would be to define bins for age for each parameter that will be age-dependent. You can do that easily with the call to make.design.data or add.design.data. If you used age in years, then you can simply add age.bins=c(0,1,10),right=FALSE to the list elements. The field right=FALSE tells the code to leave the interval open on the right except for the last interval. Thus the age intervals will be [0,1) and [1,10] which is anything less than a year and anything a year or greater. I chose 10 as an upper limit that should be ok for your case. Animals continue to age and if you had an adult (age 1) at the beginning of the study, then you can work out its maximum age or look at age in the design data field. Probably 5 will be sufficient. Now since you are using non-integer values for time.intervals and age.unit, you should look at the design data after it is created as you may need to make adjustments to your age.bins and initial.age as the value for one year may not come out exactly as 1.00...

The place where folks seem to get most confused is with initial.age and age.bins. The former is for the initial age of the animal whereas the latter are to define intervals of the age variable which is the age of the animal throughout the course of the study at each occasion.

Hopefully this clarifies the issue for you. If not feel free to respond with additional questions.

--jeff

[tasha]

Hello, I'm doing an analysis very similar to the one above.

The data is collected over 5 capture occasions (2008: Spring and Autumn, 2009: Spring and Autumn, 2010: Spring). All individuals are grouped as juveniles or adults at first capture. All juveniles are born at roughly the same time and survival is thought to vary in their first year. If juveniles were caught in both spring and autumn of the same year, they do not transition age classes until their third capture occasion. However if they are first caught in the fall, they transition age classes on their second capture occasion.

Here's how we have the model set up so far:

timeint=c(rep(c(0,0,0,0,3,0,0,0,0,9),2),0,0,0,0)
dog.process=process.data(dog,model=”HCRDMS”,strata.labels=c(‘s’,’m’,’l’,’S’,’M’,’L’), groups=c(“age”), initial.age=c(1,0), time.intervals=timeint)
dog.dll=make.design.data(dog.process)

Based on Jeff Laake’s response above, it seems that we need to “add some more structure to [the] data to split these groups of juveniles to correctly assign an initial.age.” We aren’t sure how to go about this, and would like to know if anyone can help us.

Thanks!

[jlaake]

It looks like you are close. Have you looked at the design data? That will let you know whether you have it setup correctly. The one problem I see is that it looks like you are using month as your time unit for intervals but for age you are using 1 year as a unit interval (initial.age=c(1,0)). You can either use age.unit=1/12 or change your initial ages to 12 and 0 so that age is in months. Then an adult would be anything >=12 and a juvenile would be <12. Where you may have problems is with juveniles that enter in spring versus autumn. One solution is to have 3 groups. Juv that enters in spring, juv that enters in autumn and adult. Then you can set up the initial ages so that juveniles are say 3 months old in spring and say 9 months old in autumn. Those are just examples, then each should be >12 months the following spring. Ultimately all you are trying to do here is to get the design data setup correctly with the proper age value at each time, so looking at your design data is essential. You'll want to bin age so it falls in 2 groups: juv and adult.

--jeff

[tasha]

Thanks for your reply, it was really helpful!

We changed our script so that the unit intervals for age.initial are in months as you suggested so that the unit in our capture intervals matches the unit in our age. We also recoded our data so that age at first capture reflects spring and fall entry of juveniles. Where animals age 0 are juveniles that enter in the spring, age 3 are juveniles that enter in the fall and animals that enter as adults in either spring or fall are age 12. We are now running into issues with how to bin the ages in the script so that 0 and 3 month olds are both coded as juveniles and anyone 12 and over are coded as adults. Below is a piece of the script that we tried to modify.

# create RMARK process data frame
timeint=c(rep(c(0,0,0,0,3,0,0,0,0,9),2),0,0,0,0);
#timeint=timeint[-length(timeint)] # remove last '1' from timeint
pdogf.process=process.data(pdogf,model="HCRDMS",
strata.labels=c('s','m','l','S','M','L'),
groups=c("age"),
initial.age=c(0,3,12),
time.intervals=timeint)

pdogf.ddl=make.design.data(pdogf.process, parameters=list(S=list(age.bins=c(0,12,36),right=FALSE)))

pdogf.ddl$S$agecl=pdogf.ddl$S$initial.age.class
pdogf.ddl$S$medium=0; pdogf.ddl$S$medium[pdogf.ddl$S$stratum %in% c('m','M')]=1
pdogf.ddl$S$large=0; pdogf.ddl$S$large[pdogf.ddl$S$stratum %in% c('l','L')]=1

Our Questions are:
1) Where and how should we write the definition of the age bins?
2) Do we need to define them separately for each of our parameters? (S, psi, p)?
3) How do we write age bins as a variable in the models?

Thank You very much!

[jlaake]

I couldn't quite work out what you were trying to do in your code, but it isn't initial.age that you want to use in the design data because it is static. What you want to use is age and Age which change over time. You can either use add.design.data or the cut function in R to make whatever bins you want. For example,

Code: Select all

library(RMark)
data(dipper)
dp=process.data(dipper)
ddl=make.design.data(dp)
ddl$p$age
ddl$p$Age
ddl$p$agebin=cut(ddl$p$Age,c(1,3,4,7),include.lowest=TRUE)
ddl$p$agebin


You use Age rather than age to bin because Age is numeric whereas age is a factor variable. Then you can use something like

Code: Select all

mark(dp,ddl,model.parameters=list(p=list(formula=~agebin)))


That answered question 1 and 3. Answer to question 2 is yes and note that age values will differ across parameters. S and Psi are interval parameters and age is the age at the beginning of the time interval. p is an occasion parameter and it is the age at that occasion.

[tasha]

Thank you . We tried to use add.design.data to create the bins as follows :

library(RMark)

pdogf=data.frame(ch=ch, age=as.factor(a$Age_FirstCapture), plotID=a$PitTag,stringsAsFactors=F, treatment=a$Treatment, edible8=a$Edible08, edible9=a$Edible09, pcf1=a$PCF_Ses1,pcf2=a$PCF_Ses2,pcf3=a$PCF_Ses3,pcf4=a$PCF_Ses4)
timeint=c(rep(c(0,0,0,0,3,0,0,0,0,9),2),0,0,0,0);
pdogf.process=process.data(pdogf,model="HCRDMS",
strata.labels=c('s','m','l','S','M','L'),
groups=c("age"),
initial.age=c(0,3,12),
time.intervals=timeint)
pdogf.ddl=make.design.data(pdogf.process)
pdogf.ddl=add.design.data(pdogf.process,pdogf.ddl,"S",type="age",bins=c(0,12,36),name="agecl")

However, when we look at the data using summary(pdogf.ddl$S) it does not make sense that individuals are evenly distributed among groups and there are values for age that should not be possible (e.g. 6)

group cohort age time
0 :60 1 :72 12:42 1 :18
3 :60 4 :54 3:36 4 :36
12:60 13:36 15:30 13:54
16:18 0:24 16:72
6:12
24:12
(Other):24

Furthermore, when we look at our data using View (pdogf.ddl$S) the cohorts don't appear to be entering at the right times. The 0 year olds are entering at both t=0 (spring) and t=4 (fall) when they should only enter in the spring. Similarly 3 month olds are coming in at time t=0 and t=4 when they should only be allowed to enter at time t=4 (the fall).

Do you have any suggestions as to how we can figure out what we did wrong?

Thank you.

[jlaake]

Not sure you did anything wrong. The design data are constructed for all possible combinations. It doesn't mean they are necessarily used. Other than incrementing the ages correctly, the make.design.data code doesn't know anything about what data it will get. Nowhere have you told it that certain ages only enter at certain times. It assumes that each member of a group can enter at any of the times so that is why there are the same number of entries. What I don't know is how this is going to interact with this model if it estimates abundance. I don't know a whole lot about that model. Are your critters only entering the population as young or at older ages as well from some other source. The other way you can handle it is to create a time-varying individual covariate rather than using groups but individual covariates can be problematic in cases where you are estimating abundance because you don't know the covariate values for those that are not caught.

--jeff