6

I am analyzing plant densities (n° of plants/0,36m²) across 5 different treatments (repeated on two blocks) in time (year and season). Since my outcome is plant density and repeated measures were done, I wish to analyze the data with a mixed model with a Poisson distribution.

However, on each block:plot:year:season, more than one quadrat was realized, i.e. I have more than one value of plant density for each block:plot:year:season in order to account for the within field variability. Moreover, I am interested in i) a time trend so I wish to include year as a fixed continuous predictor and 2) comparing plant densities between seasons. Hence, I was wondering if it was recommendable to include year as a fixed continuous predictor and a random categorical variable (in order to account for the fact that quadrats realized at the same sampling season on each plot are not independent).

My model (in R) with the glmer function would be the following:

mod = glmer(density ~ treatment*as.numeric(year)*season+(1|plot/as.factor(year)/season),
            family="poisson",data=dat)
  1. Could anybody specify whether this recommendable or not?
  2. Considering the quadrats are fixed within a block:plot:year:season, should I also add a quadrat ID? However, this would result in single observations per combination...
Guillaume A2
  • 311
  • 2
    (Q2) In a mixed model with Gaussian family (lmer) using an observation-level random effect does not make any sense because there is a noise term that plays its role. However, in case of Poisson GLMM like you are using here, observation-level random effect can be used to model overdispersion. See e.g. https://stats.stackexchange.com/a/83652/28666, and https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#fitting-models-with-overdispersion. In your case, (1|block:plot:year:season:quadrat) will be observation-level. – amoeba Feb 05 '18 at 14:11
  • 2
    (Q1) I will assume that year is coded as numeric (note that whenever it appears to the right of | it will automatically be treated as categorical anyway). So, treatment*year + (1|plot/year/season) is slightly unusual but does make sense. However, note that year is within-plot so you can also consider using random slopes instead: treatment*year + (year|plot) + (1|plot:year:season). – amoeba Feb 05 '18 at 14:16
  • 2
    (Q1+Q2) Assuming that quadrats have unique ids, treatment*year + (year|plot) + (1|quadrat) might be a reasonable model and a simpler way to write it. – amoeba Feb 05 '18 at 14:20
  • I really appreciate the help but this is still a little unclear. 1) For sure plot is nested in block, i.e. +(1|block/plot) 2) year (2001-2017) and season (winter, spring, summer) are common to all plots, that is to say every year and season every plot was sampled. 3) plots have unique identifiers whereas quadrats are named from A1 to H4 (32 total on all plots). ps: i edited in order to add season as fixed also. – Guillaume A2 Feb 05 '18 at 15:04
  • 2
    OK, then treatment*year*season + (year*season|plot) + (1|plot/quadrat) would make sense as as a starting point. But you will have to deal with season being categorical which might turn this model into a bit of a mess. You can't really use block as a random effect because it only has 2 levels which is way too few. – amoeba Feb 05 '18 at 15:09

1 Answers1

3

I think it's fine to include year as both fixed and random. This has been discussed a bunch, here are a few links that might help you think about justifying including the same variable as both random and fixed in the same model:

https://stats.stackexchange.com/q/155056

https://dynamicecology.wordpress.com/2015/11/04/is-it-a-fixed-or-random-effect/#comment-45407

Mixed Effects Model with Nesting

Michael Roswell
  • 351
  • 1
    Regarding / vs :, it's quite simple: (1|plot/year/season) expands to (1|plot)+(1|plot:year)+(1|plot:year:season). – amoeba Feb 05 '18 at 14:39
  • amoeba's comments above are a better answer than mine... but due to low reputation I don't know how to address that. – Michael Roswell Feb 05 '18 at 14:54