edgeR, glmLRT error
1
0
Entering edit mode
@cei-abreu-goodger-4433
Last seen 9.7 years ago
Mexico
Hi Gordon, I've been getting the following error upon running glmLRT: Error in beta[k, ] <- betaj[decr, ] NAs are not allowed in subscripted assignments which has previously been reported, for instance: https://stat.ethz.ch/pipermail/bioconductor/2011-February/038079.html I subsetted my count data (originally ~14k rows) to narrow down the culprit row, and I would like to send you this example in case it can be used to improve the glm fitting code. If it turns out to be impossible to completely avoid these kinds of errors, perhaps it would be nice to return the genes that are generating the NAs so that we can more easily remove them and proceed? I can also avoid this particular case by using only the CommonDisp or heavily increasing the prior.n, but that kind of defeats the purpose of the TagwiseDisp. Example code and sessionInfo: ################ library(edgeR) counts <- read.table("http://datos.langebio.cinvestav.mx/~cei/counts.tab") grp <- sub("_R.+", "", colnames(counts)) dge <- DGEList(counts=counts, group=grp) dge <- calcNormFactors(dge) design <- model.matrix(~0 + dge$samples$group) colnames(design) <- levels(dge$samples$group) dge <- estimateGLMCommonDisp(dge, design) dge <- estimateGLMTagwiseDisp(dge, design) fit <- glmFit(dge, design) lrt <- glmLRT(dge, fit, contrast=c(0,0,0,-1,0.5,0.5)) Error in beta[k, ] <- betaj[decr, ] : NAs are not allowed in subscripted assignments ################ R version 2.15.0 (2012-03-30) Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit) locale: [1] C/UTF-8/C/C/C/C attached base packages: [1] grDevices datasets stats graphics utils methods base other attached packages: [1] edgeR_2.6.9 limma_3.12.1 Biobase_2.16.0 BiocGenerics_0.2.0 BiocInstaller_1.4.7 ################ Many thanks, Cei -- Dr. Cei Abreu-Goodger Langebio CINVESTAV Tel: (52) 462 166 3006 cei at langebio.cinvestav.mx -- This message has been scanned for viruses and dangerous content by MailScanner, and is believed to be clean.
• 1.6k views
ADD COMMENT •
0
Entering edit mode
@gordon-smyth
Last seen 1 hour ago
WEHI, Melbourne, Australia
Dear Cei, Thanks for the nice reproducible example. We have been working on a programming solution to eliminate this occasionally occuring error in edgeR. We have a nearly complete solution in our local (non-public) version of edgeR. Running your data example with our current local code, glmLRT no longer gives an error, but the 9th row of your data generates a NA value for the likelihood ratio test statistic. That's already an improvement, but I want to eliminate the NA fits as well if possible. We will push this code out to Bioc-devel soon, but probably not before the Bioconductor conference next week. Best wishes Gordon > Date: Sun, 15 Jul 2012 18:36:16 -0500 > From: Cei Abreu-Goodger <cei at="" langebio.cinvestav.mx=""> > To: bioconductor at r-project.org > Subject: [BioC] edgeR, glmLRT error > > Hi Gordon, > > I've been getting the following error upon running glmLRT: > > Error in beta[k, ] <- betaj[decr, ] > NAs are not allowed in subscripted assignments > > which has previously been reported, for instance: > > https://stat.ethz.ch/pipermail/bioconductor/2011-February/038079.html > > I subsetted my count data (originally ~14k rows) to narrow down the > culprit row, and I would like to send you this example in case it can be > used to improve the glm fitting code. If it turns out to be impossible > to completely avoid these kinds of errors, perhaps it would be nice to > return the genes that are generating the NAs so that we can more easily > remove them and proceed? > > I can also avoid this particular case by using only the CommonDisp or > heavily increasing the prior.n, but that kind of defeats the purpose of > the TagwiseDisp. > > Example code and sessionInfo: > > ################ > library(edgeR) > > counts <- read.table("http://datos.langebio.cinvestav.mx/~cei/counts.tab") > grp <- sub("_R.+", "", colnames(counts)) > dge <- DGEList(counts=counts, group=grp) > dge <- calcNormFactors(dge) > > design <- model.matrix(~0 + dge$samples$group) > colnames(design) <- levels(dge$samples$group) > > dge <- estimateGLMCommonDisp(dge, design) > dge <- estimateGLMTagwiseDisp(dge, design) > > fit <- glmFit(dge, design) > lrt <- glmLRT(dge, fit, contrast=c(0,0,0,-1,0.5,0.5)) > > Error in beta[k, ] <- betaj[decr, ] : > NAs are not allowed in subscripted assignments > > ################ > R version 2.15.0 (2012-03-30) > Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit) > > locale: > [1] C/UTF-8/C/C/C/C > > attached base packages: > [1] grDevices datasets stats graphics utils methods base > > other attached packages: > [1] edgeR_2.6.9 limma_3.12.1 Biobase_2.16.0 > BiocGenerics_0.2.0 BiocInstaller_1.4.7 > > ################ > > > Many thanks, > > Cei > > -- > Dr. Cei Abreu-Goodger > Langebio CINVESTAV > Tel: (52) 462 166 3006 > cei at langebio.cinvestav.mx > ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:4}}
ADD COMMENT •
0
Entering edit mode
Hi guys, Many thanks for the update, I'll check the devel in a few weeks. Best, Cei On 7/16/12 11:35 PM, Gordon K Smyth wrote: > Dear Cei, > > Thanks for the nice reproducible example. > > We have been working on a programming solution to eliminate this > occasionally occuring error in edgeR. We have a nearly complete > solution in our local (non-public) version of edgeR. Running your data > example with our current local code, glmLRT no longer gives an error, > but the 9th row of your data generates a NA value for the likelihood > ratio test statistic. > > That's already an improvement, but I want to eliminate the NA fits as > well if possible. > > We will push this code out to Bioc-devel soon, but probably not before > the Bioconductor conference next week. > > Best wishes > Gordon > >> Date: Sun, 15 Jul 2012 18:36:16 -0500 >> From: Cei Abreu-Goodger <cei at="" langebio.cinvestav.mx=""> >> To: bioconductor at r-project.org >> Subject: [BioC] edgeR, glmLRT error >> >> Hi Gordon, >> >> I've been getting the following error upon running glmLRT: >> >> Error in beta[k, ] <- betaj[decr, ] >> NAs are not allowed in subscripted assignments >> >> which has previously been reported, for instance: >> >> https://stat.ethz.ch/pipermail/bioconductor/2011-February/038079.html >> >> I subsetted my count data (originally ~14k rows) to narrow down the >> culprit row, and I would like to send you this example in case it can be >> used to improve the glm fitting code. If it turns out to be impossible >> to completely avoid these kinds of errors, perhaps it would be nice to >> return the genes that are generating the NAs so that we can more easily >> remove them and proceed? >> >> I can also avoid this particular case by using only the CommonDisp or >> heavily increasing the prior.n, but that kind of defeats the purpose of >> the TagwiseDisp. >> >> Example code and sessionInfo: >> >> ################ >> library(edgeR) >> >> counts <- >> read.table("http://datos.langebio.cinvestav.mx/~cei/counts.tab") >> grp <- sub("_R.+", "", colnames(counts)) >> dge <- DGEList(counts=counts, group=grp) >> dge <- calcNormFactors(dge) >> >> design <- model.matrix(~0 + dge$samples$group) >> colnames(design) <- levels(dge$samples$group) >> >> dge <- estimateGLMCommonDisp(dge, design) >> dge <- estimateGLMTagwiseDisp(dge, design) >> >> fit <- glmFit(dge, design) >> lrt <- glmLRT(dge, fit, contrast=c(0,0,0,-1,0.5,0.5)) >> >> Error in beta[k, ] <- betaj[decr, ] : >> NAs are not allowed in subscripted assignments >> >> ################ >> R version 2.15.0 (2012-03-30) >> Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit) >> >> locale: >> [1] C/UTF-8/C/C/C/C >> >> attached base packages: >> [1] grDevices datasets stats graphics utils methods base >> >> other attached packages: >> [1] edgeR_2.6.9 limma_3.12.1 Biobase_2.16.0 >> BiocGenerics_0.2.0 BiocInstaller_1.4.7 >> >> ################ >> >> >> Many thanks, >> >> Cei >> >> -- >> Dr. Cei Abreu-Goodger >> Langebio CINVESTAV >> Tel: (52) 462 166 3006 >> cei at langebio.cinvestav.mx >> > > ______________________________________________________________________ > The information in this email is confidential and inte...{{dropped:20}}
ADD REPLY •
0
Entering edit mode
Hi Cei, This is now fixed in edgeR 2.99.3 on Bioc-devel. Your example now runs with the results given below. Note that glmLRT now operates on the fit object directly. The data object is no longer required as an argument. Best wishes Gordon lrt <- glmLRT(fit, contrast=c(0,0,0,-1,0.5,0.5)) topTags(lrt) Output: Coefficient: -1*T2_S1 0.5*T2_S2 0.5*T2_S3 logFC logCPM LR PValue FDR CA00-XX-RX1-099-H06-EB.F 16.903 19.2 37.366 9.79e-10 1.08e-08 Contig1405 2.446 18.0 10.440 1.23e-03 6.78e-03 Cc_Contig7187 2.168 17.1 7.449 6.35e-03 2.33e-02 Contig2423 -1.240 16.3 4.252 3.92e-02 1.08e-01 CA00-XX-LP1-022-G05-EB.F -1.064 12.5 2.722 9.90e-02 2.18e-01 Contig10618 0.699 13.1 0.913 3.39e-01 6.22e-01 Cc_Contig2487 -0.506 15.8 0.658 4.17e-01 6.56e-01 Contig6686 -0.387 12.6 0.337 5.62e-01 7.72e-01 Contig15077 -0.279 15.8 0.198 6.56e-01 7.82e-01 Contig7207 -0.247 12.2 0.137 7.11e-01 7.82e-01 On Tue, 17 Jul 2012, Cei Abreu-Goodger wrote: > Hi guys, > > Many thanks for the update, I'll check the devel in a few weeks. > > Best, > > Cei > > On 7/16/12 11:35 PM, Gordon K Smyth wrote: >> Dear Cei, >> >> Thanks for the nice reproducible example. >> >> We have been working on a programming solution to eliminate this >> occasionally occuring error in edgeR. We have a nearly complete >> solution in our local (non-public) version of edgeR. Running your data >> example with our current local code, glmLRT no longer gives an error, >> but the 9th row of your data generates a NA value for the likelihood >> ratio test statistic. >> >> That's already an improvement, but I want to eliminate the NA fits as >> well if possible. >> >> We will push this code out to Bioc-devel soon, but probably not before >> the Bioconductor conference next week. >> >> Best wishes >> Gordon >> >>> Date: Sun, 15 Jul 2012 18:36:16 -0500 >>> From: Cei Abreu-Goodger <cei at="" langebio.cinvestav.mx=""> >>> To: bioconductor at r-project.org >>> Subject: [BioC] edgeR, glmLRT error >>> >>> Hi Gordon, >>> >>> I've been getting the following error upon running glmLRT: >>> >>> Error in beta[k, ] <- betaj[decr, ] >>> NAs are not allowed in subscripted assignments >>> >>> which has previously been reported, for instance: >>> >>> https://stat.ethz.ch/pipermail/bioconductor/2011-February/038079.html >>> >>> I subsetted my count data (originally ~14k rows) to narrow down the >>> culprit row, and I would like to send you this example in case it can be >>> used to improve the glm fitting code. If it turns out to be impossible >>> to completely avoid these kinds of errors, perhaps it would be nice to >>> return the genes that are generating the NAs so that we can more easily >>> remove them and proceed? >>> >>> I can also avoid this particular case by using only the CommonDisp or >>> heavily increasing the prior.n, but that kind of defeats the purpose of >>> the TagwiseDisp. >>> >>> Example code and sessionInfo: >>> >>> ################ >>> library(edgeR) >>> >>> counts <- >>> read.table("http://datos.langebio.cinvestav.mx/~cei/counts.tab") >>> grp <- sub("_R.+", "", colnames(counts)) >>> dge <- DGEList(counts=counts, group=grp) >>> dge <- calcNormFactors(dge) >>> >>> design <- model.matrix(~0 + dge$samples$group) >>> colnames(design) <- levels(dge$samples$group) >>> >>> dge <- estimateGLMCommonDisp(dge, design) >>> dge <- estimateGLMTagwiseDisp(dge, design) >>> >>> fit <- glmFit(dge, design) >>> lrt <- glmLRT(dge, fit, contrast=c(0,0,0,-1,0.5,0.5)) >>> >>> Error in beta[k, ] <- betaj[decr, ] : >>> NAs are not allowed in subscripted assignments >>> >>> ################ >>> R version 2.15.0 (2012-03-30) >>> Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit) >>> >>> locale: >>> [1] C/UTF-8/C/C/C/C >>> >>> attached base packages: >>> [1] grDevices datasets stats graphics utils methods base >>> >>> other attached packages: >>> [1] edgeR_2.6.9 limma_3.12.1 Biobase_2.16.0 >>> BiocGenerics_0.2.0 BiocInstaller_1.4.7 >>> >>> ################ >>> >>> >>> Many thanks, >>> >>> Cei >>> >>> -- >>> Dr. Cei Abreu-Goodger >>> Langebio CINVESTAV >>> Tel: (52) 462 166 3006 >>> cei at langebio.cinvestav.mx ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:4}}
ADD REPLY •
0
Entering edit mode
Dear Gordon, I've tried out the new code on a larger, also problematic, dataset and it works without any problems. Many thanks, Cei On 7/27/12 4:31 AM, Gordon K Smyth wrote: > Hi Cei, > > This is now fixed in edgeR 2.99.3 on Bioc-devel. Your example now runs > with the results given below. Note that glmLRT now operates on the fit > object directly. The data object is no longer required as an argument. > > Best wishes > Gordon > > lrt <- glmLRT(fit, contrast=c(0,0,0,-1,0.5,0.5)) > topTags(lrt) > > Output: > Coefficient: -1*T2_S1 0.5*T2_S2 0.5*T2_S3 > logFC logCPM LR PValue FDR > CA00-XX-RX1-099-H06-EB.F 16.903 19.2 37.366 9.79e-10 1.08e-08 > Contig1405 2.446 18.0 10.440 1.23e-03 6.78e-03 > Cc_Contig7187 2.168 17.1 7.449 6.35e-03 2.33e-02 > Contig2423 -1.240 16.3 4.252 3.92e-02 1.08e-01 > CA00-XX-LP1-022-G05-EB.F -1.064 12.5 2.722 9.90e-02 2.18e-01 > Contig10618 0.699 13.1 0.913 3.39e-01 6.22e-01 > Cc_Contig2487 -0.506 15.8 0.658 4.17e-01 6.56e-01 > Contig6686 -0.387 12.6 0.337 5.62e-01 7.72e-01 > Contig15077 -0.279 15.8 0.198 6.56e-01 7.82e-01 > Contig7207 -0.247 12.2 0.137 7.11e-01 7.82e-01 > > > On Tue, 17 Jul 2012, Cei Abreu-Goodger wrote: > >> Hi guys, >> >> Many thanks for the update, I'll check the devel in a few weeks. >> >> Best, >> >> Cei >> >> On 7/16/12 11:35 PM, Gordon K Smyth wrote: >>> Dear Cei, >>> >>> Thanks for the nice reproducible example. >>> >>> We have been working on a programming solution to eliminate this >>> occasionally occuring error in edgeR. We have a nearly complete >>> solution in our local (non-public) version of edgeR. Running your data >>> example with our current local code, glmLRT no longer gives an error, >>> but the 9th row of your data generates a NA value for the likelihood >>> ratio test statistic. >>> >>> That's already an improvement, but I want to eliminate the NA fits as >>> well if possible. >>> >>> We will push this code out to Bioc-devel soon, but probably not before >>> the Bioconductor conference next week. >>> >>> Best wishes >>> Gordon >>> >>>> Date: Sun, 15 Jul 2012 18:36:16 -0500 >>>> From: Cei Abreu-Goodger <cei at="" langebio.cinvestav.mx=""> >>>> To: bioconductor at r-project.org >>>> Subject: [BioC] edgeR, glmLRT error >>>> >>>> Hi Gordon, >>>> >>>> I've been getting the following error upon running glmLRT: >>>> >>>> Error in beta[k, ] <- betaj[decr, ] >>>> NAs are not allowed in subscripted assignments >>>> >>>> which has previously been reported, for instance: >>>> >>>> https://stat.ethz.ch/pipermail/bioconductor/2011-February/038079.html >>>> >>>> I subsetted my count data (originally ~14k rows) to narrow down the >>>> culprit row, and I would like to send you this example in case it >>>> can be >>>> used to improve the glm fitting code. If it turns out to be impossible >>>> to completely avoid these kinds of errors, perhaps it would be nice to >>>> return the genes that are generating the NAs so that we can more easily >>>> remove them and proceed? >>>> >>>> I can also avoid this particular case by using only the CommonDisp or >>>> heavily increasing the prior.n, but that kind of defeats the purpose of >>>> the TagwiseDisp. >>>> >>>> Example code and sessionInfo: >>>> >>>> ################ >>>> library(edgeR) >>>> >>>> counts <- >>>> read.table("http://datos.langebio.cinvestav.mx/~cei/counts.tab") >>>> grp <- sub("_R.+", "", colnames(counts)) >>>> dge <- DGEList(counts=counts, group=grp) >>>> dge <- calcNormFactors(dge) >>>> >>>> design <- model.matrix(~0 + dge$samples$group) >>>> colnames(design) <- levels(dge$samples$group) >>>> >>>> dge <- estimateGLMCommonDisp(dge, design) >>>> dge <- estimateGLMTagwiseDisp(dge, design) >>>> >>>> fit <- glmFit(dge, design) >>>> lrt <- glmLRT(dge, fit, contrast=c(0,0,0,-1,0.5,0.5)) >>>> >>>> Error in beta[k, ] <- betaj[decr, ] : >>>> NAs are not allowed in subscripted assignments >>>> >>>> ################ >>>> R version 2.15.0 (2012-03-30) >>>> Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit) >>>> >>>> locale: >>>> [1] C/UTF-8/C/C/C/C >>>> >>>> attached base packages: >>>> [1] grDevices datasets stats graphics utils methods base >>>> >>>> other attached packages: >>>> [1] edgeR_2.6.9 limma_3.12.1 Biobase_2.16.0 >>>> BiocGenerics_0.2.0 BiocInstaller_1.4.7 >>>> >>>> ################ >>>> >>>> >>>> Many thanks, >>>> >>>> Cei >>>> >>>> -- >>>> Dr. Cei Abreu-Goodger >>>> Langebio CINVESTAV >>>> Tel: (52) 462 166 3006 >>>> cei at langebio.cinvestav.mx > > ______________________________________________________________________ > The information in this email is confidential and inte...{{dropped:18}}
ADD REPLY •

Login before adding your answer.

Similar Posts
Loading Similar Posts
Traffic: 785 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