Skip to content

API Reference

Model classes

scdef.scDEF

Bases: object

Single-cell Deep Exponential Families (scDEF) model.

scDEF learns hierarchical, multi-level gene expression signatures from single-cell RNA-seq data provided in an AnnData object. This model can be used for a variety of analyses including dimensionality reduction, batch correction, clustering, and visualization of cell states and gene programs.

The model fits multiple layers of latent factors ("gene signatures") to describe cellular heterogeneity at different resolutions. It supports batch correction, prior specification, and generation of corrected gene expression matrices.

Model fitting, inference routines, and additional plotting utilities are implemented as methods of this class. The stored AnnData object is updated with model results during training.

Parameters:

Name Type Description Default
adata AnnData

AnnData object containing the single-cell gene expression count matrix. Counts should be present in either adata.X or in the specified adata.layers.

required
counts_layer Optional[str]

key for adata.layers specifying which layer to use as expression counts (if not adata.X).

None
batch_key Optional[str]

key in adata.obs containing batch annotations; if provided, batch correction is performed. If None or not found, no batch correction is used.

None
batch_gene_scale bool

whether batch_key also makes gene_scale per batch. batch_key controls two independent things, and this flag decouples them:

  • The cell side (always on when batch_key resolves to >= 2 batches): batch_lib_sizes / batch_lib_ratio, i.e. the per-batch Gamma prior that cell_scale shrinks toward, plus batch_indices_onehot. This term is gene-independent, so it can only express sequencing depth, never a gene programme.
  • The gene side (this flag): gene_scale gets one row per batch, and gene_ratio / gene_ratio_init become per batch. This term is gene-specific, so it can and does absorb biology as well as technical differences.

True (default) is the historical behaviour. False keeps gene_scale a single shared row of shape (1, n_genes) with the pooled count-derived prior, exactly as when no batch_key is given, while still fitting the per-batch cell_scale prior.

True
seed Optional[int]

random seed for model initialization and stochastic routines (uses JAX's pseudo-random number generator).

42
n_factors Optional[int]

number of latent factors at the lowest layer (L0), denoted K_0 in the geometric layer-size schedule when layer_sizes is None.

100
n_layers Optional[float]

number of layers in the geometric schedule from L0 through the layer of width top_factors when layer_sizes is None. A size-1 root appended when top_factors > 1 does not count toward n_layers.

6
top_factors int

target width at the coarsest non-root layer (default 1), used with n_factors and n_layers in the geometric ladder K_l = K_0 * (K_top / K_0) ** (l / (n_layers - 1)) for l = 0, ..., n_layers - 1 (last rung fixed at K_top). When top_factors > 1, a final root layer of width 1 is appended and is not included in n_layers.

1
layer_sizes Optional[list]

explicit list of the number of factors in each scDEF layer. If None, layer sizes are set automatically.

None
layer_names Optional[list]

list of custom names for the layers. If None, layer names are enumerated as ["L0", "L1", ...].

None
logginglevel Optional[int]

verbosity level for the logger.

INFO
alpha Optional[float]

concentration parameter for the Gamma prior on z.

1.0
shrinkage_shape Optional[float]

shape parameter for shrinkage prior controlling factor usage.

1.0
shrinkage_rate Optional[float]

rate parameter for shrinkage prior controlling factor usage.

1.0
shrinkage_mean Optional[float]

target prior mean for shrinkage/factor relevance.

1.0
top_alpha Optional[float]

concentration parameter for the top layer Dirichlet prior over factor proportions.

1.0
factor_shape Optional[float]

shape of the prior distribution for factor-gene loadings matrix W.

0.1
brd_strength Optional[float]

BRD (Batch Relevance Determination) prior concentration parameter for factor relevance estimation.

1.0
brd_mean Optional[float]

mean of the BRD prior for factor relevance estimation.

1.0
use_brd Optional[bool]

if True, use BRD prior for automatic selection of active factors.

True
cell_scale_shape Optional[float]

precision/concentration parameter for cell-specific scaling priors.

1.0
gene_scale_shape Optional[float]

precision/concentration parameter for gene-specific scaling priors.

1.0
batch_cpal Optional[str]

default matplotlib color palette name used for batches.

'Dark2'
layer_cpal Optional[str]

matplotlib color palette for factors/colors at each scDEF layer.

'tab10'
lightness_mult Optional[float]

lightness multiplier to define the base color for each new scDEF layer.

0.15
set_alpha_from_cov Optional[bool]

if True, set the alpha parameter from the data coverage.

True
hierarchy_weight Optional[float]

multiplier applied to coverage-derived alpha when set_alpha_from_cov is True (ignored otherwise).

1.0
marginalize_alpha Optional[bool]

if True, infer alpha from a variational posterior during training.

False

n_gene_scale_batches property

Number of rows in gene_scale (1 when the gene side is shared).

This is n_batches under the default batch_gene_scale=True and 1 when the per-batch gene side is switched off, independently of how many batches the cell side uses.

add_batch_correction(reference_model, batch_key, **kwargs) classmethod

Warm-start a batch-corrected model from a fitted hierarchy.

See scdef.add_batch_correction for the full parameter documentation.

annotate_adata()

Write the fitted quantities onto self.adata.

Called automatically at the end of :meth:fit, so it normally does not need calling by hand — do so after changing factor_lists or the posterior means outside a fit.

Writes obs['cell_scale'] and var['gene_scale']; per layer, the hard assignment obs['<layer>'] with its score, the per-factor weights obs['<factor>_prob'], the representation obsm['X_<layer>'] with its row-normalized _probs, and the signatures uns['<layer>_signatures'] in scanpy's rank_genes_groups format.

attach_factors_to_obs(obs_key)

Attach factors to observation categories.

Parameters:

Name Type Description Default
obs_key str

key in model.adata.obs to use for attachment

required

Returns:

Type Description
List[List[str]]

list of attachment lists, one per layer

batch_elbo(rng, X, indices, local_params, global_params, num_samples, annealing_parameter, stop_gradients, stop_cell_budgets, stop_gene_budgets, alpha)

Monte Carlo estimate of the ELBO for one minibatch.

Inference internals: averaged over num_samples reparameterized draws and differentiated by JAX during training. Not needed for analysis.

clear_runtime_cache(clear_jax_cache=False)

Clear runtime-only compiled caches to reduce notebook memory pressure.

compute_factor_obs_assignment_fracs(layer_idx, factor_name, obs_key, obs_val)

Compute assignment fraction between a factor and an observation category.

Parameters:

Name Type Description Default
layer_idx int

layer index of the factor

required
factor_name str

name of the factor

required
obs_key str

key in model.adata.obs

required
obs_val str

value in obs_key to compute fraction with

required

Returns:

Type Description
float

assignment fraction value

compute_factor_obs_association_score(layer_idx, factor_name, obs_key, obs_val)

Compute association score between a factor and an observation category.

Parameters:

Name Type Description Default
layer_idx int

layer index of the factor

required
factor_name str

name of the factor

required
obs_key str

key in model.adata.obs

required
obs_val str

value in obs_key to compute association with

required

Returns:

Type Description
float

association score value

compute_factor_obs_weight_score(layer_idx, factor_name, obs_key, obs_val)

Compute weight score between a factor and an observation category.

Parameters:

Name Type Description Default
layer_idx int

layer index of the factor

required
factor_name str

name of the factor

required
obs_key str

key in model.adata.obs

required
obs_val str

value in obs_key to compute weight with

required

Returns:

Type Description
float

weight score value

compute_weight(upper_factor_name, lower_factor_name)

Compute the weight between two factors across any number of layers.

Parameters:

Name Type Description Default
upper_factor_name str

name of the upper factor

required
lower_factor_name str

name of the lower factor

required

Returns:

Type Description
float

weight value between the two factors

decompose_batch_effects(reference_model, **kwargs) classmethod

Re-learn lower layers under a frozen upper hierarchy to discover batch programs.

See scdef.decompose_batch_effects for the full parameter documentation.

elbo(rng, batch, indices, local_params, global_params, annealing_parameter, stop_gradients, stop_cell_budgets, stop_gene_budgets, alpha, min_shape=jnp.log(1e-10), max_shape=jnp.log(1000000.0), min_rate=jnp.log(1e-10), max_rate=jnp.log(10000000000.0))

Single-sample estimate of the ELBO for one minibatch.

Inference internals; :meth:batch_elbo averages this over several draws. Not needed for analysis.

filter_factors(brd_min=1.0, ard_min=0.001, clarity_min=0.5, n_eff_parents_max=1.5, brd_exceptional=None, local_l0_scores=False, batch_purity_max=None, batch_purity_soft_max=None, min_cells_upper=0.001, min_cells_lower=0.0, filter_up=True, annotate=True, upper_only=False, keep=None)

Filter irrelevant factors using BRD/ARD and hierarchy diagnostics.

Parameters:

Name Type Description Default
brd_min Optional[float]

minimum factor BRD value for layer 0 when use_brd.

1.0
ard_min Optional[float]

minimum ARD fraction of total ARD for layer 0.

0.001
clarity_min Optional[float]

minimum L0 clarity_score_01 when not using lineage avg_n_eff_parents (local_l0_scores or missing lineage columns).

0.5
n_eff_parents_max float

only when avg_n_eff_parents is present and local_l0_scores is False: keep when avg_n_eff_parents <= n_eff_parents_max (default 1.5). When brd_exceptional is set, also keep factors with BRD >= brd_exceptional. Ignored for local-L0 / clarity-only filtering.

1.5
brd_exceptional Optional[float]

if set, keep layer-0 factors with high BRD even when lineage avg_n_eff_parents exceeds n_eff_parents_max. Default None (disabled).

None
local_l0_scores bool

if True, filter by clarity_score_01 >= clarity_min; if False and avg_n_eff_parents exists, filter by n_eff_parents_max.

False
batch_purity_max Optional[float]

if set, keep layer-0 factors with hard-assignment batch_purity <= batch_purity_max (requires factor_diagnostics(..., batch_key=...)).

None
batch_purity_soft_max Optional[float]

if set, keep factors with batch_purity_soft <= batch_purity_soft_max.

None
min_cells_upper Optional[float]

minimum cells attached to upper-layer factors (fraction if <1).

0.001
min_cells_lower Optional[float]

minimum cells attached to layer-0 factors (fraction if <1).

0.0
filter_up Optional[bool]

whether to prune upper layers via inter-layer attachments.

True
annotate Optional[bool]

whether to run annotate_adata after filtering.

True
upper_only Optional[bool]

if True, only adjust upper layers (layer 0 unchanged).

False
keep Optional[Sequence[str]]

optional list of L0 factor names to retain even when they fail the BRD/ARD/hierarchy thresholds. Names are resolved like set_technical_factors (current model names or factor_obs rows via original_factor_idx).

None

fit(nmf_init=False, hierarchical_init=False, pca_key='X_pca', z_off=0.1, max_cells_init=5000, z_init_concentration=None, n_rounds=1, pretraining=False, force_decay_factor=False, root_epochs=0, collapse_l1_fraction=0.8, collapse_l0_min_factors=10, collapse_upper_min_factors=5, learn_budgets_on_refit=False, refit_top_layer=None, refit_top_factors=None, refit_rescale_relevance=True, refit_relevance_max_ratio=50.0, refit_rescale_w_by_layer_sizes=True, freeze_w=False, **kwargs)

Fit scDEF, warm-starting from a previous fit when available. Args: nmf_init: whether to initialize the model with NMF. hierarchical_init: on the first fit() call, initialize z from KMeans/Ward labels and upper-layer W containment matrices via get_hierarchical_init (L0 W uses the default prior init). Requires an existing PCA embedding in adata.obsm[pca_key]. Mutually exclusive with nmf_init. pca_key: adata.obsm key for PCA coordinates used for KMeans L0 clustering and the centroid dendrogram (default "X_pca"). z_off: inactive-cluster value in hierarchical init_z (default 0.1). Lower values reduce background score mass on non-winning factors and help low-BRD factors shrink toward zero relevance. max_cells_init: maximum number of cells to use for initialization. z_init_concentration: concentration parameter of a Gamma distribution to sample the initial z values from. Controls the initial spread of cell mass across factors: lower values explore more factors per cell. If None (default), computed as min(0.5, alpha / 10) — 10% of the prior concentration, capped at 0.5. Pass an explicit float to override. n_rounds: number of rounds to run the optimization. pretraining: whether to run pretrain before standard fit. force_decay_factor: on refit, whether to clip upper-layer sizes to the geometric template implied by n_factors, top_factors, and n_layers_schedule. If False, preserves learned per-layer dimensions and initializes all layers from previous posterior means. root_epochs: if > 0 (default 0), run a first phase with the root layer frozen, then a root-only refinement phase. In elbo, top_alpha applies to the penultimate layer z while the root is frozen, and to the root z when the root is optimized; other layers use alpha_layer (see stop_gradients on the root). Factor filter and annotate (see _learn) are skipped during the frozen-root phase and run only after the root step, using the filter / annotate values passed to fit. collapse_l1_fraction: on refit only (second and later fit() calls), scan adjacent pairs along the hierarchy (L0/L1, L1/L2, …; L0 is the bottom child layer). When the upper layer is nearly as wide as the layer below (at least this fraction of its factor count) and the lower layer exceeds collapse_l0_min_factors (L0) or collapse_upper_min_factors (L1+), drop the upper redundant layer, promote the next layer up as its replacement parent, and warm-start through composed W matrices. Ignored on the first fit(). Set to None to disable. Default 0.8. collapse_l0_min_factors: minimum L0 width required to collapse a redundant L1. collapse_upper_min_factors: minimum lower-layer width for L1/L2, L2/L3, … pairs. learn_budgets_on_refit: on refit only, keep existing cell/gene budget parameters (no re-initialization) but allow them to be optimized when True. When False (default), budgets are warm-started and frozen during the main learning phase. refit_top_layer: on refit only, truncate the hierarchy so this layer becomes the top non-root layer, preserving a final width-1 root when present and warm-starting through composed W matrices. Accepts a layer index or layer name (for example the recommended_layer_idx returned by scd.tl.find_sensible_top_layer). refit_top_factors: on refit only, rebuild a geometric hierarchy using these mixed-depth factors as the top non-root layer. Intermediate layer sizes follow the existing n_layers_schedule. refit_rescale_relevance: on refit only, rescale init_ard / init_brd warm starts so their geometric mean matches shrinkage_mean / brd_mean while preserving relative differences across factors (optionally capped by refit_relevance_max_ratio). refit_relevance_max_ratio: max/min ratio across factors after relevance rescaling on refit. Set to None to disable the ratio cap. refit_rescale_w_by_layer_sizes: on refit only, multiply each warm-started W layer by old_K / new_K so loadings match the 1 / K cold-start convention when layer widths change (for example after filter_factors). freeze_w: if True, hold every per-layer W variational parameter fixed at its current value during the entire fit. The lower-bottom gene loadings (L0W) and all parent-child layer matrices stay exactly where they are, while z, cell budgets, gene budgets, BRD, ARD, and alpha continue to learn freely. Useful for warm-starting from a fitted hierarchy (e.g. the second pass of add_batch_correction) when you want batch-specific gene scales to absorb new variance without letting the hierarchy drift. Default False. Not applied during pretraining. **kwargs: training settings, forwarded to the inference loop. The ones normally worth setting:

    * ``n_epoch`` (default ``1000``): **maximum** epochs *per round*.
      Training stops earlier once the relative improvement over the
      best loss so far, ``(best - current) / |best|``, stays below
      ``tolerance`` for ``patience`` consecutive epochs, and never
      before ``min_epochs``. With ``n_rounds > 1`` the overall
      ceiling is ``n_rounds * n_epoch``, plus any ``root_epochs``
      and pretraining epochs.
    * ``lr`` (default ``0.1``): learning rate for the global
      variational parameters; ``local_lr`` (default ``0.01``) for
      the per-cell ones.
    * ``batch_size`` (default ``256``): cells per minibatch.
    * ``num_samples`` (default ``100``): Monte Carlo samples used
      for the gradient estimate. Lower is faster and noisier.
    * ``min_epochs`` (default ``50``), ``tolerance`` (default
      ``1e-5``), ``patience`` (default ``50``): the early-stopping
      rule described above. ``tolerance`` is on the *relative*
      improvement, so it is comparable across data sets whatever the
      scale of the loss.
    * ``filter`` (default ``True``): prune the **upper** layers at
      the end of the fit, via ``filter_factors(upper_only=True)``.
      Layer 0 is deliberately left intact — filtering it is a
      separate, explicit step
      ([`filter_factors`][scdef.tl.filter_factors]), so that the thresholds
      can be chosen from the diagnostics. Has no effect when
      ``use_brd=False``.
    * ``annotate`` (default ``True``): write the cell scores and
      signatures to ``adata``, which is what puts ``X_<layer>`` in
      ``adata.obsm``. Note that under the default ``filter=True``
      this flag is not consulted: ``filter_factors`` annotates on its
      own. Passing ``annotate=False`` only suppresses annotation if
      ``filter=False`` is passed as well.

On the first call, parameters are initialized from priors (or NMF if enabled).
On subsequent calls, the model is re-initialized from the current posterior
quantities and the current `factor_lists`, enabling a fit -> filter -> fit
workflow. During refit, upper-layer sizes are clipped to the geometric
template when ``force_decay_factor`` is True before rebuilding the hierarchy.

from_hierarchy(adata, hierarchy, **kwargs) classmethod

Create a model for new data initialized from a learned hierarchy.

See scdef.from_hierarchy for the full parameter documentation.

from_reference(reference_model, adata, **kwargs) classmethod

Create a new model initialized from a fitted reference hierarchy.

See scdef.from_reference for the full parameter documentation.

get_annotations(marker_reference, gene_rankings=None)

Get annotations for factors based on marker gene reference.

Parameters:

Name Type Description Default
marker_reference Mapping[str, Sequence[str]]

dictionary mapping annotation names to gene lists

required
gene_rankings Optional[List[List[str]]]

gene rankings for each factor, if None will be computed

None

Returns:

Type Description
List[List[str]]

list of annotation lists, one per factor

get_effective_factors(brd_min=1.0, ard_min=0.001, clarity_min=0.5, n_eff_parents_max=1.5, brd_exceptional=None, local_l0_scores=False, min_cells=0.001, batch_purity_max=None, batch_purity_soft_max=None)

Indices of the layer-0 factors that pass the relevance and hierarchy criteria.

The selection :meth:filter_factors applies, exposed separately so the consequences of a threshold can be inspected before committing to it. Factors must clear brd_min and ard_min, and either clarity_min or n_eff_parents_max depending on which hierarchy diagnostics are available; brd_exceptional keeps a high-relevance factor regardless.

Returns:

Type Description
ndarray

Integer indices into layer 0.

get_hierarchical_init(pca_key='X_pca', z_on=1.0, z_off=0.1)

Build nested warm-start init_z and init_w from PCA/KMeans labels.

L0 labels come from KMeans(n_clusters=layer_sizes[0]) on adata.obsm[pca_key]. Ward linkage on the resulting L0 centroids defines nested partitions for upper layers. Only K0 centroid comparisons are needed (no O(n_cells^2) work over cells).

init_z uses soft one-hot cell assignments at every layer. init_w[0] is None (L0 W uses the default prior init). init_w[l] for l >= 1 holds parent/child containment matrices derived from the same dendrogram cuts as init_z.

Parameters:

Name Type Description Default
pca_key str

adata.obsm key for a PCA embedding (default "X_pca").

'X_pca'
z_on float

active cluster / assignment value (default 1.0).

1.0
z_off float

inactive cluster soft floor (default 0.1).

0.1

Returns:

Type Description
List[ndarray]

(init_z, init_w) lists compatible with init_var_params.

List[ndarray]

init_w[0] is None; layer l >= 1 has shape

Tuple[List[ndarray], List[ndarray]]

(layer_sizes[l], layer_sizes[l - 1]).

Raises:

Type Description
KeyError

if pca_key is missing.

ValueError

if the PCA embedding row count does not match n_cells.

get_layer_factor_orders()

Get the ordering of factors in each layer for plotting.

Returns:

Type Description
List[ndarray]

list of arrays, one per layer, containing factor indices in plotting order

get_nmf_init(max_cells=None)

NMF warm start for the layer-0 factors.

Runs scikit-learn's NMF on library-normalized counts and returns cell and gene loadings to initialize z and W. Used by fit(nmf_init=True); max_cells subsamples for speed.

get_rankings(layer_idx=0, top_genes=None, genes=True, return_scores=False, sorted_scores=True, drop_factors=None)

Get gene or factor rankings for each factor in a layer.

Parameters:

Name Type Description Default
layer_idx int

layer index to get rankings for

0
top_genes Optional[int]

number of top genes/factors to return

None
genes bool

whether to return gene rankings (True) or factor rankings (False). Gene rankings use cached confidence+mean combined scores when sorted_scores=True and drop_factors is not provided.

True
return_scores bool

whether to return scores along with rankings

False
sorted_scores bool

whether to return scores sorted by ranking

True
drop_factors Optional[List[str]]

list of factors to drop from rankings

None

Returns:

Type Description
Union[List[List[str]], Tuple[List[List[str]], List[List[float]]]]

list of rankings per factor, or tuple of (rankings, scores) if return_scores is True

get_relevances_dict()

Get dictionary of factor relevance scores.

Returns:

Type Description
Dict[str, float]

dictionary mapping factor names to relevance scores

get_signature_confidence(factor_idx, layer_idx, mc_samples=100, top_genes=10, pairwise=False)

Get confidence score for a factor signature using Monte Carlo sampling.

Parameters:

Name Type Description Default
factor_idx int

index of the factor

required
layer_idx int

layer index of the factor

required
mc_samples int

number of Monte Carlo samples to take

100
top_genes int

number of top genes to consider in each sample

10
pairwise bool

whether to compute pairwise Jaccard similarities

False

Returns:

Type Description
float

confidence score as Jaccard similarity

get_signature_sample(rng, factor_idx, layer_idx, top_genes=10, return_scores=False)

Get a single signature sample from the posterior for a factor.

Parameters:

Name Type Description Default
rng Any

JAX random number generator key

required
factor_idx int

index of the factor

required
layer_idx int

layer index of the factor

required
top_genes int

number of top genes to return

10
return_scores bool

whether to return scores along with gene names

False

Returns:

Type Description
Union[List[str], Tuple[List[str], ndarray]]

list of gene names, or tuple of (gene_names, scores) if return_scores is True

get_signatures_dict(top_genes=None, scores=False, sorted_scores=False, layer_normalize=False, drop_factors=None)

Get dictionary of gene signatures for all factors across all layers.

Parameters:

Name Type Description Default
top_genes Optional[int]

number of top genes per signature

None
scores bool

whether to return scores along with signatures

False
sorted_scores bool

whether to return scores sorted by ranking

False
layer_normalize bool

whether to normalize scores within each layer

False
drop_factors Optional[List[str]]

list of factors to exclude

None

Returns:

Type Description
Union[Dict[str, List[str]], Tuple[Dict[str, List[str]], Dict[str, ndarray]]]

dictionary mapping factor names to gene lists, or tuple of (signatures, scores) if scores is True

get_sizes_dict()

Get dictionary of factor sizes (number of cells per factor).

Returns:

Type Description
Dict[str, float]

dictionary mapping factor names to cell counts

get_summary(top_genes=10, reindex=True)

Get a text summary of the model factors and their top genes.

Parameters:

Name Type Description Default
top_genes int

number of top genes to show per factor

10
reindex bool

whether to reindex factors

True

Returns:

Type Description
str

string summary of the model

identify_mixture_factors(max_n_genes=20, thres=0.5)

Identify factors that might be better if broken apart.

Parameters:

Name Type Description Default
max_n_genes int

maximum number of genes per factor

20
thres float

threshold for identifying mixture factors

0.5

Returns:

Type Description
ndarray

array of factor indices that are mixture factors

init_var_params(init_budgets=True, init_alpha=True, init_z=None, init_w=None, init_brd=None, init_ard=None, init_gene_scale=None, nmf_init=False, z_init_concentration=0.5, **kwargs)

Initialize the variational parameters.

Called by :meth:fit on the first pass. Each init_* argument supplies a warm start for the corresponding quantity instead of drawing from the prior — this is how :func:scdef.from_reference and :func:scdef.decompose_batch_effects carry a fitted hierarchy into a new model.

load(dir_path, adata=None) classmethod

Load model from disk.

Parameters:

Name Type Description Default
dir_path Union[str, Path]

directory created by save

required
adata Optional[AnnData]

optional AnnData to attach when adata.h5ad was not saved.

None

Returns:

Type Description
scDEF

Loaded model instance.

load_adata(adata, layer=None, batch_key=None)

Attach a new AnnData to the model.

The object is copied, as at construction, so the caller's AnnData is never modified. layer names the counts layer and batch_key the batch column, with the same meaning as the constructor arguments.

make_corrected_data(layer_name='scdef_corrected')

Compute and store the low-rank reconstruction of the UMI count matrix.

The reconstructed matrix is saved to adata.layers[layer_name], providing a denoised, batch-corrected version of the expression data.

Parameters:

Name Type Description Default
layer_name str

name for the AnnData layer where the reconstructed matrix is stored

'scdef_corrected'

make_layercolors(layer_cpal='tab10', lightness_mult=0.15)

Assign a colour palette to each layer's factors.

Populates the per-layer colours the graph and UMAP plots use. layer_cpal is a matplotlib colormap name (or one per layer) and lightness_mult spreads factors within a layer by lightness.

normalize_cellscores()

Turn the per-layer cell scores into probabilities.

Row-normalizes each obsm['X_<layer>'] and writes obsm['X_<layer>_probs'] plus one obs['<factor>_prob'] column per factor. Run by :meth:annotate_adata.

pretrain(n_epoch=200, prune_alpha=None, **kwargs)

Run two-pass alpha pretraining without the final full fit pass.

Schedule: 1) Full model pass with low alpha (=1) 2) Full model pass with high alpha (pruning pressure)

