[R] Multiple tests on 2 way-ANOVA
Spencer Graves
spencer.graves at pdf.com
Sat Jul 15 17:32:42 CEST 2006
<comments in line>
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)?
> Or do I have to introduce somehow more information?
>
> 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,
> since I made already bad experience with the first example.
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
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
More information about the R-help
mailing list