How to include nested data in a mixed-effects model using R?

Question

I am analyzing data from bird foraging surveys using the lme4 package in R and I am interested in the effects of field size (area), among other variables, on swallow rate of use. The surveys took place over 2 years, and I’d like to include year as a fixed effect in the model. Each field (n = 31) was surveyed 6 times. So the data look like:

Field   Survey  Year    Area    RateOfUse
1       1       2017    3.06    0
1       2       2017    3.06    0.327
1       3       2017    3.06    0.327
1       4       2017    3.06    0.327
1       5       2017    3.06    3.92
1       6       2017    3.06    0.327
1       1       2018    3.06    0
1       2       2018    3.06    0.392
1       3       2018    3.06    2.55
1       4       2018    3.06    2.94
1       5       2018    3.06    0.588
1       6       2018    3.06    0
2
...

In order to include field as a random effect (accounting for the fact that each field has data from 6 surveys), does this make sense:

lme(RateOfUse ~ Area + Year, random = ~1|Field/Survey)

Alternatively, I could included Survey as a fixed effect to account for the fact that rate of use might have dropped off during later surveys:

lme(RateOfUse ~ Area + Year + Survey, random=~1|Field)

I also have a separate variable Field_Survey which is

Field_Survey    Year
1_1             2017
1_2             2017
1_3             2017
1_4             2017
1_5             2017
1_6             2017
1_1             2018
1_2             2018
1_3             2018
1_4             2018
1_5             2018
1_6             2018
...

So instead I could do something like:

lme(RateOfUse ~ Area + Year, random = ~1|Field_Survey)

Which of these models makes more sense? Is 31 groupings (fields) way too many? I am a bit rusty with mixed effects models and want to make sure I'm thinking about this correctly.

Answer 1

Note that random effects terms are translated to correlations. In particular, by including a random effect (e.g., a random intercept) in the Field level, you assume that RateOfUse measurements within the same Field are correlated. And by including a nested random intercept for Survey within Field, you assume that RateOfUse measurements from the same Survey within the same Field are more strongly correlated than RateOfUse measurements from different Surveys within the same Field.

Now, the question is, is this the case for your data? You can test data using a likelihood ratio test, i.e., by comparing the two models:

fm1 <- lme(RateOfUse ~ Area + Year, random = ~ 1 | Field)
fm2 <- lme(RateOfUse ~ Area + Year, random = ~ 1 | Field / Survey)
anova(fm1, fm2)

When you include the Survey as a fixed-effect, this translates to a different average for the RateOfUse per Survey level controlled for Area and Year. (Note: this will depend on whether Survey is a factor or not; if not, then Survey is included as a numeric variable, and you assume a linear relationship between the average RateOfUse and the Survey levels).