Limma: topTable without eBayes?
1
0
Entering edit mode
@gordon-smyth
Last seen 1 hour ago
WEHI, Melbourne, Australia
Why would you need or want to do that? eBayes doesn't destroy any information, so why do you need to avoid it? Gordon > Date: Tue, 31 Aug 2010 19:51:55 +0200 > From: "Hooiveld, Guido" <guido.hooiveld at="" wur.nl=""> > To: "bioconductor at stat.math.ethz.ch" <bioconductor at="" stat.math.ethz.ch=""> > Subject: [BioC] Limma: topTable without eBayes? > Content-Type: text/plain > > Dear list, > > topTable is a convenient function for displaying the top regulated genes > in an experiment. topTable is usually applied on a object that is > generated after fitting a linear model AND subsequent application of the > eBayes funtion (fit3 below). > I would like to view the results of a model WITHOUT applying the eBayes > function, preferably using topTable. That is, using topTable to > summerize the results of fit2 below. Applying topTable as such won't > work; apparently a format imposed by the eBayes function is required to > have topTable properly work. Even specifying a contrast of interest > doesn't work. > Therefore, has anyone a suggestion how to summerize the results of the > fit2 mentioned below? > > Thanks, > Guido > >> library(affy) >> library(limma) >> targets <- readTargets("targets_A23.txt") >> >> affy.data <- ReadAffy(filenames=targets$FileName) >> x.norm <- rma(affy.data) > Background correcting > Normalizing > Calculating Expression >> >> TS <- paste(targets$Strain, targets$Treatment, sep=".") >> TS <- factor(TS, levels=c("WT.Con","WT.WY","KO.Con","KO.WY")) >> design <- model.matrix(~0+TS) >> colnames(design) <- levels(TS) >> fit <- lmFit(x.norm, design) >> cont.matrix <- makeContrasts(A=WT.WY-WT.Con, B=KO.WY-KO.Con, levels=design) >> >> fit2 <- contrasts.fit(fit, cont.matrix) #<<---- how to conviently summerize results of this fit? >> >> fit3 <- eBayes(fit2) >> >> topTable(fit3) > ID A B AveExpr F P.Value adj.P.Val > 13289 1431833_a_at 3.538370 -0.091190845 7.856618 1453.7786 5.084771e-18 1.153734e-13 > 18706 1450643_s_at 2.847376 0.182291682 7.907519 1234.5101 1.739465e-17 1.620566e-13 > 17177 1449065_at 4.528581 0.125074564 6.713129 1200.7506 2.142661e-17 1.620566e-13 > 7120 1422925_s_at 3.149781 -0.007867198 8.072523 895.1004 1.945926e-16 1.103827e-12 > 7192 1422997_s_at 4.298410 0.090012125 8.311977 837.9040 3.193016e-16 1.448991e-12 > 13047 1431012_a_at 2.682180 0.228653000 8.911287 785.9448 5.159295e-16 1.698832e-12 > 6721 1422526_at 2.772829 0.123117132 8.060643 784.2985 5.240997e-16 1.698832e-12 > 16494 1448382_at 3.386068 0.193350865 9.658716 685.0981 1.442685e-15 3.703026e-12 > 296 1415965_at 3.749143 0.071163645 5.293939 683.4572 1.468807e-15 3.703026e-12 > 8911 1424716_at 2.647513 0.370083995 8.453969 665.2223 1.798165e-15 4.080035e-12 >> topTable(fit2) > Error in topTableF(fit, number = number, genelist = genelist, adjust.method = adjust.method, : > F-statistics not available. Try topTable for individual coef instead. >> topTable(fit2, coef=1) > Error in dim(data) <- dim : attempt to set an attribute on NULL >> >> sessionInfo() > R version 2.11.0 (2010-04-22) > i386-pc-mingw32 > > locale: > [1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252 LC_MONETARY=English_United States.1252 LC_NUMERIC=C LC_TIME=English_United States.1252 > > attached base packages: > [1] stats graphics grDevices utils datasets methods base > > other attached packages: > [1] moe430acdf_2.6.0 limma_3.4.0 affy_1.26.1 Biobase_2.8.0 > > loaded via a namespace (and not attached): > [1] affyio_1.16.0 preprocessCore_1.10.0 tools_2.11.0 >> > > > ------------------------------------------------ > Guido Hooiveld, PhD > Nutrition, Metabolism & Genomics Group > Division of Human Nutrition > Wageningen University > Biotechnion, Bomenweg 2 > NL-6703 HD Wageningen > the Netherlands > tel: (+)31 317 485788 > fax: (+)31 317 483342 > internet: http://nutrigene.4t.com<http: nutrigene.4t.com=""/> > email: guido.hooiveld at wur.nl ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:4}}
affy affy • 1.2k views
ADD COMMENT
0
Entering edit mode
Guido Hooiveld ★ 4.1k
@guido-hooiveld-2020
Last seen 1 day ago
Wageningen University, Wageningen, the …
Hi Gordon, Basically because I would like to retrieve the un-moderated (thus before eBayes) t-statistics. Two reasons for me to want this (a conceptual and a more practical one): I am evaluating the performance of a method that takes into account the correlation of the raw (expression) data to identify and correct for hidden variables when fitting the model (SVA methodology). AFAIK the eBayes function also takes into account the correlation structure of the data. I think it would be correct to avoid taking 'advantage' of this correlation structure twice, hence the attempt to use topTable to save the results of the fits before applying the eBayes function. Since I always use my standard lines of code that uses limma to identify DEGs, for me it would be most convenient to stick to my script (thus limma) as much as possible, hence my question. Best, Guido > -----Original Message----- > From: Gordon K Smyth [mailto:smyth at wehi.EDU.AU] > Sent: Wednesday, September 01, 2010 12:16 > To: Hooiveld, Guido > Cc: Bioconductor mailing list > Subject: [BioC] Limma: topTable without eBayes? > > Why would you need or want to do that? eBayes doesn't > destroy any information, so why do you need to avoid it? > > Gordon > > > Date: Tue, 31 Aug 2010 19:51:55 +0200 > > From: "Hooiveld, Guido" <guido.hooiveld at="" wur.nl=""> > > To: "bioconductor at stat.math.ethz.ch" > <bioconductor at="" stat.math.ethz.ch=""> > > Subject: [BioC] Limma: topTable without eBayes? > > Content-Type: text/plain > > > > Dear list, > > > > topTable is a convenient function for displaying the top regulated > > genes in an experiment. topTable is usually applied on a > object that > > is generated after fitting a linear model AND subsequent > application > > of the eBayes funtion (fit3 below). > > > I would like to view the results of a model WITHOUT applying the > > eBayes function, preferably using topTable. That is, using > topTable to > > summerize the results of fit2 below. Applying topTable as > such won't > > work; apparently a format imposed by the eBayes function is > required > > to have topTable properly work. Even specifying a contrast > of interest > > doesn't work. > > > Therefore, has anyone a suggestion how to summerize the > results of the > > fit2 mentioned below? > > > > Thanks, > > Guido > > > >> library(affy) > >> library(limma) > >> targets <- readTargets("targets_A23.txt") > >> > >> affy.data <- ReadAffy(filenames=targets$FileName) > >> x.norm <- rma(affy.data) > > Background correcting > > Normalizing > > Calculating Expression > >> > >> TS <- paste(targets$Strain, targets$Treatment, sep=".") TS <- > >> factor(TS, levels=c("WT.Con","WT.WY","KO.Con","KO.WY")) > >> design <- model.matrix(~0+TS) > >> colnames(design) <- levels(TS) > >> fit <- lmFit(x.norm, design) > >> cont.matrix <- makeContrasts(A=WT.WY-WT.Con, B=KO.WY-KO.Con, > >> levels=design) > >> > >> fit2 <- contrasts.fit(fit, cont.matrix) #<<---- how to > conviently summerize results of this fit? > >> > >> fit3 <- eBayes(fit2) > >> > >> topTable(fit3) > > ID A B AveExpr F > P.Value adj.P.Val > > 13289 1431833_a_at 3.538370 -0.091190845 7.856618 1453.7786 > > 5.084771e-18 1.153734e-13 > > 18706 1450643_s_at 2.847376 0.182291682 7.907519 1234.5101 > 1.739465e-17 1.620566e-13 > > 17177 1449065_at 4.528581 0.125074564 6.713129 1200.7506 > 2.142661e-17 1.620566e-13 > > 7120 1422925_s_at 3.149781 -0.007867198 8.072523 895.1004 > > 1.945926e-16 1.103827e-12 > > 7192 1422997_s_at 4.298410 0.090012125 8.311977 837.9040 > > 3.193016e-16 1.448991e-12 > > 13047 1431012_a_at 2.682180 0.228653000 8.911287 785.9448 > 5.159295e-16 1.698832e-12 > > 6721 1422526_at 2.772829 0.123117132 8.060643 784.2985 > 5.240997e-16 1.698832e-12 > > 16494 1448382_at 3.386068 0.193350865 9.658716 685.0981 > 1.442685e-15 3.703026e-12 > > 296 1415965_at 3.749143 0.071163645 5.293939 683.4572 > 1.468807e-15 3.703026e-12 > > 8911 1424716_at 2.647513 0.370083995 8.453969 665.2223 > 1.798165e-15 4.080035e-12 > >> topTable(fit2) > > Error in topTableF(fit, number = number, genelist = > genelist, adjust.method = adjust.method, : > > F-statistics not available. Try topTable for individual > coef instead. > >> topTable(fit2, coef=1) > > Error in dim(data) <- dim : attempt to set an attribute on NULL > >> > >> sessionInfo() > > R version 2.11.0 (2010-04-22) > > i386-pc-mingw32 > > > > locale: > > [1] LC_COLLATE=English_United States.1252 > LC_CTYPE=English_United States.1252 > LC_MONETARY=English_United States.1252 LC_NUMERIC=C > LC_TIME=English_United States.1252 > > > > attached base packages: > > [1] stats graphics grDevices utils datasets methods base > > > > other attached packages: > > [1] moe430acdf_2.6.0 limma_3.4.0 affy_1.26.1 Biobase_2.8.0 > > > > loaded via a namespace (and not attached): > > [1] affyio_1.16.0 preprocessCore_1.10.0 tools_2.11.0 > >> > > > > > > ------------------------------------------------ > > Guido Hooiveld, PhD > > Nutrition, Metabolism & Genomics Group Division of Human Nutrition > > Wageningen University Biotechnion, Bomenweg 2 > > NL-6703 HD Wageningen > > the Netherlands > > tel: (+)31 317 485788 > > fax: (+)31 317 483342 > > internet: http://nutrigene.4t.com<http: nutrigene.4t.com=""/> > > email: guido.hooiveld at wur.nl > > ______________________________________________________________________ > The information in this email is confidential and intended > solely for the addressee. > You must not disclose, forward, print or use it without the > permission of the sender. > ______________________________________________________________________ > >
ADD COMMENT
0
Entering edit mode
Dear Guido, eBayes doesn't make any use of correlation structure, so I still don't understand why it could cause you any problems. If you still want the ordinary t-statistics however, see page 54 of the limma User's Guide. There is no toptable code in limma that uses the ordinary t-statistic. I view ordinary t as a dangerous statistic to be discouraged, so I'm not willing to provide code based on it. Even simply ranking by fold change would almost always be safer and better, let alone using more sophisticated things like moderated t. Best wishes Gordon On Wed, 1 Sep 2010, Hooiveld, Guido wrote: > Hi Gordon, > > Basically because I would like to retrieve the un-moderated (thus before eBayes) t-statistics. > > Two reasons for me to want this (a conceptual and a more practical one): > I am evaluating the performance of a method that takes into account the > correlation of the raw (expression) data to identify and correct for > hidden variables when fitting the model (SVA methodology). AFAIK the > eBayes function also takes into account the correlation structure of the > data. I think it would be correct to avoid taking 'advantage' of this > correlation structure twice, hence the attempt to use topTable to save > the results of the fits before applying the eBayes function. > Since I always use my standard lines of code that uses limma to identify > DEGs, for me it would be most convenient to stick to my script (thus > limma) as much as possible, hence my question. > > Best, > Guido > > > >> -----Original Message----- >> From: Gordon K Smyth [mailto:smyth at wehi.EDU.AU] >> Sent: Wednesday, September 01, 2010 12:16 >> To: Hooiveld, Guido >> Cc: Bioconductor mailing list >> Subject: [BioC] Limma: topTable without eBayes? >> >> Why would you need or want to do that? eBayes doesn't >> destroy any information, so why do you need to avoid it? >> >> Gordon >> >>> Date: Tue, 31 Aug 2010 19:51:55 +0200 >>> From: "Hooiveld, Guido" <guido.hooiveld at="" wur.nl=""> >>> To: "bioconductor at stat.math.ethz.ch" >> <bioconductor at="" stat.math.ethz.ch=""> >>> Subject: [BioC] Limma: topTable without eBayes? >>> Content-Type: text/plain >>> >>> Dear list, >>> >>> topTable is a convenient function for displaying the top regulated >>> genes in an experiment. topTable is usually applied on a >> object that >>> is generated after fitting a linear model AND subsequent >> application >>> of the eBayes funtion (fit3 below). >> >>> I would like to view the results of a model WITHOUT applying the >>> eBayes function, preferably using topTable. That is, using >> topTable to >>> summerize the results of fit2 below. Applying topTable as >> such won't >>> work; apparently a format imposed by the eBayes function is >> required >>> to have topTable properly work. Even specifying a contrast >> of interest >>> doesn't work. >> >>> Therefore, has anyone a suggestion how to summerize the >> results of the >>> fit2 mentioned below? >>> >>> Thanks, >>> Guido >>> >>>> library(affy) >>>> library(limma) >>>> targets <- readTargets("targets_A23.txt") >>>> >>>> affy.data <- ReadAffy(filenames=targets$FileName) >>>> x.norm <- rma(affy.data) >>> Background correcting >>> Normalizing >>> Calculating Expression >>>> >>>> TS <- paste(targets$Strain, targets$Treatment, sep=".") TS <- >>>> factor(TS, levels=c("WT.Con","WT.WY","KO.Con","KO.WY")) >>>> design <- model.matrix(~0+TS) >>>> colnames(design) <- levels(TS) >>>> fit <- lmFit(x.norm, design) >>>> cont.matrix <- makeContrasts(A=WT.WY-WT.Con, B=KO.WY-KO.Con, >>>> levels=design) >>>> >>>> fit2 <- contrasts.fit(fit, cont.matrix) #<<---- how to >> conviently summerize results of this fit? >>>> >>>> fit3 <- eBayes(fit2) >>>> >>>> topTable(fit3) >>> ID A B AveExpr F >> P.Value adj.P.Val >>> 13289 1431833_a_at 3.538370 -0.091190845 7.856618 1453.7786 >>> 5.084771e-18 1.153734e-13 >>> 18706 1450643_s_at 2.847376 0.182291682 7.907519 1234.5101 >> 1.739465e-17 1.620566e-13 >>> 17177 1449065_at 4.528581 0.125074564 6.713129 1200.7506 >> 2.142661e-17 1.620566e-13 >>> 7120 1422925_s_at 3.149781 -0.007867198 8.072523 895.1004 >>> 1.945926e-16 1.103827e-12 >>> 7192 1422997_s_at 4.298410 0.090012125 8.311977 837.9040 >>> 3.193016e-16 1.448991e-12 >>> 13047 1431012_a_at 2.682180 0.228653000 8.911287 785.9448 >> 5.159295e-16 1.698832e-12 >>> 6721 1422526_at 2.772829 0.123117132 8.060643 784.2985 >> 5.240997e-16 1.698832e-12 >>> 16494 1448382_at 3.386068 0.193350865 9.658716 685.0981 >> 1.442685e-15 3.703026e-12 >>> 296 1415965_at 3.749143 0.071163645 5.293939 683.4572 >> 1.468807e-15 3.703026e-12 >>> 8911 1424716_at 2.647513 0.370083995 8.453969 665.2223 >> 1.798165e-15 4.080035e-12 >>>> topTable(fit2) >>> Error in topTableF(fit, number = number, genelist = >> genelist, adjust.method = adjust.method, : >>> F-statistics not available. Try topTable for individual >> coef instead. >>>> topTable(fit2, coef=1) >>> Error in dim(data) <- dim : attempt to set an attribute on NULL >>>> >>>> sessionInfo() >>> R version 2.11.0 (2010-04-22) >>> i386-pc-mingw32 >>> >>> locale: >>> [1] LC_COLLATE=English_United States.1252 >> LC_CTYPE=English_United States.1252 >> LC_MONETARY=English_United States.1252 LC_NUMERIC=C >> LC_TIME=English_United States.1252 >>> >>> attached base packages: >>> [1] stats graphics grDevices utils datasets methods base >>> >>> other attached packages: >>> [1] moe430acdf_2.6.0 limma_3.4.0 affy_1.26.1 Biobase_2.8.0 >>> >>> loaded via a namespace (and not attached): >>> [1] affyio_1.16.0 preprocessCore_1.10.0 tools_2.11.0 >>>> >>> >>> >>> ------------------------------------------------ >>> Guido Hooiveld, PhD >>> Nutrition, Metabolism & Genomics Group Division of Human Nutrition >>> Wageningen University Biotechnion, Bomenweg 2 >>> NL-6703 HD Wageningen >>> the Netherlands >>> tel: (+)31 317 485788 >>> fax: (+)31 317 483342 >>> internet: http://nutrigene.4t.com<http: nutrigene.4t.com=""/> >>> email: guido.hooiveld at wur.nl >> >> ______________________________________________________________________ >> The information in this email is confidential and intended >> solely for the addressee. >> You must not disclose, forward, print or use it without the >> permission of the sender. >> ______________________________________________________________________ >> >> > ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:4}}
ADD REPLY

Login before adding your answer.

Traffic: 829 users visited in the last hour
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