Diagnostics

Abacus exposes diagnostics through mmm.diagnostics.

Use this surface to check the design matrix, posterior sampling quality, and posterior predictive fit. For fitted-value plots and predictive sampling, see Posterior Predictive.

Diagnostic surfaces

mmm.diagnostics provides three groups of checks.

Area Summary method Report method What it covers
Raw input screening design_summary(X) design_report(X) Collinearity, constants, and near-constant regressors on raw input columns
MCMC mcmc_summary() mcmc_report() r_hat, ESS, divergences, BFMI, tree depth, acceptance rate
Predictive predictive_summary() predictive_report() RMSE, MAE, NRMSE, NMAE, CRPS, residual moments

The summary methods return pandas DataFrames. The report methods return typed report objects with a to_dict() method for JSON-ready export.

Raw input screening

Use design_summary(X) on the raw design matrix you want to inspect:

design = mmm.diagnostics.design_summary(X)

By default, Abacus checks:

  • all channel_columns
  • all control_columns, when present

You can limit the check to specific variables:

design = mmm.diagnostics.design_summary(
    X,
    variables=["tv", "search", "price_index"],
    vif_threshold=10.0,
    near_constant_threshold=0.99,
)

The returned table includes:

  • variable
  • mean
  • std
  • n_unique
  • dominant_share
  • is_constant
  • is_near_constant
  • vif
  • high_vif
  • max_abs_corr

design_report(X) returns a compact roll-up with matrix rank, condition number, maximum VIF, maximum absolute correlation, and lists of flagged variables.

Screening requirements

Raw input screening requires:

  • all requested columns to exist in X
  • all checked columns to be numeric

Abacus raises a ValueError if a variable is missing or non-numeric.

The method names stay design_summary() and design_report(), but the pipeline now treats them as raw input screening rather than transformed model geometry.

MCMC diagnostics

Use mcmc_summary() after fitting:

mcmc = mmm.diagnostics.mcmc_summary(
    rhat_threshold=1.01,
    ess_threshold=400.0,
)

The summary comes from arviz.summary(..., kind="diagnostics", round_to="none") and adds flag columns such as:

  • high_rhat
  • low_ess_bulk
  • low_ess_tail

mcmc_report() adds model-level diagnostics, including:

  • divergence_count
  • divergence_rate
  • divergence_status and divergence_reason
  • max_rhat
  • min_ess_bulk
  • min_ess_tail
  • bfmi_mean
  • bfmi_min
  • max_tree_depth_hits
  • max_tree_depth_observed
  • mean_acceptance_rate

MCMC summaries and reports retain unrounded diagnostics. Parameter flags use inclusive boundaries: R-hat at or above the threshold and bulk/tail ESS at or below the threshold are flagged, matching the pipeline gate comparisons. Round values only for display, after classification.

Divergence reports distinguish available evidence from unavailable evidence. Missing, empty, malformed or mismatched divergence flags produce divergence_status="unavailable", null count/rate values and a reason. Only a valid array covering the retained chain/draw coordinates can establish zero divergences. Stage 50 warns when this evidence is unavailable, preventing an all-pass diagnostic rollup. Absence alone does not establish that divergences are inapplicable to the sampler.

R-hat and ESS thresholds screen Monte Carlo exploration; passing them does not establish model validity or causal identification. Investigate retained divergences even when other diagnostics pass. For interpretation and remedies, see MCMC Diagnostics for Econometricians.

If idata is missing, Abacus raises an error and tells you to fit the model first.

Example MCMC diagnostic output:

Trace plot example Trace plot example

Predictive diagnostics

Scoring requires exactly matching observation dimensions and coordinate labels. Matching labels may appear in a different order; predictions are reordered to match the target. Missing, extra, duplicate, null or unlabelled observation coordinates raise an error. Dimensions are not implicitly broadcast. To score a subset, select the same intended observations explicitly on both arrays before calling predictive_summary_from_arrays(). Empty observations or samples are rejected.

Predictive diagnostics use the observed target and stored posterior predictive samples:

mmm.sample_posterior_predictive(
    X=X,
    random_seed=42,
    progressbar=False,
)

predictive = mmm.diagnostics.predictive_summary(original_scale=True)

The predictive summary is a one-row DataFrame with:

  • scale
  • num_observations
  • rmse
  • mae
  • nrmse
  • nmae
  • crps
  • residual_mean
  • residual_std

Abacus aligns target and prediction coordinates before flattening. That includes mixed datetime coordinate dtypes when needed.

Predictive metric definitions

Let y_i be an observed target and m_i its posterior predictive mean, with residual r_i = y_i - m_i. Abacus aligns labels and flattens all observation dimensions, including panel dimensions, into one equally weighted aggregate. num_observations counts those entries, not just unique dates.

Field Definition
rmse Square root of the mean of r_i ** 2
mae Mean of abs(r_i)
nrmse RMSE divided by max(y) - min(y) on the scored observations
nmae MAE divided by the same observed range
crps Mean continuous ranked probability score over observations, using all predictive draws
residual_mean Mean of r_i; positive means underprediction
residual_std Standard deviation of r_i with ddof=0
bias The same signed mean as residual_mean, when requested by the array helper or Stage 35

NRMSE and NMAE return NaN when the observed range is approximately zero (numpy.isclose(range, 0.0)). Range normalisation does not make arbitrary windows comparable: their ranges, composition and prediction tasks can differ. For RMSE, MAE and CRPS, lower scores are better on a comparable evaluation set; none is a causal-identification test. CRPS can be written as E|Y - y_i| - 0.5 * E|Y - Y'|, with independent predictive draws Y and Y'. It assesses a predictive distribution rather than only its mean.

Stage 35 also requests empirical coverage. For probability p, take each observation’s predictive quantiles at (1 - p) / 2 and (1 + p) / 2, then average the indicator that the observed target lies between them, including the endpoints. Entries with non-finite targets or bounds are excluded from that coverage denominator; if none remain, coverage is NaN. num_observations remains the full aligned count, not the finite coverage count. These are equal-tailed intervals, distinct from the summary facade’s HDIs. The ordinary mmm.diagnostics.predictive_summary() does not add coverage or bias columns. See holdout interpretation.

Example residual diagnostics:

Residuals over time Residuals over time

Residual histogram Residual histogram

Residuals versus fitted Residuals versus fitted

Residual autocorrelation Residual autocorrelation

Export reports

Use the report objects when you want a compact export format:

import json

report = mmm.diagnostics.mcmc_report()
payload = report.to_dict()

with open("mcmc_report.json", "w", encoding="utf-8") as handle:
    json.dump(payload, handle, indent=2)

The same pattern works for design_report(...) and predictive_report().

Pipeline outputs

The pipeline diagnostics stage uses the same retained diagnostic surfaces to write report tables and text summaries. If you run the pipeline, those stage artefacts should match the behaviour documented here.

In the structured pipeline, the raw-input screening rows in diagnostics_report.csv use the phase label raw_input_screening instead of design so the machine-readable output matches the wording here.

Common pitfalls

  • Running mcmc_summary() or mcmc_report() before fitting
  • Running predictive diagnostics before sampling posterior predictive values
  • Passing non-numeric columns into design_summary(X)
  • Treating predictive diagnostics as a substitute for design or MCMC checks