Check diagnostics#

Check sampler errors and chain behavior before interpreting posterior estimates.

import jax.numpy as jnp
import numpy as np
import tensorflow_probability.substrates.jax.bijectors as tfb
import tensorflow_probability.substrates.jax.distributions as tfd

import liesel.goose as gs
import liesel.model as lsl

Prepare the example#

These examples require the regression model and sampling results from Sample your first posterior, with sigma_sq included in the stored draws.

Read the summary#

summary = gs.Summary(results, selected=["beta", "sigma_sq"])

summary.to_dataframe().set_index("var_fqn")[
    ["mean", "sd", "mcse_mean", "ess_bulk", "rhat"]
].round(3)

	mean	sd	mcse_mean	ess_bulk	rhat
var_fqn
beta[0]	0.990	0.091	0.002	1444.217	1.004
beta[1]	1.899	0.158	0.004	1432.583	1.005
sigma_sq	1.039	0.065	0.001	2172.744	1.003

summary.aggregate_diagnostics().round(3)

	ess_bulk	ess_tail	rhat	aggregated_by
parameter
beta	1432.583	1540.084	1.005	min (ess); max (rhat)
sigma_sq	2172.744	1970.398	1.003	min (ess); max (rhat)

summary.error_df().reset_index().filter(["error_msg", "phase", "count"])

	error_msg	phase	count
0	divergent transition	warmup	50
1	divergent transition	posterior	0

aggregate_diagnostics() reports the smallest bulk/tail ESS and the largest R-hat within each parameter variable. This keeps a coefficient vector’s worst diagnostics visible without printing every element.

R-hat compares variation within and between chains.
Bulk and tail ESS describe how much information the dependent draws contain about the center and tails.
MCSE measures Monte Carlo uncertainty in an estimated summary, not the posterior uncertainty represented by its standard deviation or interval.

As a practical check, aim for R-hat below 1.01 across several chains and enough ESS for the summaries you need. Compare MCSE with the precision your analysis requires. See the Stan diagnostic guidance for thresholds and their limitations. Good diagnostics are evidence, not proof of convergence.

Inspect the paths#

The first tutorial shows trace plots of all three estimated quantities. Goose also provides plot_trace() and plot_cor() for quick trace and autocorrelation plots.

Look for drift, chains occupying different regions, and long stretches with little movement. Autocorrelation helps identify slow exploration. plot_pairs() can reveal dependence between parameters; plot_param() combines several views of one parameter.

Investigate problems#

Symptom	Next check
Non-finite starting density or acceptance probability	Inspect `model.diagnose()`, parameter support, calculations, and starting values.
Divergent HMC/NUTS transitions	Investigate posterior geometry, scaling, and parameterization; a higher target acceptance may help but is not a general repair.
Maximum NUTS tree depth	Check mixing and parameter dependence. Longer trajectories cost more, so investigate before increasing the limit.
Poor chain overlap or high R-hat	Inspect the affected parameters, starting regions, and possible separate modes.
Small ESS or large MCSE	Inspect autocorrelation and blocking. Once exploration is satisfactory, more posterior draws can improve precision.

Divergences can indicate biased exploration; a tree-depth limit is primarily an efficiency concern. Kernel error codes differ, so use their accompanying messages rather than comparing numeric codes across methods.

Check acceptance#

summary.acceptance_prob_df().reset_index()[
    ["phase", "acceptance_probability", "position_moved"]
].round(3)

	phase	acceptance_probability	position_moved
0	posterior	0.881	NaN
1	warmup	0.791	NaN

The table separates acceptance probabilities from the fraction of transitions that moved. NUTS does not report a move indicator, so its summary shows NaN. For kernels that report movement, position_moved is the fraction of transitions that changed the position.

High acceptance alone does not establish good mixing, and an exact Gibbs update does not need an accept/reject step.

To investigate adaptation, record kernel states during sampling.

Check diagnostics

Contents