Skip to content

Commit

Permalink
included compute_variance_decomposition #25
Browse files Browse the repository at this point in the history
+ cpp wrapper
+ R method
+ doc
+ test
+ updated seeAlso for estimate
  • Loading branch information
donotdespair committed Jul 16, 2024
1 parent d226173 commit db8e18b
Show file tree
Hide file tree
Showing 15 changed files with 370 additions and 20 deletions.
2 changes: 2 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
@@ -1,9 +1,11 @@
# Generated by roxygen2: do not edit by hand

S3method(compute_conditional_sd,PosteriorBSVARSIGN)
S3method(compute_fitted_values,PosteriorBSVARSIGN)
S3method(compute_historical_decompositions,PosteriorBSVARSIGN)
S3method(compute_impulse_responses,PosteriorBSVARSIGN)
S3method(compute_structural_shocks,PosteriorBSVARSIGN)
S3method(compute_variance_decompositions,PosteriorBSVARSIGN)
S3method(estimate,BSVARSIGN)
S3method(forecast,PosteriorBSVARSIGN)
export(specify_bsvarSIGN)
Expand Down
146 changes: 138 additions & 8 deletions R/compute.R
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@
#' array with attribute \code{PosteriorShocks} containing \code{S} draws of the
#' structural shocks.
#'
#' @seealso \code{\link{estimate}}, \code{\link{summary}}, \code{\link{plot}}
#' @seealso \code{\link{estimate.BSVARSIGN}}, \code{\link{summary}}, \code{\link{plot}}
#'
#' @author Xiaolei Wang \email{[email protected]} and Tomasz Woźniak \email{[email protected]}
#'
Expand All @@ -25,7 +25,7 @@
#'
#' # specify the model and set seed
#' set.seed(123)
#' sign_irf = array(matrix(c(-1, -1, 1, rep(0, 6)), nrow = 3), dim = c(3, 3, 1))
#' sign_irf = array(matrix(c(-1, -1, 1, rep(NA, 6)), nrow = 3), dim = c(3, 3, 1))
#' specification = specify_bsvarSIGN$new(oil, sign_irf = sign_irf)
#'
#' # run the burn-in
Expand Down Expand Up @@ -75,7 +75,7 @@ compute_structural_shocks.PosteriorBSVARSIGN <- function(posterior) {
#' array with attribute \code{PosteriorFitted} containing \code{S} draws from
#' the data predictive density.
#'
#' @seealso \code{\link{estimate}}, \code{\link{summary}}, \code{\link{plot}}
#' @seealso \code{\link{estimate.BSVARSIGN}}, \code{\link{summary}}, \code{\link{plot}}
#'
#' @author Xiaolei Wang \email{[email protected]} and Tomasz Woźniak \email{[email protected]}
#'
Expand All @@ -85,7 +85,7 @@ compute_structural_shocks.PosteriorBSVARSIGN <- function(posterior) {
#'
#' # specify the model and set seed
#' set.seed(123)
#' sign_irf = array(matrix(c(-1, -1, 1, rep(0, 6)), nrow = 3), dim = c(3, 3, 1))
#' sign_irf = array(matrix(c(-1, -1, 1, rep(NA, 6)), nrow = 3), dim = c(3, 3, 1))
#' specification = specify_bsvarSIGN$new(oil, sign_irf = sign_irf)
#'
#' # run the burn-in
Expand Down Expand Up @@ -142,7 +142,7 @@ compute_fitted_values.PosteriorBSVARSIGN <- function(posterior) {
#' @return An object of class PosteriorIR, that is, an \code{NxNx(horizon+1)xS} array with attribute PosteriorIR
#' containing \code{S} draws of the impulse responses.
#'
#' @seealso \code{\link{estimate}}, \code{\link{summary}}, \code{\link{plot}}
#' @seealso \code{\link{estimate.BSVARSIGN}}, \code{\link{summary}}, \code{\link{plot}}
#'
#' @author Xiaolei Wang \email{[email protected]} and Tomasz Woźniak \email{[email protected]}
#'
Expand All @@ -155,7 +155,7 @@ compute_fitted_values.PosteriorBSVARSIGN <- function(posterior) {
#'
#' # specify the model and set seed
#' set.seed(123)
#' sign_irf = array(matrix(c(-1, -1, 1, rep(0, 6)), nrow = 3), dim = c(3, 3, 1))
#' sign_irf = array(matrix(c(-1, -1, 1, rep(NA, 6)), nrow = 3), dim = c(3, 3, 1))
#' specification = specify_bsvarSIGN$new(oil, sign_irf = sign_irf)
#'
#' # run the burn-in
Expand Down Expand Up @@ -214,7 +214,7 @@ compute_impulse_responses.PosteriorBSVARSIGN <- function(posterior, horizon, sta
#' with attribute \code{PosteriorHD} containing \code{S} draws of the historical
#' decompositions.
#'
#' @seealso \code{\link{estimate}}, \code{\link{summary}}, \code{\link{plot}}
#' @seealso \code{\link{estimate.BSVARSIGN}}, \code{\link{summary}}, \code{\link{plot}}
#'
#' @author Xiaolei Wang \email{[email protected]} and Tomasz Woźniak \email{[email protected]}
#'
Expand All @@ -227,7 +227,7 @@ compute_impulse_responses.PosteriorBSVARSIGN <- function(posterior, horizon, sta
#'
#' # specify the model and set seed
#' set.seed(123)
#' sign_irf = array(matrix(c(-1, -1, 1, rep(0, 6)), nrow = 3), dim = c(3, 3, 1))
#' sign_irf = array(matrix(c(-1, -1, 1, rep(NA, 6)), nrow = 3), dim = c(3, 3, 1))
#' specification = specify_bsvarSIGN$new(oil, sign_irf = sign_irf)
#'
#' # run the burn-in
Expand Down Expand Up @@ -270,3 +270,133 @@ compute_historical_decompositions.PosteriorBSVARSIGN <- function(posterior, show

return(hd)
} # END compute_historical_decompositions.PosteriorBSVARSIGN




#' @method compute_variance_decompositions PosteriorBSVARSIGN
#'
#' @title Computes posterior draws of the forecast error variance decomposition
#'
#' @description Each of the draws from the posterior estimation of the model
#' is transformed into a draw from the posterior distribution of the forecast
#' error variance decomposition.
#'
#' @param posterior posterior estimation outcome - an object of class
#' \code{PosteriorBSVARSIGN} obtained by running the \code{estimate} function.
#' @param horizon a positive integer number denoting the forecast horizon for
#' the impulse responses computations.
#'
#' @return An object of class \code{PosteriorFEVD}, that is, an \code{NxNx(horizon+1)xS}
#' array with attribute \code{PosteriorFEVD} containing \code{S} draws of the
#' forecast error variance decomposition.
#'
#' @seealso \code{\link{compute_impulse_responses.PosteriorBSVARSIGN}}, \code{\link{estimate.BSVARSIGN}}, \code{\link{normalise_posterior}}, \code{\link{summary}}, \code{\link{plot}}
#'
#' @author Xiaolei Wang \email{[email protected]} and Tomasz Woźniak \email{[email protected]}
#'
#' @references
#' Kilian, L., & Lütkepohl, H. (2017). Structural VAR Tools, Chapter 4, In: Structural vector autoregressive analysis. Cambridge University Press.
#'
#' @examples
#' # upload data
#' data(oil)
#'
#' # specify the model and set seed
#' set.seed(123)
#' sign_irf = array(matrix(c(-1, -1, 1, rep(NA, 6)), nrow = 3), dim = c(3, 3, 1))
#' specification = specify_bsvarSIGN$new(oil, sign_irf = sign_irf)
#'
#' # run the burn-in
#' posterior = estimate(specification, 10)
#'
#' # compute forecast error variance decomposition 2 years ahead
#' fevd = compute_variance_decompositions(posterior, horizon = 8)
#'
#' # workflow with the pipe |>
#' ############################################################
#' set.seed(123)
#' oil |>
#' specify_bsvarSIGN$new(sign_irf = sign_irf) |>
#' estimate(S = 10) |>
#' compute_variance_decompositions(horizon = 8) -> fevd
#'
#' @export
compute_variance_decompositions.PosteriorBSVARSIGN <- function(posterior, horizon) {

posterior_Theta0 = posterior$posterior$Theta0
posterior_A = posterior$posterior$A
posterior_A = aperm(posterior_A, c(2, 1, 3))
N = dim(posterior_A)[2]
p = posterior$last_draw$p
S = dim(posterior_A)[3]

posterior_irf = .Call(`_bsvarSIGNs_bsvarSIGNs_ir`, posterior_A, posterior_Theta0, horizon, p, TRUE)
qqq = .Call(`_bsvarSIGNs_bsvarSIGNs_fevd`, posterior_irf)

fevd = array(NA, c(N, N, horizon + 1, S))
for (s in 1:S) fevd[,,,s] = qqq[s][[1]]
class(fevd) = "PosteriorFEVD"

return(fevd)
} # END compute_variance_decompositions.PosteriorBSVARSIGN




#' @method compute_conditional_sd PosteriorBSVARSIGN
#'
#' @title Computes posterior draws of structural shock conditional standard deviations
#'
#' @description Each of the draws from the posterior estimation of models is
#' transformed into a draw from the posterior distribution of the structural
#' shock conditional standard deviations.
#'
#' @param posterior posterior estimation outcome - an object of class
#' \code{PosteriorBSVARSIGN} obtained by running the \code{estimate} function.
#'
#' @return An object of class \code{PosteriorSigma}, that is, an \code{NxTxS}
#' array with attribute \code{PosteriorSigma} containing \code{S} draws of the
#' structural shock conditional standard deviations.
#'
#' @seealso \code{\link{estimate.BSVARSIGN}}
#'
#' @author Xiaolei Wang \email{[email protected]} and Tomasz Woźniak \email{[email protected]}
#'
#' @examples
#' # upload data
#' data(oil)
#'
#' # specify the model and set seed
#' set.seed(123)
#' sign_irf = array(matrix(c(-1, -1, 1, rep(NA, 6)), nrow = 3), dim = c(3, 3, 1))
#' specification = specify_bsvarSIGN$new(oil, sign_irf = sign_irf)
#'
#' # run the burn-in
#' posterior = estimate(specification, 10)
#'
#' # compute structural shocks' conditional standard deviations
#' sigma = compute_conditional_sd(posterior)
#'
#' # workflow with the pipe |>
#' ############################################################
#' set.seed(123)
#' oil |>
#' specify_bsvarSIGN$new(sign_irf = sign_irf) |>
#' estimate(S = 10) |>
#' compute_conditional_sd() -> csd
#'
#' @export
compute_conditional_sd.PosteriorBSVARSIGN <- function(posterior) {

Y = posterior$last_draw$data_matrices$Y
N = nrow(Y)
T = ncol(Y)
S = dim(posterior$posterior$A)[3]

posterior_sigma = array(1, c(N, T, S))
message("The model is homoskedastic. Returning an NxTxS matrix of conditional sd all equal to 1.")
class(posterior_sigma) = "PosteriorSigma"

return(posterior_sigma)
} # END compute_conditional_sd.PosteriorBSVARSIGN
4 changes: 2 additions & 2 deletions R/forecast.R
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@
#' \item{Y}{an \eqn{NxT} matrix with the data on dependent variables}
#' }
#'
#' @seealso \code{\link{estimate}}, \code{\link{summary}}, \code{\link{plot}}
#' @seealso \code{\link{estimate.BSVARSIGN}}, \code{\link{summary}}, \code{\link{plot}}
#'
#' @author Tomasz Woźniak \email{[email protected]} and Xiaolei Wang \email{[email protected]}
#'
Expand All @@ -37,7 +37,7 @@
#'
#' # specify the model and set seed
#' set.seed(123)
#' sign_irf = array(matrix(c(-1, -1, 1, rep(0, 6)), nrow = 3), dim = c(3, 3, 1))
#' sign_irf = array(matrix(c(-1, -1, 1, rep(NA, 6)), nrow = 3), dim = c(3, 3, 1))
#' specification = specify_bsvarSIGN$new(oil, sign_irf = sign_irf)
#'
#' # estimate the model
Expand Down
21 changes: 21 additions & 0 deletions inst/include/bsvarSIGNs_RcppExports.h
Original file line number Diff line number Diff line change
Expand Up @@ -172,6 +172,27 @@ namespace bsvarSIGNs {
return Rcpp::as<arma::mat >(rcpp_result_gen);
}

inline arma::field<arma::cube> bsvarSIGNs_fevd(arma::field<arma::cube>& posterior_irf) {
typedef SEXP(*Ptr_bsvarSIGNs_fevd)(SEXP);
static Ptr_bsvarSIGNs_fevd p_bsvarSIGNs_fevd = NULL;
if (p_bsvarSIGNs_fevd == NULL) {
validateSignature("arma::field<arma::cube>(*bsvarSIGNs_fevd)(arma::field<arma::cube>&)");
p_bsvarSIGNs_fevd = (Ptr_bsvarSIGNs_fevd)R_GetCCallable("bsvarSIGNs", "_bsvarSIGNs_bsvarSIGNs_fevd");
}
RObject rcpp_result_gen;
{
RNGScope RCPP_rngScope_gen;
rcpp_result_gen = p_bsvarSIGNs_fevd(Shield<SEXP>(Rcpp::wrap(posterior_irf)));
}
if (rcpp_result_gen.inherits("interrupted-error"))
throw Rcpp::internal::InterruptedException();
if (Rcpp::internal::isLongjumpSentinel(rcpp_result_gen))
throw Rcpp::LongjumpException(rcpp_result_gen);
if (rcpp_result_gen.inherits("try-error"))
throw Rcpp::exception(Rcpp::as<std::string>(rcpp_result_gen).c_str());
return Rcpp::as<arma::field<arma::cube> >(rcpp_result_gen);
}

inline arma::cube forecast_bsvarSIGNs(arma::cube& posterior_Sigma, arma::cube& posterior_A, arma::vec& X_T, arma::mat& exogenous_forecast, arma::mat& cond_forecast, const int& horizon) {
typedef SEXP(*Ptr_forecast_bsvarSIGNs)(SEXP,SEXP,SEXP,SEXP,SEXP,SEXP);
static Ptr_forecast_bsvarSIGNs p_forecast_bsvarSIGNs = NULL;
Expand Down
34 changes: 34 additions & 0 deletions inst/tinytest/test_compute_variance_decomposition.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,34 @@

data(oil)

set.seed(1)
suppressMessages(
specification_no1 <- specify_bsvarSIGN$new(oil)
)
run_no1 <- estimate(specification_no1, 3, 1, show_progress = FALSE)
fevd <- compute_variance_decompositions(run_no1, horizon = 2)

set.seed(1)
suppressMessages(
fevd2 <- oil |>
specify_bsvarSIGN$new() |>
estimate(S = 3, thin = 1, show_progress = FALSE) |>
compute_variance_decompositions(horizon = 2)
)

expect_error(
compute_variance_decompositions(run_no1),
info = "compute_variance_decompositions: specify horizon."
)

expect_equal(
sum(fevd[1,,1,1]), 100,
info = "compute_variance_decompositions: sum to 100%."
)

expect_identical(
fevd[3,3,3,3], fevd2[3,3,3,3],
info = "compute_variance_decompositions: identical for normal and pipe workflow."
)


52 changes: 52 additions & 0 deletions man/compute_conditional_sd.PosteriorBSVARSIGN.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

4 changes: 2 additions & 2 deletions man/compute_fitted_values.PosteriorBSVARSIGN.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

4 changes: 2 additions & 2 deletions man/compute_historical_decompositions.PosteriorBSVARSIGN.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

Loading

0 comments on commit db8e18b

Please sign in to comment.