[R] Tukey HSD (or other post hoc tests) following repeated measures ANOVA

Mark Difford mark_difford at yahoo.co.uk
Fri May 30 15:48:20 CEST 2008 

Hi Ullrich,

>> # what does '~1 | Subj/Cond' mean?

It is equivalent to your aov() error structure [ ... +Error(Subj/Cond) ]. 
It gives you a set of random intercepts, one for each level of your nesting
structure.

## To get some idea of what's being done it helps to have a continuous
covariate in your model.
## So look at this example:
mod1 <- lme(distance ~ age, Orthodont, random = ~ 1 | Subject)  ## random
intercepts, parallel slopes
mod2 <- lme(distance ~ age, Orthodont, random = ~ age | Subject)  ## random
intercepts, separate slopes

plot(augPred(mod1, primary=~age))  ## parallel slopes
plot(augPred(mod2, primary=~age))  ## separate slopes

In my first post I used the example I did because it's used in V&R's MASS (:
283--286), where the equivalence of the two methods is illustrated. So there
really can be no argument about whether lme() can do what aov() does.  But
lme() [and the newer, still somewhat "hard-hat" version, lmer()] can do a
whole lot more than aov() can do.  So use lme() [or lmer() if you can't get
lme() to converge].

Also, you would be much better off using the multcomp package to do
TukeyHSD. Amongst its other advantages are that you can use not only the
methods it offers for adjusting p-values, but everything offered in
p.adjust().
?p.adjust.

HTH, Mark.


