[BioC] edgeR 3.0.6: Possible (small) bug in glmLRT when passing contrasts?
Gordon K Smyth
smyth at wehi.EDU.AU
Mon Dec 17 00:09:52 CET 2012
Dear Hayssam,
Yes, thanks, I managed to introduce this bug in edgeR 3.0.6 a week ago.
Ryan Thompson alerted me to the same bug a few days ago. It was fixed in
edgeR 3.0.7 last Friday. If you update edgeR, the bug should be gone.
Best wishes
Gordon
> Date: Sun, 16 Dec 2012 04:16:15 +0100
> From: <h.soueidan at nki.nl>
> To: <bioconductor at r-project.org>
> Subject: [BioC] edgeR 3.0.6: Possible (small) bug in glmLRT when
> passing contrasts?
> Message-ID: <B3C23AFB-76D7-4C9E-B98D-22EEF45D9EF6 at nki.nl>
> Content-Type: text/plain; charset="us-ascii"
>
> Dear list,
>
> I am analyzing a multi-factor time course RNA-Seq experiment with
> batch-effects using edgeR. After I perform a glmLRT, I cannot use
> topTags to display the most significant DE tags if I passed contrasts.
> Here's what I did:
>
> [....]
> # fit is the value from a previous glmFit, design is the value from a model.matrix call
>
> # 1) We first call by using a coef
>> lrt1<-glmLRT(fit,coef=10)
>> topTags(lrt1) #Work as expected
> [....]
>
> # 2) We call with a constrast
>> cr<-makeContrasts(time48="TreatA.Time48-TreatB.Time48",levels=design)
>> lrt2 <-glmLRT(fit,contrast=cr)
>> topTags(lrt2)
>
> Error in abs(object$table$logFC) :
> Non-numeric argument to mathematical function
>
> I don't know if this is the expected behavior (I would say no).
>
>
> I browsed a bit the edgeR source code. Here's what I believe might the culprit :
>
> In glmFit.R :
>
> [....]#In glmLRT <- fun...
>
> # line 161, when using contrasts
>
> coef <- 1:ncontrasts
> if(ncontrasts < ncol(contrast)) contrast <- contrast[,qrc$pivot[coef]]
> ** logFC <- (glmfit$coefficients %*% contrast)/log(2)
> if(ncontrasts>1) {
> [....]
>
> # Line 191
> tab <- data.frame(
> *** logFC=logFC,
> logCPM=(glmfit$abundance+log(1e6))/log(2),
> [....]
>
>
> logFC is matrix, inheriting its column name from the matrix product at ** (if only one contrast was passed, will be the name of the contrast). Then, when building the value data frame at ***, the name given in the LHS of the assignment is ignored (because logFC is a matrix). The resulting tab has no logFC column.
>
> One possible fix would be to cast logFC as a vector when there's only one contrasts.
>
> A simple work-around is to rename the first column of the returned value to "logFC", that makes topTags happy.
>
> I don't think this bug is data-dependant.
>
> Best wishes,
> Hayssam.
> PS: Dear Gordon Smyth, thanks a lot for the invaluable help provided on this list! And for edgeR :)
>
>
>
>> sessionInfo()
> R version 2.15.1 (2012-06-22)
> Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit)
>
> locale:
> [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
>
> attached base packages:
> [1] grid stats graphics grDevices utils datasets methods
> [8] base
>
> other attached packages:
> [1] gplots_2.11.0 MASS_7.3-22 KernSmooth_2.23-8
> [4] caTools_1.13 bitops_1.0-4.2 gtools_2.7.0
> [7] ggplot2_0.9.3 DESeq_1.10.1 lattice_0.20-10
> [10] locfit_1.5-8 Biobase_2.18.0 BiocGenerics_0.4.0
> [13] edgeR_3.0.6 limma_3.14.3 gdata_2.12.0
> [16] stringr_0.6.2 reshape2_1.2.2 data.table_1.8.6
>
>
>
> --
> Hayssam Soueidan, Ph.D
> Bioinformatics and Statistics
> The Netherlands Cancer Institute
> Plesmanlaan 121, 1066 CX Amsterdam
> The Netherlands
>
______________________________________________________________________
The information in this email is confidential and intend...{{dropped:4}}
More information about the Bioconductor
mailing list