Skip to content

Commit

Permalink
Merge branch 'feature/spillover-mcmc' into develop
Browse files Browse the repository at this point in the history
  • Loading branch information
ygeunkim committed Nov 16, 2024
2 parents 6add121 + 6fb9e1f commit 4e7fc5b
Show file tree
Hide file tree
Showing 13 changed files with 576 additions and 214 deletions.
2 changes: 2 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -286,6 +286,7 @@ importFrom(dplyr,left_join)
importFrom(dplyr,mutate)
importFrom(dplyr,n)
importFrom(dplyr,rename)
importFrom(dplyr,select)
importFrom(foreach,"%do%")
importFrom(foreach,"%dopar%")
importFrom(foreach,foreach)
Expand Down Expand Up @@ -343,6 +344,7 @@ importFrom(tibble,add_column)
importFrom(tibble,as_tibble)
importFrom(tibble,rownames_to_column)
importFrom(tidyr,pivot_longer)
importFrom(tidyr,pivot_wider)
importFrom(tidyr,separate)
importFrom(tidyr,separate_wider_regex)
importFrom(tidyr,unite)
Expand Down
111 changes: 111 additions & 0 deletions R/misc-r.R
Original file line number Diff line number Diff line change
Expand Up @@ -171,6 +171,117 @@ get_gammaparam <- function(mode, sd) {
)
}

#' Compute Summaries from Forecast Draws
#'
#' @param draws Matrix in forms of rbind(step) x cbind(draws)
#' @param n_ahead Forecast step used
#' @param dim_data Dimension
#' @param num_draw MCMC draws
#' @param var_names Variable names
#' @param level level for lower and upper quantiles
#' @param med Get median instead of mean?
#'
#' @noRd
process_forecast_draws <- function(draws, n_ahead, dim_data, num_draw, var_names, level = .05, med = FALSE) {
mcmc_distn <-
draws %>%
unlist() %>%
array(dim = c(n_ahead, dim_data, num_draw))
if (med) {
pred_mean <- apply(mcmc_distn, c(1, 2), median)
} else {
pred_mean <- apply(mcmc_distn, c(1, 2), mean)
}
pred_se <- apply(mcmc_distn, c(1, 2), sd)
pred_lower <- apply(mcmc_distn, c(1, 2), quantile, probs = level / 2)
pred_upper <- apply(mcmc_distn, c(1, 2), quantile, probs = 1 - level / 2)
colnames(pred_mean) <- var_names
rownames(pred_mean) <- var_names
colnames(pred_se) <- var_names
rownames(pred_se) <- var_names
colnames(pred_lower) <- var_names
rownames(pred_lower) <- var_names
colnames(pred_upper) <- var_names
rownames(pred_upper) <- var_names
list(
mean = pred_mean,
sd = pred_se,
lower = pred_lower,
upper = pred_upper
)
}

#' Compute Summaries from Vector Draws
#'
#' @param dim_data Dimension
#' @param level level for lower and upper quantiles
#' @param med Get median instead of mean?
#'
#' @noRd
process_vector_draws <- function(draws, dim_data, level = .05, med = FALSE) {
mcmc_distn <- matrix(draws, ncol = dim_data)
if (med) {
pred_mean <- apply(mcmc_distn, 2, median)
} else {
pred_mean <- colMeans(mcmc_distn)
}
list(
mean = pred_mean,
sd = apply(mcmc_distn, 2, sd),
lower = apply(mcmc_distn, 2, quantile, probs = level / 2),
upper = apply(mcmc_distn, 2, quantile, probs = 1 - level / 2)
)
}

#' Compute Summaries from Dynamic Spillover
#'
#' @param dim_data Dimension
#' @param level level for lower and upper quantiles
#' @param med Get median instead of mean?
#' @param var_names Variable names
#'
#' @importFrom tibble as_tibble
#' @importFrom dplyr mutate
#' @noRd
process_dynamic_spdraws <- function(draws, dim_data, level = .05, med = FALSE, var_names) {
lapply(
draws,
function(x) {
process_vector_draws(unlist(x), dim_data = dim_data, level = level, med = med) %>%
do.call(cbind, .) %>%
as.data.frame() %>%
mutate(series = var_names)
}
) %>%
do.call(rbind, .) %>%
as_tibble()
}

#' Pivot longer spillover
#'
#' @param connect Connectedness table
#' @param col_names Column name for value
#' @noRd
gather_spillover <- function(connect, col_names = "spillover") {
connect %>%
as.data.frame() %>%
rownames_to_column(var = "series") %>%
pivot_longer(-"series", names_to = "shock", values_to = col_names)
}

