Apparently randomly ocurring error in DESeq2: in call to fitDisp, the following arguments contain NA: log_alphaSEXP, log_alpha_prior_meanSEXP
2
0
Entering edit mode
@peter-langfelder-4469
Last seen 27 days ago
United States

When running DESeq2 1.8.1, I get the following error:

in call to fitDisp, the following arguments contain NA: log_alphaSEXP, log_alpha_prior_meanSEXP

This happens when calling DESeq and varianceStabilizingTransformation.

The perplexing part is that this error is intermittent and apparently occurs at random. Running the exact same code on exact same input over and over again eventually leads to successful finish, as illustrated in the transcript below. Does some part of the calculation use a random number generator to generate starting values (although it would be great if these could always be the same), or is it a bug within compiled code (this apparently random  behaviour could be the result of using uninitialized variables)?

 

Here's a transcript in which I run DESeq on a list of 3 DESeqDataSet objects. This is not a reproducible example (the data are not defined here), it is rather just an illustration of the problem. If necessary, I can create an R session image with just the necessary information and post is somewhere for download.

 

> lapply(ds.bicov.all.1, class)
$Cortex.vs.Liver
[1] "DESeqDataSet"
attr(,"package")
[1] "DESeq2"

$Striatum.vs.Liver
[1] "DESeqDataSet"
attr(,"package")
[1] "DESeq2"

$Striatum.vs.Cortex
[1] "DESeqDataSet"
attr(,"package")
[1] "DESeq2"


# Running DESeq on the list

# Attempt 1: Error in second analysis

> deseq.all.wald.1 = lapply(ds.bicov.all.1, DESeq, minReplicatesForReplace = Inf)
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
-- standard model matrices are used for factors with two levels and an interaction,
   where the main effects are for the reference level of other factors.
   see the 'Interactions' section of the vignette for more details: vignette('DESeq2')
estimating size factors
estimating dispersions
gene-wise dispersion estimates
Error in fitBetaWrapper(ySEXP = counts(object), xSEXP = modelMatrix, nfSEXP = normalizationFactors,  :
  in call to fitBeta, the following arguments contain NA: alpha_hatSEXP
In addition: There were 50 or more warnings (use warnings() to see the first 50)

 

# Attempt 2: Error in first analysis


> deseq.all.wald.1 = lapply(ds.bicov.all.1, DESeq, minReplicatesForReplace = Inf)
estimating size factors
estimating dispersions
gene-wise dispersion estimates
Error in fitBetaWrapper(ySEXP = counts(object), xSEXP = modelMatrix, nfSEXP = normalizationFactors,  :
  in call to fitBeta, the following arguments contain NA: alpha_hatSEXP

# Attempt 3: Success


> deseq.all.wald.1 = lapply(ds.bicov.all.1, DESeq, minReplicatesForReplace = Inf)
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
-- standard model matrices are used for factors with two levels and an interaction,
   where the main effects are for the reference level of other factors.
   see the 'Interactions' section of the vignette for more details: vignette('DESeq2')
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
-- standard model matrices are used for factors with two levels and an interaction,
   where the main effects are for the reference level of other factors.
   see the 'Interactions' section of the vignette for more details: vignette('DESeq2')
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
-- standard model matrices are used for factors with two levels and an interaction,
   where the main effects are for the reference level of other factors.
   see the 'Interactions' section of the vignette for more details: vignette('DESeq2')

>

sessionInfo()
R version 3.2.2 beta (2015-08-05 r68859)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Fedora 20 (Heisenbug)

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
 [1] stats4    parallel  grid      stats     graphics  grDevices utils    
 [8] datasets  methods   base     

other attached packages:
 [1] DESeq2_1.8.1            RcppArmadillo_0.5.300.4 Rcpp_0.12.0            
 [4] GenomicRanges_1.20.5    anRichment_0.79-997     GO.db_3.1.2            
 [7] AnnotationDbi_1.30.1    GenomeInfoDb_1.4.1      IRanges_2.2.6          
[10] S4Vectors_0.6.3         Biobase_2.28.0          BiocGenerics_0.14.0    
[13] WGCNA_1.47-3            fastcluster_1.1.16      dynamicTreeCut_1.62    
[16] preprocessCore_1.30.0   sva_3.14.0              genefilter_1.50.0      
[19] mgcv_1.8-7              nlme_3.1-121            sqldf_0.4-10           
[22] RSQLite_1.0.0           DBI_0.3.1               gsubfn_0.6-6           
[25] proto_0.3-10            matrixStats_0.14.2      doParallel_1.0.8       
[28] iterators_1.0.7         foreach_1.4.2           reshape_0.8.5          
[31] Hmisc_3.16-0            ggplot2_1.0.1           Formula_1.2-1          
[34] survival_2.38-3         lattice_0.20-33         impute_1.42.0          
[37] cluster_2.0.3           class_7.3-13            MASS_7.3-43            

loaded via a namespace (and not attached):
 [1] locfit_1.5-9.1       digest_0.6.8         plyr_1.8.3          
 [4] chron_2.3-47         futile.options_1.0.0 acepack_1.3-3.3     
 [7] annotate_1.46.1      rpart_4.1-10         Matrix_1.2-2        
[10] splines_3.2.2        BiocParallel_1.2.9   geneplotter_1.46.0  
[13] stringr_1.0.0        foreign_0.8-65       munsell_0.4.2       
[16] tcltk_3.2.2          nnet_7.3-10          gridExtra_2.0.0     
[19] codetools_0.2-14     XML_3.98-1.3         xtable_1.7-4        
[22] gtable_0.1.2         magrittr_1.5         scales_0.2.5        
[25] stringi_0.5-5        XVector_0.8.0        reshape2_1.4.1      
[28] latticeExtra_0.6-26  futile.logger_1.4.1  lambda.r_1.1.7      
[31] RColorBrewer_1.1-2   tools_3.2.2          colorspace_1.2-6  

Thanks,

Peter

deseq2 • 2.3k views
ADD COMMENT
0
Entering edit mode
@mikelove
Last seen 1 hour ago
United States

hi Peter,

thanks for the report. I guess the easiest way would be if you could send a small example and code which reproduces the bug to:

maintainer("DESeq2")

There is a random number generator which is used for a procedure when the residual degrees of freedom is <= 3. But this is done in a deterministic way, by internally setting a seed, and then restoring the earlier seed when the function is finished. So this should not change through repeated calls.

ADD COMMENT
0
Entering edit mode
@mikelove
Last seen 1 hour ago
United States

hi Peter,

thanks for the code and example.

I haven't been able to raise an error on my end. I ran it through the loop 20 times.

I noticed that all the design variables are coded as numeric. Could you try switching these to factors to see if you still have the error?

Also I'd recommend setting betaPrior=FALSE because the design has an interaction term (this will be set by default in the next release).

ADD COMMENT

Login before adding your answer.

Traffic: 567 users visited in the last hour
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