197

There's a lot of discussion going on on this forum about the proper way to specify various hierarchical models using lmer.

I thought it would be great to have all the information in one place. A couple of questions to start:

  1. How to specify multiple levels, where one group is nested within the other: is it (1|group1:group2) or (1+group1|group2)?
  2. What's the difference between (~1 + ....) and (1 | ...) and (0 | ...) etc.?
  3. How to specify group-level interactions?
amoeba
  • 93,463
  • 28
  • 275
  • 317
  • 13
    The [manual](http://cran.r-project.org/web/packages/lme4/lme4.pdf) and three vignettes for the `lme4` package can be [found on CRAN](http://cran.r-project.org/web/packages/lme4/) – Henry Jul 17 '11 at 22:54
  • 5
    There are, in addition to the CRAN materials, lecture slides plus draft chapters of a book Doug is writing on (G)LMMs and R with **lme4** available from [r-forge](http://lme4.r-forge.r-project.org/) – Gavin Simpson Jul 18 '11 at 07:41
  • Direct link to the arXiv version of the JSS paper by Bates et al.: [Fitting Linear Mixed-Effects Models using lme4](https://arxiv.org/pdf/1406.5823.pdf) (in particular Section 2.2 "Understanding mixed-model formulas"). See also [the relevant section](http://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#model-specification) of Ben Bolker's FAQ. – amoeba Jan 30 '18 at 15:07
  • 8
    Arguably, the language used by `lmer` is of general statistical interest and thus is not solely a matter of programming. I am therefore voting to keep this thread open. – whuber Apr 23 '18 at 15:28
  • 2
    Another useful source which **uses Singer & Willet's nomenclature** (e.g. gammas and betas) is Magnusson's [article on lme and lmer specification](https://rpsychologist.com/r-guide-longitudinal-lme-lmer). Magnusson shows both the "classical" formulas alongside lme and lmer code. – Ben Mar 04 '19 at 23:22
  • [Here](https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#model-definition) is a comprehensive and clear summary of different formula from the maintainer of `lme4` – blueskyddd Oct 28 '21 at 15:54

3 Answers3

218

What's the difference between (~1 +....) and (1 | ...) and (0 | ...) etc.?

Say you have variable V1 predicted by categorical variable V2, which is treated as a random effect, and continuous variable V3, which is treated as a linear fixed effect. Using lmer syntax, simplest model (M1) is:

V1 ~ (1|V2) + V3

This model will estimate:

P1: A global intercept

P2: Random effect intercepts for V2 (i.e. for each level of V2, that level's intercept's deviation from the global intercept)

P3: A single global estimate for the effect (slope) of V3

The next most complex model (M2) is:

V1 ~ (1|V2) + V3 + (0+V3|V2)

This model estimates all the parameters from M1, but will additionally estimate:

P4: The effect of V3 within each level of V2 (more specifically, the degree to which the V3 effect within a given level deviates from the global effect of V3), while enforcing a zero correlation between the intercept deviations and V3 effect deviations across levels of V2.

This latter restriction is relaxed in a final most complex model (M3):

V1 ~ (1+V3|V2) + V3

In which all parameters from M2 are estimated while allowing correlation between the intercept deviations and V3 effect deviations within levels of V2. Thus, in M3, an additional parameter is estimated:

P5: The correlation between intercept deviations and V3 deviations across levels of V2

Usually model pairs like M2 and M3 are computed then compared to evaluate the evidence for correlations between fixed effects (including the global intercept).

Now consider adding another fixed effect predictor, V4. The model:

V1 ~ (1+V3*V4|V2) + V3*V4

would estimate:

P1: A global intercept

P2: A single global estimate for the effect of V3

P3: A single global estimate for the effect of V4

P4: A single global estimate for the interaction between V3 and V4

P5: Deviations of the intercept from P1 in each level of V2

P6: Deviations of the V3 effect from P2 in each level of V2

P7: Deviations of the V4 effect from P3 in each level of V2

P8: Deviations of the V3-by-V4 interaction from P4 in each level of V2

P9 Correlation between P5 and P6 across levels of V2

P10 Correlation between P5 and P7 across levels of V2

P11 Correlation between P5 and P8 across levels of V2

P12 Correlation between P6 and P7 across levels of V2

P13 Correlation between P6 and P8 across levels of V2

P14 Correlation between P7 and P8 across levels of V2

Phew, That's a lot of parameters! And I didn't even bother to list the variance parameters estimated by the model. What's more, if you have a categorical variable with more than 2 levels that you want to model as a fixed effect, instead of a single effect for that variable you will always be estimating k-1 effects (where k is the number of levels), thereby exploding the number of parameters to be estimated by the model even further.

Mike Lawrence
  • 12,691
  • 8
  • 40
  • 65
  • 1
    @Mike Lawrence Thanks for the answer! how is a 3 level model estimated then? where one grouping factor is nested within another? –  Jul 18 '11 at 23:52
  • 3
    I think DBR is referring to levels in the hierarchy. What I described is a 2-level hierarchical model, with observations nested within subjects, and DBR is asking about 3-level hierarchies, an example of which might be test items within students within schools where you want to model both students and schools as random effects, with students nested within schools. In such cases I presume that the school level deviations are first computed then the student-from-school deviations. – Mike Lawrence Jul 20 '11 at 20:07
  • Say I have one independent variable (X) at individual level and one independent variable (Z) at group level . Both are continuous variable . If the model is $$Y_{ij}=\gamma_{00}+\gamma_{10}X_{ij}+\gamma_{01}Z_{j}+\gamma_{11}X_{ij}Z_{j}+u_{1j}X_{ij}+u_{0j}+e_{ij}$$, where subscript $i$ denotes $i$th individual and $j$ denotes $j$th group . then using `lmer` syntex will the model be `Y~X+Z+(1|group)+(0+X|Z)` ,where `group` is another variable in the data frame to indicate in which group the individual belongs to ? – ABC Jul 20 '15 at 04:47
  • 1
    @MikeLawrence I've used your answer many times as everybody else, but I do have a nagging question. In the M2 model it sounds as though you should get different intercepts for every V2 level. Yet this doesn't seem to be the case when I try to reproduce it [here](http://rinterested.github.io/statistics/random_effects_simulation.html). Please note that I exchange V2 for V3 in my notation. – Antoni Parellada Jan 25 '16 at 16:42
  • @AntoniParellada Hm, I haven't looked at lmer in a while and it looks like they've changed how they access the estimates. Does ranef(m2) give you anything different from coefficients(m2)? – Mike Lawrence Jan 26 '16 at 17:26
  • Yes, it does. You can find the code [here](https://github.com/RInterested/SIMULATIONS_and_PROOFS/blob/master/athletes%20mixed%20effects) ready to copy and paste. This is confusing... The intercepts running `ranef(m2)` seem crazy...$e^{-13}$? – Antoni Parellada Jan 26 '16 at 17:34
  • @AntoniParellada Ah, I think what's happened here is that in evaluating m2 for your data, lmer has decided that the intercept sd is basically zero, so all the individual coefficients for the intercept are nearly zero (the `e-13` in `ranef(m2)` is likely attributable to floating point computation error). Note that `coef(m2)$Intercept == fixef(m2)$Intercept + ranef(m2)$Intercept` – Mike Lawrence Jan 27 '16 at 19:26
  • @MikeLawrence Great! Thank you. It makes sense, then, to expect different intercepts. I will work on tweaking my dataset. – Antoni Parellada Jan 27 '16 at 19:33
  • re: P4... 'while enforcing a zero correlation between the intercept deviations and V3 effect deviations across levels of V2'. I think this is definitionally true, but I want to cast light doubt on what lme4 does for this foruma. I remember elsewhere on this site someone speaking of simulation results that lme4 (at least for some versions) doesn't enforce 0 correlation - it just fails to estimate the correlation. I'm sorry that I can't find the entry off hand. – russellpierce Jun 28 '18 at 14:47
  • 3
    This isn't just a great answer. This represents the best of what can be found on SE: no snark, no pissing contests, no smugness, but a gentleman who took some serious time out of his day to help us Data Science students - now going on 10 years! Thank you @Mike Lawrence! – Douglas Wiley Oct 02 '21 at 02:44
86

The general trick is, as mentioned in another answer, is that the formula follows the form dependent ~ independent | grouping. The groupingis generally a random factor, you can include fixed factors without any grouping and you can have additional random factors without any fixed factor (an intercept-only model). A + between factors indicates no interaction, a * indicates interaction.

For random factors, you have three basic variants:

  1. Intercepts only by random factor: (1 | random.factor)
  2. Slopes only by random factor: (0 + fixed.factor | random.factor)
  3. Intercepts and slopes by random factor: (1 + fixed.factor | random.factor)

Note that variant 3 has the slope and the intercept calculated in the same grouping, i.e. at the same time. If we want the slope and the intercept calculated independently, i.e. without any assumed correlation between the two, we need a fourth variant:

  • Intercept and slope, separately, by random factor: (1 | random.factor) + (0 + fixed.factor | random.factor). An alternative way to write this is using the double-bar notation fixed.factor + (fixed.factor || random.factor).

There's also a nice summary in another response to this question that you should look at.

If you're up to digging into the math a bit, Barr et al. (2013) summarize the lmer syntax quite nicely in their Table 1, adapted here to meet the constraints of tableless markdown. That paper dealt with psycholinguistic data, so the two random effects are Subjectand Item.

Models and equivalent lme4 formula syntax:

    • $Y_{si} = β_0 + β_{1}X_{i} + e_{si}$
    • N/A (Not a mixed-effects model)
    • $Y_{si} = β_0 + S_{0s} + β_{1}X_{i} + e_{si} $
    • Y ∼ X+(1∣Subject)
    • $Y_{si} = β_0 + S_{0s} + (β_{1} + S_{1s})X_i + e_{si}$
    • Y ∼ X+(1 + X∣Subject)
    • $Y_{si} = β_0 + S_{0s} + I_{0i} + (β_{1} + S_{1s})X_i + e_{si}$
    • Y ∼ X+(1 + X∣Subject)+(1∣Item)
    • $Y_{si} = β_0 + S_{0s} + I_{0i} + β_{1}X_{i} + e_{si}$
    • Y ∼ X+(1∣Subject)+(1∣Item)
    • As (4), but $S_{0s}$, $S_{1s}$ independent
    • Y ∼ X+(1∣Subject)+(0 + X∣ Subject)+(1∣Item)
    • $Y_{si} = β_0 + I_{0i} + (β_{1} + S_{1s})X_i + e_{si}$
    • Y ∼ X+(0 + X∣Subject)+(1∣Item)

References:

Barr, Dale J, R. Levy, C. Scheepers und H. J. Tily (2013). Random effects structure for confirmatory hypothesis testing: Keep it maximal. Journal of Memory and Language, 68:255– 278.

Jiāgěng
  • 60
  • 5
Livius
  • 2,066
  • 14
  • 17
  • 7
    Nice. It could be better wtith information about nested '/' factors and double-bar notation '||' – skan Sep 23 '15 at 11:36
  • 2
    What about the : symbol? – eastafri Feb 01 '17 at 09:38
  • 2
    @eastafri It means that the same thing it does everywhere in R (formulas) - the interaction between two variables. – Livius Feb 01 '17 at 09:56
  • In (6), my understanding is that $S_{0s}$ and $S_{1s}$ have no correlation between them. In other words, as random variables, their covariance is $0$. To say that $S_{0s}$ and $S_{1s}$ are independent is a stronger statement and, hence, not necessarily true. Am I mistaken? – Muno Jun 24 '18 at 19:11
6

The | symbol indicates a grouping factor in mixed methods.

As per Pinheiro & Bates:

...The formula also designates a response and, when available, a primary covariate. It is given as

response ~ primary | grouping

where response is an expression for the response, primary is an expression for the primary covariate, and grouping is an expression for the grouping factor.

Depending on which method you use to perform mixed methods analysis in R, you may need to create a groupedData object to be able to use the grouping in the analysis (see the nlme package for details, lme4 doesn't seem to need this). I can't speak to the way you have specified your lmer model statements because I don't know your data. However, having multiple (1|foo) in the model line is unusual from what I have seen. What are you trying to model?

chl
  • 50,972
  • 18
  • 205
  • 364
Michelle
  • 3,640
  • 1
  • 23
  • 33