A simple criterion to design optimal non-pharmaceutical interventions for mitigating epidemic outbreaks

For mitigating the COVID-19 pandemic, much emphasis is made on implementing non-pharmaceutical interventions to keep the reproduction number below one. However, using that objective ignores that some of these interventions, like bans of public events or lockdowns, must be transitory and as short as possible because of their significant economic and societal costs. Here, we derive a simple and mathematically rigorous criterion for designing optimal transitory non-pharmaceutical interventions for mitigating epidemic outbreaks. We find that reducing the reproduction number below one is sufficient but not necessary. Instead, our criterion prescribes the required reduction in the reproduction number according to the desired maximum of disease prevalence and the maximum decrease of disease transmission that the interventions can achieve. We study the implications of our theoretical results for designing non-pharmaceutical interventions in 16 cities and regions during the COVID-19 pandemic. In particular, we estimate the minimal reduction of each region’s contact rate necessary to control the epidemic optimally. Our results contribute to establishing a rigorous methodology to design optimal non-pharmaceutical intervention policies for mitigating epidemic outbreaks.


Introduction
Since the seminal work of May & Anderson [1], the design of interventions to eradicate infectious diseases has the objective of achieving a basic (R 0 ) or effective reproduction number below one [2,3]. The underlying assumption here is that it is possible to maintain interventions for long periods, such as long-term vaccination programmes. During the COVID-19 pandemic, this same objective is guiding the design of non-pharmaceutical interventions (NPIs) [4]. However, maintaining NPIs like bans of public events or lockdowns for long periods of time is infeasible because of their substantial economic and societal costs [5,6]. Actually, instead of aiming for eradication, NPIs aim to mitigate the economic and social costs of an epidemic outbreak [7]. Nevertheless, we still lack simple guidelines to design NPIs for mitigating epidemic outbreaks, analogous to the R 0 < 1 condition for eradication.
Here, we use the classic Susceptible-Infected-Removed (SIR) epidemiological model to fully characterize the design of NPIs for mitigating epidemic outbreaks. With this aim, we consider that NPIs should achieve an optimal tradeoff between two objectives [8]. First, optimal NPIs must minimize the period in which they need to be applied, consequently minimizing their associated economic and societal costs. Second, optimal NPIs must guarantee that the disease prevalence does not exceed a specified maximum level, which for example can represent health services' capacity for that particular disease outbreak [9]. We obtain a full analytical characterization of such optimal NPIs for mitigating epidemic outbreaks, specifying the optimal intervention at each state that the epidemic can be. This characterization yields the necessary and sufficient criterion for the existence of optimal NPIs for mitigation, analogous to the R 0 < 1 condition for eradication. We find that reducing the reproduction number below one is sufficient but not necessary for their existence. Instead, for mitigation, we show that the desired maximum disease prevalence determines the necessary reduction in the reproduction number. The consequence of not reducing the reproduction number below one is that interventions must start before the disease prevalence reaches the specified maximum level. We also demonstrate numerically that the derived optimal NPIs for mitigation are robust to uncertainties in the model parameters and unmodelled epidemic dynamics (e.g. undetected infections). Finally, we explore the implications of our theoretical result by analysing the response of 16 cities and regions across the globe to the COVID-19 pandemic, finding that most regions achieved a larger-thannecessary reduction in transmission. Our results contribute to designing NPIs to optimally and robustly mitigate epidemic outbreaks.

