0
$\begingroup$

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.

$\endgroup$
3
  • $\begingroup$ Hi. I think your model is equivalent to gls(value ~ group*time , data=dat, correlation=corSymm(form= ~ 1 | Subject). Could you check ? If this is equivalent then I have a working code for multiple comparisons. $\endgroup$ Commented Aug 14, 2014 at 8:42
  • 1
    $\begingroup$ AFAIK you can create a new column int=interaction(Group, Time) and then run your model with int instead of Group*Time and also get the glht output for Int="Tukey" $\endgroup$
    – user45065
    Commented Aug 14, 2014 at 18:31
  • $\begingroup$ @StéphaneLaurent, I think the model is equivalent, if you have a working code it would be highly appreciated by me and probably by other users as well. Can you post it here ? Alternatively, if you do not want it to be publicized, you can send me privately and I will surely acknowledge it in the manuscript. $\endgroup$
    – TokyoUrban
    Commented Aug 18, 2014 at 1:04

1 Answer 1

2
$\begingroup$

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)
$\endgroup$
1
  • $\begingroup$ Excellent. Thanks a lot for the code and for the visualization tips. $\endgroup$
    – TokyoUrban
    Commented Aug 18, 2014 at 1:07

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.