spacr.ml

Workflow inputs and outputs

Regression

Match plate and well identifiers across phenotype and guide-count inputs, select the response and controls, and inspect diagnostics before interpreting hits. Direct measured responses are also supported.

Open: Home → Regression.

Inputs and outputs below include conditional alternatives. The guidance and handoff notes say which route applies.

Inputs

• Object classification scores — Saved score CSVs and, when merged, measurements/measurements.db, table png_list. Relevant tables, depending on the route: png_list. Relevant columns, depending on the route: pred, cv_predictions, ml_pred, predictions.

• Guide counts per well — Map Barcodes run folder: unique_combinations.csv and annotated_reads.h5. Well identity requires the corresponding barcode references. Relevant columns, depending on the route: count.

• Optical barcode assignments — OPS destination measurements.db: per-well geometry, phenotype alignment, nuclei and barcode tables; optional per-cycle reads. Relevant tables, depending on the route: ops_geometry, ops_phenotype, ops_objects, ops_barcodes, ops_reads.

• Measured objects — measurements/measurements.db; object tables depend on the enabled cell, nucleus, pathogen and organelle masks. Relevant tables, depending on the route: cell, nucleus, pathogen, cytoplasm. Relevant columns, depending on the route: plateID, rowID, columnID, fieldID.

Outputs

• Regression results and hits — Selected run results folder: coefficient/result CSVs, hit tables, settings and diagnostic figures.

Before this module

• Classify: Select the intended CV or ML score column and preserve plate/well identity.

• Map Barcodes: Pair guide counts with phenotype scores using consistent plate/well keys.

• OPS: Join/aggregate decoded objects to phenotype and guide inputs explicitly before Regression; this is not a direct CSV handoff.

After this module

• Run Compare: Compare compatible saved result sets and their settings.

• Hit List: Inspect ranked hits and guide agreement.

• Volcano Explorer: Inspect saved effects and adjusted significance.

• Diagnostics: Open the diagnostics actually written by this run.

• Methods & Results: Review exported prose and every traced result.

• Investigate Hit: Join compatible phenotype and object data for the chosen hit.

API reference.

Module tutorial.

Tabular Machine Learning

Train a feature-based classifier from measured objects and labels. Inspect missing-feature exclusions and grouped held-out performance before using scores.

Open: Classify → Tabular Machine Learning.

Inputs and outputs below include conditional alternatives. The guidance and handoff notes say which route applies.

Inputs

• Measured objects — measurements/measurements.db; object tables depend on the enabled cell, nucleus, pathogen and organelle masks. Relevant tables, depending on the route: cell, nucleus, pathogen, cytoplasm. Relevant columns, depending on the route: plateID, rowID, columnID, fieldID.

• Training annotations — A chosen annotation column in measurements/measurements.db, table png_list; labels belong to object identities. Relevant tables, depending on the route: png_list. Relevant columns, depending on the route: prcfo.

Outputs

• Object classification scores — Saved score CSVs and, when merged, measurements/measurements.db, table png_list. Relevant tables, depending on the route: png_list. Relevant columns, depending on the route: pred, cv_predictions, ml_pred, predictions.

• Fitted feature classifier — The fitted tabular classifier and its recorded feature list, training settings and validation results.

• Classifier evaluation bundle — Held-out predictions, labels, split metadata and calibration/leakage metrics for a saved classifier run.

API reference.

Module tutorial.

Find screen phenotypes and test which guides or genes explain them.

WHAT IT IS FOR

The Regression tile opens this module because perform_regression() is its main entry point: it joins per-well image scores to sequencing counts and estimates guide- or gene-level associations for a pooled screen. The module also contains a separate classical machine-learning workflow, generate_ml_scores(), which trains a classifier on measured single-object features and turns those predictions into the score table a regression can consume.

WHAT IT NEEDS

A regression run is configured with paired_data: ordered score/count CSV pairs whose plate and well identities agree. Older score_data and count_data lists are migrated positionally, but explicit pairs are safer. Choose the score column with dependent_variable and provide the relevant control wells, plate metadata, analysis level (guide, gene, or both), multiple-testing threshold, and either a supported regression_type or None for distribution-based selection. Family-specific settings are validated rather than silently ignored. The classical-ML path instead needs one or more measurements.db files, labelled positive and negative controls or an annotation column, and a model choice such as XGBoost, logistic regression, or random forest.

WHAT IT PRODUCES

Regression results go into a new, non-overwriting <output root>/results/<analysis kind>[_n] directory. results.csv is the primary combined table; results_grna.csv and results_gene.csv make the fitted levels explicit, and results_significant.csv records thresholded hits. The same directory holds summaries, diagnostics, volcano and plate figures, optional publication-panel packages, resource measurements, or a detailed failure report. A requested level with no fitted rows is kept as a header-only CSV so downstream tools can distinguish “tested, no rows” from a missing artifact. Classical ML writes predictions, feature-importance tables, evaluation results, and a plate heatmap beneath results.

WHAT TO DO NEXT

Read the diagnostic and failure/resource records before ranking hits, then review the significant table alongside the complete level tables and check whether guide and gene effects agree. Open the generated result panels for visual QC and retain the settings/manifests with any reported hit list. If scores do not yet exist, run generate_ml_scores(); if counts do not yet exist, create them with spacr.sequencing.generate_barecode_mapping().

Several statistical distinctions are deliberate. Guide and gene fits are separate multiple-testing families and receive separate corrections; the nominal alpha, effect-size threshold, and corrected significance cutoff are not interchangeable. A mixed model reports guide effects as shrunken BLUP predictions without guide p- or q-values, so they must not be read as a second guide significance test. Diagnostic or report-generation failures do not erase successful scientific output, while an actual regression failure is recorded and then re-raised unchanged so callers cannot mistake it for a completed run.

Classes

QuasiBinomial

Binomial GLM family scaled by a dispersion parameter (quasi-binomial).

Functions

