HeteroskedasticSVAR#

class cultivars.multivariate.structural.heteroskedacity.HeteroskedasticSVAR(result, regimes, *, shock_names=None)[source]

Bases: _IdentificationModel[SVARResult]

Identification from declared variance regimes, Rigobon (2003).

If the structural shocks’ variances shift across regimes while the impact matrix does not,

\[\Sigma_m = B \Lambda_m B', \qquad \Lambda_m \text{ diagonal}, \qquad m = 1, \dots, M,\]

then the regime covariances jointly pin down \(B\) up to column order and sign, with no zero, sign pattern or instrument. Normalizing \(\Lambda_1 = I\), write \(\Sigma_1 = L L'\); the whitened second-regime covariance \(L^{-1} \Sigma_2 L^{-\top} = Q \Lambda_2 Q'\) is symmetric, its eigendecomposition gives the rotation \(Q\), and \(B = L Q\). The rotation is unique exactly when the diagonal of \(\Lambda_2\) – the variance ratios – is distinct. With more than two regimes \(B\) is over-identified: the two-regime solution is refined by joint approximate diagonalization, and the residual off-diagonal energy is the data’s verdict on the constant-impact assumption.

Shocks are normalized to unit variance in the first regime – first in order of label appearance – and ordered by descending variance ratio in the second, so the first column is the shock whose volatility shifted most. The per-regime relative variances are reported in the summary; they are the scheme’s entire empirical content, and how far apart they sit is how strongly identified the model is.

Note

The identifying assumption is the one the scheme’s name hides in plain sight: the impact matrix is constant across the declared regimes and only the shock variances move. What the data cannot supply is meaning. The recovered shocks are labelled by their variance behaviour, and the default names are shock1 ... shockk rather than the variable names on purpose: attaching “monetary” to a column is a claim the user must argue from outside the model, and the summary says so. Two shocks whose variances shift by the same factor are indistinguishable within their span; an exact tie is refused, and a separation within sampling noise is flagged WEAK in the diagnostics.

Variables:

  • _source – The closed reduced-form result being identified.

  • _regime_labels – The regime labels in order of first appearance.

  • _assignment – The (nobs,) integer regime of each residual row.

  • _labels – The shock labels, in descending-ratio column order.

See also

  • SVARResult – the complete result returned.

  • stochastic_volatility – the same idea with the variance path estimated rather than declared.

  • non_gaussian – identification from higher moments instead of second-moment shifts.

References

Rigobon, R. (2003). Identification through heteroskedasticity. Review of Economics and Statistics, 85(4), 777-792.

Lanne, M., & Lütkepohl, H. (2008). Identifying monetary policy shocks via changes in volatility. Journal of Money, Credit and Banking, 40(6), 1131-1149.

Lewis, D. J. (2021). Identifying shocks via time-varying volatility. Review of Economic Studies, 88(6), 3086-3124.

Example

A bivariate system whose first structural shock triples its standard deviation halfway through the sample. The declared split recovers both impact columns up to scale, the shock that shifted comes first, and the first regime’s covariance factorizes exactly:

>>> import numpy as np
>>> from cultivars.multivariate.reduced_form.vector_autoregression import VAR
>>> rng = np.random.default_rng(0)
>>> B = np.array([[1.0, 0.5], [0.3, 1.0]])
>>> eps = rng.standard_normal((402, 2))
>>> eps[201:, 0] *= 3.0
>>> y = np.zeros((402, 2))
>>> for t in range(1, 402):
...     y[t] = 0.4 * y[t - 1] + B @ eps[t]
>>> res = VAR(y, order=1, names=("a", "b")).fit()
>>> labels = ["calm"] * 200 + ["volatile"] * (res.resid.shape[0] - 200)
>>> svar = HeteroskedasticSVAR(res, labels).identify()
>>> svar.is_complete, svar.shock_names, dict(svar.diagnostics)["Identification"]
(True, ('shock1', 'shock2'), 'distinct ratios')
>>> dict(svar.diagnostics)["Variances in 'volatile' vs 'calm'"]
'9.902, 0.941'
>>> unit = lambda v: v / np.linalg.norm(v)
>>> bool(np.abs(unit(svar.impact[:, 0]) - unit(B[:, 0])).max() < 0.03)
True
>>> calm = res.resid[:200]
>>> bool(np.allclose(svar.impact @ svar.impact.T, calm.T @ calm / 200))
True

Private attributes#

_assignment

_labels

_regime_labels

Private methods#

_regime_covariances()

Per-regime second moments of the residuals.