Skip to content

Commit

Permalink
included summary.PosteriorBVARPANEL #17
Browse files Browse the repository at this point in the history
+ required importing from stats
+ changed the forecast output name "Global" to "global" for consistency
  • Loading branch information
donotdespair committed Jul 25, 2024
1 parent 7980a0e commit dadfa88
Show file tree
Hide file tree
Showing 5 changed files with 253 additions and 1 deletion.
3 changes: 3 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@ S3method(estimate,PosteriorBVARPANEL)
S3method(forecast,PosteriorBVARPANEL)
S3method(plot,ForecastsPANEL)
S3method(plot,PosteriorFEVDPANEL)
S3method(summary,PosteriorBVARPANEL)
export(specify_bvarPANEL)
export(specify_panel_data_matrices)
export(specify_posterior_bvarPANEL)
Expand All @@ -19,4 +20,6 @@ importFrom(RcppTN,rtn)
importFrom(bsvars,compute_variance_decompositions)
importFrom(bsvars,estimate)
importFrom(bsvars,forecast)
importFrom(stats,quantile)
importFrom(stats,sd)
useDynLib(bvarPANELs, .registration = TRUE)
1 change: 1 addition & 0 deletions R/bvarPANELs-package.R
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,7 @@
#' @importFrom Rcpp sourceCpp
#' @importFrom R6 R6Class
#' @importFrom RcppTN rtn dtn
#' @importFrom stats sd quantile
#' @import RcppProgress
#' @note This package is currently in active development.
#' @author Tomasz Woźniak \email{wozniak.tom@pm.me}
Expand Down
2 changes: 1 addition & 1 deletion R/compute_variance_decompositions.PosteriorBVARPANEL.R
Original file line number Diff line number Diff line change
Expand Up @@ -85,7 +85,7 @@ compute_variance_decompositions.PosteriorBVARPANEL <- function(posterior, horizo
fevd[[c]] = fevd_c
}

names(fevd) = c(c_names, "Global")
names(fevd) = c(c_names, "global")

class(fevd) <- "PosteriorFEVDPANEL"
return(fevd)
Expand Down
194 changes: 194 additions & 0 deletions R/summary.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,194 @@


