fixing covid explore wip

Author: dsweber2

R/utils.R

@@ -431,9 +431,10 @@ sort_by_quantile <- function(forecasts) {
 
 
 #' Print recent targets errors.
-get_recent_targets_errors <- function(recent_minutes = 60) {
-  forecast_errors <- targets::tar_meta() %>%
-    filter(time > Sys.time() - minutes(recent_minutes), !is.na(parent), !is.na(error)) %>%
+get_recent_targets_errors <- function(time_since = minutes(60)) {
+  meta_df <- targets::tar_meta()
+  forecast_errors <- meta_df %>%
+    filter(time > Sys.time() - time_since, !is.na(parent), !is.na(error)) %>%
     arrange(desc(time)) %>%
     distinct(parent, error, .keep_all = TRUE) %>%
     select(time, parent, error) %>%
@@ -451,8 +452,8 @@ get_recent_targets_errors <- function(recent_minutes = 60) {
     }
   }
 
-  other_errors <- targets::tar_meta() %>%
-    filter(time > Sys.time() - minutes(recent_minutes), !is.na(error)) %>%
+  other_errors <- meta_df %>%
+    filter(time > Sys.time() - time_since, !is.na(error)) %>%
     arrange(desc(time)) %>%
     distinct(error, .keep_all = TRUE) %>%
     select(time, name, error)
@@ -468,4 +469,6 @@ get_recent_targets_errors <- function(recent_minutes = 60) {
       ))
     }
   }
+
+  return(invisible(meta_df %>% filter(time > Sys.time() - time_since)))
 }

reports/template.md

@@ -24,7 +24,14 @@ Simplistic/low data methods:
 - [Flu flatline](flu-notebook-flatline.html)
 - [Flu climate](flu-notebook-climate_linear.html)
 
-### Covid
+### Covid (new)
+
+- [Covid AR with population scaling](covid-notebook-scaled_pop_main.html)
+- [Covid AR with population scaling and seasonal features](covid-notebook-scaled_pop_season.html)
+- [Covid AR with population scaling, and exogenous features](covid-notebook-scaled_pop_exogenous.html)
+- [Covid Flatline](covid-notebook-flatline_forecaster.html)
+
+### Covid (old)
 
 - [Covid AR with population scaling](covid-notebook-1.html)
 - [Covid AR with population scaling and smoothed features](covid-notebook-2.html)

scripts/covid_hosp_explore.R

@@ -19,8 +19,8 @@ fetch_args <- epidatr::fetch_args_list(return_empty = FALSE, timeout_seconds = 4
 # with prototyping the pipeline.
 dummy_mode <- as.logical(Sys.getenv("DUMMY_MODE", FALSE))
 # For testing, reduce dates
-forecast_generation_dates <- forecast_generation_dates[1:10]
-forecast_dates <- forecast_dates[1:10]
+# forecast_generation_dates <- forecast_generation_dates[1:10]
+# forecast_dates <- forecast_dates[1:10]
 
 
 # ================================ FORECASTER PARAMETERS ====================
@@ -199,6 +199,7 @@ get_partially_applied_forecaster <- function(id) {
   }
 }
 
