---
title: "AR(1) stochastic volatility steady-state BVAR"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{AR(1) stochastic volatility steady-state BVAR}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---



Here we estimate a steady-state BVAR model with AR(1) stochastic volatility, see `?bvar` for details.

We will estimate the model on a quarterly US data set from Koop and Korobilis (2010) on the inflation rate $\Delta \pi_t$ (the annual percentage change in a chain-weighted GDP price index), the unemployment rate $u_t$ (seasonally adjusted civilian unemployment rate, all civilian workers aged 16 years or older) and the interest rate $r_t$ (yield on the three-month Treasury bill rate). The sample is 1953Q1-2006Q3 and we have the data vector

$$
y_t = 
\begin{pmatrix} \Delta \pi_t \\
u_t \\
r_t
\end{pmatrix}
$$

First, let's load the package, then import and plot the data.


``` r
library(SteadyStateBVAR)
data("KoopKorobilis2010")
yt <- KoopKorobilis2010
plot.ts(yt)
```

![](figure/AR(1)-1-1.png)

Let's create the bvar object which we will use throughout here.


``` r
bvar_obj <- bvar(data = yt)
```

We choose 2 lags and only a constant as the deterministic variable.


``` r
bvar_obj <- setup(bvar_obj,
                  p=2,
                  deterministic = "constant")
```

We set the overall tightness to $\lambda_1 = 0.20$, cross-equation tightness to $\lambda_2 = 0.50$ and the lag decay rate to $\lambda_3 = 1.00$. For the prior means on the first own lags, we set them to $0.6$ for $\Delta \pi_t$ and $0.9$ for $u_t$ and $r_t$. Note that the prior mean on the first own lag of inflation is set to $0.6$ instead of $0$ to reflect some degree of persistence in the series (even though it is a growth rate variable).



``` r
lambda_1 <- 0.20
lambda_2 <- 0.50
lambda_3 <- 1.00

fol_pm=c(0.6, # delta pi
         0.9,  #u
         0.9)  #R
```

Now, for the steady-state coefficients we use some toy values (let us pretend that they are expert based).
Remember that we only have a constant now, so $q=1$ and therefore $\Psi$ only has one column $\psi_1=\Psi$. Since $d_t = 1 \ \forall \ t$, we have $\Psi d_t = \mu_t$ which simplifies to $\Psi = \mu$ and as such we can directly interpret $\Psi$ as the unconditional mean, i.e. the steady-state.


``` r
theta_Psi <- 
  c(
  ppi(1.90, 2.10, interval=0.95)$mean,   #Psi: delta pi
  ppi(3.80, 4.50, interval=0.95)$mean,   #Psi: u
  ppi(2.60, 3.90, interval=0.95)$mean    #Psi: r
  )

Omega_Psi <- 
  diag(
  c(
  ppi(1.90, 2.10, interval=0.95)$var,    #Psi: delta pi
  ppi(3.80, 4.50, interval=0.95)$var,    #Psi: u
  ppi(2.60, 3.90, interval=0.95)$var     #Psi: r
  )
  )
```

Now we need to specify our stochastic volatility priors. See `?priors` for more information about the prior specification. Below I take inspiration from Carriero, Clark, and Marcellino (2024), which uses the exact same AR(1) stochastic volatility specification, but for a conventional BVAR.


``` r
k <- bvar_obj$setup$k
n_free_params_A <- bvar_obj$setup$n_free_params_A
sigma2 <- diag(bvar_obj$setup$Sigma_AR)

SV_priors_AR1 <- list(
                      theta_A            =  rep(0, n_free_params_A),
                      Omega_A            =  diag(10, n_free_params_A),
                      theta_gamma_0      =  0.1 * log(sigma2),
                      Omega_gamma_0      =  diag(2, k),
                      theta_gamma_1      =  rep(0.9, k),
                      Omega_gamma_1      =  diag(0.04, k),
                      theta_log_lambda_1 =  log(sigma2),
                      Omega_log_lambda_1 =  diag(2, k),
                      V_Phi              = (10 - k - 1) * 0.03 * diag(k),
                      m_Phi              =  10
                     )
```

Here `sigma2` contains the residual variances from AR($p$) models (the same ones we used in the Minnesota prior).

Let's put everything into the `priors()` function.


