Comparative transcriptomics aims to identify common and differing patterns (organisms, experimental conditions, tissues, and so on) between a set of conditions by using transcriptomics data The methodology to tackle the problem will vary greatly depending on which patterns we aim to compare, the amount and biological nature of them and which omics resolution. singIST is a comparative transcriptomics method that aims to compare a disease model against a human disease of interest. In drug discovery and pharmaceutical early development choosing a disease model becomes a fundamental endavour. For this reason, we present singIST as a method that compares the single-cell transcriptomics profile of a disease model, either In Vivo or In Vitro models, against a human disease that such model should mimick.
Other Bioconductor packages for comparative transcriptomics are available. One of such is HybridExpress which aims to compare the bulk transcriptomics (RNA-seq) profile for hybrids relative to their progenitor species. CoSIA is a Bioconductor package that enables organism-to-organism comparative bulk transcriptomics (RNAseq) across tissues. However, none of these packages are direct comparators to **singIST* since they are not conceived for single-cell comparative transcriptomics but for bulk transcriptomics instead.
This Vignette demonstrates how to use the singIST package. singIST quantifies the direction and mangitude of single-cell transcriptomic changes between a human reference condition and a disease model at pathway, cell type, and gene resolution. Throughout this vignette we will work with the following data:
The method is defined around three inputs:
In the following subsections we explicitly construct these inputs. The rest of the vignette uses these inputs to: (1) learn a human reference model using adaptive sparse multi-block PLS-DA, and then (2) evaluate how human-like the disease model changes are, in the learned human reference space.
singIST does not operate on a pathway gene set alone. Instead, it defines a superpathway by three things:
In many practical cases, we use the same pathway gene set for every cell type. However, singIST also supports cell-type specific gene sets. Below we create two example of superpathways used later in this vignette.
First, because singIST relies on MSigDB curated pathways, we load the
MSigDB gene set collection once and reuse it as necessary. We assume the
reference to compare against is human org = "hs" and we
extract both ENTREZID and HGNC gene symbols
id = c("SYM", "EZID").
gse <- msigdb::getMsigdb(org = "hs", id = c("SYM", "EZID"),
version = msigdb::getMsigdbVersions()[1])In this first example, we define a pathway using its MSigDB metadata, choose the cell types we want to model, and then assign the same pathway gene set to each cell type gene set. This is the simplest and most common setup for a superpathway.
# Pathway definition (curated gene set identifier)
dc_pathway <- create_pathway(standard_name = "BIOCARTA_DC_PATHWAY",
dbsource = "BIOCARTA",
collection = "c2",
subcollection = "CP")
# Superpathway = pathway + modeled cell types + per cell type gene sets
dc_superpathway <- create_superpathway(pathway_info = dc_pathway,
celltypes = c("Dendritic Cells",
"Langerhan Cells",
"T-cell"),
list())
# Assign gene sets per cell type (retrieve gene set from MSigDB)
dc_superpathway <- setGeneSetsCelltype(dc_superpathway, gse = gse)In some applications, you may want different cell types to use different gene sets. singIST supports this by allowing the user to directly provide a list of gene vectors, one per cell type.
Here we demonstrate this capability by taking the curated pathway gene set and creating different subsets for each cell type. (This is purely illustrative; in real analyses these subsets would be defined biologically rather than randomly sampled.)
# Pathway definition (curated gene set identifier)
cytokine_pathway <- create_pathway(
standard_name = "KEGG_CYTOKINE_CYTOKINE_RECEPTOR_INTERACTION",
dbsource = "KEGG",
collection = "c2",
subcollection = "CP"
)
# Superpathway = pathway + modeled cell types + per cell type gene sets
cytokine_superpathway <- create_superpathway(
pathway_info = cytokine_pathway,
celltypes = c("Dendritic Cells", "Keratinocytes", "Langerhan Cells",
"Melanocytes", "T-cell"),
list()
)
# Pull the curated gene set (retrieve gene set from MSigDB)
gene_set <- pullGeneSet(cytokine_pathway, gse = gse)
# Provide a different gene set per cell type
cytokine_superpathway$gene_sets_celltype <- list(
"Dendritic Cells" = sample(gene_set, 50),
"Keratinocytes" = sample(gene_set, 100),
"Langerhan Cells" = sample(gene_set, 10),
"Melanocytes" = sample(gene_set, 60),
"T-cell" = sample(gene_set, 160)
)The human reference data are used to train an adaptive sparse
multi-block PLS on a superpathway. In this subsection we start from the
human scRNA-seq that have already been aggregated into cell-type
pseudobulk log-normalize expression. For instance, the pseudobulk can be
performed with Seurat::AggregateExpression(). Concretely
the pseudobulk matrix:
In addition to the pseudobulk matrix, singIST requires sample level labels defining a base class and a target class (e.g., HC vs AD). The target class is the human experimental group that the disease model such mimic precisely, while the base class is the human experimental group that such differentiate from the target class.
The training procedure includes a tuning step, so the human input is not only the pseudobulk matrix and sample labels - it also includes the hyperparameter search space and the cross-validation design the user wants. In practice, we tune:
number_PLS (maximum number of
latent components to test).folds_CV and
repetition_CV.We store the former settings in a hyperparameters list
and then combine all the information into the
superpathway_input list, which is the main training input
consumed by the training procedure
fitOptimal()/ multiple_fitOptimal().
# Load the pseudobulk matrix
file_dc <- system.file("extdata", "DC_counts.rda", package = "singIST")
obj_dc <- load(file_dc)
counts_dc <- get(obj_dc[1])
# Remove the cell types that we are not going to model
counts_dc <- counts_dc[!(sub("_.*$", "", rownames(counts_dc)) %in%
c("MastC-other", "Keratinocytes", "Melanocytes")), ]Next, we define the sample identifiers and class labels.
# Load the pseudobulk matrix
sample_id_dc <- c("AD1", "AD13", "AD2", "AD3","AD4", "HC1-2", "HC3", "HC4",
"HC5")
sample_class_dc <- c(rep("AD", 5), rep("HC", 4))
base_class_dc <- "HC"
target_class_dc <- "AD"Now we define the hyperparameters.
# Initialize hyperparameter object
# Choose the array of quantile combinations that will be tested in the Cross
# Validation. Ideally one should choose a wide range of quantiles to test,
# as a narrow selection might not find the globally optimal model. The more
# quantiles one chooses the more computation time it will take, for the sake
# of simplicity we choose a small amount of quantiles to test
quantile_comb_table_dc <- as.matrix(
RcppAlgos::permuteGeneral(seq(0.1, 0.95, by = 0.3),
m = length(dc_superpathway$celltypes),
TRUE))
# Store tuning + CV design
dc_hyperparameters <- create_hyperparameters(
quantile_comb_table = quantile_comb_table_dc,
outcome_type = "binary",
number_PLS = as.integer(3),
folds_CV = as.integer(1), #LOOCV
repetition_CV = as.integer(1)
)Finally, we combine all the information into a superpathway_input list.
# Load the pseudobulk matrix
file_cytokine <- system.file("extdata", "cytokine_counts.rda", package = "singIST")
obj_cytokine <- load(file_cytokine)
counts_cytokine <- get(obj_cytokine[1])
# Remove the cell types that we are not going to model
counts_cytokine <- counts_cytokine[!(sub("_.*$", "", rownames(counts_cytokine))
%in%
c("MastC-other")), ]Next, we define the sample identifiers and class labels.
# Load the pseudobulk matrix
sample_id_cytokine <- c("AD1", "AD13", "AD2", "AD3","AD4", "HC1-2", "HC3",
"HC4", "HC5")
sample_class_cytokine <- c(rep("AD", 5), rep("HC", 4))
base_class_cytokine <- "HC"
target_class_cytokine <- "AD"Now we define the hyperparameters.
# Initialize hyperparameter object
quantile_comb_table <- as.matrix(
RcppAlgos::permuteGeneral(seq(0.05, 0.95, by = 0.2),
m = length(cytokine_superpathway$celltypes),
TRUE))
# Store tuning + CV design
cytokine_hyperparameters <- create_hyperparameters(
quantile_comb_table = quantile_comb_table,
outcome_type = "binary",
number_PLS = as.integer(2),
folds_CV = as.integer(5),
repetition_CV = as.integer(5))Finally, we combine all the information into a superpathway_input list.
cytokine_superpathway_input <- create_superpathway_input(
superpathway_info = cytokine_superpathway,
hyperparameters_info = cytokine_hyperparameters,
pseudobulk_lognorm = counts_cytokine,
sample_id = sample_id_cytokine,
sample_class = sample_class_cytokine,
base_class = base_class_cytokine,
target_class = target_class_cytokine)The disease model input is a single-cell object (either a
Seurat or SingleCellExperiment) containing
scRNA-seq measurements from the model system we want to assess. Unlike
the human input (which is provided as pseudobulk matrices), the disease
model is kept at single-cell resolution and is summarized internally by
singIST when computing recapitulation scores.
To use a disease model dataset, singIST needs:
Finally, we must provide a cell type mapping from th ehuman
superpathway cell types (names used in celltypes = ...) to
the disease model clusters. A human cell type can map to one or more
model clusters. A human cell type can map to one or more model clusters;
it is also valid to provide an empty mapping (meaning that block will
not contribute for that model).
In the following example we load the oxazolone mouse dataset shipped
in the package, standarize metadata column names to those expected by
singIST, define the mapping, and construct a
mapping organism list via
create_mapping_organism().
file_oxa <- system.file("extdata", "OXA.rda", package = "singIST")
obj_oxa <- load(file_oxa)
oxazolone <- get(obj_oxa[1])Next, we ensure the object contains the metadata fields required by
singIST. The package example dataset already contains the necessary
information, but we rename the relevant columns to the standarized names
used downstream: class, celltype_cluster, and
donor.
# The example object included in singIST uses these column positions
colnames(slot(oxazolone, "meta.data"))[c(5,11,1)] <- c("class",
"celltype_cluster",
"donor")
# If the object is a SingleCellExperiment, metadata live in colData()Now, we define the human to model cell type mapping. The names of the
list must match the human cell types used in our superpathways. Values
correspond to one or more model clusters present in
celltype_cluster.
celltype_mapping <- list(
"Dendritic Cells" = c("cDC2", "cDC1", "migratory DCs"),
"Keratinocytes" = c("Keratinocytes"),
"Langerhan Cells" = c("LC"),
"Melanocytes" = c(),
"T-cell" = c("DETC", "T"))Finally, we combine the model object and mapping into a
mapping_organism list. This list specifies the model
organism we are working with, and which class labels correspond to the
base and target conditions.
oxa_org <- create_mapping_organism(organism = "Mus musculus",
target_class = "OXA", base_class = "ETOH",
celltype_mapping = list(
"Dendritic Cells" = c("cDC2", "cDC1", "migratory DCs"),
"Keratinocytes" = c("Keratinocytes"),
"Langerhan Cells" = c("LC"),
"Melanocytes" = c(),
"T-cell" = c("DETC", "T")),
counts = oxazolone)At this stage we have constructed two key training inputs:
dc_superpathway_input and
cytokine_superpathway_input. The next step is to train an
adaptive sparse multi-block PLS-DA model (asmbPLS-DA) for each
superpathway, while selecting an optimal configuration (sparsity per
cell type block and number of PLS components) under the cross-validation
design specified in the corresponding hyperparameters
object.
singIST provides two entry points to do it:
fitOptimal() fits one superpathway_input
list.multiple_fitOptimal() is a convenience wrapper that
fits multiple superpathways in one call.If the CV design is a leave-one-out CV (LOOCV), the functions provide
with a parallelization option to make the computations faster via using
parallel = TRUE and setting the desired paralellization
scheme with BiocParallel::register().
In addition to fitting the optimal model, singIST can compute model
validation metrics, including permutation based global significance test
(npermut = 100) and resampling permutation tests
(type, e.g., jackknife or subsampling). The number of
permutation for the tests can be chosen as desired.
If the human pseudobulk input contains missing values (e.g., a cell
type block absent for some samples),
fitOptimal() / multiple_fitOptimal() automatically imputes
them during model fitting. Further details of the imputation procedure
are provided in the Appendix.
# Register a parallel backend (SnowParam works on Windows/macOS/Linux);
# on Linux you may prefer MulticoreParam.
BiocParallel::register(BiocParallel::SnowParam(workers = 2,
exportglobals = FALSE, progressbar = FALSE),
default = TRUE)
# Combine all human training inputs into a list
superpathways <- list(dc_superpathway_input,
cytokine_superpathway_input)
# Fit Optimal asmbPLS-DA models + compute validation metrics
# parallel: logical vector indicating whether the CV design should be parallel
# npermut: number of permutations used for permutation based tests
# type : validation strategy (jackknife or subsampling)
# nsubsampling: number of subsampling repeats (only used when type includes
# subsampling).
superpathways_fit <- multiple_fitOptimal(
superpathways,
parallel = c(TRUE, FALSE),
npermut = c(100),
type = c("jackknife","subsampling"),
nsubsampling = c(NULL, 5)
)
# Restore serial backend (good practice)
BiocParallel::register(BiocParallel::SerialParam(), default = TRUE) After fitting the optimal human reference models
(superpathway_fit), the final step is to quantify how well
a disease model recapitulates the human target vs base signal for each
superpathway. singIST does this by applying the disease model changes to
the human base class reference and then evaluating how “target-like” the
samples become. The resulting recapitulation outputs can be summarized
at pathway, cell type, and gene resolution.
There are two functions to do this:
singISTrecapitulations() computes recapitulations for
one fitted superpathway model.multiple_singISTrecapitulations() is a wrapper that
computes recapitulations for multiple fitted models at once.Below we compute recapitulations for both fitted superpathways
against the oxazolone mapping organism oxa_org, and then
render the results into a plotting friendly output
By default, the function assumes the organism for which the
asmbPLS-DA model was trained is “hsapiens”, it also assumes that the
gene annotation for the mapped organism is “external_gene_name”. These
parameters are by default but other cases can also be considered by
specifying the parameters from_to = "hsapiens" and
object_gene_identifiers = "external_gene_name".
# Compute singIST recapitulations for all fitted superpathways
recapitulations_multiple <- multiple_singISTrecapitulations(oxa_org,
superpathways_fit)In many workflows you will want the outputs in a standarized format
that is convenient for downstream plotting and reporting. The helper
function render_multiple_outputs() reshapes the
recapitulation objects accordingly. It makes a list that is
ready-to-plot.
During CV, singIST may encounter missing entries in the human
pseudobulk predictor (e.g., entire cell type blocks missing for some
samples). Rather than discarding those samples, singIST imputes missing
values inside each CV split in a way that preserves the multiblock
structure and avoids information leakage. This procedure is performed
automatically by fitOptimal() and
multiple_fitOptimal() and the user does not provide any
input on this, this section is merely for informative purposes.
Assume we have a multi-block dataset for which we aim to impute an asmbPLS-DA, \(D = \{\left(x_i, y_i\right) \}_{i=1}^n\) where \(x_i \in \mathbb{R}^p, y_i \in \{0,1\}\). We perform a K-fold CV, thus splitting the dataset into \(K\) folds; \(\mathcal{T}_k\) training indices, \(\mathcal{V}_k\) validation indices, \(k = 1, \dots, K\). Then we have a classic train-test dataset \(\left(X_{train}^{(k)}, Y_{train}^{(k)} \right) \in \{(x_i, y_i): i \in \mathcal{T}_k \}\) and \(\left(X_{val}^{(k)}, Y_{val}^{(k)} \right) \in \{(x_i, y_i): i \in \mathcal{V}_k \}\)
The problem lies in that \(X\) matrices we might have missing values.
We first impute \(X_{train}^{(k)}\)
using Multple Factor Analysis (MFA) imputation (via
missMDA::imputeMFA()), which explicitly accounts for the
multiblock structure (blocks correspond to modeled cell types). This
produces an imputed training matrix \(\widetilde{X}_{train}^{(k)}\).
To avoid data leakage, we do not refit the MFA model on the validation data. Instead, we compute MFA loadings and means from the imputed training fold and project each validation sample using its observed features only.
Let \(\mu \in \mathbb{R}^p\) be the vector of variable means and \(P \in \mathbb{R}^{p \times ncp}\) the matrix of MFA loadings estimated from \(\widetilde{X}_{train}^{(k)}\). For each validation sample \(x_i \in X_{val}^{(k)}\) define the set of observed feature indices:
\[O_i = \{j: (x_i)_j \neq NA\}, \quad |O_i| \geq 1\]
We estimate the latent scores \(t_i\) by solving a least-squares problem restricted to the observed features:
\[t_i = \arg \min_{u_i \in \mathbb{R}^{ncp}} ||x_{i, O_i}-\left(\mu_{O_i}+P_{O_i}u_i \right)||_2^2 = \left(P_{O_i}^TP_{O_i}\right)^{-1} P_{O_i}^T \left(x_{i,O_i}-\mu_{O_i} \right)\]
We then reconstruct the feature full vector by:
\[\widehat{x}_{ij} = \mu_{j} + P_{j\cdot} t_i\]
The imputation validation entries are defined as:
\[ \begin{align} \widehat{x}_{ij} = \begin{cases} x_{ij}, & j \in O_i\\ \widetilde{x}_{ij}, & j \notin O_i \end{cases} \quad i \in \mathcal{V}_k \end{align} \]
This preserves the same latent space used for training-fold imputation and prevents leakage from validation into the imputation model. Imputation is intended for partial missingness (e.g., one missing block while other blocks are observed). If missingness is extensive (e.g., many samples missing most blocks), model stability and interpretability may degrade; in that case consider restricting the set of modeled cell types and/or increasing sample size.
## R version 4.5.2 (2025-10-31)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.3 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so; LAPACK version 3.12.0
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 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: Etc/UTC
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] GSEABase_1.73.0 graph_1.89.1 annotate_1.89.0
## [4] XML_3.99-0.20 AnnotationDbi_1.73.0 IRanges_2.45.0
## [7] S4Vectors_0.49.0 Biobase_2.71.0 BiocGenerics_0.57.0
## [10] generics_0.1.4 msigdb_1.19.0 singIST_0.99.85
## [13] BiocStyle_2.39.0
##
## loaded via a namespace (and not attached):
## [1] matrixStats_1.5.0 spatstat.sparse_3.1-0
## [3] doParallel_1.0.17 httr_1.4.7
## [5] RColorBrewer_1.1-3 tools_4.5.2
## [7] sctransform_0.4.3 backports_1.5.0
## [9] R6_2.6.1 DT_0.34.0
## [11] lazyeval_0.2.2 uwot_0.2.4
## [13] jomo_2.7-6 withr_3.0.2
## [15] sp_2.2-0 prettyunits_1.2.0
## [17] gridExtra_2.3 progressr_0.18.0
## [19] cli_3.6.5 spatstat.explore_3.6-0
## [21] fastDummies_1.7.5 flashClust_1.01-2
## [23] sass_0.4.10 Seurat_5.4.0
## [25] mvtnorm_1.3-3 S7_0.2.1
## [27] spatstat.data_3.1-9 ggridges_0.5.7
## [29] pbapply_1.7-4 parallelly_1.46.1
## [31] limma_3.67.0 RSQLite_2.4.5
## [33] shape_1.4.6.1 ica_1.0-3
## [35] spatstat.random_3.4-3 car_3.1-3
## [37] dplyr_1.1.4 leaps_3.2
## [39] Matrix_1.7-4 abind_1.4-8
## [41] lifecycle_1.0.5 scatterplot3d_0.3-44
## [43] yaml_2.3.12 edgeR_4.9.2
## [45] carData_3.0-5 SummarizedExperiment_1.41.0
## [47] SparseArray_1.11.10 BiocFileCache_3.1.0
## [49] Rtsne_0.17 grid_4.5.2
## [51] blob_1.3.0 promises_1.5.0
## [53] dqrng_0.4.1 ExperimentHub_3.1.0
## [55] mitml_0.4-5 crayon_1.5.3
## [57] miniUI_0.1.2 lattice_0.22-7
## [59] beachmat_2.27.2 cowplot_1.2.0
## [61] KEGGREST_1.51.1 sys_3.4.3
## [63] maketools_1.3.2 pillar_1.11.1
## [65] knitr_1.51 metapod_1.19.1
## [67] GenomicRanges_1.63.1 boot_1.3-32
## [69] estimability_1.5.1 future.apply_1.20.1
## [71] codetools_0.2-20 pan_1.9
## [73] glue_1.8.0 spatstat.univar_3.1-5
## [75] data.table_1.18.0 Rdpack_2.6.4
## [77] vctrs_0.6.5 png_0.1-8
## [79] spam_2.11-3 gtable_0.3.6
## [81] cachem_1.1.0 xfun_0.55
## [83] rbibutils_2.4 S4Arrays_1.11.1
## [85] mime_0.13 Seqinfo_1.1.0
## [87] reformulas_0.4.3.1 survival_3.8-3
## [89] SingleCellExperiment_1.33.0 iterators_1.0.14
## [91] gmp_0.7-5 statmod_1.5.1
## [93] bluster_1.21.0 fitdistrplus_1.2-4
## [95] ROCR_1.0-11 nlme_3.1-168
## [97] bit64_4.6.0-1 progress_1.2.3
## [99] filelock_1.0.3 RcppAnnoy_0.0.23
## [101] bslib_0.9.0 irlba_2.3.5.1
## [103] rpart_4.1.24 KernSmooth_2.23-26
## [105] otel_0.2.0 DBI_1.2.3
## [107] nnet_7.3-20 tidyselect_1.2.1
## [109] emmeans_2.0.1 bit_4.6.0
## [111] compiler_4.5.2 curl_7.0.0
## [113] glmnet_4.1-10 httr2_1.2.2
## [115] missMDA_1.20 BiocNeighbors_2.5.0
## [117] mice_3.19.0 xml2_1.5.1
## [119] DelayedArray_0.37.0 plotly_4.11.0
## [121] checkmate_2.3.3 scales_1.4.0
## [123] lmtest_0.9-40 multcompView_0.1-10
## [125] rappdirs_0.3.3 stringr_1.6.0
## [127] digest_0.6.39 goftest_1.2-3
## [129] minqa_1.2.8 spatstat.utils_3.2-1
## [131] rmarkdown_2.30 XVector_0.51.0
## [133] htmltools_0.5.9 pkgconfig_2.0.3
## [135] lme4_1.1-38 MatrixGenerics_1.23.0
## [137] FactoMineR_2.13 dbplyr_2.5.1
## [139] fastmap_1.2.0 rlang_1.1.7
## [141] htmlwidgets_1.6.4 asmbPLS_1.0.0
## [143] shiny_1.12.1 farver_2.1.2
## [145] jquerylib_0.1.4 zoo_1.8-15
## [147] jsonlite_2.0.0 BiocParallel_1.45.0
## [149] BiocSingular_1.27.1 magrittr_2.0.4
## [151] Formula_1.2-5 scuttle_1.21.0
## [153] dotCall64_1.2 patchwork_1.3.2
## [155] Rcpp_1.1.1 reticulate_1.44.1
## [157] stringi_1.8.7 MASS_7.3-65
## [159] org.Hs.eg.db_3.22.0 AnnotationHub_4.1.0
## [161] plyr_1.8.9 parallel_4.5.2
## [163] listenv_0.10.0 ggrepel_0.9.6
## [165] deldir_2.0-4 Biostrings_2.79.4
## [167] splines_4.5.2 tensor_1.5.1
## [169] RcppAlgos_2.9.3 hms_1.1.4
## [171] locfit_1.5-9.12 igraph_2.2.1
## [173] ggpubr_0.6.2 spatstat.geom_3.6-1
## [175] ggsignif_0.6.4 RcppHNSW_0.6.0
## [177] buildtools_1.0.0 reshape2_1.4.5
## [179] biomaRt_2.67.1 ScaledMatrix_1.19.0
## [181] BiocVersion_3.23.1 evaluate_1.0.5
## [183] SeuratObject_5.3.0 scran_1.39.0
## [185] BiocManager_1.30.27 nloptr_2.2.1
## [187] foreach_1.5.2 httpuv_1.6.16
## [189] RANN_2.6.2 tidyr_1.3.2
## [191] purrr_1.2.1 polyclip_1.10-7
## [193] future_1.68.0 scattermore_1.2
## [195] ggplot2_4.0.1 rsvd_1.0.5
## [197] broom_1.0.11 xtable_1.8-4
## [199] RSpectra_0.16-2 rstatix_0.7.3
## [201] later_1.4.5 viridisLite_0.4.2
## [203] snow_0.4-4 tibble_3.3.1
## [205] memoise_2.0.1 cluster_2.1.8.1
## [207] globals_0.18.0