Contents

1 Introduction

“APL” is a package developed for computation of Association Plots, a method for visualization and analysis of single cell transcriptomics data. The main focus of “APL” is the identification of genes characteristic for individual clusters of cells from input data.

When working with APL package please cite:

Gralinska, E., Kohl, C., Fadakar, B. S., & Vingron, M. (2022). 
Visualizing Cluster-specific Genes from Single-cell Transcriptomics Data Using Association Plots. 
Journal of Molecular Biology, 434(11), 167525.

A citation can also be obtained in R by running citation("APL"). For a mathematical description of the method, please refer to the manuscript.

2 Installation

To install the APL from Bioconductor, run:

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("APL")

Alternatively the package can also be installed from GitHub:

library(devtools)
install_github("VingronLab/APL")

To additionally build the package vignette, run instead

install_github("VingronLab/APL", build_vignettes = TRUE, dependencies = TRUE)

Building the vignette will however take considerable time.

2.1 pytorch installation

In order to decrease the computation time of the singular value decomposition (SVD), we highly recommend the installation of pytorch. More information on the pytorch installation is given below. Instead of installing pytorch, users can also opt to use the R native SVD. For this, please use the argument python = FALSE wherever applicable in this vignette.

2.1.1 Install pytorch with reticulate

library(reticulate)
install_miniconda()
use_condaenv(condaenv = file.path(miniconda_path(),"envs/r-reticulate"),
             required=TRUE)
      
conda_install(envname = "r-reticulate", packages = "numpy")
conda_install(envname = "r-reticulate", packages = "pytorch")

2.1.2 Manually install pytorch with conda

To install pytorch please download the appropriate Miniconda installer for your system from the conda website. Follow the installation instructions on their website and make sure the R package reticulate is also installed before proceeding. Once installed, list all available conda environments via
conda info --envs
One of the environments should have r-reticulate in its name. Depending on where you installed it and your system, the exact path might be different. Activate the environment and install pytorch into it.

conda activate ~/.local/share/r-miniconda/envs/r-reticulate # change path accordingly.
conda install numpy
conda install pytorch

3 Preprocessing

3.1 Setup

In this vignette we will use a small data set published by Darmanis et al. (2015) consisting of 466 human adult cortical single cells sequenced on the Fluidigm platform as an example. To obtain the data necessary to follow the vignette we use the Bioconductor package scRNAseq.

Besides the package APL we will use Bioconductor packages to preprocess the data. Namely we will use SingleCellExperiment, scater and scran. However, the preprocessing could equally be performed with the single-cell RNA-seq analysis suite Seurat.

The preprocessing steps are performed according to the recommendations published in Orchestrating Single-Cell Analysis with Bioconductor by Amezquita et al. (2022). For more information about the rational behind them please refer to the book.

library(APL)
library(scRNAseq)
library(SingleCellExperiment)
library(scran)
library(scater)
set.seed(1234)

3.2 Loading the data

We start with the loading and preprocessing of the Darmanis data.

darmanis <- DarmanisBrainData()
#> see ?scRNAseq and browseVignettes('scRNAseq') for documentation
#> loading from cache
#> see ?scRNAseq and browseVignettes('scRNAseq') for documentation
#> loading from cache
darmanis
#> class: SingleCellExperiment 
#> dim: 22085 466 
#> metadata(0):
#> assays(1): counts
#> rownames(22085): 1/2-SBSRNA4 A1BG ... ZZZ3 tAKR
#> rowData names(0):
#> colnames(466): GSM1657871 GSM1657872 ... GSM1658365 GSM1658366
#> colData names(5): age c1.chip.id cell.type experiment_sample_name
#>   tissue
#> reducedDimNames(0):
#> mainExpName: NULL
#> altExpNames(0):

3.3 Normalization, PCA & Clustering

Association Plots from APL should be computed based on the normalized expression data. Therefore, we first normalize the counts from the Darmanis data and calculate both PCA and UMAP for visualizations later.

For now, APL requires the data to be clustered beforehand. The darmanis data comes already annotated, so we will use the cell types stored in the cell.type metadata column instead of performing a clustering.