reinit_factors(mixture_factors=None, init_budgets=False, exponent=1.1, **kwargs)

Re-initialize the layer-0 factors, splitting mixture factors.

Factors that mix several programmes are identified (or supplied via mixture_factors) and re-seeded, so a further fit can resolve them into separate factors instead of leaving them blended.

save(dir_path, overwrite=False, save_anndata=False)

Save model state to disk, similarly to scvi-tools.

This writes a model state pickle plus metadata to dir_path. AnnData is saved separately as adata.h5ad only when save_anndata=True.

Parameters:

Name Type Description Default
dir_path Union[str, Path]

output directory path

required
overwrite bool

whether to overwrite an existing non-empty directory

False
save_anndata bool

whether to save model.adata as adata.h5ad

False

set_factor_names()

Rebuild factor_names from the current factor_lists.

Names are <layer>_<i> numbered contiguously within each layer, so they change whenever factors are filtered — which is why anything keyed by factor name must be recomputed after filtering.

set_posterior_means()

Compute self.pmeans from the variational parameters.

The posterior means every downstream tool reads: L0z/L0W and the upper-layer equivalents, cell_scale, gene_scale, brd. Refreshed at the end of a fit.

set_posterior_variances()

Compute self.pvars, the posterior variances.

Needed by anything that reports uncertainty rather than a point estimate — the confident signatures, their per-gene confidences, and :func:scdef.plotting.factor_gene_uncertainty_boxplot. Call it after loading a saved model if those quantities were not stored.

update_model_priors(update_alpha_from_cov=True)

Recompute the prior hyperparameters from the data and layer sizes.

Called after construction and whenever the layer sizes change, so the relevance and hierarchy priors stay matched to the current architecture. With update_alpha_from_cov the layer concentration is additionally rescaled from the observed gene coverage.

update_model_size(max_n_factors=None, n_layers=None, layer_sizes=None, use_decay_factor_schedule=False)

Update latent hierarchy dimensions.

Parameters:

Name Type Description Default
max_n_factors Optional[int]

bottom-layer factor count when use_decay_factor_schedule is True (iscDEF marker layer 0 path).

None
n_layers Optional[int]

target number of geometric layers from K0 through top_factors (a final root of size 1 when top_factors > 1 is not counted), or maximum layers for the decay schedule.

None
layer_sizes Optional[List[int]]

explicit per-layer sizes. If provided, sizes are sanitized to be non-increasing and consecutive duplicates are collapsed.

None
use_decay_factor_schedule bool

if True, use decay_factor-based halving (iscDEF only). If False, use n_factors, top_factors, and n_layers_schedule on self for a geometric ladder (scDEF).

False

scdef.iscDEF

Bases: scDEF

Informed Single-cell Deep Exponential Families (iscDEF) model.

iscDEF extends scDEF with marker gene sets that shape the layer-0 loading prior W. Each typed factor is encouraged to use its own marker genes, may load other genes more weakly, and can share biology with add_other residual factors.

Tuning how much the model relies on markers vs. augments signatures

Marker reliance is controlled mainly by the Gamma prior on W at the markers layer (for markers_layer=0, that is L0). For a typed factor, each listed marker gene has prior mean loading ≈ gs_big_scale (tighter when marker_strength is high); all other genes default to ≈ gs_small_scale. Fitted signatures can still add genes beyond your list if the data and these priors allow it.

  • Stay close to the input marker lists (typing, strict gene programs): increase gs_big_scale and marker_strength; keep the default penalize_other=True so off-type markers are discouraged on the wrong factor; keep add_other small when you only have a few types to avoid ignoring the marker factors.

  • Augment markers with data-driven genes: decrease gs_big_scale and marker_strength; increase gs_small_scale or nonmarker_strength so non-marker genes are not overly suppressed on typed factors; use add_other ≥ 1 for programs not in markers_dict; set penalize_other=False if overlapping lists should not hard-reject shared genes.

List design matters as much as numeric knobs: non-overlapping marker sets per type separate factors more reliably than tuning alone. Use markers_layer=0 for one factor per type; use markers_layer>0 for coarse types at the top and finer substructure at L0.

Parameters:

Name Type Description Default
adata AnnData

AnnData object containing the gene expression count matrix. Counts must be present in either adata.X or a specified layer.

required
markers_dict Mapping[str, Sequence[str]]

dictionary mapping marker/factor names to gene lists (gene sets). These guide the formation of factors in the chosen layer.

required
add_other Optional[int]

if > 0, adds one or more other{i} residual categories. At markers_layer=0, each is a separate L0 factor. At markers_layer>0, each gets a block of L0 sub-factors (add_other * n_factors_per_marker columns at L0) and one coarse factor at the marker layer.

0
markers_layer Optional[int]

index of the layer at which gene sets are enforced as factors (0 = lowest/finest, higher = top layer). If > 0, total layers determined by this value.

0
add_root Optional[bool]

when markers_layer > 0 (default True), append a width-1 root above the marker layer. Fitting runs in two phases (main fit with frozen root, then root_epochs on the root only); see fit. Ignored when markers_layer=0.

None
cn_small_mean Optional[float]

mean prior connectivity for "small" (weakly-connected) genes between factors and gene sets.

1.0
cn_big_mean Optional[float]

mean prior connectivity for "big" (strongly-connected) genes between factors and gene sets.

10.0
cn_small_strength Optional[float]

concentration parameter for low connectivity (see scDEF prior specification).

0.1
cn_big_strength Optional[float]

concentration parameter for high connectivity.

1.0
gs_small_scale Optional[float]

prior mean scale for genes not in a factor's marker list (higher → more non-marker loading).

1.0
gs_big_scale Optional[float]

prior mean scale for genes in that factor's marker list (higher → stronger marker reliance).

10.0
marker_strength Optional[float]

Gamma prior concentration on marker-gene loadings (higher → less augmentation away from markers).

1.0
nonmarker_strength Optional[float]

prior concentration on non-marker loadings (higher → tighter; overridden to 1.0 if use_brd).

0.1
other_strength Optional[float]

prior concentration when penalizing marker genes on the wrong factor or on other rows.

0.1
penalize_other Optional[bool]

if True (default), typed factors penalize other groups' marker genes; other factors penalize all typed markers.

True
**kwargs Any

additional arguments passed to scDEF. hierarchy_weight defaults to 0.25 (scales coverage-derived alpha when set_alpha_from_cov=True).

{}

filter_factors(brd_min=1.0, ard_min=0.001, clarity_min=0.5, n_eff_parents_max=1.5, brd_exceptional=None, local_l0_scores=False, batch_purity_max=None, batch_purity_soft_max=None, min_cells_upper=0.001, min_cells_lower=0.0, filter_up=True, annotate=True, upper_only=False, keep=None)

Filter factors while preserving existing marker-based factor names.

This override keeps the base filtering behavior but restores names by subsetting the previous factor_names. This avoids marker-prefix relabeling across filter/refit workflows.

fit(nmf_init=False, max_cells_init=1024, z_init_concentration=10.0, z_init_from_score_genes=True, z_init_score_temperature=1.0, z_init_other_mass=None, z_init_other_mode='inverse_union', score_genes_layer=None, score_genes_kwargs=None, root_epochs=0, **kwargs)

Fit iscDEF, warm-starting from previous fit when available.

On refit, all layers are initialized from the previous posterior means (z and W), while BRD/ARD are initialized from layer 0. Existing marker-aware names are preserved through the refit path.

On the first fit, z can be initialized from Scanpy score_genes on typed marker sets. For markers_layer == 0 only layer 0 is set this way. For markers_layer > 0, every layer uses the same score softmax among typed markers, replicated uniformly within each marker block at that layer (and uniform other blocks when add_other is used). With z_init_other_mass=None (default), the fraction of mass on other columns at each layer equals the number of other factors at that layer divided by that layer width. With z_init_other_mode='inverse_union' (default), per-cell other init is scaled by 1 / (1 + union_marker_score) from score_genes on the union of typed marker genes; use z_init_other_mode='uniform' for constant other init. nmf_init is not used by iscDEF.

When add_root=True (markers_layer > 0 only), a width-1 root is appended and fitting runs in two phases: all layers except the root, then root_epochs on the root only (default 10 when add_root and root_epochs is 0).

Building models from a fitted model

These build a new scDEF from an existing one — or, for from_hierarchy, from a hierarchy that need not come from a fitted model. Each is also available as a classmethod on the class above (for example scdef.scDEF.decompose_batch_effects), which forwards to the function documented here.

scdef.from_reference(reference_model, adata, counts_layer=None, batch_key=None, reference_obs=None, query_obs=None, copy_cell_z=True, init_gene_scale='batch', **kwargs)

Create a new model initialized from a fitted reference hierarchy.

The new model uses adata as its data matrix and initializes global hierarchy parameters (W, BRD/ARD, alpha-related hyperparameters) from reference_model. Cell/gene budgets are initialized from the new data so modality/batch-specific scales can be learned.

