Introduction

scCS (Single-Cell Commitment Scores) is a Python package for quantifying the degree of transcriptional commitment of single cells toward specific differentiation fates, using RNA velocity as the underlying signal.

It generalizes the 2-state (homeostatic/activated) commitment score framework introduced in Kriukov et al. (2025) to arbitrary k-furcations — branching points where a progenitor population splits into k ≥ 2 terminal fates.

Three-Scorer Architecture

v0.7 introduces a three-scorer architecture, each designed for a specific experimental design:

How many experimental conditions?
       │
       ├── 1 condition  →  SingleScorer
       │                    Single-condition analysis: embedding,
       │                    scoring, driver genes, enrichment
       │
       ├── 2 conditions →  PairScorer
       │                    Pairwise comparison: permutation test,
       │                    delta-CS with bootstrap CI, mixed model
       │
       └── 3+ conditions → MultiScorer
                            Multi-comparison: omnibus tests (Kruskal-Wallis,
                            ANOVA), post-hoc (Dunn, Tukey, Conover),
                            pairwise deltas, mixed model contrasts

SingleScorer — single-condition analysis

Build a radial star embedding, compute commitment scores, per-cell fate affinities, entropy metrics, driver genes, and pathway enrichment.

PairScorer — pairwise comparison (exactly 2 conditions)

Build a shared embedding on pooled data, score each condition separately, then compare using permutation tests, delta-CS with bootstrap CI, mixed-effects models, and trajectory shift analysis.

MultiScorer — multi-condition comparison (3+ conditions)

Same shared-embedding approach as PairScorer, plus tiered statistical testing: omnibus tests (Kruskal-Wallis, ANOVA) followed by post-hoc pairwise comparisons (Dunn, Tukey HSD, Conover-Iman) with FDR/Bonferroni correction, all-pairs delta-CS, and mixed-model contrasts.

Motivation

Classical RNA velocity tools (scVelo, CellRank) describe where cells are going. scCS answers a complementary question: how strongly committed is each cell to a given fate, relative to the alternatives?

This is particularly useful when:

• Comparing commitment levels across experimental conditions (e.g. treated vs control)

• Identifying driver genes that correlate with commitment

• Quantifying reversibility of cell state transitions

• Studying multi-fate branching points (k ≥ 3)

• Testing whether commitment shifts across 3+ treatment groups

Quick Start — SingleScorer

import scCS

scorer = scCS.SingleScorer(
    adata,
    root='17',
    branches=['FateA', 'FateB', 'FateC'],
    obs_key='leiden',
)
scorer.build_embedding(ordering_metric='pseudotime')
scorer.refit_pseudotime()
scorer.fit()
result = scorer.score(n_bootstrap=500)
print(result.summary())
scorer.plot_star(result)

Quick Start — PairScorer

pscorer = scCS.PairScorer(
    adata,
    root='17',
    branches=['homeostatic', 'activated'],
    condition_obs_key='treatment',
    obs_key='leiden',
)
pscorer.build_embedding(ordering_metric='pseudotime')
pscorer.fit()
results = pscorer.score_all_conditions()
delta = pscorer.compute_delta_CS('control', 'treated')
stats = pscorer.compare_conditions(results)

Quick Start — MultiScorer

mscorer = scCS.MultiScorer(
    adata,
    root='17',
    branches=['homeostatic', 'activated'],
    condition_obs_key='treatment',
    obs_key='leiden',
)
mscorer.build_embedding(ordering_metric='pseudotime')
mscorer.fit()
results = mscorer.score_all_conditions()
omnibus = mscorer.compare_omnibus(results)
posthoc = mscorer.compare_posthoc(results, omnibus_results=omnibus)
deltas = mscorer.compute_pairwise_deltas()

Mathematical Framework

Given per-cell RNA velocity vectors projected into a radial star embedding:

1. Magnitude — magnitude_i = sqrt(vx_i² + vy_i²)

2. Angle — theta_i = atan2(vy_i, vx_i) mapped to [0°, 360°)

3. Binning — angles binned into N equal sectors of width 360°/N

4. Sector magnitude — M_sector(j) = sum of magnitudes in sector j

5. Unnormalized CS — unCS(i, j) = M_sector(i) / M_sector(j)

6. Normalized CS — nCS(i, j) = unCS(i, j) × n_cells(j) / n_cells(i)

For k fates, a full k×k pairwise matrix of unCS and nCS is computed. Per-cell scores are derived from the dot product of each cell’s velocity vector with the unit direction toward each fate centroid.

See Mathematical Framework for full derivations with LaTeX.

Entropy Metrics

scCS provides three complementary entropy metrics, each answering a different question about commitment uncertainty:

Population entropy (single scalar)

