[BioC] Loop design - biological, technical replication and contrasts
Maciej Jończyk
mjonczyk at biol.uw.edu.pl
Tue Mar 2 19:35:21 CET 2010
I also think about lower p-value cutoff (0,01).
Below I pasted information about my experimental design (it is from my
previus message).
I have 2x2 factorial experiment conducted in a loop-design (two-colour
data).
There are two maize lines (dh7, dl3) and two temperature treatments
(cold=c and control=k).
Following hybridizations were conducted : dh_c vs dh_k; dh_k vs dl_k;
dl_k vs dl_c; dl_c vs dh_c
(forming a square on a diagram)
and two (diagonal on a diagram), namely: dl_k vs dh_c and dh_k vs dl_c.
I have four biological replications of this design, including dye-swap
(i.e. two hyb. cy3-cy5 and two cy5-cy3).
In each biological replication I have also three technical replications
of each RNA source
(i.e. in dh_c vs dh_k; dl_c vs dh_c and dl_k vs dh_c, sample dh_c is
from the same RNA pool).
I don't know how analyse it in limma. I want to get following contrasts:
dh_c vs dh_k; dh_k vs dl_k; dl_k vs dl_c; dl_c vs dh_c and include
effect of temperature,
maize line and interaction.
Avehna <avhena at gmail.com> nadawca :
> > Hi Maciej:I have the same problem, I did the same procedure but
> > still I'm getting large numbers for differentially expressed
> > genes. I could reduce this number by defining p.value = 0.01 in
> > decideTests. But I'm not completely sure whether changing the
> > "method" for decideTests and/or pvalue should give better results.
> > I'm looking forward to someone else answer.Best
> > Regards,Avhena2010/3/2 Maciej Jończyk <mjonczyk at biol.uw.edu.pl> Hi
> > again, I apologise for replying to my own post, but it helps keep
> > track if someone will be interested. I analysed my data with single
> > channel analysis in limma, according to Chapter 9. of limma
> > usersguide. I changed my targets file (to make it more condensed)
> > and removed suffix which identified biological replication. So my
> > targets looks like: >nt_trg SlideNumber
> > FileName Cy3 Cy5 1 93
> > c_093_DH_K_vs_DH_CHex.gpr hk hc 2 104
> > c_104_DH_CH_vs_DH_Kex.gpr hc hk 3 116
> > c_116_DHK_vs_DHCHex.gpr hk hc 4 16
> > c_016_DH_C_vs_DH_Kex.gpr hc hk 5 94
> > c_094_DH_K_vs_DL_Kex.gpr hk lk 6 105
> > c_105_DL_K_vs_DH_Kex.gpr lk hk 7 117
> > c_117_DHK_vs_DLKex.gpr hk lk 8 139
> > c_139_DL_K_vs_DH_Kex.gpr lk hk 9 92
> > c_092_DL_CH_vs_DL_Kex.gpr lc lk 10 106
> > c_106_DL_K_vs_DL_CHex.gpr lk lc 11 118
> > c_118_DLCH_vs_DLKex.gpr lc lk 12 23
> > c_023_DL_K_vs_DL_Cex.gpr lk lc 13 95
> > c_095_DL_CH_vs_DH_CHex.gpr lc hc 14 107
> > c_107_DH_CH_vs_DL_CHex.gpr hc lc 15 119
> > c_119_DLCH_vs_DHCHex.gpr lc hc 16 136
> > c_136_DH_C_vs_DL_Cex.gpr hc lc 17 101
> > c_101_DL_K_vs_DH_CHex.gpr lk hc 18 103
> > c_103_DH_CH_vs_DL_Kex.gpr hc lk 19 121
> > c_121_DLK_vs_DHCHex.gpr lk hc 20 15
> > c_015_DH_C_vs_DL_Kex.gpr hc lk 21 100
> > c_100_DH_K_vs_DL_CHex.gpr hk lc 22 102
> > c_102_DL_CH_vs_DH_Kex.gpr lc hk 23 120
> > c_120_DHK_vs_DLCHex.gpr hk lc 24 140
> > c_140_DL_C_vs_DH_Kex.gpr lc hk I transform it to apropriate
> > form: >tgr_sc=targetsA2C(nt_trg) >tgr_sc channel.col
> > SlideNumber FileName Target 1.1
> > 1 93 c_093_DH_K_vs_DH_CHex.gpr hk 1.2
> > 2 93 c_093_DH_K_vs_DH_CHex.gpr
> > hc 2.1 1 104
> > c_104_DH_CH_vs_DH_Kex.gpr hc 2.2 2
> > 104 c_104_DH_CH_vs_DH_Kex.gpr hk 3.1 1
> > 116 c_116_DHK_vs_DHCHex.gpr hk 3.2
> > 2 116 c_116_DHK_vs_DHCHex.gpr hc 4.1
> > 1 16 c_016_DH_C_vs_DH_Kex.gpr
> > hc 4.2 2 16
> > c_016_DH_C_vs_DH_Kex.gpr hk 5.1 1
> > 94 c_094_DH_K_vs_DL_Kex.gpr hk 5.2 2
> > 94 c_094_DH_K_vs_DL_Kex.gpr lk 6.1
> > 1 105 c_105_DL_K_vs_DH_Kex.gpr lk 6.2
> > 2 105 c_105_DL_K_vs_DH_Kex.gpr hk
> > 7.1 1 117 c_117_DHK_vs_DLKex.gpr
> > hk 7.2 2 117
> > c_117_DHK_vs_DLKex.gpr lk 8.1 1
> > 139 c_139_DL_K_vs_DH_Kex.gpr lk 8.2 2
> > 139 c_139_DL_K_vs_DH_Kex.gpr hk 9.1
> > 1 92 c_092_DL_CH_vs_DL_Kex.gpr lc 9.2
> > 2 92 c_092_DL_CH_vs_DL_Kex.gpr lk
> > 10.1 1 106 c_106_DL_K_vs_DL_CHex.gpr
> > lk 10.2 2 106
> > c_106_DL_K_vs_DL_CHex.gpr lc 11.1 1
> > 118 c_118_DLCH_vs_DLKex.gpr lc 11.2 2
> > 118 c_118_DLCH_vs_DLKex.gpr lk 12.1
> > 1 23 c_023_DL_K_vs_DL_Cex.gpr lk 12.2
> > 2 23 c_023_DL_K_vs_DL_Cex.gpr
> > lc 13.1 1 95 c_095_DL_CH_vs_DH_CHex.gpr
> > lc 13.2 2 95
> > c_095_DL_CH_vs_DH_CHex.gpr hc 14.1 1
> > 107 c_107_DH_CH_vs_DL_CHex.gpr hc 14.2 2
> > 107 c_107_DH_CH_vs_DL_CHex.gpr lc 15.1
> > 1 119 c_119_DLCH_vs_DHCHex.gpr lc 15.2
> > 2 119 c_119_DLCH_vs_DHCHex.gpr hc 16.1
> > 1 136 c_136_DH_C_vs_DL_Cex.gpr hc
> > 16.2 2 136 c_136_DH_C_vs_DL_Cex.gpr
> > lc 17.1 1 101
> > c_101_DL_K_vs_DH_CHex.gpr lk 17.2 2
> > 101 c_101_DL_K_vs_DH_CHex.gpr hc 18.1 1
> > 103 c_103_DH_CH_vs_DL_Kex.gpr hc 18.2
> > 2 103 c_103_DH_CH_vs_DL_Kex.gpr lk 19.1
> > 1 121 c_121_DLK_vs_DHCHex.gpr lk
> > 19.2 2 121 c_121_DLK_vs_DHCHex.gpr
> > hc 20.1 1 15
> > c_015_DH_C_vs_DL_Kex.gpr hc 20.2 2
> > 15 c_015_DH_C_vs_DL_Kex.gpr lk 21.1 1
> > 100 c_100_DH_K_vs_DL_CHex.gpr hk 21.2
> > 2 100 c_100_DH_K_vs_DL_CHex.gpr lc 22.1
> > 1 102 c_102_DL_CH_vs_DH_Kex.gpr lc 22.2
> > 2 102 c_102_DL_CH_vs_DH_Kex.gpr
> > hk 23.1 1 120
> > c_120_DHK_vs_DLCHex.gpr hk 23.2 2
> > 120 c_120_DHK_vs_DLCHex.gpr lc 24.1 1
> > 140 c_140_DL_C_vs_DH_Kex.gpr lc 24.2
> > 2 140 c_140_DL_C_vs_DH_Kex.gpr hk Next, I
> > made design matrix >u=unique(tgr_sc$Target)
> > >f=factor(tgr_sc$Target,levels=u) >design=model.matrix(~0+f)
> > >colnames(design)=u >design hk hc lk lc 1 1 0 0 0 2
> > 0 1 0 0 3 0 1 0 0 4 1 0 0 0 5 1 0 0
> > 0 6 0 1 0 0 7 0 1 0 0 8 1 0 0 0 9 1
> > 0 0 0 10 0 0 1 0 11 0 0 1 0 12 1 0 0 0
> > 13 1 0 0 0 14 0 0 1 0 15 0 0 1 0 16 1 0
> > 0 0 17 0 0 0 1 18 0 0 1 0 19 0 0 1 0 20
> > 0 0 0 1 21 0 0 0 1 22 0 0 1 0 23 0 0 1
> > 0 24 0 0 0 1 25 0 0 0 1 26 0 1 0 0 27 0
> > 1 0 0 28 0 0 0 1 29 0 0 0 1 30 0 1 0 0
> > 31 0 1 0 0 32 0 0 0 1 33 0 0 1 0 34 0 1
> > 0 0 35 0 1 0 0 36 0 0 1 0 37 0 0 1 0 38
> > 0 1 0 0 39 0 1 0 0 40 0 0 1 0 41 1 0 0
> > 0 42 0 0 0 1 43 0 0 0 1 44 1 0 0 0 45 1
> > 0 0 0 46 0 0 0 1 47 0 0 0 1 48 1 0 0 0
> > attr(,"assign") [1] 1 1 1 1 attr(,"contrasts") attr(,"contrasts")$f
> > [1] "contr.treatment" *Is it correct form my design? I see, that it
> > simply identifies what RNA was hybridized on each array.
> > >corfit=intraspotCorrelation(nt_img_lA,design) corfit$consensus [1]
> > 0.7341876
> > >fit=lmscFit(nt_img_lAq,design,correlation=corfit$consensus) I want
> > to get contrasts "hc - hk", "lc - lk", "hc - lc", "hk - lk" and also
> > test effect of line and temperature. To do that I write this
> > command:
> >
>contr.matrix=makeContrasts(hc-hk,lc-lk,hc-lc,hk-lk,linia=(hc+hk-lc-lk)/2,temp=(hc+lc-hk-lk)/2,inter=(hc-lc)-(hk-lk),levels=design)
> > * I'm not 100% sure that it's correct.
> > >contr.fit=contrasts.fit(fit,contr.matrix)
> > >contr.fit=eBayes(contr.fit)
> >
>wynik=decideTests(contr.fit,method="global",adjust.method="BH",p.value=0.05)
> > >summary(wynik) hc - hk lc - lk hc - lc hk - lk linia temp
> > inter -1 5865 5039 3014 2685 3931 7382
> > 1113 0 30922 33433 37177 38480 35896 28364 40776 1
> > 6594 4909 3190 2216 3554 7635 1492 From
> > that it seem that there is a lot of differentially expressed genes.
> > I feel that it isn't optimal design, here technical and
> > biological replications are treated in the same manner, aren't
> > they? I've read about "duplicateCorrelation" command, is it
> > possible to combine it with single channel analysis? Or I should
> > rewrite target file (add number of replication) and rewrite
> > contrasts (e.g. hc-hk change to
> > "((hc1+hc2+hc3+hc4)-(hk1+hk2+hk3+hk4))/4 )? And if I want to
> > include a dye effect, I should only add column with 1's to my
> > design, right? Thank you for reading of my post. I'd be very
> > grateful for help. I've tried to analyse this data for a along
> > time and I think limma is the best choice. Yours sincerely, Maciej
> > Jończyk Maciej Jończyk Department of Plant Molecular
> > Ecophysiology Institute of Plant Experimental Biology Faculty of
> > Biology, University of Warsaw 02-096 Warszawa, Miecznikowa 1
> > ___________________________________ NOCC,
> > http://nocc.sourceforge.net -- This email was Anti Virus checked
> > by Astaro Security Gateway. http://www.astaro.com
> > _______________________________________________ Bioconductor mailing
> > list Bioconductor at stat.math.ethz.ch
> > https://stat.ethz.ch/mailman/listinfo/bioconductor Search the
> > archives:
> > http://news.gmane.org/gmane.science.biology.informatics.conductor
Maciej Jończyk
Department of Plant Molecular Ecophysiology
Institute of Plant Experimental Biology
Faculty of Biology, University of Warsaw
02-096 Warszawa, Miecznikowa 1
___________________________________
NOCC, http://nocc.sourceforge.net
More information about the Bioconductor
mailing list