From f98cbafdca1667404f40ac5bfb08f24e6e6af87a Mon Sep 17 00:00:00 2001 From: dsweber2 Date: Thu, 11 Jun 2026 14:04:51 -0500 Subject: [PATCH 1/6] dropping stale packages --- DESCRIPTION | 4 +--- NAMESPACE | 1 - R/model.R | 17 +++++++++++++---- R/utils.R | 12 +++--------- 4 files changed, 17 insertions(+), 17 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index db43f7e..73d6981 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -12,11 +12,9 @@ RoxygenNote: 7.3.2 Suggests: testthat (>= 3.0.0) Config/testthat/edition: 3 -Imports: +Imports: arrow, - covidcast, dplyr, - evalcast, english, jsonlite, lubridate, diff --git a/NAMESPACE b/NAMESPACE index 0c0d5a1..aa5d5b7 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -68,7 +68,6 @@ importFrom(dplyr,slice_max) importFrom(dplyr,starts_with) importFrom(dplyr,ungroup) importFrom(english,english) -importFrom(evalcast,weighted_interval_score) importFrom(jsonlite,read_json) importFrom(lubridate,days_in_month) importFrom(lubridate,make_date) diff --git a/R/model.R b/R/model.R index 1af536f..5d0e3c7 100644 --- a/R/model.R +++ b/R/model.R @@ -142,17 +142,26 @@ get_prediction <- function(test_data, taus, covariates, response, obj, return (as.data.frame(test_data)) } +#' Weighted interval score for a single observation +#' +#' Inlined from the evalcast package. +#' +#' @param taus Numeric vector of quantile levels. +#' @param residuals Numeric vector of (quantile_prediction - actual) values. +#' @param point_pred Unused; kept for interface compatibility. +#' @keywords internal +weighted_interval_score <- function(taus, residuals, point_pred) { + alpha <- 2 * pmin(taus, 1 - taus) + mean(alpha * (abs(residuals) + (residuals) * (2 * (taus >= 0.5) - 1))) +} + #' Evaluation of the test results based on WIS score -#' The WIS score calculation is based on the weighted_interval_score function -#' from the `evalcast` package from Delphi #' #' @param test_data dataframe with a column containing the prediction results of #' each requested quantile. Each row represents an update with certain #' (reference_date, report_date, location) combination. #' @template taus-template #' -#' @importFrom evalcast weighted_interval_score -#' #' @export evaluate <- function(test_data, taus, response) { n_row <- nrow(test_data) diff --git a/R/utils.R b/R/utils.R index 57f88f8..7103110 100644 --- a/R/utils.R +++ b/R/utils.R @@ -200,19 +200,13 @@ training_days_check <- function(report_date, training_days) { #' Subset list of counties to those included in the 200 most populous in the US #' +#' Requires the covidcast package, which is not installed by default. +#' #' @importFrom dplyr select %>% arrange desc pull #' @importFrom rlang .data #' @importFrom utils head get_populous_counties <- function() { - return( - covidcast::county_census %>% - dplyr::select(pop = .data$POPESTIMATE2019, fips = .data$FIPS) %>% - # Drop megacounties (states) - filter(!endsWith(.data$fips, "000")) %>% - arrange(desc(.data$pop)) %>% - pull(.data$fips) %>% - head(n=200) - ) + stop("get_populous_counties() requires the covidcast package. Install it with renv::install(\"cmu-delphi/covidcast/R-packages/covidcast\").") } #' Write a message to the console with the current time From c9068c9b3b0ef9742dbd76341eb740d398707155 Mon Sep 17 00:00:00 2001 From: dsweber2 Date: Thu, 9 Jul 2026 12:57:40 -0500 Subject: [PATCH 2/6] ref_lag not pinned to 7 days --- R/feature_engineering.R | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/R/feature_engineering.R b/R/feature_engineering.R index c253065..3d45ec8 100644 --- a/R/feature_engineering.R +++ b/R/feature_engineering.R @@ -188,8 +188,13 @@ add_lagged_terms <- function(df, value_col, refd_col, lag_col, lagged_term_list= #' @importFrom dplyr rename #' @export add_targets <- function(df, value_col, refd_col, lag_col, ref_lag, temporal_resol) { - # Add target - target_df <- df[df[[lag_col]]==ref_lag, c(refd_col, "report_date", value_col, "value_7dav")] + available_lags <- sort(unique(df[[lag_col]])) + effective_ref_lag <- min(available_lags[available_lags >= ref_lag]) + if (is.infinite(effective_ref_lag)) + stop(sprintf("ref_lag %d exceeds all available lags (max %d)", ref_lag, max(available_lags))) + if (effective_ref_lag != ref_lag) + message(sprintf("ref_lag %d not available; using next lag %d", ref_lag, effective_ref_lag)) + target_df <- df[df[[lag_col]] == effective_ref_lag, c(refd_col, "report_date", value_col, "value_7dav")] # Rename columns for clarity target_df <- target_df %>% dplyr::rename( From b6dfc4dc616abee11f0bc4a0170b097ffd19c99b Mon Sep 17 00:00:00 2001 From: dsweber2 Date: Mon, 27 Jul 2026 16:58:13 -0500 Subject: [PATCH 3/6] initial exogenous --- R/feature_engineering.R | 86 +++++++++++++++++++++++++++++++++++++++-- R/forecast.R | 6 ++- R/model.R | 9 ++--- 3 files changed, 90 insertions(+), 11 deletions(-) diff --git a/R/feature_engineering.R b/R/feature_engineering.R index 3d45ec8..445c1fa 100644 --- a/R/feature_engineering.R +++ b/R/feature_engineering.R @@ -273,6 +273,66 @@ add_params_for_dates <- function(df, refd_col, lag_col, temporal_resol="daily") return (as.data.frame(df)) } +#' Process an auxiliary reporting triangle into prefixed feature columns +#' +#' Runs a single auxiliary data frame (one geo, one signal) through the +#' fill → 7-day-average → lagged-terms → log-transform pipeline and +#' renames all value-derived columns with a `{name}_` prefix so they can +#' be safely joined to the primary preprocessed data frame. +#' +#' @param df Data frame with columns `reference_date`, `report_date`, `lag`, +#' and `value` (the signal values). +#' @param name Character scalar used as the column prefix (e.g. `"nssp"`). +#' @param lagged_term_list Numeric vector of lag values (same as used for the +#' primary signal). +#' @param temporal_resol `"daily"` or `"weekly"`. +#' @param smoothed Logical; if `FALSE` and `temporal_resol == "daily"`, a 7-day +#' moving average is computed. Otherwise `value_7dav` is set equal to +#' `value_raw`. +#' +#' @return Data frame with columns `reference_date`, `report_date`, `lag`, and +#' all value/log columns prefixed with `{name}_`. Useful predictors are +#' `{name}_log_value_7dav_lag{N}` and `{name}_log_delta_value_7dav_lag{N}`; +#' see [aux_feature_names()]. +process_aux_triangle <- function(df, name, lagged_term_list, temporal_resol, smoothed) { + filled_df <- fill_missing_updates(df, "value", "reference_date", "lag", temporal_resol) + if (!smoothed && temporal_resol == "daily") { + filled_df <- add_7davs(filled_df, "value_raw", "reference_date", "lag") + } else { + filled_df$value_7dav <- filled_df$value_raw + } + filled_df <- add_lagged_terms( + filled_df, "value_7dav", "reference_date", "lag", lagged_term_list, temporal_resol + ) + filled_df <- add_log_transformed(filled_df, lagged_term_list) + + value_cols <- grep("^(value_|log_)", colnames(filled_df), value = TRUE) + colnames(filled_df)[colnames(filled_df) %in% value_cols] <- paste0(name, "_", value_cols) + + filled_df[, c("reference_date", "report_date", "lag", paste0(name, "_", value_cols))] +} + + +#' Return the feature column names produced by an auxiliary triangle +#' +#' Gives the column names that [process_aux_triangle()] adds for a given +#' auxiliary signal, matching the log-value and log-delta features used by +#' the primary signal in [create_params_list()]. +#' +#' @param name Character scalar; the aux triangle name (must match what was +#' passed to [data_preprocessing()]). +#' @param lagged_term_list Numeric vector of lag values. +#' +#' @return Character vector of feature column names. +#' @export +aux_feature_names <- function(name, lagged_term_list) { + c( + paste0(name, "_log_value_7dav_lag", lagged_term_list), + paste0(name, "_log_delta_value_7dav_lag", lagged_term_list) + ) +} + + #' Data Preprocessing Function #' #' This function processes input data by handling missing values, computing lagged terms, @@ -292,14 +352,22 @@ add_params_for_dates <- function(df, refd_col, lag_col, temporal_resol="daily") #' @param value_type Character indicating the type of values ('count' or 'fraction'). #' @param temporal_resol Character specifying temporal resolution ('daily' or 'weekly'). #' @param smoothed Logical indicating whether smoothing should be applied. -#' -#' @importFrom dplyr full_join distinct +#' @param aux_triangles Named list of auxiliary reporting-triangle data frames. +#' Each element must have columns `reference_date`, `report_date`, `lag`, and +#' `value` (already filtered to the same single geo as `df`). Each is run +#' through the same fill → lag → log pipeline as the primary signal and +#' joined to the result; columns are prefixed with the list element name. +#' Use [aux_feature_names()] to obtain the resulting predictor column names +#' for [create_params_list()]. +#' +#' @importFrom dplyr full_join left_join distinct #' @importFrom english english #' #' @export data_preprocessing <- function(df, value_col, refd_col, lag_col, ref_lag, suffixes=c(""), lagged_term_list = NULL, value_type="count", - temporal_resol="daily", smoothed=FALSE) { + temporal_resol="daily", smoothed=FALSE, + aux_triangles = NULL) { if (value_type == "count") { if (length(value_col) > 1) warning("Multiple value column names provided; only the first one will be used.") if (length(unique(suffixes)) > 1) warning("Multiple suffixes provided; only the first one will be used.") @@ -374,6 +442,18 @@ data_preprocessing <- function(df, value_col, refd_col, lag_col, ref_lag, merged_df$inv_log_lag <- 1/(merged_df$lag + 1) merged_df <- add_params_for_dates(merged_df, "reference_date", "lag", temporal_resol) + if (!is.null(aux_triangles)) { + for (nm in names(aux_triangles)) { + aux_processed <- process_aux_triangle( + aux_triangles[[nm]], nm, lagged_term_list, temporal_resol, smoothed + ) + merged_df <- dplyr::left_join( + merged_df, aux_processed, + by = c("reference_date", "report_date", "lag") + ) + } + } + merged_df <- merged_df %>% filter(.data$lag < ref_lag) diff --git a/R/forecast.R b/R/forecast.R index 5ab1b21..f54127a 100644 --- a/R/forecast.R +++ b/R/forecast.R @@ -41,6 +41,7 @@ revision_forecast <- function(train_data, test_data, taus, smoothed_target=TRUE, lagged_term_list=NULL, params_list=NULL, + extra_params=NULL, temporal_resol="daily", lambda = 0.1, gamma = 0.1, lp_solver=LP_SOLVER, test_lag_group="", @@ -77,7 +78,7 @@ revision_forecast <- function(train_data, test_data, taus, } if (is.null(params_list)) { - params_list <- create_params_list(train_data, lagged_term_list, temporal_resol) + params_list <- create_params_list(train_data, lagged_term_list, temporal_resol, extra_params) } if (smoothed_target) { @@ -317,6 +318,7 @@ DelphiRF <- function(df, testing_start_date, taus=TAUS, smoothed_target=TRUE, lagged_term_list=NULL, params_list=NULL, + extra_params=NULL, lambda=LAMBDA, gamma=GAMMA, lag_pad=LAG_PAD, temporal_resol="daily", lp_solver=LP_SOLVER, @@ -388,7 +390,7 @@ DelphiRF <- function(df, testing_start_date, taus=TAUS, results <- revision_forecast(train_data, test_data, taus, smoothed_target, lagged_term_list, - params_list, temporal_resol, + params_list, extra_params, temporal_resol, l, g, lp_solver, test_lag_group, geo, value_type, model_save_dir, indicator, signal, geo_level, diff --git a/R/model.R b/R/model.R index 5d0e3c7..70006e3 100644 --- a/R/model.R +++ b/R/model.R @@ -341,7 +341,7 @@ generate_filename <- function(indicator, signal, #' #' @importFrom dplyr mutate select #' -create_params_list <- function(train_data, lagged_term_list, temporal_resol) { +create_params_list <- function(train_data, lagged_term_list, temporal_resol, extra_params = NULL) { params_list <- c( WEEK_ISSUES[1], Y7DAV, @@ -359,9 +359,6 @@ create_params_list <- function(train_data, lagged_term_list, temporal_resol) { paste0(dayofweek, "_issue") ) - if (temporal_resol == "daily"){ - return (c(params_list, extra_params_for_daily)) - } else { - return(params_list) - } + base_params <- if (temporal_resol == "daily") c(params_list, extra_params_for_daily) else params_list + if (!is.null(extra_params)) c(base_params, extra_params) else base_params } From b900add70d7227395a81c6236b3ba0dc08a03c94 Mon Sep 17 00:00:00 2001 From: dsweber2 Date: Tue, 28 Jul 2026 10:47:48 -0500 Subject: [PATCH 4/6] deal with timing issues --- DESCRIPTION | 2 +- NAMESPACE | 1 + R/feature_engineering.R | 11 +++-- R/preprocessing.R | 11 ++++- tests/testthat/test-feature_engineering.R | 51 +++++++++++++++++++++++ 5 files changed, 71 insertions(+), 5 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 73d6981..f929e62 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -8,7 +8,7 @@ Description: Provide real-time revision forecasts. License: MIT + file LICENSE Encoding: UTF-8 Roxygen: list(markdown = TRUE) -RoxygenNote: 7.3.2 +RoxygenNote: 7.3.3 Suggests: testthat (>= 3.0.0) Config/testthat/edition: 3 diff --git a/NAMESPACE b/NAMESPACE index aa5d5b7..310924f 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -21,6 +21,7 @@ export(WEEK_ISSUES) export(Y7DAV) export(YITL) export(add_7davs) +export(aux_feature_names) export(add_dayofweek) export(add_lagged_terms) export(add_log_transformed) diff --git a/R/feature_engineering.R b/R/feature_engineering.R index 445c1fa..03710e6 100644 --- a/R/feature_engineering.R +++ b/R/feature_engineering.R @@ -294,8 +294,10 @@ add_params_for_dates <- function(df, refd_col, lag_col, temporal_resol="daily") #' all value/log columns prefixed with `{name}_`. Useful predictors are #' `{name}_log_value_7dav_lag{N}` and `{name}_log_delta_value_7dav_lag{N}`; #' see [aux_feature_names()]. -process_aux_triangle <- function(df, name, lagged_term_list, temporal_resol, smoothed) { - filled_df <- fill_missing_updates(df, "value", "reference_date", "lag", temporal_resol) +process_aux_triangle <- function(df, name, lagged_term_list, temporal_resol, smoothed, + max_report_override = NULL) { + filled_df <- fill_missing_updates(df, "value", "reference_date", "lag", temporal_resol, + max_report_override = max_report_override) if (!smoothed && temporal_resol == "daily") { filled_df <- add_7davs(filled_df, "value_raw", "reference_date", "lag") } else { @@ -443,14 +445,17 @@ data_preprocessing <- function(df, value_col, refd_col, lag_col, ref_lag, merged_df <- add_params_for_dates(merged_df, "reference_date", "lag", temporal_resol) if (!is.null(aux_triangles)) { + primary_max_report <- max(merged_df$report_date) for (nm in names(aux_triangles)) { aux_processed <- process_aux_triangle( - aux_triangles[[nm]], nm, lagged_term_list, temporal_resol, smoothed + aux_triangles[[nm]], nm, lagged_term_list, temporal_resol, smoothed, + max_report_override = primary_max_report ) merged_df <- dplyr::left_join( merged_df, aux_processed, by = c("reference_date", "report_date", "lag") ) + merged_df <- merged_df[!is.na(merged_df[[paste0(nm, "_value_raw")]]), ] } } diff --git a/R/preprocessing.R b/R/preprocessing.R index 4371cec..20a10e2 100644 --- a/R/preprocessing.R +++ b/R/preprocessing.R @@ -56,7 +56,8 @@ fill_rows <- function(df, refd_col, lag_col, min_refd, max_refd, ref_lag) { #' @importFrom tidyr fill pivot_wider pivot_longer replace_na expand_grid #' @importFrom dplyr %>% select left_join mutate arrange distinct filter everything #' @export -fill_missing_updates <- function(df, value_col, refd_col, lag_col, temporal_resol="daily") { +fill_missing_updates <- function(df, value_col, refd_col, lag_col, temporal_resol = "daily", + max_report_override = NULL) { df <- df %>% distinct() # Remove duplicates if any if (nrow(df) == 0) { @@ -88,6 +89,14 @@ fill_missing_updates <- function(df, value_col, refd_col, lag_col, temporal_reso stop("Invalid temporal_resol. Choose either 'daily' or 'weekly'.") } + if (!is.null(max_report_override)) { + max_report_override <- as.Date(max_report_override) + if (max_report_override > max(all_report_dates)) { + extra <- seq(max(all_report_dates) + gap, max_report_override, by = gap) + all_report_dates <- c(all_report_dates, extra) + } + } + # Create a complete grid of all combinations of reference_date and report_date complete_grid <- tidyr::expand_grid( !!refd_col := all_reference_dates, diff --git a/tests/testthat/test-feature_engineering.R b/tests/testthat/test-feature_engineering.R index c4a577f..13b81c4 100644 --- a/tests/testthat/test-feature_engineering.R +++ b/tests/testthat/test-feature_engineering.R @@ -341,3 +341,54 @@ test_that("Testing add weighted related features", { expect_true(all(expected_columns %in% colnames(result))) }) + + +make_daily_tri <- function(ref_dates, lags) { + do.call(rbind, lapply(ref_dates, function(rd) { + data.frame( + reference_date = rd, + lag = lags, + value = as.numeric(as.Date(rd) - as.Date("2022-12-31")) + ) + })) +} + +test_that("data_preprocessing drops rows where aux has no coverage (earlier min)", { + ref_all <- seq(as.Date("2023-01-01"), as.Date("2023-01-15"), by = "day") + ref_late <- seq(as.Date("2023-01-08"), as.Date("2023-01-15"), by = "day") + + primary <- make_daily_tri(ref_all, 1:3) + aux <- make_daily_tri(ref_late, 1:3) + + result <- data_preprocessing( + primary, + value_col = "value", refd_col = "reference_date", lag_col = "lag", + ref_lag = 5L, temporal_resol = "daily", + aux_triangles = list(beds = aux) + ) + + expect_true(all(result$reference_date >= as.Date("2023-01-08"))) + expect_true("beds_value_raw" %in% colnames(result)) + expect_false(anyNA(result[["beds_value_raw"]])) +}) + +test_that("data_preprocessing forward-fills aux to primary max report_date (later max)", { + ref_dates <- seq(as.Date("2023-01-01"), as.Date("2023-01-10"), by = "day") + primary <- make_daily_tri(ref_dates, 1:5) + aux <- make_daily_tri(ref_dates, 1:3) + # primary max report_date = 2023-01-10 + 5 = 2023-01-15 + # aux max report_date = 2023-01-10 + 3 = 2023-01-13 + + result <- data_preprocessing( + primary, + value_col = "value", refd_col = "reference_date", lag_col = "lag", + ref_lag = 6L, temporal_resol = "daily", + aux_triangles = list(beds = aux) + ) + + # rows with lag 4 and 5 (report_dates past aux max) should be present + expect_true(any(result$lag == 4L)) + expect_true(any(result$lag == 5L)) + expect_true("beds_value_raw" %in% colnames(result)) + expect_false(anyNA(result[["beds_value_raw"]])) +}) From 787882839448f04e4d2ebef4767c2ccd1f60db29 Mon Sep 17 00:00:00 2001 From: dsweber2 Date: Tue, 28 Jul 2026 11:34:11 -0500 Subject: [PATCH 5/6] guard against either being empty --- R/feature_engineering.R | 10 +++++++ tests/testthat/test-feature_engineering.R | 36 +++++++++++++++++++++++ 2 files changed, 46 insertions(+) diff --git a/R/feature_engineering.R b/R/feature_engineering.R index 03710e6..ed35057 100644 --- a/R/feature_engineering.R +++ b/R/feature_engineering.R @@ -298,6 +298,12 @@ process_aux_triangle <- function(df, name, lagged_term_list, temporal_resol, smo max_report_override = NULL) { filled_df <- fill_missing_updates(df, "value", "reference_date", "lag", temporal_resol, max_report_override = max_report_override) + if (nrow(filled_df) == 0) { + expected_cols <- c("reference_date", "report_date", "lag", + aux_feature_names(name, lagged_term_list)) + return(setNames(data.frame(matrix(ncol = length(expected_cols), nrow = 0)), + expected_cols)) + } if (!smoothed && temporal_resol == "daily") { filled_df <- add_7davs(filled_df, "value_raw", "reference_date", "lag") } else { @@ -451,6 +457,10 @@ data_preprocessing <- function(df, value_col, refd_col, lag_col, ref_lag, aux_triangles[[nm]], nm, lagged_term_list, temporal_resol, smoothed, max_report_override = primary_max_report ) + if (nrow(aux_processed) == 0L) { + merged_df <- merged_df[0L, ] + break + } merged_df <- dplyr::left_join( merged_df, aux_processed, by = c("reference_date", "report_date", "lag") diff --git a/tests/testthat/test-feature_engineering.R b/tests/testthat/test-feature_engineering.R index 13b81c4..fc79e77 100644 --- a/tests/testthat/test-feature_engineering.R +++ b/tests/testthat/test-feature_engineering.R @@ -392,3 +392,39 @@ test_that("data_preprocessing forward-fills aux to primary max report_date (late expect_true("beds_value_raw" %in% colnames(result)) expect_false(anyNA(result[["beds_value_raw"]])) }) + +test_that("data_preprocessing returns 0 rows when aux has no data for the geo", { + ref_dates <- seq(as.Date("2023-01-01"), as.Date("2023-01-15"), by = "day") + primary <- make_daily_tri(ref_dates, 1:4) + empty_aux <- data.frame( + reference_date = as.Date(character()), + report_date = as.Date(character()), + lag = integer(), + value = numeric() + ) + + result <- data_preprocessing( + primary, + value_col = "value", refd_col = "reference_date", lag_col = "lag", + ref_lag = 5L, temporal_resol = "daily", + aux_triangles = list(beds = empty_aux) + ) + + expect_equal(nrow(result), 0L) +}) + +test_that("process_aux_triangle returns zero-row df with correct columns when aux is empty", { + empty_aux <- data.frame( + reference_date = as.Date(character()), + report_date = as.Date(character()), + lag = integer(), + value = numeric() + ) + lagged_term_list <- c(1L, 7L) + + result <- process_aux_triangle(empty_aux, "beds", lagged_term_list, "daily", TRUE) + + expect_equal(nrow(result), 0L) + expect_true(all(aux_feature_names("beds", lagged_term_list) %in% colnames(result))) + expect_true(all(c("reference_date", "report_date", "lag") %in% colnames(result))) +}) From 26beefa1a39ddad2e43efcc7e86ee9e817b678f3 Mon Sep 17 00:00:00 2001 From: dsweber2 Date: Mon, 3 Aug 2026 13:28:40 -0500 Subject: [PATCH 6/6] one-hot flexibility --- R/feature_engineering.R | 59 +++++++++++++---- R/forecast.R | 12 ++-- R/model.R | 18 +++--- tests/testthat/test-feature_engineering.R | 78 ++++++++++++++--------- 4 files changed, 109 insertions(+), 58 deletions(-) diff --git a/R/feature_engineering.R b/R/feature_engineering.R index ed35057..094753a 100644 --- a/R/feature_engineering.R +++ b/R/feature_engineering.R @@ -32,6 +32,37 @@ add_dayofweek <- function(df, time_col, suffix, wd = WEEKDAYS_ABBR) { } +#' Add grouped day-of-week one-hot columns +#' +#' Creates one binary column per group in `onehot_weekdays`. Each group is a +#' character vector of day abbreviations (from `WEEKDAYS_ABBR`); a row is 1 if +#' the date falls on any day in the group. Column names are derived from the +#' list names when present, otherwise by concatenating the day abbreviations. +#' +#' @param df A data frame containing the date column. +#' @param time_col Name of the date column. +#' @param suffix Column name suffix (e.g. `"_ref"` or `"_issue"`). +#' @param onehot_weekdays Named or unnamed list of character vectors, each +#' specifying a group of days (e.g. `list(Mon=c("Mon"), Weekends=c("Sat","Sun"))`). +#' +#' @return `df` with one additional integer column per group. +#' @export +add_grouped_dayofweek <- function(df, time_col, suffix, onehot_weekdays) { + df <- df %>% mutate({{ time_col }} := as.Date(.data[[time_col]])) + dayofweek <- as.numeric(format(df[[time_col]], format = "%u")) + group_names <- if (!is.null(names(onehot_weekdays))) { + names(onehot_weekdays) + } else { + vapply(onehot_weekdays, function(grp) paste0(grp, collapse = ""), character(1)) + } + for (ii in seq_along(onehot_weekdays)) { + day_indices <- match(onehot_weekdays[[ii]], WEEKDAYS_ABBR) + df[[paste0(group_names[ii], suffix)]] <- as.integer(dayofweek %in% day_indices) + } + df +} + + #' Add one-hot encoding for week of the month based on issue date #' #' This function calculates the week of the month for each date in the specified @@ -248,29 +279,28 @@ add_log_transformed <- function(df, lagged_term_list) { #' @param lag_col Column name representing the lag between the reference and issue date. #' @param temporal_resol A string indicating the temporal resolution ("daily" or "weekly"). #' Defaults to "daily". +#' @param onehot_weekdays Named or unnamed list of character vectors defining day groups +#' for one-hot encoding. Each group becomes one binary column per date column. +#' Names, when present, are used as column prefixes; otherwise day abbreviations are +#' concatenated. Defaults to `list(Mon=c("Mon"), Weekends=c("Sat","Sun"))`. #' #' @details -#' - If `temporal_resol` is "daily", one-hot encoded day-of-week columns are added +#' - If `temporal_resol` is "daily", one-hot encoded day-of-week group columns are added #' for both `refd_col` (reference date) and `"report_date"`. #' - One-hot encoded week-of-month columns are added for `"report_date"` in all cases. #' #' @return A modified data frame with additional date-related feature columns. #' #' @export -add_params_for_dates <- function(df, refd_col, lag_col, temporal_resol="daily") { +add_params_for_dates <- function(df, refd_col, lag_col, temporal_resol = "daily", + onehot_weekdays = list(Mon = c("Mon"), Weekends = c("Sat", "Sun"))) { df$report_date <- df[[refd_col]] + df[[lag_col]] - if (temporal_resol=="daily"){ - # Add columns for day-of-week effect - df <- add_dayofweek(df, refd_col, "_ref", WEEKDAYS_ABBR) - df <- add_dayofweek(df, "report_date", "_issue", WEEKDAYS_ABBR) - # Add columns for weekends - df$Weekends_issue <- as.integer(df$Sat_issue == 1 | df$Sun_issue == 1) - df$Weekends_ref <- as.integer(df$Sat_ref == 1 | df$Sun_ref == 1) + if (temporal_resol == "daily") { + df <- add_grouped_dayofweek(df, refd_col, "_ref", onehot_weekdays) + df <- add_grouped_dayofweek(df, "report_date", "_issue", onehot_weekdays) } - # Add columns for week-of-month effect df <- add_weekofmonth(df, "report_date", WEEK_ISSUES) - - return (as.data.frame(df)) + return(as.data.frame(df)) } #' Process an auxiliary reporting triangle into prefixed feature columns @@ -375,7 +405,8 @@ aux_feature_names <- function(name, lagged_term_list) { data_preprocessing <- function(df, value_col, refd_col, lag_col, ref_lag, suffixes=c(""), lagged_term_list = NULL, value_type="count", temporal_resol="daily", smoothed=FALSE, - aux_triangles = NULL) { + aux_triangles = NULL, + onehot_weekdays = list(Mon = c("Mon"), Weekends = c("Sat", "Sun"))) { if (value_type == "count") { if (length(value_col) > 1) warning("Multiple value column names provided; only the first one will be used.") if (length(unique(suffixes)) > 1) warning("Multiple suffixes provided; only the first one will be used.") @@ -448,7 +479,7 @@ data_preprocessing <- function(df, value_col, refd_col, lag_col, ref_lag, } merged_df$inv_log_lag <- 1/(merged_df$lag + 1) - merged_df <- add_params_for_dates(merged_df, "reference_date", "lag", temporal_resol) + merged_df <- add_params_for_dates(merged_df, "reference_date", "lag", temporal_resol, onehot_weekdays) if (!is.null(aux_triangles)) { primary_max_report <- max(merged_df$report_date) diff --git a/R/forecast.R b/R/forecast.R index f54127a..5ebb394 100644 --- a/R/forecast.R +++ b/R/forecast.R @@ -50,9 +50,10 @@ revision_forecast <- function(train_data, test_data, taus, indicator="testdata", signal="", geo_level="state", signal_suffix="", training_end_date="", - training_days =365, + training_days=365, train_models = TRUE, - make_predictions=TRUE) { + make_predictions=TRUE, + onehot_weekdays = list(Mon = c("Mon"), Weekends = c("Sat", "Sun"))) { @@ -78,7 +79,7 @@ revision_forecast <- function(train_data, test_data, taus, } if (is.null(params_list)) { - params_list <- create_params_list(train_data, lagged_term_list, temporal_resol, extra_params) + params_list <- create_params_list(train_data, lagged_term_list, temporal_resol, onehot_weekdays, extra_params) } if (smoothed_target) { @@ -329,7 +330,8 @@ DelphiRF <- function(df, testing_start_date, taus=TAUS, training_end_date="", training_days=365, train_models = TRUE, - make_predictions = TRUE) { + make_predictions = TRUE, + onehot_weekdays = list(Mon = c("Mon"), Weekends = c("Sat", "Sun"))) { testing_start_date <- as.Date(testing_start_date) @@ -396,7 +398,7 @@ DelphiRF <- function(df, testing_start_date, taus=TAUS, indicator, signal, geo_level, signal_suffix, as.character(testing_start_date), training_days, train_models, - make_predictions) + make_predictions, onehot_weekdays) test_data_list <- append(test_data_list, list(results)) } diff --git a/R/model.R b/R/model.R index 70006e3..6cddddc 100644 --- a/R/model.R +++ b/R/model.R @@ -341,23 +341,25 @@ generate_filename <- function(indicator, signal, #' #' @importFrom dplyr mutate select #' -create_params_list <- function(train_data, lagged_term_list, temporal_resol, extra_params = NULL) { +create_params_list <- function(train_data, lagged_term_list, temporal_resol, + onehot_weekdays = list(Mon = c("Mon"), Weekends = c("Sat", "Sun")), + extra_params = NULL) { params_list <- c( WEEK_ISSUES[1], Y7DAV, paste0("log_value_7dav_lag", lagged_term_list), paste0("log_delta_value_7dav_lag", lagged_term_list) ) - # Include log lag adjustments if multiple lags exist - if (length(unique(train_data$lag)) > 1){ + if (length(unique(train_data$lag)) > 1) { params_list <- c(params_list, LOG_LAG) } - dayofweek <- c("Mon", "Weekends") - extra_params_for_daily <- c( - paste0(dayofweek, "_ref"), - paste0(dayofweek, "_issue") - ) + group_names <- if (!is.null(names(onehot_weekdays))) { + names(onehot_weekdays) + } else { + vapply(onehot_weekdays, function(grp) paste0(grp, collapse = ""), character(1)) + } + extra_params_for_daily <- c(paste0(group_names, "_ref"), paste0(group_names, "_issue")) base_params <- if (temporal_resol == "daily") c(params_list, extra_params_for_daily) else params_list if (!is.null(extra_params)) c(base_params, extra_params) else base_params diff --git a/tests/testthat/test-feature_engineering.R b/tests/testthat/test-feature_engineering.R index fc79e77..8c86b22 100644 --- a/tests/testthat/test-feature_engineering.R +++ b/tests/testthat/test-feature_engineering.R @@ -189,57 +189,73 @@ test_that("add_log_transformed correctly applies log transformation", { }) -test_that("add_params_for_dates correctly adds date-related features", { - # Create a sample data frame +test_that("add_grouped_dayofweek creates one column per group with correct values", { + # Mon=1, Tue=2, ..., Sun=7 via %u format + df <- data.frame( + date = as.Date(c("2024-02-26", "2024-02-27", "2024-02-28", "2024-03-01", "2024-03-02", "2024-03-03")) + # Mon, Tue, Wed, Fri, Sat, Sun + ) + groups <- list(Mon = c("Mon"), Weekends = c("Sat", "Sun"), Other = c("Tue", "Wed", "Thurs", "Fri")) + result <- add_grouped_dayofweek(df, "date", "_ref", groups) + + expect_true(all(c("Mon_ref", "Weekends_ref", "Other_ref") %in% colnames(result))) + expect_equal(result$Mon_ref, c(1L, 0L, 0L, 0L, 0L, 0L)) + expect_equal(result$Weekends_ref, c(0L, 0L, 0L, 0L, 1L, 1L)) + expect_equal(result$Other_ref, c(0L, 1L, 1L, 1L, 0L, 0L)) + # every row sums to exactly 1 (exhaustive, non-overlapping groups) + expect_true(all(rowSums(result[, c("Mon_ref", "Weekends_ref", "Other_ref")]) == 1L)) +}) + +test_that("add_grouped_dayofweek derives column names from day abbreviations when unnamed", { + df <- data.frame(date = as.Date(c("2024-02-26", "2024-03-02"))) # Mon, Sat + result <- add_grouped_dayofweek(df, "date", "_ref", list(c("Mon"), c("Sat", "Sun"))) + expect_true("Mon_ref" %in% colnames(result)) + expect_true("SatSun_ref" %in% colnames(result)) +}) + +test_that("add_params_for_dates uses default Mon/Weekends groups in daily mode", { test_df <- data.frame( ref_date = as.Date(c("2022-01-01", "2022-01-05", "2022-01-10", "2022-02-01", "2022-02-15")), lag = c(0, 2, 5, 7, 10) ) - - # Run the function with daily resolution df_with_params <- add_params_for_dates(test_df, "ref_date", "lag", "daily") - # Check that report_date column is correctly created expect_true("report_date" %in% colnames(df_with_params)) - - # Verify day-of-week encoding is added for both reference and issue date - expect_true(all(paste0(WEEKDAYS_ABBR, "_ref") %in% colnames(df_with_params))) - expect_true(all(paste0(WEEKDAYS_ABBR, "_issue") %in% colnames(df_with_params))) - - expect_true(all(c("Weekends_issue", "Weekends_ref") %in% colnames(df_with_params))) - - # Verify that exactly one column per row is 1 for each set of one-hot encoded days - expect_true(all(rowSums(df_with_params[, paste0(WEEKDAYS_ABBR, "_ref")]) == 1)) - expect_true(all(rowSums(df_with_params[, paste0(WEEKDAYS_ABBR, "_issue")]) == 1)) - - # Verify week-of-month encoding is added for report_date + expect_true(all(c("Mon_ref", "Weekends_ref", "Mon_issue", "Weekends_issue") %in% colnames(df_with_params))) + # old per-day columns should not be present + expect_false(any(c("Tue_ref", "Wed_ref", "Thurs_ref", "Fri_ref", "Sat_ref", "Sun_ref") %in% colnames(df_with_params))) + # each row is either Mon (1,0), Weekends (0,1), or Other (0,0) — never (1,1) + expect_true(all(df_with_params$Mon_ref + df_with_params$Weekends_ref <= 1L)) + expect_true(all(df_with_params$Mon_issue + df_with_params$Weekends_issue <= 1L)) expect_true(all(WEEK_ISSUES %in% colnames(df_with_params))) - - # Check that only one week column per row has a value of 1 expect_true(all(rowSums(df_with_params[, WEEK_ISSUES]) <= 1)) }) +test_that("add_params_for_dates respects custom onehot_weekdays", { + # Dates: 2024-02-26=Mon, 2024-02-28=Wed, 2024-03-01=Fri, 2024-03-02=Sat + test_df <- data.frame( + ref_date = as.Date(c("2024-02-26", "2024-02-28", "2024-03-01", "2024-03-02")), + lag = c(0, 0, 0, 0) + ) + groups <- list(WedFri = c("Wed", "Fri"), Other = c("Mon", "Tue", "Thurs", "Sat", "Sun")) + result <- add_params_for_dates(test_df, "ref_date", "lag", "daily", onehot_weekdays = groups) + + expect_true(all(c("WedFri_ref", "Other_ref") %in% colnames(result))) + expect_equal(result$WedFri_ref, c(0L, 1L, 1L, 0L)) + expect_equal(result$Other_ref, c(1L, 0L, 0L, 1L)) +}) + test_that("add_params_for_dates correctly handles weekly resolution", { - # Create a sample data frame test_df <- data.frame( ref_date = as.Date(c("2022-03-01", "2022-03-08", "2022-03-15")), lag = c(0, 7, 14) ) - - # Run the function with weekly resolution df_with_params <- add_params_for_dates(test_df, "ref_date", "lag", "weekly") - # Check that report_date column is correctly created expect_true("report_date" %in% colnames(df_with_params)) - - # Verify that day-of-week encoding is NOT added in weekly mode - expect_false(any(paste0(WEEKDAYS_ABBR, "_ref") %in% colnames(df_with_params))) - expect_false(any(paste0(WEEKDAYS_ABBR, "_issue") %in% colnames(df_with_params))) - - # Verify week-of-month encoding is still added for report_date + # no day-of-week columns in weekly mode + expect_false(any(c("Mon_ref", "Weekends_ref", "Mon_issue", "Weekends_issue") %in% colnames(df_with_params))) expect_true(all(WEEK_ISSUES %in% colnames(df_with_params))) - - # Check that only one week column per row has a value of 1 expect_true(all(rowSums(df_with_params[, WEEK_ISSUES]) <= 1)) })