Visualization Gallery

Introduction

This vignette demonstrates various visualization techniques for scPAS analysis results. We’ll cover:

  • UMAP/t-SNE plots for risk scores
  • Cell type enrichment analysis
  • Volcano plots
  • Heatmaps
  • Violin plots
  • Spatial transcriptomics visualization

Setup and Simulated Data

library(scPAS)
library(Seurat)
library(Matrix)
library(ggplot2)
library(RColorBrewer)
library(dplyr)

# Set theme
theme_set(theme_minimal())

Create Simulated scPAS Result

set.seed(42)

# Simulate single-cell data
n_genes <- 500
n_cells <- 800

# Create count matrix
counts <- matrix(
  rpois(n_genes * n_cells, lambda = 5),
  nrow = n_genes,
  ncol = n_cells
)
rownames(counts) <- paste0("Gene", 1:n_genes)
colnames(counts) <- paste0("Cell", 1:n_cells)

# Create Seurat object
sc_obj <- CreateSeuratObject(counts = counts, project = "VisDemo")

# Add cell type labels with different proportions
cell_types <- c(
  rep("T_cells", 200),
  rep("B_cells", 150),
  rep("Monocytes", 150),
  rep("NK_cells", 100),
  rep("Fibroblasts", 100),
  rep("Endothelial", 100)
)
sc_obj$celltype <- cell_types

# Run preprocessing
sc_obj <- NormalizeData(sc_obj, verbose = FALSE)
sc_obj <- FindVariableFeatures(sc_obj, verbose = FALSE)
sc_obj <- ScaleData(sc_obj, verbose = FALSE)
sc_obj <- RunPCA(sc_obj, verbose = FALSE)
sc_obj <- RunUMAP(sc_obj, dims = 1:20, verbose = FALSE)

# Simulate scPAS results
# Make T_cells and NK_cells enriched for scPAS+
# Make B_cells and Fibroblasts enriched for scPAS-
set.seed(42)
risk_scores <- rep(NA, n_cells)

# T_cells and NK_cells: higher risk
risk_scores[cell_types %in% c("T_cells", "NK_cells")] <- 
  rnorm(sum(cell_types %in% c("T_cells", "NK_cells")), mean = 1.5, sd = 1.2)

# B_cells and Fibroblasts: lower risk
risk_scores[cell_types %in% c("B_cells", "Fibroblasts")] <- 
  rnorm(sum(cell_types %in% c("B_cells", "Fibroblasts")), mean = -1.5, sd = 1.2)

# Others: neutral
risk_scores[cell_types %in% c("Monocytes", "Endothelial")] <- 
  rnorm(sum(cell_types %in% c("Monocytes", "Endothelial")), mean = 0, sd = 0.8)

# Calculate p-values (simulated)
p_values <- 2 * pnorm(-abs(risk_scores))
fdr <- p.adjust(p_values, method = "BH")

# Classification
classification <- ifelse(
  fdr >= 0.05, "0",
  ifelse(risk_scores > 0, "scPAS+", "scPAS-")
)

# Add to metadata
sc_obj$scPAS_NRS <- risk_scores
sc_obj$scPAS_RS <- risk_scores * 0.01  # Raw scores
sc_obj$scPAS_Pvalue <- p_values
sc_obj$scPAS_FDR <- fdr
sc_obj$scPAS <- factor(classification, levels = c("scPAS-", "0", "scPAS+"))

Basic UMAP Plots

Cell Type Overview

# Define colors
celltype_colors <- c(
  "T_cells" = "#E41A1C",
  "B_cells" = "#377EB8",
  "Monocytes" = "#4DAF4A",
  "NK_cells" = "#984EA3",
  "Fibroblasts" = "#FF7F00",
  "Endothelial" = "#FFFF33"
)

DimPlot(
  sc_obj,
  group.by = "celltype",
  cols = celltype_colors,
  label = TRUE,
  label.size = 4,
  pt.size = 1
) +
  ggtitle("Cell Type Distribution") +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    legend.position = "right"
  )

