I am trying to estimate a longitudinal model in R in which there are several random intercepts that are correlated with each other, and the data are non-nested. For example, consider a simple longitudinal model in which student test scores are regressed on three random intercepts, one random school effect, one random teacher effect, and one random teacher-school match effect. The data are non-nested because teachers frequently switch schools in the longitudinal data. I am primarily interested in the random effects parameters themselves, and so I want to specify all three random effects and allow them to be correlated with each other. There are no random slopes, only random intercepts.

I tried:

model1 <- lmer(test_score ~ FEs + (1|schoolid) + (1|teacherid) + (1|matchid), data)

but this appears to force the random effects to be uncorrelated with each other. Is there a way, using lmer or any other R package to estimate a similar model in which the random intercepts are allowed to be correlated with each other (using non-nested data)?

More specifically, the above model imposes the covariance structure: $$ Cov \left[ \begin{array} \\ \theta_{s} \\ \theta_{t} \\ \theta_{m} \\ \end{array} \right] = \left[ \begin{array} \\ \sigma^{2}_{\theta_{s}}I_{s} & 0 & 0 \\ 0& \sigma^{2}_{\theta_{t}}I_{t} & 0\\ 0&0&\sigma^{2}_{\theta_{m}}I_{m}\\ \end{array} \right] $$

What I would like to estimate instead is a random effects model with the following covariance structure:

$$ Cov \left[ \begin{array} \\ \theta_{s} \\ \theta_{t} \\ \theta_{m} \\ \end{array} \right] = \left[ \begin{array} \\ \sigma^{2}_{\theta_{s}}I_{s} & \sigma^{2}_{\theta_{st}}C_{st} & \sigma^{2}_{\theta_{sm}}C_{sm} \\ \sigma^{2}_{\theta_{ts}}C_{ts} & \sigma^{2}_{\theta_{t}}I_{t} & \sigma^{2}_{\theta_{tm}}C_{tm}\\ \sigma^{2}_{\theta_{ms}}C_{ms} &\sigma^{2}_{\theta_{mt}}C_{mt} &\sigma^{2}_{\theta_{m}}I_{m}\\ \end{array} \right] $$ where $C_{st}$, for example, is a matrix of zeros and ones, with $C(i,j)=1$ where there is an observation of the pair $(school=i, teacher=j)$ in the data, and 0 if no such pair is observed.

  • $\begingroup$ If your interest is in the intercepts themselves, you may want to use fixed effects. W/ random effects you can't get estimates of the paramters, only predictions. Are you wanting to test the variances? I'm not sure this makes sense. $\endgroup$ Commented Sep 8, 2014 at 19:37
  • $\begingroup$ The model is not identifiable with fixed effects unless the matchid is assumed to be uncorrelated with schoolid and teacherid, which is precisely the assumption I am trying to avoid. The matchid would absorb all of the variation from schoolid and teacherid in this case, leaving nothing to identify those parameters. $\endgroup$ Commented Sep 8, 2014 at 19:45
  • $\begingroup$ In response to question about variances: Yes, by parameters I am referring to variance parameters. I want to know what fraction of the variation in test scores is attributable to schools vs teachers vs school-teacher pairs (where there is reason to believe that certain teachers have comparative advantages at certain schools). For further motivation on these types on problems see, for example, link $\endgroup$ Commented Sep 8, 2014 at 19:45
  • $\begingroup$ Hmmm. Do you know if there's an existing LMM framework (MLwin, HLM, etc.) that handles this structure? At a quick glance, the linked document estimates independent variances ... $\endgroup$
    – Ben Bolker
    Commented Sep 9, 2014 at 23:05


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.