# Correct specification of longitudinal model in lme4

I am trying to fit a multilevel longitudinal model and i have a question regarding how to specify it.

The data consist of about 8k observations collected from about 3k individuals at four time points. Individuals are nested in groups and there are about 200 groups. I have two different types of fixed effects: (a) repeated measures at the observation level (e.g. pred1.obs ), and (b) group level predictors that also change over time (e.g. pred2.grp). Because each group level fixed effect is also longitudinal there are 800 values (4x200 which are repeated for each member of the group at that time) but there are only 200 groups.

My question is what would be the correct specification for this model and why? e.g:

1: lmer(outcome ~ time + pred1.obs + pred2.grp + (time|id) + (time|grp))

2: lmer(outcome ~ time + pred1.obs + pred2.grp + (time|id) + (1|grp:time))

3: lmer(outcome ~ time + pred1.obs + pred2.grp + (time|id) + (time|grp) + (1|grp:time))

Thus, would lme4 correctly estimate the model if i use (time|grp) or do i need to use (1|grp:time) or the combination?

Or something else that i haven't thought of?

Many thanks, George

• it would help a lot if you wrote precisely what model are you trying to fit in mathematical notation. Commented Apr 11, 2011 at 13:49
• Hmm, I think that part of the problem is that i am not sure what exactly what the model would be in notation. For an observation level predictor X, a repeated measure group level predictor Z, and for i observation, j individuals, and g groups, Although i am not sure, I think that what i need is something like this: $\ Y_{ijg}= \beta_{0j}+\beta_{1g}+\beta_{2j}time_{ij}+ \beta_{3g}time_{ig}+ \beta_{4}X_{ij}+e_{ij}$ $\ \beta_{0j} = \gamma_{00}+u_{0j}$ $\ \beta_{1g} = \gamma_{10}+\gamma_{11}Z_{g}+u_{0g}$ $\ \beta_{2j} = \gamma_{20}+u_{2j}$ $\ \beta_{3g} = \gamma_{30}+u_{3g}$ Commented Apr 11, 2011 at 15:38

I'll imagine a concrete example, with more context, to make things easy. Assume you measure the score on test of 3k students of 200 schools and you measured each student at 4 time points (say, at each quarter). You have a covariate at student level that doesn't vary by time (like sex), that you called pred1.obs and a covariate by school that vary by time (say the number of meetings between teachers and parents until that moment in time). If this example resembles your study, than I think you have to set up a three level model (individual level, group level and time level for the groups): i = 1 ... 3000 individuals t = 1... 4 periods g = 1... 200 groups

The model would be:

y_i ~ N(a + b_[groups_g] + b.ind*pred.obs1_i, sigma^2) # 1st level
b_g = N(gamma + gamma_[time] + gamma.g[time_t]*pred2.grp, sigma.b^2) # 2nd level
gamma.g_t = N(0, sigma.gamma^2) # 3rd level


Note that you would have the slope at the second level (group level) varying by time, which makes sense, since you expect that the effect of schools on the perfomance of students may vary by time, depending of the valu of the covariate at the level of schools. I'm not that sure how to estimate this with lmer (I know how to estimate a Bayesian model using WinBugs or Jags, calling them by R). In any case, here is my suggestion.

In lme4, I'd try: First, expand pred2.grp (the covariate at group level that vary by time) to the individual level, then you would have repetead measures by individuals at the group and time level. Then:

lmer(outcome ~ pred1.obs + pred2.grp + (1|group))


You have a large number of groups, so I speculate that (depending on the setting) you may think about group as a random effect, so your (…|grp) terms are probably justified. It may also be reasonable to associate random effects with individuals ((…|id) terms). However you have time as a covariate in all your models, so I assume you are looking for a linear time effect. In the same time having the (1|grp:time)) term or any kind of (…|time)) terms in the model makes the interpretation of the time covariate difficult.

The (1|…) terms correspond to random intercepts, e.g. including (1|grp) would estimate a “mean” outcome at time=0 for each group expressed as deviation from the grand intercept; (1|id) would estimate an individual intercept.

The (time|…) terms correspond to random slopes of the time covariate, e.g. including (time|grp) would estimate a deviation from the grand slope for each group; (time|id) would estimate individual slopes. [Edit: Please note the difference between the (time|id) and (0 + time|id) expressions, see ?lmer]

In your model specifications you are using both random intercepts and random slopes. Whether you need them or which you need depends on the relationships among the variables you study. If you know the relationship, you should specify the corresponding model. Alternatively you can fit multiple models and explain the differences among the results.