---
title: "Getting Started with scSimEval: Unified Benchmarking for Single-Cell Multiomics Simulations"
author: "Kabilan S, Dr. Dwijesh Chandra Mishra, Dr. Shashi Bhushan Lal, Dr. Sudhir Srivastava, Dr. Krishna Kumar Chaturvedi, Dr. Sharanbasappa"
date: "`r Sys.Date()`"
output:
  rmarkdown::html_vignette:
    toc: true
    toc_depth: 3
    number_sections: false
    self_contained: true
vignette: >
  %\VignetteIndexEntry{Getting Started with scSimEval: Unified Benchmarking for Single-Cell Multiomics Simulations}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r knitr-setup, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.width = 7,
  fig.height = 5,
  warning = FALSE,
  message = FALSE
)
```

## 1. Introduction

Computer simulations of single-cell technologies—specifically single-cell RNA sequencing (**scRNA-seq**), single-cell ATAC sequencing (**scATAC-seq**), and paired **scRNA-seq + scATAC-seq multiomics**—are widely used in bioinformatics research. They help scientists test computational pipelines, optimize analysis workflows, and evaluate statistical models. A critical question is always: *how realistic is the simulated data compared to genuine biological experiments?*

Many existing evaluation tools require artificial "ground truth" labels that are rarely known in real experiments, such as predefined true gene networks or synthetic marker lists.

**`scSimEval`** solves this problem by providing **62 evaluation measures organized into 8 easy-to-understand categories**. It works directly with the real experimental data and simulated data you already have, without needing any artificial ground truth:

1. **Real experimental reference data** (`ref_rna`, `ref_atac`)
2. **Simulated synthetic data** (`sim_rna`, `sim_atac`)
3. **Cell type annotations** (`cell_types`)
4. **Batch or donor identifiers** (`batch_info`)
5. **Computer resources** (Runtime in seconds and peak RAM in MiB)
6. **Developmental trajectories** (automatically inferred directly from scRNA-seq counts without requiring external tools)

> **Data Modalities:** `scSimEval` supports three multiomics experimental designs:
> **Paired** (co-assay from identical cells, e.g. 10x Multiome, SHARE-seq) — full evaluation
> including all Category 7 coupling metrics; **Unpaired** (independent scRNA-seq + scATAC-seq
> from separate cells) — unimodal evaluation per modality plus population-level cross-modal
> metrics (module fidelity, network overlap, label transfer); and **Mosaic** (partial
> co-measurement, e.g. DOGMA-seq) — unimodal evaluation on full matrices plus paired coupling
> metrics applied to the co-assayed cell subset. See the dedicated
> [Multiomics vignette](demo-multiomics.html) for full details.

---

## 2. Package Architecture and Workflow

The following workflow diagram shows the step-by-step benchmarking process in `scSimEval`:

```{r workflow-diagram, echo=FALSE, fig.align='center', out.width='100%', fig.cap='**Figure 1: End-to-end benchmarking workflow of scSimEval.** The framework takes real and simulated single-cell multiomics data, evaluates 62 metrics across 8 foundational categories, and produces standardized scores, comprehensive figures, and ranking leaderboards.'}
knitr::include_graphics("figures/workflow_diagram.png")
```

The evaluation process follows four simple steps:

### Step 1: Input Real and Simulated Data
Provide your real experimental reference and simulated count matrices (`ref_rna`, `sim_rna`, `ref_atac`, `sim_atac`), along with cell type labels (`cell_types`), batch labels (`batch_info`), and computer resource records (runtime in seconds and peak RAM in MiB).

### 2.1 Input Data Formats & Metadata Requirements (Cell Types & Batch Labels)

To evaluate your datasets, `scSimEval` accepts standard single-cell matrices and optional metadata:

1. **Count Matrices (`ref`, `sim`):**
   - **Dimensions:** Features (Genes or ATAC Peaks) in **rows** $\times$ Cells in **columns**.
   - **Data Types:** R sparse matrices (`dgCMatrix`), standard numeric matrices, `data.frame`, or Bioconductor count containers (`SingleCellExperiment`, `Seurat`).

2. **Cell-Type Annotations (`cell_types`, optional):**
   - Used for computing clustering fidelity metrics (Adjusted Rand Index, Normalized Mutual Information, Silhouette width), marker gene preservation, and cell-type coloring in embedding projections.
   - **Format in R:** A character vector or factor whose length and ordering strictly match the cells (columns) in the reference count matrix:
     ```r
     # Built-in reference template in scSimEval:
     data(example_cell_types, package = "scSimEval")
     head(example_cell_types)
     # Cell_01 Cell_02 Cell_03 Cell_04 Cell_05 
     #   TypeA   TypeA   TypeA   TypeA   TypeA 
     # Levels: TypeA TypeB
     ```
   - **Format in CSV/TSV:** A 2-column table with `Cell_ID` and `Cell_Type` (or a single-column file ordered identically to the count matrix):
     ```csv
     Cell_ID,Cell_Type
     Cell_01,TypeA
     Cell_02,TypeA
     Cell_03,TypeB
     ```

3. **Batch Annotations (`batch_info`, optional):**
   - Used for Category 4 batch mixing and technical confounder alignment metrics (kBET, CMS, LISI, variance decomposition).
   - **Format in R:** A character vector or factor matching cell columns:
     ```r
     data(example_batch_info, package = "scSimEval")
     head(example_batch_info)
     # Cell_01 Cell_02 Cell_03 Cell_04 Cell_05 
     #  Batch1  Batch2  Batch1  Batch2  Batch1 
     # Levels: Batch1 Batch2
     ```
   - **Format in CSV/TSV:** A 2-column table with `Cell_ID` and `Batch`:
     ```csv
     Cell_ID,Batch
     Cell_01,Batch1
     Cell_02,Batch2
     Cell_03,Batch1
     ```

4. **Exporting Annotations from Seurat or SingleCellExperiment:**
   ```r
   # From a Seurat object:
   cell_types_df <- data.frame(Cell_ID = colnames(seurat_obj), Cell_Type = as.character(Idents(seurat_obj)))
   write.csv(cell_types_df, "cell_types.csv", row.names = FALSE)
   
   batch_df <- data.frame(Cell_ID = colnames(seurat_obj), Batch = as.character(seurat_obj$batch))
   write.csv(batch_df, "batch_annotations.csv", row.names = FALSE)

   # From a SingleCellExperiment object:
   cell_types_df <- data.frame(Cell_ID = colnames(sce), Cell_Type = as.character(colData(sce)$cell_type))
   write.csv(cell_types_df, "cell_types.csv", row.names = FALSE)
   ```

### 2.2 Metric Measurement Capacity by Data & Metadata Mode

`scSimEval` dynamically scales its evaluation depth based on the experimental data layers and metadata files provided:

| Experiment Mode | Uploaded Inputs | Metrics Measured | Active Evaluation Categories & Scope |
| :--- | :--- | :---: | :--- |
| **Single-Cell** *(scRNA / scATAC)* | Count Matrix alone *(no metadata)* | **~23 Metrics** | **Cat I** (Distributional Properties, 14 metrics), **Cat II** (Correlations & Zeroes, 6 metrics), **Cat VI** (Trajectory, 2 metrics), **Cat VIII** (Scalability, 2 metrics). Evaluates basic statistical distributions, library sizes, and sparsity curves. |
| **Single-Cell** *(scRNA / scATAC)* | Matrix + Cell Types | **~35–38 Metrics** | *Adds +12–15 metrics:* Unlocks **Cat III** (Cellular Structure & Clustering / ARI, NMI, AMI, Silhouette) and **Cat V** (Biological Signal / SimBench DE, DP, DV, top DEG Jaccard overlap). |
| **Single-Cell** *(scRNA / scATAC)* | Matrix + Batch Info | **~30 Metrics** | *Adds +7 metrics:* Unlocks **Cat IV** (Batch Effects & Confounder Mixing / Batch Silhouette, Shannon Entropy, PCR $R^2$, CMS, ISI, Seurat mixing). |
| **Single-Cell** *(scRNA / scATAC)* | **Full Metadata** *(Matrix + Cell Types + Batch)* | **Up to ~55 Metrics** | Unlocks the complete unimodal suite across **Categories I, II, III, IV, V, VI, and VIII**. |
| **Paired Multiomics** *(scRNA + scATAC)* | RNA + ATAC Matrices alone | **~48 Metrics** | Unimodal distributions for both modalities (~40 metrics) PLUS **Cat VII** (Cross-Modal Coupling: FOSCTTM cell alignment, Match@1, cross-modal generation, $R_V$ coefficient, peak-to-gene linkage). |
| **Paired Multiomics** *(scRNA + scATAC)* | **Full Metadata** *(RNA + ATAC + Cell Types + Batch)* | **All 62 Metrics** | Complete multiomics benchmark suite across all 8 canonical categories (**Categories I–VIII**). |
| **Unpaired Multiomics** *(Independent cells)* | RNA + ATAC Matrices alone | **~38–42 Metrics** | Both unimodal distributions plus population-level coupling (co-expression module fidelity, Network Jaccard, cross-modality correlation). |
| **Unpaired Multiomics** *(Independent cells)* | **RNA + ATAC + Cell Types (+ Batch)** | **Up to ~58 Metrics** | *Essential for unpaired designs:* Unlocks cluster-matched cross-modal label transfer concordance, cell-type regulatory linkage, and technical donor integration. |

### Step 2: Calculate Evaluation Measures
`scSimEval` calculates 62 evaluation measures across 8 core categories without needing any artificial ground truth.

### Step 3: Run the Complete Benchmark
Execute the complete evaluation in a single command using `evaluate_simulation_accuracy()`, `evaluate_multiomics_accuracy()`, or `evaluate_multiple_datasets()`.

### Step 4: Review Scores, Rankings, and Figures
Review clean summary tables, standardized scores ($0.00$ to $1.00$), simulator ranking leaderboards, and high-resolution figures.

---

## 3. Score Normalization and Visual Mapping Workflow

In single-cell benchmarking, different metrics have completely different units and directions. For example, runtime is in seconds, peak memory is in megabytes, statistical distances are near zero, and clustering accuracy ranges between -1 and 1. Furthermore, for some metrics, smaller values are better (such as distance, error, and runtime), whereas for others, larger values are better (such as correlation and accuracy).

To make fair and intuitive comparisons, `scSimEval` uses an automated **two-step score normalization pipeline**, followed by clear **visual mapping**:

```{r score-norm-workflow, echo=FALSE, fig.align='center', out.width='100%', fig.cap='**Figure 2: Score normalization and visual mapping workflow in scSimEval.** Box 1 shows the two-step normalization process: direction inversion for lower-is-better metrics followed by min-max scaling to a standard 0 to 1 range. Box 2 illustrates the visual mapping into bubble matrices, where bubble size reflects simulation quality and top performers are highlighted with bold squares.'}
knitr::include_graphics("figures/score_normalization_workflow.png")
```

### 3.1 Two-Step Score Normalization

1. **Step 1: Direction Inversion (Aligning Direction)**
   For all metrics where smaller values mean better simulation results (such as Kolmogorov-Smirnov distance, Wasserstein distance, root mean square error, and runtime), scores are inverted:
   $$\text{Inverted Value} = \text{Maximum} - \text{Value}$$
   After inversion, larger numbers consistently represent better simulation fidelity across every metric.

2. **Step 2: Min-Max Scaling ($0.00$ to $1.00$)**
   All values are rescaled to a standard range from $0.00$ (worst performer) to $1.00$ (best performer):
   $$\text{Standardized Score} = \frac{\text{Value} - \text{Minimum}}{\text{Maximum} - \text{Minimum}}$$

### 3.2 Visual Mapping to Bubble Matrices

* **Bubble Size:** Directly represents the standardized score ($0.00$ to $1.00$). Larger bubbles indicate higher simulation quality and closer match to real biology.
* **Top-Performer Highlighting:** Outstanding scores ($\ge 0.96$) are shown as bold squares, while standard scores are shown as filled circles.
* **Category Color Themes:** Each metric category has its own distinct color (for example, blue for distributions, green for marker genes, and orange for computer speed).
* **Missing Feature Dots:** Features or modalities not supported by a simulator appear as neutral, small grey dots.
* **Overall Leaderboard Ranking:** Simulators are arranged vertically from top to bottom by their composite average score across all criteria.

---

## 4. Installation and Getting Started

To install and load **`scSimEval`** directly from GitHub:

```r
# Install devtools if needed
if (!requireNamespace("devtools", quietly = TRUE)) install.packages("devtools")

