add conditional vax figs

exp/active-vac/fig/Makefile

@@ -68,6 +68,7 @@ ${OUTDIR}/everything_bw.png: ${INPUTS}
 ${OUTDIR}/avertvis.png: ${PLOTINS}
 ${OUTDIR}/cumulative_doses.png: ${RESDIR}/digest-doses.rds ${RESDIR}/digest-key.rds
 ${OUTDIR}/summary.png: ${RESDIR}/alt_eff.rds ${RESDIR}/digest-key.rds
+${OUTDIR}/con_summary.png: ${RESDIR}/alt_eff.rds ${RESDIR}/digest-key.rds
 
 ${OUTDIR}/%.png: %.R vis_support.rda | ${OUTDIR}
 	$(call R)
@@ -89,19 +90,15 @@ ${OUTDIR}/alt_ave_%.png: alt_ave_outcome.R vis_support.rda ${RESDIR}/alt_eff.rds
 ${RESDIR}/digest-key.rds | ${OUTDIR}
 	$(call R)
 
-${OUTDIR}/alt_eff_all.png: alt_eff_all.R vis_support.rda ${RESDIR}/alt_eff.rds \
+${OUTDIR}/%_eff_all.png: %_eff_all.R vis_support.rda ${RESDIR}/alt_eff.rds \
 ${RESDIR}/digest-key.rds ${RESDIR}/vocwindows.rds | ${OUTDIR}
 	$(call R)
 
 ${OUTDIR}/alt_ave_all.png: alt_ave_all.R vis_support.rda ${RESDIR}/alt_eff.rds \
 ${RESDIR}/digest-key.rds ${RESDIR}/vocwindows.rds | ${OUTDIR}
 	$(call R)
 
-${OUTDIR}/alt_ci_all.png: alt_ci_all.R vis_support.rda ${RESDIR}/alt_eff.rds \
-${RESDIR}/digest-key.rds ${RESDIR}/vocwindows.rds | ${OUTDIR}
-	$(call R)
-
-${OUTDIR}/alt_inc_all.png: alt_i_all.R vis_support.rda ${RESDIR}/alt_eff.rds \
+${OUTDIR}/%_ci_all.png: %_ci_all.R vis_support.rda ${RESDIR}/alt_eff.rds \
 ${RESDIR}/digest-key.rds ${RESDIR}/vocwindows.rds | ${OUTDIR}
 	$(call R)
 
@@ -121,6 +118,7 @@ altinc: $(addprefix ${OUTDIR}/,$(patsubst %,alt_inc_%.png,inf sev deaths))
 allfigs: $(patsubst %,${OUTDIR}/%.png,model_input validation_review everything)
 
 altall: $(addprefix ${OUTDIR}/,$(patsubst %,alt_%_all.png,eff ave ci inc mul) cumulative_doses.png everything.png everything_bw.png)
+conall: $(addprefix ${OUTDIR}/,$(patsubst %,conditioned_%_all.png,eff ci) con_summary.png)
 
 vis_support.rda: vis_support.R
 	$(call R)

exp/active-vac/fig/con_summary.R