apply_transformation(X, transform)	Return an sklearn FunctionTransformer for the named transform.
beta_logit(values)	A proportion on the logit scale, with the endpoints squeezed in.
binarise_response(y[, threshold, name])	Return y as a 0/1 vector for a classifier backend, refusing to guess.
calculate_p_values(X, y, model)	Return OLS-style p-values for a fitted model's coefficients.
centre_on_controls(df, dependent_variable, nc)	Subtract the negative controls' median response. Returns (df, offset).
check_and_clean_data(df, dependent_variable)	Prepare the merged count / score frame for model fitting.
check_distribution(y[, epsilon])	Check the distribution of y and recommend a regression type.
check_normality(data, variable_name[, verbose])	Check if the data is normally distributed using the Shapiro-Wilk test.
clean_controls(df, values, column)	Drop rows whose column holds one of the listed values.
create_volcano_filename(csv_path, regression_type, ...)	Build the path this run's volcano plot will be saved to.
display(*args, **kwargs)	Accept and discard display arguments when IPython is unavailable.
double_transform_warning(→ str)	Describe a response transform that is compounded by the model link.
find_optimal_threshold(y_true, y_pred_proba)	Return the probability threshold maximising F1 on the precision-recall curve.
fit_mixed_model(df, formula, dst, *[, ...])	Fit a mixed model with guides nested within genes.
fit_quality_note(→ str)	Return a one-line goodness-of-fit summary for a fitted GLM.
generate_ml_scores(settings)	Train a classical ML classifier (XGBoost / logistic / RF) on per-object features and score every well of a screen.
interpret_vision_model([settings])	Explain a spacr vision-model score using RF, permutation and SHAP importance, with per-compartment / per-channel radar plots.
label_control_condition(features, guides[, nc, pc, ...])	Label every coefficient row 'nc', 'pc', 'control' or 'other'.
load_regression_input_pairs(pairs)	Read paired inputs and resolve plate identity without filename guesses.
mcfadden_note(→ str)	Format McFadden's pseudo-R² and flag a negative value.
minimum_cell_simulation(settings[, num_repeats, ...])	Estimate the minimum number of cells per well needed for a stable well mean.
ml_analysis(df[, channel_of_interest, ...])	Train a per-object classifier on positive/negative control wells and score every row of the input DataFrame.
normalize_regression_input_pairs(settings)	Return explicit score/count rows, migrating legacy lists.
perform_mixed_model(y, X, groups[, alpha, ...])	Fit a mixed-effects linear model with groups as the random intercept.
perform_regression(settings)	Run the regression and report actionable details if it fails.
pick_glm_family_and_link(y[, name, transform])	Select the GLM family and link that suit the response.
prepare_formula(dependent_variable[, ...])	Build a fixed-effects formula for one screen-analysis level.
process_model_coefficients(model, regression_type, X, ...)	Return a DataFrame of model coefficients and p-values, one row per term.
process_reads(csv_path, fraction_threshold, plate[, ...])	Load a per-gRNA read-count CSV and return per-well normalised fractions.
process_scores(df, dependent_variable, plate[, ...])	Aggregate per-object model scores to per-well summaries, ready for regression.
regression(df, csv_path[, dependent_variable, ...])	Run the full regression pipeline: clean, fit, extract coefficients, optional volcano plot.
regression_levels(df, csv_path[, dependent_variable, ...])	Fit every level the run asked for, SEPARATELY, and return one per level.
regression_model(X, y[, regression_type, groups, ...])	Dispatch to the requested regression backend and return the fitted model.
resolve_auto_inference(data, settings, *[, ...])	Choose analysis_mode for inference='auto' from the design.
resolve_glm_transform_conflict(dependent_variable[, ...])	Resolve a transform/family-link conflict before fitting a GLM.
resolve_levels(regression_type[, level])	Which level(s) a run fits, given the backend and the level setting.
resolve_regression_src(requested, automatic)	Resolve the root directory used for regression output.
results_folder_kind(→ str)	What a run's results folder is NAMED after.
save_summary_to_file(model[, file_path])	Write model.summary().as_text() to file_path as plain text.
scale_variables(X, y)	Min-max scale the independent (X) and dependent (y) variables to [0, 1].
screen_is_blockable(→ bool)	Whether screenID can be a term in this frame's design.
select_glm_family(y)	Choose a statsmodels GLM family from the range and type of the response.
shap_analysis(model, X_train, X_test)	Build a SHAP summary beeswarm for X_test.
summary_for_console(→ str)	Return a statsmodels summary sized for terminal output.
write_plot(plot, path[, title])	Write a pyqtgraph plot out and announce it, like publish does.

Module Contents

class spacr.ml.QuasiBinomial(link=Logit(), dispersion=1.0)

Bases: statsmodels.genmod.families.Binomial

Binomial GLM family scaled by a dispersion parameter (quasi-binomial).

Parameters:

• link – statsmodels link instance. Default Logit().

• dispersion – Multiplicative variance scaling. Default 1.0.

Store the dispersion factor after delegating to Binomial.

spacr.ml.apply_transformation(X, transform)

Return an sklearn FunctionTransformer for the named transform.

Parameters:

• X – Ignored (kept for compatibility with sklearn pipeline flow).

• transform – One of 'log', 'sqrt', 'square', 'beta'. Any other value returns None.

Returns:

A FunctionTransformer or None.

spacr.ml.beta_logit(values)

A proportion on the logit scale, with the endpoints squeezed in.

transform='beta' is intended for proportional responses such as classification scores and their well aggregates, where a logarithm is not appropriate.

This is distinct from regression_type='beta', which selects a beta GLM. One transforms the response; the other selects the model family.

Parameters:

values – proportions in [0, 1], array-like; converted to a float array. Non-finite entries pass through unchanged. When any finite value is at or beyond 0 or 1 the finite values are squeezed with (y * (n - 1) + 0.5) / n before the logit.

spacr.ml.binarise_response(y, threshold=None, name='response')

Return y as a 0/1 vector for a classifier backend, refusing to guess.

The hinge backend fits a decision boundary, so it needs two classes. There are exactly two ways to get them and this function will not invent a third:

• y already holds exactly two distinct finite values (the usual case: a per-object class call aggregated to a well, or a 0/1 score). The lower value becomes 0 and the higher becomes 1, so the sign of every coefficient answers “does this gRNA push wells towards the HIGHER class”, which is the same direction the continuous models report.

• threshold is given explicitly, and y > threshold becomes 1.

A continuous response with no threshold is REFUSED. Picking a cut for the user — the mean, the median, 0.5 — would silently redefine the hypothesis being tested: on a screen whose well scores run 0.2-0.8 a median split calls half the plate positive by construction, and the resulting hit list is a plausible, unfalsifiable artefact of the split.

Parameters:

• y – Response vector (array, Series or single-column frame).

• threshold – Explicit cut; values strictly greater become 1.

• name – Name used in error messages, for a legible failure.

Returns:

numpy float array of 0.0/1.0, same length as y.

Raises:

ValueError – if y is continuous and no threshold is given, if a given threshold puts every observation in one class, or if y holds fewer than two distinct values.

Example

binarise_response([0, 1, 1, 0])            # -> [0., 1., 1., 0.]
binarise_response([2, 5, 5], )             # -> [0., 1., 1.]
binarise_response([0.2, 0.6], threshold=0.4)   # -> [0., 1.]

spacr.ml.calculate_p_values(X, y, model)

Return OLS-style p-values for a fitted model’s coefficients.

