2
$\begingroup$

I have a couple related questions about using a generalize linear mixed effects model to analyze data from an agricultural field experiment. I have found several posts that are similar to this question, but nothing that quite gets at my points of confusion.

The experiment involved 4 treatments replicated 4 times each: one treatment per field plot, for a total of 16 field plots.

On multiple dates (12) throughout the season I measured a response (counts) and a few covariates at all 16 plots. Within each plot on each date I would make multiple measurements of my response and covariates. So, measurements within a plot on a given date are non-independent (spatial pseudoreplication), and measurements through time for each plot are non-independent (temporal pseudoreplication). OK, so a linear mixed model seems appropriate.

My response data roughly fit a negative binomial distribution. OK, so a generalized linear mixed model seems appropriate.

I want to account for the variation that occurs within a plot and also the variation that occurs on different dates.

Would it be valid for me to model it this way (treating date as a factor)?

#fake data; each plot is sampled at 4 locations on 12 different dates
response=rnbinom(768,size=0.5,mu=2)
treatment=factor(rep(c("one","two","three","four"),times=192))
covariate=rnbinom(768,size=.01,mu=4)
date=factor(rep(1:12,each=64))
plotID=factor(rep(c("P1","P2","P3","P4","P5","P6","P7","P8","P9","P10","P11",
            "P12","P13","P14","P15","P16"),each=4,times=12))

#a model with random effects on the intercept
library(glmmADMB)
model1=glmmadmb(response~covariate+treatment+(1|date\plotID),family="nbinom2")

My concern is that, while plots are independent and randomly distributed, the dates are not. That is, I expect that on a given date the variance structure within a plot will be similar, but not between plots; however, I expect the variance structure to be similar across plots on a given date. Does this create problems for how I modeled the data above?

Am I even correct in thinking that the model would account for the variance that occurs on each date (1|date) and the variation that occurs within each plot on a given date (1|date:plotID)?

I could model date as a continuous fixed effect, but for reasons I won't go into I'd prefer not to do that. Regardless, I have no interest in the effects of time on my response, other than as a source of variation that I'd like to account for.

I hope my thinking about this whole thing doesn't reveal any TOO glaring gaps in understanding :)

Thanks in advance for the help, oh wise internet community!!

$\endgroup$
2
$\begingroup$

A nested random-effects formula of the form (1|date/plotID) (forward slash /, not backslash \) denotes "plot nested within date", i.e. responses vary across dates and across plots within dates. This is exactly equivalent to (and in fact is internally expanded to) (1|date)+(1|date:plotID), i.e. random effects of date and of a date-by-plot interaction. Given your crossed experimental design (each plot is measured on multiple dates, multiple plots are measured on each date, and there is more than one observation per date-plot combination), you can add a (1|plotID) term to the model.

(1|date)+(1|plotID)+(1|date:plotID), (1|date/plotID)+(1|plotID), and (1|plotID/date)+(1|date) are all equivalent; you should additionally be able to specify the model as (1|date*plotID), but at present that's not possible.

You say

My concern is that, while plots are independent and randomly distributed, the dates are not. That is, I expect that on a given date the variance structure within a plot will be similar, but not between plots; however, I expect the variance structure to be similar across plots on a given date. Does this create problems for how I modeled the data above?

I don't quite understand this, specifically what you mean by "variance structure". Most mixed models will assume in this scenario that within-plot variance is the same across plots, that different plots have the same temporal variance, and that different times have the same among-plot variance ... I don't know what the "variance structure between plots" means ... you can relax these assumptions if you absolutely need to, but (1) it can get tricky to program and (2) it's going to get hard to estimate all those separate variances.

One final thing to consider is that your experimental design in principle allows you to fit random effects of the form (treatment|plotID), i.e. looking at variation in treatment effects across plots. Schielzeth and Forstmeier 2008 comment that leaving these terms out in some cases (i.e. fitting random-intercept rather than random-slope models) may have overestimated treatment effects in some cases.

$\endgroup$
  • $\begingroup$ Thank you so much for your extended answer - especially for pointing out that my design is indeed crossed and not nested, and for pointing me in the direction of Schielzeth and Forstmeier (2008) Behavioral Ecology. If a likelihood ratio test indicated that including a (treatment|plotID) effect was not a significant improvement, do you think I should still leave it in the model to avoid Type I error? That seems to be the take-home from Schielzeth and Forstmeier. $\endgroup$ – Joe Kaser May 19 '16 at 18:17
  • 1
    $\begingroup$ I generally follow the principle of "estimate random effects if they are both theoretically and practically estimable given the qualitative (i.e., which terms in the model vary at which levels) and quantitative (i.e., how much replication/how many levels of random effects) aspects of the observational/experimental design , don't worry if they're small / insignificant" $\endgroup$ – Ben Bolker May 19 '16 at 19:48

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.