[BioC] Question on limma design matrix

Sim, Fraser Fraser_Sim at URMC.Rochester.edu
Tue Jan 13 14:49:27 CET 2009


Hi-

I have an Affy array experiment that I am analyzing that contains
matched samples from separate patients. For each patient, there are 4
arrays which comprise a 2x2 factorial design (two tissues and two
treatments). 

I wanted to incorporate both a paired/matched type approach to the limma
design matrix and also the factorial nature or the experiment. I looked
at the vignettes and came up with the following.

Note: SID refers to the sample ID

Here's my code:
	design <- model.matrix(~SID+(treatment*tissue))
	design
  (Intercept) SID1 SID2 SID3 SID4 TreatmentA tissueB TreatmentA:tissueB
1            1      0      1      0      0       0        1
0
2            1      0      1      0      0       1        1
1
3            1      0      1      0      0       0        0
0
4            1      0      1      0      0       1        0
0
5            1      0      0      0      0       0        1
0
6            1      0      0      0      0       1        1
1
7            1      0      0      0      0       0        0
0
8            1      0      0      0      0       1        0
0
9            1      1      0      0      0       0        1
0
10           1      1      0      0      0       1        1
1
11           1      1      0      0      0       0        0
0
12           1      1      0      0      0       1        0
0
13           1      0      0      1      0       0        1
0
14           1      0      0      1      0       1        1
1
15           1      0      0      1      0       0        0
0
16           1      0      0      1      0       1        0
0
17           1      0      0      0      1       1        1
1
18           1      0      0      0      1       0        1
0
19           1      0      0      0      1       1        0
0
20           1      0      0      0      1       0        0
0	

	arrayw <- arrayWeights(eset, design = design)
	fit <- lmFit(eset, design, weights=arrayw)

	contrast.matrix =
cbind(ContrastA=c(0,0,0,0,0,1,0,0),ContrastB=c(0,0,0,0,0,1,0,1),
	 ConstrastDiff=c(0,0,0,0,0,0,0,1))

	fit2 <- contrasts.fit(fit, contrast.matrix)
	fit3 <- eBayes(fit2)

Does this look correct?

ContrastA extracts the treatment difference for one tissue
ContrastB extracts the treatment difference for the other tissue
ContrastDiff extracts the interaction between the two and therefore the
genes that differ with treatment between the two tissues

Thanks,
Fraser



More information about the Bioconductor mailing list