Changes to accommodate new Figure 1.

.gitignore

@@ -11,4 +11,5 @@
 *.md
 .Rproj.user
 *.rds
-*.png
+.Rhistory
+*.png

WOR.tex

@@ -83,8 +83,8 @@ \section*{Introduction}
 
 \begin{figure}
 \centering
-\includegraphics[width=0.4\textwidth]{main.png}
-\caption{\label{mainfig}Example Benefits of Testing. A) Assuming roughly 20\% of the population (light grey region) could be vaccinated without any testing: with a highly specific test (99\%), testing can expand the percentage covered by increasing amounts depending on seroprevalence and test sensitivity (light blue: 70\% sensitive, dark blue: 90\% sensitivity). While many individuals still go unvaccinated (white background), further incremental vaccine supply would result in similarly extended coverage under pre-vaccination testing approach. B) Testing can also expand the vaccinated fraction for lower costs, again depending on the vaccine-eligible population's seroprevalence and the test's performance. For a highly specific test, and relatively high seroprevalence, \eg* 40\%, a test that is less than 27\% (if >70\% sensitive) to 35\% (if >90\% sensitive) the cost of a vaccine dose is a cheaper way to expand coverage than buying more doses of vaccine. The cost advantage of testing grows with increasing seroprevalence and test performance.}
+\includegraphics[width=\linewidth]{main.png}
+\caption{\label{mainfig}Example Benefits of Testing. A) Assuming roughly 20\% of the population could be vaccinated without any testing \cite{noauthor_gavi_nodate}: with a highly specific test (99\%), testing can expand the percentage covered by increasing amounts depending on seroprevalence and test sensitivity (light blue: 70\% sensitive, orange: 90\% sensitivity). While many individuals still go unvaccinated, further incremental vaccine supply would result in similarly extended coverage under pre-vaccination testing approach. B) Testing can also expand the vaccinated fraction for lower costs, again depending on the vaccine-eligible population's seroprevalence and the test's performance. For a highly specific test, and relatively high seroprevalence, \eg* 40\% (dashed vertical line), a test that is less than 27\% (if >70\% sensitive) to 35\% (if >90\% sensitive) the cost of a vaccine dose is a cheaper way to expand coverage than buying more doses of vaccine. The cost advantage of testing grows with increasing seroprevalence and test performance.}
 \end{figure}
 
 \section*{Methods}
@@ -149,9 +149,9 @@ \subsection*{Expanded Coverage}
 
 \subsection*{Paying for Testing versus More Doses}
 
-Similar to expanding coverage, using \eref{costlim}, we find that regions with relatively high seroprevalence could more effectively expand coverage with antibody testing than additional doses with a relatively high threshold for test cost. Per \eref{costlim}, the cost threshold is identical to the coverage expansion, and as such also increases linearly with seroprevalence and test performance. \fref{mainfig}B shows this trend for seroprevalence (moving left to right on plot) and sensitivity (shifting from lower to higher trend lines).
+Similar to expanding coverage, using \eref{costlim}, we find that regions with relatively high seroprevalence could more effectively expand coverage with antibody testing than additional doses with a relatively high threshold for test cost. Per \eref{costlim}, the cost threshold is identical to the coverage expansion, and as such also increases linearly with seroprevalence and test performance. \fref{mainfig}B shows this trend for seroprevalence (moving left to right on plot) and sensitivity (shifting from the lower to higher shaded region).
 
-Recalling the high seroprevalence (40\%), high specificity (99\%), and high sensitivity (90\%) scenario, the cost threshold would be 35\% of the cost of vaccination. Considering a per dose cost of roughly 15 USD (and ignoring the additional administration and logistics costs), a test that cost around 5 USD or less would be a better investment to expand coverage. Even at a lower sensitivity (70\%), the cost threshold remains plausible at 4 USD or less per test.
+Recalling the high seroprevalence (40\%), high specificity (99\%), and high sensitivity (90\%) scenario, the cost threshold would be 35\% of the cost of vaccination. Considering a per dose cost of roughly 15 USD (and ignoring the additional administration and logistics costs), a test that cost around 5 USD or less would be a better investment to expand coverage. Even at a lower sensitivity (70\%), the cost threshold remains plausible at 27\% of the cost, i.e. 4 USD or less per test.
 
 \section*{Discussion \&\ Conclusions}
 

fig_main.R

