How do you know the number of random effects in a mixed effects model?

Question

I am trying to fit a random slope model in R and my code is as follows:

    lmer(data=ds, Outcome ~ treatment + (0 + treatment|ID))

I get the following error message when I try running this code:

    Error: number of observations (=2035) <= number of random effects (=2035) for term (0 + treatment | ID); the random-effects parameters and the residual variance (or scale parameter) are probably unidentifiable

I have 407 unique IDs and 5 treatments, so the number of observations is 407*5 = 2035. However, I don't understand why I have 2035 random effects. Being a random slope only model, I would expect that I am estimating one random slope per unique ID and so I thought there would only be 407 random effects in my model.

I also tried fitting a random slope and intercept model as follows:

    lmer(data=ds, Outcome ~ treatment + (1 + treatment|ID))

However, I still got the same error message as above, that I have 2035 observations and the number of random effects = 2035. I get that the number of observations should be more than the number of random effects being estimated, but I still don't understand why there are 2035 random effects in this model either. I would expect that since I am estimating a random slope and intercept for each unique ID (2 random effects for each ID), I would have 407*2=814 random effects. Clearly, there is something I am missing here about how to calculate the number of random effects for these two models. Any help understanding this would be highly appreciated. Thanks.

Answer 1

This is expected behavious whenever you try to fit a model with random slopes where the variable for the random slopes is categorical and there is only one observation per treatment/group combination.

This is because the levels of a categorical variable are represented by dummy variables - essentially they are treated as different variables. So in your case, when you fit random slopes only you are asking the software to estimate 5 random slopes for each group. When you fit random intercepts and random slopes there will be 407 random intercepts, but only 4 random slopes for each group (since one level will be treated as a reference group and included in the intercept), so either way you will have 5 x 407 random effects.

The only way to solve this is by either coding the variable as numeric, if that is plausible in your study/data, or not fitting random slopes, or having more than 1 observation per treatment per group.

It may be illustrative to see this with a toy dataset:

> set.seed(1)
> dt <- expand.grid(G = LETTERS[1:4], a = LETTERS[1:2])
> dt$Y = rnorm(nrow(dt))
> dt
  G a          Y
1 A A -0.6264538
2 B A  0.1836433
3 C A -0.8356286
4 D A  1.5952808
5 A B  0.3295078
6 B B -0.8204684
7 C B  0.4874291
8 D B  0.7383247

Now we fit the models, both of which will not run for the reasons explained above.

> lmer(Y ~ a + (0 + a | G), data = dt) %>% summary()
Error: number of observations (=8) <= number of random effects (=8) for term (0 + a | G); the random-effects parameters and the residual variance (or scale parameter) are probably unidentifiable
> lmer(Y ~ a + (1 + a | G), data = dt) %>% summary()
Error: number of observations (=8) <= number of random effects (=8) for term (1 + a | G); the random-effects parameters and the residual variance (or scale parameter) are probably unidentifiable

But now we add just 1 extra row to the dataset, and they run:

> (dt <- rbind(dt, dt[1, ]))
  G a          Y
1 A A -0.6264538
2 B A  0.1836433
3 C A -0.8356286
4 D A  1.5952808
5 A B  0.3295078
6 B B -0.8204684
7 C B  0.4874291
8 D B  0.7383247
9 A A -0.6264538

> lmer(Y ~ a + (0 + a | G), data = dt) %>% summary()

Random effects:
 Groups   Name Variance  Std.Dev.  Corr 
 G        aA   1.451e+00 1.205e+00      
          aB   3.224e-01 5.678e-01 -0.04
 Residual      4.239e-15 6.511e-08     

> lmer(Y ~ a + (1 + a | G), data = dt) %>% summary()


Random effects:
 Groups   Name        Variance  Std.Dev.  Corr 
 G        (Intercept) 9.776e-01 9.887e-01      
          aB          1.222e+00 1.105e+00 -0.81
 Residual             1.159e-14 1.077e-07      
Number of obs: 9, groups:  G, 4

In the model with random slopes only we have 2 random slopes in 4 groups (8 random effects), and in the model with both random intercepts and random slopes we have 4 random intercepts and 4 random slopes.