I'm new to linear mixed-effects models and I was wondering if I could get some help in getting my model to properly work.

I have an example dataset:

data_ex <- data.frame( pnum = c(1,1,1,2,2,2,3,3,3,4,4,4,5,5,5,6,6,6,7,7,7,8,8,8,9,9,9,10,10,10),
                   group = c("1","1","1","2","2","2","3","3","3","1","1","1","2","2","2","3","3","3","1","1","1","2","2","2","3","3","3","1","1","1"),
                   day = c(1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3),
                   score_x = floor(runif(30, min=0, max=101)),
                   score_y = floor(runif(30, min=0, max=101)))

My research question: I'm interested in if there are group differences in the effect of x on y. However, each participant has been measured three times (on three days). I however do not care about this effect.

What I first did was just check the group differences on the effect of x on y like so:

lm(score_y ~ score_x * group, data_ex)

However then I realized that I am inflating my data by not accounting for the repeated measures of my participants.

I opted for aggregating my entire data across variable "day", but then I would also lose a lot of data.

Thus I wanted to try mixed-effects models. As I understand, I could account for days and participants as random slopes. I'm however not entirely sure if this would be correct for my research question.

Would I be able to answer my research question if I model my data like this?

m <- lmer(score_y  ~ score_x * group + (1 + score_x | pnum) + (1|day), data = data_ex)
m_small <- lmer(score_y ~ score_x + (1 + score_x | pnum) + (1|day), data = data_ex)
anova(m, m_small)

Thank you for your help!


2 Answers 2


I find useful to proceed as follows:

  1. Set the fixed effects. These are your predictors of interest + those that you think should be controlled for (= covariates).
  2. Among the selected fixed effects, identify those that are within-subjects and add them as by-participant random slopes. For example:
m <- lmer(y  ~ group * (within1 + within2) + (1 + within1 + within2 | pnum)

In your case, $score\_x$ is a repeatedly measured continuous variable, not an experimental factor, so you don't need to add it as a random slope.

  1. Identify other random effects, such as stimuli. Typical examples include words in a psycholinguistic experiment, or emotional pictures in psychology. You can think about it this way: just like your subjects are a small sample of the general population and you want to generalize beyond your specific sample, a set of stimuli might be a small sample of a general class (pleasant images, abstract words, etc.) and you want to draw conclusions on the general class. You can then add a random intercept for stimuli, as well as random slopes for any within-stimulus factor (which might also be within-subject, depending on your design). For example:
m <- lmer(y  ~ group * (within1 + within2) + (1 + within1 + within2 | pnum) + (1 + within2 | stimulus)

In your case, considering that you do not seem interested in the effect of $day$, an appropriate model would be:

m <- lmer(score_y  ~ group * score_x + (1 | pnum)

For a good introduction to mixed models in psychology, I can recommend Singmann & Kellen 2017.

  • $\begingroup$ The model with ... + (1 | pnum) assumes homogeneous compound symmetry, which may be too restrictive in this case. $\endgroup$ Jun 22, 2019 at 18:31

If you have a sufficient number of participants (i.e., more than 60) that have data for all three time points, then the safest option would be to fit a model with unstructured covariance matrix for the error terms. This can be done with function gls() from package nlme, e.g.,

fm <- gls(score_y ~ score_x * group, data = data_ex,
          correlation = corSymm(form = ~ 1 | pnum),
          weights = varIdent(form = ~ 1 | day))

  • 1
    $\begingroup$ (+1) I'm curious in what ways this model is superior to a random effects model? In what way is it "safest"? Or in general: When would you choose gls over lme or lmer? $\endgroup$ Jun 22, 2019 at 6:42
  • 3
    $\begingroup$ @COOLSerdash (1/2) when we have clustered data the aim is to appropriately model the covariance matrix. In case of balanced data, as in this example, and sufficient number of subjects, the best thing you can do is fit a completely unstructured matrix estimating all variances & covariances. Mixed models capture these covariances with random effects but they do impose a specific structure to them. I.e., they are more restricted than an unstructured matrix. $\endgroup$ Jun 22, 2019 at 18:26
  • 3
    $\begingroup$ @COOLSerdash (2/2) However, often in practice we have unbalanced data and/or complex designs with nested or crossed grouped data. In these cases it is practically impossible to fit an unstructured covariance matrix. And it is in these cases where random effects provide a practical and flexible compromise. $\endgroup$ Jun 22, 2019 at 18:28
  • $\begingroup$ Thanks a lot for these explanations, I appreciate it. It's very interesting as I have always assumed that mixed models with random effects are the gold standard for repeated/clustered data. $\endgroup$ Jun 22, 2019 at 18:30
  • 1
    $\begingroup$ You may find more information regarding these issues in my course notes: drizopoulos.com/courses/EMC/CE08.pdf $\endgroup$ Jun 22, 2019 at 18:37

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.