These are not valid frequentist p-values for a penalised fit, and the two callers that reach them know it in different ways. The standard error is the unpenalised rse * sqrt(diag((X'X)^-1)) while the coefficient it is divided into has been shrunk, so the test is mis-specified. The direction of the error is the one that matters here and it is the safe one: the penalty shrinks the numerator and inflates the residual in the denominator, so the statistic is too SMALL and the p-value too large. A penalised fit under-detects here; it does not manufacture hits.

lasso and elasticnet do not rely on this at all — NO_P_VALUE_TYPES routes them to a bootstrap selection frequency instead. ridge does, because it never sets a coefficient to exactly zero and so has no selection frequency to report (every feature would score 1.0), and a conservative test is a better answer than no test. tests/test_regression_orientation.py pins the null case, which is where an anticonservative version of this would show.

Parameters:

• X – Design matrix (n x p).

• y – Observed responses.

• model – Fitted estimator exposing predict and coef_.

Returns:

1D array of length p; entries are NaN when n <= p + 1.

spacr.ml.centre_on_controls(df, dependent_variable, nc)

Subtract the negative controls’ median response. Returns (df, offset).

THIS IS WHAT MAKES THE INTERCEPT MEAN SOMETHING. A fitted intercept is the response where every predictor is zero, which on a screen design is a well with no guide in it – a point that does not exist. Centred on the negative controls, the intercept is the control level, and every coefficient reads directly as “this far above or below the controls”.

The offset is returned rather than swallowed so the caller can report it: a coefficient table whose response was shifted, with nothing saying by how much, is a table nobody can compare with another run.

Parameters:

• df – the long frame the fit runs on.

• dependent_variable – the response column.

• nc – the negative-control guide or gene, as the settings name it.

Returns:

(frame, offset). The frame is a copy when it was changed and the original when it was not; offset is 0.0 when no control row could be identified, and the caller is expected to say so.

spacr.ml.check_and_clean_data(df, dependent_variable)

Prepare the merged count / score frame for model fitting.

Drops rows with a missing fraction or dependent variable, casts the identifier columns to categorical and reports (without dropping) collinear columns via VIF. The returned frame keeps only fraction, the dependent variable, gene, grna, prc, plateID, rowID, columnID, and cell_count and screenID when present, plus a computed gene_fraction column: the sum of the gene’s gRNA fractions within each well, which the regression formula regresses on.

Parameters:

• df – Merged DataFrame of counts and scores.

• dependent_variable – Name of the response column.

Returns:

The cleaned DataFrame used as the model input.

Raises:

ValueError – if a (prc, grna) pair carries more than one fraction, which makes gene_fraction ambiguous.

spacr.ml.check_distribution(y, epsilon=1e-06)

Check the distribution of y and recommend a regression type.

Parameters:

• y – Response vector.

• epsilon – How close to 0 or 1 a value may sit before it counts as a boundary case. Default 1e-6.

Returns:

One of 'logit', 'quasi_binomial', 'beta', 'ols' or 'glm', as accepted by regression()’s regression_type.

spacr.ml.check_normality(data, variable_name, verbose=False)

Check if the data is normally distributed using the Shapiro-Wilk test.

Parameters:

• data – numeric values, array-like; non-finite values are dropped and fewer than 3 remaining values returns False without testing.

• variable_name – name printed in the verbose messages only.

• verbose – print the test statistic, P value and verdict.

Returns:

True when the Shapiro-Wilk P value exceeds 0.05.

spacr.ml.clean_controls(df, values, column)

Drop rows whose column holds one of the listed values.

Parameters:

• df – Source DataFrame.

• values – List of values to remove. Anything that is not a list (a bare value included) is a no-op.

• column – Column, or list of columns, to check. None is a no-op.

Returns:

Filtered DataFrame (unchanged if column is missing or values is not a list).

spacr.ml.create_volcano_filename(csv_path, regression_type, alpha, dst)

Build the path this run’s volcano plot will be saved to.

Path construction only: nothing is read, written or created, and the .pdf in the name is not binding - spacr.plot.save_figure() rewrites the extension to whichever format the figure preference selected.

Parameters:

• csv_path – Source CSV. Only its basename with the last extension stripped becomes the <name>_volcano_plot.pdf stem, and only its directory is used, when dst is falsy. The file is never opened, so a path that does not exist is fine; a bare filename yields a bare relative result rather than a path under the working directory.

• regression_type – Prefixed to the filename, unless it is exactly 'quantile' - then alpha is prefixed instead. None is stamped literally, giving None_...: regression() calls this before check_distribution() resolves the auto-selected model, so an auto run’s plot is never named for the model it actually fitted.

• alpha – Read only on the 'quantile' branch; accepted and ignored for every other type, whatever its value. regression() passes the quantile setting here, not the penalty, so two quantiles of one screen cannot overwrite each other.

• dst – Output directory. Any falsy value, None and '' alike, falls back to the directory of csv_path. It is not created here.

Returns:

The joined path, which regression() hands to spacr.plot.volcano_plot() as save_path.

spacr.ml.display(*args, **kwargs)

Accept and discard display arguments when IPython is unavailable.

spacr.ml.double_transform_warning(name, transform, family) → str

Describe a response transform that is compounded by the model link.

For example, a log-transformed response passed to a family with a logit link fits logit(log(y)). The function returns an actionable warning before fitting; an identity link or a response without a link-like transform returns an empty string.

Parameters:

• name – Response name shown in the warning.

• transform – Transform already applied to the response.

• family – Statsmodels family whose link will be inspected.

Returns:

Warning text, or "" when the transforms do not compound.

spacr.ml.find_optimal_threshold(y_true, y_pred_proba)

Return the probability threshold maximising F1 on the precision-recall curve.

Parameters:

• y_true – Ground-truth binary labels.

• y_pred_proba – Predicted probabilities for the positive class.

Returns:

Optimal probability threshold.

spacr.ml.fit_mixed_model(df, formula, dst, *, random_row_column_effects=False, gene_column='gene', guide_column='grna', regression_backend=DEFAULT_REGRESSION_BACKEND)

Fit a mixed model with guides nested within genes.

The model treats genes as fixed effects and guides as random effects nested within genes. In statsmodels notation, groups=gene supplies the outer random intercept and vc_formula={'grna': '0 + C(grna)'} supplies the guide-within-gene variance component.

A blockable screenID supplied by prepare_formula() remains a fixed effect. With only two screen levels, a random screen variance would be estimated from one degree of freedom. The plate is not nested within the screen because plate position is already represented by the row and column structure. Single-screen data omit the constant screen term to avoid a rank-deficient design.

Parameters:

• df (pandas.DataFrame) – Model data containing the formula variables and the gene and guide grouping columns.

• formula (str) – Fixed-effects formula, normally returned by prepare_formula() with level='gene'.

• dst (path-like) – Destination for the residual histogram.

• random_row_column_effects (bool, default False) – Add row and column variance components instead of fixed terms.

• gene_column (str, default 'gene') – Column containing the outer gene groups.

• guide_column (str, default 'grna') – Column containing guides nested within each gene.

• regression_backend ({'statsmodels', 'torch'}, default 'statsmodels') – Mixed-model backend. The torch backend fits the same nested model with GPU acceleration when available.

Returns:

• mixed_model – Fitted backend-specific mixed-model result.

• coef_df (pandas.DataFrame) – Fixed effects, variance components, and guide BLUPs. Variance components and BLUPs have NaN p-values because they are not fixed-effect hypothesis tests.

Raises:

• ValueError – If required grouping columns are missing, no gene has multiple guides, or the backend cannot fit the nested design.

• MixedBackendUnavailable – If the selected mixed-model backend is unavailable.

spacr.ml.fit_quality_note(model) → str

Return a one-line goodness-of-fit summary for a fitted GLM.

McFadden’s pseudo-R-squared compares log-likelihoods, and it is the appropriate summary for a GLM with a discrete response. A Gaussian identity-link fit instead reports ordinary R-squared because its likelihood is a density and the McFadden ratio is not interpretable on the usual zero-to-one scale.

Parameters:

model – a fitted statsmodels GLM result.

Returns:

A labelled goodness-of-fit line for the console.

spacr.ml.generate_ml_scores(settings)

Train a classical ML classifier (XGBoost / logistic / RF) on per-object features and score every well of a screen.

Reads the measurements.db produced by spacr.measure.measure_crop(), merges cell/nucleus/pathogen/ cytoplasm feature tables, uses the wells marked as positive_control / negative_control (or an annotation column) as training labels, delegates fitting to ml_analysis(), and writes per-object predictions, permutation and feature-importance tables plus a plate heatmap into results/ next to the source DB.

Parameters:

settings –

Settings dict, canonicalized via spacr.settings.set_default_analyze_screen(). Key entries:

• src (str or list) — folder(s) containing measurements/measurements.db.

• channel_of_interest — 0-based channel for the recruitment ratio feature; also drives table selection.

• model_type_ml — 'xgboost', 'logistic_regression', 'random_forest'.

• positive_control / negative_control — well IDs (e.g. 'c2' / 'c1') used as training labels.

• annotation_column — override controls with a PNG-level annotation column.

• location_column — 'columnID' or 'rowID'.

• heatmap_feature — feature plotted on the plate heatmap.

• exclude, n_repeats, top_features, test_size, reg_alpha, reg_lambda, learning_rate, n_estimators, n_jobs.

• remove_low_variance_features, remove_highly_correlated_features, prune_features, cross_validation, verbose.

Returns:

The two-element list [output, plate_heatmap], where output is the 10-element result list of ml_analysis() and plate_heatmap is the plate-heatmap matplotlib figure. The CSVs and figures are written to results/ as a side effect; their paths are not returned.

Raises:

ValueError – if annotation_column is set but the png_list table lacks prcfo / that column, its object IDs do not join to the measurements, it contains fewer than two observed classes, or if heatmap_feature is not among the trained features.

Example

from spacr.ml import generate_ml_scores
settings = {
    'src': '/data/plate01',
    'channel_of_interest': 3,
    'positive_control_id': 'c2', 'negative_control_id': 'c1',
    'model_type_ml': 'xgboost', 'heatmap_feature': 'recruitment',
}
generate_ml_scores(settings)

See also

ml_analysis() — the underlying fit/evaluate routine. perform_regression() — mixed-effects regression on per-well ML scores.

spacr.ml.interpret_vision_model(settings=None)

Explain a spacr vision-model score using RF, permutation and SHAP importance, with per-compartment / per-channel radar plots.

Merges per-object measurements from measurements.db with a CSV of predicted scores, runs any combination of RF feature importance, permutation importance and SHAP over the top features, then aggregates SHAP contributions into compartment and channel radar plots so you can see which region (cell / nucleus / pathogen / cytoplasm) and which fluorescence channel drives the model.

Parameters:

settings –

Settings dict, canonicalized via spacr.settings.set_interpret_vision_model_defaults(). Key entries:

• src — folder containing measurements/measurements.db.

• scores — CSV of per-object predictions to explain.

• score_column — column of scores holding the score.

• tables — DB tables to merge (default ['cell','nucleus','pathogen','cytoplasm']).

• feature_importance / permutation_importance / shap — enable each explainer.

• top_features — cap on features shown.

• nuclei_limit / pathogen_limit — object-count caps.

• n_jobs, save.

Returns:

The merged per-object DataFrame — the measurement tables joined to the scores CSV — that the explainers were fitted on. Radar and importance plots are rendered, and with save=True importance CSVs are written alongside the DB, as side effects.

Example

from spacr.ml import interpret_vision_model
interpret_vision_model({
    'src': '/data/plate01',
    'scores': '/data/plate01/results/pred.csv',
    'score_column': 'pred',
    'shap': True, 'top_features': 30,
})

See also

spacr.submodules.interpret_vision_model() — legacy / alternative entry point returning a dict of importance DataFrames instead of the merged measurements.

spacr.ml.label_control_condition(features, guides, nc=None, pc=None, controls=None, *, strict: bool = False, verbose: bool = False)

Label every coefficient row 'nc', 'pc', 'control' or 'other'.

The condition column: what the volcano colours by, what the results panel offers in “colour by”, and – the reason this is a function rather than four lines inside process_model_coefficients() – what the EFFECT-SIZE CUT measures its spread on. A coefficient table without it is a table spacr.qt.widgets.regression_results.RegressionResultsPanel. set_threshold_method() answers “No control coefficients, so no effect-size cut” for, which is what every guide-permutation run used to get.

Precedence is nc, then pc, then the explicit controls list, so a guide named in two of them is reported once and always the same way.

Parameters:

• features – the model term per row, e.g. fraction:grna[000000_1]. nc and pc are matched as SUBSTRINGS of it, which is how a negative control given as a gene id reaches a term named for a guide.

• guides – the guide identifier per row, matched whole against controls. Both sides are compared as text, so a control list that round-tripped through a settings CSV as integers still matches.

• nc – negative-control identifier, or None for no negative control.

• pc – positive-control identifier, or None.

• controls – non-targeting guide identifiers, or None. None means “no control list” and labels nothing – it is the value perform_regression() documents for a control-free screen, and the inline version this replaced raised TypeError on it.

• strict – raise spacr.control_names.ControlNotFound when a NAMED nc or pc matches nothing. Off by default so a call on a partial frame is not an error; the run turns it on, because there a control matching nothing is a number computed against an empty set.

• verbose – print what each control resolved to and how much it matched.

Returns:

a pandas.Series of labels aligned with features.

spacr.ml.load_regression_input_pairs(pairs)

Read paired inputs and resolve plate identity without filename guesses.

Resolution order is own column, partner column, then pair-row order. Conflicting declarations are refused. Returns (count_frame, score_frame, audit_rows).

Parameters:

pairs – sequence of mappings with 'score' and 'count' table paths (either may be empty), as returned by normalize_regression_input_pairs(). Each mapping’s 'plate' is overwritten with the resolved plate label.

spacr.ml.mcfadden_note(r2) → str

Format McFadden’s pseudo-R² and flag a negative value.

A negative value means the fitted model predicts the response worse than an intercept-only model. The returned note explains that the coefficients should not be interpreted and points to a common cause: applying a response transform that duplicates the fitted family’s link.

Parameters:

r2 – Pseudo-R² value, or a value convertible to float.

Returns:

One-line diagnostic text suitable for a console or report.

spacr.ml.minimum_cell_simulation(settings, num_repeats=10, sample_size=100, tolerance=0.02, smoothing=10, increment=10, dst=None)

Estimate the minimum number of cells per well needed for a stable well mean.

For the wells with the most objects, repeatedly subsamples cells at increasing sample sizes and records the mean absolute difference from the well’s full mean. Plots the smoothed curve with a ±1 s.d. band, marks the elbow point (or settings['min_cells_per_well'] when it is set) and writes cell_min_threshold.pdf into dst.

Pass dst to keep the figure in a specific run folder. When omitted, the function uses the screen-level results folder derived from count_data for compatibility with direct notebook and script calls.

Parameters:

• settings – Requires score_data (CSV path or list of paths), dependent_variable, tolerance (int percent or float fraction) and min_cells_per_well. count_data is needed only when dst is left unset, and only to locate the figure.

• num_repeats – Subsamples drawn per sample size. Default 10.

• sample_size – Number of wells, taken largest-first by cell count, to simulate. Default 100.

• tolerance – Unused; the tolerance applied is settings['tolerance'].

• smoothing – Rolling-window width used to smooth the curve.

• increment – Step between the simulated sample sizes.

• dst – Folder for cell_min_threshold.pdf, created if missing. Default None: <folder of settings['count_data'][0]>/results.

Returns:

The elbow point’s sample size, i.e. the minimum cell count per well, for passing to process_scores().

Raises:

ValueError – if settings['tolerance'] is neither an int nor a float.

spacr.ml.ml_analysis(df, channel_of_interest=3, location_column='columnID', positive_control='c2', negative_control='c1', exclude=None, n_repeats=10, top_features=30, reg_alpha=0.1, reg_lambda=1.0, learning_rate=1e-05, n_estimators=1000, test_size=0.2, model_type='xgboost', n_jobs=-1, remove_low_variance_features=True, remove_highly_correlated_features=True, prune_features=False, cross_validation=False, verbose=False, *, split_by='well', holdout_plate=None, batch_correction='none', batch_column='plateID', batch_control_column=None, batch_control_values=None, batch_covariate_column=None, batch_combat_mean_only=False, batch_min_samples=3, batch_missing_control='error')

Train a per-object classifier on positive/negative control wells and score every row of the input DataFrame.

Called directly for one-off ML work, and internally by generate_ml_scores(). Filters features by channel, drops low-variance and highly correlated columns, splits (or CVs) train / test, fits the requested model, computes permutation and native feature importances, tunes an optimal decision threshold and writes predictions + probabilities back onto the returned DataFrame.

Parameters:

• df – Per-object feature DataFrame as produced by merging the cell/nucleus/pathogen/cytoplasm tables of a spacr.measure.measure_crop() database.

• channel_of_interest – Channel index used to select features.

• location_column – Column identifying wells / plate columns. Default 'columnID'.

• positive_control – Value(s) in location_column treated as the positive class. Default 'c2'.

• negative_control – Value(s) treated as the negative class. Default 'c1'.

• exclude – Columns to remove from feature space.

• n_repeats – Repeats for permutation importance. Default 10.

• top_features – Feature cap when prune_features=True.

• reg_alpha – XGBoost L1 penalty.

• reg_lambda – XGBoost L2 penalty.

• learning_rate – XGBoost learning rate.

• n_estimators – Tree count for tree-based models.

• test_size – Test-split fraction. Default 0.2.

• model_type – 'random_forest', 'logistic_regression', 'gradient_boosting' or 'xgboost'.

• n_jobs – Parallel job count where applicable. Default -1.

• remove_low_variance_features – Drop low-variance features.

• remove_highly_correlated_features – Drop highly correlated features.

• prune_features – If True, apply SelectKBest before training.

• cross_validation – If True, run 5-fold stratified CV.

• verbose – Log progress details.

• split_by – Independent acquisition unit for train/test splitting: 'cell', 'field', 'well' (default), or 'plate'. Legacy 'none' is an alias for 'cell'.

• batch_correction – plate correction method from spacr.batch_correction.

• batch_column – metadata column identifying plates/batches.

• batch_control_column – metadata column holding reference-control labels for control_center.

• batch_control_values – negative/reference control value(s).

• batch_min_samples – minimum rows or controls per plate.

• batch_covariate_column – Metadata column containing a biological covariate that ComBat must preserve, such as treatment, cell line, or time point. Required when batch_correction="combat"; its coefficients remain in the corrected data while estimated batch effects are removed.

• batch_combat_mean_only – If True, ComBat adjusts batch means without scaling batch variances. This can be appropriate when batches differ primarily by location or contain too few observations for stable variance estimates. Default False adjusts both means and variances.

• batch_missing_control – error or skip for missing controls.

Returns:

Tuple (output, figs) where output is a positional tuple of (scored_df, permutation_df, feature_importance_df, model, X_train, X_test, y_train, y_test, metrics_df, train_features) and figs is (permutation_fig, feature_importance_fig).

Raises:

ValueError – on unsupported model_type or when positive / negative control rows cannot be located in location_column.

Example

from spacr.ml import ml_analysis
output, figs = ml_analysis(
    df, channel_of_interest=3,
    positive_control='c2', negative_control='c1',
    model_type='xgboost',
)
scored_df = output[0]

See also

generate_ml_scores() — wraps this call with DB I/O.

spacr.ml.normalize_regression_input_pairs(settings)

Return explicit score/count rows, migrating legacy lists.

New settings store paired_data. Older files remain valid: their flat lists are zipped positionally, exactly matching the former behaviour, and the migration is reported so the invisible legacy assumption is visible.

Parameters:

settings – regression settings dictionary. paired_data is read when present, otherwise the legacy score_data and count_data lists; the dictionary is updated in place with the normalised paired_data and de-duplicated score_data/count_data lists.

Returns:

(pairs, migrated), where migrated is True when the rows came from the legacy lists.

Raises:

ValueError – when paired_data is malformed or there is not at least one score path and one count path.

spacr.ml.perform_mixed_model(y, X, groups, alpha=None, regression_backend=DEFAULT_REGRESSION_BACKEND)

Fit a mixed-effects linear model with groups as the random intercept.

Collinearity is REPORTED, never silently corrected. The previous revision reacted to any VIF above 10 by fitting

ridge = Ridge(alpha=alpha).fit(X, y)
X_ridge = ridge.coef_ * X          # "Adjust X with Ridge coefficients"
MixedLM(y, X_ridge, groups=groups)

which is not ridge regression and not a mixed model of anything. It multiplies every column by that column’s ridge coefficient, so

• a column whose ridge coefficient is 0 - which is most of them on a screen-scale one-hot design - becomes a column of zeros, and the design is singular. That is the numpy.linalg.LinAlgError: Singular matrix that regression_type='mixed' died with on real data, thrown from inside statsmodels with nothing naming the cause;

• where it did fit, every fixed effect came back multiplied by an arbitrary per-column constant, so the coefficients written to results.csv and ranked on the volcano plot were not effects on the response at all. That is the worse of the two outcomes, because it completes.

A one-hot design against an intercept ALWAYS trips VIF > 10, so this path was the normal one, not the exception.

Parameters:

• y – Response vector.

• X – Fixed-effects design matrix (DataFrame).

• groups – Cluster identifiers for the random intercept - one entry per row of X.

• alpha – Must be None. Accepted only so an old call site fails with an explanation instead of a TypeError.

• regression_backend – WHO fits it. 'statsmodels' is the default and produced every existing result; 'torch' fits the same profiled REML objective on the GPU (spacr.mixed_gpu) and returns a result object with the same attributes, so nothing downstream can tell which ran except by asking. A backend that cannot fit 'mixed' here, is not installed, or needs a GPU that is absent is REFUSED with the reason – see _require_backend().

Returns:

Fitted statsmodels MixedLMResults, or the equivalent spacr.mixed_gpu.TorchMixedResults.

Raises:

ValueError – if groups is None, if alpha is given, if groups does not align with X, or if the fixed-effects design is rank-deficient (which MixedLM would otherwise report as a bare LinAlgError from three frames deep).

spacr.ml.perform_regression(settings)

Run the regression and report actionable details if it fails.

On failure, the original exception is re-raised unchanged after a report is printed and written to the run folder. The report includes the most recent stage stored in settings['_regression_stage'], available design dimensions, and a remedy for recognized failures.

Parameters:

settings – Regression settings consumed by the fitting pipeline.

Returns:

Regression output mapping. model_data is the prepared input table, not coefficient results; fit_designs records measured design counts separately for each parametric fit level. Unrecorded counts are omitted, and permutation outputs retain their own result schema.

Raises:

Exception – Re-raises the original regression failure.

spacr.ml.pick_glm_family_and_link(y, name='', transform='')

Select the GLM family and link that suit the response.

Used by regression_type='glm' to choose a family from the data rather than from the user.

Parameters:

• y – Response vector.

• name – Response-column name printed with the selected family. The name makes clear whether the family was chosen from a derived scale.

• transform – the transform already applied, printed with the name and checked for the double transform of double_transform_warning().

Returns:

A statsmodels family instance with its link set.

Raises:

ValueError – only through _validate_poisson_response(), when the response looks like counts but cannot be one.

spacr.ml.prepare_formula(dependent_variable, random_row_column_effects=False, block_screen=False, level='grna', model_plate_position=True, intercept='fitted')

Build a fixed-effects formula for one screen-analysis level.

Parameters:

• dependent_variable (str) – Name of the response column.

• random_row_column_effects (bool, default=False) – Reserve plateID, rowID and columnID for the grouping and variance-component structure in fit_mixed_model() instead of adding them as fixed effects.

• block_screen (bool, default=False) – Add screenID as a fixed effect. Use screen_is_blockable() before enabling this for user data.

• level ({'grna', 'gene'}, default='grna') – Resolution represented by the formula. 'grna' uses fraction:grna and 'gene' uses gene_fraction:gene.

• intercept ({'fitted', 'zero', 'control', 'value'}, default='fitted') – What the intercept is. 'fitted' estimates it. 'zero' takes it out of the design, so the fit passes through the origin and a coefficient is a whole predicted score rather than a departure from a baseline. 'control' keeps the term and is completed by the caller, which centres the response on the negative controls first – the intercept is then the control level by construction. 'value' suppresses the term as well, because the caller has shifted the response by a number the user gave and the intercept is pinned at exactly that number.

• model_plate_position (bool, default=True) – Include plate position in the model. With random_row_column_effects=False it is included as fixed plate, row, and column terms; with random_row_column_effects=True plate supplies the grouping variable and row/column are variance components. New application settings default to False even though this helper retains True for API compatibility.

Returns:

str – A patsy-compatible formula for one analysis level.

Raises:

ValueError – If level is unknown or 'both', or if random plate-position effects are requested while plate position is disabled.

Notes

Guide and gene effects are fitted separately because the gene fraction is derived from its guide fractions; including both blocks in one design is rank deficient. Use regression_levels() to request both fits.

spacr.ml.process_model_coefficients(model, regression_type, X, y, nc, pc, controls, hinge_threshold=None, hinge_n_boot=200)

Return a DataFrame of model coefficients and p-values, one row per term.

Every name in REGRESSION_TYPES has a branch here. It is the same table for all of them - feature, coefficient, p_value, -log10(p_value), grna, condition - because everything downstream (the volcano plot, the hit table, the metadata merge) reads those columns and nothing else.

Parameters:

• model – The fitted object from regression_model().

• regression_type – Which backend produced it.

• X – Design matrix, used for the sklearn feature names and for the p-value approximations that need the data back.

• y – Response, likewise.

• nc – Negative-control identifier, matched against the feature name.

• pc – Positive-control identifier.

• controls – Explicit list of control gRNA identifiers.

• hinge_threshold – The binarisation cut used by the hinge fit; the bootstrap below must reproduce the SAME two classes the fit saw.

• hinge_n_boot – Bootstrap resamples used for the hinge p-values.

Returns:

Coefficient DataFrame with the row/column nuisance terms removed.

Raises:

ValueError – on an unsupported regression_type.

spacr.ml.process_reads(csv_path, fraction_threshold, plate, filter_column=None, filter_value=None, record=None, exclude_grnas=None)

Load a per-gRNA read-count CSV and return per-well normalised fractions.

Splits derived plate_row or prcfo identifiers, computes each gRNA’s fraction of the well total, applies an optional fraction-cutoff filter and returns a compact (prc, grna, fraction) frame (with gene derived from the gRNA when possible).

Parameters:

• csv_path – Path to the counts CSV, or an already-loaded DataFrame.

• fraction_threshold – Drop rows below this fraction; must be in [0, 1] or None.

• plate – Plate identifier used when no plateID column is present.

• filter_column – Column (or list of columns) to filter rows on.

• filter_value – Values (or list of values) to drop from filter_column.

• record – Optional mutable mapping that records exclusions for the persisted regression summary.

• exclude_grnas – Guide or gene identifiers to remove from the raw count table. Gene identifiers remove all associated guides. This is applied before well totals and fractions are calculated, so retained guides are normalised against the retained read count.

Returns:

DataFrame with columns prc, grna, fraction.

Raises:

ValueError – on missing required columns, invalid fraction_threshold, or when the threshold removes all rows.

spacr.ml.process_scores(df, dependent_variable, plate, min_cells_per_well=25, agg_type='mean', transform=None, regression_type='ols', invert_dependent_variable=False)

Aggregate per-object model scores to per-well summaries, ready for regression.

Ensures plateID/rowID/columnID/prc columns exist, applies an optional inversion of the raw response, aggregates by well according to agg_type (or with sum for the count models 'poisson' and 'horseshoe'), enforces min_cells_per_well and optionally transforms the aggregated response.

Parameters:

• df – Per-object score DataFrame.

• dependent_variable – Column being aggregated.

• plate – Plate identifier to stamp when the frame is single-plate; ignored (with warning) when multiple plates exist.

• min_cells_per_well – Wells with fewer objects are dropped. Default 25.

• agg_type – 'mean', 'median', 'quantile' or None.

• transform – Optional post-aggregation transform name (see apply_transformation()).

• regression_type – If 'poisson' or 'horseshoe', aggregation uses sum - both model a per-well count, not a per-well average.

• invert_dependent_variable – False/0 = no inversion; True/1 = 1 - x; -1 = 1 / x.

Returns:

(dependent_df, dependent_variable) — the per-well DataFrame and the (possibly transformed) response column name.

Raises:

ValueError – on missing identifiers, unsupported agg_type or unrecognised invert_dependent_variable.

spacr.ml.regression(df, csv_path, dependent_variable='predictions', regression_type=None, alpha=1.0, random_row_column_effects=False, nc='233460', pc='220950', controls=None, dst=None, cov_type=None, plot=False, l1_ratio=0.5, quantile=0.5, hinge_threshold=None, hinge_n_boot=200, huber_t=1.345, qc=True, spline_knots=4, spline_degree=3, legacy_volcano=False, level='grna', level_dst=None, draw_shared_panels=True, group_lasso_lambda='auto', rra_alpha=0.25, rra_permutations=10000, model_plate_position=True, regression_backend=DEFAULT_REGRESSION_BACKEND, verbose=False, transform='', intercept='fitted', intercept_value=0.0, model_data_layout='long')

Run the full regression pipeline: clean, fit, extract coefficients, optional volcano plot.

Parameters:

• df – Long-format DataFrame with gRNA/gene fractions and the dependent variable.

• csv_path – Path used to derive the volcano-plot filename.

• dependent_variable – Response column name. Default 'predictions'.

• regression_type – Model type; auto-selected via check_distribution() when None.

• regression_backend – WHO fits it, one of REGRESSION_BACKEND_ORDER. Default 'statsmodels'. It is threaded to whichever fitter this run reaches – the mixed branch and regression_model() alike – so one setting answers for the whole run, and a backend that cannot fit the chosen family is refused by name before any design is built.

• alpha – Regularisation strength for penalised models.

• random_row_column_effects – If True, fit a mixed model with random row/column effects.

• model_plate_position – Whether plateID, rowID and columnID are terms in a fixed-effects model (or the corresponding grouping/variance structure in a mixed model). Direct calls default to True for API compatibility; new application settings default to False so the terms are opt-in. See prepare_formula() for the measured costs of including or omitting them. False with random_row_column_effects=True is refused: there is nothing left to make random.

• nc – Negative-control gene identifier. Default '233460'.

• pc – Positive-control gene identifier. Default '220950'.

• controls – Explicit list of control identifiers.

• dst – Output directory for plots and summaries.

• cov_type – Optional covariance estimator for the likelihood fits.

• plot – If True, render the volcano plot after fitting.

• l1_ratio – elasticnet L1/L2 mix.

• quantile – Quantile fitted by quantile regression.

• hinge_threshold – Response cut used to binarise for hinge.

• hinge_n_boot – Bootstrap resamples behind the hinge p-values.

• huber_t – Huber tuning constant for rlm/huber.

• group_lasso_lambda – Block penalty for group_lasso.

• rra_alpha – Top fraction of the guide ranking rra aggregates.

• rra_permutations – Draws per guide count in rra’s null.

• qc – Write the regression QC suite into <dst>/regression_qc/.

• legacy_volcano – also draw the ORIGINAL matplotlib volcano. Default False. The interactive one is far faster and the house-style panel is what a run now produces; drawing both gives two volcanoes in two idioms on the same grid. Requires dst and a design matrix, so it is skipped for the mixed branch and when no destination was given.

• level – WHICH MODEL TO FIT – 'grna' (default) or 'gene'. One level, one design; 'both' is refused here because it is two fits. regression_levels() is the entry point that does both. Ignored by regression_type='mixed', which fits the gene fixed and the guide random inside it and so is already both levels.

• level_dst – Where THIS LEVEL’s figures go – the QC suite, the volcano, the publication sheet. Defaults to dst, which is what a single-level run wants. regression_levels() gives each level its own subfolder so two fits cannot overwrite each other’s regression_figure.pdf.

• draw_shared_panels – Draw the guide-fraction and response distributions, which describe the DATA and not the fit. False on the second of two fits, so the figure grid gets one copy rather than two identical ones.

• model_data_layout – 'long' preserves the historical formula with one fitted row per well-guide pair. 'wide' pivots fractions to one row per independent well before any fixed-effects estimator is fitted. Mixed models require their long nesting representation and therefore use long data even when a wide count input was supplied.

Returns:

(model, coef_df, regression_type).

spacr.ml.regression_levels(df, csv_path, dependent_variable='predictions', regression_type=None, level='both', dst=None, **kwargs)

Fit every level the run asked for, SEPARATELY, and return one per level.

THIS IS THE TWO-FIT ENTRY POINT, and the reason it exists is that the one design spaCR used to fit cannot be fitted at all. gene_fraction is the SUM of the gene’s gRNA fractions, so

y ~ fraction:grna + gene_fraction:gene + plateID + rowID + columnID

puts a block of columns and their own sums into one design. Measured on the reference TSG101 screen: 1248 parameters at rank 862 – a 386-dimensional EXACT null space – and the fit statsmodels returned had a residual sum of squares bit-identical to the one you get by adding seven times a null vector to it. See COLLINEAR_FORMULA_FRAGMENT.

Two fits, two tables, TWO CORRECTIONS. Each fit is its own multiple-testing family and is corrected within itself. Pooling them would be wrong twice over: they are not independent – same wells, and the gene regressor IS the sum of the guide regressors – and doubling the family size costs power for no protection. perform_regression() applies the correction per level and writes results_grna.csv and results_gene.csv.

regression_type='mixed' fits ONCE and returns one entry, 'gene': that model has both levels inside it already, the gene as a fixed effect and the guide as a random effect nested in the gene. Its guide output is BLUPs, which is why it cannot be split into two testing families.

Parameters:

• df – long-format DataFrame of gRNA/gene fractions and the dependent variable, passed to regression() for every level.

• csv_path – path passed to regression(), which derives the volcano-plot filename from it.

• level – 'both' (default), 'grna' or 'gene'.

• dst – the run folder. With more than one fit each level’s FIGURES go into <dst>/<level>/ so they cannot overwrite each other; the tables stay in <dst>, where every consumer looks for them.

• kwargs – passed straight through to regression().

Returns:

dict mapping level to (model, coef_df, regression_type), in fit order.

Raises:

ValueError – for a level that is not one of LEVEL_CHOICES.

spacr.ml.regression_model(X, y, regression_type='ols', groups=None, alpha=1.0, cov_type=None, weights=None, l1_ratio=0.5, quantile=0.5, hinge_threshold=None, huber_t=1.345, exposure=None, spline_knots=4, spline_degree=3, group_lasso_lambda='auto', rra_alpha=0.25, rra_permutations=10000, regression_backend=DEFAULT_REGRESSION_BACKEND, verbose=False, response_name='', transform='', glm_force_identity=False)

Dispatch to the requested regression backend and return the fitted model.

Every name in REGRESSION_TYPES is fittable here, and every one of them has a matching branch in process_model_coefficients(), so a model that fits can always be turned into a coefficient table.

The backends, and what each is for:

ols	Ordinary least squares on a continuous well response.
wls	Weighted least squares; weights is the well’s cell count, so a well of 400 cells outweighs one of 30.
rlm/huber	Robust M-estimation (Huber loss). For outlier-heavy wells: a handful of runaway wells no longer drag the fit.
glm	GLM with the family auto-selected from the response by pick_glm_family_and_link().
poisson	Poisson GLM with a log link and offset(log(exposure)), for per-well counts - so the coefficients are effects on the per-cell RATE, not on the well’s headcount.
quasi_binomial	Binomial GLM whose dispersion is estimated from the Pearson chi-square, for overdispersed fractions.
beta	Beta regression, for a fraction strictly inside (0, 1).
logit/probit	GLM-binomial on a fraction, weighted by cell count.
quantile	Quantile regression at quantile; fits the tail of the response rather than its mean.
mixed	Mixed-effects linear model with groups as the random intercept.
lasso/ridge/elasticnet	Penalised least squares.
hinge	Linear SVM (hinge loss) on a binarised response.
horseshoe	Sparse Poisson GLM with a horseshoe prior (spaCRPower’s power-analysis model), via spacr.power_model.
group_lasso	Penalised least squares with a gene’s guide columns penalised as ONE block, so a gene is selected or dropped as a set rather than one guide at a time - the penalised analogue of the mixed model’s nesting, via spacr.group_lasso.
rra	MAGeCK-style robust rank aggregation: guides ranked by their marginal effect, aggregated to the gene BY RANK with a permutation P value, via spacr.rra. It forms no joint fit, so the collinearity and the p >> n width that constrain every backend above do not reach it.

Settings a backend cannot read are REFUSED, not ignored — see REGRESSION_SETTINGS_USED.

Parameters:

• X – Design matrix (DataFrame; column names become feature names).

• y – Response variable.

• regression_type – One of REGRESSION_TYPES.

• regression_backend –

  WHO fits it – one of

  REGRESSION_BACKEND_ORDER. Default 'statsmodels', which produced every existing result. A backend that cannot fit regression_type is REFUSED here with the reason, not ignored: the two controls constrain each other in both directions

  , and a settings CSV reaches this function

  without passing a panel that could have greyed the entry out.

• groups – Cluster identifiers for the mixed model.

• alpha – Penalty weight for lasso/ridge/elasticnet and the inverse SVM margin for hinge; 'auto' / None picks it by 5-fold cross-validation for all four (mean squared error for the penalised least-squares three, balanced accuracy for hinge).

• cov_type – Covariance estimator for the likelihood fits ('HC0'..'HC3'); None for classical standard errors.

• weights – Per-observation weights - the well’s cell count. Used as var_weights by logit/probit/quasi_binomial and as the WLS weights by wls.

• l1_ratio – elasticnet mix; 1.0 is lasso, 0.0 is ridge.

• quantile – Quantile fitted by quantile regression, in (0, 1).

• hinge_threshold – Cut used to binarise a continuous response for hinge; see binarise_response().

• spline_knots – Knots per continuous covariate for spline.

• spline_degree – Polynomial degree of that basis; 3 is cubic.

• huber_t – Huber tuning constant for rlm/huber, in units of the estimated residual scale. 1.345 gives 95% efficiency under normality.

• exposure – Per-observation exposure (the well’s cell count) used as offset(log(exposure)) by horseshoe and by poisson (and by glm when it auto-selects a Poisson family).

• group_lasso_lambda – The block penalty for group_lasso. Its own key rather than alpha because it is compared against spacr.group_lasso.max_lambda(), which is a property of the design, so a value carried over from a lasso run would mean something else here.

• rra_alpha – The top fraction of the guide ranking alpha-RRA aggregates over. MAGeCK’s 0.25, which is what keeps a gene with one strong guide and three that did not cut findable.

• rra_permutations – Draws per distinct guide count in RRA’s permutation null; 10,000 puts the smallest reportable P value at 1e-4.

Returns:

Fitted statsmodels / sklearn estimator.

Raises:

ValueError – on an unsupported regression_type, or when a setting the chosen backend cannot read was set to a non-default value.

Example

import pandas as pd
X = pd.DataFrame({'Intercept': 1.0, 'fraction': [0.1, 0.5, 0.9]})
model = regression_model(X, pd.Series([0.2, 0.4, 0.7]), 'ols')
model.params['fraction']   # the recovered slope

spacr.ml.resolve_auto_inference(data, settings, *, well_column='prc', guide_column='grna')

Choose analysis_mode for inference='auto' from the design.

The simultaneous model estimates one coefficient per guide from the wells, so it needs more wells than guides – with an intercept and any plate fixed effects on top – before those coefficients are identifiable at all. Below that the design matrix is rank deficient: statsmodels still returns a number for every guide, but the numbers are one arbitrary solution out of infinitely many, and their P values describe nothing.

That is not a hypothetical. The screen this was written for has 824 guides in 587 analysed wells; the published fit had 825 parameters, rank 579 and 8 residual degrees of freedom, and refitting it did not reproduce its own coefficients.

auto therefore picks the permutation test whenever the simultaneous fit would be unidentifiable, and says so. It is deliberately conservative: it needs a real margin (_IDENTIFIABILITY_MARGIN wells per guide) rather than a bare majority, because a design that only just fits is one dropped well away from not fitting.

Anything other than inference='auto' is returned untouched, so an explicit choice is never overridden.

Parameters:

• data – the analysis table (a DataFrame); its distinct well and guide counts, and the permutation block column when present, size the design.

• settings – run settings; inference, analysis_mode, analysis_unit, agg_type and guide_permutation_block are read. It is not modified.

• well_column – column whose distinct values count the wells.

• guide_column – column whose distinct values count the guides.

Returns:

(analysis_mode, reason). reason is a sentence naming the counts, suitable for the log and for the Methods section.

spacr.ml.resolve_glm_transform_conflict(dependent_variable, transform='', available=(), regression_type='glm')

Resolve a transform/family-link conflict before fitting a GLM.

Parameters:

• dependent_variable – the response column as it stands – already the transformed one, if a transform was asked for.

• transform – the transform already applied.

• available – the column names the frame actually holds. Used to confirm the untransformed column is there before switching to it.

• regression_type – only 'glm' chooses its own family, so only 'glm' has this conflict to resolve. Everything else is returned unchanged.

Returns:

(column, transform_in_effect, force_identity, note). note explains any scale change for the run log.

A transform that is not link-like, or a regression type other than 'glm', returns the response unchanged. Resolving the column before the design matrices are built keeps the fit, coefficients, diagnostics, and goodness-of-fit summary on the same scale.

spacr.ml.resolve_levels(regression_type, level='both')

Which level(s) a run fits, given the backend and the level setting.

mixed fits ONE model that already contains both levels – the gene as a fixed effect and the guide as a random effect nested inside it – so it ignores level entirely and answers ('gene',). That is why the GUI greys the dropdown out rather than hiding it: the setting exists, but this model does not read it.

Every other backend is fixed effects only and cannot nest, so it fits one level at a time and level chooses which. 'both' is TWO FITS.

Parameters:

• regression_type – the backend name, or None (not yet chosen).

• level – 'both' (default), 'grna' or 'gene'.

Returns:

a tuple of levels to fit, in the order they are fitted.

Raises:

ValueError – for a level that is not one of LEVEL_CHOICES.

spacr.ml.resolve_regression_src(requested, automatic)

Resolve the root directory used for regression output.

A blank requested value selects automatic. An existing requested directory is used directly. If only the final path component is missing, that directory is created; missing parent directories are never created. A requested file, an unavailable parent, or a directory-creation error returns the automatic location with an explanatory message.

Parameters:

• requested – Requested output directory, or None/blank to use the automatic location.

• automatic – Existing fallback directory, normally the directory containing the first count table.

Returns:

A (path, message) tuple. message is 'automatic' when no override was requested; otherwise it describes the selected directory or the reason for falling back.

spacr.ml.results_folder_kind(settings) → str

What a run’s results folder is NAMED after.

The inference method when it decides the answer, and the regression type otherwise. Under analysis_mode='guide_permutation' the regression type is never read – ols and mixed produce byte-identical results – so a folder called ridge would name something the run did not do.

PUBLIC, AND THE ONLY COPY. A test that re-derived this rule went stale when the rule changed and reported 39 missing CSVs while every run that wrote them was fine, which is the failure the results_dir helper in tests/test_cov_ml_perform_regression.py was already written to prevent once. A suite pointing at the wrong file is worse than a silent one.

Parameters:

settings – run settings mapping, or None (treated as empty); only analysis_mode and regression_type are read.

Returns:

'guide_permutation', 'auto' when no regression type is set, or the regression type as a string.

spacr.ml.save_summary_to_file(model, file_path=SUMMARY_FILENAME)

Write model.summary().as_text() to file_path as plain text.

The content is the statsmodels text summary, never CSV – which is why the default name is SUMMARY_FILENAME and no longer summary.csv. Older runs on disk wrote mode_summary.csv; every reader in this repository accepts both, see SUMMARY_FILENAMES.

Parameters:

• model – Fitted statsmodels results object.

• file_path – Destination path. Default SUMMARY_FILENAME.

Returns:

the path written, or None if there was nothing to write.

NEVER RAISES INTO A FINISHED RUN. This is called after every table has been written; a backend whose summary() throws must not take the run down with it, and the caller is told by the None rather than by a traceback.

spacr.ml.scale_variables(X, y)

Min-max scale the independent (X) and dependent (y) variables to [0, 1].

Constant columns are passed through UNCHANGED. MinMaxScaler maps a column with zero range to all-zeros, and patsy’s intercept is exactly such a column, so scaling a design matrix used to silently delete its intercept: statsmodels then fitted a model through the origin and still printed an Intercept row, of 0.000, in the summary. Every coefficient in that fit absorbs the mean it can no longer estimate.

Parameters:

• X – Design matrix (DataFrame).

• y – Response, as a 2-D array or single-column frame.

Returns:

(X_scaled, y_scaled) - a DataFrame with X’s columns and a 2-D numpy array.

Example

X = pd.DataFrame({'Intercept': 1.0, 'a': [1.0, 2.0, 3.0]})
scale_variables(X, np.array([[0.0], [1.0], [2.0]]))[0]['Intercept']
# -> 1.0, 1.0, 1.0   (not 0.0, 0.0, 0.0)

spacr.ml.screen_is_blockable(df) → bool

Whether screenID can be a term in this frame’s design.

True only when the column exists and carries more than one distinct value. A single-screen project is the normal case and must be untouched by the design: it has no screenID at all, or one value, and either way the term would be a constant column.

The same rule spacr.measurement_scan._dummy_block() applies, stated once for the formula path so a frame cannot be blocked on by one and not the other.

Parameters:

df – the design DataFrame, or None (returns False); its screenID column is compared as strings.

spacr.ml.select_glm_family(y)

Choose a statsmodels GLM family from the range and type of the response.

A coarser rule than pick_glm_family_and_link(), which also sets the link: binary values give Binomial, any other values inside [0, 1] give QuasiBinomial, non-negative integers give Poisson and everything else Gaussian.

Parameters:

y – Response vector.

Returns:

An unfitted statsmodels family instance on its default link.

spacr.ml.shap_analysis(model, X_train, X_test)

Build a SHAP summary beeswarm for X_test.

The beeswarm is rendered with pyqtgraph so it can be embedded in the same scene-based figure workflow as other model-explanation plots.

The function returns a live FastPlot; it neither writes a file nor returns a matplotlib figure. Pass the result to write_plot() to export it in the configured figure format.

Parameters:

• model – Fitted estimator compatible with shap.Explainer.

• X_train – Training features used to seed the explainer.

• X_test – Test features to explain.

Returns:

A FastPlot holding the beeswarm, or None when Qt is unavailable or the attribution matrix cannot be plotted.

spacr.ml.summary_for_console(model, *, verbose=False, limit=CONSOLE_COEFFICIENT_LIMIT) → str

Return a statsmodels summary sized for terminal output.

When the coefficient table exceeds limit, the diagnostic header and notes are retained while the table is replaced by a pointer to the saved summary and the sortable Coefficients view. Set verbose=True to return the complete statsmodels rendering.

Parameters:

• model – Fitted model result with a summary() method.

• verbose – Return the complete summary regardless of its size.

• limit – Maximum coefficient rows printed in compact mode.

Returns:

Complete or compact plain-text model summary.

spacr.ml.write_plot(plot, path, title='')

Write a pyqtgraph plot out and announce it, like publish does.

The counterpart of spacr.figure_sink.publish() for a scene rather than a matplotlib figure: the format follows the user’s preference, the file NAME follows the format, and the written file reaches the gallery, because saved and visible are the same event.

None writes nothing, announces nothing and returns None – a plot that could not be built must not take the run down after the model has been fitted and every object scored.

Parameters:

• plot – a FastPlot, or None.

• path – where to write it; the extension may be rewritten.

• title – the name the gallery tile carries.

Returns:

the path written, or None.