set.seed(100)
clust <- quickCluster(darmanis) 
darmanis <- computeSumFactors(darmanis, cluster=clust, min.mean=0.1)
darmanis <- logNormCounts(darmanis)

dec <- modelGeneVar(darmanis)
top_darmanis <- getTopHVGs(dec, n=5000)
darmanis <- fixedPCA(darmanis, subset.row=top_darmanis) 
darmanis <- runUMAP(darmanis, dimred="PCA")
#> Found more than one class "dist" in cache; using the first, from namespace 'BiocGenerics'
#> Also defined by 'spam'
#> Found more than one class "dist" in cache; using the first, from namespace 'BiocGenerics'
#> Also defined by 'spam'

plotReducedDim(darmanis, dimred="UMAP", colour_by="cell.type")

4 Quick start

The fastest way to compute the Association Plot for a selected cluster of cells from the input data is by using a wrapper function runAPL(). runAPL() automates most of the analysis steps for ease of use.

For example, to generate an Association Plot for the oligodendrocytes we can use the following command:

runAPL(darmanis,
       assay = "logcounts",
       top = 5000,
       group = which(darmanis$cell.type == "oligodendrocytes"),
       type = "ggplot",
       python = TRUE)
#> Warning in rm_zeros(obj): Matrix contains rows with only 0s. These rows were
#> removed. If undesired set rm_zeros = FALSE.
#> 
#>  Using 74 dimensions. Subsetting.

#> 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |======================================================================| 100%
#> 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |=======================                                               |  33%
  |                                                                            
  |===============================================                       |  67%
  |                                                                            
  |======================================================================| 100%

The generated Association Plot is computed based on the log-normalized count matrix. By default runAPL uses the top 5,000 most variable genes in the data, but the data can be subset to any number of genes by changing the value for the argument top. The dimensionality of the CA is determined automatically by the elbow rule described below (see here). This default behavior can be overriden by setting the dimensions manually (parameter dims). The cluster-specificity score (\(S_\alpha\)) for each gene is also calculated (score = TRUE). In order to better explore the data, type can be set to "plotly" to obtain an interactive plot. runAPL has many arguments to further customize the output and fine tune the calculations. Please refer to the documentation (?runAPL) for more information. The following sections in this vignette will discuss the choice of dimensionality and the \(S_\alpha\)-score.

5 Step-by-step way of computing Association Plots

Alternatively, Association Plots can be computed step-by-step. This allows to adjust the Association Plots to user’s needs. Below we explain each step of the process of generating Association Plots.

5.1 Correspondence Analysis

The first step of Association Plot computations is correspondence analysis (CA). CA is a data dimensionality reduction method similar to PCA, however it allows for a simultaneous embedding of both cells and genes from the input data in the same space. In this example we perform CA on the log-normalized count matrix of the darmanis brain data.

# Computing CA on logcounts
logcounts <- logcounts(darmanis)
ca <- cacomp(obj = logcounts,
             top = 5000,
             python = TRUE)
#> Warning in rm_zeros(obj): Matrix contains rows with only 0s. These rows were
#> removed. If undesired set rm_zeros = FALSE.

# The above is equivalent to:
# ca <- cacomp(obj = darmanis,
#              assay = "logcounts",
#              top = 5000,
#              python = TRUE)

The function cacomp accepts as an input any matrix with non-negative entries, be it a single-cell RNA-seq, bulk RNA-seq or other data. For ease of use, cacomp accepts also SingleCellExperiment and Seurat objects, however for these we additionally have to specify via the assay and/or slot (for Seurat) parameter from where to extract the data. Importantly, in order to ensure the interpretability of the results cacomp (and related functions such as runAPL) requires that the input matrix contains both row and column names.

When performing a feature selection before CA, we can set the argument top to the desired number of genes with the highest variance across cells from the input data to retain for further analysis. By default, only the top 5,000 most variable genes are kept as a good compromise between computational time and keeping the most relevant genes. If we want to ensure however that even marker genes of smaller clusters are kept, we can increase the number of genes.

