library(bsvars)
data(us_fiscal_lsuw)
spec = specify_bsvar$new(us_fiscal_lsuw)
burn = estimate(spec, S = 1000)
post = estimate(burn, S = 10000)
\[ \]
\[ \]
bsvars package features
bsvars models and identification
bsvars modeling of monetary policy
\[ \]
\[ \]
Lecture Slides as a Website
GitHub repo to reproduce the slides and results
bsvars package features
bsvars package features
\[ \]
- Bayesian estimation of Structural VARs
- seven heteroskedastic processes
- identification through exclusion restrictions, heteroskedasticity, and non-normality
- efficient and fast Gibbs samplers
bsvars package features
\[ \]
- excellent computational speed
- frontier econometric techniques
- compiled code using cpp via Rcpp and RcppArmadillo
- data analysis in R
bsvars package features
\[ \]
- package and loading
- simple model setup
- simple estimation
bsvars package features
\[ \]
- structural analyses
- predictive analyses
bsvars package features
\[ \]
- workflow with the pipe
library(bsvars)
data(us_fiscal_lsuw)
us_fiscal_lsuw |>
specify_bsvar$new() |>
estimate(S = 1000) |>
estimate(S = 10000) -> post
post |> compute_impulse_responses(horizon = 12) -> irfs
post |> compute_variance_decompositions(horizon = 12) -> fevd
post |> compute_historical_decompositions() -> hds
post |> compute_structural_shocks() -> ss
post |> compute_conditional_sd() -> csds
post |> forecast(horizon = 4) -> forebsvars package features
- progress bar
bsvars models and identification
bsvars models
\[ \]
Structural VAR
\[\begin{align} \text{reduced form:}&&\mathbf{y}_t &= \mathbf{A}\mathbf{x}_t + \boldsymbol{\varepsilon}_t \\ \text{structural form:}&&\mathbf{B}_0\boldsymbol{\varepsilon}_t &= \mathbf{u}_t \\ \text{structural shocks:}&&\mathbf{u}_t\mid\mathbf{x}_t &\sim N\left( \mathbf{0}_N, \text{diag}\left(\boldsymbol{\sigma}_t^2\right) \right) \end{align}\]
- interpretable structural specification
- identification through exclusion restrictions, heteroskedasticity, and/or non-normality
- facilitates application of frontier numerical techniques
bsvars models
Reduced form hierarchical prior
\[\begin{align} \text{autoregressive slopes:}&& [\mathbf{A}]_{n\cdot}'\mid\gamma_{A.n} &\sim N_{Np+1}\left( \underline{\mathbf{m}}_{n.A}, \gamma_{A.n}\underline{\Omega}_A \right) \\ \text{autoregressive shrinkage:}&&\gamma_{A.n} | s_{A.n} &\sim IG2\left(s_{A.n}, \underline{\nu}_A\right)\\ \text{local scale:}&&s_{A.n} | s_{A} &\sim G\left(s_{A}, \underline{a}_A\right)\\ \text{global scale:}&&s_{A} &\sim IG2\left(\underline{s}_{s_A}, \underline{\nu}_{s_A}\right) \end{align}\]
- Minnesota prior mean and shrinkage decay with increasing lags
- Flexibility in shrinkage and scale hyper-parameters
- 3-level equation-specific local-global hierarchical prior
bsvars models
Structural form hierarchical prior
\[\begin{align} \text{exclusion restrictions:}&& [\mathbf{B}_0]_{n\cdot} = \mathbf{b}_n\mathbf{V}_n\\ \text{structural relations:}&& \mathbf{B}_0\mid\gamma_{B.1},\dots,\gamma_{B.N} &\sim |\det(\mathbf{B}_0)|^{\underline{\nu}_B - N}\exp\left\{-\frac{1}{2} \sum_{n=1}^{N} \gamma_{B.n}^{-1} \mathbf{b}_n\mathbf{b}_n' \right\} \\ \text{structural shrinkage:}&&\gamma_{B.n} | s_{B.n} &\sim IG2\left(s_{B.n}, \underline{\nu}_b\right)\\ \text{local scale:}&&s_{B.n} | s_{B} &\sim G\left(s_{B}, \underline{a}_B\right)\\ \text{global scale:}&&s_{B} &\sim IG2\left(\underline{s}_{s_B}, \underline{\nu}_{s_B}\right) \end{align}\]
- Highly adaptive equation-by-equation exclusion restrictions
- Likelihood-shape preserving prior
- Flexibility in shrinkage and scale hyper-parameters
- 3-level equation-specific local-global hierarchical prior
bsvars models
\[ \]
Volatility models
- Homoskedastic \(\sigma_{n.t}^2 = 1\)
- Stochastic Volatility: non-centred and centred
- Stationary Markov-switching heteroskedastisity
- Sparse Markov-switching heteroskedastisity
Non-normal models
- Finite mixture of normal components
- Sparse mixture of normal components
bsvars models
Non-centred Stochastic Volatility
\[\begin{align} \text{conditional variance:}&&\sigma_{n.t}^2 &= \exp\left\{\omega_n h_{n.t}\right\}\\ \text{log-volatility:}&&h_{n.t} &= \rho_n h_{n.t-1} + v_{n.t}\\ \text{volatility innovation:}&&\tilde{v}_{n.t}&\sim N\left(0,1\right)\\ \text{autoregression:}&&\rho_n &\sim U(-1,1)\\ \text{volatility of the log-volatility:}&&\omega_n\mid \sigma_{\omega.n}^2 &\sim N\left(0, \sigma_{\omega.n}^2\right)\\ \text{hierarchy:}&&\sigma_{\omega.n}^2 \mid s_\sigma&\sim G(s_\sigma, \underline{a}_\sigma) \end{align}\]
- excellent volatility forecasting performance
- standardization around \(\sigma_{n.t}^2 = 1\)
- homoskedasticity verification by testing \(\omega_n = 0\)
bsvars models
Centred Stochastic Volatility
\[\begin{align} \text{conditional variance:}&&\sigma_{n.t}^2 &= \exp\left\{\tilde{h}_{n.t}\right\}\\ \text{log-volatility:}&&\tilde{h}_{n.t} &= \rho_n \tilde{h}_{n.t-1} + \tilde{v}_{n.t}\\ \text{volatility innovation:}&&\tilde{v}_{n.t}&\sim N\left(0,\sigma_v^2\right)\\ \text{autoregression:}&&\rho_n &\sim U(-1,1)\\ \text{volatility of the log-volatility:}&&\sigma_v^2 \mid s_v &\sim IG2(s_v, \underline{a}_v)\\ \text{hierarchy:}&&s_v \mid s_\sigma &\sim G(s_\sigma, \underline{a}_\sigma) \end{align}\]
- excellent volatility forecasting performance
- weak standardisation
- no homoskedasticity verification available
bsvars models
Stochastic Volatility: conditional variances
bsvars identification (simplified)
\[\begin{align} &\\ \mathbf\Sigma &= \mathbf{B}_0^{-1}\mathbf{B}_0^{-1\prime}\\[1ex] \end{align}\]
- \(\mathbf\Sigma\) can be estimated using data easily
- The relationship presents a system of equations to be solved for \(\mathbf{B}_0\)
- \(\mathbf\Sigma\) is a symmetric \(N\times N\) matrix
- \(\mathbf\Sigma\) has \(N(N+1)/2\) unique elements (equations)
- \(\mathbf{B}_0\) is an \(N\times N\) matrix with \(N^2\) unique elements to estimate
- We cannot estimate all elements of \(\mathbf{B}_0\) using \(N(N+1)/2\) equations
- \(\mathbf{B}_0\) is
not identified
bsvars identification (simplified)
\[\begin{align} &\\ \mathbf\Sigma &= \mathbf{B}_0^{-1}\mathbf{B}_0^{-1\prime}\\[1ex] \end{align}\]
Identification.
- Only \(N(N+1)/2\) elements in \(\mathbf{B}_0\) can be estimated
- Impose \(N(N-1)/2\) restrictions on \(\mathbf{B}_0\) to solve the equation
- This identifies the rows of \(\mathbf{B}_0\) up to a sign
- Change the sign of any number of \(\mathbf{B}_0\) rows and \(\mathbf\Sigma\) will not change
- Often \(\mathbf{B}_0\) is made lower-triangular
bsvars identification (simplified)
Let \(N=2\)
\[\begin{align} \begin{bmatrix}\sigma_1^2&\sigma_{12}\\ \sigma_{12}&\sigma_2^2\end{bmatrix} &\qquad \begin{bmatrix}B_{0.11}&B_{0.12}\\ B_{0.21}&B_{0.22}\end{bmatrix}\\[1ex] \end{align}\]
- 3 unique elements in \(\mathbf\Sigma\) - 3 equations in the system
- 4 elements in \(\mathbf{B}_0\) cannot be estimated
Identification.
\[\begin{align} \begin{bmatrix}\sigma_1^2&\sigma_{12}\\ \sigma_{12}&\sigma_2^2\end{bmatrix} &\qquad \begin{bmatrix}B_{0.11}& 0\\ B_{0.21}&B_{0.22}\end{bmatrix}\\[1ex] \end{align}\]
- 3 equations identify 3 elements in \(\mathbf{B}_0\)
bsvars identification (simplified)
Identification through Heteroskedasticity.
Suppose that:
- there are two covariances, \(\mathbf\Sigma_1\) and \(\mathbf\Sigma_2\), associated with the sample
- matrix \(\mathbf{B}_0\) does not change over time
- structural shocks are heteroskedastic with covariances \(\text{diag}\left(\boldsymbol\sigma_1^2\right)\) and \(\text{diag}\left(\boldsymbol\sigma_2^2\right)\)
\[\begin{align} \mathbf\Sigma_1 &= \mathbf{B}_0^{-1}\text{diag}\left(\boldsymbol\sigma_1^2\right)\mathbf{B}_0^{-1\prime}\\[1ex] \mathbf\Sigma_2 &= \mathbf{B}_0^{-1}\text{diag}\left(\boldsymbol\sigma_2^2\right)\mathbf{B}_0^{-1\prime} \end{align}\]
bsvars identification (simplified)
Identification through Heteroskedasticity.
\[\begin{align} \mathbf\Sigma_1 &= \mathbf{B}_0^{-1}\text{diag}\left(\boldsymbol\sigma_1^2\right)\mathbf{B}_0^{-1\prime}\\[1ex] \mathbf\Sigma_2 &= \mathbf{B}_0^{-1}\text{diag}\left(\boldsymbol\sigma_2^2\right)\mathbf{B}_0^{-1\prime} \end{align}\]
Identification.
- \(\mathbf\Sigma_1\) and \(\mathbf\Sigma_2\) contain \(N^2+N\) unique elements
- All \(N^2\) elements of \(\mathbf{B}_0\) can be estimated
- Both \(N\)-vectors \(\boldsymbol\sigma_1^2\) and \(\boldsymbol\sigma_2^2\) can be estimated due to additional restriction: \(E\left[\text{diag}\left(\boldsymbol\sigma_i^2\right)\right] = \mathbf{I}_N\)
bsvars identification (simplified)
The setup can be generalised to conditional heteroskedasticity of structural shocks
\[\begin{align} u_t |Y_{t-1} &\sim N_N\left(\mathbf{0}_N, \text{diag}\left(\boldsymbol\sigma_t^2\right)\right)\\[1ex] \mathbf\Sigma_t &= \mathbf{B}_0^{-1}\text{diag}\left(\boldsymbol\sigma_t^2\right)\mathbf{B}_0^{-1\prime} \end{align}\]
Identification.
- Matrix \(\mathbf{B}_0\) is identified up to its rows’ sign change and equations’ reordering
- shocks are identified if changes in their conditional variances are non-proportional
- Structural shocks’ conditional variances \(\boldsymbol\sigma_t^2\) can be estimated
Heteroskedasticity Modeling.
Choose any (conditional) variance model for \(\boldsymbol\sigma_t^2\) that fits the data well.
bsvars identification verification
Non-centred Stochastic Volatility.
A structural shock is homoskedastic if \[ \omega_n = 0 \]
\[ \]
- if \(\omega_n = 0\) the shock is homoskedastic with \(\sigma_{nt}^2 = 1\)
- the only homoskedastic shock in the system is identified
- two or more homoskedastic shocks have not identified through heteroskedasticity
- heteroskedastic shocks are identified with probability 1
bsvars identification verification
Savage-Dickey Density Ratio.
Verify the restriction through the posterior odds ratio using the SDDR: \[ SDDR = \frac{\Pr[\omega_n = 0 | data]}{\Pr[\omega_n \neq 0 | data]}= \frac{p(\omega_n = 0 | data)}{p(\omega_n = 0 )} \]
Estimate the marginal prior and posterior ordinate using analytical or numerical integration: \[\begin{align} \hat{p}(\omega_n = 0 | data) &= \frac{1}{S}\sum_{s=1}^S p(\omega_n = 0 | data, \mathbf{A}^{(s)}, \mathbf{B}_0^{(s)}, \dots)\\ \lim\limits_{\omega_n\rightarrow 0}\hat{p}(\omega_n = 0) &= \Gamma\left(\underline{a}_{\sigma}+\frac{3}{2}\right)\left[\sqrt{2\pi s_{\sigma}}\left(\underline{a}_{\sigma}^2-\frac{1}{4}\right)\Gamma\left(\underline{a}_{\sigma}\right)\right]^{-1} \end{align}\]
bsvars modeling of monetary policy
bsvars modeling of monetary policy
Consider a system of four variables:
\[\begin{align} y_t = \begin{bmatrix} \Delta rgdp_t & \pi_t & CR_t & EX_t \end{bmatrix}' \end{align}\]
- \(\Delta rgdp_t\) - real Gross Domestic Product growth
- \(\pi_t\) - Consumer Price Index inflation rate
- \(CR_t\) - Cash Rate Target - Australian nominal interest rate
- \(EX_t\) - USD/AUD exchange rate
\[ \]
- monthly data from July 1969 to September 2023
- quarterly variables interpolated to monthly frequency
bsvars modeling of monetary policy
bsvars modeling of monetary policy
Foreign sector.
- \(tot_t^{(AUD)}\) - Australian terms of trade
- \(\Delta rgdp_t^{(USA)}\) - U.S. real Gross Domestic Product growth
- \(\pi_t^{(USA)}\) - U.S. CPI inflation rate
- \(FFR_t\) - Federal Funds Rate
\[ \]
- quarterly variables interpolated to monthly frequency
- contemporaneous values and the first lag of the foreign sector variables are included as exogenous variables
bsvars modeling of monetary policy
bsvars modeling of monetary policy
Lower-triangular system.
\[\begin{align} \begin{bmatrix} B_{0.11}&0&0&0\\ B_{0.21}&B_{0.22}&0&0\\ B_{0.31}&B_{0.32}&B_{0.33}&0\\ B_{0.41}&B_{0.42}&B_{0.43}&B_{0.44} \end{bmatrix} \begin{bmatrix} \Delta rgdp_t \\ \pi_t \\ CR_t \\ EX_t \end{bmatrix} &= \dots + \begin{bmatrix} u_t^{ad} \\ u_t^{as} \\ u_t^{mps} \\ u_t^{ex} \end{bmatrix} \end{align}\]
Identified shocks.
- \(u_t^{mps}\) - monetary policy shock identified via Taylor’s Rule
- \(u_t^{ex}\) - currency shock
bsvars modeling of monetary policy
Extended system.
\[\begin{align} \begin{bmatrix} B_{0.11}&0&0&0\\ B_{0.21}&B_{0.22}&0&0\\ B_{0.31}&B_{0.32}&B_{0.33}&B_{0.34}\\ B_{0.41}&B_{0.42}&B_{0.43}&B_{0.44} \end{bmatrix} \begin{bmatrix} \Delta rgdp_t \\ \pi_t \\ CR_t \\ EX_t \end{bmatrix} &= \dots + \begin{bmatrix} u_t^{ad} \\ u_t^{as} \\ u_t^{mps} \\ u_t^{ex} \end{bmatrix} \end{align}\]
Identified shocks.
- \(u_t^{mps}\) - identified via Taylor’s Rule extended by exchange rate
- \(u_t^{ex}\) - currency shock indistinguishable from the monetary policy shock
- identification through heteroskedasticity helps identifying the shocks
bsvars modeling of monetary policy
Unrestricted system.
\[\begin{align} \begin{bmatrix} B_{0.11}&B_{0.12}&B_{0.13}&B_{0.14}\\ B_{0.21}&B_{0.22}&B_{0.23}&B_{0.24}\\ B_{0.31}&B_{0.32}&B_{0.33}&B_{0.34}\\ B_{0.41}&B_{0.42}&B_{0.43}&B_{0.44} \end{bmatrix} \begin{bmatrix} \Delta rgdp_t \\ \pi_t \\ CR_t \\ EX_t \end{bmatrix} &= \dots + \begin{bmatrix} u_{1.t} \\ u_{2.t} \\ u_{3.t} \\ {4.t} \end{bmatrix} \end{align}\]
Identified shocks.
- identification through heteroskedasticity required
- can shocks be labelled?
bsvars modeling of monetary policy
General setup.
bsvars modeling of monetary policy
Specify the lower-triangular SVAR-SV.
spec_lt = specify_bsvar_sv$new( # specify an SVAR-SV
as.matrix(y), # endogenous variables
p = p, # lag order
exogenous = x # exogenous variables
)
A_ols = t(solve( # compute A_ols
tcrossprod(spec_lt$data_matrices$X),
tcrossprod(spec_lt$data_matrices$X, spec_lt$data_matrices$Y)
))
spec_lt$prior$A = A_ols # prior mean of A set to A_ols
spec_lt<BSVARSV>
Public:
centred_sv: FALSE
clone: function (deep = FALSE)
data_matrices: DataMatricesBSVAR, R6
get_data_matrices: function ()
get_identification: function ()
get_prior: function ()
get_starting_values: function ()
identification: IdentificationBSVARs, R6
initialize: function (data, p = 1L, B, exogenous = NULL, centred_sv = FALSE,
p: 12
prior: PriorBSVARSV, PriorBSVAR, R6
starting_values: StartingValuesBSVARSV, StartingValuesBSVAR, R6bsvars modeling of monetary policy
Specify the extended SVAR-SV.
[,1] [,2] [,3] [,4]
[1,] TRUE FALSE FALSE FALSE
[2,] TRUE TRUE FALSE FALSE
[3,] TRUE TRUE TRUE TRUE
[4,] TRUE TRUE TRUE TRUEbsvars modeling of monetary policy
Specify the unrestricted SVAR-SV.
bsvars modeling of monetary policy
Estimate the lower-triangular SVAR-SV.
**************************************************|
bsvars: Bayesian Structural Vector Autoregressions|
**************************************************|
Gibbs sampler for the SVAR-SV model |
Non-centred SV model is estimated |
**************************************************|
Progress of the MCMC simulation for 50 draws
Every 2nd draw is saved via MCMC thinning
Press Esc to interrupt the computations
**************************************************|bsvars modeling of monetary policy
Estimate the extended SVAR-SV.
**************************************************|
bsvars: Bayesian Structural Vector Autoregressions|
**************************************************|
Gibbs sampler for the SVAR-SV model |
Non-centred SV model is estimated |
**************************************************|
Progress of the MCMC simulation for 50 draws
Every 2nd draw is saved via MCMC thinning
Press Esc to interrupt the computations
**************************************************|bsvars modeling of monetary policy
Estimate the unrestricted SVAR-SV.
**************************************************|
bsvars: Bayesian Structural Vector Autoregressions|
**************************************************|
Gibbs sampler for the SVAR-SV model |
Non-centred SV model is estimated |
**************************************************|
Progress of the MCMC simulation for 50 draws
Every 2nd draw is saved via MCMC thinning
Press Esc to interrupt the computations
**************************************************|bsvars modeling of monetary policy
Extended SVAR-SV results.
bsvars modeling of monetary policy
Extended SVAR-SV results.
bsvars modeling of monetary policy
Extended SVAR-SV results.
bsvars modeling of monetary policy
Extended SVAR-SV results.
bsvars modeling of monetary policy
Extended SVAR-SV results.
Verify if shocks are homoskedastic.
\[ SDDR = \frac{p(\omega_n = 0\mid data)}{p(\omega_n = 0)} \]
post_ex |> verify_volatility() -> sddr
vv = cbind(sddr$logSDDR, 1 / (1 + exp(sddr$logSDDR)))
colnames(vv) = c("log(SDDR)", "Pr[w_n != 0|data]")
rownames(vv) = paste0("shock ", 1:N)
round(vv, 3) log(SDDR) Pr[w_n != 0|data]
shock 1 -242.729 1.000
shock 2 -218.592 1.000
shock 3 -3.973 0.982
shock 4 -244.784 1.000bsvars modeling of monetary policy
Extended SVAR-SV results.
bsvars modeling of monetary policy
Extended SVAR-SV results.
bsvars modeling of monetary policy
Extended SVAR-SV results.
Verify if exogenous variables are needed with SDDR.
\[ SDDR = \frac{p(\mathbf{A}_d = \mathbf{0}_{N \times d}\mid data)}{p(\mathbf{A}_d = \mathbf{0}_{N \times d})} \]
A0 = array(NA, dim(A_ols))
A0[,50:57] = 0
post_ex |> verify_autoregression(hypothesis = A0) -> sddr_Ad
vv = cbind(sddr_Ad$logSDDR, 1 / (1 + exp(sddr_Ad$logSDDR)))
colnames(vv) = c("log(SDDR)", "Pr[A_d != 0|data]")
rownames(vv) = paste0("H0: A_d = 0")
round(vv, 3) log(SDDR) Pr[A_d != 0|data]
H0: A_d = 0 -1.346985e+15 1bsvars modeling of monetary policy
Extended SVAR-SV results.
Verify if \(EX_t\) Granger causes \(\Delta rgdp_t\).
\[ SDDR = \frac{p(\mathbf{A}_{1.14} = \dots = \mathbf{A}_{12.14} = 0\mid data)}{p(\mathbf{A}_{1.14} = \dots = \mathbf{A}_{12.14} = 0)} \]
A0 = array(NA, dim(A_ols))
A0[1,seq(from = 4, to = 48, by = 4)] = 0
post_ex |> verify_autoregression(hypothesis = A0) -> sddr_Gc
vv = cbind(sddr_Gc$logSDDR, 1 / (1 + exp(sddr_Gc$logSDDR)))
colnames(vv) = c("log(SDDR)", "Pr[A_{i.14} != 0|data] for all i")
rownames(vv) = paste0("H0: A_{i.14} = 0")
round(vv, 3) log(SDDR) Pr[A_{i.14} != 0|data] for all i
H0: A_{i.14} = 0 -2.835044e+16 1bsvars modeling of monetary policy
Lower-triangular SVAR-SV results.
bsvars modeling of monetary policy
Unrestricted SVAR-SV results.
