R: lmer - correct syntax?

Question

I'm wondering if a treatment has an effect on my mite population. Therefore I've got a dataset with repeated measurements, some data is missing.

data:

ID   Treatment    Mites   Time    Location    StartPopulation    X1bib
ID1  Control      7       t1      Loc1        5                  10000
ID1  Control      8       t2      Loc1        5                  10000
ID1  Control      10      t3      Loc1        5                  10000
ID1  Control      11      t7      Loc1        5                  10000
ID2  Control      12      t1      Loc2        11                 13000
ID2  Control              t2      Loc2        11                 13000
ID2  Control      14      t3      Loc2        11                 13000
ID3  Treatment    20      t1      Loc1        20                 12000
ID3  Treatment    22      t2      Loc1        20                 12000
ID3  Treatment            t3      Loc1        20                 12000
ID4  Treatment    20      t1      Loc11       18                 11500
and so on..
totally: 110 IDs; 7 different measurements (Time)

variables:

ID:              factor, unique ID for each population
Treatment:       factor ("Treatment" or "Control")
Mites:           numeric, the variable I'm interested in
Time:            factor with total 7 levels
Location:        factor with total 11 levels
StartPopulation: numeric (mean of Mites for t=-3, -2, -1 -> before Treatment started)
X1bib:           numeric

I'm interested if my Treatment changed the Mites count - and if yes if there's an increase in it's effect over time. StartPopulation sure had an influence on Mites, otherFactor and Location could've had also.

As I use a mixed model I'd like to use glmer in R. My syntax looks like this: (changed it, thank you for your answers so far)

PPP <- glmer(Mites ~ Treatment * Time + StartPopulation + X1bib + (1|ID) + (1|Location), data=vat_database, family=poisson)

which outputs:

Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: poisson  ( log )
Formula: Mites ~ Treatment * Time + StartPopulation + X1bib + (1 | ID) +      (1 | Location)
   Data: vat_database
     AIC      BIC   logLik deviance df.resid 
     Inf      Inf     -Inf      Inf      349 
Random effects:
 Groups   Name        Std.Dev.
 ID       (Intercept) 1       
 Location (Intercept) 1       
Number of obs: 367, groups:  ID, 78; Location, 9
Fixed Effects:
                  (Intercept)             TreatmentTreatment                     Timevmf_A2                     Timevmf_A3                     Timevmf_K1                     Timevmf_K2  
                    2.418e-01                      5.342e-01                      3.252e-01                      5.389e-01                      5.725e-01                      1.102e+00  
                   Timevmf_K3                     Timevmf_K4                StartPopulation                          X1bib  TreatmentTreatment:Timevmf_A2  TreatmentTreatment:Timevmf_A3  
                    1.079e+00                      7.893e-01                      1.486e-01                     -1.331e-06                     -4.664e-01                     -5.453e-01  
TreatmentTreatment:Timevmf_K1  TreatmentTreatment:Timevmf_K2  TreatmentTreatment:Timevmf_K3  TreatmentTreatment:Timevmf_K4  
                   -4.513e-01                     -5.476e-01                     -4.477e-01                     -6.858e-01  
fit warnings:
Some predictor variables are on very different scales: consider rescaling
convergence code 0; 1 optimizer warnings; 81500 lme4 warnings

Am I right considering that on Time="vmf_K1" my Treatment Mite Population was -4.513e-01 smaller than my Control Mite Population? How about significances?

Answer 1

As gung pointed out, you need to use a Poisson distribution because you have count data as the dependent variable, and this means you need to use glmer() instead of (lmer). Luckily, this doesn't require much change in the syntax.

You're using a mixed model because you have repeated measures from the same 'subjects' (whatever that is in this study?). So basically, by adding (1|ID), you're allowing each 'subject' to have its own intercept, or rather, a random value added to the intercept. These random effects are assumed to follow a normal distribution and in many cases this is an acceptable assumption.

By adding (1|Time), you're allowing each time point (1-3) to have its own random effect, or intercept. This makes little sense because firstly, there are too few levels (a general recommendation often repeated is at least 5-6 levels, in this case you have 3), and secondly, if there is an effect of time on Mites, it probably follow some kind of pattern and is not random. Time should not be used as a random effect, but you could allow each individual to have its own effect of time on mites, and in this case you would use (Time|ID). With this syntax, each 'subject' will have their own random slope of the effect of Time on Mites, but I don't think this is what you're asking for.

So perhaps it's best to stick with (1|ID) as the random effects component. Let's now discuss the fixed effects.

You're interested in treatment so obviously it should be included as a fixed effect. You're also interested in the effect of Time, and the effect of treatment over time. This is called an interaction effect and it can be added to the model by typing Treatment * Time in the formula.

When adding the interaction term you will get a new line in the output:

Treatment
Time
Treatment:Time

The first two simply denotes the effect of Treatment and Time on Mites, separately. The last one is the interaction, or the specific effect of Time in the treatment group (assuming that control is the reference). If the interaction term Treatment:Time has low standard errors (or gives a low p-value if you're calculating it) it means that the effect of time is different for the treatment group than for the control group.

So I think your model might look like this:

PPP <- glmer(Mites ~ Treatment * Time + StartPopulation + Location + otherFactor + (1|ID), data=vat_database, family=poisson)

I hope this helps!

Answer 2

Besides suggestions made by others, you may wish to use "Location" as a random effect (eg 1|Location). Its not clear how many time points you have; if you have only 3 then you should not use time as a random effect. If you have 7 then it is ok to use it as a random effect. If the amount of variation over time varies between locations, you might want include a random effect that allows for that (1|Location:Time). With Poisson regression you should check for overdispersion and may want to include an observation-level random effect to account for that. Create a numeric index column like this

df$indix <- 1:dim(df)

Fabulous ecology-based information is on Poisson regression is available at Ben Bolkers RPubs site. He provides code in one of his entries for a function to check for overdispersion. There are also R packages that have functions.