Solved – $\textit{Post hoc}$ test after repeated measures ANOVA (LME + Multcomp)

anovapost-hocrepeated measures

I have a dataset where the samples are dependent, belong to different groups, and measurements were taken over time. It looks like this:

Subject Group   Time    Value
S1  G1  12h 5.55
S1  G1  24h 7.63
S1  G1  36h 9.88
S2  G2  12h 3.26
S2  G2  24h 4.57
S2  G2  36h 6.44
S3  G3  12h 3.23
S3  G3  24h 4.10
S3  G3  36h 5.57
S4  G1  12h 5.65
S4  G1  24h 7.89
S4  G1  36h 10.43
S5  G2  12h 4.18
S5  G2  24h 4.93
S5  G2  36h 6.70
S6  G3  12h 3.53
S6  G3  24h 4.52
S6  G3  36h 6.25
S7  G1  12h 6.38
S7  G1  24h 8.68
S7  G1  36h 11.30
S8  G2  12h 4.73
S8  G2  24h 5.16
S8  G2  36h 6.75
S9  G3  12h 3.92
S9  G3  24h 4.58
S9  G3  36h 6.91

I would like to compare the groups using repeated measures ANOVA, and after that check where the differences are using a post hoc test.

So far, I have used the NLME package,

lme_data <- lme(Value~Group*Time, data=data, random = ~1| Subject)

Followed by:

summary(glht(lme_data, linfct=mcp(Group = "Tukey"), test = adjusted(type = "bonferroni")))

or

summary(glht(lme_data, linfct=mcp(Time = "Tukey"), test = adjusted(type = "bonferroni")))

that I saw in a previous post.

The question is, while the glht works fine for "Time" and "Group" separately, I would like to check the Time:Group interaction, to know which groups are different AND where those differences are.

Doing it with ANOVA and TukeyHSD is straightforward for independent samples, but somehow for repeated measures ANOVA I have been struggling.

Best Answer

I would like to check the Time:Group interaction, to know which groups are different AND where those differences are.

If I understand your question correctly, you would like to tease apart the interaction effect between Group and Time. One possible approach is to perform various tests for all the combinations of the two factors (and maybe plot them out with the effect estimates and their standard errors). For example,

library(phia)    
testInteractions(lme_data, custom=list(Group=c(1,0,0), Time=c(1,0,0)))
testInteractions(lme_data, custom=list(Group=c(0,1,0), Time=c(1,0,0)))
testInteractions(lme_data, custom=list(Group=c(0,0,1), Time=c(1,0,0)))

show the effects for each group at 12h of Time and their significance. Similarly, with

testInteractions(lme_data, custom=list(Group=c(1,0,0), Time=c(1,0,0)))
testInteractions(lme_data, custom=list(Group=c(1,0,0), Time=c(0,1,0)))
testInteractions(lme_data, custom=list(Group=c(1,0,0), Time=c(0,0,1)))

you obtain the effects at each time point for group G1. Furthermore, you can test all the pairwise comparisons,

testInteractions(lme_data, custom=list(Group=c(1,-1,0), Time=c(1,0,0)))
testInteractions(lme_data, custom=list(Group=c(1,0,-1), Time=c(1,0,0)))
testInteractions(lme_data, custom=list(Group=c(0,1,-1), Time=c(1,0,0)))
...

testInteractions(lme_data, custom=list(Group=c(1,0,0), Time=c(1,-1,0)))
testInteractions(lme_data, custom=list(Group=c(1,0,0), Time=c(1,0,-1)))
testInteractions(lme_data, custom=list(Group=c(1,0,0), Time=c(0,1,-1)))
...

With all these effects combined, you should be able to have a detailed picture about the interaction.

To visualize these effects, plot them out with:

library(effects)
plot(allEffects(lme_data))

and

library(lsmeans)
lsmip(lme_data, Group~Time)
lsmip(lme_data, Time~Group)
Related Question