Doubly-robust semiparametric inference#

LFC() estimates per-gene log-fold changes using a doubly-robust AIPW estimator. It requires the augmented covariate matrix W = [X | U] where U are the latent factors from fit_gcate(), and produces a DataFrame with effect estimates, standard errors, and BH-adjusted p-values.

For screens with hundreds of perturbations use gcate_lfc_batch(), which runs GCATE and LFC in batches to keep peak memory bounded.

Effect-size columns#

The returned DataFrame reports the treatment-versus-control effect on both natural-log and base-2 scales:

Column

Definition

tau

Natural logarithm of the treated/control mean ratio.

std

Standard error of tau on the natural-log scale.

log2fc

Base-2 log fold change, exactly tau / log(2).

log2fc_se

Standard error of log2fc, exactly std / log(2).

log2fc_se is the standard error of the estimated log2 fold change, not the sample standard deviation of expression. The rescaling leaves stat, p-values, adjusted p-values, and discoveries unchanged. tau and std remain available for backward compatibility. The same columns are returned by gcate_lfc_batch; compatible caches made by older versions are upgraded in memory when loaded.

Choosing the variance estimator#

LFC defaults to usevar='unequal' (Welch inference), which estimates the treatment and control variances separately. Prefer this default when the treatment and control sample sizes, or their propensity-weighted effective sample sizes, are meaningfully unbalanced. Also retain it when arm-specific pseudo-outcome variances may differ and for case-control, bulk, and donor-level pseudo-bulk analyses. Independence of the rows does not imply equal treatment-arm variances: disease severity, biological response, residual composition, library size, treatment imbalance, and heterogeneous expression can all make pooled inference anti-conservative.

For a small, approximately balanced perturbation comparison, usevar='pooled' may provide better power when the independent sampling units and arm-specific pseudo-outcome variances are reasonably comparable. Treat it as an opt-in, empirically justified analysis rather than an automatic small-sample choice. Balanced counts alone do not justify pooling in a case-control study. Pooled inference can produce much smaller standard errors and substantially more discoveries. For a deliberately justified batched analysis, pass lfc_kwargs=dict(usevar='pooled').

There is no universal sample-size ratio at which the recommendation changes. Compare nominal arm sizes, propensity-weighted effective sample sizes, arm-specific pseudo-outcome variability, and the stability of discoveries under both estimators. Retain usevar='unequal' when these diagnostics do not support pooling.

Donor-level independence alone does not establish equal arm variances. For example, the SEA-AD tutorial uses usevar='unequal' because disease severity, inter-individual response, residual cell composition, and library-size variation can produce different gene-wise variability between disease groups.

Welch inference does not itself model within-subject correlation. Repeated cells from the same donor or experimental unit should still be pseudo-bulked or handled with cluster-aware inference; usevar='unequal' only protects against unequal arm variances.

Propensity diagnostics#

estimate_propensity_scores() estimates propensity scores without fitting the outcome model. Use K=5 to obtain out-of-fold scores for positivity and overfitting diagnostics. summarize_propensity_scores() reports overlap, tail mass, and inverse-weight effective sample size, while plot_propensity_scores() compares treatment and control distributions.

Both the standalone estimator and LFC use class-balanced logistic propensity fitting by default. This preserves historical causarray behavior and ensures that standalone overlap diagnostics describe the same nuisance model used for effect estimation. Pass class_weight=None to estimate_propensity_scores or ps_class_weight=None to LFC for a calibrated-probability sensitivity analysis.

LFC uses the standard AIPW pseudo-outcome, which may be negative for individual cells even though its counterfactual mean is positive. Individual pseudo-outcomes are never clipped because doing so can bias the arm means, particularly when a large shared control group is compared with much smaller treatment groups.

For a calibrated-propensity sensitivity analysis, use LFC(..., ps_class_weight=None) and diagnose the matching scores with estimate_propensity_scores(..., class_weight=None).

Treatment-specific covariate diagnostics#

