"missing value where TRUE/FALSE needed" error occurring with certain seeds (for nb_fast and poisson both with theta.ml)

[diegoalexespi]

Hi. I am receiving the following error on my counts data when attempting to use the vst function.

> z <- sctransform::vst(my_counts, verbosity = 2, method = "poisson", theta_estimation_fun = "theta.ml")
Calculating cell attributes from input UMI matrix: log_umi
Variance stabilizing transformation of count matrix of size 15887 by 2691
Model formula is y ~ log_umi
Get Negative Binomial regression parameters per gene
Using 2000 genes, 2691 cells
  |==============================================================                                                             |  50%Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 't': missing value where TRUE/FALSE needed
In addition: There were 50 or more warnings (use warnings() to see the first 50)

The error only appears when using either method = "nb_fast" or method = "poisson" with theta_estimation_fun = "theta.ml". The error disappears when using method = "glmGamPoi" or method = "qpoisson". The error also disappears for method = "nb_fast" and method = "poisson" if I use theta_estimation_fun = "theta.mm". Moreover, the error appearance for "nb_fast" and "poisson" depends on the seed:

> set.seed(90835)
> z <- sctransform::vst(my_counts, verbosity = 2, method = "poisson", theta_estimation_fun = "theta.ml")
Calculating cell attributes from input UMI matrix: log_umi
Variance stabilizing transformation of count matrix of size 15887 by 2691
Model formula is y ~ log_umi
Get Negative Binomial regression parameters per gene
Using 2000 genes, 2691 cells
  |                                                                                                                           |   0%Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 't': missing value where TRUE/FALSE needed
In addition: There were 50 or more warnings (use warnings() to see the first 50)
> set.seed(1234)
> z <- sctransform::vst(my_counts, verbosity = 2, method = "poisson", theta_estimation_fun = "theta.ml")
Calculating cell attributes from input UMI matrix: log_umi
Variance stabilizing transformation of count matrix of size 15887 by 2691
Model formula is y ~ log_umi
Get Negative Binomial regression parameters per gene
Using 2000 genes, 2691 cells
  |==============================================================                                                             |  50%Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 't': missing value where TRUE/FALSE needed
In addition: There were 50 or more warnings (use warnings() to see the first 50)
> set.seed(123)
> z <- sctransform::vst(my_counts, verbosity = 2, method = "poisson", theta_estimation_fun = "theta.ml")
Calculating cell attributes from input UMI matrix: log_umi
Variance stabilizing transformation of count matrix of size 15887 by 2691
Model formula is y ~ log_umi
Get Negative Binomial regression parameters per gene
Using 2000 genes, 2691 cells
  |===========================================================================================================================| 100%
Found 83 outliers - those will be ignored in fitting/regularization step

Second step: Get residuals using fitted parameters for 15887 genes
  |===========================================================================================================================| 100%
Calculating gene attributes
Wall clock passed: Time difference of 21.56692 secs
There were 50 or more warnings (use warnings() to see the first 50)
>

It seems that the error also ends at similar places for nb_fast and poisson when the seed is the same for both:

