[BioC] design matrix Limma design for paired t-test
Ingrid Mercier
Ingrid.Mercier at ipbs.fr
Wed Jun 13 14:46:20 CEST 2012
Thanks a lot Belinda !!
I mistaked so I replaced Time=Treat by Time only, and it's good.
So, I have a last question : I 'm confused with the differents coef in
topTable.
I get genes but I tested several coef without understanding their
significance.
Somebody can explain me what mean coef="TreatT", or coef= "Time18",coef=
" Donor5",coef= " Donor6", coef= "Donor7",coef= " Donor8".
My main objective is to identidy the differential expressed genes
between the Control donors and Treated Donors at 4 hours or 18 hours.
I have no idea, which coef I have to use it.
Cheers,
Ingrid
Ingrid MERCIER
Mycobacterial Interactions with Host Cells Team
Institute of Pharmacology& Structural Biology
CNRS - University of Toulouse
BP 64182
F-31077 Toulouse Cedex France
Tel +33 (0)5 61 17 54 63
Le 13/06/2012 08:45, Belinda Phipson a écrit :
> Hi Ingrid
>
> The problem with your code is the following line:
>> Time=Treat=factor(Targets$Time)
> Where you essentially set the time factor equal to the treat factor.
>
> Cheers,
> Belinda
>
>
> -----Original Message-----
> From: bioconductor-bounces at r-project.org
> [mailto:bioconductor-bounces at r-project.org] On Behalf Of Ingrid Mercier
> Sent: Wednesday, 13 June 2012 1:02 AM
> To: bioconductor at r-project.org; smyth at wehi.edu.au
> Subject: [BioC] design matrix Limma design for paired t-test
>
> Dear list and Gordon,
>
> I have some troubles to computed a moderated paired t-test in the linear
> model.
> Here is my experimental plan :
>
> I used a single channel Agilent microarray.
> 2 types of cells : Control (S) and Treated (T)
> Fives human donors : 4-5-6-7-8
> Two times of treatment : 4 hours and 18 hours
>
> I want to compare teh differential expresed genes between my C versus T at 4
> hours and then at 18 hours.
>
> Here is my design :
>
>
> My targets frame is :
>> Targets
> X FileName Treatment
> Donor Time
> 1 DC_4_4 US10463851_252665214446_S01_GE1_1010_Sep10_1_2.txt T
> 4 4
> 2 SC_4_4 US10463851_252665214448_S01_GE1_1010_Sep10_1_2.txt C
> 4 4
> 3 DC_18_4 US10463851_252665214447_S01_GE1_1010_Sep10_1_2.txt T
> 4 18
> 4 SC_18_4 US10463851_252665214444_S01_GE1_1010_Sep10_1_3.txt C
> 4 18
> 5 DC_4_5 US10463851_252665214448_S01_GE1_1010_Sep10_1_4.txt T
> 5 4
> 6 SC_4_5 US10463851_252665214444_S01_GE1_1010_Sep10_1_1.txt C
> 5 4
> 7 DC_18_5 US10463851_252665214446_S01_GE1_1010_Sep10_1_3.txt T
> 5 18
> 8 SC_18_5 US10463851_252665214447_S01_GE1_1010_Sep10_1_4.txt C
> 5 18
> 9 DC_4_6 US10463851_252665214445_S01_GE1_1010_Sep10_1_4.txt T
> 6 4
> 10 SC_4_6 US10463851_252665214447_S01_GE1_1010_Sep10_1_3.txt C
> 6 4
> 11 DC_18_6 US10463851_252665214448_S01_GE1_1010_Sep10_1_3.txt T
> 6 18
> 12 SC_18_6 US10463851_252665214445_S01_GE1_1010_Sep10_1_3.txt C
> 6 18
> 13 DC_4_7 US10463851_252665214444_S01_GE1_1010_Sep10_1_4.txt T
> 7 4
> 14 SC_4_7 US10463851_252665214445_S01_GE1_1010_Sep10_1_2.txt C
> 7 4
> 15 DC_18_7 US10463851_252665214447_S01_GE1_1010_Sep10_1_1.txt T
> 7 18
> 16 SC_18_7 US10463851_252665214446_S01_GE1_1010_Sep10_1_1.txt C
> 7 18
> 17 DC_4_8 US10463851_252665214444_S01_GE1_1010_Sep10_1_2.txt T
> 8 4
> 18 SC_4_8 US10463851_252665214446_S01_GE1_1010_Sep10_1_4.txt C
> 8 4
> 19 DC_18_8 US10463851_252665214445_S01_GE1_1010_Sep10_1_1.txt T
> 8 18
> 20 SC_18_8 US10463851_252665214448_S01_GE1_1010_Sep10_1_1.txt C
> 8 18
>
>
> then I create my design matrix :
>
>> Donor
> [1] 4 4 4 4 5 5 5 5 6 6 6 6 7 7 7 7 8 8 8 8
> Levels: 4 5 6 7 8
>> Treat=factor(Targets$Treatment,levels=c("C","T"))
>> Treat
> [1] T C T C T C T C T C T C T C T C T C T C
> Levels: C T
>> Time=Treat=factor(Targets$Time)
>> Time
> [1] 4 4 18 18 4 4 18 18 4 4 18 18 4 4 18 18 4 4 18 18
> Levels: 4 18
>
>> design=model.matrix(~Donor+Treat+Time)
>> design
> (Intercept) Donor5 Donor6 Donor7 Donor8 Treat18 Time18
> 1 1 0 0 0 0 0 0
> 2 1 0 0 0 0 0 0
> 3 1 0 0 0 0 1 1
> 4 1 0 0 0 0 1 1
> 5 1 1 0 0 0 0 0
> 6 1 1 0 0 0 0 0
> 7 1 1 0 0 0 1 1
> 8 1 1 0 0 0 1 1
> 9 1 0 1 0 0 0 0
> 10 1 0 1 0 0 0 0
> 11 1 0 1 0 0 1 1
> 12 1 0 1 0 0 1 1
> 13 1 0 0 1 0 0 0
> 14 1 0 0 1 0 0 0
> 15 1 0 0 1 0 1 1
> 16 1 0 0 1 0 1 1
> 17 1 0 0 0 1 0 0
> 18 1 0 0 0 1 0 0
> 19 1 0 0 0 1 1 1
> 20 1 0 0 0 1 1 1
> attr(,"assign")
> [1] 0 1 1 1 1 2 3
> attr(,"contrasts")
> attr(,"contrasts")$Donor
> [1] "contr.treatment"
>
> attr(,"contrasts")$Treat
> [1] "contr.treatment"
>
> attr(,"contrasts")$Time
> [1] "contr.treatment"
>
>
> In this design matrix I think something is wrong, because of the column
> Treat18 is the same as Time18.
> I don't understand why.
> So, the following code failed, and the differential expressed genes are
> odds.
>
> Somebody can help me !!! Thanks all.
>
>
>> fit=lmFit(test_norm,design)
> Coefficients not estimable: Time18
> Message d'avis :
> Partial NA coefficients for 34183 probe(s)
>> fit2=eBayes(fit)
> Message d'avis :
> In ebayes(fit = fit, proportion = proportion, stdev.coef.lim =
> stdev.coef.lim, :
> Estimation of var.prior failed - set to default value
>
>
>> table = topTable(fit2,1, number=5000,
> p.value=0.05,adjust.method="BH",sort.by="logFC",lfc=2)
>> head(table)
> ID logFC AveExpr t P.Value adj.P.Val
> B
> 6509 A_33_P3396434 18.44159 18.41239 245.14490 1.308161e-31 2.353520e-28
> 53.41519
> 22398 A_33_P3223592 18.25824 18.24591 242.75647 1.545005e-31 2.514901e-28
> 53.36821
> 10771 A_33_P3244165 18.21029 18.02229 90.76191 2.796577e-24 2.467615e-23
> 44.59915
> 6149 A_33_P3346552 18.14780 18.12098 207.18556 2.282464e-30 1.147374e-27
> 52.50960
> 23554 A_33_P3210160 18.08158 18.21026 239.64192 1.924175e-31 2.560908e-28
> 53.30521
> 20924 A_33_P3286278 18.04425 18.07312 179.72121 2.558128e-29 5.025546e-27
> 51.56876
>
>
> Best,
>
> Ingrid
>
>
More information about the Bioconductor
mailing list