Oddly large beta and real CIs

[kquail]

I ran a Burnham model with 4 grouping factors: site, sex, feed treatment, and age (juvenile or adult). I used code from this example https://www.rdocumentation.org/packages/RMark/versions/3.0.0/topics/Burnham to bin the design data for S and p such that juveniles become adults in their second year.

Code: Select all

#binning for the S parameter
band.ddl=add.design.data(data = band.process,
                       ddl = band.ddl,
                       parameter="S",
                       type = "age",       
                       bins = c(0,1,11),    #0 and 1+
                       right = FALSE,      #bin intervals are open on the right and open on the left
                       name = "age",
                       replace = TRUE)     
# do the same to the p parameter
band.ddl=add.design.data(data = band.process,
                       ddl = band.ddl,
                       parameter="p",
                       type = "age",
                       bins = c(0,1,11),
                       right = FALSE,
                       name = "age",
                       replace = TRUE)


I then ran a model { S(~time * fed + sex + age + site)p(~age)r(~1)F(~1) } where both S and p could vary with age. I found that the beta estimate CIs for p [-859, 860] and F [-1912,1948] are huge and there is only 1 real estimate for p for adults, whereas I would expect another also for juveniles. When trying different models, by varying which parameters are allowed to vary, it seems like F has huge CIs for every one, but p has huge CIs when it is allowed to vary by age. The example on the website where I found the code to change binning allowed S and F to vary by age and the results look normal, but when I instead allowed p to vary by age, the same thing occurred (large beta CIs and only 1 real estimate for p, although F looked fine). For F maybe the reason is that F essentially equals 1? For p though maybe this has something to do with p occurring at the time of the sampling occasion, while the other parameters technically occur between sampling occasions? Maybe there is some issue with exactly when the individuals change over from juvenile to adult? Or perhaps I am missing something else?

I’m hoping you can recreate my issue with the LD1 example dataset. I apologize if this is not an RMark-specific question. Thanks in advance!

[jlaake]

Not certain how I could recreate it with a different set of data. I suggest you do the following.

1. Look at ddl$p to see if your age binning worked as you thought. Always look at ddl for the parameter to identify problems.

2. Start with a simple model (~1) for each parameter and build the model upwards. Is this a harvested population? If harvest rate is varying then I'd expect recovery to vary as well. Beta CI values are wide when parameters are at boundaries and parameters are often at boundaries with sparse data and overly complex models in relation to the amount of data. Do not try to interpret results of complex models until you have evaluated simpler models as well and determine that the complex models are the best.

Jeff

[kquail]

Hi Jeff, thanks for the suggestions.
When processing the data, I set initial ages to adult = 1 and juvenile = 0. In ddl$p I see that there are no rows where Age =0. 'age' is always [1,11] --the bin for adult ages. Similarly, there are no rows where cohort = time. These such rows do exist in ddl$S.

I suppose this is because we can't calculate p for the very first occasion where an animal is encountered, only after? So perhaps this is actually a MARK question, but does that mean that I can't calculate p separately for juveniles if they become an adult by the next occasion after they were first encountered? Any ideas how to get around that?

In response to your second suggestion: Yes, this is a harvested population, and harvest might vary slightly across years. That is a good point that r possibly also varies by time and site. I had not thought to include that in my model set.
I have over 8000 individuals in the dataset, but sample size is definitely not equally distributed across all the groups. For example, there are 3-4x more individuals marked as juveniles than as adults and 2-3x more individuals from Site1 than Site2.
I did not specify in the first post, but I did, indeed create a full set of models of varying complexity. It just turned out that the best models based on AIC were the most complex (i.e. all 5 parameters plus the interaction included for S, and either p.dot (top model) or p.age (second model)). The p.dot model suggests that p is about 0.48, so it doesn't seem to be at a boundary. On the other hand, F =1. Site, Sex, and treatment are all important parameters, but age might be uninformative for both S and p. The main reason I want a separate p for adults and juveniles is so that I have a value for a follow-up analysis on only juveniles where I hoped to use the juvenile p from this analysis to help more accurately estimate S (which I expect to be rather low in that analysis).

[jlaake]

When processing the data, I set initial ages to adult = 1 and juvenile = 0. In ddl$p I see that there are no rows where Age =0. 'age' is always [1,11] --the bin for adult ages. Similarly, there are no rows where cohort = time. These such rows do exist in ddl$S.

I suppose this is because we can't calculate p for the very first occasion where an animal is encountered, only after? So perhaps this is actually a MARK question, but does that mean that I can't calculate p separately for juveniles if they become an adult by the next occasion after they were first encountered? Any ideas how to get around that?

Correct. Phi/S and F are interval parameters which cover from the start occasion to the next occasion (e.g. survival of a juvenile to being an adult one year later). Whereas p is an occasion parameter which is relevant to a specific occasion. The live-dead models are based on releases of animals so it is only time after released that is modeled. This is the same as CJS models of survival.

Only way around it is to use a Jolly-Seber model which includes first capture in the likelihood - like robust models and closed models where one of the goals is to estimate abundance. You could use the live capture data with POPAN if you want. I don't know of any live-dead models that use the POPAN structure in MARK but I have not kept up with all of the MARK developments since I retired. Maybe someone else can chime in on that point.

With regard to you other data set, consider combining the two sets of data in a POPAN formulation if that is appropriate. But if your goal is to estimate survival, not sure why you will need p for juveniles for the other weaker data sets if you are using CJS type release model.

[kquail]

That makes sense. Thanks for suggesting to try those other models instead.

The thought process behind wanting p (and actually now that I think about it, probably also F) for juveniles from the full dataset is because the second analysis I want to run is with individuals marked as chicks. Presumably chick survival will be much lower and this datset has a smaller sample size, so we were thinking that fixing p at at the known value for juveniles (the age at which chicks are resighted) should make the calculation of S more accurate.

[jlaake]

You can estimate F for juveniles with the existing data set. It is an interval parameter like S. But rather than fixing p for juveniles like you are thinking I'd recommend analyzing the data sets together and sharing parameters that are supported. Fixing parameters at a value from a separate analysis ignores the uncertainty in that parameter. You would just be fooling yourself.

Jeff