@@ -14,10 +14,10 @@ $$\\[.4in]$$
1414
1515``` {r echo=FALSE}
1616knitr::opts_chunk$set(
17- fig.align = "center",
18- message = FALSE,
19- warning = FALSE,
20- cache = FALSE
17+ fig.align = "center",
18+ message = FALSE,
19+ warning = FALSE,
20+ cache = FALSE
2121)
2222knitr::opts_knit$set(root.dir = here::here())
2323ggplot2::theme_set(ggplot2::theme_bw())
@@ -40,24 +40,24 @@ library(ggplot2)
4040library(plotly)
4141
4242if (params$disease == "flu") {
43- epi_data <- readr::read_csv("https://data.cdc.gov/resource/ua7e-t2fy.csv?$limit=20000&$select=weekendingdate,jurisdiction,totalconfflunewadm", show_col_types = FALSE) %>%
44- rename(time_value = weekendingdate, value = totalconfflunewadm, geo_value = jurisdiction)
43+ epi_data <- readr::read_csv("https://data.cdc.gov/resource/ua7e-t2fy.csv?$limit=20000&$select=weekendingdate,jurisdiction,totalconfflunewadm", show_col_types = FALSE) %>%
44+ rename(time_value = weekendingdate, value = totalconfflunewadm, geo_value = jurisdiction)
4545} else if (params$disease == "covid") {
46- epi_data <- readr::read_csv("https://data.cdc.gov/resource/ua7e-t2fy.csv?$limit=20000&$select=weekendingdate,jurisdiction,totalconfc19newadm", show_col_types = FALSE) %>%
47- rename(time_value = weekendingdate, value = totalconfc19newadm, geo_value = jurisdiction)
46+ epi_data <- readr::read_csv("https://data.cdc.gov/resource/ua7e-t2fy.csv?$limit=20000&$select=weekendingdate,jurisdiction,totalconfc19newadm", show_col_types = FALSE) %>%
47+ rename(time_value = weekendingdate, value = totalconfc19newadm, geo_value = jurisdiction)
4848} else {
49- stop("Invalid disease")
49+ stop("Invalid disease")
5050}
5151today <- Sys.Date()
5252epi_data <- epi_data %>%
53- filter(geo_value %nin% c("AS", "USA", "VI", "PR", "MP", "GU")) %>%
54- mutate(time_value = as.Date(time_value), geo_value = tolower(geo_value)) %>%
55- mutate(
56- target_end_date = time_value,
57- value = ifelse(value == 0, NA, value)
58- ) %>%
59- arrange(geo_value, time_value) %>%
60- as_epi_df(as_of = today)
53+ filter(geo_value %nin% c("AS", "USA", "VI", "PR", "MP", "GU")) %>%
54+ mutate(time_value = as.Date(time_value), geo_value = tolower(geo_value)) %>%
55+ mutate(
56+ target_end_date = time_value,
57+ value = ifelse(value == 0, NA, value)
58+ ) %>%
59+ arrange(geo_value, time_value) %>%
60+ as_epi_df(as_of = today)
6161
6262aheads <- -1:3
6363forecast_date <- today
@@ -69,18 +69,18 @@ subset_geos <- c("ca", "ny", "fl", "tx", "pa")
6969
7070``` {r}
7171quantile_forecast <- map(aheads, ~ forecaster_baseline_linear(epi_data, .x, log = FALSE)) %>%
72- bind_rows() %>%
73- filter(geo_value %in% subset_geos) %>%
74- mutate(forecaster = "baseline_linear")
72+ bind_rows() %>%
73+ filter(geo_value %in% subset_geos) %>%
74+ mutate(forecaster = "baseline_linear")
7575truth_data <- epi_data %>%
76- filter(geo_value %in% subset_geos, time_value > forecast_date - 365 * 2)
76+ filter(geo_value %in% subset_geos, time_value > forecast_date - 365 * 2)
7777train_fit <- quantile_forecast %>%
78- filter(is.na(quantile)) %>%
79- mutate(time_value = target_end_date)
78+ filter(is.na(quantile)) %>%
79+ mutate(time_value = target_end_date)
8080quantile_forecast <- quantile_forecast %>% filter(!is.na(quantile))
8181train_data <- epi_data %>%
82- filter(geo_value %in% subset_geos, time_value > forecast_date - 30) %>%
83- summarize(start = min(time_value), stop = max(time_value))
82+ filter(geo_value %in% subset_geos, time_value > forecast_date - 30) %>%
83+ summarize(start = min(time_value), stop = max(time_value))
8484
8585# Bad forecast filters
8686filter_geos <- filter_forecast_geos(quantile_forecast, truth_data)
@@ -92,7 +92,7 @@ Filter geo recommendations: `r paste(filter_geos, collapse = ", ")`.
9292
9393``` {r}
9494p <- plot_forecasts(quantile_forecast, forecast_date, truth_data = truth_data, relevant_period = train_data) +
95- geom_point(data = train_fit, mapping = aes(x = time_value, y = value), size = 0.25)
95+ geom_point(data = train_fit, mapping = aes(x = time_value, y = value), size = 0.25)
9696
9797ggplotly(p, tooltip = "text", height = 5000, width = 1700)
9898```
@@ -101,18 +101,18 @@ ggplotly(p, tooltip = "text", height = 5000, width = 1700)
101101
102102``` {r}
103103quantile_forecast <- map(aheads, ~ forecaster_baseline_linear(epi_data, .x, log = TRUE)) %>%
104- bind_rows() %>%
105- filter(geo_value %in% subset_geos) %>%
106- mutate(forecaster = "baseline_linear")
104+ bind_rows() %>%
105+ filter(geo_value %in% subset_geos) %>%
106+ mutate(forecaster = "baseline_linear")
107107truth_data <- epi_data %>%
108- filter(geo_value %in% subset_geos, time_value > forecast_date - 365 * 2)
108+ filter(geo_value %in% subset_geos, time_value > forecast_date - 365 * 2)
109109train_fit <- quantile_forecast %>%
110- filter(is.na(quantile)) %>%
111- mutate(time_value = target_end_date)
110+ filter(is.na(quantile)) %>%
111+ mutate(time_value = target_end_date)
112112quantile_forecast <- quantile_forecast %>% filter(!is.na(quantile))
113113train_data <- epi_data %>%
114- filter(geo_value %in% subset_geos, time_value > forecast_date - 30) %>%
115- summarize(start = min(time_value), stop = max(time_value))
114+ filter(geo_value %in% subset_geos, time_value > forecast_date - 30) %>%
115+ summarize(start = min(time_value), stop = max(time_value))
116116
117117filter_geos <- filter_forecast_geos(quantile_forecast, truth_data)
118118```
@@ -123,7 +123,7 @@ Filter geo recommendations: `r paste(filter_geos, collapse = ", ")`.
123123
124124``` {r}
125125p <- plot_forecasts(quantile_forecast, forecast_date, truth_data = truth_data, relevant_period = train_data) +
126- geom_point(data = train_fit, mapping = aes(x = time_value, y = value), size = 0.25)
126+ geom_point(data = train_fit, mapping = aes(x = time_value, y = value), size = 0.25)
127127
128128ggplotly(p, tooltip = "text", height = 5000, width = 1700)
129129```
0 commit comments