Bayesian Inference of Structural Vector Autoregressions (SVAR) with the `bvartools` package
with tags r var svar vector-autoregression bvartoolsThe bvartools allows to perform Bayesian inference of Vector autoregressive (VAR) models, including structural VARs. This post guides through the Bayesian inference of SVAR models in R using the bvartools package.
Data
For this illustration we generate an artificial data set with three endogenous variables, which follows the data generating process
\[y_t = A_1 y_{t - 1} + B \epsilon_t,\]
where
\[ A_1 = \begin{bmatrix} 0.3 & 0.12 & 0.69 \\ 0 & 0.3 & 0.48 \\ 0.24 & 0.24 & 0.3 \end{bmatrix} \text{, } B = \begin{bmatrix} 1 & 0 & 0 \\ -0.14 & 1 & 0 \\ -0.06 & 0.39 & 1 \end{bmatrix} \text{ and } \epsilon_t \sim N(0, I_3). \]
# Reset random number generator for reproducibility
set.seed(24579)
tt <- 500 # Number of time series observations
# Coefficient matrix
A_1 <- matrix(c(0.3, 0, 0.24,
0.12, 0.3, 0.24,
0.69, 0.48, 0.3), 3)
# Structural coefficients
B <- diag(1, 3)
B[lower.tri(B)] <- c(-0.14, -0.06, 0.39)
# Generate series
sd_sigma <- 1
series <- matrix(rnorm(3, 0, sd_sigma), 3, tt + 1) # Raw series with zeros
for (i in 2:(tt + 1)){
series[, i] <- A_1 %*% series[, i - 1] + B %*% rnorm(3, 0, sd_sigma)
}
series <- ts(t(series)) # Convert to time series object
dimnames(series)[[2]] <- c("S1", "S2", "S3") # Rename variables
# Plot the series
plot.ts(series, main = "Artificial time series")
Model generation
library(bvartools)
# Generate basic model
mod <- create_bvarmodel(series,
p = 1,
deterministic = "none",
structural = TRUE,
error = "gamma",
iterations = 5000,
burnin = 2000)Prior specification
mod_w_prior <- add_priors(mod,
coef = list(v_i = 0),
sigma = list(shape = 3, rate = .0001))Draw posteriors
mod_w_prior <- add_initial_values(mod_w_prior, method = "ols")
post <- add_posterior_coefficients(mod_w_prior)Evaluation
Summary statistics
summary(post)##
## Bayesian SVAR model with p = 1
##
## Endogenous variables: S1, S2, S3
##
## Variable: S1
##
## Mean SD Naive SD Time-series SD 2.5% 50% 97.5%
## S1.l1 0.3293 0.02931 0.0004146 0.0004462 0.27127 0.3289 0.3875 *
## S2.l1 0.1125 0.04041 0.0005714 0.0005714 0.03299 0.1129 0.1922 *
## S3.l1 0.6807 0.04105 0.0005806 0.0005806 0.60108 0.6804 0.7626 *
## S1 1.0000 0.00000 0.0000000 0.0000000 1.00000 1.0000 1.0000 *
## S2 0.0000 0.00000 0.0000000 0.0000000 0.00000 0.0000 0.0000
## S3 0.0000 0.00000 0.0000000 0.0000000 0.00000 0.0000 0.0000
##
## Variable: S2
##
## Mean SD Naive SD Time-series SD 2.5% 50% 97.5%
## S1.l1 0.05834 0.03334 0.0004715 0.0004715 -0.007253 0.05793 0.1257
## S2.l1 0.31065 0.04071 0.0005757 0.0005757 0.233159 0.30901 0.3926 *
## S3.l1 0.63546 0.05282 0.0007470 0.0007695 0.533515 0.63427 0.7389 *
## S1 0.18132 0.04546 0.0006429 0.0006623 0.093218 0.18069 0.2716 *
## S2 1.00000 0.00000 0.0000000 0.0000000 1.000000 1.00000 1.0000 *
## S3 0.00000 0.00000 0.0000000 0.0000000 0.000000 0.00000 0.0000
##
## Variable: S3
##
## Mean SD Naive SD Time-series SD 2.5% 50% 97.5%
## S1.l1 0.25716 0.03311 0.0004682 0.0004582 0.19305 0.25673 0.3228 *
## S2.l1 0.16797 0.04150 0.0005869 0.0005869 0.08465 0.16852 0.2463 *
## S3.l1 0.19279 0.05879 0.0008314 0.0008314 0.07859 0.19242 0.3108 *
## S1 0.05137 0.04562 0.0006452 0.0006452 -0.03864 0.05162 0.1399
## S2 -0.31421 0.04384 0.0006200 0.0007002 -0.39941 -0.31465 -0.2264 *
## S3 1.00000 0.00000 0.0000000 0.0000000 1.00000 1.00000 1.0000 *
##
## Variance-covariance matrix:
##
## Mean SD Naive SD Time-series SD 2.5% 50% 97.5%
## S1_S1 0.9862 0.06224 0.0008802 0.0008630 0.8735 0.9838 1.119 *
## S1_S2 0.0000 0.00000 0.0000000 0.0000000 0.0000 0.0000 0.000
## S1_S3 0.0000 0.00000 0.0000000 0.0000000 0.0000 0.0000 0.000
## S2_S2 1.0372 0.06554 0.0009269 0.0009269 0.9195 1.0333 1.175 *
## S2_S3 0.0000 0.00000 0.0000000 0.0000000 0.0000 0.0000 0.000
## S3_S3 0.9864 0.06228 0.0008807 0.0009067 0.8711 0.9842 1.116 *Compare with actual values
solve(B) %*% A_1## [,1] [,2] [,3]
## [1,] 0.30000 0.120000 0.690000
## [2,] 0.04200 0.316800 0.576600
## [3,] 0.24162 0.123648 0.116526solve(B)## [,1] [,2] [,3]
## [1,] 1.0000 0.00 0
## [2,] 0.1400 1.00 0
## [3,] 0.0054 -0.39 1Plotting posterior distributions
t
plot(post, type = "boxplot")
Forecasts
post <- add_forecast_input(post, n_ahead = 10)
post <- add_posterior_forecasts(post)
pred_svar <- predict(post, n_ahead = 10)
plot(pred_svar)
Impulse response functions
irf_svar <- irf(post, n_ahead = 10, impulse = "S1", response = "S3", type = "sir")
plot(irf_svar)
Variance decomposition
fevd_svar <- fevd(post, n_ahead = 10, response = "S3", type = "sir")
plot(fevd_svar)
Compare with frequentist results
Estimate basic VAR
library(vars)
freq_var <- VAR(series, p = 1, type = "none")Estimate structural coefficients
# Specify A0
a <- diag(1, 3)
a[lower.tri(a)] <- NA
freq_var_a <- SVAR(freq_var, Amat = a, max.iter = 1000)
freq_var_a##
## SVAR Estimation Results:
## ========================
##
##
## Estimated A matrix:
## S1 S2 S3
## S1 1.00000 0.0000 0
## S2 0.18177 1.0000 0
## S3 0.05078 -0.3132 1Impulse response function
freq_var_irf <- irf(freq_var_a, n.ahead = 10, impulse = "S1", response = "S3")
plot(freq_var_irf)
Literature
Chan, J., Koop, G., Poirier, D. J., & Tobias, J. L. (2019). Bayesian Econometric Methods (2nd ed.). Cambridge: University Press.
Lütkepohl, H. (2007). New Introduction to Multiple Time Series Analyis (2nd ed.). Berlin: Springer.
Pfaff, B. (2008). VAR, SVAR and SVEC models: Implementation within R package vars. Journal of Statistical Software, 27(4).
