[R-sig-ME] mcmcpvalue and contrasts

Hank Stevens HStevens at MUOhio.edu
Mon Feb 25 23:42:48 CET 2008


Hi Folks,
I wanted to double check that my intuition makes sense.

Examples of mcmcpvalue that I have seen use treatment "contrast" coding.
However, in more complex designs, testing overall effects of a factor  
might be better done with other contrasts, such as sum or Helmert  
contrasts.

My Contention:
Different contrasts test different hypothesis, and therefore result in  
different P-values. This consequence of contrasts differs from  
analysis of variance, as in anova( lm(Y ~ X1*X2) ).

*** This is right, isn't it? ***

I provide a self-contained example.

library(lme4) # v. 0.99875-9
library(coda)
library(SASmixed)

mcmcpvalue <- function(samp) {
    ##From Bates pers comm. September 14, 2006 2:59:23 PM EDT
    ## elementary version that creates an empirical p-value for the
    ## hypothesis that the columns of samp have mean zero versus a
    ## general multivariate distribution with elliptical contours.
# samp <- rnorm(10000, m=3)
    ## differences from the mean standardized by the observed
    ## variance-covariance factor
    std <- backsolve(chol(var(samp)),
                     cbind(0, t(samp)) - colMeans(samp),
                     transpose = TRUE)
    sqdist <- colSums(std * std)
    sum(sqdist[-1] > sqdist[1])/nrow(samp)
}


options(contrasts=c("contr.treatment","contr.poly"))
print(fm1 <- lmer(resistance ~ ET * position + (1|Grp), Semiconductor),
       corr=F)
c1 <- mcmcsamp(fm1, 50000)

options(contrasts=c("contr.helmert","contr.poly"))
print(fm2 <- lmer(resistance ~ ET * position + (1|Grp), Semiconductor),
       corr=F)
c2 <- mcmcsamp(fm2, 50000)

  mcmcpvalue(c1[, 2:4 ])
[1] 0.32458
 > mcmcpvalue(c2[, 2:4 ])
[1] 0.0249
 > mcmcpvalue(c1[, 5:7 ])
[1] 0.45376
 > mcmcpvalue(c2[, 5:7 ])
[1] 0.16302
 > mcmcpvalue(c1[, 8:16 ])
[1] 0.6373
 > mcmcpvalue(c2[, 8:16 ])
[1] 0.88808
 >
 > sessionInfo()
R version 2.6.2 (2008-02-08)
i386-apple-darwin8.10.1

locale:
C

attached base packages:
[1] stats     graphics  utils     datasets  grDevices methods
[7] base

other attached packages:
[1] SASmixed_0.4-2    xtable_1.5-2      reshape_0.8.0
[4] coda_0.13-1       lme4_0.99875-9    Matrix_0.999375-4
[7] lattice_0.17-4    CarbonEL_0.1-4

loaded via a namespace (and not attached):
[1] grid_2.6.2
 >

Dr. Hank Stevens, Assoicate Professor
338 Pearson Hall
Botany Department
Miami University
Oxford, OH 45056

Office: (513) 529-4206
Lab: (513) 529-4262
FAX: (513) 529-4243
http://www.cas.muohio.edu/~stevenmh/
http://www.cas.muohio.edu/ecology
http://www.muohio.edu/botany/

"If the stars should appear one night in a thousand years, how would men
believe and adore." -Ralph Waldo Emerson, writer and philosopher  
(1803-1882)




More information about the R-sig-mixed-models mailing list