summarize_treatment_associations compares every observed covariate or latent factor with each treatment using shared all-zero controls. It reports pairwise Spearman correlations, standardized mean differences, and Benjamini–Hochberg adjusted p-values. plot_treatment_associations displays either effect-size measure as a heatmap. These are descriptive diagnostics: there is deliberately no automatic threshold or drop decision.

By default the adjustment pools every treatment-by-covariate test into one family, which with many perturbations is large and correspondingly conservative. Pass bh_scope='per_treatment' to adjust within each treatment’s own block; the n_tests_in_family column records how many tests entered each row’s correction either way. In the heatmap, spearman_rho uses a fixed (-1, 1) colour range so that panels stay comparable across subsets, while the unbounded standardized mean difference scales to the data; pass vmax to set a symmetric limit explicitly.

When a scientifically justified sensitivity analysis uses a different propensity design for each treatment, refit_propensity_scores accepts a treatment-specific mapping such as {'Satb2': ['U9']}. Supplying existing raw scores refits only the named treatment columns and carries the others over unchanged, up to the shared clip, which is applied to the whole returned matrix so a single consistent bound reaches LFC; pass clip=None to leave carried-over scores exactly as supplied. For logistic L2 propensity models, penalty_factors_by_treatment can instead retain a covariate while shrinking its coefficient more strongly. For example, {'Satb2': {'log_library_size': 10}} applies ten times the ordinary ridge penalty to standardized log-library size for Satb2 only. This weighted penalty is implemented by rescaling that feature during both fitting and prediction; it is intentionally unavailable for tree and ensemble propensity models. The updated scores can then be passed to LFC together with cached outcome predictions:

associations = summarize_treatment_associations(
    A, W_A, covariate_names=covariate_names,
    covariate_types=covariate_types,
)
pi_filtered, audit = refit_propensity_scores(
    A, W_A,
    pi_hat=estimation['pi_hat_raw'],
    covariate_names=covariate_names,
    penalty_factors_by_treatment={
        'Satb2': {'log_library_size': 10},
    },
)
filtered_results, _ = LFC(
    Y, W, A, W_A,
    Y_hat=estimation['Y_hat'], pi_hat=pi_filtered,
)

Removing a treatment predictor does not create overlap in the underlying population and can omit a genuine measured confounder. Treat filtered fits as sensitivity analyses, distinguish pre-treatment covariates from possible post-treatment variables, and compare propensity overlap, effective sample sizes, and effect estimates before and after filtering.

Filtering can go too far. If it leaves a constant design, the propensity model degenerates to a covariate-free constant and AIPW reduces to an unweighted contrast; refit_propensity_scores raises a RuntimeWarning and sets degenerate_design in the returned audit table. A very large penalty factor shrinks a coefficient without making the design constant, so check score_std in the same table: a value near zero means the scores carry almost no covariate information even though nothing was formally dropped.

Post-hoc diagnostic masks#

align_test_mask aligns a Boolean treatment-by-gene diagnostic mask to an existing causarray result table by (trt, gene_names) labels. It accepts a named DataFrame, a two-level indexed Series, or an array with explicit treatment and gene names. The returned one-dimensional mask follows the original row order of the result table, so it remains correct for reordered or batched results:

keep = align_test_mask(
    df_res, support_mask,
    treatment_names=perturbation_names,
    gene_names=gene_names,
)
df_res_flagged = df_res.assign(support_keep=keep)

Mask alignment does not refit effects, mutate df_res, or change standard errors and p-values. A mask chosen after inspecting outcomes should be reported as a sensitivity diagnostic; retain the original BH-adjusted p-values rather than redefining the testing family post hoc. The Replogle tutorial demonstrates this workflow for several expression-support rules.

The result includes mean_control, mean_treated, and estimable. These are computed from the unclipped pseudo-outcomes. For numerical stability, valid aggregate arm means are floored at thres_diff only when constructing the log ratio and delta-method denominator. A pair with a nonfinite or nonpositive raw aggregate mean remains non-estimable, so the aggregate floor cannot create an extreme discovery from an invalid estimate.

