I’m having a hard time figuring out what type of model I should use to analyse data with a rather complexe design.
Experimental goal and protocol: The goal is to study the effect of beehive color and inner temperature on bee activity. The experimental layout is: The number of bees at the entrance to 48 beehives was sampled over several months. Beehives were of different colors (4 colors) and the outside temperature of the beehives was controlled (3 temperatures). Because it was not possible to set this up in a random manner, the beehives were 1/ first grouped into three lots (each lot had 16 beehives), with three different outside temperatures (TempOut: 10, 15, 20) ; and 2/ in each lot, the 16 beehives were grouped into 4 Groups (4 beehives in each group), painted in different colors (Color: factor noted 1 to 4). I would call this a nested plot design. In addition, there was a wasp infestation close by, during month 6 of the experiment, causing the number of bees to drop. We consider here that this variable has three modalities (Wasp: "Before" infestation: 5 months of data, "During" infestation: one month of data and "After" infestation: 6 months of data). Observations (number of bees at entrance, inner temperature of beehive noted "TempIn") was done every week during 12 months.
Here’s an illustation, which shows the experimental design.
We would like to analyse the effect of (1) beehive color "Color" (categorical, 4 levels), (2) beehive inner temperature "TempIn" (continuous), and (3) wasp infestation "Wasp" (categorical, 3 levels) on bee activity "BeeActivity" (continuous).
Questions : (1) What type of statistical model should I test?
There is clearly a temporal autocorrelation in the data, as beehives were sampled multiple times. Also, the set-up is clearly not randomized.
I’m thinking of three solutions:
i. applying a mixed model effect, with a temporal autocorrelation, with three fixed effects (ColorTempInWasp), one random effect (Group) and specifying that data from each beehive is temporally correlated.
The problem I see is that the Groups (12 of them, with 4 beehives in each group) received different treatments (Color*TempIn). This is different from random effects such as « subject » (different patients in a test for example) which all receive the same treatments.
ii. applying a mixed model effect, with a temporal autocorrelation, with three fixed effects (ColorTempInWasp), one random effect (TempOut) and specifying that data from each beehive is temporally correlated. The problem is that TempIn and TempOut will be well correlated.
iii. The other thing I was thinking of was to average bee activity for the 4 beehives in each group, so that I would have 12 average measures of bee activity per sampling date. I could then test a linear model of the type: gls(BeeActivity ∼ ColorTempInWasp, na.action = na.omit, data = Beehive, correlation = corAR1(form =∼ DoY)) The assumption here, I think, would be that these average measures of bee activity randomly received the different treatments. However, I don’t think this works because the hives were still grouped into « lots » with constant temperatures within the « lots ». These kind of look like blocs in a randomized block design, but the blocs are not repeated.
(2) How would the model(s) be correctly written in R?
For the mixed model, I was thinking of something like:
M5 <- lme(BeeActivity ∼ Color*TempIn*Wasp, random= ∼1 | Group, method = "REML", correlation = corAR1(form =∼ DoY), data = Beehive)
Can I set both a random structure and an autocorrelation structure in the same model?
Help much appreciated! Thanks.
Example dataset attached here (partially re-constructed)