Risk Score Visualization

FeaturePlot(
  sc_obj,
  features = "scPAS_NRS",
  pt.size = 1,
  order = TRUE
) +
  scale_color_gradient2(
    low = "royalblue",
    mid = "white",
    high = "indianred",
    midpoint = 0,
    name = "NRS"
  ) +
  ggtitle("Normalized Risk Score") +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14)
  )

Cell Classification

# Classification colors
class_colors <- c("scPAS-" = "royalblue", "0" = "gray80", "scPAS+" = "indianred")

DimPlot(
  sc_obj,
  group.by = "scPAS",
  cols = class_colors,
  pt.size = 1,
  order = c("0", "scPAS-", "scPAS+")  # Draw significant cells on top
) +
  ggtitle("scPAS Classification") +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14)
  )

Combined Multi-Panel Plot

library(patchwork)

p1 <- DimPlot(sc_obj, group.by = "celltype", cols = celltype_colors, label = TRUE) +
  ggtitle("Cell Types") + NoLegend()

p2 <- FeaturePlot(sc_obj, features = "scPAS_NRS", pt.size = 0.8) +
  scale_color_gradient2(low = "royalblue", mid = "white", high = "indianred", midpoint = 0) +
  ggtitle("Risk Score")

p3 <- DimPlot(sc_obj, group.by = "scPAS", cols = class_colors, 
              order = c("0", "scPAS-", "scPAS+")) +
  ggtitle("Classification")

p1 | p2 | p3

Cell Type Enrichment Analysis

Proportion Bar Plot

# Calculate proportions
prop_data <- as.data.frame(table(sc_obj$celltype, sc_obj$scPAS))
colnames(prop_data) <- c("CellType", "scPAS", "Count")

# Calculate proportions within each cell type
prop_data <- prop_data %>%
  dplyr::group_by(CellType) %>%
  dplyr::mutate(
    Total = sum(Count),
    Proportion = Count / Total
  ) %>%
  dplyr::ungroup()

ggplot(prop_data, aes(x = CellType, y = Proportion, fill = scPAS)) +
  geom_bar(stat = "identity", position = "stack", width = 0.7) +
  scale_fill_manual(values = class_colors) +
  labs(
    x = "Cell Type",
    y = "Proportion",
    title = "scPAS Classification by Cell Type",
    fill = "Classification"
  ) +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(hjust = 0.5, face = "bold"),
    legend.position = "right"
  ) +
  coord_flip()

Enrichment Heatmap

# Calculate enrichment (observed / expected)
enrichment_matrix <- reshape2::dcast(prop_data, CellType ~ scPAS, value.var = "Count")
rownames(enrichment_matrix) <- enrichment_matrix$CellType
enrichment_matrix <- as.matrix(enrichment_matrix[, -1])

# Calculate expected counts
total_per_class <- colSums(enrichment_matrix)
total_per_type <- rowSums(enrichment_matrix)
grand_total <- sum(enrichment_matrix)

expected <- outer(total_per_type, total_per_class) / grand_total
enrichment <- log2((enrichment_matrix + 1) / (expected + 1))

# Plot heatmap
enrichment_df <- reshape2::melt(enrichment)
colnames(enrichment_df) <- c("CellType", "scPAS", "Enrichment")

ggplot(enrichment_df, aes(x = scPAS, y = CellType, fill = Enrichment)) +
  geom_tile(color = "white", size = 0.5) +
  geom_text(aes(label = round(Enrichment, 2)), color = "black", size = 4) +
  scale_fill_gradient2(
    low = "royalblue",
    mid = "white",
    high = "indianred",
    midpoint = 0,
    name = "log2(O/E)"
  ) +
  labs(
    x = "scPAS Classification",
    y = "Cell Type",
    title = "Cell Type Enrichment in scPAS Classes"
  ) +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    axis.text.x = element_text(angle = 0)
  )

Volcano-Style Plot