@@ -0,0 +1,100 @@
+
+.pkgs <- c("data.table", "ggplot2", "cabputils")
+.pkgs |> sapply(require, character.only = TRUE) |> all() |> stopifnot()
+
+.args <- commandArgs(args = c(
+  file.path("fig", "vis_support.rda"),
+  file.path("fig", "process", "alt_eff.rds"),
+  file.path("fig", "process", "digest-key.rds"),
+  file.path("fig", "output", "con_summary.png")
+))
+
+load(.args[1])
+
+overdates <- as.Date(c("2021-05-27", "2021-11-26", "2022-03-07"))
+
+scn.dt <- readRDS(.args[3])[inf_con == TRUE][, .(
+  scenario, quar, alloc = fifelse(pas_vac, pas_alloc, act_alloc),
+  act_vac
+)]
+
+dt <- readRDS(.args[2])[
+  (outcome %in% c("inf", "deaths")) & (scenario %in% scn.dt$scenario)
+]
+
+q.dt <- dt |>
+  DT(,c("c.value", "c.averted") := .(cumsum(value), cumsum(averted)), by=setdiff(key(dt), "date")) |>
+  DT(date %in% overdates) |>
+  quantile(j=.(c.value, c.averted, c.effectiveness), sampleby = "realization", probs = qprobs(c(`90`=0.9)))
+
+
+plt.dt <- q.dt[scn.dt, on = .(scenario), nomatch = 0]
+ref.dt <- CJ(
+  outcome = c("inf", "deaths"),
+  date = overdates,
+  measure = "c.effectiveness",
+  qmed = 0, quar = FALSE, alloc = c("LIC", "MIC", "HIC", "USA"),
+  act_vac = "none"
+)
+plt.dt <- rbind(plt.dt, ref.dt, fill = TRUE)
+plt.dt[date == overdates[1], variant := "alpha"]
+plt.dt[date == overdates[2], variant := "delta"]
+plt.dt[date == overdates[3], variant := "omicron"]
+
+plt.dt$act_vac <- factor(plt.dt$act_vac, levels = c("ring", "none", "age", "risk"))
+
+scale_shape_quar <- rejig(
+  scale_shape_manual,
+  name = "Extra NPI", labels = c(nonpi="None", wquar = "Quarantine Contacts"),
+  values = c(wquar=21, nonpi=19),
+  guide = guide_legend(title.position = "top", title.hjust = 0.5, order = 1)
+)
+
+p <- ggplot(plt.dt[measure == "c.effectiveness"]) + aes(
+  x=variant, color = act_vac,
+  shape = c("nonpi","wquar")[quar+1]
+) +
+  geom_pointrange(
+    aes(y=qmed, ymin=q90l, ymax=q90h),
+    data = \(dt) dt[quar == FALSE],
+    position = position_dodge(width = 0.5),
+    size = 0.5, stroke = 0
+  ) +
+  geom_pointrange(
+    aes(y=qmed, ymin=q90l, ymax=q90h),
+    data = \(dt) dt[quar == TRUE],
+    position = position_dodge(width = 0.5), fill = "white",
+    size = 0.4, stroke = 0.4
+  ) +
+  facet_grid(
+    outcome ~ alloc, scales = "free_y", switch = "y",
+    labeller = labeller(outcome = c(inf = "Infections", deaths = "Deaths"))
+  ) +
+  coord_cartesian(clip = "off") +
+  theme_minimal() +
+  theme(
+    strip.placement = "outside", legend.position = "bottom",
+    #panel.spacing.y = unit(1.5, "line"), panel.spacing.x = unit(1, "line"),
+    legend.text = element_text(size = rel(.75)),
+    panel.grid.major.x = element_blank(),
+    panel.border = element_rect(fill = NA, color = "grey")
+  ) +
+  # scale_linetype_quar(
+  #   guide = guide_legend(title.position = "top", title.hjust = 0.5, order = 1)
+  # ) +
+  scale_shape_quar() +
+  scale_color_strategy(
+    breaks = c("ring", "none", "age", "risk"),
+    values = c(none = "black", ring = "#fb6502", risk = "#3b90db", age = "#209033")
+  ) +
+  scale_x_discrete(name = NULL) +
+  scale_y_continuous(name = "Cumulative Effectiveness After Each Variant")
+
+g <- ggplotGrob(p)
+id <- which(g$layout$name == "guide-box")
+g$layout[id, c("l","r")] <- c(1, ncol(g))
+grDevices::png(tail(.args,1), width = 9, height = 6, bg = "white", units = "in", res = 300)
+grid::grid.draw(g)
+grDevices::dev.off()
+
+#store(obj = g, args = .args, width = 9, height = 6, bg = "white")

exp/active-vac/fig/conditioned_ci_all.R

