Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

27 coherent bsvars notation #29

Merged
merged 4 commits into from
Jul 15, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
13 changes: 10 additions & 3 deletions R/estimate.BSVARSIGN.R
Original file line number Diff line number Diff line change
Expand Up @@ -96,16 +96,23 @@ estimate.BSVARSIGN = function(specification, S, thin = 1, show_progress = TRUE)

# get the inputs to estimation
# prior = specification$last_draw$prior$get_prior()
prior = specification$prior
prior = specification$prior$get_prior()
starting_values = specification$starting_values$get_starting_values()
identification = specification$identification$get_identification()
max_tries = identification$max_tries
data_matrices = specification$data_matrices$get_data_matrices()
p = specification$p

prior$B = t(prior$A)
prior$Ysoc = t(prior$Ysoc)
prior$Xsoc = t(prior$Xsoc)
prior$Ysur = t(prior$Ysur)
prior$Xsur = t(prior$Xsur)
Y = t(data_matrices$Y)
X = t(data_matrices$X)

# estimation
qqq = .Call(`_bsvarSIGNs_bsvar_sign_cpp`, S, p,
data_matrices$Y, data_matrices$X, identification$VB,
qqq = .Call(`_bsvarSIGNs_bsvar_sign_cpp`, S, p, Y, X, identification$VB,
identification$sign_irf, identification$sign_narrative,
identification$sign_B, identification$zero_irf,
prior, starting_values, show_progress, thin, max_tries)
Expand Down
66 changes: 32 additions & 34 deletions R/specify_bsvarSIGN.R
Original file line number Diff line number Diff line change
Expand Up @@ -139,23 +139,20 @@ specify_prior_bsvarSIGN = R6::R6Class(
"PriorBSVARSIGN",

public = list(
#' @field p a positive integer - the number of lags.
p = 1,

#' @field hyper a \code{(N+3)xS} matrix of hyper-parameters \eqn{\mu, \delta, \lambda, \psi}.
hyper = matrix(),

#' @field B a \code{KxN} normal prior mean matrix for the autoregressive
#' @field A a \code{NxK} normal prior mean matrix for the autoregressive
#' parameters.
B = matrix(),
A = matrix(),

#' @field V a \code{KxK} matrix determining the normal prior column-specific
#' covariance for the autoregressive parameters.
V = matrix(),

#' @field Vp a \code{px1} vector of lag shrinkage level.
Vp = matrix(),

#' @field Vd a \code{(d+1)x1} vector of (large) variances for constants.
Vd = matrix(),

#' @field S an \code{NxN} matrix determining the inverted-Wishart prior scale
#' of error terms covariance matrix.
S = matrix(),
Expand All @@ -167,22 +164,22 @@ specify_prior_bsvarSIGN = R6::R6Class(
#' @field data an \code{TxN} matrix of observations.
data = matrix(),

#' @field Y an \code{TxN} matrix of dependent variables.
#' @field Y an \code{NxT} matrix of dependent variables.
Y = matrix(),

#' @field X an \code{TxK} matrix of independent variables.
#' @field X an \code{KxT} matrix of independent variables.
X = matrix(),

#' @field Ysoc an \code{NxN} matrix with the sum-of-coefficients dummy observations.
Ysoc = matrix(),

#' @field Xsoc an \code{NxK} matrix with the sum-of-coefficients dummy observations.
#' @field Xsoc an \code{KxN} matrix with the sum-of-coefficients dummy observations.
Xsoc = matrix(),

#' @field Ysur an \code{NxN} matrix with the single-unit-root dummy observations.
Ysur = matrix(),

#' @field Xsur an \code{NxK} matrix with the single-unit-root dummy observations.
#' @field Xsur an \code{KxN} matrix with the single-unit-root dummy observations.
Xsur = matrix(),

#' @field mu.scale a positive scalar - the shape of the gamma prior for \eqn{\mu}.
Expand Down Expand Up @@ -224,7 +221,7 @@ specify_prior_bsvarSIGN = R6::R6Class(
#' prior = specify_prior_bsvarSIGN$new(oil, p = 1)
#' prior$B # show autoregressive prior mean
#'
initialize = function(data, p, exogenous = NULL, stationary = rep(FALSE, dim(data)[2])) {
initialize = function(data, p, exogenous = NULL, stationary = rep(FALSE, ncol(data))) {

stopifnot("Argument p must be a positive integer number." = p > 0 & p %% 1 == 0)

Expand All @@ -245,9 +242,7 @@ specify_prior_bsvarSIGN = R6::R6Class(
B = matrix(0, K, N)
B[1:N,] = diag(!stationary)

Vp = (1:p)^-2
Vd = rep(100, 1 + d)
V = diag(c(kronecker(Vp, rep(1, N)), Vd))
V = diag(c(kronecker((1:p)^-2, rep(1, N)), rep(100, 1 + d)))

s2.ols = rep(NA, N)
for (n in 1:N) {
Expand Down Expand Up @@ -282,20 +277,18 @@ specify_prior_bsvarSIGN = R6::R6Class(
# }
# Xstar = cbind(Xstar, c(rep(0, N), 1), matrix(0, N + 1, d))

self$Y = Y
self$X = X
self$Vp = Vp
self$Vd = Vd

self$p = p
self$hyper = hyper
self$B = B
self$A = t(B)
self$V = V
self$S = diag(N)
self$nu = N + 2
self$Ysoc = Ysoc
self$Xsoc = Xsoc
self$Ysur = Ysur
self$Xsur = Xsur
self$Y = t(Y)
self$X = t(X)
self$Ysoc = t(Ysoc)
self$Xsoc = t(Xsoc)
self$Ysur = t(Ysur)
self$Xsur = t(Xsur)
self$mu.scale = scale
self$mu.shape = shape
self$delta.scale = scale
Expand All @@ -316,11 +309,10 @@ specify_prior_bsvarSIGN = R6::R6Class(
#'
get_prior = function(){
list(
p = self$p,
hyper = self$hyper,
B = self$B,
A = self$A,
V = self$V,
Vp = self$Vp,
Vd = self$Vd,
S = self$S,
nu = self$nu,
Ysoc = self$Ysoc,
Expand Down Expand Up @@ -373,10 +365,16 @@ specify_prior_bsvarSIGN = R6::R6Class(
init = narrow_hyper(model, hyper)
prior = self$get_prior()

prior$B = t(prior$A)
prior$Ysoc = t(prior$Ysoc)
prior$Xsoc = t(prior$Xsoc)
prior$Ysur = t(prior$Ysur)
prior$Xsur = t(prior$Xsur)

result = stats::optim(
init,
\(x) -log_posterior_hyper(extend_hyper(hyper, model, matrix(x)),
model, self$Y, self$X, prior),
model, t(self$Y), t(self$X), prior),
method = 'L-BFGS-B',
lower = rep(0, length(init)),
upper = init * 100,
Expand All @@ -393,7 +391,8 @@ specify_prior_bsvarSIGN = R6::R6Class(
variance = e$vectors %*% diag(as.vector(1 / abs(e$values))) %*% t(e$vectors)
}

self$hyper = sample_hyper(S, burn_in, mode, model, self$Y, self$X, variance, prior)
self$hyper = sample_hyper(S, burn_in, mode, model,
t(self$Y), t(self$X), variance, prior)
self$hyper = self$hyper[, -(1:burn_in)]
} # END estimate_hyper

Expand Down Expand Up @@ -716,15 +715,14 @@ specify_bsvarSIGN = R6::R6Class(
B[lower.tri(B, diag = TRUE)] = TRUE

self$data_matrices = bsvars::specify_data_matrices$new(data, p, exogenous)
self$data_matrices$Y = t(self$data_matrices$Y)
self$data_matrices$X = t(self$data_matrices$X)
self$identification = specify_identification_bsvarSIGN$new(N,
sign_irf,
sign_narrative,
sign_B,
zero_irf,
max_tries)
self$prior = specify_prior_bsvarSIGN$new(data, p, exogenous, stationary)
self$prior = specify_prior_bsvarSIGN$new(data, p, exogenous,
stationary)
self$starting_values = bsvars::specify_starting_values_bsvar$new(N, self$p, d)
}, # END initialize

Expand Down
Binary file modified data/oil.rda
Binary file not shown.
4 changes: 2 additions & 2 deletions inst/include/bsvarSIGNs_RcppExports.h
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,7 @@ namespace bsvarSIGNs {
return Rcpp::as<arma::field<arma::cube> >(rcpp_result_gen);
}

inline Rcpp::List bsvar_sign_cpp(const int& S, const int& lags, const arma::mat& Y, const arma::mat& X, const arma::field<arma::mat>& VB, const arma::cube& sign_irf, const arma::mat& sign_narrative, const arma::mat& sign_B, const arma::field<arma::mat>& Z, const Rcpp::List& prior, const Rcpp::List& starting_values, const bool show_progress = true, const int thin = 100, const int& max_tries = 10000) {
inline Rcpp::List bsvar_sign_cpp(const int& S, const int& p, const arma::mat& Y, const arma::mat& X, const arma::field<arma::mat>& VB, const arma::cube& sign_irf, const arma::mat& sign_narrative, const arma::mat& sign_B, const arma::field<arma::mat>& Z, const Rcpp::List& prior, const Rcpp::List& starting_values, const bool show_progress = true, const int thin = 100, const int& max_tries = 10000) {
typedef SEXP(*Ptr_bsvar_sign_cpp)(SEXP,SEXP,SEXP,SEXP,SEXP,SEXP,SEXP,SEXP,SEXP,SEXP,SEXP,SEXP,SEXP,SEXP);
static Ptr_bsvar_sign_cpp p_bsvar_sign_cpp = NULL;
if (p_bsvar_sign_cpp == NULL) {
Expand All @@ -98,7 +98,7 @@ namespace bsvarSIGNs {
RObject rcpp_result_gen;
{
RNGScope RCPP_rngScope_gen;
rcpp_result_gen = p_bsvar_sign_cpp(Shield<SEXP>(Rcpp::wrap(S)), Shield<SEXP>(Rcpp::wrap(lags)), Shield<SEXP>(Rcpp::wrap(Y)), Shield<SEXP>(Rcpp::wrap(X)), Shield<SEXP>(Rcpp::wrap(VB)), Shield<SEXP>(Rcpp::wrap(sign_irf)), Shield<SEXP>(Rcpp::wrap(sign_narrative)), Shield<SEXP>(Rcpp::wrap(sign_B)), Shield<SEXP>(Rcpp::wrap(Z)), Shield<SEXP>(Rcpp::wrap(prior)), Shield<SEXP>(Rcpp::wrap(starting_values)), Shield<SEXP>(Rcpp::wrap(show_progress)), Shield<SEXP>(Rcpp::wrap(thin)), Shield<SEXP>(Rcpp::wrap(max_tries)));
rcpp_result_gen = p_bsvar_sign_cpp(Shield<SEXP>(Rcpp::wrap(S)), Shield<SEXP>(Rcpp::wrap(p)), Shield<SEXP>(Rcpp::wrap(Y)), Shield<SEXP>(Rcpp::wrap(X)), Shield<SEXP>(Rcpp::wrap(VB)), Shield<SEXP>(Rcpp::wrap(sign_irf)), Shield<SEXP>(Rcpp::wrap(sign_narrative)), Shield<SEXP>(Rcpp::wrap(sign_B)), Shield<SEXP>(Rcpp::wrap(Z)), Shield<SEXP>(Rcpp::wrap(prior)), Shield<SEXP>(Rcpp::wrap(starting_values)), Shield<SEXP>(Rcpp::wrap(show_progress)), Shield<SEXP>(Rcpp::wrap(thin)), Shield<SEXP>(Rcpp::wrap(max_tries)));
}
if (rcpp_result_gen.inherits("interrupted-error"))
throw Rcpp::internal::InterruptedException();
Expand Down
12 changes: 6 additions & 6 deletions inst/tinytest/test_specify.R
Original file line number Diff line number Diff line change
Expand Up @@ -8,13 +8,13 @@ spec = specify_bsvarSIGN$new(oil)
expect_identical(class(spec)[1],
"BSVARSIGN")

expect_identical(dim(spec$data_matrices$Y)[1],
dim(spec$data_matrices$X)[1])
expect_identical(dim(spec$data_matrices$Y)[2],
dim(spec$data_matrices$X)[2])

expect_identical(length(spec$data_matrices$get_data_matrices()),
2L)

expect_identical(dim(spec$prior$B)[1], dim(spec$data_matrices$X)[2])
expect_identical(dim(spec$prior$A)[2], dim(spec$data_matrices$X)[1])

expect_identical(class(spec$prior$get_prior()),
"list")
Expand Down Expand Up @@ -51,13 +51,13 @@ spec = specify_bsvarSIGN$new(oil,
expect_identical(class(spec)[1],
"BSVARSIGN")

expect_identical(dim(spec$data_matrices$Y)[1],
dim(spec$data_matrices$X)[1])
expect_identical(dim(spec$data_matrices$Y)[2],
dim(spec$data_matrices$X)[2])

expect_identical(length(spec$data_matrices$get_data_matrices()),
2L)

expect_identical(dim(spec$prior$B)[1], dim(spec$data_matrices$X)[2])
expect_identical(dim(spec$prior$A)[2], dim(spec$data_matrices$X)[1])

expect_identical(class(spec$prior$get_prior()),
"list")
Expand Down
9 changes: 5 additions & 4 deletions inst/varia/nicetry.R
Original file line number Diff line number Diff line change
Expand Up @@ -4,17 +4,18 @@ name = sapply(1:7, \(i) data$ShortDescr[[i]][[1]])
data = data$y
colnames(data) = name

data = optimism
# data = optimism

spec = specify_bsvarSIGN$new(data, p = 4)

start_time = Sys.time()

spec$prior$estimate_hyper(S = 10000, burn_in = 5000,
mu = TRUE, delta = TRUE, lambda = TRUE, psi = TRUE)
# post = estimate(spec, S = 1000)
# irf = compute_impulse_responses(post, horizon = 40)
# plot(irf, probability = 0.68)
post = estimate(spec, S = 100)
class(post) = "PosteriorBSVAR"
irf = compute_impulse_responses(post, horizon = 40)
plot(irf, probability = 0.68)

end_time = Sys.time()
end_time - start_time
Expand Down
18 changes: 8 additions & 10 deletions man/specify_prior_bsvarSIGN.Rd

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

12 changes: 6 additions & 6 deletions src/RcppExports.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -127,12 +127,12 @@ RcppExport SEXP _bsvarSIGNs_bsvarSIGNs_ir(SEXP posterior_BSEXP, SEXP posterior_T
return rcpp_result_gen;
}
// bsvar_sign_cpp
Rcpp::List bsvar_sign_cpp(const int& S, const int& lags, const arma::mat& Y, const arma::mat& X, const arma::field<arma::mat>& VB, const arma::cube& sign_irf, const arma::mat& sign_narrative, const arma::mat& sign_B, const arma::field<arma::mat>& Z, const Rcpp::List& prior, const Rcpp::List& starting_values, const bool show_progress, const int thin, const int& max_tries);
static SEXP _bsvarSIGNs_bsvar_sign_cpp_try(SEXP SSEXP, SEXP lagsSEXP, SEXP YSEXP, SEXP XSEXP, SEXP VBSEXP, SEXP sign_irfSEXP, SEXP sign_narrativeSEXP, SEXP sign_BSEXP, SEXP ZSEXP, SEXP priorSEXP, SEXP starting_valuesSEXP, SEXP show_progressSEXP, SEXP thinSEXP, SEXP max_triesSEXP) {
Rcpp::List bsvar_sign_cpp(const int& S, const int& p, const arma::mat& Y, const arma::mat& X, const arma::field<arma::mat>& VB, const arma::cube& sign_irf, const arma::mat& sign_narrative, const arma::mat& sign_B, const arma::field<arma::mat>& Z, const Rcpp::List& prior, const Rcpp::List& starting_values, const bool show_progress, const int thin, const int& max_tries);
static SEXP _bsvarSIGNs_bsvar_sign_cpp_try(SEXP SSEXP, SEXP pSEXP, SEXP YSEXP, SEXP XSEXP, SEXP VBSEXP, SEXP sign_irfSEXP, SEXP sign_narrativeSEXP, SEXP sign_BSEXP, SEXP ZSEXP, SEXP priorSEXP, SEXP starting_valuesSEXP, SEXP show_progressSEXP, SEXP thinSEXP, SEXP max_triesSEXP) {
BEGIN_RCPP
Rcpp::RObject rcpp_result_gen;
Rcpp::traits::input_parameter< const int& >::type S(SSEXP);
Rcpp::traits::input_parameter< const int& >::type lags(lagsSEXP);
Rcpp::traits::input_parameter< const int& >::type p(pSEXP);
Rcpp::traits::input_parameter< const arma::mat& >::type Y(YSEXP);
Rcpp::traits::input_parameter< const arma::mat& >::type X(XSEXP);
Rcpp::traits::input_parameter< const arma::field<arma::mat>& >::type VB(VBSEXP);
Expand All @@ -145,15 +145,15 @@ BEGIN_RCPP
Rcpp::traits::input_parameter< const bool >::type show_progress(show_progressSEXP);
Rcpp::traits::input_parameter< const int >::type thin(thinSEXP);
Rcpp::traits::input_parameter< const int& >::type max_tries(max_triesSEXP);
rcpp_result_gen = Rcpp::wrap(bsvar_sign_cpp(S, lags, Y, X, VB, sign_irf, sign_narrative, sign_B, Z, prior, starting_values, show_progress, thin, max_tries));
rcpp_result_gen = Rcpp::wrap(bsvar_sign_cpp(S, p, Y, X, VB, sign_irf, sign_narrative, sign_B, Z, prior, starting_values, show_progress, thin, max_tries));
return rcpp_result_gen;
END_RCPP_RETURN_ERROR
}
RcppExport SEXP _bsvarSIGNs_bsvar_sign_cpp(SEXP SSEXP, SEXP lagsSEXP, SEXP YSEXP, SEXP XSEXP, SEXP VBSEXP, SEXP sign_irfSEXP, SEXP sign_narrativeSEXP, SEXP sign_BSEXP, SEXP ZSEXP, SEXP priorSEXP, SEXP starting_valuesSEXP, SEXP show_progressSEXP, SEXP thinSEXP, SEXP max_triesSEXP) {
RcppExport SEXP _bsvarSIGNs_bsvar_sign_cpp(SEXP SSEXP, SEXP pSEXP, SEXP YSEXP, SEXP XSEXP, SEXP VBSEXP, SEXP sign_irfSEXP, SEXP sign_narrativeSEXP, SEXP sign_BSEXP, SEXP ZSEXP, SEXP priorSEXP, SEXP starting_valuesSEXP, SEXP show_progressSEXP, SEXP thinSEXP, SEXP max_triesSEXP) {
SEXP rcpp_result_gen;
{
Rcpp::RNGScope rcpp_rngScope_gen;
rcpp_result_gen = PROTECT(_bsvarSIGNs_bsvar_sign_cpp_try(SSEXP, lagsSEXP, YSEXP, XSEXP, VBSEXP, sign_irfSEXP, sign_narrativeSEXP, sign_BSEXP, ZSEXP, priorSEXP, starting_valuesSEXP, show_progressSEXP, thinSEXP, max_triesSEXP));
rcpp_result_gen = PROTECT(_bsvarSIGNs_bsvar_sign_cpp_try(SSEXP, pSEXP, YSEXP, XSEXP, VBSEXP, sign_irfSEXP, sign_narrativeSEXP, sign_BSEXP, ZSEXP, priorSEXP, starting_valuesSEXP, show_progressSEXP, thinSEXP, max_triesSEXP));
}
Rboolean rcpp_isInterrupt_gen = Rf_inherits(rcpp_result_gen, "interrupted-error");
if (rcpp_isInterrupt_gen) {
Expand Down
Loading
Loading