# Create volcano-style data
volcano_data <- data.frame(
  NRS = sc_obj$scPAS_NRS,
  neg_log_FDR = -log10(sc_obj$scPAS_FDR + 1e-10),
  Classification = sc_obj$scPAS
)

ggplot(volcano_data, aes(x = NRS, y = neg_log_FDR, color = Classification)) +
  geom_point(alpha = 0.6, size = 1.5) +
  scale_color_manual(values = class_colors) +
  geom_hline(yintercept = -log10(0.05), linetype = "dashed", color = "gray40") +
  geom_vline(xintercept = 0, linetype = "dashed", color = "gray40") +
  annotate("text", x = max(volcano_data$NRS) * 0.7, y = -log10(0.05) + 0.5, 
           label = "FDR = 0.05", color = "gray40", size = 3.5) +
  labs(
    x = "Normalized Risk Score",
    y = "-log10(FDR)",
    title = "Cell-Level Volcano Plot"
  ) +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    legend.position = "bottom"
  )

Violin Plots

Risk Score by Cell Type

VlnPlot(
  sc_obj,
  features = "scPAS_NRS",
  group.by = "celltype",
  cols = celltype_colors,
  pt.size = 0
) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray40") +
  stat_summary(fun = median, geom = "point", size = 3, color = "black") +
  labs(
    x = "Cell Type",
    y = "Normalized Risk Score",
    title = "Risk Score Distribution by Cell Type"
  ) +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(hjust = 0.5, face = "bold"),
    legend.position = "none"
  )

Split Violin by Classification

# Custom split violin
ggplot(sc_obj@meta.data, aes(x = celltype, y = scPAS_NRS, fill = scPAS)) +
  geom_violin(scale = "width", trim = TRUE, alpha = 0.7) +
  scale_fill_manual(values = class_colors) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray40") +
  labs(
    x = "Cell Type",
    y = "Normalized Risk Score",
    title = "Risk Score by Cell Type and Classification",
    fill = "scPAS"
  ) +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(hjust = 0.5, face = "bold")
  )

Box Plots with Statistical Tests

# Calculate summary statistics
library(dplyr)
summary_stats <- sc_obj@meta.data %>%
  group_by(celltype) %>%
  summarise(
    median_NRS = median(scPAS_NRS),
    mean_NRS = mean(scPAS_NRS),
    n = n(),
    .groups = "drop"
  ) %>%
  arrange(desc(median_NRS))

# Reorder cell types
sc_obj@meta.data$celltype_ordered <- factor(
  sc_obj$celltype,
  levels = summary_stats$celltype
)

ggplot(sc_obj@meta.data, aes(x = celltype_ordered, y = scPAS_NRS, fill = celltype_ordered)) +
  geom_boxplot(outlier.size = 0.5, alpha = 0.8) +
  scale_fill_manual(values = celltype_colors[summary_stats$celltype]) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray40") +
  labs(
    x = "Cell Type (ordered by median NRS)",
    y = "Normalized Risk Score",
    title = "Risk Score Distribution (Ordered)"
  ) +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(hjust = 0.5, face = "bold"),
    legend.position = "none"
  ) +
  annotate("text", x = 0.5, y = max(sc_obj$scPAS_NRS) * 0.9, 
           label = paste("n =", nrow(sc_obj@meta.data)), hjust = 0, size = 3.5)

Density Plots

ggplot(sc_obj@meta.data, aes(x = scPAS_NRS, fill = celltype)) +
  geom_density(alpha = 0.5) +
  scale_fill_manual(values = celltype_colors) +
  geom_vline(xintercept = 0, linetype = "dashed", color = "gray40") +
  labs(
    x = "Normalized Risk Score",
    y = "Density",
    title = "Risk Score Density by Cell Type",
    fill = "Cell Type"
  ) +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    legend.position = "right"
  )

Pie Chart Summary

