[BioC] edgeR, glmLRT error

[Cei Abreu-Goodger]

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