@@ -3,6 +3,11 @@ suppressPackageStartupMessages({
   require(data.table)
   require(ggplot2)
   require(patchwork)
+  require(magrittr)
+  require(dplyr)
+  require(scales)
+  require(forcats)
+  require(tibble)
 })
 
 .args <- if (interactive()) c(
@@ -19,148 +24,143 @@ refspec <- .99
 
 ref[, notest_sp := covaxcoverage*seropos ]
 ref[, notest_sn := covaxcoverage*(1-seropos) ]
-ref[, test70_sp := notest_sp*(PPDmul(refspec, .7, seropos) - 1) ]
-ref[, test70_sn := notest_sn*(PPDmul(refspec, .7, seropos) - 1) ]
-ref[, test90_sp := notest_sp*(PPDmul(refspec, .9, seropos)-1) - test70_sp ]
-ref[, test90_sn := notest_sn*(PPDmul(refspec, .9, seropos)-1) - test70_sn ]
+ref[, test70_sp := notest_sp*(PPDmul(refspec, .7, seropos))]
+ref[, test70_sn := notest_sn*(PPDmul(refspec, .7, seropos))]
+ref[, test90_sp := notest_sp*(PPDmul(refspec, .9, seropos))]
+ref[, test90_sn := notest_sn*(PPDmul(refspec, .9, seropos))]
 
 ref.mlt <- melt.data.table(
   ref,
   id.vars = "seropos"
 )[, c("scenario", "serostatus") := tstrsplit(variable, "_")]
 
 ref.mlt[,
-  ord := 2*c("notest"=2,"test70"=1,"test90"=0)[scenario]+c(sp=0, sn=1)[serostatus]
-]
+        ord := 2*c("notest"=2,"test70"=1,"test90"=0)[scenario]+c(sp=0, sn=1)[serostatus]
+        ]
 
-ref[, costlim70 := PPDmul(refspec, .7, seropos) - 1 ]
-ref[, costlim90 := PPDmul(refspec, .9, seropos) - 1 ]
 
-lblsize <- 3
+lblsize <- 2
 lnsize <- 1
 
