Comparing the CEA's "cubic model" to quadratic and quartic models

comparing-polynomial-models-covid.R

library(tidyverse)
library(broom)

US <- read_csv("https://raw.githubusercontent.com/nytimes/covid-19-data/master/us.csv") %>%
  mutate(new_deaths = deaths - lag(deaths)) %>%
  filter(date >= "2020-02-26")

models <- tibble(degrees = 2:4) %>%
  mutate(model = map(degrees, ~ lm(log(new_deaths + 1) ~ poly(date, .), data = US)))

# Predict the future
US_future <- US %>%
  bind_rows(tibble(date = Sys.Date() + 1:20))

models %>%
  mutate(augmented = map(model, augment, newdata = US_future)) %>%
  unnest(augmented) %>%
  mutate(type = c("Linear", "Quadratic", "Cubic", "Quartic")[degrees]) %>%
  mutate(type = fct_reorder(type, degrees)) %>%
  ggplot(aes(date, new_deaths)) +
  geom_line() +
  geom_line(aes(y = exp(.fitted) - 1, color = type), lty = 2) +
  labs(title = "How a polynomial model's predictions change based on # of degrees",
       subtitle = "Fit to log(deaths + 1). Source of death counts: NY Times",
       y = "Daily US deaths",
       x = "") +
  facet_wrap(~ type, scales = "free", ncol = 1) +
  scale_y_continuous(labels = comma_format()) +
  theme(legend.position = "none")