1 Model Validation

1.1 Comparison of Hazards

### Hazards based on Model
pred_dat <- readRDS("../../Data/Intermediate/coxdata.RDS")
pred_dat <- pred_dat %>% filter(!(never_tested == "T" & age > 45))
covardf <- with(pred_dat, expand.grid(
  age.young = unique(age.young),
  race2 = unique(race2),
  region.ewa = unique(region.ewa),
  snap.grp3 = unique(snap.grp3)))
# define the survival object in weeks
sdfwk <- with(pred_dat, Surv(mlt*4, tested))
# fit model
fit_gomp_scalewk <- flexsurvreg(sdfwk ~
                                  age.young +
                                  race2 +
                                  region.ewa +
                                  snap.grp3,
                                anc = list(shape =~ age.young),
                                data = pred_dat,
                                weights = ego.wawt,
                                dist = "gompertz")

haz_gompwk <- summary(fit_gomp_scalewk, newdata = covardf,
                      type = "hazard", tidy = TRUE)



pred_haz <- haz_gompwk %>%
  select(age.young, race = race2, snap.grp3, region.ewa, wlt = time,
         prop = est) %>%
  mutate(num = NA, num_tested = NA, type = "param target")

### Hazards based on the simulated data
sim_dat <- readRDS("../../EpiModel/AE/sim_epimodel3/sim_on_2021-03-11_at_2033.rds")
haz_check <- sim_dat$temp[[1]]$check_haz %>% mutate(prop = num_tested / num,
                                           type = "sim")

### Combing Data
all_dat <- bind_rows(pred_haz, haz_check)

1.2 Observed hazards within subgroups

1.2.1 Age.Young

all_dat %>% ggplot(aes(x = wlt, y = prop, color = age.young,
                       linetype = type)) + geom_smooth(se = FALSE) +
  xlab("Weeks Since Last Test") + 
  facet_grid(race ~ snap.grp3) + coord_cartesian(xlim = c(0, 500)) +
  guides(linetype = guide_legend(override.aes = list(color = "black")))

1.2.2 Race

all_dat %>% ggplot(aes(x = wlt, y = prop, color = race,
                       linetype = type)) + geom_smooth(se = FALSE) +
  xlab("Weeks Since Last Test") + 
  facet_grid(age.young ~ snap.grp3) + coord_cartesian(xlim = c(0, 500)) +
  guides(linetype = guide_legend(override.aes = list(color = "black")))

1.2.3 Region

all_dat %>% ggplot(aes(x = wlt, y = prop, color = region.ewa,
                       linetype = type)) + geom_smooth(se = FALSE) +
  xlab("Weeks Since Last Test") + 
  facet_grid(race ~ snap.grp3) + coord_cartesian(xlim = c(0, 500)) +
  guides(linetype = guide_legend(override.aes = list(color = "black")))

1.2.4 Snap.grp3

all_dat %>% ggplot(aes(x = wlt, y = prop, color = snap.grp3,
                       linetype = type)) + geom_smooth(se = FALSE) +
  xlab("Weeks Since Last Test") + 
  facet_grid(race ~ age.young) + coord_cartesian(xlim = c(0, 500)) +
  guides(linetype = guide_legend(override.aes = list(color = "black")))

1.3 Time since last test (in weeks) for individuals on PrEP:

We expect that all individuals on PrEP will have tested last, no more than 13 weeks after their last test.

all_elig <- which(sim_dat$attr[[1]]$ever.part == 1 &
                      (sim_dat$attr[[1]]$diag.status == 0 |
                       is.na(sim_dat$attr[[1]]$diag.status)))

prep_elig <- all_elig[which(sim_dat$attr[[1]]$prepStat[all_elig] == 1)]
time_since <- sim_dat$control$nsteps - sim_dat$attr[[1]]$last.neg.test

hist(time_since[prep_elig], main = "Time Since Last Test", xlab = "Weeks")

rolling_year_sum <- function(x, years = 1, avg = FALSE){
  val <- data.table::frollsum(fill(data.frame(x), "x")$x, years * 52)
  if (avg) {val <- val / (years * 52)}
  return(val)
}