Parameters:

Name Type Description Default
reference_model 'scDEF'

a fitted scDEF providing the hierarchy.

required
adata AnnData

AnnData for the new model.

required
counts_layer Optional[str]

counts layer key in adata.

None
batch_key Optional[str]

batch annotation column in adata.obs.

None
reference_obs Optional[str]

reference batch label (for gene-scale init).

None
query_obs Optional[str]

query batch label (for gene-scale init).

None
copy_cell_z bool

copy per-cell z warm starts for shared cells.

True
init_gene_scale Union[str, ndarray]

how to initialize per-batch gene_scale variational means before the first fit.

  • 'batch' (default): use per-batch count means from load_adata (1 / gene_ratio_init).
  • 'reference': broadcast the reference model's fitted pmeans['gene_scale'] to every batch.
  • array: explicit (n_genes,) or (n_batches, n_genes) means.
'batch'
**kwargs Any

additional keyword arguments forwarded to the model constructor.

{}

Returns:

Type Description
'scDEF'

A new (unfitted) scDEF model with hierarchy warm-started from the reference.

scdef.add_batch_correction(reference_model, batch_key, *, adata=None, counts_layer=None, copy_cell_z=True, freeze_w=False, learn_budgets=True, n_epoch=400, lr=0.05, tolerance=0.0001, from_reference_kwargs=None, **fit_kwargs)

Warm-start a batch-corrected model from a fitted hierarchy.

Designed for the workflow: 1. Fit reference_model without a batch_key to learn the factor hierarchy on the unbatched signal (optionally followed by filter_factors). 2. Call this function to construct a new model that shares the same hierarchy (factor_lists, layer sizes, W, BRD, ARD) and re-fits it under per-batch gene-scale priors so batch effects are absorbed by gene_scale, not by the hierarchy.

Parameters:

Name Type Description Default
reference_model 'scDEF'

a fitted scDEF providing the hierarchy.

required
batch_key str

column in adata.obs to use as the new batch annotation.

required
adata Optional[AnnData]

AnnData for the second pass. Defaults to reference_model.adata.

None
counts_layer Optional[str]

counts layer for the new adata.

None
copy_cell_z bool

whether to copy per-cell z warm starts for shared cells.

True
freeze_w bool

hold every per-layer W fixed during the second fit.

False
learn_budgets bool

allow per-batch gene-scale and per-cell budgets to move.

True
n_epoch int

epochs for the second-pass fit.

400
lr float

learning rate for the second-pass fit.

0.05
tolerance float

early-stopping tolerance for the second-pass fit.

0.0001
from_reference_kwargs Optional[Mapping[str, Any]]

extra kwargs forwarded to from_reference.

None
**fit_kwargs Any

additional kwargs forwarded to model.fit().

{}

Returns:

Type Description
'scDEF'

The new fitted model with batch correction applied.

scdef.decompose_batch_effects(reference_model, *, adata=None, counts_layer=None, batch_cell_scale=True, top_layer=1, n_epoch=400, lr=0.05, tolerance=0.0001, nmf_init=False, init_gene_scale='reference', **fit_kwargs)

Re-learn lower layers under a frozen upper hierarchy to discover batch programs.

Two-stage workflow:

  1. reference_model was fitted with a batch_key, producing a hierarchy where per-batch gene_scale absorbed between-batch variance.
  2. This function creates a new model with the per-batch gene_scale switched off, warm-starts all W from the reference, and re-learns all layers up to top_layer. At the boundary (top_layer), only W is re-learned while z stays fixed — preserving the cell-to-group assignments as the structural constraint. Layers below top_layer are fully re-learned (both W and z). Layers above top_layer remain completely fixed.

Which half of batch_key is discarded. batch_key in scDEF controls two independent quantities:

  • the gene sidegene_scale with one row per batch, i.e. a gene-specific per-batch multiplier;
  • the cell sidebatch_lib_sizes / batch_lib_ratio, the per-batch Gamma prior that cell_scale shrinks toward, i.e. a gene-independent per-batch multiplier.

Only the gene side must go. A gene-specific term can express a gene programme, so leaving it in place would let it re-absorb exactly the structure this function is trying to surface in the L0 factors. A gene-independent term cannot express a programme at all — it is one number per batch, so the most it can represent is sequencing depth. Historically both were discarded together; batch_cell_scale=True (the default) now keeps the cell side.

Keeping the cell side is expected to reduce the depth component that would otherwise have to land in z: with batch_key=None the model shrinks every cell_scale toward a single Gamma fitted to the pooled library-size mean and variance, which is misspecified for both batches when they differ in depth. This expectation has not been validated by refitting; no claim is made here about the effect on the resulting factors.

L0 factor BRD and ARD are re-initialized from model priors rather than copied from the reference, so factor relevance can be re-estimated during decomposition.

How the shared gene_scale is warm-started. The count-derived prior alone leaves a large reconstruction gap here, because the batch-key model often learns per-batch scales orders of magnitude above it. The default (init_gene_scale='reference') instead solves for the shared scale directly: with the likelihood mean (z @ W) * cell_scale * gene_scale, summing over cells at this model's own initialization gives U_g = ((cell_scale^T z_0) W_0)_g, and s_g = X.sum(0)_g / U_g is the maximum-likelihood shared scale for the pooled marginal of gene g. It is computed from the target data and this model's warm starts, so it needs no transplant of the reference's fitted scale, and it applies whether the reference had one gene_scale row or several.

This is exact for the pooled marginal only. A single shared row still cannot fit two batches whose levels for a gene differ — no shared scale can. That residual is left deliberately unabsorbed: pushing it into z @ W is what makes per-batch structure visible in the re-learned lower layers, which is the point of the decomposition. What the pooled MLE removes is the part that is not structure, namely a systematic level offset in the warm start.

With top_layer=1 (default): - L0: W warm-started and re-learned, z re-learned - L1: W warm-started and re-learned, z frozen - L2+: fully frozen

With top_layer=2: - L0: W warm-started and re-learned, z re-learned - L1: W warm-started and re-learned, z re-learned - L2: W warm-started and re-learned, z frozen - L3+: fully frozen

Parameters:

Name Type Description Default
reference_model 'scDEF'

a fitted scDEF that was trained with batch_key.

required
adata Optional[AnnData]

AnnData for the second stage. Defaults to reference_model.adata.

None
counts_layer Optional[str]

counts layer for adata.

None
batch_cell_scale bool

if True (default), carry the reference model's batch_key into the decomposed model with batch_gene_scale=False, so the per-batch cell_scale prior is kept while gene_scale stays a single shared row. Degrades gracefully (with a log message) to no batch key when the reference has no batch_key, when that key is absent from the target adata.obs, or when fewer than two batches are present in the cells being fitted. Set to False to reproduce the historical construction, which discarded both sides.

True
top_layer int

the highest layer whose W is re-learned. Its z remains frozen as the structural anchor. Default 1.

1
n_epoch int

training epochs for the re-learning phase.

400
lr float

learning rate for the re-learning phase.

0.05
tolerance float

early-stopping tolerance.

0.0001
nmf_init bool

if True, initialize L0 W via NMF on the data instead of warm-starting from the reference. Default False.

False
init_gene_scale Union[str, ndarray]

warm start for the shared gene_scale in the decomposed model.

  • "reference" (default): the pooled-marginal MLE X.sum(0)_g / ((cell_scale^T z_0) W_0)_g, which reproduces the observed pooled counts exactly at initialization (see above). Falls back to the geometric mean of reference_model.pmeans['gene_scale'] across batches when nmf_init=True — layer-0 W is then set by NMF, so a profile derived from the reference W would not describe the model's actual initialization — and, per gene, for genes with no factor support.
  • "prior": only the count-derived prior mean (1 / gene_ratio_init).
  • an explicit (n_genes,) array.
'reference'
**fit_kwargs Any

additional keyword arguments forwarded to _learn.

{}

Returns:

Type Description
'scDEF'

A new fitted model whose lower-layer factors reveal batch-specific

'scDEF'

and shared gene programs under the frozen upper-layer cell assignments.

scdef.from_hierarchy(adata, hierarchy, counts_layer=None, batch_key=None, init_brd=None, init_ard=None, init_z=None, **kwargs)

Create a model for new data initialized from a learned hierarchy.

hierarchy can be either a fitted scDEF model (preferred) or an explicit sequence of W matrices. When a model is passed, current factor_lists are respected and the corresponding W submatrices, BRD/ARD, and hyperparameters are copied.

Parameters:

Name Type Description Default
adata AnnData

AnnData for the new model.

required
hierarchy Union['scDEF', Sequence[ndarray]]

a fitted scDEF model or a sequence of W matrices.

required
counts_layer Optional[str]

counts layer key in adata.

None
batch_key Optional[str]

batch annotation column in adata.obs.

None
init_brd Optional[ndarray]

explicit BRD initialization (overrides reference).

None
init_ard Optional[ndarray]

explicit ARD initialization (overrides reference).

None
init_z Optional[Sequence[ndarray]]

explicit per-layer z initialization.

None
**kwargs Any

additional keyword arguments forwarded to the model constructor.

{}

Returns:

Type Description
'scDEF'

A new (unfitted) scDEF model initialized from the hierarchy.

Tools

scdef.tl

Tooling utilities for scDEF.

add_l0_lineage_aggregate_scores(per_factor, layer_names)

Add lineage-averaged clarity and effective parents for layer-0 factors only.

For each L0 factor, follows best_parent upward through layers L1 … L{n-2}, collecting clarity_score_01 and n_eff_parents from each factor along the path (same definitions as compute_hierarchy_scores). Stores the mean over those positions in avg_clarity and avg_n_eff_parents on the L0 row only; other rows get NaN.

This is called automatically from compute_hierarchy_scores; it remains public for advanced use on a pre-built per_factor frame.

This captures cases where an L0 factor maps cleanly to one L1 parent (low local n_eff_parents) while that parent is ambiguous relative to L2, by letting lineage averages reflect uncertainty higher in the hierarchy.

Parameters:

Name Type Description Default
per_factor DataFrame

per-factor scores from compute_hierarchy_scores. Factor identity for lookups uses row index labels (not the child_factor column when present).

required
layer_names Sequence[str]

ordered model layer names (model.layer_names). The walk length is len(layer_names) - 1 (one score per child layer, L0 through L{n-2}); layer order is not inferred from strings in the frame so non-lexicographic names stay correct.

required

Returns:

Type Description
DataFrame

Copy of per_factor with two additional float columns.

annotate_factors(model, annotations)

Attach descriptive annotations to factors in adata.uns['factor_obs'].

Annotations are stored in the annotation column of factor_obs, keyed by the resolved factor rows (see _resolve_factor_obs_names). Factor names may be current model names (e.g. L0_4) even after filtering.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance.

required
annotations Mapping[str, str]

mapping {factor_name: description}.

required

Returns:

Type Description
DataFrame

Updated factor_obs dataframe.

assign_confident(model, n_samples=500, tau=0.3, credible_level=0.9, key_added='confident', rng_key=None, exclude_technical=False, exclude_batch_technical=False, batch_technical_top_layer=None)

Pick the finest scDEF layer at which each cell is confidently assigned.

For each cell c and layer k (restricted to filtered factors model.factor_lists[k]), draws n_samples reparameterized samples z^(s) ~ q(z_{c,k}) from the log-normal variational posterior and normalizes each sample: ẑ^(s) = z^(s) / sum(z^(s)). exclude_technical / exclude_batch_technical narrow that per-layer candidate set further, and all scores below are then computed among the remaining factors only.

