Skip to content

Commit

Permalink
add warning messages, fixed things, add tests, wrote the doc about wa…
Browse files Browse the repository at this point in the history
…rnings
  • Loading branch information
nicorbtt committed Nov 30, 2023
1 parent 92c7427 commit dd23119
Show file tree
Hide file tree
Showing 3 changed files with 62 additions and 33 deletions.
51 changes: 31 additions & 20 deletions R/reconc.R
Original file line number Diff line number Diff line change
Expand Up @@ -53,23 +53,23 @@
warning = TRUE
warning_code = c(warning_code, 1)
warning_msg = c(warning_msg,
paste("<add warning message>", n_eff))
paste0("Importance Sampling: all the weights are zeros. This is probably caused by a strong incoherence between bottom and upper base forecasts. An uniform distribution for the weights is used. effective_sample_size= ", round(n_eff,2),"."))
}

# 2. n_eff < threshold
if (n_eff < n_eff_min) {
warning = TRUE
warning_code = c(warning_code, 2)
warning_msg = c(warning_msg,
paste("<add warning message>", n_eff))
paste0("Importance Sampling: effective_sample_size= ", round(n_eff,2), " (< ", n_eff_min,")."))
}

# 3. n_eff < p*n, e.g. p = 0.05
if (n_eff < p_n_eff*n) {
warning = TRUE
warning_code = c(warning_code, 3)
warning_msg = c(warning_msg,
paste("<add warning message>", n_eff))
paste0("Importance Sampling: effective_sample_size= ", round(n_eff,2), " (< ", round(p_n_eff * 100, 2),"%)."))
}

