Can this be used on scRNA-seq dataset?

[Jiayi-Zheng]

Hello, I have some scRNA-seq dataset and was trying to adopt variancePartition, with considering the expression data of each single cell as a sample, with two other categories of labels to determine variance partition. However, I've been encountering various of errors when running variancePartition on the normalized count matrix of scRNA-seq, from first telling me singular fit of model to now Error in stats::as.formula(formula) : object 'form' not found Calls: fitExtractVarPartModel ... fitExtractVarPartModel -> .fitExtractVarPartModel -> <Anonymous> I'm not sure if this model can be fitted for scRNA-seq dataset, or that some other processing is needed, like shortlisting genes, or pseudocell construction?

Thank you!

[GabrielHoffman]

Can you show sessionInfo()

[Jiayi-Zheng]

@GabrielHoffman Thank you so much for your reply! Attached is it:

> sessionInfo()
R version 4.1.0 (2021-05-18)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: CentOS Linux 8

Matrix products: default
BLAS/LAPACK: /usr/lib64/libopenblasp-r0.3.12.so

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

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

other attached packages:
 [1] Matrix_1.5-3             variancePartition_1.22.0 Biobase_2.52.0
 [4] BiocGenerics_0.38.0      scales_1.2.1             BiocParallel_1.26.2
 [7] limma_3.48.3             ggplot2_3.4.3            MuDataSeurat_0.0.0.9000
[10] Seurat_4.9.9.9060        SeuratObject_4.9.9.9091  sp_2.0-0

loaded via a namespace (and not attached):
  [1] minqa_1.2.5            Rtsne_0.16             colorspace_2.1-0
  [4] deldir_1.0-9           colorRamps_2.3.1       ellipsis_0.3.2
  [7] ggridges_0.5.4         RcppHNSW_0.5.0         spatstat.data_3.0-1
 [10] leiden_0.4.3           listenv_0.9.0          ggrepel_0.9.3
 [13] bit64_4.0.5            RSpectra_0.16-1        fansi_1.0.4
 [16] codetools_0.2-18       splines_4.1.0          doParallel_1.0.16
 [19] polyclip_1.10-4        spam_2.9-1             jsonlite_1.8.7
 [22] nloptr_2.0.3           pbkrtest_0.5.2         broom_1.0.3
 [25] ica_1.0-3              cluster_2.1.2          png_0.1-8
 [28] uwot_0.1.16            shiny_1.7.5            sctransform_0.4.0
 [31] spatstat.sparse_3.0-2  compiler_4.1.0         httr_1.4.7
 [34] backports_1.2.1        fastmap_1.1.1          lazyeval_0.2.2
 [37] cli_3.6.1              later_1.3.1            prettyunits_1.1.1
 [40] htmltools_0.5.6        tools_4.1.0            igraph_1.3.5
 [43] dotCall64_1.0-2        gtable_0.3.4           glue_1.6.2
 [46] RANN_2.6.1             reshape2_1.4.4         dplyr_1.1.3
 [49] Rcpp_1.0.11            scattermore_1.2        vctrs_0.6.3
 [52] spatstat.explore_3.2-3 nlme_3.1-152           progressr_0.14.0
 [55] iterators_1.0.13       lmtest_0.9-40          spatstat.random_3.1-6
 [58] stringr_1.5.0          globals_0.16.2         lme4_1.1-31
 [61] mime_0.12              miniUI_0.1.1.1         lifecycle_1.0.3
 [64] irlba_2.3.5.1          gtools_3.9.4           goftest_1.2-3
 [67] future_1.33.0          MASS_7.3-54            zoo_1.8-12
 [70] hms_1.1.2              promises_1.2.1         spatstat.utils_3.0-3
 [73] RColorBrewer_1.1-3     reticulate_1.32.0      pbapply_1.7-2
 [76] gridExtra_2.3          stringi_1.7.12         fastDummies_1.7.3
 [79] foreach_1.5.1          caTools_1.18.2         boot_1.3-28
 [82] bitops_1.0-7           rlang_1.1.1            pkgconfig_2.0.3
 [85] matrixStats_1.0.0      lattice_0.20-44        ROCR_1.0-11
 [88] purrr_1.0.2            tensor_1.5             patchwork_1.1.3
 [91] htmlwidgets_1.6.2      cowplot_1.1.1          bit_4.0.4
 [94] tidyselect_1.2.0       parallelly_1.36.0      RcppAnnoy_0.0.21
 [97] plyr_1.8.8             magrittr_2.0.3         R6_2.5.1
[100] gplots_3.1.3           generics_0.1.3         pillar_1.9.0
[103] withr_2.5.0            fitdistrplus_1.1-11    survival_3.2-12
[106] abind_1.4-5            tibble_3.2.1           future.apply_1.11.0
[109] crayon_1.5.2           hdf5r_1.3.8            KernSmooth_2.23-20
[112] utf8_1.2.3             spatstat.geom_3.2-5    plotly_4.10.2
[115] progress_1.2.2         grid_4.1.0             data.table_1.14.8
[118] digest_0.6.33          xtable_1.8-4           tidyr_1.3.0
[121] httpuv_1.6.11          munsell_0.5.0          viridisLite_0.4.2

The first error message I encountered while trying variancepartition on scRNA-seq dataset (log normalized counts, 16k protein coding genes, categorized by cell type and protocol) was

boundary (singular) fit: see help('isSingular')
Calls: fitExtractVarPartModel ... .fitExtractVarPartModel -> checkModelStatus -> checkModelStatus

I then did some correlation test between expression of genes, and removed the genes with colinearity (threshold by pearson correlation of normalized expression > 0.8 or <-0.8), and then I got

Error in stats::as.formula(formula) : object 'form' not found Calls: fitExtractVarPartModel ... fitExtractVarPartModel -> .fitExtractVarPartModel -> <Anonymous>

I'm not sure if I should continue to fix things based on this error message or that the single cell format is not so compatible with the package... Meanwhile, I was trying to do pseudocell (averaging normalized expression of ~30-40 single cells and aggregate to one cell, as a made-up bulk seq expression, I trialed it with the 16k genes and again got boundary (singular) fit: see help('isSingular')

[GabrielHoffman]

That version is pretty old. This is resolved in recent versions. Can you install R 4.3 and then install variancePartition

[Jiayi-Zheng]

Thank you so much! Now it's working!