causarray.DR_learner.LFC(Y, W, A, W_A=None, family='nb', offset=False, Y_hat=None, pi_hat=None, cross_est=False, K=None, mask=None, usevar: Literal['unequal', 'pooled'] = 'unequal', thres_min=0.01, thres_diff=0.01, eps_var=0.0001, fdx=False, fdx_alpha=0.05, fdx_c=0.1, verbose=False, backend: str = 'auto', ps_clip=(0.01, 0.99), ps_class_weight='balanced', **kwargs)#

Estimate log-fold changes of treatment effects (LFCs) using AIPW.

Fits a doubly-robust AIPW estimator for the log-ratio of counterfactual means E[Y(1)] / E[Y(0)]. Call this after fit_gcate() to incorporate estimated latent factors into the covariate matrix W.

Parameters:
Yarray, shape (n, p)

Count matrix of outcomes.

Warray, shape (n, d)

Covariate matrix, typically [X | U] where U are the latent factors from GCATE.

Aarray, shape (n, a)

Binary treatment indicator matrix.

W_Aarray or None, shape (n, d_A)

Covariate matrix for the propensity model. If None, W is used.

familystr

GLM family for the outcome model: 'nb' (default) or 'poisson'.

offsetbool or array-like

Log-scale offset for the outcome model. True computes size factors automatically; False or None disables the offset.

Y_hatarray or None, shape (n, p, a, 2)

Pre-computed counterfactual predictions. When provided, cross-fitting is skipped.

pi_hatarray or None, shape (n, a)

Pre-computed propensity scores. When provided, propensity fitting is skipped.

cross_estbool

Whether to use two-fold cross-estimation for nuisance parameters. An explicit K takes precedence.

Kint or None

Number of nuisance-estimation folds. None uses 2 when cross_est=True and 1 otherwise.

maskarray or None, shape (n, a)

Boolean mask indicating eligible cells for each treatment. It limits propensity-model fitting and final estimand computation.

usevarstr

Variance estimator for the AIPW pseudo-outcomes:

  • 'unequal' (default, v0.0.6+): Welch variance s₀²/n₀ + s₁²/n₁ with Welch-Satterthwaite degrees of freedom; p-values use the t-distribution. Prefer this estimator when treatment and control sample sizes or effective sample sizes are meaningfully unbalanced, when arm-specific pseudo-outcome variances may differ, and for case-control, bulk, and donor-level pseudo-bulk analyses. Independence of the rows does not imply equal treatment-arm variances, and biological heterogeneity or imbalance can make pooled inference anti-conservative.

  • 'pooled': pooled-variance estimator (s² + eps_var) / n. For a small, approximately balanced perturbation comparison, this estimator may provide better power when the independent sampling units and arm-specific pseudo-outcome variances are reasonably comparable. Treat it as an opt-in, empirically justified analysis, not as an automatic choice for every small study. Balanced sample counts alone are insufficient for case-control data, where biological heterogeneity commonly favors 'unequal'. Pooled inference can produce substantially smaller standard errors and many more discoveries.

There is no universal arm-size ratio at which the choice should switch. Inspect nominal and propensity-weighted effective sample sizes, arm-specific pseudo-outcome variability, and sensitivity of the discoveries. When those diagnostics are uncertain, retain 'unequal'.

'unequal' accommodates arm-specific variance but does not model within-donor or within-subject correlation. Repeated cells from the same biological unit should still be pseudo-bulked or analyzed with a cluster-aware method; changing usevar alone does not remove pseudoreplication.

Changed in version 0.0.6: Default changed from 'pooled' to 'unequal'. The 'unequal' formula was also corrected from (s₀²/n₀ + s₁²/n₁)/2 to the standard Welch form, which shrinks t-statistics by ≈ √2 relative to v0.0.5. Pass usevar='pooled' to recover pre-v0.0.6 behaviour.

