[BioC] Help expanding limmaUsersGuide section 8.7 to 3x2 factorial

[Jenny Drnevich]

(My previous reply got held up for moderator approval because somehow 
it was too big - too much extra text at the bottom. This is the 
shortened version)

Hi Gordon,

Thank you!! Now I can either use

 > solve(crossprod(design),t(design))

as the contrast matrix to turn the main effects & interaction into 
the group means, and then do a second contrast matrix to get the 
pairwise comparisons I want, or I could first do the cell means model, then use

 > t(solve(crossprod(design),t(design)))

as the contrast matrix to get the main effects & interactions. I 
never quite understood matrix algebra, but at least now I can go back 
and forth either way that's convenient. I assume this also extends to 
n-way ANOVAs?

 > Pop <- gl(2,4,labels=c("BMA","BVO"))
 > Trt <- gl(2,2,8,labels=c("c","e"))
 > Time <- gl(2,1,8,labels=c("0","1"))
 > contrasts(Pop) <- contr.sum(2)
 > contrasts(Trt) <- contr.sum(2)
 > contrasts(Time) <- contr.sum(2)
 > design <- model.matrix(~Pop*Trt*Time)
 > rownames(design) <- paste(Pop,Trt,Time,sep=".")
 > solve(crossprod(design),t(design))
                 BMA.c.0 BMA.c.1 BMA.e.0 BMA.e.1 BVO.c.0 BVO.c.1 
BVO.e.0 BVO.e.1
(Intercept)       0.125   0.125   0.125   0.125   0.125   0.125   0.125   0.125
Pop1              0.125   0.125   0.125   0.125  -0.125  -0.125  -0.125  -0.125
Trt1              0.125   0.125  -0.125  -0.125   0.125   0.125  -0.125  -0.125
Time1             0.125  -0.125   0.125  -0.125   0.125  -0.125   0.125  -0.125
Pop1:Trt1         0.125   0.125  -0.125  -0.125  -0.125  -0.125   0.125   0.125
Pop1:Time1        0.125  -0.125   0.125  -0.125  -0.125   0.125  -0.125   0.125
Trt1:Time1        0.125  -0.125  -0.125   0.125   0.125  -0.125  -0.125   0.125
Pop1:Trt1:Time1   0.125  -0.125  -0.125   0.125  -0.125   0.125   0.125  -0.125

Thanks again,
Jenny


At 02:35 AM 1/13/2011, Gordon K Smyth wrote:
>Hi Jenny,
>
>Here is a way to get the equivalent coefficients to put into your 
>contrast  matrix for a 4x2
>factorial:
>
> > Pop <- gl(4,2,labels=c("BMA","BVO","BWF","CAR"))
> > Trt <- gl(2,1,8,labels=c("c","e"))
> > contrasts(Pop) <- contr.sum(4)
> > contrasts(Trt) <- contr.sum(2)
> > design <- model.matrix(~Pop*Trt)
> > rownames(design) <- paste(Pop,Trt,sep=".")
> > solve(crossprod(design),t(design))
>               BMA.c  BMA.e  BVO.c  BVO.e  BWF.c  BWF.e  CAR.c  CAR.e
>(Intercept)  0.125  0.125  0.125  0.125  0.125  0.125  0.125  0.125 
>Pop1         0.375  0.375
>-0.125 -0.125 -0.125 -0.125 -0.125 -0.125 Pop2        -0.125 
>-0.125  0.375  0.375 -0.125 -0.125
>-0.125 -0.125 Pop3        -0.125 -0.125 -0.125 -0.125  0.375  0.375 
>-0.125 -0.125 Trt1
>0.125 -0.125  0.125 -0.125  0.125 -0.125  0.125 -0.125 
>Pop1:Trt1    0.375 -0.375 -0.125  0.125
>-0.125  0.125 -0.125  0.125 Pop2:Trt1   -0.125  0.125  0.375 -0.375 
>-0.125  0.125 -0.125  0.125
>Pop3:Trt1   -0.125  0.125 -0.125  0.125  0.375 -0.375 -0.125  0.125
>
>Best wishes
>Gordon
>
>PS. Sending again, because the email didn't seem to get to BioC the 
>first time.