Merge branch 'ndefries/backfill/speed2' into ndefries/backfill/speed-join-v-merge-order-matters

Author: nmdefries

backfill_corrections/delphiBackfillCorrection/DESCRIPTION

@@ -23,9 +23,11 @@ Imports:
     jsonlite,
     lubridate,
     tidyr,
-    zoo,
     utils,
-    parallel
+    parallel,
+    purrr,
+    vctrs,
+    RcppRoll
 Suggests:
     knitr (>= 1.15),
     rmarkdown (>= 1.4),

backfill_corrections/delphiBackfillCorrection/NAMESPACE

@@ -21,6 +21,7 @@ export(read_data)
 export(read_params)
 export(run_backfill)
 import(covidcast)
+importFrom(RcppRoll,roll_mean)
 importFrom(arrow,read_parquet)
 importFrom(dplyr,"%>%")
 importFrom(dplyr,across)
@@ -48,6 +49,7 @@ importFrom(lubridate,make_date)
 importFrom(lubridate,month)
 importFrom(lubridate,year)
 importFrom(parallel,detectCores)
+importFrom(purrr,map_dfc)
 importFrom(quantgen,quantile_lasso)
 importFrom(readr,write_csv)
 importFrom(stats,coef)
@@ -60,8 +62,7 @@ importFrom(stringr,str_split)
 importFrom(tibble,tribble)
 importFrom(tidyr,crossing)
 importFrom(tidyr,drop_na)
-importFrom(tidyr,fill)
 importFrom(tidyr,pivot_longer)
 importFrom(tidyr,pivot_wider)
 importFrom(utils,head)
-importFrom(zoo,rollmeanr)
+importFrom(vctrs,vec_fill_missing)

backfill_corrections/delphiBackfillCorrection/R/constants.R

@@ -19,10 +19,12 @@ SQRTSCALE <-c('sqrty0', 'sqrty1', "sqrty2")
 LOG_LAG <-"inv_log_lag"
 
 # Dates
+DATE_FORMAT <- "%Y-%m-%d"
 WEEKDAYS_ABBR <- c("Mon", "Tue", "Wed", "Thurs", "Fri", "Sat") # wd
 WEEK_ISSUES <- c("W1_issue", "W2_issue", "W3_issue") # wm
 TODAY <- Sys.Date()
 