thres_minfloat

Genes whose maximum counterfactual mean is below this threshold are excluded (reported as tau=0, padj=NaN).

thres_difffloat

Genes whose counterfactual means differ by less than this value are excluded.

eps_varfloat

Small constant added to per-arm variances to prevent division by zero.

fdxbool

Whether to apply FDX control (P(FDP > fdx_c) < fdx_alpha).

fdx_alphafloat

Significance level for FDX control.

fdx_cfloat

FDP threshold for FDX control.

verbosebool

Print progress information.

backendstr

GLM backend: "auto" (default), "fast" (force crispyx), or "original" (force statsmodels).

ps_cliptuple(float, float)

Bounds applied to propensity scores used by AIPW. Raw, unclipped scores remain available as estimation['pi_hat_raw'].

ps_class_weightstr, dict or None

Class weighting for the propensity model. 'balanced' remains the default to limit nuisance-model drift; pass None for calibrated probabilities.

**kwargs

Additional arguments forwarded to the GLM fitting functions.

Returns:
df_resDataFrame

Test results with natural-log effect estimate tau and standard error std, together with their base-2 equivalents log2fc and log2fc_se. The fold change is treatment relative to control, and log2fc_se is a standard error (not a sample standard deviation). The result also contains inference columns, raw mean_control and mean_treated counterfactual means, and an estimable flag (plus trt for multiple treatments).

New in version 0.0.8: Added the log2fc and log2fc_se convenience columns. The original natural-log tau and std columns remain unchanged.

causarray.DR_learner.LFC_batch(*args, **kwargs)#

Deprecated alias for gcate_lfc_batch().

Deprecated since version Use: gcate_lfc_batch instead. LFC_batch will be removed in a future release.

causarray.DR_learner.VIM(eta_est, X, id_covs, **kwargs)#

Estimate variable importance measures (VIM) for heterogeneous treatment effects.

Decomposes treatment effect variance into components explained by each covariate using conditional average treatment effect (CATE) regression.

Parameters:
eta_estarray, shape (n, p)

Influence function values from LFC() or compute_causal_estimand().

Xarray, shape (n, d)

Covariate matrix.

id_covsint or array-like of int

Column indices of X to compute VIM for. An integer k is treated as range(k).

Returns:
estimationdict

Dictionary with keys:

'CATE', 'CATE_lower', 'CATE_upper'array, shape (n_covs, n, p)

Conditional average treatment effect and pointwise confidence band.

'VTE'array, shape (p,)

Total variance of the treatment effect (marginal).

'CVTE'array, shape (n_covs, p)

Conditional variance of the treatment effect given each covariate.

'VIM_mean'array, shape (n_covs, p)

VIM point estimate CVTE / VTE - 1 for each covariate and gene.

'VIM_sd'array, shape (n_covs, p)

Standard deviation of the VIM estimate.

causarray.DR_learner.compute_causal_estimand(estimand, Y, W, A, W_A=None, family='nb', offset=False, Y_hat=None, pi_hat=None, mask=None, fdx=False, fdx_B=1000, fdx_alpha=0.05, fdx_c=0.1, verbose=False, random_state=0, backend: str = 'auto', K=1, ps_clip=(0.01, 0.99), ps_class_weight='balanced', **kwargs)#

Estimate causal treatment effects using AIPW with a user-supplied estimand.

Parameters:
estimandcallable

Function that maps influence function values (etas, A) to (eta_est, tau_est, var_est[, df_eff]). See LFC() for an example implementation.

Yarray, shape (n, p)

Count matrix of outcomes.

Warray, shape (n, d)

Covariate matrix (including latent factors from GCATE).

Aarray, shape (n, a)

Binary treatment indicator matrix.

W_Aarray or None, shape (n, d_A)

Covariate matrix for the propensity model. If None, W is used.

familystr

