Problem in limma-voom
1
0
Entering edit mode
Yunshun Chen ▴ 840
@yunshun-chen-5451
Last seen 26 days ago
Australia
Dear Ming, Regarding to your first question, a much easier way is to create your design matrix without the intercept. > design <- model.matrix(~0+RasTum) Then if you intend to pool both RasP class and RasP1Hit class together to compare to RasPNot class, you can use the following contrast: > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ RasP.Tumor)-RasPNot.Tumor, levels=design ) Similarly, for your second question, you can avoid confusion by creating your design matrix without the intercept. Suppose the six columns of your design matrix correspond to: Ras(Tum), Ras1H(Tum), RasNot(Tum), Ras(Nor), Ras1H(Nor), RasNot(Nor). Then you can make contrasts for any comparisons you are interested. For example: > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) Hope that helps. Yunshun Chen ------------------------------------------------------------------- Hi, Dr. Gordon and all: I run into some limma design issue not sure what they are: I had three tumor types want to compare with each other: RasP, RasP1Hit, RasPNot here are the critical commands: >targets<-targets[targets$RasType!="RasOnly" & targets$Type=="Tumor",] > show(head(targets)); Subject SampleName Type RasPType RasType 1 TCGA_38_4625 T4625_01A Tumor RasP RasP 5 TCGA_38_4627 T4627_01A Tumor RasP RasP 7 TCGA_38_4632 T4632_01A Tumor RasP RasP 9 TCGA_44_2655 T2655_01A Tumor RasP RasP 11 TCGA_44_2657 T2657_01A Tumor RasP RasP 13 TCGA_44_2661 T2661_01A Tumor RasP RasP > RasTum<-paste(targets$RasType, targets$Type,sep=".") > RasTum<-factor(RasTum, levels=c("RasP1Hit.Tumor","RasP.Tumor","RasPNot.Tumor")); > show(RasTum); [1] RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor [6] RasP.Tumor RasP.Tumor RasP.Tumor RasP1Hit.Tumor RasP1Hit.Tumor [11] RasP.Tumor RasPNot.Tumor RasP1Hit.Tumor RasP.Tumor RasP1Hit.Tumor [16] RasP1Hit.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor RasP1Hit.Tumor [21] RasP1Hit.Tumor RasP1Hit.Tumor RasPNot.Tumor RasP.Tumor RasPNot.Tumor [26] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor [31] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor [36] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasP.Tumor [41] RasP.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor RasPNot.Tumor Levels: RasP1Hit.Tumor RasP.Tumor RasPNot.Tumor > design<-model.matrix(~RasTum); > colnames(design)<-sub("RasTum","",colnames(design)); > colnames(design)[1] <- "Intercept" > rownames(design)<-targets$Sample > show(design[1:5,]) Intercept RasP.Tumor RasPNot.Tumor T4625_01A 1 1 0 T4627_01A 1 1 0 T4632_01A 1 1 0 T2655_01A 1 1 0 T2657_01A 1 1 0 ... > con.matrix<-makeContrasts( + RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot .Tumor=RasP.Tumor-RasPNot.Tumor, + RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor, + levels=design) > show(con.matrix); Contrasts Levels RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor Intercept 0 0 RasP.Tumor 1 1 RasPNot.Tumor -1 -1 Contrasts Levels RasP1Hit.Tumor_RasPNot.Tumor Intercept 0 RasP.Tumor 0 RasPNot.Tumor -1 My first question is: as shown in show(RasTum), there are three subtype of tumors I want to compare between: Levels: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor shown in show(design[1:5,]), everything is relative to RasP1Hit.Tumor. So in my con.matrix<-makeContrasts command: RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor in particular, RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor since everything is relative to RasP1Hit.Tumor. However, for contrast: RasPRasP1Hit.Tumor_RasPNot.Tumor, I intended to pool both RasP class and RasP1Hit class together to compare to RasPNot class, I end up with: RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor the right side is the same as RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor because everything is relative to RasP1Hit.Tumor. So what is the right way to set up contrast RasPRasP1Hit.Tumor_RasPNot.Tumor in makeContrasts command? kind of confused here? Because of the setting, the results are as below: > show(summary(de <- decideTests(fit))); RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor -1 0 0 0 17764 17764 1 0 0 RasP1Hit.Tumor_RasPNot.Tumor -1 0 0 17764 1 0 My 2nd question is: I did a similar design with including all normal samples, since at the same model I want to get tumor vs normal contrasts as well. > con.matrix<-makeContrasts( + RasP.Tumor_RasP.Normal = RasP.Tumor-RasP.Normal,RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal = RasP.Tumor-RasP.Normal-RasP1Hit.Normal, + RasPNot.Tumor_RasPNot.Normal=RasPNot.Tumor-RasPNot.Normal, + RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot .Tumor=RasP.Tumor-RasPNot.Tumor, + RasP.Normal_RasPNot.Normal=RasP.Normal- RasPNot.Normal,RasPRasP1Hit.Normal_Ra sPNot.Normal=RasP1Hit.Normal+RasP.Normal-RasPNot.Normal, + Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+-RasP.Normal- RasPNot.Normal, + Tumor_Normal_2=RasP.Tumor+ RasPNot.Tumor-RasP.Normal-RasPNot.Normal-RasP1Hit.Normal, + levels=design) here everything is relative to RasP1Hit.Tumor. and the result as below: > show(summary(de <- decideTests(fit))); RasP.Tumor_RasP.Normal RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal -1 3947 4522 0 9912 8830 1 4274 4781 RasPNot.Tumor_RasPNot.Normal RasP.Tumor_RasPNot.Tumor -1 4809 24 0 8461 18102 1 4863 7 RasPRasP1Hit.Tumor_RasPNot.Tumor RasP.Normal_RasPNot.Normal -1 24 608 0 18102 17158 1 7 367 RasPRasP1Hit.Normal_RasPNot.Normal Tumor_Normal Tumor_Normal_2 -1 2001 5611 5580 0 13934 6728 6651 1 2198 5794 5902 RasP.Tumor_RasPNot.Tumor contrast has small number of DEGs (24+7) at FDR 5% compared to using tumor data alone for the model shown earlier above, this might be cuased by filtering as well and I am fine with that small difference. What is bothering me is RasP.Normal_RasPNot.Normal contrast has about 608+367=~1k DEGs, which is way more than RasP.Tumor_RasPNot.Tumor contrast. I did plot plotMDS and see those normals are mixed together (I did see two subgroups seems having batch effects, but for each batch (if there is), there are RasP.Normal and RasPNot.Normal mixed there), all seem not supposed to having 1k DEGs either. Plus, the tumor contrast RasP.Tumor_RasPNot.Tumor has such small set of DEGs. Biologically, it might be true (population, ethics group, batch effect etc), but I just want to make sure what I did is reasonable. Any advice would be highly appreciated! Thanks a lot for your help! Best Ming ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:4}}
limma limma • 1.3k views
ADD COMMENT •
0
Entering edit mode
Ming ▴ 380
@ming-yi-6363
Last seen 19 months ago
United States
Dear Yunshun: Thanks so much for your advice, especially design <- model.matrix(~0+RasTum), which is new to me without intercept, and seems a great idea. But I am still not quite following what exactly the following: > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ RasP.Tumor)-RasPNot.Tumor, levels=design ) Intend to pool RasP1Hit.Tumor and RasP.Tumor as one group to compare to RasPNot class, why 0.5*? why not directly RasP_RasPNot=(RasP1Hit.Tumor+ RasP.Tumor)-RasPNot.Tumor? As extension to a more combined case, such as when I build model with all of these tumor types and their corresponding normals, to get contrast of overall tumors vs normals, I need pool all of the tumor types: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor and also need to pool all of the normal types: RasP1Hit.Normal, RasP.Normal and RasPNot.Normal and so to make the contrast, shall I do the similar thing: > con.matrix <- makeContrasts( Tumor_Normal=(1/3)*(RasP1Hit.Tumor+ RasP.Tumor+RasPNot.Tumor)-(1/3)*( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) Or > con.matrix <- makeContrasts( Tumor_Normal=(RasP1Hit.Tumor+ RasP.Tumor+RasPNot.Tumor)-( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) Which one is the right one? Bottom line is I am confused by what 0.5* means in the makeContrasts command, what is the meaning of that? Similarly, the way you just advised: con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) is very nice with design matrix column as Ras.Tum, Ras1H.Tum, RasNot.Tum, Ras.Nor, Ras1H.Nor, RasNot.Nor. Again, why Tum_Nor=1/3*c(1,1,1,-1,-1,-1) rather than Tum_Nor=c(1,1,1,-1,-1,-1)? In other word, if I use full column names to set up the contrast: I have to do: Tum_Nor=1/3 *( Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor- Ras1H.Nor-RasNot.Nor) Rather than directly: Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor Not sure what would be difference and impact of the first one with 1/3*. I thought usually we used the second way (straightforward): Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor Does this mean the second way is incorrect? If yes, why incorrect? Thanks again for your help! I never had seen such design way in the limma or edgeR documentation (I mean 0.5* or 1/3* etc) Best Ming > From: yuchen@wehi.EDU.AU > To: yi02@hotmail.com > CC: bioconductor@r-project.org > Subject: RE: [BioC] Problem in limma-voom > Date: Wed, 26 Feb 2014 11:19:56 +1100 > > Dear Ming, > > Regarding to your first question, a much easier way is to create your design > matrix without the intercept. > > design <- model.matrix(~0+RasTum) > > Then if you intend to pool both RasP class and RasP1Hit class together to > compare to RasPNot class, you can use the following contrast: > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > Similarly, for your second question, you can avoid confusion by creating > your design matrix without the intercept. > Suppose the six columns of your design matrix correspond to: Ras(Tum), > Ras1H(Tum), RasNot(Tum), Ras(Nor), Ras1H(Nor), RasNot(Nor). > Then you can make contrasts for any comparisons you are interested. > For example: > > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , > Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > Hope that helps. > > Yunshun Chen > > > ------------------------------------------------------------------- > > Hi, Dr. Gordon and all: > > > > I run into some limma design issue not sure what they are: > > > > I had three tumor types want to compare with each other: RasP, RasP1Hit, > RasPNot > > > > here are the critical commands: > > >targets<-targets[targets$RasType!="RasOnly" & targets$Type=="Tumor",] > > > show(head(targets)); > Subject SampleName Type RasPType RasType > 1 TCGA_38_4625 T4625_01A Tumor RasP RasP > 5 TCGA_38_4627 T4627_01A Tumor RasP RasP > 7 TCGA_38_4632 T4632_01A Tumor RasP RasP > 9 TCGA_44_2655 T2655_01A Tumor RasP RasP > 11 TCGA_44_2657 T2657_01A Tumor RasP RasP > 13 TCGA_44_2661 T2661_01A Tumor RasP RasP > > > > > RasTum<-paste(targets$RasType, targets$Type,sep=".") > > RasTum<-factor(RasTum, > levels=c("RasP1Hit.Tumor","RasP.Tumor","RasPNot.Tumor")); > > show(RasTum); > [1] RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor > [6] RasP.Tumor RasP.Tumor RasP.Tumor RasP1Hit.Tumor > RasP1Hit.Tumor > [11] RasP.Tumor RasPNot.Tumor RasP1Hit.Tumor RasP.Tumor > RasP1Hit.Tumor > [16] RasP1Hit.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > RasP1Hit.Tumor > [21] RasP1Hit.Tumor RasP1Hit.Tumor RasPNot.Tumor RasP.Tumor > RasPNot.Tumor > [26] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > RasPNot.Tumor > [31] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > RasPNot.Tumor > [36] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasP.Tumor > [41] RasP.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > RasPNot.Tumor > Levels: RasP1Hit.Tumor RasP.Tumor RasPNot.Tumor > > > > design<-model.matrix(~RasTum); > > colnames(design)<-sub("RasTum","",colnames(design)); > > colnames(design)[1] <- "Intercept" > > rownames(design)<-targets$Sample > > show(design[1:5,]) > Intercept RasP.Tumor RasPNot.Tumor > T4625_01A 1 1 0 > T4627_01A 1 1 0 > T4632_01A 1 1 0 > T2655_01A 1 1 0 > T2657_01A 1 1 0 > ... > > > con.matrix<-makeContrasts( > + > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > .Tumor=RasP.Tumor-RasPNot.Tumor, > + RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor, > + levels=design) > > > > show(con.matrix); > Contrasts > Levels RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > Intercept 0 0 > RasP.Tumor 1 1 > RasPNot.Tumor -1 -1 > Contrasts > Levels RasP1Hit.Tumor_RasPNot.Tumor > Intercept 0 > RasP.Tumor 0 > RasPNot.Tumor -1 > > > > > My first question is: as shown in show(RasTum), there are three subtype of > tumors I want to compare between: > > Levels: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor > > shown in show(design[1:5,]), everything is relative to RasP1Hit.Tumor. So in > my con.matrix<-makeContrasts command: > > RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > in particular, RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor since everything > is relative to RasP1Hit.Tumor. > > However, for contrast: RasPRasP1Hit.Tumor_RasPNot.Tumor, I intended to pool > both RasP class and RasP1Hit class together to compare to RasPNot class, I > end up with: > > RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > the right side is the same as > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > because everything is relative to RasP1Hit.Tumor. So what is the right way > to set up contrast RasPRasP1Hit.Tumor_RasPNot.Tumor in makeContrasts > command? kind of confused here? > > > > Because of the setting, the results are as below: > > > show(summary(de <- decideTests(fit))); > RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > -1 0 0 > 0 17764 17764 > 1 0 0 > RasP1Hit.Tumor_RasPNot.Tumor > -1 0 > 0 17764 > 1 0 > > > > > > My 2nd question is: I did a similar design with including all normal > samples, since at the same model I want to get tumor vs normal contrasts as > well. > > > con.matrix<-makeContrasts( > + RasP.Tumor_RasP.Normal = > RasP.Tumor-RasP.Normal,RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal = > RasP.Tumor-RasP.Normal-RasP1Hit.Normal, > + RasPNot.Tumor_RasPNot.Normal=RasPNot.Tumor-RasPNot.Normal, > + > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > .Tumor=RasP.Tumor-RasPNot.Tumor, > + > RasP.Normal_RasPNot.Normal=RasP.Normal- RasPNot.Normal,RasPRasP1Hit.Normal_Ra > sPNot.Normal=RasP1Hit.Normal+RasP.Normal-RasPNot.Normal, > + Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+-RasP.Normal- RasPNot.Normal, > + Tumor_Normal_2=RasP.Tumor+ > RasPNot.Tumor-RasP.Normal-RasPNot.Normal-RasP1Hit.Normal, > + levels=design) > > > > here everything is relative to RasP1Hit.Tumor. and the result as below: > > > > show(summary(de <- decideTests(fit))); > RasP.Tumor_RasP.Normal RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal > -1 3947 4522 > 0 9912 8830 > 1 4274 4781 > RasPNot.Tumor_RasPNot.Normal RasP.Tumor_RasPNot.Tumor > -1 4809 24 > 0 8461 18102 > 1 4863 7 > RasPRasP1Hit.Tumor_RasPNot.Tumor RasP.Normal_RasPNot.Normal > -1 24 608 > 0 18102 17158 > 1 7 367 > RasPRasP1Hit.Normal_RasPNot.Normal Tumor_Normal Tumor_Normal_2 > -1 2001 5611 5580 > 0 13934 6728 6651 > 1 2198 5794 5902 > > > > RasP.Tumor_RasPNot.Tumor contrast has small number of DEGs (24+7) at FDR 5% > compared to using tumor data alone for the model shown earlier above, this > might be cuased by filtering as well and I am fine with that small > difference. What is bothering me is RasP.Normal_RasPNot.Normal contrast has > about 608+367=~1k DEGs, which is way more than RasP.Tumor_RasPNot.Tumor > contrast. I did plot plotMDS and see those normals are mixed together (I did > see two subgroups seems having batch effects, but for each batch (if there > is), there are RasP.Normal and RasPNot.Normal mixed there), all seem not > supposed to having 1k DEGs either. Plus, the tumor contrast > RasP.Tumor_RasPNot.Tumor has such small set of DEGs. Biologically, it might > be true (population, ethics group, batch effect etc), but I just want to > make sure what I did is reasonable. > > > > Any advice would be highly appreciated! > > Thanks a lot for your help! > > > > Best > > Ming > > > ______________________________________________________________________ > The information in this email is confidential and inte...{{dropped:9}}
ADD COMMENT •
0
Entering edit mode
Dear Yunshun: I quickly tested the difference, based on numbes of DEGs, apparently, with or without 0.5* (or 1/3*) seem not have any impact on the final results of DEGs. > show(summary(de <- decideTests(fit))); ... Tumor_Normal -1 6365 0 5941 1 5944 Tumor_Normal_2 Tumor_Normal_3 -1 6365 6365 0 5941 5941 1 5944 5944 the makeContrasts command is as below: con.matrix<-makeContrasts( ... Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly .Tumor-RasP1Hit.Normal-RasP.Normal-RasPNot.Normal-RasOnly.Normal, Tumor_Normal_2=0.25*(RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOn ly.Tumor)-0.25*(RasP1Hit.Normal+RasP.Normal+RasPNot.Normal+RasOnly.Nor mal), Tumor_Normal_3=0.25*c(1,-1,1,-1,1,-1,1,-1),levels=design) it seems that all of the settting work exactly the same way, does not matter whether there is 0.5* or not (here is 0.25*). Please confirm what I found. or there is a reason for 0.5* etc. Thx a lot! Best Ming > From: yi02@hotmail.com > To: yuchen@wehi.edu.au > Date: Wed, 26 Feb 2014 15:06:14 +0000 > CC: bioconductor@r-project.org > Subject: Re: [BioC] Problem in limma-voom > > Dear Yunshun: > > Thanks so much for your advice, especially design <- model.matrix(~0+RasTum), which is new to me without intercept, and seems a great idea. But I am still not quite following what exactly the following: > > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > Intend to pool RasP1Hit.Tumor and RasP.Tumor as one group to compare to RasPNot class, why 0.5*? why not directly RasP_RasPNot=(RasP1Hit.Tumor+ RasP.Tumor)-RasPNot.Tumor? > > As extension to a more combined case, such as when I build model with all of these tumor types and their corresponding normals, to get contrast of overall tumors vs normals, I need pool all of the tumor types: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor and also need to pool all of the normal types: RasP1Hit.Normal, RasP.Normal and RasPNot.Normal and so to make the contrast, shall I do the similar thing: > > > con.matrix <- makeContrasts( Tumor_Normal=(1/3)*(RasP1Hit.Tumor+ > RasP.Tumor+RasPNot.Tumor)-(1/3)*( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) > > Or > > con.matrix <- makeContrasts( Tumor_Normal=(RasP1Hit.Tumor+ > RasP.Tumor+RasPNot.Tumor)-( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) > > Which one is the right one? Bottom line is I am confused by what 0.5* means in the makeContrasts command, what is the meaning of that? > > Similarly, the way you just advised: > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > is very nice with design matrix column as Ras.Tum, > Ras1H.Tum, RasNot.Tum, Ras.Nor, Ras1H.Nor, RasNot.Nor. > > Again, why Tum_Nor=1/3*c(1,1,1,-1,-1,-1) rather than Tum_Nor=c(1,1,1,-1,-1,-1)? In other word, if I use full column names to set up the contrast: > I have to do: Tum_Nor=1/3 *( Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor- Ras1H.Nor-RasNot.Nor) > Rather than directly: > Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor > > Not sure what would be difference and impact of the first one with 1/3*. I thought usually we used the second way (straightforward): Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor > > Does this mean the second way is incorrect? If yes, why incorrect? > > Thanks again for your help! I never had seen such design way in the limma or edgeR documentation (I mean 0.5* or 1/3* etc) > > Best > > Ming > > > > From: yuchen@wehi.EDU.AU > > To: yi02@hotmail.com > > CC: bioconductor@r-project.org > > Subject: RE: [BioC] Problem in limma-voom > > Date: Wed, 26 Feb 2014 11:19:56 +1100 > > > > Dear Ming, > > > > Regarding to your first question, a much easier way is to create your design > > matrix without the intercept. > > > design <- model.matrix(~0+RasTum) > > > > Then if you intend to pool both RasP class and RasP1Hit class together to > > compare to RasPNot class, you can use the following contrast: > > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > > > Similarly, for your second question, you can avoid confusion by creating > > your design matrix without the intercept. > > Suppose the six columns of your design matrix correspond to: Ras(Tum), > > Ras1H(Tum), RasNot(Tum), Ras(Nor), Ras1H(Nor), RasNot(Nor). > > Then you can make contrasts for any comparisons you are interested. > > For example: > > > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , > > Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > > > Hope that helps. > > > > Yunshun Chen > > > > > > ------------------------------------------------------------------- > > > > Hi, Dr. Gordon and all: > > > > > > > > I run into some limma design issue not sure what they are: > > > > > > > > I had three tumor types want to compare with each other: RasP, RasP1Hit, > > RasPNot > > > > > > > > here are the critical commands: > > > > >targets<-targets[targets$RasType!="RasOnly" & targets$Type=="Tumor",] > > > > > show(head(targets)); > > Subject SampleName Type RasPType RasType > > 1 TCGA_38_4625 T4625_01A Tumor RasP RasP > > 5 TCGA_38_4627 T4627_01A Tumor RasP RasP > > 7 TCGA_38_4632 T4632_01A Tumor RasP RasP > > 9 TCGA_44_2655 T2655_01A Tumor RasP RasP > > 11 TCGA_44_2657 T2657_01A Tumor RasP RasP > > 13 TCGA_44_2661 T2661_01A Tumor RasP RasP > > > > > > > > > RasTum<-paste(targets$RasType, targets$Type,sep=".") > > > RasTum<-factor(RasTum, > > levels=c("RasP1Hit.Tumor","RasP.Tumor","RasPNot.Tumor")); > > > show(RasTum); > > [1] RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor > > [6] RasP.Tumor RasP.Tumor RasP.Tumor RasP1Hit.Tumor > > RasP1Hit.Tumor > > [11] RasP.Tumor RasPNot.Tumor RasP1Hit.Tumor RasP.Tumor > > RasP1Hit.Tumor > > [16] RasP1Hit.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > > RasP1Hit.Tumor > > [21] RasP1Hit.Tumor RasP1Hit.Tumor RasPNot.Tumor RasP.Tumor > > RasPNot.Tumor > > [26] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > > RasPNot.Tumor > > [31] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > > RasPNot.Tumor > > [36] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasP.Tumor > > [41] RasP.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > > RasPNot.Tumor > > Levels: RasP1Hit.Tumor RasP.Tumor RasPNot.Tumor > > > > > > > design<-model.matrix(~RasTum); > > > colnames(design)<-sub("RasTum","",colnames(design)); > > > colnames(design)[1] <- "Intercept" > > > rownames(design)<-targets$Sample > > > show(design[1:5,]) > > Intercept RasP.Tumor RasPNot.Tumor > > T4625_01A 1 1 0 > > T4627_01A 1 1 0 > > T4632_01A 1 1 0 > > T2655_01A 1 1 0 > > T2657_01A 1 1 0 > > ... > > > > > con.matrix<-makeContrasts( > > + > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > > .Tumor=RasP.Tumor-RasPNot.Tumor, > > + RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor, > > + levels=design) > > > > > > > show(con.matrix); > > Contrasts > > Levels RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > > Intercept 0 0 > > RasP.Tumor 1 1 > > RasPNot.Tumor -1 -1 > > Contrasts > > Levels RasP1Hit.Tumor_RasPNot.Tumor > > Intercept 0 > > RasP.Tumor 0 > > RasPNot.Tumor -1 > > > > > > > > > > My first question is: as shown in show(RasTum), there are three subtype of > > tumors I want to compare between: > > > > Levels: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor > > > > shown in show(design[1:5,]), everything is relative to RasP1Hit.Tumor. So in > > my con.matrix<-makeContrasts command: > > > > RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor > > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > > > in particular, RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor since everything > > is relative to RasP1Hit.Tumor. > > > > However, for contrast: RasPRasP1Hit.Tumor_RasPNot.Tumor, I intended to pool > > both RasP class and RasP1Hit class together to compare to RasPNot class, I > > end up with: > > > > RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > the right side is the same as > > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > because everything is relative to RasP1Hit.Tumor. So what is the right way > > to set up contrast RasPRasP1Hit.Tumor_RasPNot.Tumor in makeContrasts > > command? kind of confused here? > > > > > > > > Because of the setting, the results are as below: > > > > > show(summary(de <- decideTests(fit))); > > RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > > -1 0 0 > > 0 17764 17764 > > 1 0 0 > > RasP1Hit.Tumor_RasPNot.Tumor > > -1 0 > > 0 17764 > > 1 0 > > > > > > > > > > > > My 2nd question is: I did a similar design with including all normal > > samples, since at the same model I want to get tumor vs normal contrasts as > > well. > > > > > con.matrix<-makeContrasts( > > + RasP.Tumor_RasP.Normal = > > RasP.Tumor-RasP.Normal,RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal = > > RasP.Tumor-RasP.Normal-RasP1Hit.Normal, > > + RasPNot.Tumor_RasPNot.Normal=RasPNot.Tumor-RasPNot.Normal, > > + > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > > .Tumor=RasP.Tumor-RasPNot.Tumor, > > + > > RasP.Normal_RasPNot.Normal=RasP.Normal- RasPNot.Normal,RasPRasP1Hit.Normal_Ra > > sPNot.Normal=RasP1Hit.Normal+RasP.Normal-RasPNot.Normal, > > + Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+-RasP.Normal- RasPNot.Normal, > > + Tumor_Normal_2=RasP.Tumor+ > > RasPNot.Tumor-RasP.Normal-RasPNot.Normal-RasP1Hit.Normal, > > + levels=design) > > > > > > > > here everything is relative to RasP1Hit.Tumor. and the result as below: > > > > > > > show(summary(de <- decideTests(fit))); > > RasP.Tumor_RasP.Normal RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal > > -1 3947 4522 > > 0 9912 8830 > > 1 4274 4781 > > RasPNot.Tumor_RasPNot.Normal RasP.Tumor_RasPNot.Tumor > > -1 4809 24 > > 0 8461 18102 > > 1 4863 7 > > RasPRasP1Hit.Tumor_RasPNot.Tumor RasP.Normal_RasPNot.Normal > > -1 24 608 > > 0 18102 17158 > > 1 7 367 > > RasPRasP1Hit.Normal_RasPNot.Normal Tumor_Normal Tumor_Normal_2 > > -1 2001 5611 5580 > > 0 13934 6728 6651 > > 1 2198 5794 5902 > > > > > > > > RasP.Tumor_RasPNot.Tumor contrast has small number of DEGs (24+7) at FDR 5% > > compared to using tumor data alone for the model shown earlier above, this > > might be cuased by filtering as well and I am fine with that small > > difference. What is bothering me is RasP.Normal_RasPNot.Normal contrast has > > about 608+367=~1k DEGs, which is way more than RasP.Tumor_RasPNot.Tumor > > contrast. I did plot plotMDS and see those normals are mixed together (I did > > see two subgroups seems having batch effects, but for each batch (if there > > is), there are RasP.Normal and RasPNot.Normal mixed there), all seem not > > supposed to having 1k DEGs either. Plus, the tumor contrast > > RasP.Tumor_RasPNot.Tumor has such small set of DEGs. Biologically, it might > > be true (population, ethics group, batch effect etc), but I just want to > > make sure what I did is reasonable. > > > > > > > > Any advice would be highly appreciated! > > > > Thanks a lot for your help! > > > > > > > > Best > > > > Ming > > > > > > ______________________________________________________________________ > > The information in this email is confidential and inte...{{dropped:9}} > > _______________________________________________ > Bioconductor mailing list > Bioconductor@r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor [[alternative HTML version deleted]]
ADD REPLY •
0
Entering edit mode
Hi, YunShun: To follow up with my last email, I did test in other contrasts: con.matrix<-makeContrasts( ... RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor+RasP1Hit.Tumor-RasPNot.Tum or,RasPRasP1Hit.Tumor_RasPNot.Tumor_2=0.5*(RasP.Tumor+RasP1Hit.Tumor)- RasPNot.Tumor, Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly.Tumor- RasP1Hit.Normal-RasP.Normal-RasPNot.Normal-RasOnly.Normal, Tumor_Normal_2=0.25*(RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOn ly.Tumor)-0.25*(RasP1Hit.Normal+RasP.Normal+RasPNot.Normal+RasOnly.Nor mal), Tumor_Normal_3=0.25*c(1,-1,1,-1,1,-1,1,-1),levels=design) here is what I got: > show(summary(de <- decideTests(fit))); ... RasPRasP1Hit.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor_2 -1 2657 0 0 2316 18250 1 13277 0 RasPRasP1Hit.Normal_RasPNot.Normal_2 Tumor_Normal Tumor_Normal_2 -1 457 6365 6365 0 17731 5941 5941 1 62 5944 5944 Tumor_Normal_3 -1 6365 0 5941 1 5944 No impact on Tumor vs Normal contrast, but certinlay has big impact on RasPRasP1Hit.Tumor_RasPNot.Tumor contrast. Now I am really confused. Any advice? Thanks again and best Ming > From: yi02@hotmail.com > To: yuchen@wehi.edu.au > Date: Wed, 26 Feb 2014 17:32:28 +0000 > CC: bioconductor@r-project.org > Subject: Re: [BioC] Problem in limma-voom > > Dear Yunshun: > > > > I quickly tested the difference, based on numbes of DEGs, apparently, with or without 0.5* (or 1/3*) seem not have any impact on the final results of DEGs. > > show(summary(de <- decideTests(fit))); > > ... > > Tumor_Normal > -1 6365 > 0 5941 > 1 5944 > > Tumor_Normal_2 Tumor_Normal_3 > -1 6365 6365 > 0 5941 5941 > 1 5944 5944 > > > the makeContrasts command is as below: > > con.matrix<-makeContrasts( > ... Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly.Tumor-RasP1Hit.Normal-RasP .Normal-RasPNot.Normal-RasOnly.Normal, > Tumor_Normal_2=0.25*(RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+Ras Only.Tumor)-0.25*(RasP1Hit.Normal+RasP.Normal+RasPNot.Normal+RasOnly.N ormal), > Tumor_Normal_3=0.25*c(1,-1,1,-1,1,-1,1,-1),levels=design) > > > > it seems that all of the settting work exactly the same way, does not matter whether there is 0.5* or not (here is 0.25*). Please confirm what I found. or there is a reason for 0.5* etc. Thx a lot! > > > > Best > > > > Ming > > > > > From: yi02@hotmail.com > > To: yuchen@wehi.edu.au > > Date: Wed, 26 Feb 2014 15:06:14 +0000 > > CC: bioconductor@r-project.org > > Subject: Re: [BioC] Problem in limma-voom > > > > Dear Yunshun: > > > > Thanks so much for your advice, especially design <- model.matrix(~0+RasTum), which is new to me without intercept, and seems a great idea. But I am still not quite following what exactly the following: > > > > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > > > Intend to pool RasP1Hit.Tumor and RasP.Tumor as one group to compare to RasPNot class, why 0.5*? why not directly RasP_RasPNot=(RasP1Hit.Tumor+ RasP.Tumor)-RasPNot.Tumor? > > > > As extension to a more combined case, such as when I build model with all of these tumor types and their corresponding normals, to get contrast of overall tumors vs normals, I need pool all of the tumor types: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor and also need to pool all of the normal types: RasP1Hit.Normal, RasP.Normal and RasPNot.Normal and so to make the contrast, shall I do the similar thing: > > > > > con.matrix <- makeContrasts( Tumor_Normal=(1/3)*(RasP1Hit.Tumor+ > > RasP.Tumor+RasPNot.Tumor)-(1/3)*( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) > > > > Or > > > con.matrix <- makeContrasts( Tumor_Normal=(RasP1Hit.Tumor+ > > RasP.Tumor+RasPNot.Tumor)-( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) > > > > Which one is the right one? Bottom line is I am confused by what 0.5* means in the makeContrasts command, what is the meaning of that? > > > > Similarly, the way you just advised: > > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > > > is very nice with design matrix column as Ras.Tum, > > Ras1H.Tum, RasNot.Tum, Ras.Nor, Ras1H.Nor, RasNot.Nor. > > > > Again, why Tum_Nor=1/3*c(1,1,1,-1,-1,-1) rather than Tum_Nor=c(1,1,1,-1,-1,-1)? In other word, if I use full column names to set up the contrast: > > I have to do: Tum_Nor=1/3 *( Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor- Ras1H.Nor-RasNot.Nor) > > Rather than directly: > > Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor > > > > Not sure what would be difference and impact of the first one with 1/3*. I thought usually we used the second way (straightforward): Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor > > > > Does this mean the second way is incorrect? If yes, why incorrect? > > > > Thanks again for your help! I never had seen such design way in the limma or edgeR documentation (I mean 0.5* or 1/3* etc) > > > > Best > > > > Ming > > > > > > > From: yuchen@wehi.EDU.AU > > > To: yi02@hotmail.com > > > CC: bioconductor@r-project.org > > > Subject: RE: [BioC] Problem in limma-voom > > > Date: Wed, 26 Feb 2014 11:19:56 +1100 > > > > > > Dear Ming, > > > > > > Regarding to your first question, a much easier way is to create your design > > > matrix without the intercept. > > > > design <- model.matrix(~0+RasTum) > > > > > > Then if you intend to pool both RasP class and RasP1Hit class together to > > > compare to RasPNot class, you can use the following contrast: > > > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > > > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > > > > > Similarly, for your second question, you can avoid confusion by creating > > > your design matrix without the intercept. > > > Suppose the six columns of your design matrix correspond to: Ras(Tum), > > > Ras1H(Tum), RasNot(Tum), Ras(Nor), Ras1H(Nor), RasNot(Nor). > > > Then you can make contrasts for any comparisons you are interested. > > > For example: > > > > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , > > > Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > > > > > Hope that helps. > > > > > > Yunshun Chen > > > > > > > > > ------------------------------------------------------------------- > > > > > > Hi, Dr. Gordon and all: > > > > > > > > > > > > I run into some limma design issue not sure what they are: > > > > > > > > > > > > I had three tumor types want to compare with each other: RasP, RasP1Hit, > > > RasPNot > > > > > > > > > > > > here are the critical commands: > > > > > > >targets<-targets[targets$RasType!="RasOnly" & targets$Type=="Tumor",] > > > > > > > show(head(targets)); > > > Subject SampleName Type RasPType RasType > > > 1 TCGA_38_4625 T4625_01A Tumor RasP RasP > > > 5 TCGA_38_4627 T4627_01A Tumor RasP RasP > > > 7 TCGA_38_4632 T4632_01A Tumor RasP RasP > > > 9 TCGA_44_2655 T2655_01A Tumor RasP RasP > > > 11 TCGA_44_2657 T2657_01A Tumor RasP RasP > > > 13 TCGA_44_2661 T2661_01A Tumor RasP RasP > > > > > > > > > > > > > RasTum<-paste(targets$RasType, targets$Type,sep=".") > > > > RasTum<-factor(RasTum, > > > levels=c("RasP1Hit.Tumor","RasP.Tumor","RasPNot.Tumor")); > > > > show(RasTum); > > > [1] RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor > > > [6] RasP.Tumor RasP.Tumor RasP.Tumor RasP1Hit.Tumor > > > RasP1Hit.Tumor > > > [11] RasP.Tumor RasPNot.Tumor RasP1Hit.Tumor RasP.Tumor > > > RasP1Hit.Tumor > > > [16] RasP1Hit.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > > > RasP1Hit.Tumor > > > [21] RasP1Hit.Tumor RasP1Hit.Tumor RasPNot.Tumor RasP.Tumor > > > RasPNot.Tumor > > > [26] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > [31] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > [36] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasP.Tumor > > > [41] RasP.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > Levels: RasP1Hit.Tumor RasP.Tumor RasPNot.Tumor > > > > > > > > > > design<-model.matrix(~RasTum); > > > > colnames(design)<-sub("RasTum","",colnames(design)); > > > > colnames(design)[1] <- "Intercept" > > > > rownames(design)<-targets$Sample > > > > show(design[1:5,]) > > > Intercept RasP.Tumor RasPNot.Tumor > > > T4625_01A 1 1 0 > > > T4627_01A 1 1 0 > > > T4632_01A 1 1 0 > > > T2655_01A 1 1 0 > > > T2657_01A 1 1 0 > > > ... > > > > > > > con.matrix<-makeContrasts( > > > + > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > > > .Tumor=RasP.Tumor-RasPNot.Tumor, > > > + RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor, > > > + levels=design) > > > > > > > > > > show(con.matrix); > > > Contrasts > > > Levels RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > > > Intercept 0 0 > > > RasP.Tumor 1 1 > > > RasPNot.Tumor -1 -1 > > > Contrasts > > > Levels RasP1Hit.Tumor_RasPNot.Tumor > > > Intercept 0 > > > RasP.Tumor 0 > > > RasPNot.Tumor -1 > > > > > > > > > > > > > > > My first question is: as shown in show(RasTum), there are three subtype of > > > tumors I want to compare between: > > > > > > Levels: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor > > > > > > shown in show(design[1:5,]), everything is relative to RasP1Hit.Tumor. So in > > > my con.matrix<-makeContrasts command: > > > > > > RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor > > > > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > > > > > > > in particular, RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor since everything > > > is relative to RasP1Hit.Tumor. > > > > > > However, for contrast: RasPRasP1Hit.Tumor_RasPNot.Tumor, I intended to pool > > > both RasP class and RasP1Hit class together to compare to RasPNot class, I > > > end up with: > > > > > > RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > the right side is the same as > > > > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > because everything is relative to RasP1Hit.Tumor. So what is the right way > > > to set up contrast RasPRasP1Hit.Tumor_RasPNot.Tumor in makeContrasts > > > command? kind of confused here? > > > > > > > > > > > > Because of the setting, the results are as below: > > > > > > > show(summary(de <- decideTests(fit))); > > > RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > > > -1 0 0 > > > 0 17764 17764 > > > 1 0 0 > > > RasP1Hit.Tumor_RasPNot.Tumor > > > -1 0 > > > 0 17764 > > > 1 0 > > > > > > > > > > > > > > > > > > My 2nd question is: I did a similar design with including all normal > > > samples, since at the same model I want to get tumor vs normal contrasts as > > > well. > > > > > > > con.matrix<-makeContrasts( > > > + RasP.Tumor_RasP.Normal = > > > RasP.Tumor-RasP.Normal,RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal = > > > RasP.Tumor-RasP.Normal-RasP1Hit.Normal, > > > + RasPNot.Tumor_RasPNot.Normal=RasPNot.Tumor-RasPNot.Normal, > > > + > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > > > .Tumor=RasP.Tumor-RasPNot.Tumor, > > > + > > > RasP.Normal_RasPNot.Normal=RasP.Normal- RasPNot.Normal,RasPRasP1Hit.Normal_Ra > > > sPNot.Normal=RasP1Hit.Normal+RasP.Normal-RasPNot.Normal, > > > + Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+-RasP.Normal- RasPNot.Normal, > > > + Tumor_Normal_2=RasP.Tumor+ > > > RasPNot.Tumor-RasP.Normal-RasPNot.Normal-RasP1Hit.Normal, > > > + levels=design) > > > > > > > > > > > > here everything is relative to RasP1Hit.Tumor. and the result as below: > > > > > > > > > > show(summary(de <- decideTests(fit))); > > > RasP.Tumor_RasP.Normal RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal > > > -1 3947 4522 > > > 0 9912 8830 > > > 1 4274 4781 > > > RasPNot.Tumor_RasPNot.Normal RasP.Tumor_RasPNot.Tumor > > > -1 4809 24 > > > 0 8461 18102 > > > 1 4863 7 > > > RasPRasP1Hit.Tumor_RasPNot.Tumor RasP.Normal_RasPNot.Normal > > > -1 24 608 > > > 0 18102 17158 > > > 1 7 367 > > > RasPRasP1Hit.Normal_RasPNot.Normal Tumor_Normal Tumor_Normal_2 > > > -1 2001 5611 5580 > > > 0 13934 6728 6651 > > > 1 2198 5794 5902 > > > > > > > > > > > > RasP.Tumor_RasPNot.Tumor contrast has small number of DEGs (24+7) at FDR 5% > > > compared to using tumor data alone for the model shown earlier above, this > > > might be cuased by filtering as well and I am fine with that small > > > difference. What is bothering me is RasP.Normal_RasPNot.Normal contrast has > > > about 608+367=~1k DEGs, which is way more than RasP.Tumor_RasPNot.Tumor > > > contrast. I did plot plotMDS and see those normals are mixed together (I did > > > see two subgroups seems having batch effects, but for each batch (if there > > > is), there are RasP.Normal and RasPNot.Normal mixed there), all seem not > > > supposed to having 1k DEGs either. Plus, the tumor contrast > > > RasP.Tumor_RasPNot.Tumor has such small set of DEGs. Biologically, it might > > > be true (population, ethics group, batch effect etc), but I just want to > > > make sure what I did is reasonable. > > > > > > > > > > > > Any advice would be highly appreciated! > > > > > > Thanks a lot for your help! > > > > > > > > > > > > Best > > > > > > Ming > > > > > > > > > ______________________________________________________________________ > > > The information in this email is confidential and inte...{{dropped:9}} > > > > _______________________________________________ > > Bioconductor mailing list > > Bioconductor@r-project.org > > https://stat.ethz.ch/mailman/listinfo/bioconductor > > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor > > [[alternative HTML version deleted]] > > _______________________________________________ > Bioconductor mailing list > Bioconductor@r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor [[alternative HTML version deleted]]
ADD REPLY •
0
Entering edit mode
Hi Ming, Suppose we have three groups A, B and C. When you pool A and B together and compare with C, you can think of it as comparing the average of A and B with C. That is: 0/5*(A+B) - C. You will have the same testing results in terms of the number of DEGs if you use: (A+B)-2C instead of the above. However, the coefficients (or the logFCs in your case) will be doubled and cannot be interpreted as the logFCs for your comparisons. Regards, Yunshun Chen From: Ming Yi [mailto:yi02@hotmail.com] Sent: Thursday, 27 February 2014 7:01 AM To: yuchen@wehi.EDU.AU Cc: Bioconductor mailing list Subject: RE: [BioC] Problem in limma-voom Hi, YunShun: To follow up with my last email, I did test in other contrasts: con.matrix<-makeContrasts( ... RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor+RasP1Hit.Tumor- RasPNot.Tumor,Ras PRasP1Hit.Tumor_RasPNot.Tumor_2=0.5*(RasP.Tumor+RasP1Hit.Tumor)-RasPNo t.Tumo r, Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly.Tumor-RasP1Hit.Normal-RasP .Normal-RasPN ot.Normal-RasOnly.Normal, Tumor_Normal_2=0.25*(RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly.Tumor)-0.25*(RasP1Hit.Normal+RasP .Norma l+RasPNot.Normal+RasOnly.Normal), Tumor_Normal_3=0.25*c(1,-1,1,-1,1,-1,1,-1),levels=design) here is what I got: > show(summary(de <- decideTests(fit))); ... RasPRasP1Hit.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor_2 -1 2657 0 0 2316 18250 1 13277 0 RasPRasP1Hit.Normal_RasPNot.Normal_2 Tumor_Normal Tumor_Normal_2 -1 457 6365 6365 0 17731 5941 5941 1 62 5944 5944 Tumor_Normal_3 -1 6365 0 5941 1 5944 No impact on Tumor vs Normal contrast, but certinlay has big impact on RasPRasP1Hit.Tumor_RasPNot.Tumor contrast. Now I am really confused. Any advice? Thanks again and best Ming > From: yi02@hotmail.com > To: yuchen@wehi.edu.au > Date: Wed, 26 Feb 2014 17:32:28 +0000 > CC: bioconductor@r-project.org > Subject: Re: [BioC] Problem in limma-voom > > Dear Yunshun: > > > > I quickly tested the difference, based on numbes of DEGs, apparently, with or without 0.5* (or 1/3*) seem not have any impact on the final results of DEGs. > > show(summary(de <- decideTests(fit))); > > ... > > Tumor_Normal > -1 6365 > 0 5941 > 1 5944 > > Tumor_Normal_2 Tumor_Normal_3 > -1 6365 6365 > 0 5941 5941 > 1 5944 5944 > > > the makeContrasts command is as below: > > con.matrix<-makeContrasts( > ... Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly.Tumor-RasP1Hit.Normal-RasP .Normal-RasPN ot.Normal-RasOnly.Normal, > Tumor_Normal_2=0.25*(RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly.Tumor)-0.25*(RasP1Hit.Normal+RasP .Norma l+RasPNot.Normal+RasOnly.Normal), > Tumor_Normal_3=0.25*c(1,-1,1,-1,1,-1,1,-1),levels=design) > > > > it seems that all of the settting work exactly the same way, does not matter whether there is 0.5* or not (here is 0.25*). Please confirm what I found. or there is a reason for 0.5* etc. Thx a lot! > > > > Best > > > > Ming > > > > > From: yi02@hotmail.com > > To: yuchen@wehi.edu.au > > Date: Wed, 26 Feb 2014 15:06:14 +0000 > > CC: bioconductor@r-project.org > > Subject: Re: [BioC] Problem in limma-voom > > > > Dear Yunshun: > > > > Thanks so much for your advice, especially design <- model.matrix(~0+RasTum), which is new to me without intercept, and seems a great idea. But I am still not quite following what exactly the following: > > > > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > > > Intend to pool RasP1Hit.Tumor and RasP.Tumor as one group to compare to RasPNot class, why 0.5*? why not directly RasP_RasPNot=(RasP1Hit.Tumor+ RasP.Tumor)-RasPNot.Tumor? > > > > As extension to a more combined case, such as when I build model with all of these tumor types and their corresponding normals, to get contrast of overall tumors vs normals, I need pool all of the tumor types: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor and also need to pool all of the normal types: RasP1Hit.Normal, RasP.Normal and RasPNot.Normal and so to make the contrast, shall I do the similar thing: > > > > > con.matrix <- makeContrasts( Tumor_Normal=(1/3)*(RasP1Hit.Tumor+ > > RasP.Tumor+RasPNot.Tumor)-(1/3)*( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) > > > > Or > > > con.matrix <- makeContrasts( Tumor_Normal=(RasP1Hit.Tumor+ > > RasP.Tumor+RasPNot.Tumor)-( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) > > > > Which one is the right one? Bottom line is I am confused by what 0.5* means in the makeContrasts command, what is the meaning of that? > > > > Similarly, the way you just advised: > > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > > > is very nice with design matrix column as Ras.Tum, > > Ras1H.Tum, RasNot.Tum, Ras.Nor, Ras1H.Nor, RasNot.Nor. > > > > Again, why Tum_Nor=1/3*c(1,1,1,-1,-1,-1) rather than Tum_Nor=c(1,1,1,-1,-1,-1)? In other word, if I use full column names to set up the contrast: > > I have to do: Tum_Nor=1/3 *( Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor) > > Rather than directly: > > Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor > > > > Not sure what would be difference and impact of the first one with 1/3*. I thought usually we used the second way (straightforward): Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor > > > > Does this mean the second way is incorrect? If yes, why incorrect? > > > > Thanks again for your help! I never had seen such design way in the limma or edgeR documentation (I mean 0.5* or 1/3* etc) > > > > Best > > > > Ming > > > > > > > From: yuchen@wehi.EDU.AU > > > To: yi02@hotmail.com > > > CC: bioconductor@r-project.org > > > Subject: RE: [BioC] Problem in limma-voom > > > Date: Wed, 26 Feb 2014 11:19:56 +1100 > > > > > > Dear Ming, > > > > > > Regarding to your first question, a much easier way is to create your design > > > matrix without the intercept. > > > > design <- model.matrix(~0+RasTum) > > > > > > Then if you intend to pool both RasP class and RasP1Hit class together to > > > compare to RasPNot class, you can use the following contrast: > > > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > > > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > > > > > Similarly, for your second question, you can avoid confusion by creating > > > your design matrix without the intercept. > > > Suppose the six columns of your design matrix correspond to: Ras(Tum), > > > Ras1H(Tum), RasNot(Tum), Ras(Nor), Ras1H(Nor), RasNot(Nor). > > > Then you can make contrasts for any comparisons you are interested. > > > For example: > > > > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , > > > Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor-RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > > > > > Hope that helps. > > > > > > Yunshun Chen > > > > > > > > > ------------------------------------------------------------------- > > > > > > Hi, Dr. Gordon and all: > > > > > > > > > > > > I run into some limma design issue not sure what they are: > > > > > > > > > > > > I had three tumor types want to compare with each other: RasP, RasP1Hit, > > > RasPNot > > > > > > > > > > > > here are the critical commands: > > > > > > >targets<-targets[targets$RasType!="RasOnly" & targets$Type=="Tumor",] > > > > > > > show(head(targets)); > > > Subject SampleName Type RasPType RasType > > > 1 TCGA_38_4625 T4625_01A Tumor RasP RasP > > > 5 TCGA_38_4627 T4627_01A Tumor RasP RasP > > > 7 TCGA_38_4632 T4632_01A Tumor RasP RasP > > > 9 TCGA_44_2655 T2655_01A Tumor RasP RasP > > > 11 TCGA_44_2657 T2657_01A Tumor RasP RasP > > > 13 TCGA_44_2661 T2661_01A Tumor RasP RasP > > > > > > > > > > > > > RasTum<-paste(targets$RasType, targets$Type,sep=".") > > > > RasTum<-factor(RasTum, > > > levels=c("RasP1Hit.Tumor","RasP.Tumor","RasPNot.Tumor")); > > > > show(RasTum); > > > [1] RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor > > > [6] RasP.Tumor RasP.Tumor RasP.Tumor RasP1Hit.Tumor > > > RasP1Hit.Tumor > > > [11] RasP.Tumor RasPNot.Tumor RasP1Hit.Tumor RasP.Tumor > > > RasP1Hit.Tumor > > > [16] RasP1Hit.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > > > RasP1Hit.Tumor > > > [21] RasP1Hit.Tumor RasP1Hit.Tumor RasPNot.Tumor RasP.Tumor > > > RasPNot.Tumor > > > [26] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > [31] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > [36] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasP.Tumor > > > [41] RasP.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > Levels: RasP1Hit.Tumor RasP.Tumor RasPNot.Tumor > > > > > > > > > > design<-model.matrix(~RasTum); > > > > colnames(design)<-sub("RasTum","",colnames(design)); > > > > colnames(design)[1] <- "Intercept" > > > > rownames(design)<-targets$Sample > > > > show(design[1:5,]) > > > Intercept RasP.Tumor RasPNot.Tumor > > > T4625_01A 1 1 0 > > > T4627_01A 1 1 0 > > > T4632_01A 1 1 0 > > > T2655_01A 1 1 0 > > > T2657_01A 1 1 0 > > > ... > > > > > > > con.matrix<-makeContrasts( > > > + > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > > > .Tumor=RasP.Tumor-RasPNot.Tumor, > > > + RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor, > > > + levels=design) > > > > > > > > > > show(con.matrix); > > > Contrasts > > > Levels RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > > > Intercept 0 0 > > > RasP.Tumor 1 1 > > > RasPNot.Tumor -1 -1 > > > Contrasts > > > Levels RasP1Hit.Tumor_RasPNot.Tumor > > > Intercept 0 > > > RasP.Tumor 0 > > > RasPNot.Tumor -1 > > > > > > > > > > > > > > > My first question is: as shown in show(RasTum), there are three subtype of > > > tumors I want to compare between: > > > > > > Levels: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor > > > > > > shown in show(design[1:5,]), everything is relative to RasP1Hit.Tumor. So in > > > my con.matrix<-makeContrasts command: > > > > > > RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor > > > > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > > > > > > > in particular, RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor since everything > > > is relative to RasP1Hit.Tumor. > > > > > > However, for contrast: RasPRasP1Hit.Tumor_RasPNot.Tumor, I intended to pool > > > both RasP class and RasP1Hit class together to compare to RasPNot class, I > > > end up with: > > > > > > RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > the right side is the same as > > > > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > because everything is relative to RasP1Hit.Tumor. So what is the right way > > > to set up contrast RasPRasP1Hit.Tumor_RasPNot.Tumor in makeContrasts > > > command? kind of confused here? > > > > > > > > > > > > Because of the setting, the results are as below: > > > > > > > show(summary(de <- decideTests(fit))); > > > RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > > > -1 0 0 > > > 0 17764 17764 > > > 1 0 0 > > > RasP1Hit.Tumor_RasPNot.Tumor > > > -1 0 > > > 0 17764 > > > 1 0 > > > > > > > > > > > > > > > > > > My 2nd question is: I did a similar design with including all normal > > > samples, since at the same model I want to get tumor vs normal contrasts as > > > well. > > > > > > > con.matrix<-makeContrasts( > > > + RasP.Tumor_RasP.Normal = > > > RasP.Tumor-RasP.Normal,RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal = > > > RasP.Tumor-RasP.Normal-RasP1Hit.Normal, > > > + RasPNot.Tumor_RasPNot.Normal=RasPNot.Tumor-RasPNot.Normal, > > > + > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > > > .Tumor=RasP.Tumor-RasPNot.Tumor, > > > + > > > RasP.Normal_RasPNot.Normal=RasP.Normal- RasPNot.Normal,RasPRasP1Hit.Normal_Ra > > > sPNot.Normal=RasP1Hit.Normal+RasP.Normal-RasPNot.Normal, > > > + Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+-RasP.Normal- RasPNot.Normal, > > > + Tumor_Normal_2=RasP.Tumor+ > > > RasPNot.Tumor-RasP.Normal-RasPNot.Normal-RasP1Hit.Normal, > > > + levels=design) > > > > > > > > > > > > here everything is relative to RasP1Hit.Tumor. and the result as below: > > > > > > > > > > show(summary(de <- decideTests(fit))); > > > RasP.Tumor_RasP.Normal RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal > > > -1 3947 4522 > > > 0 9912 8830 > > > 1 4274 4781 > > > RasPNot.Tumor_RasPNot.Normal RasP.Tumor_RasPNot.Tumor > > > -1 4809 24 > > > 0 8461 18102 > > > 1 4863 7 > > > RasPRasP1Hit.Tumor_RasPNot.Tumor RasP.Normal_RasPNot.Normal > > > -1 24 608 > > > 0 18102 17158 > > > 1 7 367 > > > RasPRasP1Hit.Normal_RasPNot.Normal Tumor_Normal Tumor_Normal_2 > > > -1 2001 5611 5580 > > > 0 13934 6728 6651 > > > 1 2198 5794 5902 > > > > > > > > > > > > RasP.Tumor_RasPNot.Tumor contrast has small number of DEGs (24+7) at FDR 5% > > > compared to using tumor data alone for the model shown earlier above, this > > > might be cuased by filtering as well and I am fine with that small > > > difference. What is bothering me is RasP.Normal_RasPNot.Normal contrast has > > > about 608+367=~1k DEGs, which is way more than RasP.Tumor_RasPNot.Tumor > > > contrast. I did plot plotMDS and see those normals are mixed together (I did > > > see two subgroups seems having batch effects, but for each batch (if there > > > is), there are RasP.Normal and RasPNot.Normal mixed there), all seem not > > > supposed to having 1k DEGs either. Plus, the tumor contrast > > > RasP.Tumor_RasPNot.Tumor has such small set of DEGs. Biologically, it might > > > be true (population, ethics group, batch effect etc), but I just want to > > > make sure what I did is reasonable. > > > > > > > > > > > > Any advice would be highly appreciated! > > > > > > Thanks a lot for your help! > > > > > > > > > > > > Best > > > > > > Ming > > > > > > > > > ______________________________________________________________________ > > > The information in this email is confidential and inte...{{dropped:9}} > > > > _______________________________________________ > > Bioconductor mailing list > > Bioconductor@r-project.org > > https://stat.ethz.ch/mailman/listinfo/bioconductor > > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor > > [[alternative HTML version deleted]] > > _______________________________________________ > Bioconductor mailing list > Bioconductor@r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:6}}
ADD REPLY •
0
Entering edit mode
Thx, Yunshun! I noticed that you are one of the authors of edgeR. Appreciated very much your help! It would be nice to include these details in documentation of both edgeR and limma. Thx again, Ming From: yuchen@wehi.EDU.AU To: yi02@hotmail.com CC: bioconductor@r-project.org Subject: RE: [BioC] Problem in limma-voom Date: Thu, 27 Feb 2014 10:33:47 +1100 Hi Ming, Suppose we have three groups A, B and C. When you pool A and B together and compare with C, you can think of it as comparing the average of A and B with C.That is: 0/5*(A+B) – C. You will have the same testing results in terms of the number of DEGs if you use: (A+B)-2C instead of the above.However, the coefficients (or the logFCs in your case) will be doubled and cannot be interpreted as the logFCs for your comparisons. Regards,Yunshun Chen From: Ming Yi [mailto:yi02@hotmail.com] Sent: Thursday, 27 February 2014 7:01 AM To: yuchen@wehi.EDU.AU Cc: Bioconductor mailing list Subject: RE: [BioC] Problem in limma-voom Hi, YunShun: To follow up with my last email, I did test in other contrasts: con.matrix<-makeContrasts( ... RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor+RasP1Hit.Tumor-RasPNot.Tum or,RasPRasP1Hit.Tumor_RasPNot.Tumor_2=0.5*(RasP.Tumor+RasP1Hit.Tumor)- RasPNot.Tumor, Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly.Tumor- RasP1Hit.Normal-RasP.Normal-RasPNot.Normal-RasOnly.Normal, Tumor_Normal_2=0.25*(RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOn ly.Tumor)-0.25*(RasP1Hit.Normal+RasP.Normal+RasPNot.Normal+RasOnly.Nor mal), Tumor_Normal_3=0.25*c(1,-1,1,-1,1,-1,1,-1),levels=design) here is what I got: > show(summary(de <- decideTests(fit))); ... RasPRasP1Hit.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor_2 -1 2657 0 0 2316 18250 1 13277 0 RasPRasP1Hit.Normal_RasPNot.Normal_2 Tumor_Normal Tumor_Normal_2 -1 457 6365 6365 0 17731 5941 5941 1 62 5944 5944 Tumor_Normal_3 -1 6365 0 5941 1 5944 No impact on Tumor vs Normal contrast, but certinlay has big impact on RasPRasP1Hit.Tumor_RasPNot.Tumor contrast. Now I am really confused. Any advice? Thanks again and best Ming > From: yi02@hotmail.com > To: yuchen@wehi.edu.au > Date: Wed, 26 Feb 2014 17:32:28 +0000 > CC: bioconductor@r-project.org > Subject: Re: [BioC] Problem in limma-voom > > Dear Yunshun: > > > > I quickly tested the difference, based on numbes of DEGs, apparently, with or without 0.5* (or 1/3*) seem not have any impact on the final results of DEGs. > > show(summary(de <- decideTests(fit))); > > ... > > Tumor_Normal > -1 6365 > 0 5941 > 1 5944 > > Tumor_Normal_2 Tumor_Normal_3 > -1 6365 6365 > 0 5941 5941 > 1 5944 5944 > > > the makeContrasts command is as below: > > con.matrix<-makeContrasts( > ... Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnly .Tumor-RasP1Hit.Normal-RasP.Normal-RasPNot.Normal-RasOnly.Normal, > Tumor_Normal_2=0.25*(RasP.Tumor+ RasPNot.Tumor+RasP1Hit.Tumor+RasOnl y.Tumor)-0.25*(RasP1Hit.Normal+RasP.Normal+RasPNot.Normal+RasOnly.Norm al), > Tumor_Normal_3=0.25*c(1,-1,1,-1,1,-1,1,-1),levels=design) > > > > it seems that all of the settting work exactly the same way, does not matter whether there is 0.5* or not (here is 0.25*). Please confirm what I found. or there is a reason for 0.5* etc. Thx a lot! > > > > Best > > > > Ming > > > > > From: yi02@hotmail.com > > To: yuchen@wehi.edu.au > > Date: Wed, 26 Feb 2014 15:06:14 +0000 > > CC: bioconductor@r-project.org > > Subject: Re: [BioC] Problem in limma-voom > > > > Dear Yunshun: > > > > Thanks so much for your advice, especially design <- model.matrix(~0+RasTum), which is new to me without intercept, and seems a great idea. But I am still not quite following what exactly the following: > > > > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > > > Intend to pool RasP1Hit.Tumor and RasP.Tumor as one group to compare to RasPNot class, why 0.5*? why not directly RasP_RasPNot=(RasP1Hit.Tumor+ RasP.Tumor)-RasPNot.Tumor? > > > > As extension to a more combined case, such as when I build model with all of these tumor types and their corresponding normals, to get contrast of overall tumors vs normals, I need pool all of the tumor types: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor and also need to pool all of the normal types: RasP1Hit.Normal, RasP.Normal and RasPNot.Normal and so to make the contrast, shall I do the similar thing: > > > > > con.matrix <- makeContrasts( Tumor_Normal=(1/3)*(RasP1Hit.Tumor+ > > RasP.Tumor+RasPNot.Tumor)-(1/3)*( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) > > > > Or > > > con.matrix <- makeContrasts( Tumor_Normal=(RasP1Hit.Tumor+ > > RasP.Tumor+RasPNot.Tumor)-( RasP1Hit.Normal+RasP.Normal +RasPNot.Normal), levels=design) > > > > Which one is the right one? Bottom line is I am confused by what 0.5* means in the makeContrasts command, what is the meaning of that? > > > > Similarly, the way you just advised: > > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > > > is very nice with design matrix column as Ras.Tum, > > Ras1H.Tum, RasNot.Tum, Ras.Nor, Ras1H.Nor, RasNot.Nor. > > > > Again, why Tum_Nor=1/3*c(1,1,1,-1,-1,-1) rather than Tum_Nor=c(1,1,1,-1,-1,-1)? In other word, if I use full column names to set up the contrast: > > I have to do: Tum_Nor=1/3 *( Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor- Ras1H.Nor-RasNot.Nor) > > Rather than directly: > > Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor > > > > Not sure what would be difference and impact of the first one with 1/3*. I thought usually we used the second way (straightforward): Tum_Nor= Ras.Tum+Ras1H.Tum+RasNot.Tum-Ras.Nor-Ras1H.Nor-RasNot.Nor > > > > Does this mean the second way is incorrect? If yes, why incorrect? > > > > Thanks again for your help! I never had seen such design way in the limma or edgeR documentation (I mean 0.5* or 1/3* etc) > > > > Best > > > > Ming > > > > > > > From: yuchen@wehi.EDU.AU > > > To: yi02@hotmail.com > > > CC: bioconductor@r-project.org > > > Subject: RE: [BioC] Problem in limma-voom > > > Date: Wed, 26 Feb 2014 11:19:56 +1100 > > > > > > Dear Ming, > > > > > > Regarding to your first question, a much easier way is to create your design > > > matrix without the intercept. > > > > design <- model.matrix(~0+RasTum) > > > > > > Then if you intend to pool both RasP class and RasP1Hit class together to > > > compare to RasPNot class, you can use the following contrast: > > > > con.matrix <- makeContrasts( RasP_RasPNot=0.5*(RasP1Hit.Tumor+ > > > RasP.Tumor)-RasPNot.Tumor, levels=design ) > > > > > > Similarly, for your second question, you can avoid confusion by creating > > > your design matrix without the intercept. > > > Suppose the six columns of your design matrix correspond to: Ras(Tum), > > > Ras1H(Tum), RasNot(Tum), Ras(Nor), Ras1H(Nor), RasNot(Nor). > > > Then you can make contrasts for any comparisons you are interested. > > > For example: > > > > con.matrix <- cbind( Tum_Nor=1/3*c(1,1,1,-1,-1,-1) , > > > Ras1H_RasNot=0.5*c(0,1,-1,0,1,-1), RasAllNor- RasNotNor=c(0,0,0,0.5,0.5,-1) ) > > > > > > Hope that helps. > > > > > > Yunshun Chen > > > > > > > > > ------------------------------------------------------------------- > > > > > > Hi, Dr. Gordon and all: > > > > > > > > > > > > I run into some limma design issue not sure what they are: > > > > > > > > > > > > I had three tumor types want to compare with each other: RasP, RasP1Hit, > > > RasPNot > > > > > > > > > > > > here are the critical commands: > > > > > > >targets<-targets[targets$RasType!="RasOnly" & targets$Type=="Tumor",] > > > > > > > show(head(targets)); > > > Subject SampleName Type RasPType RasType > > > 1 TCGA_38_4625 T4625_01A Tumor RasP RasP > > > 5 TCGA_38_4627 T4627_01A Tumor RasP RasP > > > 7 TCGA_38_4632 T4632_01A Tumor RasP RasP > > > 9 TCGA_44_2655 T2655_01A Tumor RasP RasP > > > 11 TCGA_44_2657 T2657_01A Tumor RasP RasP > > > 13 TCGA_44_2661 T2661_01A Tumor RasP RasP > > > > > > > > > > > > > RasTum<-paste(targets$RasType, targets$Type,sep=".") > > > > RasTum<-factor(RasTum, > > > levels=c("RasP1Hit.Tumor","RasP.Tumor","RasPNot.Tumor")); > > > > show(RasTum); > > > [1] RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor RasP.Tumor > > > [6] RasP.Tumor RasP.Tumor RasP.Tumor RasP1Hit.Tumor > > > RasP1Hit.Tumor > > > [11] RasP.Tumor RasPNot.Tumor RasP1Hit.Tumor RasP.Tumor > > > RasP1Hit.Tumor > > > [16] RasP1Hit.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > > > RasP1Hit.Tumor > > > [21] RasP1Hit.Tumor RasP1Hit.Tumor RasPNot.Tumor RasP.Tumor > > > RasPNot.Tumor > > > [26] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > [31] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > [36] RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasPNot.Tumor RasP.Tumor > > > [41] RasP.Tumor RasP.Tumor RasP.Tumor RasPNot.Tumor > > > RasPNot.Tumor > > > Levels: RasP1Hit.Tumor RasP.Tumor RasPNot.Tumor > > > > > > > > > > design<-model.matrix(~RasTum); > > > > colnames(design)<-sub("RasTum","",colnames(design)); > > > > colnames(design)[1] <- "Intercept" > > > > rownames(design)<-targets$Sample > > > > show(design[1:5,]) > > > Intercept RasP.Tumor RasPNot.Tumor > > > T4625_01A 1 1 0 > > > T4627_01A 1 1 0 > > > T4632_01A 1 1 0 > > > T2655_01A 1 1 0 > > > T2657_01A 1 1 0 > > > ... > > > > > > > con.matrix<-makeContrasts( > > > + > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > > > .Tumor=RasP.Tumor-RasPNot.Tumor, > > > + RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor, > > > + levels=design) > > > > > > > > > > show(con.matrix); > > > Contrasts > > > Levels RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > > > Intercept 0 0 > > > RasP.Tumor 1 1 > > > RasPNot.Tumor -1 -1 > > > Contrasts > > > Levels RasP1Hit.Tumor_RasPNot.Tumor > > > Intercept 0 > > > RasP.Tumor 0 > > > RasPNot.Tumor -1 > > > > > > > > > > > > > > > My first question is: as shown in show(RasTum), there are three subtype of > > > tumors I want to compare between: > > > > > > Levels: RasP1Hit.Tumor, RasP.Tumor and RasPNot.Tumor > > > > > > shown in show(design[1:5,]), everything is relative to RasP1Hit.Tumor. So in > > > my con.matrix<-makeContrasts command: > > > > > > RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor > > > > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > > > > > > > in particular, RasP1Hit.Tumor_RasPNot.Tumor=-RasPNot.Tumor since everything > > > is relative to RasP1Hit.Tumor. > > > > > > However, for contrast: RasPRasP1Hit.Tumor_RasPNot.Tumor, I intended to pool > > > both RasP class and RasP1Hit class together to compare to RasPNot class, I > > > end up with: > > > > > > RasPRasP1Hit.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > the right side is the same as > > > > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor-RasPNot.Tumor > > > > > > because everything is relative to RasP1Hit.Tumor. So what is the right way > > > to set up contrast RasPRasP1Hit.Tumor_RasPNot.Tumor in makeContrasts > > > command? kind of confused here? > > > > > > > > > > > > Because of the setting, the results are as below: > > > > > > > show(summary(de <- decideTests(fit))); > > > RasP.Tumor_RasPNot.Tumor RasPRasP1Hit.Tumor_RasPNot.Tumor > > > -1 0 0 > > > 0 17764 17764 > > > 1 0 0 > > > RasP1Hit.Tumor_RasPNot.Tumor > > > -1 0 > > > 0 17764 > > > 1 0 > > > > > > > > > > > > > > > > > > My 2nd question is: I did a similar design with including all normal > > > samples, since at the same model I want to get tumor vs normal contrasts as > > > well. > > > > > > > con.matrix<-makeContrasts( > > > + RasP.Tumor_RasP.Normal = > > > RasP.Tumor-RasP.Normal,RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal = > > > RasP.Tumor-RasP.Normal-RasP1Hit.Normal, > > > + RasPNot.Tumor_RasPNot.Normal=RasPNot.Tumor-RasPNot.Normal, > > > + > > > RasP.Tumor_RasPNot.Tumor=RasP.Tumor- RasPNot.Tumor,RasPRasP1Hit.Tumor_RasPNot > > > .Tumor=RasP.Tumor-RasPNot.Tumor, > > > + > > > RasP.Normal_RasPNot.Normal=RasP.Normal- RasPNot.Normal,RasPRasP1Hit.Normal_Ra > > > sPNot.Normal=RasP1Hit.Normal+RasP.Normal-RasPNot.Normal, > > > + Tumor_Normal=RasP.Tumor+ RasPNot.Tumor+-RasP.Normal- RasPNot.Normal, > > > + Tumor_Normal_2=RasP.Tumor+ > > > RasPNot.Tumor-RasP.Normal-RasPNot.Normal-RasP1Hit.Normal, > > > + levels=design) > > > > > > > > > > > > here everything is relative to RasP1Hit.Tumor. and the result as below: > > > > > > > > > > show(summary(de <- decideTests(fit))); > > > RasP.Tumor_RasP.Normal RasPRasP1Hit.Tumor_RasPRasP1Hit.Normal > > > -1 3947 4522 > > > 0 9912 8830 > > > 1 4274 4781 > > > RasPNot.Tumor_RasPNot.Normal RasP.Tumor_RasPNot.Tumor > > > -1 4809 24 > > > 0 8461 18102 > > > 1 4863 7 > > > RasPRasP1Hit.Tumor_RasPNot.Tumor RasP.Normal_RasPNot.Normal > > > -1 24 608 > > > 0 18102 17158 > > > 1 7 367 > > > RasPRasP1Hit.Normal_RasPNot.Normal Tumor_Normal Tumor_Normal_2 > > > -1 2001 5611 5580 > > > 0 13934 6728 6651 > > > 1 2198 5794 5902 > > > > > > > > > > > > RasP.Tumor_RasPNot.Tumor contrast has small number of DEGs (24+7) at FDR 5% > > > compared to using tumor data alone for the model shown earlier above, this > > > might be cuased by filtering as well and I am fine with that small > > > difference. What is bothering me is RasP.Normal_RasPNot.Normal contrast has > > > about 608+367=~1k DEGs, which is way more than RasP.Tumor_RasPNot.Tumor > > > contrast. I did plot plotMDS and see those normals are mixed together (I did > > > see two subgroups seems having batch effects, but for each batch (if there > > > is), there are RasP.Normal and RasPNot.Normal mixed there), all seem not > > > supposed to having 1k DEGs either. Plus, the tumor contrast > > > RasP.Tumor_RasPNot.Tumor has such small set of DEGs. Biologically, it might > > > be true (population, ethics group, batch effect etc), but I just want to > > > make sure what I did is reasonable. > > > > > > > > > > > > Any advice would be highly appreciated! > > > > > > Thanks a lot for your help! > > > > > > > > > > > > Best > > > > > > Ming > > > > > > > > > ______________________________________________________________________ > > > The information in this email is confidential and inte...{{dropped:9}} > > > > _______________________________________________ > > Bioconductor mailing list > > Bioconductor@r-project.org > > https://stat.ethz.ch/mailman/listinfo/bioconductor > > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor > > [[alternative HTML version deleted]] > > _______________________________________________ > Bioconductor mailing list > Bioconductor@r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:11}}
ADD REPLY •

Login before adding your answer.

Similar Posts
Loading Similar Posts
Traffic: 732 users visited in the last hour
Content Search
Users
Tags
Badges
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6