ess_bulk#

cultivars.bayes.chains.ess_bulk(*chains)[source]#

Bulk effective sample size, for central posterior quantities.

The number of independent draws that would carry the same information about the centre of the posterior as the correlated draws supplied. With \(S\) pooled draws and integrated autocorrelation time \(\tau\),

\[\widehat{S}_{\text{eff}} = \frac{S}{\hat\tau}, \qquad \hat\tau = -1 + 2 \sum_{t=0}^{K} \hat P_t, \qquad \hat P_t = \hat\rho_{2t} + \hat\rho_{2t+1},\]

where \(\hat\rho_t\) is the autocorrelation at lag \(t\) pooled across the split chains through the between-chain variance, so that a chain which has not mixed shows as a slowly decaying correlation rather than as a small within-chain variance. The truncation \(K\) is Geyer’s (1992) initial monotone sequence: the paired sums are cut at the first non-positive pair and then forced non-increasing. The computation is on the rank-normalized split chains of rhat(), and the estimate is capped at \(S \log_{10} S\), since an antithetic chain can otherwise report more effective draws than the estimator’s own precision supports.

Parameters:

*chains (ArrayLike) – One or more (S,) or (S, d) chains of the same shape, with S at least eight and finite throughout.

Return type:

float | NDArray[double]

Returns:

A float for (S,) chains, else (d,); the total across chains, not per chain. nan for a quantity that never moved.

Raises:

Note

Vehtari et al. (2021) recommend at least 100 effective draws per chain before any other diagnostic is read, because \(\widehat{R}\) and the standard errors are themselves noisy below that. The threshold convergence() applies is min_ess times the number of chains.

See also

ess_tail(): The same estimator on quantile indicators, for

interval endpoints.

mcse(): The standard error that the untransformed version of

this size implies for the posterior mean.

References

Geyer, C. J. (1992). Practical Markov chain Monte Carlo.

Statistical Science, 7(4), 473-483.

Vehtari, A., Gelman, A., Simpson, D., Carpenter, B., & Bürkner,

P.-C. (2021). Rank-normalization, folding, and localization: An improved \(\widehat{R}\) for assessing convergence of MCMC. Bayesian Analysis, 16(2), 667-718.

Example

Independent draws have an effective size close to their count:

>>> import numpy as np
>>> rng = np.random.default_rng(2)
>>> 1400 < float(ess_bulk(rng.standard_normal(2000))) < 2600
True

A first-order autoregression with coefficient \(\phi\) has \(\tau = (1 + \phi) / (1 - \phi)\), so at \(\phi = 0.9\) the 2000 draws are worth about a hundred:

>>> x = np.empty(2000)
>>> x[0] = 0.0
>>> for t in range(1, 2000):
...     x[t] = 0.9 * x[t - 1] + rng.standard_normal()
>>> 40 < float(ess_bulk(x)) < 250
True
On this page