GLM family for the outcome model: 'nb' (default) or 'poisson'.

offsetbool or array-like

Log-scale offset for the outcome model. True computes size factors automatically; False or None disables the offset.

Y_hatarray or None, shape (n, p, a, 2)

Pre-computed counterfactual predictions. When provided, cross-fitting is skipped.

pi_hatarray or None, shape (n, a)

Pre-computed propensity scores. When provided, propensity fitting is skipped.

Kint

Number of folds used for nuisance estimation. The default 1 preserves in-sample fitting.

ps_cliptuple(float, float)

Bounds applied to propensity scores used by AIPW.

ps_class_weightstr, dict or None

Class weighting for the propensity model. 'balanced' preserves the established nuisance fit; pass None for calibrated probabilities.

maskarray or None, shape (n, a)

Boolean mask indicating eligible cells for each treatment. It limits propensity-model fitting and final estimand computation.

fdxbool

Whether to apply FDX control (P(FDP > fdx_c) < fdx_alpha).

fdx_Bint

Number of bootstrap samples for FDX control.

fdx_alphafloat

Significance level for FDX control.

fdx_cfloat

FDP threshold for FDX control.

backendstr

GLM backend: "auto" (default), "fast" (force crispyx), or "original" (force statsmodels).

verbosebool

Print progress information.

**kwargs

Additional arguments forwarded to the GLM fitting functions.

Returns:
df_resDataFrame

Test results produced by estimand. An estimand may optionally return a fifth dictionary whose arrays are added as diagnostic columns.

causarray.DR_learner.gcate_lfc_batch(Y, X, A, r, W_A=None, batch_size=10, n_batches=None, max_cells=2000, n_ctrl=2000, family='nb', offset=True, warm_start_U=False, cache_path=None, random_state=0, verbose=False, gcate_kwargs=None, lfc_kwargs=None, **kwargs)#

Batch-wise GCATE + doubly-robust LFC estimation.

Partitions perturbations into chunks of batch_size, runs fit_gcate_batch() to estimate per-batch latent confounders, then calls LFC() on each batch independently. All large intermediate arrays (res_1, res_2, Y_hat, pi_hat) are freed immediately after each batch so that peak memory is bounded by one batch’s worth of data regardless of the total number of perturbations.

Results can optionally be cached to an HDF5 file (cache_path) so that interrupted runs can be resumed without re-processing completed batches.

Parameters:
Yarray-like or DataFrame, shape (n, p)

Count matrix.

Xarray, shape (n, d)

Covariate matrix (intercept column should be included).

Aarray-like or DataFrame, shape (n, a)

Binary treatment indicator matrix; control cells have all-zero rows.

rint

Number of latent factors.

W_Aarray or None, shape (n, d_A)

Propensity-score covariate matrix. If None, X is used.

batch_sizeint

Perturbations per batch (default 10). Ignored when n_batches is set. Batches are sized evenly with numpy.array_split() so the last batch is never drastically smaller than the others.

n_batchesint or None

Total number of batches. When set, overrides batch_size and perturbations are split as evenly as possible across exactly n_batches batches (e.g. n_batches=2 on a 29-pert dataset gives two batches of 15 and 14).

max_cellsint or None

Maximum pert cells per batch (default 2 000). None means no cap. Ctrl cells are added on top so the actual batch size is at most n_ctrl + max_cells. The cap is rarely active because typical Perturb-seq datasets have only a few hundred cells per perturbation.

n_ctrlint

Number of ctrl cells in the fixed subsample (default 2 000).

familystr

GLM family (default 'nb').

offsetbool or array-like

Offset specification passed to fit_gcate_batch().

warm_start_Ubool

Passed to fit_gcate_batch().

cache_pathstr or None

Path to an HDF5 file used for incremental caching. When set:

  • On entry, any already-computed batches are loaded from the store and their indices are skipped by fit_gcate_batch().

  • After each new batch, the result DataFrame is appended to the store under key /batch_{i:04d}.

  • On exit, all batches (cached + newly computed) are concatenated and returned.

