[R] Multiple tests on 2 way-ANOVA

Spencer Graves spencer.graves at pdf.com
Sat Jul 15 17:32:42 CEST 2006

Grathwohl, Dominik, LAUSANNE, NRC-BAS wrote:
> Dear r-helpers,
>
> I have a question about multiple testing.
> Here an example that puzzles me:
> All matrixes and contrast vectors are presented in treatment contrasts.
>
> 1. example:
> library(multcomp)
> n<-60; sigma<-20
> # n = sample size per group
> # sigma standard deviation of the residuals
>
> cov1 <- matrix(c(3/4,-1/2,-1/2,-1/2,1,0,-1/2,0,1), nrow = 3, ncol=3, byrow=TRUE,
> 	dimnames = list(c("A", "B", "C"), c("C.1", "C.2", "C.3")))
> # cov1 = variance covariance matrix of the beta coefficients of a
> # 2x2 factorial design (see Piantadosi 2005, p. 509)
>
> cm1 <- matrix(c(0, 1, 0, 0, 0, 1), nrow = 2, ncol=3, byrow=TRUE,
> 	dimnames = list(c("A", "B"), c("C.1", "C.2", "C.3")))
> # cm1 = contrast matrix for main effects
>
> v1 <- csimint(estpar=c(100, 6, 5), df=4*n-3, covm=cov1*sigma^2/n, cmatrix=cm1, conf.level=0.95)
> summary(v1)
>
> The adjusted p-values are almost the Bonferroni p-values.
> If I understood right: You need not to adjust for multiple testing
> on main effects in a 2x2 factorial design
> assuming the absence of interaction.

SG:  Where did you get this idea?  A p value of 0.05 says that if the
null hypothesis of no effect is true, a result at least as extreme as
that observed work actually occur with probability 0.05.  Thus, with 2
independent tests, the probability of getting a result at least that
extreme in one or both of the tests is 1-(1-0.05)^2 = 0.0975, which is
almost 2*0.05.  Thus, if I were to consider only main effects in a 2x2
factorial design, this is what I would get from Bonferroni.

> I do not think that there is a bug,
> I want to understand, why multcomp does adjust for multiple tests
> having all information about the design of the trial (variance covariance matrix)?
>
> 2. example:
> And I have second question: How do I proper correct for multiple testing
> if I want to estimate in the presence of interaction the two average main effects.
> Can some one point me to some literature where I can learn these things?
> Here the example, 2x2 factorial with interaction, estimation of average main effects:
>
> cov2 <- matrix(
> c(1,-1,-1, 1,
>  -1, 2, 1,-2,
>  -1, 1, 2,-2,
>   1,-2,-2, 4)
> , nrow=4, ncol=4, byrow=TRUE)
> cm2 <- matrix(c(0, 1, 0, 1/2, 0, 0, 1, 1/2), nrow = 2, ncol=4, byrow=TRUE,
> 	dimnames = list(c("A", "B"), c("C.1", "C.2", "C.3", "C.4")))
> v2 <- csimint(estpar=c(100, 6, 5, 2), df=4*n-4, covm=cov2*sigma^2/n, cmatrix=cm2, conf.level=0.95)
> summary(v2)

SG:  The Bonferroni p-value is the observed times the number of rows in
the contrast matrix.  The number of columns is irrelevant to Bonferroni.

SG:  I'm not sure, but I believe that the the adjusted p value would
likely be close to (if not exactly) the rank of the contrast matrix;
given the rank, it is (I think) independent of the number of rows and
columns.

SG:  These two assertions are consistent with the following example,
where I increase the number of dimensions by a factor of 4 without
changing the rank.  The Bonferroni p value increased by a factor of 4
while the adjusted p value did not change, as predicted.

cm2.4 <- rbind(cm2, cm2, cm2, cm2)
v2.4 <- csimint(estpar=c(100, 6, 5, 2), df=4*n-4,
covm=cov2*sigma^2/n,
cmatrix=cm2.4, conf.level=0.95)
summary(v2.4)
>
> I do not believe that this is the most efficient way for doing this,

SG:  I hope this reply converts the "bad experience" to "good".  As for
efficiency, you did very well by including such simple but elegant
examples.  Your post might have more efficiently elicited more and more
elegant responses sooner with a more carefully chosen Subject, perhaps
like "Multiple Comparisons Questions".  However, the selection of a
possible better subject might rely on information you didn't have.

Hope this helps.
Spencer Graves
>
> My R.version:
>
> platform i386-pc-mingw32
> arch     i386
> os       mingw32
> system   i386, mingw32
> status
> major    2
> minor    2.1
> year     2005
> month    12
> day      20
> svn rev  36812
> language R
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help