diff --git a/analysis/figures/figure2f_h.py b/analysis/figures/figure2f_h.py
index e2170d7e..db6c7c49 100644
--- a/analysis/figures/figure2f_h.py
+++ b/analysis/figures/figure2f_h.py
@@ -35,7 +35,7 @@ def makeFigure():
     return f
 
 
-def groupDrugs(X: anndata, label_name: str):
+def groupDrugs(X: anndata.AnnData, label_name: str):
     """Groups drugs of similar category"""
     labels = X.obs[label_name]
     names = np.unique(labels)
diff --git a/analysis/figures/figure3.py b/analysis/figures/figure3.py
index 1afcb24a..e3c400cf 100644
--- a/analysis/figures/figure3.py
+++ b/analysis/figures/figure3.py
@@ -104,20 +104,6 @@ def makeFigure():
     return f
 
 
-def plot_cell_perc_corr(cellDF: pd.DataFrame, pop1: str, pop2: str, ax: Axes):
-    """Plots correlation of cell percentages against each other"""
-    newDF = pd.DataFrame()
-    newDF2 = pd.DataFrame()
-    newDF[[pop1, "Condition"]] = cellDF.loc[cellDF["Cell Type"] == pop1][
-        ["Cell Type Percentage", "Condition"]
-    ]
-    newDF2[[pop2, "Condition"]] = cellDF.loc[cellDF["Cell Type"] == pop2][
-        ["Cell Type Percentage", "Condition"]
-    ]
-    newDF = newDF.merge(newDF2, on="Condition")
-    sns.scatterplot(newDF, x=pop1, y=pop2, hue="Condition", ax=ax)
-
-
 def plot_cell_perc_comp_corr(
     X: anndata.AnnData, cellDF: pd.DataFrame, pop: str, comp: int, ax: Axes, unique=None
 ):
diff --git a/analysis/figures/figure4e_k.py b/analysis/figures/figure4e_k.py
index f8ad3cf5..ddc15a8a 100644
--- a/analysis/figures/figure4e_k.py
+++ b/analysis/figures/figure4e_k.py
@@ -53,7 +53,7 @@ def makeFigure():
 
 
 def plot_correlation_cmp_cell_count_perc(
-    X: anndata, cmp: int, cellcountDF: pd.DataFrame, ax: Axes, cellPerc=True
+    X: anndata.AnnData, cmp: int, cellcountDF: pd.DataFrame, ax: Axes, cellPerc=True
 ):
     """Plot component weights by cell type count or percentage for a cell type"""
     yt = np.unique(X.obs["Condition"])
diff --git a/analysis/figures/figureS9e_f.py b/analysis/figures/figureS9e_f.py
index 9d36dad4..338dcd92 100644
--- a/analysis/figures/figureS9e_f.py
+++ b/analysis/figures/figureS9e_f.py
@@ -48,8 +48,8 @@ def makeFigure():
 
 def logreg_weights_scores(
     X: anndata.AnnData, y: pd.Series, prediction: str
-) -> pd.DataFrame:
-    """Fit logistic regression model, return coefficients of that model"""
+) -> tuple[pd.DataFrame, float]:
+    """Fit logistic regression, returning its coefficients and accuracy score."""
     status = y[prediction]
     cond_factors = np.array(X.uns["Pf2_A"])
     lr = logistic_regression("accuracy").fit(cond_factors, status)
diff --git a/analysis/gating.py b/analysis/gating.py
index 3b58113d..ed3cd3c7 100644
--- a/analysis/gating.py
+++ b/analysis/gating.py
@@ -5,7 +5,11 @@
 
 
 def gateThomsonCellsLeiden(X) -> npt.ArrayLike:
-    """Manually gates cell types for Thomson PaCMAP"""
+    """Manually gates cell types for Thomson PaCMAP.
+
+    Run by hand to (re)build ThomsonCellTypes.parquet, which
+    :func:`gateThomsonCells` then reads; deliberately has no in-repo caller.
+    """
     sc.pp.neighbors(X, n_neighbors=15, use_rep="projections", random_state=0)
     sc.tl.leiden(X, resolution=3, random_state=0)
     X.obs["Cell Type"] = X.obs.leiden.replace(thomson_layer1).astype(str)
@@ -26,7 +30,11 @@ def gateThomsonCells(X) -> npt.ArrayLike:
 
 
 def Thomson_Doublet():
-    """Detects doublets in scRNA-seq"""
+    """Detects doublets in scRNA-seq.
+
+    Run by hand to (re)build ThomsonDoublets.csv; deliberately has no in-repo
+    caller.
+    """
     import hdf5plugin  # noqa: F401
     import vsparse
 
@@ -49,38 +57,6 @@ def Thomson_Doublet():
     X.obs["doublet"].to_csv("analysis/data/Thomson/ThomsonDoublets.csv")
 
 
