ashr icon indicating copy to clipboard operation
ashr copied to clipboard

Strange error ashr

Open william-denault opened this issue 10 months ago • 3 comments

Hello,

I was using ashr on simingly non problematic data set (attached).

I got the following error

df <- as.data.frame(read_csv("strange_ash_error.csv")) 

ash(df$x, df$sd)
 

[strange_ash_error.csv](https://github.com/user-attachments/files/18967160/strange_ash_error.csv)

Error in if (!all(nonzero_cols)) { : 
  missing value where TRUE/FALSE needed
In addition: Warning message:
In log(1 - exp(lpb - lpa)) : NaNs produced

I am a bit surprised that this cause problem given that hist( df$x/df$sd, nclass=100) seems not to pathologic

here the session Info


sessionInfo()
R version 4.3.3 (2024-02-29 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 11 x64 (build 26120)

Matrix products: default


locale:
[1] LC_COLLATE=English_United States.utf8 
[2] LC_CTYPE=English_United States.utf8   
[3] LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C                          
[5] LC_TIME=English_United States.utf8    

time zone: Etc/GMT-1
tzcode source: internal

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

other attached packages:
 [1] readr_2.1.5                             
 [2] fsusieR_0.2.94                          
 [3] testthat_3.2.3                          
 [4] TxDb.Hsapiens.UCSC.hg38.knownGene_3.18.0
 [5] GenomicFeatures_1.54.4                  
 [6] Gviz_1.46.1                             
 [7] GenomicRanges_1.54.1                    
 [8] GenomeInfoDb_1.38.8                     
 [9] org.Hs.eg.db_3.18.0                     
[10] AnnotationDbi_1.64.1                    
[11] IRanges_2.36.0                          
[12] S4Vectors_0.40.2                        
[13] Biobase_2.62.0                          
[14] AnnotationHub_3.10.1                    
[15] BiocFileCache_2.10.2                    
[16] dbplyr_2.5.0                            
[17] BiocGenerics_0.48.1                     

loaded via a namespace (and not attached):
  [1] later_1.4.1                   BiocIO_1.12.0                
  [3] bitops_1.0-9                  filelock_1.0.3               
  [5] tibble_3.2.1                  XML_3.99-0.18                
  [7] rpart_4.1.23                  lifecycle_1.0.4              
  [9] mixsqp_0.3-54                 rprojroot_2.0.4              
 [11] vroom_1.6.5                   lattice_0.22-5               
 [13] ensembldb_2.26.1              MASS_7.3-60.0.1              
 [15] backports_1.5.0               magrittr_2.0.3               
 [17] Hmisc_5.2-2                   rmarkdown_2.29               
 [19] yaml_2.3.10                   remotes_2.5.0                
 [21] httpuv_1.6.15                 sessioninfo_1.2.3            
 [23] pkgbuild_1.4.6                cowplot_1.1.3                
 [25] DBI_1.2.3                     RColorBrewer_1.1-3           
 [27] abind_1.4-8                   pkgload_1.4.0                
 [29] zlibbioc_1.48.2               purrr_1.0.2                  
 [31] AnnotationFilter_1.26.0       biovizBase_1.50.0            
 [33] RCurl_1.98-1.16               nnet_7.3-20                  
 [35] VariantAnnotation_1.48.1      rappdirs_0.3.3               
 [37] GenomeInfoDbData_1.2.11       irlba_2.3.5.1                
 [39] wavethresh_4.7.3              codetools_0.2-19             
 [41] DelayedArray_0.28.0           xml2_1.3.6                   
 [43] tidyselect_1.2.1              matrixStats_1.5.0            
 [45] base64enc_0.1-3               GenomicAlignments_1.38.2     
 [47] ellipsis_0.3.2                Formula_1.2-5                
 [49] tools_4.3.3                   progress_1.2.3               
 [51] Rcpp_1.0.14                   glue_1.8.0                   
 [53] gridExtra_2.3                 SparseArray_1.2.4            
 [55] xfun_0.50                     MatrixGenerics_1.14.0        
 [57] usethis_3.1.0                 dplyr_1.1.4                  
 [59] withr_3.0.2                   BiocManager_1.30.25          
 [61] fastmap_1.2.0                 latticeExtra_0.6-30          
 [63] caTools_1.18.3                digest_0.6.37                
 [65] truncnorm_1.0-9               R6_2.6.1                     
 [67] mime_0.12                     colorspace_2.1-1             
 [69] jpeg_0.1-10                   dichromat_2.0-0.1            
 [71] biomaRt_2.58.2                RSQLite_2.3.9                
 [73] generics_0.1.3                data.table_1.16.4            
 [75] rtracklayer_1.62.0            prettyunits_1.2.0            
 [77] httr_1.4.7                    htmlwidgets_1.6.4            
 [79] S4Arrays_1.2.1                pkgconfig_2.0.3              
 [81] gtable_0.3.6                  blob_1.2.4                   
 [83] XVector_0.42.0                brio_1.1.5                   
 [85] htmltools_0.5.8.1             profvis_0.4.0                
 [87] smashr_1.2-9                  ProtGenerics_1.34.0          
 [89] scales_1.3.0                  png_0.1-8                    
 [91] ashr_2.2-63                   knitr_1.49                   
 [93] rstudioapi_0.17.1             tzdb_0.4.0                   
 [95] rjson_0.2.23                  checkmate_2.3.2              
 [97] curl_6.2.0                    cachem_1.1.0                 
 [99] stringr_1.5.1                 BiocVersion_3.18.1           
[101] parallel_4.3.3                miniUI_0.1.1.1               
[103] foreign_0.8-86                RcppZiggurat_0.1.6           
[105] restfulr_0.0.15               desc_1.4.3                   
[107] pillar_1.10.1                 vctrs_0.6.5                  
[109] urlchecker_1.0.1              promises_1.3.2               
[111] xtable_1.8-4                  cluster_2.1.6                
[113] htmlTable_2.4.3               evaluate_1.0.3               
[115] invgamma_1.1                  cli_3.6.3                    
[117] compiler_4.3.3                Rsamtools_2.18.0             
[119] rlang_1.1.5                   crayon_1.5.3                 
[121] SQUAREM_2021.1                LaplacesDemon_16.1.6         
[123] interp_1.1-6                  fs_1.6.5                     
[125] stringi_1.8.4                 deldir_2.0-4                 
[127] BiocParallel_1.36.0           munsell_0.5.1                
[129] Biostrings_2.70.3             lazyeval_0.2.2               
[131] devtools_2.4.5                Matrix_1.6-5                 
[133] BSgenome_1.70.2               hms_1.1.3                    
[135] bit64_4.6.0-1                 ggplot2_3.5.1                
[137] KEGGREST_1.42.0               shiny_1.10.0                 
[139] SummarizedExperiment_1.32.0   interactiveDisplayBase_1.40.0
[141] Rfast_2.1.4                   memoise_2.0.1                
[143] RcppParallel_5.1.10           bit_4.5.0.1            

william-denault avatar Feb 25 '25 15:02 william-denault

@william-denault I'm not quite sure, but I'll note that

fit <- ash(df$x,df$sd,mixcompdist = "normal")

worked for me.

pcarbo avatar Feb 25 '25 15:02 pcarbo

@william-denault @pcarbo There is a floating point issue when evaluating log_comp_dens_conv.unimix on this data set:

Browse[1]> mean(abs(lpa - lpb) < .Machine$double.eps)
[1] 0.5059023

As a result, the log likelihood matrix is filled with NaN and Inf.

aksarkar avatar Feb 25 '25 23:02 aksarkar

@william-denault The data are weird: you have many data points that are very close to zero and also have very small s.e.'s. There could be something wrong with these data.

Image

pcarbo avatar Feb 26 '25 03:02 pcarbo