This lets you resume an interrupted run by re-calling the function with the same cache_path — completed batches are not re-run.

random_stateint

RNG seed.

verbosebool

Print per-batch timing.

gcate_kwargsdict or None

Extra keyword arguments forwarded to fit_gcate_batch() (and ultimately fit_gcate()). E.g.:

gcate_kwargs=dict(backend='fast',
                  kwargs_es_1=dict(max_iters=10, rel_tol=2e-4),
                  kwargs_es_2=dict(max_iters=10, rel_tol=2e-4))
lfc_kwargsdict or None

Extra keyword arguments forwarded to LFC() (e.g. usevar, fdx, thres_min). Retain the default usevar='unequal' when arm sizes or effective sample sizes are meaningfully unbalanced, when arm-specific variability may differ, and for case-control, bulk, or pseudo-bulk analyses. For a small, approximately balanced perturbation comparison with comparable pseudo-outcome variability, lfc_kwargs=dict(usevar='pooled') may improve power. There is no universal balance threshold; compare the relevant diagnostics and retain 'unequal' when uncertain.

**kwargs

Additional arguments forwarded to both fit_gcate_batch() and LFC(). When a key collides with gcate_kwargs / lfc_kwargs, the stage-specific dict wins — this lets you scope a kwarg to one stage (e.g. gcate_kwargs=dict( backend='fast') paired with a top-level backend='original' targeting LFC).

Returns:
df_resDataFrame

Concatenated result from all batches. Includes natural-log tau and std columns, base-2 log2fc and log2fc_se columns, and a 'batch' column with the 0-based batch index. Older compatible caches containing only tau and std are upgraded in memory.

causarray.DR_estimation.estimate_propensity_scores(A, X_A, K=1, ps_model='logistic', mask=None, clip=None, random_state=0, verbose=False, class_weight='balanced', **kwargs)#

Estimate per-treatment propensity scores.

Each treatment is compared with the shared all-zero control group. With K > 1, every returned score is predicted by a model that did not train on that cell. Logistic models use class_weight='balanced' by default, matching LFC() and historical causarray fits. Pass class_weight=None for calibrated treatment probabilities.

Parameters:
Aarray-like, shape (n,) or (n, a)

Binary treatment indicators. Rows containing only zeros are controls.

X_Aarray-like, shape (n, d_A)

Covariates used by the propensity model, including an intercept column when fit_intercept=False.

Kint, optional

Number of folds. 1 fits and predicts on all eligible cells; values greater than one produce out-of-fold predictions.

ps_model{‘logistic’, ‘random_forest_cv’, ‘ensemble’}, optional

Propensity model.

maskarray-like or None, shape (n,) or (n, a)

Optional per-treatment eligibility mask for model fitting.

cliptuple(float, float) or None, optional

Bounds applied after prediction. None returns raw probabilities.

random_stateint, optional

Random seed used for fold construction and supported estimators.

class_weightstr, dict or None, optional

Class weighting for logistic propensity estimation. The default 'balanced' matches LFC(); pass None for calibrated probabilities.

Returns:
pi_hatndarray, shape (n, a)

Estimated probabilities P(A_j=1 | X_A).

causarray.DR_estimation.refit_propensity_scores(A, X_A, drop_by_treatment=None, pi_hat=None, treatment_names=None, covariate_names=None, penalty_factors_by_treatment=None, K=1, ps_model='logistic', mask=None, clip=None, random_state=0, verbose=False, class_weight='balanced', **kwargs)#

Refit propensity scores with treatment-specific covariate filtering.

This helper supports sensitivity analyses in which different treatments omit different observed covariates or latent factors, or apply stronger L2 regularization to selected covariates. When existing scores are supplied, only treatments named in either treatment-specific mapping are refitted; the remaining columns are carried over from pi_hat unchanged, up to the shared clip described below. Outcome models are not fit by this function and cached Y_hat values can be reused in LFC().