``` r
bvar_obj <- priors(bvar_obj,
                   lambda_1 = lambda_1,
                   lambda_2 = lambda_2,
                   lambda_3 = lambda_3,
                   first_own_lag_prior_mean =fol_pm,
                   theta_Psi = theta_Psi,
                   Omega_Psi = Omega_Psi,
                   SV = TRUE,
                   SV_type = "AR1",
                   SV_priors = SV_priors_AR1)
```

Now we can fit the model. Note that we can use arguments from `rstan::sampling()` such as `control` where we can tweak `max_treedepth` and `adapt_delta`.


``` r
bvar_obj <- fit(bvar_obj,
                H = 40,
                d_pred = matrix(rep(1, 40)),
                iter = 4000,
                warmup = 1000,
                chains = 2,
                cores = 2,
                control = list(max_treedepth = 12, adapt_delta = 0.999))
```

Now lets see the posterior means


``` r
summary(bvar_obj, stat="mean", t = 215) #t = 215 for covariance matrix
#> Posterior mean estimates
#> ------------------------
#> 
#> 
#> beta
#> --------------------------------------------------------------------------------             
#>               delta pi     u     r
#>   delta pi.l1     1.27  0.02  0.17
#>   u.l1           -0.09  1.17 -0.16
#>   r.l1            0.00 -0.01  1.04
#>   delta pi.l2    -0.28  0.01 -0.11
#>   u.l2            0.07 -0.23  0.17
#>   r.l2           -0.01  0.02 -0.11
#> --------------------------------------------------------------------------------
#> 
#> 
#> Psi
#> --------------------------------------------------------------------------------          
#>            [,1]
#>   delta pi 2.00
#>   u        4.28
#>   r        3.51
#> --------------------------------------------------------------------------------
#> 
#> 
#> Sigma_u,t (t = 215)
#> --------------------------------------------------------------------------------
#>          delta pi     u     r
#> delta pi     0.07 -0.01  0.02
#> u           -0.01  0.03 -0.02
#> r            0.02 -0.02  0.17
#> --------------------------------------------------------------------------------
#> 
#> 
#> A
#> --------------------------------------------------------------------------------          
#>            delta pi    u r
#>   delta pi     1.00 0.00 0
#>   u            0.13 1.00 0
#>   r           -0.24 0.45 1
#> --------------------------------------------------------------------------------
#> 
#> 
#> gamma_0
#> --------------------------------------------------------------------------------
#> delta pi        u        r 
#>    -0.14    -0.15    -0.09 
#> --------------------------------------------------------------------------------
#> 
#> 
#> gamma_1
#> --------------------------------------------------------------------------------
#> delta pi        u        r 
#>     0.95     0.95     0.94 
#> --------------------------------------------------------------------------------
#> 
#> 
#> Phi
#> --------------------------------------------------------------------------------          
#>            delta pi    u    r
#>   delta pi     0.06 0.04 0.06
#>   u            0.04 0.07 0.07
#>   r            0.06 0.07 0.14
#> --------------------------------------------------------------------------------
```

You can always look at the `stanfit` object `bvar_obj$fit$stan` directly if you want. Note that
the `z`'s below are not parameters per se, they are simply used in a reparameterization trick to sample
the log volatilities more efficiently.