# Install scSimEval from GitHub
devtools::install_github("kabilanbio/scSimEval")

# Load the package
library(scSimEval)
```

```{r setup, include=FALSE}
library(scSimEval)
```

`scSimEval` includes built-in example datasets so you can test all functions immediately:
* **`example_scrna`**: Paired real and simulated scRNA-seq matrices ($60\text{ genes} \times 80\text{ cells}$) with cell type labels and batches.
* **`example_scatac`**: Paired real and simulated scATAC-seq chromatin accessibility matrices ($60\text{ peaks} \times 80\text{ cells}$).
* **`example_multiomics`**: A paired bundle containing both RNA and ATAC modalities, cell labels, and computer resource logs.

```{r load-data}
# Load example datasets
data(example_scrna)
data(example_scatac)
data(example_multiomics)

# Check dimensions
cat("Real scRNA-seq dimensions:", dim(example_scrna$ref), "\n")
cat("Simulated scRNA-seq dimensions:", dim(example_scrna$sim), "\n")
cat("Cell types present:", levels(example_scrna$cell_types), "\n")
```

> **Package Example Datasets vs. Shiny App Demo Data:**
> The R package datasets above (`example_scrna`, `example_scatac`, `example_multiomics`) are small, runnable raw count matrices (60 features x 80 cells) engineered for rapid demonstration and unit testing in standard R terminal sessions.
> In contrast, the **Shiny App Studio** (`launch_scSimEval_app()`) bundles a rich pre-computed multi-dataset benchmark evaluation suite (`demo_benchmark_data.rds`, 6 simulators evaluated across all 62 metrics on 2 datasets) allowing interactive exploration of leaderboards, bubble matrices, and diagnostic plots without waiting for long evaluation runs.

---

## 5. Evaluating Data Properties and Distributions

The function `evaluate_simulation_accuracy()` evaluates 1D and 2D statistical distributions and zero-count patterns between real and simulated matrices:

```{r unimodal-eval}
# Evaluate unimodal accuracy on scRNA-seq
unimodal_res <- evaluate_simulation_accuracy(
  ref_data = example_scrna$ref,
  sim_data = example_scrna$sim,
  compute_bivariate = FALSE,
  verbose = FALSE
)