-def getHiResOldLupus(X) -> npt.ArrayLike:
-    """Manually gates cell types for SLE PaCMAP"""
-    X.obs["Cell Type Old2"] = X.obs["Cell Type Old"].astype(str)
-    X.obs.cell_type_lympho = X.obs.cell_type_lympho.astype(str)
-    X.obs.loc[X.obs["cell_type_lympho"] != "nan", "Cell Type Old2"] = X.obs.loc[
-        X.obs["cell_type_lympho"] != "nan"
-    ].cell_type_lympho.values
-    X.obs = X.obs.replace({"Cell Type Old2": cell_type_conv})
-
-    return X
-
-
-cell_type_conv = {
-    "T4_naive": "T4 Naive",
-    "B_naive": "B Naive",
-    "CytoT_GZMH+": "T8 GZMH",
-    "T4_em": "T4 EM",
-    "NK_dim": "NK Dim",
-    "CytoT_GZMK+": "T8 GZMK",
-    "T8_naive": "T8 Naive",
-    "ncM": "nCM",
-    "B_mem": "B mem",
-    "T4_reg": "T4 Reg",
-    "cDC": "cDC2",
-    "T_mait": "T mait",
-    "B_plasma": "B plasma",
-    "NK_bright": "NK Bright",
-    "B_atypical": "B Atypical",
-    "cM": "CM",
-}
-
-
 thomson_layer1 = {
     "0": "NK Cells",
     "1": "T Cells",
@@ -179,7 +155,8 @@ def getHiResOldLupus(X) -> npt.ArrayLike:
 }
 
 
