diff --git a/.covrignore b/.covrignore
new file mode 100644
index 00000000..ddb6d3dd
--- /dev/null
+++ b/.covrignore
@@ -0,0 +1,4 @@
+# This does the API call to packit - not 
+# easy to test without mocking all the network.
+
+R/packit.R
diff --git a/DESCRIPTION b/DESCRIPTION
index 2e25561a..8eea4ac0 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: stoner
 Title: Support for Building VIMC Montagu Touchstones, using Dettl
-Version: 0.1.19
+Version: 0.1.20
 Authors@R:
     c(person("Wes", "Hinsley",role = c("aut", "cre", "cst", "dnc", "elg", "itr", "sng", "ard"),
            email = "w.hinsley@imperial.ac.uk"),
@@ -20,6 +20,7 @@ Imports:
     cli,
     DBI,
     data.table,
+    httr2,
     jsonlite,
     lgr,
     dplyr,
@@ -27,6 +28,8 @@ Imports:
     prettyunits,
     readr,
     rlang,
+    shiny,
+    shinyjs,
     testthat,
     utils,
     withr
@@ -35,6 +38,7 @@ RoxygenNote: 7.3.3
 Roxygen: list(markdown = TRUE)
 Suggests:
     knitr,
+    mockery,
     rcmdcheck,
     rmarkdown,
     RPostgres
diff --git a/NAMESPACE b/NAMESPACE
index d71f8c5b..0fa89610 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,5 +1,7 @@
 # Generated by roxygen2: do not edit by hand
 
+export(fetch_packit)
+export(stochastic_explorer)
 export(stone_dump)
 export(stone_extract)
 export(stone_load)
@@ -16,17 +18,22 @@ export(stoner_calculate_dalys)
 export(stoner_dalys_for_db)
 import(arrow)
 import(dplyr)
+import(httr2)
 import(readr)
 importFrom(arrow,write_parquet)
 importFrom(data.table,as.data.table)
 importFrom(data.table,rbindlist)
+importFrom(grDevices,recordPlot)
 importFrom(graphics,lines)
+importFrom(graphics,par)
 importFrom(magrittr,"%>%")
 importFrom(rlang,":=")
+importFrom(shiny,runApp)
 importFrom(stats,median)
 importFrom(stats,quantile)
 importFrom(testthat,expect_equal)
 importFrom(testthat,expect_true)
+importFrom(tools,file_ext)
 importFrom(utils,capture.output)
 importFrom(utils,read.csv)
 importFrom(utils,write.csv)
diff --git a/NEWS.md b/NEWS.md
index 0639fd81..86600b47 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,8 +1,13 @@
-# Stoner 0.1.19
+# Stoner 0.1.20
 
-* add stoner_stochastic_central to create central parquet from standardised stochastics.
+* Allow missing yll for processing older stochastics
+
+# Changes in 0.1.19
+
+* Add stoner_stochastic_central to create central parquet from standardised stochastics.
 
 # Changes in 0.1.18
+
 * Add stoner_stochastic_standardise for converting (drop)box stochastics to standard form
 * Add stoner_stochastic_graphs for quick plotting of stochastics.
 
diff --git a/R/packit.R b/R/packit.R
new file mode 100644
index 00000000..087b0981
--- /dev/null
+++ b/R/packit.R
@@ -0,0 +1,77 @@
+##' Download a file from a packit, for example an artefact containing
+##' central burden estimates, which we could then plot on top of
+##' stochastics. But more generally, we can use this for fetching any
+##' file from packet (ie, from the Montagu Reporting Portal).
+##'
+##' @export
+##' @title Fetch packet from a packit server
+##' @import httr2
+##' @importFrom tools file_ext
+##' @param packet_id The id of the packet containing the artefact.
+##' @param filename The filename of the file within the packet.
+##' @param server By default, the URL to the packit API on Montagu,
+##' but this can be set to other packit API's if we want.
+##' @returns The filename of the temporary file which has been downloaded.
+fetch_packit <- function(packet_id, filename,
+  server = "https://montagu.vaccineimpact.org/packit/api/") {
+
+  # First we have to create a client, and a flow...
+
+  client <- httr2::oauth_client(
+    id = "orderly",
+    token_url = sprintf("%s/deviceAuth/token", server),
+    name = "orderly"
+  )
+
+  # You will be asked to type a code at this point...
+
+  flow <- httr2::oauth_flow_device(
+    client = client,
+    auth_url = sprintf("%s/deviceAuth", server),
+    pkce = FALSE,
+    scope = NULL,
+    open_browser = FALSE,
+    auth_params = list(),
+    token_params = list()
+  )
+
+  # Now we have a flow$access_token we can use, to get
+  # a one-time token to download the file.
+
+  req <- httr2::request(sprintf("%s/packets/%s/files/token", server, packet_id))
+  req <- req |>
+    httr2::req_headers("Accept" = "application/json",
+                       "Content-Type" = "application/json",
+                       "Authorization" = paste("Bearer", flow$access_token)) |>
+         httr2::req_method("POST") |>
+         httr2::req_body_json(list(paths = list(filename)))
+
+  x <- httr2::req_perform(req)
+  ott <- httr2::resp_body_json(x)$id
+
+  # This is the URL to the file we want, including the token.
+
+  url <- sprintf("%s/packets/%s/file?path=%s&token=%s&filename=%s&inline=false",
+                 server, packet_id, filename, ott, filename)
+
+  req <- httr2::request(url)
+  req <- req |>
+    httr2::req_headers("Authorization" = paste("Bearer", flow$access_token))
+
+
+  # And finally, we download the file in chunks.
+
+  out <- tempfile(fileext = paste0(".", tools::file_ext(filename)))
+  outfile <- file(out, open = "wb")
+  con <- httr2::req_perform_connection(req, blocking = TRUE)
+  while (!httr2::resp_stream_is_complete(con)) {
+    chunk <- httr2::resp_stream_raw(con, kb = 32)
+    if (length(chunk) == 0) break
+    writeBin(chunk, outfile)
+  }
+  close(outfile)
+  out
+}
+
+
+
diff --git a/R/stochastic_files.R b/R/stochastic_files.R
index 0f3f0dbe..a185813c 100644
--- a/R/stochastic_files.R
+++ b/R/stochastic_files.R
@@ -31,14 +31,30 @@
 ##' these to the simpler names. Processing Rubella stochastic files without
 ##' this set to TRUE will fail - so while we should always do this, keeping
 ##' the parameter makes it more clear in the code what we're doing and why.
+##' @param hepb2019_fix In 2019 (and 2017), HepB deaths and cases were subdivided into
+##' number of different causes. This flag combines those into the single
+##' appropriate burden outcome.
+##' @param hib2019_fix In 2019 (and 2017), Hib deaths and cases were subdivided into
+##' number of different causes. This flag combines those into the single
+##' appropriate burden outcome.
 ##' @param missing_run_id_fix Some groups in the past have omitted run_id
 ##' from the files, but included them in the filenames. This fix inserts
 ##' that into the files if the index parameter indicates we have 200 runs to
 ##' process.
-
+##' @param allow_missing_yll yll was introduced in 2023? This flag allows
+##' it to be missing for processing older stochastics.
+##' @param allow_missing_dalys Some early groups did not provide dalys; this
+##' flag allows dalys to be skipped.
+##' @param allow_missing_indexes In some early runs, different groups
+##' provided different numbers of files for each scenario, because some
+##' countries did not implement particular coverage campaigns. This
+##' flag needs to be TRUE for those groups, but the default is FALSE,
+##' since it's rare, and we generally want errors for missing files.
 stone_stochastic_standardise <- function(
     group, in_path, out_path, scenarios, files, index = 1,
-    rubella_fix = TRUE, missing_run_id_fix = TRUE) {
+    rubella_fix = TRUE, hepb2019_fix = TRUE, hib2019_fix = TRUE,
+    missing_run_id_fix = TRUE, allow_missing_yll = TRUE,
+    allow_missing_dalys = TRUE, allow_missing_indexes = FALSE) {
 
   dir.create(out_path, showWarnings = FALSE, recursive = TRUE)
   if ((length(files) == 1) && (grepl(":scenario", files))) {
@@ -47,36 +63,115 @@ stone_stochastic_standardise <- function(
       files[j] <- gsub(":scenario", scenarios[j], files[j])
     }
   }
-
   for (i in seq_along(scenarios)) {
     message(scenarios[i])
     all_data <- as.data.frame(data.table::rbindlist(lapply(index, function(j) {
       cat("\r", j)
       file <- gsub(":index", j, files[i])
-      d <- read.csv(file.path(in_path, file))
+      filepath <- file.path(in_path, file)
+      if (!file.exists(filepath) && (allow_missing_indexes)) {
+        return(NULL)
+      }
+      d <- read.csv(filepath)
       d$country_name <- NULL
 
-      # Fixes needed to standardise Rubella
+      # Various accumulated fixes for e/non-standard stochastics
+      # The flags are all true by default - Stoner will break untidily
+      # if the flags are turned off, but the problem occurs.
+
       if (rubella_fix) {
-        names(d)[names(d) == "rubella_deaths_congenital"] <- "deaths"
-        names(d)[names(d) == "rubella_cases_congenital"] <- "cases"
-        d$rubella_infections <- NULL
+        if ("rubella_deaths_congenital" %in% names(d)) {
+          message("Converting rubella_deaths_congenital to deaths")
+          names(d)[names(d) == "rubella_deaths_congenital"] <- "deaths"
+        }
+        if ("rubella_cases_congenital" %in% names(d)) {
+          message("Converting rubella_cases_congenital to cases")
+          names(d)[names(d) == "rubella_cases_congenital"] <- "cases"
+        }
+        if ("rubella_infections" %in% names(d)) {
+          message("Ignoring rubella_infections")
+          d$rubella_infections <- NULL
+        }
+      }
+
+      if (hepb2019_fix) {
+        if (!"cases" %in% names(d)) {
+          d$cases <- 0
+          cda_cases <- c("hepb_cases_acute_severe", "hepb_cases_dec_cirrh" ,
+                         "hepb_cases_hcc")
+          li_cases <- c("hepb_cases_acute_symp", "hepb_cases_fulminant",
+                        "hepb_cases_chronic",
+                        "hepb_chronic_symptomatic_in_acute_phase")
+          ic_cases <- c("hepb_cases_acute_severe", "hepb_cases_comp_cirrh",
+                        "hepb_cases_hcc_no_cirrh")
+          for (i in unique(c(cda_cases, li_cases, ic_cases))) {
+            if (i %in% names(d)) {
+              message(sprintf("Including %s in cases", i))
+              d$cases <- d$cases + d[[i]]
+              d[[i]] <- NULL
+            }
+          }
+        }
+
+        if (!"deaths" %in% names(d)) {
+          d$deaths <- 0
+          cda_deaths <- c("hepb_deaths_acute", "hepb_deaths_dec_cirrh",
+                          "hepb_deaths_hcc")
+          li_deaths <- c("hepb_deaths_acute", "hepb_deaths_total_cirrh",
+                         "hepb_deaths_hcc")
+          ic_deaths <- c("hepb_deaths_acute", "hepb_deaths_comp_cirrh",
+                         "hepb_deaths_dec_cirrh", "hepb_deaths_hcc")
+
+          for (i in unique(c(cda_deaths, li_deaths, ic_deaths))) {
+            if (i %in% names(d)) {
+              message(sprintf("Including %s in deaths", i))
+              d$deaths <- d$deaths + d[[i]]
+              d[[i]] <- NULL
+            }
+          }
+        }
       }
 
-      # Detect where run_id is missing, but in filenames
+      if (hib2019_fix) {
+        if (("cases_pneumo" %in% names(d)) &&
+            ("cases_men" %in% names(d))) {
+
+          message("cases = cases_men + cases_pneumo")
+          d$cases <- d$cases_pneumo + d$cases_men
+          d$cases_pneumo <- NULL
+          d$cases_men <- NULL
+        }
+        if (("deaths_pneumo" %in% names(d)) &&
+            ("deaths_men" %in% names(d))) {
+          message("deaths = deaths_men + deaths_pneumo")
+          d$deaths <- d$deaths_pneumo + d$deaths_men
+          d$deaths_pneumo <- NULL
+          d$deaths_men <- NULL
+        }
+      }
 
       if (missing_run_id_fix) {
         if ((!"run_id" %in% names(d)) && (length(index) == 200)) d$run_id <- j
       }
 
+
       # Round to integer, as per guidance. (Not using as.integer, as that
       # has limits on how large numbers can be, so we are just truncating
       # digits here)
 
-      d$dalys <- round(d$dalys)
+      if (("dalys" %in% names(d)) || (!allow_missing_dalys)) {
+        d$dalys <- round(d$dalys)
+      } else {
+        message("Dalys missing. (Ignored)")
+      }
+
       d$deaths <- round(d$deaths)
       d$cases <- round(d$cases)
-      d$yll <- round(d$yll)
+      if (("yll" %in% names(d)) || (!allow_missing_yll)) {
+        d$yll <- round(d$yll)
+      } else {
+        message("yll missing. (Ignored)")
+      }
       d$cohort_size <- round(d$cohort_size)
 
       d
diff --git a/R/stochastic_graphs.R b/R/stochastic_graphs.R
index eee71841..172369f0 100644
--- a/R/stochastic_graphs.R
+++ b/R/stochastic_graphs.R
@@ -1,3 +1,11 @@
+age_string <- function(ages) {
+  if (is.null(ages)) {
+    "all ages"
+  } else if (identical(unique(sort(ages)), as.numeric(min(ages):max(ages)))) {
+    sprintf("age %d..%d", min(ages), max(ages))
+  } else "selected ages"
+}
+
 ##' Draw a stochastic plot showing all the different runs, with the mean,
 ##' median, 5% and 95% quantiles shown.
 ##'
@@ -6,8 +14,9 @@
 ##' @import dplyr
 ##' @importFrom rlang :=
 ##' @import arrow
-##' @importFrom graphics lines
+##' @importFrom graphics lines par
 ##' @importFrom stats quantile median
+##' @importFrom grDevices recordPlot
 ##' @param base The folder in which the standardised stochastic files are found.
 ##' @param touchstone The touchstone name (for the graph title)
 ##' @param disease The disease, used for building the filename and graph title.
@@ -21,23 +30,61 @@
 ##' @param by_cohort If TRUE, then age is subtracted from year to convert it to
 ##' year of birth before aggregating.
 ##' @param log If TRUE, then use a logged y-axis.
+##' @param packit_id If set, then read central burden estimates from a file
+##' within a packit on the Montagu packit server.
+##' @param packit_file Used with packit_id to specify the filename of an RDS
+##' file providing burden estimates. We expect to find scenario, year, age,
+##' country, burden_outcome and value fields in the table.
+##' @param include_quantiles Default TRUE, select whether to plot the
+##' 5% and 95% quantile lines.
+##' @param include_mean Default TRUE, select whether to plot the mean.
+##' @param include_median Default TRUE, select whether to plot the median.
+##' @param scenario2 Default NULL; if set, then the burdens from this
+##' scenario will be subtracted from those in `scenario` - ie, this plots
+##' an impact graph of applying the second scenario. For many graphs that
+##' use this, the result will be positive numbers, representing cases
+##' or deaths averted.
 
 stone_stochastic_graph <- function(base, touchstone, disease, group, country,
                                    scenario, outcome, ages = NULL,
-                                   by_cohort = FALSE, log = FALSE) {
+                                   by_cohort = FALSE, log = FALSE,
+                                   packit_id = NULL, packit_file = NULL,
+                                   include_quantiles = TRUE,
+                                   include_mean = TRUE,
+                                   include_median = TRUE,
+                                   scenario2 = NULL) {
 
   d <- prepare_graph_data(base, touchstone, disease, group, country,
                          scenario, outcome, ages, by_cohort)
+  age <- age_string(ages)
 
-  title <- sprintf("%s, %s, %s\n%s, %s\n", touchstone, disease, group,
+  title <- sprintf("%s, %s, %s, %s\n%s, %s\n", touchstone, disease, group, age,
                    scenario, country)
+  outcome_ylab <- outcome
+
+  if (!is.null(scenario2)) {
+    d2 <- prepare_graph_data(base, touchstone, disease, group, country,
+                            scenario2, outcome, ages, by_cohort)
+    title <- sprintf("%s, %s, %s, %s\n%s, %s\n", touchstone, disease, group, age,
+                     sprintf("Impact of %s ->\n%s", scenario, scenario2), country)
+    d <- d[order(d$year, d$run_id), ]
+    d2 <- d2[order(d2$year, d2$run_id), ]
+    d[[outcome]] <- d[[outcome]] - d2[[outcome]]
+    outcome_ylab <- paste(outcome_ylab, "averted")
+  }
 
   runs <- max(d$run_id)
   miny <- max(1, min(d[[outcome]]))
   maxy <- max(d[[outcome]])
   log <- if (log) "y" else ""
 
-  plot(ylab = outcome, xlab = if (by_cohort) "Birth Cohort" else "year",
+  if (!is.null(packit_id)) {
+    central <- prepare_central_data(packit_id, packit_file,
+                                    country, scenario, outcome, ages, by_cohort)
+    # To be continued...
+  }
+  par(mar = c(5, 4, 5, 2))
+  plot(ylab = outcome_ylab, xlab = if (by_cohort) "Birth Cohort" else "year",
        x = d$year[d$run_id == 1], y = d[[outcome]][d$run_id == 1], type="l",
        col = "#b0b0b0", ylim = c(miny, maxy), main = title, log = log)
 
@@ -54,10 +101,17 @@ stone_stochastic_graph <- function(base, touchstone, disease, group, country,
       q95 = quantile(.data[[outcome]], 0.95),
       .groups = "drop"
     )
-  lines(x = avgs$year, y = avgs$mean, col = "#ff4040", lwd = 2)
-  lines(x = avgs$year, y = avgs$median, col = "#00ff00", lwd = 2)
-  lines(x = avgs$year, y = avgs$q05, col = "#202020", lwd = 2)
-  lines(x = avgs$year, y = avgs$q95, col = "#202020", lwd = 2)
+  if (include_mean) {
+    lines(x = avgs$year, y = avgs$mean, col = "#ff4040", lwd = 2)
+  }
+  if (include_median) {
+    lines(x = avgs$year, y = avgs$median, col = "#00ff00", lwd = 2)
+  }
+  if (include_quantiles) {
+    lines(x = avgs$year, y = avgs$q05, col = "#202020", lwd = 2)
+    lines(x = avgs$year, y = avgs$q95, col = "#202020", lwd = 2)
+  }
+  recordPlot()
 }
 
 
@@ -82,3 +136,53 @@ prepare_graph_data <- function(base, touchstone, disease, group, country,
   d
 }
 
+prepare_central_data <- function(packit_id, packit_file,
+  country, scenario, outcome, ages, by_cohort) {
+
+  central <- readRDS(fetch_packit(packit_id, packit_file))
+  central <- central[central$country == country, ]
+  central <- central[central$scenario == scenario, ]
+  central <- central[central$burden_outcome == outcome, ]
+  if (!is.null(ages)) {
+    central <- central[central$age %in% ages, ]
+  }
+  if (by_cohort) {
+    central$year <- central$year - central$age
+  }
+  names(central)[names(central) == "value"] <- outcome
+  central <- central %>% group_by(.data$year) %>%
+    summarise(
+      !!outcome := sum(.data[[outcome]], na.rm = TRUE),
+      .groups = "drop")
+  central$run_id <- 0
+  central
+}
+
+##' Launch a Shiny app to allow interactive plotting of
+##' standardised stochastic data, burden estimates,
+##' impacts, comparisons between touchstones, and
+##' comparisons between modelling groups.
+##'
+##' @export
+##' @importFrom shiny runApp
+##' @title Stochastic plot
+##' @param data_dir The location of the standardised stochastic folder
+##' hierarchy; this can be a local path, a fully-qualified network path on
+##' windows, or a mount point on linux or Mac.
+stochastic_explorer <- function(
+  data_dir = "//wpia-hn2.hpc.dide.ic.ac.uk/vimc_stochastics") {
+  if (!dir.exists(data_dir)) {
+    cli::cli_abort(c(
+      "x" = "Cannot access the path/mount: {.path {data_dir}}",
+      "i" = "Please check you can see this path normally. If not:",
+      "*" = "You need ZScaler on and connected if you're not within DIDE",
+      "*" = "(But strongly advise remote desktop into DIDE - large files)",
+      "*" = "On linux, ensure the mount is correctly set up.",
+      "*" = "Check with DIDE IT that you have access to VIMC files",
+      "*" = "Check your general internet access."
+    ))
+  }
+
+  assign("data_dir", data_dir, envir = .GlobalEnv)
+  runApp(system.file("app", package = "stoner"))
+}
diff --git a/README.md b/README.md
index c115f9f7..163a6e3b 100644
--- a/README.md
+++ b/README.md
@@ -4,7 +4,7 @@
 [![R build status](https://github.com/vimc/stoner/workflows/R-CMD-check/badge.svg)](https://github.com/vimc/stoner/actions)
 <!-- badges: end -->
 
-<pre>
+<pre style="line-height: 1;">
      .-'''-. ,---------.    ,-----.    ,---.   .--.    .-''-.  .-------.
     / _     \\          \ .'  .-,  '.  |    \  |  |  .'_ _   \ |  _ _   \
    (`' )/`--' `--.  ,---'/ ,-.|  \ _ \ |  ,  \ |  | / ( ` )   '| ( ' )  |
diff --git a/inst/app/app.R b/inst/app/app.R
new file mode 100644
index 00000000..bef62580
--- /dev/null
+++ b/inst/app/app.R
@@ -0,0 +1,342 @@
+library(shiny)
+library(shinyjs)
+library(arrow)
+
+source("utils.R")
+
+app_ui <- function() {
+  useShinyjs()
+
+  make_multi_input <- function(prefix, name, count, choices,
+                               n1 = name, n2 = NULL) {
+    if (count == 1) {
+      return(list(selectInput(sprintf("%s_%s", prefix, name),
+                              paste0(name, ":"),
+                              choices = choices, selected = NULL)))
+    } else {
+      return(list(selectInput(sprintf("%s_%s1", prefix, name), n1,
+                              choices = choices, selected = NULL),
+                  selectInput(sprintf("%s_%s2", prefix, name), n2,
+                              choices = choices, selected = NULL)))
+    }
+  }
+
+  make_panel <- function(title, prefix, n_touchstones, n_scenarios, n_groups) {
+    c0 <- character(0)
+    ts <- get_touchstones()
+
+    touchstones <- make_multi_input(prefix, "touchstone", n_touchstones, ts,
+                                    "Touchstone 1:", "Touchstone 2:")
+    scenarios <- make_multi_input(prefix, "scenario", n_scenarios, c0,
+                                    "Base Scenario:", "Compare Scenario:")
+    groups <- make_multi_input(prefix, "group", n_groups, c0,
+                               "Group 1:", "Group 2:")
+
+    if ((n_touchstones > 1) || (n_groups > 1)) {
+      graphs <- list(
+        plotOutput(sprintf("%s_main_plot1", prefix)),
+        plotOutput(sprintf("%s_main_plot2", prefix)))
+    } else {
+      graphs <- list(
+        plotOutput(sprintf("%s_main_plot", prefix)))
+    }
+
+    sidebar <- append(touchstones, list(
+      selectInput(sprintf("%s_disease", prefix), "Disease:", choices = c0)))
+    sidebar <- append(sidebar, groups)
+    sidebar <- append(sidebar, scenarios)
+    sidebar <- append(sidebar, list(
+      selectInput(sprintf("%s_country", prefix), "Country:", choices = c0),
+      selectInput(sprintf("%s_outcome", prefix), "Y-Axis:", choices = c0),
+      radioButtons(sprintf("%s_year", prefix), "X-Axis:", inline = TRUE,
+                   choices = c("Calendar", "Cohort")),
+      radioButtons(sprintf("%s_ages", prefix), "Age:", inline = TRUE,
+                   choices = c("All", "Under 5")),
+      checkboxGroupInput(
+        sprintf("%s_options", prefix), "Plot Options:",
+        choices = c("Quantiles", "Median", "Mean", "Log-Y"),
+        selected = c("Quantiles", "Median", "Mean", "Log-Y"),
+        inline = TRUE),
+      actionButton(sprintf("%s_plot_btn", prefix), "Plot")
+    ))
+
+    tabPanel(title, sidebarLayout(sidebarPanel(sidebar), mainPanel(graphs)))
+  }
+
+  fluidPage(
+    tags$style(HTML("
+      .shiny-busy {
+        cursor: wait !important;
+      }
+      .selectize-input {
+        white-space: nowrap;
+        overflow: hidden;
+        text-overflow: ellipsis;
+      }
+
+      .shiny-options-group {
+        display: grid;
+        grid-template-columns: repeat(2, 1fr);  /* 2 per row */
+        gap: 4px 12px;  /* row gap, column gap */
+      }
+
+      .shiny-options-group .checkbox label {
+        display: flex;
+        align-items: centre;
+        margin: 0;
+        padding-left: 0;
+        gap: 6px;
+      }
+
+      .control-label {
+        margin-bottom: 2px;
+        font-size: 12px;
+      }
+
+      .selectize-input {
+        min-height: 30px;
+        padding: 4px 8px;
+        font-size: 13px;
+      }
+
+      .form-group {
+        margin-bottom: 8px;
+      }
+
+      .shiny-plot-output {
+        padding-top: 5px;
+      }
+    ")),
+    titlePanel("Stochastic Explorer"),
+    tabsetPanel(
+      make_panel("Burden", "b", 1, 1, 1),
+      make_panel("Burden/TS", "bts", 2, 1, 1),
+      make_panel("Burden/MG", "bmg", 1, 1, 2),
+      make_panel("Impact", "i", 1, 2, 1),
+      make_panel("Impact/TS", "its", 2, 2, 1),
+      make_panel("Impact/MG", "img", 1, 2, 2)
+    )
+  )
+}
+
+app_server <- function(input, output, session) {
+
+  add_observers <- function(prefix, n_touchstone = 1, n_group = 1,
+                            n_scenario = 1) {
+    if (n_touchstone == 1) {
+      it1 <- sprintf("%s_touchstone", prefix)
+      it2 <- NULL
+      observeEvent(input[[it1]], {
+                   update_touchstone(session,
+                                     input[[it1]],
+                                     NULL,
+                                     input,
+                                     prefix)})
+
+    } else {
+      it1 <- sprintf("%s_touchstone1", prefix)
+      it2 <- sprintf("%s_touchstone2", prefix)
+
+      observeEvent(input[[it1]], {
+                   update_touchstone(session,
+                                      input[[it1]],
+                                      input[[it2]],
+                                      input,
+                                      prefix)})
+
+      observeEvent(input[[it2]], {
+                   update_touchstone(session,
+                                     input[[it1]],
+                                     input[[it2]],
+                                     input,
+                                     prefix)})
+    }
+
+    id <- sprintf("%s_disease", prefix)
+    observeEvent(input[[id]], {
+                 update_disease(session,
+                                input[[it1]],
+                                if (is.null(it2)) NULL else input[[it2]],
+                                input[[id]],
+                                input,
+                                prefix)})
+
+    if (n_group == 1) {
+      ig1 <- sprintf("%s_group", prefix)
+      ig2 <- NULL
+      observeEvent(input[[ig1]], {
+                   update_group(session,
+                                input[[it1]],
+                                if (is.null(it2)) NULL else input[[it2]],
+                                input[[id]],
+                                input[[ig1]],
+                                NULL,
+                                input,
+                                prefix)})
+
+    } else {
+      ig1 <- sprintf("%s_group1", prefix)
+      ig2 <- sprintf("%s_group2", prefix)
+      observeEvent(input[[ig1]], {
+                   update_group(session,
+                                input[[it1]],
+                                if (is.null(it2)) NULL else input[[it2]],
+                                input[[id]],
+                                input[[ig1]],
+                                input[[ig2]],
+                                input,
+                                prefix)})
+
+      observeEvent(input[[ig2]], {
+                   update_group(session,
+                                input[[it1]],
+                                if (is.null(it2)) NULL else input[[it2]],
+                                input[[id]],
+                                input[[ig1]],
+                                input[[ig2]],
+                                input,
+                                prefix)})
+    }
+
+    if (n_scenario == 1) {
+      is1 <- sprintf("%s_scenario", prefix)
+      is2 <- NULL
+      observeEvent(input[[is1]], {
+                   update_scenario(session,
+                                   input[[it1]],
+                                   if (is.null(it2)) NULL else input[[it2]],
+                                   input[[id]],
+                                   input[[ig1]],
+                                   if (is.null(ig2)) NULL else input[[ig2]],
+                                   input[[is1]],
+                                   NULL,
+                                   input, prefix)})
+    } else {
+      is1 <- sprintf("%s_scenario1", prefix)
+      is2 <- sprintf("%s_scenario2", prefix)
+      observeEvent(input[[is1]], {
+                   update_scenario(session,
+                                   input[[it1]],
+                                   if (is.null(it2)) NULL else input[[it2]],
+                                   input[[id]],
+                                   input[[ig1]],
+                                   if (is.null(ig2)) NULL else input[[ig2]],
+                                   input[[is1]],
+                                   input[[is2]],
+                                   input,
+                                   prefix)})
+
+      observeEvent(input[[is2]], {
+                   update_scenario(session,
+                                   input[[it1]],
+                                   if (is.null(it2)) NULL else input[[it2]],
+                                   input[[id]],
+                                   input[[ig1]],
+                                   if (is.null(ig2)) NULL else input[[ig2]],
+                                   input[[is1]],
+                                   input[[is2]],
+                                   input,
+                                   prefix)})
+    }
+
+    ic <- sprintf("%s_country", prefix)
+    observeEvent(input[[ic]], {
+                 update_country(session,
+                                input[[it1]],
+                                if (is.null(it2)) NULL else input[[it2]],
+                                input[[id]],
+                                input[[ig1]],
+                                if (is.null(ig2)) NULL else input[[ig2]],
+                                input[[is1]],
+                                if (is.null(is2)) NULL else input[[is2]],
+                                input[[ic]],
+                                input,
+                                prefix)})
+
+    n_graphs <- 1 + ((n_touchstone * n_group) > 1)
+    if (n_graphs == 1) {
+      graphs <- sprintf("%s_main_plot", prefix)
+    } else {
+      graphs <- sprintf("%s_main_plot%d", prefix, 1:2)
+    }
+    button <- sprintf("%s_plot_btn", prefix)
+
+    plot_reactive <- eventReactive(input[[button]], {
+      ages <- NULL
+      if (input[[sprintf("%s_ages", prefix)]] == "Under 5") {
+        ages <- 0:4
+      }
+
+      check <- function(pre, x1, x2, type) {
+        if (is.null(x2)) return(TRUE)
+        if ((grepl(pre, prefix)) && (input[[x1]] == input[[x2]])) {
+          showModal(modalDialog(
+            title = "Warning",
+            sprintf("Same %s selected. Dull comparison.", type),
+            easyClose = TRUE,
+            footer = modalButton("OK")
+          ))
+          return(FALSE)
+        }
+        TRUE
+      }
+
+      if (!check("mg", ig1, ig2, "group")) return(NULL)
+      if (!check("ts", it1, it2, "touchstone")) return(NULL)
+      if (!check("i", is1, is2, "scenario")) return(NULL)
+
+      list(
+        ages = ages,
+        opts = input[[sprintf("%s_options", prefix)]],
+
+        touchstone = c(input[[it1]], if (!is.null(it2)) input[[it2]] else NULL),
+        group      = c(input[[ig1]], if (!is.null(ig2)) input[[ig2]] else NULL),
+        scenario   = c(input[[is1]], if (!is.null(is2)) input[[is2]] else NULL),
+
+        disease = input[[id]],
+        country = input[[ic]],
+        outcome = input[[sprintf("%s_outcome", prefix)]],
+        year    = input[[sprintf("%s_year", prefix)]]
+      )
+    })
+
+    for (g in seq_along(graphs)) {
+      local({
+        gg <- g
+
+        output[[graphs[gg]]] <- renderPlot({
+          pr <- plot_reactive()
+          req(pr)
+
+          it <- pr$touchstone[min(gg, length(pr$touchstone))]
+          ig <- pr$group[min(gg, length(pr$group))]
+
+          stoner::stone_stochastic_graph(
+            base = data_dir,
+            touchstone = it,
+            disease = pr$disease,
+            group = ig,
+            country = pr$country,
+            scenario = pr$scenario[1],
+            scenario2 = if (length(pr$scenario) > 1) pr$scenario[2] else NULL,
+            outcome = pr$outcome,
+            by_cohort = pr$year == "Cohort",
+            ages = pr$ages,
+            log = "Log-Y" %in% pr$opts,
+            include_median = "Median" %in% pr$opts,
+            include_quantiles = "Quantiles" %in% pr$opts,
+            include_mean = "Mean" %in% pr$opts
+          )
+        })
+      })
+    }
+  }
+  add_observers("b")
+  add_observers("bts", 2, 1, 1)
+  add_observers("bmg", 1, 2, 1)
+  add_observers("i", 1, 1, 2)
+  add_observers("its", 2, 1, 2)
+  add_observers("img", 1, 2, 2)
+}
+
+shinyApp(app_ui, app_server)
diff --git a/inst/app/utils.R b/inst/app/utils.R
new file mode 100644
index 00000000..daaa9cc3
--- /dev/null
+++ b/inst/app/utils.R
@@ -0,0 +1,180 @@
+# Helper functions to dig information from the
+# shared file system.
+
+# Functions ending in _ts will find things common
+# between two touchstones, for comparisons.
+
+# Functions ending in _mg will find things common
+# between to modelling groups, for comparisons.
+
+get_touchstones <- function() {
+  sort(unique(basename(list.dirs(data_dir, recursive = FALSE))))
+}
+
+get_diseases <- function(touchstone1, touchstone2) {
+  lookup <- function(touchstone, res = NULL) {
+    if (is.null(touchstone)) return(res)
+    path <- file.path(data_dir, touchstone)
+    dirs <- basename(list.dirs(path, recursive = FALSE))
+    unique(unlist(lapply(strsplit(dirs, "_"), `[[`, 1)))
+  }
+  res <- lookup(touchstone1)
+  sort(unique(res[res %in% lookup(touchstone2, res)]))
+}
+
+get_groups <- function(touchstone1, touchstone2, disease) {
+  lookup <- function(touchstone, disease, res = NULL) {
+    if (is.null(touchstone)) return(res)
+    path <- file.path(data_dir, touchstone)
+    dirs <- basename(list.dirs(path, recursive = FALSE))
+    dirs <- dirs[substr(dirs, 1, nchar(disease) + 1) == paste0(disease, "_")]
+    unique(substring(dirs, nchar(disease) + 2))
+  }
+  res <- lookup(touchstone1, disease)
+  sort(unique(res[res %in% lookup(touchstone2, disease, res)]))
+}
+
+get_scenarios <- function(touchstone1, touchstone2, disease, group1, group2) {
+  lookup <- function(touchstone, disease, group, res = NULL) {
+    if ((is.null(touchstone)) || (is.null(group))) return(res)
+    path <- file.path(data_dir, touchstone, paste(disease, group, sep = "_"))
+    files <- basename(list.files(path, recursive = FALSE))
+    files <- gsub(".pq", "", files)
+    files <- substring(files, nchar(group) + 2)
+    ends <- unlist(lapply(gregexpr("_", files), `[[`, 1)) - 1
+    unique(substring(files, 1, ends))
+  }
+
+  res <- lookup(touchstone1, disease, group1)
+  res <- res[res %in% lookup(touchstone2, disease, group1, res)]
+  res <- res[res %in% lookup(touchstone1, disease, group2, res)]
+  sort(unique(res[res %in% lookup(touchstone2, disease, group2, res)]))
+}
+
+get_countries <- function(touchstone1, touchstone2, disease, group1, group2,
+                          scenario1, scenario2) {
+  lookup <- function(touchstone, disease, group, scenario, res = NULL) {
+    if ((is.null(touchstone)) || (is.null(group)) || (is.null(scenario))) {
+      return(res)
+    }
+
+    path <- file.path(data_dir, touchstone, paste(disease, group, sep = "_"))
+    files <- basename(list.files(path, recursive = FALSE))
+    files <- gsub(".pq", "", files)
+    files <- substring(files, nchar(group) + 2)
+    files <- files[substring(files, 1, nchar(scenario) + 1) == paste0(scenario, "_")]
+    ends <- unlist(lapply(gregexpr("_", files), `[[`, 1)) + 1
+    unique(substring(files, ends))
+  }
+  res <- lookup(touchstone1, disease, group1, scenario1)
+  res <- res[res %in% lookup(touchstone2, disease, group1, scenario1, res)]
+  res <- res[res %in% lookup(touchstone1, disease, group2, scenario1, res)]
+  res <- res[res %in% lookup(touchstone2, disease, group2, scenario1, res)]
+  res <- res[res %in% lookup(touchstone1, disease, group1, scenario2, res)]
+  res <- res[res %in% lookup(touchstone2, disease, group1, scenario2, res)]
+  res <- res[res %in% lookup(touchstone1, disease, group2, scenario2, res)]
+  sort(unique(res[res %in% lookup(touchstone2, disease, group2, scenario2, res)]))
+}
+
+get_outcomes <- function(touchstone1, touchstone2, disease, group1, group2,
+                         scenario1, scenario2, country) {
+  lookup <- function(touchstone, disease, group, scenario, country, res = NULL) {
+    if ((is.null(touchstone)) || (is.null(group)) || (is.null(scenario))) {
+      return(res)
+    }
+    path <- file.path(data_dir, touchstone, paste(disease, group, sep = "_"))
+    thefile <- sprintf("%s_%s_%s.pq", group, scenario, country)
+    tbl <- arrow::read_parquet(file.path(path, thefile), col_select = NULL)
+    cols <- names(tbl)
+    sort(cols[!cols %in% c("disease", "run_id", "year", "age",
+                           "country", "cohort_size")])
+  }
+  res <- lookup(touchstone1, disease, group1, scenario1, country)
+  res <- res[res %in% lookup(touchstone2, disease, group1, scenario1, country, res)]
+  res <- res[res %in% lookup(touchstone1, disease, group2, scenario1, country, res)]
+  res <- res[res %in% lookup(touchstone2, disease, group2, scenario1, country, res)]
+  res <- res[res %in% lookup(touchstone1, disease, group1, scenario2, country, res)]
+  res <- res[res %in% lookup(touchstone2, disease, group1, scenario2, country, res)]
+  res <- res[res %in% lookup(touchstone1, disease, group2, scenario2, country, res)]
+  sort(unique(res[res %in% lookup(touchstone2, disease, group2, scenario2, country, res)]))
+}
+
+# GUI helpers.
+
+# Update a dropdown. Keep the previously selected value if possible,
+# otherwise select the first one.
+
+update_dropdown_keep <- function(session, id, choices, selected) {
+  if ((is.null(selected)) ||
+      (!selected %in% choices)) selected <- choices[[1]]
+  updateSelectInput(session, id, choices = choices, selected = selected)
+  selected
+}
+
+# The next events recurse through the dropdowns updating available
+# options. In a messy attempt to avoid duplication, the prefix is
+# "b" for the simple burden panel, "i" for impact,
+# "bts" for multi-touchstone burden, "its" for multi-touchstone impact
+# "bmg" for multi-group burden, and "img" for multi-group impact.
+
+# It's messy as sometimes we have two touchstones, groups, or scenarios.
+
+update_touchstone <- function(session, ts1, ts2, input, prefix) {
+  d <- paste(prefix, "disease", sep = "_")
+  diseases <- get_diseases(ts1, ts2)
+  disease <- update_dropdown_keep(session, d, diseases, input[[d]])
+  update_disease(session, ts1, ts2, disease, input, prefix)
+}
+
+update_disease <- function(session, ts1, ts2, disease, input, prefix) {
+  if (is.null(disease) || (disease == "")) return()
+  groups <- get_groups(ts1, ts2, disease)
+  if (grepl("mg", prefix)) {
+    g1 <- paste(prefix, "group1", sep = "_")
+    g2 <- paste(prefix, "group2", sep = "_")
+    group1 <- update_dropdown_keep(session, g1, groups, input[[g1]])
+    group2 <- update_dropdown_keep(session, g2, groups, input[[g2]])
+    update_group(session, ts1, ts2, disease, group1, group2, input, prefix)
+  } else {
+    g <- paste(prefix, "group", sep = "_")
+    group1 <- update_dropdown_keep(session, g, groups, input[[g]])
+    group2 <- NULL
+  }
+  update_group(session, ts1, ts2, disease, group1, group2, input, prefix)
+}
+
+update_group <- function(session, ts1, ts2, disease, g1, g2, input, prefix) {
+  if (is.null(g1) || (g1 == "")) return()
+  scenarios <- get_scenarios(ts1, ts2, disease, g1, g2)
+  if (grepl("i", prefix)) {
+    sc <- sprintf("%s_scenario1", prefix)
+    s1 <- update_dropdown_keep(session, sc, scenarios, input[[sc]])
+    sc <- sprintf("%s_scenario2", prefix)
+    s2 <- update_dropdown_keep(session, sc, scenarios, input[[sc]])
+  } else {
+    sc <- sprintf("%s_scenario", prefix)
+    s1 <- update_dropdown_keep(session, sc, scenarios, input[[sc]])
+    s2 <- NULL
+  }
+  update_scenario(session, ts1, ts2, disease, g1, g2, s1, s2, input, prefix)
+}
+
+update_scenario <- function(session, ts1, ts2, disease, g1, g2, s1, s2,
+                            input, prefix) {
+  if (is.null(s1) || (s1 == "")) return()
+
+  cc <- paste(prefix, "country", sep = "_")
+  countries <- get_countries(ts1, ts2, disease, g1, g2, s1, s2)
+  country <- update_dropdown_keep(session, cc, countries, input[[cc]])
+  update_country(session, ts1, ts2, disease, g1, g2, s1, s2, country,
+                 input, prefix)
+}
+
+update_country <- function(session, ts1, ts2, disease, g1, g2, s1, s2, country,
+                           input, prefix) {
+  if (is.null(country) || (country == "")) return()
+  if (is.null(g1) || (g1 == "")) return()
+  o <- paste(prefix, "outcome", sep = "_")
+  outcomes <- get_outcomes(ts1, ts2, disease, g1, g2, s1, s2, country)
+  update_dropdown_keep(session, o, outcomes, input[[o]])
+}
diff --git a/man/fetch_packit.Rd b/man/fetch_packit.Rd
new file mode 100644
index 00000000..e77e94f2
--- /dev/null
+++ b/man/fetch_packit.Rd
@@ -0,0 +1,29 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/packit.R
+\name{fetch_packit}
+\alias{fetch_packit}
+\title{Fetch packet from a packit server}
+\usage{
+fetch_packit(
+  packet_id,
+  filename,
+  server = "https://montagu.vaccineimpact.org/packit/api/"
+)
+}
+\arguments{
+\item{packet_id}{The id of the packet containing the artefact.}
+
+\item{filename}{The filename of the file within the packet.}
+
+\item{server}{By default, the URL to the packit API on Montagu,
+but this can be set to other packit API's if we want.}
+}
+\value{
+The filename of the temporary file which has been downloaded.
+}
+\description{
+Download a file from a packit, for example an artefact containing
+central burden estimates, which we could then plot on top of
+stochastics. But more generally, we can use this for fetching any
+file from packet (ie, from the Montagu Reporting Portal).
+}
diff --git a/man/stochastic_explorer.Rd b/man/stochastic_explorer.Rd
new file mode 100644
index 00000000..f3683ca2
--- /dev/null
+++ b/man/stochastic_explorer.Rd
@@ -0,0 +1,19 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/stochastic_graphs.R
+\name{stochastic_explorer}
+\alias{stochastic_explorer}
+\title{Stochastic plot}
+\usage{
+stochastic_explorer(data_dir = "//wpia-hn2.hpc.dide.ic.ac.uk/vimc_stochastics")
+}
+\arguments{
+\item{data_dir}{The location of the standardised stochastic folder
+hierarchy; this can be a local path, a fully-qualified network path on
+windows, or a mount point on linux or Mac.}
+}
+\description{
+Launch a Shiny app to allow interactive plotting of
+standardised stochastic data, burden estimates,
+impacts, comparisons between touchstones, and
+comparisons between modelling groups.
+}
diff --git a/man/stone_stochastic_graph.Rd b/man/stone_stochastic_graph.Rd
index 9fa15c93..5e40693b 100644
--- a/man/stone_stochastic_graph.Rd
+++ b/man/stone_stochastic_graph.Rd
@@ -14,7 +14,13 @@ stone_stochastic_graph(
   outcome,
   ages = NULL,
   by_cohort = FALSE,
-  log = FALSE
+  log = FALSE,
+  packit_id = NULL,
+  packit_file = NULL,
+  include_quantiles = TRUE,
+  include_mean = TRUE,
+  include_median = TRUE,
+  scenario2 = NULL
 )
 }
 \arguments{
@@ -40,6 +46,26 @@ if left as NULL, then all ages are used and aggregated.}
 year of birth before aggregating.}
 
 \item{log}{If TRUE, then use a logged y-axis.}
+
+\item{packit_id}{If set, then read central burden estimates from a file
+within a packit on the Montagu packit server.}
+
+\item{packit_file}{Used with packit_id to specify the filename of an RDS
+file providing burden estimates. We expect to find scenario, year, age,
+country, burden_outcome and value fields in the table.}
+
+\item{include_quantiles}{Default TRUE, select whether to plot the
+5\% and 95\% quantile lines.}
+
+\item{include_mean}{Default TRUE, select whether to plot the mean.}
+
+\item{include_median}{Default TRUE, select whether to plot the median.}
+
+\item{scenario2}{Default NULL; if set, then the burdens from this
+scenario will be subtracted from those in \code{scenario} - ie, this plots
+an impact graph of applying the second scenario. For many graphs that
+use this, the result will be positive numbers, representing cases
+or deaths averted.}
 }
 \description{
 Draw a stochastic plot showing all the different runs, with the mean,
diff --git a/man/stone_stochastic_standardise.Rd b/man/stone_stochastic_standardise.Rd
index 343275fd..d010e207 100644
--- a/man/stone_stochastic_standardise.Rd
+++ b/man/stone_stochastic_standardise.Rd
@@ -12,7 +12,12 @@ stone_stochastic_standardise(
   files,
   index = 1,
   rubella_fix = TRUE,
-  missing_run_id_fix = TRUE
+  hepb2019_fix = TRUE,
+  hib2019_fix = TRUE,
+  missing_run_id_fix = TRUE,
+  allow_missing_yll = TRUE,
+  allow_missing_dalys = TRUE,
+  allow_missing_indexes = FALSE
 )
 }
 \arguments{
@@ -47,10 +52,30 @@ these to the simpler names. Processing Rubella stochastic files without
 this set to TRUE will fail - so while we should always do this, keeping
 the parameter makes it more clear in the code what we're doing and why.}
 
+\item{hepb2019_fix}{In 2019 (and 2017), HepB deaths and cases were subdivided into
+number of different causes. This flag combines those into the single
+appropriate burden outcome.}
+
+\item{hib2019_fix}{In 2019 (and 2017), Hib deaths and cases were subdivided into
+number of different causes. This flag combines those into the single
+appropriate burden outcome.}
+
 \item{missing_run_id_fix}{Some groups in the past have omitted run_id
 from the files, but included them in the filenames. This fix inserts
 that into the files if the index parameter indicates we have 200 runs to
 process.}
+
+\item{allow_missing_yll}{yll was introduced in 2023? This flag allows
+it to be missing for processing older stochastics.}
+
+\item{allow_missing_dalys}{Some early groups did not provide dalys; this
+flag allows dalys to be skipped.}
+
+\item{allow_missing_indexes}{In some early runs, different groups
+provided different numbers of files for each scenario, because some
+countries did not implement particular coverage campaigns. This
+flag needs to be TRUE for those groups, but the default is FALSE,
+since it's rare, and we generally want errors for missing files.}
 }
 \description{
 Convert a modelling group's stochastic files into an intermediate
diff --git a/pkgdown/extra.css b/pkgdown/extra.css
index a5229e31..4a0ae85a 100644
--- a/pkgdown/extra.css
+++ b/pkgdown/extra.css
@@ -1,3 +1,7 @@
 .row > main {
   hyphens: none;
 }
+
+pre {
+  line-height: 1;
+}
diff --git a/scripts/process_stoch_201710gavi.R b/scripts/process_stoch_201710gavi.R
new file mode 100644
index 00000000..bcb89f47
--- /dev/null
+++ b/scripts/process_stoch_201710gavi.R
@@ -0,0 +1,298 @@
+base_in_path <- "//wpia-hn2.hpc.dide.ic.ac.uk/vimc_stochastics_dropbox/latest/201710gavi"
+base_out_path <- "//wpia-hn2.hpc.dide.ic.ac.uk/vimc_stochastics/201710gavi"
+
+# In case useful for looking up scenarios
+
+vault <- vaultr::vault_client(login = "github")
+password <- vault$read("/secret/vimc/database/production/users/readonly")$password
+con <- DBI::dbConnect(RPostgres::Postgres(),
+                      dbname = "montagu",
+                      host = "montagu.vaccineimpact.org",
+                      port = 5432, password = password,
+                      user = "readonly")
+
+fetch_scenarios <- function(disease, touchstone='201710gavi-5') {
+  sort(DBI::dbGetQuery(con, "
+    SELECT scenario_description FROM
+           scenario JOIN scenario_description
+           ON scenario.scenario_description = scenario_description.id
+           WHERE disease = $1 AND touchstone = $2", list(disease, touchstone))$scenario_description)
+}
+
+# Let's unleash the cow
+
+setwd("Q:/testcow")
+writeLines("vimc/stoner@VIMC-9230", "pkgdepends.txt")
+hipercow::hipercow_init(driver = "dide-windows")
+hipercow::hipercow_provision()
+# Network/memory might be too much for more than a job per node.
+hres <- hipercow::hipercow_resources(cores = 1L, exclusive = TRUE)
+
+
+###############
+# HepB
+
+scenarios <- c("hepb-no-vaccination",
+               "hepb-hepb-routine-with",
+               "hepb-bd-routine-with",
+               "hepb-bd-routine-with-hepb-routine-with",
+               "hepb-bd-routine-best-hepb-routine-with")
+
+stoner::stone_stochastic_standardise(
+  group = "CDA-Razavi",
+  in_path = file.path(base_in_path, "CDA-Razavi"),
+  out_path = file.path(base_out_path, "HepB_CDA-Razavi"),
+  scenarios = scenarios,
+  files = "Devin Razavi-Shearer - stochastic_burden_template_HepB-CDA-Razavi_:scenario_:index.csv.xz",
+  index = 1:98,
+  allow_missing_indexes = TRUE)
+
+hipercow::task_create_expr(resources = hres, expr =
+  stoner::stone_stochastic_standardise(
+    group = "Li",
+    in_path = file.path(base_in_path, "Li"),
+    out_path = file.path(base_out_path, "HepB_Li"),
+    scenarios = scenarios,
+    files = ":scenario:index.csv.xz",
+    index = 1:200))
+
+hipercow::task_create_expr(resources = hres, expr =
+  stoner::stone_stochastic_standardise(
+    group = "IC-Hallett",
+    in_path = file.path(base_in_path, "IC-Hallett"),
+    out_path = file.path(base_out_path, "HepB_IC-Hallett"),
+    scenarios = scenarios,
+    files = "stochastic_burden_est_:scenario_:index.csv.xz",
+    index = 1:200))
+
+
+###############
+# HIB
+
+scenarios <- c("hib-no-vaccination", "hib-routine-gavi")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-Hib"),
+  out_path = file.path(base_out_path, "Hib_JHU-Tam"),
+  scenarios = scenarios,
+  files = c("No_Hib_all.csv.xz", "With_Hib_all.csv.xz"))
+
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Clark",
+  in_path = file.path(base_in_path, "LSHTM-Clark-Hib"),
+  out_path = file.path(base_out_path, "Hib_LSHTM-Clark"),
+  scenarios = scenarios,
+  files = c("Kaja Abbas - stochastic-burden-estimate.201710gavi-5.Hib_LSHTM-CLark_standard-Hib-no-vaccination.csv.xz",
+            "Kaja Abbas - stochastic-burden-estimate.201710gavi-5.Hib_LSHTM-CLark_standard-Hib-routine-gavi.csv.xz"))
+
+
+###############
+# HPV
+
+scenarios <- c("hpv-no-vaccination", "hpv-routine-gavi",
+               "hpv-campaign-gavi")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "Harvard-Sweet",
+  in_path = file.path(base_in_path, "Harvard-Sweet"),
+  out_path = file.path(base_out_path, "HPV_Harvard-Sweet"),
+  scenarios = scenarios,
+  files = c("1. central_burden_template_HPV-Harvard-Sweet No Vaccine_PSA-:index_08.22.19.csv.xz",
+            "2. central_burden_template_HPV-Harvard-Sweet Routine_PSA-:index_09.09.19.csv.xz",
+            "3. central_burden_template_HPV-Harvard-Sweet Campaign_PSA-:index_08.22.19.csv.xz"),
+  index = 1:200))
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Jit",
+  in_path = file.path(base_in_path, "LSHTM-Jit-HPV"),
+  out_path = file.path(base_out_path, "HPV_LSHTM-Jit"),
+  scenarios =  scenarios,
+  files = "Kaja Abbas - stochastic-burden-estimate.201710gavi-5.HPV_LSHTM-Jit_standard-:scenario.csv.xz"))
+
+###############
+# JE
+
+scenarios <- c("je-routine-no-vaccination", "je-campaign-gavi",
+               "je-routine-gavi")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "OUCRU-Clapham",
+  in_path = file.path(base_in_path, "OUCRU-Clapham"),
+  out_path = file.path(base_out_path, "JE_OUCRU-Clapham"),
+  scenarios =  scenarios,
+  files = "Duy Nguyen - stochastic-burden-template.201710gavi-6.JE_OUCRU-Clapham_standard_:scenario_:index.csv.xz",
+  index = 1:200))
+
+stub <- "Sean Moore - stochastic_burden_est_JE_UND-Moore"
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "UND-Moore",
+  in_path = file.path(base_in_path, "UND-Moore"),
+  out_path = file.path(base_out_path, "JE_UND-Moore"),
+  scenarios =  scenarios,
+  files = c("Sean Moore - stochastic_burden_est_JE_UND-Moore_je-no-vaccination.csv.xz",
+            "Sean Moore - stochastic_burden_est_JE_UND-Moore_je-campaign-gavi.csv.xz",
+            "Sean Moore - stochastic_burden_est_JE_UND-Moore_je-routine-gavi.csv.xz")))
+
+
+
+###############
+# Measles
+scenarios <- c("measles-no-vaccination", "measles-mcv1-gavi",
+               "measles-mcv2-gavi", "measles-campaign-gavi")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "PSU-Ferrari",
+  in_path = file.path(base_in_path, "PSU-Ferrari"),
+  out_path = file.path(base_out_path, "Measles_PSU-Ferrari"),
+  scenarios = scenarios,
+  files = c("Matthew Ferrari - stochastic_burden_novax_Measles-PSU-Ferrari-5.csv.xz",
+            "Matthew Ferrari - stochastic_burden_mcv1_Measles-PSU-Ferrari-5.csv.xz",
+            "Matthew Ferrari - stochastic_burden_mcv2_Measles-PSU-Ferrari-5.csv.xz",
+            "Matthew Ferrari - stochastic_burden_sia_Measles-PSU-Ferrari-5.csv.xz")))
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Jit",
+  in_path = file.path(base_in_path, "LSHTM-Jit-Measles"),
+  out_path = file.path(base_out_path, "Measles_LSHTM-Jit"),
+  scenarios = scenarios,
+  files = c("stochastic_burden_template_Measles-LSHTM-Jit-no_vaxx.csv.xz",
+            "Petra Klepac - stochastic_burden_template_Measles-LSHTM-Jit-mcv1-only.csv.xz",
+            "stochastic_burden_template_Measles-LSHTM-Jit-mcv1_mcv2.csv.xz",
+            "stochastic_burden_template_Measles-LSHTM-Jit-mcv1-mcv2-campaigns.csv.xz")))
+
+###############
+# MenA
+
+scenarios = c("mena-no-vaccination", "mena-campaign-gavi",
+              "mena-routine-gavi")
+
+stoner::stone_stochastic_standardise(
+  group = "Cambridge-Trotter",
+  in_path = file.path(base_in_path, "Cambridge-Trotter"),
+  out_path = file.path(base_out_path, "MenA_Cambridge-Trotter"),
+  scenarios = scenarios,
+  files = c("Cambridge-Trotter-mena-no-vacc.csv.xz",
+            "Cambridge-Trotter-mena-campaign.csv.xz",
+            "Cambridge-Trotter-mena-routine.csv.xz"))
+
+stoner::stone_stochastic_standardise(
+  group = "KPW-Jackson",
+  in_path = file.path(base_in_path, "KPW-Jackson"),
+  out_path = file.path(base_out_path, "MenA_KPW-Jackson"),
+  scenarios = scenarios,
+  files = c("Mike Jackson - stochastic_burden_est_MenA-KPW-Jackson_mena-no-vaccination_:index.csv.xz",
+            "Mike Jackson - stochastic_burden_est_MenA-KPW-Jackson_mena-campaign-gavi_:index.csv.xz",
+            "Mike Jackson - stochastic_burden_est_MenA-KPW-Jackson_mena-routine-gavi_:index.csv.xz"),
+  index = 1:2,
+  allow_missing_indexes = TRUE)
+
+###############
+# PCV
+
+scenarios <- c("pcv-no-vaccination", "pcv-routine-gavi")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-PCV"),
+  out_path = file.path(base_out_path, "PCV_JHU-Tam"),
+  scenarios = scenarios,
+  files = c("No_PCV_all.csv.xz", "With_PCV_all.csv.xz"))
+
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Clark",
+  in_path = file.path(base_in_path, "LSHTM-Clark-PCV"),
+  out_path = file.path(base_out_path, "PCV_LSHTM-Clark"),
+  scenarios = scenarios,
+  files = c("Kaja Abbas - stochastic-burden-estimate.201710gavi-5.PCV_LSHTM-CLark_standard-PCV-no-vaccination.csv.xz",
+            "Kaja Abbas - stochastic-burden-estimate.201710gavi-5.PCV_LSHTM-CLark_standard-PCV-routine-gavi.csv.xz"))
+
+###############
+# Rota
+
+scenarios <- c("rota-no-vaccination", "rota-routine-gavi")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-Rota"),
+  out_path = file.path(base_out_path, "Rota_JHU-Tam"),
+  scenarios = scenarios,
+  files = c("No_Rota_all.csv.xz", "With_Rota_all.csv.xz"))
+
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Clark",
+  in_path = file.path(base_in_path, "LSHTM-Clark-Rota"),
+  out_path = file.path(base_out_path, "Rota_LSHTM-Clark"),
+  scenarios = scenarios,
+  files = c("Kaja Abbas - stochastic-burden-estimate.201710gavi-5.Rota_LSHTM-CLark_standard-Rota-no-vaccination.csv.xz",
+            "Kaja Abbas - stochastic-burden-estimate.201710gavi-5.Rota_LSHTM-CLark_standard-Rota-routine-gavi.csv.xz"))
+
+###############
+# Rubella
+
+scenarios <- c(
+  "rubella-routine-no-vaccination", "rubella-campaign-gavi",
+  "rubella-rcv1-gavi", "rubella-rcv2-gavi")
+
+#hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "JHU-Lessler",
+  in_path = file.path(base_in_path, "JHU-Lessler"),
+  out_path = file.path(base_out_path, "Rubella_JHU-Lessler"),
+  scenarios = scenarios,
+  files = c("Amy Winter - stochastic_burden_est-rubella-no-vaccination_:index.csv.xz",
+            "Amy Winter - stochastic_burden_est-rubella-campaign-gavi_:index.csv.xz",
+            "Amy Winter - stochastic_burden_est-rubella-routine-gavi_:index.csv.xz",
+            "Amy Winter - stochastic_burden_est-rubella-rcv2-gavi_:index.csv.xz"),
+  index = 1:10)
+
+scenarios <- c(
+  "rubella-routine-no-vaccination", "rubella-campaign-gavi",
+  "rubella-rcv2-gavi")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "PHE-Vynnycky",
+  in_path = file.path(base_in_path, "PHE-Vynnycky"),
+  out_path = file.path(base_out_path, "Rubella_PHE-Vynnycky"),
+  scenarios = scenarios,
+  files = c("stochastic_burden_estimate_Rubella-PHE-Vynnycky_rubella-routine-no-vaccination_:index.csv.xz",
+            "stochastic_burden_estimate_Rubella-PHE-Vynnycky_rubella-campaign-gavi_:index.csv.xz",
+            "stochastic_burden_estimate_Rubella-PHE-Vynnycky_rubella-routine-gavi_:index.csv.xz"),
+  index = 1:200))
+
+
+###############
+# YF
+
+scenarios <- c("yf-no-vaccination", "yf-preventive-gavi",
+               "yf-routine-gavi")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "IC-Garske",
+  in_path = file.path(base_in_path, "IC-Garske"),
+  out_path = file.path(base_out_path, "YF_IC-Garske"),
+  scenarios = scenarios,
+  files = c("Katy Gaythorpe - NEW_NEW_stochastic_burden_estimate_YF_no-vaccination_:index.csv.xz",
+            "Katy Gaythorpe - NEW_NEW_stochastic_burden_estimate_YF_preventive-gavi_:index.csv.xz",
+            "Katy Gaythorpe - NEW_NEW_stochastic_burden_estimate_YF_routine-gavi_:index.csv.xz"),
+  index = 1:200))
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "UND-Perkins",
+  in_path = file.path(base_in_path, "UND-Perkins"),
+  out_path = file.path(base_out_path, "YF_UND-Perkins"),
+  scenarios = scenarios,
+  files = c("John Huber - stochastic_burden_est_YF_UND-Perkins_yf-no-vaccination-gavi_:index.csv.xz",
+            "John Huber - stochastic_burden_est_YF_UND-Perkins_yf-preventive-gavi_:index.csv.xz",
+            "John Huber - stochastic_burden_est_YF_UND-Perkins_yf-routine-gavi_:index.csv.xz"),
+  index = 1:200))
diff --git a/scripts/process_stoch_201910gavi.R b/scripts/process_stoch_201910gavi.R
new file mode 100644
index 00000000..1ed13183
--- /dev/null
+++ b/scripts/process_stoch_201910gavi.R
@@ -0,0 +1,419 @@
+base_in_path <- "//wpia-hn2.hpc.dide.ic.ac.uk/vimc_stochastics_dropbox/latest/201910gavi"
+base_out_path <- "//wpia-hn2.hpc.dide.ic.ac.uk/vimc_stochastics/201910gavi"
+
+# In case useful for looking up scenarios
+
+vault <- vaultr::vault_client(login = "github")
+password <- vault$read("/secret/vimc/database/production/users/readonly")$password
+con <- DBI::dbConnect(RPostgres::Postgres(),
+                      dbname = "montagu",
+                      host = "montagu.vaccineimpact.org",
+                      port = 5432, password = password,
+                      user = "readonly")
+
+fetch_scenarios <- function(disease) {
+  sort(DBI::dbGetQuery(con, "
+    SELECT scenario_description FROM
+           scenario JOIN scenario_description
+           ON scenario.scenario_description = scenario_description.id
+           WHERE disease = $1 AND touchstone='201910gavi-5'", disease)$scenario_description)
+}
+
+# Let's unleash the cow
+
+setwd("Q:/testcow")
+writeLines("vimc/stoner@VIMC-9230", "pkgdepends.txt")
+hipercow::hipercow_init(driver = "dide-windows")
+hipercow::hipercow_provision()
+# Network/memory might be too much for more than a job per node.
+hres <- hipercow::hipercow_resources(cores = 1L, exclusive = TRUE)
+
+###############
+# Cholera
+
+scenarios = c("cholera-no-vaccination", "cholera-campaign-default")
+
+stoner::stone_stochastic_standardise(
+  group = "IVI-Kim",
+  in_path = file.path(base_in_path, "IVI-Kim-Cholera"),
+  out_path = file.path(base_out_path, "Cholera_IVI-Kim"),
+  scenarios = scenarios,
+  files = c("Jong-Hoon Kim - stoch_output_Cholera_novacc_20210902.csv.xz",
+            "Jong-Hoon Kim - stoch_output_Cholera_campaign_20210902.csv.xz")
+  )
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Lee",
+  in_path = file.path(base_in_path, "JHU-Lee"),
+  out_path = file.path(base_out_path, "Cholera_JHU-Lee"),
+  scenarios = scenarios,
+  files = c("Kaiyue Zou - stochastic-burden-template.201910gavi-5.Cholera_no-vaccination.csv.xz",
+            "Kaiyue Zou - stochastic-burden-template.201910gavi-5.Cholera_campaign-default.csv.xz")
+)
+
+###############
+# HepB
+
+scenarios <- c("hepb-no-vaccination",
+               "hepb-bd-default-hepb-routine-default",
+               "hepb-bd-routine-bestcase-hepb-routine-bestcase",
+               "hepb-bd-routine-bestcase",
+               "hepb-bd-routine-default",
+               "hepb-hepb-routine-bestcase",
+               "hepb-hepb-routine-default",
+               "hepb-stop")
+
+stub <- "Ivane Gamkrelidze - stochastic-burden-template.201910gavi-4.HepB_CDA-Razavi"
+
+stoner::stone_stochastic_standardise(
+  group = "CDA-Razavi",
+  in_path = file.path(base_in_path, "CDA-Razavi"),
+  out_path = file.path(base_out_path, "HepB_CDA-Razavi"),
+  scenarios = scenarios,
+  files =
+    c(sprintf("%s_all_hepb-no-vaccination.csv.xz", stub),
+      sprintf("%s_all_hepb-bd-default-hepb-routine-default.csv.xz", stub),
+      sprintf("%s_all_hepb-bd-routine-bestcase-hepb-routine-bestcase.csv.xz", stub),
+      sprintf("%s_bd_hepb-bd-routine-bestcase.csv.xz", stub),
+      sprintf("%s_bd_hepb-bd-routine-default.csv.xz", stub),
+      sprintf("%s_non_bd_hepb-hepb-routine-bestcase.csv.xz", stub),
+      sprintf("%s_non_bd_hepb-hepb-routine-default.csv.xz", stub),
+      sprintf("%s_all_hepb-stop.csv.xz", stub))
+)
+
+hipercow::task_create_expr(resources = hres, expr =
+  stoner::stone_stochastic_standardise(
+    group = "Li",
+    in_path = file.path(base_in_path, "Li"),
+    out_path = file.path(base_out_path, "HepB_Li"),
+    scenarios = scenarios,
+    files = ":scenario:index.csv.xz",
+    index = 1:200))
+
+hipercow::task_create_expr(resources = hres, expr =
+  stoner::stone_stochastic_standardise(
+    group = "IC-Hallett",
+    in_path = file.path(base_in_path, "IC-Hallett"),
+    out_path = file.path(base_out_path, "HepB_IC-Hallett"),
+    scenarios = scenarios,
+    files = "stochastic_burden_est_HepB-IC-Hallett_:scenario_:index.csv.xz",
+    index = 1:200))
+
+
+###############
+# HIB
+
+scenarios <- c("hib-no-vaccination", "hib-routine-default",
+               "hib-routine-bestcase")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-Hib"),
+  out_path = file.path(base_out_path, "Hib_JHU-Tam"),
+  scenarios = scenarios,
+  files = c("novac:index.csv.xz", "default:index.csv.xz", "best:index.csv.xz"),
+  index = 1:14)
+
+hipercow::task_create_expr(resources = hres, expr =
+  stoner::stone_stochastic_standardise(
+  group = "LSHTM-Clark",
+  in_path = file.path(base_in_path, "LSHTM-Clark_Hib"),
+  out_path = file.path(base_out_path, "Hib_LSHTM-Clark"),
+  scenarios = scenarios,
+  files = c("VIMC_Hib_PSA_NoVax.csv.xz", "VIMC_Hib_PSA_Default.csv.xz",
+            "VIMC_Hib_PSA_Best.csv.xz")))
+
+###############
+# HPV
+
+scenarios <- c("hpv-no-vaccination", "hpv-routine-bestcase",
+               "hpv-campaign-bestcase", "hpv-routine-default",
+               "hpv-campaign-default")
+
+stub <- "stochastic-burden-est.201910gavi-5.HPV_Harvard-Sweet"
+
+stoner::stone_stochastic_standardise(
+  group = "Harvard-Sweet",
+  in_path = file.path(base_in_path, "Harvard-Sweet"),
+  out_path = file.path(base_out_path, "HPV_Harvard-Sweet"),
+  scenarios = scenarios,
+  files = c(sprintf("%s_novacc_run_:index.csv.xz", stub),
+            sprintf("%s_routine-bestcase_run_:index.csv.xz", stub),
+            sprintf("%s_campaign-bestcase_run_:index.csv.xz", stub),
+            sprintf("%s_routine-default_run_:index.csv.xz", stub),
+            sprintf("%s_campaign-default_run_:index.csv.xz", stub)),
+  index = 1:200)
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Jit",
+  in_path = file.path(base_in_path, "LSHTM-Jit_HPV"),
+  out_path = file.path(base_out_path, "HPV_LSHTM-Jit"),
+  scenarios =  scenarios,
+  files =
+    c("stochastic-burden-novaccination_201910gavi-4_hpv-no-vaccination.csv.xz",
+      "stochastic-burden-vaccination_201910gavi-4_hpv-campaign-bestcase.csv.xz",
+      "stochastic-burden-vaccination_201910gavi-4_hpv-campaign-default.csv.xz",
+      "stochastic-burden-vaccination_201910gavi-4_hpv-routine-bestcase.csv.xz",
+      "stochastic-burden-vaccination_201910gavi-4_hpv-routine-default.csv.xz")))
+
+###############
+# JE
+
+scenarios <- c("je-routine-no-vaccination", "je-campaign-default",
+               "je-campaign-bestcase", "je-routine-default",
+               "je-routine-bestcase")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "OUCRU-Clapham",
+  in_path = file.path(base_in_path, "OUCRU-Clapham"),
+  out_path = file.path(base_out_path, "JE_OUCRU-Clapham"),
+  scenarios =  scenarios,
+  files = c(
+    "Template_Stochastic_Naive4_correcting_:index.csv.xz",
+    "Template_Stochastic_Campaign_Default4_correcting_:index.csv.xz",
+    "Template_Stochastic_Campaign_Best4_correcting_:index.csv.xz",
+    "Template_Stochastic_Routine_Default4_correcting_:index.csv.xz",
+    "Template_Stochastic_Routine_Best4_correcting_:index.csv.xz"),
+  index = 1:200))
+
+stub <- "Sean Moore - stochastic_burden_est_JE_UND-Moore"
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "UND-Moore",
+  in_path = file.path(base_in_path, "UND-Moore"),
+  out_path = file.path(base_out_path, "JE_UND-Moore"),
+  scenarios =  scenarios,
+  files = c(
+    sprintf("%s_je-campaign-bestcase.csv.xz", stub),
+    sprintf("%s_je-campaign-default.csv.xz", stub),
+    sprintf("%s_je-no-vaccination.csv.xz", stub),
+    sprintf("%s_je-routine-bestcase.csv.xz", stub),
+    sprintf("%s_je-routine-default.csv.xz", stub))))
+
+
+
+###############
+# Measles
+scenarios <- c("measles-no-vaccination", "measles-mcv1-default",
+               "measles-mcv2-default", "measles-mcv1-bestcase",
+               "measles-mcv2-bestcase", "measles-campaign-default",
+               "measles-campaign-only-default", "measles-campaign-bestcase",
+               "measles-campaign-only-bestcase", "measles-stop")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "PSU-Ferrari",
+  in_path = file.path(base_in_path, "PSU-Ferrari"),
+  out_path = file.path(base_out_path, "Measles_PSU-Ferrari"),
+  scenarios = scenarios,
+  files = c(
+    "Heather Santos - novax_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "Heather Santos - default_mcv1_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "Heather Santos - default_mcv2_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "Heather Santos - bestcase_mcv1_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "Heather Santos - bestcase_mcv2_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "Heather Santos - default_campaign_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "Heather Santos - default_campaign_only_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "Heather Santos - bestcase_campaign_only_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz",
+    "Heather Santos - stop_stochastic:index_burden_Measles-PSU-Ferrari.csv.xz"),
+  index = 1:8))
+
+stub <- "stochastic_burden_estimate_measles-LSHTM-Jit"
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Jit",
+  in_path = file.path(base_in_path, "LSHTM-Jit_Measles"),
+  out_path = file.path(base_out_path, "Measles_LSHTM-Jit"),
+  scenarios = scenarios,
+  files = c(sprintf("%s-no-vaccination_Portnoy.csv.xz", stub),
+            sprintf("%s-mcv1-default_Portnoy.csv.xz", stub),
+            sprintf("%s-mcv2-default_Portnoy.csv.xz", stub),
+            sprintf("%s-mcv1-bestcase_Portnoy.csv.xz", stub),
+            sprintf("%s-mcv2-bestcase_Portnoy.csv.xz", stub),
+            sprintf("%s-campaign-default_Portnoy.csv.xz", stub),
+            sprintf("%s-campaign-only-default_Portnoy.csv.xz", stub),
+            sprintf("%s-campaign-bestcase_Portnoy.csv.xz", stub),
+            sprintf("%s-campaign-only-bestcase_Portnoy.csv.xz", stub),
+            sprintf("%s-stop_Portnoy.csv.xz", stub))))
+
+###############
+# MenA
+
+scenarios = c("mena-no-vaccination", "mena-campaign-default",
+              "mena-campaign-bestcase", "mena-routine-default",
+              "mena-routine-bestcase")
+stub <- "Andromachi Karachaliou - stochastic-burden.201910gavi-4.MenA_Cambridge-Trotter"
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "Cambridge-Trotter",
+  in_path = file.path(base_in_path, "Cambridge-Trotter"),
+  out_path = file.path(base_out_path, "MenA_Cambridge-Trotter"),
+  scenarios = scenarios,
+  files = c(sprintf("%s_no-vaccination_:index.csv.xz", stub),
+            sprintf("%s_campaign-default_:index.csv.xz", stub),
+            sprintf("%s_campaign-bestcase_:index.csv.xz", stub),
+            sprintf("%s_routine-default_:index.csv.xz", stub),
+            sprintf("%s_routine-bestcase_:index.csv.xz", stub)),
+  index = 1:52))
+
+stub <- "Michael Jackson - stochastic_burden_est_MenA_KPWA"
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "KPW-Jackson",
+  in_path = file.path(base_in_path, "KPW-Jackson"),
+  out_path = file.path(base_out_path, "MenA_KPW-Jackson"),
+  scenarios = scenarios,
+  files = c(sprintf("%s_both_bestcase_:index.csv.xz", stub),
+            sprintf("%s_both_default_:index.csv.xz", stub),
+            sprintf("%s_campaign_bestcase_:index.csv.xz", stub),
+            sprintf("%s_campaign_default_:index.csv.xz", stub),
+            sprintf("%s_none_default_:index.csv.xz", stub)),
+  index = 1:26))
+
+###############
+# PCV
+
+scenarios <- c("pcv-no-vaccination", "pcv-routine-default",
+               "pcv-routine-bestcase")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-PCV"),
+  out_path = file.path(base_out_path, "PCV_JHU-Tam"),
+  scenarios = scenarios,
+  files = c("novac:index.csv.xz", "default:index.csv.xz", "best:index.csv.xz"),
+  index = 1:14)
+
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Clark",
+  in_path = file.path(base_in_path, "LSHTM-Clark_PCV"),
+  out_path = file.path(base_out_path, "PCV_LSHTM-Clark"),
+  scenarios = scenarios,
+  files = c("VIMC_Sp_PSA_NoVax.csv.xz", "VIMC_Sp_PSA_Default.csv.xz",
+            "VIMC_Sp_PSA_Best.csv.xz"))
+
+###############
+# Rota
+
+scenarios <- c("rota-no-vaccination", "rota-routine-default",
+               "rota-routine-bestcase")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-Rota"),
+  out_path = file.path(base_out_path, "Rota_JHU-Tam"),
+  scenarios = scenarios,
+  files = c("novac:index.csv.xz", "default:index.csv.xz", "best:index.csv.xz"),
+  index = 1:14)
+
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Clark",
+  in_path = file.path(base_in_path, "LSHTM-Clark_Rota"),
+  out_path = file.path(base_out_path, "Rota_LSHTM-Clark"),
+  scenarios = scenarios,
+  files = c("Hira Tanvir - VIMC_Rota_PSA_NoVax.csv.xz",
+            "Hira Tanvir - VIMC_Rota_PSA_Default.csv.xz",
+            "Hira Tanvir - VIMC_Rota_PSA_Best.csv.xz"))
+
+stoner::stone_stochastic_standardise(
+  group = "Emory-Lopman",
+  in_path = file.path(base_in_path, "Emory-Lopman"),
+  out_path = file.path(base_out_path, "Rota_Emory-Lopman"),
+  scenarios = scenarios,
+  files = paste0("Molly Steele - stochastic-burden.201910gavi-4.",
+                "Rota_Emory-Lopman_:scenario.csv.xz"))
+
+###############
+# Rubella
+
+scenarios <- c(
+  "rubella-routine-no-vaccination", "rubella-campaign-default",
+  "rubella-campaign-bestcase", "rubella-rcv1-default",
+  "rubella-rcv1-bestcase", "rubella-rcv1-rcv2-default",
+  "rubella-rcv1-rcv2-bestcase", "rubella-rcv2-default",
+  "rubella-rcv2-bestcase", "rubella-stop")
+
+stub <- "Amy Winter - stochastic_burden_est-rubella"
+
+id10 <- hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "JHU-Lessler",
+  in_path = file.path(base_in_path, "JHU-Lessler"),
+  out_path = file.path(base_out_path, "Rubella_JHU-Lessler"),
+  scenarios = scenarios,
+  files = c(sprintf("%s-no-vaccination_:index.csv.xz", stub),
+            sprintf("%s-campaign-default_:index.csv.xz", stub),
+            sprintf("%s-campaign-bestcase_:index.csv.xz", stub),
+            sprintf("%s-rcv1-default_:index.csv.xz", stub),
+            sprintf("%s-rcv1-bestcase_:index.csv.xz", stub),
+            sprintf("%s-rcv1-rcv2-default_:index.csv.xz", stub),
+            sprintf("%s-rcv1-rcv2-bestcase_:index.csv.xz", stub),
+            sprintf("%s-rcv2-default_:index.csv.xz", stub),
+            sprintf("%s-rcv2-bestcase_:index.csv.xz", stub),
+            sprintf("%s-stop_:index.csv.xz", stub)),
+  index = 1:12))
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "PHE-Vynnycky",
+  in_path = file.path(base_in_path, "PHE-Vynnycky"),
+  out_path = file.path(base_out_path, "Rubella_PHE-Vynnycky"),
+  scenarios = scenarios,
+  files = "stochastic_burden_est_:scenario_country:index.csv.xz",
+  index = 1:112))
+
+
+###############
+# Typhoid
+
+scenarios = c("typhoid-no-vaccination",
+                "typhoid-campaign-default", "typhoid-routine-default")
+
+stoner::stone_stochastic_standardise(
+  group = "IVI-Kim",
+  in_path = file.path(base_in_path, "IVI-Kim-Typhoid"),
+  out_path = file.path(base_out_path, "Typhoid_IVI-Kim"),
+  scenarios = scenarios,
+  files = c("Jong-Hoon Kim - stoch_Typhoid_novacc.csv.xz",
+           "Jong-Hoon Kim - stoch_Typhoid_campaign.csv.xz",
+           "Jong-Hoon Kim - stoch_Typhoid_routine.csv.xz"))
+
+stoner::stone_stochastic_standardise(
+  group = "Yale-Pitzer",
+  in_path = file.path(base_in_path, "Yale-Pitzer"),
+  out_path = file.path(base_out_path, "Typhoid_Yale-Pitzer"),
+  scenarios = scenarios,
+  files = c("Virginia Pitzer - 2021-02-18 17.00.26 - stochastic_output_TF-Yale-Pitzer_novacc.csv.xz",
+            "Virginia Pitzer - 2021-02-18 16.58.03 - stochastic_output_TF-Yale-Pitzer_campaign.csv.xz",
+            "Virginia Pitzer - 2021-02-18 16.59.14 - stochastic_output_TF-Yale-Pitzer_camproutine.csv.xz"))
+
+###############
+# YF
+
+scenarios <- c("yf-no-vaccination", "yf-preventive-bestcase",
+               "yf-preventive-default", "yf-routine-bestcase",
+               "yf-routine-default", "yf-stop")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "IC-Garske",
+  in_path = file.path(base_in_path, "IC-Garske"),
+  out_path = file.path(base_out_path, "YF_IC-Garske"),
+  scenarios = scenarios,
+  files = "stochastic-burden-estimates.201910gavi-4_YF_IC-Garske_:scenario_:index.csv.xz",
+  index = 1:200))
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "UND-Perkins",
+  in_path = file.path(base_in_path, "UND-Perkins"),
+  out_path = file.path(base_out_path, "YF_UND-Perkins"),
+  scenarios = scenarios,
+  files = "stochastic_burden_est_YF_UND-Perkins_:scenario_:index.csv.xz",
+  index = 1:200))
diff --git a/scripts/process_stoch_202110gavi.R b/scripts/process_stoch_202110gavi.R
new file mode 100644
index 00000000..210d016c
--- /dev/null
+++ b/scripts/process_stoch_202110gavi.R
@@ -0,0 +1,386 @@
+base_in_path <- "//wpia-hn2.hpc.dide.ic.ac.uk/vimc_stochastics_dropbox/latest/202110gavi"
+base_out_path <- "//wpia-hn2.hpc.dide.ic.ac.uk/vimc_stochastics/202110gavi"
+
+vault <- vaultr::vault_client(login = "github")
+password <- vault$read("/secret/vimc/database/production/users/readonly")$password
+con <- DBI::dbConnect(RPostgres::Postgres(),
+                      dbname = "montagu",
+                      host = "montagu.vaccineimpact.org",
+                      port = 5432, password = password,
+                      user = "readonly")
+
+fetch_scenarios <- function(disease) {
+  sort(DBI::dbGetQuery(con, "
+    SELECT scenario_description FROM
+           scenario JOIN scenario_description
+           ON scenario.scenario_description = scenario_description.id
+           WHERE disease = $1 AND touchstone='202110gavi-3'", disease)$scenario_description)
+}
+
+# Let's unleash the cow
+
+setwd("Q:/testcow")
+writeLines("vimc/stoner@VIMC-9230", "pkgdepends.txt")
+hipercow::hipercow_init(driver = "dide-windows")
+hipercow::hipercow_provision()
+# Network/memory might be too much for more than a job per node.
+hres <- hipercow::hipercow_resources(cores = 1L, exclusive = TRUE)
+
+###############
+# Cholera
+
+scenarios = c("cholera-no-vaccination", "cholera-campaign-default")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "IVI-Kim",
+  in_path = file.path(base_in_path, "IVI-Kim-Cholera"),
+  out_path = file.path(base_out_path, "Cholera_IVI-Kim"),
+  scenarios = scenarios,
+  files = c("Jong-Hoon Kim - stoch_Cholera_campaign_20211221T00.csv.xz",
+            "Jong-Hoon Kim - stoch_Cholera_campaign_20211222T212131.csv.xz")
+  ))
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Lee",
+  in_path = file.path(base_in_path, "JHU-Lee-Cholera"),
+  out_path = file.path(base_out_path, "Cholera_JHU-Lee"),
+  scenarios = scenarios,
+  files = c("Kaiyue Zou - no-vaccination.csv.xz",
+            "Kaiyue Zou - campaign-default.csv.xz"))
+
+###############
+# HepB
+
+scenarios <- c("hepb-no-vaccination",
+               "hepb-bd-routine-default",
+               "hepb-bd-routine-ia2030_target",
+               "hepb-bd-routine-ia2030_target-hepb-routine-ia2030_target",
+               "hepb-bd-default-hepb-routine-default",
+               "hepb-hepb-routine-default",
+               "hepb-hepb-routine-ia2030_target")
+
+hipercow::task_create_expr(resources = hres, expr =
+  stoner::stone_stochastic_standardise(
+    group = "Li",
+    in_path = file.path(base_in_path, "Li"),
+    out_path = file.path(base_out_path, "HepB_Li"),
+    scenarios = scenarios,
+    files = ":scenario:index.csv.xz",
+    index = 1:200))
+
+hipercow::task_create_expr(resources = hres, expr =
+  stoner::stone_stochastic_standardise(
+    group = "IC-Hallett",
+    in_path = file.path(base_in_path, "IC-Hallett"),
+    out_path = file.path(base_out_path, "HepB_IC-Hallett"),
+    scenarios = scenarios,
+    files = "stochastic_burden_est_HepB-IC-Hallett_:scenario_:index.csv.xz",
+    index = 1:200))
+
+###############
+# HIB
+
+scenarios <- c("hib-no-vaccination", "hib-routine-default",
+               "hib-routine-ia2030_target")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-Carter-Hib"),
+  out_path = file.path(base_out_path, "Hib_JHU-Tam"),
+  scenarios = scenarios,
+  files = c("hib-no-vaccination-LiST.csv.xz",
+            "hib-routine-default-LiST.csv.xz",
+            "hib-routine-ia2030_target-LiST.csv.xz"))
+
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Clark",
+  in_path = file.path(base_in_path, "LSHTM-Clark_Hib"),
+  out_path = file.path(base_out_path, "Hib_LSHTM-Clark"),
+  scenarios = scenarios,
+  files = c("Kaja Abbas - PSA_202110gavi-3_hib-no-vaccination.csv.xz",
+            "Kaja Abbas - PSA_202110gavi-3_hib-routine-default.csv.xz",
+            "Kaja Abbas - PSA_202110gavi-3_hib-routine-ia2030_target.csv.xz"))
+
+###############
+# HPV
+
+scenarios <- c("hpv-no-vaccination", "hpv-campaign-default",
+               "hpv-routine-default","hpv-campaign-ia2030_target",
+               "hpv-routine-ia2030_target")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "Harvard-Sweet",
+  in_path = file.path(base_in_path, "Harvard-Sweet"),
+  out_path = file.path(base_out_path, "HPV_Harvard-Sweet"),
+  scenarios = scenarios,
+  files = c("Allison Portnoy - stochastic-burden-est.novacc_run_200.csv.xz",
+            "Allison Portnoy - stochastic-burden-est.coverage_202110gavi-3_hpv-campaign-default_run_:index.csv.xz",
+            "Allison Portnoy - stochastic-burden-est.coverage_202110gavi-3_hpv-routine-default_run_:index.csv.xz",
+            "Allison Portnoy - stochastic-burden-est.coverage_202110gavi-3_hpv-campaign-ia2030_target_run_:index.csv.xz",
+            "Allison Portnoy - stochastic-burden-est.coverage_202110gavi-3_hpv-routine-ia2030_target_run_:index.csv.xz"),
+  index = 1:200))
+
+id <- hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Jit",
+  in_path = file.path(base_in_path, "LSHTM-Jit_HPV"),
+  out_path = file.path(base_out_path, "HPV_LSHTM-Jit"),
+  scenarios =  scenarios,
+  files =
+    c("stochastic-burden-novaccination_all_202110gavi-3_hpv-no-vaccination.csv.xz",
+      "stochastic-burden-vaccination_all_202110gavi-3_hpv-campaign-default.csv.xz",
+      "stochastic-burden-vaccination_all_202110gavi-3_hpv-routine-default.csv.xz",
+      "stochastic-burden-vaccination_all_202110gavi-3_hpv-campaign-ia2030_target.csv.xz",
+      "stochastic-burden-vaccination_all_202110gavi-3_hpv-routine-ia2030_target.csv.xz")))
+
+###############
+# JE
+
+scenarios <- c("je-routine-no-vaccination", "je-campaign-default",
+               "je-campaign-ia2030_target", "je-routine-default",
+               "je-routine-ia2030_target")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "NUS-Clapham",
+  in_path = file.path(base_in_path, "NUS-Clapham-JE"),
+  out_path = file.path(base_out_path, "JE_NUS-Clapham"),
+  scenarios =  scenarios,
+  files = c(
+    "Naive_Stochastic__:index.csv.xz",
+    "Campaign_Stochastic_:index.csv.xz",
+    "Campaign_Target_Stochastic_:index.csv.xz",
+    "Routine_Stochastic_:index.csv.xz",
+    "Routine_Target_Stochastic_:index.csv.xz"),
+  index = 1:200))
+
+stub <- "Sean Moore - stochastic_burden_est_JE_UND-Moore"
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "UND-Moore",
+  in_path = file.path(base_in_path, "UND-Moore-JE"),
+  out_path = file.path(base_out_path, "JE_UND-Moore"),
+  scenarios =  scenarios,
+  files = sprintf("%s_:scenario.csv.xz", stub)))
+
+
+###############
+# Measles
+scenarios <- c("measles-no-vaccination",
+               "measles-campaign-default",
+               "measles-campaign-ia2030_target",
+               "measles-campaign-only-default",
+               "measles-campaign-only-ia2030_target",
+               "measles-mcv1-default",
+               "measles-mcv1-ia2030_target",
+               "measles-mcv2-default",
+               "measles-mcv2-ia2030_target")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "PSU-Ferrari",
+  in_path = file.path(base_in_path, "PSU-Ferrari-Measles"),
+  out_path = file.path(base_out_path, "Measles_PSU-Ferrari"),
+  scenarios = scenarios,
+  files = "coverage_202110gavi-3_:scenario.csv.xz"))
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Jit",
+  in_path = file.path(base_in_path, "LSHTM-Jit_Measles"),
+  out_path = file.path(base_out_path, "Measles_LSHTM-Jit"),
+  scenarios = scenarios,
+  files = c("Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-no-vaccination.csv.xz",
+            "Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-campaign-default.csv.xz",
+            "Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-campaign-ia2030_target.csv.xz",
+            "Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-campaign-only-default.csv.xz",
+            "Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-campaign-only-ia2030_target.csv.xz",
+            "Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-mcv1-default.csv.xz",
+            "Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-mcv1-ia2030_target.csv.xz",
+            "Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-mcv2-default.csv.xz",
+            "Han Fu - stochastic_burden_estimate_measles-LSHTM-Jit-mcv2-ia2030_target.csv.xz"
+            )))
+
+###############
+# MenA
+
+scenarios = c("mena-no-vaccination", "mena-campaign-default",
+              "mena-campaign-ia2030_target", "mena-routine-default",
+              "mena-routine-ia2030_target", "mena-booster-default")
+stub <- "Andromachi Karachaliou - stochastic-burden.202110gavi-2.MenA_Cambridge-Trotter"
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "Cambridge-Trotter",
+  in_path = file.path(base_in_path, "Cambridge-Trotter"),
+  out_path = file.path(base_out_path, "MenA_Cambridge-Trotter"),
+  scenarios = scenarios,
+  files = c(sprintf("%s_no_vaccination_:index.csv.xz", stub),
+            sprintf("%s_campaign_default_:index.csv.xz", stub),
+            sprintf("%s_campaign-ia2030_target_:index.csv.xz", stub),
+            sprintf("%s_routine_default_:index.csv.xz", stub),
+            sprintf("%s_routine-ia2030_target_:index.csv.xz", stub),
+            sprintf("%s_booster_:index.csv.xz", stub)),
+  index = 1:26))
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "KPW-Jackson",
+  in_path = file.path(base_in_path, "KPW-Jackson-MenA"),
+  out_path = file.path(base_out_path, "MenA_KPW-Jackson"),
+  scenarios = scenarios,
+  files = c("stochastic_burden_est_MenA_KPWA_none_default_:index.csv.xz",
+            "stochastic_burden_est_MenA_KPWA_campaign_default_:index.csv.xz",
+            "stochastic_burden_est_MenA_KPWA_campaign_ia2030_target_:index.csv.xz",
+            "stochastic_burden_est_MenA_KPWA_routine_default_:index.csv.xz",
+            "stochastic_burden_est_MenA_KPWA_routine_ia2030_target_:index.csv.xz",
+            "stochastic_burden_est_MenA_KPWA_booster_default_:index.csv.xz"),
+  index = 1:26))
+
+###############
+# PCV
+
+scenarios <- c("pcv-no-vaccination", "pcv-routine-default",
+               "pcv-routine-ia2030_target")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-Carter-PCV"),
+  out_path = file.path(base_out_path, "PCV_JHU-Tam"),
+  scenarios = scenarios,
+  files = ":scenario-LiST.csv.xz")
+
+stoner::stone_stochastic_standardise(
+  group = "NUS-Chen",
+  in_path = file.path(base_in_path, "LSHTM-NUS-Chen_PCV"),
+  out_path = file.path(base_out_path, "PCV_NUS-Chen"),
+  scenarios = scenarios,
+  files = "stochastic_burden_est_:scenario.csv.xz")
+
+###############
+# Rota
+
+scenarios <- c("rota-no-vaccination", "rota-routine-default",
+               "rota-routine-ia2030_target")
+
+stoner::stone_stochastic_standardise(
+  group = "JHU-Tam",
+  in_path = file.path(base_in_path, "JHU-Tam-Carter-Rota"),
+  out_path = file.path(base_out_path, "Rota_JHU-Tam"),
+  scenarios = scenarios,
+  files = ":scenario-LiST.csv.xz")
+
+stoner::stone_stochastic_standardise(
+  group = "LSHTM-Clark",
+  in_path = file.path(base_in_path, "LSHTM-Clark_Rota"),
+  out_path = file.path(base_out_path, "Rota_LSHTM-Clark"),
+  scenarios = scenarios,
+  files = "Kaja Abbas - PSA_202110gavi-3_:scenario.csv.xz")
+
+stoner::stone_stochastic_standardise(
+  group = "Emory-Lopman",
+  in_path = file.path(base_in_path, "Emory-Lopman"),
+  out_path = file.path(base_out_path, "Rota_Emory-Lopman"),
+  scenarios = scenarios,
+  files = c("Aniruddha Deshpande - stochastic_burden_est_lopman_no_vaccination_2022_01_31.csv.xz",
+            "Aniruddha Deshpande - stochastic_burden_est_lopman_routine_2022_01_31.csv.xz",
+            "Aniruddha Deshpande - stochastic_burden_est_lopman_ia2030_target_2022_01_31.csv.xz"))
+
+###############
+# Rubella
+
+scenarios <- c("rubella-routine-no-vaccination",
+  "rubella-campaign-default", "rubella-campaign-ia2030_target",
+  "rubella-rcv1-default", "rubella-rcv1-ia2030_target",
+  "rubella-rcv1-rcv2-default", "rubella-rcv1-rcv2-ia2030_target",
+  "rubella-rcv2-default", "rubella-rcv2-ia2030_target")
+
+id10 <- hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "JHU-Lessler",
+  in_path = file.path(base_in_path, "JHU-UGA-Winter-Rubella"),
+  out_path = file.path(base_out_path, "Rubella_JHU-Lessler"),
+  scenarios = scenarios,
+  files = "Amy Winter - stochastic_burden_est-:scenario_:index.csv.xz",
+  index = 1:11))
+
+id11<-hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "PHE-Vynnycky",
+  in_path = file.path(base_in_path, "PHE-Vynnycky_Rubella"),
+  out_path = file.path(base_out_path, "Rubella_PHE-Vynnycky"),
+  scenarios = scenarios,
+  files = c("VIMC_NV_RCV1RCV2Camp_country:index.csv.xz",
+            "VIMC_DF_Camp_country:index.csv.xz",
+            "VIMC_IA_Camp_country:index.csv.xz",
+            "VIMC_DF_RCV1Camp_country:index.csv.xz",
+            "VIMC_IA_RCV1Camp_country:index.csv.xz",
+            "VIMC_DF_RCV1RCV2_country:index.csv.xz",
+            "VIMC_IA_RCV1RCV2_country:index.csv.xz",
+            "VIMC_DF_RCV1RCV2Camp_country:index.csv.xz",
+            "VIMC_IA_RCV1RCV2Camp_country:index.csv.xz"),
+  index = 1:112))
+
+
+
+###############
+# Typhoid
+
+scenarios = c("typhoid-no-vaccination",
+              "typhoid-campaign-default",
+              "typhoid-campaign-ia2030_target",
+              "typhoid-routine-default",
+              "typhoid-routine-ia2030_target")
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "IVI-Kim",
+  in_path = file.path(base_in_path, "IVI-Kim-Typhoid"),
+  out_path = file.path(base_out_path, "Typhoid_IVI-Kim"),
+  scenarios = scenarios,
+  files = c("Jong-Hoon Kim - stoch_Typhoid_novacc_20220509T233408.csv.xz",
+            "Jong-Hoon Kim - stoch_Typhoid_campaign-default_20220510T103619.csv.xz",
+            "Jong-Hoon Kim - stoch_Typhoid_campaign-ia2030_20220510T105521.csv.xz",
+            "Jong-Hoon Kim - stoch_Typhoid_routine-default_20220510T103357.csv.xz",
+            "Jong-Hoon Kim - stoch_Typhoid_routine-ia2030_20220514.csv.xz")))
+
+hipercow::task_create_expr(resources = hres, expr =
+stoner::stone_stochastic_standardise(
+  group = "Yale-Pitzer",
+  in_path = file.path(base_in_path, "Yale-Pitzer-Typhoid"),
+  out_path = file.path(base_out_path, "Typhoid_Yale-Pitzer"),
+  scenarios = scenarios,
+  files = c("Holly Burrows - stochastic_burden_est_TF-Yale-Burrows-novacc_202110.csv.xz",
+            "Holly Burrows - stochastic_burden_est_TF-Yale-Burrows_campaign-default_202110.csv.xz",
+            "Holly Burrows - stochastic_burden_est_TF-Yale-Burrows_campaign-IA2030_202110.csv.xz",
+            "Holly Burrows - stochastic_burden_est_TF-Yale-Burrows_routine-default_202110.csv.xz",
+            "Holly Burrows - stochastic_burden_est_TF-Yale-Burrows_routine-IA2030_202110.csv.xz")))
+
+
+
+###############
+# YF
+
+scenarios <- c("yf-no-vaccination", "yf-preventive-ia2030_target",
+               "yf-preventive-default", "yf-routine-ia2030_target",
+               "yf-routine-default")
+
+stub <- "Keith Fraser - stochastic-burden-estimates.202110gavi-3_YF_IC-Garske"
+stoner::stone_stochastic_standardise(
+  group = "IC-Garske",
+  in_path = file.path(base_in_path, "IC-Garske"),
+  out_path = file.path(base_out_path, "YF_IC-Garske"),
+  scenarios = scenarios,
+  files = sprintf("%s_:scenario_:index.csv.xz", stub),
+  index = 1:200)
+
+stoner::stone_stochastic_standardise(
+  group = "UND-Perkins",
+  in_path = file.path(base_in_path, "UND-Perkins-YF"),
+  out_path = file.path(base_out_path, "YF_UND-Perkins"),
+  scenarios = scenarios,
+  files = "stochastic_burden_est_YF_UND-Perkins_:scenario_:index.csv.xz",
+  index = 1:200)
diff --git a/scripts/process_stoch_202310gavi.R b/scripts/process_stoch_202310gavi.R
index 0c4915b4..b2eaaefa 100644
--- a/scripts/process_stoch_202310gavi.R
+++ b/scripts/process_stoch_202310gavi.R
@@ -20,7 +20,7 @@ stoner::stone_stochastic_standardise(
 stoner::stone_stochastic_standardise(
   group = "JHU-Lee",
   in_path = file.path(base_in_path, "Cholera-JHU-Lee"),
-  out_path = file.path(base_out_path, "Cholera-JHU-Lee"),
+  out_path = file.path(base_out_path, "Cholera_JHU-Lee"),
   scenarios = scenarios,
   files = c("stochastic-burden-template.202310gavi-4.Cholera_standard_template.529.no-vaccination Christina Alam.csv.xz",
             "stochastic-burden-template.202310gavi-4.Cholera_standard_template.529.ocv1-default_one Christina Alam.csv.xz",
diff --git a/tests/testthat/test_stochastic_files.R b/tests/testthat/test_stochastic_files.R
index d1a668ab..ac19da8e 100644
--- a/tests/testthat/test_stochastic_files.R
+++ b/tests/testthat/test_stochastic_files.R
@@ -147,3 +147,186 @@ test_that("Create central works", {
 
   expect_true(p1$cases[p1$country == "NPL"] == 2 * p1$cases[p1$country == "LAP"])
 })
+
+test_that("Rubella fix works", {
+  df <- fake_data()
+  df$rubella_deaths_congenital <- df$deaths
+  df$rubella_cases_congenital <- df$cases
+  df$deaths <- NULL
+  df$cases <- NULL
+  df$rubella_infections <- sample(nrow(df))
+
+  tmpin <- tempdir()
+  tmpout <- tempdir()
+  tmpfile <- tempfile(tmpdir = tmpin)
+  write.csv(df, paste0(tmpfile, "_opt"), row.names = FALSE)
+
+  stone_stochastic_standardise(
+    group = "north_pole_rub", in_path = tmpin, out_path = tmpout,
+    scenarios = "opt",
+    files = paste0(basename(tmpfile), "_opt")
+  )
+
+  files <- list.files(path = tmpout)
+  expect_true("north_pole_rub_opt_LAP.pq" %in% files)
+  pq <- arrow::read_parquet(file.path(tmpout, "north_pole_rub_opt_LAP.pq"))
+  expect_true("cases" %in% names(pq))
+  expect_true("deaths" %in% names(pq))
+  expect_false("rubella_deaths_congenital" %in% names(pq))
+  expect_false("rubella_cases_congenital" %in% names(pq))
+  expect_false("rubella_infections" %in% names(pq))
+  expect_true(all(pq$cases == df$rubella_cases_congenital))
+  expect_true(all(pq$deaths == df$rubella_deaths_congenital))
+})
+
+test_that("Hib/PCV fix works", {
+  df <- fake_data()
+  df$cases_pneumo <- df$cases + 1
+  df$cases_men <- df$cases + 2
+  df$cases <- NULL
+  df$deaths_pneumo <- df$deaths + 3
+  df$deaths_men <- df$deaths + 4
+  df$deaths <- NULL
+
+  tmpin <- tempdir()
+  tmpout <- tempdir()
+  tmpfile <- tempfile(tmpdir = tmpin)
+  write.csv(df, paste0(tmpfile, "_opt"), row.names = FALSE)
+
+  stone_stochastic_standardise(
+    group = "north_pole_hib", in_path = tmpin, out_path = tmpout,
+    scenarios = "opt",
+    files = paste0(basename(tmpfile), "_opt")
+  )
+
+  files <- list.files(path = tmpout)
+  expect_true("north_pole_hib_opt_LAP.pq" %in% files)
+  pq <- arrow::read_parquet(file.path(tmpout, "north_pole_hib_opt_LAP.pq"))
+  expect_true("cases" %in% names(pq))
+  expect_true("deaths" %in% names(pq))
+  expect_false("cases_men" %in% names(pq))
+  expect_false("cases_pneumo" %in% names(pq))
+  expect_false("deaths_men" %in% names(pq))
+  expect_false("deaths_pneumo" %in% names(pq))
+
+  expect_true(all(pq$cases == df$cases_pneumo + df$cases_men))
+  expect_true(all(pq$deaths == df$deaths_pneumo + df$deaths_men))
+
+})
+
+test_that("HepB fix works", {
+  df <- fake_data()
+  df$hepb_cases_acute_severe <- df$cases + 1
+  df$hepb_cases_dec_cirrh <- df$cases + 2
+  df$hepb_cases_hcc <- df$cases + 3
+  df$hepb_cases_acute_symp <- df$cases + 4
+  df$hepb_cases_fulminant <- df$cases + 5
+  df$hepb_cases_chronic <- df$cases + 6
+  df$hepb_chronic_symptomatic_in_acute_phase <- df$cases + 7
+  df$hepb_cases_comp_cirrh <- df$cases + 8
+  df$hepb_cases_hcc_no_cirrh <- df$cases + 9
+  df$cases <- NULL
+  df$hepb_deaths_acute <- df$deaths + 1
+  df$hepb_deaths_dec_cirrh <- df$deaths + 2
+  df$hepb_deaths_hcc <- df$deaths + 3
+  df$hepb_deaths_total_cirrh <- df$deaths + 4
+  df$hepb_deaths_comp_cirrh <- df$deaths + 5
+  df$deaths <- NULL
+
+  tmpin <- tempdir()
+  tmpout <- tempdir()
+  tmpfile <- tempfile(tmpdir = tmpin)
+  write.csv(df, paste0(tmpfile, "_opt"), row.names = FALSE)
+
+  stone_stochastic_standardise(
+    group = "north_pole_hepb", in_path = tmpin, out_path = tmpout,
+    scenarios = "opt",
+    files = paste0(basename(tmpfile), "_opt")
+  )
+
+  files <- list.files(path = tmpout)
+  expect_true("north_pole_hepb_opt_LAP.pq" %in% files)
+  pq <- arrow::read_parquet(file.path(tmpout, "north_pole_hepb_opt_LAP.pq"))
+
+  expect_true("cases" %in% names(pq))
+  expect_true("deaths" %in% names(pq))
+  expect_false("hepb_cases_acute_severe" %in% names(pq))
+  expect_false("hepb_cases_dec_cirrh" %in% names(pq))
+  expect_false("hepb_cases_hcc" %in% names(pq))
+  expect_false("hepb_cases_acute_symp" %in% names(pq))
+  expect_false("hepb_cases_fulminant" %in% names(pq))
+  expect_false("hepb_cases_chronic" %in% names(pq))
+  expect_false("hepb_chronic_symptomatic_in_acute_phase" %in% names(pq))
+  expect_false("hepb_cases_comp_cirrh" %in% names(pq))
+  expect_false("hepb_cases_comp_hcc_no_cirrh" %in% names(pq))
+
+  expect_false("hepb_deaths_acute" %in% names(pq))
+  expect_false("hepb_deaths_dec_cirrh" %in% names(pq))
+  expect_false("hepb_deaths_hcc" %in% names(pq))
+  expect_false("hepb_deaths_total_cirrh" %in% names(pq))
+  expect_false("hepb_deaths_comp_cirrh" %in% names(pq))
+
+  expect_true(all(pq$cases == df$hepb_cases_acute_severe + df$hepb_cases_dec_cirrh +
+                          df$hepb_Cases_hcc + df$hepb_cases_acute_symp +
+                          df$hepb_cases_fulminant + df$hepb_cases_chronic +
+                          df$hepb_chronic_symptomatic_in_acute_phase +
+                          df$hepb_cases_comp_cirrh + df$hepb_cases_comp_hcc_no_cirrh))
+
+  expect_true(all(pq$deaths == df$hepb_deaths_acute + df$hepb_deaths_dec_cirrh +
+                           df$hep_deaths_hcc + df$hepb_deaths_total_cirrh +
+                           df$hepb_deaths_comp_cirrh))
+})
+
+test_that("Missing yll/dalys works", {
+  df <- fake_data()
+  df$dalys <- NULL
+  df$yll <- NULL
+  tmpin <- tempdir()
+  tmpout <- tempdir()
+  tmpfile <- tempfile(tmpdir = tmpin)
+  write.csv(df, paste0(tmpfile, "_opt"), row.names = FALSE)
+
+  stone_stochastic_standardise(
+    group = "north_pole_flu", in_path = tmpin, out_path = tmpout,
+    scenarios = "opt",
+    files = paste0(basename(tmpfile), "_opt")
+  )
+  files <- list.files(path = tmpout)
+  expect_true("north_pole_flu_opt_LAP.pq" %in% files)
+  pq <- arrow::read_parquet(file.path(tmpout, "north_pole_flu_opt_LAP.pq"))
+  expect_false("dalys" %in% names(pq))
+  expect_false("yll" %in% names(pq))
+})
+
+test_that("Different file count per scenario is handled", {
+  fake <- fake_data()
+  tmpin <- tempdir()
+  tmpout <- tempdir()
+  tmpfile <- tempfile(tmpdir = tmpin)
+
+  fake$country <- "POL"
+  write.csv(fake, paste0(tmpfile, "_optimistic_1"), row.names = FALSE)
+  write.csv(fake, paste0(tmpfile, "_fatalistic_1"), row.names = FALSE)
+  fake$country <- "NOR"
+  write.csv(fake, paste0(tmpfile, "_optimistic_2"), row.names = FALSE)
+  write.csv(fake, paste0(tmpfile, "_fatalistic_2"), row.names = FALSE)
+  fake$country <- "LAP"
+  write.csv(fake, paste0(tmpfile, "_fatalistic_3"), row.names = FALSE)
+
+  stone_stochastic_standardise(
+    group = "north_pole_lurgy",
+    in_path = tmpin,
+    out_path = tmpout,
+    scenarios = c("optimistic", "fatalistic"),
+    files = paste0(basename(tmpfile), "_:scenario_:index"),
+    index = 1:3,
+    allow_missing_indexes = TRUE
+  )
+  files <- list.files(path = tmpout)
+  expect_true("north_pole_lurgy_optimistic_POL.pq" %in% files)
+  expect_true("north_pole_lurgy_optimistic_NOR.pq" %in% files)
+  expect_false("north_pole_lurgy_optimistic_LAP.pq" %in% files)
+  expect_true("north_pole_lurgy_fatalistic_POL.pq" %in% files)
+  expect_true("north_pole_lurgy_fatalistic_NOR.pq" %in% files)
+  expect_true("north_pole_lurgy_fatalistic_LAP.pq" %in% files)
+})
diff --git a/tests/testthat/test_stochastic_graphs.R b/tests/testthat/test_stochastic_graphs.R
index 243d3d3a..58a27358 100644
--- a/tests/testthat/test_stochastic_graphs.R
+++ b/tests/testthat/test_stochastic_graphs.R
@@ -66,4 +66,92 @@ test_that("stochastic_graph data transforms", {
     base, touchstone, disease, group, country,
     scenario, "deaths", log = TRUE))
 
+  expect_no_error(stone_stochastic_graph(
+    base, touchstone, disease, group, country,
+    scenario, "deaths", scenario2 = scenario))
+
+  # Packit gets called if needed
+
+  fake_result <- mockery::mock("fake_result")
+  mockery::stub(stone_stochastic_graph, "prepare_central_data", fake_result)
+
+  expect_no_error(stone_stochastic_graph(
+    base, touchstone, disease, group, country,
+    scenario, "deaths", packit_id = "123",
+    packit_file = "file.csv"))
+
+  mockery::expect_called(fake_result, 1)
+  mockery::expect_args(fake_result, 1, "123", "file.csv", country,
+                       scenario, "deaths", NULL, FALSE)
+
+})
+
+test_that("stochastic_explorer data_dir handling", {
+  expect_error(stochastic_explorer(file.path(tempdir(), "potato", "salad")),
+                      "Cannot access the path/mount")
+})
+
+test_that("Can launch shiny app", {
+  fake_path <- tempdir()
+  runApp_called <- FALSE
+  runApp_arg <- NULL
+
+  local_mocked_bindings(
+    runApp = function(app_dir) {
+      runApp_called <<- TRUE
+      runApp_arg <<- app_dir
+      invisible(NULL)
+    },
+    .env = environment(stochastic_explorer)
+  )
+
+  withr::with_envvar(c(), {
+    if (exists("data_dir", envir = .GlobalEnv)) {
+      rm(data_dir, envir = .GlobalEnv)
+    }
+    stochastic_explorer(data_dir = fake_path)
+    expect_true(exists("data_dir", envir = .GlobalEnv))
+    expect_equal(get("data_dir", envir = .GlobalEnv), fake_path)
+    expect_true(runApp_called)
+  })
+})
+
+test_that("Age formats are reasonable", {
+  expect_equal(age_string(NULL), "all ages")
+  expect_equal(age_string(c(5,4,3,2,1,5,4,3,2,1)), "age 1..5")
+  expect_equal(age_string(c(2,4,6,8)), "selected ages")
+})
+
+test_that("Parsing central from packit works", {
+  # Packit gets called if needed
+
+  fake <- data.frame(
+    scenario_type = "RSV-rout", scenario = "RSV-rout",
+    year = c(rep(2000, 4), rep(2001, 4), rep(2000, 4), rep(2001, 4)),
+    age = c(rep(0, 8), rep(1, 8)),
+    country = "RFP",
+    burden_outcome = rep(c("cases", "dalys", "deaths", "yll"), 2),
+    value = 1:16)
+
+  rds <- tempfile(fileext = ".rds")
+  saveRDS(fake, rds)
+
+  fetch_fake <- function(id, file) rds
+  mockery::stub(prepare_central_data, "fetch_packit", fetch_fake)
+
+  res <- prepare_central_data("123", "file.csv",
+    "RFP", "RSV-rout", "deaths", 0:5, TRUE)
+
+  # Data in for death is: (year, age, deaths)
+  # 2000, 0, 3
+  # 2000, 1, 11
+  # 2001, 0, 7
+  # 2001, 1, 15
+  # For cohort - this should become...
+  # 1999, 11   (2000 year 1, were born in 1999)
+  # 2000, 18   (2000 year 0, and 2001 year 1 born in 2000)
+  # 2001, 7    (2001 year 0)
+
+  expect_true(all.equal(res$year, c(1999, 2000, 2001)))
+  expect_true(all.equal(res$deaths, c(11, 18, 7)))
 })
diff --git a/vignettes/Stochastics.Rmd b/vignettes/Stochastics.Rmd
index 10bafeed..9be11c9d 100644
--- a/vignettes/Stochastics.Rmd
+++ b/vignettes/Stochastics.Rmd
@@ -4,7 +4,7 @@ author: "Wes Hinsley"
 date: "`r Sys.Date()`"
 output: rmarkdown::html_vignette
 vignette: >
-  %\VignetteIndexEntry{stoner}
+  %\VignetteIndexEntry{Stochastics}
   %\VignetteEngine{knitr::rmarkdown}
   %\VignetteEncoding{UTF-8}
 ---
diff --git a/vignettes/Usage.Rmd b/vignettes/Usage.Rmd
index 89d8b17e..52681084 100644
--- a/vignettes/Usage.Rmd
+++ b/vignettes/Usage.Rmd
@@ -4,7 +4,7 @@ author: "Wes Hinsley"
 date: "`r Sys.Date()`"
 output: rmarkdown::html_vignette
 vignette: >
-  %\VignetteIndexEntry{stoner}
+  %\VignetteIndexEntry{Usage}
   %\VignetteEngine{knitr::rmarkdown}
   %\VignetteEncoding{UTF-8}
 ---