I am modeling change over time in group psychotherapy subjects using R and lme4. My data have the following structure:

  • subject (id)
  • time (code 1-10 for equally spaced repeated measures)
  • outcome (for every repeated measure)
  • treatment (0/1 for psychotherapy/waiting list control)

My first two-level model with random slopes and intercepts works well and is simple:

lmer(outcome ~ time * treatment + (time | subject), data=data, REML=FALSE)

Now I was wondering if I should use a three-level partially nested model because the group psychotherapy subjects are nested within therapist (there where several therapist providing the treatment), but controls were non-nested. I am guessing I should account at least for main effects of therapists as argued by De Jong, Moerbeek & Van der Leeden (2010).

De Jong, K., Moerbeek, M., & Van der Leeden, R. (2010). A priori power analysis in longitudinal three-level multilevel models: an example with therapist effects. Psychotherapy Research, 20(3), 273-284.

I have found a very useful resource in the following link under "Partially nested models": http://rpsychologist.com/r-guide-longitudinal-lme-lmer

The author gives the following code for a model virtually identical to the one I would like to test:

lmer(outcome ~ time * treatment + (1 | group:subject) + (0 + time | therapist:subject) + (0 + time:treatment | group) + (0 + treatment | group), data=data) 

The data he presents is virtually identical but he adds in the model the "group" variable. I don't understand why this is because treatment/control grouping is the same as treatment/nontreatment grouping. If a subject received treatment than he is in the experimental group, if not he in the control group. How would you write this three-level partially nested model? I feel confused. Thank you.


1 Answer 1


Question was eloquently answered by Thierry on Stackoverflow.

I post his response, maybe it will be helpful for some one.

Keep it simple. Just add a dummy therapist none to the subjects without therapist. Then fit the model below.

lmer(outcome ~ time * treatment + (time | therapist/subject), data=data)

therapist none is confounded with treatment waiting list. Therapist is a random effect and thus penalised. Treatment is a fixed effect and not penalised. Hence all information will go to the treatment effect for waiting list and the therapist effect for none will be zero.


Your Answer

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

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