1 Back to Outline

This section includes information on drug assistance program diagnostics. You can go back to table of contents if you would like.

Other related documents are:

if (is.null(params$sim_data)) {
  sim_loc <- "../../EpiModel/AE/sim_epimodel3/sim_on_2021-06-13_at_2033.rds"
}else{
  sim_loc <- params$sim_data
}
attr_names <- EpiModelWHAMPDX::attr_names
attr_names$insurance <- NULL
suppressPackageStartupMessages(library(tidyverse))
suppressPackageStartupMessages(library(magrittr))
suppressPackageStartupMessages(library(cowplot))
suppressPackageStartupMessages(library(ggrepel))
library(here)
options(dplyr.summarise.inform = FALSE)
library(EpiModelWHAMPDX)
sim_dat <- readRDS(sim_loc)
knitr::opts_chunk$set(echo = TRUE, message = FALSE,
                      warning = FALSE, fig.width = 8)
cur_targs <- readRDS(here("Data", "EpiModelSims", "WHAMP.dx.targs.rds"))

2 DAP Outcomes

prop.dap <- list(
  adap.prop = list(name = "prop.adap",
                   sec_title = "Proportion of ART users in ADAP",
                   plot_name = "Proportion in ADAP",
                   plot_cap = "Among ART users.  Target source: WA DOH",
                   plot_ylab = "Proportion",
                   plt_type = "line",
                   vars = c("num.adap.", "tx.i.num."),
                   sum_fun = function(x, y) { x / y }
  ),
  pdap.prop = list(name = "prop.pdap",
                   sec_title = "Proportion of those on PrEP",
                   plot_name = "Proportion on PDAP",
                   plot_cap =
                     "Among individuals on PrEP, Target source: WA DOH",
                   plot_ylab = "Proportion",
                   plt_type = "line",
                   vars = c("num.pdap.", "neg.prep.num."),
                   sum_fun = function(x, y) { x / y }
  )
)
print_many_plots(prop.dap, sim_dat, targ_df = cur_targs, num_hash = 2)

2.1 Proportion of ART users in ADAP

2.1.1 ovr

2.1.2 race

2.1.3 region

2.1.4 age.grp

2.1.5 snap5

2.2 Proportion of those on PrEP

2.2.1 ovr

2.2.2 race

2.2.3 region

2.2.4 age.grp

2.2.5 snap5

2.3 Yearly Costs

if ("num" %in% names(sim_dat$epi)) {
  epi.dat <- sim_dat$epi
  num.sims <- ncol(epi.dat$num)
  if (num.sims == 1) {
    num.steps <- 1
  }else{
    num.steps <- length(epi.dat$num)
  }
}else{
  epi.dat <- sim_dat$epi[[1]]
  num.sims <- 1
  num.steps <- length(epi.dat$num)
}

if (num.sims > 1) {
  pop.sizes <- apply(epi.dat$num, 1, mean)
}else{
  pop.sizes <- epi.dat$num
  if (class(pop.sizes) %in% "data.frame") {
    pop.sizes <- pop.sizes[, 1, drop = TRUE]
  }
}

sim_pop_sizes <- data.frame(pop.size = pop.sizes)
sim_start_year <- round(
  sim_dat$control$year.start - sim_dat$control$start / 52)
sim_pop_sizes$year <- (1:nrow(sim_pop_sizes) + 1) / 52 + sim_start_year
sim_pop_sizes <- sim_pop_sizes %>% filter(year %% 1 == 0.0000000)
sim_pop_sizes$tr.pop <- NA_real_
sim_pop_sizes$year <- round(sim_pop_sizes$year)

for (yrv in 2014:2019) {
  sim_pop_sizes[which(sim_pop_sizes$year == yrv), "tr.pop"]  <-
    eval(parse(text = paste0("sum(WApopdata::msm.pop.totals_", 
               yrv, "$pop.all$num)")))
}

sim_pop_sizes$tr.pop[sim_pop_sizes$year < 2014] <- 
  sim_pop_sizes$tr.pop[sim_pop_sizes$year == 2014]

sim_pop_sizes$tr.pop[sim_pop_sizes$year > 2019] <- 
  sim_pop_sizes$tr.pop[sim_pop_sizes$year == 2019] * 
  WApopdata::growth_rate ** 
  (sim_pop_sizes$year[sim_pop_sizes$year > 2019] - 2019)

mult_df <- sim_pop_sizes %>% mutate(
  multipl = tr.pop / pop.size
)

## TODO: Make this function with multiple simulations
yrl_cost <- map(sim_dat$whamp, function(x) {
  x$adap_cost_annual_cate %>% 
  mutate(year = year + sim_start_year) %>% 
  group_by(year) %>% summarise(cost = sum(cost))
})
all_costs <- reduce(yrl_cost, left_join, by = "year") 
fin_cost <- data.frame(
  "year" = all_costs$year, 
  "cost" = apply(all_costs %>% select(-year), 1, mean, na.rm = TRUE)
  )
fin_cost <- left_join(fin_cost, mult_df, by = "year") %>%
  mutate(true_cost = round(cost * multipl / 1000000, 1))

yrl_cost <- map(sim_dat$whamp, function(x) {
  x$adap_cost_annual_cate %>% 
  mutate(year = year + sim_start_year) %>% 
  group_by(year) %>% summarise(cost = sum(cost))
})

yrl_cost_pdap <- map(sim_dat$whamp, function(x) {
  x$pdap_cost_annual_cate %>% 
  mutate(year = year + sim_start_year) %>% 
  group_by(year) %>% summarise(cost = sum(cost))})
all_costs_pdap <- reduce(yrl_cost_pdap, left_join, by = "year") 
fin_cost_pdap <- data.frame(
  "year" = all_costs_pdap$year, 
  "cost" = apply(all_costs_pdap %>% select(-year), 1, mean, na.rm = TRUE)
  )
fin_cost_pdap <- left_join(fin_cost_pdap, mult_df, by = "year") %>%
  mutate(true_cost = round(cost * multipl / 1000000, 1))

bind_rows(bind_cols(fin_cost, program = "ADAP"),
          bind_cols(fin_cost_pdap, program = "PDAP")) %>%
  select(year, true_cost, program) %>% 
  pivot_wider(names_from = "program", values_from = true_cost) %>%
  filter(!(is.na(ADAP) & is.na(PDAP))) %>%
  gt::gt() %>% gt::tab_header(
    title = "Cost of drug assitance programs in millions",
    subtitle = gt::html("Based on simulated costs scaled up to <br> match the size of Washington state")
  ) %>% gt::cols_label(year = "Year")
Cost of drug assitance programs in millions
Based on simulated costs scaled up to
match the size of Washington state
Year ADAP PDAP
2015 8.3 NA
2016 19.4 NA
2017 29.4 NA
2018 47.7 2.2
2019 51.9 1.4
2020 53.5 2.0
2021 50.9 2.1
2022 50.6 2.4
2023 54.9 2.6
2024 52.1 2.9
2025 51.4 2.9
2026 54.8 3.1
2027 54.1 3.3
2028 53.9 3.6
2029 52.9 3.8
2030 55.1 4.0
2031 54.4 4.2
2032 51.0 4.2
2033 50.4 4.3