The output of cacomp is an object of class cacomp:

ca
#> cacomp object with 466 columns, 5000 rows and 465 dimensions.
#> Calc. standard coord.:  std_coords_rows, std_coords_cols
#> Calc. principal coord.: prin_coords_rows, prin_coords_cols
#> Calc. APL coord.:       
#> Explained inertia:      1.1% Dim1, 0.7% Dim2

As can be seen in the summarized output, by default both types of coordinates in the CA space (principal and standardized) are calculated. Once the coordinates for the Association Plot are calculated, they will also be shown in the output of cacomp. Slots are accessed through an accessor function:

cacomp_slot(ca, "std_coords_cols")[1:5, 1:5]
#>                  Dim1       Dim2       Dim3        Dim4       Dim5
#> GSM1657871 -0.1067086 -0.1253297  0.5046485 -0.32461780  2.4262520
#> GSM1657872  1.5101357 -0.4299418 -0.1219273  0.07942823  0.2362747
#> GSM1657873 -0.2237619 -0.2610148  0.1955599  0.16963578  2.1477040
#> GSM1657874 -0.5680872 -0.1138251  0.4725071  0.04409559 -0.6708562
#> GSM1657875 -0.4739344 -0.2505648  0.4384626  0.10316896  2.8910899

In the case of SingleCellExperiment and Seurat objects, we can alternatively set return_input = TRUE to get the input object back, with the CA results computed by “APL” and stored in the appropriate slot for dimension reduction. This also allows for using the plotting functions that come with these packages:

darmanis <- cacomp(obj = darmanis,
               assay = "logcounts",
               top = 5000,
               return_input = TRUE,
               python = TRUE)
#> Warning in rm_zeros(obj): Matrix contains rows with only 0s. These rows were
#> removed. If undesired set rm_zeros = FALSE.

plotReducedDim(darmanis, dimred="CA",
               ncomponents = c(1,2),
               colour_by="cell.type")

plotReducedDim(darmanis, dimred="CA",
               ncomponents = c(3,4),
               colour_by="cell.type")

However, some functions such as apl_coords() require information that cannot be stored in the single-cell container objects. It is therefore often easier to work with a cacomp object instead. We can convert Seurat or SingleCellExperiment objects which have CA results stored to a cacomp object using the function as.cacomp():

# Converting the object pbmc to cacomp
ca <- as.cacomp(darmanis)

5.2 Reducing the number of CA dimensions

When working with high-dimensional data, after singular value decomposition there will often be many dimensions which are representing the noise in the data. In order to minimize the noise, it is generally recommended to reduce the dimensionsionality of the data before generating Association Plots.

The number of dimensions to retain can be computed using the function pick_dims. This function offers three standard methods which we implemented:

  • elbow rule (method = "elbow_rule") - the number of dimensions to retain is calculated based on scree plots generated for randomized data, and corresponds to a point in the plot where the band of randomized singular values enters the band of the original singular values,

  • 80% rule (method = "maj_inertia") - only those first dimensions are retained which in total account for >= 80% of total inertia,

  • average rule (method = "avg_inertia") - only those dimensions are retained which account for more inertia than a single dimension on average.

Additionally, the user can compute a scree plot to choose the number of dimensions by themselves:

pick_dims(ca, method = "scree_plot") +
  xlim(c(0,20))
#> Warning: Removed 445 rows containing missing values (`position_stack()`).
#> Warning: Removed 1 rows containing missing values (`geom_col()`).
#> Warning: Removed 445 rows containing missing values (`geom_line()`).

In the scree plot above we can see that the first dimension explains only ~1% of the total inertia and we observe the “jump” in the scree plot at roughly 5 dimensions. The first dimensions however explain only a small amount of the total inertia.

Here we compute the number of dimensions using the elbow rule. For demonstration, only three data permutations are computed:

pd <- pick_dims(ca,
                mat = logcounts(darmanis),
                method = "elbow_rule",
                reps = 3, 
                python = TRUE)
pd
#> [1] 74

In this case the elbow rule leads to a much higher number of dimensions.

