2
$\begingroup$

I want to fit a random slope and intercept model to my data (see below), with rating being regressed on position. The plot below shows the results for all subjects. The blue lines are separate linear regressions fitted to each individual. The red lines show the results from the mixed effects model. I am a bit surprised, that the difference between the lines is quite big in some cases. For example, the lines for subject 4 or 7 are very different. In the full dataset, the difference is sometimes even more pronounced.

Fitting one mixed effects model or many single models is a different thing, of course. Still I do not quite understand why such big differences arise. I am not very familiar with mixed models, so my naive expectation was that the two lines would not deviate too much, given that slopes and intercepts can vary freely.

Am I wrong to expect smaller deviations and why do they arise? (or did I go wrong at some point)?

enter image description here

The R code for the model and plot follows. The data is found below.

library(nlme)
library(ggplot2)

m <- lme(rating ~ 1+ position, random= ~ 1 + position | subject, data=x)

ggplot(x, aes(position, rating)) +
  geom_point(color="grey") + 
  geom_smooth(method="lm", se = FALSE) +
  geom_line(aes(y=predict(m)), color="red") +
  facet_wrap( ~ subject) 

Data:

x <- structure(list(position = c(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 
1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 
2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 
3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 
4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 0, 1, 2, 3, 4, 
5, 6, 7, 8, 9), subject = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 
10L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 
6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 8L, 8L, 
8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 
9L, 9L), .Label = c("1", "2", "3", "4", "5", "6", "7", "8", "9", 
"10"), class = "factor"), rating = c(2, 3, 4, 3, 2, 4, 4, 3, 
2, 1, 4, 2, 4, 3, 2, 3, 4, 3, 5, 4, 3, 2, 2, 2, 4, 2, 4, 3, 2, 
2, 3, 5, 3, 4, 4, 4, 3, 2, 3, 5, 4, 5, 2, 3, 4, 2, 4, 4, 1, 2, 
4, 5, 4, 2, 3, 4, 3, 2, 2, 2, 4, 5, 4, 4, 5, 2, 3, 4, 3, 2, 4, 
3, 4, 4, 4, 3, 4, 5, 4, 5, 4, 3, 4, 5, 4, 5, 4, 5, 5, 5, 1, 3, 
3, 4, 3, 3, 5, 3, 5, 3)), .Names = c("position", "subject", "rating"
), row.names = c(NA, -100L), class = "data.frame")
$\endgroup$
  • $\begingroup$ "given that slopes and intercepts can vary freely" That's not correct for mixed models. They can vary, but not freely. A distribution (IIRC multivariate Gaussian) is assumed. You might want to study this and other literature. Your lm model has more degrees of freedom and can thus vary more "freely". $\endgroup$ – Roland Apr 4 '16 at 8:24

Your Answer

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

Browse other questions tagged or ask your own question.