[BioC] edgeR, glmLRT error

[Gordon K Smyth]

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}}