covid_animate.R

remotes::install_github("kjhealy/covdata")
library(covdata)
library(dplyr)
library(ggplot2)
library(ggrepel)
library(tidyr)

covus_wide <- covus %>% 
  select(date, state, measure, count) %>% 
  pivot_wider(id_cols = c(date, state), 
              names_from = measure, 
              values_from = count)

##### define df_covus #####
excluding <- c("MP","VI", "AS", "GU")

data("uspop")
uspop$state_abbr[uspop$state == "District of Columbia"] <- "DC"

df_covus <- covus_wide %>% 
  filter(!state %in% excluding) %>% 
  left_join(y = uspop %>% 
              filter(sex == "Both Sexes", 
                     hisp_label == "Total") %>% 
              select(state_abbr, pop) %>% 
              # add missing PR population
              rbind(data.frame(state_abbr = "PR",
                               pop=3193694)), # census.gov
            by = c("state" = "state_abbr")) %>% 
  group_by(state) %>% 
  mutate(daily_cases = positive - lead(positive)) %>% 
  mutate(daily_cases_capita = daily_cases / pop)

library(zoo)
library(statebins)
library(gganimate)

df_covus2 <- df_covus %>% 
  select(date, state, daily_cases_capita) %>% 
  group_by(state) %>% 
  mutate(weekly_avg_percapita = 
           rollapply(daily_cases_capita,
                     7, mean, align = "left", fill = NA),
         weekly_avg_permille = weekly_avg_percapita * 1000)

df_covus3 <- df_covus2
df_covus3$month <-  format(df_covus3$date,"%B")

df_covus3 <- df_covus3 %>% 
  group_by(state, month) %>% 
  summarise(monthly_avg_permille = 
              mean(weekly_avg_permille, 
                   na.rm = TRUE))

# add missing rows
missing_df <- 
  data.frame(state = df_covus3$state %>% 
               unique(), 
             month = c(rep("January", 52),
                     rep("February", 52)), 
             monthly_avg_permille = 0) %>% 
  filter(!state == "WA")

df_covus3 <- rbind(df_covus3, missing_df)

df_covus3$month <-  df_covus3$month %>% 
  factor(levels = c("January","February",
                    "March","April",
                    "May","June",
                    "July"))

plot <- df_covus3 %>% 
  ggplot(aes(state = state, 
             fill = monthly_avg_permille)) + 
  geom_statebins(na.rm = FALSE) + 
  coord_equal() + 
  scale_fill_viridis_c(option = "magma") + 
  theme_void() +
  theme(legend.title = element_blank(),
        legend.position = c(0.12, 0.8), 
        legend.justification = c(0, 0),
        legend.direction = "horizontal") + 
  transition_states(month,
                    transition_length = 10,
                    state_length = 2,
                    wrap = FALSE) +
  labs(title = '{closest_state}, 2020',
       subtitle = "average daily new COVID-19 cases per thousand residents",
       caption = "gist.github.com/jmclawson") +
  theme(plot.title = element_text(hjust = 0.5,
                                  size=14, 
                                  face="bold"),
        plot.subtitle = element_text(hjust = 0.5))

animate(plot, start_pause = 5, end_pause = 5)

anim_save("COVID-19_by_month.gif")

##### Louisiana Specific, by Parish #####

covid_la <- nytcovcounty %>% 
  filter(state == "Louisiana") %>% 
  group_by(county) %>% 
  mutate(daily_cases = cases - lag(cases, 
                                   order_by = date, 
                                   default = 0),
         daily_deaths = deaths - lag(deaths, 
                                     order_by = date,
                                     default = 0)) %>% 
  mutate(daily_cases = 
           rollapply(daily_cases,
                     5, mean, align = "right", fill = NA))

covid_la$daily_cases[covid_la$daily_cases<0] <- 0

covid_la$subregion <- louisiana$county %>% 
  tolower() %>% 
  gsub(pattern = "\\.", replacement = "", x = .) %>% 
  gsub(pattern = "lasalle", replacement = "la salle", x = .)

map_louisiana <- 
  map_data("county") %>% 
  filter(region=="louisiana")

library(maps)

map_louisiana$polyname <- 
  paste(map_louisiana$region, 
        map_louisiana$subregion, 
        sep=",")

map_louisiana <- 
  map_louisiana %>% 
  left_join(maps::county.fips, "polyname")

map_louisiana$fips[map_louisiana$subregion=="st martin"] <- 22099

map_louisiana$fips <- as.character(map_louisiana$fips)
covid_la$fips <- as.character(covid_la$fips)

covid_la <- covid_la %>% 
  filter(date >= as.Date("2020-06-01"),
         !county == "Unknown")

# map_louisiana <- map_louisiana %>% 
#   select(long, lat, fips)
# 
# colnames(map_louisiana) <- 
#   c("x", "y", "id")

covid_la1 <- left_join(map_louisiana, covid_la, by="fips")

library(transformr)

la_plot <- ggplot(covid_la1, 
       aes(x = long, y=lat,
           group=group,
           fill=daily_cases)) +
  geom_polygon(color="gray") +
  coord_map(projection = "albers", lat0=39, lat1=45) + 
  theme_void() +
  scale_fill_viridis_c(option = "magma") + 
  theme(legend.title = element_blank(),
        legend.position = c(0.53, 0.52), 
        legend.justification = c(0, 0),
        legend.direction = "horizontal") +
  transition_states(date,
                    transition_length = 2,
                    state_length = 2,
                    wrap = FALSE) +
  labs(title = '{closest_state}',
       subtitle = "5-day average of daily new cases per parish",
       caption = "gist.github.com/jmclawson  ") +
  theme(plot.title = element_text(hjust = 0.5,
                                  size=14, 
                                  face="bold"),
        plot.subtitle = element_text(hjust = 0.5))

animate(la_plot, start_pause = 5, end_pause = 15)

anim_save("COVID-19_louisiana.gif")