Crossed fixed effects model specification including nesting and repeated measures using glmm in R

Question

Background:

I am interested in looking at the effects that Culture, Treatment and Time have on my response variable Y as well as the interactions between them whilst accounting for the nested variation between tanks and the variation caused by measuring the same individuals through time.

• Culture and Treatment are crossed based on a figure found here: http://www.nature.com/nmeth/journal/v11/n10/full/nmeth.3137.html

• As “Tanks” are exclusive to each Treatment, Tank is nested within Treatment.

• The same 50 individuals per tank are measured through time, so Subjects are repeated measures.

• As both Culture and Treatment have less than 5 levels each, I have treated them as fixed effects and I would like to extract fixed effects coefficients.

• The Gamma distribution with link="log" was used because my response variable Y are non-integers (i.e. 0.3, 0.4 etc.).

• I have checked the random effects observation numbers as mentioned by Dr. Ben Bolker here: (Have I correctly specified my model in lmer?), and the models above are able to capture the nested random nature (N=18).

• I have chosen the lme4 package because it is able to handle crossed designs.

• I think that either (1|Treatment:Tank) or (1|Tank) are equivalent since the nesting of Tank within Treatment can be “discovered” from the way the data is structured, but I would like to keep it as (1|Treatment:Tank) to remind myself. Inferred from Dr. Doug Bates here: https://stat.ethz.ch/pipermail/r-sig-mixed-models/2009q3/002790.html

Potential models:

gamma_3 <- glmer(Y~Treatment*Culture*Time+(1|Treatment:Tank)+(1|Subject),
                 data=raw_data, family=Gamma(link="log")) 

gamma_4 <- glmer(Y~Treatment*Culture*Time+(1|Treatment:Tank)+(1|Treatment)+
                 (1|Culture)+(1|Subject), data=raw_data, family=Gamma(link="log"))

Questions:

1. Have I correctly captured the crossed nature of Culture and Treatment in model gamma_3? I understand that if they were treated as Random Effects it would be correct to specify them as (1|Culture) + (1|Treatment) as suggested by Dr. Doug Bates and Dr. Ben Bolker, respectively: https://stat.ethz.ch/pipermail/r-sig-mixed-models/2009q3/002790.html and http://glmm.wikidot.com/faq

   However, I have found less resources on how to specify crossed fixed effects.

2. In model gamma_4 I have added the above structure mentioned in question 1 [(1|Culture) + (1|Treatment)] in addition to keeping Culture and Treatment as Fixed effects to capture their crossed nature. However, this seems to be advised against by Dr. Ben Bolker here: Have I correctly specified my model in lmer?

3. Is (1|Subject) sufficient to specify repeated measures in the glmm() context?

4. Finally, I am getting the following error message although the model does produced an output:

   a) fixed-effect model matrix is rank deficient so dropping 5 columns / coefficients

   This was surprising to me because I have heaps of data. I have 50 individuals per tank, 3 tanks per treatment etc.

   b) Warning message:In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model failed to converge with max|grad| = 0.101624 (tol = 0.001, component 1)

   Why would the model fail to converge but still produced an output? Any clues on how I can fix this or glean some more information on the nature of the error?

---

Note: Also, a tag for crossed designs does not seem to exist. Potentially someone with a high enough reputation could create a new tag if they find it useful or correct for this community?

Answer 1

1. Yes. * in (R's version of) Wilkinson-Rogers notation specifies "interaction plus all lower level terms", so (ignoring Time for the moment) Treatment*Culture is equivalent to 1 + Treatment + Culture + Treatment:Culture. If you wanted Treatment nested within Culture you would use Culture/Treatment = 1 + Culture + Culture:Treatment (i.e. omit main effect of Treatment. Note, however, that for fixed effects specifying nesting vs crossing changes the model parameterization, but the overall fit (e.g. number of parameters, predictions of the model, log-likelihood, etc.) is the same for nesting vs. crossing.

2. This is a bad idea, not just because it's almost impossible to fit a random effect when the grouping variable has only 2 levels, but because your model is now overfitted - you've included Culture as both a random effect and a fixed effect.

3. Depends what you mean by "repeated measures". It does specify that samples within an individual are not independent. However, you might want to consider variation in slopes among individuals (Time|Subject) (if Time is being considered as a continuous covariate), or consider temporal autocorrelation (not possible in lme4, but in nlme using e.g. correlation=corAR1()

4. (a) You're probably treating Time as a factor (categorical predictor) and have some times missing for some culture/treatment combinations? (That's technically a "message", not an "error message".) attr(getME(fitted_model,"X"),"col.dropped") should tell you the column numbers and names of the variables that were dropped. (b) have you looked at ?convergence ?

You could actually use nlme::lme for this problem if you wanted, by (1) using log-Normal instead of Gamma responses, (2) fitting Subject nested within Tank ...

lme(Y~Treatment*Culture+Time,random=~1|Tank/Subject),data=raw_data)