[BioC] [limma] Strange results in contrasts with dye-swap
mjonczyk at biol.uw.edu.pl
mjonczyk at biol.uw.edu.pl
Thu Aug 4 11:27:51 CEST 2011
Dear List Members,
I've spent last day searching archive to find out how to
account for dye effect in dye swap common-reference
design.
I've found a lot of useful information but there are
few things I still don't understand.
I have time-course experiment with 7 time points, 2 conditions (c, k),
4 biological replications (with dye-swap).
Here is my targets frame (powt stands for biol. replication):
> trg
FileName Cy3 Cy5 powt
1 _024_1_vs_7ex.gpr k0 k2 1
2 _228_k0_vs_16ex.gpr k0 c12 1
3 _229_k0_vs_8ex.gpr k0 c4 1
4 _230_k0_vs_15ex.gpr k0 k10 1
5 _235_k0_vs_19ex.gpr k0 k14 1
6 _237_k0_vs_14ex.gpr k0 c10 1
7 _238_k0_vs_10ex.gpr k0 c6 1
8 _239_k0_vs_12ex.gpr k0 c8 1
9 _240_k0_vs_13ex.gpr k0 k8 1
10 _248_k0_vs_17ex.gpr k0 k12 1
11 _249_k0_vs_9ex.gpr k0 k4 1
12 _253_k0_vs_11ex.gpr k0 k6 1
13 _254_k0_vs_6ex.gpr k0 c2 1
14 _256_k0_vs_18ex.gpr k0 c14 1
15 _189_K8_vs_K0ex.gpr k8 k0 2
16 _190_K4_vs_K0ex.gpr k4 k0 2
17 _191_K6_vs_K0ex.gpr k6 k0 2
18 _192_C8_vs_K0ex.gpr c8 k0 2
19 _193_C10_vs_K0ex.gpr c10 k0 2
20 _194_C2_vs_K0ex.gpr c2 k0 2
21 _198_C14_vs_K0ex.gpr c14 k0 2
22 _200_K12_vs_K0ex.gpr k12 k0 2
23 _202_K10_vs_K0ex.gpr k10 k0 2
24 _206_C6_vs_K0ex.gpr c6 k0 2
25 _207_C12_vs_K0ex.gpr c12 k0 2
26 _208_K2_vs_K0ex.gpr k2 k0 2
27 _209_K14_vs_K0ex.gpr k14 k0 2
28 _210_C4_vs_K0ex.gpr c4 k0 2
29 _201_k0_vs_C2ex.gpr k0 c2 3
30 _203_k0_vs_K6ex.gpr k0 k6 3
31 _205_k0_vs_C12ex.gpr k0 c12 3
32 _217_k0_vs_K14ex.gpr k0 k14 3
33 _219_k0_vs_C4ex.gpr k0 c4 3
34 _231_k0_vs_C10ex.gpr k0 c10 3
35 _232_k0_vs_K8ex.gpr k0 k8 3
36 _233_k0_vs_C8ex.gpr k0 c8 3
37 _234_k0_vs_C14ex.gpr k0 c14 3
38 _235_k0_vs_C6ex.gpr k0 c6 3
39 _236_k0_vs_K10ex.gpr k0 k10 3
40 _237_k0_vs_K2ex.gpr k0 k2 3
41 _238_k0_vs_K4ex.gpr k0 k4 3
42 _239_k0_vs_K12ex.gpr k0 k12 3
43 _004_8_vs_1ex.gpr c8 k0 4
44 _005_10_vs_1ex.gpr c10 k0 4
45 _006_2_vs_1ex.gpr c2 k0 4
46 _007_15_vs_1ex.gpr k14 k0 4
47 _008_4_vs_1ex.gpr c4 k0 4
48 _010_3_vs_1ex.gpr k2 k0 4
49 _011_9_vs_1ex.gpr k8 k0 4
50 _012_14_vs_1ex.gpr c14 k0 4
51 _014_12_vs_1ex.gpr c12 k0 4
52 _026_13_vs_1ex.gpr k12 k0 4
53 _036_11_vs_1ex.gpr k10 k0 4
54 _037_6_vs_1ex.gpr c6 k0 4
55 _043_7_vs_1ex.gpr k6 k0 4
56 _044_5_vs_1ex.gpr k4 k0 4
I'm interested in contrasts between "c" and "k" in each separate time point.
I also want to account for dye-effect.
PROBLEM: I compared results from analysis with, without de-effect and also with
contrast
for dye-effect and it looks strange to me. Here are details of my analysis.
I've made design as specified in manual (here without dye-effect):
> design=modelMatrix(trg,ref="k0")
design
c10 c12 c14 c2 c4 c6 c8 k10 k12 k14 k2 k4 k6 k8
[1,] 0 0 0 0 0 0 0 0 0 0 1 0 0 0
[2,] 0 1 0 0 0 0 0 0 0 0 0 0 0 0
[3,] 0 0 0 0 1 0 0 0 0 0 0 0 0 0
[4,] 0 0 0 0 0 0 0 1 0 0 0 0 0 0
[5,] 0 0 0 0 0 0 0 0 0 1 0 0 0 0
[6,] 1 0 0 0 0 0 0 0 0 0 0 0 0 0
[7,] 0 0 0 0 0 1 0 0 0 0 0 0 0 0
[8,] 0 0 0 0 0 0 1 0 0 0 0 0 0 0
[9,] 0 0 0 0 0 0 0 0 0 0 0 0 0 1
[10,] 0 0 0 0 0 0 0 0 1 0 0 0 0 0
[11,] 0 0 0 0 0 0 0 0 0 0 0 1 0 0
[12,] 0 0 0 0 0 0 0 0 0 0 0 0 1 0
[13,] 0 0 0 1 0 0 0 0 0 0 0 0 0 0
[14,] 0 0 1 0 0 0 0 0 0 0 0 0 0 0
[15,] 0 0 0 0 0 0 0 0 0 0 0 0 0 -1
[16,] 0 0 0 0 0 0 0 0 0 0 0 -1 0 0
[17,] 0 0 0 0 0 0 0 0 0 0 0 0 -1 0
[18,] 0 0 0 0 0 0 -1 0 0 0 0 0 0 0
[19,] -1 0 0 0 0 0 0 0 0 0 0 0 0 0
[20,] 0 0 0 -1 0 0 0 0 0 0 0 0 0 0
[21,] 0 0 -1 0 0 0 0 0 0 0 0 0 0 0
[22,] 0 0 0 0 0 0 0 0 -1 0 0 0 0 0
[23,] 0 0 0 0 0 0 0 -1 0 0 0 0 0 0
[24,] 0 0 0 0 0 -1 0 0 0 0 0 0 0 0
[25,] 0 -1 0 0 0 0 0 0 0 0 0 0 0 0
[26,] 0 0 0 0 0 0 0 0 0 0 -1 0 0 0
[27,] 0 0 0 0 0 0 0 0 0 -1 0 0 0 0
[28,] 0 0 0 0 -1 0 0 0 0 0 0 0 0 0
[29,] 0 0 0 1 0 0 0 0 0 0 0 0 0 0
[30,] 0 0 0 0 0 0 0 0 0 0 0 0 1 0
[31,] 0 1 0 0 0 0 0 0 0 0 0 0 0 0
[32,] 0 0 0 0 0 0 0 0 0 1 0 0 0 0
[33,] 0 0 0 0 1 0 0 0 0 0 0 0 0 0
[34,] 1 0 0 0 0 0 0 0 0 0 0 0 0 0
[35,] 0 0 0 0 0 0 0 0 0 0 0 0 0 1
[36,] 0 0 0 0 0 0 1 0 0 0 0 0 0 0
[37,] 0 0 1 0 0 0 0 0 0 0 0 0 0 0
[38,] 0 0 0 0 0 1 0 0 0 0 0 0 0 0
[39,] 0 0 0 0 0 0 0 1 0 0 0 0 0 0
[40,] 0 0 0 0 0 0 0 0 0 0 1 0 0 0
[41,] 0 0 0 0 0 0 0 0 0 0 0 1 0 0
[42,] 0 0 0 0 0 0 0 0 1 0 0 0 0 0
[43,] 0 0 0 0 0 0 -1 0 0 0 0 0 0 0
[44,] -1 0 0 0 0 0 0 0 0 0 0 0 0 0
[45,] 0 0 0 -1 0 0 0 0 0 0 0 0 0 0
[46,] 0 0 0 0 0 0 0 0 0 -1 0 0 0 0
[47,] 0 0 0 0 -1 0 0 0 0 0 0 0 0 0
[48,] 0 0 0 0 0 0 0 0 0 0 -1 0 0 0
[49,] 0 0 0 0 0 0 0 0 0 0 0 0 0 -1
[50,] 0 0 -1 0 0 0 0 0 0 0 0 0 0 0
[51,] 0 -1 0 0 0 0 0 0 0 0 0 0 0 0
[52,] 0 0 0 0 0 0 0 0 -1 0 0 0 0 0
[53,] 0 0 0 0 0 0 0 -1 0 0 0 0 0 0
[54,] 0 0 0 0 0 -1 0 0 0 0 0 0 0 0
[55,] 0 0 0 0 0 0 0 0 0 0 0 0 -1 0
[56,] 0 0 0 0 0 0 0 0 0 0 0 -1 0 0
I've constructed contrast matrix
>
contrast.matrix=makeContrasts(c2-k2,c4-k4,c6-k6,c8-k8,c10-k10,c12-k12,c14-k14,levels=design)
> contrast.matrix
Contrasts
Levels c2 - k2 c4 - k4 c6 - k6 c8 - k8 c10 - k10 c12 - k12 c14 - k14
c10 0 0 0 0 1 0 0
c12 0 0 0 0 0 1 0
c14 0 0 0 0 0 0 1
c2 1 0 0 0 0 0 0
c4 0 1 0 0 0 0 0
c6 0 0 1 0 0 0 0
c8 0 0 0 1 0 0 0
k10 0 0 0 0 -1 0 0
k12 0 0 0 0 0 -1 0
k14 0 0 0 0 0 0 -1
k2 -1 0 0 0 0 0 0
k4 0 -1 0 0 0 0 0
k6 0 0 -1 0 0 0 0
k8 0 0 0 -1 0 0 0
> fit.userguide=lmFit(img_lA_ncav,design)
> fit2.ug=contrasts.fit(fit.userguide,contrast.matrix)
> fit2.ug=eBayes(fit2.ug)
> test.ug=decideTests(fit2.ug,method="global",adjust.method="fdr",p.value=0.05)
After that I have done the same analysis with dye-effect:
> design_dye=cbind(DyeEffect=1,design)
> design_dye
DyeEffect c10 c12 c14 c2 c4 c6 c8 k10 k12 k14 k2 k4 k6 k8
[1,] 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0
[2,] 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0
[3,] 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0
[4,] 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0
[5,] 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0
[6,] 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0
[7,] 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0
[8,] 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0
[9,] 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1
[10,] 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0
[11,] 1 0 0 0 0 0 0 0 0 0 0 0 1 0 0
[12,] 1 0 0 0 0 0 0 0 0 0 0 0 0 1 0
[13,] 1 0 0 0 1 0 0 0 0 0 0 0 0 0 0
[14,] 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0
[15,] 1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1
[16,] 1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0
[17,] 1 0 0 0 0 0 0 0 0 0 0 0 0 -1 0
[18,] 1 0 0 0 0 0 0 -1 0 0 0 0 0 0 0
[19,] 1 -1 0 0 0 0 0 0 0 0 0 0 0 0 0
[20,] 1 0 0 0 -1 0 0 0 0 0 0 0 0 0 0
[21,] 1 0 0 -1 0 0 0 0 0 0 0 0 0 0 0
[22,] 1 0 0 0 0 0 0 0 0 -1 0 0 0 0 0
[23,] 1 0 0 0 0 0 0 0 -1 0 0 0 0 0 0
[24,] 1 0 0 0 0 0 -1 0 0 0 0 0 0 0 0
[25,] 1 0 -1 0 0 0 0 0 0 0 0 0 0 0 0
[26,] 1 0 0 0 0 0 0 0 0 0 0 -1 0 0 0
[27,] 1 0 0 0 0 0 0 0 0 0 -1 0 0 0 0
[28,] 1 0 0 0 0 -1 0 0 0 0 0 0 0 0 0
[29,] 1 0 0 0 1 0 0 0 0 0 0 0 0 0 0
[30,] 1 0 0 0 0 0 0 0 0 0 0 0 0 1 0
[31,] 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0
[32,] 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0
[33,] 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0
[34,] 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0
[35,] 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1
[36,] 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0
[37,] 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0
[38,] 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0
[39,] 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0
[40,] 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0
[41,] 1 0 0 0 0 0 0 0 0 0 0 0 1 0 0
[42,] 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0
[43,] 1 0 0 0 0 0 0 -1 0 0 0 0 0 0 0
[44,] 1 -1 0 0 0 0 0 0 0 0 0 0 0 0 0
[45,] 1 0 0 0 -1 0 0 0 0 0 0 0 0 0 0
[46,] 1 0 0 0 0 0 0 0 0 0 -1 0 0 0 0
[47,] 1 0 0 0 0 -1 0 0 0 0 0 0 0 0 0
[48,] 1 0 0 0 0 0 0 0 0 0 0 -1 0 0 0
[49,] 1 0 0 0 0 0 0 0 0 0 0 0 0 0 -1
[50,] 1 0 0 -1 0 0 0 0 0 0 0 0 0 0 0
[51,] 1 0 -1 0 0 0 0 0 0 0 0 0 0 0 0
[52,] 1 0 0 0 0 0 0 0 0 -1 0 0 0 0 0
[53,] 1 0 0 0 0 0 0 0 -1 0 0 0 0 0 0
[54,] 1 0 0 0 0 0 -1 0 0 0 0 0 0 0 0
[55,] 1 0 0 0 0 0 0 0 0 0 0 0 0 -1 0
[56,] 1 0 0 0 0 0 0 0 0 0 0 0 -1 0 0
> fit.userguide.dye=lmFit(img_lA_ncav,design_dye)
>
contrast.matrix.d=makeContrasts(p2="c2-k2",p4="c4-k4",p6="c6-k6",p8="c8-k8",p10="c10-k10",p12="c12-k12",p14="c14-k14",levels=design_dye)
> fit2.ug.d=contrasts.fit(fit.userguide.dye,contrast.matrix.d)
> fit2.ug.d=eBayes(fit2.ug.d)
> test.ug.d=decideTests(fit2.ug.d,method="global",adjust.method="BH",p.value=0.05)
Then I included dye effect in contrasts:
>
contrast.matrix.d2=makeContrasts(DyeEffect,p2="c2-k2",p4="c4-k4",p6="c6-k6",p8="c8-k8",p10="c10-k10",p12="c12-k12",p14="c14-k14",levels=design_dye)
> contrast.matrix.d2
Contrasts
Levels DyeEffect p2 p4 p6 p8 p10 p12 p14
DyeEffect 1 0 0 0 0 0 0 0
c10 0 0 0 0 0 1 0 0
c12 0 0 0 0 0 0 1 0
c14 0 0 0 0 0 0 0 1
c2 0 1 0 0 0 0 0 0
c4 0 0 1 0 0 0 0 0
c6 0 0 0 1 0 0 0 0
c8 0 0 0 0 1 0 0 0
k10 0 0 0 0 0 -1 0 0
k12 0 0 0 0 0 0 -1 0
k14 0 0 0 0 0 0 0 -1
k2 0 -1 0 0 0 0 0 0
k4 0 0 -1 0 0 0 0 0
k6 0 0 0 -1 0 0 0 0
k8 0 0 0 0 -1 0 0 0
> fit2.ug.d2=contrasts.fit(fit.userguide.dye,contrast.matrix.d2)
> fit2.ug.d2=eBayes(fit2.ug.d2)
> test.ug.d2=decideTests(fit2.ug.d2,method="global",adjust.method="BH",p.value=0.05)
Results:
TEST WITHOUT DYE EFFECT
> summary(test.ug)
c2 - k2 c4 - k4 c6 - k6 c8 - k8 c10 - k10 c12 - k12 c14 - k14
-1 115 96 377 141 263 175 265
0 43082 43048 42326 42973 42694 42842 42752
1 196 249 690 279 436 376 376
TEST WITH DYE EFFECT
> summary(test.ug.d)
p2 p4 p6 p8 p10 p12 p14
-1 246 217 686 317 530 316 526
0 42797 42755 41594 42530 42103 42446 42239
1 350 421 1113 546 760 631 628
TEST WITH DYE EFFECT AND CONTRAST FOR DYE EFFECT
> summary(test.ug.d2)
DyeEffect p2 p4 p6 p8 p10 p12 p14
-1 6660 452 402 1103 593 858 539 883
0 30803 42353 42261 40659 41853 41359 41888 41564
1 5930 588 730 1631 947 1176 966 946
*QUESTION*
While I understand that second test (dye-effect included in model) can give more
significant genes (adjustment for dye-effect) than model without dye-effect, I
can't understand why the third test (with contrast for dye-effect) gives even
more significant genes. In the third analysis more contrasts have been performed
so there should be lower p-value cut-off and consequently one could suppose that
there will be less significant genes.
*OTHER QUESTIONS*
1. Is second model (WITH DYE EFFECT) correct?
2. So my anova model is y=Treatment x + DyeEffect, right?
3. How Array effect is taken into account in analysis?
4. Should I include biol-replication effect in analysis (as block)?
Thanks for any help,
Best Regards,
Maciej
________________________________________________________
Maciej Jończyk, MSc
Department of Plant Molecular Ecophysiology
Institute of Plant Experimental Biology
Faculty of Biology, University of Warsaw
02-096 Warszawa, Miecznikowa
1
More information about the Bioconductor
mailing list