-# Taken directly from PopAlign and used to annotate level 1
+# Taken directly from PopAlign and used to annotate level 1. Reference data
+# consulted when curating the layer maps above by hand -- no code reads it.
 marker_genes_1 = {
     "Monocytes": [
         "CD14",
@@ -211,7 +188,8 @@ def getHiResOldLupus(X) -> npt.ArrayLike:
 }
 
 
-# Take from various resources to annotate level 2
+# Taken from various resources to annotate level 2. Reference data consulted
+# when curating the layer maps above by hand -- no code reads it.
 marker_genes_2 = {
     "B cells": ["PXK", "MS4A1", "CD19", "CD74", "CD79A", "BANK1", "PTPRC", "CR2"],
     "B Memory": ["NPIB15", "BACH2", "IL7", "NMBR", "MS4A1", "MBL2", "LY86", "CD27"],
diff --git a/scrise/annotation_alignment.py b/scrise/annotation_alignment.py
index a54f4cee..7a6ff9f0 100644
--- a/scrise/annotation_alignment.py
+++ b/scrise/annotation_alignment.py
@@ -294,6 +294,40 @@ def cell_type_alignment(
     )
 
 
+def _component_metrics(
+    loadings_matrix: np.ndarray,
+    codes: np.ndarray,
+    n_types: int,
+    signed: bool,
+    n_permutations: int,
+    rng: np.random.Generator,
+) -> tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray, np.ndarray]:
+    """Per-component AUROCs, permutation p-values, tau, eta^2 and epsilon^2."""
+    n_cells, n_comps = loadings_matrix.shape
+    aurocs_mat = np.zeros((n_comps, n_types), dtype=float)
+    p_vals_mat = np.ones((n_comps, n_types), dtype=float)
+    tau_vec = np.zeros(n_comps, dtype=float)
+    eta2_vec = np.zeros(n_comps, dtype=float)
+    eps2_vec = np.zeros(n_comps, dtype=float)
+
+    for comp_idx in range(n_comps):
+        y = loadings_matrix[:, comp_idx].astype(float)
+        if signed:
+            y = np.abs(y)
+
+        auc = compute_auroc_per_cell_type(y, codes, n_types)
+        aurocs_mat[comp_idx] = auc
+        tau_vec[comp_idx] = compute_tau(auc, baseline=0.0)
+        eta2_vec[comp_idx] = compute_eta_squared(y, codes, n_types)
+        eps2_vec[comp_idx] = compute_kruskal_epsilon_squared(y, codes, n_types)
+
+        p_vals_mat[comp_idx] = _component_p_values(
+            y, codes, n_types, n_cells, auc, n_permutations, rng
+        )
+
+    return aurocs_mat, p_vals_mat, tau_vec, eta2_vec, eps2_vec
+
+
 def score_cell_type_alignment(
     data: anndata.AnnData | np.ndarray | pd.DataFrame,
     cell_types: pd.Series | np.ndarray | str | None = None,
@@ -353,26 +387,9 @@ def score_cell_type_alignment(
     rng = _as_generator(random_state)
 
     component_labels = [i + 1 for i in range(n_comps)]
-    aurocs_mat = np.zeros((n_comps, n_types), dtype=float)
-    p_vals_mat = np.ones((n_comps, n_types), dtype=float)
-    tau_vec = np.zeros(n_comps, dtype=float)
-    eta2_vec = np.zeros(n_comps, dtype=float)
-    eps2_vec = np.zeros(n_comps, dtype=float)
-
-    for comp_idx in range(n_comps):
-        y = loadings_matrix[:, comp_idx].astype(float)
-        if signed:
-            y = np.abs(y)
-
-        auc = compute_auroc_per_cell_type(y, codes, n_types)
-        aurocs_mat[comp_idx] = auc
-        tau_vec[comp_idx] = compute_tau(auc, baseline=0.0)
-        eta2_vec[comp_idx] = compute_eta_squared(y, codes, n_types)
-        eps2_vec[comp_idx] = compute_kruskal_epsilon_squared(y, codes, n_types)
-
-        p_vals_mat[comp_idx] = _component_p_values(
-            y, codes, n_types, n_cells, auc, n_permutations, rng
-        )
+    aurocs_mat, p_vals_mat, tau_vec, eta2_vec, eps2_vec = _component_metrics(
+        loadings_matrix, codes, n_types, signed, n_permutations, rng
+    )
 
     # Joint BH FDR correction across all (component x cell_type) tests
     if aurocs_mat.size > 1:
diff --git a/scrise/factor_io.py b/scrise/factor_io.py
index 2816af08..6043636f 100644
--- a/scrise/factor_io.py
+++ b/scrise/factor_io.py
@@ -111,28 +111,28 @@ def _read_raw_dataset(raw_path: str) -> anndata.AnnData:
         return anndata.read_h5ad(raw_path)
 
 
-def _attach_raw_data(adata: anndata.AnnData, raw_path: str) -> None:
-    """Match the factors' cells and genes against raw data and attach ``X``."""
-    raw = _read_raw_dataset(raw_path)
-
+def _subset_raw_cells(adata: anndata.AnnData, raw: anndata.AnnData) -> anndata.AnnData:
+    """Align the raw cells to the factors' barcodes, retrying via ``cell_barcode``."""
     if not isinstance(raw.obs, pd.DataFrame):
         raise TypeError("raw.obs must be an in-memory pandas DataFrame.")
 
-    # Match cells by index or cell_barcode column
     if not np.all(adata.obs_names.isin(raw.obs_names)) and "cell_barcode" in raw.obs:
         raw.obs.index = pd.Index(raw.obs["cell_barcode"].astype(str))
 
-    common_cells = adata.obs_names[adata.obs_names.isin(raw.obs_names)]
-    if len(common_cells) == 0:
+    if not adata.obs_names.isin(raw.obs_names).any():
         raise ValueError(
             "No matching cell barcodes found between factors and raw data."
         )
-    raw_sub = raw[adata.obs_names, :].copy()
+    return raw[adata.obs_names, :].copy()
 
+
+def _subset_raw_genes(
+    adata: anndata.AnnData, raw_sub: anndata.AnnData
+) -> anndata.AnnData:
+    """Align the raw genes to the factors' names, retrying via ``gene_ids``."""
     if not isinstance(raw_sub.var, pd.DataFrame):
         raise TypeError("raw_sub.var must be an in-memory pandas DataFrame.")
 
-    # Match genes
     if (
         not np.all(adata.var_names.isin(raw_sub.var_names))
         and "gene_ids" in raw_sub.var
@@ -140,10 +140,16 @@ def _attach_raw_data(adata: anndata.AnnData, raw_path: str) -> None:
     ):
         raw_sub.var.index = pd.Index(raw_sub.var["gene_ids"].astype(str))
 
-    common_genes = adata.var_names[adata.var_names.isin(raw_sub.var_names)]
-    if len(common_genes) == 0:
+    if not adata.var_names.isin(raw_sub.var_names).any():
         raise ValueError("No matching gene names found between factors and raw data.")
-    raw_sub = raw_sub[:, adata.var_names].copy()
+    return raw_sub[:, adata.var_names].copy()
+
+
+def _attach_raw_data(adata: anndata.AnnData, raw_path: str) -> None:
+    """Match the factors' cells and genes against raw data and attach ``X``."""
+    raw_sub = _subset_raw_genes(
+        adata, _subset_raw_cells(adata, _read_raw_dataset(raw_path))
+    )
 
     from parafac2.normalize import prepare_dataset
 
diff --git a/scrise/factorization.py b/scrise/factorization.py
index ea17061f..a18fe01f 100644
--- a/scrise/factorization.py
+++ b/scrise/factorization.py
@@ -269,6 +269,12 @@ def pf2(
         Additional keyword arguments forwarded to ``parafac2_nd`` (e.g.
         ``n_inner``, ``callback``), for options not otherwise exposed here.
 
+    condition_key : str, optional (default: None)
+        Column in ``X.obs`` holding the condition labels, used to derive
+        ``condition_unique_idxs`` when that column is not already present.
+    adata : anndata.AnnData, optional (default: None)
+        Alias for ``X``; supply either one, not both.
+
     Returns
     -------
     anndata.AnnData
diff --git a/scrise/plotting/factors.py b/scrise/plotting/factors.py
index 7ec7c813..e71a98a1 100644
--- a/scrise/plotting/factors.py
+++ b/scrise/plotting/factors.py
@@ -111,6 +111,12 @@ def plot_condition_factors(
         Custom colors for condition group labels. If None, uses default palette.
     group_cond : bool, optional (default: False)
         If True and cond_group_labels provided, sorts conditions by group.
+    control_pattern : str, optional (default: None)
+        Substring identifying the control conditions whose median and spread
+        set the normalization. Overrides ``ThomsonNorm``'s implied 'CTRL'.
+    control_conditions : Sequence[str], optional (default: None)
+        Explicit control condition names, taking precedence over
+        ``control_pattern``. With neither given, every condition contributes.
     """
     yt = pd.Series(np.unique(data.obs[cond]))
     X = np.array(data.uns["Pf2_A"])
@@ -152,6 +158,22 @@ def plot_condition_factors(
     ax.set(xlabel="Component")
 
 
+def _draw_unit_scaled_heatmap(X: np.ndarray, yt, ax: Axes) -> None:
+    """Draw a factor heatmap rescaled so the largest magnitude sits at +/-1."""
+    X = X / np.max(np.abs(np.asarray(X)))
+    sns.heatmap(
+        data=X,
+        xticklabels=np.arange(1, X.shape[1] + 1),
+        yticklabels=yt,
+        ax=ax,
+        center=0,
+        cmap=cmap,
+        vmin=-1,
+        vmax=1,
+    )
+    ax.set(xlabel="Component")
+
+
 def plot_eigenstate_factors(data: anndata.AnnData, ax: Axes):
     """Plot eigen-state factors as a heatmap showing cell state patterns.
 
@@ -168,23 +190,8 @@ def plot_eigenstate_factors(data: anndata.AnnData, ax: Axes):
     ax : matplotlib.axes.Axes
         Matplotlib axes object to plot on.
     """
-    rank = data.uns["Pf2_B"].shape[1]
-    xticks = np.arange(1, rank + 1)
-    X = data.uns["Pf2_B"]
-    X = X / np.max(np.abs(np.array(X)))
-    yt = np.arange(1, rank + 1)
-
-    sns.heatmap(
-        data=X,
-        xticklabels=xticks,
-        yticklabels=yt,
-        ax=ax,
-        center=0,
-        cmap=cmap,
-        vmin=-1,
-        vmax=1,
-    )
-    ax.set(xlabel="Component")
+    X = np.asarray(data.uns["Pf2_B"])
+    _draw_unit_scaled_heatmap(X, np.arange(1, X.shape[1] + 1), ax)
 
 
 def plot_gene_factors(data: anndata.AnnData, ax: Axes, weight=0.08, trim=True):
@@ -209,7 +216,6 @@ def plot_gene_factors(data: anndata.AnnData, ax: Axes, weight=0.08, trim=True):
     trim : bool, optional (default: True)
         If True, filters genes based on the weight parameter. If False, shows all genes.
     """
-    rank = data.varm["Pf2_C"].shape[1]
     X = np.array(data.varm["Pf2_C"])
     yt = data.var.index.values
 
@@ -224,22 +230,7 @@ def plot_gene_factors(data: anndata.AnnData, ax: Axes, weight=0.08, trim=True):
         yt = yt[kept_idxs]
 
     ind = reorder_table(X)
-    X = X[ind]
-    X = X / np.max(np.abs(X))
-    yt = yt[ind]
-    xticks = np.arange(1, rank + 1)
-
-    sns.heatmap(
-        data=X,
-        xticklabels=xticks,
-        yticklabels=yt,
-        ax=ax,
-        center=0,
-        cmap=cmap,
-        vmin=-1,
-        vmax=1,
-    )
-    ax.set(xlabel="Component")
+    _draw_unit_scaled_heatmap(X[ind], yt[ind], ax)
 
 
 def reorder_table(projs: np.ndarray) -> np.ndarray:
diff --git a/scrise/plotting/general.py b/scrise/plotting/general.py
index e636c2ca..efc9abf8 100644
--- a/scrise/plotting/general.py
+++ b/scrise/plotting/general.py
@@ -139,6 +139,17 @@ def avegene_per_status(X: anndata.AnnData, gene: str, cellType="Cell Type"):
     return df
 
 
+def _two_gene_frame(
+    X: anndata.AnnData, hue: str, cellType: str
+) -> tuple[pd.Index, pd.DataFrame]:
+    """Mean-centered expression of exactly two genes, tagged by hue and cell type."""
+    assert X.shape[1] == 2
+    dataDF = X.to_df().subtract(X.var["means"].values)
+    dataDF[hue] = X.obs[hue].values
+    dataDF["Cell Type"] = X.obs[cellType].values
+    return X.var_names, dataDF
+
+
 def gene_plot_cells(
     X: anndata.AnnData,
     hue: str,
@@ -149,12 +160,7 @@ def gene_plot_cells(
     cellType="Cell Type",
 ):
     """Plots two genes on either a per cell or per cell type basis"""
-    assert X.shape[1] == 2
-    genes = X.var_names
-    dataDF = X.to_df()
-    dataDF = dataDF.subtract(X.var["means"].values)
-    dataDF[hue] = X.obs[hue].values
-    dataDF["Cell Type"] = X.obs[cellType].values
+    genes, dataDF = _two_gene_frame(X, hue, cellType)
     alpha = 1
 
     if average:
@@ -189,12 +195,7 @@ def plot_cell_gene_corr(
     cellType="Cell Type",
 ):
     """Plots two genes on either a per cell or per cell type basis"""
-    assert X.shape[1] == 2
-    genes = X.var_names
-    dataDF = X.to_df()
-    dataDF = dataDF.subtract(X.var["means"].values)
-    dataDF[hue] = X.obs[hue].values
-    dataDF["Cell Type"] = X.obs[cellType].values
+    genes, dataDF = _two_gene_frame(X, hue, cellType)
 
     meanDF = dataDF.groupby([hue, "Cell Type"], observed=True).mean().reset_index()
     pivoted = meanDF.pivot(index=hue, columns="Cell Type", values=genes)
diff --git a/scrise/plotting/pacmap.py b/scrise/plotting/pacmap.py
index e1005fd2..c2050d8c 100644
--- a/scrise/plotting/pacmap.py
+++ b/scrise/plotting/pacmap.py
@@ -70,6 +70,33 @@ def _render_pacmap(
     return assign_labels(ax)
 
 
+def _plot_continuous_pacmap(
+    points: np.ndarray,
+    values: np.ndarray,
+    ax: Axes,
+    cmap,
+    span: tuple[float, float],
+    title: str,
+) -> None:
+    """Shade the embedding by a continuous per-cell value, with a colorbar."""
+    data = pd.DataFrame(points, columns=["x", "y"])
+    data["val_cat"] = values
+
+    _render_pacmap(
+        points=points,
+        data=data,
+        agg_expr=ds.mean("val_cat"),
+        ax=ax,
+        cmap=cmap,
+        span=span,
+    )
+
+    lo, hi = span
+    psm = plt.pcolormesh([[lo, hi], [lo, hi]], cmap=cmap)
+    plt.colorbar(psm, ax=ax)
+    ax.set(title=title)
+
+
 def plot_gene_pacmap(gene: str, X: anndata.AnnData, ax: Axes, clip_outliers=0.9995):
     """Plot PaCMAP embedding colored by gene expression levels.
 
@@ -100,21 +127,7 @@ def plot_gene_pacmap(gene: str, X: anndata.AnnData, ax: Axes, clip_outliers=0.99
     values /= np.max(values)
 
     points = np.array(X.obsm["X_pf2_PaCMAP"])
-    data = pd.DataFrame(points, columns=["x", "y"])
-    data["val_cat"] = values
-
-    _render_pacmap(
-        points=points,
-        data=data,
-        agg_expr=ds.mean("val_cat"),
-        ax=ax,
-        cmap=cmap,
-        span=(0, 1),
-    )
-
-    psm = plt.pcolormesh([[0, 1], [0, 1]], cmap=cmap)
-    plt.colorbar(psm, ax=ax)
-    ax.set(title=f"{gene}")
+    _plot_continuous_pacmap(points, values, ax, cmap, (0.0, 1.0), f"{gene}")
 
 
 def plot_wp_pacmap(X: anndata.AnnData, cmp: int, ax: Axes, cbarMax: float = 1.0):
@@ -145,22 +158,10 @@ def plot_wp_pacmap(X: anndata.AnnData, cmp: int, ax: Axes, cbarMax: float = 1.0)
     cmap = sns.diverging_palette(240, 10, as_cmap=True)
 
     values /= np.max(np.abs(values))
-    data = pd.DataFrame(points, columns=["x", "y"])
-    data["val_cat"] = values
-
-    _render_pacmap(
-        points=points,
-        data=data,
-        agg_expr=ds.mean("val_cat"),
-        ax=ax,
-        cmap=cmap,
-        span=(-cbarMax, cbarMax),
+    _plot_continuous_pacmap(
+        points, values, ax, cmap, (-cbarMax, cbarMax), f"Cmp. {cmp}"
     )
 
-    psm = plt.pcolormesh([[-cbarMax, cbarMax], [-cbarMax, cbarMax]], cmap=cmap)
-    plt.colorbar(psm, ax=ax)
-    ax.set(title=f"Cmp. {cmp}")
-
 
 def plot_labels_pacmap(
     X: anndata.AnnData,
diff --git a/scrise/rank_selection.py b/scrise/rank_selection.py
index 2621d8c4..a45beca7 100644
--- a/scrise/rank_selection.py
+++ b/scrise/rank_selection.py
@@ -427,6 +427,40 @@ def _bicv_trial(
     return results
 
 
+def _resolve_dataset_alias(
+    X: anndata.AnnData | None, adata: anndata.AnnData | None
+) -> anndata.AnnData:
+    """Return whichever of the ``X``/``adata`` aliases was actually passed."""
+    if X is None:
+        X = adata
+    if X is None:
+        raise ValueError("Either X or adata must be provided.")
+    return X
+
+
+def _validate_split_params(
+    n_repeats: int, held_out_cell_frac: float, held_out_gene_frac: float
+) -> None:
+    """Reject repeat counts and held-out fractions that cannot form a split."""
+    if not (0 < held_out_cell_frac < 1) or not (0 < held_out_gene_frac < 1):
+        raise ValueError(
+            "held_out_cell_frac and held_out_gene_frac must both be between 0 and 1."
+        )
+    if n_repeats < 1:
+        raise ValueError("n_repeats must be at least 1.")
+
+
+def _ensure_condition_idxs(X: anndata.AnnData, condition_key: str | None) -> None:
+    """Fill in ``condition_unique_idxs`` from ``condition_key`` when absent."""
+    if "condition_unique_idxs" in X.obs:
+        return
+    if condition_key is None or condition_key not in X.obs:
+        raise KeyError(
+            "X.obs must contain 'condition_unique_idxs', or provide 'condition_key' pointing to a valid column in X.obs."
+        )
+    X.obs["condition_unique_idxs"] = pd.Categorical(X.obs[condition_key]).codes
+
+
 def _resolve_bicv_inputs(
     X: anndata.AnnData | None,
     adata: anndata.AnnData | None,
@@ -443,27 +477,11 @@ def _resolve_bicv_inputs(
     rejects rank requests that cannot yield a well-posed trial at these
     held-out fractions.
     """
-    if X is None and adata is not None:
-        X = adata
-    if X is None:
-        raise ValueError("Either X or adata must be provided.")
+    X = _resolve_dataset_alias(X, adata)
     if ranks is None:
         raise ValueError("ranks must be provided.")
-
-    if not (0 < held_out_cell_frac < 1) or not (0 < held_out_gene_frac < 1):
-        raise ValueError(
-            "held_out_cell_frac and held_out_gene_frac must both be between 0 and 1."
-        )
-    if n_repeats < 1:
-        raise ValueError("n_repeats must be at least 1.")
-
-    if "condition_unique_idxs" not in X.obs:
-        if condition_key is not None and condition_key in X.obs:
-            X.obs["condition_unique_idxs"] = pd.Categorical(X.obs[condition_key]).codes
-        else:
-            raise KeyError(
-                "X.obs must contain 'condition_unique_idxs', or provide 'condition_key' pointing to a valid column in X.obs."
-            )
+    _validate_split_params(n_repeats, held_out_cell_frac, held_out_gene_frac)
+    _ensure_condition_idxs(X, condition_key)
 
     X = X.to_memory() if hasattr(X, "to_memory") else X
 
@@ -554,6 +572,12 @@ def bicv(
         here. See :func:`scrise.factorization.pf2`'s ``normalize_slices``
         for how it can help with unequal cell counts across conditions.
 
+    condition_key : str, optional (default: None)
+        Column in ``X.obs`` holding the condition labels, used to derive
+        ``condition_unique_idxs`` when that column is not already present.
+    adata : anndata.AnnData, optional (default: None)
+        Alias for ``X``; supply either one, not both.
+
     Returns
     -------
     pandas.DataFrame
diff --git a/scrise/tests/conftest.py b/scrise/tests/conftest.py
index 7deb65e0..474dec94 100644
--- a/scrise/tests/conftest.py
+++ b/scrise/tests/conftest.py
@@ -13,7 +13,7 @@
 
 matplotlib.use("Agg")
 
-import matplotlib.pyplot as plt  # noqa: E402
+import matplotlib.pyplot as plt
 
 
 @pytest.fixture(autouse=True)
diff --git a/scrise/tests/test_plotting_factors.py b/scrise/tests/test_plotting_factors.py
index 7f48ed89..728a790a 100644
--- a/scrise/tests/test_plotting_factors.py
+++ b/scrise/tests/test_plotting_factors.py
@@ -33,7 +33,7 @@ def test_plot_condition_factors_smoke(rank, n_cond):
     # hit log(negative).
     adata.uns["Pf2_A"] = np.abs(adata.uns["Pf2_A"]) + 0.1
 
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_condition_factors(adata, ax)
     assert ax.get_xlabel() == "Component"
     assert len(ax.get_yticklabels()) == n_cond
@@ -49,7 +49,7 @@ def test_plot_condition_factors_with_group_labels_and_legend():
         ["groupA" if i % 2 == 0 else "groupB" for i in range(len(yt))]
     )
 
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_condition_factors(adata, ax, cond_group_labels=group_labels, group_cond=True)
     # A legend patch per unique group should have been added.
     legend = ax.get_legend()
@@ -59,7 +59,7 @@ def test_plot_condition_factors_with_group_labels_and_legend():
 
 def test_plot_eigenstate_factors_smoke():
     adata = make_mock_factored_adata(rank=4)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_eigenstate_factors(adata, ax)
     assert ax.get_xlabel() == "Component"
     assert len(ax.get_xticklabels()) == 4
@@ -68,18 +68,18 @@ def test_plot_eigenstate_factors_smoke():
 @pytest.mark.parametrize("trim", [True, False])
 def test_plot_gene_factors_smoke(trim):
     adata = make_mock_factored_adata(n_genes=30, rank=3)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_gene_factors(adata, ax, trim=trim)
     assert ax.get_xlabel() == "Component"
 
 
 def test_plot_gene_factors_trim_reduces_or_keeps_gene_count():
     adata = make_mock_factored_adata(n_genes=30, rank=3)
-    fig1, ax1 = plt.subplots()
+    _fig1, ax1 = plt.subplots()
     plot_gene_factors(adata, ax1, trim=False, weight=0.08)
     n_all = len(ax1.get_yticklabels())
 
-    fig2, ax2 = plt.subplots()
+    _fig2, ax2 = plt.subplots()
     plot_gene_factors(adata, ax2, trim=True, weight=0.08)
     n_trimmed = len(ax2.get_yticklabels())
 
@@ -106,7 +106,7 @@ def test_plot_bicv_r2x_smoke_and_axis_labels():
     X = make_synthetic_pf2_data(n_cond=4, n_genes=20, rank=2, seed=0)
     results = bicv(X, [2, 4], n_repeats=1, random_state=0, max_iter=30)
 
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_bicv_r2x(results, ax)
 
     assert ax.get_xlabel() == "Rank"
@@ -122,15 +122,15 @@ def test_plot_condition_factors_control_pattern_and_conditions():
     cond_names = [f"CTRL_{c % 2}" if c < 2 else f"TRT_{c}" for c in cond_codes]
     adata.obs["Condition"] = pd.Categorical(cond_names)
 
-    fig1, ax1 = plt.subplots()
+    _fig1, ax1 = plt.subplots()
     plot_condition_factors(adata, ax1, control_pattern="CTRL")
     assert len(ax1.get_yticklabels()) > 0
 
-    fig2, ax2 = plt.subplots()
+    _fig2, ax2 = plt.subplots()
     plot_condition_factors(adata, ax2, control_conditions=["CTRL_0", "CTRL_1"])
     assert len(ax2.get_yticklabels()) > 0
 
-    fig3, ax3 = plt.subplots()
+    _fig3, ax3 = plt.subplots()
     plot_condition_factors(adata, ax3, ThomsonNorm=True)
     assert len(ax3.get_yticklabels()) > 0
 
@@ -138,6 +138,6 @@ def test_plot_condition_factors_control_pattern_and_conditions():
 def test_plot_gene_factors_no_genes_pass_weight_raises():
     adata = make_mock_factored_adata(n_genes=10, rank=2)
     adata.varm["Pf2_C"] = np.full_like(adata.varm["Pf2_C"], 0.01)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     with pytest.raises(ValueError, match="No genes exceeded the weight threshold"):
         plot_gene_factors(adata, ax, weight=0.5, trim=True)
diff --git a/scrise/tests/test_plotting_general.py b/scrise/tests/test_plotting_general.py
index 360c8758..bfe13964 100644
--- a/scrise/tests/test_plotting_general.py
+++ b/scrise/tests/test_plotting_general.py
@@ -40,7 +40,7 @@ def _two_gene_adata(n_cells=40, seed=0):
 
 
 def test_rotate_xaxis_and_yaxis_set_tick_rotation():
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     rotate_xaxis(ax, rotation=45)
     rotate_yaxis(ax, rotation=30)
     # matplotlib doesn't expose rotation as a simple getter pre-draw; the
@@ -107,13 +107,13 @@ def test_plot_avegene_per_celltype_smoke():
     X = _two_gene_adata(n_cells)
     X.obs["Condition"] = np.random.choice(["cond_0", "cond_1"], size=n_cells)
 
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_avegene_per_celltype(X, ["gene_a", "gene_b"], ax)
 
 
 def test_gene_plot_cells_smoke_and_shape_assertion():
     X = _two_gene_adata(20)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     gene_plot_cells(X, hue="hue", ax=ax)
 
     # gene_plot_cells asserts exactly 2 genes/columns are present.
@@ -122,7 +122,7 @@ def test_gene_plot_cells_smoke_and_shape_assertion():
         var=pd.DataFrame({"means": [0.0, 0.0, 0.0]}, index=["a", "b", "c"]),
         obs=pd.DataFrame({"hue": ["x"] * 5, "Cell Type": ["t"] * 5}),
     )
-    fig2, ax2 = plt.subplots()
+    _fig2, ax2 = plt.subplots()
     with pytest.raises(AssertionError):
         gene_plot_cells(X_three_genes, hue="hue", ax=ax2)
 
@@ -132,7 +132,7 @@ def test_plot_cell_gene_corr_smoke_with_missing_pivot_columns():
     pivoting, plot_cell_gene_corr should fall back to an empty frame rather
     than raising a KeyError."""
     X = _two_gene_adata(20)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_cell_gene_corr(
         X, hue="hue", cells=["type_0", "type_1"], ax=ax, unique=["grp_0"]
     )
@@ -140,7 +140,7 @@ def test_plot_cell_gene_corr_smoke_with_missing_pivot_columns():
 
 def test_plot_r2x_smoke_and_axes_labels():
     X = make_synthetic_pf2_data(n_cond=4, n_genes=15, rank=2, seed=0)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_r2x(X, np.array([1, 2, 3]), ax, compress=None)
 
     assert ax.get_xlabel() == "Number of Components"
diff --git a/scrise/tests/test_plotting_pacmap.py b/scrise/tests/test_plotting_pacmap.py
index 3199f745..4c76f90e 100644
--- a/scrise/tests/test_plotting_pacmap.py
+++ b/scrise/tests/test_plotting_pacmap.py
@@ -37,7 +37,7 @@ def test_to_hex_returns_hex_strings():
 
 
 def test_assign_labels_sets_pacmap_axis_labels():
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     out = assign_labels(ax)
     assert out is ax
     assert ax.get_xlabel() == "PaCMAP1"
@@ -49,7 +49,7 @@ def test_assign_labels_sets_pacmap_axis_labels():
 def test_plot_gene_pacmap_smoke():
     adata = make_mock_factored_adata(n_cells=50, n_genes=10, with_embedding=True)
     gene = adata.var_names[0]
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_gene_pacmap(gene, adata, ax)
     assert ax.get_title() == gene
 
@@ -57,14 +57,14 @@ def test_plot_gene_pacmap_smoke():
 @pytest.mark.parametrize("cmp", [1, 3])
 def test_plot_wp_pacmap_smoke(cmp):
     adata = make_mock_factored_adata(n_cells=50, rank=3, with_embedding=True)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_wp_pacmap(adata, cmp=cmp, ax=ax)
     assert ax.get_title() == f"Cmp. {cmp}"
 
 
 def test_plot_labels_pacmap_smoke_and_legend_matches_categories():
     adata = make_mock_factored_adata(n_cells=60, with_embedding=True)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_labels_pacmap(adata, labelType="Cell Type", ax=ax)
 
     legend = ax.get_legend()
@@ -78,7 +78,7 @@ def test_plot_labels_pacmap_condition_filter_collapses_to_other():
     present = adata.obs["Cell Type"].unique().tolist()
     keep = [present[0]]
 
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_labels_pacmap(adata, labelType="Cell Type", ax=ax, condition=keep)
 
     legend = ax.get_legend()
diff --git a/scrise/tests/test_plotting_stability.py b/scrise/tests/test_plotting_stability.py
index ed376cb4..7908902d 100644
--- a/scrise/tests/test_plotting_stability.py
+++ b/scrise/tests/test_plotting_stability.py
@@ -96,7 +96,7 @@ def test_resample_preserves_shape_and_samples_with_replacement():
 
 def test_plot_fms_diff_ranks_smoke():
     X = make_synthetic_pf2_data(n_cond=4, n_genes=15, rank=2, seed=0)
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_fms_diff_ranks(X, ax, ranksList=[2], runs=1, compress=None)
 
     assert ax.get_xlabel() == "Component"
@@ -108,7 +108,7 @@ def test_plot_fms_percent_drop_smoke():
     X = make_synthetic_pf2_data(
         n_cond=3, n_genes=10, rank=2, seed=0, cells_per_cond=(20, 30)
     )
-    fig, ax = plt.subplots()
+    _fig, ax = plt.subplots()
     plot_fms_percent_drop(X, ax, percentList=np.array([0, 10]), runs=1, rank=2)
     assert ax.get_ylabel() == "FMS"
     assert ax.get_ylim() == (0.0, 1.0)