+
 # ================================ TARGETS =================================
 # ================================ PARAMETERS TARGETS ======================
 parameter_targets <- list2(
@@ -232,7 +233,7 @@ data_targets <- list2(
     }
   ),
   tar_target(
-    name = hhs_latest_data,
+    name = hhs_evaluation_data,
     command = {
       epidatr::pub_covidcast(
         source = "hhs",
@@ -242,20 +243,20 @@ data_targets <- list2(
         geo_values = "*",
         time_values = "*",
         fetch_args = fetch_args
-      )
+      ) %>%
+        select(signal, geo_value, time_value, value) %>%
+        daily_to_weekly(keys = c("geo_value", "signal")) %>%
+        select(signal, geo_value, target_end_date = time_value, true_value = value) %>%
+        # Correction for timing offsets
+        mutate(target_end_date = target_end_date + 3)
     }
   ),
   tar_target(
-    name = hhs_evaluation_data,
+    name = state_geo_values,
     command = {
-      hhs_latest_data %>%
-        select(signal, geo_value, time_value, value) %>%
-        daily_to_weekly(keys = c("geo_value", "signal")) %>%
-        rename(
-          true_value = value,
-          target_end_date = time_value
-        ) %>%
-        select(signal, geo_value, target_end_date, true_value)
+      hhs_evaluation_data %>%
+        pull(geo_value) %>%
+        unique()
     }
   ),
   tar_target(
@@ -490,13 +491,7 @@ forecasts_and_scores <- tar_map(
         forecast_scaled <- forecast
         actual_eval_data <- hhs_evaluation_data
       }
-      # Fix for timing offsets
-      actual_eval_data <- actual_eval_data %>% mutate(target_end_date = target_end_date + 3)
-      state_geo_values <- actual_eval_data %>%
-        pull(geo_value) %>%
-        unique()
       forecast_scaled <- forecast_scaled %>%
-        filter(geo_value %in% state_geo_values) %>%
         mutate(forecast_date = forecast_date + 3, target_end_date = target_end_date + 3) %>%
         rename("model" = "id")
 
@@ -512,34 +507,52 @@ combined_forecasts_and_scores <- rlang::list2(
   tar_combine(
     delphi_forecasts,
     forecasts_and_scores[["forecast"]],
-    command = dplyr::bind_rows(!!!.x) %>% rename(forecaster = id)
+    command = {
+      dplyr::bind_rows(!!!.x) %>%
+        rename(forecaster = id) %>%
+        filter(geo_value %in% state_geo_values) %>%
+        mutate(forecast_date = forecast_date + 3, target_end_date = target_end_date + 3)
+    }
   ),
   tar_combine(
     delphi_scores,
     forecasts_and_scores[["score"]],
-    command = dplyr::bind_rows(!!!.x) %>% rename(forecaster = id)
+    command = {
+      dplyr::bind_rows(!!!.x) %>%
+        rename(forecaster = id) %>%
+        filter(geo_value %in% state_geo_values)
+    }
   ),
 )
 external_forecasts_and_scores <- rlang::list2(
   tar_target(
-    external_forecasts,
+    outside_forecaster_subset,
+    command = c("COVIDhub-baseline", "COVIDhub-trained_ensemble", "COVIDhub_CDC-ensemble")
+  ),
+  tar_target(
+    external_forecasts_file,
     command = {
       s3load("covid19_forecast_hub_2023.rds", bucket = "forecasting-team-data", verbose = FALSE)
       full_results
     }
   ),
   tar_target(
-    external_scores,
+    external_forecasts,
     command = {
-      actual_eval_data <- hhs_evaluation_data %>%
-        mutate(target_end_date = target_end_date + 3)
-      cmu_forecast_dates <- ref_time_values + 3
-      filtered_forecasts <- external_forecasts %>%
+      external_forecasts_file %>%
+        filter(geo_value %in% state_geo_values, forecaster %in% outside_forecaster_subset) %>%
         mutate(forecast_date = forecast_date + 5, target_end_date = target_end_date + 5) %>%
-        filter(forecast_date %in% cmu_forecast_dates) %>%
-        rename(model = forecaster) %>%
+        filter(forecast_date %in% (ref_time_values + 3)) %>%
         rename(prediction = value) %>%
-        filter(!is.na(geo_value))
+        mutate(prediction = prediction * 7)
+    }
+  ),
+  tar_target(
+    external_scores,
+    command = {
+      actual_eval_data <- hhs_evaluation_data
+      filtered_forecasts <- external_forecasts %>%
+        rename(model = forecaster)
       evaluate_predictions(forecasts = filtered_forecasts, truth_data = actual_eval_data) %>%
         rename(forecaster = model)
     }
@@ -551,14 +564,28 @@ joined_forecasts_and_scores <- rlang::list2(
   tar_target(
     family_notebooks,
     command = {
-      actual_eval_data <- hhs_evaluation_data %>%
-        mutate(target_end_date = target_end_date + 3)
+      actual_eval_data <- hhs_evaluation_data
       delphi_forecaster_subset <- forecaster_parameter_combinations[[forecaster_families]]$id
-      outside_forecaster_subset <- c("COVIDhub-baseline", "COVIDhub-ensemble")
+
       filtered_forecasts <- joined_forecasts %>%
         filter(forecaster %in% c(delphi_forecaster_subset, outside_forecaster_subset))
       filtered_scores <- joined_scores %>%
         filter(forecaster %in% c(delphi_forecaster_subset, outside_forecaster_subset))
+
+      # TODO: Write an assert to make sure that these dates are similar. It's a bit tricky.
+      # actual_eval_data %>%
+      #   filter(target_end_date > "2023-09-01") %>%
+      #   distinct(target_end_date) %>%
+      #   pull(target_end_date) %>%
+      #   sort()
+      # filtered_forecasts %>%
+      #   distinct(target_end_date) %>%
+      #   pull(target_end_date) %>%
+      #   sort()
+      # filtered_scores %>%
+      #   distinct(target_end_date) %>%
+      #   pull(target_end_date) %>%
+      #   sort()
       forecaster_parameters <- forecaster_parameter_combinations[[forecaster_families]]
       rmarkdown::render(
         "scripts/reports/comparison-notebook.Rmd",

scripts/reports/comparison-notebook.Rmd

@@ -40,23 +40,23 @@ library(purrr)
 ```
 
 ```{r}
+# outside_forecaster_subset <- c("COVIDhub-baseline", "COVIDhub-trained_ensemble", "COVIDhub_CDC-ensemble")
+# i <- 1
 # params <- list(
-#     forecaster_parameter_combinations = tar_read(forecaster_parameter_combinations),
-#     forecaster_family = "scaled_pop_season",
-#     joined_forecasts = tar_read(joined_forecasts),
-#     joined_scores = tar_read(joined_scores),
-#     truth_data = tar_read(hhs_evaluation_data) %>%
-#         select(-population) %>%
-#         mutate(target_end_date = target_end_date + 3),
-#     disease = "flu"
+#   forecaster_family = forecaster_families[[i]],
+#   forecaster_parameters = forecaster_parameter_combinations[[i]],
+#   forecasts = tar_read(joined_forecasts) %>% filter(forecaster %in% c(forecaster_parameter_combinations[[i]]$id, outside_forecaster_subset)),
+#   scores = tar_read(joined_scores) %>% filter(forecaster %in% c(forecaster_parameter_combinations[[i]]$id, outside_forecaster_subset)),
+#   truth_data = tar_read(hhs_evaluation_data),
+#   disease = "covid"
 # )
 
 if (params$disease == "flu") {
   base_forecaster_name <- "FluSight-baseline"
   ensemble_forecaster_name <- "FluSight-ensemble"
 } else {
   base_forecaster_name <- "COVIDhub-baseline"
-  ensemble_forecaster_name <- "COVIDhub-ensemble"
+  ensemble_forecaster_name <- "COVIDhub_CDC-ensemble"
 }
 
 # Load scores and filter them, get global variables
@@ -96,17 +96,23 @@ The table is sorted by ascending WIS and contains all the forecasters in this no
 if (params$disease == "flu") {
   ignore_keys <- c("forecaster", "keys_to_ignore", "pop_scaling")
 } else {
-  ignore_keys <- c("forecaster", "keys_to_ignore")
+  ignore_keys <- c("forecaster", "keys_to_ignore", "outcome")
 }
 
 param_table <- params$forecaster_parameters %>%
   select(-any_of(ignore_keys)) %>%
   {
     if ("n_training" %in% colnames(.)) {
-      (.) %>% mutate(n_training = as.character(n_training))
+      . <- (.) %>% mutate(n_training = as.character(n_training))
     } else {
       .
     }
+    if ("trainer" %in% colnames(.)) {
+      . <- (.) %>% mutate(trainer = trainer[[1]])
+    } else {
+      .
+    }
+    .
   } %>%
   full_join(
     scores %>%

scripts/targets-common.R

@@ -4,10 +4,10 @@ suppressPackageStartupMessages({
 
 # On tanka, we have 64 cores, but we leave some free to try to reduce thrashing
 # and to allow for other users.
-if (parallel::detectCores() < 30) {
-  num_workers <- parallel::detectCores() - 1L
+if (parallel::detectCores() == 64) {
+  num_workers <- 30L
 } else {
-  num_workers <- parallel::detectCores() - 20L
+  num_workers <- parallel::detectCores() - 1L
 }
 
 main_controller <- crew_controller_local(