In cmhoove14/LEMMAABMv4: COVID ABM

knitr::opts_chunk$set(echo = TRUE)

knitr::opts_chunk$set(echo = FALSE)
knitr::opts_knit$set(root.dir = rprojroot::find_rstudio_root_file())

devtools::load_all()
library(mapview)

# Get sims file-----------------
sim_output <- readRDS("data/outputs/ABMv4_bta_calibrate_Dec_Start_n5_bta0.1-0.3_2021-01-31.rds")
  bta_sweeps <- rep(seq(0.1,0.3,length.out = 6), each = 5)
  nsims <- length(sim_output)

# Get observed datasets for comparison --------------
source("data/get/COVID_CA_get_latest.R")

# Census tracts and agents
SF_cts <- readRDS("data/processed/SF_cts_sf.rds")
sf_pop <- readRDS("data/processed/SF_agents_processed.rds")

Observed vs simulated hospitalizations

sim_hosps <- rbindlist(lapply(1:nsims, function(s){
  df <- as.data.frame(sim_output[[i]][["epi_curve"]]) %>% 
    filter(state == "Ih") %>% 
    mutate(ndays = row_number(),
           Date = as.Date(ref_date+ndays),
           iter = i,
           bta = bta_sweeps[i])

  return(df)

}))

sf_hosp %>% 
  filter(Date > ref_date) %>% 
  ggplot() +
    geom_col(aes(x = Date, y = HOSP_tot), 
             col = "darkblue", fill = "blue",
             alpha = 0.4) +
    geom_line(data = sim_hosps,
              aes(x = Date, y = N, col = as.factor(bta), group = iter),
              alpha = 0.7) +
    theme_classic() +
    labs(x = "Date", y = "Hospitalizations",
         title = "Hospitalizations sim compared to observed",
         col = expression(beta))

Observed vs simulated detected cases

sim_cases <- rbindlist(lapply(1:nsims, function(s){
  df <- as.data.frame(sim_output[[i]][["linelist_tests"]]) %>% 
    mutate(iter = i,
           bta = bta_sweeps[i])

  return(df)

}))

tests_sum_by_date <- 
  sim_cases %>% 
  group_by(iter, Date, bta) %>% 
  summarise(n_tests = n(),
            n_pos = sum(test_pos),
            n_Ih = sum(state == "Ih"),
            n_Im = sum(state %in% c("Im", "Imh")),
            n_Ipa = sum(state %in% c("Ip", "Ia")),
            per_pos = n_pos/n_tests)


# Sums of testing data
tests_sum_by_date %>% 
  ggplot() +
  geom_col(data = sf_test,
           aes(x = Date, y = pos),
           fill = "grey50", alpha = 0.5) +
  geom_line(aes(x = as.Date(Date), y = n_pos, col = factor(bta), group = iter)) +
  theme_classic() +
  labs(y = "Positive Tests",
       title = "Confirmed cases compared to observed")

tests_sum_by_date %>% 
  ggplot() +
  geom_col(data = sf_test,
           aes(x = Date, y = pct),
           fill = "grey50", alpha = 0.5) +
  geom_line(aes(x = as.Date(Date), y = per_pos, col = factor(bta), group = iter)) +
  theme_classic() +
  labs(title = "Percent positive compared to observed")

Census tract level distribution of cases

sf_pop_sum <- sf_pop %>% 
  group_by(geoid) %>% 
  summarise(Pop = n())

SF_pop_cts <- SF_cts %>% 
  left_join(sf_pop_sum, by = c("GEOID" = "geoid")) %>%