``` r
print(bvar_obj$fit$stan)
#> Inference for Stan model: steady_state_bvar_AR1_stochastic_volatility.
#> 2 chains, each with iter=4000; warmup=1000; thin=1; 
#> post-warmup draws per chain=3000, total post-warmup draws=6000.
#> 
#>                        mean se_mean    sd  2.5%   25%   50%   75%  97.5% n_eff Rhat
#> beta[1,1]              1.27    0.00  0.06  1.16  1.23  1.27  1.31   1.38  5160    1
#> beta[1,2]              0.02    0.00  0.04 -0.05  0.00  0.02  0.05   0.09  5901    1
#> beta[1,3]              0.17    0.00  0.08  0.00  0.11  0.17  0.22   0.33  6117    1
#> beta[2,1]             -0.09    0.00  0.04 -0.16 -0.11 -0.09 -0.07  -0.02  5951    1
#> beta[2,2]              1.17    0.00  0.05  1.06  1.13  1.17  1.21   1.27  5008    1
#> beta[2,3]             -0.16    0.00  0.08 -0.30 -0.21 -0.16 -0.11   0.00  5132    1
#> beta[3,1]              0.00    0.00  0.02 -0.03 -0.01  0.00  0.01   0.04  5985    1
#> beta[3,2]             -0.01    0.00  0.02 -0.04 -0.02 -0.01  0.00   0.02  6203    1
#> beta[3,3]              1.04    0.00  0.06  0.93  1.00  1.04  1.08   1.16  5446    1
#> beta[4,1]             -0.28    0.00  0.06 -0.39 -0.32 -0.28 -0.24  -0.16  5177    1
#> beta[4,2]              0.01    0.00  0.04 -0.07 -0.02  0.01  0.03   0.08  6059    1
#> beta[4,3]             -0.11    0.00  0.08 -0.27 -0.17 -0.11 -0.05   0.05  6217    1
#> beta[5,1]              0.07    0.00  0.03  0.00  0.05  0.07  0.09   0.14  6188    1
#> beta[5,2]             -0.23    0.00  0.05 -0.33 -0.26 -0.23 -0.19  -0.13  4843    1
#> beta[5,3]              0.17    0.00  0.07  0.03  0.12  0.17  0.22   0.31  5318    1
#> beta[6,1]             -0.01    0.00  0.02 -0.04 -0.02 -0.01  0.01   0.03  6261    1
#> beta[6,2]              0.02    0.00  0.02 -0.01  0.01  0.02  0.04   0.06  6238    1
#> beta[6,3]             -0.11    0.00  0.06 -0.22 -0.14 -0.11 -0.07   0.00  5555    1
#> Psi[1,1]               2.00    0.00  0.05  1.90  1.96  2.00  2.03   2.09 12709    1
#> Psi[2,1]               4.28    0.00  0.18  3.93  4.16  4.28  4.40   4.62 10680    1
#> Psi[3,1]               3.51    0.00  0.33  2.83  3.29  3.52  3.74   4.16 10833    1
#> z[1,1]                 0.01    0.00  0.41 -0.76 -0.28 -0.01  0.27   0.84  7481    1
#> z[1,2]                 1.25    0.00  0.43  0.43  0.96  1.23  1.54   2.15  7952    1
#> z[1,3]                -1.02    0.01  0.53 -2.01 -1.36 -1.04 -0.68   0.06  7495    1
#> z[2,1]                 0.00    0.01  0.94 -1.85 -0.63 -0.01  0.63   1.82 11636    1
#> z[2,2]                 0.24    0.01  0.97 -1.70 -0.42  0.25  0.89   2.11 10038    1
#> z[2,3]                -0.08    0.01  0.98 -1.98 -0.76 -0.08  0.58   1.85 12678    1
#> z[3,1]                 0.12    0.01  0.96 -1.75 -0.52  0.13  0.77   2.01 11131    1
#> z[3,2]                 0.26    0.01  1.01 -1.75 -0.42  0.28  0.95   2.22 11982    1
#> z[3,3]                -0.11    0.01  1.00 -2.09 -0.79 -0.10  0.57   1.83 14312    1
#> z[4,1]                -0.38    0.01  0.95 -2.24 -1.04 -0.38  0.29   1.42 11655    1
#> z[4,2]                -0.22    0.01  0.97 -2.13 -0.86 -0.22  0.43   1.67 13269    1
#> z[4,3]                -0.24    0.01  0.99 -2.18 -0.92 -0.22  0.43   1.71 13366    1
#> z[5,1]                -0.25    0.01  0.95 -2.11 -0.90 -0.26  0.39   1.60 11303    1
#> z[5,2]                -0.24    0.01  0.99 -2.19 -0.91 -0.24  0.44   1.68 11485    1
#> z[5,3]                -0.23    0.01  0.99 -2.15 -0.89 -0.22  0.44   1.65 12966    1
#> z[6,1]                -0.15    0.01  0.96 -2.00 -0.81 -0.15  0.49   1.74 14683    1
#> z[6,2]                -0.11    0.01  0.97 -1.97 -0.76 -0.11  0.52   1.84 13910    1
#> z[6,3]                -0.19    0.01  0.98 -2.15 -0.84 -0.19  0.47   1.74 11176    1
#> z[7,1]                 0.08    0.01  0.97 -1.80 -0.57  0.08  0.74   2.00 13398    1
#> z[7,2]                 0.03    0.01  0.97 -1.85 -0.62  0.01  0.67   1.92 12425    1
#> z[7,3]                -0.14    0.01  0.98 -2.05 -0.80 -0.16  0.53   1.86 13863    1
#> z[8,1]                 0.27    0.01  0.92 -1.58 -0.33  0.28  0.89   2.05 11512    1
#> z[8,2]                 0.04    0.01  0.99 -1.92 -0.64  0.04  0.73   1.96 12908    1
#> z[8,3]                -0.08    0.01  0.99 -2.05 -0.74 -0.08  0.57   1.85 13953    1
#> z[9,1]                 0.43    0.01  0.93 -1.38 -0.21  0.43  1.06   2.26 13284    1
#> z[9,2]                 0.23    0.01  0.98 -1.67 -0.43  0.22  0.87   2.15 12659    1
#> z[9,3]                 0.01    0.01  0.96 -1.80 -0.63  0.02  0.66   1.92 11587    1
#> z[10,1]               -0.15    0.01  0.92 -1.96 -0.77 -0.15  0.47   1.69 11896    1
#> z[10,2]                0.13    0.01  0.96 -1.78 -0.53  0.13  0.78   2.03 12824    1
#> z[10,3]               -0.14    0.01  0.99 -2.03 -0.80 -0.15  0.52   1.77 12533    1
#> z[11,1]               -0.21    0.01  0.93 -2.06 -0.83 -0.21  0.42   1.59 13594    1
#> z[11,2]                0.10    0.01  0.96 -1.75 -0.54  0.09  0.73   2.05 11469    1
#> z[11,3]               -0.21    0.01  0.97 -2.17 -0.86 -0.20  0.45   1.66 11544    1
#> z[12,1]               -0.32    0.01  0.96 -2.19 -0.97 -0.32  0.32   1.59 13902    1
#> z[12,2]                0.21    0.01  0.96 -1.71 -0.43  0.22  0.85   2.08 13973    1
#> z[12,3]               -0.19    0.01  1.00 -2.17 -0.85 -0.18  0.49   1.80 11531    1
#> z[13,1]               -0.11    0.01  0.92 -1.94 -0.75 -0.10  0.51   1.69 12618    1
#> z[13,2]                0.37    0.01  0.98 -1.53 -0.25  0.37  1.02   2.33 11793    1
#> z[13,3]               -0.10    0.01  0.97 -1.97 -0.77 -0.10  0.55   1.84 11747    1
#> z[14,1]               -0.19    0.01  0.97 -2.05 -0.87 -0.18  0.48   1.66 12023    1
#> z[14,2]                0.40    0.01  0.98 -1.50 -0.26  0.40  1.05   2.32 11069    1
#> z[14,3]               -0.14    0.01  0.97 -2.06 -0.80 -0.13  0.49   1.81 12887    1
#> z[15,1]               -0.14    0.01  0.95 -2.04 -0.80 -0.15  0.50   1.73 12416    1
#> z[15,2]                0.35    0.01  0.97 -1.55 -0.29  0.36  1.00   2.28 13583    1
#> z[15,3]               -0.06    0.01  1.00 -2.03 -0.73 -0.05  0.62   1.89 13108    1
#> z[16,1]                0.06    0.01  0.96 -1.83 -0.58  0.05  0.71   1.91 13162    1
#> z[16,2]                0.48    0.01  0.96 -1.43 -0.16  0.48  1.11   2.37 12278    1
#> z[16,3]                0.01    0.01  0.98 -1.86 -0.66  0.02  0.69   1.91 12646    1
#> z[17,1]                0.11    0.01  0.95 -1.78 -0.52  0.11  0.76   1.95 11898    1
#> z[17,2]                0.55    0.01  0.97 -1.34 -0.12  0.54  1.21   2.45 12485    1
#> z[17,3]                0.07    0.01  0.97 -1.82 -0.58  0.07  0.73   1.95 13700    1
#> z[18,1]                0.17    0.01  0.94 -1.67 -0.46  0.16  0.80   2.00  9939    1
#> z[18,2]                0.72    0.01  1.00 -1.30  0.05  0.73  1.37   2.69 10622    1
#> z[18,3]                0.16    0.01  0.99 -1.77 -0.50  0.16  0.83   2.10 12037    1
#> z[19,1]                0.37    0.01  0.96 -1.58 -0.26  0.39  1.03   2.22  9302    1
#> z[19,2]                0.84    0.01  0.97 -1.07  0.20  0.85  1.49   2.70 10167    1
#> z[19,3]                0.16    0.01  0.99 -1.78 -0.51  0.15  0.84   2.08 12696    1
#> z[20,1]                0.13    0.01  0.95 -1.73 -0.50  0.13  0.77   1.98 13397    1
#> z[20,2]                0.51    0.01  0.96 -1.36 -0.14  0.51  1.17   2.41 10028    1
#> z[20,3]                0.13    0.01  0.96 -1.73 -0.51  0.12  0.77   2.02 12549    1
#> z[21,1]               -0.07    0.01  0.95 -1.94 -0.71 -0.06  0.58   1.78 11022    1
#> z[21,2]                0.09    0.01  0.98 -1.84 -0.57  0.07  0.77   2.02 13856    1
#> z[21,3]                0.17    0.01  0.97 -1.73 -0.50  0.17  0.84   2.05 15289    1
#> z[22,1]                0.16    0.01  0.94 -1.67 -0.46  0.15  0.79   2.01 11092    1
#> z[22,2]                0.26    0.01  0.97 -1.61 -0.39  0.27  0.93   2.15 11579    1
#> z[22,3]                0.28    0.01  0.99 -1.67 -0.38  0.28  0.94   2.22 12314    1
#> z[23,1]               -0.24    0.01  0.96 -2.15 -0.87 -0.24  0.39   1.62 13287    1
#> z[23,2]               -0.03    0.01  1.00 -1.96 -0.71 -0.04  0.64   1.92 13724    1
#> z[23,3]               -0.11    0.01  0.99 -2.05 -0.76 -0.11  0.55   1.83 13234    1
#> z[24,1]               -0.22    0.01  0.95 -2.09 -0.86 -0.22  0.41   1.68 13098    1
#> z[24,2]                0.11    0.01  0.95 -1.75 -0.52  0.10  0.75   1.98 11951    1
#> z[24,3]               -0.04    0.01  0.97 -1.95 -0.70 -0.03  0.61   1.87 11868    1
#> z[25,1]               -0.11    0.01  0.94 -1.93 -0.73 -0.12  0.52   1.74 12240    1
#> z[25,2]                0.14    0.01  0.97 -1.76 -0.51  0.14  0.79   2.03 13125    1
#> z[25,3]                0.06    0.01  0.97 -1.83 -0.59  0.07  0.69   1.99 11364    1
#> z[26,1]                0.21    0.01  0.95 -1.67 -0.44  0.21  0.86   2.06 11865    1
#> z[26,2]                0.31    0.01  0.95 -1.59 -0.30  0.30  0.94   2.19 12368    1
#> z[26,3]                0.16    0.01  0.98 -1.74 -0.50  0.16  0.82   2.10 13009    1
#> z[27,1]               -0.17    0.01  0.92 -1.92 -0.80 -0.17  0.46   1.62 10960    1
#>  [ reached 'max' / getOption("max.print") -- omitted 3783 rows ]
#> 
#> Samples were drawn using NUTS(diag_e) at Tue Jul 28 05:20:52 2026.
#> For each parameter, n_eff is a crude measure of effective sample size,
#> and Rhat is the potential scale reduction factor on split chains (at 
#> convergence, Rhat=1).
```