# Preview library size distribution metrics
lib_summary <- subset(unimodal_res$metrics_summary_table, Property == "library_size")
knitr::kable(head(lib_summary, 8), digits = 4, caption = "Library Size Statistical Distance Measures")
```

### 5.1 Comparative Distribution QC Plot (`plot_distribution_qc`)

To visually inspect whether simulated counts match real biological distributions, `plot_distribution_qc()` generates a clean 14-panel layout comparing properties such as library size, mean expression, gene variance, zero proportions, and cell-cell correlations:

```{r plot-dist-qc-fig, echo=FALSE, fig.align='center', out.width='100%', fig.cap='**Figure 3: 14-panel comparative distribution QC layout.** Compares real biological data (warm brick red) and simulated data (steel blue) across single-cell properties, bivariate relationships, and biological signal retention.'}
knitr::include_graphics("figures/comparative_distribution_qc.png")
```

```r
# Generate the 14-panel comparative QC plot
plot_distribution_qc(
  ref_data   = example_scrna$ref,
  sim_data   = example_scrna$sim,
  title      = "Data Properties Quality Control",
  ref_name   = "Real Data",
  sim_name   = "Simulated Data",
  cell_types = example_scrna$cell_types
)
```

---

## 6. Evaluating Cell Clustering and Group Structure

To test whether simulated data preserve distinct cell populations and cluster boundaries without requiring external truth, `evaluate_clustering_metrics()` calculates cluster separation and concordance:

```{r clustering-eval}
clust_res <- evaluate_clustering_metrics(
  data = example_scrna$sim,
  cell_types = example_scrna$cell_types,
  ref_data = example_scrna$ref
)