Parameters:
Aarray-like, shape (n, a)

Binary treatment indicator matrix.

X_Aarray-like, shape (n, d_A)

Full propensity-model design before treatment-specific filtering.

drop_by_treatmentmapping or None

Treatment names or indices mapped to covariate names or indices to remove for that treatment. Defaults to no removals.

pi_hatarray-like or None, shape (n, a)

Existing raw propensity scores. If supplied, treatments absent from drop_by_treatment are not refitted.

treatment_names, covariate_namessequence, optional

Column labels, inferred from DataFrames when possible.

penalty_factors_by_treatmentmapping or None

Treatment names or indices mapped to {covariate: factor} mappings. A factor greater than one applies that multiple of the ordinary L2 penalty to the named coefficient. This is implemented by dividing the covariate by sqrt(factor) during both fitting and prediction and is available only for ps_model='logistic' with an L2 penalty. A factor of one leaves the covariate unchanged. Interpret relative penalties on a common scale; prep_causarray_data standardizes log-library size.

K, ps_model, mask, random_state, verbose, class_weight, **kwargs

Passed to estimate_propensity_scores() for each refitted model.

cliptuple(float, float) or None, optional

Bounds applied to the whole returned matrix, refitted columns and carried-over columns alike, so that a single consistent bound reaches LFC(). Pass None to leave carried-over scores exactly as supplied.

Returns:
pi_updatedndarray, shape (n, a)

Updated propensity scores.

reportDataFrame

Audit table of refitted treatments, retained/dropped covariates, and the resulting score spread. degenerate_design flags a treatment whose retained design is constant, and score_std reports the standard deviation of its refitted scores on eligible rows; both make a collapsed propensity model visible without reading the warning stream.

Warns:
RuntimeWarning

If filtering leaves a constant design for some treatment, in which case its scores carry no covariate information.

Notes

A very large penalty factor shrinks a coefficient towards zero without making the design constant. That case leaves degenerate_design False but drives score_std towards zero.

New in version 0.0.9.

Diagnostics for treatment overlap, propensity scores, and result masks.

causarray.diagnostics.align_test_mask(results: DataFrame, test_mask: ndarray | DataFrame | Series, treatment_names: Sequence[Hashable] | None = None, gene_names: Sequence[Hashable] | None = None, treatment_col: str = 'trt', gene_col: str = 'gene_names') ndarray#

Align a treatment-by-gene diagnostic mask to an LFC result table.

Alignment uses treatment and gene labels rather than row position. This makes it possible to annotate or subset existing causarray results with a post-hoc diagnostic rule without refitting effects or changing their standard errors and p-values. Labels are compared after conversion to strings, matching the labeling convention of batch LFC results.

Parameters:
resultsDataFrame

Result table containing one row per treatment-gene test.

test_maskndarray, DataFrame, or Series

Boolean diagnostic mask. A two-dimensional array requires treatment_names and gene_names. A DataFrame uses its index as treatment names and columns as gene names unless names are supplied explicitly. A Series must have a two-level (treatment, gene) MultiIndex.

treatment_names, gene_namessequence, optional

Labels for the rows and columns of a two-dimensional mask.

treatment_col, gene_colstr, optional

Columns containing treatment and gene labels in results.

Returns:
alignedndarray of bool, shape (n_tests,)

Mask values in the original row order of results.

Raises:
ValueError

If labels are missing or duplicated, the mask is malformed or non-Boolean, or a result test is absent from the mask.

Notes

This function only aligns a diagnostic mask. It does not modify the multiple-testing family or recompute adjusted p-values.

New in version 0.0.9.

causarray.diagnostics.plot_propensity_scores(A, pi_hat, treatments=None, treatment_names=None, overlap_bounds=(0.05, 0.95), bins=40, max_panels=4, axes=None, clip_bounds=(0.01, 0.99))#

