Question about application of GLMM with Poisson

Question

I'm working on revising stats for a manuscript involving male reproductive success of deer. We measured three variables (body size, antler size, and age) of male deer in a captive population on an annual basis over a 6 year period, and goal was to determine the relative influence of each variable on annual reproductive success (i.e., number of fawns produced each year / not lifetime reproductive success). I understand that my predictors are collinear and the issues created; however, there really is no way around including collinear predictors in our model as PCA's, etc. would essentially destroy the core question of our research. We plan to use model averaging to evaluate predictors in the end so ran global model first. Here is current code:

global.model = glmer(Fawn ~ Age + I(Age^2) + Age*AvgAge + BodySize + I(BodySize^2) + 
                          BodySize*AvgAge + SSCM + I(SSCM^2) + SSCM*AvgAge + AvgAge + 
                          (1|Sire) + (1|Year), 
                     data=datum, family=poisson, na.action="na.fail")

Quadratic effects were included for predictors due to expected non-linearities. Two random terms were included to account for the individual potential sires being sampled multiple times during the study and year (input as a factor) effects. AvgAge is a term related to population demographics, and interactions with predictors are included.

So here are my questions:

1. A reviewer suggested that there is temporal correlation in my data (e.g. male age at year one is correlated with male age at year two) that needs to be addressed. He suggested including a temporal autocorrelation structure, or including Year as a numeric predictor to deal with this. Am I missing something here or isn't this correlation the whole purpose of including the random effect for each male, coded as (1|Sire) in this case? Also, including year as a numeric predictor really mucks things up because:

   1. I'm not interested in the specific effect of Year,
   2. there was a good bit of variability in the number of males we measured each year, and
   3. it makes an already complex model more complex.

2. Is the following code sufficient to use to screen for overdispersion? So I assume here that if the ratio is < 1 then you likely do not have issues with overdispersion?

   overdisp.glmer(test)
   # Residual deviance: 84.153 on 105 degrees of freedom (ratio: 0.801)
   

Answer 1

It is possible that your data has temporal autocorrelation that the model you describe would not account for. Specifically, temporal autocorrelation in the model residuals with respect to sire age. For a given sire, might look something like this...

The Sire intercept would not account for this type of autocorrelation. I don't see that including Year as a numeric predictor would address this either, but nesting the Year intercept within the Sire intercept would, e.g. (1|Sire/Year). However, this approach would effectively give every observation a unique random intercept, so might be overkill. An alternate option would be a generalized additive mixed model (e.g. gamm function in {mgcv}, which can incorporate any of the {nlme} correlation structures).

Regardless, you might first consider visually inspecting for within-sire temporal autocorrelation in the residuals (e.g. see the acf function). Since you only have 6 years for each sire, identifying any clear pattern might be challenging. (As an aside, the ecologist in me doesn't see any a priori reason for strong temporal autocorrelation in male deer reproductive success, above and beyond what your model would capture with the Year intercept).

As for overdispersion, I believe that method is ok, though the webpage referenced in the overdisp.glmer function's documentation notes "Do PLEASE note the usual, and extra, caveats noted here: this is an APPROXIMATE estimate of an overdispersion parameter".

Another way to assess over-dispersion would be to compare the model you describe to a quasi-Poisson model which you can fit with various libraries (can even fit with glmer by adding a random intercept for every observation; see here for more information). You could compare the Poisson and quasi-Poisson models with AICc or examine their prediction intervals overlaid with the raw data.