cat("Average Silhouette Width:", round(clust_res$silhouette, 4), "\n")
cat("Davies-Bouldin Index:", round(clust_res$davies_bouldin, 4), "\n")
cat("Adjusted Rand Index (ARI):", round(clust_res$ARI, 4), "\n")
cat("Normalized Mutual Information (NMI):", round(clust_res$NMI, 4), "\n")
```

---

## 7. Evaluating Batch Effects and Technical Confounders

When evaluating multi-batch or multi-donor simulations, `evaluate_batch_metrics()` checks whether technical batch differences are properly represented without overshadowing biological variation:

```{r batch-eval}
batch_res <- evaluate_batch_metrics(
  data = example_scrna$sim,
  batch_info = example_scrna$batch_info,
  cell_types = example_scrna$cell_types,
  verbose = FALSE
)

cat("Batch Shannon Entropy:", round(batch_res$shannon_entropy, 4), "\n")
cat("PC Regression R2 (PCR):", round(batch_res$pcr_r2, 4), "\n")
cat("Cross-Batch Transfer Accuracy:", round(batch_res$cross_batch_accuracy, 4), "\n")
```

---

## 8. Evaluating Cell Lineages and Developmental Trajectories

Instead of requiring external pseudotime labels, `scSimEval` automatically infers differentiation paths directly from scRNA-seq expression counts and compares the real and simulated trajectories:

```{r trajectory-eval}
# Automatically infer trajectories and evaluate fidelity
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
)

cat("Pseudotime Spearman Correlation:", round(traj_res$pseudotime_correlation, 4), "\n")
cat("Lineage Tree Branch Height RMSE:", round(as.numeric(traj_res$tree_height_rmse), 4), "\n")
```

---

## 9. Evaluating Biological Signals and Marker Genes (DEGs)

Simulated datasets must preserve the marker genes and expression fold-changes found in genuine biological experiments. `scSimEval` integrates 15 differentially expressed gene (DEG) fidelity measures from three major benchmarking frameworks:
* **Simpipe (Duo et al., 2024):** True DEG ratio, $p$-value uniformity on non-DEGs, and cell-identity classification accuracy and $F_1$ score.
* **SimBench (Cao et al., 2021):** Symmetric error on DEG proportions, fold-change Pearson and Spearman correlations, and top-DEG Jaccard overlap.
* **Shaky Foundations (Crowell et al., 2023):** Group silhouette width, group separation discrepancy, and percent variance explained (PVE).

```r
# Evaluate marker gene and biological signal retention
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,
  fdr_cutoff            = 0.05,
  logfc_cutoff          = 0.5,
  top_n_de              = 50
)

# View the first 10 metrics
knitr::kable(head(deg_res$deg_summary_table, 10), digits = 4)
```

---

## 10. Evaluating Computational Resource Usage (Scalability)

Simulating thousands or millions of single cells requires fast and memory-efficient software. `plot_scalability_benchmark()` generates a comprehensive 6-panel dashboard comparing runtime, memory usage, Pareto tradeoff frontiers, multithreading speedup, and cell throughput:

```{r plot-scalability-fig, echo=FALSE, fig.align='center', out.width='100%', fig.cap='**Figure 4: Computational scalability benchmark suite.** 6-panel dashboard comparing execution runtime, peak RAM, runtime-memory tradeoff, CPU efficiency, cost footprint, and throughput across simulation methods.'}
knitr::include_graphics("figures/scalability_benchmark.png")
```

```r
# Generate the 6-panel computational scalability dashboard
plot_scalability_benchmark(
  scalability_data = sim_comparison_list,
  mode = "composite",
  panels = "6panel"
)
```

---

## 11. Master Evaluation Pipeline: `evaluate_multiomics_accuracy`

The flagship function `evaluate_multiomics_accuracy()` runs all evaluation categories in one command and returns a tidy summary table:

```{r master-pipeline}
master_eval <- evaluate_multiomics_accuracy(
  ref_multi = example_multiomics$ref_multi,
  sim_multi = example_multiomics$sim_multi,
  cell_types = example_multiomics$cell_types,
  batch_info = example_multiomics$batch_info,
  memory_mb = example_multiomics$resource_stats$memory_mb,
  elapsed_time = example_multiomics$resource_stats$elapsed_time,
  verbose = FALSE
)