#' @title Provides posterior estimation summary for Bayesian Hierarchical Panel
#' Vector Autoregressions
#'
#' @description Provides posterior mean, standard deviations, as well as 5 and 95
#' percentiles of the parameters for all \code{C} countries.
#'
#' @param object an object of class \code{PosteriorBVARPANEL} obtained using the
#' \code{estimate()} function applied to
#' Vector Autoregressions containing draws from the posterior distribution of
#' the parameters.
#' @param ... additional arguments affecting the summary produced.
#'
#' @return A list reporting the posterior mean, standard deviations, as well as 5 and 95
#' percentiles of the country-specific and global parameters.
#'
#' @method summary PosteriorBVARPANEL
#'
#' @seealso \code{\link{estimate.BVARPANEL}}, \code{\link{specify_bvarPANEL}}
#'
#' @author Tomasz Woźniak \email{wozniak.tom@pm.me}
#'
#' @examples
#' # upload data
#' data(ilo_cubic_panel) # load the data
#' data(ilo_exogenous_variables) # load the exogenous variables
#'
#' set.seed(123)
#'
#' # specify the model
#' specification = specify_bvarPANEL$new(ilo_cubic_panel, exogenous = ilo_exogenous_variables)
#' burn_in = estimate(specification, 10) # run the burn-in
#' posterior = estimate(burn_in, 10) # estimate the model
#' summary(posterior)
#'
#' # workflow with the pipe |>
#' ############################################################
#' set.seed(123)
#' ilo_cubic_panel |>
#' specify_bvarPANEL$new(exogenous = ilo_exogenous_variables) |>
#' estimate(S = 10) |>
#' estimate(S = 10) |>
#' summary()
#'
#' @export
summary.PosteriorBVARPANEL = function(
object,
...
) {

S = dim(object$posterior$A_c)[4]
C = dim(object$posterior$A_c)[3]
N = dim(object$posterior$A_c)[2]
K = dim(object$posterior$A_c)[1]
p = object$last_draw$p
d = K - N * p

out = list()
param = c("A", "Sigma")
country_names = names(object$last_draw$data_matrices$Y)

# country-specific parameter summary
for (c in 1:C) {

out_c = list()
out_c$Sigma = list()
out_c$A = list()

for (n in 1:N) {

Sigma_c = matrix(object$posterior$Sigma_c[n,1:n,c,], nrow = n)
out_c$Sigma[[n]] = matrix(
cbind(
apply(Sigma_c, 1, mean),
apply(Sigma_c, 1, sd),
apply(Sigma_c, 1, quantile, probs = 0.05),
apply(Sigma_c, 1, quantile, probs = 0.95)
),
ncol = 4
)
colnames(out_c$Sigma[[n]]) = c("mean", "sd", "5% quantile", "95% quantile")
rownames(out_c$Sigma[[n]]) = paste0("Sigma[", n, ",", 1:n, "]")

A_c = object$posterior$A_c[,n,c,]
out_c$A[[n]] = cbind(
apply(A_c, 1, mean),
apply(A_c, 1, sd),
apply(A_c, 1, quantile, probs = 0.05),
apply(A_c, 1, quantile, probs = 0.95)
)
colnames(out_c$A[[n]]) = c("mean", "sd", "5% quantile", "95% quantile")

Anames = c(
paste0(
rep("lag", p * N),
kronecker((1:p), rep(1, N)),
rep("_var", p * N),
kronecker((1:N), rep(1, p))
),
"const"
)
if (d > 1) {
Anames = c(Anames, paste0("exo", 1:(d - 1)))
}
rownames(out_c$A[[n]]) = Anames
} # END n loop

names(out_c$Sigma) = paste0("equation", 1:N)
names(out_c$A) = paste0("equation", 1:N)

out[[c]] = out_c
} # END c loop

names(out) = country_names


# global parameter summary
out_g = list()
out_g$A = list()
out_g$Sigma = list()
out_g$V = list()
out_g$hyper = list()

for (n in 1:N) {

Sigma = matrix(object$posterior$Sigma[n,1:n,], nrow = n)
out_g$Sigma[[n]] = matrix(
cbind(
apply(Sigma, 1, mean),
apply(Sigma, 1, sd),
apply(Sigma, 1, quantile, probs = 0.05),
apply(Sigma, 1, quantile, probs = 0.95)
),
ncol = 4
)
colnames(out_g$Sigma[[n]]) = c("mean", "sd", "5% quantile", "95% quantile")
rownames(out_g$Sigma[[n]]) = paste0("Sigma[", n, ",", 1:n, "]")

A = object$posterior$A[,n,]
out_g$A[[n]] = cbind(
apply(A, 1, mean),
apply(A, 1, sd),
apply(A, 1, quantile, probs = 0.05),
apply(A, 1, quantile, probs = 0.95)
)
colnames(out_g$A[[n]]) = c("mean", "sd", "5% quantile", "95% quantile")
rownames(out_g$A[[n]]) = Anames

} # END n loop

names(out_g$Sigma) = paste0("equation", 1:N)
names(out_g$A) = paste0("equation", 1:N)


for (k in 1:K) {

V = matrix(object$posterior$V[k,1:k,], nrow = k)
out_g$V[[k]] = matrix(
cbind(
apply(V, 1, mean),
apply(V, 1, sd),
apply(V, 1, quantile, probs = 0.05),
apply(V, 1, quantile, probs = 0.95)
),
ncol = 4
)
colnames(out_g$V[[k]]) = c("mean", "sd", "5% quantile", "95% quantile")
rownames(out_g$V[[k]]) = paste0("V[", k, ",", 1:k, "]")

} # END k loop

hyper = t(cbind(
object$posterior$nu,
object$posterior$m,
object$posterior$w,
object$posterior$s
))
out_g$hyper = matrix(
cbind(
apply(hyper, 1, mean),
apply(hyper, 1, sd),
apply(hyper, 1, quantile, probs = 0.05),
apply(hyper, 1, quantile, probs = 0.95)
),
ncol = 4
)
colnames(out_g$hyper) = c("mean", "sd", "5% quantile", "95% quantile")
rownames(out_g$hyper) = c("nu", "m", "w", "s")

out$global = out_g

return(out)
} # END summary.PosteriorBVARPANEL
54 changes: 54 additions & 0 deletions man/summary.PosteriorBVARPANEL.Rd

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

0 comments on commit dadfa88

Please sign in to comment.