Characterizing optimal non-pharmaceutical interventions 2.1. Optimal epidemic mitigation using non-pharmaceutical interventions
Our objective is to characterize the reduction in the disease transmission that is optimal for each state in which the epidemic outbreak can be. For this, we leverage on the mathematical tractability of the SIR model [10], where the state can be characterized by the pair (S, I) Here, S is the proportion of the population that is susceptible to the disease, and I is the disease prevalence (i.e. the proportion of the population that is infected); see figure 1a. We discuss later other more detailed epidemic models. The epidemic state changes with time t as the disease is transmitted, producing the trajectory (S(t), I(t)) for t ≥ 0. For epidemic mitigation, we consider that the goal is keeping the disease prevalence below a specified level I max ∈ (0, 1]. A main factor determining this constant is health services' capacity in the sense that a prevalence above I max causes higher mortality due to hospital saturation [9]. In general, the selection of I max could depend on other social and economic factors of the specific population where the outbreak occurs. To keep I(t) ≤ I max , we assume we can apply one or several NPIs that reduce disease transmission by the factor (1 − u), for some u ∈ [0, 1]; see figure 1a. The NPIs achieve no reduction when u = 0, and they completely stop transmission when u = 1. Different NPIs correspond to particular values of u. For instance, a study of NPIs during the COVID-19 pandemic [11] found that closing most non-essential business corresponds to u ≈ 0.25, closing schools and universities to u ≈ 0.37, and limiting gatherings to at most 10 people to u ≈ 0.42. In practice, it can be unfeasible to fully stop the disease transmission (i.e. u < 1) because of inherent challenges like asymptomatic transmission [12], or because of the necessity of maintaining a working economy [13]. Therefore, we upper-bound the reduction by u max ∈ (0, 1). We say that u is admissible if u ∈ [0, u max ]. Different admissible NPIs can keep the disease prevalence below I max . For instance, 'intervention 1' in the example of figure 1b,c keeps this restriction and it has an 'effective duration' of 120 days. Here, the effective duration of an intervention is the time interval between the start of the outbreak and the last time that a non-zero intervention is applied (figure 1d). 'Intervention 2' of figure 1b,c also keeps the restriction For the optimal NPI design problem, the objective is to design the intervention u*(t) with minimal effective duration such that u*(t) ∈ [0, u max ] and I(t) ≤ I max for all t ≥ 0. (b,c) The response of the SIR model for two interventions (parameters are β = 0.52, γ = 1/7, I 0 = 8.855 × 10 −7 and S 0 = 1 − I 0 ). Both interventions 1 and 2 satisfy u(t) ≤ u max and guarantee that I(t) ≤ I max . Actually, intervention 2 is the optimal one derived using our analysis: it is the intervention with minimal effective duration satisfying I(t) ≤ I max . (d) The effective duration of an intervention measures the interval between the start of the outbreak and the last time that a non-zero intervention is applied. In this example, the effective duration of intervention 1 is 120 days, while the effective duration of intervention 2 is 69 days.
royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20200803 I(t) ≤ I max , but its effective duration is only 69 days. To design the optimal NPI for mitigating epidemic outbreaks, we ask for the intervention with minimal effective duration. Specifically, we ask for the admissible intervention u*(S(t), I(t)) required now (i.e. at the current state) such that: (1) it minimizes the effective duration of the intervention and (2) it ensures that the prevalence can be maintained below I max for all future time by using some admissible intervention. If the optimal NPI problem has a solution u*, then u*(S, I) characterizes the optimal reduction in the disease transmission that the NPIs should achieve if the epidemic state is (S, I ). In particular, u* gives the optimal way to start and stop the NPIs.