summary_df <- master_eval$benchmark_summary_table
cat("Total evaluated metric instances:", nrow(summary_df), "\n")
table(summary_df$Category)
```

---

## 12. Visualization Suite

`scSimEval` provides built-in plotting functions rendering high-resolution figures (supporting 600 DPI for publication and reporting needs).

### 12.1 Multi-Dimensional Benchmarking Bubble Matrix (`plot_benchmark_bubble_matrix`)

The flagship visualization **`plot_benchmark_bubble_matrix()`** summarizes multi-simulator performance across the 8 evaluation categories. Bubble sizes indicate normalized simulation fidelity (larger bubbles = better performance):

```{r plot-bubble-matrix-fig, echo=FALSE, fig.align='center', out.width='100%', fig.cap='**Figure 5: Flagship multi-dimensional benchmarking bubble matrix.** Ranks simulation methods top-to-bottom across the 8 evaluation categories. Bubble size reflects standardized fidelity score ($0.00$ to $1.00$), and top performers are highlighted with bold squares.'}
knitr::include_graphics("figures/benchmark_bubble_matrix.png")
```

```r
# Generate benchmark bubble matrix
plot_benchmark_bubble_matrix(
  data = sim_comparison_list,
  bubble_size_range = c(2, 7.5),
  title = "Benchmarking Single-Cell Simulation Methods",
  subtitle = "Comparative performance across 8 evaluation categories"
)
```

### 12.2 Evaluation Summary Horizontal Bar Matrix (`plot_evaluation_summary`)

To present an executive overview of overall simulation accuracy, `plot_evaluation_summary()` ranks candidate simulators from best to worst and displays average scores across categories:

```{r eval-summary-bars, echo=FALSE, fig.align='center', out.width='100%', fig.cap='**Figure 6: Executive evaluation summary horizontal bar matrix.** Ranks simulators across the 8 categories and displays overall composite performance with exact score annotations.'}
knitr::include_graphics("figures/evaluation_summary_bars.png")
```

```r
# Generate horizontal summary bar chart
plot_evaluation_summary(
  data = sim_comparison_list,
  title = "Single-Cell Simulator Evaluation Summary",
  subtitle = "Standardized benchmark scores across categories and overall ranking"
)
```

### 12.2b Simulator Performance Leaderboard (`compute_method_leaderboard`)

To obtain exact numeric rankings and composite fidelity scores across all 62 criteria, `compute_method_leaderboard()` extracts the definitive simulator performance leaderboard:

```r
# Generate the method performance leaderboard
leaderboard <- compute_method_leaderboard(demo$benchmark_summary_table)
knitr::kable(leaderboard, caption = "Overall Simulator Performance Leaderboard (Ranked by Standardized Fidelity)")
```

Every metric is first converted into a direction-aware standardized score in $[0.00, 1.00]$ (where $1.00$ universally represents best observed performance; distance and error measures are inverted). The **Fidelity_Score** is the arithmetic mean across all standardized metric scores, and **Average_Fidelity** expresses this as an overall fidelity percentage ($0\%$ to $100\%$).

### 12.3 Metric Score Distributions Across Categories (`plot_metric_boxplots`)

To inspect variability across individual metrics and simulators, `plot_metric_boxplots()` displays boxplots with jittered data points across all 8 categories:

```{r plot-boxplots-fig, echo=FALSE, fig.align='center', out.width='100%', fig.cap='**Figure 7: Metric score distributions across the 8 categories.** Displays standardized scores ($[0, 1]$, higher is better) with individual data points and boxplots across candidate simulators.'}
knitr::include_graphics("figures/metric_boxplots.png")
```

```r
# Generate 8-category standardized boxplots
plot_metric_boxplots(
  benchmark_data = sim_list,
  score_type     = "normalized",
  ncol           = 4
)
```

### 12.4 Comprehensive Benchmark Heatmap (`plot_metric_heatmap`)

To examine exact numeric values for every metric without any loss of detail, `plot_metric_heatmap()` displays all simulators on the horizontal axis and all 62 measures on the vertical axis, printing **exact raw values inside each cell**:

```{r plot-heatmap-fig, echo=FALSE, fig.align='center', out.width='88%', fig.cap='**Figure 8: Complete benchmark metric heatmap.** Displays exact unnormalized raw scores in bold text inside every cell, grouped cleanly across the 8 evaluation categories.'}
knitr::include_graphics("figures/metric_heatmap.png")
```

```r
# Render full benchmark heatmap with exact numbers
plot_metric_heatmap(
  benchmark_data    = sim_comparison_list,
  scale_fill        = "relative",
  facet_by_category = TRUE
)
```

### 12.5 Dimension Reduction and Ordination Plots (PCA & MDS)

To discover which simulation methods behave similarly across all metrics, `scSimEval` includes Principal Component Analysis (`plot_metric_pca`) and Multidimensional Scaling (`plot_metric_mds`):

```{r plot-pca-cat1, echo=FALSE, fig.align='center', out.width='85%', fig.cap='**Figure 9: Principal Component Analysis (PCA) ordination of simulation methods.** Projects simulators based on their metric profiles, showing global affinities and key discriminating metric vectors.'}
knitr::include_graphics("figures/individual_category_plots/pca_cat1_distribution.png")
```

```r
# Run PCA biplot across benchmark metric profiles
plot_metric_pca(
  benchmark_data = sim_comparison_list,
  top_n_loadings = 8
)

