R - Model specification for repeated measures in linear mixed-effect model (lmer)

Question

I am quite new in linear mixed-effects models, and trying to figure out how to properly design my model according to what I want to show (especially with the nesting and random effect). I am reading several tutorials and examples in forums, but it’s still a bit confusing. My several attempts ended up obtaining very different results from a simple repeated ANOVA. Therefore, I suspect that I am doing something wrong.

The experimental design is quite simple:

• I have two groups of subjects (experimental group and control group), each subject receiving one specific treatment.
• In each group, I have males and females.
• I measure the treatment response after 1, 2, and 3 months (three time points).

Ultimately, I would like to estimate:

• The effect of treatment on my measure (fixed)
• The effect of sex on my measure (fixed)
• Whether my measure changes over time (over the three months) (random?)
• The interaction between treatment, sex and change over time

|| ID |  Treat  | Sex | Month | Measure ||
||----|---------|-----|-------|---------||
|| s1 | treat1  |  M  |   1   |   324   ||
|| s1 | treat1  |  M  |   2   |   328   ||
|| s1 | treat1  |  M  |   3   |   379   ||
|| s2 | treat1  |  F  |   1   |   327   ||
|| s2 | treat1  |  F  |   2   |   360   ||
|| s2 | treat1  |  F  |   3   |   399   ||
|| s3 | treat2  |  M  |   1   |   336   ||
|| s3 | treat2  |  M  |   2   |   431   ||
|| s3 | treat2  |  M  |   3   |   351   ||
|| s4 | treat2  |  F  |   1   |   320   ||
|| s4 | treat2  |  F  |   2   |   305   ||
|| s4 | treat2  |  F  |   3   |   376   ||

I would like to avoid having the three levels of months in my model, but rather a global effect of time over my measurements. I guess the my variable Month is nested within treatment, but my different trials to compute the model failed.

model <- lmer(Measure ~ Treat * Sex * Month + (1 + Treat * Sex | Month), data)

Happy to provide more details...

Thank you very much for you help and insight!

EDIT 1:

Many thanks for the help. I understand much better what I need to do.

I tried to model my data using each individual as random effect.

model <- lmer(Score ~ Treat * Sex * Month + (1 | ID), data)

This is what I obtain:

Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Score ~ Treat * Sex * Month + (1 | ID)
   Data: data

REML criterion at convergence: 271.7

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-1.78970 -0.52428  0.03797  0.44426  2.62502 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 1.068    1.033   
 Residual             1.166    1.080   
Number of obs: 84, groups:  ID, 31

Fixed effects:
                      Estimate Std. Error       df t value Pr(>|t|)    
(Intercept)           15.08581    0.58668 61.55813  25.714   <2e-16 ***
TreatEXP              -0.83115    0.86683 62.39386  -0.959   0.3413    
SexM                   0.26652    0.82943 61.58985   0.321   0.7490    
Month8                 0.04437    0.65994 50.90799   0.067   0.9467    
Month12                1.23650    0.59705 47.52705   2.071   0.0438 *  
TreatEXP:SexM         -0.02797    1.17231 60.15952  -0.024   0.9810    
TreatEXP:Month8       -1.52711    0.92560 49.28606  -1.650   0.1053    
TreatEXP:Month12      -1.79569    0.88186 47.58390  -2.036   0.0473 *  
SexM:Month8            1.24144    0.88969 49.31780   1.395   0.1692    
SexM:Month12           1.14472    0.86380 48.26876   1.325   0.1913    
TreatEXP:SexM:Month8  -0.20806    1.22646 48.05375  -0.170   0.8660    
TreatEXP:SexM:Month12 -1.37308    1.20781 47.48363  -1.137   0.2613

I’m not quite sure why I have Month8 and Month12 instead of just Month. I rather would like to see an overall effect of time on my dependent variable Score, a bit what I obtain when I perform an ANOVA.

fit_all <- aov_ez('ID', 'Score', data, between = c('Treat', 'Sex'), within = c('Month'))

And this is what I obtain with the ANOVA:

Univariate Type III Repeated-Measures ANOVA Assuming Sphericity

                 Sum Sq num Df Error SS den Df   F value    Pr(>F)    
(Intercept)     13625.4      1   82.352     18 2978.1552 < 2.2e-16 ***
Treat              80.5      1   82.352     18   17.5897 0.0005454 ***
Sex                14.9      1   82.352     18    3.2582 0.0878224 .  
Treat:Sex           0.6      1   82.352     18    0.1407 0.7119409    
Month               8.9      2   31.631     36    5.0489 0.0116734 *  
Treat:Month        27.5      2   31.631     36   15.6423 1.291e-05 ***
Sex:Month           2.3      2   31.631     36    1.2852 0.2889756    
Treat:Sex:Month     0.6      2   31.631     36    0.3276 0.7228025

As you can see there is a massive difference between LME model and ANOVA...

Again, thank you very much for your comments, they are extremely helpful.

Answer 1

I think what you want is:

model <- lmer(Measure ~ Treat * Sex * Month + (1|ID), data)

This is what the above model fits:

1. Fixed effect for treatment, sex, time and the 3-way interaction between these terms (It's unclear from your question if you really want a three way interaction or if you're actually looking for two 2-way interactions). You can use time as a continuous covariate (not categorical). This will allow you to model the effect of time on your outcome in a linear manner. The estimated coefficient for Time will also tell you if the outcome is changing over time (although you have to take some care in interpretation given the many interaction terms).

2. A random intercept for each ID via 1|ID. As your unit of analysis is a person, this is the most relevant random effect to include. One way of thinking about that is that your experimental sample is randomly drawn from a larger population and thus, each individual represents a realized random sample from this population. The random intercept allows the intercept to vary for each patient and gives each patient their own intercept. Scientifically, we can think of it as representing inter-individual variability.

The model you have currently specified fits a random intercept for "Month". You cannot fit a fixed effect and a random intercept for the same variable as they are competing for the same information [1]. Moreover, that model isn't specifying nested random effects and at any rate, thinking of month as being nested within treatment is not correct. Nesting occurs when a level of a categorical variable belongs to one (and only one) upper level factor. The standard example here is students within a classroom. Any student you pick could belong to only one classroom. So you can think of students as being "nested" within classrooms. [2] For your experimental design you can think of the nesting structure as: Treatment -> Individual. But we generally don't model treatments as random effects because in most cases the treatments don't represent a random sample of possible treatments.