# Summary counts
pie_data <- as.data.frame(table(sc_obj$scPAS))
colnames(pie_data) <- c("Classification", "Count")
pie_data$Percentage <- pie_data$Count / sum(pie_data$Count) * 100
pie_data$Label <- paste0(pie_data$Classification, "\n", round(pie_data$Percentage, 1), "%")

ggplot(pie_data, aes(x = "", y = Count, fill = Classification)) +
  geom_bar(stat = "identity", width = 1, color = "white") +
  coord_polar(theta = "y") +
  scale_fill_manual(values = class_colors) +
  geom_text(aes(label = Label), position = position_stack(vjust = 0.5), size = 4) +
  labs(title = "scPAS Classification Summary") +
  theme_void() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    legend.position = "none"
  )

Publication-Ready Figure

# Create a publication-ready multi-panel figure
p1 <- DimPlot(sc_obj, group.by = "celltype", cols = celltype_colors, label = TRUE, pt.size = 0.5) +
  ggtitle("A. Cell Types") + theme(legend.position = "none")

p2 <- FeaturePlot(sc_obj, features = "scPAS_NRS", pt.size = 0.5) +
  scale_color_gradient2(low = "royalblue", mid = "white", high = "indianred", midpoint = 0, name = "NRS") +
  ggtitle("B. Risk Score")

p3 <- DimPlot(sc_obj, group.by = "scPAS", cols = class_colors, pt.size = 0.5,
              order = c("0", "scPAS-", "scPAS+")) +
  ggtitle("C. Classification")

p4 <- VlnPlot(sc_obj, features = "scPAS_NRS", group.by = "celltype", 
              cols = celltype_colors, pt.size = 0) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  ggtitle("D. Risk by Cell Type") + NoLegend() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

p5 <- ggplot(prop_data, aes(x = CellType, y = Proportion, fill = scPAS)) +
  geom_bar(stat = "identity", position = "stack", width = 0.7) +
  scale_fill_manual(values = class_colors) +
  coord_flip() +
  labs(x = "", y = "Proportion", fill = "") +
  ggtitle("E. Proportion by Cell Type")

p6 <- ggplot(pie_data, aes(x = "", y = Count, fill = Classification)) +
  geom_bar(stat = "identity", width = 1, color = "white") +
  coord_polar(theta = "y") +
  scale_fill_manual(values = class_colors) +
  geom_text(aes(label = Label), position = position_stack(vjust = 0.5), size = 3) +
  theme_void() +
  ggtitle("F. Overall Summary") +
  theme(legend.position = "none")

# Combine
(p1 | p2 | p3) / (p4 | p5 | p6) +
  plot_annotation(
    title = "scPAS Analysis Results",
    theme = theme(plot.title = element_text(hjust = 0.5, face = "bold", size = 16))
  )

Saving Plots

# Save individual plots
ggsave("umap_celltype.pdf", plot = p1, width = 8, height = 6, dpi = 300)
ggsave("umap_risk_score.pdf", plot = p2, width = 8, height = 6, dpi = 300)
ggsave("umap_classification.pdf", plot = p3, width = 8, height = 6, dpi = 300)

# Save combined figure
combined_fig <- (p1 | p2 | p3) / (p4 | p5 | p6)
ggsave("scPAS_combined_figure.pdf", plot = combined_fig, width = 14, height = 10, dpi = 300)
ggsave("scPAS_combined_figure.png", plot = combined_fig, width = 14, height = 10, dpi = 300)

Session Information

sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04 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.32.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=en_US.UTF-8    
#>  [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] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#>  [1] future_1.75.0      survminer_0.5.2    ggpubr_1.0.0       survival_3.8-9    
#>  [5] dplyr_1.2.1        patchwork_1.3.2    RColorBrewer_1.1-3 Seurat_5.5.1      
#>  [9] SeuratObject_5.4.0 sp_2.2-3           scPAS_1.0.4        Matrix_1.7-5      
#> [13] ggplot2_4.0.3      rmarkdown_2.31    
#> 
#> loaded via a namespace (and not attached):
#>   [1] sys_3.4.3              jsonlite_2.0.0         magrittr_2.0.5        
#>   [4] spatstat.utils_3.2-4   farver_2.1.2           vctrs_0.7.3           
#>   [7] ROCR_1.0-12            spatstat.explore_3.8-1 rstatix_1.0.0         
#>  [10] htmltools_0.5.9        broom_1.0.13           Formula_1.2-5         
#>  [13] sass_0.4.10            sctransform_0.4.3      parallelly_1.48.0     
#>  [16] KernSmooth_2.23-26     bslib_0.11.0           htmlwidgets_1.6.4     
#>  [19] ica_1.0-3              plyr_1.8.9             plotly_4.12.0         
#>  [22] zoo_1.8-15             cachem_1.1.0           commonmark_2.0.0      
#>  [25] buildtools_1.0.0       igraph_2.3.3           mime_0.13             
#>  [28] lifecycle_1.0.5        pkgconfig_2.0.3        R6_2.6.1              
#>  [31] fastmap_1.2.0          fitdistrplus_1.2-6     shiny_1.14.0          
#>  [34] digest_0.6.39          tensor_1.5.1           RSpectra_0.16-2       
#>  [37] irlba_2.3.7            labeling_0.4.3         progressr_1.0.0       
#>  [40] spatstat.sparse_3.2-0  httr_1.4.8             polyclip_1.10-7       
#>  [43] abind_1.4-8            compiler_4.6.1         withr_3.0.3           
#>  [46] S7_0.2.2               backports_1.5.1        carData_3.0-6         
#>  [49] fastDummies_1.7.6      ggsignif_0.6.4         MASS_7.3-66           
#>  [52] tools_4.6.1            lmtest_0.9-40          otel_0.2.0            
#>  [55] httpuv_1.6.17          future.apply_1.20.2    goftest_1.2-3         
#>  [58] glue_1.8.1             nlme_3.1-170           promises_1.5.0        
#>  [61] gridtext_0.1.6         grid_4.6.1             Rtsne_0.17            
#>  [64] cluster_2.1.8.2        reshape2_1.4.5         generics_0.1.4        
#>  [67] gtable_0.3.6           spatstat.data_3.1-9    preprocessCore_1.75.0 
#>  [70] tidyr_1.3.2            data.table_1.18.4      xml2_1.6.0            
#>  [73] car_3.1-5              spatstat.geom_3.8-1    RcppAnnoy_0.0.23      
#>  [76] markdown_2.0           ggrepel_0.9.8          RANN_2.6.2            
#>  [79] pillar_1.11.1          stringr_1.6.0          spam_2.11-4           
#>  [82] RcppHNSW_0.7.0         later_1.4.8            splines_4.6.1         
#>  [85] ggtext_0.1.2           lattice_0.22-9         deldir_2.0-4          
#>  [88] tidyselect_1.2.1       maketools_1.3.2        miniUI_0.1.2          
#>  [91] pbapply_1.7-4          knitr_1.51             gridExtra_2.3.1       
#>  [94] litedown_0.10          scattermore_1.2        xfun_0.60             
#>  [97] matrixStats_1.5.0      stringi_1.8.7          lazyeval_0.2.3        
#> [100] yaml_2.3.12            evaluate_1.0.5         codetools_0.2-20      
#> [103] tibble_3.3.1           cli_3.6.6              uwot_0.2.4            
#> [106] xtable_1.8-8           reticulate_1.46.0      jquerylib_0.1.4       
#> [109] Rcpp_1.1.2             globals_0.19.1         spatstat.random_3.5-0 
#> [112] png_0.1-9              spatstat.univar_3.2-0  parallel_4.6.1        
#> [115] dotCall64_1.2          listenv_1.0.0          viridisLite_0.4.3     
#> [118] scales_1.4.0           ggridges_0.5.7         purrr_1.2.2           
#> [121] rlang_1.3.0            cowplot_1.2.0