The confidence score is defined on the gap between the cell-level winner and its nearest competitor, so the score is invariant to layer size K_k (only the top two factors enter, dormant factors don't inflate it):

  • Cell-level winner f* = argmax_f E_s[ẑ_f].
  • Per-sample gap gap^(s) = ẑ^(s)_{f*} - max_{g != f*} ẑ^(s)_g. Equals 0 when the winner and runner-up are tied, 1 when the winner holds all the mass, and goes negative in samples where some other factor out-competes f* (so identity flipping directly penalizes the score).
  • Effect size = E_s[gap^(s)] — posterior-mean margin by which the winner beats its nearest competitor.
  • Posterior SD = SD_s[gap^(s)] — uncertainty of that margin (diagnostic).
  • Confidence = quantile_{1 - credible_level}({gap^(s)}) — the empirical (1 - credible_level)-quantile of the gap across posterior samples. Reads as:

    "In at least credible_level of posterior samples, the winning factor f* had at least confidence more normalized mass than any competing factor."

Layer-size invariant: a value of e.g. 0.3 means "winner leads runner-up by 30 percentage points of normalized mass" regardless of the number of factors at the layer.

Two auxiliary diagnostic scores are also computed:

  • winner_probability[c, k] = max_f P_s[argmax = f] — how often the winner is the argmax across samples. Identity-stability only; ignores magnitude (a tight 0.51/0.49 scores 1.0 here).
  • entropy_confidence[c, k] = 1 - H(p) / log(K_k) — normalized entropy of the argmax distribution.

Selection rule (finest-that-clears). For each cell, the "best" layer is the finest (lowest index) multi-factor layer (K_k >= 2) whose confidence clears tau. If no multi-factor layer clears tau, the cell is assigned to the top-most single-factor layer (typically the root) as a "stem-cell-like" catch-all. This encodes the biological intuition:

  • Terminally differentiated cells have a clear dominant factor at L0 with tight posterior → assigned at L0.
  • Partially differentiated cells are ambiguous between siblings at L0 but clear at their parent layer L1 → assigned at L1.
  • Stem-like cells are ambiguous everywhere → assigned at the root.

Writes:

  • adata.obsm[f"{key_added}_effect_size"](n_cells, n_layers) float. Posterior-mean gap between winner and nearest competitor.
  • adata.obsm[f"{key_added}_posterior_sd"](n_cells, n_layers) float. Posterior SD of the gap.
  • adata.obsm[f"{key_added}_confidence"](n_cells, n_layers) float. Lower empirical quantile of the gap — the layer-size-invariant score that tau gates on.
  • adata.obsm[f"{key_added}_winner_mass"](n_cells, n_layers) float. Diagnostic: posterior mean of the winner's normalized mass (E[ẑ_{f*}]). Not K-invariant.
  • adata.obsm[f"{key_added}_winner_probability"](n_cells, n_layers) float. Diagnostic: posterior argmax-identity probability.
  • adata.obsm[f"{key_added}_entropy_confidence"](n_cells, n_layers) float. Diagnostic: normalized-entropy score.
  • adata.obsm[f"{key_added}_argmax_factor"](n_cells, n_layers) int, slot indices into the effective candidate list of each layer — the filtered factor list minus any factors excluded by exclude_technical / exclude_batch_technical (-1 if the layer has no candidates left).
  • adata.obs[f"{key_added}_confidence_{layer_name}"] per layer.
  • adata.obs[f"{key_added}_argmax_{layer_name}"] per layer (factor name).
  • adata.obs[f"{key_added}_best_layer"] — layer name of the chosen layer.
  • adata.obs[f"{key_added}_best_layer_idx"] — integer index of that layer in model.layer_names (0 = finest). -1 if no layer was chosen.
  • adata.obs[f"{key_added}_best_factor_idx"] — slot within the best layer's filtered factor list.
  • adata.obs[f"{key_added}_factor"] — factor name at the best layer.
  • adata.obs[f"{key_added}_best_effect_size"] — effect size at the best layer.
  • adata.obs[f"{key_added}_best_posterior_sd"] — posterior SD at the best layer.
  • adata.obs[f"{key_added}_best_confidence"] — combined confidence at the best layer.
  • adata.obs[f"{key_added}_depth_score"] — assignment-centric depth in [0, 1]: best_layer_index / (n_layers - 1), so 0 at the finest layer (L0) and 1 at the coarsest layer (root). If n_layers == 1, the score is 0 for assigned cells. Cells with no valid layer (best_layer_index < 0) get NaN.
  • adata.uns[key_added] — metadata (layer names, sizes, tau, n_samples, credible_level, metric name, fit revision).

Parameters:

Name Type Description Default
model scDEF

scDEF model instance (must already be fitted).

required
n_samples int

number of Monte Carlo samples S drawn per cell/layer.

500
tau float

minimum confidence (lower quantile of the winner-runner-up gap) for a multi-factor layer to be eligible as the cell's best layer (in [0, 1]). Reads as "in at least credible_level of posterior samples, the winner must lead the runner-up by at least tau of the normalized mass".

0.3
credible_level float

posterior credibility of the lower bound, in (0, 1). The confidence is the empirical (1 - credible_level)-quantile of the per-sample gap. Default 0.9 → "in 90% of posterior samples, the winner led by at least confidence". Set higher (e.g. 0.95) for a stricter score.

0.9
key_added str

prefix used for all written keys in adata.

'confident'
rng_key Optional[Any]

optional jax.random key; if None, derived from model.seed.

None
exclude_technical bool

if True, factors marked technical in factor_obs (see set_technical_factors) are removed from the candidate set at every layer, so no cell is ever assigned to an ambient/stress program. Confidence is then computed among the remaining factors only.

False
exclude_batch_technical bool

if True, factors marked batch_technical (see set_batch_technical_factors) are removed from the candidate set at every layer below batch_technical_top_layer, and any cell whose layer-0 winner is one of them is rolled up to that layer. Per-batch splits only exist below the roll-up layer, so the flag is deliberately not applied at or above it.

False
batch_technical_top_layer Optional[int]

layer index the batch-technical roll-up targets. Defaults to adata.uns['batch_technical_top_layer'], recorded by decompose_batch_effects, else 1.

None
Example

After flagging per-batch splits, their cells attach to the

batch-corrected L1 parent instead of overshooting to the root.

scdef.tl.set_batch_technical_factors(model, splits) scdef.tl.assign_confident(model, exclude_batch_technical=True) model.adata.obs["confident_best_layer"].value_counts()

batch_structure_report(model, batch_key=None, group_layer=1, min_group_cells=80, random_seed=0, reference=None)

Describe how batch structure appears among the layer-0 factors.

Summarises the shape of that structure — how batch-skewed each L0 factor is, whether it has an opposite-batch sibling under the same parent, whether it is confined to one branch or overlaid across several, and how separable the batches are inside each branch — so that the analyst can decide what to filter, correct (factor_batch_correction) or keep.

Any fitted model with two layers works, provided batch_key is a column of adata.obs. The model need not have been fitted with that key, and need not have come from decompose_batch_effects. What the report means does depend on which it is:

  • a plain fit with no batch_key — the batch was never corrected, so this is the raw batch structure as the factorization happened to capture it. The natural first look, before deciding whether to use a batch key at all.
  • a decomposed model — the upper layers are batch-corrected and frozen while L0 was re-learned without per-batch gene scales, so batch structure is pushed into L0 and read against a corrected hierarchy. This is the configuration the shape buckets were designed around.
  • a fit made with batch_key — what is left is the residual structure the per-batch gene_scale did not absorb, which is a smaller and different quantity. Useful for auditing a correction, but do not read its magnitudes as the batch effect in the data.

It deliberately emits no verdict. Nothing here distinguishes a technical per-batch duplication of one cell type from a genuine condition-specific biological program: in paired reference data the two are indistinguishable on every column below (same eff_parents ~1, same opp_batch_sibling, same near-ceiling branch_auc). Only the experimental design can settle that, which is why shape is a geometric label and never a cause.

This function supersedes the removed verdict-style "which factors look technical" suggester, whose flag columns did not hold up: its split rule fired on the IFN data's biological monocyte factors, and the batch_split_corr it ranked on saturates at 0.79-0.99 across factor pairs sharing no genes, placing the genuine pbmcs2b Cytotoxic-T split below unrelated pairs. The same geometry is reported here, descriptively.

No cell-type or other biological annotation is read. The inputs are model.pmeans, model.factor_lists / factor_names, adata.obs[batch_key], and — for the two optional gene-side columns — the reference fit's per-batch gene_scale.

Cell side and gene side. Every column below except the last two is cell-side: it asks which cells score on a factor and how they are spread over batches. gene_scale_affinity is gene-side: it asks whether the factor's gene programme is the one the reference fit's per-batch gene_scale absorbed, and never looks at a cell. A factor can be perfectly mixed across batches and still be built from exactly those genes, or the reverse, so the two are worth reading together.

.. warning::

gene_scale_affinity is not a technical-vs-biological score and must not be sorted on as if it were. The per-batch gene_scale absorbs whatever differs between batches at the gene level, and what that is depends entirely on what the batch key encodes. On Kang CTRL/STIM the top-scoring factor is the interferon-response factor — the biology the experiment is about, which must be kept. On pbmcs2b (two runs of one donor) the top-scoring factor is a stress/ambient programme — an artefact to remove. Identical statistic, opposite verdicts; only the experimental design settles it. This is the same trap that the removed "which factors look technical" suggester fell into.

Notation: kept0 = model.factor_lists[0], Z = pmeans['L0z'][:, kept0], W the connection weights from L0 up to the kept group_layer factors, a0(c) = argmax_k Z[c, k] the hard L0 assignment, S_k the cells assigned to factor k, and Znorm the row-normalized Z (each cell's loading distributed over the kept L0 factors).

Columns, one row per kept L0 factor:

  • n_cells: |S_k|. Sums to n_obs over the frame. This is the only column that counts cells with a missing batch_key value; every batch-derived column below is computed on the labelled cells alone.
  • dom_batch: modal batch of S_k (empty string if the factor has no cells with a batch label).
  • frac_dom_batch: fraction of the batch-labelled cells of S_k that lie in dom_batch. Hard purity, with a floor at the largest batch prior and a ceiling of 1.
  • batch_purity_soft: max_b mass_b / sum_b mass_b where mass_b = sum_{c in batch b} Znorm[c, k] over all cells. Floored near 1/n_batches and strongly compressed relative to frac_dom_batch, because normalized loading has a wide background across cells that the factor does not explain.
  • loading_ratio: median Znorm[:, k] in dom_batch over the median in the other batches, taken over all cells. See the note on estimator sensitivity below — read it as a rough direction check, not as a ranking.
  • eff_parents: exp(H(p)) for p the factor's W column normalized to a distribution, i.e. the effective number of group_layer parents it loads on. Range [1, n_kept_group]. Values below ~1.1 are indistinguishable from the numerical noise in tiny W entries; read the column as a coarse "one parent" vs "several parents" flag.
  • parent: current name of the argmax group_layer parent.
  • opp_batch_sibling: True when another kept L0 factor shares this parent and has a different dom_batch — the geometry of a per-batch split of one branch, whatever its cause.
  • branch_auc: cross-validated AUC (3-fold stratified logistic regression on log1p(Z), standardized in-fold) for predicting batch from the L0 scores of the cells in parent's branch. A branch-level quantity: it is identical for all children of a parent and cannot separate siblings. NaN when the branch has fewer than min_group_cells labelled cells, fewer than two batches, or a batch with fewer than 15 cells.
  • shape: descriptive bucket, geometry only, never cause:

================== ========================================================== branch_split frac_dom_batch >= 0.7, eff_parents < 1.5, has an opposite-batch sibling. One branch appearing as two batch-skewed halves. branch_skewed same but with no opposite-batch sibling. A batch-skewed factor whose branch has no counterpart half. overlaid frac_dom_batch >= 0.7 and eff_parents >= 1.5. A batch-skewed program spread over several branches rather than confined to one. balanced everything else, i.e. frac_dom_batch < 0.7. ================== ==========================================================

  • gene_scale_affinity_<batch>, one per batch, plus gene_scale_affinity_max and gene_scale_affinity_batch: present only when the reference fit's per-batch gene_scale is available (see reference). Each per-batch column is the Spearman correlation, over all genes, between the factor's gene loadings and the per-gene log-ratio of that batch's gene_scale against the other batches. _max is the largest of them and _batch names which batch attains it — the top_score and top_batch of get_factor_batch_gene_scale_affinity. With exactly two batches the two contrasts are mirror images, so the per-batch columns are exact negatives of each other, _max is non-negative and _batch only says which side; the magnitude is the information. With many batches the frame gets correspondingly wide — attrs['gene_scale_affinity'] holds the same numbers. Read the warning above before using any of them.

Rows are sorted by frac_dom_batch descending.

result.attrs carries batch_key, group_layer, min_group_cells, random_seed, the shape cut points (frac_dom_min, eff_parents_max) and branch_summary: a DataFrame indexed by kept group_layer factor with n_cells (cells whose hard assignment at group_layer is this factor), n_l0_children (kept L0 factors whose argmax parent is this factor), batch_auc and max_child_frac_dom. When the gene-side columns are present it also carries gene_scale_affinity: the full factors-by-batches frame, with one column per batch rather than only the best one.

Parameters:

Name Type Description Default
model scDEF

any fitted scDEF model with at least two layers; see above for how the reading changes with how it was fitted.

required
batch_key Optional[str]

key in adata.obs holding the batch labels, needing at least two observed values. Defaults to the model's own batch_key when it has one — from being constructed with it, or from decompose_batch_effects, which carries it — and raises if there is neither. The model itself need not have been fitted with the key you pass: any obs column works.

None
group_layer int

layer whose factors define the branches (default 1, the layer usually frozen by the decomposition). Above 1 the per-layer W slices are chained, as in factor_batch_correction.

1
min_group_cells int

minimum labelled cells in a branch before its batch_auc is estimated; smaller branches report NaN.

80
random_seed int

seed for the cross-validation splits.

0
reference Union[scDEF, DataFrame, ndarray, None]

where the per-batch gene_scale contrast for the two gene-side columns comes from — a fitted reference scDEF model, a genes-by-batches DataFrame of log-ratios, or a (n_batches, n_genes) array, as accepted by get_factor_batch_gene_scale_affinity. None (default) uses the profile decompose_batch_effects stored on this model. If neither is available the two columns are simply omitted and the rest of the report is unaffected — a model decomposed before that record existed needs the reference passed explicitly.

None

Returns:

Type Description
DataFrame

DataFrame indexed by current L0 factor name with columns n_cells,

DataFrame

dom_batch, frac_dom_batch, batch_purity_soft,

DataFrame

loading_ratio, eff_parents, parent, opp_batch_sibling,

DataFrame

branch_auc and shape, sorted by frac_dom_batch descending.

Raises:

Type Description
KeyError

batch_key is not in adata.obs.

ValueError

no batch_key was given and the model carries none; batch_key has fewer than two observed values; or group_layer is out of range.

Note

Three caveats worth carrying into any reading of the frame.

loading_ratio is estimator-sensitive. Taking the median over all cells measures the factor's background loading level rather than its loading where it is active, and it can point at the opposite batch from frac_dom_batch. Restricting the median to S_k compresses every factor into a narrow band; using means tracks frac_dom_batch much more closely. The all-cell median is reported for continuity, but the column should not be used to rank factors.

parent and branch_auc use different definitions of the hierarchy: parent is the argmax over the W column, while the branch is defined by the argmax over the group_layer cell scores. These disagree for some factors, and where they do, branch_auc is measured on a cell population that largely excludes the factor's own cells.

Small n_cells rows (a few tens of cells) have a binomial standard error on frac_dom_batch of 0.05-0.10 and should not be compared with rows of several hundred cells.

Example

ref = scdef.scDEF(adata, counts_layer="counts", batch_key="stim") ref.fit() model = scdef.scDEF.decompose_batch_effects(ref, top_layer=1) report = scdef.tl.batch_structure_report(model, batch_key="stim") report.loc[report["shape"] == "branch_split", ["parent", "dom_batch"]] report.attrs["branch_summary"].head()

build_differentiation_paths(model, rel_parent_weight=0.25, abs_parent_weight=0.0, max_paths_per_leaf=8, key_added='differentiation_paths')

Build hierarchy-consistent differentiation paths (top->leaf).

Paths are generated by walking from each layer-0 factor upward through adjacent-layer parent weights, keeping all parent links that pass adaptive filtering:

w >= max(abs_parent_weight, rel_parent_weight * max_parent_weight_for_child).

build_transition_paths(model, rel_parent_weight=0.2, abs_parent_weight=0.0, max_path_len=5, max_paths_per_pair=5, terminal_layer_idx=0, sources=None, targets=None, key_added='transition_paths')

Build transition paths on a soft inter-layer factor graph.

De novo mode: sources=None and targets=None -> all terminal pairs. Targeted mode: provide both sources and targets factor-name lists.

compute_hierarchy_scores(model, use_filtered=False, filter_upper_layers=True, factor_weight='uniform', eps=1e-12)

Compute per-factor and global hierarchy scores from learned W matrices.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
use_filtered bool

whether to use model.factor_lists / model.factor_names

False
filter_upper_layers bool

when use_filtered is False, whether to still use filtered factors for layers > 0 (both as parents and children)

True
factor_weight str

weighting scheme for factors, either "uniform" or "usage"

'uniform'
eps float

small epsilon value for numerical stability

1e-12

Returns:

Type Description
Dict[str, Any]

dict containing per_factor (index = child factor name, child_factor

Dict[str, Any]

column retained; includes avg_clarity and avg_n_eff_parents for L0),

Dict[str, Any]

per_transition, global_score, and global_ambiguity.

compute_within_group_pairwise_dissimilarity(model, layer, obs_key, metric='jsd', eps=1e-12)

Compute within-group pairwise cell dissimilarity for one layer.

Cells are represented by normalized factor memberships from model.adata.obsm[f"X_{layer_name}"]. Pairwise distances are computed within each category of obs_key and summarized per group.

Results are cached in model.adata.uns['within_group_pairwise_dissimilarity'].

drop_technical(model)

Remove factors marked technical from factor_lists and re-annotate.

Technical flags in factor_obs are cleared for the remaining factors.

Parameters:

Name Type Description Default
model scDEF

fitted scDEF model instance

required

factor_batch_correction(model, reduce='sum', key_added='X_L0_batch_corrected', labels_key_added='batch_corrected', top_layer=None)

Apply the batch-technical correction to the scores and to the labels.

Factors flagged by set_batch_technical_factors are per-batch views of a program that the batch-corrected top_layer already represents once. This removes them from the layer-0 representation and writes the same correction as cell-level labels, so an embedding, a heatmap and a UMAP colouring all describe one corrected view.

Flagged factors are grouped by their top_layer parent and handled by group size:

  • Two or more flagged siblings under one parent -- the per-batch halves of a single program -- collapse into one merged column, labelled by joining the members with + in layer order (e.g. "L0_7+L0_14").
  • A lone flagged factor under a parent, with no flagged counterpart -- a batch-skewed program with no opposite-batch half -- has nothing to merge with, so its column is dropped. The cells it claimed are then described by the factors they still score on, and their labels roll up to the parent.

Anything not flagged is untouched and keeps its own column, including a batch-restricted factor that is genuine biology (an ISG program, say) and any non-flagged sibling of a flagged one. An embedding built on the result therefore mixes batches where the split was judged technical and keeps them apart where it was not.

Two adata.obs columns are always written, from the same grouping: f"L0_{labels_key_added}" labels merged siblings by their joined name and a lone flagged factor by its parent, while labels_key_added labels every flagged cell by its parent. All other cells keep their layer-0 label in both.

This addresses one shape of batch structure only: the one that appears as an extra column per batch. It does nothing about a factor that both batches use and merely score differently in magnitude -- on pbmcs2b, correcting the branch_split columns moves the median within-branch batch AUC only 0.970 -> 0.930, because what is left is not confined to a clean pair of sibling columns.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance with cell scores annotated (X_L0 in adata.obsm) and factors flagged by set_batch_technical_factors.

required
reduce Literal['sum', 'max']

how to merge a group of two or more flagged siblings. "sum" (default) adds them, which is exact when the halves are disjoint -- the usual case, since each cell is dominated by the half from its own batch. "max" takes the per-cell peak instead. The two differ only for cells carrying real mass on both halves, where "sum" reports the combined program and "max" the stronger half alone. Irrelevant to a lone flagged factor, which is dropped either way.

'sum'
key_added str

adata.obsm key for the corrected matrix. The column labels go to adata.uns[key_added + '_factors'], the per-column member names to ..._members, the flagged factors to ..._batch_technical and the dropped ones to ..._dropped, so a cached matrix can be checked against the current flags.

'X_L0_batch_corrected'
labels_key_added str

base name for the two adata.obs columns.

'batch_corrected'
top_layer Optional[int]

the batch-corrected parent layer flagged factors group and roll up to. Defaults to adata.uns['batch_technical_top_layer'], recorded by decompose_batch_effects as the layer whose z it froze, else 1. Override only to inspect a different grouping.

None

Returns:

Type Description
None

None. Everything is written to model.adata: the corrected

None

(n_cells, n_corrected_factors) score matrix to

None

obsm[key_added] — with fewer columns than X_L0 whenever

None

anything was flagged — its column labels to

None

uns[key_added + '_factors'], and the two label columns to obs.

Raises:

Type Description
KeyError

X_L0 is missing, the model has no parent layer to group by, or the connection weights are unavailable.

ValueError

reduce is not "sum" or "max", top_layer is out of range, the stored scores are stale, or every column was dropped.

Example

model = scdef.scDEF.decompose_batch_effects(ref, top_layer=1) rep = scdef.tl.batch_structure_report(model) flagged = rep.index[rep["shape"].isin(["branch_split", "branch_skewed"])] scdef.tl.set_batch_technical_factors(model, flagged) scdef.tl.factor_batch_correction(model) scdef.pl.umap(model, color=["L0_batch_corrected", "batch_corrected"])

factor_diagnostics(model, recompute=False, batch_key=None, sensible_top_n_eff_parents_max=1.5, sensible_top_min_best_parent_prob=None, sensible_top_min_clear_children=2, sensible_top_ignore_root=True, sensible_top_use_filtered=True, confidence_threshold=0.9, tau_quantile=0.99, min_effect=None, mc_samples=100, random_seed=0, batch_split_min_batch_frac=0.7, gene_scale_reference=None)

Compute/store factor diagnostics in model.adata.uns['factor_obs'].

Populates per-factor hierarchy scores plus ARD, BRD, n_cells (hard argmax posterior z assignments among kept factors — the same rule as annotate_adata / make_graph(..., assignments=True)), and optional batch metrics when batch_key is set.

Also runs set_confident_signatures so plotting helpers (make_graph, pl.factor_diagnostics(color='signature_confidence'), etc.) can use cached signatures without a separate call.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
recompute bool

if True, force recomputation of the cached fixed upper-layer factor subset used for clarity scores, even if the fit revision did not change.

False
batch_key Optional[str]

key in model.adata.obs used to compute per-factor batch metrics. Defaults to the model's own batch_key when it has one — from being constructed with it, or from decompose_batch_effects, which carries it — so a batch-aware model gets its batch metrics without being told twice. A model with no batch key, or one whose key has fewer than two observed values, gets no batch columns. batch_purity uses hard winner cells (argmax variational z). batch_purity_soft uses the batch distribution of per-cell memberships from X_<layer>_probs (or row-normalized X_<layer> / posterior z if probs are missing). Both are 1 - entropy / log(n_batches). Also writes dom_batch / frac_dom_batch (the factor's dominant batch and its share of the factor's cells) and, for layer 0, the batch_split_* columns from _compute_l0_batch_split. All are plottable via scdef.pl.factor_diagnostics.

None
sensible_top_n_eff_parents_max float

threshold used to classify sensible-top factors on the hierarchy walk.

1.5
sensible_top_min_best_parent_prob Optional[float]

optional best-parent probability threshold used alongside n_eff_parents_max when deciding whether a child merge is clear.

None
sensible_top_min_clear_children int

parent escape-hatch count; a parent is accepted if it owns at least this many clear best-parent children, even when weighted ambiguity is high.

2
sensible_top_ignore_root bool

whether to ignore a width-1 final root layer when classifying sensible-top factors.

True
sensible_top_use_filtered bool

whether to use model.factor_lists / model.factor_names for hierarchy transitions.

True
confidence_threshold float 0.9
tau_quantile float 0.99
min_effect Optional[float] None
mc_samples int 100
random_seed int 0
batch_split_min_batch_frac float

center of the batch-skew ramp that weights candidate partners in _compute_l0_batch_split (default 0.7: a candidate at 0.6 contributes nothing, one at 0.8 contributes fully). Only used when batch_key is set.

0.7
gene_scale_reference Optional[Any]

optional gene-side batch diagnostic, off by default. Pass the REFERENCE fit (the model made with batch_key and per-batch gene_scale) — or anything else get_factor_batch_gene_scale_affinity accepts — to add two layer-0 columns: gene_scale_affinity, the Spearman correlation between the factor's gene loadings and the per-batch gene_scale log-ratio of the batch it matches best, and gene_scale_affinity_batch, that batch. High values say the factor is built from the genes the reference fit's per-batch term absorbed; they do not say the factor is technical (on a condition-like batch key that programme is the biology of the experiment). See that function's docstring.

None

filter_factors(model, batch_key=None, diagnostics_kwargs=None, **filter_kwargs)

Filter factors and refresh the diagnostics and signatures in one step.

model.filter_factors renames factors, which invalidates everything keyed by those names — the stored signatures, the hierarchies, and the frozen upper-layer subset used by factor_diagnostics. Calling the two separately leaves the model in that in-between state, where scd.pl.make_graph(show_signatures=True) raises until diagnostics are re-run. This wrapper does both, so the model is immediately usable::

scd.tl.filter(model, batch_key="Experiment", brd_min=1.0)

Equivalent to::

model.filter_factors(brd_min=1.0)
scd.tl.factor_diagnostics(model, batch_key="Experiment")

Parameters:

Name Type Description Default
model scDEF

scDEF model instance.

required
batch_key Optional[str]

passed to factor_diagnostics, which computes the batch metrics (batch_purity, frac_dom_batch, batch_split_corr). Defaults to the model's own batch_key when it has one, so a batch-aware model does not need it restated here. It plays no part in the filtering itself.

None
diagnostics_kwargs Optional[Mapping[str, Any]]

optional extra keyword arguments for factor_diagnostics (e.g. {"mc_samples": 200}).

None
**filter_kwargs Any

forwarded to filter_factors (brd_min, ard_min, n_eff_parents_max, keep, ...).

{}
Example

scdef.tl.filter(model, batch_key="Experiment", brd_min=1.0) scdef.pl.make_graph(model, show_signatures=True) # works right away

find_sensible_top_factors(model, n_eff_parents_max=1.5, min_best_parent_prob=None, min_clear_children=2, ignore_root=True, use_filtered=True, recompute=False, store=True)

Return factors marked as sensible top factors in factor_obs.

The annotation is computed in scdef.tools.factor.factor_diagnostics and materialized on model.adata.uns['factor_obs'] as is_sensible_top_factor.

find_sensible_top_layer(model, n_eff_parents_max=1.5, min_best_parent_prob=None, min_clear_fraction=0.8, min_clear_children=2, ignore_root=True, use_filtered=True, store=True)

Find the coarsest hierarchy layer supported by confident merges.

Each candidate parent is scored by both the W-row-weighted average of its children's n_eff_parents and the number of clear best-parent children it owns. A parent counts as a clear merge target when either the weighted average is at most n_eff_parents_max or it has at least min_clear_children clear best-parent children (children that name it as their best_parent and themselves pass the threshold). The transition is accepted only when the fraction of clear parents is at least min_clear_fraction.

get_batch_specific_genes_from_gene_scale(model, *, eps=1e-12, log_base=2.0, reference='mean_other_batches')

Per-gene log-ratios of batch-specific gene_scale vs a reference profile.

After fitting (and annotate_adata), scDEF stores inferred positive gene scaling factors in model.pmeans["gene_scale"] with shape (n_batches, n_genes) when batch_key was set with at least two batches. Higher scale for a gene in a batch means the model explains more variance / signal for that gene in that batch (relative to the Gamma prior mean encoded in gene_ratio).

For each batch b, this computes::

log_ratio[g, b] = log( scale[b, g] + eps ) - log( ref[g] + eps )

with log at the chosen base (default log2), and ref either the mean of the other batches at gene g (default) or the global mean across all batches at g.

Parameters:

Name Type Description Default
model 'scDEF'

Fitted scDEF model whose pmeans['gene_scale'] has shape (n_batches, n_genes) with n_batches >= 2.

required
eps float

Small constant for numerical stability.

1e-12
log_base float

Logarithm base (2 for log2, np.e for natural log).

2.0
reference Literal['mean_other_batches', 'global_mean']

mean_other_batches compares each batch to the mean of the remaining batches per gene. global_mean compares to the mean across all batches (same reference column for every batch).

'mean_other_batches'

Returns:

Type Description
DataFrame

DataFrame indexed like model.adata.var_names, one column per batch

DataFrame

label in model.batches (length must match gene_scale batch

DataFrame

dimension). Positive entries indicate relatively higher gene_scale in

DataFrame

that batch vs the chosen reference.

Raises:

Type Description
ValueError

If gene_scale has fewer than two batch rows or is missing.

get_batch_technical_factors(model)

Current model names of the factors marked batch_technical.

The counterpart of set_batch_technical_factors, and the batch-side analogue of get_technical_factors. Names are translated to the current model.factor_names entries the same way, so they can be compared directly against the live model, and flagged factors the model no longer keeps are omitted.

The two flags mean different things and are not interchangeable. A technical factor is a candidate for deletion by drop_technical, and the flag propagates up the tree. A batch_technical factor is a layer-0 per-batch view of a program the corrected parent layer already represents: nothing is deleted and nothing propagates — factor_batch_correction merges or drops it in the corrected representation and leaves the model itself untouched.

Returns:

Type Description
List[str]

Current layer-0 factor names flagged batch-technical, or [] if none

List[str]

are flagged or diagnostics have not been computed.

get_biological_signature(model, top_genes=10)

Gene signature of the top-layer factor — the programme every cell shares.

Looks up the cached confident signature of f"{top_layer_name}_0". Unlike get_technical_signature and get_global_signature, which pool several factors into a relevance-weighted consensus, this reads a single factor's list.

Factors flagged technical are dropped before the lookup, which in practice only matters if the top-layer factor is itself flagged — in that case the result is empty.

With the default top_factors=1 the top layer has exactly one factor. If the model was built with a wider top layer, only its first factor is read.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance. Requires factor_diagnostics (for the technical column) and set_confident_signatures.

required
top_genes int

maximum number of genes to return.

10

Returns:

Type Description
List[str]

The top gene names for the top-layer factor, or [] if it has no

List[str]

cached signature or is itself flagged technical.

get_confident_signatures(model, layer_idx=0, confidence_threshold=0.9, tau_quantile=0.99, min_effect=None, max_genes=None, mc_samples=100, random_seed=0, return_confidences=False)

Get confidence-based signatures per factor using posterior mean/variance.

For each factor independently, this computes a per-factor threshold tau = quantile(E[W_k,:], tau_quantile) and keeps genes that satisfy P(W_k,g > tau) >= confidence_threshold under a normal approximation using the posterior mean and variance of W.

For layer_idx > 0, confidences are estimated with Monte Carlo sampling from the variational posterior via model.get_signature_sample.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
layer_idx int

layer index to use

0
confidence_threshold float

minimum posterior confidence to keep a gene

0.9
tau_quantile float

quantile of factor mean loadings used as threshold tau

0.99
min_effect Optional[float]

optional minimum posterior mean loading E[W_k,g]

None
max_genes Optional[int]

optional maximum number of genes to keep per factor

None
mc_samples int

number of Monte Carlo samples used for layer_idx > 0 confidence estimation

100
random_seed int

random seed for Monte Carlo sampling in upper layers

0
return_confidences bool

whether to also return per-gene confidence arrays

False

Genes are ranked by a combined DE-style score that uses both confidence and posterior mean loading: score = E[W_k,g] * -log10(1 - confidence_k,g).

Returns:

Type Description
Union[Dict[str, List[str]], Tuple[Dict[str, List[str]], Dict[str, ndarray]]]

Dictionary mapping factor names to confident gene lists. If

Union[Dict[str, List[str]], Tuple[Dict[str, List[str]], Dict[str, ndarray]]]

return_confidences is True, also returns a dictionary mapping

Union[Dict[str, List[str]], Tuple[Dict[str, List[str]], Dict[str, ndarray]]]

factor names to confidence arrays aligned with each gene list.

get_factor_annotations(model, factor_names)

Look up factor_obs['annotation'] values for factor names.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance.

required
factor_names Sequence[str]

factor names as in model.factor_names[layer].

required

Returns:

Type Description
List[Optional[str]]

List parallel to factor_names with annotation strings or None.

get_factor_batch_gene_scale_affinity(model, reference=None, *, reference_batches=None, top_n_genes=None, eps=1e-12, log_base=2.0, reference_mode='mean_other_batches')

Score each layer-0 factor by its affinity for the reference fit's per-batch gene_scale contrast.

Named to match get_batch_specific_genes_from_gene_scale, which produces the contrast this reads: that function answers which genes the reference fit's per-batch gene_scale term singled out, this one answers which factors of the decomposed fit look like those genes — hence "factor ... gene_scale affinity".

What it measures. A reference fit made with batch_key (and batch_gene_scale=True) carries a per-batch, per-gene multiplier, pmeans['gene_scale'] of shape (n_batches, n_genes), which soaks up the gene-level differences between batches so the hierarchy does not have to. decompose_batch_effects then freezes the batch-corrected upper hierarchy and re-learns layer 0 without per-batch gene scales, so whatever that term was holding has to reappear in the layer-0 factors. This function measures how much of it landed in each factor: for every batch b it takes the per-gene log-ratio log(gene_scale[b, g] / ref[g]) from get_batch_specific_genes_from_gene_scale, and for every kept layer-0 factor k computes the Spearman rank correlation between that vector and the factor's gene loadings pmeans['L0W'][k, :], over all genes.

Rank correlation over all genes is deliberate. The loadings and the log-ratios live on incomparable scales and both have heavy tails, so a Pearson correlation would be decided by a handful of genes; and restricting to a factor's top genes throws away the (informative) fact that a factor down-weights the batch-associated genes. top_n_genes is offered for exploration but changes the estimand — it conditions on the factor's own top genes, so scores are no longer comparable across factors with different loading profiles.

This is gene-side evidence. Every other batch diagnostic in scdef is cell-side — batch_purity / batch_purity_soft and the batch_split_* columns of factor_diagnostics. Those ask which cells score on a factor and how they are distributed over batches. This one never looks at a cell: it asks whether the factor's gene programme is the one the reference fit's per-batch term absorbed. A factor can be perfectly mixed across batches at the cell level and still be built out of exactly those genes, and vice versa.

.. warning::

This does not separate technical from biological factors, and must not be used as if it did. The per-batch gene_scale absorbs whatever differs between batches at the gene level, and what that is depends entirely on what the batch key encodes. Two worked cases, both scored with this function:

  • Kang CTRL/STIM (batch_key='stim') — the reference gene_scale learned the interferon response: a 16-gene ISG panel sits at median rank 48 of 4000 in the STIM log-ratio (chance is 2000), with IFIT1 #8, ISG15 #16 and CXCL10 #18 at a 62.3x per-batch ratio. So here the top-scoring factor is the ISG factor — which is biology, the entire point of the experiment, and must be kept.
  • pbmcs2b (batch_key='Experiment', two runs of one donor) — the same term captured platelet/ambient genes, lncRNA and a dissociation/immediate-early block, so the top-scoring factor is a stress-response factor — an artefact to remove.

Identical statistic, opposite verdicts. Read the score as "how strongly this factor matches the gene programme the reference fit's per-batch scale absorbed" — a pointer to the batch-associated programme. Inspect the flagged factor's genes (model.adata.uns['L0_signatures'], get_confident_signatures) and decide from the experimental design. A high score on a condition-like batch key is expected and is not grounds for dropping the factor.

Parameters:

Name Type Description Default
model 'scDEF'

the decomposed scDEF model (after decompose_batch_effects) whose layer-0 factors are scored. Its pmeans['L0W'] is full width; only the factors in model.factor_lists[0] are scored, and they are labelled with model.factor_names[0].

required
reference Union['scDEF', DataFrame, ndarray, None]

where the per-batch gene_scale contrast comes from. One of

  • a fitted reference scDEF model (batch_key set, gene_scale of shape (n_batches, n_genes)) — the log-ratios are computed from it with get_batch_specific_genes_from_gene_scale;
  • a DataFrame of precomputed log-ratios, genes-by-batches, indexed by var_names — exactly what that function returns;
  • a raw array of gene_scale, shape (n_batches, n_genes), matching pmeans['gene_scale'].

None (default) looks for adata.uns['reference_gene_scale_log_ratios'] or adata.uns['reference_gene_scale'] on the decomposed model, which decompose_batch_effects records from the fit it decomposed, and raises if neither is present. Supplying it explicitly is only needed for a model decomposed before that record existed, or to score against a different reference.

None
reference_batches Optional[Sequence[str]]

batch labels for the raw-array form of reference (ignored otherwise). Defaults to batch_0, batch_1, ...

None
top_n_genes Optional[int]

if given, restrict each factor's correlation to that factor's own top-N genes by loading. Default None uses every gene, which is the intended statistic; see above.

None
eps float

numerical-stability constant passed through to the log-ratios.

1e-12
log_base float

log base for the log-ratios. Irrelevant to the result — a rank correlation is invariant to any monotone rescaling — and kept only so the underlying table matches what you would get from get_batch_specific_genes_from_gene_scale.

2.0
reference_mode Literal['mean_other_batches', 'global_mean']

mean_other_batches (default) or global_mean, passed straight through. With exactly two batches the two columns are mirror images either way, so the two scores of a factor are exact negatives of each other.

'mean_other_batches'

Returns:

Type Description
DataFrame

A DataFrame indexed by model.factor_names[0] (one row per kept

DataFrame

layer-0 factor), with one column per batch label holding the Spearman

DataFrame

correlation of that factor's loadings with that batch's gene-scale

DataFrame

log-ratio, plus

DataFrame
  • top_batch — the batch whose column is largest for that factor;
DataFrame
  • top_score — that correlation, negative only when the factor is anti-correlated with every batch's contrast;
DataFrame
  • abs_top_score — its magnitude, for ranking regardless of sign.
DataFrame

attrs['n_genes_used'], attrs['reference_mode'] and

DataFrame

attrs['top_n_genes'] record the settings. Rows are in

DataFrame

factor_names[0] order, not sorted by score.

Raises:

Type Description
ValueError

the reference and the model disagree on the gene axis (both shapes are named in the message); the reference gene_scale has a single batch row, so no contrast exists; model.pmeans['L0W'] is missing or its gene dimension does not match adata.n_vars; or top_n_genes is not positive.

TypeError

reference is not a model, DataFrame or array.

Example
decomposed = scdef.scDEF.load("ifn_batch_model")
ref = scdef.scDEF.load("ifn_model")          # fitted with batch_key='stim'
aff = scdef.tl.get_factor_batch_gene_scale_affinity(decomposed, ref)
aff.sort_values("abs_top_score", ascending=False).head()

# then LOOK at the genes before concluding anything:
top = aff["abs_top_score"].idxmax()
decomposed.adata.uns["L0_signatures"][top][:10]

get_global_factors(model, layer_idx=0, n_eff_parents_min=1.5, exclude_technical=True)

Return L{layer_idx} factors shared across lineages (high effective parents).

Uses avg_n_eff_parents for layer 0 when available; otherwise local n_eff_parents. Factors with score >= n_eff_parents_min are returned.

Parameters:

Name Type Description Default
model 'scDEF'

fitted scDEF model with factor_obs stored.

required
layer_idx int

child layer to query (default 0 = L0).

0
n_eff_parents_min float

minimum effective-parent score for a global factor.

1.5
exclude_technical bool

drop factors marked technical in factor_obs.

True

Returns:

Type Description
List[str]

List of factor names in the current model view at layer_idx.

get_global_signature(model, top_genes=10, return_scores=False)

Consensus gene signature over global layer-0 factors.

Requires make_global_hierarchy (or make_hierarchies) to have been run so model.adata.uns['global_hierarchy'] exists.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
top_genes int

number of top genes to return

10
return_scores bool

if True, also return consensus scores

False

Returns:

Type Description
Union[List[str], Tuple[List[str], ndarray]]

Gene list, or (genes, scores) when return_scores=True.

get_hierarchy(model, simplified=True, drop_factors=None)

Get a dictionary containing the polytree contained in the scDEF graph.

Parameters:

Name Type Description Default
simplified Optional[bool]

whether to collapse single-child nodes

True
drop_factors Optional[Sequence[str]]

factors to drop from the hierarchy

None

Returns: hierarchy: the dictionary containing the hierarchy

get_lineage_factors(model, top_factor_label, layer_idx=0, n_eff_parents_max=1.5, prob_min=0.5, exclude_technical=True)

Return factors at layer_idx that are clear descendants of top_factor_label.

top_factor_label may name a factor at any layer (including the top). Factors at layer_idx are included when the upward best_parent chain reaches that ancestor and each step on the path is unambiguous (n_eff_parents <= n_eff_parents_max, best_parent_prob >= prob_min).

When layer_idx == ancestor_layer - 1, these are the direct children of the ancestor; for lower layer_idx, they are deeper descendants with the same clarity constraints on every step up to the ancestor.

Parameters:

Name Type Description Default
model 'scDEF'

fitted scDEF / iscDEF / sscDEF model with factor_obs stored.

required
top_factor_label str

ancestor factor name (any layer; current model.factor_names or factor_obs index).

required
layer_idx int

child layer to query (default 0 = L0).

0
n_eff_parents_max float

maximum effective parents on each step toward the ancestor.

1.5
prob_min float

minimum best_parent_prob on each step toward the ancestor.

0.5
exclude_technical bool

drop factors marked technical in factor_obs.

True

Returns:

Type Description
List[str]

List of factor names in the current model view at layer_idx.

get_obs_score_rankings(model, layer, obs_key, obs_values, mode='fracs', ascending=False, recompute=False)

Return per-obs-value factor rankings by observation association score.

This reads cached matrices from model.adata.uns['obs_scores'] (written by scd.pl.obs_scores). If cache is missing/stale for the requested key/model, it is recomputed on demand for the requested obs_key and mode.

get_obs_value_specific_factors(model, layer, obs_key, obs_values, mode='fracs', min_specificity=0.0, top_n=None, recompute=False, return_scores=False)

Get factors specific to each obs value in a layer.

Specificity is defined within the provided obs_values as: specificity = score(obs_value) - max(score(other_obs_values)). Higher values indicate stronger specificity for that obs category.

get_stored_confident_signatures(model, layer_idx=0, max_genes=None, return_confidences=False, return_combined_scores=False, return_signature_confidences=False)

Load precomputed confident signatures (and optional scores) from cache.

get_technical_factors(model)

Current model names of the factors marked technical in factor_obs.

factor_obs rows are keyed by the factor names in place when diagnostics ran; after filter_factors the model renames factors contiguously. The returned names are always the current model.factor_names entries, so they can be compared directly against the live model. Technical factors that are no longer kept by the model are omitted.

get_technical_signature(model, top_genes=10, return_scores=False)

Consensus gene signature over the factors flagged as technical.

Pools the layer-0 gene rankings of the factors in the technical hierarchy into one ranked list, so the variation drop_technical would remove can be inspected before removing it.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance. Requires make_technical_hierarchy (or make_hierarchies) to have run, so adata.uns['technical_hierarchy'] exists.

required
top_genes int

maximum number of genes to return.

10
return_scores bool

also return the consensus score of each returned gene.

False

Returns:

Type Description
Union[List[str], Tuple[List[str], ndarray]]

The top gene names, or (genes, scores) when return_scores is True.

gsea(model, libs=('KEGG_2019_Human',), custom_gene_sets=None, organism='Human', background_genes=None, layers=None, top_genes=None, cutoff=0.05, outdir=None)

Run Enrichr pathway enrichment for cached signatures across layers.

This utility uses signatures from scd.tl.get_stored_confident_signatures and does not rely on model-level ranking by raw W. Online libraries in libs are fetched to local dicts and merged with custom_gene_sets so each factor is tested against one combined gene-set universe using a single gp.enrich call. By default, runs for all layers and stores per-layer results in adata.uns['factor_enrichments'].

make_biological_hierarchy(model)

Make the biological hierarchy of the model.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required

Returns:

Name Type Description
biological_hierarchy Dict[str, Sequence[str]]

dictionary containing the biological hierarchy

make_global_hierarchy(model)

Make the global (shared-across-lineages) hierarchy of the model.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required

Returns:

Name Type Description
global_hierarchy Dict[str, Sequence[str]]

dictionary with synthetic root global_top and global layer-0 factors as direct children.

make_hierarchies(model)

Store the biological, technical, and global hierarchies of the model.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required

make_technical_hierarchy(model)

Make the technical hierarchy of the model.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required

Returns:

Name Type Description
technical_hierarchy Dict[str, Sequence[str]]

dictionary containing the technical hierarchy

multilayer_umap(model, layers=None, weights=None, normalize_per_layer=False, use_log=False, metric='euclidean', key_added='multilayer', eps=1e-08, neighbors_kwargs=None, umap_kwargs=None)

Compute a single UMAP from the concatenation of all scDEF layers.

Each cell's representation is the concatenation of its soft factor assignments X_{layer_name} across the selected layers, producing a multi-resolution signature that encodes identity at every level of the hierarchy simultaneously.

This tends to produce a lineage-aware embedding that smooths through differentiation trajectories: terminally differentiated cells form tight clusters by shared fine-layer factor, progenitor cells bridge siblings through their shared parent-layer mass, and stem-like cells with diffuse scores sit in their own regime. Contrast with umap (per-layer), which only captures a single scale.

By default it uses raw per-layer scores as stored in X_{layer}. Set normalize_per_layer=True to convert each layer's row to proportions before optional log/weighting; this makes geometry reflect relative composition at each layer rather than absolute score scale.

Use weights to bias the embedding toward fine or coarse resolution.

Single-factor layers (K_k == 1, typically the root) are skipped automatically — every cell has mass 1 there, so they add no discriminative signal.

Parameters:

Name Type Description Default
model scDEF

fitted scDEF model.

required
layers Optional[Sequence[int]]

layer indices to include. If None, uses all layers with K_k >= 2, in ascending order.

None
weights Optional[Sequence[float]]

optional per-layer multiplicative weights applied to each layer's sub-vector before concatenation. Must have the same length as layers. Larger weight → that layer has more influence on the embedding geometry. Default: uniform.

None
normalize_per_layer bool

if True, row-normalize each selected layer block so rows sum to 1 before optional log and weighting. Useful when you want distances to depend on relative factor composition rather than total per-layer score magnitude.

False
use_log bool

if True, replace each sub-vector by log(X + eps) before applying weights/concatenation. Helpful when fine resolution is dominated by a single factor and small proportions get crushed by Euclidean distance.

False
metric str

distance metric for sc.pp.neighbors.

'euclidean'
key_added str

suffix used to store results. Writes adata.obsm[f"X_{key_added}"] (the concatenated representation) and adata.obsm[f"X_umap_{key_added}"] (the UMAP embedding).

'multilayer'
eps float

floor for log when use_log=True.

1e-08
neighbors_kwargs Optional[Dict[str, object]]

extra kwargs forwarded to sc.pp.neighbors.

None
umap_kwargs Optional[Dict[str, object]]

extra kwargs forwarded to sc.tl.umap.

None

Returns:

Type Description
ndarray

The concatenated representation of shape (n_cells, sum K_k).

multilevel_paga(model, neighbors_rep='X_L0', layers=None, reuse_pos=True, layout='fa', random_seed=0, **paga_kwargs)

Compute and cache multilevel PAGA results for plotting.

score_paths(model, paths_key='transition_paths', key_added=None, normalize_per_layer=True, min_affinity=0.0, eps=1e-12)

Score per-cell position and affinity on previously built paths.

Writes
  • adata.obsm[f"{prefix}_positions"] (n_cells, n_paths)
  • adata.obsm[f"{prefix}_affinities"] (n_cells, n_paths)
  • per-path summary table in adata.uns[paths_key]["path_stats"].

set_batch_technical_factors(model, factors)

Mark layer-0 factors as batch-technical (per-batch splits of one type).

Unlike set_technical_factors, this flag does not propagate up the tree and no factor is deleted. Batch-technical factors are layer-0 per-batch views of a program that the parent layer already represents once, so that parent is the roll-up target and propagating the flag into it would throw away the very signal the roll-up depends on.

This works on any fitted model with a batch column in adata.obs, not only on one from decompose_batch_effects — a plain fit made without a batch_key leaves batch structure in the factors directly, and can be flagged and corrected the same way. See batch_structure_report for how the reading differs between the two.

This records which factors are batch-technical and nothing else. The roll-up target comes from adata.uns['batch_technical_top_layer'], written by decompose_batch_effects as the layer whose z it froze — the only place that knows it. Consumers (factor_batch_correction, assign_confident) read it by default and take an override argument, so there is no second copy here to fall out of step with the decomposition.

Sets factor_obs['batch_technical'].

Parameters:

Name Type Description Default
model scDEF

scDEF model instance.

required
factors Sequence[str]

factor names to flag. Names are resolved against the current model.factor_names (see _resolve_factor_obs_names), so names taken from a filtered model are safe.

required

Raises:

Type Description
ValueError

if any name in factors cannot be resolved.

Example

import scdef

Stage 1: batch-corrected reference fit.

ref = scdef.scDEF(adata, counts_layer="counts", batch_key="Experiment") ref.fit()

Stage 2: re-learn L0/L1 without the batch key to expose batch programs.

model = scdef.scDEF.decompose_batch_effects(ref, top_layer=1)

Describe the batch geometry, then decide from the design which

branch splits are technical (here both batches are one donor).

rep = scdef.tl.batch_structure_report(model, batch_key="Experiment") splits = rep.index[rep["shape"] == "branch_split"] scdef.tl.set_batch_technical_factors(model, splits) scdef.tl.factor_batch_correction(model, top_layer=1)

Cells owned by a split are assigned at their L1 parent instead.

scdef.tl.assign_confident(model, exclude_batch_technical=True)

set_cell_entropies(model, layers=None, key_suffix='entropy', effective_suffix='effective_n_factors', normalize=True, eps=1e-12)

Compute per-cell assignment entropy and store one column per layer.

For each selected layer, uses model.adata.obsm[f"X_{layer_name}"] to build per-cell membership probabilities and computes Shannon entropy.

If normalize=True, entropy is divided by log(n_factors_layer) so values are approximately in [0, 1] (for layers with >1 factors).

Also stores an effective number of factors per cell, defined as exp(H) where H is the non-normalized Shannon entropy.

Returns:

Type Description
List[str]

List of created/updated entropy column names.

set_confident_signatures(model, confidence_threshold=0.9, tau_quantile=0.99, min_effect=None, mc_samples=100, random_seed=0)

Precompute and cache confident signatures/scores for all layers.

Stores signatures, per-gene confidences, combined scores, and per-factor weighted signature Jaccard confidences (posterior stability of each confident gene list, weighted by combined_scores) in model.adata.uns['confident_signatures'] for reuse by plotting.

set_factor_signatures(model, signatures=None, top_genes=10)

Store a signature per factor in adata.uns['factor_signatures'].

With signatures=None the confident signatures of every layer are pooled into one {factor_name: genes} mapping and stored, truncated to top_genes; pass a mapping instead to store curated lists.

Note that set_confident_signatures already writes this key, so the signatures=None path only re-writes it at a different length. Nothing inside scdef reads uns['factor_signatures'] — the plots take their gene lists from uns['confident_signatures'] instead — so this is for downstream use and for overriding the stored lists by hand.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance. With signatures=None, requires set_confident_signatures to have been run.

required
signatures Optional[Dict[str, List[str]]]

mapping of factor name to gene list. None builds it from the cached confident signatures.

None
top_genes int

genes per factor when building from the cache.

10

Returns:

Type Description
Dict[str, List[str]]

The stored mapping.

set_global_factors(model, factors=None, layer_idx=0, n_eff_parents_min=1.5, exclude_technical=True)

Mark global (shared-across-lineages) factors in factor_obs.

Global factors are identified from hierarchy diagnostics (high effective parents). They are excluded from make_biological_hierarchy. Use drop_technical to remove technical factors from the active model.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
factors Optional[Sequence[str]]

explicit factor names to mark as global (resolved like set_technical_factors). When None, uses get_global_factors.

None
layer_idx int

child layer for automatic selection (default 0).

0
n_eff_parents_min float

minimum effective-parent score when factors is None.

1.5
exclude_technical bool

do not mark technical factors as global.

True

set_technical_factors(model, factors=None, brd_min=1.0, ard_min=0.001, clarity_min=0.5, n_eff_parents_max=1.5, brd_exceptional=None, local_l0_scores=False, batch_purity_max=None, batch_purity_soft_max=None, min_cells_lower=0.0)

Set the technical factors of the model.

Technical factors must be layer 0 factors.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
factors Optional[Sequence[str]]

list of factor names to mark as technical. Names are resolved against the current model.factor_names (and translated to the corresponding factor_obs rows via original_factor_idx), so it is safe to pass names from the model after filter_factors(). When provided, criteria-based selection is skipped.

None
brd_min Optional[float]

minimum BRD threshold for keeping biological layer-0 factors when factors is None.

1.0
ard_min Optional[float]

minimum ARD fraction threshold for keeping biological layer-0 factors when factors is None.

0.001
clarity_min Optional[float]

minimum L0 clarity when not using lineage avg_n_eff_parents.

0.5
n_eff_parents_max float

only for lineage diagnostics: ceiling on avg_n_eff_parents (default 1.5; matches scd.pl.factor_diagnostics). When brd_exceptional is set, factors with BRD >= brd_exceptional are kept regardless.

1.5
brd_exceptional Optional[float]

if set, high-BRD escape hatch when lineage effective parents exceed n_eff_parents_max. Default None (disabled).

None
local_l0_scores bool

if True, biological factors are chosen using n_eff_parents and n_eff_parents_max instead of lineage averages / clarity_min.

False
batch_purity_max Optional[float]

if set, layer-0 factors with hard batch_purity above this value are not biological (requires factor_diagnostics(..., batch_key=...)).

None
batch_purity_soft_max Optional[float]

if set, same for soft batch_purity_soft. Same semantics as filter_factors / factor_diagnostics plot.

None
min_cells_lower Optional[float]

minimum cell-count criterion for keeping biological layer-0 factors when factors is None. Same semantics as scDEF.filter_factors(..., min_cells_lower=...).

0.0
Notes

When factors is None, the candidate pool is restricted to the layer-0 factors currently kept in model.factor_lists[0]. Already filtered-out factors are never re-introduced as technical.

umap(model, layers=None, use_log=False, metric='euclidean')

Compute UMAP embeddings for each scDEF layer.

The resulting embeddings are stored in model.adata.obsm[f"X_umap_{layer_name}"] for each layer. Any pre-existing adata.obsm['X_umap'], adata.uns['neighbors'], adata.uns['umap'], and adata.obsp neighbor graphs are restored afterward (removed if absent), so generic scanpy calls keep using the original embedding and graph.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
layers Optional[List[int]]

which layers to compute UMAPs for, in processing order. If None, all layers with more than one factor are used, coarse-to-fine (descending layer index).

None
use_log bool

whether to use log-transformed cell-factor weights for the neighbor graph computation.

False
metric str

distance metric for neighbors computation.

'euclidean'

When a corrected layer-0 representation is already present in adata.obsm['X_L0_batch_corrected'] (written by factor_batch_correction), it is embedded as well and stored as adata.obsm['X_umap_L0_corrected'], in addition to the per-layer embeddings. This function never builds that representation.

Plotting

scdef.pl

Plotting utilities for scDEF.

biological_hierarchy(model, **kwargs)

Plot the biological hierarchy of the model.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
**kwargs Any

keyword arguments passed to make_graph

{}

Returns:

Type Description
Graph

Graphviz Graph object

cell_entropies(model, thres=0.9, entropy_suffix='entropy', effective_suffix='effective_n_factors', show=True)

Plot cell entropies and factor numbers across layers.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
thres float

Threshold for cumulative sum calculation

0.9
entropy_suffix str

suffix used by per-layer entropy columns in adata.obs

'entropy'
effective_suffix str

suffix used by per-layer effective-factor columns

'effective_n_factors'
show bool

Whether to show the plot

True

continuous_obs_scores(model, obs_keys, mode='correlations', vmax=None, vmin=None, **kwargs)

Plot the correlations between a set of cell annotations and factors.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
obs_keys Sequence[str]

the keys in model.adata.obs to use

required
mode Literal['correlations']

how to compute scores

'correlations'
**kwargs Any

plotting keyword arguments

{}

factor_diagnostics(model, brd_min=1.0, ard_min=0.001, clarity_min=0.5, batch_purity_max=None, batch_purity_soft_max=None, n_eff_parents_max=1.5, brd_exceptional=None, figsize=(6, 4), ax=None, annotate_factors=False, annotation_fontsize=8, annotation_alpha=0.8, all_factors=False, local_l0_scores=False, x='BRD', y=None, color=None, size='ARD', batch_panel=None, show=True, _draw_colorbar=True, _draw_size_legend=True)

Diagnostic scatter plot of layer-0 factors with flexible axis/color/size mapping.

By default plots BRD vs effective parents with marker size scaled by ARD. When batch_purity_max or batch_purity_soft_max is set and color is not overridden, points are colored by the corresponding batch purity (sizes still default to ARD).

Any per-factor column of factor_obs can be mapped to an axis, color or size. That includes the batch-split diagnostics written by scdef.tools.factor_diagnostics(..., batch_key=...), so a within-cell-type per-batch split can be inspected directly, e.g.::

scdef.pl.factor_diagnostics(
    model, x="avg_n_eff_parents", y="batch_split_corr", color="frac_dom_batch"
)

When those diagnostics exist and the axes are left at their defaults, a second panel showing them is drawn automatically (see batch_panel). Both batch_split_corr and frac_dom_batch are dense — every kept factor has a value — so nothing is silently dropped from either panel. A high batch_split_corr at a balanced frac_dom_batch (top-left of that panel) is the batch-corrected cell-type factor, not a split half.

This plot renders diagnostics; it issues no verdict, and neither does any tool it draws from. To describe the shape of the batch structure — which factors are batch-skewed, which have an opposite-batch sibling under the same parent, and how separable the batches are inside each branch — use batch_structure_report::

rep = scdef.tl.batch_structure_report(model, batch_key="Experiment")
splits = rep.index[rep["shape"] == "branch_split"]

shape is geometry, never cause: only the experimental design can say whether a per-batch split is a technical duplication or a genuine condition-specific program.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
brd_min float

minimum BRD filter threshold

1.0
ard_min float

minimum ARD filter threshold (fraction of total ARD)

0.001
clarity_min float

used for the horizontal cutoff when filtering on local n_eff_parents (not lineage avg_n_eff_parents): cutoff is effective_parents_from_clarity(clarity_min, K_parents).

0.5
batch_purity_max Optional[float]

optional upper bound on hard-assignment batch purity. Factors pass when batch_purity <= batch_purity_max. If provided and color is not set, the scatter is colored by batch_purity with the colorbar threshold at this value.

None
batch_purity_soft_max Optional[float]

optional upper bound on soft batch purity (from X_<layer>_probs). Factors pass when batch_purity_soft <= batch_purity_soft_max. If provided and color is not set (and batch_purity_max is None), color defaults to batch_purity_soft. Requires scdef.tools.factor_diagnostics(..., batch_key=...).

None
n_eff_parents_max float

used when filtering on lineage avg_n_eff_parents (local_l0_scores=False and column present): dashed line at this value and pass rule y <= n_eff_parents_max (default 1.5). When brd_exceptional is set, factors with BRD >= brd_exceptional also pass.

1.5
brd_exceptional Optional[float]

if set, dashed vertical line on BRD and pass rule allowing high-BRD factors to be kept even when effective parents exceed n_eff_parents_max. Default None (disabled).

None
figsize tuple

Figure size (if ax is None)

(6, 4)
ax Optional[Axes]

matplotlib Axes to plot on

None
annotate_factors bool

whether to annotate each point with its factor label

False
annotation_fontsize int

fontsize for factor text annotations

8
annotation_alpha float

alpha value for factor text annotations

0.8
all_factors bool

if True, plot diagnostics for all layer-0 factors from the complete snapshot model.adata.uns['factor_obs_full'] (including factors that were filtered out). Default (False) plots the current view model.adata.uns['factor_obs'], which after model.filter_factors() contains only kept factors.

False
local_l0_scores bool

when y is not set, use layer-0 n_eff_parents on the y-axis instead of lineage avg_n_eff_parents.

False
x FactorDiagQuantity

quantity for the x-axis (ARD, BRD, n_eff_parents, avg_n_eff_parents, batch_purity, batch_purity_soft, batch_split_corr, frac_dom_batch, signature_confidence, n_cells). batch_split_corr and frac_dom_batch require scdef.tools.factor_diagnostics(..., batch_key=...).

'BRD'
y Optional[FactorDiagQuantity]

quantity for the y-axis; default follows local_l0_scores / avg_n_eff_parents availability.

None
color Optional[FactorDiagQuantity]

quantity for marker color. Default None: use batch_purity when batch_purity_max is set, else batch_purity_soft when batch_purity_soft_max is set, else uncolored markers. In the two-panel batch layout the default is batch_purity for both panels, so the same factor reads the same way in each; pass color explicitly to override the left panel.

None
size Optional[FactorDiagQuantity]

quantity for marker size. Default ARD. Pass None for fixed marker size.

'ARD'
batch_panel Optional[bool]

draw a second panel with the per-batch split diagnostics (x='frac_dom_batch', y='batch_split_corr', size='n_cells', color='batch_purity') beside the usual one. None (default) enables it automatically when the batch diagnostics exist and the axes were left at their defaults; passing x/y or an ax keeps the single panel. True forces both panels (requires ax=None and the batch diagnostics); False forces one. In the two-panel layout the left panel keeps whatever x/y/color/size you pass.

None
show bool

whether to show the plot

True

Returns:

Type Description
Optional[Union[Axes, ndarray]]

Axes object if show is False, None otherwise. With two panels, an array

Optional[Union[Axes, ndarray]]

of the two Axes.

factor_gene_uncertainty_boxplot(model, factor, layer_idx=0, max_genes=50, mc_samples=100, random_seed=0, whisker_quantiles=(0.05, 0.95), sort_by='mean', color_by_confidence=False, confidence_tau_quantile=0.99, confidence_cmap='viridis', confidence_vmin=0.0, confidence_vmax=1.0, add_confidence_colorbar=True, show_confidence_cutoff_line=False, confidence_include_threshold=0.9, confidence_cutoff_line_kwargs=None, show_tau_quantile_line=False, tau_quantile_line_kwargs=None, xtick_rotation=90.0, figsize=(12, 4), ax=None, show=True)

Plot per-gene uncertainty boxes for a factor using posterior mean/variance.

Genes are sorted by posterior mean loading or confidence for the selected factor. For each gene, a box is drawn from an approximated posterior distribution of W using the stored posterior mean and variance: - box: 25th to 75th percentile - median: 50th percentile - whiskers: configurable quantiles (default 5th/95th)

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
factor Union[int, str]

factor index (within kept factors of the layer) or factor name

required
layer_idx int

layer index to visualize

0
max_genes Optional[int]

maximum number of genes to plot; if None, plot all genes

50
mc_samples int

number of posterior samples for layer_idx > 0

100
random_seed int

random seed for upper-layer posterior sampling

0
whisker_quantiles Tuple[float, float]

lower/upper quantiles for whiskers

(0.05, 0.95)
sort_by Literal['mean', 'confidence']

sorting criterion for genes; "mean" (default) sorts by posterior mean loading, "confidence" sorts by posterior confidence P(score > tau)

'mean'
color_by_confidence bool

whether to color each box by posterior confidence P(W > tau)

False
confidence_tau_quantile float

factor-wise quantile used to define tau for confidence coloring

0.99
confidence_cmap str

matplotlib colormap name for confidence coloring

'viridis'
confidence_vmin float

lower bound of confidence colormap normalization

0.0
confidence_vmax float

upper bound of confidence colormap normalization

1.0
add_confidence_colorbar bool

whether to add a confidence colorbar when color_by_confidence is True

True
show_confidence_cutoff_line bool

whether to draw a vertical line marking the last plotted gene with confidence >= confidence_include_threshold

False
confidence_include_threshold float

confidence threshold used for the cutoff line (default 0.9)

0.9
confidence_cutoff_line_kwargs Optional[Dict[str, Any]]

optional kwargs passed to ax.axvline for styling the cutoff line

None
show_tau_quantile_line bool

whether to draw a horizontal line at tau, where tau is the quantile threshold set by confidence_tau_quantile

False
tau_quantile_line_kwargs Optional[Dict[str, Any]]

optional kwargs passed to ax.axhline for styling the tau quantile line

None
xtick_rotation float

x tick label rotation in degrees

90.0
figsize Tuple[float, float]

figure size if ax is None

(12, 4)
ax Optional[Axes]

matplotlib axis to draw on

None
show bool

whether to show the plot

True

Returns:

Type Description
Optional[Axes]

Axes object if show is False, None otherwise.

factor_genes(model, thres=0.9, show=True)

Plot number of genes in factors across layers.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
thres float

threshold for cumulative sum calculation

0.9
show bool

whether to show the plot

True

Returns:

Type Description
Optional[Figure]

Figure object if show is False, None otherwise

factor_gini(model, idx, thres=0.9, show=True)

Plot Gini coefficient for a specific factor.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
idx int

Factor index to plot

required
thres float

Threshold for cumulative sum calculation

0.9
show bool

Whether to show the plot

True

factors_bars(model, obs_keys, sort_layer_factors=True, orders=None, sharey=True, layers=None, vmax=None, vmin=None, fontsize=12, title_fontsize=12, legend_fontsize=8, figsize=(10, 4), total=False, show=True)

Plot factor scores as bar charts.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
obs_keys Union[str, List[str]]

observation keys to plot

required
sort_layer_factors bool

whether to sort factors by layer

True
orders Optional[List[ndarray]]

custom factor orders

None
sharey bool

whether to share y-axis across subplots

True
layers Optional[List[int]]

which layers to plot

None
vmax Optional[float]

maximum value for y-axis

None
vmin Optional[float]

minimum value for y-axis

None
fontsize int

font size for labels

12
title_fontsize int

title font size

12
legend_fontsize int

legend font size

8
figsize Tuple[float, float]

figure size

(10, 4)
total bool

whether to plot total scores

False
show bool

whether to show the plot

True

gini_brd(model, normalize=False, figsize=(4, 4), alpha=0.6, fontsize=12, legend_fontsize=10, show=True, ax=None)

Plot Gini coefficient vs BRD scores.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
normalize bool

whether to normalize BRD scores

False
figsize Tuple[float, float]

figure size

(4, 4)
alpha float

transparency level

0.6
fontsize int

font size for labels

12
legend_fontsize int

font size for legend

10
show bool

whether to show the plot

True
ax Optional[Axes]

matplotlib axes to plot on

None

Returns:

Type Description
Optional[Axes]

Axes object if show is False, None otherwise

global_hierarchy(model, show_signatures=True, **kwargs)

Plot the global (shared-across-lineages) factor hierarchy.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
show_signatures bool

whether to show gene signatures

True
**kwargs Any

keyword arguments passed to make_technical_hierarchy_graph

{}

Returns:

Type Description
Graph

Graphviz Graph object

layers_obs(model, obs_keys, obs_mats, obs_clusters, obs_vals_dict, sort_layer_factors=True, orders=None, layers=None, vmax=None, vmin=None, cb_title='', cb_title_fontsize=10, fontsize=12, title_fontsize=12, pad=0.1, shrink=0.7, figsize=(10, 4), xticks_rotation=90.0, cmap=None, show=True, rasterized=False, **kwargs)

Plot observation matrices across layers.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
obs_keys Union[str, List[str]]

observation keys to plot

required
obs_mats Dict[str, Dict[int, ndarray]]

observation matrices dictionary

required
sort_layer_factors bool

whether to sort factors by layer

True
orders Optional[List[ndarray]]

custom factor orders

None
layers Optional[List[int]]

which layers to plot

None
vmax Optional[float]

maximum value for colormap

None
vmin Optional[float]

minimum value for colormap

None
cb_title str

colorbar title

''
cb_title_fontsize int

colorbar title font size

10
fontsize int

font size for labels

12
title_fontsize int

title font size

12
pad float

padding for colorbar

0.1
shrink float

shrink factor for colorbar

0.7
figsize Tuple[float, float]

figure size

(10, 4)
xticks_rotation float

rotation angle for x-axis ticks

90.0
cmap Optional[str]

colormap name

None
show bool

whether to show the plot

True
rasterized bool

whether to rasterize the plot

False
**kwargs Any

additional plotting keyword arguments

{}

loss(model, figsize=(4, 4), fontsize=12, ax=None, show=True)

Plot training loss over epochs.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
figsize Tuple[float, float]

figure size

(4, 4)
fontsize int

font size for labels

12
ax Optional[Axes]

matplotlib axes to plot on

None
show bool

whether to show the plot

True

Returns:

Type Description
Optional[Axes]

Axes object if show is False, None otherwise

make_graph(model, hierarchy=None, show_all=False, factor_annotations=None, top_factor=None, show_signatures=True, drop_factors=None, root_signature=None, root_ranking=None, enrichments=None, show_enrichments=False, top_genes=None, show_batch_counts=False, filled=None, wedged=None, assignments=True, color_edges=True, show_confidences=False, n_cells_label=False, n_cells=False, node_size_max=2.0, node_size_min=0.05, scale_level=False, show_label=True, gene_score=None, gene_cmap='viridis', path=None, path_color='red', path_node_penwidth=2.5, confident_assignments=False, confident_key='confident', shell=False, r=2.0, r_decay=0.8, root_shape=None, bottom_layer=0, **fontsize_kwargs)

Make Graphviz-formatted scDEF graph.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
hierarchy Optional[Dict[str, Sequence[str]]]

dictionary containing the polytree to draw instead of the whole graph

None
show_all Optional[bool]

whether to show all factors even post filtering

False
factor_annotations Optional[Dict[str, str]]

factor annotations to include in the node labels

None
top_factor Optional[str]

only include factors below this factor

None
show_signatures Optional[bool]

whether to show the ranked gene signatures in the node labels

True
drop_factors Optional[List[str]]

list of factors to drop from the graph

None
root_signature Optional[List[str]]

root signature to display

None
root_ranking Optional[List[str]]

root ranking to display

None
enrichments Optional[DataFrame]

enrichment results dataframe to include in node labels

None
show_enrichments Optional[bool]

whether to show enrichment terms in node labels

False
top_genes Optional[Union[int, List[int]]]

number of genes from each signature to be shown in the node labels

None
show_batch_counts Optional[bool]

whether to show the number of cells from each batch that attach to each factor

False
filled Optional[Union[str, Dict[str, float]]]

key from model.adata.obs to use to fill the nodes with, or dictionary of factor scores

None
wedged Optional[str]

key from model.adata.obs to use to wedge the nodes with

None
assignments Optional[bool]

whether to use the assignments of cells to factors to wedge the nodes, rather than the scores

True
color_edges Optional[bool]

whether to color the graph edges according to the upper factors

True
show_confidences Optional[bool]

whether to show precomputed signature Jaccard confidence (from scd.tl.set_confident_signatures) for each factor

False
n_cells_label Optional[bool]

whether to show the number of cells that attach to the factor

False
n_cells Optional[bool]

whether to scale the node sizes by the number of cells that attach to the factor

False
node_size_max Optional[float]

maximum node size when scaled by cell numbers

2.0
node_size_min Optional[float]

minimum node size when scaled by cell numbers

0.05
scale_level Optional[bool]

whether to scale node sizes per level instead of across all levels

False
show_label Optional[bool]

whether to show labels on nodes

True
gene_score Optional[str]

color the nodes by the score they attribute to a gene, normalized by layer. Overrides filled and wedged

None
gene_cmap Optional[str]

colormap to use for gene_score

'viridis'
path Optional[Union[Sequence[str], Dict[str, Any]]]

ordered factor names along a path, or a dict with key "nodes" (e.g. differentiation or transition path nodes). Names outside the hierarchy view are dropped when hierarchy is set. Edge highlights map each consecutive pair to the drawn weight edge (coarser layer → adjacent finer layer), so coarse→fine chains match, and transition zigzags (e.g. L0→L1→L0→…) work when each hop is between adjacent layers. If hierarchy is set, a pair is highlighted when it is a tree edge in either direction; otherwise the same adjacent-layer rule applies. Sets Graphviz color (node border and edge stroke) to path_color; fillcolor is unchanged. With gene_score and color_edges, path edges use path_color instead of the parent gene-based edge color.

None
path_color str

stroke color for highlighted path nodes (border) and edges

'red'
path_node_penwidth float

Graphviz penwidth for nodes on path (border thickness). Ignored when path is not set. Default 2.5; use 1.0 for the usual thin border (no emphasis).

2.5
confident_assignments bool

if True, attach cells using adata.obs[f"{confident_key}_factor"] (cross-layer assignment from assign_confident) instead of per-layer adata.obs[layer_name]. Requires that column to exist.

False
confident_key str

key_added prefix used with assign_confident (default "confident").

'confident'
shell Optional[bool]

whether to use shell layout. If the model has a final width-1 root layer, shell layout omits that root and plots the remaining hierarchy.

False
r Optional[float]

radius parameter for shell layout

2.0
r_decay Optional[float]

radius decay parameter for shell layout

0.8
root_shape Optional[str]

Graphviz node shape for root-layer factors (e.g. "diamond", "box", "hexagon"). When None the default ellipse is used.

None
bottom_layer int

lowest layer index to include (default 0). Set to 1 to omit the finest layer (L0) when it has too many factors.

0
**fontsize_kwargs Any

keyword arguments to adjust the fontsizes according to the gene scores

{}

Returns:

Type Description
Graph

Graphviz Graph object

multilevel_paga(model, neighbors_rep='X_L0', layers=None, figsize=(16, 4), reuse_pos=True, recompute=False, fontsize=12, show=True, **paga_kwargs)

Plot cached multilevel PAGA graphs across scDEF layers.

obs_cell_factor_heatmap(model, subset_obs_key, subset_obs, group_obs_key, layer_idx=0, factors=None, values='score', cluster_cells=True, sort_layer_factors=True, merge_batch_technical=False, group_order=None, show_group_track=True, figsize=None, figwidth=10.0, row_height=0.004, group_track_width=0.05, cmap='viridis', norm=None, vmin=None, vmax=None, colorbar_label=None, group_separator_linewidth=1.0, group_separator_color='white', group_separator_alpha=0.9, xlabel='Factors', factor_label_rotation=90.0, factor_fontsize=8, label_fontsize=10, colorbar_label_fontsize=9, colorbar_tick_fontsize=8, show_annotations=False, annotation_fontsize=8, annotation_rotation=45.0, save=None, show=True)

Heatmap of per-cell factor scores or probabilities for one or more obs subsets.

Restricts to cells matching subset_obs in subset_obs_key (for example one patient or coarse cell type), orders rows by group_obs_key (for example treatment or subtype), and within each group optionally sorts cells by Ward hierarchical clustering on the displayed factor vectors. Columns are factors from layer_idx. By default columns follow the hierarchy graph (model.get_layer_factor_orders()); pass factors to select and order columns explicitly.

When subset_obs contains multiple values, one heatmap panel is drawn per value as a vertical stack of subplots with a shared x-axis (factor columns). The group_obs_key color strip is drawn to the left of each heatmap.

Requires model.annotate() (or equivalent) so layer scores and/or probabilities are stored in adata.obsm as X_<layer> and X_<layer>_probs.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance.

required
subset_obs_key str

adata.obs column selecting the outer group (e.g. patient).

required
subset_obs Union[str, Sequence[str]]

Value or values of subset_obs_key to include.

required
group_obs_key str

adata.obs column defining row blocks (e.g. treatment).

required
layer_idx int

Layer whose factors form the heatmap columns.

0
factors Optional[Sequence[Union[str, int]]]

Optional list of factor names (or indices within the layer) to plot. Columns appear in this order; overrides sort_layer_factors.

None
values Literal['score', 'prob']

"score" for raw cell factor scores; "prob" for normalized probabilities from X_<layer>_probs. With merge_batch_technical, probabilities are instead re-derived by row-normalizing the merged matrix, so a merged column holds the probability mass of the whole group rather than of one batch half (identical to summing the members' probabilities when reduce="sum").

'score'
cluster_cells bool

If True, hierarchically cluster cells within each group_obs_key block (Ward linkage on Euclidean distance).

True
sort_layer_factors bool

If True and factors is None, order columns with model.get_layer_factor_orders(). This applies with merge_batch_technical too: a merged column's members share a parent, so they are contiguous in the hierarchy ordering and the merged column simply takes their place.

True
merge_batch_technical bool

If True, replace each group of batch-technical splits with one merged column via factor_batch_correction, so the program reads as a single column firing across both batch blocks instead of two columns each firing in one. Merged columns are labelled "L0_7+L0_14" so the corrected split stays visible. Layer 0 only. factors may reference merged labels. With values="prob", probabilities are re-derived from the merged scores (see values).

False
group_order Optional[Sequence[str]]

Optional explicit order of group_obs_key categories; default uses categorical order when available, else sorted unique values.

None
show_group_track bool

Draw a narrow color strip for group_obs_key to the left of each heatmap.

True
figsize Optional[Tuple[float, float]]

Figure size (width, height) in inches. If None, inferred from figwidth, row_height, and the number of cells per panel.

None
figwidth float

Figure width in inches when figsize is None.

10.0
row_height float

Height in inches per cell row when figsize is None.

0.004
group_track_width float

Gridspec width ratio for the group strip (default matches the colorbar column width).

0.05
cmap Union[str, Colormap]

Matplotlib colormap name or Colormap instance for the heatmap.

'viridis'
norm Optional[Normalize]

Optional normalization for the heatmap colormap (e.g. matplotlib.colors.Normalize). If set, overrides vmin/vmax.

None
vmin Optional[float]

Optional lower bound for the color scale (ignored when norm is set).

None
vmax Optional[float]

Optional upper bound for the color scale (ignored when norm is set).

None
colorbar_label Optional[str]

Optional colorbar title; inferred from values if None.

None
group_separator_linewidth float

Line width of horizontal separators between group_obs_key blocks.

1.0
group_separator_color str

Color of subgroup separator lines.

'white'
group_separator_alpha float

Alpha of subgroup separator lines.

0.9
xlabel str

X-axis label.

'Factors'
factor_label_rotation float

Rotation of factor names on the x-axis.

90.0
factor_fontsize int

Font size for factor name tick labels.

8
label_fontsize int

Font size for the x-axis label.

10
colorbar_label_fontsize int

Font size for the colorbar title.

9
colorbar_tick_fontsize int

Font size for colorbar tick labels.

8
show_annotations bool

If True, show factor_obs['annotation'] labels on a secondary x-axis at the top of the first (top-row) heatmap panel (45° by default).

False
annotation_fontsize int

Font size for top annotation labels.

8
annotation_rotation float

Rotation of top annotation labels in degrees.

45.0
save Optional[str]

If set, save the figure to this path (in addition to returning it when show=False).

None
show bool

Whether to call plt.show() and return None. If False, returns the Figure so it can be saved or further customized.

True

Returns:

Type Description
Optional[Figure]

Figure if show is False, else None.

obs_factor_dotplot(model, obs_key, layer_idx, cluster_rows=True, cluster_cols=True, figsize=(8, 2), s_min=100, s_max=500, titlesize=12, labelsize=12, legend_fontsize=12, legend_titlesize=12, cmap='viridis', logged=False, width_ratios=[5, 1, 1], show_ylabel=True, show=True)

Plot dotplot showing factor assignments for observations.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
obs_key str

key in model.adata.obs to use for grouping

required
layer_idx int

layer index to plot

required
cluster_rows bool

whether to cluster rows

True
cluster_cols bool

whether to cluster columns

True
figsize Tuple[float, float]

figure size

(8, 2)
s_min int

minimum circle size

100
s_max int

maximum circle size

500
titlesize int

title font size

12
labelsize int

label font size

12
legend_fontsize int

legend font size

12
legend_titlesize int

legend title font size

12
cmap str

colormap name

'viridis'
logged bool

whether to log transform colors

False
width_ratios List[float]

width ratios for subplots

[5, 1, 1]
show_ylabel bool

whether to show y-axis label

True
show bool

whether to show the plot

True

Returns:

Type Description
Optional[Figure]

Figure object if show is False, None otherwise

obs_scores(model, obs_keys, hierarchy=None, mode='fracs', vmax=None, vmin=None, **kwargs)

Plot the association between a set of cell annotations and factors.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
obs_keys Sequence[str]

the keys in model.adata.obs to use

required
hierarchy Optional[Dict[str, Sequence[str]]]

the polytree to restrict the associations to

None
mode Literal['f1', 'fracs', 'weights', 'prob', 'soft_prec', 'score']

scoring method:

  • "f1": hard F1 score (argmax assignments).
  • "fracs": hard precision / factor purity (P(obs=V | assigned to F)).
  • "weights": soft F1 score (harmonic mean of soft precision and soft recall, using X_<layer>_probs).
  • "prob": soft mean membership / soft recall (E[p(F) | obs=V]).
  • "soft_prec": soft precision (fraction of each factor's total soft mass in the obs category).
  • "score": mean raw cell-factor score (X_<layer>) per obs category — the obs-averaged analogue of obs_cell_factor_heatmap(values="score"). Values are not bounded to [0, 1].
'fracs'
**kwargs Any

plotting keyword arguments

{}

path_embedding(model, path_id='auto', paths_key='transition_paths', score_key=None, basis='umap_multilayer', min_affinity=0.0, affinity_alpha_range=(0.15, 1.0), cmap='viridis', point_size=10.0, show_background=True, background_color='lightgray', background_alpha=0.2, obs_key=None, obs_order=None, ncols=3, ax=None, show=True)

Plot cells on an embedding colored by position along one path.

This visualization uses outputs from scd.tl.score_paths: - adata.obsm[f"{score_key}_positions"] - adata.obsm[f"{score_key}_affinities"].

Color encodes path position (0->1), and point alpha scales with path affinity. If obs_key is provided, draws one facet per category value.

path_trajectory_heatmap(model, path_id, paths_key='differentiation_paths', score_key=None, min_sort_affinity=0.05, path_mass_eps=1e-12, genes_per_factor=3, smoothing=50, figwidth=8, gene_height=0.28, block_spacing=1, ytick_prefix_layer=True, genes=None, annotation_obs_key=None, subset_obs_key=None, subset_obs=None, heatmap_cmap='RdYlBu_r', factor_heatmap_cmap='viridis', colorbar_gap=0.16, xlabel='Cells', normalize=True, save=None, show=True)

Trajectory heatmap along a stored multi-layer path (differentiation or transition).

Uses adata.uns[paths_key]['paths'][path_id] factor nodes (root→leaf for differentiation paths) only for cell ordering along the path.

Default mode: per-node factor scores from adata.obs['<factor>_score'] (annotate_adata / fit) plus confident genes from scd.tl.set_confident_signatures(model).

Custom genes (genes=[...]): one heatmap block, one row per gene; no factor score rows and no confident-signature cache required.

Cell order

Prefer scd.tl.score_paths matrices {score_key}_positions and {score_key}_affinities (same column order as paths). Cells must meet min_sort_affinity and have finite positions. If no such cells remain, falls back to the same logic as trajectory_heatmap: anchor cells on the L0 terminus (differentiation) or transition source/target on L0, then sort by probability-weighted index along nodes using X_<layer>_probs.

If genes is set to a non-empty list of gene names, skips per-factor blocks (no <factor>_score rows and no confident signatures). Renders a single heatmap block with one row per gene (same path-based cell order). genes_per_factor is ignored in that mode.

Parameters:

Name Type Description Default
model scDEF

fitted scDEF model.

required
path_id int

path_id field stored with the path, or a list index.

required
paths_key str

adata.uns key from build_differentiation_paths or build_transition_paths.

'differentiation_paths'
score_key Optional[str]

prefix for obsm position/affinity matrices from score_paths; defaults to paths_key.

None
min_sort_affinity float

when using score_paths output, minimum affinity along the path for a cell to be included and ordered.

0.05
path_mass_eps float

minimum summed path-node probability mass in fallback mode.

1e-12
ytick_prefix_layer bool

if True, factor row labels are '<layer> <name>'.

True
genes Optional[Sequence[str]]

optional list of adata.var_names entries to plot as rows only.

None
normalize bool

min–max scale each gene expression row after smoothing (default True); factor score rows are always scaled. Same as trajectory_heatmap.

True

pathway_scores(model, pathways, top_genes=20, **kwargs)

Plot the association between a set of cell annotations and a set of gene signatures.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
pathways DataFrame

a pandas DataFrame containing PROGENy pathways

required
top_genes Optional[int]

number of top genes to consider

20
**kwargs Any

plotting keyword arguments

{}

qc(model, figsize=(8, 12), show=True)

Plot QC metrics for scDEF run.

Plots include: loss over epochs, BRD vs Gini coefficient, learned vs observed cell scales, learned vs observed gene scales, and biological relevance determination. If trace diagnostics are available (e.g. n_eff_parents_trace), a trace-oriented layout is used.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
figsize Tuple[float, float]

figure size in inches

(8, 12)
show bool

whether to show the plot

True

Returns: Figure object if show is False, None otherwise

relevance(model, mode='brd', thres=None, iqr_mult=None, show_yticks=False, scale='linear', normalize=False, fontsize=14, legend_fontsize=12, xlabel='Factor', ylabel='Relevance', color=False, show=True, ax=None, **kwargs)

Plot relevance determination scores.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
mode Literal['brd', 'ard']

mode to plot, either "brd" or "ard"

'brd'
thres Optional[float]

threshold value for relevance cutoff

None
iqr_mult Optional[float]

multiplier for IQR-based threshold

None
show_yticks bool

whether to show y-axis ticks

False
scale Literal['linear', 'log']

scale for y-axis, either "linear" or "log"

'linear'
normalize bool

whether to normalize relevance scores

False
fontsize int

font size for labels

14
legend_fontsize int

font size for legend

12
xlabel str

label for x-axis

'Factor'
ylabel str

label for y-axis

'Relevance'
color bool

whether to color bars by factor type

False
show bool

whether to show the plot

True
ax Optional[Axes]

matplotlib axes to plot on

None
**kwargs Any

additional plotting keyword arguments

{}

Returns:

Type Description
Optional[Axes]

Axes object if show is False, None otherwise

scale(model, scale_type, figsize=(4, 4), alpha=0.6, fontsize=12, legend_fontsize=10, ax=None, show=True)

Plot learned scale factors vs observed scales.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
scale_type Literal['cell', 'gene']

type of scale to plot, either "cell" or "gene"

required
figsize Tuple[float, float]

figure size

(4, 4)
alpha float

transparency level

0.6
fontsize int

font size for labels

12
legend_fontsize int

font size for legend

10
ax Optional[Axes]

matplotlib axes to plot on

None
show bool

whether to show the plot

True

Returns:

Type Description
Optional[Axes]

Axes object if show is False, None otherwise

scales(model, figsize=(8, 4), alpha=0.6, fontsize=12, legend_fontsize=10, show=True)

Plot both cell and gene scales.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
figsize Tuple[float, float]

figure size

(8, 4)
alpha float

transparency level

0.6
fontsize int

font size for labels

12
legend_fontsize int

font size for legend

10
show bool

whether to show the plot

True

Returns:

Type Description
Optional[Figure]

Figure object if show is False, None otherwise

signatures_scores(model, obs_keys, markers, top_genes=10, hierarchy=None, **kwargs)

Plot the association between a set of cell annotations and a set of gene signatures.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
obs_keys Sequence[str]

the keys in model.adata.obs to use

required
markers Mapping[str, Sequence[str]]

a dictionary with keys corresponding to model.adata.obs[obs_keys] and values to gene lists

required
top_genes Optional[int]

number of genes to consider in the score computations

10
hierarchy Optional[Dict[str, Sequence[str]]]

the polytree to restrict the associations to

None
**kwargs Any

plotting keyword arguments

{}

technical_hierarchy(model, show_signatures=True, **kwargs)

Plot the technical hierarchy of the model.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
show_signatures bool

whether to show gene signatures

True
**kwargs Any

keyword arguments passed to make_graph

{}

Returns:

Type Description
Graph

Graphviz Graph object

trajectory_heatmap(model, factor_path, layer_idx=0, l0_path=None, genes_per_factor=3, smoothing=50, figwidth=8, gene_height=0.28, block_spacing=1, annotation_obs_key=None, subset_obs_key=None, subset_obs=None, heatmap_cmap='RdYlBu_r', factor_heatmap_cmap='viridis', colorbar_gap=0.16, xlabel='Cells', normalize=True, save=None, show=True)

Plot stacked trajectory heatmap from precomputed confident signatures.

Requires scd.tl.set_confident_signatures(model) to be run beforehand.

Cell sorting

Cells are restricted to those assigned to factor_path (and optional subset_obs filter), then ordered by a soft trajectory progress score. For each selected cell, the function takes its factor probabilities (from X_<layer>_probs), keeps only the columns for the path factors, and computes a probability-weighted average of path positions 0..len(path)-1. This gives a continuous progress coordinate used to sort cells from early to late along the path.

If layer_idx > 0 and l0_path is provided, the sorting weights are computed on layer 0 probabilities instead. If that yields near-zero total path weight for all selected cells, sorting falls back to the layer_idx path probabilities.

normalize (default True): min–max scale each gene row to [0, 1] after smoothing; if False, use smoothed raw expression. Factor score rows are always min–max scaled.

umap(model, color=[], layers=None, figsize=(16, 4), fontsize=12, legend_fontsize=10, rasterized=True, n_legend_cols=1, factor_subset=None, show=True)

Plot pre-computed UMAPs for different layers.

UMAP embeddings must have been computed first via scdef.tl.umap. Each panel temporarily points adata.obsm['X_umap'] at X_umap_{layer_name}; the original X_umap is restored afterward. When layers is None, panels follow ascending layer index (L0, L1, ...); otherwise the given layers order is used.

Parameters:

Name Type Description Default
model scDEF

scDEF model instance

required
color Union[str, List[str]]

color key(s) to use for coloring

[]
layers Optional[List[Union[int, str]]]

which layers to plot, in panel order. Entries are layer indices; a string entry is used as an embedding name directly, so layers=['L0_corrected'] plots adata.obsm['X_umap_L0_corrected']. scdef.tl.umap writes that embedding whenever a corrected representation is already present in adata.obsm['X_L0_batch_corrected']; it never builds the correction itself, so run the correction first. Indices and names can be mixed, e.g. layers=[0, 'L0_corrected'] for a before/after pair.

None
figsize Tuple[float, float]

figure size

(16, 4)
fontsize int

font size for labels

12
legend_fontsize int

legend font size

10
rasterized bool

whether to rasterize the plot

True
n_legend_cols int

number of columns in legend

1
factor_subset Optional[List[str]]

subset of factors to plot

None
show bool

whether to show the plot

True

within_group_pairwise_dissimilarity(model, layer, obs_key, metric='jsd', kind='box', showfliers=False, figsize=(8, 4), show=True)

Plot within-group pairwise dissimilarity distributions.

Uses cached results from model.adata.uns['within_group_pairwise_dissimilarity'] when available, otherwise computes them with scd.tl.compute_within_group_pairwise_dissimilarity.