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 |
|---|---|
|
Natural logarithm of the treated/control mean ratio. |
|
Standard error of |
|
Base-2 log fold change, exactly |
|
Standard error of |
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 matrixW.- Parameters:
- Yarray, shape (n, p)
Count matrix of outcomes.
- Warray, shape (n, d)
Covariate matrix, typically
[X | U]whereUare 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,Wis used.- familystr
GLM family for the outcome model:
'nb'(default) or'poisson'.- offsetbool or array-like
Log-scale offset for the outcome model.
Truecomputes size factors automatically;FalseorNonedisables 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
Ktakes precedence.- Kint or None
Number of nuisance-estimation folds.
Noneuses 2 whencross_est=Trueand 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 variances₀²/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; changingusevaralone 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₁)/2to the standard Welch form, which shrinks t-statistics by ≈ √2 relative to v0.0.5. Passusevar='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; passNonefor calibrated probabilities.- **kwargs
Additional arguments forwarded to the GLM fitting functions.
- Returns:
- df_resDataFrame
Test results with natural-log effect estimate
tauand standard errorstd, together with their base-2 equivalentslog2fcandlog2fc_se. The fold change is treatment relative to control, andlog2fc_seis a standard error (not a sample standard deviation). The result also contains inference columns, rawmean_controlandmean_treatedcounterfactual means, and anestimableflag (plustrtfor multiple treatments).New in version 0.0.8: Added the
log2fcandlog2fc_seconvenience columns. The original natural-logtauandstdcolumns remain unchanged.
- causarray.DR_learner.LFC_batch(*args, **kwargs)#
Deprecated alias for
gcate_lfc_batch().Deprecated since version Use:
gcate_lfc_batchinstead.LFC_batchwill 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()orcompute_causal_estimand().- Xarray, shape (n, d)
Covariate matrix.
- id_covsint or array-like of int
Column indices of
Xto compute VIM for. An integerkis treated asrange(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 - 1for 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]). SeeLFC()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,Wis used.- familystr
GLM family for the outcome model:
'nb'(default) or'poisson'.- offsetbool or array-like
Log-scale offset for the outcome model.
Truecomputes size factors automatically;FalseorNonedisables 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
1preserves 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; passNonefor 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, runsfit_gcate_batch()to estimate per-batch latent confounders, then callsLFC()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,Xis used.- batch_sizeint
Perturbations per batch (default 10). Ignored when
n_batchesis set. Batches are sized evenly withnumpy.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_sizeand perturbations are split as evenly as possible across exactlyn_batchesbatches (e.g.n_batches=2on a 29-pert dataset gives two batches of 15 and 14).- max_cellsint or None
Maximum pert cells per batch (default 2 000).
Nonemeans no cap. Ctrl cells are added on top so the actual batch size is at mostn_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 ultimatelyfit_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 defaultusevar='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()andLFC(). When a key collides withgcate_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-levelbackend='original'targeting LFC).
- Returns:
- df_resDataFrame
Concatenated result from all batches. Includes natural-log
tauandstdcolumns, base-2log2fcandlog2fc_secolumns, and a'batch'column with the 0-based batch index. Older compatible caches containing onlytauandstdare 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 useclass_weight='balanced'by default, matchingLFC()and historical causarray fits. Passclass_weight=Nonefor 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.
1fits 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.
Nonereturns 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'matchesLFC(); passNonefor 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_hatunchanged, up to the sharedclipdescribed below. Outcome models are not fit by this function and cachedY_hatvalues can be reused inLFC().- 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_treatmentare 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 bysqrt(factor)during both fitting and prediction and is available only forps_model='logistic'with an L2 penalty. A factor of one leaves the covariate unchanged. Interpret relative penalties on a common scale;prep_causarray_datastandardizes 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(). PassNoneto 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_designflags a treatment whose retained design is constant, andscore_stdreports 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_designFalse but drivesscore_stdtowards 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_namesandgene_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_boundsis forwarded tosummarize_propensity_scores()for the returned summary; passNonewhenpi_hatholds raw, unclipped scores.Changed in version 0.0.9:
clip_boundsis 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 theLFC(ps_clip=(0.01, 0.99))default.clipped_fractioncounts scores sitting on those bounds. PassNonefor raw, unclipped scores;clipped_fractionis thenNaNrather than a misleading0.0, since clipping cannot be inferred from the values alone.- .. versionchanged:: 0.0.9
clip_bounds=NonereportsNaNinstead of0.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-adjustedpadj,n_tests_in_family, andstandardized_mean_differencefor every treatment-by-covariate pair.n_tests_in_familyrecords how many tests entered that row’s adjustment, so the correction is self-describing.
New in version 0.0.9: ..