+# Signals we want to make predictions for
 INDICATORS_AND_SIGNALS <- tibble::tribble(
   ~indicator, ~signal, ~name_suffix, ~sub_dir,
   "changehc", "covid", "", "chng",

backfill_corrections/delphiBackfillCorrection/R/io.R

@@ -244,7 +244,7 @@ get_training_date_range <- function(params) {
   if (params$train_models) {
     if (params_element_exists_and_valid(params, "training_end_date")) {
       # Use user-provided end date.
-      training_end_date <- as.Date(params$training_end_date)
+      training_end_date <- as.Date(params$training_end_date, DATE_FORMAT)
     } else {
       # Default end date is today.
       training_end_date <- default_end_date

backfill_corrections/delphiBackfillCorrection/R/main.R

@@ -121,7 +121,7 @@ run_backfill <- function(df, params,
             ) %>%
             drop_na()
           geo_test_data <- filter(combined_df,
-              issue_date %in% params$test_dates
+              issue_date %in% as.character(params$test_dates)
             ) %>%
             drop_na()
 
@@ -300,7 +300,7 @@ main <- function(params,
       input_group$indicator, input_group$signal, params, input_group$sub_dir
     )
     if (length(files_list) == 0) {
-      warning("No files found for indicator indicator ", input_group$indicator,
+      warning("No files found for indicator ", input_group$indicator,
               " signal ", input_group$signal, ", skipping")
       next
     }
@@ -331,9 +331,6 @@ main <- function(params,
       refd_col = refd_col, lag_col = lag_col, issued_col = issued_col
     )
 
-    input_data[[refd_col]] <- as.Date(input_data[[refd_col]], "%Y-%m-%d")
-    input_data[[issued_col]] <- as.Date(input_data[[issued_col]], "%Y-%m-%d")
-    
     # Check available training days
     training_days_check(input_data[[issued_col]], params$training_days)
 

backfill_corrections/delphiBackfillCorrection/R/model.R

@@ -183,11 +183,8 @@ evaluate <- function(test_data, taus) {
 exponentiate_preds <- function(test_data, taus) {
   pred_cols = paste0("predicted_tau", taus)
 
-  # Drop original predictions and join on exponentiated versions
-  test_data = bind_cols(
-    select(test_data, -starts_with("predicted")),
-    exp(test_data[, pred_cols])
-  )
+  # Replace original predictions with exponentiated versions
+  test_data[, pred_cols] <- exp(test_data[, pred_cols])
 
   return(test_data)
 }
@@ -211,17 +208,22 @@ get_model <- function(model_path, train_data, covariates, tau,
         " does not exist; training new model")
     }
     # Quantile regression
+    ## TODO: how does the speed compare using GLPK? Apparently it's faster on smaller
+    # models.
     obj <- quantile_lasso(as.matrix(train_data[covariates]),
                          train_data$log_value_target, tau = tau,
                          lambda = lambda, standardize = FALSE, lp_solver = lp_solver)
 
     # Save model to cache.
     create_dir_not_exist(dirname(model_path))
+    ## TODO: save() is fairly slow. Since we're not sharing the model files, can we
+    # use saveRDS() instead?
     save(obj, file=model_path)
   } else {
     # Load model from cache invisibly. Object has the same name as the original
     # model object, `obj`.
     msg_ts("Loading from ", model_path)
+    ## TODO: readRDS()
     load(model_path)
   }
 

backfill_corrections/delphiBackfillCorrection/R/preprocessing.R

@@ -25,8 +25,11 @@
 fill_rows <- function(df, refd_col, lag_col, min_refd, max_refd, ref_lag) {
   # Full list of lags
   # +30 to have values for calculating 7-day averages
-  lags <- min(df[[lag_col]]): (ref_lag + 30) 
-  refds <- seq(min_refd, max_refd, by="day") # Full list reference date
+  lags <- min(df[[lag_col]]): (ref_lag + 30)
+  # Full list reference dates
+  refds <- as.character(
+    seq(as.Date(min_refd, DATE_FORMAT), as.Date(max_refd, DATE_FORMAT), by="day")
+  )
   row_inds_df <- setNames(
     as.data.frame(crossing(refds, lags)),
     c(refd_col, lag_col)
@@ -45,49 +48,59 @@ fill_rows <- function(df, refd_col, lag_col, min_refd, max_refd, ref_lag) {
 #' @template refd_col-template
 #' @template lag_col-template
 #' 
-#' @importFrom tidyr fill pivot_wider pivot_longer
-#' @importFrom dplyr %>% everything select
+#' @importFrom tidyr pivot_wider pivot_longer
+#' @importFrom purrr map_dfc
+#' @importFrom vctrs vec_fill_missing
 #' 
 #' @export
 fill_missing_updates <- function(df, value_col, refd_col, lag_col) {
-  pivot_df <- df[order(df[[lag_col]], decreasing=FALSE), ] %>%
-    pivot_wider(id_cols=lag_col, names_from=refd_col, values_from=value_col)
+  pivot_df <- pivot_wider(
+    df[order(df[[lag_col]], decreasing=FALSE), ],
+    id_cols=lag_col, names_from=refd_col, values_from=value_col
+  )
 
   if (any(diff(pivot_df[[lag_col]]) != 1)) {
     stop("Risk exists in forward filling")
   }
-  pivot_df <- fill(pivot_df, everything(), .direction="down")
+
+  pivot_df <- map_dfc(pivot_df, function(col) {
+    vec_fill_missing(col, direction="down")
+  })
 
   # Fill NAs with 0s
   pivot_df[is.na(pivot_df)] <- 0
 
   backfill_df <- pivot_longer(pivot_df,
     -lag_col, values_to="value_raw", names_to=refd_col
   )
-  backfill_df[[refd_col]] = as.Date(backfill_df[[refd_col]])
 
   return (as.data.frame(backfill_df))
 }
 
 #' Calculate 7 day moving average for each issue date
+#'
 #' The 7dav for date D reported on issue date D_i is the average from D-7 to D-1
+#'
 #' @param pivot_df Data Frame where the columns are issue dates and the rows are 
 #'    reference dates
 #' @template refd_col-template
 #' 
-#' @importFrom zoo rollmeanr
-#' 
+#' @importFrom RcppRoll roll_mean
+#'
 #' @export
 get_7dav <- function(pivot_df, refd_col) {
-  for (col in colnames(pivot_df)) {
-    if (col == refd_col) next
-    pivot_df[, col] <- rollmeanr(pivot_df[, col], 7, align="right", fill=NA)
-  }
+  pivot_df <- cbind(
+    # Keep time values at the front
+    pivot_df[, refd_col],
+    # Compute moving average of all non-refd columns
+    RcppRoll::roll_mean(
+      as.matrix(pivot_df[, names(pivot_df)[names(pivot_df) != refd_col]]),
+      7L, align = "right", fill = NA
+    )
+  )
   backfill_df <- pivot_longer(pivot_df,
     -refd_col, values_to="value_raw", names_to="issue_date"
   )
-  backfill_df[[refd_col]] = as.Date(backfill_df[[refd_col]])
-  backfill_df[["issue_date"]] = as.Date(backfill_df[["issue_date"]])
   return (as.data.frame(backfill_df))
 }
 
@@ -99,7 +112,7 @@ get_7dav <- function(pivot_df, refd_col) {
 #' 
 #' @export
 add_shift <- function(df, n_day, refd_col) {
-  df[, refd_col] <- as.Date(df[, refd_col]) + n_day
+  df[[refd_col]] <- as.character(as.Date(df[[refd_col]], DATE_FORMAT) + n_day)
   return (df)
 }
 
@@ -113,7 +126,7 @@ add_shift <- function(df, n_day, refd_col) {
 #' 
 #' @export
 add_dayofweek <- function(df, time_col, suffix, wd = WEEKDAYS_ABBR) {
-  dayofweek <- as.numeric(format(df[[time_col]], format="%u"))
+  dayofweek <- as.numeric(format(as.Date(df[[time_col]], DATE_FORMAT), format="%u"))
   for (i in seq_along(wd)) {
     df[, paste0(wd[i], suffix)] <- as.numeric(dayofweek == i)
   }
@@ -159,7 +172,7 @@ get_weekofmonth <- function(date) {
 #' 
 #' @export
 add_weekofmonth <- function(df, time_col, wm = WEEK_ISSUES) {
-  weekofmonth <- get_weekofmonth(df[[time_col]])
+  weekofmonth <- get_weekofmonth(as.Date(df[[time_col]], DATE_FORMAT))
   for (i in seq_along(wm)) {
     df[, paste0(wm[i])] <- as.numeric(weekofmonth == i)
   }
@@ -174,16 +187,17 @@ add_weekofmonth <- function(df, time_col, wm = WEEK_ISSUES) {
 #' @template lag_col-template
 #' @template ref_lag-template
 #' 
-#' @importFrom dplyr %>% full_join left_join
+#' @importFrom dplyr full_join left_join
 #' @importFrom tidyr pivot_wider drop_na
 #' 
 #' @export
 add_7davs_and_target <- function(df, value_col, refd_col, lag_col, ref_lag) {
-  df$issue_date <- df[[refd_col]] + df[[lag_col]]
-  pivot_df <- df[order(df$issue_date, decreasing=FALSE), ] %>%
-    pivot_wider(id_cols=refd_col, names_from="issue_date", 
-                values_from=value_col)
-  
+  df$issue_date <- as.character(as.Date(df[[refd_col]], DATE_FORMAT) + df[[lag_col]])
+  pivot_df <- pivot_wider(
+    df[order(df$issue_date, decreasing=FALSE), ],
+    id_cols=refd_col, names_from="issue_date", values_from=value_col
+  )
+
   # Add 7dav avg
   avg_df <- get_7dav(pivot_df, refd_col)
   avg_df <- add_shift(avg_df, 1, refd_col) # 7dav until yesterday

backfill_corrections/delphiBackfillCorrection/R/utils.R

@@ -91,13 +91,13 @@ read_params <- function(path = "params.json", template_path = "params.json.templ
       stop("`test_dates` setting in params must be a length-2 list of dates")
     }
     params$test_dates <- seq(
-      as.Date(params$test_dates[1]),
-      as.Date(params$test_dates[2]),
+      as.Date(params$test_dates[1], DATE_FORMAT),
+      as.Date(params$test_dates[2], DATE_FORMAT),
       by="days"
     )
   }
   if (params_element_exists_and_valid(params, "training_end_date")) {
-    if (as.Date(params$training_end_date) > TODAY) {
+    if (as.Date(params$training_end_date, DATE_FORMAT) > TODAY) {
       stop("training_end_date can't be in the future")
     }
   }
@@ -189,7 +189,9 @@ validity_checks <- function(df, value_types, num_col, denom_col, signal_suffixes
 #' @param issue_date contents of input data's `issue_date` column
 #' @template training_days-template
 training_days_check <- function(issue_date, training_days) {
-  valid_training_days = as.integer(max(issue_date) - min(issue_date)) + 1
+  valid_training_days = as.integer(
+    as.Date(max(issue_date), DATE_FORMAT) - as.Date(min(issue_date), DATE_FORMAT)
+  ) + 1
   if (training_days > valid_training_days) {
     warning(sprintf("Only %d days are available at most for training.", valid_training_days))
   }