233
$\begingroup$

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?
$\endgroup$
8
  • 15
    $\begingroup$ The manual and three vignettes for the lme4 package can be found on CRAN $\endgroup$
    – Henry
    Commented Jul 17, 2011 at 22:54
  • 5
    $\begingroup$ 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 $\endgroup$ Commented Jul 18, 2011 at 7:41
  • 12
    $\begingroup$ 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. $\endgroup$
    – whuber
    Commented Apr 23, 2018 at 15:28
  • 3
    $\begingroup$ Another useful source which uses Singer & Willet's nomenclature (e.g. gammas and betas) is Magnusson's article on lme and lmer specification. Magnusson shows both the "classical" formulas alongside lme and lmer code. $\endgroup$
    – Ben
    Commented Mar 4, 2019 at 23:22
  • 2
    $\begingroup$ Here is a comprehensive and clear summary of different formula from the maintainer of lme4 $\endgroup$
    – blueskyddd
    Commented Oct 28, 2021 at 15:54

3 Answers 3

251
$\begingroup$

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.

$\endgroup$
10
  • 2
    $\begingroup$ @Mike Lawrence Thanks for the answer! how is a 3 level model estimated then? where one grouping factor is nested within another? $\endgroup$
    – DBR
    Commented Jul 18, 2011 at 23:52
  • 3
    $\begingroup$ 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. $\endgroup$ Commented Jul 20, 2011 at 20:07
  • 3
    $\begingroup$ @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. Please note that I exchange V2 for V3 in my notation. $\endgroup$ Commented Jan 25, 2016 at 16:42
  • 1
    $\begingroup$ @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 $\endgroup$ Commented Jan 27, 2016 at 19:26
  • 10
    $\begingroup$ 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! $\endgroup$ Commented Oct 2, 2021 at 2:44
110
$\begingroup$

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.

$\endgroup$
6
  • 7
    $\begingroup$ Nice. It could be better wtith information about nested '/' factors and double-bar notation '||' $\endgroup$
    – skan
    Commented Sep 23, 2015 at 11:36
  • 3
    $\begingroup$ What about the : symbol? $\endgroup$
    – eastafri
    Commented Feb 1, 2017 at 9:38
  • 2
    $\begingroup$ @eastafri It means that the same thing it does everywhere in R (formulas) - the interaction between two variables. $\endgroup$
    – Livius
    Commented Feb 1, 2017 at 9:56
  • $\begingroup$ 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? $\endgroup$
    – Muno
    Commented Jun 24, 2018 at 19:11
  • 2
    $\begingroup$ @shadowtalker (X|Subject) is the same as (1+X|Subject) -- like most formulae bits in R, the intercept is there unless you explicitly suppress it. $\endgroup$
    – Livius
    Commented Nov 20, 2023 at 23:04
9
$\begingroup$

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?

$\endgroup$