cultivars.multivariate.large_dim.gibbs

The Gibbs-sampled BVAR: every prior the conjugate model must refuse.

The conjugate model buys its exact posterior by tying the coefficient prior to the innovation covariance through a Kronecker product, and two whole families of priors cannot pay that price. Litterman’s original prior – cross-variable coefficients shrunk harder than own lags – breaks the factorization by design. And the adaptive shrinkage hierarchies – the horseshoe, spike-and-slab selection, Dirichlet-Laplace, Normal-Gamma – have no fixed variance to factor at all: each coefficient’s variance is a latent state with its own full conditional. This model is where both families estimate. The prior is stated independently of the covariance, the independent Normal-Wishart pairing of Koop and Korobilis (2010),

\[\beta \sim N(\underline\beta, \underline V), \qquad \Sigma \sim IW(\underline S, \underline\nu),\]

and the posterior is reached by Gibbs sampling: one joint generalized-least-squares draw of all coefficients given the covariance,

\[\beta \mid \Sigma, y \sim N(\bar\beta, \bar V), \qquad \bar V^{-1} = \underline V^{-1} + \Sigma^{-1} \otimes X'X,\]

exact and ordering-free; an inverse-Wishart draw of the covariance given the coefficients from the residual scatter; and, for an adaptive prior, one sweep of its scale hierarchy’s own exact conditionals, which rewrites \(\underline V\) before the next coefficient draw.

Two commitments shape the surface. First, what is given up is stated rather than hidden. With the prior independent of the covariance the evidence has no closed form; the marginal likelihood offered is Chib’s (1995) estimate from the Gibbs output, which the two-block sampler makes exact for a static moment prior, and which is refused – not approximated – under an adaptive prior, whose scale block would need reduced runs, and under a dummy-observation prior, whose artificial rows the sampler treats as data. Draws are a Markov chain, not independent, so burn-in and thinning are real choices, and the result’s convergence says whether they were enough. And the joint coefficient draw factorizes nothing, so each sweep costs a Cholesky of a \((k^2 p + k n_{det})\)-square precision: exact at any moderate dimension, and deliberately not pretending to Banbura-scale systems, which belong to the conjugate path. Second, the adaptive prior’s latent diagnostic is a first-class output, Rao-Blackwellized over sweeps: for a selection prior GibbsBVARResult.inclusion_probabilities() – the posterior probability that each lag coefficient is in the slab, which is the George-Sun-Ni restriction search read as a result rather than a procedure – and for a global-local prior GibbsBVARResult.shrinkage_scales(), how much freedom each standardized coefficient retained.

Layout. GibbsBVAR validates the specification on _GibbsBayesianVectorAutoRegressionModel in _internals, whose _fit_gibbs builds the design, lays the prior’s moments and any dummy rows out through _gibbs_static_inputs, runs the sweep – the GLS coefficient draw inline, _draw_inverse_wishart from _core._samplers for the covariance, and the prior’s own _draw_scales and _scale_variance for an adaptive hierarchy – and packs a _VectorGibbsFit; its _marginal_likelihood evaluates the Chib identity. GibbsBVARResult extends _VectorPosteriorDrawsResult, which supplies the credible intervals, the posterior-band impulse responses, the posterior predictive, the stability share, the chain diagnostics and the posterior replications. The priors are in priors; the conjugate model is bayesian; the same sampler with stochastic volatility on the innovations is volatility.

References

Koop, G., & Korobilis, D. (2010). Bayesian multivariate time series methods for empirical macroeconomics. Foundations and Trends in Econometrics, 3(4), 267-358.

George, E. I., Sun, D., & Ni, S. (2008). Bayesian stochastic search for VAR model restrictions. Journal of Econometrics, 142(1), 553-580.

Carvalho, C. M., Polson, N. G., & Scott, J. G. (2010). The horseshoe estimator for sparse signals. Biometrika, 97(2), 465-480.

Griffin, J. E., & Brown, P. J. (2010). Inference with normal-gamma prior distributions in regression problems. Bayesian Analysis, 5(1), 171-188.

Chib, S. (1995). Marginal likelihood from the Gibbs output. Journal of the American Statistical Association, 90(432), 1313-1321.

Example

A sparse VAR(1) – own lags only – under three priors: the selection prior finds the sparsity pattern, the horseshoe shrinks the cross lags, and Chib’s evidence ranks the static priors:

>>> import numpy as np
>>> from cultivars.bayes.priors import (
...     HorseshoePrior,
...     IndependentNormalWishartPrior,
...     SpikeAndSlabPrior,
... )
>>> rng = np.random.default_rng(0)
>>> y = np.zeros((200, 3))
>>> for t in range(1, 200):
...     y[t] = 0.6 * y[t - 1] + rng.standard_normal(3)
>>> selected = GibbsBVAR(y, order=1, prior=SpikeAndSlabPrior()).fit(
...     n_draws=300, n_burn=100, seed=0
... )
>>> included = selected.inclusion_probabilities()[0] > 0.5
>>> bool(np.array_equal(included, np.eye(3, dtype=bool)))
True
>>> shrunk = GibbsBVAR(y, order=1, prior=HorseshoePrior()).fit(
...     n_draws=300, n_burn=100, seed=0
... )
>>> cross = ~np.eye(3, dtype=bool)
>>> bool(np.abs(shrunk.coefficients[0][cross]).max() < 0.15)
True
>>> loose = GibbsBVAR(y, order=1).fit(n_draws=300, n_burn=100, seed=0)
>>> tight = GibbsBVAR(
...     y, order=1, prior=IndependentNormalWishartPrior(tightness=0.05)
... ).fit(n_draws=300, n_burn=100, seed=0)
>>> bool(loose.marginal_likelihood().log_value > tight.marginal_likelihood().log_value)
True

Classes

GibbsBVAR	Non-conjugate Bayesian VAR: independent Normal-Wishart, sampled by Gibbs.
GibbsBVARResult	A fitted Gibbs BVAR: a sampled posterior, and the prior's own diagnostics.