[BioC] Limma: topTable without eBayes?
Jenny Drnevich
drnevich at illinois.edu
Tue Aug 31 23:36:16 CEST 2010
Hi Guido,
topTable for a single coef needs the $t, $p.value and $lods that are
created from eBayes. You could manually compute these to be regular
t-statistics instead of those from eBayes. See section 10 of the
limmaUsersGuide(). Try this (I didn't test it):
fit2$t <- fit2$coef/fit2$stdev.unscaled/fit2$sigma
fit2$p.value <- 2 * pt(-abs(fit2$t), df = fit2$df.residual)
The $lods is harder; if you care about it you'll have to figure out
how to duplicate it from the ebayes code by using modified inputs,
which is what I did for the $t and $p-value; if you don't care about
it but need something for topTable to work, I'd just do:
fit2$lods <- fit3$lods
HTH,
Jenny
At 12:51 PM 8/31/2010, Hooiveld, Guido wrote:
>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
>
>
>
> [[alternative HTML version deleted]]
>
>_______________________________________________
>Bioconductor mailing list
>Bioconductor at stat.math.ethz.ch
>https://stat.ethz.ch/mailman/listinfo/bioconductor
>Search the archives:
>http://news.gmane.org/gmane.science.biology.informatics.conductor
More information about the Bioconductor
mailing list