[BioC] What's the difference of the two approach for two groups in limma?

Wed Apr 6 23:09:55 CEST 2011

Hi Jim,

Yes, they are exactly same now. Thanks a lots.

Yours sincerely,

Jianhong Ou

jianhong.ou at umassmed.edu

On Apr 6, 2011, at 4:57 PM, James W. MacDonald wrote:

> Hi Jianhong Ou,
> 
> On 4/6/2011 3:07 PM, Ou, Jianhong wrote:
>> Hi all,
>> 
>> In the usersguide.pdf of limma, there are two approach for two groups (8.4, Two Groups: Common Reference)
>> 
>> At first, I thought they are same, but...
>> 
>> [code]
>>> sd<- 0.3*sqrt(4/rchisq(100,df=4))
>>> y<- matrix(rnorm(100*6,sd=sd),100,6)
>>> rownames(y)<- paste("Gene",1:100)
>>> y[1:2,4:6]<- y[1:2,4:6] + 2
>>> design<- cbind(Grp1=1,Grp2vs1=c(0,0,0,1,1,1))
>>> fit<- lmFit(y,design)
>>> fit<- eBayes(fit)
>>> topTable(fit,coef=1,n=20,sort.by="logFC")
> 
> You made a mistake here. You want coef = 2, which is the contrast 
> between groups 1 and 2. What your code is doing is testing to see if the 
> average of the first group is equal to zero or not.
> 
> Best,
> 
> Jim
> 
>>         ID      logFC          t     P.Value adj.P.Val          B
>> 33 Gene 33  1.2800739  4.9383702 0.001156486 0.1156486 -0.4485569
>> 23 Gene 23 -0.9290832 -1.6732940 0.133021448 0.9913757 -4.9814041
>> 27 Gene 27  0.8616516  3.3965657 0.009487765 0.4743882 -2.4838934
>> 76 Gene 76  0.7509497  1.7891908 0.111590341 0.9913757 -4.8250883
>> 54 Gene 54 -0.5448046 -1.5837617 0.152121135 0.9913757 -5.0988743
>> 16 Gene 16 -0.4970755 -1.7635064 0.116040512 0.9913757 -4.8601096
>> 94 Gene 94 -0.4921761 -2.2254793 0.056881214 0.9913757 -4.2056166
>> 93 Gene 93 -0.4574475 -1.3155906 0.224961686 0.9913757 -5.4296133
>> 87 Gene 87 -0.4030807 -2.1728549 0.061734408 0.9913757 -4.2822231
>> 92 Gene 92 -0.3708210 -1.3967269 0.200232145 0.9913757 -5.3332685
>> 38 Gene 38 -0.3693720 -1.6087183 0.146559311 0.9913757 -5.0664427
>> 41 Gene 41  0.3296712  0.3589144 0.728998939 0.9913757 -6.2080737
>> 51 Gene 51 -0.3236529 -1.0426366 0.327764650 0.9913757 -5.7250738
>> 40 Gene 40  0.3003745  1.3450728 0.215687305 0.9913757 -5.3950122
>> 83 Gene 83  0.2948265  1.5802816 0.152911783 0.9913757 -5.1033766
>> 47 Gene 47  0.2914631  1.4141449 0.195242647 0.9913757 -5.3121375
>> 39 Gene 39 -0.2828721 -1.5449138 0.161161348 0.9913757 -5.1488494
>> 69 Gene 69  0.2757269  1.2502012 0.246750193 0.9913757 -5.5046099
>> 19 Gene 19  0.2684281  1.2757686 0.238027563 0.9913757 -5.4755794
>> 3   Gene 3 -0.2642721 -1.6381776 0.140233597 0.9913757 -5.0278441
>>> design<- cbind(Grp1=c(1,1,1,0,0,0),Grp2=c(0,0,0,1,1,1))
>>> cont.matrix<- makeContrasts(G1vsG2=Grp1-Grp2, levels=design)
>>> fit1<- lmFit(y,design)
>>> fit1<- contrasts.fit(fit1, cont.matrix)
>>> fit1<- eBayes(fit1)
>>> topTable(fit1,coef=1,n=20,sort.by="logFC")
>>         ID      logFC          t      P.Value   adj.P.Val          B
>> 2   Gene 2 -2.0639600 -7.2886732 8.728373e-05 0.008728373  1.9349063
>> 1   Gene 1 -1.9752877 -5.0131716 1.053468e-03 0.035115593 -0.6626835
>> 33 Gene 33  1.9059513  5.1993066 8.380943e-04 0.035115593 -0.4227207
>> 23 Gene 23 -1.6921631 -2.1549881 6.347344e-02 0.793418043 -4.8748108
>> 93 Gene 93 -1.5367983 -3.1252258 1.421720e-02 0.313735106 -3.3751312
>> 27 Gene 27  1.0978401  3.0600772 1.568676e-02 0.313735106 -3.4759186
>> 54 Gene 54 -0.7961582 -1.6365658 1.405731e-01 0.949255469 -5.6234924
>> 37 Gene 37 -0.6926892 -2.2573983 5.412387e-02 0.773198083 -4.7195118
>> 42 Gene 42 -0.6895975 -0.9536791 3.683231e-01 0.952227637 -6.4283974
>> 76 Gene 76  0.6887321  1.1603291 2.795521e-01 0.952227637 -6.2154077
>> 87 Gene 87 -0.6277337 -2.3927605 4.383747e-02 0.730624439 -4.5120531
>> 51 Gene 51  0.6097573  1.3889781 2.024874e-01 0.952227637 -5.9464611
>> 41 Gene 41 -0.5810686 -0.4473240 6.665707e-01 0.952227637 -6.7979655
>> 15 Gene 15  0.5801261  1.8016224 1.094946e-01 0.898763677 -5.3938078
>> 7   Gene 7 -0.5487280 -0.4418824 6.703406e-01 0.952227637 -6.8005963
>> 80 Gene 80 -0.5459118 -1.5641098 1.566352e-01 0.949255469 -5.7209503
>> 39 Gene 39 -0.5139565 -1.9848398 8.263317e-02 0.875815090 -5.1285443
>> 16 Gene 16 -0.5085620 -1.2758028 2.380161e-01 0.952227637 -6.0835464
>> 95 Gene 95  0.5040813  0.9693600 3.609117e-01 0.952227637 -6.4133578
>> 69 Gene 69  0.4815805  1.5440263 1.613734e-01 0.949255469 -5.7475589
>> [/code]
>> 
>> I want to know why the output are different. Which one should I believe?
>> 
>> Thanks a lot.
>> 
>> Yours sincerely,
>> 
>> Jianhong Ou
>> 
>> jianhong.ou at umassmed.edu
>> 
>> _______________________________________________
>> Bioconductor mailing list
>> Bioconductor at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/bioconductor
>> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor
> 
> -- 
> James W. MacDonald, M.S.
> Biostatistician
> Douglas Lab
> University of Michigan
> Department of Human Genetics
> 5912 Buhl
> 1241 E. Catherine St.
> Ann Arbor MI 48109-5618
> 734-615-7826
> **********************************************************
> Electronic Mail is not secure, may not be read every day, and should not be used for urgent or sensitive issues 
> 
> 