# Run classical MDS ordination
plot_metric_mds(
  benchmark_data = sim_comparison_list,
  point_size = 4.5
)
```

### 12.6 Cell Embeddings (t-SNE, UMAP, & PCA) and Low-Dimensional Quality Metrics

Beyond numerical metric profiles, directly projecting real biological cells and simulated cells into low-dimensional embedding spaces (such as UMAP, t-SNE, or PCA) enables researchers to visually inspect whether synthetic datasets faithfully capture phenotypic clusters, biological manifolds, and technical depth variations.

`scSimEval` provides a high-level suite of functions for embedding computation, comparative plotting, and quantitative low-dimensional quality assessment:

```{r cell-embeddings-demo, eval=TRUE}
# 1. Compute low-dimensional embeddings across reference and simulations
emb_data <- compute_dataset_embeddings(
  reference   = example_scrna$ref,
  simulated   = list("Splatter" = example_scrna$sim),
  reduction   = "pca",
  n_pcs       = 10,
  cell_types  = example_scrna$cell_types,
  batch       = example_scrna$batch_info
)

# Preview embedding data structure
head(emb_data, 3)

# 2. Plot comparative cell embeddings (faceted across datasets)
plot_dataset_embeddings(
  embedding_data = emb_data,
  reduction      = "pca",
  layout         = "facet",
  color_by       = "cell_type",
  pt_size        = 1.0,
  alpha          = 0.8
)

# 3. Compute quantitative low-dimensional quality metrics
quality_table <- compute_embedding_quality_metrics(emb_data)
knitr::kable(quality_table, digits = 4, caption = "Quantitative Low-Dimensional Quality Metrics")
```

The quantitative embedding metrics evaluated include:
- **Mean Silhouette Score**: Quantifies cluster separability in the low-dimensional embedding space.
- **Adjusted Rand Index (ARI)**: Evaluates unsupervised cluster recovery fidelity against known biological cell types using $k$-means on the principal components.
- **Mean Library Size & Mean Detected Features**: Verifies whether cellular sequencing depth and gene detection rates are aligned.
- **Library Size Discrepancy (%)**: Quantifies the percentage deviation in total cellular counts relative to the reference.

---

## 13. Multi-Dataset and Multi-Method Benchmarking

When conducting large-scale benchmarking studies across multiple patient samples, tissues, or experimental cohorts, `evaluate_multiple_datasets()` processes all datasets simultaneously:

```{r multi-dataset-eval}
# Set up a collection of datasets
dataset_collection <- list(
  "Method 1-scRNA-seq" = list(ref = example_scrna$ref, sim = example_scrna$sim),
  "Method 1-scATAC-seq" = list(ref = example_scatac$ref, sim = example_scatac$sim),
  "Method 2-scRNA-seq" = list(ref = example_scrna$ref, sim = example_scrna$sim),
  "Method 2-scATAC-seq" = list(ref = example_scatac$ref, sim = example_scatac$sim)
)

# Run consolidated benchmarking in one step
consolidated_results <- evaluate_multiple_datasets(
  datasets = dataset_collection,
  pair_by_prefix = TRUE,
  compute_bivariate = FALSE,
  verbose = FALSE
)