-p.coverage <- ggplot(ref.mlt) +
-  aes(seropos, value, group = ord, alpha = serostatus, fill = scenario) +
-  geom_area() +
-  coord_cartesian(xlim = c(0.2, 0.8), ylim = c(0, 0.4), expand = FALSE) +
-  scale_alpha_manual(
-    values = c("sn"=.5, "sp"=.75),
-    guide = "none"
-  ) +
-  scale_fill_manual(
-    values = c(notest="grey5", test70=blues9[5], test90=blues9[9]),
-    guide = "none"
-  ) +
-  geom_line(aes(seropos, y=covaxcoverage, linetype = "notest"), ref, inherit.aes = FALSE, show.legend = F, size = lnsize) +
-  geom_line(aes(seropos, y=(costlim70+1)*covaxcoverage, linetype = "test70"), ref, inherit.aes = FALSE, show.legend = F, size = lnsize) +
-  geom_line(aes(seropos, y=(costlim90+1)*covaxcoverage, linetype = "test90"), ref, inherit.aes = FALSE, show.legend = F, size = lnsize) +
-  annotate(
-    "label",
-    x=0.5, y= covaxcoverage * .90,
-    label = expression(symbol('\257')*"  "*symbol('\257')*'  Coverage without Testing  '*symbol('\257')*"  "*symbol('\257')),
-    label.size = 0, fill = alpha(c("white"), 0.5),
-    size = lblsize
-  ) +
-  annotate(
-    "label",
-    x = (0.8-0.2)*.95 + 0.2, y = (PPDmul(refspec, .7, 0.8)*covaxcoverage + covaxcoverage)/2,
-    label = '... >70% sensitive test',
-    hjust = "right",
-    label.size = 0, fill = alpha(c("white"), 0.5),
-    size = lblsize
-  ) +
-  annotate(
-    "label",
-    x = (0.8-0.2)*.95 + 0.2,
-    y = (PPDmul(refspec, .9, 0.8)+PPDmul(refspec, .7, 0.8))*covaxcoverage/2,
-    label = 'additional coverage with\n >90% sensitive test',
-    hjust = "right",
-    label.size = 0, fill = alpha(c("white"), 0.5),
-    size = lblsize
-  ) +
-  annotate(
-    "text",
-    x = (0.8-0.2)*.05 + .2,
-    y = covaxcoverage/2,
-    label = 'Seronegative\nIndividuals',
-    hjust = "left",
-    color = "black",
-    size = lblsize
-  ) +
-  annotate(
-    "text",
-    x = (0.8-0.2)*.95 + .2,
-    y = covaxcoverage/2,
-    label = 'Seropositive\nIndividuals',
-    hjust = "right",
-    color = "white",
-    size = lblsize
-  ) +
-  scale_y_continuous("Total Immunized %", breaks = seq(0,0.5,0.1), minor_breaks = NULL, labels = function(b) sprintf("%i%%", as.integer(b*100))) +
-  scale_x_continuous("% Seropositive", minor_breaks = NULL, labels = function(b) sprintf("%i%%", as.integer(b*100))) +
-  scale_linetype_manual(
-    values = c(notest="solid", test70="longdash", test90="dotted")
-  ) +
+cols <- c('seropos', 'value')
+ref.mlt2 <- data.table::copy(ref.mlt)[ ,
+                                       c("value", "seropos") := list(floor((100)*value),
+                                                                     round(100*seropos))
+                                       ][(seropos %% 10) == 0]
+
+
+scale_x <-     scale_x_continuous(labels = scales::percent,
+                                  breaks = seq(0,1,by=0.20),
+                                  limits = c(0.2, 0.8), expand = expansion(mult = 0.1)) 
+
+scale_palette <- data.frame(color = c(notest = "black", 
+                                      test70 = "lightskyblue", 
+                                      test90 = "orange"),
+                            labels = c("No test",
+                                       "70% ≤ Sensitivity < 90%", 
+                                       "Sensitivity ≥ 90%")) %>% 
+  tibble::rownames_to_column("breaks")
+
+
+scale_fill <- scale_fill_manual(
+  values = setNames(scale_palette$color, scale_palette$breaks),
+  labels = scale_palette$labels,
+  name   = "Scenario") 
+
+scale_color <- scale_color_manual(
+  values = setNames(scale_palette$color, scale_palette$breaks),
+  labels = scale_palette$labels,
+  name   = "Scenario")
+
+fig1 <- ref.mlt %>%
+  dplyr::mutate(Serostatus = forcats::fct_recode(serostatus,
+                                                 `Sero-negative` = 'sn',
+                                                 `Sero-positive` = 'sp'),
+                scenario = factor(scenario)) %>%
+  dplyr::bind_rows(.,
+                   dplyr::group_by_at(.,
+                                      .vars = vars(-c(Serostatus, 
+                                                      serostatus, ord, value,
+                                                      variable))
+                                      ) %>%
+                     summarise(value = sum(value)) %>%
+                     mutate(Serostatus = "All")) %>%
+  ggplot(data = ., aes(y = value, x = seropos)) +
+  geom_line(aes(color = scenario)) +
   theme_minimal() +
-  theme(
-    plot.margin = margin(r = 1, unit = "line")
-  )
+  ylab("Total immunized, %") + 
+  xlab("Percent seropositive") +
+  scale_y_continuous(labels = scales::percent, limits = c(0,NA), expand = c(0,0)) +
+  scale_x +
+  theme(legend.position  = 'bottom',
+        legend.direction = 'horizontal',
+        legend.spacing.y = unit(0, 'mm')) +
+  scale_color +
+  facet_grid(cols = vars(Serostatus))#, margins = 'Serostatus')
+
+annotations <- list(data.frame(x = 0.6,
+                               y = (PPDmul(refspec, .7, .6) - 1)/2,
+                               label = 'prefer tests\nif sensitivity ≥ 70%',
+                               scenario = 'test70',
+                               hjust = 1),
+                    data.frame(x = 0.6,
+                               y = (PPDmul(refspec, .7, .8) + PPDmul(refspec, .9, .8) - 2)/3,
+                               label = 'prefer tests\nif sensitivity ≥ 90%',
+                               scenario = 'test90',
+                               hjust = 0.5),
+                    data.frame(x = .6,
+                               y = .75,
+                               label = 'prefer additional\nvaccine doses',
+                               scenario = 'notest',
+                               hjust = 0)) %>%
+  dplyr::bind_rows(.)
 
