scPAS : Single-Cell Phenotype-Associated Subpopulation identifier
Aimin Xie
2024-08-25
Introduction
This scPAS R package is a bioinformatics tool designed to identify phenotype-associated cell subpopulations in single-cell data by integrating bulk data. This approach provides both quantitative and qualitative estimates of the association between cells and phenotypes. The core function of this package is scPAS. It requires three input data: single-cell RNA-seq data, bulk-seq data, and the phenotype of patients in bulk samples.
Installation
The scPAS is developed under R (version >= 4.0.5). scPAS R package can be installed from Github using devtools:
devtools::install_github("aiminXie/scPAS")If it cannot be installed automatically, install the following dependencies:
install.packages('Seurat')
install.packages('Rcpp')
install.packages('Matrix')
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("preprocessCore")library(scPAS)
library(Matrix)
library(Seurat)Apply scPAS with Cox regression
In this example, we used TCGA-LIHC bulk samples with survival information to identify survival-related cell subpopulations within the LIHC scRNA-seq data.
Load data
The data used in this tutorial is stored at: https://drive.google.com/drive/folders/1Qgdx7YWSkJJ-ZjYHFt9tD33ZC85vQk8o?usp=sharing
load(file = './data/data_survival.rda')
dim(sc_dataset)
#> [1] 21287 8853
head(bulk_dataset[,1:10])
#> TCGA-2Y-A9GS-01 TCGA-2Y-A9GT-01 TCGA-2Y-A9GU-01 TCGA-2Y-A9GV-01
#> ARHGEF10L 11.5199 12.2544 12.6434 11.9114
#> HIF3A 2.4425 4.6559 1.8083 3.7751
#> RNF17 3.7881 0.0000 0.0000 0.0000
#> RNF10 11.5109 11.7863 12.3187 11.7253
#> RNF11 10.6574 11.0133 10.5722 10.9586
#> RNF13 10.2623 10.4574 9.9284 10.5863
#> TCGA-2Y-A9GW-01 TCGA-2Y-A9GX-01 TCGA-2Y-A9GY-01 TCGA-2Y-A9GZ-01
#> ARHGEF10L 11.8108 10.6616 10.6321 12.4218
#> HIF3A 2.8340 3.2456 7.1429 5.7972
#> RNF17 0.0000 0.0000 7.9222 0.7507
#> RNF10 11.8382 11.5872 11.5770 11.2699
#> RNF11 10.3071 10.7500 9.9313 10.7820
#> RNF13 10.1785 10.4880 10.0815 10.5819
#> TCGA-2Y-A9H0-01 TCGA-2Y-A9H1-01
#> ARHGEF10L 10.8164 10.2715
#> HIF3A 3.7561 4.3687
#> RNF17 0.0000 0.0000
#> RNF10 12.7789 11.2257
#> RNF11 10.0653 10.1549
#> RNF13 10.6751 9.9593
head(phenotype)
#> time status
#> TCGA-2Y-A9GS-01 724 1
#> TCGA-2Y-A9GT-01 1624 1
#> TCGA-2Y-A9GU-01 1939 0
#> TCGA-2Y-A9GV-01 2532 1
#> TCGA-2Y-A9GW-01 1271 1
#> TCGA-2Y-A9GX-01 2442 0Prepare the scRNA-seq data
Single-cell data needs to be processed using the Seurat pipeline to create a Seurat object suitable for visualization. To avoid the impact of extremely sparsely expressed genes on the model, they need to be removed beforehand.
keep_genes <- rownames(sc_dataset)[rowSums(sc_dataset@assays$RNA@counts>1) >= 20]
sc_dataset <- subset(sc_dataset,features = keep_genes)
sc_dataset <- run_Seurat(sc_dataset,verbose = FALSE)
library(RColorBrewer)
qual_col_pals = brewer.pal.info[brewer.pal.info$category == 'qual',]
col_vector = unlist(mapply(brewer.pal, qual_col_pals$maxcolors, rownames(qual_col_pals)))
set.seed(seed = 1234)
cellType_col <- sample(col_vector, 7)
names(cellType_col) <- levels(sc_dataset$celltype)
UMAP_celltype <- DimPlot(sc_dataset, reduction ="umap",group.by="celltype",label =T, cols = cellType_col)
UMAP_celltype
Prepare the bulk data and phenotype
The bulk data required by scPAS is a matrix. Similarly, genes that are not expressed in the majority of samples in bulk data will also be removed.
bulk_dataset <- as.matrix(bulk_dataset)
class(bulk_dataset)
#> [1] "matrix" "array"
dim(bulk_dataset)
#> [1] 18807 370
bulk_dataset = bulk_dataset[apply(bulk_dataset, 1, function(x)sum(x==0) < 0.25 *ncol(bulk_dataset)),]
dim(bulk_dataset)
#> [1] 15874 370The clinical survival phenotype requires a two-column matrix with columns named ‘time’ and ‘status’. The first column contains the survival time. The second column is a binary variable, with ‘1’ indicating event (e.g.recurrence of cancer or death), and ‘0’ indicating right-censored.
head(phenotype)
#> time status
#> TCGA-2Y-A9GS-01 724 1
#> TCGA-2Y-A9GT-01 1624 1
#> TCGA-2Y-A9GU-01 1939 0
#> TCGA-2Y-A9GV-01 2532 1
#> TCGA-2Y-A9GW-01 1271 1
#> TCGA-2Y-A9GX-01 2442 0The phenotypic data must match the samples of the bulk expression profile.
all(rownames(phenotype)==colnames(bulk_dataset))
#> [1] TRUERun scPAS without imputation
The input provided above is survival data; therefore, scPAS needs to fit a Cox regression model (family = ‘cox’).
imputation = F , without imputation.
nfeature = 3000 indicates that the top 3000 highly variable genes are selected for model training.
network_class = ‘SC’ indicates gene-gene similarity networks derived from single-cell data.
time.start <- Sys.time()
sc_dataset <- scPAS(bulk_dataset = bulk_dataset,sc_dataset = sc_dataset,assay = 'RNA',phenotype = phenotype,imputation = F, nfeature = 3000,alpha = 0.01,network_class = 'SC',family = 'cox')
#> [1] "Step 1:Quantile normalization of bulk data."
#> [1] "Step 2: Extracting single-cell expression profiles...."
#> [1] "Step 3: Constructing a gene-gene similarity by single cell data...."
#> [1] "Step 4: Optimizing the network-regularized sparse regression model...."
#> [1] "Perform cox regression on the given clinical outcomes:"
#> [1] "alpha = 0.01"
#> [1] "lambda = 3.57200645154233"
#> [1] "scPAS identified 237 rick+ features and 228 rick- features."
#> [1] "The percentage of selected feature is: 17.521%"
#>
#> [1] "|**************************************************|"
#> [1] "Step 5: calculating quantified risk scores...."
#> [1] "Step 6: qualitative identification by permutation test program with 2000 times random perturbations"
#> [1] "Finished."
time.end <- Sys.time()
time.end-time.start
#> Time difference of 5.911283 minsscPAS will return a Seurat object, and the predicted results are added to the metadata.
head(sc_dataset@meta.data)
#> orig.ident nCount_RNA nFeature_RNA celltype
#> single_set1_AAACCTGAGGCGTACA-1 single_set1 2193 867 CAF
#> single_set1_AAACGGGAGATCGATA-1 single_set1 3501 1213 CAF
#> single_set1_AAAGCAAAGATCGGGT-1 single_set1 7609 2360 CAF
#> single_set1_AAATGCCGTCTCAACA-1 single_set1 3116 1052 CAF
#> single_set1_AACACGTCACGGCTAC-1 single_set1 1980 1036 TEC
#> single_set1_AACCGCGAGACGCTTT-1 single_set1 9948 2691 CAF
#> percent.mt RNA_snn_res.0.6 seurat_clusters
#> single_set1_AAACCTGAGGCGTACA-1 5.861298 5 5
#> single_set1_AAACGGGAGATCGATA-1 3.769339 5 5
#> single_set1_AAAGCAAAGATCGGGT-1 5.173526 12 12
#> single_set1_AAATGCCGTCTCAACA-1 4.671717 5 5
#> single_set1_AACACGTCACGGCTAC-1 1.813725 4 4
#> single_set1_AACCGCGAGACGCTTT-1 2.696030 5 5
#> scPAS_RS scPAS_NRS scPAS_Pvalue scPAS_FDR
#> single_set1_AAACCTGAGGCGTACA-1 -0.2325927711 -2.4165031 0.007835197 0.06021268
#> single_set1_AAACGGGAGATCGATA-1 -0.1546594066 -1.7643536 0.038836229 0.14452342
#> single_set1_AAAGCAAAGATCGGGT-1 0.2166658461 1.3943235 0.081609976 0.21313383
#> single_set1_AAATGCCGTCTCAACA-1 0.0823817980 0.7093568 0.239051551 0.36730204
#> single_set1_AACACGTCACGGCTAC-1 -0.0004370191 -0.2546144 0.399510460 0.45672341
#> single_set1_AACCGCGAGACGCTTT-1 0.0611285679 0.2426248 0.404148050 0.45947383
#> scPAS
#> single_set1_AAACCTGAGGCGTACA-1 0
#> single_set1_AAACGGGAGATCGATA-1 0
#> single_set1_AAAGCAAAGATCGGGT-1 0
#> single_set1_AAATGCCGTCTCAACA-1 0
#> single_set1_AACACGTCACGGCTAC-1 0
#> single_set1_AACCGCGAGACGCTTT-1 0The scPAS provides both quantitative (scPAS_NRS) and qualitative (scPAS) prediction results:
library(ggplot2)
UMAP_scPAS <- DimPlot(sc_dataset, reduction = 'umap',
group.by = 'scPAS',raster=FALSE,
cols = c("0"='grey',"scPAS+"='indianred1',"scPAS-"='royalblue'),
pt.size = 0.001, order = c("scPAS-","scPAS+")) + ggtitle(label = NULL) + labs(col = "scPAS")
UMAP_RS <- FeaturePlot(object = sc_dataset,raster=FALSE,features = 'scPAS_NRS',max.cutoff = 2,min.cutoff = -2) + scale_colour_gradientn(colours = rev(brewer.pal(n = 10, name = "RdBu")))+
ggtitle(label = NULL) + labs(col = "Risk score")
UMAP_celltype | UMAP_RS | UMAP_scPAS
Run scPAS with imputation
imputation = T , implement single-cell data imputation using the default KNN method.
time.start <- Sys.time()
sc_dataset <- scPAS(bulk_dataset = bulk_dataset,sc_dataset = sc_dataset,assay = 'RNA',phenotype = phenotype,imputation = T, nfeature = 3000,alpha = 0.01,network_class = 'SC',family = 'cox')
#> [1] "Step 1:Quantile normalization of bulk data."
#> [1] "Step2: Imputation of missing values in single cell RNA-sequencing data with KNN...."
#> [1] "Step 2: Extracting single-cell expression profiles...."
#> [1] "Step 3: Constructing a gene-gene similarity by single cell data...."
#> [1] "Step 4: Optimizing the network-regularized sparse regression model...."
#> [1] "Perform cox regression on the given clinical outcomes:"
#> [1] "alpha = 0.01"
#> [1] "lambda = 3.57200645154233"
#> [1] "scPAS identified 232 rick+ features and 233 rick- features."
#> [1] "The percentage of selected feature is: 17.521%"
#>
#> [1] "|**************************************************|"
#> [1] "Step 5: calculating quantified risk scores...."
#> [1] "Step 6: qualitative identification by permutation test program with 2000 times random perturbations"
#> [1] "Finished."
time.end <- Sys.time()
time.end-time.start
#> Time difference of 3.014233 minsA new assay (imputation) has been created to store the imputed expression profiles.
sc_dataset
#> An object of class Seurat
#> 21550 features across 8853 samples within 2 assays
#> Active assay: imputation (10775 features, 0 variable features)
#> 1 layer present: data
#> 1 other assay present: RNA
#> 3 dimensional reductions calculated: pca, tsne, umapUsing imputed expression profiles for prediction, scPAS will capture more cells.
library(ggplot2)
UMAP_scPAS <- DimPlot(sc_dataset, reduction = 'umap',
group.by = 'scPAS',raster=FALSE,
cols = c("0"='grey',"scPAS+"='indianred1',"scPAS-"='royalblue'),
pt.size = 0.001, order = c("scPAS-","scPAS+")) + ggtitle(label = NULL) + labs(col = "scPAS")
UMAP_RS <- FeaturePlot(object = sc_dataset,raster=FALSE,features = 'scPAS_NRS',max.cutoff = 2,min.cutoff = -2) + scale_colour_gradientn(colours = rev(brewer.pal(n = 10, name = "RdBu")))+
ggtitle(label = NULL) + labs(col = "Risk score")
UMAP_celltype | UMAP_RS | UMAP_scPAS
Other imputation methods can also be used. Users can manually create a new assay (e.g., imputation) to store the imputed expression profiles. By changing the assay parameter in the scPAS function to the corresponding assay name (e.g., assay = ‘imputation’), scPAS will use the imputed expression profiles for prediction.
Apply the trained model to independent bulk data
The returned Seurat object not only contains the prediction results for each single cell but also provides the trained model.
names(sc_dataset@misc$scPAS_para)
#> [1] "alpha" "lambda" "family" "Coefs" "bulk" "phenotype"
#> [7] "Network"
head(sort(sc_dataset@misc$scPAS_para$Coefs))
#> SRP14 OAZ2 MGMT SOCS2 MTHFS C16orf54
#> -0.018268366 -0.013174652 -0.011565087 -0.011088680 -0.010492243 -0.009823842
tail(sort(sc_dataset@misc$scPAS_para$Coefs))
#> NET1 NR1H3 ADAMTS5 ASF1A RARS GTPBP4
#> 0.01357959 0.01625871 0.01657540 0.01908916 0.01976357 0.02352664Applying this model to other independent bulk data can enable the prediction of prognostic risk for bulk samples.
load(file = './data/LIRI_data.Rdata')
LIRI_exp[1:6,1:6]
#> SP106993 SP106998 SP107004 SP107007 SP107010 SP107012
#> 1/2-SBSRNA4 0.1441744 0.2660846 0.2031348 0.3873066 0.3810149 0.285765
#> A1BG 3.9890176 6.9508735 6.1227629 5.5053963 5.3647062 5.786084
#> A1BG-AS1 2.6502620 5.1971995 4.7794116 3.9968752 4.1382148 4.502823
#> A1CF 3.9264937 2.9610495 3.6155205 2.6217911 3.7064656 3.691380
#> A2LD1 1.8809964 2.3581755 2.2436378 1.7367255 2.0600142 1.949026
#> A2M 4.8589420 6.5323214 3.5344285 5.1471874 5.6751418 6.860621
head(LIRI_sur)
#> sample donor_id status time
#> 1 SP106993 DO48672 0 450
#> 2 SP106998 DO48674 0 1260
#> 3 SP107004 DO48677 0 1140
#> 4 SP107007 DO48679 0 1200
#> 5 SP107010 DO48681 0 900
#> 6 SP107012 DO48682 1 120Use the scPAS.prediction function to perform predictions on independent data.
LIRI_pre_res <- scPAS.prediction(model = sc_dataset,test.data = LIRI_exp)
head(LIRI_pre_res)
#> sample scPAS_RS scPAS_NRS scPAS_Pvalue scPAS_FDR scPAS
#> SP106993 SP106993 -0.10727771 -1.728757 4.192632e-02 5.590176e-02 0
#> SP106998 SP106998 -0.60846291 -7.421802 5.776874e-14 6.302044e-13 scPAS-
#> SP107004 SP107004 -0.50019297 -5.537983 1.529879e-08 6.799460e-08 scPAS-
#> SP107007 SP107007 0.09447085 1.091426 1.375427e-01 1.642301e-01 0
#> SP107010 SP107010 0.42192908 4.367016 6.297773e-06 1.821043e-05 scPAS+
#> SP107012 SP107012 0.55258674 6.780482 5.988799e-12 4.791039e-11 scPAS+Survival analysis shows that applying this model to the ICGC-LIRI dataset can effectively predict patient survival outcomes.
library("survival")
library("survminer")
LIRI_pre_res$status <- LIRI_sur$status
LIRI_pre_res$time <- LIRI_sur$time
fit <- survfit(Surv(time, status) ~ scPAS, data = LIRI_pre_res)
survplot_LIRI <- ggsurvplot(
fit, title='LIRI',
data = LIRI_pre_res,
size = 1,
conf.int = F,
pval.method = TRUE,
pval = TRUE,
risk.table = F,
risk.table.height = 0.25,
xlab ='Days',
ylab='Survival probability',
legend=c(0.7,0.2),
legend.title=''
)
survplot_LIRI$plot
Apply the trained model to spatial transcriptomics data
ST_data <- readRDS(file = './data/HCC_ST_L2.RDS')
SP1 <- SpatialDimPlot(ST_data, label = F, repel = T, label.size = 4,group.by = 'seurat_clusters',alpha = 0,
pt.size.factor = 1.6) + theme(legend.position = 'None')
SP2 <- SpatialDimPlot(ST_data, label = T,stroke = NA, repel = T, label.size = 4,group.by = 'seurat_clusters',alpha = 0.9,
pt.size.factor = 1.35) + theme(legend.position = "bottom")
SP1 | SP2
Directly applying the trained model to spatial transcriptomics data can reveal the spatial localization of prognostic-related cells
ST_data <- scPAS.prediction(model = sc_dataset,test.data = ST_data,assay = 'Spatial',imputation = T,independent = T)
#> [1] "Step2: Imputation of missing values in single cell RNA-sequencing data with KNN...."
head(ST_data)
#> orig.ident nCount_Spatial nFeature_Spatial seurat_cluster
#> AAACAACGAATAGTTC-1 SeuratProject 2587 1681 0
#> AAACAAGTATCTCCCA-1 SeuratProject 1549 1067 0
#> AAACAATCTACTAGCA-1 SeuratProject 3180 1911 3
#> AAACACCAATAACTGC-1 SeuratProject 981 625 6
#> AAACAGAGCGACTCCT-1 SeuratProject 2630 1698 0
#> AAACAGCTTTCAGAAG-1 SeuratProject 454 378 4
#> AAACAGGGTCTATATT-1 SeuratProject 442 371 4
#> AAACAGTGTTCCTGGG-1 SeuratProject 1104 736 6
#> AAACATGGTGAGAGGA-1 SeuratProject 1627 967 6
#> AAACATTTCCCGGATT-1 SeuratProject 2083 1427 0
#> seurat_clusters histology scPAS_RS scPAS_NRS scPAS_Pvalue
#> AAACAACGAATAGTTC-1 0 malignant 0.1121530 1.866797 3.096497e-02
#> AAACAAGTATCTCCCA-1 0 malignant 0.2814544 3.931696 4.217440e-05
#> AAACAATCTACTAGCA-1 3 malignant 0.1433746 2.812488 2.457994e-03
#> AAACACCAATAACTGC-1 6 unnamed -0.4778185 -4.289498 8.953882e-06
#> AAACAGAGCGACTCCT-1 0 malignant 0.2989895 5.044205 2.277051e-07
#> AAACAGCTTTCAGAAG-1 4 unnamed -0.2359508 -1.860167 3.143094e-02
#> AAACAGGGTCTATATT-1 4 unnamed -0.3812576 -2.855067 2.151388e-03
#> AAACAGTGTTCCTGGG-1 6 unnamed -0.3141847 -2.820325 2.398749e-03
#> AAACATGGTGAGAGGA-1 6 unnamed -0.4743717 -4.407568 5.226890e-06
#> AAACATTTCCCGGATT-1 0 malignant 0.1691941 2.842201 2.240164e-03
#> scPAS_FDR scPAS
#> AAACAACGAATAGTTC-1 4.414625e-02 scPAS+
#> AAACAAGTATCTCCCA-1 1.959965e-04 scPAS+
#> AAACAATCTACTAGCA-1 5.213861e-03 scPAS+
#> AAACACCAATAACTGC-1 6.241689e-05 scPAS-
#> AAACAGAGCGACTCCT-1 4.168220e-06 scPAS+
#> AAACAGCTTTCAGAAG-1 4.473080e-02 scPAS-
#> AAACAGGGTCTATATT-1 4.670918e-03 scPAS-
#> AAACAGTGTTCCTGGG-1 5.113284e-03 scPAS-
#> AAACATGGTGAGAGGA-1 4.157241e-05 scPAS-
#> AAACATTTCCCGGATT-1 4.835206e-03 scPAS+The results show that cells associated with poor survival are primarily located in regions enriched with malignant liver cancer cells.
p3 <- SpatialFeaturePlot(ST_data, features = "scPAS_NRS",stroke=NA,pt.size.factor = 1.27,min.cutoff = -2,max.cutoff = 2) +
ggtitle(label = NULL) + theme(legend.position = "top") + scale_fill_gradientn(name='Risk score',colours = Seurat:::SpatialColors(n = 100))
p4 <- SpatialDimPlot(ST_data,group.by = 'histology',cols = c("malignant"='indianred1',"non-malignant"= NA), stroke=NA,
alpha = c(1,0.1),pt.size.factor = 1.27) + ggtitle(label = NULL) + labs(col = "pathologists") + theme(legend.position = "top")
p5 <- SpatialDimPlot(ST_data,group.by = 'scPAS',cols = c("0"='grey',"scPAS+"='indianred1',"scPAS-"='royalblue'), stroke=NA,
alpha = c(1,0.1),pt.size.factor = 1.27) + ggtitle(label = NULL) + labs(col = "scPAS") + theme(legend.position = "top")
p4 | p3 | p5
Information about the current R session
sessionInfo()
#> R version 4.3.1 (2023-06-16 ucrt)
#> Platform: x86_64-w64-mingw32/x64 (64-bit)
#> Running under: Windows 11 x64 (build 22631)
#>
#> Matrix products: default
#>
#>
#> locale:
#> [1] LC_COLLATE=Chinese (Simplified)_China.utf8
#> [2] LC_CTYPE=Chinese (Simplified)_China.utf8
#> [3] LC_MONETARY=Chinese (Simplified)_China.utf8
#> [4] LC_NUMERIC=C
#> [5] LC_TIME=Chinese (Simplified)_China.utf8
#>
#> time zone: Asia/Shanghai
#> tzcode source: internal
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] survminer_0.4.9 ggpubr_0.6.0 survival_3.5-5
#> [4] ggplot2_3.4.3 preprocessCore_1.64.0 RColorBrewer_1.1-3
#> [7] Seurat_5.0.3 SeuratObject_5.0.2 sp_2.0-0
#> [10] Matrix_1.7-0 scPAS_0.2.0
#>
#> loaded via a namespace (and not attached):
#> [1] rstudioapi_0.15.0 jsonlite_1.8.7 magrittr_2.0.3
#> [4] spatstat.utils_3.0-3 farver_2.1.1 rmarkdown_2.24
#> [7] vctrs_0.6.3 ROCR_1.0-11 spatstat.explore_3.2-1
#> [10] rstatix_0.7.2 htmltools_0.5.6 broom_1.0.6
#> [13] sass_0.4.7 sctransform_0.4.1 parallelly_1.36.0
#> [16] KernSmooth_2.23-21 bslib_0.5.1 htmlwidgets_1.6.2
#> [19] ica_1.0-3 plyr_1.8.8 plotly_4.10.2
#> [22] zoo_1.8-12 cachem_1.0.8 igraph_1.5.1
#> [25] mime_0.12 lifecycle_1.0.3 pkgconfig_2.0.3
#> [28] R6_2.5.1 fastmap_1.1.1 fitdistrplus_1.1-11
#> [31] future_1.33.0 shiny_1.7.5 digest_0.6.33
#> [34] colorspace_2.1-0 patchwork_1.1.3 tensor_1.5
#> [37] RSpectra_0.16-1 irlba_2.3.5.1 labeling_0.4.2
#> [40] progressr_0.14.0 km.ci_0.5-6 fansi_1.0.4
#> [43] spatstat.sparse_3.0-2 httr_1.4.7 polyclip_1.10-4
#> [46] abind_1.4-5 compiler_4.3.1 withr_2.5.0
#> [49] backports_1.4.1 carData_3.0-5 fastDummies_1.7.3
#> [52] highr_0.10 ggsignif_0.6.4 MASS_7.3-60
#> [55] tools_4.3.1 lmtest_0.9-40 httpuv_1.6.11
#> [58] future.apply_1.11.0 goftest_1.2-3 glue_1.6.2
#> [61] nlme_3.1-162 promises_1.2.1 grid_4.3.1
#> [64] Rtsne_0.16 cluster_2.1.4 reshape2_1.4.4
#> [67] generics_0.1.3 gtable_0.3.3 spatstat.data_3.0-1
#> [70] KMsurv_0.1-5 tidyr_1.3.0 data.table_1.14.8
#> [73] car_3.1-2 utf8_1.2.3 spatstat.geom_3.2-4
#> [76] RcppAnnoy_0.0.21 ggrepel_0.9.3 RANN_2.6.1
#> [79] pillar_1.9.0 stringr_1.5.0 spam_2.10-0
#> [82] RcppHNSW_0.6.0 later_1.3.1 splines_4.3.1
#> [85] dplyr_1.1.2 lattice_0.21-8 deldir_1.0-9
#> [88] tidyselect_1.2.1 miniUI_0.1.1.1 pbapply_1.7-2
#> [91] knitr_1.43 gridExtra_2.3 scattermore_1.2
#> [94] xfun_0.40 matrixStats_1.0.0 stringi_1.7.12
#> [97] lazyeval_0.2.2 yaml_2.3.7 evaluate_0.21
#> [100] codetools_0.2-19 tibble_3.2.1 cli_3.6.1
#> [103] uwot_0.1.16 xtable_1.8-4 reticulate_1.31
#> [106] munsell_0.5.0 jquerylib_0.1.4 survMisc_0.5.6
#> [109] Rcpp_1.0.11 globals_0.16.2 spatstat.random_3.1-5
#> [112] png_0.1-8 parallel_4.3.1 ellipsis_0.3.2
#> [115] dotCall64_1.1-1 listenv_0.9.0 viridisLite_0.4.2
#> [118] scales_1.2.1 ggridges_0.5.4 leiden_0.4.3
#> [121] purrr_1.0.2 rlang_1.1.4 cowplot_1.1.1