# Preview consolidated results overview
knitr::kable(consolidated_results$dataset_overview, digits = 4, caption = "Consolidated Benchmark Overview")
```

---

## 14. Interactive Web Application (Shiny Studio)

For researchers who prefer an interactive graphical interface rather than writing R code, `scSimEval` includes an embedded Shiny web application (`scSimEvalApp`). You can explore benchmark results, filter metrics, view bubble plots, and export high-resolution (600 DPI) figures directly in your web browser:

```r
# Launch the interactive web app in your default browser
launch_scSimEval_app()
```

### 14.1 The 7 Modular Application Tabs

The Shiny app provides a unified, production-ready graphical workspace structured into seven primary navigation tabs:

1. **Tab 1: Home** — Executive overview of the 62 evaluation measures, biological KPI dimensions, and an interactive operational tour across the four ingestion modes (with all-time pageview tracking in the cloud app).
2. **Tab 2: Data Hub** — Comprehensive data ingestion centre featuring 4 operational modes:
   - *Mode 1: Pre-Computed Benchmark (Demo Data)*: Instant 1-click loading of a multi-dataset benchmark study across 6 simulators (744 metric observations) for zero-wait exploration.
   - *Mode 2: Single-Cell RNA/ATAC Evaluation*: Evaluates empirical reference against simulated count matrices with optional cell-type and batch metadata.
   - *Mode 3: Multi-Omics Co-Assays*: Evaluates paired or unpaired multiomics simulations assessing cross-modal coupling and feature alignment.
   - *Mode 4: Upload Pre-Computed Benchmark Table*: Instant loading of previously calculated benchmark CSV or RDS files.
   *Includes an interactive matrix inspector (Subtab 2), a live **Dataset Properties Summary** QC table (cells, features, sparsity %, library sizes), and an automated **HPC Batch Script Generator** (SLURM / PBS / Bash) for cluster execution.*
3. **Tab 3: Comparative Bubble Matrix & Leaderboard** — Flagship 62-metric comparative matrix with a multi-dataset selector bar (`All Datasets (Average)` vs. individual datasets `Dataset 1`, `Dataset 2`, ...), Minimalist vs. Publication layout toggles, confirmation plot buttons, 600 DPI figure downloads, and **Download Standalone R Script** for reproducible analysis workflows outside Shiny. Subtab 2 provides visual simulator ranking podiums, composite fidelity rankings, and custom category weight sliders.
4. **Tab 4: Visualizations** — 8 dedicated diagnostic panels equipped with global multi-dataset selectors, active dataset badges, confirmation controls, and workflow R script downloads:
   - *1. Evaluation Summary*: Horizontal bar chart ranking simulators across categories and composite fidelity.
   - *2. Distribution QC*: Comparative ECDF curves, KS distances, and gene/cell mean-variance trends.
   - *3. Scalability Benchmark*: Runtime (seconds) vs. Peak RAM (**MiB**), resource cost index, and cell throughput rate.
   - *4. Metric Plots*: Category and metric score boxplots/barplots with Raw vs. Normalized view toggles.
   - *5. Metric Heatmap*: Unsupervised hierarchical clustering of simulators and metrics with exact numbers in bold.
   - *6. PCA Ordination*: Principal component biplots projecting simulators in metric space.
   - *7. MDS Metric Space*: Multidimensional scaling projection of simulator algorithmic proximity.
   - *8. Cell Embeddings*: Dimension reduction projections (UMAP, t-SNE, PCA) comparing reference cells to simulated cells, with quantitative cluster fidelity tables.
5. **Tab 5: Benchmark Results & Export** — Master Benchmark Scores Table featuring per-column dropdown filtering, explicit dual-value reporting (`Score_Raw` in native statistical units, runtime in seconds, peak RAM in **MiB**, and normalized `Score` in [0, 1]), alongside exports for Excel (`.xlsx`), CSV, TSV, RDS, and full ZIP archives.
6. **Tab 6: Help & Getting Started** — Comprehensive user manual, formula glossary, metric measurement capacity reference table, hardware sizing guidelines, and input file templates.
7. **Tab 7: Contact** — ICAR-IASRI institutional details, author profiles, and citation information.

---

## 15. Supported Multiomics Data Types and Metric Compatibility

Not all multiomics evaluation metrics are valid for every experimental design. Understanding
which design your simulated data belongs to is critical for selecting the correct functions.

### 15.1 Three Experimental Designs

| Data Type | Definition | Example Technologies |
|---|---|---|
| **Paired** | RNA and ATAC measured from **exactly the same cells** | 10x Multiome, SHARE-seq, SNARE-seq |
| **Unpaired** | RNA and ATAC from **separate cells** of the same tissue | Independent 10x scRNA + scATAC experiments |
| **Mosaic** | **Some cells** have both modalities; others have only one | DOGMA-seq, mosaic pooling designs |

### 15.2 Metric Compatibility by Data Type

| Metric Group | Paired | Unpaired | Mosaic |
|---|:---:|:---:|:---:|
| Unimodal evaluation (Cat 1–6, 8) per modality | ✅ | ✅ | ✅ |
| Co-expression module fidelity | ✅ | ✅ | ✅ |
| Peak co-accessibility fidelity | ✅ | ✅ | ✅ |
| Cross-modal label transfer | ✅ | ✅ | ✅ |
| Network Jaccard similarity | ✅ | ✅ | ✅ |
| **FOSCTTM / Match@1** | ✅ | ❌ | ✅ paired subset |
| **Cross-modal generation fidelity** | ✅ | ❌ | ✅ paired subset |
| **Peak-to-gene regulatory coupling** | ✅ | ❌ | ✅ paired subset |
| **`evaluate_multiomics_accuracy()`** | ✅ direct | ❌ | ✅ paired subset |

### 15.3 Recommended Functions by Data Type

**Paired data** — use the master pipeline directly:
```r
evaluate_multiomics_accuracy(ref_multi, sim_multi, cell_types, ...)
```

**Unpaired data** — evaluate modalities separately, then consolidate:
```r
# Per-modality unimodal evaluation
rna_eval  <- evaluate_simulation_accuracy(ref_rna,  sim_rna)
atac_eval <- evaluate_simulation_accuracy(ref_atac, sim_atac)

# Population-level cross-modal metrics (no cell pairing required)
calc_coexpression_module_fidelity(ref_rna, sim_rna)
calc_peak_coaccessibility_fidelity(ref_atac, sim_atac)
evaluate_cross_modal_prediction(sim_atac, cell_types)

# Consolidate side-by-side
evaluate_multiple_datasets(list(RNA = list(ref=ref_rna, sim=sim_rna),
                                ATAC = list(ref=ref_atac, sim=sim_atac)))
```

**Mosaic data** — subset co-assayed cells first, then combine:
```r
# Step 1: Extract co-assayed subset
sim_rna_paired  <- sim_rna_full[, paired_cell_idx]
sim_atac_paired <- sim_atac_full[, paired_cell_idx]

# Step 2: Paired coupling metrics on subset
evaluate_multiomics_accuracy(
  ref_multi = list(rna = ref_rna_paired, atac = ref_atac_paired),
  sim_multi = list(rna = sim_rna_paired, atac = sim_atac_paired)
)

