# What is a valid post-hoc analysis for a three-way repeated measures ANOVA?

I've performed a three-way repeated measures ANOVA; what post-hoc analyses are valid?

This is a fully balanced design (2x2x2) with one of the factors having a within-subjects repeated measure. I'm aware of multivariate approaches to repeated measures ANOVA in R, but my first instinct is to proceed with a simple aov() style of ANOVA:

aov.repeated <- aov(DV ~ IV1 * IV2 * Time + Error(Subject/Time), data=data)

DV = response variable

IV1 = independent variable 1 (2 levels, A or B)

IV2 = independent variable 2 (2 levels, Yes or No)

IV3 = Time (2 levels, Before or After)

Subject = Subject ID (40 total subjects, 20 for each level of IV1: nA = 20, nB = 20)

summary(aov.repeated)

Error: Subject
Df Sum Sq Mean Sq F value   Pr(>F)
IV1       1   5969  5968.5  4.1302 0.049553 *
IV2       1   3445  3445.3  2.3842 0.131318
IV1:IV2   1  11400 11400.3  7.8890 0.007987 **
Residuals 36  52023  1445.1
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Error: Subject:Time
Df Sum Sq Mean Sq F value   Pr(>F)
Time            1    149   148.5  0.1489 0.701906
IV1:Time        1    865   864.6  0.8666 0.358103
IV2:Time        1  10013 10012.8 10.0357 0.003125 **
IV1:IV2:Time    1    852   851.5  0.8535 0.361728
Residuals      36  35918   997.7
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Alternatively, I was thinking about using the nlme package for a lme style ANOVA:

aov.repeated2 <- lme(DV ~ IV1 * IV2 * Time, random = ~1|Subject/Time, data=data)
summary(aov.repeated2)

Fixed effects: DV ~ IV1 * IV2 * Time
Value Std.Error DF   t-value p-value
(Intercept)                      99.2  11.05173 36  8.975972  0.0000
IV1                              19.7  15.62950 36  1.260437  0.2156
IV2                              65.9  15.62950 36  4.216385  0.0002 ***
Time                             38.2  14.12603 36  2.704228  0.0104 *
IV1:IV2                         -60.8  22.10346 36 -2.750701  0.0092 **
IV1:Time                        -26.2  19.97722 36 -1.311494  0.1980
IV2:Time                        -57.8  19.97722 36 -2.893295  0.0064 **
IV1:IV2:Time                     26.1  28.25206 36  0.923826  0.3617

My first instinct post-hoc of significant 2-way interactions with Tukey contrasts using glht() from multcomp package:

data$IV1IV2int <- interaction(data$IV1, data$IV2) data$IV2Timeint <- interaction(data$IV2, data$Time)

aov.IV1IV2int <- lme(DV ~ IV1IV2int, random = ~1|Subject/Time, data=data)
aov.IV2Timeint <- lme(DV ~ IV2Timeint, random = ~1|Subject/Time, data=data)

IV1IV2int.posthoc <- summary(glht(aov.IV1IV2int, linfct = mcp(IV1IV2int = "Tukey")))
IV2Timeint.posthoc <- summary(glht(aov.IV2Timeint, linfct = mcp(IV2Timeint = "Tukey")))

IV1IV2int.posthoc
#A.Yes - B.Yes == 0        0.94684
#B.No - B.Yes == 0         0.01095 *
#A.No - A.Yes == 0         0.80785
#A.No - B.No == 0          0.00346 **

IV2Timeint.posthoc
#No.After - Yes.After == 0           0.0142 *
#Yes.Before - Yes.After == 0         0.0558 .
#No.Before - No.After == 0           0.1941
#No.Before - Yes.Before == 0         0.8616

The main problem I see with these post-hoc analyses are some comparisons that aren't useful for my hypotheses.

Any suggestions for an appropriate post-hoc analysis are greatly appreciated, thanks.

-
Your random-effect model looks strange: / is used to denote nesting (as typically seen in a split-plot experiment), unlike its use in the Error term of aov() where it mainly indicates how to build error strata. –  chl May 14 '12 at 20:09
@chl I formatted the Error term of aov() this way to specify that Time is the within-groups factor. From Baron Error(subj/(color + shape)) seems to be used in the same way. –  RobJackson28 May 14 '12 at 20:24
@chl Thank you for bringing up the lme model, I'm unclear on the proper usage of /. How would you specify Time as the within-groups factor as in Error() with aov()? –  RobJackson28 May 14 '12 at 20:46

I think statisticians will tell you that there is always a problem with any post hoc analysis because seeing the data may influence what you look at and you could be biased becuase you are hunting for significant results. The FDA in clinical trial studies requires that the statistical plan be completely spelled out in the protocol. in a linear model you certainly could prespecify the contrasts that you would like to look at in the event that the ANOVA or ANCOVA finds an overall difference. Such prespecified contrasts would be fine to look at as long as the usual treatment for multiplicity is also part of it.

-
This is pretty much the problem I'm having, given the analysis was handed over to me (without any a priori statistical planning other than "let's run a bunch of t-tests"). I've managed to distil the basic hypotheses that were originally intended, but I'm having a little trouble with the post-hoc syntax. Justifying all these steps to the experimenter is my main goal, so as to avoid the t-test dogma. The ultimate goal: make statistical planning a necessity for future experimental designs. For the time being, I've got to work with what I have. –  RobJackson28 May 14 '12 at 20:39
Then all I would add is that if you want to go ahead with post-hoc analyses, I see no problem as long as you do proper multiplicity adjustment. –  Michael Chernick May 14 '12 at 20:55
Am I correct in assuming that multiciplicity adjustments are analogous to family-wise error adjustments? E.g., Tukey's, Bonferroni, etc.? –  RobJackson28 May 14 '12 at 21:00
Exactly right. bootstrap and permutation methods are also available in SAS for example with Proc MULTTEST. See the work of Westfall and Young. –  Michael Chernick May 14 '12 at 21:11
Thank you for the assistance @Michael, I appreciate it. However, I'm still unclear on the syntax to use in R. Specifically, I'm unsure if it's most appropriate to manually specify contrast matrices for the relevant Tukey contrasts using glht(), or to perform all the comparisons by default. Additionally, I'm not sure how to properly deal with the repeated measure of Time in terms of post-hoc. –  RobJackson28 May 15 '12 at 1:18
show 1 more comment

If you have a software package like SAS you would probably use proc mixed to do the the repeated measures mixed model and if you specify which contrast you want to use SAS will handle it properly for you. You may also be able to do it with the repeated option in PROC GLM but be careful because they behave differently and make different assumptions. The repeated observations are usually correalted because they have something common. I often have repeated measures on the same patient at different time points. So in computing the contrasts the covariance terms enter into the problem.

-