res = list(warning = warning,
Expand Down Expand Up @@ -127,32 +127,40 @@
#'
#' If `in_type[[i]]`='params', then `base_forecast[[i]]` is a vector containing the estimated:
#'
#' * mean and sd for the Gaussian base forecast if `distr[[i]]`='gaussian', see \link[stats]{Normal},;
#' * mean and sd for the Gaussian base forecast if `distr[[i]]`='gaussian', see \link[stats]{Normal};
#' * lambda for the Poisson base forecast if `distr[[i]]`='poisson', see \link[stats]{Poisson};
#' * mu and size for the negative binomial base forecast if `distr[[i]]`='nbinom', see \link[stats]{NegBinomial}.
#'
#' See the description of the parameters `in_type` and `distr` for more details.
#'
#' The order of the `base_forecast` list is given by the order of the time series in the summing matrix.
#'
#' Warnings are triggered from the Importance Sampling step if:
#'
#' * weights are all zeros;
#' * the effective sample size is < 200;
#' * the effective sample size is < 5% of the sample size (`num_samples` if `in_type` is 'params' or the size of the base forecast if if `in_type` is 'samples').
#'
#' @param S summing matrix (n x n_bottom).
#' @param base_forecasts a list containing the base_forecasts, see details.
#' @param in_type a string or a list of length n. If it is a list the i-th element is a string with two possible values:
#' @param S Summing matrix (n x n_bottom).
#' @param base_forecasts A list containing the base_forecasts, see details.
#' @param in_type A string or a list of length n. If it is a list the i-th element is a string with two possible values:
#'
#' * 'samples' if the i-th base forecasts are in the form of samples;
#' * 'params' if the i-th base forecasts are in the form of estimated parameters.
#'
#' If it `in_type` is a string it is assumed that all base forecasts are of the same type.
#'
#' @param distr a string or a list of length n describing the type of base forecasts. If it is a list the i-th element is a string with two possible values:
#' @param distr A string or a list of length n describing the type of base forecasts. If it is a list the i-th element is a string with two possible values:
#'
#' * 'continuous' or 'discrete' if `in_type[[i]]`='samples';
#' * 'gaussian', 'poisson' or 'nbinom' if `in_type[[i]]`='params'.
#'
#' If `distr` is a string it is assumed that all distributions are of the same type.
#'
#' @param num_samples number of samples drawn from the reconciled distribution.
#' @param seed seed for reproducibility.
#' @param num_samples Number of samples drawn from the reconciled distribution.
#' @param suppress_warnings Logical. If \code{TRUE}, no warnings about effective sample size
#' are triggered. If \code{FALSE}, warnings are generated. Default is \code{FALSE}. See Details.
#' @param seed Seed for reproducibility.
#'
#' @return A list containing the reconciled forecasts. The list has the following named elements:
#'
Expand Down Expand Up @@ -239,7 +247,7 @@ reconc_BUIS <- function(S,
in_type,
distr,
num_samples = 2e4,
suppressWarnings = FALSE, #TODO add to doc
suppress_warnings = FALSE,
seed = NULL) {
set.seed(seed)

Expand Down Expand Up @@ -314,12 +322,14 @@ reconc_BUIS <- function(S,
distr_ = distr_H[[hi]]
)
check_weights.res = .check_weigths(weights)
if (check_weights.res$warning & !suppressWarnings) {
if (check_weights.res$warning & !suppress_warnings) {
warning_msg = check_weights.res$warning_msg
# TODO add information to the string
# add information to the warning message
upper_fromS_i = which(lapply(seq_len(nrow(S)), function(i) sum(abs(S[i,] - c))) == 0)
warning_msg = paste(warning_msg, "Check upper ts at S-row-index:", upper_fromS_i)
warning(warning_msg)
for (wmsg in warning_msg) {
wmsg = paste(wmsg, paste0("Check the upper forecast at index: ", upper_fromS_i,"."))
warning(wmsg)
}
}
B[, b_mask] = .resample(B[, b_mask], weights)
}
Expand All @@ -337,18 +347,19 @@ reconc_BUIS <- function(S,
)
}
check_weights.res = .check_weigths(weights)
if (check_weights.res$warning & !suppressWarnings) {
if (check_weights.res$warning & !suppress_warnings) {
warning_msg = check_weights.res$warning_msg
# TODO add information to the string
# add information to the warning message
upper_fromS_i = c()
for (gi in 1:nrow(G)) {
c = G[gi, ]
upper_fromS_i = c(upper_fromS_i,
which(lapply(seq_len(nrow(S)), function(i) sum(abs(S[i,] - c))) == 0))
}
warning_msg = paste(warning_msg, "Check upper ts at S-row-index:",
paste0("{",paste(upper_fromS_i, collapse = ","), "}"))
warning(warning_msg)
for (wmsg in warning_msg) {
wmsg = paste(wmsg, paste0("Check the upper forecasts at index: ", paste0("{",paste(upper_fromS_i, collapse = ","), "}.")))
warning(wmsg)
}
}
B = .resample(B, weights)
}
Expand Down
25 changes: 18 additions & 7 deletions man/reconc_BUIS.Rd

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

19 changes: 13 additions & 6 deletions tests/testthat/test-weightswarn.R
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
test_that("Test shrinkage estimator", {
test_that("Test effective sample size", {
S = matrix(data = c(1,1,1,0,0,1), nrow=3, byrow = TRUE)

# -----------
Expand All @@ -17,6 +17,10 @@ test_that("Test shrinkage estimator", {
expect_equal(check_w$warning_code, 1)
expect_equal(check_w$n_eff, n)

# Try the warning message
# base_forecast = list(u,b1,b2)
# a = bayesRecon::reconc_BUIS(S, base_forecast, in_type = "samples", distr = list("continuous","discrete","discrete"), seed=42)

# -----------
n = 199
b1 = rpois(n=n, lambda = 3)
Expand All @@ -33,6 +37,10 @@ test_that("Test shrinkage estimator", {
expect_equal(check_w$warning_code, c(1,2))
expect_equal(check_w$n_eff, n)

# Try the warning message
# base_forecast = list(u,b1,b2)
# a = bayesRecon::reconc_BUIS(S, base_forecast, in_type = "samples", distr = list("continuous","discrete","discrete"), seed=42)

# -----------
n = 2000
b1 = rpois(n=n, lambda = 3)
Expand All @@ -48,9 +56,8 @@ test_that("Test shrinkage estimator", {
expect_equal(check_w$warning, TRUE)
expect_equal(check_w$warning_code, c(2,3))
expect_equal(check_w$n_eff < 200, TRUE)

# Try the warning message
# base_forecast = list(u,b1,b2)
# a = bayesRecon::reconc_BUIS(S, base_forecast, in_type = "samples", distr = list("continuous","discrete","discrete"), seed=42)
})


# base_forecast = list(u,b1,b2)
#
# a = bayesRecon::reconc_BUIS(S, base_forecast, in_type = "samples", distr = list("continuous","discrete","discrete"), seed=42)

0 comments on commit dd23119

Please sign in to comment.