Introduction
Single-cell assay for transposase-accessible chromatin sequencing (scATAC-seq) measures chromatin accessibility at single-cell resolution. Unlike scRNA-seq expression counts, scATAC-seq data have unique properties: * High Sparsity: Over 90% to 98% of values in peak-by-cell matrices are zeroes because each locus has at most 2 copies per cell. * Near-Binary Nature: Values primarily represent open (accessible) vs. closed chromatin states. * Peak Co-Accessibility: Genomic loci regulated together exhibit correlated accessibility across cells.
This tutorial shows how to evaluate synthetic scATAC-seq data using
scSimEval.
1. Loading the Packaged scATAC-seq Dataset
scSimEval includes an empirical reference and simulated
scATAC-seq peak accessibility dataset (example_scatac):
library(scSimEval)
# Load packaged scATAC-seq example data
data(example_scatac)
cat("Reference scATAC-seq dimensions:", dim(example_scatac$ref), "\n")
#> Reference scATAC-seq dimensions: 60 80
cat("Simulated scATAC-seq dimensions:", dim(example_scatac$sim), "\n")
#> Simulated scATAC-seq dimensions: 60 80
cat("Cell types:", levels(example_scatac$cell_types), "\n")
#> Cell types: TypeA TypeB
cat("Batches:", levels(example_scatac$batch_info), "\n")
#> Batches: Batch1 Batch2The matrices contain integer insertion counts across genomic peaks ().
2. Epigenomic Distribution & Summary Properties
We evaluate whether simulated peak insertion frequencies and cell
library sizes match real distributions using
evaluate_simulation_accuracy():
atac_res <- evaluate_simulation_accuracy(
ref_data = example_scatac$ref,
sim_data = example_scatac$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 peak accessibility and library size statistical distances
head(atac_res$metrics_summary_table, 8)
#> Category Property Metric Value
#> 1 Distributional Properties library_size MAD 2.0000000000
#> 2 Distributional Properties library_size KS 0.1750000000
#> 3 Distributional Properties library_size MAE 2.3375000000
#> 4 Distributional Properties library_size RMSE 2.4874685928
#> 5 Distributional Properties library_size OV 0.8548143838
#> 6 Distributional Properties library_size Bhattacharyya 0.0001055641
#> 7 Distributional Properties library_size Wasserstein 2.3375000000
#> 8 Distributional Properties library_size ECDF_DiffArea 0.08634259263. Epigenomic Cluster & Cell-Type Concordance
In scATAC-seq, identifying cell types requires clustering cells based
on chromatin accessibility profiles. We test whether the simulator
preserves genuine biological separation using
evaluate_clustering_metrics():
atac_clust <- evaluate_clustering_metrics(
data = example_scatac$sim,
cell_types = example_scatac$cell_types,
ref_data = example_scatac$ref
)
cat("Average Silhouette Width (Simulated):", round(atac_clust$silhouette, 4), "\n")
#> Average Silhouette Width (Simulated): 0.0558
cat("Davies-Bouldin Index:", round(atac_clust$davies_bouldin, 4), "\n")
#> Davies-Bouldin Index: 3.6911
cat("Adjusted Rand Index (ARI):", round(atac_clust$ARI, 4), "\n")
#> Adjusted Rand Index (ARI): 0.7627
cat("Normalized Mutual Information (NMI):", round(atac_clust$NMI, 4), "\n")
#> Normalized Mutual Information (NMI): 0.7214. Epigenomic Batch Integration & Confounder Handling
When simulating multi-sample experiments, batch variations can
confound biological signals. We evaluate batch mixing metrics using
evaluate_batch_metrics():
atac_batch <- evaluate_batch_metrics(
data = example_scatac$sim,
batch_info = example_scatac$batch_info,
cell_types = example_scatac$cell_types,
verbose = FALSE
)
cat("Batch Shannon Entropy:", round(atac_batch$shannon_entropy, 4), "\n")
#> Batch Shannon Entropy: 0.9671
cat("Principal Component Regression (PCR R2):", round(atac_batch$pcr_r2, 4), "\n")
#> Principal Component Regression (PCR R2): 0.0125
cat("Cross-Batch Transfer Accuracy:", round(atac_batch$cross_batch_accuracy, 4), "\n")
#> Cross-Batch Transfer Accuracy: 0.8255. Chromatin Peak Co-Accessibility & Regulatory Coupling
In real biological cells, peaks that are close to each other or regulated by the same transcription factor complexes show correlated accessibility:
# Evaluate peak co-accessibility fidelity
coacc_res <- calc_peak_coaccessibility_fidelity(
ref_atac = example_scatac$ref,
sim_atac = example_scatac$sim
)
cat("Co-Accessibility RV Coefficient:", round(coacc_res$rv_coefficient, 4), "\n")
#> Co-Accessibility RV Coefficient: 0.6014
cat("Peak Co-Accessibility Pearson Correlation:", round(coacc_res$coaccessibility_pearson, 4), "\n")
#> Peak Co-Accessibility Pearson Correlation: 0.1128
cat("Peak Co-Accessibility Spearman Correlation:", round(coacc_res$coaccessibility_spearman, 4), "\n")
#> Peak Co-Accessibility Spearman Correlation: 0.0902
cat("Frobenius Matrix Distance:", round(coacc_res$frobenius_distance, 4), "\n")
#> Frobenius Matrix Distance: 0.00536. Master Unimodal scATAC-seq Pipeline
Execute the comprehensive evaluation workflow in a single coordinated command:
master_atac <- evaluate_simulation_accuracy(
ref_data = example_scatac$ref,
sim_data = example_scatac$sim,
elapsed_time = 18.2,
memory_mb = 512.4,
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 tidy summary table
head(master_atac$metrics_summary_table)
#> Category Property Metric Value
#> 1 Distributional Properties library_size MAD 2.0000000000
#> 2 Distributional Properties library_size KS 0.1750000000
#> 3 Distributional Properties library_size MAE 2.3375000000
#> 4 Distributional Properties library_size RMSE 2.4874685928
#> 5 Distributional Properties library_size OV 0.8548143838
#> 6 Distributional Properties library_size Bhattacharyya 0.00010556417. Conclusion
scSimEval effectively benchmarks scATAC-seq simulations,
validating peak insertion distributions, extreme sparsity, cell type
clustering, batch mixing, and chromatin peak co-accessibility without
requiring synthetic ground-truth annotations.
sessionInfo()
#> R version 4.6.1 (2026-06-24 ucrt)
#> Platform: x86_64-w64-mingw32/x64
#> Running under: Windows 11 x64 (build 26200)
#>
#> Matrix products: default
#> LAPACK version 3.12.1
#>
#> locale:
#> [1] LC_COLLATE=English_India.utf8 LC_CTYPE=English_India.utf8
#> [3] LC_MONETARY=English_India.utf8 LC_NUMERIC=C
#> [5] LC_TIME=English_India.utf8
#>
#> time zone: Asia/Calcutta
#> tzcode source: internal
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] scSimEval_0.99.7
#>
#> loaded via a namespace (and not attached):
#> [1] RColorBrewer_1.1-3 jsonlite_2.0.0
#> [3] magrittr_2.0.5 ggbeeswarm_0.7.3
#> [5] SuppDists_1.1-9.9 modeltools_0.2-24
#> [7] farver_2.1.2 nloptr_2.2.1
#> [9] rmarkdown_2.32 provenance_4.4
#> [11] fs_2.1.0 ragg_1.5.2
#> [13] vctrs_0.7.3 memoise_2.0.1
#> [15] minqa_1.2.8 htmltools_0.5.9
#> [17] S4Arrays_1.13.0 curl_8.0.0
#> [19] BiocNeighbors_2.7.3 SparseArray_1.13.2
#> [21] overlapping_2.5 sass_0.4.10
#> [23] bslib_0.12.0 htmlwidgets_1.6.4
#> [25] desc_1.4.3 testthat_3.3.2
#> [27] httr2_1.3.0 cachem_1.1.0
#> [29] igraph_2.3.3 lifecycle_1.0.5
#> [31] pkgconfig_2.0.3 rsvd_1.0.5
#> [33] Matrix_1.7-6 R6_2.6.1
#> [35] fastmap_1.2.0 clue_0.3-68
#> [37] rbibutils_2.4.1 MatrixGenerics_1.25.0
#> [39] digest_0.6.39 S4Vectors_0.51.10
#> [41] scater_1.41.2 irlba_2.3.7
#> [43] IsoplotR_7.0 textshaping_1.0.5
#> [45] GenomicRanges_1.65.4 RSQLite_3.53.3
#> [47] beachmat_2.29.3 philentropy_0.10.0
#> [49] filelock_1.0.3 abind_1.4-8
#> [51] compiler_4.6.1 proxy_0.4-29
#> [53] bit64_4.8.6 withr_3.0.3
#> [55] CellMixS_1.29.0 S7_0.2.2
#> [57] BiocParallel_1.47.0 viridis_0.6.5
#> [59] DBI_1.3.0 MASS_7.3-66
#> [61] DelayedArray_0.39.8 bluster_1.23.1
#> [63] tools_4.6.1 vipor_0.4.7
#> [65] otel_0.2.0 beeswarm_0.4.0
#> [67] prabclus_2.3-5 nnet_7.3-21
#> [69] glue_1.8.1 nlme_3.1-171
#> [71] grid_4.6.1 waddR_1.27.0
#> [73] cluster_2.1.8.3 ade4_1.7-24
#> [75] generics_0.1.4 gtable_0.3.6
#> [77] class_7.3-24 tidyr_1.3.2
#> [79] ScaledMatrix_1.21.0 BiocSingular_1.29.1
#> [81] XVector_0.52.0 flexmix_2.3-21
#> [83] BiocGenerics_0.59.12 ggrepel_0.9.8
#> [85] RANN_2.6.3 pillar_1.11.1
#> [87] limma_3.99.0 clValid_0.7
#> [89] robustbase_0.99-7 splines_4.6.1
#> [91] dplyr_1.2.1 BiocFileCache_3.3.0
#> [93] lattice_0.22-9 bit_4.6.0
#> [95] tidyselect_1.2.1 SingleCellExperiment_1.35.2
#> [97] locfit_1.5-9.12 scuttle_1.23.2
#> [99] clusterSim_0.51-6 knitr_1.52
#> [101] gridExtra_2.3.1 reformulas_0.4.4
#> [103] IRanges_2.47.5 Seqinfo_1.3.2
#> [105] edgeR_4.99.6 SummarizedExperiment_1.43.0
#> [107] stats4_4.6.1 xfun_0.60
#> [109] Biobase_2.73.2 statmod_1.5.2
#> [111] diptest_0.77-2 randtests_1.0.2
#> [113] brio_1.1.5 arm_1.15-3
#> [115] matrixStats_1.5.0 DEoptimR_1.2-2
#> [117] yaml_2.3.12 boot_1.3-32
#> [119] kSamples_1.2-12 eva_0.2.7
#> [121] evaluate_1.0.5 codetools_0.2-20
#> [123] MLmetrics_1.1.3 kernlab_0.9-33
#> [125] tibble_3.3.1 cli_3.6.6
#> [127] RcppParallel_6.2.1 systemfonts_1.3.2
#> [129] Rdpack_2.6.6 jquerylib_0.1.4
#> [131] dichromat_2.0-1 Rcpp_1.1.2
#> [133] dbplyr_2.6.0 coda_0.19-4.1
#> [135] parallel_4.6.1 pkgdown_2.2.1
#> [137] ggplot2_4.0.3 blob_1.3.0
#> [139] mclust_6.1.3 lme4_2.0-6
#> [141] viridisLite_0.4.3 ggridges_0.5.7
#> [143] scales_1.4.0 e1071_1.7-17
#> [145] purrr_1.2.2 fpc_2.2-15
#> [147] rlang_1.3.0 cowplot_1.2.0