We can forecast


``` r
forecast(bvar_obj, pi = 0.68, show_all = TRUE)
```

![](figure/AR(1)-2-1.png)![](figure/AR(1)-2-2.png)![](figure/AR(1)-2-3.png)

Let us plot the log volatility estimates and predictions


``` r
stochastic_volatility_plot(bvar_obj, ci = 0.95, vol = "log_lambda")
```

![](figure/AR(1)-3-1.png)![](figure/AR(1)-3-2.png)![](figure/AR(1)-3-3.png)

Let us plot the estimates and predictions of the implied innovation standard deviations


``` r
stochastic_volatility_plot(bvar_obj, vol = "sd")
```

![](figure/AR(1)-4-1.png)![](figure/AR(1)-4-2.png)![](figure/AR(1)-4-3.png)

We can also produce orthogonalized IRFs


``` r
IRF(bvar_obj, method = "OIRF", t=215, ci=0.68) #latest t
```

![](figure/AR(1)-5 -1.png)


## References

Carriero, A., Clark, T. E., and Marcellino, M. (2024).
Capturing macro-economic tail risks with Bayesian vector autoregressions. 
*Journal of Money, Credit and Banking*, 56(5), pp. 1099–1127.

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

Villani, M. (2009). Steady-state priors for vector autoregressions. *Journal of Applied Econometrics*, 24(4), pp. 630–650.