# Step 3: Full unimodal evaluation on complete matrices
evaluate_simulation_accuracy(ref_rna_full, sim_rna_full)
evaluate_simulation_accuracy(ref_atac_full, sim_atac_full)
```

> See the dedicated **[Multiomics Vignette](demo-multiomics.html)** for complete
> step-by-step workflows and code examples for all three data types.

## 16. Large-Scale Multi-Dataset Benchmarking on HPC Clusters

In comprehensive simulation benchmark studies, researchers often evaluate simulation tools across multiple diverse biological tissues, cell lines, or conditions (e.g., evaluating 6 simulators across 16 distinct single-cell datasets, producing $16 \times 6 = 96$ simulation comparisons and 5,952 metric observations across all 62 measures).

To enable seamless, beginner-friendly execution of such studies without crashing computer memory, `scSimEval` provides `benchmark_simulation_suite()`:

```r
# Define reference datasets (supports matrices or RDS file paths on disk)
refs <- list(
  "PBMC_10k"       = pbmc_ref,
  "Pancreas_Baron" = pancreas_ref,
  "Brain_Cortex"   = cortex_ref
  # ... up to 16 datasets
)

# Define simulator outputs per dataset
sims <- list(
  "Splatter"  = list("PBMC_10k" = splat_pbmc, "Pancreas_Baron" = splat_panc, "Brain_Cortex" = splat_cortex),
  "scDesign3" = list("PBMC_10k" = scd3_pbmc,  "Pancreas_Baron" = scd3_panc,  "Brain_Cortex" = scd3_cortex),
  "SCRIP"     = list("PBMC_10k" = scrip_pbmc, "Pancreas_Baron" = scrip_panc, "Brain_Cortex" = scrip_cortex),
  "SymSim"    = list("PBMC_10k" = sym_pbmc,   "Pancreas_Baron" = sym_panc,   "Brain_Cortex" = sym_cortex),
  "dyngen"    = list("PBMC_10k" = dyn_pbmc,   "Pancreas_Baron" = dyn_panc,   "Brain_Cortex" = dyn_cortex),
  "simATAC"   = list("PBMC_10k" = sima_pbmc,  "Pancreas_Baron" = sima_panc,  "Brain_Cortex" = sima_cortex)
)

# Run high-performance evaluation
mset_bench <- benchmark_simulation_suite(
  datasets    = refs,
  simulators  = sims,
  threads     = 8,        # multi-core HPC acceleration
  auto_sparse = TRUE,     # coerces dense matrices to dgCMatrix to minimize RAM
  run_gc      = TRUE,     # garbage collection after each run prevents memory bloat
  save_path   = "my_multidataset_benchmark.rds" # direct Shiny-ready bundle
)

# Render comparative bubble matrix:
# Automatically averages metrics across all 16 datasets under each simulator.
# Simulators are shown on the vertical axis; no cluttered dataset axes.
plot_benchmark_bubble_matrix(mset_bench)
```

In the interactive **scSimEval Studio** (`launch_scSimEval_app()`), users can navigate to **Data Hub**, choose **Mode 1 (Pre-Computed Benchmark)** or **Mode 4 (Upload Saved Results)** to instantly load multi-dataset studies or upload HPC-generated RDS files.

---

## 17. Extension to Other Single-Cell Data Types

While `scSimEval` is primarily optimized and demonstrated for **scRNA-seq**, **scATAC-seq**, and **paired scRNA-seq + scATAC-seq multiomics**, all 62 evaluation measures operate on standard numeric count or continuous matrices. 

Consequently, the package can be readily extended to benchmark simulations across other single-cell data modalities, including:
* **CITE-seq / REAP-seq:** Passing RNA expression counts as Modality 1 and antibody-derived tag (ADT) surface protein counts as Modality 2.
* **Single-Cell DNA Methylation:** Passing cell-level methylation levels or $\beta$-values as numeric input matrices.
* **Spatial Omics:** Passing spatial gene expression counts alongside complementary spatial coordinates or morphological features.

---

## 18. Summary

`scSimEval` provides a unified, objective, and ground-truth-free evaluation toolkit for benchmarking single-cell and multiomics simulation methods. With 62 metrics organized into 8 clear categories, automated score normalization, and intuitive visualizations, it enables researchers to select and validate simulation methods with confidence.

---

## 19. Authors and Maintainers

`scSimEval` is developed and maintained by:

* **Kabilan S** (Ph.D Student and Maintainer, ICAR-IASRI) — `kabilan151414@gmail.com`
* **Dr. Dwijesh Chandra Mishra** (Author & Thesis Guide, ICAR-IASRI) — `dwij.mishra@gmail.com`
* **Dr. Shashi Bhushan Lal** (Author, ICAR-IASRI) — `sblall16@gmail.com`
* **Dr. Sudhir Srivastava** (Author, ICAR-IASRI) — `sudhir0401bm@gmail.com`
* **Dr. Krishna Kumar Chaturvedi** (Author, ICAR-IASRI) — `kkcchaturvedi@gmail.com`
* **Dr. Sharanbasappa** (Contributor, ICAR-IASRI) — `smadival509@gmail.com`

* **GitHub Repository:** [https://github.com/kabilanbio/scSimEval](https://github.com/kabilanbio/scSimEval)
* **Bug Reports & Issues:** [https://github.com/kabilanbio/scSimEval/issues](https://github.com/kabilanbio/scSimEval/issues)

---

## 18. Session Information

```{r session-info}
sessionInfo()
```