Probabilistic Forecasts¶

Conventional VARs produce point forecasts. A Bayesian VAR produces a full posterior predictive distribution over future paths. This means every forecast comes with calibrated uncertainty — wide bands when the model is unsure, narrow when the data are informative.

import numpy as np
import pandas as pd

from impulso import VAR, VARData
from impulso.samplers import NUTSSampler

Setup¶

We repeat the data-generating process from the quickstart tutorial. The DGP is a VAR(1) with three macro variables — GDP growth, inflation, and an interest rate. If you’ve already worked through that notebook, the setup code below will be familiar.

rng = np.random.default_rng(42)
T = 200
n_vars = 3

A_true = np.array([
    [0.6, 0.0, -0.1],
    [0.2, 0.5, 0.0],
    [0.0, 0.15, 0.4],
])

y = np.zeros((T, n_vars))
for t in range(1, T):
    y[t] = A_true @ y[t - 1] + rng.standard_normal(n_vars) * 0.1

index = pd.date_range("2000-01-01", periods=T, freq="QS")
data = VARData(endog=y, endog_names=["gdp_growth", "inflation", "rate"], index=index)

sampler = NUTSSampler(draws=500, tune=500, chains=2, cores=1, random_seed=42)
fitted = VAR(lags=1, prior="minnesota").fit(data, sampler=sampler)
fitted

FittedVAR(n_lags=1, data=VARData(endog_names=['gdp_growth', 'inflation', 'rate'], exog_names=None), var_names=['gdp_growth', 'inflation', 'rate'], has_exog=False)

Point forecasts¶

Call .forecast(steps=8) to produce an 8-step-ahead forecast. The result is a ForecastResult object that holds the full posterior predictive draws. The .median() method extracts the central tendency — the posterior median at each horizon.

fcast = fitted.forecast(steps=8)
fcast.median()

Each row is a forecast horizon (1 through 8 quarters ahead). The values converge toward the unconditional mean of the process as the horizon increases — a hallmark of stationary VARs.

Credible intervals¶

The .hdi() method computes the highest density interval at a given probability level. An 89% HDI means 89% of the posterior forecast mass falls within these bounds. We use 89% rather than 95% following the ArviZ convention — it avoids the false precision of round numbers.

hdi = fcast.hdi(prob=0.89)

print("Lower bounds:")
print(hdi.lower)
print("\nUpper bounds:")
print(hdi.upper)

Lower bounds:
   gdp_growth  inflation      rate
0   -0.029950   0.053070 -0.071551
1   -0.035203   0.015770 -0.049894
2   -0.034322  -0.007591 -0.037466
3   -0.034431  -0.022282 -0.030866
4   -0.037699  -0.032357 -0.027206
5   -0.037708  -0.038754 -0.028479
6   -0.037572  -0.038684 -0.025764
7   -0.037351  -0.041122 -0.024178

Upper bounds:
   gdp_growth  inflation      rate
0    0.003275   0.089627 -0.037677
1    0.013563   0.067530 -0.006803
2    0.023137   0.051853  0.008579
3    0.027775   0.041415  0.014705
4    0.027247   0.033236  0.018755
5    0.028361   0.028338  0.017704
6    0.029564   0.028623  0.020326
7    0.030857   0.026834  0.021439

The intervals widen at longer horizons. This is expected: uncertainty compounds over time because each forecast step propagates parameter uncertainty forward.

Visualise the forecast¶

The .plot() method produces a fan chart showing the median forecast with shaded credible bands for each variable.

fig = fcast.plot()

The fan chart shows the posterior median (line) and 89% HDI (shaded region) for each variable. The bands widen at longer horizons, reflecting compounding uncertainty. GDP growth and the interest rate show the widest bands, consistent with their stronger cross-variable dependencies in the DGP.

Tidy export¶

For downstream analysis or dashboarding, .to_dataframe() returns the median forecast in a tidy DataFrame format.

fcast.to_dataframe()

Summary¶

Bayesian VAR forecasts provide more than point predictions. The full posterior predictive distribution lets you quantify and communicate forecast uncertainty honestly. For structural questions — what happens to inflation when the central bank raises rates? — see the Structural Analysis tutorial.