working example

Author: dajmcdon

.Rbuildignore

@@ -8,3 +8,4 @@
 ^_pkgdown\.yml$
 ^docs$
 ^pkgdown$
+^musings$

.gitignore

@@ -4,3 +4,4 @@
 .Ruserdata
 docs
 inst/doc
+.DS_Store

R/check-train_window.R

@@ -0,0 +1,118 @@
+#' Check Training Window Length
+#'
+#' `check_train_window` creates a *specification* of a recipe
+#'  check that will check if there is insufficient training data
+#'
+#' @inheritParams check_missing
+#' @min_train_window Positive integer. The minimum amount of training
+#'   data time points required
+#'   to fit a predictive model. Using less results causes downstream
+#'   fit calls to return minimal objects rather than crashing.
+#' @param warn If `TRUE` the check will throw a warning instead
+#'   of an error when failing.
+#' @param train_window The number of days of training data.
+#'   This is `NULL` until computed by [prep()].
+#' @template check-return
+#' @family checks
+#' @export
+#'
+check_train_window <-
+  function(recipe,
+           ...,
+           role = NA,
+           skip = FALSE,
+           trained = FALSE,
+           min_train_window = 20,
+           warn = TRUE,
+           train_length,
+           id = rand_id("train_window_check_")) {
+    add_check(
+      recipe,
+      check_train_window_new(
+        terms = dplyr::enquos(...),
+        role = role,
+        skip = skip,
+        trained = trained,
+        min_train_window = min_train_window,
+        warn = warn,
+        train_length = train_length,
+        id = id
+      )
+    )
+  }
+
+## Initializes a new object
+check_train_window_new <-
+  function(terms, role, skip, trained, min_train_window, warn,
+           train_length, id) {
+    check(
+      subclass = "train_window",
+      terms = terms,
+      role = role,
+      skip = skip,
+      trained = trained,
+      min_train_window = min_train_window,
+      warn = warn,
+      train_length = train_length,
+      id = id
+    )
+  }
+
+
+prep.check_train_window <- function(x,
+                                    training,
+                                    info = NULL,
+                                    ...) {
+
+  train_length <- nrow(training)
+
+
+  check_train_window_new(
+    terms      = x$terms,
+    role       = x$role,
+    trained    = TRUE,
+    skip       = x$skip,
+    warn       = x$warn,
+    min_train_window = min_train_window,
+    warn = warn,
+    train_length = train_length,
+    id         = x$id
+  )
+}
+
+bake.check_range <- function(object,
+                             new_data,
+                             ...) {
+
+  mtw <- object$min_train_window
+  stopifnot(is.numeric(mtw), length(mtw) == 1L, mtw == as.integer(mtw))
+
+  n <- nrow(new_data)
+  n.complete <- sum(complete.cases(new_data))
+
+  msg <- NULL
+  if (n < mtw) {
+    msg <- paste0(msg, "Total available rows of data is ", n,
+                  "\n < min_train_window ", mtw, ".\n")
+  }
+  if (n.complete < mtw) {
+    msg <- paste0(msg, "Total complete rows of data is ", n.complete,
+                  "\n < min_train_window ", mtw, ".\n")
+  }
+
+  if (object$warn & !is.null(msg)) {
+    rlang::warn(msg)
+  } else if (!is.null(msg)) {
+    rlang::abort(msg)
+  }
+
+  as_tibble(new_data)
+}
+
+print.check_train_window <-
+  function(x, width = max(20, options()$width - 30), ...) {
+    title <- "Checking number of training observations"
+    invisible(x)
+  }
+
+

R/epi_ahead.R

@@ -121,8 +121,11 @@ bake.step_epi_ahead <- function(object, new_data, ...) {
   }
 
   grid <- tidyr::expand_grid(
-    col = object$columns, lag_val = -object$ahead, ahead_val = object$ahead) %>%
-    dplyr::mutate(newname = glue::glue("{object$prefix}{ahead_val}_{col}")) %>%
+    col = object$columns, lag_val = -object$ahead) %>%
+    dplyr::mutate(
+      ahead_val = -lag_val,
+      newname = glue::glue("{object$prefix}{ahead_val}_{col}")
+    ) %>%
     dplyr::select(-ahead_val)
 
   ## ensure no name clashes