# Compute the amount of inertia explained by each of the dimensions
D <- cacomp_slot(ca, "D")
expl_inertia <- (D^2/sum(D^2))*100

# Compute the amount of intertia explained 
# by the number of dimensions defined by elbow rule
sum(expl_inertia[seq_len(pd)])
#> [1] 28.00469

In this example the elbow rule suggests to keep 74 dimensions that explain 28% of the total inertia from the data.

Finally, we can reduce the dimensionality of the data to the desired number of dimensions:

ca <- subset_dims(ca, dims = pd)

5.3 Association Plots

When working with single-cell transcriptomics data we are often interested in which genes are associated to a cluster of cells. To reveal such genes we can compute an Association Plot for a selected cluster of cells. In the following example we want to generate an Association Plot for the cluster of endothelial cells:

# Specifying a cell cluster of interest
endo <- which(darmanis$cell.type == "endothelial")

# Calculate Association Plot coordinates for endothelial cells
ca <- apl_coords(ca, group = endo)

After computing the coordinates of genes and cells in the Association Plot we are able to plot the results using the apl function.


# endothelial marker genes
marker_genes <- c("APOLD1", "TM4SF1", "SULT1B1", "ESM1", "SELE")

# Plot APL
apl(ca,
    row_labs = TRUE,
    rows_idx = marker_genes, 
    type = "ggplot") # type = "plotly" for an interactive plot

In the Association Plot all genes are represented by blue circles. The further to the right a gene is located the more associated it is with the chosen cluster of cells and the lower the y-axis value, the more specific it is for the selected cluster. Additionally, it is possible to highlight in the Association Plot any set of genes. In the example above we highlighted five genes (APOLD1, TM4SF1, SULT1B1, ESM1, SELE) which are known to be marker genes for endothelial cells. As we can see in the plot, they are located in the right part of the plot, which confirms their specificity for endothelial cells.

By default we plot only the genes in the Association Plot. To also display the cells in the Association Plot, use the argument show_cols = TRUE. This way we can identify other cells which show similar expression profiles to the cells of interest. Cells that belong to the cluster of interest will be colored in red, and all remaining cells will be colored in violet. Furthermore, an interactive plot in which you can hover over genes to see their name can be created by setting type = "plotly".

5.4 Association Plots with the \(S_\alpha\)-scores

The \(S_\alpha\)-score allows us to rank genes by their specificity for a selected cell cluster, and is computed for each gene from the Association Plot separately. The higher the \(S_\alpha\)-score of a gene, the more characteristic its expression for the investigated cell cluster. The \(S_\alpha\)-scores can be computed using the apl_score function. To display the \(S_\alpha\)-scores in the Association Plot, we can use the argument show_score = TRUE in the apl function:

# Compute S-alpha score
# For the calculation the input matrix is also required.
ca <- apl_score(ca,
                mat = logcounts(darmanis),
                reps = 5,
                python = TRUE)
apl(ca,
    show_score = TRUE,
    type = "ggplot") 

By default, only genes that have a \(S_\alpha\)-score larger than 0 are colored as these tend to be genes of interest and we consider them as cluster-specific genes. This cutoff can be easily changed through the score_cutoff argument to apl().

The \(S_\alpha\)-scores are stored in the "APL_score" slot and can be accessed as follows:

head(cacomp_slot(ca, "APL_score"))
#>             Rowname    Score Row_num Rank
#> PRKCDBP     PRKCDBP 2.816535    2390    1
#> IFITM1       IFITM1 2.725801    2081    2
#> LOC646268 LOC646268 2.633763      81    3
#> OR2A20P     OR2A20P 2.633763    4888    4
#> OR4C46       OR4C46 2.633763    4889    5
#> FOXC1         FOXC1 2.493272    2937    6

To see the expression of genes with the highest \(S_\alpha\)-scores (or any selected genes) across all cell types from the data we can use plotting functions provided by scater:

scores <- cacomp_slot(ca, "APL_score")

plotExpression(darmanis, 
               features = head(scores$Rowname,3),
               x = "cell.type", 
               colour_by = "cell.type")


