unexpected results for age structured models

[rmaluw]

Hey folks,

I am having some trouble understanding why I am getting some particular results. I have found various iterations of what I am trying to do on the forum before, but not the issue that I am running into. Apologies if there is a relevant post that answers this that I did not find.

I have a 27 year dataset for a birds that we survey (see snippet of data below). If the bird was seen/banded as a nestling, it is marked as a juvenile, if banded as an adult, then it is marked as older. Juveniles that are seen the following year are then considered adults/older. My main questions of interest are pretty standard: does survival vary between first years and adults, and does survival vary over time?

Here is where I am confused. When I run the age model (with detection varying with time), I get first year survival as constant (as expected). However, while the survival for 1+ years (adults) is overall higher as expected, it seems to change with the age of the adult, so clearly I have set something up incorrectly. I suspect it is something very simple, but I cannot for the life of me figure out what it is. Any help would be appreciated. Thank you!

Code: Select all

su.pr <- process.data(surv, model = "CJS",
                      age.var = 1,
                      initial.ages = c(0,1),
                      groups = ("age"))


su.ddl = make.design.data(su.pr,
                          parameters = list(Phi = list(age.bins = c(0,1,27)),
                                            p = list(age.bins = c(1,2,27))),
                          right = FALSE)

run.su = function(){
 
  Phi.dot = list(formula = ~ 1)
  Phi.age = list(formula = ~ age)
  Phi.time = list(formula = ~ -1 + time)
  Phi.age.time = list(formula = ~ age + time)
 
  p.dot = list(formula = ~ 1)
  p.time = list(formula = ~ -1 + time)
 
  su.model.list = create.model.list("CJS")
 
  su.results = mark.wrapper(su.model.list,
                            data = su.pr, ddl = su.ddl)
 
 
  return(su.results)
}

su.res <- run.su()

su.res


[jlaake]

Look at the design data for Phi. My guess is that age is not being binned. Possibly confusion between variable named age in data which is the initial age and the design data field also named age which is the age of the animal throughout the sampling. Always look at your design data before you start fitting models or if there are problems.

[rmaluw]

This ended up being exactly the issue. I ended up adding a new field in the design data based off the existing Age field, similar to the kestrel example posted a few years ago, and that solved the issue. Thanks!

Code: Select all

su.ddl$Phi$adult=0
su.ddl$Phi$adult[su.ddl$Phi$Age>=1]=1


[rmaluw]

I have another question about this - I ran an age+time model (survival varying for both adults and juveniles with time), and got a rather odd result. Juvenile survival is consistently lower, but adult and juvenile survival is correlated to a degree that I don't believe. The ddl and PIM look correct to me, so I am unsure where I went wrong. Note that age is assigned according to the 'adult' group I added, which is the last column in the ddl.

Any ideas? I am going to run this in MARK to compare, since I think this is probably an issue arising from something I did incorrectly in my code.

Thank you!



PIM:



ddl:

[cooch]

I have another question about this - I ran an age+time model (survival varying for both adults and juveniles with time), and got a rather odd result. Juvenile survival is consistently lower, but adult and juvenile survival is correlated to a degree that I don't believe.

It would seem you don't understand additive models. They enforce a parallelism (on the logit scale, say) which, if you think about it, makes the estimates perfectly correlated. As discussed at length in chapter 6, estimates for one group in an additive model are simply the set of estimates for the other group plus an additive constant.

As such, such perfect parallelism almost never consitutes a 'meaningful' biological model. Its simply a structural model that is intermediate between a model with full 'group' x time interaction, and a model without group. It may be more parsimonious - since it has fewer parameters than the model with interactions. Biological meaning? Generally zilch. Applied utility? As one of a model set that you're averaging over, potentially worth exploring.

[jlaake]

I'm not trying to pile on here but only to state that RMark was designed to help save time and prevent errors in construction of PIMs and design matrices. It is absolutely essential that you understand what the R formula is specifying in terms of the design matrix and model for your data. There is ample documentation in R on formulae and I have stated in the RMark documentation that you should use model.matrix with your formula and design data so you can understand the design matrix you are creating. I also recommend using and learning the MARK interface before using MARK. Without a good understanding of MARK (refer to Cooch and White ebook) it is easy to misunderstand your results.

[rmaluw]

Sorry about the late reply - I wanted to reread Ch 6. It seems I left out the interaction term for the full model, which in hindsight seems obvious. I am going through and running in MARK as well, so I can compare between the two since the that interface is much more intuitive for me.

These are fair criticisms, and a good motivation to go back through and relearn MARK. I have spent considerable time learning and using MARK (not RMark) in the past, but have not used it since. I was hoping to get through some exploratory analysis by brushing up on relevant chapters and not having to learn it in its entirety from scratch, but obviously that limits my understanding of what I am doing. I appreciate the responses.