[R] Fisher's LSD multiple comparisons in a two-way ANOVA

Jinsong Zhao jszhao at yeah.net
Fri Apr 6 11:01:13 CEST 2012


On 2012-04-05 10:49, Richard M. Heiberger wrote:
> Here is your example.  The table you displayed in gigawiz ignored the
> two-way factor structure
> and interpreted the data as a single factor with 6 levels.  I created
> the interaction of
> a and b to get that behavior.
> ## your example, with data stored in a data.frame
> tmp <- data.frame(x=c(76, 84, 78, 80, 82, 70, 62, 72,
>                      71, 69, 72, 74, 66, 74, 68, 66,
>                      69, 72, 72, 78, 74, 71, 73, 67,
>                      86, 67, 72, 85, 87, 74, 83, 86,
>                      66, 68, 70, 76, 78, 76, 69, 74,
>                      72, 72, 76, 69, 69, 82, 79, 81),
>                    a=factor(rep(c("A1", "A2"), each = 24)),
>                    b=factor(rep(c("B1", "B2", "B3"), each=8, times=2)))
> x.aov <- aov(x ~ a*b, data=tmp)
> summary(x.aov)
> ## your request
> require(multcomp)
> tmp$ab <- with(tmp, interaction(a, b))
> xi.aov <- aov(x ~ ab, data=tmp)
> summary(xi.aov)
> xi.glht <- glht(xi.aov, linfct=mcp(ab="Tukey"))
> confint(xi.glht)
>
> ## graphs
> ## boxplots
> require(lattice)
> bwplot(x ~ ab, data=tmp)
> ## interaction plot
> ## install.packages("HH")  ## if you don't have HH yet
> require(HH)
> interaction2wt(x ~ a*b, data=tmp)
>

Thank you very much for the demonstration.

There is still a small difference between the results of glht() and the 
the table displayed in gigawiz. I try my best to figure out, but fail...

By the way, I have a similar question. I built a ANOVA model:

activity ~ pH * I * f

             Df Sum Sq Mean Sq F value   Pr(>F)
pH           1   1330    1330  59.752 2.15e-10 ***
I            1    137     137   6.131   0.0163 *
f            6  23054    3842 172.585  < 2e-16 ***
pH:I         1    152     152   6.809   0.0116 *
pH:f         6    274      46   2.049   0.0741 .
I:f          6   5015     836  37.544  < 2e-16 ***
pH:I:f       6    849     142   6.356 3.82e-05 ***
Residuals   56   1247      22
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.

Now, how can I do a multi-comparison on `pH:I'?

Do I need to do separate ANOVA for each `pH' or `I', just as that in 
demo("MMC.WoodEnergy", "HH")? And then do multi-comparison on `I' or 
`pH' in each separate ANOVA?

Thanks again.

Regards,
Jinsong



More information about the R-help mailing list