From 0e7d90f072a6770c209f53d9bcf5ee5b2aef4720 Mon Sep 17 00:00:00 2001 From: Dmitry Shemetov Date: Wed, 25 Sep 2024 19:01:42 -0700 Subject: [PATCH 1/3] fix: update for compatibility with epiprocess==0.9.0 Co-authored-by: Daniel McDonald --- .Rbuildignore | 3 +- DESCRIPTION | 5 +- R/autoplot.R | 2 +- R/cdc_baseline_forecaster.R | 6 +- R/epi_recipe.R | 12 +- R/epi_workflow.R | 3 +- R/flusight_hub_formatter.R | 4 +- R/key_colnames.R | 15 ++- R/step_epi_slide.R | 121 +++++++++--------- R/utils-misc.R | 16 +-- data-raw/grad_employ_subset.R | 2 +- data/grad_employ_subset.rda | Bin 8491 -> 9546 bytes man/autoplot-epipred.Rd | 2 + man/cdc_baseline_forecaster.Rd | 6 +- man/epi_slide_wrapper.Rd | 4 +- man/flusight_hub_formatter.Rd | 4 +- man/step_epi_slide.Rd | 26 ++-- tests/testthat/test-epi_recipe.R | 2 +- tests/testthat/test-key_colnames.R | 14 +- tests/testthat/test-layer_add_forecast_date.R | 2 + tests/testthat/test-layer_add_target_date.R | 2 + tests/testthat/test-pad_to_end.R | 2 +- tests/testthat/test-step_epi_slide.R | 50 ++++---- tests/testthat/test-step_training_window.R | 2 +- vignettes/articles/sliding.Rmd | 55 ++++---- vignettes/arx-classifier.Rmd | 20 ++- vignettes/epipredict.Rmd | 9 +- 27 files changed, 198 insertions(+), 191 deletions(-) diff --git a/.Rbuildignore b/.Rbuildignore index f1a8c3636..510725267 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -19,4 +19,5 @@ ^DEVELOPMENT\.md$ ^doc$ ^Meta$ -^.lintr$ \ No newline at end of file +^.lintr$ +^.venv$ \ No newline at end of file diff --git a/DESCRIPTION b/DESCRIPTION index a7df5ba5a..37134c160 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: epipredict Title: Basic epidemiology forecasting methods -Version: 0.0.20 +Version: 0.0.21 Authors@R: c( person("Daniel", "McDonald", , "daniel@stat.ubc.ca", role = c("aut", "cre")), person("Ryan", "Tibshirani", , "ryantibs@cmu.edu", role = "aut"), @@ -23,8 +23,7 @@ URL: https://github.com/cmu-delphi/epipredict/, https://cmu-delphi.github.io/epipredict BugReports: https://github.com/cmu-delphi/epipredict/issues/ Depends: - epiprocess (>= 0.8.0), - epiprocess (< 0.9.0), + epiprocess (>= 0.9.0), parsnip (>= 1.0.0), R (>= 3.5.0) Imports: diff --git a/R/autoplot.R b/R/autoplot.R index d35850fd6..648c74e33 100644 --- a/R/autoplot.R +++ b/R/autoplot.R @@ -131,7 +131,7 @@ autoplot.epi_workflow <- function( if (length(extra_keys) == 0L) extra_keys <- NULL edf <- as_epi_df(edf, as_of = object$fit$meta$as_of, - additional_metadata = list(other_keys = extra_keys) + other_keys = extra_keys %||% character() ) if (is.null(predictions)) { return(autoplot( diff --git a/R/cdc_baseline_forecaster.R b/R/cdc_baseline_forecaster.R index 74af5e443..b2e7434e2 100644 --- a/R/cdc_baseline_forecaster.R +++ b/R/cdc_baseline_forecaster.R @@ -29,11 +29,11 @@ #' mutate(deaths = pmax(death_rate / 1e5 * pop * 7, 0)) %>% #' select(-pop, -death_rate) %>% #' group_by(geo_value) %>% -#' epi_slide(~ sum(.$deaths), before = 6, new_col_name = "deaths") %>% +#' epi_slide(~ sum(.$deaths), .window_size = 7, .new_col_name = "deaths_7dsum") %>% #' ungroup() %>% #' filter(weekdays(time_value) == "Saturday") #' -#' cdc <- cdc_baseline_forecaster(weekly_deaths, "deaths") +#' cdc <- cdc_baseline_forecaster(weekly_deaths, "deaths_7dsum") #' preds <- pivot_quantiles_wider(cdc$predictions, .pred_distn) #' #' if (require(ggplot2)) { @@ -47,7 +47,7 @@ #' geom_line(aes(y = .pred), color = "orange") + #' geom_line( #' data = weekly_deaths %>% filter(geo_value %in% four_states), -#' aes(x = time_value, y = deaths) +#' aes(x = time_value, y = deaths_7dsum) #' ) + #' scale_x_date(limits = c(forecast_date - 90, forecast_date + 30)) + #' labs(x = "Date", y = "Weekly deaths") + diff --git a/R/epi_recipe.R b/R/epi_recipe.R index 88ba605cd..0ff9efee8 100644 --- a/R/epi_recipe.R +++ b/R/epi_recipe.R @@ -95,7 +95,7 @@ epi_recipe.epi_df <- keys <- key_colnames(x) # we know x is an epi_df var_info <- tibble(variable = vars) - key_roles <- c("geo_value", "time_value", rep("key", length(keys) - 2)) + key_roles <- c("geo_value", rep("key", length(keys) - 2), "time_value") ## Check and add roles when available if (!is.null(roles)) { @@ -499,8 +499,11 @@ prep.epi_recipe <- function( if (!is_epi_df(training)) { # tidymodels killed our class # for now, we only allow step_epi_* to alter the metadata - training <- dplyr::dplyr_reconstruct( - as_epi_df(training), before_template + metadata <- attr(before_template, "metadata") + training <- as_epi_df( + training, + as_of = metadata$as_of, + other_keys = metadata$other_keys %||% character() ) } training <- dplyr::relocate(training, all_of(key_colnames(training))) @@ -579,8 +582,7 @@ bake.epi_recipe <- function(object, new_data, ..., composition = "epi_df") { new_data <- as_epi_df( new_data, as_of = meta$as_of, - # avoid NULL if meta is from saved older epi_df: - additional_metadata = meta$additional_metadata %||% list() + other_keys = meta$other_keys %||% character() ) } new_data diff --git a/R/epi_workflow.R b/R/epi_workflow.R index b059a81d0..f448f4aff 100644 --- a/R/epi_workflow.R +++ b/R/epi_workflow.R @@ -98,7 +98,8 @@ is_epi_workflow <- function(x) { fit.epi_workflow <- function(object, data, ..., control = workflows::control_workflow()) { object$fit$meta <- list( max_time_value = max(data$time_value), - as_of = attributes(data)$metadata$as_of + as_of = attr(data, "metadata")$as_of, + other_keys = attr(data, "metadata")$other_keys ) object$original_data <- data diff --git a/R/flusight_hub_formatter.R b/R/flusight_hub_formatter.R index c91f738ae..3e0eb1aaa 100644 --- a/R/flusight_hub_formatter.R +++ b/R/flusight_hub_formatter.R @@ -67,11 +67,11 @@ abbr_to_location <- function(abbr) { #' mutate(deaths = pmax(death_rate / 1e5 * pop * 7, 0)) %>% #' select(-pop, -death_rate) %>% #' group_by(geo_value) %>% -#' epi_slide(~ sum(.$deaths), before = 6, new_col_name = "deaths") %>% +#' epi_slide(~ sum(.$deaths), .window_size = 7, .new_col_name = "deaths_7dsum") %>% #' ungroup() %>% #' filter(weekdays(time_value) == "Saturday") #' -#' cdc <- cdc_baseline_forecaster(weekly_deaths, "deaths") +#' cdc <- cdc_baseline_forecaster(weekly_deaths, "deaths_7dsum") #' flusight_hub_formatter(cdc) #' flusight_hub_formatter(cdc, target = "wk inc covid deaths") #' flusight_hub_formatter(cdc, target = paste(horizon, "wk inc covid deaths")) diff --git a/R/key_colnames.R b/R/key_colnames.R index c69d1a628..b9ebde5dc 100644 --- a/R/key_colnames.R +++ b/R/key_colnames.R @@ -1,19 +1,20 @@ #' @export key_colnames.recipe <- function(x, ...) { - possible_keys <- c("geo_value", "time_value", "key") - keys <- x$var_info$variable[x$var_info$role %in% possible_keys] - keys[order(match(keys, possible_keys))] %||% character(0L) + geo_key <- x$var_info$variable[x$var_info$role %in% "geo_value"] + time_key <- x$var_info$variable[x$var_info$role %in% "time_value"] + keys <- x$var_info$variable[x$var_info$role %in% "key"] + c(geo_key, keys, time_key) %||% character(0L) } #' @export key_colnames.epi_workflow <- function(x, ...) { # safer to look at the mold than the preprocessor mold <- hardhat::extract_mold(x) - possible_keys <- c("geo_value", "time_value", "key") molded_names <- names(mold$extras$roles) - keys <- map(mold$extras$roles[molded_names %in% possible_keys], names) - keys <- unname(unlist(keys)) - keys[order(match(keys, possible_keys))] %||% character(0L) + geo_key <- names(mold$extras$roles[molded_names %in% "geo_value"]$geo_value) + time_key <- names(mold$extras$roles[molded_names %in% "time_value"]$time_value) + keys <- names(mold$extras$roles[molded_names %in% "key"]$key) + c(geo_key, keys, time_key) %||% character(0L) } kill_time_value <- function(v) { diff --git a/R/step_epi_slide.R b/R/step_epi_slide.R index 9714971fa..c7d3f9fbd 100644 --- a/R/step_epi_slide.R +++ b/R/step_epi_slide.R @@ -19,13 +19,18 @@ #' argument must be named `.x`. A common, though very difficult to debug #' error is using something like `function(x) mean`. This will not work #' because it returns the function mean, rather than `mean(x)` -#' @param before,after the size of the sliding window on the left and the right -#' of the center. Usually non-negative integers for data indexed by date, but -#' more restrictive in other cases (see [epiprocess::epi_slide()] for details). -#' @param f_name a character string of at most 20 characters that describes -#' the function. This will be combined with `prefix` and the columns in `...` -#' to name the result using `{prefix}{f_name}_{column}`. By default it will be determined -#' automatically using `clean_f_name()`. +#' @param .window_size the size of the sliding window, required. Usually a +#' non-negative integer will suffice (e.g. for data indexed by date, but more +#' restrictive in other time_type cases (see [epiprocess::epi_slide()] for +#' details). For example, set to 7 for a 7-day window. +#' @param .align a character string indicating how the window should be aligned. +#' By default, this is "right", meaning the slide_window will be anchored with +#' its right end point on the reference date. (see [epiprocess::epi_slide()] +#' for details). +#' @param f_name a character string of at most 20 characters that describes the +#' function. This will be combined with `prefix` and the columns in `...` to +#' name the result using `{prefix}{f_name}_{column}`. By default it will be +#' determined automatically using `clean_f_name()`. #' #' @template step-return #' @@ -37,53 +42,55 @@ #' rec <- epi_recipe(jhu) %>% #' step_epi_slide(case_rate, death_rate, #' .f = \(x) mean(x, na.rm = TRUE), -#' before = 6L +#' .window_size = 7L #' ) #' bake(prep(rec, jhu), new_data = NULL) -step_epi_slide <- - function(recipe, - ..., - .f, - before = 0L, - after = 0L, - role = "predictor", - prefix = "epi_slide_", - f_name = clean_f_name(.f), - skip = FALSE, - id = rand_id("epi_slide")) { - if (!is_epi_recipe(recipe)) { - cli_abort("This recipe step can only operate on an {.cls epi_recipe}.") - } - .f <- validate_slide_fun(.f) - epiprocess:::validate_slide_window_arg(before, attributes(recipe$template)$metadata$time_type) - epiprocess:::validate_slide_window_arg(after, attributes(recipe$template)$metadata$time_type) - arg_is_chr_scalar(role, prefix, id) - arg_is_lgl_scalar(skip) +step_epi_slide <- function(recipe, + ..., + .f, + .window_size = NULL, + .align = c("right", "center", "left"), + role = "predictor", + prefix = "epi_slide_", + f_name = clean_f_name(.f), + skip = FALSE, + id = rand_id("epi_slide")) { + if (!is_epi_recipe(recipe)) { + cli_abort("This recipe step can only operate on an {.cls epi_recipe}.") + } + .f <- validate_slide_fun(.f) + if (is.null(.window_size)) { + cli_abort("step_epi_slide: `.window_size` must be specified.") + } + epiprocess:::validate_slide_window_arg(.window_size, attributes(recipe$template)$metadata$time_type) + .align <- rlang::arg_match(.align) + arg_is_chr_scalar(role, prefix, id) + arg_is_lgl_scalar(skip) - recipes::add_step( - recipe, - step_epi_slide_new( - terms = enquos(...), - before = before, - after = after, - .f = .f, - f_name = f_name, - role = role, - trained = FALSE, - prefix = prefix, - keys = key_colnames(recipe), - columns = NULL, - skip = skip, - id = id - ) + recipes::add_step( + recipe, + step_epi_slide_new( + terms = enquos(...), + .window_size = .window_size, + .align = .align, + .f = .f, + f_name = f_name, + role = role, + trained = FALSE, + prefix = prefix, + keys = key_colnames(recipe), + columns = NULL, + skip = skip, + id = id ) - } + ) +} step_epi_slide_new <- function(terms, - before, - after, + .window_size, + .align, .f, f_name, role, @@ -96,8 +103,8 @@ step_epi_slide_new <- recipes::step( subclass = "epi_slide", terms = terms, - before = before, - after = after, + .window_size = .window_size, + .align = .align, .f = .f, f_name = f_name, role = role, @@ -119,8 +126,8 @@ prep.step_epi_slide <- function(x, training, info = NULL, ...) { step_epi_slide_new( terms = x$terms, - before = x$before, - after = x$after, + .window_size = x$.window_size, + .align = x$.align, .f = x$.f, f_name = x$f_name, role = x$role, @@ -165,8 +172,8 @@ bake.step_epi_slide <- function(object, new_data, ...) { # } epi_slide_wrapper( new_data, - object$before, - object$after, + object$.window_size, + object$.align, object$columns, c(object$.f), object$f_name, @@ -190,7 +197,7 @@ bake.step_epi_slide <- function(object, new_data, ...) { #' @importFrom dplyr bind_cols group_by ungroup #' @importFrom epiprocess epi_slide #' @keywords internal -epi_slide_wrapper <- function(new_data, before, after, columns, fns, fn_names, group_keys, name_prefix) { +epi_slide_wrapper <- function(new_data, .window_size, .align, columns, fns, fn_names, group_keys, name_prefix) { cols_fns <- tidyr::crossing(col_name = columns, fn_name = fn_names, fn = fns) # Iterate over the rows of cols_fns. For each row number, we will output a # transformed column. The first result returns all the original columns along @@ -204,10 +211,10 @@ epi_slide_wrapper <- function(new_data, before, after, columns, fns, fn_names, g result <- new_data %>% group_by(across(all_of(group_keys))) %>% epi_slide( - before = before, - after = after, - new_col_name = result_name, - f = function(slice, geo_key, ref_time_value) { + .window_size = .window_size, + .align = .align, + .new_col_name = result_name, + .f = function(slice, geo_key, ref_time_value) { fn(slice[[col_name]]) } ) %>% diff --git a/R/utils-misc.R b/R/utils-misc.R index af064b37c..b4d1c28b7 100644 --- a/R/utils-misc.R +++ b/R/utils-misc.R @@ -33,15 +33,14 @@ check_pname <- function(res, preds, object, newname = NULL) { grab_forged_keys <- function(forged, workflow, new_data) { - keys <- c("geo_value", "time_value", "key") forged_roles <- names(forged$extras$roles) - extras <- dplyr::bind_cols(forged$extras$roles[forged_roles %in% keys]) + extras <- dplyr::bind_cols(forged$extras$roles[forged_roles %in% c("geo_value", "time_value", "key")]) # 1. these are the keys in the test data after prep/bake new_keys <- names(extras) # 2. these are the keys in the training data old_keys <- key_colnames(workflow) # 3. these are the keys in the test data as input - new_df_keys <- key_colnames(new_data, extra_keys = setdiff(new_keys, keys[1:2])) + new_df_keys <- key_colnames(new_data, extra_keys = setdiff(new_keys, c("geo_value", "time_value"))) if (!(setequal(old_keys, new_df_keys) && setequal(new_keys, new_df_keys))) { cli::cli_warn(c( "Not all epi keys that were present in the training data are available", @@ -49,12 +48,11 @@ grab_forged_keys <- function(forged, workflow, new_data) { )) } if (is_epi_df(new_data)) { - extras <- as_epi_df(extras) - attr(extras, "metadata") <- attr(new_data, "metadata") - } else if (all(keys[1:2] %in% new_keys)) { - l <- list() - if (length(new_keys) > 2) l <- list(other_keys = new_keys[-c(1:2)]) - extras <- as_epi_df(extras, additional_metadata = l) + meta <- attr(new_data, "metadata") + extras <- as_epi_df(extras, as_of = meta$as_of, other_keys = meta$other_keys %||% character()) + } else if (all(c("geo_value", "time_value") %in% new_keys)) { + if (length(new_keys) > 2) other_keys <- new_keys[!new_keys %in% c("geo_value", "time_value")] + extras <- as_epi_df(extras, other_keys = other_keys %||% character()) } extras } diff --git a/data-raw/grad_employ_subset.R b/data-raw/grad_employ_subset.R index ae063d22f..38719a02e 100644 --- a/data-raw/grad_employ_subset.R +++ b/data-raw/grad_employ_subset.R @@ -101,6 +101,6 @@ ncol(gemploy) grad_employ_subset <- gemploy %>% as_epi_df( as_of = "2022-07-19", - additional_metadata = list(other_keys = c("age_group", "edu_qual")) + other_keys = c("age_group", "edu_qual") ) usethis::use_data(grad_employ_subset, overwrite = TRUE) diff --git a/data/grad_employ_subset.rda b/data/grad_employ_subset.rda index 3d74741cbd65cfcd114a7a224baa610478ee49e0..9380b43b5072674686e79ae494d9be9af2afc14a 100644 GIT binary patch literal 9546 zcmV-QCAHc@T4*^jL0KkKS^D6;*#L9I|NsC0|NsC0|NsAg|NsC0|NsC0|NsC0|NsB* z|9}7g|Nr1CyZ`_IH~@2`000aCIrIu5A|vIh?qHx6-$<{8Iq(1g00000008)X3QAEx z6beu&UH||O9CE9D+%)hWd)^0K096jEf1?0s0Wg^~(WWLQ0T_V7X`zrA8X5+MnhCTC zgu_4rV5fjgqZ1}WMrxVmJxrs~gvqJtrkON)MkYq4(?KurN(Ai)|kXu!k(mlz=kFxCQOWs zCYc&wAi^Gy#9}mIFa%&GfM5hli73&vJw`GzO!XxFQ%_=Qn0Zm-Pf(}nOrEEx+iIt% z-l??Oo~fEPq->|8dY*{djF~YS4LqiRdNPNidW@Q8k*BprH3LVcr>W{RX`|HjHl{{~ zO&+F(Xh*4#5h##NAOy{$CV?iIQ}ZKJXlPT(=pfY}pr@4bJdrdQ6V&vEMxIDC4Wejh z{ZZ;3NHlsQOw&+k(9_a1^&X>X4^gI#G(A9hL)0Fl5&DzVH29kP7WLkB7b}r+<#QUj zMT*9;)Q)1V>DM5QMqnKwbGK@tD5C(ZFh3++m^Y5Mm9~( zS9ej&ZOzW>w=7)gt;@T$T<+@RT@9IwxbB)Zxw%q1jdODD>D^rAlXDi=E2M6lk#iX7 z%UvYIa~rwTbPnsQo0ZP)h})X(F6Oyi3!I~K@LR>>qUF-(_g99+c$W?{7ah)BZpod= z+lq5@F6)`y&Pm+s?i<71*F!>HMR{I1<+{9&4b8gl?;bnEQ+;Z?&Fu0 zqsb#~_+&gxNjaosv@x+KSS-FF>1iMgHKcXwWr zuDmYgOy?InyA95}h0eQg&6mcRDf4ySd!muID=IySsNsbGg$e zcPDppc67UNuI_Pnb=}uH9FFeiM|Ui4yP8GL^4>gnyq-DZvE*IiUCX;SNDkq5G>3D$ zhZj4YJB~X9QArXIWD^uo3Wy{KmSE{Q5K)4}mNecK4dD=iA_zl%slI9~CAJIJgs*S9 z${_J_#J40<8=WrX6Gov~H$Yv~=hb5%Xf>%qt|ma)tn5P?PTysQs2FFhubP=a_%QO%wwWix>UlOut`K-jBCsxBrEYK=L zGicKq=LVbOi&JLTz8JcTwB1-sRfbL(j8|b=vzYWu+8cgxKf|Hg@_aU?(NT-G%SM^m zDYMueTU&8D+Y&XgSkcY6UqINc-?a3MwPRV=(XIAh9KBxSbDNshuH5Kh7rckuXOW~F zG}biZINUR^+#z}B@m*8Dn8R6eFVRjtD?bm>v~F*@BWtqDZ>9!v2iDm`QQCBEbxm}Q zq*@9g5wpv4o+mYN*{!zb^3Q9QWtiJ2k|bmjkb;FWX-k@uG0r6-OTD?rw8BhU&{cLt zCXl~M&7W@FM1`_pI0o~utd2@Tq0(3f>nAZ%twnat!04}QPLaQ}n{>Pz4t$5+r~TXb zoj355Nvt%Yr&fJzw3z{60***K+X7+@c(8@UB?{DNVB-jbT%hRBtXO31>=WSH5S&oHcv)_2DP)hY$_!5d$Nx!^$tTQJ4^%Az|jHvY#NTOD7K*EE8gQBn~ z#qgx%kHXP~1PG|uFTQWov(n%&GFE@a?bpxQ2f+O#{PIwGJdch);U6_sm)TZY2k=jq zQTq5EeCN*kK4nJ_><7KO^?PGKe=lpka6kI-5Thg0w@0jfo(&?N-yQjYRBj~p&ThfG zAPSL+h>ImYpi9jsJqt`gX4)0lNrP9 zCNaDwP9lr~!2(8Oep94UNfdIcp~C6KxK!X6GX9X5;+Z-2lxW5Cka!R*1~@BG+7g?y zQi@fPxwSJ_89!SuklRRh|E25590f#eFuj;u(6!T{L0kU>sIdlPdC=EL)$F7Hf&6Z8_S z#~PASX48LGCPiBM-^)A`FXSjJ6WBC^SH_;L97FHqSuy2OqmK@Ss@A#z)+8f+ZTwK& z7eG%k)`+7pzy}1NrqDPdrmmL2W^Niq=g7PDn1Fn=&kTGuVHOlH_p`F?fr+XQ3V)#C zxM-Briz~#4(64&I#&yn^xQrLY#g|q%lJ!{jtY*7#4v)L3_ld8Z)7vI*M*s3NO)VXZ<3bfhS(@|TGQlu+kST0i9BVczKz2^pu zqBaF>^~waE4XIAQC3`EVTzguaf^SE*H+Kj(i;af?Qn9sdxG}JUAVJZOi!w^;_9{J2 z1}*%uTbX^X=+X>EpxC>bcRjI{ z371e)T06`>R^CfaI(JGQfIOst)#(r2Gqt%2wppxAScp{m56RfwXzUuudnrR(vW}Ka zCkhA#;+Tl4e>(2!4UA5DI)v`2oDq0pWVpbquYiWwtbEW~^h86Se1WJyAzZ&#pv=M43vP8eiwjb`zjwUvp+ z?&u^~S3T`HX2snM)>k-u{WYps!qlqHl5&}q(r2h_8(z<0nljR!`0CAsxF|s$Wafbl zLZCWwQy^LJaiUReRdiE|wFwECw6ie+$pXTy6$+5`VJ}k3A|cr^LSdM`yfrG~+UbdU zl8u<0Av_?#ff2QUpr-BD5qZ5Rs4~5TtHgwyGgv z;>4mT7WQ&TK)BvdD{ikbldj%g-corUDCcl!2Ot3`C?&ZW2%v|VP&R-;JG&JIfYb=aN!tYq6arF;QUrTj z8kK=4ZF$#x&vW0_!ek#y6L1BaJIV^!HUY#+5UG|kD))_D66V)=KW(=ac#~yQOBGFH zt*Y>rwSLc4b8Km^p=z;Ps=QlP^ZTg{y0&)}Hi7Y~4FHEwStzxr6dF>CK`dB>Dp~^B zuA(?a5EeUKrs9&tR0Ti2Afgzes)JixTG&7mxO$#7mGSf?Za^KYu-Y_;3LG}sKutP` znk7XzLJD{QDs>2|VTh@hMN)AHqN3SAK?OE21Qk5Zr;2*f?Szs*5h=Wlu2^Ldh(%W3 zRq*yrv{cHfyrtn2BB>%MdV=Mb?{4c!u*aQwjR;8y5-0)aIf)y2y;WU=3(hvx8%C)j zR}M8ca==AvdQ$3xnE;T3O0-l+G-OazJ9owcL=wFrW&o3uCLr+C14#mcf_Sk6QE})1 zBt#-X0t8YD2>^&HBt$?8BA_WE8kmBBv^?K0UxdgZWi#_+RkUij3weim6q9uLVDF5$vl0dCY z_@DegRo*@1g@7%-&?OK7@x&t@SOXAj&$ILe1tQM z?tY-qwucuwfQV>KK;ciF{e1ma*eGCh@HG(wC^zm5uMhNImrq+z&2f6J?9Gf!*MzyS z+G2u(s2#Tns3!z|bv3tz1Q$ii3SoD~S%r7mdXc`}4G_tg*#V-L;S&ABAm1H#D`|r_ zpW{!rl%*;EMH+qODL;&*DLu6Q4@bdH<~uA!aiz-KCd6|5l(~JJ(p)c6^E&ZzwaN z4yYnV*aFg3GG zx0ohaVda^4415PHk(ggi<#Sv{_V^vCdym@akJ5M?Vf7(KA)v1obq-V5apNGn zs^{IxEVx|M%gpegL#jGTEK+Bw5F{9cBbkYG9>2ODkN-}h;lzn@vRJtnlQe@BGT-I# z z85yoWT8f6(HJ6xjFt*z zyxDFMCX@qg6jVt?NwEba8Yv`Agv6OF%w$3e0~Scdu`^_>8Zu}`AW9&O{UyoH5e{Q? zSaV!UyK}nkYo>^d1xjMuiEE@;nIb5{LX1}~ms|m2cU>k(N(hLYx&(nG2#Qji=pcfE zl)@GJ#>G zC`v>lQ$bClWXMT|Bus>iQV{|K0+S^q4Jgr=ymwfH0Fp@&CYS_Nn25zOn*|U87GP-8 zB*baqbW)N938Mxy%4tkvAjZKF-E?gbk`p5m6DW{2GE9iUnGCW}Af$M3iZK!a8JJY? z>IoW(WZ1D9n!Qhj$P}l1W5JM3j_c45LO!Ws#9GCW#c0jDVIT#6^@U0!tvN zV8&otNs=^#l1+dN0Y){56s8E4QjHN1Ln2^NjTDG#HZ*Kf=V}UKLYSmT#8VQ9OBBOF zvV#i-jT$0I(<~%xR8ff1G>n>L<-kZ8F&Zj1G|8sL0i;VMNEDMPCQ}he(J?3p#F{cj zreKL7ATN@=8^AtEuP+nVm~8z#YkN+`)jj8nUEqa}pYX(1yN z5T=PKD2!3DsM2j3PxET3{;#%ww$--w>Os%D6Jm$Yd4aoexSvPP;zSkN1>BEd9kv;= zB=_Bd0&JcC=kPqesAyv_cBW<8IRrrU;k1PA|F-K@B6J=p#>EI_xn~%XIbr#nXEPj$#os1!VpesA^RD`TqMMDBE8}O7KjKyAQpra60l?;BP0=MLMKT<3i-%`Xg`?T zSpXONhy~n!Yv|e-BdPqxC!Qh4!slzh&(RP*qhU8L^Z;3Q3vMWLvjXzmB+vFdV_RA zNigfcv;_;u`*HN$|EK#-uH24+ae(0OP8X@?TJ z0$QSVVhZdKB>?474pjsdtxhT56*#wUfK?O=X~E#oPB&+HP%gIs5G@Don&`wqv<9#c z3rHYWU=R>~A}E$!Q)AsauMz(a+<2|jQ98F&+kj)0K?5(@uxEF{op-y4rpdO{!IOB|=br3BgiRXwbyHKw90~qLCP`Jyf1VFk7H4T!~;O}XQYoM@+ z_W|uL-A8O4MreY&M7rycy>%21l86;pk3QqUaAP zieqJ&zfn^mK!Ub5XoGQR2sAdKwe)XV9b>T$h#JI)7a2hegA!!)M7o?a=65==C4+Df zakXfQ1@u6OgOUh7iQ$)5uI8x+$k^WH1C4B|wH_-NK7k9CFR{8ymBB8hzgyXSt;k-` zG&#W&oEANKp`7*RdoBL+@?Z3okVe1A4i z7f1WsS&sXJy-YX_&+hFE63sf}wFAwm81jL1Z&9(u4pJ99v=b zgca0Jh5|b#M51NrkR$qjpu+OZUv{Ru19dcsuedgy#QJ^LpEx_)d z%c=lyYE40q=zb8koRHD;Vb)9w>kKr5c`&$8Kg+;v1}1j_uZ39I{}gHvR7Fa3QnZ4r zpLJbc3}WFjK&T->)C1ph^|%yJam&0c_bw)I1zowGiu3W#E4V1u__PDke$KFQj&L8t1_0KBjp_lT1WkURsOB#fyW?j(Cn5yeb z`XMrjl2oKZ%6+$yOZZ2JHKEF|rzIgHe77{Pq=^vElJc=D@CyP?AvJAESyPlQ#W>8d z@{AyLM6OXzoAf_)kvrIF7-vtiMrEf@0|L2_n4C+5PLxN6<62ppEe`Y@WpgVCwUAL+ z5Vk)@0c#9r(wbAoU#Psh$FCevDctKAok}(n*}+FEG)3US*cNWUhr*{_ePPcK`k25w z3+iKI8)b)?zPy&osY#Vm=Q=TZtO?6^JzR!~MaBojaMgu7O0d|@?_a>Hq02-!_IamX z^KhQhTy>>c>mvf^2s~v`mb9TNwL+djEat*Sz9~*rcDsa+@k3G=7r=i{T(-3>W;Hlk zy&P~7wLh67pKCK!0H(gTBZc97iDU;inmiqS9I6=2G&oJO!j3PcJ;8M<3dbpxqk58r z3qNtc@!2$S3#E0UC72TF;Wi1?1Hn$+&@9>qC0;xl(|=E8ER4r07~Ppf_~WQBR z)pnk0^wJiYozRKzs7$>6K;aW}B zyH=^qoD3rPS7PyBJX%5Rj9y@JY{u5)dDT8{t|pfgYQ6(eS@XONJi36kP23aVQffz@ zJnTf+leB_imZBP3|BE(+N%Z5;wl$XjS5PgP6@+w3C4Hwp{|0+Kz5zsZ+%W=8aLZ*F zp+xbysBruporTjCGqRaI>8~*4BTH%AL_?^wIJ@;;MZ6wr$S5w95GXK!%cdaAC&DQM zlyx=UGc=`cCEE!W4KB|BhNSPfHfR%p3E18v0VG2Z*ED_%V|T4$Yl|Rx*>`p!5=K@Q zVPMM~W<^XMUOiM?$s(C#fO6-=d$4oR`zEA?3zs!TgrurR?>_*b#UOUBD&mxnNAF=I zJv~2ke!T~`K1BxXnU)6T09rBX%#5B^MJjwZsNWxMrK-->+DRxaewT8&G^iHL9mA6- zZI}j~V)tth>}SI{v=eiZ7uRwZ!qDzgCg*V(K--Yqr%xc;jfYZh7J;`T-S)roQl51R z{mg}*u{>XQ+keHDt!5kfC5k$|_8}=#b$H1Llp=lzE9jV;>wD5C0;8U@Ni+h6d!A_w zWpjf?&ug4O(KFymWKT#ON)px_nr7icc%DwLdPtQXsxD$TwNogxm$^ttKr)i1m|$=f z)4>8t22f0{07jWa>YNlZqq?U$wm5*O2j3!}&n{QB@*LRj2h=aU%X29+q#X6>T&aEq zCLsRZ8jQ<6H_QVxgW-Q*2nc>gz$gX`Bp?tnU=$7Eq-okgn==%bS+nR%xlq1^@c_Tu z&9<{VmI7$f1YjH}nz*F{6f*TfCK0Z#T`-I~9v4LD)TVEsXJ`5MOyq|UGs+pv&| zCwIdRA-^+Zyl>6;6xQ}M4g%o}1y9!ayAKc8@xER^SEK3i%Bkx++T`QcQI**nGs_E_ zCWxNq(ed#cq3;8DPCexDa)9&UM)OCAEsfRPjSEhPW~mO%{!{*@{K}95^c2|uC`Dz2 zWS~J;ke>hYJFT(E_J5L6Kp;wqAPtN5N%D^&;+ELFW$hQ{FZ~ymFjQ@kbnND~y7odN zm7+XpzRX2Bw!n`DX#ts)hqm>L#zpe|*wdS(5Q=gATF?ykFiBiAo#t^i>N=C|*K>}T zK@OpOjLq1pfYPc@%Fr|b0)oxV9W)lBw;0GMURGgS7aMZ75Yn8qSf1pgN&bCz99KP< zmh|+LM;B$UWuIg-qVyQ46v>|zCTR41GadesQ^m z5uu?&(6*0dCR}*sN~%YjqPl2thU9#OlS%*+o1UWSzpz0oM2(_Qn=3_{Els*uwM^tW zDvjSE#r%$)#(2)3oAMrbq)9y3NHJo~*;f$I3v-giF_#Z2h(s>kR&`pC$KbB*fWx0y z5L=)vT6+F8({E2NPE2n=;RUJ%hH%-1yP9hz(;q}C?pjhh2YxT9d38-|-Q%p@xc zez8?R#fHXem3?Si^^*%%1cxf0T4c!403(BcC>Tk+&o9W|zn>!{G6N&p5Hbvm^xQe_ z;z=jz=sPU-cVAD^p7u{e<5RWx9V}$PhVu-GTm--wPVK=Vly1xrCCm?0V1Wvq^7dP@ zmE9-M{4WP~7L&i=5cjm9y>dd;_j(yQErvXQy~o|Dap~2 zFrd5Si+=Jl95||!i){eUj8N}1cj-tDc{!R94A_Vq?zVubqW6?R-vkijx-?LFet|Sv zsAPl%Bl^9AQWQ0NEgb+ETLGAr4KAKW-Yualiv{#Z4PZ}zoh5`raMFGkP zfyv4W3z9%i^RU*|gp$P!3`D`A5<-AsVkV1*5tK&>1aw>_5jW6NeD*NKP#@HS@4EUR zJo8AD5MF1ili&E%rHTrB2r3wj{Q*ov1P|ZcX7Cu60GODdZ^JSYh>j)-@65^BN0$x%>gnDu8B@N0rI=J*UCnZn#7VE{bT}G%DVN1#d7!JW@h3inu}oUWHi#gbb_^y=t^8l46DALAWBF?c&i#euMeUOS9%`{C_Gb*$IAGiA8ia-L2 zNCof!1BT!Ld(~ch-h&XJ8C5^pGyu?}g^fU~n21%yTl*BMh z8kscpPevn5nrdo1lM(7;Q)#K>^+ukaM9>;-1vZS6)Wm6%Bg%RZ5RpwMihD^t7)?x? zd8iLm*qKME`lq3$r=b}%Xk^4>XaEfZKxohbqd;Ud(9xgL6RH{;gK^OoN00LkD5t9G_0GI#(01=SHFaRb1000000003n z000DFKq8Vrg)ym#klKwj@TQtHZ6g{8XwxT&JxvV(qX8N;13(5rq3RlFG}A*PO*Cj| z0Lh>RfDJMRO#o;B02*irl!Y`VM9Jvb2#RJhHldI-n3VMWNv4{fQ}mgMPg5h*(?*&b zQ`Buh4X7TF!g(@{sQpvadYdHK3F(y6M#4=D${SIkk?J(`jXbBR>KjwkdY-4L>UxHs z1zwb+7D_(tsKzyOT)QJtXg7CU?&>Nn-E>EHUEQ}`+eXx9bsKJ#1#WH1Hlu7rMRRS; zZ8^r~Hs;-!&gj=FEpvvE*F{e1jm~a&b|*Ty!(7)3yMm5sb;YF7$k#DPE1cY`oW&Jh zVd2=k5o3kQxm>PxZtJAj>2Tb+u1z>QuI|h)2JmX|S9ccne-8(TIlH=ZtX?N|9tGDN zyPcs!mUE=KcXkNox^nLBMUG>;i>23dlIy#+n+T;0bd$8)-Ob4Pb> z@owbs<8nKlySqnl?Cy+rbaR=>*Ctm^?n|!e-LCGgxhF2@ySXLO>!$8HM~p7cO@*i@9~x!=&r8sUy3&+1XvkaYU2~p$R0CB&3o- zB!M5HG#P~^j{yd7q`$6k#IY!P$Z*zM$y6ktnE^4RvvtR-AmuYKZTUo_f4uoPy2v1h zUAcoD!9-!9Ftt-6GNyDoSB3|}kg3hA^F3wMrFC~^#6N7)aU||mkpYOE*9Pu18V$Lo zT(Q{SOf(UNb*1!Oaf>tv)(OEJol#_x;h1dSGe%cz96v1Wi-Lm*M_@GAXSTdoQ%<RAJ(wDYIgTSMU57F~0Etyj$~6D%FBG~_#Doxkgg)F}r})fdvwEoNu9_4aLE zrf{mU7$0`{@E5}^ce>ho+Az;nxgoCkGe(jvT~l#U{v#faV}$6onf39{`LpF0*~8=< z42^C4sm{RGR-x8@3kTkPB^tP#+bx`B)ZDhqHm=VIuyPk#yDhDtCR(_r#U5mjUmOePuatA`AeWpEnOlcRiGSv-sc`)&yld8j6AQfRWB+EJ&xFDB= z)Cq9gl5#V2YLQv7fC|a=y8Y{Z#H#JCen5_w+0OFytp#t2ec@M^8|(ThYx%Y;MhhER zau-II@###w>5e9P#qRK9)u`zbbm4cELAq!gpf>P@|s) zjr1w2#t>$_@4I9`3!H3YSpP^z6;9O%)2KEzO?m5ggV$DRTm!YgDL4bD*$4( zMySsnGn&>#gFghJ)ki@7C)_TQ*mkEp;uhh z-6X3n4)p(6MJ@dY*CKN4iDElmP2w$#9jij}P(FZ^pvP#QeFEobY0-J3ah$&^<*cGo zR;DOhx~U#rbAQ4RJ~{#@PR`THqWRCR&+@MV6&Vc6KU`fgC}p!Cu49O znKGwYM3`-)-vk#Jycu@)3zf{4H-{>sg{x~YYGl(G9tzA2+Krm-J9osuKbBgEQyB;I zehRr-Z6+X{_nie78w%?w1rcB#SWy7eKQ7GMOTofeMWQhPwP|5zFw)LX{&G``3;^Vl1)Z z;`>f_hjL<51oBXw@gURFUkBK3qpd zmFgzi?asl#b)Ag|S3u-dBCD@QGuC@!#+7lWjTEm8zc%Ldj7V6`at5rA7=%M~qTL(I zML8OpW%X90c}Pn!pV6k?smCk_5{L;$u(iz*$>lb-KB9>}2sRg~=-pMaCxJA|RHT!G zF?m8Ppg@?WsT|R)44cNGPm2m*?v#1l3s_*&8f}(8kVjdmu^V?u5V~$^rO6spddxW( zaBxygE3SODvE}khvv2P6T07k!Yd?fG+c9}j?KP;zTSanb4>-$-aNfBox}GZxZAU6< ztg9{*cm?gQ(W#)KAn|b1XlbNM=_^sI`P2o{1Q7dxiFgPg1P53VC8-?ij!G*o85L({ z=8eBqWv7{!7o0BuKi62+=C+8iyiOK>achDk!Z6Sak}|60t%1u{67Dgx%UPv^y0hB1 zuyS1xNFwir=LIE}1Y;QsHJbx_JH#`h(?tpHUa%ptSlJ64B{HgI4wtHo&g9Livb;1B zSf*J`-NCV(*I6;Wp12=F4`gyEHN~Kt=2JEZXJ-u4Tg?V|t-4$@VQ9kl#)kl)q=VoyD?x(EnQ3NV(F)`;ItnK<%jY9```)r5x{xr05e5@b!CG*E z6_^f~Gdq`Tj2@-Hb+B~cx{)FSD_BLUSwf3#;CPLXpI)pRTivhs}@Su#)v zN%E$;rGip2!q!Eq55Xcy8meg2k}KB$)Pf{RphT;HLa0%aBtlXJt60c@6`v+01Xmkl zB!q1RgBdaEJZrcKM!MnBQql4=GhEFC z1<7Wqr$k0_c$<7hbT&;Hn|qes66LQG$?f7BcqNw2G+Qe>TWch^j}M=C+iTqAH>N@IhA$ zKt(%R_?1>bL&d}3aI3Rj&RVR@b{yMku|WeHq)-z~q9=f;r}7}Cqy-A-DY=cDR_q{1NFX*H7|4Ml#sV!QBmpD@4rWF^KF=?6 zRg$fBUaiQI5tC}$%4Qg^9je=*G5(8fjYf!-U(2S(Rw#}>x>Fz&&Wj-ul{5tlAyp)# zmut(gr6iQ`%#luNedc z?qYFuf$zV@N*v2Tlb#KP-|m6Q5flI!*+NJne@=N@*!b2eD0UAXQ4mx-9}~~!RtDQ9 zf-BdMP*FXkQ2DR`E1BNnkKpHqaDVXUg@B-E1Qi2*XF!CQL*;ro{wWJU=aEHn7i)%U z#`qKkE67nx;FwbQGbp&rtqGcVU7x`Ir_lN~DvQ4zsh$U6za6uT3$=9zY#UueY)DAs zBALN8Gv;xqyogu_2JHj@--k@AJzG2pVp(>KfX|7h7ZXPbBCWRnYfGx%kq)!&>Ph2$wa6>Rc`Qm~~eV zgV3$uy6idgp^%QT>GUsCGocOD;4(OmQ0~vL_!^Ew=q+jOPRwXzu|~sRh$y2bID<>3W^D`hPh587{SZnpy2zqB}~LCE40l zc)pW_x>cO1LJPlO(6057spMU0b+xh{;NqU)qwV5*e=J+GB4ya!;rDtw9OJIYqWC?Y zw8(*h^3tdFRrEG({%MXE)VyAuyC$j731P{Q*)|u%f_Sg9Jl^(RepX<%*wSxW6r@?? z+f~w_tC3b-JcDmw+@i;4XB1?aUyoF?_e9pi8|I!y%b&My;m$yf_M-% z6}pb@JxH%jg;fIvkGFV2Lw$b*8fH&)L!?=YU6uGEcYxgWcQpKfSf*0?!%tR z-v_Pdx7);$l$t_iq_SayGYrf^2x2mrB$$d>H5sLnB+$boWLPu?rjSXc5sM~8BQd0a z(gPA%7|BBt0g3>d1|b_LNFz~znS?YNDKT0CD2TCyWMf8yAu@>>Xo#?4q}C}hg9)G+ zqM$Jp%0nec7Bp)jU?E6}h9tpnTXfJw%>-j?ltzjz6__SWhy+O)B8bEwWHo|~8VpoL z3`j;{MI|MA^Qg}5w>c=xQ**BCoZT#J6(JS7x-}a!2}+7fcXiwWib*-SpaOHUP2IX_ zk}!*#r4ls=g52OpkeLF+P*Nt5BNS&hW{j9*rddNIVrG;hP%=$-ZqkVo3Q3ldW+^5i zxpXNI2xzgnxFZsYicKU5Bt~GOY(Nzr9lD5&3>ZlmlO{?rWRWb;WYjQb z8Yoi&z!=byQHV4|n-M9aW>O;{P>^JyNkln#cv|d4VP;5z8Z00NQjrp1Vwn*oiV_S` zXu?oIi~=T1hD?BwiI`1<#$;%ug2q5HCSjDNG-OId!jh4qjTA;SV`3>N!4N>rfee-` zP!TX;69_W|Ofn>r#)Q;N2$KUbWWf=%#I(w3Xv~n%N-45N2q=*;5GA650VOg~kfRZh zsHBP_G=`)ogfk?g6lelFxrA((jRci4NXCX3DKr)$!4)MK(+nFVVx)vDSfZu{NYaoD zhCrhPRLq2wjFOTpVv@;{0!j$X)LTVK2G}rYlLj!FCTtYSB8rU~Brz6)7BW$w3{j#P zj7mmA2%s$nzXq#s!>wz4j@sA6*TnCUZ6Y_8mm!_|seuvAm?#L$xpy)hI-q7K;=U|Z(s3A{3{icZw!ihN5j84_(<&zds=ExOT2+<^h+imcNY1_!-`(;eReBnc)1frrK5b4qE4h;Na|kE`#uyVnSF z|6_1IkcKgqV@euwg>OmPR@Gkev{F=f6tA@_YZkz_9_e=vylYNK%#V#T;Y5Gsiw zBxHgu2t@qVRjbZaOCKD)TtFWrfMeg1Ej#(-fkZK?Dr~@OionJ`dW!8qTSY{LM2~<` zXwC(3)bRlo2#8K=0?S^awDk1t#wZ%)6{td_b~q{+rXq1B)^TNG^`K@Vtj0b`P&Ciy^2a5!NLMTm9=Wfjt0@KBkxDVAq@@G?W20+tV$ zu0Fg0a}c$!y&#lsbYkY>-%il0c1!hX3V}V5rVd3pyz+o)uDy{WVEvE3Q zv2DZ`_WT-MC#fiwnp1WWHg+w za1c_l%NhX{RXC|rdE(EJrNCfA~nj*cvhIi(pHNU zf)G}cqD&PiOz)O;J6D zye(sfNI=M-0?M)NA|b-)Y(PUS2xFwoIN)P|0gi@B8M>B@kRmxO6kMDGMQQ?|ro^k- zn$j{!vR*p!dH3J)n*)q`3==MG4*5V$odFTo?G+Pony8dEZqKULAV-a`Iec@1xeBXv z+NJ(i2#!MV4pLEZl@lH)u4pS265672M#Hk%b?+{*Jiga^bS@hI9|L2s$OCY(+;Rdt zg1Znz`?6y2=M`BdQ7PqFg&7$wsgt~?D}_E-XJJwv6K}utHud^-UFqdh{Hsd3u< z)}1~iY~EO8og~@${jJPDN~VzVZ{xxIwi9R8`<=LbC(uIMk6p&8(#eYx9$J6T&E z-D2WGK2)BeH$cKP#le0KYQNQGmCx-lY%KCAiwQAG9yZf)-_az0yTabfV4iMKaWhlo z3>Iu}*kZ(`0Xb?P*GL&cPA?WXej@vV&z8&{>QYI7u21&&9}Zt>vuscgN)?pNXqG0? zrPe8AZ|o8|u=kXEe|LPnQpw}#4%)u&tjLpMifq+?Mmg&rUWH=j59iZ|=_Ql#Ar~(0 z{=&U%2j~&w6gS$!9Bd`(z8=N$FqY?k0n|I32G5K{4CP&J&!EnbWJv{L@fY zR?%oo^j=%QKJUu%i|R=w2s73P@Hgf10Bd39xa`p=&(J8Fm}X5YfwHl(z^H|&X6us_@QxGKF>DeS&Eg2O zc#3wkD-qfd%%qa1w%=;U)N@*V)n|v#Eh3aY+LZjP<;z+Kmm;Ij-VT#t=1?u1bj?KK>1EB~b1|(IsK3`j zJjx3MMimERgN9X9gHAfE2%3PzjMcSfaV9S#gp8LZYf5@juxez1CeIIJHsg+R#s+ap zVA%m^XxdTXPDdU{D;otzMky!g_m)2%b{4^#6@-LteiyI zt=1i>HkY(Tf}6$ch^b*~7qk@#6>7&JF$!wUR~Jg&$fha<2{Q9eWr;Eoi~rLsBzT%cQHpj?td+fwKy!rfm_ z>CWtwr9L0ueJ@ULM;r`65b+J;Kn|16qAk#eN9=U)Nn56A1``5SVA3W9t4PBVN?<0D zhEdXCx(;O%Z@W#=SI#F*1l!If=G#{2xE8{zX6Tf6O}aM8@i%S4tC8E>l-1l3jXdaX z5+SEGKhvb_z8pJa(yFeql%%dx=aRcMa$MsRp44J9zP6OBIhePMDCe8Q>Pn62Fu=FH z4eUN8eH?ee2Z!c#oAa!&Qs$WwfxgotOKj6F5Q%Ln4jyO zJ=wx@_$7S@1)?3theGoVL5liSq_`RzZ6U!X&E33QXGz*o-LUC7{l_87>-0q}zS6E6^g4N4GDDZmO~zAphkJ6z^_cHcV#rFzma ze3p_$nCcPV`r>kv4^NBLNUZz%ayf#&?55!;ofnBDb4#0EA$QTl=QFs;Mp!c<)Tm09 zZ`xyt#90@P4{GW$G!sL#40*m+TAwbFafECO^1m;>v8&a;ks_T3&77IkZ!5z zVf##GA|%wzB;YHMGEt<{^Q9Nq!lGd;u=B5=oh+v2PFBv@y|I1528yCx5|pusr5edW zG6@^v6D#@;NSziUY23471G_b)Xyi#$1Y}H!F=It_0+SwMTFI@WGCrzF%0-#B-G5`% zySJ6vbQ@lSCeysb@qBMdzR-F%}8CI(fN%q?u9B5N4o`mVkXmUZbM`kPS6?QA@xdCl6TRR(e%D%+y{kwFz$9oWSzz9zjLIo)T8lccp3Bfs{1Qh0qf(m7e z46}{=DPuy8jldU5*tQAUz+yfJEqilia_+9t4{d$wWhG3s8 z{yvP<2%`+9#`BYLGwEona3owrRH*>1U<&}QB!y}~mH`0oP^I_?0KY&7Y*-RhFRuF~ z5Nhj4wMmd`88TM(wWN_N0uXO?qeNWTiRo%BvY|}UBIfE*UNiyUkWlIQM?pj6^*#2e z+7;ljK>CmnIZ;6`1Q&u|C@FM;oc{=Kg-YIF4vTz}f)_$4cG@@U5HH0I7X1_+i$y4a zZbB;VeNzTV3D+90Gb=F>DNE;`BCX`AStu7pfdfC$VoWjxVz<}}iNQfj;{NKaMi_<+axn_U0^ diff --git a/man/autoplot-epipred.Rd b/man/autoplot-epipred.Rd index 27bfdf5f7..10236eb98 100644 --- a/man/autoplot-epipred.Rd +++ b/man/autoplot-epipred.Rd @@ -121,4 +121,6 @@ arx <- arx_forecaster(jhu, "death_rate", c("case_rate", "death_rate"), args_list = arx_args_list(ahead = 14L) ) autoplot(arx, .max_facets = 6) +NULL + } diff --git a/man/cdc_baseline_forecaster.Rd b/man/cdc_baseline_forecaster.Rd index cd3c4ed67..0c7f1e436 100644 --- a/man/cdc_baseline_forecaster.Rd +++ b/man/cdc_baseline_forecaster.Rd @@ -44,11 +44,11 @@ weekly_deaths <- case_death_rate_subset \%>\% mutate(deaths = pmax(death_rate / 1e5 * pop * 7, 0)) \%>\% select(-pop, -death_rate) \%>\% group_by(geo_value) \%>\% - epi_slide(~ sum(.$deaths), before = 6, new_col_name = "deaths") \%>\% + epi_slide(~ sum(.$deaths), .window_size = 7, .new_col_name = "deaths_7dsum") \%>\% ungroup() \%>\% filter(weekdays(time_value) == "Saturday") -cdc <- cdc_baseline_forecaster(weekly_deaths, "deaths") +cdc <- cdc_baseline_forecaster(weekly_deaths, "deaths_7dsum") preds <- pivot_quantiles_wider(cdc$predictions, .pred_distn) if (require(ggplot2)) { @@ -62,7 +62,7 @@ if (require(ggplot2)) { geom_line(aes(y = .pred), color = "orange") + geom_line( data = weekly_deaths \%>\% filter(geo_value \%in\% four_states), - aes(x = time_value, y = deaths) + aes(x = time_value, y = deaths_7dsum) ) + scale_x_date(limits = c(forecast_date - 90, forecast_date + 30)) + labs(x = "Date", y = "Weekly deaths") + diff --git a/man/epi_slide_wrapper.Rd b/man/epi_slide_wrapper.Rd index 0c05b7650..d67db1c88 100644 --- a/man/epi_slide_wrapper.Rd +++ b/man/epi_slide_wrapper.Rd @@ -6,8 +6,8 @@ \usage{ epi_slide_wrapper( new_data, - before, - after, + .window_size, + .align, columns, fns, fn_names, diff --git a/man/flusight_hub_formatter.Rd b/man/flusight_hub_formatter.Rd index b43bc0ac2..b2be9b4fe 100644 --- a/man/flusight_hub_formatter.Rd +++ b/man/flusight_hub_formatter.Rd @@ -52,11 +52,11 @@ weekly_deaths <- case_death_rate_subset \%>\% mutate(deaths = pmax(death_rate / 1e5 * pop * 7, 0)) \%>\% select(-pop, -death_rate) \%>\% group_by(geo_value) \%>\% - epi_slide(~ sum(.$deaths), before = 6, new_col_name = "deaths") \%>\% + epi_slide(~ sum(.$deaths), .window_size = 7, .new_col_name = "deaths_7dsum") \%>\% ungroup() \%>\% filter(weekdays(time_value) == "Saturday") -cdc <- cdc_baseline_forecaster(weekly_deaths, "deaths") +cdc <- cdc_baseline_forecaster(weekly_deaths, "deaths_7dsum") flusight_hub_formatter(cdc) flusight_hub_formatter(cdc, target = "wk inc covid deaths") flusight_hub_formatter(cdc, target = paste(horizon, "wk inc covid deaths")) diff --git a/man/step_epi_slide.Rd b/man/step_epi_slide.Rd index 46bb386ad..242f8e312 100644 --- a/man/step_epi_slide.Rd +++ b/man/step_epi_slide.Rd @@ -8,8 +8,8 @@ step_epi_slide( recipe, ..., .f, - before = 0L, - after = 0L, + .window_size = NULL, + .align = c("right", "center", "left"), role = "predictor", prefix = "epi_slide_", f_name = clean_f_name(.f), @@ -41,19 +41,25 @@ argument must be named \code{.x}. A common, though very difficult to debug error is using something like \code{function(x) mean}. This will not work because it returns the function mean, rather than \code{mean(x)}} -\item{before, after}{the size of the sliding window on the left and the right -of the center. Usually non-negative integers for data indexed by date, but -more restrictive in other cases (see \code{\link[epiprocess:epi_slide]{epiprocess::epi_slide()}} for details).} +\item{.window_size}{the size of the sliding window, required. Usually a +non-negative integer will suffice (e.g. for data indexed by date, but more +restrictive in other time_type cases (see \code{\link[epiprocess:epi_slide]{epiprocess::epi_slide()}} for +details). For example, set to 7 for a 7-day window.} + +\item{.align}{a character string indicating how the window should be aligned. +By default, this is "right", meaning the slide_window will be anchored with +its right end point on the reference date. (see \code{\link[epiprocess:epi_slide]{epiprocess::epi_slide()}} +for details).} \item{role}{For model terms created by this step, what analysis role should they be assigned? \code{lag} is default a predictor while \code{ahead} is an outcome.} \item{prefix}{A character string that will be prefixed to the new column.} -\item{f_name}{a character string of at most 20 characters that describes -the function. This will be combined with \code{prefix} and the columns in \code{...} -to name the result using \verb{\{prefix\}\{f_name\}_\{column\}}. By default it will be determined -automatically using \code{clean_f_name()}.} +\item{f_name}{a character string of at most 20 characters that describes the +function. This will be combined with \code{prefix} and the columns in \code{...} to +name the result using \verb{\{prefix\}\{f_name\}_\{column\}}. By default it will be +determined automatically using \code{clean_f_name()}.} \item{skip}{A logical. Should the step be skipped when the recipe is baked by \code{\link[=bake]{bake()}}? While all operations are baked @@ -80,7 +86,7 @@ jhu <- case_death_rate_subset \%>\% rec <- epi_recipe(jhu) \%>\% step_epi_slide(case_rate, death_rate, .f = \(x) mean(x, na.rm = TRUE), - before = 6L + .window_size = 7L ) bake(prep(rec, jhu), new_data = NULL) } diff --git a/tests/testthat/test-epi_recipe.R b/tests/testthat/test-epi_recipe.R index ed27d88c0..a4cbb00b4 100644 --- a/tests/testthat/test-epi_recipe.R +++ b/tests/testthat/test-epi_recipe.R @@ -59,7 +59,7 @@ test_that("epi_recipe formula works", { time_value = seq(as.Date("2020-01-01"), by = 1, length.out = 5), geo_value = "ca", z = "dummy_key" - ) %>% epiprocess::as_epi_df(additional_metadata = list(other_keys = "z")) + ) %>% epiprocess::as_epi_df(other_keys = "z") # with an additional key r <- epi_recipe(y ~ x + geo_value, tib) diff --git a/tests/testthat/test-key_colnames.R b/tests/testthat/test-key_colnames.R index d55a515ca..3b3118740 100644 --- a/tests/testthat/test-key_colnames.R +++ b/tests/testthat/test-key_colnames.R @@ -30,12 +30,12 @@ test_that("key_colnames extracts additional keys when they are present", { value = 1:length(geo_value) + 0.01 * rnorm(length(geo_value)) ) %>% as_epi_df( - additional_metadata = list(other_keys = c("state", "pol")) + other_keys = c("state", "pol") ) expect_identical( key_colnames(my_data), - c("geo_value", "time_value", "state", "pol") + c("geo_value", "state", "pol", "time_value") ) my_recipe <- epi_recipe(my_data) %>% @@ -43,16 +43,10 @@ test_that("key_colnames extracts additional keys when they are present", { step_epi_naomit() # order of the additional keys may be different - expect_setequal( - key_colnames(my_recipe), - c("geo_value", "time_value", "state", "pol") - ) + expect_equal(key_colnames(my_recipe), c("geo_value", "state", "pol", "time_value")) my_workflow <- epi_workflow(my_recipe, linear_reg()) %>% fit(my_data) # order of the additional keys may be different - expect_setequal( - key_colnames(my_workflow), - c("geo_value", "time_value", "state", "pol") - ) + expect_equal(key_colnames(my_workflow), c("geo_value", "state", "pol", "time_value")) }) diff --git a/tests/testthat/test-layer_add_forecast_date.R b/tests/testthat/test-layer_add_forecast_date.R index 9595b47b6..915804d6b 100644 --- a/tests/testthat/test-layer_add_forecast_date.R +++ b/tests/testthat/test-layer_add_forecast_date.R @@ -93,6 +93,8 @@ test_that("forecast date works for daily", { unclass() %>% as.data.frame() %>% mutate(time_value = as.POSIXlt(time_value)$year + 1900L) %>% + group_by(geo_value, time_value) %>% + summarize(case_rate = mean(case_rate), death_rate = mean(death_rate), .groups="drop") %>% as_epi_df() expect_error(predict(wf1, latest_yearly)) diff --git a/tests/testthat/test-layer_add_target_date.R b/tests/testthat/test-layer_add_target_date.R index e5349839b..f1fa3f217 100644 --- a/tests/testthat/test-layer_add_target_date.R +++ b/tests/testthat/test-layer_add_target_date.R @@ -104,6 +104,8 @@ test_that("target date works for daily and yearly", { unclass() %>% as.data.frame() %>% mutate(time_value = as.POSIXlt(time_value)$year + 1900L) %>% + group_by(geo_value, time_value) %>% + summarize(case_rate = mean(case_rate), death_rate = mean(death_rate), .groups = "drop") %>% as_epi_df() expect_error(predict(wf1, latest_bad)) diff --git a/tests/testthat/test-pad_to_end.R b/tests/testthat/test-pad_to_end.R index 0ea6244b0..6949f06ac 100644 --- a/tests/testthat/test-pad_to_end.R +++ b/tests/testthat/test-pad_to_end.R @@ -32,6 +32,6 @@ test_that("test set padding works", { # make sure it maintains the epi_df dat <- dat %>% dplyr::rename(geo_value = gr1) %>% - as_epi_df() + as_epi_df(other_keys = "gr2") expect_s3_class(pad_to_end(dat, "geo_value", 2), "epi_df") }) diff --git a/tests/testthat/test-step_epi_slide.R b/tests/testthat/test-step_epi_slide.R index 29e046eae..e90c317f8 100644 --- a/tests/testthat/test-step_epi_slide.R +++ b/tests/testthat/test-step_epi_slide.R @@ -7,32 +7,23 @@ edf <- data.frame( value = c(2:21, 3:22) ) %>% as_epi_df() - r <- epi_recipe(edf) -rolled_before <- edf %>% - group_by(geo_value) %>% - epi_slide(value = mean(value), before = 3L) %>% - pull(value) -rolled_after <- edf %>% - group_by(geo_value) %>% - epi_slide(value = mean(value), after = 3L) %>% - pull(value) test_that("epi_slide errors when needed", { # not an epi_recipe - expect_error(recipe(edf) %>% step_epi_slide(value, .f = mean, before = 6L)) + expect_error(recipe(edf) %>% step_epi_slide(value, .f = mean, .window_size = 7L)) # non-scalar args - expect_error(r %>% step_epi_slide(value, .f = mean, before = c(3L, 6L))) - expect_error(r %>% step_epi_slide(value, .f = mean, after = c(3L, 6L))) + expect_error(r %>% step_epi_slide(value, .f = mean, .window_size = c(3L, 6L))) + expect_error(r %>% step_epi_slide(value, .f = mean, .align = c("right", "left"))) expect_error(r %>% step_epi_slide(value, .f = mean, skip = c(TRUE, FALSE))) expect_error(r %>% step_epi_slide(value, .f = mean, role = letters[1:2])) expect_error(r %>% step_epi_slide(value, .f = mean, prefix = letters[1:2])) expect_error(r %>% step_epi_slide(value, .f = mean, id = letters[1:2])) # wrong types - expect_error(r %>% step_epi_slide(value, .f = mean, before = 1.5)) - expect_error(r %>% step_epi_slide(value, .f = mean, after = 1.5)) + expect_error(r %>% step_epi_slide(value, .f = mean, .window_size = 1.5)) + expect_error(r %>% step_epi_slide(value, .f = mean, .align = 1.5)) expect_error(r %>% step_epi_slide(value, .f = mean, skip = "a")) expect_error(r %>% step_epi_slide(value, .f = mean, role = 1)) expect_error(r %>% step_epi_slide(value, .f = mean, prefix = 1)) @@ -45,31 +36,40 @@ test_that("epi_slide errors when needed", { test_that("epi_slide handles different function specs", { cfun <- r %>% - step_epi_slide(value, .f = "mean", before = 3L) %>% + step_epi_slide(value, .f = "mean", .window_size = 4L) %>% prep(edf) %>% bake(new_data = NULL) + expected_out <- edf %>% + group_by(geo_value) %>% + epi_slide(~ mean(.x$value), .window_size = 4L) %>% + ungroup() %>% + rename(epi_slide__.f_value = slide_value) + expect_equal(cfun, expected_out) ffun <- r %>% - step_epi_slide(value, .f = mean, before = 3L) %>% + step_epi_slide(value, .f = mean, .window_size = 4L) %>% prep(edf) %>% bake(new_data = NULL) + expect_equal(ffun, expected_out) # formula NOT currently supported expect_error( lfun <- r %>% - step_epi_slide(value, .f = ~ mean(.x, na.rm = TRUE), before = 3L), + step_epi_slide(value, .f = ~ mean(.x, na.rm = TRUE), .window_size = 4L), regexp = "cannot be a formula." ) + # expect_equal(lfun, rolled_before) blfun <- r %>% - step_epi_slide(value, .f = function(x) mean(x, na.rm = TRUE), before = 3L) %>% + step_epi_slide(value, .f = function(x) mean(x, na.rm = TRUE), .window_size = 4L) %>% prep(edf) %>% bake(new_data = NULL) + expected_out <- edf %>% + group_by(geo_value) %>% + epi_slide(~ mean(.x$value, na.rm = TRUE), .window_size = 4L) %>% + ungroup() %>% + rename(epi_slide__.f_value = slide_value) + expect_equal(blfun, expected_out) nblfun <- r %>% - step_epi_slide(value, .f = \(x) mean(x, na.rm = TRUE), before = 3L) %>% + step_epi_slide(value, .f = \(x) mean(x, na.rm = TRUE), .window_size = 4L) %>% prep(edf) %>% bake(new_data = NULL) - - expect_equal(cfun[[4]], rolled_before) - expect_equal(ffun[[4]], rolled_before) - # expect_equal(lfun[[4]], rolled_before) - expect_equal(blfun[[4]], rolled_before) - expect_equal(nblfun[[4]], rolled_before) + expect_equal(nblfun, expected_out) }) diff --git a/tests/testthat/test-step_training_window.R b/tests/testthat/test-step_training_window.R index f49668a40..a9f2170d3 100644 --- a/tests/testthat/test-step_training_window.R +++ b/tests/testthat/test-step_training_window.R @@ -84,7 +84,7 @@ test_that("step_training_window works with multiple keys", { expect_equal(nrow(p4), 12L) expect_equal(ncol(p4), 5L) expect_s3_class(p4, "epi_df") - expect_named(p4, c("geo_value", "time_value", "additional_key", "x", "y")) + expect_named(p4, c("geo_value", "additional_key", "time_value", "x", "y")) expect_equal( p4$time_value, rep(c( diff --git a/vignettes/articles/sliding.Rmd b/vignettes/articles/sliding.Rmd index ec6f67359..31cc7b9b0 100644 --- a/vignettes/articles/sliding.Rmd +++ b/vignettes/articles/sliding.Rmd @@ -54,7 +54,7 @@ applying `epi_slide()` to the latest snapshot of the data. First, we download the version history (ie. archive) of the percentage of doctor's visits with CLI (COVID-like illness) computed from medical insurance claims and the number of new confirmed COVID-19 cases per 100,000 population -(daily) for all 50 states from the COVIDcast API. +(daily) for all 50 states from the COVIDcast API.
@@ -69,7 +69,6 @@ versions for the less up-to-date input archive. theme_set(theme_bw()) y <- readRDS("all_states_covidcast_signals.rds") - y <- purrr::map(y, ~ select(.x, geo_value, time_value, version = issue, value)) x <- epix_merge( @@ -93,15 +92,11 @@ output. ```{r arx-kweek-preliminaries, warning = FALSE} # Latest snapshot of data, and forecast dates -x_latest <- epix_as_of(x, max_version = max(x$versions_end)) -fc_time_values <- seq( - from = as.Date("2020-08-01"), - to = as.Date("2021-11-01"), - by = "1 month" -) +x_latest <- epix_as_of(x, version = max(x$versions_end)) +fc_time_values <- seq(from = as.Date("2021-11-01"), to = as.Date("2021-11-01"), by = "1 month") aheads <- c(7, 14, 21, 28) -k_week_ahead <- function(epi_df, outcome, predictors, ahead = 7, engine) { +forecast_k_week_ahead <- function(epi_df, outcome, predictors, ahead = 7, engine) { epi_slide( epi_df, ~ arx_forecaster( @@ -109,9 +104,9 @@ k_week_ahead <- function(epi_df, outcome, predictors, ahead = 7, engine) { args_list = arx_args_list(ahead = ahead) )$predictions %>% select(-geo_value), - before = 120 - 1, - ref_time_values = fc_time_values, - new_col_name = "fc" + .window_size = 120, + .ref_time_values = fc_time_values, + .new_col_name = "fc" ) %>% select(geo_value, time_value, starts_with("fc")) %>% mutate(engine_type = engine$engine) @@ -121,22 +116,22 @@ k_week_ahead <- function(epi_df, outcome, predictors, ahead = 7, engine) { ```{r make-arx-kweek} # Generate the forecasts and bind them together fc <- bind_rows( - map( - aheads, - ~ k_week_ahead( - x_latest, "case_rate", c("case_rate", "percent_cli"), .x, - engine = linear_reg() - ) - ) %>% list_rbind(), - map( - aheads, - ~ k_week_ahead( - x_latest, "case_rate", c("case_rate", "percent_cli"), .x, - engine = rand_forest(mode = "regression") - ) - ) %>% list_rbind() -) %>% - pivot_quantiles_wider(fc_.pred_distn) + map(aheads, ~ forecast_k_week_ahead( + x_latest, + outcome = "case_rate", + predictors = c("case_rate", "percent_cli"), + ahead = .x, + engine = linear_reg() + )), + map(aheads, ~ forecast_k_week_ahead( + x_latest, + outcome = "case_rate", + predictors = c("case_rate", "percent_cli"), + ahead = .x, + engine = rand_forest(mode = "regression") + )) +) +pivot_quantiles_wider(fc_.pred_distn) ``` Here, `arx_forecaster()` does all the heavy lifting. It creates leads of the @@ -235,11 +230,11 @@ can_latest <- epix_as_of(can, max_version = max(can$DT$version)) can_fc <- bind_rows( map( aheads, - ~ k_week_ahead(can_latest, "cr_7dav", "cr_7dav", .x, linear_reg()) + ~ forecast_k_week_ahead(can_latest, "cr_7dav", "cr_7dav", .x, linear_reg()) ) %>% list_rbind(), map( aheads, - ~ k_week_ahead( + ~ forecast_k_week_ahead( can_latest, "cr_7dav", "cr_7dav", .x, boost_tree(mode = "regression", trees = 20) ) diff --git a/vignettes/arx-classifier.Rmd b/vignettes/arx-classifier.Rmd index ae1641cce..b2a2bbf8e 100644 --- a/vignettes/arx-classifier.Rmd +++ b/vignettes/arx-classifier.Rmd @@ -50,10 +50,7 @@ jhu <- case_death_rate_subset %>% geo_value %in% c("ca", "fl", "tx", "ny", "nj") ) -out <- arx_classifier(jhu, - outcome = "case_rate", - predictors = "case_rate" -) +out <- arx_classifier(jhu, outcome = "case_rate", predictors = "case_rate") out$predictions ``` @@ -93,7 +90,8 @@ relying on the default of 0.25. We can do this by passing 0.5 to the `breaks` argument in `arx_class_args_list()` as follows: ```{r} -out_break_0.5 <- arx_classifier(jhu, +out_break_0.5 <- arx_classifier( + jhu, outcome = "case_rate", predictors = "case_rate", args_list = arx_class_args_list( @@ -142,8 +140,8 @@ the present? To answer this question, we can create a predictive model for upswings and downswings of case rates rather than the raw case rates themselves. For this situation, our target is to predict whether there is an increase in case rates -or not. Following -[McDonald, Bien, Green, Hu, et al.(2021)](https://www.pnas.org/doi/full/10.1073/pnas.2111453118), +or not. Following +[McDonald, Bien, Green, Hu, et al.(2021)](https://www.pnas.org/doi/full/10.1073/pnas.2111453118), we look at the relative change between $Y_{l,t}$ and $Y_{l, t+a}$, where the former is the case rate at location $l$ at time $t$ and the latter is the rate for that location at @@ -152,7 +150,7 @@ with two classes $$\begin{align} Z_{l,t} = \left\{\begin{matrix} -\text{up,} & \text{if } Y_{l,t}^\Delta > 0.25\\ +\text{up,} & \text{if } Y_{l,t}^\Delta > 0.25\\ \text{not up,} & \text{otherwise} \end{matrix}\right. \end{align}$$ @@ -166,7 +164,7 @@ $$\begin{align} \pi_{\text{up}}(x) &= Pr(Z_{l, t} = \text{up}|x) = \frac{e^{g_{\text{up}}(x)}}{1 + e^{g_{\text{up}}(x)}}, \\ \pi_{\text{not up}}(x)&= Pr(Z_{l, t} = \text{not up}|x) = 1 - Pr(Z_{l, t} = \text{up}|x) = \frac{1}{1 + e^{g_{\text{up}}(x)}} \end{align}$$ -where +where $$ g_{\text{up}}(x) = \log\left ( \frac{\Pr(Z_{l, t} = \text{up} \vert x)}{\Pr(Z_{l, t} = \text{not up} \vert x)} \right ) = \beta_{10} + \beta_{11}Y_{l,t}^\Delta + \beta_{12}Y_{l,t-7}^\Delta + \beta_{13}Y_{l,t-14}^\Delta. @@ -223,7 +221,7 @@ require access to the training data. The other optional arguments for controlling the growth rate calculation (that can be inputted as `additional_gr_args`) can be found in the documentation for -`epiprocess::growth_rate()` and the related +`epiprocess::growth_rate()` and the related `vignette("growth_rate", package = "epiprocess")`. ### Visualizing the results @@ -280,4 +278,4 @@ to start with using the built-in classifier for ostensibly simple projects and begin to implement your own when the modelling project takes a complicated turn. To get some practice on coding up a classifier by hand, consider translating this binary classification model example to an `epi_workflow`, akin to that in -`vignette("preprocessing-and-models")`. +`vignette("preprocessing-and-models")`. diff --git a/vignettes/epipredict.Rmd b/vignettes/epipredict.Rmd index 7e24b04c6..1925de2fb 100644 --- a/vignettes/epipredict.Rmd +++ b/vignettes/epipredict.Rmd @@ -110,8 +110,6 @@ the *same set* of `geo_value`'s and `time_value`'s could actually be different. For more details, see [`{epiprocess}`](https://cmu-delphi.github.io/epiprocess/articles/epiprocess.html). - - ## Why doesn't this package already exist? As described above: @@ -121,7 +119,7 @@ preprocessing, training, and prediction, bound together, through a package calle `{workflows}`. We built `{epipredict}` on top of that setup. In this way, you CAN use almost everything they provide. -* However, `{workflows}` doesn't do postprocessing. And nothing in the -verse +* However, `{workflows}` doesn't do postprocessing. And nothing in the -verse handles _panel data_. * The tidy-team doesn't have plans to do either of these things. (We checked). @@ -131,7 +129,7 @@ handles _panel data_. etc.[^2] Our group has not prioritized these sorts of models for epidemic forecasting, but one could also integrate these methods into our framework. -[^2]: These are [`{timetk}`](https://business-science.github.io/timetk/index.html) +[^2]: These are [`{timetk}`](https://business-science.github.io/timetk/index.html) and [`{modeltime}`](https://business-science.github.io/timetk/index.html). There are *lots* of useful methods there than can be used to do fairly complex machine learning methodology, though not directly for panel data and not for direct @@ -327,6 +325,7 @@ the `time_value`, `geo_value`, and any additional keys so that these are availab when necessary. The `epi_recipe` from `out_gb` can be extracted from the result: + ```{r} extract_recipe(out_gb$epi_workflow) ``` @@ -441,7 +440,7 @@ But ideally, a user could create their own forecasters by building up the components we provide. In other vignettes, we try to walk through some of these customizations. -To illustrate everything above, here is (roughly) the code for the +To illustrate everything above, here is (roughly) the code for the `flatline_forecaster()` applied to the `case_rate`. ```{r} From b3e3189da970ac341fa0294b02cf49ee3899b0ac Mon Sep 17 00:00:00 2001 From: Dmitry Shemetov Date: Fri, 27 Sep 2024 14:16:28 -0700 Subject: [PATCH 2/3] styler: style --- tests/testthat/test-layer_add_forecast_date.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/testthat/test-layer_add_forecast_date.R b/tests/testthat/test-layer_add_forecast_date.R index 915804d6b..6d0e637c8 100644 --- a/tests/testthat/test-layer_add_forecast_date.R +++ b/tests/testthat/test-layer_add_forecast_date.R @@ -94,7 +94,7 @@ test_that("forecast date works for daily", { as.data.frame() %>% mutate(time_value = as.POSIXlt(time_value)$year + 1900L) %>% group_by(geo_value, time_value) %>% - summarize(case_rate = mean(case_rate), death_rate = mean(death_rate), .groups="drop") %>% + summarize(case_rate = mean(case_rate), death_rate = mean(death_rate), .groups = "drop") %>% as_epi_df() expect_error(predict(wf1, latest_yearly)) From 374cb2f1673529e13571cb6e1fe11a920f385d02 Mon Sep 17 00:00:00 2001 From: Dmitry Shemetov Date: Fri, 27 Sep 2024 15:52:18 -0700 Subject: [PATCH 3/3] doc: fix vignettes --- R/layer_add_forecast_date.R | 1 + R/layer_add_target_date.R | 1 + 2 files changed, 2 insertions(+) diff --git a/R/layer_add_forecast_date.R b/R/layer_add_forecast_date.R index 02395f960..3d5ea010b 100644 --- a/R/layer_add_forecast_date.R +++ b/R/layer_add_forecast_date.R @@ -104,6 +104,7 @@ slather.layer_add_forecast_date <- function(object, components, workflow, workflows::extract_preprocessor(workflow)$template, "metadata" )$time_type if (expected_time_type == "week") expected_time_type <- "day" + if (expected_time_type == "integer") expected_time_type <- "year" validate_date( forecast_date, expected_time_type, call = rlang::expr(layer_add_forecast_date()) diff --git a/R/layer_add_target_date.R b/R/layer_add_target_date.R index 9176fb593..094ec8501 100644 --- a/R/layer_add_target_date.R +++ b/R/layer_add_target_date.R @@ -90,6 +90,7 @@ slather.layer_add_target_date <- function(object, components, workflow, workflows::extract_preprocessor(workflow)$template, "metadata" )$time_type if (expected_time_type == "week") expected_time_type <- "day" + if (expected_time_type == "integer") expected_time_type <- "year" if (!is.null(object$target_date)) { target_date <- object$target_date