2
$\begingroup$

I have read this question and answers:

Crossed vs nested random effects: how do they differ and how are they specified correctly in lme4?

however, I am struggling to understand why, provided that the factors are coded correctly, the nesting is equivalent as a random intercept for the interaction between the two factors, along with random intercepts for the upper level factor. The latter makes sense to me, but the former (the interaction) does not. I am looking for an intuitive explanation about why this is so.

As a specific example, using the answer above, suppose that we have classes nested within schools:

> library(lme4)
> dt <- read.table("http://bayes.acs.unt.edu:8083/BayesContent/class/Jon/R_SC/Module9/lmm.data.txt", header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE)

> dt$classID <- paste(dt$school, dt$class, sep=".")

> m0 <- lmer(extro ~ open + agree + social + (1 | school) + (1 | classID), data = dt)
> summary(m0)

> m1 <- lmer(extro ~ open + agree + social + (1 | school) + (1|classID:school), data = dt)
> summary(m1)

These 2 models are equivalent, and the output for the random effects sections are:

Random effects:
Groups   Name        Variance Std.Dev.
 classID  (Intercept)  8.2043  2.8643  
 school   (Intercept) 93.8409  9.6872  
 Residual              0.9684  0.9841  
Number of obs: 1200, groups:  classID, 24; school, 6

and

Random effects:
 Groups         Name        Variance Std.Dev.
 classID:school (Intercept)  8.2043  2.8643  
 school         (Intercept) 93.8409  9.6872  
 Residual                    0.9684  0.9841  
Number of obs: 1200, groups:  classID:school, 24; school, 6

Why is the random effect for classID in the first model identical to the random effect for the interaction classID:school in the 2nd ?

$\endgroup$
3
$\begingroup$

Why is the random effect for classID in the first model identical to the random effect for the interaction classID:school in the 2nd ?

That is because classID and classID:school are indeed the same, in this dataset. That is, the nesting structure is explicit in the way that the factor classID is coded - each level of classID "belongs" to one and only one level of school. That is precisely why the answer to that linked question has this line of code:

 > dt$classID <- paste(dt$school, dt$class, sep=".")

which creates a new factor for which the levels are unique to each level of school. As explained in that answer, the original dataset does not have explicit nesting - it could be fully cross-classified (or just "crossed" for short) or it could be fully nested.

By specifying random intercepts for the interaction between the two factors, we also make exactly the same structure as the line of code above, but without having to create a new variable to do so.

To see this more clearly:

> dt <- read.table("http://bayes.acs.unt.edu:8083/BayesContent/class/Jon/R_SC/Module9/lmm.data.txt", header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE)

> xtabs(~ school + class, dt)

        class
school  a  b  c  d
   I   50 50 50 50
   II  50 50 50 50
   III 50 50 50 50
   IV  50 50 50 50
   V   50 50 50 50
   VI  50 50 50 50

there is an ambiguity about whether these are (fully) nested, or (fully) cross-classified (crossed). If we know that they are nested then we specify the random intercepts in lme4 as:

(1|School/Class) or equivalently (1|School) + (1|Class:School)

whereas if they are crossed then we specify the random intercepts as

(1|School) + (1|Class)

Notice that, in the answer given to the linked question, this problem is resolved, if the data are nested, by introducing a new variable, classID such that each class ID is unique to the school to which it "belongs" (which is the definition of nesting):

> dt$classID <- paste(dt$school, dt$class, sep=".")
> xtabs(~ school + classID, dt)

      classID
school I.a I.b I.c I.d II.a II.b II.c II.d III.a III.b III.c III.d IV.a IV.b
   I    50  50  50  50    0    0    0    0     0     0     0     0    0    0
   II    0   0   0   0   50   50   50   50     0     0     0     0    0    0
   III   0   0   0   0    0    0    0    0    50    50    50    50    0    0
   IV    0   0   0   0    0    0    0    0     0     0     0     0   50   50
   V     0   0   0   0    0    0    0    0     0     0     0     0    0    0
   VI    0   0   0   0    0    0    0    0     0     0     0     0    0    0

      classID
school IV.c IV.d V.a V.b V.c V.d VI.a VI.b VI.c VI.d
   I      0    0   0   0   0   0    0    0    0    0
   II     0    0   0   0   0   0    0    0    0    0
   III    0    0   0   0   0   0    0    0    0    0
   IV    50   50   0   0   0   0    0    0    0    0
   V      0    0  50  50  50  50    0    0    0    0
   VI     0    0   0   0   0   0   50   50   50   50

which shows that each level of class occurs only in one level of school.

However, if instead, we form the interaction between class and school manually:

> dt$int.school.class <- interaction(dt$school, dt$class)

> xtabs(~ school + int.school.class, dt)

then we obtain exactly the same cross-tabulation:

> xtabs(~ school + int.school.class, dt)
      int.school.class
school I.a II.a III.a IV.a V.a VI.a I.b II.b III.b IV.b V.b VI.b I.c II.c III.c
   I    50    0     0    0   0    0  50    0     0    0   0    0  50    0     0
   II    0   50     0    0   0    0   0   50     0    0   0    0   0   50     0
   III   0    0    50    0   0    0   0    0    50    0   0    0   0    0    50
   IV    0    0     0   50   0    0   0    0     0   50   0    0   0    0     0
   V     0    0     0    0  50    0   0    0     0    0  50    0   0    0     0
   VI    0    0     0    0   0   50   0    0     0    0   0   50   0    0     0
      int.school.class
school IV.c V.c VI.c I.d II.d III.d IV.d V.d VI.d
   I      0   0    0  50    0     0    0   0    0
   II     0   0    0   0   50     0    0   0    0
   III    0   0    0   0    0    50    0   0    0
   IV    50   0    0   0    0     0   50   0    0
   V      0  50    0   0    0     0    0  50    0
   VI     0   0   50   0    0     0    0   0   50

[The cross tabulations might look different at first glance but this is solely because the ordering of the columns is different]

When we specify the interaction in the random part of the model formula (1|Class:School) that is what lme4 is doing internally - without the need to creating a new factor.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.