-p.cost <- ggplot(ref) +
+ref[, costlim70 := PPDmul(refspec, .7, seropos) - 1 ]
+ref[, costlim90 := PPDmul(refspec, .9, seropos) - 1 ]
+
+
+p.cost <- 
+  ggplot(ref) +
+  geom_vline(lty = 2, xintercept = 0.4) +
   aes(seropos) +
-  geom_blank(aes(linetype = "notest")) +
-  geom_ribbon(aes(ymin=0, ymax=costlim70), fill = alpha(blues9[5], .8), show.legend = F) +
-  geom_ribbon(aes(ymin=costlim70, ymax=costlim90), fill = alpha(blues9[9], .8), show.legend = F) +
-  geom_line(aes(y=costlim70, linetype = "test70"), size = lnsize) +
-  geom_line(aes(y=costlim90, linetype = "test90"), size = lnsize) +
-  annotate(
-    "label",
-    x = (0.8-0.2)*.975 + .2,
-    y = (PPDmul(refspec, .7, .6) - 1)/3*2,
-    label = '...if sensitivity >70%',
-    hjust = "right",
-    label.size = 0, fill = alpha(c("white"), 0.5),
-    vjust = "top",
-    size = lblsize
-  ) +
-  annotate(
-    "label",
-    x = (0.8-0.2)*.975 + .2,
-    y = (PPDmul(refspec, .7, .8) + PPDmul(refspec, .9, .8) - 2)/2,
-    label = 'prefer tests\nif sensitivity > 90%',
-    hjust = "right",
-    label.size = 0, fill = alpha(c("white"), 0.5),
-    size = lblsize
-  ) +
-  annotate(
-    "text",
-    x = .3,
-    y = .75+.125,
-    label = 'prefer vaccine doses',
-    hjust = "left",
-    size = lblsize
-  ) +
-  coord_cartesian(xlim = c(0.2, 0.8), ylim = c(0, 1), expand = FALSE) +
-  scale_x_continuous("% Seropositive", minor_breaks = NULL, labels = function(b) sprintf("%i%%", as.integer(b*100))) +
-  scale_y_continuous("Testing Cost % of Vaccination Dose Cost", labels = function(b) sprintf("%i%%", as.integer(b*100))) +
-  scale_linetype_manual(
-    "Test Performance",
-    breaks = c("notest", "test70", "test90"),
-    labels = c(
-      "notest"="(No Test)",
-      "test70"="70% Sensitivity",
-      "test90"="90% Sensitivity"
-    ),
-    values = c("solid", "longdash", "dotted"),
-    guide = "none"
-  ) +
+  #geom_blank(aes(linetype = "notest")) +
+  geom_ribbon(aes(ymin=0, ymax=costlim70), 
+              fill = 'lightskyblue',
+              alpha = 0.5,
+              show.legend = FALSE) +
+  geom_ribbon(aes(ymin=costlim70, ymax=costlim90), 
+              fill = "orange",
+              alpha = 0.5,
+              show.legend = FALSE) +
+  scale_color + 
+  scale_fill +
+  scale_y_continuous(limits = c(0, 1), expand = c(0,0),
+                     "Testing Cost as % of\nVaccination Dose Cost", 
+                     labels = function(b) sprintf("%i%%", as.integer(b*100))) +
+  scale_x +
+  xlab("Percent seropositive") +
   theme_minimal() +
   theme(
-    plot.margin = margin(r = 1, unit = "line")
+    plot.margin = margin(r = 1, unit = "line"),
+    legend.position = 'none'
   )
 
-res.p <- p.coverage / p.cost +
-  plot_annotation(tag_levels = "A") +
-  plot_layout(guides = "collect")
+for (i in 1:nrow(annotations)){
+  p.cost <- p.cost + 
+    geom_label(
+      x = annotations[i,"x"],
+      y = annotations[i,"y"],
+      label = annotations[i,"label"],
+      fill = dplyr::filter(scale_palette, breaks == annotations[i,"scenario"]) %>%
+        dplyr::pull(color) %>%
+        sub(pattern = 'black', replacement = 'white', x = .) ,
+      size = lblsize)
+}
+
+res.p <- (fig1 + p.cost +
+            plot_annotation(tag_levels = "A") +
+            plot_layout(guides = "collect", widths = c(3,1))) &
+  theme(legend.position = 'bottom')
+
+ggsave(tail(.args, 1), res.p, height = 3, width = 9, dpi = 600)
 
-ggsave(tail(.args, 1), res.p, height = 7.5, width = 4, dpi = 600)