Ullrich Ecker wrote:
> 
> Great, thanks a lot, guys!
> 
> Only thing is I now have two working versions, that yield *slightly* 
> different results.
> 
> 
> ACCaov <- aov(Acc ~ Cond + Error(Subj/Cond), WMU3C)
> ACClme <- lme(Acc ~ Cond, random = ~1 | Subj/Cond, WMU3C) # what does 
> '~1 | Subj/Cond' mean?
> summary(glht(ACClme, linfct=mcp(Cond="Tukey")))
> 
> yielding
> 
> Linear Hypotheses:
>             Estimate Std. Error z value p value
> 2 - 1 == 0 0.392560   0.027210  14.427  <1e-05 ***
> 3 - 1 == 0 0.400372   0.027210  14.714  <1e-05 ***
> 3 - 2 == 0 0.007812   0.025442   0.307    0.95
> 
> 
> and
> 
> 
> ACCaov <- aov(Acc ~ Cond + Error(Subj/Cond), WMU3C)
> ACCaov2 <- aov(terms(Acc ~ Subj + Cond, WMU3C))  # gives same result 
> as 1st aov, but yields aov not aovlist
> ACCtukey <- TukeyHSD(ACCaov2, "Cond"); ACCtukey
> 
> yielding
> 
> $Cond
>           diff         lwr        upr     p adj
> 2-1 0.3637756  0.29889404 0.42865721 0.0000000
> 3-1 0.3715881  0.30670654 0.43646971 0.0000000
> 3-2 0.0078125 -0.05329192 0.06891692 0.9504659
> 
> 
> I am trying to work my way through this (so I'm one of the good guys, 
> at least I'm trying to understand my stats... ;-)), but any hint to 
> what solution may be more appropriate to my problem would be much
> appreciated.
> 
> Cheers,
> 
> Ullrich
> 
> 
> 
> At 03:16 PM 30/05/2008, you wrote:
> 
>>Hi Ullrich,
>>
>> >> The model is
>> >> RT.aov <- aov(RT~Cond + Error(Subj/Cond), WMU3C)
>> >> I understand that TukeyHSD only works with an aov object, but that
>> >> RT.aov is an aovlist object.
>>
>>You want to use lme() in package nlme, then glht() in the multcomp
package.
>>This will give you multiplicity adjusted p-values and confidence
intervals.
>>
>>## Example
>>require(MASS)         ## for oats data set
>>require(nlme)         ## for lme()
>>require(multcomp)  ## for multiple comparison stuff
>>
>>Aov.mod <- aov(Y ~ N + V + Error(B/V), data = oats)
>>Lme.mod <- lme(Y ~ N + V, random = ~1 | B/V, data = oats)
>>
>>summary(Aov.mod)
>>anova(Lme.mod)
>>
>>summary(Lme.mod)
>>summary(glht(Lme.mod, linfct=mcp(V="Tukey")))
>>
>>HTH, Mark.
>>
>>
>>Ullrich Ecker wrote:
>> >
>> > Hi everyone,
>> >
>> > I am fairly new to R, and I am aware that others have had this
>> > problem before, but I have failed to solve the problem from previous
>> > replies I found in the archives.
>> >
>> > As this is such a standard procedure in psychological science, there
>> > must be an elegant solution to this...I think.
>> >
>> > I would much appreciate a solution that even I could understand... ;-)
>> >
>> >
>> > Now, I want to calculate a post-hoc test following up a within-subjects
>> > ANOVA.
>> >
>> > The dv is reaction time (RT), there is a 3-level Condition factor
>> > (Cond; within-subjects), a number of subjects (Subj), and the
>> > dataframe is called WMU3C.
>> >
>> > The model is
>> >
>> >  > RT.aov <- aov(RT~Cond + Error(Subj/Cond), WMU3C)
>> >
>> > I understand that TukeyHSD only works with an aov object, but that
>> > RT.aov is an aovlist object.
>> >
>> >  > class(RT.aov)
>> > [1] "aovlist" "listof"
>> >
>> > I've tried to work around it using the "maiz" example in the MMC
>> > documentation of the HH package (a solution previously recommended),
>> > but I couldn't get it to work: My best shot was to calculate another
>> > aov avoiding the error term (I don't see how this could be a feasible
>> > approach, but that's how I understood the MMC example) and a contrast
>> > vector (contrasting conditions 2 and 3):
>> >
>> > I have to admit that I don't quite understand what I'm doing here
>> > (not that you couldn't tell)
>> >
>> >  > RT2.aov <- aov(terms(RT~Subj*Cond, WMU3C))
>> >  > Cond.lmat <- c(0,1,-1)
>> >  > Tukey <- glht.mmc(RT2.aov, focus = "Cond", focus.lmat = Cond.lmat)
>> >
>> > yielding
>> >
>> > Error in mvt(lower = carg$lower, upper = carg$upper, df = df, corr =
>> > carg$corr,  :
>> >    NA/NaN/Inf in foreign function call (arg 6)
>> > In addition: Warning message:
>> > In cov2cor(covm) : diagonal has non-finite entries
>> >
>> >  > Tukey
>> >        height
>> >
>> >
>> >
>> > Thank you very much for your help!
>> >
>> > Ullrich
>> >
>> >
>> > Dr Ullrich Ecker
>> > Postdoctoral Research Associate
>> > Cognitive Science Laboratories
>> > School of Psychology (Mailbag M304)
>> > Room 211 Sanders Building
>> > University of Western Australia
>> > 35 Stirling Hwy
>> > Crawley WA 6009
>> > Australia
>> > Office: +61 8 6488 3266
>> > Mobile: +61 4 5822 0072
>> > Fax: +61 8 6488 1006
>> > E-mail: ullrich.ecker at uwa.edu.au
>> > Web: www.cogsciwa.com
>> >       [[alternative HTML version deleted]]
>> >
>> > ______________________________________________
>> > R-help at r-project.org mailing list
>> > https://stat.ethz.ch/mailman/listinfo/r-help
>> > PLEASE do read the posting guide
>> > http://www.R-project.org/posting-guide.html
>> > and provide commented, minimal, self-contained, reproducible code.
>> >
>> >
>>
>>--
>>View this message in context: 
>>http://www.nabble.com/Tukey-HSD-%28or-other-post-hoc-tests%29-following-repeated-measures-ANOVA-tp17508294p17553029.html
>>Sent from the R help mailing list archive at Nabble.com.
>>
>>______________________________________________
>>R-help at r-project.org mailing list
>>https://stat.ethz.ch/mailman/listinfo/r-help
>>PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
>>and provide commented, minimal, self-contained, reproducible code.
> 
> 	[[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 

-- 
View this message in context: http://www.nabble.com/Tukey-HSD-%28or-other-post-hoc-tests%29-following-repeated-measures-ANOVA-tp17508294p17559307.html
Sent from the R help mailing list archive at Nabble.com.

More information about the R-help mailing list