plotReducedDim(darmanis, 
               dimred="UMAP",
               colour_by= scores$Rowname[2])

As expected, the 3 most highly scored genes are over-expressed in the endothelial cells. Due to the small size of the data set and number of cells in the cluster (only 20 out of 466 cells are endothelial cells) some cluster specific genes are only expressed in a few cells. Most data sets nowadays are significantly larger so this should not be a major concern and it can further be mitigated by performing a more stringent feature selection before CA.

5.5 Visualization of CA

In addition to Association Plots “APL” produces also other forms of the output. For instance, we can use “APL” to generate a two- and three-dimensional correspondence analysis projection of the data. The so-called biplot visualizes both cells and genes from the input data and can be created using the function ca_biplot. Alternatively, a three-dimensional data projection plot can be generated using the function ca_3Dplot. To generate such biplots a cacomp object is required.


# Specifying a cell cluster of interest
endo <- which(darmanis$cell.type == "endothelial")

# Creating a static plot
ca_biplot(ca, col_labels = endo, type = "ggplot")


# Creating an interactive plot
# ca_biplot(ca, type = "plotly", col_labels = platelets)

# 3D plot
# ca_3Dplot(ca, col_labels = platelets)

The above described plots give us a quick overview of the first 2 dimensions of the data (more dimensions can be plotted). As shown in the commented-out code, to interactively explore the projection of the data type = "plotly" can be set.

6 APL and GO enrichment analysis

After computing an Association Plot and identifying a set of genes specific for a selected cluster of cells we might be interested in conducting a Gene Ontology (GO) enrichment analysis of the identified gene set. To conduct a GO enrichment analysis of microglia specific genes as idenitfied using an Association Plot, we first need to compute the coordinates of the genes in the Association Plot for microglia cells, as well as the \(S_\alpha\)-score for each gene:

# Get indices of microglia cells
microglia <- which(darmanis$cell.type == "microglia")

# Calculate Association Plot coordinates of the genes and the $S_\alpha$-scores
ca <- apl_coords(ca, group = microglia)

ca <- apl_score(ca,
                mat = logcounts(darmanis),
                reps = 5,
                python = TRUE)

Now we can conduct GO enrichment analysis as implemented in the package topGO using the most cluster-specific genes from the Association Plot. By default we use all genes with an \(S_\alpha\)-score higher than 0, but the cutoff may have to be adjusted depending on the dataset. In the example below we restrict it to genes with a \(S_\alpha\)-score higher than 1 to restrict it to truly significant genes. In case that no \(S_\alpha\)-scores were calculated, one can also choose to use the ngenes (by default 1000) genes with the highest x-coordinates by setting use_coords = TRUE.

enr <- apl_topGO(ca,
                 ontology = "BP",
                 organism = "hs",
                 score_cutoff = 1)
head(enr)
#>        GO.ID
#> 1 GO:0045087
#> 2 GO:0019886
#> 3 GO:0032755
#> 4 GO:0032760
#> 5 GO:0002250
#> 6 GO:0002503
#>                                                                                Term
#> 1                                                            innate immune response
#> 2 antigen processing and presentation of exogenous peptide antigen via MHC class II
#> 3                                   positive regulation of interleukin-6 production
#> 4                           positive regulation of tumor necrosis factor production
#> 5                                                          adaptive immune response
#> 6                        peptide antigen assembly with MHC class II protein complex
#>   Annotated Significant Expected raw.p.value
#> 1       278          57    17.67     6.5e-12
#> 2        15          11     0.95     5.8e-11
#> 3        39          15     2.48     4.4e-09
#> 4        45          16     2.86     5.0e-09
#> 5       166          41    10.55     2.5e-08
#> 6        11           8     0.70     3.3e-08

The function plot_enrichment() was implemented to visualize the topGO results in form of a dotplot.

plot_enrichment(enr)

Microglia cells are innate immune cells of the brain and as such the most highly scored genes are enriched in gene sets related to the immune response and microglia specific gene sets as one would expect.

Session info