#' Pivot longer spillover summaries
#'
#' @param distn Connectedness table distribution
#' @param prefix Column names prefix
#'
#' @noRd
join_long_spillover <- function(connect, prefix = "spillover") {
gather_spillover(connect$mean, col_names = prefix) %>%
left_join(gather_spillover(connect$lower, col_names = paste(prefix, "lower", sep = "_")), by = c("series", "shock")) %>%
left_join(gather_spillover(connect$upper, col_names = paste(prefix, "upper", sep = "_")), by = c("series", "shock")) %>%
left_join(gather_spillover(connect$sd, col_names = paste(prefix, "sd", sep = "_")), by = c("series", "shock"))
}

#' Define Minnesota Group Matrix
#'
#' This function creates a matrix with group index
Expand Down
64 changes: 57 additions & 7 deletions R/print-spillover.R
Original file line number Diff line number Diff line change
Expand Up @@ -8,24 +8,47 @@ print.bvharspillover <- function(x, digits = max(3L, getOption("digits") - 3L),
cat("Directional spillovers:\n")
cat("variables (i) <- shocks (j)\n")
cat("========================\n")
if (is.list(x$connect)) {
connect <- x$connect$mean
to <- x$to$mean
from <- x$from$mean
tot <- x$tot$mean
} else {
connect <- x$connect
to <- x$to
from <- x$from
tot <- x$tot
}
connect <- rbind(connect, "to_spillovers" = to)
connect <- cbind(connect, "from_spillovers" = c(from, tot))
print(
x$connect,
connect,
digits = digits,
print.gap = 2L,
quote = FALSE
)
cat("\n*Lower right corner: Total spillover\n")
cat("------------------------\n")
cat("Net spillovers:\n")
if (is.list(x$net)) {
net <- x$net$mean
} else {
net <- x$net
}
print(
x$net,
net,
digits = digits,
print.gap = 2L,
quote = FALSE
)
cat("\nNet pairwise spillovers:\n")
if (is.list(x$net_pairwise)) {
net_pairwise <- x$net_pairwise$mean
} else {
net_pairwise <- x$net_pairwise
}
print(
x$net_pairwise,
net_pairwise,
digits = digits,
print.gap = 2L,
quote = FALSE
Expand All @@ -43,15 +66,20 @@ knit_print.bvharspillover <- function(x, ...) {
#' @param digits digit option to print
#' @param ... not used
#' @importFrom utils head
#' @importFrom dplyr mutate n select
#' @importFrom tidyr pivot_wider
#' @order 2
#' @export
print.bvhardynsp <- function(x, digits = max(3L, getOption("digits") - 3L), ...) {
cat("Dynamics of spillover effect:\n")
cat(sprintf("Forecast using %s\n", x$process))
cat(sprintf("Forecast step: %d\n", x$ahead))
cat("========================\n")
is_mcmc <- !is.vector(x$tot)
cat("Total spillovers:\n")
cat(sprintf("# A vector: %d\n", length(x$tot)))
if (!is_mcmc) {
cat(sprintf("# A vector: %d\n", length(x$tot)))
}
print(
head(x$tot),
digits = digits,
Expand All @@ -66,25 +94,47 @@ print.bvhardynsp <- function(x, digits = max(3L, getOption("digits") - 3L), ...)
# print.gap = 2L,
# quote = FALSE
# )
if (is_mcmc) {
dim_data <- nrow(x$to) / length(x$index)
to_distn <-
x$to %>%
select("series", "mean") %>%
mutate(id = rep(x$index, each = dim_data)) %>%
pivot_wider(names_from = "series", values_from = "mean")
from_distn <-
x$from %>%
select("series", "mean") %>%
mutate(id = rep(x$index, each = dim_data)) %>%
pivot_wider(names_from = "series", values_from = "mean")
net_distn <-
x$net %>%
select("series", "mean") %>%
mutate(id = rep(x$index, each = dim_data)) %>%
pivot_wider(names_from = "series", values_from = "mean")
} else {
to_distn <- x$to
from_distn <- x$from
net_distn <- x$net
}
cat("To spillovers:\n")
print(
x$to,
to_distn,
digits = digits,
print.gap = 2L,
quote = FALSE
)
cat("------------------------\n")
cat("From spillovers:\n")
print(
x$from,
from_distn,
digits = digits,
print.gap = 2L,
quote = FALSE
)
cat("------------------------\n")
cat("Net spillovers:\n")
print(
x$net,
net_distn,
digits = digits,
print.gap = 2L,
quote = FALSE
Expand Down
Loading

0 comments on commit 4e7fc5b

Please sign in to comment.