Optimal non-pharmaceutical interventions for mitigation exist without reducing the reproduction number below one
Our first main result is a complete analytical characterization of the optimal NPIs for outbreak mitigation in the SIR model (see box 1 for a summary and electronic supplementary material, note S1, for details). To understand how these optimal NPIs work, note that the SIR model predicts a safe zone of states (S, I) where, without any further interventions, the disease prevalence will not exceed I max (blue zone in figure 2a-c). The safe zone is characterized by the inequality I F R0 (S), where R 0 is the basic reproduction number of the outbreak in the population, and the function F R is defined in equation (2.1) of box 1. The goal of the optimal NPIs is thus to reach this safe zone as fast as possible without violating the restriction I(t) ≤ I max . The ability to achieve this goal depends on the epidemic state. That is, we can partition the plane (S, I ) in two regions: those states from which it is possible to reach the safe zone without exceeding I max ( feasible states), and those where it is impossible (unfeasible states). We find these two regions are characterized by the separating curve F Rc (S), where we call R c := (1 − u max )R 0 the controlled reproduction number (figure 2a-c). Note that R c describes the maximum reduction in the basic reproduction number that (constant) admissible interventions can achieve. Therefore, R c < 1 is the necessary and sufficient condition that a constant and permanent admissible intervention (i.e. u(t) ≡ const. for all t ≥ 0) needs to satisfy to eradicate a disease outbreak in the SIR model. However, for outbreak mitigation, our analysis shows that feasible states exists without achieving disease eradication (white regions in figure 2b,c). This result is important because it Box 1. Optimal NPIs for the Susceptible-Infected-Removed (SIR) model.
The SIR model with interventions u(t) ∈ [0, u max ] reducing disease transmission takes the form Here, S(t) and I(t) are the proportion of the population that is susceptible or infected at time t ≥ 0, respectively. We denote by (S 0 , I 0 ) the initial state at t = 0. The parameters of the SIR model are the (effective) contact rate β ≥ 0, and the mean residence time of infected individuals γ ≥ 0 (in units of day −1 ). By assuming S 0 ≈ 1, these two parameters yield the basic reproduction number R 0 = β/γ. We are interested in reaching the safe zone The safe zone is the largest set with the following property: if, for any given time t 1 , the state (S 1 , I 1 ) belongs to S, we can set u = 0 henceforth and still have I(t) ≤ I max for all t ≥ t 1 . That is, when S is reached, we can terminate the intervention with the assurance that a possible rebound in the disease prevalence will not exceed I max .
Our goal is to steer an arbitrary initial state (S 0 , I 0 ) to the safe zone S in minimal time without violating the constraint I(t) ≤ I max . We say that an intervention achieving this goal is an optimal intervention.
In electronic supplementary material, note S1, we prove that the existence of an optimal intervention is characterized by the separating curve F Rc as follows: 1. An optimal intervention exists if and only if the initial state (S 0 , I 0 ) lies below this separating curve (i.e. I 0 F Rc (S 0 )).
Above, R c : = (1 − u max )R 0 is the controlled reproduction number. Moreover: 2. If it exists, the optimal intervention u* at the state (S, I) is Above, the curve S ¼ C(I) is defined in electronic supplementary material, note S1, while S* denotes the intersection of S ¼ C(I) and I ¼ F Rc (S).
Code to calculate the optimal interventions is provided as electronic supplementary material.
royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20200803 proves that optimal NPIs for epidemic mitigation do not require reducing the basic reproduction number below one (i.e. without achieving R c < 1).

A design criterion for optimal non-pharmaceutical interventions
We demonstrated above that optimal NPIs for outbreak mitigation exist even when R c > 1. However, how large can R c be before NPIs keeping I(t) ≤ I max do not exist? When S(0) → 1, our characterization shows that such NPIs exists if and only if The above inequality is our second main result, connecting the specified maximum disease prevalence I max with the outbreak's controlled reproduction number R c = (1 − u max )R 0 (electronic supplementary material, note S2). The inequality (2.3) governs the existence of NPIs for mitigating epidemic outbreaks, in analogy to how the condition R c < 1 works for disease eradication. Note that R c < 1 is a sufficient condition for the existence for NPIs, but the inequality (2.3) shows that this condition is far from necessary. If I max > 0, there exists R c > 1 for which NPIs exist (figure 2d). Note also that the maximum feasible R c increases with I max .
We can use (2.3) to design NPIs for outbreak mitigation as follows. Consider an infectious disease outbreak with a given R 0 and that the specified maximum prevalence is I max . Then, the inequality (2.3) gives the criterion to design NPIs by providing the range of disease transmission reduction u max that the NPIs should attain. In particular, it provides the minimal reduction u Ã max in the disease transmission required for the existence of NPIs.   royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20200803 is the maximum admissible controlled reproduction number (orange point in figure 2d ). Therefore, if an outbreak in the population has R 0 = 3, then the minimal required reduction

