From b2e07fff496550a120784f8a1425a8266bb4f587 Mon Sep 17 00:00:00 2001 From: Stefan Haan Date: Tue, 3 Dec 2024 21:21:18 +0100 Subject: [PATCH 1/4] shocks can only be applied to factors --- R/RcppExports.R | 4 ++++ R/irf.R | 20 ++++++++++++++++++ R/utility_functions.R | 9 ++++----- src/RcppExports.cpp | 15 ++++++++++++++ src/irf.cpp | 47 +++++++++++++++++++++++++++++++++++++++++++ 5 files changed, 90 insertions(+), 5 deletions(-) create mode 100644 R/irf.R create mode 100644 src/irf.cpp diff --git a/R/RcppExports.R b/R/RcppExports.R index cf9c280..daf777d 100644 --- a/R/RcppExports.R +++ b/R/RcppExports.R @@ -23,6 +23,10 @@ my_gig <- function(n, lambda, chi, psi) { .Call(`_bayesianVARs_my_gig`, n, lambda, chi, psi) } +irf_cpp <- function(coefficients, factor_loadings, f_shock, ahead) { + .Call(`_bayesianVARs_irf_cpp`, coefficients, factor_loadings, f_shock, ahead) +} + sample_PHI_cholesky <- function(PHI, PHI_prior, Y, X, U, d_sqrt, V_prior) { .Call(`_bayesianVARs_sample_PHI_cholesky`, PHI, PHI_prior, Y, X, U, d_sqrt, V_prior) } diff --git a/R/irf.R b/R/irf.R new file mode 100644 index 0000000..f02979d --- /dev/null +++ b/R/irf.R @@ -0,0 +1,20 @@ +# x should be of class bayesianVARs_bvar +irf <- function(x, factor_shock, ahead=8) { + if (x$sigma_type != "factor") { + stop("impulse reponse functions are only available for factor models") + } + + number_of_factors <- dim(x$facload)[2] + if (length(factor_shock) != number_of_factors) { + stop("the shock vector to the factors must have dimension equals to the number of factors") + } + + ret <- irf_cpp( + x$PHI, + x$facload, + factor_shock, + ahead + ) + dimnames(ret)[[2]] <- dimnames(x$Y)[[2]] + ret +} diff --git a/R/utility_functions.R b/R/utility_functions.R index f370a3d..83e03ab 100644 --- a/R/utility_functions.R +++ b/R/utility_functions.R @@ -1733,7 +1733,6 @@ predict.bayesianVARs_bvar <- function(object, ahead = 1L, each = 1L, stable = TR # relevant mod settings sv_indicator <- which(object$heteroscedastic==TRUE) - intercept <- object$intercept # data preparation variables <- colnames(object$Y) @@ -2038,10 +2037,10 @@ predict_old <- function(object, n.ahead, stable = TRUE, LPL = FALSE, Y_obs = NA, #' predictions <- predict(mod, ahead = 1L) #' print(predictions) print.bayesianVARs_predict <- function(x, ...){ - cat(paste("\nGeneric functions for bayesianVARs_predict objects:\n", - " - summary.bayesianVARs_predict(),\n", - " - pairs.bayesianVARs_predict(),\n", - " - plot.bayesianVARs_predict() (alias for pairs.bayesianVARs_predict()).\n")) + cat("\nGeneric functions for bayesianVARs_predict objects:\n", + " - summary.bayesianVARs_predict(),\n", + " - pairs.bayesianVARs_predict(),\n", + " - plot.bayesianVARs_predict() (alias for pairs.bayesianVARs_predict()).\n") invisible(x) } diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index 98fffac..7c309b9 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -56,6 +56,20 @@ BEGIN_RCPP return rcpp_result_gen; END_RCPP } +// irf_cpp +arma::cube irf_cpp(const arma::cube& coefficients, const arma::cube& factor_loadings, const arma::colvec& f_shock, const arma::uword ahead); +RcppExport SEXP _bayesianVARs_irf_cpp(SEXP coefficientsSEXP, SEXP factor_loadingsSEXP, SEXP f_shockSEXP, SEXP aheadSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< const arma::cube& >::type coefficients(coefficientsSEXP); + Rcpp::traits::input_parameter< const arma::cube& >::type factor_loadings(factor_loadingsSEXP); + Rcpp::traits::input_parameter< const arma::colvec& >::type f_shock(f_shockSEXP); + Rcpp::traits::input_parameter< const arma::uword >::type ahead(aheadSEXP); + rcpp_result_gen = Rcpp::wrap(irf_cpp(coefficients, factor_loadings, f_shock, ahead)); + return rcpp_result_gen; +END_RCPP +} // sample_PHI_cholesky arma::mat sample_PHI_cholesky(const arma::mat PHI, const arma::mat& PHI_prior, const arma::mat& Y, const arma::mat& X, const arma::mat& U, const arma::mat& d_sqrt, const arma::mat& V_prior); RcppExport SEXP _bayesianVARs_sample_PHI_cholesky(SEXP PHISEXP, SEXP PHI_priorSEXP, SEXP YSEXP, SEXP XSEXP, SEXP USEXP, SEXP d_sqrtSEXP, SEXP V_priorSEXP) { @@ -167,6 +181,7 @@ END_RCPP static const R_CallMethodDef CallEntries[] = { {"_bayesianVARs_bvar_cpp", (DL_FUNC) &_bayesianVARs_bvar_cpp, 21}, {"_bayesianVARs_my_gig", (DL_FUNC) &_bayesianVARs_my_gig, 4}, + {"_bayesianVARs_irf_cpp", (DL_FUNC) &_bayesianVARs_irf_cpp, 4}, {"_bayesianVARs_sample_PHI_cholesky", (DL_FUNC) &_bayesianVARs_sample_PHI_cholesky, 7}, {"_bayesianVARs_dmvnrm_arma_fast", (DL_FUNC) &_bayesianVARs_dmvnrm_arma_fast, 4}, {"_bayesianVARs_predh", (DL_FUNC) &_bayesianVARs_predh, 6}, diff --git a/src/irf.cpp b/src/irf.cpp new file mode 100644 index 0000000..4bb4531 --- /dev/null +++ b/src/irf.cpp @@ -0,0 +1,47 @@ +#include + +using namespace Rcpp; +using namespace arma; + +inline void shift_and_insert( + arma::mat& X, //the columns of X should be y1,y2,y3, y1.l1,y2.l1,y3.l1,...,1 + const arma::mat& new_y //what to insert in y1,y2,y3 +) { + for (uword i = X.n_cols-2; new_y.n_cols <= i; i--) { + X.col(i) = X.col(i-new_y.n_cols); + } + X.head_cols(new_y.n_cols) = new_y; +} + +// [[Rcpp::export]] +arma::cube irf_cpp( + const arma::cube& coefficients, //rows: lagged variables + intercept, columns: variables, slices: draws + const arma::cube& factor_loadings, //rows: variables, columns: factors, slices: draws + const arma::colvec& f_shock, //columns: how much each of the factors is shocked + const arma::uword ahead //how far to predict ahead +) { + const uword n_variables = coefficients.n_cols; + const uword n_posterior_draws = coefficients.n_slices; + arma::cube irf(ahead+1, n_variables, n_posterior_draws, arma::fill::none); + + // trace out n_posterior_draws paths + arma::mat current_predictors(n_posterior_draws, coefficients.n_rows, arma::fill::zeros); + // all paths start with the some shock at t=0 to the variables + // the shock is to the factors, so the uncertainty of the factor loadings + // translates to uncertainty to the shock to the variables + for (uword r = 0; r < n_posterior_draws; r++) { + const rowvec y_shock = (factor_loadings.slice(r) * f_shock).t(); + irf.slice(r).row(0) = y_shock; + current_predictors.head_cols(n_variables).row(r) = y_shock; + } + + for (uword t = 1; t <= ahead; t++) { + for(uword r = 0; r < n_posterior_draws; r++) { + irf.slice(r).row(t) = current_predictors.row(r) * coefficients.slice(r); + } + // shift everything and make predictions the new predictors at lag zero + shift_and_insert(current_predictors, irf.row_as_mat(t)); + } + + return irf; +} From 3f4587fcef60efee6b9d801fbacd80ad06342aa0 Mon Sep 17 00:00:00 2001 From: Stefan Haan Date: Wed, 4 Dec 2024 20:00:23 +0100 Subject: [PATCH 2/4] plotting of bayesianVARs_irf objects --- DESCRIPTION | 2 +- NAMESPACE | 2 ++ R/irf.R | 12 +++++++-- R/plots.R | 75 +++++++++++++++++++++++++++++++++++++++++++++++++++++ man/irf.Rd | 14 ++++++++++ 5 files changed, 102 insertions(+), 3 deletions(-) create mode 100644 man/irf.Rd diff --git a/DESCRIPTION b/DESCRIPTION index 3f400b3..1f18a94 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -49,4 +49,4 @@ Config/testthat/edition: 3 Encoding: UTF-8 LazyData: true Roxygen: list(markdown = TRUE) -RoxygenNote: 7.3.1 +RoxygenNote: 7.3.2 diff --git a/NAMESPACE b/NAMESPACE index e711666..2031efb 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -7,6 +7,7 @@ S3method(fitted,bayesianVARs_bvar) S3method(pairs,bayesianVARs_predict) S3method(plot,bayesianVARs_bvar) S3method(plot,bayesianVARs_fitted) +S3method(plot,bayesianVARs_irf) S3method(plot,bayesianVARs_predict) S3method(predict,bayesianVARs_bvar) S3method(print,bayesianVARs_bvar) @@ -18,6 +19,7 @@ S3method(summary,bayesianVARs_draws) S3method(summary,bayesianVARs_predict) S3method(vcov,bayesianVARs_bvar) export(bvar) +export(irf) export(my_gig) export(posterior_heatmap) export(specify_prior_phi) diff --git a/R/irf.R b/R/irf.R index f02979d..9f84098 100644 --- a/R/irf.R +++ b/R/irf.R @@ -1,4 +1,10 @@ -# x should be of class bayesianVARs_bvar +#' Impulse response functions +#' +#' Effect of a shock on the variables over time +#' +#' @param x An object of type `bayesianVARs_bvar` +#' +#' @export irf <- function(x, factor_shock, ahead=8) { if (x$sigma_type != "factor") { stop("impulse reponse functions are only available for factor models") @@ -15,6 +21,8 @@ irf <- function(x, factor_shock, ahead=8) { factor_shock, ahead ) - dimnames(ret)[[2]] <- dimnames(x$Y)[[2]] + + colnames(ret) <- colnames(x$Y) + class(ret) <- "bayesianVARs_irf" ret } diff --git a/R/plots.R b/R/plots.R index 102a445..862ef66 100644 --- a/R/plots.R +++ b/R/plots.R @@ -623,3 +623,78 @@ plot.bayesianVARs_predict <- function(x, dates = NULL, vars = "all", ahead = NUL invisible(x) } + +#' @export +plot.bayesianVARs_irf <- function(x, vars = "all", + quantiles = c(0.05,0.25,0.5,0.75,0.95), + n_col = 1L, + ...){ + + n_ahead <- nrow(x) + + if(length(vars)==1L & any(vars == "all")){ + vars <- 1:ncol(x) + }else if(is.character(vars)){ + if(any(base::isFALSE(vars %in% colnames(x)))){ + stop("Elements of 'vars' must coincide with 'colnames(x)'!") + }else{ + vars <- which(colnames(x) %in% vars) + } + }else if(is.numeric(vars)){ + vars <- as.integer(vars) + if(any(vars > ncol(x))){ + stop("'max(vars)' must be at most 'ncol(x)'!") + } + } + var_names <- colnames(x) + + t <- seq_len(n_ahead) + dates <- t-1 + + quantiles <- sort(quantiles) + nr_quantiles <- length(quantiles) + nr_intervals <- floor(nr_quantiles/2) + even <- nr_quantiles%%2 == 0 + + pred_quants <- apply(x, 1:2, quantile, quantiles) + + oldpar <- par(no.readonly = TRUE) + on.exit(par(oldpar), add = TRUE) + M <- length(vars) + par(mfrow=c(min(5,ceiling(M/n_col)),n_col), mar = c(2,2,2,1), mgp = c(2,.5,0)) + for(j in seq_along(vars)){ + plot(t, rep(0, n_ahead), type = "n", xlab="", ylab = "", + xaxt="n", ylim = range(pred_quants[,,vars[j]])) + + axis(side=1, at = t, labels = dates[t]) + mtext(var_names[j], side = 3) + + if(nr_intervals>0){ + for(r in seq.int(nr_intervals)){ + alpha_upper <- 0.4 + alpha_lower <- 0.2 + alphas <- seq(alpha_lower, alpha_upper, length.out = nr_intervals) + if(nr_intervals==1){ + alphas <- alpha_upper + } + polygon(x = c(t, rev(t)), + y = c(pred_quants[r,,vars[j]], rev(pred_quants[r+1,,vars[j]])), + col = scales::alpha("red", alphas[r]), + border = NA) + if(length(quantiles)>2){ + polygon(x = c(t, rev(t)), + y = c(pred_quants[nrow(pred_quants)+1-r,,vars[j]], + rev(pred_quants[nrow(pred_quants)-r,,vars[j]])), + col = scales::alpha("red", alphas[r]), + border = NA) + } + } + if(!even){ + lines(t, pred_quants[ceiling(length(quantiles)/2),,vars[j]], + col = "red", lwd = 2) + } + } + } + + invisible(x) +} diff --git a/man/irf.Rd b/man/irf.Rd new file mode 100644 index 0000000..48b84b5 --- /dev/null +++ b/man/irf.Rd @@ -0,0 +1,14 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/irf.R +\name{irf} +\alias{irf} +\title{Impulse response functions} +\usage{ +irf(x, factor_shock, ahead = 8) +} +\arguments{ +\item{x}{An object of type \code{bayesianVARs_bvar}} +} +\description{ +Effect of a shock on the variables over time +} From b192f1f94daec647af5995fbffc6ba0a0629e852 Mon Sep 17 00:00:00 2001 From: Stefan Haan Date: Wed, 4 Dec 2024 21:45:37 +0100 Subject: [PATCH 3/4] implement irf for cholesky model --- R/RcppExports.R | 4 ++-- R/irf.R | 23 +++++++++++++++-------- src/RcppExports.cpp | 11 ++++++----- src/irf.cpp | 22 +++++++++++++++++----- 4 files changed, 40 insertions(+), 20 deletions(-) diff --git a/R/RcppExports.R b/R/RcppExports.R index daf777d..313348a 100644 --- a/R/RcppExports.R +++ b/R/RcppExports.R @@ -23,8 +23,8 @@ my_gig <- function(n, lambda, chi, psi) { .Call(`_bayesianVARs_my_gig`, n, lambda, chi, psi) } -irf_cpp <- function(coefficients, factor_loadings, f_shock, ahead) { - .Call(`_bayesianVARs_irf_cpp`, coefficients, factor_loadings, f_shock, ahead) +irf_cpp <- function(coefficients, factor_loadings, U_vecs, shock, ahead) { + .Call(`_bayesianVARs_irf_cpp`, coefficients, factor_loadings, U_vecs, shock, ahead) } sample_PHI_cholesky <- function(PHI, PHI_prior, Y, X, U, d_sqrt, V_prior) { diff --git a/R/irf.R b/R/irf.R index 9f84098..c85c817 100644 --- a/R/irf.R +++ b/R/irf.R @@ -5,20 +5,27 @@ #' @param x An object of type `bayesianVARs_bvar` #' #' @export -irf <- function(x, factor_shock, ahead=8) { - if (x$sigma_type != "factor") { - stop("impulse reponse functions are only available for factor models") +irf <- function(x, shock, ahead=8) { + if (x$sigma_type == "factor") { + number_of_factors <- dim(x$facload)[2] + if (length(shock) != number_of_factors) { + stop("the shock vector must have dimension equals to the number of factors") + } } - - number_of_factors <- dim(x$facload)[2] - if (length(factor_shock) != number_of_factors) { - stop("the shock vector to the factors must have dimension equals to the number of factors") + else if (x$sigma_type == "cholesky") { + if (length(shock) != ncol(mod$PHI)) { + stop("the shock vector must have dimension equals to the number of variables") + } + } + else { + stop("unknown model") } ret <- irf_cpp( x$PHI, x$facload, - factor_shock, + x$U, + shock, ahead ) diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index 7c309b9..68cc4c5 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -57,16 +57,17 @@ BEGIN_RCPP END_RCPP } // irf_cpp -arma::cube irf_cpp(const arma::cube& coefficients, const arma::cube& factor_loadings, const arma::colvec& f_shock, const arma::uword ahead); -RcppExport SEXP _bayesianVARs_irf_cpp(SEXP coefficientsSEXP, SEXP factor_loadingsSEXP, SEXP f_shockSEXP, SEXP aheadSEXP) { +arma::cube irf_cpp(const arma::cube& coefficients, const arma::cube& factor_loadings, const arma::mat& U_vecs, const arma::colvec& shock, const arma::uword ahead); +RcppExport SEXP _bayesianVARs_irf_cpp(SEXP coefficientsSEXP, SEXP factor_loadingsSEXP, SEXP U_vecsSEXP, SEXP shockSEXP, SEXP aheadSEXP) { BEGIN_RCPP Rcpp::RObject rcpp_result_gen; Rcpp::RNGScope rcpp_rngScope_gen; Rcpp::traits::input_parameter< const arma::cube& >::type coefficients(coefficientsSEXP); Rcpp::traits::input_parameter< const arma::cube& >::type factor_loadings(factor_loadingsSEXP); - Rcpp::traits::input_parameter< const arma::colvec& >::type f_shock(f_shockSEXP); + Rcpp::traits::input_parameter< const arma::mat& >::type U_vecs(U_vecsSEXP); + Rcpp::traits::input_parameter< const arma::colvec& >::type shock(shockSEXP); Rcpp::traits::input_parameter< const arma::uword >::type ahead(aheadSEXP); - rcpp_result_gen = Rcpp::wrap(irf_cpp(coefficients, factor_loadings, f_shock, ahead)); + rcpp_result_gen = Rcpp::wrap(irf_cpp(coefficients, factor_loadings, U_vecs, shock, ahead)); return rcpp_result_gen; END_RCPP } @@ -181,7 +182,7 @@ END_RCPP static const R_CallMethodDef CallEntries[] = { {"_bayesianVARs_bvar_cpp", (DL_FUNC) &_bayesianVARs_bvar_cpp, 21}, {"_bayesianVARs_my_gig", (DL_FUNC) &_bayesianVARs_my_gig, 4}, - {"_bayesianVARs_irf_cpp", (DL_FUNC) &_bayesianVARs_irf_cpp, 4}, + {"_bayesianVARs_irf_cpp", (DL_FUNC) &_bayesianVARs_irf_cpp, 5}, {"_bayesianVARs_sample_PHI_cholesky", (DL_FUNC) &_bayesianVARs_sample_PHI_cholesky, 7}, {"_bayesianVARs_dmvnrm_arma_fast", (DL_FUNC) &_bayesianVARs_dmvnrm_arma_fast, 4}, {"_bayesianVARs_predh", (DL_FUNC) &_bayesianVARs_predh, 6}, diff --git a/src/irf.cpp b/src/irf.cpp index 4bb4531..ff1a3b6 100644 --- a/src/irf.cpp +++ b/src/irf.cpp @@ -16,21 +16,33 @@ inline void shift_and_insert( // [[Rcpp::export]] arma::cube irf_cpp( const arma::cube& coefficients, //rows: lagged variables + intercept, columns: variables, slices: draws - const arma::cube& factor_loadings, //rows: variables, columns: factors, slices: draws - const arma::colvec& f_shock, //columns: how much each of the factors is shocked + const arma::cube& factor_loadings, //rows: variables, columns: factors, slices: draws, + const arma::mat& U_vecs, //rows: entries of a (variables x variables)-upper triagonal matrix with ones on the diagonals, cols: draws + const arma::colvec& shock, //columns: how much each of the factors is shocked const arma::uword ahead //how far to predict ahead ) { const uword n_variables = coefficients.n_cols; const uword n_posterior_draws = coefficients.n_slices; + const bool is_factor_model = factor_loadings.n_cols > 0; + arma::cube irf(ahead+1, n_variables, n_posterior_draws, arma::fill::none); // trace out n_posterior_draws paths arma::mat current_predictors(n_posterior_draws, coefficients.n_rows, arma::fill::zeros); // all paths start with the some shock at t=0 to the variables - // the shock is to the factors, so the uncertainty of the factor loadings - // translates to uncertainty to the shock to the variables + const arma::uvec upper_indices = arma::trimatu_ind(arma::size(n_variables, n_variables), 1); for (uword r = 0; r < n_posterior_draws; r++) { - const rowvec y_shock = (factor_loadings.slice(r) * f_shock).t(); + rowvec y_shock; + if (is_factor_model) { + y_shock = (factor_loadings.slice(r) * shock).t(); + } + else { + arma::mat U(n_variables, n_variables, arma::fill::eye); + U(upper_indices) = U_vecs.col(r); + U = U.t(); + y_shock = solve(arma::trimatl(U), shock).t(); + } + irf.slice(r).row(0) = y_shock; current_predictors.head_cols(n_variables).row(r) = y_shock; } From 7a0f9a3f9623a60a76b2d0f7d2b44e3873a21041 Mon Sep 17 00:00:00 2001 From: Stefan Haan Date: Thu, 5 Dec 2024 21:46:33 +0100 Subject: [PATCH 4/4] write some documentation for irf --- R/irf.R | 19 +++++++++++++++++-- man/irf.Rd | 25 ++++++++++++++++++++++--- 2 files changed, 39 insertions(+), 5 deletions(-) diff --git a/R/irf.R b/R/irf.R index c85c817..291944f 100644 --- a/R/irf.R +++ b/R/irf.R @@ -1,8 +1,23 @@ #' Impulse response functions #' -#' Effect of a shock on the variables over time +#' Effect of a shock on the factors or variables over time #' -#' @param x An object of type `bayesianVARs_bvar` +#' @param x An object of type `bayesianVARs_bvar`. +#' @param shock A vector of shocks. If a factor model was used, then `shock` is +#' applied to the factors, and its dimension must be the number of factors. +#' +#' If the Cholesky model was used, then the dimension of `shock` is expected +#' to be equals the number of equations. Note that the contemporaneous effects +#' of a shock in the Cholesky model depend on the ordering of equations. +#' +#' @return Returns a `bayesianVARs_irf` object +#' +#' @examples +#' train_data <- 100 * usmacro_growth[,c("GDPC1", "PCECC96", "GPDIC1", "AWHMAN", "GDPCTPI", "CES2000000008x", "FEDFUNDS", "GS10", "EXUSUKx", "S&P 500")] +#' prior_sigma <- specify_prior_sigma(data=train_data, type="cholesky") +#' mod <- bvar(train_data, lags=2L, draws=2000, prior_sigma=prior_sigma) +#' ir <- irf(mod, shock=c(0,0,0,0,0,0,1,0,0,0)) +#' plot(ir, n_col=2) #' #' @export irf <- function(x, shock, ahead=8) { diff --git a/man/irf.Rd b/man/irf.Rd index 48b84b5..dada9fc 100644 --- a/man/irf.Rd +++ b/man/irf.Rd @@ -4,11 +4,30 @@ \alias{irf} \title{Impulse response functions} \usage{ -irf(x, factor_shock, ahead = 8) +irf(x, shock, ahead = 8) } \arguments{ -\item{x}{An object of type \code{bayesianVARs_bvar}} +\item{x}{An object of type \code{bayesianVARs_bvar}.} + +\item{shock}{A vector of shocks. If a factor model was used, then \code{shock} is +applied to the factors, and its dimension must be the number of factors. + +\if{html}{\out{
}}\preformatted{If the Cholesky model was used, then the dimension of `shock` is expected +to be equals the number of equations. Note that the contemporaneous effects +of a shock in the Cholesky model depend on the ordering of equations. +}\if{html}{\out{
}}} +} +\value{ +Returns a \code{bayesianVARs_irf} object } \description{ -Effect of a shock on the variables over time +Effect of a shock on the factors or variables over time +} +\examples{ +train_data <- 100 * usmacro_growth[,c("GDPC1", "PCECC96", "GPDIC1", "AWHMAN", "GDPCTPI", "CES2000000008x", "FEDFUNDS", "GS10", "EXUSUKx", "S&P 500")] +prior_sigma <- specify_prior_sigma(data=train_data, type="cholesky") +mod <- bvar(train_data, lags=2L, draws=2000, prior_sigma=prior_sigma) +ir <- irf(mod, shock=c(0,0,0,0,0,0,1,0,0,0)) +plot(ir, n_col=2) + }