Bio Single Cell Differential Abundance

Version Compatibility

Reference examples tested with: miloR 2.0+, scCODA 0.1.9+, sccomp 1.8+, speckle 1.0+

Before using code patterns, verify installed versions match. If versions differ:

• Python: pip show <package> then help(module.function) to check signatures
• R: packageVersion('<pkg>') then ?function_name to verify parameters

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

Differential Abundance Testing

"Did cell-type proportions change between conditions?" -> Test whether populations expanded or contracted between groups, accounting for the fact that proportions are not independent.

• R (cluster-free): miloR - build kNN graph, define neighborhoods, testNhoods() with a GLM and SpatialFDR
• Python/R (cluster-based): scCODA (Bayesian Dirichlet-multinomial), sccomp (Bayesian, outlier-robust), propeller (speckle, arcsin-sqrt + limma)

Governing principle

Composition data live on a SIMPLEX: proportions sum to 1, so they are NOT independent - when one population expands, every other proportion is mechanically forced down even if its absolute count never changed. Running a per-cluster t-test (or Wilcoxon) on proportions across samples is therefore invalid: it ignores the negative correlation the constraint imposes, treats each cell type as a free measurement, and produces correlated false positives (one true expansion drags down the rest, which then test as spurious "depletions"). Valid methods model the joint composition: either a Dirichlet-multinomial / log-ratio model with a reference (scCODA, sccomp) or a variance-stabilizing transform plus a linear model (propeller), or they sidestep hard clusters entirely by testing abundance on the kNN graph (Milo). Replicates are samples, not cells. The unit of replication for a composition claim is the biological sample/donor; thousands of cells from one donor are one draw. Differential abundance needs biological replicates per condition (Milo, scCODA, sccomp, propeller all model sample-level counts), and few replicates (n<3-4/group) leave abundance shifts underpowered and unstable - more donors help, more cells per donor barely do. With n=1 per condition the donor is perfectly confounded with condition: the effect is unidentifiable, not merely underpowered, yet scCODA and sccomp will still emit confident credible_effects() that are pure donor idiosyncrasy - require >=2 (ideally 3-4) biological replicates per group before believing any abundance call. Differential abundance and differential expression are different questions and confound each other. A pseudobulk or cluster-level "DE" signal between conditions can be pure composition: if a cluster mixes substates and treatment shifts their ratio, the aggregated profile changes although no gene changed expression in any cell - differential abundance masquerading as differential expression, invisible if only DE is run. Always pair a condition-DE analysis (single-cell/markers-annotation, differential-expression/deseq2-basics) with a differential-abundance test and interpret them jointly.

Choosing a differential-abundance method

scCODA and sccomp are cluster-based and Bayesian and report credible/FDR-controlled effects; Milo is cluster-free and catches shifts within a cell type that clustering averages away; propeller is the fast frequentist option. Run a cluster-based method and Milo when feasible and reconcile. When methods compete, verify current best practice against installed docs.

The reference-cell-type choice in scCODA

Compositional analysis is always relative to something. scCODA fixes one cell type as the reference assumed unchanged by the covariates, and reports every other type's change relative to it; the verdict can flip with a different reference. Choose a cell type that is biologically stable and abundant across all samples, or use reference_cell_type='automatic' (scCODA picks a type with low dispersion present in all samples). A reference that actually changes will bias all other calls. sccomp avoids a hard reference by modeling all groups jointly; Milo avoids it via the graph.

Adjusting for nuisance covariates and confounded designs

Goal: Adjust the abundance model for technical or biological nuisances (sequencing batch, timing, sex, age) and recognize when adjustment cannot help.

Approach: Add the nuisance as an extra additive term in the model formula with the condition of interest last; the test then reports the condition effect holding the nuisance constant. The nuisance column must vary within each condition - if a batch is perfectly confounded with condition (e.g. all controls sequenced in batch 1, all treated in batch 2), the term is unidentifiable and the test is invalid; the fix is experimental (multiplex conditions across batches), not statistical.

# Milo: batch added before condition; batch column lives in design.df
design <- distinct(as.data.frame(colData(milo))[, c('sample', 'batch', 'condition')])
rownames(design) <- design$sample
da <- testNhoods(milo, design = ~ batch + condition, design.df = design, reduced.dim = 'PCA')

# scCODA: additive patsy formula; covariate columns must be in the count table
data = dat.from_pandas(counts, covariate_columns=['sample', 'batch', 'condition'])
model = mod.CompositionalAnalysis(data, formula='batch + condition', reference_cell_type='automatic')

# sccomp: nuisance added to formula_composition (and optionally formula_variability)
res <- sccomp_estimate(counts_tbl, formula_composition = ~ batch + condition, .sample = sample, .cell_group = cell_type, .count = count, cores = 1)

Diagnose confounding before modeling: cross-tabulate batch x condition; if a batch maps to a single condition, no covariate term recovers the effect. Build Milo's kNN graph on a batch-corrected embedding, but keep batch in the GLM design as well, since integration and design adjustment address different residual structure.