@@ -143,7 +146,7 @@ bake.step_epi_ahead <- function(object, new_data, ...) {
     by = ok
   )
 
-  dplyr::full_join(new_data, lagged, by = object$keys) %>%
+  dplyr::full_join(new_data, lagged, by = ok) %>%
     dplyr::group_by(dplyr::across(dplyr::all_of(ok[-1]))) %>%
     dplyr::arrange(time_value) %>%
     dplyr::ungroup()

R/epi_lag.R

@@ -106,7 +106,7 @@ bake.step_epi_lag <- function(object, new_data, ...) {
     by = ok
   )
 
-  dplyr::full_join(new_data, lagged, by = object$keys) %>%
+  dplyr::full_join(new_data, lagged, by = ok) %>%
     dplyr::group_by(dplyr::across(dplyr::all_of(ok[-1]))) %>%
     dplyr::arrange(time_value) %>%
     dplyr::ungroup()

R/epi_recipe.R

@@ -107,7 +107,7 @@ epi_recipe.epi_df <-
       var_info$role <- roles
     } else {
       var_info <- var_info %>% dplyr::filter(!(variable %in% keys))
-      var_info$role <- NA
+      var_info$role <- "raw"
     }
     ## Now we add the keys when necessary
     var_info <- dplyr::union(

musings/example-recipe.R

@@ -0,0 +1,65 @@
+library(tidyverse)
+library(covidcast)
+library(delphi.epidata)
+library(epiprocess)
+library(tidymodels)
+x <- covidcast(
+  data_source = "jhu-csse",
+  signals = "confirmed_7dav_incidence_prop",
+  time_type = "day",
+  geo_type = "state",
+  time_values = epirange(20200301, 20211231),
+  geo_values = "*"
+) %>%
+  fetch_tbl() %>%
+  select(geo_value, time_value, case_rate = value)
+
+y <- covidcast(
+  data_source = "jhu-csse",
+  signals = "deaths_7dav_incidence_prop",
+  time_type = "day",
+  geo_type = "state",
+  time_values = epirange(20200301, 20211231),
+  geo_values = "*"
+) %>%
+  fetch_tbl() %>%
+  select(geo_value, time_value, death_rate = value)
+
+x <- x %>%
+  full_join(y, by = c("geo_value", "time_value")) %>%
+  as_epi_df()
+rm(y)
+
+xx <- x %>% filter(time_value > "2021-12-01")
+
+
+# Baseline AR3
+r <- epi_recipe(x) %>% # if we add this as a class, maybe we get better
+                       # behaviour downstream?
+  step_epi_ahead(death_rate, ahead = 7) %>%
+  step_epi_lag(death_rate, lag = c(0, 7, 14)) %>%
+  step_epi_lag(case_rate, lag = c(0, 7, 14)) %>%
+  step_naomit(all_predictors()) %>%
+  # below, `skip` means we don't do this at predict time
+  # we should probably do something useful here to avoid user error
+  step_naomit(all_outcomes(), skip = TRUE)
+
+
+slm <- linear_reg()
+
+slm_fit <- workflow() %>%
+  add_recipe(r) %>%
+  add_model(slm) %>%
+  fit(data = x)
+
+x_latest <- x %>%
+  filter(!is.na(case_rate), !is.na(death_rate)) %>%
+  group_by(geo_value) %>%
+  slice_tail(n = 15) # have lag 0,...,14, so need 15 for a complete case
+
+pp <- predict(slm_fit, new_data = x_latest) # drops the keys...
+
+xl <- x %>%
+  filter(!is.na(case_rate), !is.na(death_rate)) %>%
+  group_by(geo_value) %>%
+  slice_tail(n = 14)