In single-cell RNA sequencing (scRNA-seq), evaluating the quality of simulated count data is important before testing analysis tools such as clustering, batch correction, or differential expression algorithms.
This tutorial shows how to evaluate synthetic scRNA-seq data
step-by-step using scSimEval without
needing any artificial ground truth.
scSimEval includes pre-packaged real and simulated
scRNA-seq count matrices (example_scrna):
library(scSimEval)
# Load packaged scRNA-seq example data
data(example_scrna)
# Examine dimensions and structure
cat("Reference scRNA-seq dimensions:", dim(example_scrna$ref), "\n")
#> Reference scRNA-seq dimensions: 60 80
cat("Simulated scRNA-seq dimensions:", dim(example_scrna$sim), "\n")
#> Simulated scRNA-seq dimensions: 60 80
cat("Cell types represented:", levels(example_scrna$cell_types), "\n")
#> Cell types represented: TypeA TypeB
cat("Batches represented:", levels(example_scrna$batch_info), "\n")
#> Batches represented: Batch1 Batch2The dataset consists of: * ref: Feature-by-cell matrix
of real experimental scRNA-seq read counts (60 genes \(\times\) 80 cells). * sim:
Feature-by-cell matrix of simulated counts generated by a simulation
method. * cell_types: Factor of cell type labels. *
batch_info: Factor of technical batch identifiers.
We evaluate how closely the simulated gene expression distributions
match the real reference data using
evaluate_simulation_accuracy():
unimodal_res <- evaluate_simulation_accuracy(
ref_data = example_scrna$ref,
sim_data = example_scrna$sim,
compute_bivariate = FALSE,
verbose = FALSE
)
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
# Inspect statistical distance metrics for library size
lib_summary <- subset(unimodal_res$metrics_summary_table, Property == "library_size")
head(lib_summary, 6)
#> Category Property Metric Value
#> 1 Distributional Properties library_size MAD 1.500000e+01
#> 2 Distributional Properties library_size KS 2.125000e-01
#> 3 Distributional Properties library_size MAE 1.506250e+01
#> 4 Distributional Properties library_size RMSE 1.571504e+01
#> 5 Distributional Properties library_size OV 8.370258e-01
#> 6 Distributional Properties library_size Bhattacharyya 1.188839e-04We evaluate whether the simulated data preserves distinct cell types
in unsupervised and supervised clustering using
evaluate_clustering_metrics():
clust_res <- evaluate_clustering_metrics(
data = example_scrna$sim,
cell_types = example_scrna$cell_types,
ref_data = example_scrna$ref
)
# Preview clustering concordance metrics
cat("Average Silhouette Width (Simulated):", round(clust_res$silhouette, 4), "\n")
#> Average Silhouette Width (Simulated): 0.0337
cat("Davies-Bouldin Index:", round(clust_res$davies_bouldin, 4), "\n")
#> Davies-Bouldin Index: 4.9116
cat("Adjusted Rand Index (ARI):", round(clust_res$ARI, 4), "\n")
#> Adjusted Rand Index (ARI): 0.0671
cat("Normalized Mutual Information (NMI):", round(clust_res$NMI, 4), "\n")
#> Normalized Mutual Information (NMI): 0.0884We test whether technical batch differences are properly mixed
without masking biological differences using
evaluate_batch_metrics():
batch_res <- evaluate_batch_metrics(
data = example_scrna$sim,
batch_info = example_scrna$batch_info,
cell_types = example_scrna$cell_types,
verbose = FALSE
)
# Preview batch mixing metrics
cat("Batch Shannon Entropy:", round(batch_res$shannon_entropy, 4), "\n")
#> Batch Shannon Entropy: 0.9796
cat("Principal Component Regression (PCR R2):", round(batch_res$pcr_r2, 4), "\n")
#> Principal Component Regression (PCR R2): 0.0121
cat("Cross-Batch Transfer Accuracy:", round(batch_res$cross_batch_accuracy, 4), "\n")
#> Cross-Batch Transfer Accuracy: 0.6125We assess marker gene effect sizes and top marker overlap using
evaluate_deg_fidelity():
deg_res <- evaluate_deg_fidelity(
ref_data = example_scrna$ref,
sim_data = example_scrna$sim,
ref_celltypes = example_scrna$cell_types,
sim_celltypes = example_scrna$cell_types
)
# Display DEG summary metrics
head(deg_res$deg_summary_table)
#> Framework Metric Value
#> 1 Simpipe (Duo et al., 2024) DEG_Ratio 0.0000000
#> 2 Simpipe (Duo et al., 2024) PValue_Uniformity_Chisq 8.0000000
#> 3 Simpipe (Duo et al., 2024) Distribution_Score 1.0000000
#> 4 Simpipe (Duo et al., 2024) Classifier_Accuracy 0.3750000
#> 5 Simpipe (Duo et al., 2024) Classifier_Macro_F1 0.3650794
#> 6 Simpipe (Duo et al., 2024) Classifier_Macro_Recall 0.3968254
#> Direction
#> 1 Target = 1.0 (Optimal balance)
#> 2 Lower is better (closer to Uniform)
#> 3 Higher is better (1.0 = Uniform null)
#> 4 Higher is better (Group predictability)
#> 5 Higher is better (Balanced F1)
#> 6 Higher is better (Sensitivity)We infer developmental pseudotime and lineage branching directly from
scRNA-seq expression counts using
evaluate_trajectory_metrics():
traj_res <- evaluate_trajectory_metrics(
ref_data = example_scrna$ref,
sim_data = example_scrna$sim,
cell_types_ref = example_scrna$cell_types,
cell_types_sim = example_scrna$cell_types
)
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
cat("Pseudotime Spearman Correlation:", round(traj_res$pseudotime_correlation, 4), "\n")
#> Pseudotime Spearman Correlation: 1
cat("Lineage Tree Branch Height RMSE:", round(as.numeric(traj_res$tree_height_rmse), 4), "\n")
#> Lineage Tree Branch Height RMSE: 1.8672evaluate_simulation_accuracyInstead of running individual steps, execute the entire unimodal evaluation suite in a single coordinated command:
eval_res <- evaluate_simulation_accuracy(
ref_data = example_scrna$ref,
sim_data = example_scrna$sim,
elapsed_time = 12.5,
memory_mb = 340.2,
compute_bivariate = FALSE,
verbose = FALSE
)
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
#> Metric: 'bhattacharyya' with unit: 'log'; comparing: 2 vectors
# Preview tidy evaluation summary
head(eval_res$metrics_summary_table)
#> Category Property Metric Value
#> 1 Distributional Properties library_size MAD 1.500000e+01
#> 2 Distributional Properties library_size KS 2.125000e-01
#> 3 Distributional Properties library_size MAE 1.506250e+01
#> 4 Distributional Properties library_size RMSE 1.571504e+01
#> 5 Distributional Properties library_size OV 8.370258e-01
#> 6 Distributional Properties library_size Bhattacharyya 1.188839e-04scSimEval includes built-in functions for
high-resolution diagnostic plotting:
# 14-Panel Comparative Distribution QC Plot
plot_distribution_qc(
ref_data = example_scrna$ref,
sim_data = example_scrna$sim,
title = "scRNA-seq Simulation Quality Control"
)This demonstration illustrates how scSimEval
comprehensively assesses scRNA-seq simulations across distributions,
zero patterns, clustering, batch integration, marker genes, and
trajectory dynamics without requiring artificial ground-truth
labels.
sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04.1 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] scSimEval_0.99.7 rmarkdown_2.32
#>
#> loaded via a namespace (and not attached):
#> [1] RColorBrewer_1.1-3 sys_3.4.3
#> [3] jsonlite_2.0.0 magrittr_2.0.5
#> [5] ggbeeswarm_0.7.3 SuppDists_1.1-9.9
#> [7] modeltools_0.2-25 farver_2.1.2
#> [9] nloptr_2.2.1 provenance_4.4
#> [11] vctrs_0.7.3 memoise_2.0.1
#> [13] minqa_1.2.8 htmltools_0.5.9
#> [15] S4Arrays_1.13.2 curl_8.0.0
#> [17] BiocNeighbors_2.7.3 SparseArray_1.13.4
#> [19] overlapping_2.5 sass_0.4.10
#> [21] bslib_0.12.0 testthat_3.3.2
#> [23] httr2_1.3.0 cachem_1.1.0
#> [25] buildtools_1.0.0 igraph_2.3.4
#> [27] lifecycle_1.0.5 pkgconfig_2.0.3
#> [29] rsvd_1.0.5 Matrix_1.7-6
#> [31] R6_2.6.1 fastmap_1.2.0
#> [33] rbibutils_2.4.1 MatrixGenerics_1.25.0
#> [35] clue_0.3-68 digest_0.6.39
#> [37] S4Vectors_0.51.10 scater_1.41.2
#> [39] irlba_2.4.1 IsoplotR_7.0
#> [41] GenomicRanges_1.65.4 RSQLite_3.53.3
#> [43] beachmat_2.29.3 philentropy_0.10.0
#> [45] filelock_1.0.3 labeling_0.4.3
#> [47] abind_1.4-8 compiler_4.6.1
#> [49] proxy_0.4-29 bit64_4.8.6
#> [51] withr_3.0.3 CellMixS_1.29.0
#> [53] S7_0.2.2 BiocParallel_1.47.0
#> [55] viridis_0.6.5 DBI_1.3.0
#> [57] MASS_7.3-66 DelayedArray_0.39.8
#> [59] bluster_1.23.1 tools_4.6.1
#> [61] vipor_0.4.7 otel_0.2.0
#> [63] beeswarm_0.4.0 prabclus_2.3-5
#> [65] nnet_7.3-21 glue_1.8.1
#> [67] nlme_3.1-171 grid_4.6.1
#> [69] waddR_1.27.0 cluster_2.1.8.3
#> [71] ade4_1.7-24 generics_0.1.4
#> [73] gtable_0.3.6 class_7.3-24
#> [75] tidyr_1.3.2 ScaledMatrix_1.21.0
#> [77] BiocSingular_1.29.2 XVector_0.53.0
#> [79] flexmix_2.3-21 BiocGenerics_0.59.12
#> [81] ggrepel_0.9.8 RANN_2.6.3
#> [83] pillar_1.11.1 limma_3.99.1
#> [85] clValid_0.7 robustbase_0.99-7
#> [87] splines_4.6.1 dplyr_1.2.1
#> [89] BiocFileCache_3.3.0 lattice_0.23-1
#> [91] bit_4.6.0 tidyselect_1.2.1
#> [93] SingleCellExperiment_1.35.2 locfit_1.5-9.12
#> [95] scuttle_1.23.2 maketools_1.3.2
#> [97] clusterSim_0.51-6 knitr_1.52
#> [99] gridExtra_2.3.1 reformulas_0.4.4
#> [101] IRanges_2.47.5 Seqinfo_1.3.2
#> [103] edgeR_4.99.6 SummarizedExperiment_1.43.0
#> [105] stats4_4.6.1 xfun_0.61
#> [107] Biobase_2.73.2 statmod_1.5.2
#> [109] diptest_0.77-2 randtests_1.0.2
#> [111] brio_1.1.5 arm_1.15-3
#> [113] matrixStats_1.5.0 DEoptimR_1.2-2
#> [115] yaml_2.3.12 boot_1.3-32
#> [117] eva_0.2.7 evaluate_1.0.5
#> [119] codetools_0.2-20 MLmetrics_1.1.3
#> [121] kSamples_1.2-12 kernlab_0.9-33
#> [123] tibble_3.3.1 cli_3.6.6
#> [125] RcppParallel_6.2.1 Rdpack_2.6.6
#> [127] jquerylib_0.1.4 Rcpp_1.1.2
#> [129] dbplyr_2.6.0 coda_0.19-4.1
#> [131] parallel_4.6.1 ggplot2_4.0.3
#> [133] blob_1.3.0 mclust_6.1.3
#> [135] lme4_2.0-6 viridisLite_0.4.3
#> [137] ggridges_0.5.7 scales_1.4.0
#> [139] e1071_1.7-17 purrr_1.2.2
#> [141] fpc_2.2-15 rlang_1.3.0
#> [143] cowplot_1.2.0