Optimal non-pharmaceutical interventions for outbreak mitigation are simple
For any epidemic state, the optimal transmission reduction takes a simple form which can be described by colouring the (S, I) plane; see top row of figure 3. Here, for all states in the white region the optimal intervention is no intervention; for all states in the yellow region the optimal intervention is u*(S, I ) = u max . There are regions (specifically lines) where the optimal intervention switches frequently between u* = 0 and u* = u max producing a so-called 'singular arc' that slides along the two regions, leading to an 'average' intervention u* ∈ [0, u max ]. In general, we find that the optimal NPIs have four phases: a first one where no intervention is needed, a second phase where interventions start with maximum strength, a third phase of gradual decrease of interventions, and a 'final push' where the maximum interventions are re-applied for a short period to reach the safe zone faster. We illustrate the above behaviour in three qualitatively different cases. The first case is when the optimal intervention starts just when the disease prevalence reaches I max (figure 3a). This case occurs when the interventions are strong enough to stop the rise in prevalence at I max regardless of the remaining fraction the population that is still susceptible to the disease. Our analysis shows that this case occurs if and only if u max is large enough to render R c = (1 − u max )R 0 ≤ 1. When the initial susceptible population is close to 1 (pink trajectory in  In this case, the optimal intervention starts when the disease prevalence reaches I max . Afterwards, the intervention decreases in a hyperbolic arc until reaching the point S = S*. At that time, the intervention becomes maximum in the 'final push' to reach the safe zone. (b) For u max = 0.58, the controlled reproduction number is R c = (1 − u max )R 0 = 1.52 > 1. Here, F R c (1) . 0, implying that the epidemic still can be mitigated for initial states with S 0 ≈ 1 and I 0 ≈ 0 ( pink trajectory). In this case, the optimal intervention starts when the initial condition hits the separating curve below I max at t = 35. At that instant, the intervention starts with the maximum value u max , and continues in that form until the trajectory reaches I max . (c) Choosing u max = 0.4 yields R c = 2.184 > 1. In this case, the optimal intervention problem does not have a solution for all initial states S 0 > 0.85. This is illustrated by pink trajectory: even when applying the maximum intervention from the start, I(t) will grow beyond I max .
royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20200803 figure 3a), the optimal intervention first waits until the disease prevalence reaches I max . At that time, the optimal NPI stops the disease prevalence exactly at I max , and then it gradually decreases its magnitude to ensure that the disease prevalence slides along I max as the susceptible population decreases. When the susceptible population reaches the threshold S*, the optimal intervention is again the maximum one ( figure 3a). This 'final push' allows reaching the safe zone faster, releasing the interventions sooner. The middle and bottom panels of figure 3a show the resulting disease prevalence and optimal interventions as a function of time.
Note that a smaller initial susceptible population yields other trajectories (green and purple in figure 3a). The second case is when an 'early' intervention is necessary before the disease prevalence reaches I max (figure 3b). This case happens when the admissible reduction in the contact rate cannot immediately stop the disease prevalence at I max if the susceptible population is large at that time. We find this case occurs if and only if u max is small in the sense that R c = (1 − u max )R 0 > 1. Here, a trajectory may hit the yellow region before reaching I max (pink trajectory in figure 3b). When that happens, the optimal intervention starts with the maximum reduction u* = u max . Then it maintains this maximum reduction to 'slide' the trajectory between the yellow and white regions. Once the trajectory reaches I max , the magnitude of the optimal intervention decreases to slide the trajectory along I max . Again, the final push occurs when the susceptible population reaches the point S*.
The third case is when the initial state (S 0 , I 0 ) lies in the unfeasible region (figure 3c). This case occurs when u max is so small that, even if the maximum admissible intervention u = u max is applied from the start of the outbreak, the disease prevalence will exceed I max (pink trajectory in figure 3c). In this case, the optimal intervention problem is unfeasible because it is impossible to achieve I(t) ≤ I max . However, note that using u* = u max yields the smallest prevalence peak. Other trajectories that start with a smaller proportion of susceptible individuals remain feasible (green and purple in figure 3c). In particular, note that the threshold S* decreases as R c increases.