aids_testing <- list(
  gen.tests = list(name = "aids.test.rate",
                   sec_title = "Rate of aids testing",
                   plot_name = "Testing Rate for HIV",
                   plot_cap = "Among individuals with AIDS",
                   # plot_cap = "Four year average of AIDS tests / Number eligable for AIDS testing",
                   plot_ylab = "Rate",
                   plt_type = "line",
                   vars = c("aidsrel.tests.", "aids.tst.elig."),
                   sum_fun = function(x, y) {
                     rolling_year_sum(x, 4) / rolling_year_sum(y, 4)
                   }
  )
)
aids.test.int <- sim_dat$param$vl.aids.int/2
ratesAIDS <- 1/aids.test.int
whamp_targs <- readRDS("../../Data/EpiModelSims/WHAMP.dx.targs.rds")
targ_df <- bind_rows(whamp_targs,
                     data.frame(subcat = "ovr", low = ratesAIDS * 0.97,
                                high = ratesAIDS * 1.03,
                                targ = ratesAIDS,
                                measure = "aids.test.rate", cat = "ovr",
                                sub_cat = "Overall"))

print_many_plots(aids_testing, 2, targ_df = targ_df)

1.4 Rate of aids testing

1.4.1 ovr

1.4.2 race

1.4.3 region

1.4.4 age.grp

2 Data Validation

2.1 Non-Prep tests (in last two years)

testing_try <- list(
  nonprep.tests = list(name = "gen.tests",
                   plot_name = "Number of tests", 
                   plot_cap = "Non-PrEP related tests", 
                   plot_ylab = "Two year rolling total",
                   sec_title = "Number", 
                   plt_type = "line", 
                   vars = c("tot.tests.", "prep.tests."),
                   sum_fun = function(x, y) {
                     rolling_year_sum(x, years = 2) - 
                       rolling_year_sum(y, years = 2)
                   } 
  ),
  nonprep.tests.prop = list(
    name = "tests2yr.noprep",
    # plot_name = "Two year rolling total of non-prep related tests / \n number negative and not on prep",
    plot_name = "Rate of Testing",
    plot_ylab = "Two year rolling average",
    plot_ylab = "Non-PrEP related tests",
    sec_title = "Number per 2 person years", 
    plt_type = "line", 
    vars = c("tot.tests.", "prep.tests.", 
             "num.", "neg.prep.num."),
    sum_fun = function(x, y, z1, z2) {
      z3 <- z2
      z3[is.na(z3)] <- 0
      (rolling_year_sum(x, years = 2) - 
          rolling_year_sum(y, years = 2)) / 
        rolling_year_sum(z1 - z3, years = 2, avg = TRUE)
    } 
  )
)

print_many_plots(testing_try, targ_df = targ_df)

2.1.1 Number

2.1.1.1 ovr

2.1.1.2 race

2.1.1.3 region

2.1.1.4 age.grp

2.1.2 Number per 2 person years

2.1.2.1 ovr

2.1.2.2 race

2.1.2.3 region

2.1.2.4 age.grp

2.2 Time since last test among non-prep users (in months)

time_since_noprep <- list(
  mean.time = list(name = "mlt.noprep",
                   plot_name = "Mean Time", 
                   plot_ylab = "Months", 
                   plot_cap = "Non-PrEP related tests.", 
                   plt_type = "line", 
                   vars = "mean.t.since.last.dx.noprep.",
                   sum_fun = function(x){x/4}
  ),
  median.time = list(name = "mlt.noprep.med",
                   plot_name = "Median Time",
                   plot_ylab = "Months", 
                   plot_cap = "Non-PrEP related tests.", 
                   plt_type = "line", 
                   vars = "median.t.since.last.dx.noprep.",
                   sum_fun = function(x){x/4})
)

print_many_plots(time_since_noprep, targ_df = targ_df)

2.2.1 Mean Time

2.2.1.1 ovr

2.2.1.2 race

2.2.1.3 region

2.2.1.4 age.grp

2.2.2 Median Time

2.2.2.1 ovr

2.2.2.2 race

2.2.2.3 region

2.2.2.4 age.grp

2.3 New diagnoses (number in last year)

2.3.1 New (non-late) diagnoses

2.3.1.1 ovr

2.3.1.2 race

2.3.1.3 region

2.3.1.4 age.grp

2.3.2 Number Undiagnosed

2.3.2.1 ovr

2.3.2.2 race

2.3.2.3 region

2.3.2.4 age.grp

2.3.3 Late diagnoses

2.3.3.1 ovr

2.3.3.2 race

2.3.3.3 region

2.3.3.4 age.grp

2.4 New Diagnoses (rate)

2.4.1 New Dx Rate

2.4.1.1 ovr

2.4.1.2 race

2.4.1.3 region

2.4.1.4 age.grp

2.4.2 Undiagnosed fraction

2.4.2.1 ovr

2.4.2.2 race

2.4.2.3 region

2.4.2.4 age.grp

2.4.3 Late diagnosis fraction

2.4.3.1 ovr

2.4.3.2 race

2.4.3.3 region

2.4.3.4 age.grp