CITE-seq setup for scviR, 2025
May 12, 2025
Contents
1 Overview
This vignette concerns a 2025 update of scviR.
The objective is to allow exploration of CITE-seq data with TOTALVI as illustrated in notebooks provided in the scvi-tools project. Ultimately it would be desirable to compare the analyses of the OSCA book to those produced with TOTALVI, but at this time we are focused on tool interoperability.
As of May 2025, use BiocManager to install scviR in R 4.5 or above:
BiocManager::install("scviR")Note that this package uses basilisk primarily to pin versions of associated software. We expose python objects in the global environment. When the required API comes into focus, more isolation of python operations and objects will be established.
2 Acquire CITE-seq data from 10x human discovery platform (15k cells, 349 proteins)
library(scviR)
HDP.h5 = cacheCiteseqHDPdata()
mdata1 = muonR()$read_10x_h5(HDP.h5)
mdata1$mod["rna"]$var_names_make_unique()
reticulate::py_run_string('r.mdata1.mod["rna"].layers["counts"] = r.mdata1.mod["rna"].X.copy()')
mdata1## MuData object with n_obs × n_vars = 15530 × 18431
## var: 'gene_ids', 'feature_types', 'ENS.GENE', 'IsoCtrl', 'genome', 'map_rna', 'pattern', 'read', 'sequence', 'uniprot'
## 2 modalities
## rna: 15530 x 18082
## var: 'gene_ids', 'feature_types', 'ENS.GENE', 'IsoCtrl', 'genome', 'map_rna', 'pattern', 'read', 'sequence', 'uniprot'
## layers: 'counts'
## prot: 15530 x 349
## var: 'gene_ids', 'feature_types', 'ENS.GENE', 'IsoCtrl', 'genome', 'map_rna', 'pattern', 'read', 'sequence', 'uniprot'3 Use the scvi-tools approach to preprocessing
Filter genes using scanpy.
scr = scanpyR()
scr$pp$normalize_total(mdata1$mod["rna"])
scr$pp$log1p(mdata1$mod["rna"])
#
scr$pp$highly_variable_genes(
mdata1$mod["rna"],
n_top_genes=4000L,
flavor="seurat_v3",
layer="counts",
)Add the filtered data to the MuData instance.
py_run_string('r.mdata1.mod["rna_subset"] = r.mdata1.mod["rna"][:, r.mdata1.mod["rna"].var["highly_variable"]].copy()')
mdata1 = MuDataR()$MuData(mdata1$mod)Produce dense versions of quantification matrices, and “update”.
4 Prepare for TOTALVI
Text of notebook:
Now we run `setup_mudata`, which is the MuData analog to `setup_anndata`.
The caveat of this workflow is that we need to provide this function which
modality of the `mdata` object contains each piece of data. So for example,
the batch information is in `mdata.mod["rna"].obs["batch"]`. Therefore, in the `modalities`
argument below we specify that the `batch_key` can be
found in the `"rna_subset"` modality of the MuData object.
Notably, we provide `protein_layer=None`. This means scvi-tools will pull
information from `.X` from the modality specified in `modalities` (`"protein"`
in this case). In the case of RNA, we want to use the counts,
which we stored in `mdata.mod["rna"].layers["counts"]`.5 Train TOTALVI
Here’s the model:
Use model$module to see the complete architecture.
Perform truncated training:
n_epochs = 50L
n_epochs.cpu = 5L
acc = "cpu"
tchk = try(reticulate::import("torch"))
if (!inherits(tchk, "try-error") && tchk$backends$mps$is_available()) acc = "mps"
if (!inherits(tchk, "try-error") && tchk$cuda$is_available()) acc = "gpu"
runtim.cpu = system.time(model$train(max_epochs=n_epochs.cpu, accelerator = "cpu"))
runtim = system.time(model$train(max_epochs=n_epochs, accelerator = acc))The gpu processor was used leading to an average clock time per epoch of 0.99 seconds.
By contrast, the average epoch runtime for a very short run with CPU only is 10.84.
Extract the ELBO criteria. The plot below includes both the CPU and GPU (if available) results.
total_ep = n_epochs + n_epochs.cpu
val_elbo = unlist(model$history$elbo_validation)
tr_elbo = model$history$elbo_train$elbo_trai
plot(1:total_ep, tr_elbo, type="l", xlab="epoch", ylab="ELBO", ylim=c(0,7000))
lines(1:total_ep, val_elbo, col="blue")
legend(3, 6000, col=c("black", "blue"), lty=1, legend=c("train", "validate"))
6 Session info
## R version 4.5.0 Patched (2025-05-01 r88189)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.2 LTS
##
## Matrix products: default
## BLAS: /media/volume/biocgpu/biocbuild/bbs-3.22-bioc-gpu/R/lib/libRblas.so
## LAPACK: /media/volume/biocgpu/biocbuild/bbs-3.22-bioc-gpu/R/lib/libRlapack.so; LAPACK version 3.12.1
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_GB LC_COLLATE=C
## [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
##
## time zone: America/New_York
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] scviR_1.9.11 shiny_1.10.0
## [3] basilisk_1.21.0 reticulate_1.42.0
## [5] scater_1.37.0 ggplot2_3.5.2
## [7] scuttle_1.19.0 SingleCellExperiment_1.31.0
## [9] SummarizedExperiment_1.39.0 Biobase_2.69.0
## [11] GenomicRanges_1.61.0 GenomeInfoDb_1.45.3
## [13] IRanges_2.43.0 S4Vectors_0.47.0
## [15] BiocGenerics_0.55.0 generics_0.1.3
## [17] MatrixGenerics_1.21.0 matrixStats_1.5.0
## [19] BiocStyle_2.37.0
##
## loaded via a namespace (and not attached):
## [1] DBI_1.2.3 gridExtra_2.3 httr2_1.1.2
## [4] rlang_1.1.6 magrittr_2.0.3 compiler_4.5.0
## [7] RSQLite_2.3.11 mgcv_1.9-3 dir.expiry_1.17.0
## [10] png_0.1-8 vctrs_0.6.5 pkgconfig_2.0.3
## [13] crayon_1.5.3 fastmap_1.2.0 dbplyr_2.5.0
## [16] XVector_0.49.0 labeling_0.4.3 promises_1.3.2
## [19] rmarkdown_2.29 UCSC.utils_1.5.0 ggbeeswarm_0.7.2
## [22] tinytex_0.57 purrr_1.0.4 bit_4.6.0
## [25] xfun_0.52 cachem_1.1.0 beachmat_2.25.0
## [28] jsonlite_2.0.0 blob_1.2.4 later_1.4.2
## [31] DelayedArray_0.35.1 BiocParallel_1.43.0 irlba_2.3.5.1
## [34] parallel_4.5.0 R6_2.6.1 bslib_0.9.0
## [37] RColorBrewer_1.1-3 limma_3.65.0 jquerylib_0.1.4
## [40] Rcpp_1.0.14 bookdown_0.43 knitr_1.50
## [43] splines_4.5.0 httpuv_1.6.16 Matrix_1.7-3
## [46] tidyselect_1.2.1 abind_1.4-8 yaml_2.3.10
## [49] viridis_0.6.5 codetools_0.2-20 curl_6.2.2
## [52] lattice_0.22-7 tibble_3.2.1 basilisk.utils_1.21.0
## [55] withr_3.0.2 evaluate_1.0.3 BiocFileCache_2.99.0
## [58] pillar_1.10.2 BiocManager_1.30.25 filelock_1.0.3
## [61] scales_1.4.0 xtable_1.8-4 glue_1.8.0
## [64] pheatmap_1.0.12 tools_4.5.0 BiocNeighbors_2.3.0
## [67] ScaledMatrix_1.17.0 grid_4.5.0 nlme_3.1-168
## [70] beeswarm_0.4.0 BiocSingular_1.25.0 vipor_0.4.7
## [73] cli_3.6.5 rsvd_1.0.5 rappdirs_0.3.3
## [76] S4Arrays_1.9.0 viridisLite_0.4.2 dplyr_1.1.4
## [79] gtable_0.3.6 sass_0.4.10 digest_0.6.37
## [82] SparseArray_1.9.0 ggrepel_0.9.6 farver_2.1.2
## [85] memoise_2.0.1 htmltools_0.5.8.1 lifecycle_1.0.4
## [88] httr_1.4.7 statmod_1.5.0 mime_0.13
## [91] bit64_4.6.0-1