Optimal non-pharmaceutical interventions for outbreak mitigation are robust
To evaluate the optimal NPIs for outbreak mitigation in more realistic scenarios, we numerically analysed their performance in three epidemic models with uncertain epidemic parameters and more detailed epidemic dynamics (see details in electronic supplementary material, note S3). In all cases, we consider that the basic reproduction number has been estimated asR 0 using an SIR model, and that the optimal NPIs for mitigation are designed using this estimate. Then, these optimal NPIs are applied to an outbreak with possibly different epidemic dynamics and possibly different R 0 . Note that estimation errors in R 0 will affect the correct start of the NPIs and the 'final push' for reaching the safe zone.
In the first scenario, we consider an outbreak with SIR dynamics where the strength of the NPIs is uncertain. We model this uncertainty replacing u by ku in the model equations, where k ∈ (0, 1). Then, for example, k = 0.9 (resp. k = 1.1) represents a 10% underestimation (resp. overestimation) of the NPI strength. Across outbreaks with different R 0 's and an uncertainty of 10% in the intervention's strength, we find that the disease prevalence is maintained below I max as long as R 0 is not underestimated (figure 4a). In the second scenario, we consider an SEIR outbreak with an incubation period for the disease. For an incubation period of 7 days, as typical for a COVID-19 infection, the optimal NPIs maintain the disease prevalence below I max if R 0 < 2.5 and its value is estimated with an error of below 30% (solid yellow and orange in figure 4b). For larger R 0 or a larger incubation period, the disease prevalence may exceed I max (red in figure 4b).
For the final scenario, we consider an SEIIR model with an incubation period of 7 days and with a fraction p ∈ [0, 1] of infected individuals that are asymptomatic and thus remain hidden to the epidemic surveillance system. The goal is to maintain the prevalence of symptomatic individuals below I max , without knowing the fraction of asymptomatic individuals. This situation occurs during the COVID-19 pandemic, where between p = 0.55 and p = 0.8 of infections are asymptomatic [14]. For p < 0.7 and R 0 < 3.64, the optimal NPIs maintain the disease prevalence of symptomatic individuals below or very close to I max if the estimation error for R 0 is below 30% (dotted and solid lines in figure 4c). An outbreak with low R 0 produces a maximum disease prevalence of symptomatic individuals below I max , which may result in a larger effective duration of the interventions. Overall, these numerical results show that the optimal NPIs are robust against a wide range of parameter uncertainty and unmodelled dynamics, provided that the estimation error in the outbreak's basic reproduction number does not exceed 30%.

Designing optimal non-pharmaceutical interventions for mitigating the COVID-19 pandemic
To explore the implications of our simple criterion for designing NPIs for outbreak mitigation, we analysed how 16 cities and regions implemented NPIs during the COVID-19 pandemic. For each region or city, we constructed I max using the number of available intensive care beds during the first months of the pandemic, considering that a fraction of the infected individuals will require them (electronic supplementary material, note S4). The values for I max that we obtain range from 2.87 × 10 −3 for Lima (Peru) to 109.78 × 10 −3 for Boston (USA), reflecting the large heterogeneity of the available health services across the globe (figure 5a). With this information, we calculated the maximum feasible R Ã c for each region using our design criterion of inequality (2.3). Since R Ã c is a monotone function of I max , we find that R Ã c follows the same trend as I max (figure 5b). The smallest R Ã c ¼ 1:08 occurs for Lima and the largest R Ã c ¼ 1:75 for Boston. Note that in both cases R Ã c . 1. This result implies that, for the R 0 of a region's disease outbreak, a successful mitigation of the outbreak requires NPI policies that achieve at least a reduction u Ã max such that (1 À u Ã max )R 0 R Ã c . Next, we investigated the minimal reduction u Ã max in transmission required to achieve those upper bounds for the COVID-19 pandemic. For this, we first collected information for the R 0 in each region calculated at the start of the pandemic and when the NPIs were inactive (electronic supplementary material, note S3). We find a median nominal R 0 of 2.2, with Tokyo having the smallest one (R 0 = 1.3) and Madrid having the largest one (R 0 = 3.11); see figure 5c. From these values of R 0 , we calculated the minimal required reduction u Ã max per region or city (blue in figure 5d).   Figure 4. Optimal non-pharmaceutical interventions are robust. For all panels, the estimated parameters used for constructing the optimal NPIs areĝ ¼ 1=7, b ¼ 0:52, I max = 0.1, u max = 0.6. We consider a population of N = 8.855 × 10 6 as in Mexico City, and the initial conditions I(0) = 1/N and S(0) = 1 − 1/N. If the models contain other state variables, they were initialized at zero. The optimal NPIs are constructed assumingR 0 ¼b =ĝ, while the actual epidemic dynamics has a possibly different R 0 = β/γ. Panels show results for outbreaks with three values of R 0 : low (yellow), medium (orange) and large (red). (a) SIR model where the reduction in the disease transmission by the NPIs is uncertain. We model this case replacing u by ku in the model equations. Panel shows the results for k = 1.1 (dotted), k = 1 (solid) and k = 0.9 (dashed). (b) SEIR model where exposed individuals do not transmit the infection, with λ > 0 the incubation period. Panel shows the results for λ = 1/5 (dotted), λ = 1/7 (solid) and λ = 1/11 (dashed). (c) SEIIR model with λ = 1/7 and two classes of infected individuals (symptomatic and asymptomatic). Here, p ∈ [0, 1] is the proportion of exposed individuals that become asymptomatic. The vertical axis denotes the disease prevalence for symptomatic individuals. The panel shows the results for p = 0.55 (dotted), p = 0.7 (solid) and p = 0.8 (dashed). royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20200803 the city with large I max ends up requiring a smaller u Ã max (e.g. Boston with u Ã max ¼ 0:26 and Lima with u Ã max ¼ 0:50). To evaluate the feasibility of achieving the minimal reduction predicted by our analysis, we collected data for the average mobility reduction in each region during the NPIs (grey in figure 5d; electronic supplementary material, note S4). Considering this average mobility reduction as a proxy for the reduction in disease transmission, we find that all regions achieved a greater than necessary reduction. For example, Delhi attained a mobility reduction of 0.84, while the minimal necessary reduction in transmission according to our analysis is u Ã max ¼ 0:42. Other regions are in the boundary. For example, New South Wales attained a mobility reduction of 0.48, while the minimal necessary reduction in transmission was u Ã max ¼ 0:44. Overall, across regions, we find a median excess of 0.22 in the reduction of mobility compared to the minimal reduction in transmission u Ã max predicted by our analysis.

