irfCompute#
Purpose#
Compute orthogonalized impulse response functions using Cholesky identification.
Format#
- irf = irfCompute(result, n_ahead)#
- Parameters:
result (struct) – an instance of a
varResultorbvarResultstructure.n_ahead (scalar) – number of horizons to compute (e.g., 20).
var_names (Mx1 string array) – Optional keyword, override variable names from the estimation result.
quiet (scalar) – Optional keyword, set to 1 to suppress printed output. Default = 0.
- Returns:
irf (struct) –
An instance of an
irfResultstructure containing:irf.irf
Array of (n_ahead+1) mxm matrices.
irf.irf[h+1]is the mxm impulse response matrix at horizon h. Element [i, j] is the response of variable i to a one-standard-deviation shock to variable j. Index 1 is the impact response (h=0).irf.n_ahead
Scalar, number of horizons computed.
irf.m
Scalar, number of variables.
irf.var_names
Mx1 string array, variable names.
irf.ident
String, identification method:
"cholesky"or"generalized".
Model#
An impulse response function (IRF) traces the dynamic effect of a one-standard-deviation structural shock to variable \(j\) on variable \(i\) over \(h\) periods.
For a VAR(p) in companion form \(Y_t = F Y_{t-1} + G \varepsilon_t\), the reduced-form IRF at horizon \(h\) is:
where \(F\) is the \(mp \times mp\) companion matrix and \(J = [I_m \; 0 \; \cdots \; 0]\) selects the first \(m\) rows.
Cholesky identification: To give shocks a structural interpretation, the reduced-form innovations are orthogonalized via the Cholesky factorization \(\Sigma = P P'\) where \(P\) is lower triangular. The structural IRF is:
Element \(\Theta_h[i, j]\) is the response of variable \(i\) at horizon \(h\) to a one-standard-deviation shock to variable \(j\).
Identification assumption: Cholesky identification imposes a recursive causal ordering. Variable 1 can affect all others contemporaneously; variable \(m\) is affected by all others but affects none contemporaneously. This assumption is appropriate when there is a natural fast-to-slow ordering (e.g., financial variables respond faster than real activity). Algorithm ———
Extract companion matrix \(F\) and Cholesky factor \(P = \text{chol}(\Sigma)'\) from the VAR estimates.
Iterate: For \(h = 0, 1, \ldots, n\_ahead\):
\[\Theta_h = J \, F^h \, J' \, P\]The companion power \(F^h\) is computed iteratively (matrix multiplication, not matrix exponentiation) for numerical stability.
Store \(\Theta_0, \Theta_1, \ldots, \Theta_{n\_ahead}\) as an array of \(m \times m\) matrices.
Complexity: \(O(n\_ahead \cdot m^2 p^2)\) — dominated by the \(mp \times mp\) matrix multiplications. Sub-millisecond for typical systems. Examples ——–
Monetary Policy Shock#
Trace the effect of a federal funds rate shock on GDP and CPI:
new;
library timeseries;
data = loadd(getGAUSSHome("pkgs/timeseries/examples/macro.dat"));
// Variable ordering: GDP (slow), CPI (medium), FFR (fast policy instrument)
// This ordering means: FFR shocks can affect GDP and CPI contemporaneously,
// but GDP shocks take one period to reach FFR.
result = varFit(data, 4);
struct irfResult irf;
irf = irfCompute(result, 20);
Output:
================================================================================
Impulse Response Functions (cholesky)
Horizons: 0-20
================================================================================
Shock to: GDP
h GDP CPI FFR
--------------------------------------------------------------------------------
0 0.5280 0.0456 0.0919
1 0.1859 0.0810 0.2753
2 0.1600 0.0612 0.4442
...
================================================================================
The impact response (h=0) shows that a 1-SD GDP shock raises GDP by 0.528, CPI by 0.046, and FFR by 0.092 — consistent with the central bank responding to output movements within the quarter.
IRF from BVAR with Shrinkage#
new;
library timeseries;
data = loadd(getGAUSSHome("pkgs/timeseries/examples/macro.dat"));
ctl = bvarControlCreate();
ctl.p = 4;
struct bvarResult br;
br = bvarFit(data, ctl, quiet=1);
// IRF at the posterior mean of B and Sigma
irf = irfCompute(br, 20);
For posterior IRF bands (credible intervals), use irfSvCompute() with
an SV-BVAR result.
Plotting IRFs#
Reshape IRF results into a plot-ready dataframe:
new;
library timeseries;
data = loadd(getGAUSSHome("pkgs/timeseries/examples/macro.dat"));
result = varFit(data, 4, quiet=1);
irf = irfCompute(result, 20, quiet=1);
// Get plot data: (n_ahead+1) x (m*m) matrix with column names
struct DataFrame plot_data;
plot_data = irfPlotData(irf);
// Plot GDP response to FFR shock
plotXY(seqa(0, 1, 21), plot_data[., "GDP<-FFR"]);
Troubleshooting#
IRFs don’t decay to zero:
If the VAR is near-nonstationary (max eigenvalue close to 1), IRFs can be very
persistent. This is not a numerical error — it reflects the model’s dynamics.
Check result.max_eigenvalue. For near-unit-root systems, consider:
Using longer horizons (40-60 periods instead of 20).
Differencing the data.
Adding sum-of-coefficients priors (
bvarFit()with lambda6 > 0).
IRFs are sensitive to variable ordering:
This is inherent to Cholesky identification — different orderings produce different
structural shocks. If the ordering is uncertain, use girfCompute() (generalized
IRF, ordering-invariant) or svarIdentify() (sign restrictions).
Impact response has wrong sign: Check the variable ordering. In Cholesky identification, the first variable’s shock is unrestricted; later variables’ shocks are residualized. A monetary policy variable (FFR) should typically be ordered last so its shock is “purged” of contemporaneous output and price movements. Remarks ——-
Identification: The Cholesky decomposition of \(\Sigma\) is used to orthogonalize the innovations. The ordering of variables in the data determines the recursive causal structure: variable 1 can affect all others contemporaneously, variable 2 can affect variables 3, …, m but not 1, and so on.
To change the identification ordering, reorder the columns of the data
before calling varFit() or bvarFit().
For ordering-invariant responses, use girfCompute() (generalized IRF,
Pesaran & Shin 1998). For theory-based identification with sign/zero restrictions,
see svarIdentify().
For BVAR, the IRF is computed at the posterior mean of B and \(\Sigma\).
For posterior IRF bands, use irfSvCompute() with an bvarSvResult.
Indexing convention:
irf.irf[1, ., .] is the impact response (h=0). irf.irf[h+1, ., .] is the response
at horizon h. Element irf.irf[h+1, i, j] is the response of variable i to
a shock to variable j.
Verification
————
Verified against R vars::irf() with boot=FALSE at \(10^{-6}\) tolerance
on a 2-variable VAR(1) with known DGP. Tests cover impact values, decay patterns
at h=1 and h=2, and the Cholesky lower-triangularity constraint (zero upper-off-diagonal
at h=0). See gausslib-var/tests/r_benchmark.rs.
Additionally verified against ECB BEAR Cholesky IRFs on matched-prior BVAR(4),
covering all 9 shock-response pairs at horizons 0, 10, and 20 (17 tests).
See crossval/bear_matched_irf.e.
References
———-
Lutkepohl, H. (2005). New Introduction to Multiple Time Series Analysis. Springer. Chapter 2.3 (IRF computation), Chapter 9 (structural identification).
Pesaran, M.H. and Y. Shin (1998). “Generalized impulse response analysis in linear multivariate models.” Economics Letters, 58(1), 17-29.
Sims, C.A. (1980). “Macroeconomics and reality.” Econometrica, 48(1), 1-48.
Library#
timeseries
Source#
irf.src
See also
Functions irfSvCompute(), girfCompute(), fevdCompute(), hdCompute(), irfPlotData(), svarIdentify()
See also
Guides Choosing a VAR Model, Verification and Cross-Validation