@@ -0,0 +1,76 @@
+
+.pkgs <- c("data.table", "ggplot2", "scales", "ggh4x", "cabputils", "geomtextpath")
+.pkgs |> sapply(require, character.only = TRUE) |> all() |> stopifnot()
+
+#' assumes R project at the experiment root level
+.args <- commandArgs(args = c(
+  file.path("fig", "vis_support.rda"),
+  file.path("fig", "process", c("alt_eff.rds", "digest-key.rds", "vocwindows.rds")),
+  file.path("fig", "output", "conditioned_ci_all.png")
+))
+
+load(.args[1])
+
+intfilter <- if (interactive()) expression(realization < 10) else expression(realization >= 0)
+
+#' comes key'd
+inc.dt <- readRDS(.args[2])[
+  eval(intfilter) & eval(outfilter)
+][
+  eval(datefilter)
+][, .(
+  scenario, realization, date, outcome, value, averted
+)]
+
+inc.dt[order(date),
+  c("c.value", "c.averted") := .(
+    cumsum(value), cumsum(averted)
+  ), by=.(scenario, realization, outcome)
+][, i.c.value := c.value + c.averted ]
+
+intscns <- inc.dt[, unique(scenario)]
+
+scn.dt <- readRDS(.args[3])
+
+takeover.wins <- readRDS(.args[4])
+
+intscn.dt <- scn.dt[scenario %in% intscns]
+# reconstructing reference scenarios
+refscn.dt <- scn.dt[quar == FALSE & pas_vac == TRUE & act_vac == "none"]
+
+incref.dt <- copy(inc.dt)[
+  intscn.dt, on=.(scenario)
+][(act_vac == "ring") & (quar == FALSE)][ # only need to go from one reference
+  refscn.dt, on =.(act_alloc = pas_alloc, inf_con, quar)
+][,.(
+  c.value = i.c.value
+), by=.(i.scenario, realization, date, outcome)
+][refscn.dt, on=.(i.scenario = scenario)]
+setnames(incref.dt, "i.scenario", "scenario")
+
+plt.dt <- setkeyv(
+  rbind(
+    inc.dt[, c(key(inc.dt), "c.value"), with = FALSE][intscn.dt, on=.(scenario)],
+    incref.dt
+  ),
+  union(key(inc.dt), colnames(scn.dt))
+)[inf_con == TRUE]
+
+rm(inc.dt)
+gc()
+
+plt.qs <- plt.prep(plt.dt, j = .(c.value))
+
+mm.ref <- plt.qs[,.(ymin = min(q90l), ymax = max(q90h)),by=.(outcome)]
+mm.ref[, ymin := ymin - .15*(ymax-ymin) ]
+mm.ref[, yspan := ymax - ymin ]
+tw <- takeover.wins[q == 0.5][CJ(measure, outcome = mm.ref$outcome), on=.(measure)][mm.ref, on=.(outcome)]
+tw[, end := pmin(end, plt.qs[, max(date)])]
+tw[, mid := start + (end-start)/2 ]
+
+p <- allplot(
+  plt.qs, yl = "Per 10k, Cumulative\nIncidence of ...",
+  withRef = FALSE, ins = voc.band(tw)
+)
+
+ggsave(tail(.args, 1), p, height = 6, width = 10, bg = "white")

exp/active-vac/fig/conditioned_eff_all.R

@@ -0,0 +1,50 @@
+
+.pkgs <- c("data.table", "ggplot2", "scales", "ggh4x", "cabputils", "geomtextpath")
+
+stopifnot(all(sapply(.pkgs, require, character.only = TRUE)))
+
+#' assumes R project at the experiment root level
+.args <- commandArgs(args = c(
+  file.path("fig", "vis_support.rda"),
+  file.path("fig", "process", c("alt_eff.rds", "digest-key.rds", "vocwindows.rds")),
+  file.path("fig", "output", "conditioned_eff_all.png")
+))
+
+load(.args[1])
+
+#' comes key'd
+eff.dt <- readRDS(.args[2])[
+  eval(datefilter) & eval(outfilter)
+][, .(
+  scenario, realization, date, outcome, c.effectiveness
+)]
+
+intscns <- eff.dt[, unique(scenario)]
+
+scn.dt <- readRDS(.args[3])[scenario %in% intscns]
+
+takeover.wins <- readRDS(.args[4])
+
+plt.dt <- setkeyv(
+  eff.dt[scn.dt, on=.(scenario)],
+  union(key(eff.dt), colnames(scn.dt))
+)[inf_con == TRUE]
+
+rm(eff.dt)
+gc()
+
+plt.qs <- plt.prep(plt.dt, j = .(c.effectiveness))
+
+mm.ref <- plt.qs[,.(ymin = min(q90l), ymax = max(q90h)),by=.(outcome)]
+mm.ref[, ymin := ymin - .15*(ymax-ymin) ]
+mm.ref[, yspan := ymax - ymin ]
+tw <- takeover.wins[q == 0.5][CJ(measure, outcome = mm.ref$outcome), on=.(measure)][mm.ref, on=.(outcome)]
+tw[, end := pmin(end, plt.qs[, max(date)])]
+tw[, mid := start + (end-start)/2 ]
+
+p <- allplot(
+  plt.qs, yl = "Cumulative Effectiveness\nAgainst Incidence of ...",
+  withRef = TRUE, ins = voc.band(tw)
+)
+
+ggsave(tail(.args, 1), p, height = 6, width = 10, bg = "white")