> set.seed(90835)
> z <- sctransform::vst(my_counts, verbosity = 2, method = "nb_fast", theta_estimation_fun = "theta.ml")
Calculating cell attributes from input UMI matrix: log_umi
Variance stabilizing transformation of count matrix of size 15887 by 2691
Model formula is y ~ log_umi
Get Negative Binomial regression parameters per gene
Using 2000 genes, 2691 cells
  |                                                                                                                           |   0%Error in while ((it <- it + 1) < limit && abs(del) > eps) { : 
  missing value where TRUE/FALSE needed
In addition: There were 50 or more warnings (use warnings() to see the first 50)
> set.seed(90835)
> z <- sctransform::vst(my_counts, verbosity = 2, method = "poisson", theta_estimation_fun = "theta.ml")
Calculating cell attributes from input UMI matrix: log_umi
Variance stabilizing transformation of count matrix of size 15887 by 2691
Model formula is y ~ log_umi
Get Negative Binomial regression parameters per gene
Using 2000 genes, 2691 cells
  |                                                                                                                           |   0%Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 't': missing value where TRUE/FALSE needed
In addition: There were 50 or more warnings (use warnings() to see the first 50)
>

Any help is appreciated! Happy to provide more information to help as best as I can. I tried to reproduce the error with pbmc_small and pbmc3k, but I have not been able to. Here is my sessionInfo():

R version 4.0.3 (2020-10-10)
Platform: x86_64-apple-darwin17.0 (64-bit)
Running under: macOS Catalina 10.15.7

Matrix products: default
BLAS:   /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

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

other attached packages:
 [1] sctransform_0.3.1       pbmc3k.SeuratData_3.1.4 Nebulosa_0.99.94        patchwork_1.0.1         readxl_1.3.1           
 [6] readr_1.4.0             kableExtra_1.2.1        cowplot_1.1.0           Seurat_3.9.9.9002       ggplot2_3.3.2          
[11] tibble_3.0.4            tidyr_1.1.2             dplyr_1.0.2             knitr_1.30              magrittr_1.5           

loaded via a namespace (and not attached):
  [1] backports_1.1.10            plyr_1.8.6                  igraph_1.2.6                lazyeval_0.2.2             
  [5] splines_4.0.3               listenv_0.8.0               usethis_1.6.3               GenomeInfoDb_1.25.11       
  [9] digest_0.6.25               htmltools_0.5.0             fansi_0.4.1                 memoise_1.1.0              
 [13] tensor_1.5                  cluster_2.1.0               ks_1.11.7                   ROCR_1.0-11                
 [17] remotes_2.2.0               globals_0.13.1              matrixStats_0.57.0          prettyunits_1.1.1          
 [21] colorspace_1.4-1            rappdirs_0.3.1              rvest_0.3.6                 ggrepel_0.8.2              
 [25] xfun_0.18                   callr_3.5.1                 crayon_1.3.4                RCurl_1.98-1.2             
 [29] jsonlite_1.7.1              spatstat_1.64-1             spatstat.data_1.4-3         survival_3.2-7             
 [33] zoo_1.8-8                   glue_1.4.2                  polyclip_1.10-0             gtable_0.3.0               
 [37] zlibbioc_1.35.0             XVector_0.29.3              webshot_0.5.2               leiden_0.3.3               
 [41] DelayedArray_0.15.16        pkgbuild_1.1.0              future.apply_1.6.0          SingleCellExperiment_1.11.8
 [45] BiocGenerics_0.35.4         abind_1.4-5                 scales_1.1.1                pheatmap_1.0.12            
 [49] mvtnorm_1.1-1               miniUI_0.1.1.1              Rcpp_1.0.5                  viridisLite_0.3.0          
 [53] xtable_1.8-4                reticulate_1.16             bit_4.0.4                   rsvd_1.0.3                 
 [57] mclust_5.4.6                stats4_4.0.3                htmlwidgets_1.5.2           httr_1.4.2                 
 [61] RColorBrewer_1.1-2          ellipsis_0.3.1              ica_1.0-2                   pkgconfig_2.0.3            
 [65] uwot_0.1.8.9001             deldir_0.1-29               tidyselect_1.1.0            rlang_0.4.8                
 [69] reshape2_1.4.4              later_1.1.0.1               munsell_0.5.0               cellranger_1.1.0           
 [73] tools_4.0.3                 cli_2.1.0                   generics_0.0.2              devtools_2.3.2             
 [77] ggridges_0.5.2              evaluate_0.14               stringr_1.4.0               fastmap_1.0.1              
 [81] yaml_2.2.1                  goftest_1.2-2               fs_1.5.0                    processx_3.4.4             
 [85] bit64_4.0.5                 fitdistrplus_1.1-1          SeuratData_0.2.1            purrr_0.3.4                
 [89] RANN_2.6.1                  pbapply_1.4-3               future_1.19.1               nlme_3.1-149               
 [93] mime_0.9                    xml2_1.3.2                  hdf5r_1.3.3                 compiler_4.0.3             
 [97] rstudioapi_0.11             curl_4.3                    plotly_4.9.2.1              png_0.1-7                  
[101] testthat_2.3.2              spatstat.utils_1.17-0       glmGamPoi_1.1.13            stringi_1.5.3              
[105] ps_1.4.0                    desc_1.2.0                  lattice_0.20-41             Matrix_1.2-18              
[109] SeuratDisk_0.0.0.9013       vctrs_0.3.4                 pillar_1.4.6                lifecycle_0.2.0            
[113] BiocManager_1.30.10         lmtest_0.9-38               RcppAnnoy_0.0.16            data.table_1.13.0          
[117] bitops_1.0-6                irlba_2.3.3                 httpuv_1.5.4                GenomicRanges_1.41.6       
[121] R6_2.4.1                    promises_1.1.1              KernSmooth_2.23-17          gridExtra_2.3              
[125] IRanges_2.23.10             sessioninfo_1.1.1           codetools_0.2-16            pkgload_1.1.0              
[129] MASS_7.3-53                 assertthat_0.2.1            SummarizedExperiment_1.19.9 rprojroot_1.3-2            
[133] withr_2.3.0                 S4Vectors_0.27.14           GenomeInfoDbData_1.2.4      mgcv_1.8-33                
[137] parallel_4.0.3              hms_0.5.3                   grid_4.0.3                  rpart_4.1-15               
[141] DelayedMatrixStats_1.11.1   rmarkdown_2.4               MatrixGenerics_1.1.6        Rtsne_0.15                 
[145] Biobase_2.49.1              shiny_1.5.0   

[ChristophH]

Hi,
Thank you for this detailed error report. Before going into the problem, could you try the current version of sctransform in the develop branch?

remotes::install_github("ChristophH/sctransform@develop")

[diegoalexespi]

Installing from the current version in the develop branch appears to have worked, thank you so much! There are no more errors of this kind apparent. This seems related to issue #65 in case that helps others.

[fuxins]

Hi,
I'm getting the same error.
I’ve installed the current version of sctransform from the develop branch, but it seems not work for me
Here's my commands and error messages

>sob<-CreateSeuratObject(
   countMetrix,
   project = projectid,
   assay = "RNA",
   names.field = 1,
   names.delim = "_",
   meta.data = NULL,
   min.cells = 5,
   min.features = 500
 )
Warning: Feature names cannot have underscores ('_'), replacing with dashes ('-')
>sob_sct <- SCTransform(sob,verbose=TRUE)
Calculating cell attributes from input UMI matrix: log_umi
Variance stabilizing transformation of count matrix of size 17032 by 428
Model formula is y ~ log_umi
Get Negative Binomial regression parameters per gene
Using 2000 genes, 428 cells
  |                                                                                                                                          |   0%Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 't': missing value where TRUE/FALSE needed
In addition: Warning messages:
1: In trigamma(th + y) : NaNs produced
2: In trigamma(th) : NaNs produced

countMetrix_partial.txt

Any advice would be appreciated!

[ChristophH]

The maximum number of features (i.e. genes) detected in the matrix you provided is 198. When I create the Seurat object with min.features = 500 it removes all cells from the downstream analysis and fails. When I set the limit to 5 all cells pass and sctransform works.
The full matrix from your example looks bigger, but I don't know exactly what is causing the error. If you have updated sctransform and Seurat, make sure to restart your R session.

[fuxins]

Thanks a lot! I restart R session then it works.