Discussion
Our results provide a complete analytical characterization of the optimal NPIs for mitigating epidemic outbreaks in the SIR model. We also show that these optimal NPIs are robust as they can 'work' in epidemic models with more complicated dynamics. The SIR model is a minimal strategic model of the general population dynamics of a disease. Although this model ignores critical epidemiological phenomena, using the SIR model allows us to leverage on its mathematical tractability to obtain a complete characterization of the optimal NPIs for outbreak mitigation. The feedback form u*(S, I) of the optimal intervention reflects such complete characterization, prescribing the optimal action to perform if the epidemic is at any state (S, I). This feedback strategy should be contrasted to most other studies applying optimal control to epidemic outbreaks, where the derived optimal intervention u*(t, S 0 , I 0 ) is an open-loop function of time [15][16][17][18]. The open-loop intervention gives the optimal action at any time given a particular initial state (S 0 , I 0 ). However, it does not tell us the optimal action if the epidemic is not in the exact state predicted by the model. Understanding the optimal action to perform at any state has the crucial advantage of allowing us to apply this knowledge to any model, and to reality. Feedback can give control strategies the required robustness to work on real systems despite large uncertainties and unknown dynamics [19,20], and we numerically confirmed that our optimal NPIs for outbreak mitigation have such robustness. Indeed, other works have also found that interventions derived from the SIR model can work in detailed agent-based models of epidemic outbreaks [21]. Future work could analyse the robustness of optimal interventions when the epidemic state is not entirely known. This situation may happen when significant delays exist in reporting new infections, or when tests for identifying infected individuals are limited. For example, control-theoretical techniques like the construction of observers and predictors allow applying our optimal interventions when the only available information is the disease prevalence, and when this information is obtained with a significant delay [22].
Our framework could also guide the complete characterization of optimal NPIs for mitigating epidemic outbreaks using more detailed models or more detailed optimization objectives, but this is likely very challenging. Indeed, deriving such complete characterization for very detailed models can be unreasonable, considering the tradeoff between how detailed is a model and how much we can trust its predictions [23]. Note also that our approach could be applied to calculate the optimal NPIs in the presence of a constant vaccination rate by modifying the SIR model accordingly (e.g. [24]). royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20200803