Plot propensity distributions for treatment and control cells.

clip_bounds is forwarded to summarize_propensity_scores() for the returned summary; pass None when pi_hat holds raw, unclipped scores.

Changed in version 0.0.9: clip_bounds is no longer fixed at its default.

causarray.diagnostics.plot_treatment_associations(summary, value='spearman_rho', treatments=None, covariates=None, alpha=None, vmax=None, ax=None)#

Plot a heatmap of treatment-by-covariate associations.

Parameters:
summaryDataFrame

Output from summarize_treatment_associations().

value{‘spearman_rho’, ‘standardized_mean_difference’}, optional

Statistic displayed in the heatmap.

treatments, covariatessequence, optional

Ordered subsets to display.

alphafloat or None, optional

If provided, mark cells whose adjusted p-value is at most alpha. No significance threshold is applied by default.

vmaxfloat or None, optional

Symmetric colour limit, so the heatmap spans (-vmax, vmax). When omitted, 'spearman_rho' uses the fixed range (-1, 1) given by its natural bounds, which keeps panels comparable across subsets, while the unbounded 'standardized_mean_difference' scales to the largest finite magnitude on display.

axmatplotlib Axes, optional

Axes on which to draw the heatmap.

Returns:
fig, axmatplotlib Figure and Axes

The populated figure and axes.

New in version 0.0.9: ..

causarray.diagnostics.summarize_propensity_scores(A, pi_hat, treatment_names=None, overlap_bounds=(0.05, 0.95), clip_bounds=(0.01, 0.99), bins=40)#

Summarize overlap and inverse-weight stability for each treatment.

Other perturbations are excluded from a treatment’s diagnostic comparison; each row compares that treatment with shared all-zero controls.

Parameters:
A, pi_hat, treatment_names, overlap_bounds, bins

Treatments, scores, labels, the reference overlap window, and the histogram resolution used for overlap_ratio.

clip_boundstuple(float, float) or None, optional

The bounds the caller already applied to pi_hat, defaulting to the LFC(ps_clip=(0.01, 0.99)) default. clipped_fraction counts scores sitting on those bounds. Pass None for raw, unclipped scores; clipped_fraction is then NaN rather than a misleading 0.0, since clipping cannot be inferred from the values alone.

.. versionchanged:: 0.0.9

clip_bounds=None reports NaN instead of 0.0.

causarray.diagnostics.summarize_treatment_associations(A, Z, treatment_names=None, covariate_names=None, covariate_types=None, bh_scope='global')#

Summarize covariate association with each treatment.

Each treatment is compared with shared all-zero controls; rows assigned to other treatments are excluded. Spearman p-values are adjusted with the Benjamini-Hochberg procedure over the family selected by bh_scope. The returned statistics are diagnostics and do not automatically identify variables that should be removed from an adjustment set.

Parameters:
Aarray-like, shape (n, a)

Binary treatment indicator matrix. DataFrame column names are retained.

Zarray-like, shape (n, q)

Observed covariates and/or estimated latent factors to diagnose.

treatment_namessequence, optional

Treatment labels. Inferred from a DataFrame or generated when omitted.

covariate_namessequence, optional

Covariate labels. Inferred from a DataFrame or generated when omitted.

covariate_typessequence, optional

Labels such as 'observed' and 'latent'. Defaults to 'observed' for every column.

bh_scope{‘global’, ‘per_treatment’}, optional

Multiple-testing family. 'global' (the default) adjusts across every finite treatment-by-covariate test at once, which with many treatments is a large and correspondingly conservative family. 'per_treatment' adjusts within each treatment’s own block, which is the right choice when each treatment is screened on its own.

Returns:
summaryDataFrame

Long-form table containing spearman_rho, pvalue, BH-adjusted padj, n_tests_in_family, and standardized_mean_difference for every treatment-by-covariate pair. n_tests_in_family records how many tests entered that row’s adjustment, so the correction is self-describing.

New in version 0.0.9: ..