Revision 7492ec296dab70e3eda5479fedc56f8310cc4417 authored by Peter Ashcroft on 04 March 2021, 08:12:50 UTC, committed by GitHub on 04 March 2021, 08:12:50 UTC
1 parent 2fa6a7f
quarantine.Rmd
---
title: "Quantifying the impact of quarantine duration on COVID-19 transmission"
author:
- "Peter Ashcroft\\textsuperscript{1}, Sonja Lehtinen\\textsuperscript{1}, Daniel C. Angst\\textsuperscript{1}, Nicola Low\\textsuperscript{2}, and Sebastian Bonhoeffer\\textsuperscript{1}"
- "\\textsuperscript{1}*Institute of Integrative Biology, ETH Zurich, Switzerland*,"
- "\\textsuperscript{2}*Institute of Social and Preventive Medicine, University of Bern, Bern, Switzerland*"
output:
pdf_document:
citation_package: natbib
fig_caption: true
latex_engine: pdflatex
header-includes:
\setcounter{secnumdepth}{2}
\usepackage{multirow}
bibliography: 2020-quarantine.bib
biblio-style: apalike
---
```{r setup, include = FALSE}
knitr::opts_chunk$set(
echo = FALSE,
cache = TRUE,
warning = FALSE,
fig.width = 8,
fig.height = 4,
fig.pos = "ht",
eval.after = "fig.cap",
out.extra = ""
)
library(tidyverse)
library(ggplot2)
library(cowplot)
library(RColorBrewer)
library(magick)
library(parallel)
library(stats)
#' Some plot functions
plotTheme <- theme_bw() +
theme(
plot.background = element_rect(fill = "white", colour = "white"),
panel.grid = element_blank(),
panel.grid.major.y = element_line(size = 0.1, colour = "grey90"),
strip.background = element_blank()
)
#' Reduce text size
plotTheme$text$size <- 12
thinLine <- 0.8
thickLine <- 1.2
dodge <- 0.4
my_geom_line <- function(...) {
geom_line(size = thinLine, position = position_dodge(width = dodge), ...)
}
my_geom_point <- function(...) {
geom_point(size = 2*thickLine, position = position_dodge(width = dodge), ...)
}
my_geom_errorbar <- function(...) {
geom_errorbar(aes(ymin = lower, ymax = upper),
width = 0, size = 0.4, alpha = 0.5,
position = position_dodge(width = dodge), ...)
}
dayLabels <- function(x) {
labs <- paste(x, ifelse(x == 1, "day", "days"))
names(labs) <- x
return(labs)
}
#' Parallel loop function
ParLapply <- function(X, FUN, ..., PARALLEL = TRUE, SEED = NULL) {
if (PARALLEL) {
if (is.null(SEED)) SEED <- 111
num.cores <- min(c(length(X), max(detectCores() - 1, 1)))
cl <- makeCluster(num.cores, type = "FORK")
clusterSetRNGStream(cl, SEED)
out <- parLapply(cl, X, FUN, ...)
stopCluster(cl)
return(out)
} else {
out <- lapply(X, FUN, ...)
return(out)
}
}
```
# Abstract
The large number of individuals placed into quarantine because of possible SARS-CoV-2 exposure has high societal and economic costs.
There is ongoing debate about the appropriate duration of quarantine, particularly since the fraction of individuals who eventually test positive is perceived as being low.
We use empirically-determined distributions of incubation period, infectivity, and generation time to quantify how the duration of quarantine affects onward transmission from traced contacts of confirmed SARS-CoV-2 cases and from returning travellers.
We also consider the roles of testing followed by release if negative (test-and-release), reinforced hygiene, adherence, and symptoms in calculating quarantine efficacy.
We show that there are quarantine strategies based on a test-and-release protocol that, from an epidemiological viewpoint, perform almost as well as a 10 day quarantine, but with fewer person days spent in quarantine.
The findings apply to both travellers and contacts, but the specifics depend on the context.
# Introduction
Quarantining individuals with high risk of recent infection is one of the pillars of the non-pharmaceutical interventions to control the ongoing severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) pandemic \citep{kucharski:TheLancetInfectiousDiseases:2020}.
Owing to the large fraction of transmission of SARS-CoV-2 that is pre-symptomatic or asymptomatic
\citep{ashcroft:SMW:2020,buitrago-garcia:PLOSMedicine:2020,ferretti:medRxiv:2020,he:NatMed:2020}, quarantine can prevent a substantial fraction of onward transmission that would not be detected otherwise.
In mathematical modelling studies, it was estimated that thermal screening at airports would allow more than 50\% of infected travellers to enter the general population \citep{quilty:Eurosurveillance:2020,gostic:eLife:2020}, which could have been prevented by mandatory quarantine.
Quarantine is also a fundamental part of the test-trace-isolate-quarantine (TTIQ) intervention strategy to break chains of transmission within a community \citep{salathe:SwissMed.Wkly.:2020}.
With the high or increasing case numbers that are observed in many places around the globe, however, more and more people are being placed into quarantine.
There is ongoing public debate about the appropriateness of quarantine and its duration.
Quarantine lowers onward transmission in two ways: first, preventing transmission prior to symptom onset (with the assumption that symptomatic individuals will isolate) and decreasing overall transmission from persistently asymptomatic individuals.
The appropriate length of quarantine thus depends on both incubation period and the temporal profile of infectiousness.
In theory, quarantine periods could be avoided altogether through widespread and regular testing programmes, but the low sensitivity of reverse transcriptase PCR (RT-PCR) tests, particularly in early infection \citep{kucirka:AnnalsofInternalMedicine:2020}, and limitations on testing capacity in most countries, preclude this approach.
Quarantine has high economic, societal and psychological costs \citep{nicola:InternationalJournalofSurgery:2020,brooks:TheLancet:2020}.
It restricts individual freedoms \citep{parmet:NEJM:2020}, although the level of restriction imposed is generally judged to be proportionate, given the severity of coronavirus disease 2019 (COVID-19).
The low number of individuals placed in quarantine that turn out to be infected is another argument that is given against quarantine.
Individuals are generally placed into quarantine for one of two reasons: either they have been identified as a recent close contact of a confirmed SARS-CoV-2 case by contact tracing, or they have returned from recent travel to an area with community transmission that has been assessed to pose significant epidemiological risk \citep{WHO:quarantine}.
These groups of quarantined individuals differ in two important ways: compared with traced contacts, travel returners may have lower probability of being infected and have less precise information about the likely time of exposure.
This raises the question whether the duration of quarantine should be the same for these two groups of individuals.
To our knowledge there are no published analyses of surveillance data that directly assess the impact of duration of quarantine on SARS-CoV-2 transmission \citep{nussbaumer-streit:Cochrane:2020}.
In this study, we present a mathematical model that allows quantification of the effects of changing quarantine duration.
We use the distributions of incubation time (time from infection to onset of symptoms), infectivity (infectiousness as a function of days since symptom onset), and generation time (difference of time points of infection between infector and infectee).
These distributions have been estimated by \citet{ferretti:medRxiv:2020}, combining multiple empirical studies of documented SARS-CoV-2 transmission pairs \citep{ferretti:Science:2020,xia:medRxiv:2020,cheng:JAMAInternMed:2020,he:NatMed:2020}.
Using the model, we explore the effect of duration of quarantine for both traced contacts of confirmed SARS-CoV-2 cases and for returning travellers on the fraction of prevented onward transmission.
We assess the effects of test-and-release strategies and the time delay between test and result.
These considerations are particularly important given that multiple testing has been shown to be of little benefit \citep{clifford:medRxiv:2020}.
We also address the role of pre-symptomatic patients becoming symptomatic and therefore being isolated independent of quarantine.
Furthermore, as one of the arguments for shortening the duration of quarantine is to increase the number of people complying with the recommendation, we investigate by how much adherence needs to increase to offset the increased transmission due to earlier release from quarantine.
Finally, we assess the role of reinforced individual-level prevention measures, such as mask wearing, for those released early from quarantine.
Making policy decisions about the duration of quarantine fundamentally requires specifying how the effectiveness of quarantine relates to its costs.
The effectiveness can be measured in terms of the overall reduction of transmission, while economic, societal, and individual costs are likely a function of the number of days spent in quarantine.
In addition to the epidemiological outcome, which considers only the reduction in transmission, we also present results based on the ratio of transmission prevented to the average number of days spent in quarantine.
# Results
```{r incubation-distribution-definition-Li, include = F, eval = F}
#' Parameters
incParams <- data.frame(meanlog = 1.42, sdlog = 0.661)
incParams$mean <- exp(incParams$meanlog + (incParams$sdlog^2)/2)
#' Function
getIncubationPeriod <- function(times, params, CDF = FALSE) {
if (CDF) return(plnorm(q = times, meanlog = params$meanlog, sdlog = params$sdlog))
else return(dlnorm(x = times, meanlog = params$meanlog, sdlog = params$sdlog))
}
```
```{r incubation-distribution-definition-Ferretti, include = F}
# #' Code from Ferretti et al.
# inc<-function(x){
# (
# dlnorm(x,meanlog = 1.621, sdlog = 0.418)+ # Lauer
# dlnorm(x,meanlog = 1.425, sdlog = 0.669)+ # Li
# dlnorm(x,meanlog = 1.57, sdlog = 0.65)+ # Bi
# dlnorm(x,meanlog = 1.53, sdlog = 0.464)+ # Jiang???
# dlnorm(x,meanlog = 1.611, sdlog = 0.472)+ # Linton
# dlnorm(x,meanlog = 1.54, sdlog = 0.47)+ # Zhang
# dlnorm(x,meanlog = 1.857, sdlog = 0.547) # Ma
# )/7
# }
# incCum<-function(x){
# (
# plnorm(x,meanlog = 1.621, sdlog = 0.418)+
# plnorm(x,meanlog = 1.425, sdlog = 0.669)+
# plnorm(x,meanlog = 1.57, sdlog = 0.65)+
# plnorm(x,meanlog = 1.53, sdlog = 0.464)+
# plnorm(x,meanlog = 1.611, sdlog = 0.472)+
# plnorm(x,meanlog = 1.54, sdlog = 0.47)+
# plnorm(x,meanlog = 1.857, sdlog = 0.547)
# )/7
# }
# Incubation period distribution is averaged across 7 reported distributions:
incParams <- rbind(
data.frame(study = "Bi", n = 183, meanlog = 1.570, sdlog = 0.650), # Reported
data.frame(study = "Lauer", n = 181, meanlog = 1.621, sdlog = 0.418), # Reported
data.frame(study = "Li", n = 10, meanlog = 1.434, sdlog = 0.661), # From source code of He et al.
data.frame(study = "Linton", n = 158, meanlog = 1.611, sdlog = 0.472), # with(list(mu = 5.6, sigma = 2.8), c(mean = log(mu^2 / sqrt(mu^2 + sigma^2)), sd = sqrt(log(1 + sigma^2 / mu^2))))
data.frame(study = "Ma", n = 587, meanlog = 1.857, sdlog = 0.547), # with(list(mu = 7.44, sigma = 4.39), c(mean = log(mu^2 / sqrt(mu^2 + sigma^2)), sd = sqrt(log(1 + sigma^2 / mu^2))))
data.frame(study = "Zhang", n = 49, meanlog = 1.540, sdlog = 0.470), # Reported
data.frame(study = "Jiang", n = 2015, meanlog = 1.530, sdlog = 0.464) # Unknown...
)
incParams$study <- factor(incParams$study)
#' Mean and SD of this incubation time
incParams$mean <- mean(exp(incParams$meanlog + (incParams$sdlog^2)/2))
#sqrt(mean((exp(incParams$sdlog^2)-1) * exp(2*incParams$meanlog + incParams$sdlog^2))) # SD
#' PDF and CDF of the incubation period
getIncubationPeriod <- function(times, params, CDF = F) {
y <- sapply(levels(params$study), function(study) {
if (CDF) { # Cumulative density function
return(plnorm(q = times,
meanlog = params[params$study == study, "meanlog"],
sdlog = params[params$study == study, "sdlog"]))
} else { # Probability density function
return(dlnorm(x = times,
meanlog = params[params$study == study, "meanlog"],
sdlog = params[params$study == study, "sdlog"]))
}
})
#' Average over the studies
if (length(times) == 1) return(mean(y))
else return(apply(y,1,mean))
}
```
```{r incubation-periods-plot, eval = F}
times <- seq(0,20,0.1)
#' Define the incubation distribution from meta-distribution
incDist <- data.frame(
t = times,
pdf = getIncubationPeriod(times = times, params = incParams),
CDF = getIncubationPeriod(times = times, params = incParams, CDF = TRUE),
study = "combined"
)
#' Calculate incubations distributions for each study individually
incDistSep <- lapply(levels(incParams$study), function(study) {
params <- incParams[incParams$study == study, ]
params$study <- factor(study, levels = study)
data.frame(
t = times,
pdf = getIncubationPeriod(times = times, params = params),
CDF = getIncubationPeriod(times = times, params = params, CDF = TRUE),
study = study
)
})
df <- rbind(incDist, do.call(rbind,incDistSep))
df$study <- factor(df$study, levels = c("combined", levels(incParams$study)))
colours <- c("black", brewer.pal(length(levels(incParams$study)),"Set1"))
pdfPlot <- ggplot(df, aes(x = t, y = pdf, colour = study)) +
geom_line(aes(size = ifelse(study == "combined", 1.5, 1.2))) +
scale_colour_manual(values = colours) +
guides(size = F) +
coord_cartesian(xlim = c(0,20), ylim = c(0,0.23), expand = F) +
labs(x = "incubation period (days)", y = "probability density") +
ggtitle("incubation period distribution") +
plotTheme +
theme(plot.title = element_text(size = plotTheme$text$size, hjust = 0.5),
legend.position = c(0.7,0.8))
cdfPlot <- ggplot(df, aes(x = t, y = CDF, colour = study)) +
geom_line(aes(size = ifelse(study == "combined", 1.5, 1.2))) +
scale_colour_manual(values = colours) +
guides(size = F) +
coord_cartesian(xlim = c(0,20), ylim = c(0,1), expand = F) +
labs(x = "incubation period (days)", y = "cumulative probability") +
ggtitle("cumulative incubation period distribution") +
plotTheme +
theme(plot.title = element_text(size = plotTheme$text$size, hjust = 0.5),
legend.position = "none")
plot_grid(pdfPlot, cdfPlot, align = "hv", axis = "tb", nrow = 1)
```
```{r serial-interval-data, include = F, eval = F}
#' Data from Ferretti (40), Xia (32), Zhang (35), He (66)
info_data <- read.csv("data/Ferretti-TableTransmissionPairs_addendum.tsv", sep = "\t", stringsAsFactors = F)
location <- info_data$Country
info_data[,"Exponential.growth"] <- 0 + (info_data[,"Exponential.growth"] == "YES")
info_data <- cbind(
apply(info_data[,c("T1L","T1R","s1","T2L","T2R","s2","Tr","Exponential.growth")], 2, as.numeric),
info_data[,"Source.manuscript", drop = F]
)
info_data <- as.data.frame(info_data)
#' Exposure intervals for non-tested cases in Cheng
dnont <- read.csv("data/Ferretti-Taiwan_transmission_pair_contact_exposure_v0514.csv")
dnont <- dnont[is.na(dnont$test),]
dnont <- data.frame(
left_exposure = as.numeric(difftime(dnont[,"time2exposure_start"], dnont[,"onset_source"], units = "d")),
right_exposure = as.numeric(difftime(dnont[,"time2exposure_end"], dnont[,"onset_source"], units = "d"))
)
dnont <- dnont[dnont[,"left_exposure"] <= dnont[,"right_exposure"],]
dnont <- dnont[dnont$left_exposure > -100,]
#' Data of symptomatic cases from Cheng (18)
dt <- read.csv("data/Ferretti-Taiwan_transmission_pair_v0514.csv")
dt <- dt[!is.na(dt$onset),]
#
dtrans <- data.frame(
left_exposure = as.numeric(difftime(dt[,"time2exposure_start"], dt[,"onset_source"], units = "d")) - 0.5,
right_exposure = as.numeric(difftime(dt[,"time2exposure_end"], dt[,"onset_source"],units = "d")) + 0.5,
case.noncase = "case",
serial_interval = as.numeric(difftime(dt[,"onset"], dt[,"onset_source"], units = "d"))
)
dtrans$left_exposure[dtrans$left_exposure < -100] <- NA
# Joining Ferretti (40), Xia (32), Zhang (35), He (66), and Cheng (18) datasets
info_data_taiwan <- data.frame(
T1L = NA,
T1R = NA,
s1 = 0,
T2L = dtrans$left_exposure + 0.5,
T2R = dtrans$right_exposure - 0.5,
s2 = dtrans$serial_interval,
Tr = 100,
Exponential.growth = 0,
Source.manuscript = "Cheng"
)
serial_interval_data <- rbind(info_data, info_data_taiwan)
save(serial_interval_data, file = "data/Ferretti-serial_interval_data.RData")
```
```{r serial-interval-data-plot, eval = F}
load("data/Ferretti-serial_interval_data.RData")
data <- serial_interval_data
#' Compute serial interval from symptom onsets
data$SI <- data$s2 - data$s1
data$Source.manuscript <- factor(data$Source.manuscript,
levels = c("Zhang", "He", "Ferretti", "Xia", "Cheng"))
minSI <- min(data$SI) - 1.5
maxSI <- max(data$SI) + 1.5
ggplot() +
geom_histogram(data = data,
aes(x = SI, fill = Source.manuscript),
breaks = seq(from = minSI, to = maxSI)) +
scale_fill_brewer(palette = "Set1") +
coord_cartesian(expand = F) + # no gaps between bars and axes
labs(fill = "Dataset:",
x = "Serial interval (days)",
y = "Count") +
plotTheme +
theme(legend.position = c(0.8, 0.8))
```
```{r fit-infectivity-profile, include = F, eval = F}
load("data/Ferretti-serial_interval_data.RData")
#' Maximum duration of infection (days)
M <- 30
mintW <- 2 * M + 2
#' Growth rate of the epidemic (doubling time of 5 days)
r <- log(2)/5
#' Discrete incubation time distribution (integrated probability per day)
Inc <- getIncubationPeriod(seq(0,2*M) + 0.5, params = incParams, CDF = T) -
getIncubationPeriod(pmax(seq(0,2*M) - 0.5, 0), params = incParams, CDF = T)
#' Student t-distribution for infectivity profile
infectivityProfile <- function(x, pars, CDF = F) {
if (CDF) {
return(pt((x - pars["shift"])/exp(pars["log_scale"]), df = exp(pars["log_df"])))
} else {
return(dt((x - pars["shift"])/exp(pars["log_scale"]), df = exp(pars["log_df"])) / exp(pars["log_scale"]))
}
}
#' Limits for the exposure time integrals (this is evaluated per transmission pair)
getLimits <- function(pair_data) {
v <- list()
v$T1L <- max(pair_data$T1L, pair_data$s1 - 2*M, na.rm = T)
v$T1R <- min(pair_data$T1R, pair_data$s1, pair_data$s2, pair_data$T2R, na.rm = T)
v$T2L <- max(pair_data$T2L, pair_data$s2 - 2*M, v$T1L, na.rm = T)
v$T2R <- min(pair_data$T2R, pair_data$s2, na.rm = T)
return(v)
}
#' Perform the double integral (or sum) per transmission pair
#' to get probability of observing serial interval
getProbSI <- function(pair_data, W, limits) {
sum(
sapply(seq(limits$T1L, limits$T1R), function(t1) {
exp(r * pair_data$Exponential.growth * t1) *
Inc[pair_data$s1 - t1 + 1] *
sum(
sapply(seq(max(t1, limits$T2L), limits$T2R), function(t2) {
W[t2 - pair_data$s1 + mintW] * Inc[pair_data$s2 - t2 + 1]
})
)
})
)
}
#' log-likelihood of a single transmission pair (1 row of data)
getLLH <- function(pair_data, W) {
limits <- getLimits(pair_data)
log(getProbSI(pair_data, W, limits))
}
#' Evaluate the new profile as a function of the optimisation parameters,
#' and then compute the negative log-likelihood for each transmission pair
getTotalLLH <- function(pars, data, w) {
x <- seq(0, (4*M + mintW)) - (mintW - 1)
W <- w(x + 0.5, pars, CDF = T) - w(x - 0.5, pars, CDF = T)
return(sum(apply(data, 1, function(data_row) getLLH(pair_data = as.list(data_row), W))))
}
#' Set up and run optimiser
data <- subset(
serial_interval_data,
select = -c(Source.manuscript),
subset = (serial_interval_data$Source.manuscript %in% c("Ferretti", "Xia", "He", "Cheng"))
)
optim_out <- optim(par = c(shift = 0, log_scale = 0, log_df = 0),
fn = getTotalLLH, data = data, w = infectivityProfile,
control = list(fnscale = -1))
optim_out <- optim(par = optim_out$par,
fn = getTotalLLH, data = data, w = infectivityProfile,
control = list(fnscale = -1))
#' Confidence interval
times <- seq(-10, 10, 0.1)
upper <- rep(NA, length(times))
lower <- rep(NA, length(times))
#' Threshold likelihood value for 95%
minLLH <- optim_out$value - qchisq(0.95, df = length(optim_out$par))/2
#' Parameters of the best fitting model
initpars <- optim_out$par
pars <- initpars
df_param <- data.frame(shift = pars["shift"], scale = exp(pars["log_scale"]),
df = exp(pars["log_df"]), llh = optim_out$value)
#' Create a random vector on a sphere in the 3D parameter space
set.seed(42)
dir <- qnorm(runif(length(optim_out$par)))
dir <- dir/sqrt(sum(dir^2))
i <- 0
while (i < 1000) {
#' Perturb parameters along the vector (Monte Carlo)
pars <- pars + dir*runif(1)*0.2
llh <- getTotalLLH(pars, data, infectivityProfile)
if(llh > minLLH) {
#' Which points make up the confidence interval
infProf <- infectivityProfile(times, pars)
lower <- pmin(infProf, lower, na.rm = T)
upper <- pmax(infProf, upper, na.rm = T)
df_param <- rbind(df_param, data.frame(shift = pars["shift"], scale = exp(pars["log_scale"]), df = exp(pars["log_df"]), llh = llh))
i <- i + 1
} else {
#' Change the direction of perturbation
dir <- qnorm(runif(length(optim_out$par)))
dir <- dir/sqrt(sum(dir^2))
pars <- initpars
}
}
#' Save parameters and likelihood value
infParamsLLH <- df_param
save(infParamsLLH, file = "data/infParamsLLH.RData")
```
```{r plot-infectivity-profile, eval = F}
df_plot <- data_frame(
days = times,
infectivityProfile = infectivityProfile(times, optim_out$par),
incubation = getIncubationPeriod(times, params = incParams),
lower = lower,
upper = upper
)
# new_df <- lapply(seq_len(nrow(df_param)), function(i) {
# data.frame(
# id = factor(i),
# llh = df_param[[i,"llh"]],
# days = times,
# infProf = infectivityProfile(times, c(shift = df_param[[i,"shift"]], log_scale = log(df_param[[i,"scale"]]), log_df = log(df_param[[i,"df"]]))))
# })
# new_df <- do.call(rbind, new_df)
ggplot(data = df_plot, aes(x = days, y = infectivityProfile)) +
geom_ribbon(aes(ymin = lower, ymax = upper), color = "transparent", fill = "grey") +
geom_line() +
geom_vline(xintercept = 0, linetype = "dotted") +
#geom_line(data = new_df, aes(x = days, y = infProf, colour = llh, group = id)) +
#scale_colour_viridis_c() +
labs(x = "TOST (days)", y = "Probability density (per day)") +
coord_cartesian(xlim = c(-7.5, 10), expand = F) +
plotTheme
```
```{r infectivity-profile-definition, include = F}
#' Parameters
#infParams <- data.frame(shift = -0.0747, scale = 1.8567, df = 3.3454)
infParams <- data.frame(shift = -0.0776, scale = 1.857, df = 3.345)
infParams$mean <- infParams$shift/infParams$scale
#' Function
getInfectivityProfile <- function(times, params, CDF = FALSE) {
if (CDF) return(pt(q = (times - params$shift)/params$scale, df = params$df))
else return(dt(x = (times - params$shift)/params$scale, df = params$df)/params$scale)
}
```
```{r fit-generation-time, include = F, eval = F}
load("data/Ferretti-serial_interval_data.RData")
#' Maximum duration of infection (days)
M <- 30
#' Growth rate of the epidemic (doubling time of 5 days)
r <- log(2)/5
#' Discrete incubation time distribution (integrated probability per day)
Inc <- getIncubationPeriod(seq(0,2*M) + 0.5, params = incParams, CDF = T) -
getIncubationPeriod(pmax(seq(0,2*M) - 0.5, 0), params = incParams, CDF = T)
#' Weibull distribution for generation time
generationTime <- function(x, pars, CDF = F) {
if (CDF) {
return(pweibull(x, shape = exp(pars["log_shape"]), scale = exp(pars["log_scale"])))
} else {
return(dweibull(x, shape = exp(pars["log_shape"]), scale = exp(pars["log_scale"])))
}
}
#' Limits for the exposure time integrals (this is evaluated per transmission pair)
getLimits <- function(pair_data) {
v <- list()
v$T1L <- max(pair_data$T1L, pair_data$s1 - 2*M, na.rm = T)
v$T1R <- min(pair_data$T1R, pair_data$s1, pair_data$s2, pair_data$T2R, na.rm = T)
v$T2L <- max(pair_data$T2L, pair_data$s2 - 2*M, v$T1L, na.rm = T)
v$T2R <- min(pair_data$T2R, pair_data$s2, na.rm = T)
return(v)
}
#' Perform the double integral (or sum) per transmission pair
#' to get probability of observing serial interval
getProbSI <- function(pair_data, W, limits) {
sum(
sapply(seq(limits$T1L, limits$T1R), function(t1) {
exp(r * pair_data$Exponential.growth * t1) *
Inc[pair_data$s1 - t1 + 1] *
sum(
sapply(seq(max(t1, limits$T2L), limits$T2R), function(t2) {
W[t2 - t1 + 1] * Inc[pair_data$s2 - t2 + 1]
})
)
})
)
}
#' log-likelihood of a single transmission pair (1 row of data)
getLLH <- function(pair_data, W) {
limits <- getLimits(pair_data)
log(getProbSI(pair_data, W, limits))
}
#' Evaluate the new profile as a function of the optimisation parameters,
#' and then compute the negative log-likelihood for each transmission pair
getTotalLLH <- function(pars, data, w) {
x <- seq(0, 4*M)
W <- w(x + 0.5, pars, CDF = T) - w(pmax(x - 0.5,0), pars, CDF = T)
return(sum(apply(data, 1, function(data_row) getLLH(pair_data = as.list(data_row), W))))
}
#' Set up and run optimiser
data <- subset(
serial_interval_data,
select = -c(Source.manuscript),
subset = (serial_interval_data$Source.manuscript %in% c("Ferretti", "Xia", "He", "Cheng"))
)
optim_out <- optim(par = c(log_shape = 0, log_scale = 0),
fn = getTotalLLH, data = data, w = generationTime,
control = list(fnscale = -1))
optim_out <- optim(par = optim_out$par,
fn = getTotalLLH, data = data, w = generationTime,
control = list(fnscale = -1))
#' Confidence interval
times <- seq(0, 20, 0.1)
upper <- rep(NA, length(times))
lower <- rep(NA, length(times))
#' Threshold likelihood value for 95%
minLLH <- optim_out$value - qchisq(0.95, df = length(optim_out$par))/2
#' Parameters of the best fitting model
initpars <- optim_out$par
pars <- initpars
df_param <- data.frame(shape = exp(pars["log_shape"]), scale = exp(pars["log_scale"]),
llh = optim_out$value)
#' Create a random vector on a sphere in the 3D parameter space
set.seed(422)
dir <- qnorm(runif(length(optim_out$par)))
dir <- dir/sqrt(sum(dir^2))
i <- 0
while (i < 1000) {
#' Perturb parameters along the vector (Monte Carlo)
pars <- pars + dir*runif(1)*0.2
llh <- getTotalLLH(pars, data, generationTime)
if(llh > minLLH) {
#' Which points make up the confidence interval
genTime <- generationTime(times, pars)
lower <- pmin(genTime, lower, na.rm = T)
upper <- pmax(genTime, upper, na.rm = T)
df_param <- rbind(df_param, data.frame(shape = exp(pars["log_shape"]), scale = exp(pars["log_scale"]), llh = llh))
i <- i + 1
} else {
#' Change the direction of perturbation
dir <- qnorm(runif(length(optim_out$par)))
dir <- dir/sqrt(sum(dir^2))
pars <- initpars
}
}
#' Save parameters and likelihood value
genParamsLLH <- df_param
save(genParamsLLH, file = "data/genParamsLLH.RData")
```
```{r plot-generation-time, eval = F}
df_plot <- data_frame(
days = times,
generationTime = generationTime(times, optim_out$par),
incubation = getIncubationPeriod(times, params = incParams),
lower = lower,
upper = upper
)
# new_df <- lapply(seq_len(nrow(df_param)), function(i) {
# data.frame(
# id = factor(i),
# llh = df_param[[i,"llh"]],
# days = times,
# genTime = generationTime(times, c(log_shape = log(df_param[[i,"shape"]]), log_scale = log(df_param[[i,"scale"]]))))
# })
# new_df <- do.call(rbind, new_df)
ggplot(data = df_plot, aes(x = days, y = generationTime)) +
geom_ribbon(aes(ymin = lower, ymax = upper), color = "transparent", fill = "grey") +
geom_line() +
geom_vline(xintercept = 0, linetype = "dotted") +
#geom_line(data = new_df, aes(x = days, y = genTime, colour = llh, group = id)) +
#scale_colour_viridis_c() +
labs(x = "generation time (days)", y = "Probability density (per day)") +
coord_cartesian(xlim = c(0,15), expand = F) +
plotTheme
```
```{r generation-time-distribution-definition, include = F}
#' Parameters
#genParams <- data.frame(shape = 3.2862, scale = 6.1244)
genParams <- data.frame(shape = 3.277, scale = 6.127)
genParams$mean <- genParams$scale * gamma(1 + 1/genParams$shape)
#' Function
getGenDist <- function(times, params, CDF = FALSE) {
if (CDF) return(pweibull(q = times, shape = params$shape, scale = params$scale))
else dweibull(x = times, shape = params$shape, scale = params$scale)
}
```
```{r save-distributions, include = F}
stepSize <- 1e-2
times <- seq(-25, 25, stepSize)
i0 <- which(times == 0)
iMax <- length(times)
#' Define the incubation distribution, which is independent of the generation time parameters
incDist <- data.frame(
t = times,
pdf = getIncubationPeriod(times = times, params = incParams),
CDF = getIncubationPeriod(times = times, params = incParams, CDF = TRUE)
)
#' Infectivity profile
infProfMLE <- data.frame(
t = times,
pdf = getInfectivityProfile(times = times, params = infParams),
CDF = getInfectivityProfile(times = times, params = infParams, CDF = T)
)
#' Now compute the confidence interval
load("data/infParamsLLH.RData")
upper <- rep(NA, length(times))
lower <- rep(NA, length(times))
for (id in seq_len(nrow(infParamsLLH))) {
ip <- getInfectivityProfile(times = times, params = infParamsLLH[id,])
lower <- pmin(ip, lower, na.rm = T)
upper <- pmax(ip, upper, na.rm = T)
}
infProfMLE$lower <- lower
infProfMLE$upper <- upper
#' Generation time distribution
genDistMLE <- data.frame(
t = times,
pdf = getGenDist(times = times, params = genParams),
CDF = getGenDist(times = times, params = genParams, CDF = T)
)
#' Now compute the confidence interval
load("data/genParamsLLH.RData")
upper <- rep(NA, length(times))
lower <- rep(NA, length(times))
for (id in seq_len(nrow(genParamsLLH))) {
gt <- getGenDist(times = times, params = genParamsLLH[id,])
lower <- pmin(gt, lower, na.rm = T)
upper <- pmax(gt, upper, na.rm = T)
}
genDistMLE$lower <- lower
genDistMLE$upper <- upper
#' Integral of g(t')Q(t'+t) from 0 to \infty with 0 \le t < \infty
integralMLE <- data.frame(
t = times,
J = sapply(seq_along(times), function(i) {
g <- incDist[seq(i0, iMax), "pdf"]
Qt <- genDistMLE[seq(i, iMax - i0 + i), "CDF"]
Qt[is.na(Qt)] <- 1
sum(g * Qt) * stepSize
}, USE.NAMES = F)
)
#' Save results
save(times, i0, iMax, stepSize, incDist, infProfMLE, genDistMLE, integralMLE, file = "data/savedDistributions.RData")
```
```{r save-distribution-samples, include = F, eval = F}
#' Load the MLE distributions
load("data/savedDistributions.RData")
#' Infectivity profile
load("data/infParamsLLH.RData")
infProfLLH <- ParLapply(seq_len(nrow(infParamsLLH)), function(id) {
data.frame(
t = times,
pdf = getInfectivityProfile(times = times, params = infParamsLLH[id,]),
CDF = getInfectivityProfile(times = times, params = infParamsLLH[id,], CDF = T)
)
})
#' Generation time distribution
load("data/genParamsLLH.RData")
genDistLLH <- ParLapply(seq_len(nrow(infParamsLLH)), function(id) {
data.frame(
t = times,
pdf = getGenDist(times = times, params = genParamsLLH[id,]),
CDF = getGenDist(times = times, params = genParamsLLH[id,], CDF = T)
)
})
#' Integral of g(t')Q(t'+t) from 0 to \infty with 0 \le t < \infty
integralLLH <- ParLapply(seq_len(nrow(infParamsLLH)), function(id) {
data.frame(
t = times,
J = sapply(seq_along(times), function(i) {
g <- incDist[seq(i0, iMax), "pdf"]
Qt <- genDistLLH[[id]][seq(i, iMax - i0 + i), "CDF"]
Qt[is.na(Qt)] <- 1
sum(g * Qt) * stepSize
}, USE.NAMES = F)
)
})
save(infProfLLH, genDistLLH, integralLLH, file = "data/savedDistributionsLLH.RData")
#' Create a zip file with the data files:
#' - data/savedDistributions.RData
#' - data/savedDistributionsLLH.RData
#' - data/infParamsLLH.RData
#' - data/genParamsLLH.RData
tar(tarfile = "archivedDistributions.tar.gz", files = c("data/savedDistributions.RData", "data/savedDistributionsLLH.RData",
"data/infParamsLLH.RData", "data/genParamsLLH.RData"), compression = "gzip", tar = "tar")
```
```{r getIntegral, include = F}
#' Function for integrals of q
getIntegral <- function(upper, lower, tE, params) {
getGenDist(times = upper - tE, params = params, CDF = T) -
getGenDist(times = lower - tE, params = params, CDF = T)
}
```
```{r utility-function, include = F}
#' Utility function definition
getUtility <- function(s = 1, efficacy, time) {
s * efficacy / time
}
```
```{r falseNeg, include = F}
#' False negative probabilities
falseNeg <- approxfun(x = c(0,1,4,5,6,7,8,9,21), y = c(1,1,0.67,0.38,0.25,0.21,0.20,0.21,0.66))
```
```{r get-CI, include = F}
#' Function for computing confidence intervals
getCI <- function(df.list, variable) {
lower <- df.list[[1]][[variable]]
upper <- df.list[[1]][[variable]]
for (df in df.list) {
lower <- pmin(lower, df[[variable]])
upper <- pmax(upper, df[[variable]])
}
#' Add to df
df <- df.list[[1]]
df$lower <- lower
df$upper <- upper
return(df[,names(df) != variable])
}
```
## Model description
In the absence of quarantine, individuals that are infected with SARS-CoV-2 can infect further individuals in the population.
In the model, the timing of onward transmission from an infected individual is determined by the generation time distribution, which describes the time interval between the infection of an infector and infectee (see Figure \ref{fig:distributions}).
To quantify how much transmission is prevented by quarantining individuals who have been infected with SARS-CoV-2, we need to know the time at which the individual was exposed ($t_E$), as well as when they enter ($t_Q$) and are released from ($t_R$) quarantine.
The fraction of transmission that is prevented by quarantine is then the total transmission probability (i.e. the area under the curve) that lies between $t_Q$ and $t_R$ (Figure \ref{fig:model}).
We refer to this fraction of prevented transmission as quarantine efficacy, and is defined in Eq. \eqref{eq:F} in Materials and methods.
Unless otherwise stated, we assume that adherence to quarantine is 100\%.
```{r model, fig.cap = caption}
caption <- "Quantifying the impact of quarantine using a mathematical model.
Here the y-axis represents the probability of transmission.
These infectivity curves are a schematic representation of the generation time distribution shown in Figure \\ref{fig:distributions}.
A) Traced contacts are exposed to an infector at a known time $t_E=0$ and then enter quarantine at time $t_Q$.
Some transmission can occur prior to quarantine.
Under the standard quarantine protocol, the contact is quarantined until time $t_R$, and no transmission is assumed to occur during this time.
The area under the infectivity curve between $t_Q$ and $t_R$ (blue) is the fraction of transmission that is prevented by quarantine.
Transmission can occur after the individual leaves quarantine.
Under the test-and-release protocol, quarantined individuals are tested at time $t_T$ and released at time $t_R$ if they receive a negative test result.
Otherwise the individual is isolated until they are no longer infectious.
The probability that an infected individual returns a false-negative test result, and therefore is prematurely released, depends on the timing of the test relative to infection ($t_T-t_E$) \\citep{kucirka:AnnalsofInternalMedicine:2020}.
B) For returning travellers, the time of exposure $t_E$ is unknown and we assume infection could have occurred on any day of the trip.
The travellers enter quarantine immediately upon return at time $t_Q=0$, and then leave quarantine at time $t_R$ under the standard quarantine protocol.
Test-and-release quarantine proceeds as in panel A."
ggdraw() + draw_image(image_read_pdf("timeline3.pdf", density = 600), scale = 1.0)
```
Under the standard quarantine strategy, all potentially-exposed individuals are quarantined for the same duration.
An alternative approach is the test-and-release strategy, which uses virological testing during quarantine to release individuals with a negative test result earlier.
Individuals with a positive test result are isolated until they are no longer infectious.
The timing of the test ($t_T$) is important due to the substantial false-negative rate of the RT-PCR test in the early stages of infection \citep{kucirka:AnnalsofInternalMedicine:2020}.
A false-negative test result would release an infected individual into the community prematurely, leading to further transmission (Figure \ref{fig:model}A).
In this case, quarantine efficacy is defined as the expected fraction of transmission that is prevented by quarantine across false-negative and positive testing individuals, as defined in Eq. \eqref{eq:Ftest} in Materials and methods.
As well as the epidemiological benefit of quarantine (i.e. the fraction of transmission prevented by quarantining an infected individual), we can also quantify the economic and societal costs in terms of the expected number of person days spent in quarantine.
We can then define the utility of a quarantine strategy as the ratio between the quarantine efficacy and the average time spent in quarantine, i.e. the transmission prevented per day spent in quarantine, as defined in Eq. \eqref{eq:utility} in Materials and methods.
This utility measure is dependent on the fraction of individuals in quarantine that are infected.
This definition of utility should be considered as an example of such a utility function, but this may not be the best way to quantify quarantine utility.
Details of the calculations used can be found in Materials and methods.
Further extensions to the model, including the role of reinforced hygiene measures, asymptomatic infections, and adherence to quarantine are described in Appendix 1.
## Quarantining traced contacts of confirmed SARS-CoV-2 cases
Traced contacts have a known (last) time of exposure to a confirmed case.
There is usually a delay between this exposure time and the start of quarantine.
Under the standard quarantine protocol, traced contacts are released from quarantine once a number of days have passed after the last exposure time.
In Switzerland, for example, quarantine lasts until 10 days after the last exposure.
Any shortening of a traced contact's quarantine duration will lead to an increase in transmission from that individual if they are infected, but the degree of increase depends on the extent of the shortening.
The expected onward transmission that is prevented by quarantine shows the diminishing return of increasing the quarantine duration (black line in Figure \ref{fig:contacts}A).
Increasing quarantine duration beyond 10 days shows almost no additional benefit (Figure \ref{fig:contacts-quarantine}A): the standard quarantine protocol (here with a three-day delay between exposure and the start of quarantine) can maximally prevent 90.8\% [95\% confidence interval (CI): 79.6\%,97.6\%] of onward transmission from an infected traced contact, while release on day 10 prevents 90.1\% [CI: 76.0\%,97.5\%].
```{r contacts-params, include = F}
#' Parameters to scan over
tE <- 0
tQ <- 3
tR.vals <- seq(tQ,10)
DeltaT.vals <- seq(3,0)
s <- 0.1
tR.compare <- 10
tEnd <- 10
#' Plot colours
cols <- c("#000000",colorRampPalette(brewer.pal(9, "Oranges"))(12)[c(12,10,8,5)])
names(cols) <- c("noTest", DeltaT.vals)
colours <- scale_colour_manual(
name = expression("delay between test and release" ~ (t[R]-t[T])),
values = cols,
labels = c(noTest = "no test", dayLabels(DeltaT.vals)), aesthetics = c("colour","fill"),
guide = guide_legend(title.position = "left", title.hjust = 0.5, nrow = 1)
)
```
```{r contacts-frac, dependson = c("contacts-params"), include = F}
#' Fraction of transmission prevented without testing
getFraction.noTest <- function(params) {
fraction <- ifelse(tR.vals >= tQ,
getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params),
NA)
data.frame(
tE = tE,
tQ = tQ,
tT = NA,
tR = tR.vals,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR.vals - tE,
fraction = fraction
)
}
#' Compute fraction for MLE parameters
frac <- getFraction.noTest(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Fraction of transmission prevented based on test-and-release
getFraction.test <- function(params, r = 0) {
#' Fraction of transmission before and after release
fraction.before <- getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params)
fraction.after <- getIntegral(upper = tEnd, lower = tR.vals, tE = tE, params = params)
#' Now loop over different testing delays
lapply(DeltaT.vals, function(DeltaT) {
tT.vals <- tR.vals - DeltaT
fraction <- ifelse(tT.vals >= tQ,
fraction.before +
(1 - falseNeg(tT.vals - tE) + r * falseNeg(tT.vals - tE)) * fraction.after,
NA)
data.frame(
tE = tE,
tQ = tQ,
tT = tT.vals,
tR = tR.vals,
DeltaT = factor(DeltaT, levels = c("noTest", DeltaT.vals)),
release = tR.vals - tE,
fraction = fraction
)
}) %>% bind_rows()
}
#' Compute fraction for MLE parameters
fracTest <- getFraction.test(params = genParams, r = 0)
#' Evaluate fraction from test-and-release for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction.test(params = genParamsLLH[id,], r = 0)
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
fracTest <- merge(fracTest, CI, sort = F)
#' Combine standard and test-and-release
combined <- rbind(frac,fracTest)
#' Maximum preventable fraction of transmission
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE, params = genParams)
maxPreventableLLH <- sapply(seq_len(nrow(genParamsLLH)), function(id) {
getIntegral(upper = Inf, lower = tQ, tE = tE, params = genParamsLLH[id,])
})
maxPreventableCI <- c(lower = min(maxPreventableLLH), upper = max(maxPreventableLLH))
#' Plot
fracPlot <- ggplot(combined, aes(x = release, y = fraction, colour = DeltaT, fill = DeltaT)) +
geom_hline(yintercept = maxPreventable, colour = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,1)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(3,10), ylim = c(0,1)) +
labs(x = expression("days until release after exposure" ~ (t[R]-t[E])), y = "fraction of transmission\nprevented by quarantine") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
#' Extract and remove legend
legend <- get_legend(fracPlot)
fracPlot <- fracPlot + theme(legend.position = "None")
print(maxPreventable)
print(maxPreventableCI)
print(frac[frac$tR == 10, c("fraction","lower","upper")])
print(fracTest[fracTest$tT == 6 & fracTest$tR == 6, c("fraction","lower","upper")])
print(fracTest[fracTest$tT == 5 & fracTest$tR == 7, c("fraction","lower","upper")])
#' Save source data
write.csv(combined, file = "data/Figure 2-source data 1.csv", row.names = F)
```
```{r contacts-utility, dependson = c("contacts-params"), include = F}
#' Function to calculate relative utility of standard quarantine
getRelUtility.noTest <- function(params) {
fraction.compare <- getIntegral(upper = tR.compare, lower = tQ, tE = tE, params = params)
utility.compare <- getUtility(s = s, efficacy = fraction.compare, time = tR.compare - tQ)
fraction <- getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params)
utility <- ifelse(tR.vals >= tQ,
getUtility(s = s, efficacy = fraction, time = tR.vals - tQ),
NA)
data.frame(
s = s,
tE = tE,
tQ = tQ,
tT = NA,
tR = tR.vals,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR.vals - tE,
relUtility = utility / utility.compare
)
}
#' Relative utility for MLE parameters
relUtility <- getRelUtility.noTest(params = genParams)
#' Evaluate relative utility for all parameter sets \theta
relUtilityLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relUtilityLLH, "relUtility")
#' Add CI to results
relUtility <- merge(relUtility, CI, sort = F)
#' Relative utility of test-and-release strategy compared to standard
getRelUtility.test <- function(params, r = 0) {
#' Standard quarantine strategy results
fraction.compare <- getIntegral(upper = tR.compare, lower = tQ, tE = tE, params = params)
utility.compare <- getUtility(s = s, efficacy = fraction.compare, time = tR.compare - tQ)
#' Test-and-release results
frac.before <- getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params)
frac.after <- getIntegral(upper = tEnd, lower = tR.vals, tE = tE, params = params)
time.before <- tR.vals - tQ
time.after <- tEnd - tR.vals
lapply(DeltaT.vals, function(DeltaT) {
tT.vals <- tR.vals - DeltaT
fraction <- ifelse(tT.vals >= tQ,
frac.before +
(1 - falseNeg(tT.vals - tE) + r * falseNeg(tT.vals - tE)) * frac.after,
NA)
utility <- getUtility(s = s, efficacy = fraction,
time = time.before + s * (1 - falseNeg(tT.vals - tE)) * time.after)
data.frame(
s = s,
tE = tE,
tQ = tQ,
tT = tT.vals,
tR = tR.vals,
DeltaT = factor(DeltaT, levels = c("noTest",DeltaT.vals)),
release = tR.vals - tE,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}
#' Relative utility using the MLE parameters
relUtilityTest <- getRelUtility.test(params = genParams, r = 0)
#' Evaluate relative utility for all parameter sets \theta
relUtilityTestLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility.test(params = genParamsLLH[id,], r = 0)
})
#' Compute CI
CI <- getCI(relUtilityTestLLH, "relUtility")
#' Add CI to results
relUtilityTest <- merge(relUtilityTest, CI, sort = F)
#' Combine standard and test-and-release
relUtilityCombined <- rbind(relUtility, relUtilityTest)
#' Plot
utilityPlot <- ggplot(relUtilityCombined, aes(x = release, y = relUtility, colour = DeltaT, fill = DeltaT)) +
#geom_vline(xintercept = tR.compare - tE, color = "darkgrey", size = thickLine) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,1)) +
coord_cartesian(xlim = c(3,10), ylim = c(0,3)) +
labs(x = expression("days until release after exposure" ~ (t[R]-t[E])), y = "\nrelative utility of quarantine") +
plotTheme + theme(legend.position = "None")
print(relUtilityTest[relUtilityTest$tT == 3 & relUtilityTest$DeltaT == 0, c("relUtility","lower","upper")])
print(relUtilityTest[relUtilityTest$tT == 5 & relUtilityTest$tR == 7, c("relUtility","lower","upper")])
print(relUtilityTest[relUtilityTest$tT == 6 & relUtilityTest$tR == 6, c("relUtility","lower","upper")])
#' Save source data
write.csv(relUtilityCombined, file = "data/Figure 2-source data 2.csv", row.names = F)
```
```{r contacts, dependson = c(-1,-2,-3), fig.cap = caption}
caption <- "Quantifying the impact of quarantine for traced contacts.
A) The fraction of transmission that is prevented by quarantining an infected contact.
Quarantine begins at time $t_Q=3$ after exposure at time $t_E=0$, i.e. there is a three day delay between exposure and the start of quarantine.
Under the standard quarantine protocol (black), individuals are released without being tested [Eq. \\eqref{eq:F}].
The test-and-release protocol (colours) requires a negative test result before early release, otherwise individuals remain isolated until they are no longer infectious (day 10) [Eq. \\eqref{eq:Ftest}].
Colour intensity represents the delay between test and release (from 0 to 3 days).
The grey line represents the maximum attainable prevention by increasing the time of release, while keeping $t_Q=3$ fixed.
B) The relative utility of the quarantine scenarios in panel A compared to the standard protocol 10 day quarantine [Eq. \\eqref{eq:relUtility}].
Utility is defined as the fraction of transmission prevented per day spent in quarantine.
The grey line represents equal utilities (relative utility of 1).
We assume that the fraction of individuals in quarantine that are infected is $10\\%$, and that there are no false-positive test results.
Error bars reflect the uncertainty in the generation time distribution."
plot_grid(
plot_grid(fracPlot, utilityPlot,
nrow = 1, axis = "tb", align = "hv", labels = "AUTO"),
legend,
ncol = 1, rel_heights = c(1,0.2)
)
```
The maximum attainable prevention also applies to the test-and-release strategy: the onward transmission prevented under a test-and-release strategy will always be below this level (coloured lines in Figure \ref{fig:contacts}A).
This is because of the chance of prematurely releasing an infectious individual who received a false-negative test result.
On the other hand, it is always better to test a person prior to release from quarantine, so that individuals with asymptomatic and pre-symptomatic infections can be detected and prevented from being released.
Hence, these scenarios provide upper and lower bounds for the efficacy of the test-and-release strategy.
The fraction of transmission that is prevented increases if we test later in quarantine, because we not only increase the duration of quarantine but also reduce the false-negative probability.
The delay between testing and release from quarantine can have a substantial effect on the efficacy.
Current laboratory-based RT-PCR tests have a typical turnaround of 24-48 hours \citep{quilty:medRxiv:2020}.
This delay is primarily operational, and so could be reduced by increasing test throughput.
There are also rapid antigen-detection tests, which can provide same day results, but with lower sensitivity and specificity than RT-PCR tests \citep{guglielmi:Nature:2020}.
Here we assume that tests have the same sensitivity and specificity regardless of the delay to result.
Comparing to a test with 2 days delay until result, we observe that using a rapid test with same day release can reduce the quarantine duration of individuals with a negative test result by one day while maintaining the same efficacy (Figure \ref{fig:contacts}A): the extra accuracy gained by waiting one extra day until testing balances the increased transmission caused by reducing the duration.
E.g., a rapid test on day six has roughly the same efficacy (80.5\% [CI: 67.9\%,88.7\%]) as testing on day five and releasing on day seven (82.3\% [CI: 68.2\%,93.4\%]), while shortening the quarantine duration of individuals with a negative test result from seven to six days.
In Figure \ref{fig:contacts} we have assumed a fixed delay of three days between exposure and the start of quarantine.
Shortening this delay increases the maximum efficacy of quarantine because pre-quarantine transmission is reduced (Figure \ref{fig:contacts-quarantine}A).
If the duration of quarantine is longer than 10 days, then little can be gained in terms of prevention by quarantining for longer, but reducing the delay between exposure and quarantine does lead to increased efficacy.
Note that we have assumed that the contact was infected at the last time of exposure.
If there have been multiple contacts between them and the index case, then transmission may have occurred earlier and we would overestimate the efficacy of quarantine.
For the standard quarantine strategy, the duration of quarantine is independent of whether individuals in quarantine are infected.
Therefore, the utility of the standard quarantine strategy (i.e. the ratio of efficacy to duration) is directly proportional to the fraction of individuals in quarantine that are infected.
By comparing two different standard quarantine strategies through their relative utility (i.e. the ratio of the utilities), we can eliminate the dependence on the fraction of infecteds in quarantine (see Materials and methods).
Therefore, the argument that we should shorten quarantine because of the low probability of quarantined individuals being infected is misguided in this situation.
By calculating the relative utility for the standard quarantine strategy compared to the baseline 10 day quarantine, we observe that there is a quarantine strategy (release after seven days) which maximises the ratio between the fraction of transmission prevented and the number of days spent in quarantine (black line in Figure \ref{fig:contacts}B).
The optimal strategy lies between six to eight days if we vary the delay between exposure and the start of quarantine (Figure \ref{fig:contacts-quarantine}B).
Under the test-and-release quarantine protocol the average time spent in quarantine is dependent on the fraction of infecteds in quarantine; only the infected individuals can test positive and face a longer period of isolation (i.e. we assume there are no false-positive test results).
Hence the utility of the test-and-release strategy, as well as the relative utility of test-and-release compared to the standard quarantine protocol, is dependent on the fraction of individuals in quarantine that are infected.
In Figure \ref{fig:contacts}B we fix the fraction of infecteds in quarantine to 10\%.
By calculating the relative utility for the test-and-release quarantine strategies shown in Figure \ref{fig:contacts}A compared to the baseline 10 day quarantine, we see that testing-and-releasing before day 10 always increases the utility (Figure \ref{fig:contacts}B).
Testing on day five and releasing test-negative individuals on day seven has a relative utility of 1.53 [CI: 1.45,1.62] compared to a standard 10 day quarantine.
Reducing the delay between test and result leads to a corresponding increase in utility: a rapid test (zero delay between test and result) on day six has a relative utility of 1.90 [CI: 1.83,1.98] for an almost equivalent efficacy.
In Figure \ref{fig:contacts} we have made the following assumptions:
i) individuals released from quarantine have -- in the post-quarantine phase -- the same transmission probability as individuals who were not quarantined;
ii) adherence to quarantine is 100\%;
and iii) the transmission prevented by quarantine for cases who develop symptoms is attributed to quarantine.
We now relax these assumptions to assess their impact on quarantine efficacy.
Reinforced prevention measures post-quarantine, where individuals who are released from quarantine must adhere to strict hygiene and social distancing protocols until 10 days after exposure have passed, can reduce post-quarantine transmission.
Considering a 50\% reduction of post-quarantine transmission leads to large increases in both efficacy and utility for early testing strategies, but with diminishing returns as the release date is increased towards day 10 (Figure \ref{fig:contacts-reduced}; see Appendix 1: Reinforced prevention measures after early release).
Note that we assume no contribution to the number of days spent in quarantine in the utility function due to mask wearing and social distancing in the post-release phase.
Adherence to quarantine is unlikely to be 100\% and could depend on the proposed duration of quarantine.
For simplicity we treat adherence to quarantine as a binary variable: a fraction of individuals adhere to quarantine completely for the proposed duration, while the remaining fraction do not undergo any quarantine.
We now ask: by how much would the fraction of those who adhere to quarantine have to increase to maintain the efficacy of quarantine if the duration is shortened?
In the absence of testing during quarantine, shortening from 10 to five days would require almost three times as many individuals to adhere to the quarantine guidelines in order to maintain the same overall efficacy (relative adherence 2.90 [CI: 2.15,4.36]; black line in Figure \ref{fig:contacts-further}A).
This three-fold increase is not possible if adherence to the 10 day strategy is already above 33\%, as the maximum adherence cannot exceed 100\%; the required increase in adherence grows rapidly as quarantine is shortened and soon becomes infeasible.
Hence the argument of shortening quarantine to increase adherence is of limited use.
Shortening to seven days (without testing) may be effective provided that adherence can increase by 30\% (relative adherence 1.30 [CI: 1.08,1.55]).
Under the test-and-release strategy, however, the efficacy of the standard 10 day quarantine can be matched with release on day five or six if adherence is also increased by 30\%.
Releasing earlier than day five would seemingly be infeasible given the sharp increase in adherence required.
```{r contacts-further-adherence, include = F}
#' Parameters to scan over
tE <- 0
tQ <- 3
tR.vals <- seq(tQ,10)
DeltaT.vals <- seq(3,0)
tR.compare <- 10
tEnd <- 10
#' Plot colours
cols <- c("#000000",colorRampPalette(brewer.pal(9, "Oranges"))(12)[c(12,10,8,5)])
names(cols) <- c("noTest", DeltaT.vals)
colours <- scale_colour_manual(
name = expression("delay between test and release" ~ (t[R]-t[T])),
values = cols,
labels = c(noTest = "no test", dayLabels(DeltaT.vals)), aesthetics = c("colour","fill"),
guide = guide_legend(title.position = "top", title.hjust = 0.5,
nrow = 2, byrow = T)
)
#' Function to calculate relative adherence of standard quarantine
getRelAdherence.noTest <- function(params) {
fraction.compare <- getIntegral(upper = tR.compare, lower = tQ, tE = tE, params = params)
fraction <- ifelse(tR.vals > tQ,
getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params),
NA)
data.frame(
tE = tE,
tQ = tQ,
tT = NA,
tR = tR.vals,
tR.compare = tR.compare,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR.vals - tE,
relAdherence = fraction.compare / fraction
)
}
#' Relative adherence for MLE parameters
relAdherence <- getRelAdherence.noTest(params = genParams)
#' Evaluate relative utility for all parameter sets \theta
relAdherenceLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelAdherence.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relAdherenceLLH, "relAdherence")
#' Add CI to results
relAdherence <- merge(relAdherence, CI, sort = F)
#' Relative adherence of test-and-release strategy compared to standard
getRelAdherence.test <- function(params, r = 0) {
#' Standard quarantine strategy results
fraction.compare <- getIntegral(upper = tR.compare, lower = tQ, tE = tE, params = params)
#' Test-and-release results
frac.before <- getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params)
frac.after <- getIntegral(upper = tEnd, lower = tR.vals, tE = tE, params = params)
lapply(DeltaT.vals, function(DeltaT) {
tT.vals <- tR.vals - DeltaT
fraction <- ifelse(tT.vals >= tQ,
frac.before +
(1 - falseNeg(tT.vals - tE) + r * falseNeg(tT.vals - tE)) * frac.after,
NA)
data.frame(
tE = tE,
tQ = tQ,
tT = tT.vals,
tR = tR.vals,
tR.compare = tR.compare,
DeltaT = factor(DeltaT, levels = c("noTest",DeltaT.vals)),
release = tR.vals - tE,
relAdherence = fraction.compare / fraction
)
}) %>% bind_rows()
}
#' Relative adherence using the MLE parameters
relAdherenceTest <- getRelAdherence.test(params = genParams, r = 0)
#' Evaluate relative utility for all parameter sets \theta
relAdherenceTestLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelAdherence.test(params = genParamsLLH[id,], r = 0)
})
#' Compute CI
CI <- getCI(relAdherenceTestLLH, "relAdherence")
#' Add CI to results
relAdherenceTest <- merge(relAdherenceTest, CI, sort = F)
#' Combine standard and test-and-release
relAdherenceCombined <- rbind(relAdherence, relAdherenceTest)
#' Plot
adherencePlot <- ggplot(relAdherenceCombined, aes(x = release, y = relAdherence, colour = DeltaT, fill = DeltaT)) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
#geom_vline(xintercept = nCompare, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10)) +
coord_cartesian(xlim = c(3,10), ylim = c(0,4)) +
labs(x = expression("days until release after exposure" ~ (t[R]-t[E])), y = "relative adherence required\nto maintain quarantine efficacy") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
print(relAdherenceCombined[relAdherenceCombined$tR == 5, c("tT","relAdherence","lower","upper")])
print(relAdherenceCombined[relAdherenceCombined$tR == 6, c("tT","relAdherence","lower","upper")])
print(relAdherenceCombined[relAdherenceCombined$tR == 7, c("tT","relAdherence","lower","upper")])
#' Save source data
write.csv(relAdherenceCombined, file = "data/Figure 3-source data 1.csv", row.names = F)
```
```{r contacts-further-asymptomatic, include = F}
tE <- 0
tQ <- 3
tR.vals <- seq(0,15)
a.vals <- seq(0,1,0.25)
delay <- 0
#' Fraction of contacts that occur during quarantine for a given level of asymptomatic vs symptomatic
getFraction <- function(params) {
fraction.asymptomatic <- ifelse(tR.vals >= tQ,
getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params),
NA)
fraction.symptomatic <- sapply(tR.vals, function(tR) {
if(tR < tQ) return(NA)
tS.vals <- seq(tE, tR - delay, stepSize)
symptoms.before <- stepSize * sum(
getIncubationPeriod(times = tS.vals - tE, params = incParams, CDF = F) *
getIntegral(upper = pmax(tS.vals + delay, tQ), lower = tQ, tE = tE, params = params)
)
symptoms.after <- (1 - getIncubationPeriod(times = tR - delay - tE, params = incParams, CDF = T)) *
getIntegral(upper = tR, lower = tQ, tE = tE, params = params)
return(symptoms.before + symptoms.after)
})
lapply(a.vals, function(a) {
data.frame(
tE = tE,
tQ = tQ,
tR = tR.vals,
delay = delay,
release = tR.vals - tE,
a = factor(a, levels = a.vals),
fraction = a * fraction.asymptomatic + (1 - a) * fraction.symptomatic
)
}) %>% bind_rows()
}
#' Compute fraction for MLE parameters
frac <- getFraction(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Plot
labs <- scales::percent(as.numeric(levels(frac$a)))
names(labs) <- levels(frac$a)
asymptomaticPlot <- ggplot(frac, aes(x = release, y = fraction, colour = a, fill = a)) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
scale_colour_viridis_d(option = "viridis", direction = -1, end = 0.9, aesthetics = c("colour","fill"),
name = expression("fraction asymptomatic" ~ (a)),
labels = labs,
guide = guide_legend(title.position = "top", title.hjust = 0.5, nrow = 2, byrow = T)) +
scale_x_continuous(breaks = seq(0,10)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(3,10), ylim = c(0,1)) +
labs(x = expression("days until release after exposure" ~ (t[R]-t[E])), y = "fraction of transmission\nprevented by quarantine") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
print(frac[frac$a == "0.25",])
#' Save source data
write.csv(frac, file = "data/Figure 3-source data 2.csv", row.names = F)
```
```{r contacts-further, dependson = c(-1,-2), fig.height = 4.5, fig.cap = caption}
caption <- "How adherence and symptoms affect quarantine efficacy for traced contacts.
A) The fold-change in adherence to a new quarantine strategy that is required to maintain efficacy of the baseline 10 day standard strategy.
Quarantine strategies are the same as in Figure \\ref{fig:contacts} (standard = black, test-and-release = colours).
The grey line represents equal adherence (relative adherence of 1).
B) The impact of symptomatic cases on the fraction of total onward transmission per infected traced contact that is prevented by standard (no test) quarantine [Appendix 1-Eq. \\eqref{eq:Fasymptomatic}].
We assume that symptomatic individuals will immediately self-isolate at symptom onset.
The time of symptom onset is determined by the incubation period distribution (see Figure \\ref{fig:distributions}D).
The curve for $100\\%$ asymptomatic cases corresponds to the black curve in Figure \\ref{fig:contacts}A.
As in Figure \\ref{fig:contacts}, we fix the time of exposure at $t_E=0$ and the time of entering quarantine at $t_Q=3$ days.
Error bars reflect the uncertainty in the generation time distribution."
plot_grid(adherencePlot, asymptomaticPlot,
nrow = 1, axis = "tb", align = "hv", labels = "AUTO"
)
```
As a final consideration, we note that our quantification of the fraction of transmission prevented by quarantine is more relevant to individuals with persistently asymptomatic SARS-CoV-2 infection than to those who develop symptoms during quarantine and are subsequently isolated.
If symptomatic cases go into isolation once symptoms appear then quarantine has no further impact on transmission after symptom onset, as these cases would anyway be isolated.
To account for this, we can modify the model such that cases are removed from the infectious pool upon symptom onset (see Appendix 1).
For example, in a fully-asymptomatic population a 10 day quarantine can prevent 90.1\% [CI: 76.0\%,97.5\%] of transmission.
However, if 25\% of cases are asymptomatic then only 50.8\% [CI: 42.8\%,56.5\%] of transmission is prevented by quarantine, while 39.3\% is prevented by the self-isolation of symptomatic cases (Figure \ref{fig:contacts-further}B).
We assume self-isolation occurs immediately after symptom onset, but any delay between symptom onset and self-isolation would mean that more transmission is prevented by quarantine (Figure \ref{fig:contacts-asymptomatic-delay}).
The fraction of transmission prevented by quarantine is an increasing function of the fraction of asymptomatic cases (Figure \ref{fig:contacts-further}B).
This means that we likely overestimate the efficacy of quarantine as we are also counting transmission that could be prevented by isolation following symptom onset.
Furthermore, we have assumed that the false-negative rate is the same between symptomatic and asymptomatic cases.
If the test is less sensitive (higher false-negative probability) for asymptomatic cases, then quarantine efficacy would be further reduced.
## Quarantining returning travellers
The rules for whether travellers returning from abroad are quarantined are frequently changed according to the epidemiological scenario in the travel destination and/or in the home country.
A high risk of infection while abroad due to high prevalence, or the possibility of returning with a new virological variant, can lead to the imposition or reinstatement of quarantine measures \citep{russell:TheLancetPublicHealth:2020}.
Countries that have already eliminated the infection may be even stricter in their quarantine approach to prevent new community-transmission clusters from being seeded.
Here we do not discuss these scenarios or the concept of relative risk, we simply quantify how effective quarantine strategies would be at preventing transmission if the returning traveller was infected while abroad.
Should quarantine rules be instated or modified, these results can help determine the optimal quarantine duration and/or testing strategy.
The timing of infection of a traveller during a trip abroad is generally unknown.
We assume that infection could have happened on each day of the trip with equal probability.
Quarantine begins immediately upon return, which we refer to as day zero, and lasts for a number of days (e.g. currently 10 days in Switzerland) from this timepoint (Figure \ref{fig:model}B).
We consider the fraction of \emph{local} transmission that is prevented by quarantine.
That is, the fraction of the transmission that could occur in the local country that is prevented by quarantine [Eq. \eqref{eq:Flocal}].
For a seven day trip, as in Figure \ref{fig:travellers}, the maximum transmission that could occur in the local country is 73.3\% [CI: 65.7\%,80.3\%].
The remaining infectivity potential was already used up before arrival.
A standard (no test) 10 day quarantine will prevent 99.9\% [CI: 98.0\%,100.0\%] of local transmission if the individual was infected during a seven day trip (Figure \ref{fig:travellers}A).
There is little benefit to gain by increasing the duration of quarantine beyond 10 days.
On the other hand, standard quarantine efficacy decreases quickly as the duration is shortened.
The test-and-release strategy can improve the efficacy of shorter-duration quarantines.
Testing on day five and releasing on day seven (to account for test processing delays) performs similarly to 10 day quarantine, preventing 98.5\% [CI: 95.5\%,99.6\%] of local transmission (Figure \ref{fig:travellers}A).
Testing and releasing on day six (i.e. no delay between test and result) still prevents 97.8\% [CI: 94.4\%,99.0\%] of local transmission.
Hence, if the rapid test has the same sensitivity and specificity as the laboratory-based RT-PCR test, then the duration of quarantine of individuals with a negative test result can be shortened by one day with minimal loss in efficacy compared to a test with a 48h turnaround.
The timing of the test can have significant impact on prevented transmission; late testing reduces the false-negative probability but increases the stay in quarantine.
An important consequence of this is that testing on arrival is a poor strategy for limiting transmission:
testing and releasing on day zero would prevent only 35.2\% [CI: 35.1\%,35.3\%] of local transmission, while testing on arrival and releasing after two days prevents only 54.1\% [CI: 49.5\%,59.4\%].
As was the case for the traced contacts, the fraction of local transmission prevented by standard quarantine bounds the efficacy of the test-and-release quarantine strategy from below (Figure \ref{fig:travellers}A).
```{r travellers-params, include = F}
#' Parameters to scan over
y <- 7
tE.vals <- seq(-y,0)
tQ <- 0
tR.vals <- seq(0,10)
DeltaT.vals <- seq(3,0)
s <- 0.1
tR.compare <- 10
tEnd <- 10
#' Plot colours
cols <- c("#000000",colorRampPalette(brewer.pal(9, "Oranges"))(12)[c(12,10,8,5)])
names(cols) <- c("noTest", DeltaT.vals)
colours <- scale_colour_manual(
name = expression("delay between test and release" ~ (t[R]-t[T])),
values = cols,
labels = c(noTest = "no test", dayLabels(DeltaT.vals)), aesthetics = c("colour","fill"),
guide = guide_legend(title.position = "left", title.hjust = 0.5, nrow = 1)
)
```
```{r travellers-frac, dependson = c("travellers-params"), include = F}
#' Fraction of transmission prevented without testing
getFraction.noTest <- function(params) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
lapply(tR.vals, function(tR) {
fraction <- mean(ifelse(maxPreventable < 1e-10, 1, getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
data.frame(
y = y,
tQ = tQ,
tT = NA,
tR = tR,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR - tQ,
fraction = fraction
)
}) %>% bind_rows()
}
#' Compute fraction for MLE parameters
frac <- getFraction.noTest(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Fraction of transmission prevented based on test-and-release
getFraction.test <- function(params, r = 0) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
lapply(DeltaT.vals, function(DeltaT) {
lapply(tR.vals, function(tR) {
tT <- tR - DeltaT
fraction.before <- getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)
fraction.after <- getIntegral(upper = tEnd, lower = tR, tE = tE.vals, params = params)
fraction <- mean(ifelse(maxPreventable < 1e-10, 1, (fraction.before + (1 - falseNeg(tT-tE.vals) + r * falseNeg(tT-tE.vals)) * fraction.after)/maxPreventable))
data.frame(
y = y,
tQ = tQ,
tT = tT,
tR = tR,
DeltaT = factor(DeltaT, levels = c("noTest", DeltaT.vals)),
release = tR - tQ,
fraction = fraction
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Compute fraction for MLE parameters
fracTest <- getFraction.test(params = genParams, r = 0)
#' Evaluate fraction from test-and-release for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction.test(params = genParamsLLH[id,], r = 0)
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
fracTest <- merge(fracTest, CI, sort = F)
#' Combine standard and test-and-release
combined <- rbind(frac, fracTest)
#' Maximum preventable fraction of transmission
maxPreventable <- mean(getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = genParams))
maxPreventableLLH <- sapply(seq_len(nrow(genParamsLLH)), function(id) {
mean(getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = genParamsLLH[id,]))
})
maxPreventableCI <- c(fraction = maxPreventable, lower = min(maxPreventableLLH), upper = max(maxPreventableLLH))
#' Plot
fracPlot <- ggplot(combined, aes(x = release, y = fraction, colour = DeltaT, fill = DeltaT)) +
geom_hline(yintercept = 1, colour = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(0,10), ylim = c(0,1)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "fraction of local transmission\nprevented by quarantine") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
#' Extract and remove legend
legend <- get_legend(fracPlot)
fracPlot <- fracPlot + theme(legend.position = "None")
print(maxPreventableCI)
print(frac[frac$tR == 10, c("tT","tR","fraction","lower","upper")])
print(fracTest[fracTest$tT == 5 & fracTest$tR == 7, c("tT","tR","fraction","lower","upper")])
print(fracTest[fracTest$tT == 6 & fracTest$tR == 6, c("tT","tR","fraction","lower","upper")])
print(fracTest[fracTest$tT == 0 & fracTest$tR == 0, c("tT","tR","fraction","lower","upper")])
print(fracTest[fracTest$tT == 0 & fracTest$tR == 2, c("tT","tR","fraction","lower","upper")])
#' Save source data
write.csv(combined, file = "data/Figure 4-source data 1.csv", row.names = F)
```
```{r travellers-utility, dependson = c("travellers-params"), include = F}
#' Function to calculate relative utility of standard quarantine
getRelUtility.noTest <- function(params) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
fraction.compare <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR.compare, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
utility.compare <- getUtility(s = s, efficacy = fraction.compare, time = tR.compare - tQ)
lapply(tR.vals, function(tR) {
fraction <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
utility <- ifelse(tR > tQ,
getUtility(s = s, efficacy = fraction, time = tR - tQ),
NA)
data.frame(
y = y,
tQ = tQ,
tT = NA,
tR = tR,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR - tQ,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}
#' Relative utility for MLE parameters
relUtility <- getRelUtility.noTest(params = genParams)
#' Evaluate relative utility for all parameter sets \theta
relUtilityLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relUtilityLLH, "relUtility")
#' Add CI to results
relUtility <- merge(relUtility, CI, sort = F)
#' Relative utility of test-and-release strategy compared to standard
getRelUtility.test <- function(params, r = 0) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
#' Standard quarantine strategy results
fraction.compare <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR.compare, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
utility.compare <- getUtility(s = s, efficacy = fraction.compare, time = tR.compare - tQ)
#' Test-and-release results
lapply(DeltaT.vals, function(DeltaT) {
lapply(tR.vals, function(tR) {
tT <- tR - DeltaT
frac.before <- getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)
frac.after <- getIntegral(upper = Inf, lower = tR, tE = tE.vals, params = params)
fraction <- mean(ifelse(maxPreventable < 1e-10, 1,
(frac.before + (1 - falseNeg(tT-tE.vals) + r * falseNeg(tT-tE.vals)) * frac.after)/maxPreventable))
time.before <- tR - tQ
time.after <- tEnd - tR
time <- mean(time.before + s * (1 - falseNeg(tT-tE.vals)) * time.after)
utility <- ifelse(tT >= tQ,
getUtility(s = s, efficacy = fraction,
time = time),
NA)
data.frame(
y = y,
tQ = tQ,
tT = tT,
tR = tR,
DeltaT = factor(DeltaT, levels = c("noTest",DeltaT.vals)),
release = tR - tQ,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Relative utility using the MLE parameters
relUtilityTest <- getRelUtility.test(params = genParams, r = 0)
#' Evaluate relative utility for all parameter sets \theta
relUtilityTestLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility.test(params = genParamsLLH[id,], r = 0)
})
#' Compute CI
CI <- getCI(relUtilityTestLLH, "relUtility")
#' Add CI to results
relUtilityTest <- merge(relUtilityTest, CI, sort = F)
#' Combine standard and test-and-release
relUtilityCombined <- rbind(relUtility, relUtilityTest)
#' Plot
utilityPlot <- ggplot(relUtilityCombined, aes(x = release, y = relUtility, colour = DeltaT, fill = DeltaT)) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
coord_cartesian(xlim = c(0,10), ylim = c(0,4)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "\nrelative utility of quarantine") +
plotTheme + theme(legend.position = "None")
print(relUtilityTest[relUtilityTest$tT == 0 & relUtilityTest$tR == 0, c("tT","tR","relUtility","lower","upper")])
#' Save source data
write.csv(relUtilityCombined, file = "data/Figure 4-source data 2.csv", row.names = F)
```
```{r travellers, dependson = c(-1,-2,-3), fig.cap = caption}
caption <- "Quantifying the impact of quarantine for returning travellers.
A) The fraction of \\emph{local} transmission that is prevented by quarantining an infected traveller returning from a seven day trip.
Quarantine begins upon return at time $t_Q=0$, and we assume exposure could have occurred at any time during the trip, i.e. $-7 \\le t_E \\le 0$.
Under the standard quarantine protocol (black), individuals are released without being tested [Eq. \\eqref{eq:FlocalAve}].
The test-and-release protocol (colours) requires a negative test result before early release, otherwise individuals remain isolated until they are no longer infectious (day 10).
Colour intensity represents the delay between test and release (from 0 to 3 days).
While extended quarantine can prevent 100\\% of \\emph{local} transmission (grey line), this represents 73.3\\% [CI: 65.7\\%,80.3\\%] of the \\emph{total} transmission potential (see Figure \\ref{fig:travellers-quarantine}A).
The remaining transmission occurred before arrival.
B) The relative utility of the quarantine scenarios in panel A compared to the standard protocol 10 day quarantine [Eq. \\eqref{eq:relUtility}].
Utility is defined as the \\emph{local} fraction of transmission that is prevented per day spent in quarantine.
The grey line represents equal utilities (relative utility of 1).
We assume that the fraction of individuals in quarantine that are infected is $10\\%$, and that there are no false-positive test results.
Error bars reflect the uncertainty in the generation time distribution."
plot_grid(
plot_grid(fracPlot, utilityPlot,
nrow = 1, axis = "tb", align = "hv", labels = "AUTO"),
legend,
ncol = 1, rel_heights = c(1,0.2)
)
```
We again measure the utility of quarantine by calculating the efficacy (local transmission prevented across all individuals in quarantine, assuming 100% adherence) per day spent in quarantine, and then comparing these utilities for different quarantine strategies to the utility of the standard 10 day quarantine through the relative utility (Figure \ref{fig:travellers}B).
In the absence of testing the duration of quarantine, and hence the relative utility, is independent of the fraction of individuals in quarantine that are infected.
For travellers returning from a seven-day trip, the relative utility is a decreasing function of quarantine duration (black line in Figure \ref{fig:travellers}B).
The maximum utility strategy would then be to shorten quarantine as much as possible.
As was the case for traced contacts, under the test-and-release quarantine protocol the average time spent in quarantine, the utility, and the relative utility compared to the standard 10 day quarantine will depend on the fraction of individuals in quarantine that are infected.
This fraction may change depending on disease prevalence at the travel destination and the duration of travel.
E.g., the infected fraction of travellers returning from a long stay in a high-risk country is likely to be higher than the infected fraction of travellers returning from a short stay to a low risk country.
In Figure \ref{fig:travellers}B we keep this fraction fixed at 10\%.
Early testing greatly reduces the average duration of quarantine and hence leads to increased utility, despite the low fraction of transmission that is prevented (coloured lines in Figure \ref{fig:travellers}B).
The average quarantine duration increases linearly with the fraction of infecteds in quarantine [Eq. \eqref{eq:duration} in Materials and methods].
The ratio of quarantine efficacy to the average quarantine duration will also increase, such that quarantine is of higher utility if the fraction of infecteds is higher.
However, the relative utility of test-and-release quarantine compared to the standard 10 day protocol will decrease and approach one as the fraction of infecteds increases.
Hence, if the disease prevalence among those returning from travel abroad is high then test-and-release may not bring substantial benefits over the standard 10 day protocol.
Our assumption that infection occurs with uniform probability across each day of a trip leads to some interesting results.
Returning travellers that have been infected on a short journey will have, on average, used up less of their infectivity potential by the time they return than a traveller who was infected on a long journey.
Hence, the \emph{total} transmission that can be prevented by a long quarantine period (e.g. 10 days) upon arrival, is greater for short trips (Figure \ref{fig:travellers-quarantine}A).
When considering the fraction of \emph{local} transmission that can be prevented by quarantine, then shorter quarantine durations have a greater impact on long than short trips (Figure \ref{fig:travellers-quarantine}C).
Again, this is because, on average, the traveller on a long trip would have been exposed earlier and they will be infectious for a shorter time period after arrival.
If an individual traveller is to be quarantined then the optimum duration of quarantine, based on our metric of utility, would depend on the duration of their travel, with shorter journeys requiring longer quarantine (Figure \ref{fig:travellers-quarantine}B,D).
This might be counter-intuitive because individuals who have been on longer journeys to high risk countries have a higher probability of being infected.
The absolute utility (transmission prevented by quarantine across all individuals in quarantine divided by the average quarantine duration) of quarantining such individuals could indeed be higher than for individuals returning from shorter journeys.
However, here, we are not considering the question of whether to quarantine or not, but we are assuming that the individual is quarantined and are trying to optimise the duration of quarantine in response to the expected infection dynamics.
We observe an almost-linear response between quarantine duration and the relative utility of the standard (no test) quarantine: for every day that quarantine is shortened, we see the same additive increase in relative utility (black line in Figure \ref{fig:travellers}B).
This almost-linear response is coincidental to the seven day trip duration: longer or shorter trips show non-linear responses (Figure \ref{fig:travellers-quarantine}D).
Trips shorter than seven days have a maximum relative utility of between four and seven days, while trips longer than seven days have maximum utility for maximally-shortened quarantine durations.
Enforcing additional hygiene and social distancing guidelines post-quarantine can increase both efficacy and utility, but with diminishing returns as the release date is increased (Figure \ref{fig:travellers-reduced}).
As discussed for traced contacts, the loss of efficacy due to shortening quarantine could be offset by increasing quarantine adherence.
Shortening from 10 to five days would require adherence to increase by 20\% (relative adherence 1.20 [CI: 1.12,1.35]) in order to maintain the same overall efficacy (Figure \ref{fig:travellers-further}A).
With test-and-release this required increase in adherence is even smaller.
We note that the change in adherence required to balance a change in efficacy for shortened quarantine durations is dependent on the travel duration, with short travel durations requiring a greater increase in adherence compared with longer travel durations.
# Discussion
Quarantine is one of the most important measures in controlling the ongoing SARS-CoV-2 epidemic due to the large fraction of presymptomatic and asymptomatic transmission.
A quarantine period of 10 days from exposure, as currently implemented in Switzerland, is long enough to prevent almost all onward transmission from infected contacts of confirmed cases from the point of entering quarantine: increasing the duration of quarantine beyond 10 days has no extra benefit.
Reducing the delay to quarantining individuals increases the fraction of total transmission that is preventable.
The same 10 day quarantine duration will prevent almost all local onward transmission from infected travel returners from the time of arrival, independent of their travel duration.
Any decrease in the duration of quarantine of an infected individual will result in increased onward transmission.
Furthermore, our analyses suggest that this increase in transmission cannot realistically be compensated by increased adherence for significantly shortened quarantine (fewer than five days).
However, there are diminishing returns for each day that we add to quarantine: increasing the duration from 10 days has a negligible effect in terms of reduced transmission.
One therefore has to assess how much human cost, measured in terms of days spent in quarantine, we are willing to spend to prevent disease transmission.
By comparing the ratio of prevented transmission to quarantine duration, we have shown that maximal utility strategies can exist.
This ratio is maximised for quarantine durations of six to eight days after exposure for traced contacts, and potentially less for returning travellers depending on their duration of travel.
Importantly, under this metric the fraction of individuals in quarantine that are infected does not affect the optimal duration of quarantine.
Therefore, the argument that we should shorten quarantine because of the low probability of being infected is misguided under our definition of utility and in the absence of testing during quarantine.
A test-and-release strategy will lead to a lower average quarantine duration across infected and non-infected individuals.
However, due to the considerable false-negative probability of the RT-PCR test \citep{kucirka:AnnalsofInternalMedicine:2020}, this strategy also leads to increased transmission as infectious individuals are prematurely released.
Nevertheless, test-and-release strategies prevent more transmission than releasing without testing, and hence test-and-release increases the utility of quarantine.
Reducing the delay between test and result leads to further reduced transmission and increased utility, and reinforcing individual prevention measures after release is also effective for short quarantine periods.
The ratio of transmission prevented versus days spent in quarantine is only one possible definition of utility.
Defining the appropriate function is ultimately a policy question: the economic, societal, and individual costs are likely a function of the number of days spent in quarantine, but we cannot determine the shape of this function.
Furthermore, the local epidemiological situation will dictate which metric of quarantine efficacy is to be optimised.
In situations where the goal is to prevent the (re)introduction of SARS-CoV-2, one should focus on maximising the reduction of transmission (and hence minimising the transmission risk).
If the virus is already endemic then considering the trade-off between transmission reduction and quarantine duration could determine the optimum strategy.
Another perspective is that the utility of preventing transmission is crucially dependent on whether it brings the effective reproductive number under one.
Ultimately, bringing the reproductive number below one requires a combination of effective measures including isolation, physical distancing, hygiene, contact tracing, and quarantine \citep{kucharski:TheLancetInfectiousDiseases:2020}.
Effective quarantine is only possible in the presence of efficient contact tracing to find the potentially exposed individuals in a short time, as well as surveillance of disease prevalence to identify high-risk travel.
Further reducing the time taken to quarantine a contact after exposure and reducing the delay between test and result will allow average quarantine durations to be shorter, which increases the benefit to cost ratio of quarantine.
The scenarios of returning travellers and traced contacts of confirmed SARS-CoV-2 cases differ in the probability of having been exposed and infected and on the information available about the likely window of exposure.
The impact of quarantining returning travellers depends on the duration of travel, and whether we consider the local prevention of transmission or the total transmission prevented by quarantine.
However, a single test done immediately after return can only prevent a small fraction of the transmission from a returning traveller because of the false-negative rate of the RT-PCR test early in infection.
Therefore testing should be postponed until as late as possible, and utilising rapid tests could be crucial if their performance characteristics are acceptable.
This same principle also applies to traced contacts.
Our findings are aligned with those of two recent simulation studies which estimate the role that quarantine plays in limiting transmission from returning travellers \citep{clifford:medRxiv:2020} and from traced contacts \citep{quilty:medRxiv:2020}.
Our results are based on the latest estimates of the generation time distribution of COVID-19 \citep{ferretti:medRxiv:2020}.
Potential limitations to our approach could be that these distributions may change throughout the epidemic, particularly depending on how people respond to symptoms \citep{ali:Science:2020}.
Furthermore, these distributions, and also the test sensitivity profile, could be different between persistently asymptomatic and symptomatic individuals \citep{buitrago-garcia:PLOSMedicine:2020}, which ultimately lead to an overestimation of how much transmission is prevented by quarantine.
In addition, we have primarily assumed that symptom onset during quarantine has no impact on quarantine efficacy.
However, this symptomatic transmission should not be counted towards the efficacy of quarantine as the infected individual should already self-isolate after symptom onset.
We have quantified this effect and have shown that this assumption leads us to overestimate quarantine efficacy.
For travellers, another consideration is that lengthy quarantine is seen as a deterrent to travel to high risk areas \citep{IATA:2020}.
Any shortening of quarantine may lead to an increase in travel volume, potentially leading to a compounded increase in disease transmission.
In the absence of empirical data about the effectiveness of different durations of quarantine, mathematical modelling can be used objectively to explore the fraction of onward transmission by infected contacts or returning travellers that can be prevented.
However, determining the optimal quarantine strategy to implement depends on the impact that shortening the duration of quarantine has on individuals, society, and the economy, versus how much weight is assigned to a consequential increase in transmission.
Both the individual, societal and economic impact, as well as the weight of transmission increase, will have to be considered based on the current epidemiological situation.
We have shown that there are quarantine strategies based on a test-and-release protocol that, from an epidemiological viewpoint, perform almost as well as the standard 10 day quarantine, but with a lower cost in terms of person days spent in quarantine.
This applies to both travellers and contacts, but the specifics depend on the context.
# Materials and methods
## Quantifying the benefit of quarantine
For an infected individual who was exposed at time $t_E$, the fraction of transmission that is prevented by the standard quarantine strategy is given by the area under the generation time distribution, $q(t)$ (Figure \ref{fig:distributions}B), between the times at which the individual enters ($t_Q$) and leaves ($t_R$) quarantine \citep{grantz:medRxiv:2020}, i.e.
\begin{equation}
F_{\rm qs}(t_E,t_Q,t_R) = \int_{t_Q}^{t_R} {\rm d}t \, q(t-t_E).
\label{eq:F}
\end{equation}
The duration of time that the individual spends in quarantine is then $D_{\rm qs}=t_R-t_Q$.
The test-and-release strategy uses virological testing during quarantine to release individuals with a negative test result and to place those with a positive test result into isolation.
As illustrated in Figure \ref{fig:model}A, a test is issued at time $t_T \ge t_Q$.
If the test is negative, the individual is released when the test result arrives at time $t_R$.
Otherwise, the individual is isolated until they are no longer infectious.
One challenge with this strategy is the high probability of a false-negative RT-PCR test result (i.e. an infectious individual is prematurely released into the community).
As reported by \citet{kucirka:AnnalsofInternalMedicine:2020}, the false-negative rate is 100\% on days 0 and 1 post infection, falling to 67\% (day 4), 38\% (day 5), 25\% (day 6), 21\% (day 7), 20\% (day 8), and 21\% (day 9), before rising to 66\% on day 21.
We use linear interpolation and label this function $f(t)$, the false-negative probability on day $t$ after infection.
The fraction of transmission prevented by quarantining an infected individual under the test-and-release strategy is:
\begin{equation}
F_{\rm qtr}(t_E,t_Q,t_T,t_R) = f(t_T-t_E) \int_{t_Q}^{t_R} {\rm d}t \, q(t-t_E) + \bigl[1-f(t_T-t_E)\bigr] \int_{t_Q}^{t_{\rm end}} {\rm d}t \, q(t-t_E),
%= \int_{t_Q}^{t_R} {\rm d}t \, q(t-t_E) + \bigl[1-f(t_T-t_E)\bigr] \int_{t_R}^{t_{\rm end}} {\rm d}t \, q(t-t_E),
\label{eq:Ftest}
\end{equation}
where the first term captures the fraction of individuals who receive a false-negative test result and are released at time $t_R$, and the second term captures individuals who return a positive test and are subsequently isolated until they are no longer infectious at time $t_{\rm end}$.
A further challenge with this false-negative rate is that it was calculated by \citet{kucirka:AnnalsofInternalMedicine:2020} from symptomatic cases only.
Here we assume that this test sensitivity profile also applies to asymptomatic cases.
Quarantine is applied pre-emptively, such that we do not know the infection status of individuals when they enter quarantine.
If only a fraction $s$ of the individuals that are quarantined are infected, then the average reduction in transmission across all individuals in quarantine is $sF$, where $F$ is the fraction of transmission prevented when an infected individual is quarantined [i.e. Eq. \eqref{eq:F} or Eq. \eqref{eq:Ftest}].
For the standard quarantine protocol, the average number of days spent in quarantine is independent of $s$: all individuals are quarantined for the same duration.
However, under the test-and-release protocol, only the individuals who are actually infected can test positive and remain isolated after $t_R$.
All non-infected individuals ($1-s$) will receive a negative test result and are released at time $t_R$.
Among the infected individuals in quarantine ($s$), a fraction $f(t_T-t_E)$ will receive a false-negative test result and will be released at time $t_R$, while the remaining fraction [$1-f(t_T-t_E)$] will receive a positive test result and are isolated until they are no longer infectious.
Hence the average number of days spent in quarantine for test-and-release is
\begin{equation}
\begin{aligned}
D_{\rm qtr} &= (1-s)(t_R-t_Q) + s\left[f(t_T-t_E)(t_R-t_Q) + \bigl[1-f(t_T-t_E)\bigr](t_{\rm end}-t_Q)\right] \\
&= (t_R-t_Q) + s[1-f(t_T-t_E)] (t_{\rm end}-t_R),
\end{aligned}
\label{eq:duration}
\end{equation}
where $s[1-f(t_T-t_E)]$ is the fraction of quarantined individuals who return a positive test result.
We see that the average test-and-release quarantine duration increases linearly with the fraction of individuals in quarantine that are infected ($s$).
Model parameters and timepoints are summarised in Table \ref{tab:glossary}.
\begin{table}[h]
\begin{tabular}{p{1cm} p{6cm} p{6cm}}
Value & Definition & Notes \\ \hline
$q(t)$ & Generation time distribution & Weibull distribution: shape~=~3.277, scale = 6.127 \\ \hline
$t_E$ & Time of exposure & $t_E=0$ for traced contacts \\ \hline
$t_Q$ & Time at which quarantine begins & $t_Q=0$ for returning travellers \\ \hline
$t_R$ & Time of release from quarantine & - \\ \hline
$t_T$ & Time of test & - \\ \hline
$t_{\rm end}$ & End of infectiousness & $t_{\rm end}=t_E+10$ days \\ \hline
$g(t)$ & Incubation period distribution & Meta-log-normal distribution (Appendix 1: Distribution parameters) \\ \hline
$t_S$ & Time of symptom onset & $t_S=t_E+$incubation period \\ \hline
$D_{\rm qs}$ & Realised average duration of standard quarantine & $D_{\rm qs}=t_R-t_Q$ \\ \hline
+$D_{\rm qtr}$ & Realised average duration of test-and-release quarantine & See Eq. \eqref{eq:duration} \\ \hline
+$F_{\rm qs}(\cdot)$, $F_{\rm qtr}(\cdot)$ & Quarantine efficacy; the fraction of transmission prevented by quarantining an infected individual & See Eqs. \eqref{eq:F} and \eqref{eq:Ftest} \\ \hline
$y$ & Duration of travel journey (days) & - \\ \hline
$s$ & Fraction of individuals in quarantine that are infected & - \\ \hline
$f(t)$ & Probability of returning a false-negative test result if tested $t$ days after exposure & From \citet{kucirka:AnnalsofInternalMedicine:2020} \\ \hline
$r$ & Reduction of transmission under reinforced prevention measures post-quarantine & - \\ \hline
$\alpha(D)$ & Probability to adhere to quarantine of duration $D$ & - \\ \hline
$a$ & Fraction of persistently asymptomatic cases & - \\ \hline
$\Delta$ & Delay between symptom onset and isolation (days) & See Appendix 1: Persistently asymptomatic infections and the role of self-isolation\\ \hline
\end{tabular}
\caption{Summary of terms used in the mathematical model.}
\label{tab:glossary}
\end{table}
## Transmission reduction versus days spent in quarantine
One possible metric to relate the effectiveness of quarantine to its negative impact on society is to consider the ratio between the amount of overall transmission prevented and the number of person days spent in quarantine.
We refer to this ratio as the utility of quarantine.
Concretely, for an efficacy $F$ [$F_{\rm qs}$ or $F_{\rm qtr}$ as defined by Eqs. \eqref{eq:F} or \eqref{eq:Ftest}, respectively], fraction of individuals in quarantine that are infected $s$, and average time spent in quarantine $D$ ($D_{\rm qs}$ or $D_{\rm qtr}$), we define the utility as
\begin{equation}
U(s,F,D) = \frac{sF}{D}.
\label{eq:utility}
\end{equation}
We can then compare the utility of two quarantine strategies by calculating the relative utility, i.e. the ratio between the two utilities:
\begin{equation}
{\rm RU}(s,F,D,F^*,D^*) = \dfrac{s F/D}{s F^*/D^*} = \frac{F/D}{F^*/D^*},
\label{eq:relUtility}
\end{equation}
where $F$ and $D$ are the efficacy and duration of quarantine of a new strategy, and $F^*$ and $D^*$ are the efficacy and duration of the baseline quarantine strategy to which we compare.
The efficacies $F$ and $F^*$ in Eq. \eqref{eq:relUtility} are independent of fraction of individuals in quarantine that are infected $s$.
For the standard quarantine strategy, the durations $D=D_{\rm qs}$ and $D^*=D_{\rm qs}^*$ are also independent of $s$, and hence the relative utility of the standard quarantine strategy is independent of $s$.
For the test-and-release strategy, however, the duration is a linearly-increasing function of $s$ [$D=D_{\rm qtr}(s)$; Eq. \eqref{eq:duration}].
Hence the relative utility of the test-and-release strategy is dependent on $s$:
\begin{equation}
{\rm RU}[s,F_{\rm qtr},D_{\rm qtr}(s),F_{\rm qs}^*,D_{\rm qs}^*]
= \frac{F_{\rm qtr}/D_{\rm qtr}(s)}{F_{\rm qs}^*/D_{\rm qs}^*}.
\label{eq:relUtilityTest}
\end{equation}
In Appendix 1 we show that the relative utility of the test-and-release quarantine strategy is a decreasing function of $s$.
## Traced contacts versus returning travellers
We consider the scenarios of a traced contact and a returning traveller differently, because the values of $t_E$, $t_Q$, and $t_R$ are implemented differently in each case.
### Traced contacts
Following a positive test result, a confirmed index case has their recent close contacts traced.
From contact tracing interviews, we know the date of last exposure between index case and a contact ($t_E$), which we assume is the time of infection of the contact.
The contacts begin quarantine at time $t_Q \ge t_E$.
The delay between exposure and the start of quarantine represents the sum of the delay to the index case receiving a positive test result and the further delay to tracing the contacts.
Under the standard quarantine protocol the traced contacts are quarantined for a number of days after their last exposure.
For example, in Switzerland quarantine lasts until $t_R=t_E+10$ days, but may be longer or shorter depending on individual states' regulations.
Note that the actual time spent in quarantine is $D_{\rm qs} = t_R-t_Q$ days, which is typically shorter than 10 days due to the delay between exposure and the start of quarantine.
For convenience, we set $t_E=0$ for the traced contacts, without loss of generality.
### Returning travellers
Unlike traced contacts, we generally do not know when travellers were (potentially) exposed.
This means that quarantine starts from the date that they return ($t_Q=0$) and lasts until time $t_R$ (Figure \ref{fig:model}B).
For simplicity, we assume that a traveller was infected with uniform probability at some time over a travel period of duration $y$ days.
For each possible exposure time $-y \le t_E \le 0$ during the trip, we can compute the fraction of transmission prevented using Eq. \eqref{eq:F}, and then take the average over the possible exposure times.
This represents the average fraction of the \emph{total} transmission potential that is prevented by quarantining this traveller:
\begin{equation}
\overline{F}_{\rm qs}^{\rm (total)}(y,t_R) = \frac{1}{y+1}\sum_{t_E=-y}^0 \int_0^{t_R} {\rm d}t \, q(t-t_E),
\label{eq:FtotalAve}
\end{equation}
where we have used $t_Q=0$.
For each exposure time $-y \le t_E \le 0$ we can also compute the \emph{local} fraction of transmission prevented by quarantine, which is the fraction of transmission prevented by quarantine divided by the maximum amount of transmission that could occur in the local country, i.e.
\begin{equation}
F_{\rm qs}^{\rm (local)}(t_E,t_R) = \dfrac{\int_0^{t_R} {\rm d}t \, q(t-t_E)}{\int_0^\infty {\rm d}t \, q(t-t_E)},
\label{eq:Flocal}
\end{equation}
where we have again used $t_Q=0$.
Taking the average over the possible exposure times $-y \le t_E \le 0$ we have
\begin{equation}
\overline{F}_{\rm qs}^{\rm (local)}(y,t_R) = \frac{1}{y+1}\sum_{t_E=-y}^0 F_{\rm qs}^{\rm (local)}(t_E,t_R).
\label{eq:FlocalAve}
\end{equation}
## Interactive app
To complement the results in this manuscript, and to allow readers to investigate different quarantine scenarios, we have developed an online interactive application.
This can be found at \url{https://ibz-shiny.ethz.ch/covidDashboard/quarantine}.
\clearpage
# Appendix 1 {-}
\setcounter{equation}{0}
\setcounter{figure}{0}
\setcounter{table}{0}
\renewcommand\theHequation{appendix.\arabic{equation}}
# Utility and relative utility of test-and-release quarantine {-}
From Eq. \eqref{eq:utility} of Materials and methods, we can write the utility of the test-and-release strategy as
\begin{equation}
U(s,F_{\rm qtr},D_{\rm qtr}(s)) = \frac{sF_{\rm qtr}}{D_{\rm qtr}(s)},
\label{eq:utilityTest}
\end{equation}
where $F_{\rm qtr}$ is the quarantine efficacy [Eq. \eqref{eq:Ftest} in Materials and methods], $s$ is the fraction of individuals in quarantine that are infected $s$, and $D_{\rm qtr}(s)$ is the average time spent in quarantine [Eq. \eqref{eq:duration} in Materials and methods].
The duration $D_{\rm qtr}(s)$ is a linear function of $s$, which we can write simply as $D_{\rm qtr}(s) = ms+b$ [from Eq. \eqref{eq:duration} we have $m = [1-f(t_T-t_E)] (t_{\rm end}-t_R)$ and $b = (t_R-t_Q)$].
We now ask, how does this utility change if we increase $s$?
Taking the derivative of Eq. \eqref{eq:utilityTest} with respect to $s$, we recover
\begin{equation}
\frac{{\rm d}U}{{\rm d}s} = \frac{{\rm d}}{{\rm d}s} \left(\frac{sF_{\rm qtr}}{ms+b}\right) = \frac{bF_{\rm qtr}}{(ms+b)^2} > 0.
\end{equation}
Hence, any increase in $s$ leads to an increase in utility.
The relative utility of the test-and-release quarantine strategy to the standard quarantine strategy is defined in Eq. \eqref{eq:relUtilityTest} in Materials and methods.
Again taking the derivative with respect to $s$, we recover
\begin{equation}
\frac{{\rm d}\,{\rm RU}}{{\rm d}s} = \frac{{\rm d}}{{\rm d}s} \left(\frac{F_{\rm qtr}/(ms+b)}{F_{\rm qs}^*/D_{\rm qs}^*}\right) = \frac{-m D_{\rm qs}^* F_{\rm qtr}}{F_{\rm qs}^*(ms+b)^2} < 0.
\end{equation}
Hence, any increase in $s$ leads to a decrease in the relative utility of the test-and-release strategy compared to the standard quarantine strategy.
# Reinforced prevention measures after early release {-}
We further consider the possibility that individuals who are released early from quarantine are asked to strengthen hygiene, mask wearing, and social distancing protocols until the end of the infectious period.
We assume that these practices reduce transmission by a fraction $r$, such that the onward transmission prevented by quarantining an infected individual and reinforcing hygiene measures is
\begin{equation}
F_{\rm qs}^{\rm (reduced)}(t_E,t_Q,t_R) = F_{\rm qs}(t_E,t_Q,t_R) + r \int_{t_R}^{t_{\rm end}} {\rm d}t \, q(t-t_E),
\label{eq:Freduced1}
\end{equation}
for the standard quarantine protocol, and
\begin{equation}
F_{\rm qtr}^{\rm (reduced)}(t_E,t_Q,t_T,t_R) = F_{\rm qtr}(t_E,t_Q,t_T,t_R) + r f(t_T-t_E) \int_{t_R}^{t_{\rm end}} {\rm d}t \, q(t-t_E),
\label{eq:Freduced2}
\end{equation}
for the test-and-release protocol.
# Adherence to quarantine {-}
For the majority of our results we have assumed that quarantine is completely adhered too.
Because of this assumption we will overestimate the efficacy of quarantine at the level of the population, as it is likely that adherence will be less than 100\%.
Adherence could be included as a time-varying property of an individual, such that the probability that an individual follows the quarantine guidelines is high at the beginning of quarantine, but is waning as the duration spent in quarantine increases.
However, for simplicity, we assume adherence is binary; either an individual completes the full duration of quarantine, or they do not enter quarantine at all.
We denote the probability to adhere to quarantine as $\alpha(D)$, which depends on the quarantine duration $D$.
The average fraction of transmission prevented by (standard) quarantine is then $s \alpha(t_R-t_Q) F_{\rm qs}(t_E,t_Q,t_R)$, where $s$ is the fraction of individuals in quarantine that are infected, which we assume is independent of quarantine duration, and we have used $D=D_{\rm qs}=t_R-t_Q$.
We do not know the adherence probabilities $\alpha(D)$.
However, for two quarantine strategies with release dates $t_R$ and $t_R^*$ to have the same efficacy they must satisfy
\begin{equation}
\begin{aligned}
s \alpha(t_R-t_Q) F_{\rm qs}(t_E,t_Q,t_R) &= s \alpha(t_R^*-t_Q) F_{\rm qs}(t_E,t_Q,t_R^*)\\
\Rightarrow \quad \frac{\alpha(t_R-t_Q)}{\alpha(t_R^*-t_Q)} &= \frac{F_{\rm qs}(t_E,t_Q,t_R^*)}{F_{\rm qs}(t_E,t_Q,t_R)}.
\end{aligned}
\end{equation}
That is to say, the change in the fraction of transmission prevented by quarantine must be compensated by an inverse change in the adherence: A strategy which prevents half as much transmission as another would require adherence to be doubled to be equally effective.
We therefore define the required relative adherence as
\begin{equation}
{\rm RA}_{\rm qs}(t_E,t_Q,t_R,t_R^*) = \frac{F_{\rm qs}(t_E,t_Q,t_R^*)}{F_{\rm qs}(t_E,t_Q,t_R)}.
\label{eq:adherence}
\end{equation}
This definition of relative adherence is directly extended to the test-and-release strategy, which we compare to the baseline standard protocol:
\begin{equation}
{\rm RA}_{\rm qtr}(t_E,t_Q,t_T,t_R,t_R^*) = \frac{F_{\rm qs}(t_E,t_Q,t_R^*)}{F_{\rm qtr}(t_E,t_Q,t_T,t_R)}.
\end{equation}
# Persistently asymptomatic infections and the role of self-isolation {-}
If an individual develops symptoms, is tested, and ultimately tests positive while in quarantine, we can remove them from the infectious pool as they would have to isolate themselves.
Importantly, this symptomatic individual would be removed from the infectious pool whether they have been placed in quarantine or not.
Therefore, this symptomatic transmission should not be counted towards the efficacy of quarantine.
Let $a$ be the fraction of asymptomatic cases, who will be quarantined using the standard strategy from time $t_Q$ until $t_R$.
The symptomatic cases (which make up a fraction $1-a$ of cases) will develop symptoms at time $t_S$, as described by the incubation period distribution shown in Figure \ref{fig:distributions}D.
We assume that the symptomatic cases would be isolated shortly after they develop symptoms at time $t_S+\Delta$, so these individuals are effectively quarantined until time $\min(t_R,t_S+\Delta)$.
We further assume equal transmissibility of persistently asymptomatic and symptomatic infections and that both are described by the same generation time distribution.
This assumption might be an overestimate as onward transmission from persistently asymptomatic cases is less than onward transmission from symptomatic cases \citep{buitrago-garcia:PLOSMedicine:2020}.
For each traced contact who is put into quarantine, the fraction of infections that would be prevented by quarantine is
\begin{equation}
\begin{aligned}
F_{\rm qs}^{\rm (asymp)}(t_E,t_Q,t_R,a,\Delta) &=
a \int_{t_Q}^{t_R} {\rm d}t \, q(t-t_E) \\
&+ (1-a) \int_{t_E}^\infty {\rm d}t_S \, g(t_S-t_E) \int_{t_Q}^{\min(t_R,t_S+\Delta)} {\rm d}t \, q(t-t_E),
\end{aligned}
\label{eq:Fasymptomatic}
\end{equation}
where $g(t)$ is the incubation period distribution, and the outer integral over $t_S$ is the averaging over the possible times of symptom onset.
Note that this formulation assumes that the timing of onward transmission is independent of the incubation period (see \citet{ferretti:medRxiv:2020} and \cite{Lehtinen:Interface:2021} for further discussion of this assumption).
Unless otherwise stated, we assume $a=1$ and $\Delta=0$.
# Confidence intervals {-}
The primary source of uncertainty in the outcomes of this model come from the generation time distribution, which is inferred from the empirical serial interval distribution combined with the incubation period distribution \citep{ferretti:medRxiv:2020}.
Following \citet{ferretti:medRxiv:2020}, we use a likelihood ratio test to extract sample parameter sets for the generation time distribution that lie within the 95\% confidence interval.
Concretely, we first identify the parameter set $\hat{\theta}$ for the generation time distribution which maximises the likelihood of observing the empirical serial interval distribution.
The likelihood function and the fitting process is described in detail by \citet{ferretti:medRxiv:2020}.
The generation time distribution is described by a Weibull distribution (with $n=2$ parameters).
We then randomly sample the parameter space of the generation time distribution, and keep 1,000 parameter sets whose likelihood satisfies $\ln \mathcal{L}(\theta) > \ln \mathcal{L}(\hat{\theta}) - \lambda_2/2$, where $\lambda_2$ is the 95\% quantile of a $\chi^2$ distribution with $n=2$ degrees of freedom.
These parameter sets, as shown in Appendix 1-figure \ref{fig:gen-dist-params}, define the 95\% confidence interval for the generation time distribution (Figure \ref{fig:distributions}B).
We then use these sampled parameter sets to calculate quarantine efficacy, and the extrema of these efficacies across all of these parameter sets determines the 95\% confidence interval of the efficacy.
```{r gen-dist-params, fig.width = 5, fig.height = 3, fig.cap = caption}
caption <- "Log-likelihood values [$\\ln\\mathcal{L}(\\theta)-\\ln\\mathcal{L}(\\hat{\\theta})$] for random parameter samples of the generation time distribution.
These samples define the 95\\% confidence interval for the generation time distribution parameters.
Red dot is the maximum likelihood parameter combination as shown in Appendix 1-table \\ref{tab:distribution-parameters}."
ggplot(genParamsLLH, aes(x = shape, y = scale, colour = llh - max(llh))) +
geom_point(size = 0.8) +
geom_point(data = genParamsLLH[genParamsLLH$llh == max(genParamsLLH$llh),], colour = "red", size = 2) +
scale_colour_viridis_c() +
coord_cartesian(xlim = c(2,5), ylim = c(5,7)) +
labs(x = "shape parameter", y = "scale parameter", colour = expression(Delta*"log-likelihood")) +
plotTheme + theme(panel.grid.major.y = element_blank())
#' Save source data
write.csv(genParamsLLH, file = "data/Appendix 1-figure 1-source data 1.csv", row.names = F)
```
# Distribution parameters {-}
The parameters that define the incubation period distribution, generation time distribution, and infectivity profile, are shown in Appendix 1-table \ref{tab:distribution-parameters}.
```{r, distribution-properties, include = F}
#' Incubation period properties:
t <- seq(0,100,stepSize)
list(
# mean = sum(t * getIncubationPeriod(times = t, params = incParams) * stepSize),
mean = incParams$mean[[1]],
sd = sqrt(sum(t^2 * getIncubationPeriod(times = t, params = incParams) * stepSize) - incParams$mean[[1]]^2),
median = approx(x = getIncubationPeriod(times = t, params = incParams, CDF = T), y = t, xout = 0.5)$y
)
#' Generation time properties
t <- seq(0,100,stepSize)
list(
mean = genParams$mean,
sd = with(list(a = genParams$shape, b = genParams$scale),
sqrt(b^2 * (gamma(1 + 2/a) - (gamma(1 + 1/a))^2))),
median = approx(x = getGenDist(times = t, params = genParams, CDF = T), y = t, xout = 0.5)$y
)
#' Infectivity profile properties
t <- seq(-50,50,stepSize)
list(
mean = infParams$mean,
sd = sqrt(sum(t^2 * getInfectivityProfile(times = t, params = infParams) * stepSize) - infParams$mean^2),
median = approx(x = getInfectivityProfile(times = t, params = infParams, CDF = T), y = t, xout = 0.5)$y
)
```
\begin{table}[h]
\begin{tabular}{p{1.8cm} p{1.8cm} l r}
Distribution & Shape & Parameters & Properties \\ \hline
\multirow{7}{=}{Incubation period} & \multirow{7}{=}{Meta-log-normal} &
meanlog = 1.570, sdlog = 0.650 (Bi)& \\
& & meanlog = 1.621, sdlog = 0.418 (Lauer)& \\
& & meanlog = 1.434, sdlog = 0.661 (Li)& mean = 5.723\\
& & meanlog = 1.611, sdlog = 0.472 (Linton)& sd = 3.450 \\
& & meanlog = 1.857, sdlog = 0.547 (Ma)& median = 4.936 \\
& & meanlog = 1.540, sdlog = 0.470 (Zhang)& \\
& & meanlog = 1.530, sdlog = 0.464 (Jiang)& \\ \hline
\multirow{3}{=}{Generation time} & \multirow{3}{=}{Weibull} & & mean = 5.494 \\
& & shape = 3.277, scale = 6.127 & sd = 1.845 \\
& & & median = 5.479 \\ \hline
\multirow{3}{=}{Infectivity profile} & \multirow{3}{=}{Shifted Student's \emph{t}} & & mean = -0.042 \\
& & shift = -0.078, scale = 1.857, df = 3.345 & sd = 2.876 \\
& & & median = -0.078 \\ \hline
\end{tabular}
\caption{Parameters of the distributions used in this work.
The meta-log-normal distribution is the average of seven reported log-normal distributions \citep{ferretti:medRxiv:2020}.
The shifted Student's \emph{t} distribution for the infectivity profile is defined in \texttt{R} by \texttt{dt((x-shift)/scale, df)/scale} \citep{ferretti:medRxiv:2020}.}
\label{tab:distribution-parameters}
\end{table}
\clearpage
# Figure supplements {-}
\setcounter{figure}{0}
\renewcommand\thefigure{1-figure~supplement~\arabic{figure}}
```{r distributions, fig.height = 4.5, fig.cap = caption}
caption <- "A) The timeline of infection for an infector--infectee transmission pair.
The generation time is defined as the time interval between subsequent infections ($t_2-t_1$), while the serial interval is defined as the time between symptoms onsets in this transmission pair ($t_{S_2}-t_{S_1}$).
The incubation period describes the time between infection and symptom onset in a single individual (e.g. $t_{S_1}-t_1$).
B) The generation time distribution follows a Weibull distribution, and is inferred from the serial interval distribution \\citep{ferretti:medRxiv:2020}.
C) The infectivity profile (describing the time between symptom onset in the infector and infection of the infectee) follows a shifted Student's \\emph{t}-distribution, and is also inferred from the serial interval distribution \\citep{ferretti:medRxiv:2020}.
D) The distribution of incubation times follows a meta-distribution constructed from the average of seven reported log-normal distributions \\citep{bi:TheLancetInfectiousDiseases:2020,jiang:medRxiv:2020,lauer:AnnInternMed:2020,li:NEJM:2020,linton:J.Clin.Med.:2020,ma:medRxiv:2020,zhang:TheLancetInfectiousDiseases:2020}, as described in \\citet{ferretti:medRxiv:2020}."
#' Load the distributions
load("data/savedDistributions.RData")
#' Y-limits for consistent plotting
yLim <- c(0,0.35)
labY <- 0.212
#' Plot distributions together
distPlots <- plot_grid(
#' Generation time
ggplot(genDistMLE, aes(x = t, y = pdf)) +
geom_ribbon(aes(ymin = lower, ymax = upper), colour = "transparent", fill = "grey") +
geom_vline(xintercept = genParams$mean, linetype = "dashed", colour = "darkgrey") +
annotate("text", label = paste0(" mean = ", format(round(genParams$mean, 1), nsmall = 1), " days"), size = plotTheme$text$size/3, x = genParams$mean, y = labY, hjust = 0, vjust = 0) +
geom_line() +
coord_cartesian(xlim = c(0,15), ylim = yLim, expand = F) +
labs(x = "generation time", y = "probability density") +
ggtitle("generation time dist.") +
plotTheme + theme(plot.title = element_text(size = plotTheme$text$size, hjust = 0.5)),
#' Infectivity profile
ggplot(infProfMLE, aes(x = t, y = pdf)) +
geom_ribbon(aes(ymin = lower, ymax = upper), colour = "transparent", fill = "grey") +
geom_vline(xintercept = infParams$mean, linetype = "dashed", colour = "darkgrey") +
annotate("text", label = paste0(" mean = ", format(round(infParams$mean, 1), nsmall = 1), " days"), size = plotTheme$text$size/3, x = infParams$mean, y = labY, hjust = 0, vjust = 0) +
geom_line() +
coord_cartesian(xlim = c(-10,15), ylim = yLim, expand = F) +
labs(x = "days post symptoms", y = "probability density") +
ggtitle("infectivity profile") +
plotTheme + theme(plot.title = element_text(size = plotTheme$text$size, hjust = 0.5)),
#' Incubation period
ggplot(incDist, aes(x = t, y = pdf)) +
geom_vline(xintercept = incParams$mean, linetype = "dashed", colour = "darkgrey") +
annotate("text", label = paste0(" mean = ", format(round(incParams$mean, 1), nsmall = 1), " days"), size = plotTheme$text$size/3, x = incParams$mean, y = labY, hjust = 0, vjust = 0) +
geom_line() +
coord_cartesian(xlim = c(0,20), ylim = yLim, expand = F) +
labs(x = "incubation period (days)", y = "probability density") +
ggtitle("incubation period") +
plotTheme + theme(plot.title = element_text(size = plotTheme$text$size, hjust = 0.5)),
nrow = 1, align = "hv", axis = "tb",
#labels = "AUTO"
labels = c("B","C","D")
)
plot_grid(
ggdraw() + draw_image(magick::image_read_pdf("TTIQ-timeIntervals.pdf", density = 600), scale = 0.9),
distPlots, ncol = 1, labels = c("A",""), rel_heights = c(1,1.4)
)
#' Save source data
combined <- rbind(
cbind(genDistMLE, distribution = "generation"),
cbind(infProfMLE, distribution = "infectivity"),
cbind(incDist, lower = NA, upper = NA, distribution = "incubation")
)
write.csv(combined, file = "data/Figure 1-figure supplement 1-source data 1.csv", row.names = F)
```
\setcounter{figure}{0}
\renewcommand\thefigure{2-figure~supplement~\arabic{figure}}
```{r contacts-quarantine-params, include = F}
#' Parameters to scan over
tE <- 0
tQ.vals <- seq(0,4)
tR.vals <- seq(0,15)
tR.compare <- 10
#' Plot colours
colours <- scale_colour_viridis_d(
option = "inferno", begin = 0.1, end = 0.9, aesthetics = c("colour","fill"),
#name = expression("delay to quarantine" ~ (Delta[Q])),
name = expression(atop("delay between exposure and", "the start of quarantine" ~ (t[Q]-t[E]))),
labels = dayLabels(tQ.vals - tE),
guide = guide_legend(title.position = "left", title.hjust = 0.5, nrow = 1)
)
```
```{r contacts-quarantine-frac, dependson = c("contacts-quarantine-params"), include = F}
#' Function for computing fraction
getFraction <- function(params) {
lapply(tQ.vals, function(tQ) {
fraction <- ifelse(tR.vals > tQ,
getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params),
NA)
data.frame(
tE = tE,
tQ = tQ,
tR = tR.vals,
delay = factor(tQ - tE, levels = tQ.vals - tE),
release = tR.vals - tE,
fraction = fraction
)
}) %>% bind_rows()
}
#' Compute fraction for MLE parameters
frac <- getFraction(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Plot
fracPlot <- ggplot(frac, aes(x = release, y = fraction, colour = delay, fill = delay)) +
#geom_ribbon(aes(ymin = lower, ymax = upper), alpha = 0.4, colour = "transparent") +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(0,15), ylim = c(0,1)) +
labs(x = expression("days until release after exposure" ~ (t[R]-t[E])), y = "fraction of transmission\nprevented by quarantine") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
#' Extract and remove legend
legend <- get_legend(fracPlot)
fracPlot <- fracPlot + theme(legend.position = "none")
#' Save source data
write.csv(frac, file = "data/Figure 2-figure supplement 1-source data 1.csv", row.names = F)
```
```{r contacts-quarantine-utility, dependson = c("contacts-quarantine-params"), include = F}
#' Function to calculate relative utility
getRelUtility <- function(params) {
relUtility <- lapply(tQ.vals, function(tQ) {
fraction.compare <- getIntegral(upper = tR.compare, lower = tQ, tE = tE, params = params)
utility.compare <- getUtility(efficacy = fraction.compare, time = tR.compare - tQ)
fraction <- getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params)
utility <- ifelse(tR.vals > tQ,
getUtility(efficacy = fraction, time = tR.vals - tQ),
NA)
data.frame(
tE = tE,
tQ = tQ,
tR = tR.vals,
delay = factor(tQ - tE, levels = tQ.vals - tE),
release = tR.vals - tE,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}
#' Relative utility for MLE parameters
relUtility <- getRelUtility(params = genParams)
#' Evaluate relative utility for all parameter sets \theta
relUtilityLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relUtilityLLH, "relUtility")
#' Add CI to results
relUtility <- merge(relUtility, CI, sort = F)
#' Plot
utilityPlot <- ggplot(relUtility, aes(x = release, y = relUtility, colour = delay, fill = delay)) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
#geom_vline(xintercept = tR.compare - tE, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
coord_cartesian(xlim = c(0,15), ylim = c(0,2)) +
labs(x = expression("days until release after exposure" ~ (t[R]-t[E])), y = "\nrelative utility of quarantine") +
plotTheme + theme(legend.position = "None")
#' Save source data
write.csv(relUtility, file = "data/Figure 2-figure supplement 1-source data 2.csv", row.names = F)
```
```{r contacts-quarantine, dependson = c(-1,-2,-3), fig.cap = caption}
caption <- "Quantifying the effect of duration and delay for the standard quarantine protocol (no test) for traced contacts.
A) The fraction of transmission that is prevented by quarantining an infected contact [Eq. \\eqref{eq:F}].
We fix the time of exposure to $t_E=0$, and quarantine begins after a delay of 0-4 days (colour).
B) The relative utility of different quarantine durations compared to release on day 10 [Eq. \\eqref{eq:relUtility}].
Error bars reflect the uncertainty in the generation time distribution."
plot_grid(
plot_grid(fracPlot, utilityPlot, nrow = 1, axis = "tb", align = "hv", labels = "AUTO"),
legend, ncol = 1, rel_heights = c(1,0.2)
)
```
```{r contacts-reduced-params, include = F}
#' Parameters to scan over
tE <- 0
tQ <- 3
tR.vals <- seq(tQ,10)
DeltaT.vals <- seq(3,0)
s <- 0.1
r <- 0.5
tR.compare <- 10
tEnd <- 10
#' Plot colours
cols <- c("#000000",colorRampPalette(brewer.pal(9, "Oranges"))(12)[c(12,10,8,5)])
names(cols) <- c("noTest", DeltaT.vals)
colours <- scale_colour_manual(
name = expression("delay between test and release" ~ (t[R]-t[T])),
values = cols,
labels = c(noTest = "no test", dayLabels(DeltaT.vals)), aesthetics = c("colour","fill"),
guide = guide_legend(title.position = "left", title.hjust = 0.5, nrow = 1)
)
```
```{r contacts-reduced-frac, dependson = c("contacts-reduced-params"), include = F}
#' Fraction of transmission prevented without testing
getFraction.noTest <- function(params) {
fraction <- ifelse(tR.vals >= tQ,
getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params) + r * getIntegral(upper = tEnd, lower = tR.vals, tE = tE, params = params),
NA)
data.frame(
r = r,
tE = tE,
tQ = tQ,
tT = NA,
tR = tR.vals,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR.vals - tE,
fraction = fraction
)
}
#' Compute fraction for MLE parameters
frac <- getFraction.noTest(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Fraction of transmission prevented based on test-and-release
getFraction.test <- function(params, r = 0) {
#' Fraction of transmission before and after release
fraction.before <- getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params)
fraction.after <- getIntegral(upper = tEnd, lower = tR.vals, tE = tE, params = params)
#' Now loop over different testing delays
lapply(DeltaT.vals, function(DeltaT) {
tT.vals <- tR.vals - DeltaT
fraction <- ifelse(tT.vals >= tQ,
fraction.before +
(1 - falseNeg(tT.vals - tE) + r * falseNeg(tT.vals - tE)) * fraction.after,
NA)
data.frame(
r = r,
tE = tE,
tQ = tQ,
tT = tT.vals,
tR = tR.vals,
DeltaT = factor(DeltaT, levels = c("noTest", DeltaT.vals)),
release = tR.vals - tE,
fraction = fraction
)
}) %>% bind_rows()
}
fracTest <- getFraction.test(params = genParams, r = r)
#' Evaluate fraction from test-and-release for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction.test(params = genParamsLLH[id,], r = r)
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
fracTest <- merge(fracTest, CI, sort = F)
#' Combine standard and test-and-release
combined <- rbind(frac, fracTest)
#' Maximum preventable fraction of transmission
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE, params = genParams)
maxPreventableLLH <- sapply(seq_len(nrow(genParamsLLH)), function(id) {
getIntegral(upper = Inf, lower = tQ, tE = tE, params = genParamsLLH[id,])
})
maxPreventableCI <- c(min(maxPreventableLLH), max(maxPreventableLLH))
#' Plot
fracPlot <- ggplot(combined, aes(x = release, y = fraction, colour = DeltaT, fill = DeltaT)) +
geom_hline(yintercept = maxPreventable, colour = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line(linetype = "longdash") +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,1)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(3,10), ylim = c(0,1)) +
labs(x = expression("days until release after exposure" ~ (t[R]-t[E])), y = "fraction of transmission prevented \nby quarantine and reinforced hygiene") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
#' Extract and remove legend
legend <- get_legend(fracPlot)
fracPlot <- fracPlot + theme(legend.position = "None")
print(maxPreventable)
print(maxPreventableCI)
print(frac[frac$tR == 10, c("fraction","lower","upper")])
#' Save source data
write.csv(combined, file = "data/Figure 2-figure supplement 2-source data 1.csv", row.names = F)
```
```{r contacts-reduced-utility, dependson = c("contacts-reduced-params"), include = F}
#' Function to calculate relative utility of standard quarantine
getRelUtility.noTest <- function(params) {
fraction.compare <- getIntegral(upper = tR.compare, lower = tQ, tE = tE, params = params)
utility.compare <- getUtility(s = s, efficacy = fraction.compare, time = tR.compare - tQ)
fraction <- getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params) + r * getIntegral(upper = tEnd, lower = tR.vals, tE = tE, params = params)
utility <- ifelse(tR.vals > tQ,
getUtility(s = s, efficacy = fraction, time = tR.vals - tQ),
NA)
data.frame(
r = r,
tE = tE,
tQ = tQ,
tT = NA,
tR = tR.vals,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR.vals - tE,
relUtility = utility / utility.compare
)
}
#' Relative utility for MLE parameters
relUtility <- getRelUtility.noTest(params = genParams)
#' Evaluate relative utility for all parameter sets \theta
relUtilityLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relUtilityLLH, "relUtility")
#' Add CI to results
relUtility <- merge(relUtility, CI, sort = F)
#' Relative utility of test-and-release strategy compared to standard
getRelUtility.test <- function(params, r = 0) {
#' Standard quarantine strategy results
fraction.compare <- getIntegral(upper = tR.compare, lower = tQ, tE = tE, params = params)
utility.compare <- getUtility(s = s, efficacy = fraction.compare, time = tR.compare - tQ)
#' Test-and-release results
frac.before <- getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params)
frac.after <- getIntegral(upper = tEnd, lower = tR.vals, tE = tE, params = params)
time.before <- tR.vals - tQ
time.after <- tEnd - tR.vals
lapply(DeltaT.vals, function(DeltaT) {
tT.vals <- tR.vals - DeltaT
fraction <- ifelse(tT.vals >= tQ,
frac.before +
(1 - falseNeg(tT.vals - tE) + r * falseNeg(tT.vals - tE)) * frac.after,
NA)
utility <- getUtility(s = s, efficacy = fraction,
time = time.before + s * (1 - falseNeg(tT.vals - tE)) * time.after)
data.frame(
r = r,
tE = tE,
tQ = tQ,
tT = tT.vals,
tR = tR.vals,
DeltaT = factor(DeltaT, levels = c("noTest", DeltaT.vals)),
release = tR.vals - tE,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}
#' Relative utility using the MLE parameters
relUtilityTest <- getRelUtility.test(params = genParams, r = r)
#' Evaluate relative utility for all parameter sets \theta
relUtilityTestLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility.test(params = genParamsLLH[id,], r = r)
})
#' Compute CI
CI <- getCI(relUtilityTestLLH, "relUtility")
#' Add CI to results
relUtilityTest <- merge(relUtilityTest, CI, sort = F)
#' Combine standard and test-and-release
relUtilityCombined <- rbind(relUtility, relUtilityTest)
#' Plot
utilityPlot <- ggplot(relUtilityCombined, aes(x = release, y = relUtility, colour = DeltaT, fill = DeltaT)) +
#geom_vline(xintercept = tR.compare - tE, color = "darkgrey", size = thickLine) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line(linetype = "longdash") +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,1)) +
coord_cartesian(xlim = c(3,10), ylim = c(0,4.5)) +
labs(x = expression("days until release after exposure" ~ (t[{R^"-"}]-t[E])), y = "relative utility of quarantine\n and reinforced hygiene") +
plotTheme + theme(legend.position = "None")
print(relUtilityTest[relUtilityTest$tT == 3 & relUtilityTest$DeltaT == 0,c("relUtility","lower","upper")])
print(relUtilityTest[relUtilityTest$tT == 5 & relUtilityTest$DeltaT == 2,c("relUtility","lower","upper")])
print(relUtilityTest[relUtilityTest$tT == 6 & relUtilityTest$DeltaT == 0,c("relUtility","lower","upper")])
#' Save source data
write.csv(relUtilityCombined, file = "data/Figure 2-figure supplement 2-source data 2.csv", row.names = F)
```
```{r contacts-reduced, dependson = c(-1,-2,-3), fig.cap = caption}
caption <- "Quantifying the impact of quarantine and reinforced hygiene measures for traced contacts.
A) The fraction of transmission that is prevented by quarantining an infected contact and enforcing strict hygiene measures after release (see Appendix 1: Reinforced prevention measures after early release for details).
The scenarios are the same as in Figure \\ref{fig:contacts} (i.e. exposure at time $t_E=0$ and quarantine entry at time $t_Q=3$), but we reduce post-quarantine transmission by $r=50\\%$ until day 10, after which further transmission is unlikely.
The grey line represents the maximum attainable prevention by increasing release time, but keeping $t_Q=3$ fixed.
B) The relative utility of the quarantine and hygiene scenarios in panel A compared to the standard protocol 10 day quarantine [Eq. \\eqref{eq:relUtility}].
The grey line represents equal utilities (relative utility of 1).
We assume that the fraction of individuals in quarantine that are infected is $s=10\\%$, and that there are no false-positive test results.
Error bars reflect the uncertainty in the generation time distribution."
plot_grid(
plot_grid(fracPlot, utilityPlot, nrow = 1, axis = "tb", align = "hv", labels = "AUTO"),
legend, ncol = 1, rel_heights = c(1,0.2)
)
```
\setcounter{figure}{0}
\renewcommand\thefigure{3-figure~supplement~\arabic{figure}}
```{r contacts-asymptomatic-delay, fig.height = 4.5, fig.cap = caption}
caption <- "How the delay between symptom onset and self-isolation affects quarantine efficacy for traced contacts.
The y-axis is the fraction of total onward transmission per infected traced contact that is prevented by standard (no test) quarantine [Appendix 1-Eq. \\eqref{eq:Fasymptomatic}].
Each panel corresponds to a different delay between symptom onset and self-isolation in symptomatic individuals.
The left-most panel (zero delay) corresponds to Figure \\ref{fig:contacts-further}B.
The time of symptom onset is determined by the incubation period distribution (see Figure \\ref{fig:distributions}D).
As in Figures \\ref{fig:contacts} and \\ref{fig:contacts-further}, we fix the time of exposure at $t_E=0$ and the time of entering quarantine at $t_Q=3$ days.
Error bars reflect the uncertainty in the generation time distribution."
tE <- 0
tQ <- 3
tR.vals <- seq(0,15)
a.vals <- seq(0,1,0.25)
delay.vals <- seq(0,3)
#' Fraction of contacts that occur during quarantine for a given level of asymptomatic vs symptomatic
getFraction <- function(params) {
lapply(delay.vals, function(delay) {
fraction.asymptomatic <- ifelse(tR.vals >= tQ,
getIntegral(upper = tR.vals, lower = tQ, tE = tE, params = params),
NA)
fraction.symptomatic <- sapply(tR.vals, function(tR) {
if(tR < tQ) return(NA)
tS.vals <- seq(tE, tR - delay, stepSize)
symptoms.before <- stepSize * sum(
getIncubationPeriod(times = tS.vals - tE, params = incParams, CDF = F) *
getIntegral(upper = pmax(tS.vals + delay, tQ), lower = tQ, tE = tE, params = params)
)
symptoms.after <- (1 - getIncubationPeriod(times = tR - delay - tE, params = incParams, CDF = T)) *
getIntegral(upper = tR, lower = tQ, tE = tE, params = params)
return(symptoms.before + symptoms.after)
})
lapply(a.vals, function(a) {
data.frame(
tE = tE,
tQ = tQ,
tR = tR.vals,
delay = factor(dayLabels(delay), levels = dayLabels(delay.vals)),
release = tR.vals - tE,
a = factor(a, levels = a.vals),
fraction = a * fraction.asymptomatic + (1 - a) * fraction.symptomatic
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Compute fraction for MLE parameters
frac <- getFraction(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Plot
labs <- scales::percent(as.numeric(levels(frac$a)))
names(labs) <- levels(frac$a)
ggplot(frac, aes(x = release, y = fraction, colour = a, fill = a)) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
facet_grid(~ delay, labeller = label_both) +
scale_colour_viridis_d(option = "viridis", direction = -1, end = 0.9, aesthetics = c("colour","fill"),
name = expression("fraction asymptomatic" ~ (a)),
labels = labs,
guide = guide_legend(title.position = "left", title.hjust = 0.5, nrow = 1, byrow = T)) +
scale_x_continuous(breaks = seq(0,10)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(3,10), ylim = c(0,1)) +
labs(x = expression("days until release after exposure" ~ (t[R]-t[E])), y = "fraction of transmission\nprevented by quarantine") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
#' Save source data
write.csv(frac, file = "data/Figure 3-figure supplement 1-source data 1.csv", row.names = F)
```
\setcounter{figure}{0}
\renewcommand\thefigure{4-figure~supplement~\arabic{figure}}
```{r travellers-quarantine-params, include = F}
#' Parameters to sweep over
y.vals <- c(1,2,3,5,7,10,14)
tQ <- 0
tR.vals <- seq(0,15)
tR.compare <- 10
#' Consistent colours
colours <- scale_colour_viridis_d(
option = "inferno", begin = 0.1, end = 0.9, aesthetics = c("colour","fill"),
name = expression("duration of travel" ~ (y)), labels = dayLabels(y.vals),
guide = guide_legend(title.position = "left", title.hjust = 0.5, nrow = 1)
)
```
```{r travellers-quarantine-frac-total, dependson = c("travellers-quarantine-params"), include = F}
#' Fraction of total transmission prevented without testing
getFraction <- function(params) {
lapply(y.vals, function(y) {
tE.vals <- seq(-y,0)
lapply(tR.vals, function(tR) {
fraction <- mean(getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params))
data.frame(
y = factor(y, levels = y.vals),
tQ = tQ,
tR = tR,
release = tR - tQ,
fraction = fraction
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Evaluate for MLE parameters
frac <- getFraction(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Plot
fracTotalPlot <- ggplot(frac, aes(x = release, y = fraction, colour = y, fill = y)) +
#geom_hline(yintercept = 1, colour = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(0,10), ylim = c(0,1)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "fraction of total transmission\nprevented by quarantine") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
#' Extract and remove legend
legend <- get_legend(fracTotalPlot)
fracTotalPlot <- fracTotalPlot + theme(legend.position = "None")
#' Save source data
write.csv(frac, file = "data/Figure 4-figure supplement 1-source data 1.csv", row.names = F)
```
```{r travellers-quarantine-utility-total, dependson = c("travellers-quarantine-total-params"), include = F}
#' Relative utility compared to n=10 days
getRelUtility <- function(params) {
lapply(y.vals, function(y) {
tE.vals <- seq(-y,0)
fraction.compare <- mean(getIntegral(upper = tR.compare, lower = tQ, tE = tE.vals, params = params))
utility.compare <- getUtility(efficacy = fraction.compare, time = tR.compare - tQ)
lapply(tR.vals, function(tR) {
fraction <- mean(getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params))
utility <- ifelse(tR > tQ,
getUtility(efficacy = fraction, time = tR - tQ),
NA)
data.frame(
y = factor(y, levels = y.vals),
tQ = tQ,
tR = tR,
release = tR - tQ,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Evaluate for MLE parameters
relUtility <- getRelUtility(params = genParams)
#' Evaluate utility for all parameter sets \theta
relUtilityLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relUtilityLLH, "relUtility")
#' Add CI to results
relUtility <- merge(relUtility, CI, sort = F)
#' Plot
utilityTotalPlot <- ggplot(relUtility, aes(x = release, y = relUtility, colour = y, fill = y)) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
coord_cartesian(xlim = c(0,10), ylim = c(0,2.5)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "\nrelative utility of quarantine") +
plotTheme + theme(legend.position = "None")
#' Save source data
write.csv(relUtility, file = "data/Figure 4-figure supplement 1-source data 2.csv", row.names = F)
```
```{r travellers-quarantine-frac-local, dependson = c("travellers-quarantine-params"), include = F}
#' Fraction of local transmission prevented without testing
getFraction <- function(params) {
lapply(y.vals, function(y) {
tE.vals <- seq(-y,0)
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
lapply(tR.vals, function(tR) {
fraction <- mean(ifelse(maxPreventable < 1e-10, 1, getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
data.frame(
y = factor(y, levels = y.vals),
tQ = tQ,
tR = tR,
release = tR - tQ,
fraction = fraction
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Compute for MLE parameters
frac <- getFraction(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Plot
fracLocalPlot <- ggplot(frac, aes(x = release, y = fraction, colour = y, fill = y)) +
geom_hline(yintercept = 1, colour = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(0,10), ylim = c(0,1)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "fraction of local transmission\nprevented by quarantine") +
plotTheme + theme(legend.position = "none")
#' Save source data
write.csv(frac, file = "data/Figure 4-figure supplement 1-source data 3.csv", row.names = F)
```
```{r travellers-quarantine-utility-local, dependson = c("travellers-quarantine-params"), include = F}
#' Relative utility compared to n=10 days
getRelUtility <- function(params) {
lapply(y.vals, function(y) {
tE.vals <- seq(-y,0)
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
fraction.compare <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR.compare, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
utility.compare <- getUtility(efficacy = fraction.compare, time = tR.compare - tQ)
lapply(tR.vals, function(tR) {
fraction <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
utility <- ifelse(tR > tQ,
getUtility(efficacy = fraction, time = tR - tQ),
NA)
data.frame(
y = factor(y, levels = y.vals),
tQ = tQ,
tR = tR,
release = tR - tQ,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Calculate for MLE parameters
relUtility <- getRelUtility(params = genParams)
#' Evaluate utility for all parameter sets \theta
relUtilityLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relUtilityLLH, "relUtility")
#' Add CI to results
relUtility <- merge(relUtility, CI, sort = F)
#' Plot
utilityLocalPlot <- ggplot(relUtility, aes(x = release, y = relUtility, colour = y, fill = y)) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
coord_cartesian(xlim = c(0,10), ylim = c(0,2.5)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "\nrelative utility of quarantine") +
plotTheme + theme(legend.position = "None")
#' Save source data
write.csv(relUtility, file = "data/Figure 4-figure supplement 1-source data 4.csv", row.names = F)
```
```{r travellers-quarantine, dependson = c(-1,-2,-3,-4,-5), fig.height = 8, fig.cap = caption}
caption <- "Quantifying the effect of travel duration and quarantine duration for the standard quarantine protocol (no test) for returning travellers.
A) The fraction of \\emph{total} transmission that is prevented by quarantining an infected traveller [Eq. \\eqref{eq:FtotalAve}].
B) The relative utility of the different quarantine durations in panel A compared to release on day 10, based on the total fraction of transmission prevented.
C) The fraction of \\emph{local} transmission that is prevented by quarantining an infected traveller [Eq. \\eqref{eq:FlocalAve}].
D) The relative utility of the different quarantine durations in panel C compared to release on day 10, based on the local fraction of transmission prevented.
Colours represent the duration of travel $y$ and we assume infection can occur with equal probability on each day $t_E$ which satisfies $-y \\le t_E \\le 0$.
Quarantine begins at time $t_Q=0$, which is the time of arrival.
Error bars reflect the uncertainty in the generation time distribution."
plot_grid(
plot_grid(fracTotalPlot, utilityTotalPlot, fracLocalPlot, utilityLocalPlot,
nrow = 2, axis = "tb", align = "hv", labels = "AUTO"),
legend,
ncol = 1, rel_heights = c(1,0.1)
)
```
```{r travellers-reduced-params, include = F}
#' Parameters to scan over
y <- 7
tE.vals <- seq(-y,0)
tQ <- 0
tR.vals <- seq(0,10)
DeltaT.vals <- seq(3,0)
s <- 0.1
r <- 0.5
tR.compare <- 10
tEnd <- 10
#' Plot colours
cols <- c("#000000",colorRampPalette(brewer.pal(9, "Oranges"))(12)[c(12,10,8,5)])
names(cols) <- c("noTest", DeltaT.vals)
colours <- scale_colour_manual(
name = expression("delay between test and release" ~ (t[R]-t[T])),
values = cols,
labels = c(noTest = "no test", dayLabels(DeltaT.vals)), aesthetics = c("colour","fill"),
guide = guide_legend(title.position = "left", title.hjust = 0.5, nrow = 1)
)
```
```{r travellers-reduced-frac, dependson = c("travellers-reduced-params"), include = F}
#' Fraction of transmission prevented without testing
getFraction.noTest <- function(params) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
lapply(tR.vals, function(tR) {
fraction <- mean(ifelse(maxPreventable < 1e-10, 1, (getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params) + r * getIntegral(upper = tEnd, lower = tR, tE = tE.vals, params = params))/maxPreventable))
data.frame(
r = r,
y = y,
tQ = tQ,
tT = NA,
tR = tR,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR - tQ,
fraction = fraction
)
}) %>% bind_rows()
}
#' Compute fraction for MLE parameters
frac <- getFraction.noTest(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Fraction of transmission prevented based on test-and-release
getFraction.test <- function(params, r = 0) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
lapply(DeltaT.vals, function(DeltaT) {
lapply(tR.vals, function(tR) {
tT <- tR - DeltaT
fraction.before <- getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)
fraction.after <- getIntegral(upper = tEnd, lower = tR, tE = tE.vals, params = params)
fraction <- mean(ifelse(maxPreventable < 1e-10, 1, (fraction.before + (1 - falseNeg(tT-tE.vals) + r * falseNeg(tT-tE.vals)) * fraction.after)/maxPreventable))
data.frame(
r = r,
y = y,
tQ = tQ,
tT = tT,
tR = tR,
DeltaT = factor(DeltaT, levels = c("noTest", DeltaT.vals)),
release = tR - tQ,
fraction = fraction
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Compute fraction for MLE parameters
fracTest <- getFraction.test(params = genParams, r = r)
#' Evaluate fraction from test-and-release for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction.test(params = genParamsLLH[id,], r = r)
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
fracTest <- merge(fracTest, CI, sort = F)
#' Combine standard and test-and-release
combined <- rbind(frac, fracTest)
#' Maximum preventable fraction of transmission
maxPreventable <- mean(getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = genParams))
maxPreventableLLH <- sapply(seq_len(nrow(genParamsLLH)), function(id) {
mean(getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = genParamsLLH[id,]))
})
maxPreventableCI <- c(lower = min(maxPreventableLLH), upper = max(maxPreventableLLH))
#' Plot
fracPlot <- ggplot(combined, aes(x = release, y = fraction, colour = DeltaT, fill = DeltaT)) +
geom_hline(yintercept = 1, colour = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line(linetype = "longdash") +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(0,10), ylim = c(0,1)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "fraction of local trans. prevented \nby quarantine and reinforced hygiene") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
#' Extract and remove legend
legend <- get_legend(fracPlot)
fracPlot <- fracPlot + theme(legend.position = "None")
#' Save source data
write.csv(combined, file = "data/Figure 4-figure supplement 2-source data 1.csv", row.names = F)
```
```{r travellers-reduced-utility, dependson = c("travellers-reduced-params"), include = F}
#' Function to calculate relative utility of standard quarantine
getRelUtility.noTest <- function(params) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
fraction.compare <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR.compare, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
utility.compare <- getUtility(s = s, efficacy = fraction.compare, time = tR.compare - tQ)
lapply(tR.vals, function(tR) {
fraction <- mean(ifelse(maxPreventable < 1e-10, 1,
(getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params) + r * getIntegral(upper = tEnd, lower = tR, tE = tE.vals, params = params))/maxPreventable))
utility <- ifelse(tR > tQ,
getUtility(s = s, efficacy = fraction, time = tR - tQ),
NA)
data.frame(
r = r,
y = y,
tQ = tQ,
tT = NA,
tR = tR,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR - tQ,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}
#' Relative utility for MLE parameters
relUtility <- getRelUtility.noTest(params = genParams)
#' Evaluate relative utility for all parameter sets \theta
relUtilityLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relUtilityLLH, "relUtility")
#' Add CI to results
relUtility <- merge(relUtility, CI, sort = F)
#' Relative utility of test-and-release strategy compared to standard
getRelUtility.test <- function(params, r = 0) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
#' Standard quarantine strategy results
fraction.compare <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR.compare, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
utility.compare <- getUtility(s = s, efficacy = fraction.compare, time = tR.compare - tQ)
#' Test-and-release results
lapply(DeltaT.vals, function(DeltaT) {
lapply(tR.vals, function(tR) {
tT <- tR - DeltaT
frac.before <- getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)
frac.after <- getIntegral(upper = tEnd, lower = tR, tE = tE.vals, params = params)
fraction <- mean(ifelse(maxPreventable < 1e-10, 1,
(frac.before + (1 - falseNeg(tT-tE.vals) + r * falseNeg(tT-tE.vals)) * frac.after)/maxPreventable))
time.before <- tR - tQ
time.after <- tEnd - tR
time <- mean(time.before + s * (1 - falseNeg(tT-tE.vals)) * time.after)
utility <- ifelse(tT >= tQ,
getUtility(s = s, efficacy = fraction,
time = time),
NA)
data.frame(
r = r,
y = y,
tQ = tQ,
tT = tT,
tR = tR,
DeltaT = factor(DeltaT, levels = c("noTest",DeltaT.vals)),
release = tR - tQ,
relUtility = utility / utility.compare
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Relative utility using the MLE parameters
relUtilityTest <- getRelUtility.test(params = genParams, r = r)
#' Evaluate relative utility for all parameter sets \theta
relUtilityTestLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelUtility.test(params = genParamsLLH[id,], r = r)
})
#' Compute CI
CI <- getCI(relUtilityTestLLH, "relUtility")
#' Add CI to results
relUtilityTest <- merge(relUtilityTest, CI, sort = F)
#' Combine standard and test-and-release
relUtilityCombined <- rbind(relUtility, relUtilityTest)
#' Plot
utilityPlot <- ggplot(relUtilityCombined, aes(x = release, y = relUtility, colour = DeltaT, fill = DeltaT)) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line(linetype = "longdash") +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
coord_cartesian(xlim = c(0,10), ylim = c(0,4)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "relative utility of quarantine\n and reinforced hygiene") +
plotTheme + theme(legend.position = "None")
#' Save source data
write.csv(relUtilityCombined, file = "data/Figure 4-figure supplement 2-source data 2.csv", row.names = F)
```
```{r travellers-reduced, dependson = c(-1,-2,-3), fig.cap = caption}
caption <- "Quantifying the impact of quarantine and reinforced hygiene measures for returning travellers.
A) The fraction of \\emph{local} transmission that is prevented by quarantining an infected traveller and enforcing strict hygiene measures after release (see Appendix 1:Reinforced prevention measures after early release for details).
The scenarios are the same as in Figure \\ref{fig:travellers} (i.e. exposure occurs with equal probability between day minus seven and return at day zero, $-y \\le t_E \\le 0$, and quarantine starts at time $t_Q=0$), but we reduce post-quarantine transmission by $r=50\\%$ until day 10, after which further transmission is unlikely.
B) The relative utility of the quarantine and hygiene scenarios in panel A compared to the standard protocol 10 day quarantine [Eq. \\eqref{eq:relUtility}].
We assume that the fraction of individuals in quarantine that are infected is 10\\%, and that there are no false-positive test results.
Error bars reflect the uncertainty in the generation time distribution."
plot_grid(
plot_grid(fracPlot, utilityPlot,
nrow = 1, axis = "tb", align = "hv", labels = "AUTO"),
legend,
ncol = 1, rel_heights = c(1,0.2)
)
```
```{r travellers-further-adherence, include = F}
#' Parameters to consider
y <- 7
tE.vals <- seq(-y,0)
tQ <- 0
tR.vals <- seq(tQ,10)
DeltaT.vals <- seq(3,0)
tR.compare <- 10
tEnd <- 10
cols <- c("#000000",colorRampPalette(brewer.pal(9, "Oranges"))(12)[c(12,10,8,5)])
names(cols) <- c("noTest", DeltaT.vals)
colours <- scale_colour_manual(
name = expression("delay between test and release" ~ (t[R]-t[T])),
values = cols,
labels = c(noTest = "no test", dayLabels(DeltaT.vals)), aesthetics = c("colour","fill"),
guide = guide_legend(title.position = "top", title.hjust = 0.5,
nrow = 2, byrow = T)
)
#' Relative adherence to maintain quarantine efficacy without testing
getRelAdherence.noTest <- function(params) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
fraction.compare <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR.compare, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
lapply(tR.vals, function(tR) {
fraction <- mean(ifelse(maxPreventable < 1e-10, 1, getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
data.frame(
y = y,
tQ = tQ,
tT = NA,
tR = tR,
DeltaT = factor("noTest", levels = c("noTest", DeltaT.vals)),
release = tR - tQ,
relAdherence = ifelse(tR > tQ,
fraction.compare / fraction,
NA)
)
}) %>% bind_rows()
}
#' Evaluate for MLE parameters
relAdherence <- getRelAdherence.noTest(params = genParams)
#' Evaluate adherence for all parameter sets \theta
relAdherenceLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelAdherence.noTest(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(relAdherenceLLH, "relAdherence")
#' Add CI to results
relAdherence <- merge(relAdherence, CI, sort = F)
#' Fraction of transmission prevented based on test-and-release
getRelAdherence.test <- function(params, r = 0) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
fraction.compare <- mean(ifelse(maxPreventable < 1e-10, 1,
getIntegral(upper = tR.compare, lower = tQ, tE = tE.vals, params = params)/maxPreventable))
lapply(DeltaT.vals, function(DeltaT) {
lapply(tR.vals, function(tR) {
tT <- tR - DeltaT
fraction.before <- getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)
fraction.after <- getIntegral(upper = tEnd, lower = tR, tE = tE.vals, params = params)
fraction <- mean(ifelse(maxPreventable < 1e-10, 1, (fraction.before + (1 - falseNeg(tT-tE.vals) + r * falseNeg(tT - tE.vals)) * fraction.after)/maxPreventable))
data.frame(
y = y,
tQ = tQ,
tT = tT,
tR = tR,
DeltaT = factor(DeltaT, levels = c("noTest", DeltaT.vals)),
release = tR - tQ,
relAdherence = ifelse(tR >= tQ,
fraction.compare / fraction,
NA)
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Relative adherence using the MLE parameters
relAdherenceTest <- getRelAdherence.test(params = genParams, r = 0)
#' Evaluate relative utility for all parameter sets \theta
relAdherenceTestLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getRelAdherence.test(params = genParamsLLH[id,], r = 0)
})
#' Compute CI
CI <- getCI(relAdherenceTestLLH, "relAdherence")
#' Add CI to results
relAdherenceTest <- merge(relAdherenceTest, CI, sort = F)
#' Combine standard and test-and-release
relAdherenceCombined <- rbind(relAdherence, relAdherenceTest)
#' Plot
adherencePlot <- ggplot(relAdherenceCombined, aes(x = release, y = relAdherence, colour = DeltaT, fill = DeltaT)) +
geom_hline(yintercept = 1, color = "darkgrey", size = thickLine) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
colours +
scale_x_continuous(breaks = seq(0,10,2)) +
coord_cartesian(xlim = c(0,10), ylim = c(0,4.5)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "relative adherence required\nto maintain quarantine efficacy") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
print(relAdherence[relAdherence$tR == 5, c("tR","relAdherence","lower","upper")])
#' Save source data
write.csv(relAdherenceCombined, file = "data/Figure 4-figure supplement 3-source data 1.csv", row.names = F)
```
```{r travellers-further-asymptomatic, include = F}
#' Parameters to consider
y <- 7
tE.vals <- seq(-y,0)
tQ <- 0
tR.vals <- seq(0,10)
a.vals <- seq(0,1,0.25)
delay <- 0
#' Fraction of contacts that occur during quarantine for a given level of asymptomatic vs symptomatic
getFraction <- function(params) {
maxPreventable <- getIntegral(upper = Inf, lower = tQ, tE = tE.vals, params = params)
lapply(tR.vals, function(tR) {
fraction.asymptomatic <- getIntegral(upper = tR, lower = tQ, tE = tE.vals, params = params)
fraction.symptomatic <- sapply(tE.vals, function(tE) {
tS.vals <- seq(tE, tR - delay, stepSize)
symptoms.before <- stepSize * sum(
getIncubationPeriod(times = tS.vals - tE, params = incParams, CDF = F) *
getIntegral(upper = pmax(tS.vals + delay, tQ), lower = tQ, tE = tE, params = params)
)
symptoms.after <- (1 - getIncubationPeriod(times = tR - delay - tE, params = incParams, CDF = T)) *
getIntegral(upper = tR, lower = tQ, tE = tE, params = params)
return(symptoms.before + symptoms.after)
})
lapply(a.vals, function(a) {
fraction <- mean(ifelse(maxPreventable < 1e-10, 1,
(a * fraction.asymptomatic + (1 - a) * fraction.symptomatic)/maxPreventable))
data.frame(
y = y,
tQ = tQ,
tR = tR,
release = tR - tQ,
delay = delay,
a = factor(a, levels = a.vals),
fraction = fraction
)
}) %>% bind_rows()
}) %>% bind_rows()
}
#' Evaluate for MLE parameters
frac <- getFraction(params = genParams)
#' Evaluate fraction for all parameter sets \theta
fracLLH <- lapply(seq_len(nrow(genParamsLLH)), function(id) {
getFraction(params = genParamsLLH[id,])
})
#' Compute CI
CI <- getCI(fracLLH, "fraction")
#' Add CI to results
frac <- merge(frac, CI, sort = F)
#' Plot
labs <- scales::percent(as.numeric(levels(frac$a)))
names(labs) <- levels(frac$a)
asymptomaticPlot <- ggplot(frac, aes(x = release, y = fraction, colour = a, fill = a)) +
my_geom_errorbar() +
my_geom_line() +
my_geom_point() +
scale_colour_viridis_d(option = "viridis", direction = -1, end = 0.9, aesthetics = c("colour","fill"),
name = expression("fraction asymptomatic" ~ (a)),
labels = labs,
guide = guide_legend(title.position = "top", title.hjust = 0.5, nrow = 2, byrow = T)) +
scale_x_continuous(breaks = seq(0,10,2)) +
scale_y_continuous(labels = scales::percent) +
coord_cartesian(xlim = c(0,10), ylim = c(0,1)) +
labs(x = expression("day of release after arrival" ~ (t[R]-t[Q])), y = "fraction of local transmission\nprevented by quarantine") +
plotTheme + theme(legend.position = "bottom", legend.background = element_rect(colour = "black"), legend.key.height = unit(0.5,"cm"))
#' Save source data
write.csv(frac, file = "data/Figure 4-figure supplement 3-source data 2.csv", row.names = F)
```
```{r travellers-further, dependson = c(-1,-2), fig.height = 4.5, fig.cap = caption}
caption <- "How adherence and symptoms affect quarantine efficacy for returning travellers.
A) The fold-change in adherence to a new quarantine strategy that is required to maintain efficacy (local fraction of transmission prevented) of the baseline 10 day standard strategy.
Quarantine strategies are the same as in Figure \\ref{fig:travellers} (standard = black, test-and-release = colours).
The grey line represents equal adherence (relative adherence of 1).
B) The impact of symptomatic cases on the fraction of local transmission per infected traveller that is prevented by standard (no test) quarantine [Appendix 1-Eq. \\eqref{eq:Fasymptomatic}].
We assume that symptomatic individuals will immediately self-isolate at symptom onset.
The time of symptom onset is determined by the incubation period distribution (see Figure \\ref{fig:distributions}D).
The curve for $a=100\\%$ corresponds to the black curve in Figure \\ref{fig:travellers}A.
For both panels, as in Figure \\ref{fig:travellers}, we fix the trip duration to seven days and assume exposure can occur at any time $-y \\le t_E \\le 0$.
Quarantine begins at time $t_Q=0$.
Error bars reflect the uncertainty in the generation time distribution."
plot_grid(adherencePlot, asymptomaticPlot,
nrow = 1, axis = "tb", align = "hv", labels = "AUTO")
```
\clearpage
# References
Computing file changes ...