#> R version 4.3.2 Patched (2023-11-13 r85521)
#> Platform: x86_64-pc-linux-gnu (64-bit)
#> Running under: Ubuntu 22.04.3 LTS
#> 
#> Matrix products: default
#> BLAS:   /home/biocbuild/bbs-3.18-bioc/R/lib/libRblas.so 
#> LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0
#> 
#> 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] SparseM_1.81                org.Hs.eg.db_3.18.0        
#>  [3] AnnotationDbi_1.64.1        scater_1.30.1              
#>  [5] ggplot2_3.4.4               scran_1.30.2               
#>  [7] scuttle_1.12.0              scRNAseq_2.16.0            
#>  [9] SingleCellExperiment_1.24.0 SummarizedExperiment_1.32.0
#> [11] Biobase_2.62.0              GenomicRanges_1.54.1       
#> [13] GenomeInfoDb_1.38.6         IRanges_2.36.0             
#> [15] S4Vectors_0.40.2            BiocGenerics_0.48.1        
#> [17] MatrixGenerics_1.14.0       matrixStats_1.2.0          
#> [19] APL_1.6.1                   BiocStyle_2.30.0           
#> 
#> loaded via a namespace (and not attached):
#>   [1] ProtGenerics_1.34.0           spatstat.sparse_3.0-3        
#>   [3] bitops_1.0-7                  httr_1.4.7                   
#>   [5] RColorBrewer_1.1-3            tools_4.3.2                  
#>   [7] sctransform_0.4.1             utf8_1.2.4                   
#>   [9] R6_2.5.1                      lazyeval_0.2.2               
#>  [11] uwot_0.1.16                   withr_3.0.0                  
#>  [13] sp_2.1-3                      prettyunits_1.2.0            
#>  [15] gridExtra_2.3                 progressr_0.14.0             
#>  [17] cli_3.6.2                     spatstat.explore_3.2-6       
#>  [19] fastDummies_1.7.3             labeling_0.4.3               
#>  [21] sass_0.4.8                    topGO_2.54.0                 
#>  [23] Seurat_5.0.1                  spatstat.data_3.0-4          
#>  [25] ggridges_0.5.6                pbapply_1.7-2                
#>  [27] Rsamtools_2.18.0              parallelly_1.36.0            
#>  [29] limma_3.58.1                  RSQLite_2.3.5                
#>  [31] FNN_1.1.4                     generics_0.1.3               
#>  [33] BiocIO_1.12.0                 ica_1.0-3                    
#>  [35] spatstat.random_3.2-2         dplyr_1.1.4                  
#>  [37] GO.db_3.18.0                  Matrix_1.6-5                 
#>  [39] ggbeeswarm_0.7.2              fansi_1.0.6                  
#>  [41] abind_1.4-5                   lifecycle_1.0.4              
#>  [43] edgeR_4.0.14                  yaml_2.3.8                   
#>  [45] SparseArray_1.2.3             BiocFileCache_2.10.1         
#>  [47] Rtsne_0.17                    grid_4.3.2                   
#>  [49] blob_1.2.4                    dqrng_0.3.2                  
#>  [51] promises_1.2.1                ExperimentHub_2.10.0         
#>  [53] crayon_1.5.2                  miniUI_0.1.1.1               
#>  [55] lattice_0.22-5                beachmat_2.18.0              
#>  [57] cowplot_1.1.3                 GenomicFeatures_1.54.3       
#>  [59] KEGGREST_1.42.0               magick_2.8.2                 
#>  [61] metapod_1.10.1                pillar_1.9.0                 
#>  [63] knitr_1.45                    rjson_0.2.21                 
#>  [65] future.apply_1.11.1           codetools_0.2-19             
#>  [67] leiden_0.4.3.1                glue_1.7.0                   
#>  [69] data.table_1.15.0             vctrs_0.6.5                  
#>  [71] png_0.1-8                     spam_2.10-0                  
#>  [73] org.Mm.eg.db_3.18.0           gtable_0.3.4                 
#>  [75] cachem_1.0.8                  xfun_0.42                    
#>  [77] S4Arrays_1.2.0                mime_0.12                    
#>  [79] survival_3.5-7                statmod_1.5.0                
#>  [81] bluster_1.12.0                interactiveDisplayBase_1.40.0
#>  [83] ellipsis_0.3.2                fitdistrplus_1.1-11          
#>  [85] ROCR_1.0-11                   nlme_3.1-164                 
#>  [87] bit64_4.0.5                   progress_1.2.3               
#>  [89] filelock_1.0.3                RcppAnnoy_0.0.22             
#>  [91] bslib_0.6.1                   irlba_2.3.5.1                
#>  [93] vipor_0.4.7                   KernSmooth_2.23-22           
#>  [95] colorspace_2.1-0              DBI_1.2.1                    
#>  [97] tidyselect_1.2.0              bit_4.0.5                    
#>  [99] compiler_4.3.2                curl_5.2.0                   
#> [101] graph_1.80.0                  BiocNeighbors_1.20.2         
#> [103] xml2_1.3.6                    DelayedArray_0.28.0          
#> [105] plotly_4.10.4                 bookdown_0.37                
#> [107] rtracklayer_1.62.0            scales_1.3.0                 
#> [109] lmtest_0.9-40                 rappdirs_0.3.3               
#> [111] stringr_1.5.1                 digest_0.6.34                
#> [113] goftest_1.2-3                 spatstat.utils_3.0-4         
#> [115] rmarkdown_2.25                XVector_0.42.0               
#> [117] htmltools_0.5.7               pkgconfig_2.0.3              
#> [119] sparseMatrixStats_1.14.0      highr_0.10                   
#> [121] dbplyr_2.4.0                  fastmap_1.1.1                
#> [123] ensembldb_2.26.0              rlang_1.1.3                  
#> [125] htmlwidgets_1.6.4             shiny_1.8.0                  
#> [127] DelayedMatrixStats_1.24.0     farver_2.1.1                 
#> [129] jquerylib_0.1.4               zoo_1.8-12                   
#> [131] jsonlite_1.8.8                BiocParallel_1.36.0          
#> [133] BiocSingular_1.18.0           RCurl_1.98-1.14              
#> [135] magrittr_2.0.3                GenomeInfoDbData_1.2.11      
#> [137] dotCall64_1.1-1               patchwork_1.2.0              
#> [139] munsell_0.5.0                 Rcpp_1.0.12                  
#> [141] viridis_0.6.5                 reticulate_1.35.0            
#> [143] stringi_1.8.3                 zlibbioc_1.48.0              
#> [145] MASS_7.3-60.0.1               AnnotationHub_3.10.0         
#> [147] plyr_1.8.9                    parallel_4.3.2               
#> [149] listenv_0.9.1                 ggrepel_0.9.5                
#> [151] deldir_2.0-2                  Biostrings_2.70.2            
#> [153] splines_4.3.2                 tensor_1.5                   
#> [155] hms_1.1.3                     locfit_1.5-9.8               
#> [157] igraph_2.0.1.1                spatstat.geom_3.2-8          
#> [159] RcppHNSW_0.6.0                reshape2_1.4.4               
#> [161] biomaRt_2.58.2                ScaledMatrix_1.10.0          
#> [163] BiocVersion_3.18.1            XML_3.99-0.16.1              
#> [165] evaluate_0.23                 SeuratObject_5.0.1           
#> [167] BiocManager_1.30.22           httpuv_1.6.14                
#> [169] RANN_2.6.1                    tidyr_1.3.1                  
#> [171] purrr_1.0.2                   polyclip_1.10-6              
#> [173] future_1.33.1                 scattermore_1.2              
#> [175] rsvd_1.0.5                    xtable_1.8-4                 
#> [177] restfulr_0.0.15               AnnotationFilter_1.26.0      
#> [179] RSpectra_0.16-1               later_1.3.2                  
#> [181] viridisLite_0.4.2             tibble_3.2.1                 
#> [183] beeswarm_0.4.0                memoise_2.0.1                
#> [185] GenomicAlignments_1.38.2      cluster_2.1.6                
#> [187] globals_0.16.2