6
$\begingroup$

I'm considering a mixed-effects model to try to understand factors that influence the number of ticks sampled on wild rodents. My data is nested so that I have one tick count per rodent, multiple rodents per site and multiple sites per year (sites are repeated across years but not every site is present in all years).

So far, without fixed effects, my model looks something like this:

glmer(Ticks~(1|Year)+(1|Site), family=poisson, data=tickdata)

I understand that this will account for random variation between sites and between years. My main concern is that I think the effect of specific sites will change between years, e.g. some sites may recieve more rainfall than others, but the identity of the sites with the highest rainfall differs between years.

So do I need some sort of interaction term? Maybe something like this:

glmer(Ticks~(1|Year)+(1|Site)+(1|Year:Site), family=poisson, data=tickdata) 
$\endgroup$
4
  • $\begingroup$ I wonder if it makes more sense to nest year at the bottom of the hierarchy: observation in year nested in ticks, nested in rodents, nested in sites? $\endgroup$
    – Alexis
    May 6, 2014 at 14:28
  • $\begingroup$ @Alexis, you don't need to nest in this case (and probably shouldn't) $\endgroup$
    – Ben Bolker
    May 6, 2014 at 21:35
  • $\begingroup$ Because...? Are there no year-specific determinants of tick prevalence? Are there no site-specific determinants of tick prevalence? Rodent-specific? Ignoring such structural heterogeneity of variation surely guarantees biased estimates of the standard errors? $\endgroup$
    – Alexis
    May 6, 2014 at 21:49
  • 1
    $\begingroup$ I mean that you can cross random effects rather than nesting: exactly as the OP said, there are site-specific effects, year-specific effects, and site-by-year effects. (You could also consider rodent-specific effects, which would be observation-level effects/characterizing overdispersion if each rodent is measured once. Otherwise (if rodents are recaptured/recounted) you could also include observation-level effects). $\endgroup$
    – Ben Bolker
    May 7, 2014 at 1:43

1 Answer 1

9
$\begingroup$

Have you tried it? That sounds like it should be fine.

set.seed(101)
## generate fully crossed design:
d <- expand.grid(Year=2000:2010,Site=1:30)
## sample 70% of the site/year comb to induce lack of balance
d <- d[sample(1:nrow(d),size=round(0.7*nrow(d))),]
## now get Poisson-distributed number of obs per site/year
library(plyr)
d <- ddply(d,c("Site","Year"),transform,rep=seq(rpois(1,lambda=10)))
library(lme4)
d$ticks <- simulate(~1+(1|Year)+(1|Site)+(1|Year:Site),
                    family=poisson,newdata=d,
                    newparams=list(beta=2, ## mean(log(ticks))=2
                               theta=c(1,1,1)))[[1]]
mm <- glmer(ticks~1+(1|Year)+(1|Site)+(1|Year:Site),
                    family=poisson,data=d)

We get out approximately what we put in -- equal variances at each level:

## Generalized linear mixed model fit by maximum likelihood (Laplace
##   Approximation) [glmerMod]
##  Family: poisson  ( log )
## Formula: ticks ~ 1 + (1 | Year) + (1 | Site) + (1 | Year:Site)
##    Data: d
## 
##      AIC      BIC   logLik deviance df.resid 
##  12487.3  12510.2  -6239.7  12479.3     2267 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9944 -0.6842 -0.0726  0.6010  3.8532 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Year:Site (Intercept) 1.0818   1.0401  
##  Site      (Intercept) 1.0490   1.0242  
##  Year      (Intercept) 0.9787   0.9893  
## Number of obs: 2271, groups:  Year:Site, 231 Site, 30 Year, 11
## 
## Fixed effects:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   2.1952     0.3593   6.109    1e-09 ***
## ---
## Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

You may want to include an observation-level random effect to allow for overdispersion (see the "grouse ticks" example in http://rpubs.com/bbolker/glmmchapter)

$\endgroup$
5
  • $\begingroup$ Trying your code here, I get an error: Error in UseMethod("simulate") : no applicable method for 'simulate' applied to an object of class "formula". Any thoughts on what that may be about? $\endgroup$ Jan 13, 2015 at 17:57
  • $\begingroup$ are you sure (1) you loaded lme4 (2) you have a recent version of lme4 ? $\endgroup$
    – Ben Bolker
    Jan 13, 2015 at 18:11
  • $\begingroup$ Seems I had an older version of lme4. Sorry to bother you. $\endgroup$ Jan 13, 2015 at 18:19
  • $\begingroup$ Why did we need to sample 70% of the site/year comb to induce lack of balance here? $\endgroup$ Apr 8, 2017 at 13:36
  • $\begingroup$ I no longer remember why I did that. $\endgroup$
    – Ben Bolker
    Apr 8, 2017 at 14:15

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.