optimHess problem in apeglm/DESeq2
1
0
Entering edit mode
ATpoint ★ 4.6k
@atpoint-13662
Last seen 4 hours ago
Germany

Hi Mike, based on https://www.biostars.org/p/9559740/#9559740 I made this MRE which errors in R 4.3 - Bioc 3.17 and latest DESeq2. I just checked that it runs fine in the Bioc 3.16 Docker container but errors in the 3.17 one. Hope that helps tracking it down.

suppressMessages(library(DESeq2))
set.seed(1)
dds <- makeExampleDESeqDataSet(n=50000, m=20)
dds <- DESeq(dds, quiet=TRUE)
res <- lfcShrink(dds, coef=2, type="apeglm", quiet=TRUE)
#> Warning in nbinomGLM(x = x, Y = YNZ, size = size, weights = weightsNZ, offset =
#> offsetNZ, : the line search routine failed, possibly due to insufficient
#> numeric precision
#> Error in optimHess(par = init, fn = nbinomFn, gr = nbinomGr, x = x, y = y, : nicht endlicher Wert von optim angegeben
sessionInfo()
#> R version 4.3.0 (2023-04-21 ucrt)
#> Platform: x86_64-w64-mingw32/x64 (64-bit)
#> Running under: Windows 11 x64 (build 22000)
#> 
#> Matrix products: default
#> 
#> 
#> locale:
#> [1] LC_COLLATE=German_Germany.utf8  LC_CTYPE=German_Germany.utf8   
#> [3] LC_MONETARY=German_Germany.utf8 LC_NUMERIC=C                   
#> [5] LC_TIME=German_Germany.utf8    
#> 
#> time zone: Etc/GMT-1
#> tzcode source: internal
#> 
#> attached base packages:
#> [1] stats4    stats     graphics  grDevices utils     datasets  methods  
#> [8] base     
#> 
#> other attached packages:
#>  [1] DESeq2_1.40.1               SummarizedExperiment_1.30.1
#>  [3] Biobase_2.60.0              MatrixGenerics_1.12.0      
#>  [5] matrixStats_0.63.0          GenomicRanges_1.52.0       
#>  [7] GenomeInfoDb_1.36.0         IRanges_2.34.0             
#>  [9] S4Vectors_0.38.1            BiocGenerics_0.46.0        
#> 
#> loaded via a namespace (and not attached):
#>  [1] utf8_1.2.3              bitops_1.0-7            lattice_0.21-8         
#>  [4] digest_0.6.31           magrittr_2.0.3          evaluate_0.21          
#>  [7] grid_4.3.0              mvtnorm_1.1-3           fastmap_1.1.1          
#> [10] plyr_1.8.8              Matrix_1.5-4            fansi_1.0.4            
#> [13] scales_1.2.1            numDeriv_2016.8-1.1     codetools_0.2-19       
#> [16] emdbook_1.3.12          cli_3.6.1               bbmle_1.0.25           
#> [19] rlang_1.1.1             crayon_1.5.2            XVector_0.40.0         
#> [22] munsell_0.5.0           reprex_2.0.2            withr_2.5.0            
#> [25] DelayedArray_0.26.2     yaml_2.3.7              S4Arrays_1.0.1         
#> [28] tools_4.3.0             parallel_4.3.0          BiocParallel_1.34.1    
#> [31] coda_0.19-4             bdsmatrix_1.3-6         colorspace_2.1-0       
#> [34] ggplot2_3.4.2           locfit_1.5-9.7          GenomeInfoDbData_1.2.10
#> [37] vctrs_0.6.2             R6_2.5.1                lifecycle_1.0.3        
#> [40] zlibbioc_1.46.0         fs_1.6.2                MASS_7.3-58.4          
#> [43] pkgconfig_2.0.3         pillar_1.9.0            gtable_0.3.3           
#> [46] glue_1.6.2              Rcpp_1.0.10             xfun_0.39              
#> [49] tibble_3.2.1            rstudioapi_0.14         knitr_1.42             
#> [52] htmltools_0.5.5         apeglm_1.22.0           rmarkdown_2.21         
#> [55] compiler_4.3.0          RCurl_1.98-1.12
Created on 2023-05-10 with reprex v2.0.2
DESeq2 apeglm • 2.4k views
ADD COMMENT
1
Entering edit mode
@mikelove
Last seen 2 days ago
United States

Thanks ATpoint !

Yeah weird, i still have the previous release and it's fine, I needed to update anyway, so will check this once I have it set up.

set.seed(1)
dds <- makeExampleDESeqDataSet(n=50000, m=20)
dds <- DESeq(dds, quiet=TRUE)
res <- lfcShrink(dds, coef=2, type="apeglm", quiet=TRUE)
packageVersion("DESeq2")
[1] 1.38.3
ADD COMMENT
0
Entering edit mode

I'm having trouble getting 1.40 binary on my Mac and it doesn't appear in 1.39.8. Will let you know once I've got a setup that I can debug.

ADD REPLY
0
Entering edit mode

More to the general audience than to you, but I really like Docker to get different environments going reproducibly, especially for projects that run over years when the used Bioc version is long deprecated. You could conveniently use docker run -d -p 8787:8787 -e PASSWORD=bioc bioconductor/bioconductor_docker:RELEASE_3_17 to spin up the current (or any) Bioconductor image and then log into the RStudio-Server session via localhost/8787. Nice thing is that the environment is constant (Ubuntu) and any possible external dependency for any Bioc package is already installed.

ADD REPLY
0
Entering edit mode

Ok I can now recreate this with 1.40.1 on a Mac. Will work on it.

ADD REPLY
1
Entering edit mode

I believe this should be resolved here:

https://github.com/azhu513/apeglm/issues/5

If someone can check with devtools::install_github("azhu513/apeglm") I can port the fix to release (apeglm).

Note that I was able to recreate the bug in release with fewer genes:

set.seed(1)
dds <- makeExampleDESeqDataSet(n=10000, m=20)
ADD REPLY
0
Entering edit mode

Yes, works now, thank you!

ADD REPLY
0
Entering edit mode

Ok pushed just now to 1.22.1 (apeglm) should be online by Sunday.

Thanks for your reprex.

ADD REPLY

Login before adding your answer.

Traffic: 621 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