How evenly is total velocity mass distributed across fate sectors? H_pop = 0 → all mass in one sector; H_pop = 1 → uniform. Note: can be misleadingly high for split-committed populations where cells are individually decisive but split between fates.

Per-fate cell entropy (shape (k,))

For each fate j: how individually decisive are cells about fate j? Computed as the mean binary entropy of each cell’s affinity score s_ij treated as [s_ij, 1−s_ij], averaged over all cells. Low → cells are sharply committed (or sharply not committed) to fate j. High → cells are ambiguous about fate j (scores cluster near 0.5).

NN-smoothed per-cell entropy (shape (n_cells,))

For each cell: average cell_scores over its k nearest neighbors in the scCS embedding (X_sccs), then compute k-way Shannon entropy on the smoothed scores. Removes single-cell velocity noise while preserving local commitment structure. Use plot_nn_entropy_elbow() to choose the optimal k.

Workflow

AnnData (with scVelo velocity)
       │
       ▼
SingleScorer(adata, root=..., branches=[...], obs_key=...)
       │
       ├── build_embedding(ordering_metric='pseudotime')
       │                          → radial star layout in obsm['X_sccs']
       ├── refit_pseudotime()     → fix arm coverage with subset pseudotime
       │
       ├── fit()                  → builds FateMap, projects velocity
       │
       ├── score(k_nn=15)         → CommitmentScoreResult
       │     ├── M_sector         (per-fate velocity magnitude)
       │     ├── pairwise_unCS    (k×k unnormalized scores)
       │     ├── pairwise_nCS     (k×k cell-count corrected scores)
       │     ├── cell_scores      (per-cell fate affinities)
       │     ├── population_entropy   (scalar)
       │     ├── per_fate_entropy     (shape k)
       │     └── nn_cell_entropy      (shape n_cells)
       │
       ├── plot_star()                → radial embedding visualization
       ├── plot_commitment_bar()      → unCS/nCS bar chart
       ├── plot_rose()                → polar velocity magnitude plot
       ├── plot_pairwise_cs()         → k×k heatmap
       ├── plot_nn_entropy_elbow()    → choose optimal k_nn
       ├── plot_expression_trends()   → gene expression vs CS axis
       ├── get_velocity_drivers()     → ranked driver genes per fate
       ├── get_deg_drivers()          → DEG analysis per fate arm
       ├── get_velocity_fate_drivers() → velocity-fate Spearman correlation
       └── get_enrichment()           → pathway enrichment per fate

PairScorer(adata, root=..., branches=[...],
           condition_obs_key=..., obs_key=...)
       │
       ├── build_embedding() / fit() / refit_pseudotime()
       │
       ├── score_all_conditions()     → dict[condition → CommitmentScoreResult]
       ├── compute_delta_CS()         → ΔnCS matrix with bootstrap CI
       ├── compare_conditions()       → permutation test (k=2)
       ├── fit_mixed_model()          → LMM with replicate random effect
       ├── trajectory_shift()         → KS + Wasserstein pseudotime shift
       │
       ├── plot_star_grid()           → side-by-side star panels per condition
       ├── plot_rose_grid()           → side-by-side rose plots per condition
       ├── plot_affinity_distributions() → violin/box/strip per fate
       ├── plot_delta_cs_heatmap()    → diverging ΔCS heatmap with CI
       ├── plot_compare_conditions_bar() → grouped nCS bar chart
       ├── plot_commitment_vector_radar() → radar chart of commitment vectors
       └── plot_trajectory_shift()    → KDE pseudotime distributions

MultiScorer(adata, root=..., branches=[...],
            condition_obs_key=..., obs_key=...)
       │
       ├── build_embedding() / fit()
       │
       ├── score_all_conditions()     → dict[condition → CommitmentScoreResult]
       │
       ├── Tier 2: Omnibus + post-hoc
       │   ├── compare_omnibus()      → Kruskal-Wallis / ANOVA per fate
       │   ├── compare_posthoc()      → Dunn / Tukey / Conover pairwise
       │   └── compute_pairwise_deltas() → ΔCS for ALL condition pairs
       │
       ├── Tier 3: Mixed models + trajectory
       │   ├── fit_mixed_model()      → LMM with replicate random effect
       │   ├── fit_mixed_model_contrasts() → custom condition contrasts
       │   └── trajectory_shift()     → KS + Wasserstein pseudotime shift
       │
       ├── plot_omnibus_summary()     → fates × conditions heatmap + significance
       ├── plot_posthoc_heatmap()     → condition × condition p-value heatmap
       ├── plot_pairwise_delta_grid() → grid of ΔCS heatmaps for all pairs
       └── [all PairScorer plots also available]

Citation

If you use scCS in your research, please cite:

Kriukov et al. (2025) Single-cell transcriptome of myeloid cells in response to transplantation of human retinal neurons reveals reversibility of microglial activation. DOI: 10.XXXX