Analytical estimation of maximum fraction of infected individuals with one-shot non-pharmaceutical intervention in a hybrid epidemic model

NF\fnmNaoya Fujiwara TO\fnmTomokatsu Onaga TW\fnmTakayuki Wada ST\fnmShouhei Takeuchi JS\fnmJunji Seto TN\fnmTomoki Nakaya KA\fnmKazuyuki Aihara \orgnameGraduate School of Information Sciences, Tohoku University, \street6-3-09 Aoba, Aramaki-aza Aoba-ku, \postcode980-8579 \citySendai, Miyagi, \cnyJapan \orgnameInstitute of Industrial Science, The University of Tokyo, \street4-6-1 Komaba, \postcode153-8505 \cityMeguro-ku, Tokyo, \cnyJapan \orgnameCenter for Spatial Information Science, The University of Tokyo, \street5-1-5 Kashiwanoha, \postcode277-8508 \cityKashiwa, Chiba, \cnyJapan \orgnameFrontier Research Institute for Interdisciplinary Sciences, Tohoku University, \streetAramaki aza Aoba 6-3, Aoba-ku, \postcode980-8578 \citySendai, Miyagi, \cnyJapan \orgnameGraduate School of Human Life Science Human Life Science Course, Osaka City University, \street3-3-138, Sugimoto, Sumiyoshi-ku, \postcode558-8585 \cityOsaka, Osaka, \cnyJapan \orgnameFaculty of Nursing and Nutrition, University of Nagasaki, \street1-1-1 Manabino, Nagayo-cho, Nishi-Sonogi-gun, \postcode851-2195 \cityNagasaki, \cnyJapan \orgnameDepartment of Microbiology, Yamagata Prefectural Institute of Public Health, \street1-6-6 Toka-machi, \postcode990-0031 \cityYamagata, Yamagata, \cnyJapan \orgnameGraduate School of Environmental Studies, Tohoku University, \streetAoba, 468-1, Aramaki, Aoba-ku, \postcode980-8572 \citySendai, Miyagi, \cnyJapan \orgnameInternational Research Center for Neurointelligence, The University of Tokyo, \street7-3-1 Hongo, \postcode113-0033 \cityBunkyo-ku, Tokyo, \cnyJapan

Abstract

\parttitle

Background Facing a global epidemic of new infectious diseases such as COVID-19, non-pharmaceutical interventions (NPIs), which reduce transmission rates without medical actions, are being implemented around the world to mitigate spreads. One of the problems in assessing the effects of NPIs is that different NPIs have been implemented at different times based on the situation of each country; therefore, few assumptions can be shared about how the introduction of policies affects the patient population. Mathematical models can contribute to further understanding these phenomena by obtaining analytical solutions as well as numerical simulations.

\parttitle

Methods and Results In this study, an NPI was introduced into the SIR model for a conceptual study of infectious diseases under the condition that the transmission rate was reduced to a fixed value only once within a finite time duration, and its effect was analyzed numerically and theoretically. It was analytically shown that the maximum fraction of infected individuals and the final size could be larger if the intervention starts too early. The analytical results also suggested that more individuals may be infected at the peak of the second wave with a stronger intervention.

\parttitle

Conclusions This study provides quantitative relationship between the strength of a one-shot intervention and the reduction in the number of patients with no approximation. This suggests the importance of the strength and time of NPIs, although detailed studies are necessary for the implementation of NPIs in complicated real-world environments as the model used in this study is based on various simplifications.

infectious diseases,

pandemics,

non-pharmaceutical interventions,

hybrid dynamical systems,

keywords:

\startlocaldefs\endlocaldefs

{fmbox}\dochead

Research

{artnotes}

{abstractbox}

Introduction

Because of the global spread of COVID-19, our human society is facing a major public health crisis. The COVID-19 pandemic is caused by an emerging pathogen, SARS-CoV-2, for which there is no immunized population, causing an overshooting increase in the number of infected patients and depleting medical resources in many countries. Medical institutions are facing a difficult situation in which they must control second transmissions while treating critically ill patients, and as the number of patients increases, the medical system becomes swiftly tighter. When the number of patients exceeds the capacity, the quality of medical care deteriorates drastically, and the number of medical devices required for life support reaches its limit. This situation further increases the fatality rate of this infectious disease and causes serious damage to our society.

No effective treatment for COVID-19 has been established yet as of September, 2020, and only a public health approach can function as a control measure for the epidemic. To mitigate the spread of COVID-19, each country is implementing non-pharmaceutical interventions (NPIs) [1] to regulate social activities. NPIs comprise policies such as case isolation, voluntary home quarantine, closure of schools and universities, social distancing, stopping mass gatherings, and border closure. In major European countries, these NPIs were implemented, depending on the epidemic situation, during the first part of spring in 2020, with different timings and intensities [2, 3]. Although the first phase of the epidemic appeared to be suppressed by these mitigation measures, the re-epidemic became clearer in almost all countries because of deregulation after the first epidemic.

Because the control of epidemics by NPIs has caused a situation involving the imposition of strong restrictions on human socioeconomic activities, it would be desirable to study in advance the optimal timing, intensity, and duration of interventions that could bring about more promising results with minimal damage to society. Regarding COVID-19, for the purpose of ex-post verification, the influence of NPIs implemented in each country on the effective reproduction number was estimated. In European countries, the correlation between the decay of the effective reproduction number and the implementation of various NPIs has been verified [3]. The impact of travel limitations in China on the spread of infection has been discussed [4, 5]. However, a reliable estimation of the effects of NPIs is difficult, because the differences in NPI strategies employed by each country are strongly related to various background factors, such as the epidemic situations, social structure, legal systems, and culture [2].

Mathematical modelling is an important method for estimating the effect of NPIs. In particular, the global pandemic of COVID-19 revealed that a situation in which only NPIs are effective against an emerging infectious disease is possible in societies in the 2020s. Such a situation had been “neglected” as a practical research target, which enhances the importance of theoretical approaches. Recently, compartmental models such as the SIR model [6, 7, 8, 9] have been extended to estimate the effects of NPIs on the number of patients [10, 11, 12, 13]. In other settings, optimal policies where intervention intensity can change continuously over time have been discussed in the context of minimizing objective functions [14, 15, 16, 17]. These solutions are very useful if we can estimate the effect of NPIs on the change in transmission parameters precisely.

The impacts of NPIs should be evaluated as a discontinuous change of the transmission rate in a model to represent the temporal discontinuity of intervention in the real world. The dynamics of the number of infected patients in continuous time in a system that includes discrete parameters, state spaces, and continuous-time dynamics can be modeled as a hybrid dynamical system [18, 19, 20, 21, 22]. This framework has been applied to mathematical models of infectious diseases [23, 24, 25]. In the simplest case, the one-shot intervention model, in which the intervention is implemented only once during the epidemic, can be used to discuss the theoretical dependence of the effects of NPIs on the timing and intensity [26]. As the COVID-19 epidemic continues, the accumulation of theoretical research on the effects brought about by NPIs has become even more significant. Recently, compartmental models with intervention have been studied both numerically [27, 28] and using some analytical methods [29, 30] in line with the COVID-19 epidemic.

In this study, we provide exact solutions of a simple SIR model with one-shot intervention, represented by a single discrete reduction in the transmission parameter during an epidemic. These solutions describe the dependence of the peak number of infected patients on the reproduction number under consideration of the implementation of NPIs and intervention timing. Theoretical and numerical analyses revealed non-trivial relations among the intensity of suppression of pandemics via NPIs, the number of infected individuals at the peaks, and the final size of infection cases. The methods and results shown in this study provide basic theoretical understanding in the context of the evaluation of NPIs.

Materials and methods

In this study, we focus on the dependence of the maximum fraction of infected individuals on the timing of the NPIs. Note that this analysis is motivated by the COVID-19 epidemic, but we consider a hypothetical epidemic, whose properties are not necessarily the same as that of COVID-19.

Model

We here introduce the SIR model, where the time evolution of the fraction of susceptible ( $s(t)$ ), infected ( $i(t)$ ), and removed ( $r(t)$ ) individuals is given by the following ordinary differential equation (ODE):

$\displaystyle\frac{ds}{dt}$	$\displaystyle=$	$\displaystyle-\beta si,$	(1)
$\displaystyle\frac{di}{dt}$	$\displaystyle=$	$\displaystyle\beta si-\gamma i,$	(2)
$\displaystyle\frac{dr}{dt}$	$\displaystyle=$	$\displaystyle\gamma i,$	(3)

taking into account the one-shot intervention and the second wave after the intervention. For simplicity, the total population is assumed to be unity, that is, $s(t)+i(t)+r(t)=1$ holds. The summation of the right-hand sides of Eqs. (1)–(3) vanishes, which guarantees conservation of the total population. Therefore, as seen below, only Eqs. (1) and (2) are numerically integrated to obtain the time evolutions of $s(t)$ , $i(t)$ , and $r(t)$ . In this model, the basic reproduction number and the effective reproduction number at time $t$ are given by $R_{0}=\beta/\gamma$ and $R_{t}=\beta s(t)/\gamma$ , respectively.

Kermack and McKendrick [6] derived the equation for the maximum fraction of infected individuals and showed that the final size is obtained by solving a transcendental equation with a given $R_{0}$ in the SIR model without any intervention [7, 8, 9]. By applying the technique they employed in the derivation, we here provide the relationship between the fractions of susceptible and removed individuals at arbitrary times $t_{0}$ and $t_{1}$ , that is, $s(t_{0})$ , $r(t_{0})$ , and $s(t_{1})$ , $r(t_{1})$ as

\displaystyle s(t_{1})

\displaystyle=

\displaystyle s(t_{0})\exp\bigg{\{}-\frac{\beta}{\gamma}[r(t_{1})-r(t_{0})]\bigg{\}}.

(4)

See Supplementary Section 1, Additional File 1 for details of the derivation. All analytical results presented in this paper are derived based on this equality.

Non-pharmaceutical intervention

In the present framework, an NPI in an isolated population is represented by a change in the transmission rate $\beta$ . Let $\beta=\beta_{\rm off}(>\gamma)$ be the transmission rate without the intervention, and it is switched to $\beta_{\rm on}(<\beta_{\rm off})$ when the intervention starts at $t=t_{\rm on}$ , and restored to $\beta_{\rm off}$ at $t=t_{\rm off}=t_{\rm on}+\Delta t$ . Let the corresponding basic reproduction numbers be $R_{0,\rm off}=\beta_{\rm off}/\gamma$ and $R_{0,\rm on}=\beta_{\rm on}/\gamma$ , respectively. Here, $\Delta t$ denotes the duration of the intervention. In this paper, we study the effect of one-shot intervention, one of the simplest implementation schemes, where the intervention is implemented only once. The mathematical methods employed in this study can be applied to more complex cases, for example, where multiple interventions are implemented intermittently. If the fraction of susceptible individuals remains large enough and the herd immunity is not achieved after the intervention, a second wave occurs. It is thus necessary to consider this second wave to evaluate the effect of the intervention.

The one-shot intervention setting has been numerically studied by Bootsma and Ferguson [26], showing that the final size, which depends on the timing $t_{\rm on}$ and $t_{\rm off}$ , can be smaller for weaker intervention with larger $R_{0,\rm on}$ . Di Lauro et al. [28] numerically studied the dependence of the final size and the peak fraction of infected individuals on the timing of the intervention, where the intervention was assumed to start when the fraction of infected and recovered individuals exceeds a threshold value. In particular, they concluded that the onset timing should be chosen so that the two peaks during and after the intervention are comparable. Morris et al. [29] showed that the maximum fraction of infected individuals with one-shot interventions can approach that achieved by the optimal intervention, which requires an unrealistic intervention, such as $R_{0,\rm on}=0$ . Sadeghi et al. [30] also suggested the existence of the optimal timing of the intervention based on discussion using a linearized equation and numerical simulation.

Analysis

Equations (1) and (2) were numerically integrated using ODEPACK [31], which consists of nine Fortran solvers for the initial value problem for ODE systems. The time step was set to $0.001$ . The nonintervention transmission rate and the recovery rate are fixed as $\beta_{\rm off}=2/7\ {\rm days}^{-1}$ and $\gamma=1/7\ {\rm day}^{-1}$ , that is, $R_{0,\rm off}=2$ . The initial condition $s(0)=1-\epsilon$ , $i(0)=\epsilon$ , and $r(0)=0$ , where $\epsilon=0.001$ is employed in this study. As discussed in detail below, the system behaves qualitatively differently for different values of $R_{0,\rm on}$ . Here, we primarily focus on the cases for $R_{0,\rm on}=1.4$ (Fig. 1) and $R_{0,\rm on}=0.7$ (Fig. 2), corresponding to relatively weaker and stronger intervention intensities, respectively. Note that $R_{0,\rm on}>1$ in the former case implies that the infection may spread even in the presence of the intervention. The time series without the intervention is shown in Figs. 1(A) and 2(A). They are identical, because $R_{0,\rm on}$ does not affect the dynamics without the intervention. Figures 1(B)–(E) and 2(B)–(E) show the time series with an intervention with a constant intervention duration $\Delta t=60$ days. The second wave can occur if the herd immunity is not achieved when the intervention ends (Figs. 1(C) and 2(B), (C)).

Refer to caption — Figure 1: \csentenceWeak intervention with large $R_{0,\rm on}=1.4>1$ . Parameters are $\beta_{\rm off}=2/7\ {\rm days}^{-1}$ and $\gamma=1/7\ {\rm days}^{-1}$ . (A) Time series of the fraction of infected and removed individuals, $i(t)$ in blue and $r(t)$ in orange, without the intervention. Time series with the intervention with $\beta_{\rm on}=1.4/7\ {\rm days}^{-1}$ and the intervention duration, $\Delta t=60$ days, with onset times (B) $t_{\rm on}=10$ days, (C) $t_{\rm on}=19.2$ days, (D) $t_{\rm on}=33.4$ days, and (E) $t_{\rm on}=65$ days, are depicted. The intervention is implemented for $t_{\rm on}\leq t\leq t_{\rm off}=t_{\rm on}+\Delta t$ and are represented by grey intervals in these panels. Conditions for the peaks of the fraction of infected individuals are given by $s(t)=1/R_{0,\rm off}$ , when the intervention is not implemented, and by $s(t)=1/R_{0,\rm on}$ during intervention. The final size of the outbreak for each case is represented by $r(\infty)$ . Panels (F) and (G) represent the maximum fraction of infected individuals $i_{\rm max}$ normalized by that without the intervention $i_{\rm max}^{0}$ , plotted in terms of $r(t_{\rm on})$ and $\Delta r=r(t_{\rm off})-r(t_{\rm on})$ , and $t_{\rm on}$ and $\Delta t$ , respectively. Symbols (A)-(E) in panels (F) and (G) denote the intervention timings of the time series of the corresponding panels in (A)-(E). Symbols (i)-(iv) in panel (F) denote the timing that the maximum infected fraction is observed, as described in the main text. The boundaries between regions shown in panel (F) are obtained analytically. See the main text for the details.

Maximum fraction of infected individuals

In addition to the numerically obtained time series, we analytically show that peaks of the fraction of infected individuals can appear at the following four timings by applying Eq. (4): after the intervention, during the intervention, at the onset of the intervention, and before the intervention. Note that there can be two peaks of the fraction of infected individuals if a second wave occurs. In such a case, the peak with a larger fraction gives the global maximum. Let us describe the four cases with respect to the timing of the maximum in detail.

•

(i) The maximum appears after the intervention (Figs. 1(B), (C), 2(B), and (C)). If the intervention ends before achieving herd immunity, a peak is observed during the second wave after the intervention. The fraction of infected individuals for this peak gives the global maximum if it is higher than the first peak before or during the intervention.
•

(ii) The maximum appears during the intervention (Fig. 1(D)). If the effective reproduction number $R_{t}$ declines and crosses unity during the intervention, there is a peak in this timing. The condition $R_{0,\rm on}>1$ is necessary for the existence of this peak, because the effective reproduction number has to be larger than unity at the onset of the intervention. Therefore, this peak does not appear for $R_{0,\rm on}<1$ . This peak is the global maximum if it is larger than the second peak. The peak also appears at this timing in Fig. 1(C), but the second peak is slightly higher than this peak.
•

(iii) The maximum appears at the onset of the intervention at $t=t_{\rm on}$ (Fig. 2(D)). If $R_{t_{\rm on}}<1$ holds when the intervention starts, the fraction of infected individuals decreases during the intervention, and we observe a local peak at $t=t_{\rm on}$ . There may be another peak if herd immunity is not achieved during the intervention. If the fraction of infected individuals at the second peak is less than at this peak, the timing of the global maximum is given at $t_{\rm on}$ . Note that this region also appears for $R_{0,\rm on}>1$ , although the time series is not shown in Figs. 1(B)–(E).
•

(iv) The maximum appears before the onset of the intervention (Figs. 1(E) and 2(E)). The fraction of infected individuals reaches its maximum before the intervention. This case implies that the intervention starts too late and fails to mitigate outbreaks in terms of the maximum fraction of infected individuals.

The fraction of infected individuals at the peaks can be calculated analytically. Conditions for the peaks are given in terms of the fraction of susceptible individuals as $s(t)=1/R_{0,{\rm off}}=\gamma/\beta_{\rm off}$ without the intervention and $s(t)=1/R_{0,{\rm on}}=\gamma/\beta_{\rm on}$ with the intervention. The maximum fraction of infected individuals for cases (ii), (iii), and (iv) does not depend on the fraction of removed individuals at the offset of the intervention $r(t_{\rm off})$ , because these peaks appear before $t_{\rm off}$ . For case (i), the maximum fraction of infected individuals $i_{\rm max}$ is given as

\displaystyle i_{\rm max}

\displaystyle=

\displaystyle 1-\left(1-\frac{R_{0,\rm on}}{R_{0,\rm off}}\right)\Delta r-\frac{1}{R_{0,\rm off}}\bigg{[}1+\log\left(R_{0,\rm off}\right)\bigg{]},

(5)

which depends on $\Delta r:=r(t_{\rm off})-r(t_{\rm on})$ , the difference in the fraction of removed individuals between the onset and offset of the intervention. See Supplementary Section 2, Additional File 1 for the explicit form of $i_{\rm max}$ for all cases and its derivation.

The boundaries between regions (i) and (ii), (i) and (iii), (ii) and (iii), and (iii) and (iv) can be obtained analytically with respect to $r(t_{\rm on})$ and $\Delta r$ (Figs. 1(F) and 2(F)). On the boundaries, the fraction of infected individuals at two peaks is comparable (Figs. 1(C) and 2(C)). See Supplementary Section 3, Additional File 1 for details of the derivation.

Final size with intervention

The final size also reflects the effect of the intervention. The final size of removed individuals $r(\infty)$ with the intervention is obtained by solving the equation

\displaystyle r(\infty)

\displaystyle=

\displaystyle 1-\exp\bigg{[}\big{(}R_{0,\rm off}-R_{0,\rm on}\big{)}\Delta r\bigg{]}\exp\bigg{[}-R_{0,\rm off}\ r(\infty)\bigg{]},

(6)

in a self-consistent manner [32]. Specifically, as $r(\infty)$ appears on both sides, this equation can be solved numerically or using the Lambert $W$ function, except for some special cases. This equation implies that the final size depends on $\Delta r=r(t_{\rm off})-r(t_{\rm on})$ . Another important implication of this equation is that $r(t_{\rm off})$ has the upper bound $\tilde{r}$ depending on $r(t_{\rm on})$ , which is given by

\displaystyle\tilde{r}(r(t_{\rm on}))

\displaystyle=

\displaystyle 1-\exp\bigg{[}-(R_{0,\rm off}-R_{0,\rm on})r(t_{\rm on})\bigg{]}\exp\bigg{[}-R_{0,\rm on}\ \tilde{r}\bigg{]}.

(7)

See Supplementary Section 4, Additional File 1 for the details of the derivation.

Using this equality, one can show that the final size in the presence of the intervention is always smaller than that without the intervention. For

\displaystyle R_{0,\rm on}

\displaystyle<

\displaystyle\frac{R_{0,\rm off}}{R_{0,\rm off}-1}\log\left(R_{0,\rm off}\right),

(8)

one can achieve $r(\infty)\approx 1-1/R_{0,\rm off}$ by setting $t_{\rm on}$ properly, which is the smallest prevalence to achieve herd immunity, with an intervention duration $\Delta t$ that is large enough [32]. See Supplementary Section 5, Additional File 1 for details of the derivation. Numerical results regarding the final size are summarized in Supplementary Section 6, Additional File 1.

Results

We report the numerical and analytical results, when the reproduction number under the intervention $R_{0,\rm on}$ is large (Fig. 1) and small (Fig. 2), showing qualitatively different behaviors.

In Figs. 1(F), (G), 2(F), and (G), the dependence of the maximum fraction of infected individuals on the timing of the intervention is plotted. Here, $i_{\rm max}$ is normalized by that in the absence of the intervention $i_{\rm max}^{0}$ (Figs. 1(A) and 2(A)). As the maximum infected fraction drops in the presence of the intervention, $i_{\rm max}/i_{\rm max}^{0}$ is less than unity and quantifies the effectiveness of the intervention in terms of the maximum fraction of infected individuals. As this ratio decreases, the intervention shows more success in reducing the maximum fraction.

It is difficult to obtain the time series of the SIR model analytically without any approximations, for example, linearization, or the method presented in [6]. Therefore, $i_{\rm max}$ is numerically computed with different onset and offset times for the intervention, $t_{\rm on}$ and $t_{\rm off}$ (Figs. 1(G) and 2(G)). However, $i_{\rm max}$ can be analytically calculated with respect to the fraction of the recovered individuals at the onset, $r(t_{\rm on})$ , and the offset, $r(t_{\rm off})$ , of the intervention (Figs. 1(F) and 2(F)). Note that there exists a one-to-one correspondence between $\{t_{\rm on},t_{\rm off}\}$ (panels (F)) and $\{r(t_{\rm on}),r(t_{\rm off})\}$ (panels (G)) in Figs. 1 and 2.

Weak intervention (large $R_{0,\rm on}=1.4$ )

In this case, the reproduction number is larger than unity even in the presence of the intervention. Therefore, the fraction of infected individuals may increase during the intervention period. The peaks of infected individuals can be observed during the intervention, and the maximum infected fraction can appear at any of the four timings (i)–(iv) classified above. Figures 1(B)–(E) show that $i_{\rm max}$ is minimized in the intermediate onset time $t_{\rm on}$ near (C), where the peaks during and after the intervention are comparable. This non-monotonic dependence on $t_{\rm on}$ is clearly visualized in Fig. 1(G). The peak of infected individuals during the intervention is smaller than that without intervention. This intervention mitigates the second wave.

As shown in Eq. (5) and Fig. 1(F), the maximum fraction of infected individuals is linear in $\Delta r$ if the peak of the second wave is the maximum, that is, case (i). This is verified by the fact that the contours lie horizontally in case (i). To clarify this point, the contours are explicitly shown in Supplementary Section 7, Additional File 1. If the fraction of the infected individuals reaches the maximum during or at the onset of the intervention (cases (ii) and (iii), respectively), the maximum fraction depends only on $t_{\rm on}$ , which has one-to-one correspondence to $r(t_{\rm on})$ . This is verified in Fig. 1(F), (G), and Supplementary Section 7, Additional File 1. If the maximum appears before the intervention, that is, case (iv), the maximum is independent of both $t_{\rm on}$ and $t_{\rm off}$ .

Onset and offset times for the intervention with constant $\Delta t$ corresponding to Figs. 1(A)–(E) are plotted in Figs. 1(F) and (G). As suggested by the time series, the maximum infected fraction is smallest in the intermediate intervention onset $t_{\rm on}$ , near (C). It is clear from panel (F) that this point is located close to the boundary between cases (i) and (ii), where peaks during and after the intervention are comparable. In case (ii), the maximum does not depend on $\Delta r$ , which implies that a longer intervention does not reduce the maximum. Note that the relationship between $(t_{\rm on},t_{\rm off})$ and $(r(t_{\rm on}),r(t_{\rm off}))$ is non-monotonous.

Strong intervention (small $R_{0,\rm on}=0.7$ )

When the transmission rate is small during the intervention, the maximum infected fraction is minimized in the intermediate starting time of the intervention $t_{\rm on}$ , namely, early implementation of the intervention does not necessarily minimize the infected fraction (Figs. 2(F), (G)). This result is intuitively understood as follows. If the intervention starts too early, the infection does not spread because of the small intervention transmission rate. Therefore, the second wave after the intervention is large, thus the early intervention is not effective. If the timing of the intervention is characterized in terms of $r(t_{\rm on})$ and $r(t_{\rm off})$ (Fig. 2(F)), the maximum of the second wave depends only on $\Delta r$ . If the maximum fraction is found during the intervention, its value is independent of $t_{\rm off}$ and $r(t_{\rm off})$ . The peak does not appear during the intervention, that is, case (ii) does not appear because $R_{t,\rm on}<1$ holds.

Maximum infected fraction $\bar{i}_{\rm max}$ versus reproduction number under the intervention $R_{0,\rm on}$

For each $R_{0,\rm on}$ , there exist onset and offset timings of the intervention that minimize the maximum fraction of infected individuals $i_{\rm max}$ (Figs. 1(F), (G) and 2(F), (G)). Let this value be $\bar{i}_{\rm max}(R_{0,\rm on})$ . Figure 3 plots the dependence of $\bar{i}_{\rm max}$ on $R_{0,\rm on}$ . As seen in the figure, it is minimized at a non-trivial intermediate value of $R_{0,\rm on}=R_{0,\rm on}^{\ast}\approx 1.23>1$ . This implies that the maximum fraction of infected individuals is minimized for a weak intervention under the one-shot condition. For $R_{0,\rm on}\geq R_{0,\rm on}^{\ast}$ , $\bar{i}_{\rm max}$ is achieved at the boundary between regions (i) and (ii) at $t_{\rm on}=0$ , where the peaks during and after the intervention are comparable (Fig. 1(C)). For $R_{0,\rm on}^{\ast}\geq R_{0,\rm on}\geq R_{0,\rm on}^{(1)}$ , the boundary between regions (i) and (ii) with $t_{\rm on}>0$ gives $\bar{i}_{\rm max}$ , where $R_{0,\rm on}^{(1)}\approx 1.08$ is the parameter value below which region (ii) does not exist. For strong intervention $R_{0,\rm on}\leq R_{0,\rm on}^{(1)}$ , $\bar{i}_{\rm max}$ is found at the boundary between regions (i) and (iii), where the peaks of the onset and after the intervention are comparable (Fig. 2(C)), at $t_{\rm on}>0$ . Conditions for $\bar{i}_{\rm max}$ are analytically derived in all cases, shown by the solid line in Fig. 3. The conditions for $\bar{i}_{\rm max}$ for $R_{0,\rm on}\geq R_{0,\rm on}^{(1)}$ can be explicitly solved. The condition for small $R_{0,\rm on}\leq R_{0,\rm on}^{(1)}$ cannot be solved explicitly, and the parametric equations for $\bar{i}_{\rm max}$ and $R_{0,\rm on}$ are used to plot the theoretical curve in Fig. 3. See Supplementary Section 8, Additional File 1 for details of the derivation.

Intervention strategies

It is possile to optimize the onset and offset timings of the NPI by minimizing an objective function under certain constraints. In the present framework, this can be formulated as an optimization problem in a hybrid nonlinear dynamical system. The optimal intervention strategy depends on the objective function. Here, we discuss the following two simple scenarios to minimize $i_{\rm max}$ . In general, more complex objective functions can be used for optimization. The following discussion may provide some intuition for considering such cases, taking into account the second wave.

Minimizing $i_{\rm max}$ with a constraint in the intervention duration

Let us start with a case where the intervention duration $\Delta t$ is less than a certain value. As $i_{\rm max}$ is monotonically decreasing in $\Delta t$ for a fixed $t_{\rm on}$ , we can assume that $\Delta t$ is a constant, and $t_{\rm on}$ is varied. This case has been studied in Figs. 1(B)–(E) and 2(B)–(E). This scheme fixes the intervention duration, so it is easier to anticipate the economic impact of the intervention, which depends on the duration of the intervention, than in the next scenario based on the fraction of recovered individuals.

As discussed above, $i_{\rm max}$ is minimized in the intermediate onset time $t_{\rm on}$ , when the first peak during the intervention and the second peak after the intervention are comparable. For larger $R_{0,\rm on}$ , $t_{\rm on}$ giving the minimum $i_{\rm max}$ converges to zero for large $\Delta t$ (Fig. 1(G)), but converges to non-zero for large $\Delta t$ for smaller $R_{0,\rm on}$ (Fig. 2(G)). These results are understood as follows: for a strong intervention with small $R_{0,\rm on}$ , the intervention immediately suppresses the fraction of infected individuals. Therefore, the early onset of the intervention prevents more individuals from achieving immunity during the intervention and eventually increases the maximum fraction of infected individuals.

Minimizing intervention duration with a constraint in $i_{\rm max}$

Another possible constraint is to minimize the intervention duration $\Delta t$ , keeping the maximum fraction of infected individuals $i_{\rm max}$ constant. This corresponds to choosing $t_{\rm on}$ along a contour of $i_{\rm max}$ , which prevents the overcapacity of medical support. In Figs. 1(G) and 2(G), $\Delta t$ is minimized for an intermediate $t_{\rm on}$ . Namely, a contour that crosses point (C) ( $\Delta t=60$ ) reaches $\Delta t\approx 20$ if the onset time is later than $t_{\rm on}=19.2$ (Fig. 1(C)). As we have already discussed, the maximum fraction depends on $\Delta r$ in region (i). It takes a shorter time to achieve the same $\Delta r$ if the intervention starts later in this case. As the maximum fraction of infected individuals depends on $r(t_{\rm on})$ but is independent of $\Delta r$ in cases (ii) and (iii), the optimal $r(t_{\rm on})$ is determined by the boundary between regions (i) and (ii) or that between (i) and (iii) in Figs. 1(F) and 2(F). Evidently this timing is an intermediate value of $t_{\rm on}>0$ . It may be easier to set a plan for the intervention with respect to $r(t_{\rm on})$ and $r(t_{\rm off})$ , rather than the timing $t_{\rm on}$ and $t_{\rm off}$ .

Similar optimization problems can be considered using the final size as the objective functions. These problems are discussed in Supplementary Section 9, Additional File 1.

Conclusion

In the present study, in consideration of the actual COVID-19 situation, we studied the situation in which the reproduction number of an infectious disease is temporarily reduced by implementing an NPI once during the epidemic, using a simple mathematical model. The results provide theoretical implications as to how strong NPIs should be introduced during an epidemic of an emerging infectious disease. If the effective reproduction number during the intervention is too small, the fraction of infected individuals at the peak in the second wave may be higher than the first peak. It was also shown numerically and analytically that the fraction of infected individuals can also increase if the intervention is started too early. The upper limit of medical capacity is an essential practical constraint. In particular, for infectious diseases such as COVID-19, which is too emergent to expect effective treatments, it is more important to avoid exhausting the medical system. This study suggests that it will be necessary to be alert for a larger second wave that may occur after strong intervention in such cases.

We analytically derived the peak fractions of infected individuals in the SIR model with the one-shot intervention. These analytical results suggest that the peak fractions can be smaller with non-trivial intervention timings. The maximum fraction is smallest for $R_{0,\rm on}>1$ , that is, the intervention reproduction number is not in the disease-free regime. In the literature, Bootsma and Ferguson [26] showed that the final size can be minimized for non-trivial intervention timing and $R_{0,\rm on}$ numerically. Di Lauro et al. [28] numerically showed that the peak fraction of infected individuals also depends on the timing of the intervention. We obtained analytical expressions for these quantities in this study. The analytical results for the final size are presented in Supplementary Section 6, Additional File 1. It should be noted that the formulae for the peak fraction depend on the timing of the peak, resulting in various cases compared with the final size. Some analytical results are available for this system. Sadeghi et al. [30] explained these non-trivial effects based on solutions of the linearized equation, which exponentially grows and decays without and with the intervention, respectively. Linearization is one of the simplest approximations and is applicable in this case; in particular, this approximation is useful in discussions regarding timing. Conversely, linearization cannot be used to discuss the important case where $R_{0,\rm on}>1$ , as the linearized equation cannot explain the declining number of infected individuals. The proposed method in this study provides a unified framework, including the cases where linearization is not feasible. Morris et al. derived [29] an equation for the peaks of the fraction of infected individuals in terms of $s(t_{\rm on})$ , $s(t_{\rm off})$ , $i(t_{\rm on})$ , and $i(t_{\rm off})$ . In this work, we further show the peak fractions in terms of the two parameters $r(t_{\rm on})$ and $r(t_{\rm off})$ for a general initial condition.

Although concrete measures and guidelines for COVID-19 are required, it is emphasized again that the results in this study were derived using a simplified model with many assumptions. First, the results presented in this study are based on the simplified SIR model, which takes into account neither the realistic pathology of COVID-19 nor societal response. This model assumes uniform and random contact within a group and does not consider interactions between different subgroups in the population. Recovered patients are assumed to have complete immunity in this model. These assumptions are not applicable to COVID-19, where there are still many unclear factors regarding heterogeneity in contact networks and the immune response of patients after recovery. If a society develops other public health measures during the intervention, the basic reproduction numbers before and after the NPI may be different. Next, this study is a one-shot intervention model, in which the transmission coefficient returns to the original value after a single intervention. Practically, each re-pandemic requires multiple intermittent interventions [33], making the intervention process much more complex [34]. Furthermore, the model is based on a deterministic dynamical system of an isolated population. Therefore, important NPI measures such as border control cannot be estimated in the present framework. The deterministic nature of the model assumes that the infectious disease cannot be eradicated, as the number of infected individuals remains non-zero for a finite time. If the population is small enough and no imported cases are assumed, strong and early intervention, which is not necessarily recommended in this study, may eradicate the disease, and a second wave does not occur. Strong intervention would be necessary in other cases, for example, when the number of infected patients approaches the capacity of the medical system. These complex situations may be analyzed in detail by extending the methods presented in this study, which may lead to different conclusions from those reached in this study.

In this study, the effect of NPIs was modeled as a hybrid dynamical system, which may further enable us to approach more refined models in future investigations. The influence of NPIs with respect to political decisions and behavioral changes of people can be expressed more accurately by introducing a hybrid dynamical system. In recent years, dynamical systems theory and control theory have been developed, and phenomena specific to hybrid systems such as Zeno solutions and sliding motions have been discussed [18, 19, 20, 21, 22]. Some studies, such as [35] and [36], have proposed control of infectious diseases using sliding mode control. The optimal policy under NPIs can be discussed by modeling the effect of economic damages associated with execution and then minimizing the cost function. Discrete changes in parameters such as the transmission coefficient in NPI implementation are strongly linked to the intensity of measures, suppression of economic activity, and changes in human mobility. For example, the correlation between the decrease in human mobility with NPIs and the effective reproduction number for the COVID-19 pandemic has been studied [37, 38]. Such studies can contribute to modeling the costs of NPIs. The mathematical model of the epidemic suppression effect can be constructed using a hybrid dynamical system, taking into account the negative socioeconomic impact of NPIs.

Competing interests

The authors declare that they have no competing interests.

Author’s contributions

N.F., T.W., and K.A. designed the study. N.F. derived the analytical results. T.O., S.T., and K.A. provided comments from theoretical viewpoints. T.O. conducted the numerical simulations. T.W., J.S., and T.N. discussed the implications of the study with respect to public health.

Acknowledgements

K.A. and N.F. are grateful to Prof. H. Inaba for their valuable comments. T.N. and N.F. are supported by Starting Grants for Research toward Resilient Society (SGRRS), Tohoku University. N.F. is supported by JSPS KAKENHI Grant Number JP18K11462. K.A. is partially supported by Moonshot R&D Grant Number JPMJMS2021.

References

[1] Aledort, J.E., Lurie, N., Wasserman, J., Bozzette, S.A.: Non-pharmaceutical public health interventions for pandemic influenza: an evaluation of the evidence base. BMC public health 7(1), 208 (2007)
[2] Cohen, J., Kupferschmidt, K.: Countries test tactics in ‘war’ against covid-19. Science 367(6484), 1287–1288 (2020)
[3] Flaxman, S., Mishra, S., Gandy, A., Unwin, H., Coupland, H., Mellan, T., Zhu, H., Berah, T., Eaton, J., Perez Guzman, P., et al.: Report 13: Estimating the number of infections and the impact of non-pharmaceutical interventions on covid-19 in 11 european countries (2020)
[4] Chinazzi, M., Davis, J.T., Ajelli, M., Gioannini, C., Litvinova, M., Merler, S., Pastore y Piontti, A., Mu, K., Rossi, L., Sun, K., Viboud, C., Xiong, X., Yu, H., Halloran, M.E., Longini, I.M., Vespignani, A.: The effect of travel restrictions on the spread of the 2019 novel coronavirus (covid-19) outbreak. Science 368(6489), 395–400 (2020)
[5] Kraemer, M.U., Yang, C.-H., Gutierrez, B., Wu, C.-H., Klein, B., Pigott, D.M., Du Plessis, L., Faria, N.R., Li, R., Hanage, W.P., et al.: The effect of human mobility and control measures on the covid-19 epidemic in china. Science 368(6490), 493–497 (2020)
[6] Kermack, W.O., McKendrick, A.G.: A contribution to the mathematical theory of epidemics. Proceedings of the royal society of london. Series A, Containing papers of a mathematical and physical character 115(772), 700–721 (1927)
[7] Anderson, R.M., May, R.M.: Infectious Diseases of Humans: Dynamics and Control. Oxford university press, Oxford (1992)
[8] Chowell, G., Brauer, F.: The basic reproduction number of infectious diseases: computation and estimation using compartmental epidemic models. In: Mathematical and Statistical Estimation Approaches in Epidemiology, pp. 1–30. Springer, Dordrecht (2009)
[9] Inaba, H.: Age-structured Population Dynamics in Demography and Epidemiology. Springer, Singapore (2017)
[10] Funk, S., Salathé, M., Jansen, V.A.A.: Modelling the influence of human behaviour on the spread of infectious diseases: a review. Journal of The Royal Society Interface 7(50), 1247–1256 (2010)
[11] Merler, S., Ajelli, M., Pugliese, A., Ferguson, N.M.: Determinants of the spatiotemporal dynamics of the 2009 h1n1 pandemic in europe: Implications for real-time modelling. PLOS Computational Biology 7(9), 1–13 (2011)
[12] He, D., Dushoff, J., Day, T., Ma, J., Earn, D.J.: Inferring the causes of the three waves of the 1918 influenza pandemic in england and wales. Proceedings of the Royal Society B: Biological Sciences 280(1766), 20131345 (2013)
[13] Wearing, H.J., Rohani, P., Keeling, M.J.: Appropriate models for the management of infectious diseases. PLOS Medicine 2(8), 320 (2005)
[14] Lin, F., Muthuraman, K., Lawley, M.: An optimal control theory approach to non-pharmaceutical interventions. BMC infectious diseases 10(1), 32 (2010)
[15] Perrings, C., Castillo-Chavez, C., Chowell, G., Daszak, P., Fenichel, E.P., Finnoff, D., Horan, R.D., Kilpatrick, A.M., Kinzig, A.P., Kuminoff, N.V., et al.: Merging economics and epidemiology to improve the prediction and management of infectious disease. EcoHealth 11(4), 464–475 (2014)
[16] Zhou, Y., Yang, K., Zhou, K., Liang, Y.: Optimal vaccination policies for an sir model with limited resources. Acta biotheoretica 62(2), 171–181 (2014)
[17] Ketcheson, D.I.: Optimal control of an sir epidemic through finite-time non-pharmaceutical intervention. arXiv preprint arXiv:2004.08848 (2020)
[18] Filippov, A.F.: Differential Equations with Discontinuous Righthand Sides: Control Systems. Kluwer, Dordrecht (1988)
[19] Lunze, J., Lamnabhi-Lagarrigue, F.: Handbook of Hybrid Systems Control: Theory, Tools, Applications. Cambridge University Press, Cambridge (2009)
[20] Aihara, K., Suzuki, H.: Theory of hybrid dynamical systems and its applications to biological and medical systems. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 368(1930), 4893–4914 (2010)
[21] di Bernardo, M., Hogan, S.J.: Discontinuity-induced bifurcations of piecewise smooth dynamical systems. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 368(1930), 4915–4935 (2010)
[22] Heemels, W.P.M.H., De Schutter, B., Lunze, J., Lazar, M.: Stability analysis and controller synthesis for hybrid dynamical systems. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 368(1930), 4937–4960 (2010)
[23] Agur, Z., Cojocaru, L., Mazor, G., Anderson, R.M., Danon, Y.L.: Pulse mass measles vaccination across age cohorts. Proceedings of the National Academy of Sciences 90(24), 11698–11702 (1993)
[24] Wang, A., Xiao, Y., Cheke, R.A.: Global dynamics of a piece-wise epidemic model with switching vaccination strategy. Discrete and Continuous Dynamical Systems-Series B (DCDS-B) 19(9), 2915–2940 (2014)
[25] Chladná, Z., Kopfová, J., Rachinskii, D., Rouf, S.C.: Global dynamics of sir model with switched transmission rate. Journal of Mathematical Biology 80(4), 1209–1233 (2020)
[26] Bootsma, M.C., Ferguson, N.M.: The effect of public health measures on the 1918 influenza pandemic in us cities. Proceedings of the National Academy of Sciences 104(18), 7588–7593 (2007)
[27] Anderson, R.M., Heesterbeek, H., Klinkenberg, D., Hollingsworth, T.D.: How will country-based mitigation measures influence the course of the covid-19 epidemic? The Lancet 395(10228), 931–934 (2020)
[28] Di Lauro, F., Kiss, I.Z., Miller, J.: The timing of one-shot interventions for epidemic control. medRxiv 2020.03.02.20030007 (2020)
[29] Morris, D.H., Rossine, F.W., Plotkin, J.B., Levin, S.A.: Optimal, near-optimal, and robust epidemic control. arXiv:2004.02209 (2020)
[30] Sadeghi, M., Greene, J., Sontag, E.: Universal features of epidemic models under social distancing guidelines. bioRxiv 2020.06.21.163931 (2020)
[31] Hindmarsh, A.C.: Odepack, a systematized collection of ode solvers. Scientific computing, 55–64 (1983)
[32] Fujiwara, N., Onaga, T., Wada, T., Aihara, K.: Effects of infection control policies on the final size. Seisan Kenkyu (in Japanese) 72(2), 141–143 (2020)
[33] Ferguson, N., Laydon, D., Nedjati Gilani, G., Imai, N., Ainslie, K., Baguelin, M., Bhatia, S., Boonyasiri, A., Cucunuba Perez, Z., Cuomo-Dannenburg, G., et al.: Report 9: Impact of non-pharmaceutical interventions (npis) to reduce covid19 mortality and healthcare demand (2020)
[34] Kissler, S.M., Tedijanto, C., Goldstein, E., Grad, Y.H., Lipsitch, M.: Projecting the transmission dynamics of sars-cov-2 through the postpandemic period. Science 368(6493), 860–868 (2020)
[35] Xiao, Y., Xu, X., Tang, S.: Sliding mode control of outbreaks of emerging infectious diseases. Bulletin of mathematical biology 74(10), 2403–2422 (2012)
[36] Khalili Amirabadi, R., Heydari, A., Zarrabi, M.: Analysis and control of seir epedemic model via sliding mode control. Advanced Modeling and Optimization 18(1), 153–162 (2016)
[37] Badr, H.S., Du, H., Marshall, M., Dong, E., Squire, M.M., Gardner, L.M.: Association between mobility patterns and covid-19 transmission in the usa: a mathematical modelling study. The Lancet Infectious Diseases (2020)
[38] Yabe, T., Tsubouchi, K., Fujiwara, N., Wada, T., Sekimoto, Y., Ukkusuri, S.V.: Non-compulsory measures sufficiently reduced human mobility in japan during the covid-19 epidemic. Scientific Reports 10, 18053 (2020)

Additional Files

Additional file 1 — Supplementary texts

Detailed discussions carried out in the main text are given.

Supplementary materials: Analytical estimation of maximum fraction of infected individuals with one-shot non-pharmaceutical intervention in a hybrid epidemic model

S1 Method: Derivation of Eq. (4) [6, 7, 8, 9]

We want to derive Eq. (4) in the main text from the evolution equation of the SIR model

$\displaystyle\frac{ds}{dt}$	$\displaystyle=$	$\displaystyle-\beta si,$	(S9)
$\displaystyle\frac{di}{dt}$	$\displaystyle=$	$\displaystyle\beta si-\gamma i,$	(S10)
$\displaystyle\frac{dr}{dt}$	$\displaystyle=$	$\displaystyle\gamma i.$	(S11)

We integrate Eq. (1) from time $t_{0}$ to $t_{1}$ :

\displaystyle s(t_{1})

\displaystyle=

\displaystyle s(t_{0})\exp\left[-\beta\int_{t_{0}}^{t_{1}}i(t^{\prime})dt^{\prime}\right].

(S13)

The integral on the right-hand side can be replaced with the integration of Eq. (3),

\displaystyle r(t_{1})-r(t_{0})

\displaystyle=

\displaystyle\gamma\int_{t_{0}}^{t_{1}}i(t^{\prime})dt^{\prime}.

(S14)

Then, we obtain

\displaystyle s(t_{1})

\displaystyle=

\displaystyle s(t_{0})\exp\bigg{\{}-\frac{\beta}{\gamma}[r(t_{1})-r(t_{0})]\bigg{\}},

(S15)

which is Eq. (4) in the main text. We apply this equation to various $t_{0}$ and $t_{1}$ .

S2 Method: Analytical expressions for the maximum fraction of infected individuals

As described in the main text, there are four possible timings where the maximum fraction of infected individuals $i_{\rm max}$ is observed. We summarize the analytical forms of the maximum fraction in all cases with respect to the fractions of removed individuals at the onset $r(t_{\rm on})$ and the offset $r(t_{\rm off})$ of the intervention, that is, we assume that the timing of the onset and offset of the intervention is determined by the fraction of removed individuals.

•

(i) The maximum appears after the intervention:

\displaystyle i_{\rm max}

\displaystyle=

\displaystyle 1-\left(1-\frac{R_{0,\rm on}}{R_{0,\rm off}}\right)\Delta r-\frac{1}{R_{0,\rm off}}\bigg{[}1+\log\left(R_{0,\rm off}\right)\bigg{]},

(S17)

where $\Delta r=r(t_{\rm off})-r(t_{\rm on})$ .

•

(ii) The maximum appears during the intervention:

\displaystyle i_{\rm max}

\displaystyle=

\displaystyle 1+\left(\frac{R_{0,\rm off}}{R_{0,\rm on}}-1\right)r(t_{\rm on})-\frac{1}{R_{0,\rm on}}\bigg{[}1+\log\left(R_{0,\rm on}\right)\bigg{]}.

(S18)

•

(iii) The maximum appears at the onset of the intervention:

$\displaystyle i_{\rm max}$ $\displaystyle=$ $\displaystyle 1-\exp\left[-R_{0,\rm off}\ r(t_{\rm on})\right]-r(t_{\rm on}).$ (S19)
•

(iv) The maximum appears before the onset of the intervention:

$\displaystyle i_{\rm max}$ $\displaystyle=$ $\displaystyle 1-\frac{1}{R_{0,\rm off}}\bigg{[}1+\log\left(R_{0,\rm off}\right)\bigg{]}.$ (S20)

Derivation

Here, we derive Eqs. (S17)–(S20). The conditions for the fraction of infected individuals to be the local peak are given by $s(t)=\gamma/\beta_{\rm off}=1/R_{0,\rm off}$ without the intervention and $s(t)=\gamma/\beta_{\rm on}=1/R_{0,\rm on}$ with the intervention.

First, let us apply Eq. (4) until the onset time of the intervention, with $\beta=\beta_{\rm off}$ , $t_{1}=t_{\rm on}$ , and $t_{0}=0$ . Using the initial condition $s(0)=1-\epsilon$ , $i(0)=\epsilon$ , $r(0)=0$ , we obtain

\displaystyle s(t_{\rm on})

\displaystyle=

\displaystyle(1-\epsilon)\exp\bigg{[}-\frac{\beta_{\rm off}}{\gamma}r(t_{\rm on})\bigg{]}.

(S21)

In the following, we take the limit $\epsilon\rightarrow 0$ .

Substituting Eq. (S21) into Eq. (4), one can derive the fraction of susceptible and removed individuals during and after the intervention.

(i) Maximum fraction after the intervention

We derive the fraction of infected individuals after the intervention. Let $t_{p}^{\rm after}>t_{\rm off}$ be the time at which the infected fraction reaches the peak after the intervention. The relationship between $s(t_{\rm off})$ and $r(t_{\rm off})$ is derived by substituting Eq. (S21) with $\epsilon\rightarrow 0$ into Eq. (4) with $t_{0}=t_{\rm on}$ and $t_{1}=t_{\rm off}$ as follows:

\displaystyle s(t_{\rm off})

\displaystyle=

\displaystyle\exp\bigg{[}-\frac{\beta_{\rm on}}{\gamma}r(t_{\rm off})\bigg{]}\exp\bigg{[}-\frac{\beta_{\rm off}-\beta_{\rm on}}{\gamma}r(t_{\rm on})\bigg{]}.

(S22)

Then, by applying Eq. (4) setting $t_{0}=t_{\rm off}$ and $t_{1}=t_{p}^{\rm after}$ , we have

\displaystyle s(t_{p}^{\rm after})

\displaystyle=

\displaystyle\exp\bigg{[}-\frac{\beta_{\rm off}}{\gamma}r(t_{p}^{\rm after})\bigg{]}\exp\bigg{\{}\frac{\beta_{\rm off}-\beta_{\rm on}}{\gamma}\bigg{[}r(t_{\rm off})-r(t_{\rm on})\bigg{]}\bigg{\}}.

(S23)

The condition for the peak after the intervention is $s(t_{p}^{\rm after})=\gamma/\beta_{\rm off}$ . Substituting this condition into Eq. (S23), the fractions of removed and infected individuals are given by:

$\displaystyle r(t_{p}^{\rm after})$	$\displaystyle=$	$\displaystyle\frac{\beta_{\rm off}-\beta_{\rm on}}{\beta_{\rm off}}[r(t_{\rm off})-r(t_{\rm on})]+\frac{\gamma}{\beta_{\rm off}}\log\left(\frac{\beta_{\rm off}}{\gamma}\right),$	(S24)
$\displaystyle i(t_{p}^{\rm after})$	$\displaystyle=$	$\displaystyle 1-s(t_{p}^{\rm after})-r(t_{p}^{\rm after})$	(S25)
	$\displaystyle=$	$\displaystyle 1-\frac{\beta_{\rm off}-\beta_{\rm on}}{\beta_{\rm off}}[r(t_{\rm off})-r(t_{\rm on})]-\frac{\gamma}{\beta_{\rm off}}\bigg{[}1+\log\left(\frac{\beta_{\rm off}}{\gamma}\right)\bigg{]},$	(S26)

which is equivalent to Eq. (S17) with $R_{0,\rm off}=\beta_{\rm off}/\gamma$ and $R_{0,\rm on}=\beta_{\rm on}/\gamma$ . Equation (S26) gives $i_{\rm max}$ if the maximum appears after the intervention. Note that the fraction of infected individuals depends only on the onset and offset timings as $r(t_{\rm off})-r(t_{\rm on})$ . The later the onset or sooner the offset is, the higher the peak is.

(ii) Maximum fraction during the intervention

If the effective reproduction number at the onset of the intervention is greater than unity, $\beta_{\rm on}s(t_{\rm on})/\gamma>1$ , and that at the offset is less than unity, $\beta_{\rm on}s(t_{\rm off})/\gamma<1$ , then there is a peak during the intervention. Let $t_{p}^{\rm during}$ be the time at which the infected fraction reaches the peak, where $t_{\rm on}<t_{p}^{\rm during}<t_{\rm off}$ . We obtain the relationship between susceptible and removed individuals during the intervention by applying Eq. (4) with $t_{0}=t_{\rm on}$ and $t_{1}=t_{p}^{\rm during}$ and using Eq. (S21) as the relationship between $s(t_{\rm on})$ and $r(t_{\rm on})$ :

\displaystyle s(t_{p}^{\rm during})

\displaystyle=

\displaystyle\exp\bigg{[}-\frac{\beta_{\rm on}}{\gamma}r(t_{p}^{\rm during})\bigg{]}\exp\bigg{[}-\frac{\beta_{\rm off}-\beta_{\rm on}}{\gamma}r(t_{\rm on})\bigg{]}.

(S27)

The condition for the peak during the intervention is $s(t_{p}^{\rm during})=\gamma/\beta_{\rm on}$ . By substituting this condition into Eq. (S27), we obtain

\displaystyle r(t_{p}^{\rm during})

\displaystyle=

\displaystyle\frac{\gamma}{\beta_{\rm on}}\bigg{[}-\frac{1}{\gamma}(\beta_{\rm off}-\beta_{\rm on})r(t_{\rm on})+\log\left(\frac{\beta_{\rm on}}{\gamma}\right)\bigg{]}.

(S28)

Then, the fraction of infected individuals at this peak is

	$\displaystyle i(t_{p}^{\rm during})$	$\displaystyle=$	$\displaystyle 1-s(t_{p}^{\rm during})-r(t_{p}^{\rm during})$		(S29)
		$\displaystyle=$	$\displaystyle 1+\frac{\beta_{\rm off}-\beta_{\rm on}}{\beta_{\rm on}}r(t_{\rm on})-\frac{\gamma}{\beta_{\rm on}}\bigg{[}1+\log\left(\frac{\beta_{\rm on}}{\gamma}\right)\bigg{]},$		(S30)

which is equivalent to Eq. (S18). As $\beta_{\rm off}>\beta_{\rm on}$ and $\frac{\beta_{\rm off}-\beta_{\rm on}}{\beta_{\rm on}}>0$ , the fraction of infected individuals at this peak linearly increases compared with the fraction of removed individuals at the onset of the intervention $r(t_{\rm on})$ . Equation (S30) is the condition for the maximum if this peak is higher than another peak.

The conditions for the existence of this peak are $\beta_{\rm on}s(t_{\rm on})>\gamma$ and $\beta_{\rm on}s(t_{\rm off})<\gamma$ , which are interpreted as

	$\displaystyle r(t_{\rm on})$	$\displaystyle<$	$\displaystyle\frac{\gamma}{\beta_{\rm off}}\log\left(\frac{\beta_{\rm on}}{\gamma}\right),$		(S31)
	$\displaystyle r(t_{\rm off})$	$\displaystyle>$	$\displaystyle\frac{\gamma}{\beta_{\rm on}}\bigg{[}-\frac{1}{\gamma}(\beta_{\rm off}-\beta_{\rm on})r(t_{\rm on})+\log\left(\frac{\beta_{\rm on}}{\gamma}\right)\bigg{]}.$		(S32)

Therefore, $\beta_{\rm on}>\gamma$ is required for the existence of this peak.

(iii) Maximum fraction at the onset of the intervention

At the onset of the intervention, the fraction of infected individuals is given as

	$\displaystyle i(t_{\rm on})$	$\displaystyle=$	$\displaystyle 1-s(t_{\rm on})-r(t_{\rm on})$		(S33)
		$\displaystyle=$	$\displaystyle 1-\exp\bigg{[}-\frac{\beta_{\rm off}}{\gamma}r(t_{\rm on})\bigg{]}-r(t_{\rm on}),$		(S34)

by substituting Eq. (S21).

This peak appears if the effective reproduction number before the intervention is larger than unity $\beta_{\rm off}s(t_{\rm on})>\gamma$ and that after the onset of the intervention is less than unity $\beta_{\rm on}s(t_{\rm on})<\gamma$ . These conditions are summarized in terms of $r(t_{\rm on})$ as follows:

\displaystyle\max\left\{0,\frac{\gamma}{\beta_{\rm off}}\log\left(\frac{\beta_{\rm on}}{\gamma}\right)\right\}<r(t_{\rm on})<\frac{\gamma}{\beta_{\rm off}}\log\left(\frac{\beta_{\rm off}}{\gamma}\right),

(S35)

by substituting the above conditions into Eq. (S21).

(iv) Maximum fraction before the intervention

Let $t_{1}=t_{p}^{\rm before}$ be the time at which the peak appears before the onset of the intervention. We apply Eq. (4) with $\beta=\beta_{\rm off}$ , $t_{0}=0$ , and $t_{1}=t_{p}^{\rm before}$ . Then, we obtain

\displaystyle r(t_{p}^{\rm before})

\displaystyle=

\displaystyle\frac{\gamma}{\beta_{\rm off}}\log\left[\frac{1}{s(t_{p}^{\rm before})}\right].

(S36)

The peak condition $s(t_{p}^{\rm before})=\gamma/\beta_{\rm off}$ leads to

\displaystyle r(t_{p}^{\rm before})

\displaystyle=

\displaystyle\frac{\gamma}{\beta_{\rm off}}\log\left(\frac{\beta_{\rm off}}{\gamma}\right).

(S37)

Substituting this equation into the conservation of total population $i(t)=1-s(t)-r(t)$ , we obtain

\displaystyle i(t_{p}^{\rm before})

\displaystyle=

\displaystyle 1-\frac{\gamma}{\beta_{\rm off}}\left[1+\log\left(\frac{\beta_{\rm off}}{\gamma}\right)\right].

(S38)

This peak appears for $\displaystyle r(t_{\rm on})>\frac{\gamma}{\beta_{\rm off}}\log\left(\frac{\beta_{\rm off}}{\gamma}\right)$ , meaning a late onset of the intervention.

S3 Method: Derivation of boundaries between regions of different maxima

As shown in Figs. 1(F) and 2(F) in the main text, the timing giving the maximum fraction of infected individuals switches in the $(r(t_{\rm on}),\Delta r)$ plane, where $\Delta r=r(t_{\rm off})-r(t_{\rm on})$ . These figures suggest possible transitions between regions (i) and (ii), between (i) and (iii), between (ii) and (iii), and between (iii) and (iv) listed above. At the boundaries, the peaks of the different timings are expected to be equal. Based on this idea, we sketch the derivations of the equations for the boundaries between these timings.

Boundaries between regions (i) and (ii)

In regions where peaks of infected individuals during and after the intervention coexist, the global maximum is given by

$\displaystyle i_{\rm max}$	$\displaystyle=$	$\displaystyle\max[i(t_{p}^{\rm after}),i(t_{p}^{\rm during})]$	(S40)
	$\displaystyle=$	$\displaystyle\max\bigg{\{}1-\frac{\beta_{\rm off}-\beta_{\rm on}}{\beta_{\rm off}}[r(t_{\rm off})-r(t_{\rm on})]-\frac{\gamma}{\beta_{\rm off}}\bigg{[}1+\log\left(\frac{\beta_{\rm off}}{\gamma}\right)\bigg{]},$
		$\displaystyle 1+\frac{\beta_{\rm off}-\beta_{\rm on}}{\beta_{\rm on}}r(t_{\rm on})-\frac{\gamma}{\beta_{\rm on}}\bigg{[}1+\log\left(\frac{\beta_{\rm on}}{\gamma}\right)\bigg{]}\bigg{\}}.$

There occurs a transition of the timing giving the global maximum, and the condition at the boundary is $i(t_{p}^{\rm after})=i(t_{p}^{\rm during})$ , that is,

\displaystyle\Delta r

\displaystyle=

\displaystyle-\frac{\beta_{\rm off}}{\beta_{\rm on}}r(t_{\rm on})+\frac{\gamma}{\beta_{\rm on}}\left[1+\frac{\beta_{\rm off}\log\left(\frac{\beta_{\rm on}}{\gamma}\right)-\beta_{\rm on}\log\left(\frac{\beta_{\rm off}}{\gamma}\right)}{\beta_{\rm off}-\beta_{\rm on}}\right].

(S41)

This condition gives the boundary between regions (i) and (ii) in Fig. 1(F) in the main text. For $\Delta r$ smaller than this condition, the peak of the second wave is larger than that during the intervention.

Boundaries between (i) and (iii)

Two peaks at the onset of the intervention and after the intervention can coexist. The maximum fraction of infected individuals is

	$\displaystyle i_{\rm max}$	$\displaystyle=$	$\displaystyle\max[i(t_{\rm on}),i(t_{p}^{\rm after})]$
		$\displaystyle=$	$\displaystyle\max\bigg{\{}e^{-\frac{\beta_{\rm off}}{\gamma}r(t_{\rm on})}+r(t_{\rm on}),\frac{\beta_{\rm off}-\beta_{\rm on}}{\beta_{\rm off}}[r(t_{\rm off})-r(t_{\rm on})]+\frac{\gamma}{\beta_{\rm off}}\bigg{[}1+\log\left(\frac{\beta_{\rm off}}{\gamma}\right)\bigg{]}\bigg{\}}.$

The condition for the boundary $i(t_{\rm on})=i(t_{p}^{\rm after})$ can be rewritten as

\displaystyle\Delta r

\displaystyle=

\displaystyle\frac{\beta_{\rm off}}{\beta_{\rm off}-\beta_{\rm on}}\left[e^{-\frac{\beta_{\rm off}}{\gamma}r(t_{\rm on})}+r(t_{\rm on})\right]-\frac{\gamma}{\beta_{\rm off}-\beta_{\rm on}}\bigg{[}1+\log\left(\frac{\beta_{\rm off}}{\gamma}\right)\bigg{]}.

This boundary is shown in Figs. 1(F) and 2(F) in the main text.

Boundaries between (ii) and (iii)

Peaks at the onset and during the intervention cannot coexist. As shown in Eqs. (S31) and (S35),

\displaystyle r(t_{\rm on})

\displaystyle=

\displaystyle\frac{\gamma}{\beta_{\rm off}}\log\left(\frac{\beta_{\rm on}}{\gamma}\right),

(S45)

gives the boundary between regions (ii) and (iii). This condition does not depend on $\Delta r$ (Fig. 1(F) in the main text).

Boundaries between (iii) and (iv)

As shown in Eqs. (S35) and (S37),

\displaystyle r(t_{\rm on})

\displaystyle=

\displaystyle\frac{\gamma}{\beta_{\rm off}}\log\left(\frac{\beta_{\rm off}}{\gamma}\right),

(S46)

gives the condition for the boundary between regions (iii) and (iv). As in the case of the boundary between (ii) and (iii), this condition does not depend on $\Delta r$ (Figs. 1(F) and 2(F) in the main text).

S4 Method: Derivation of the final size equation with the intervention

As discussed in Eq. (S23), the relationship between $s(t_{1})$ and $r(t_{1})$ after the intervention at $t_{1}>t_{\rm off}$ is

\displaystyle s(t_{1})

\displaystyle=

\displaystyle\exp\bigg{[}-\frac{\beta_{\rm off}}{\gamma}r(t_{1})\bigg{]}\exp\bigg{\{}\frac{\beta_{\rm off}-\beta_{\rm on}}{\gamma}\bigg{[}r(t_{\rm off})-r(t_{\rm on})\bigg{]}\bigg{\}}.

(S47)

We obtain the final size equation for $t\rightarrow\infty$ with the intervention by substituting this into the condition $s(\infty)+r(\infty)=1$ as

\displaystyle r(\infty)

\displaystyle=

\displaystyle 1-\exp\bigg{[}-\frac{\beta_{\rm off}}{\gamma}r(\infty)\bigg{]}\exp\bigg{\{}\frac{\beta_{\rm off}-\beta_{\rm on}}{\gamma}\bigg{[}r(t_{\rm off})-r(t_{\rm on})\bigg{]}\bigg{\}},

(S48)

which is Eq. (6) in the main text.

S5 Appendix: Features of the final size equation

Let us study the features of Eq. (6). Note that we assume that $R_{0,\rm on}$ is fixed in this section.

First, we show that the final size monotonically decreases with respect to $r(t_{\rm off})$ for a fixed $r(t_{\rm on})$ . Regareing Eq. (6) as an implicit function of $r(\infty)$ and $\Delta r$ , we obtain

\displaystyle\bigg{\{}1-\frac{\beta_{\rm off}}{\gamma}[1-r(\infty)]\bigg{\}}\frac{\partial r(\infty)}{\partial(\Delta r)}

\displaystyle=

\displaystyle-\frac{\beta_{\rm off}-\beta_{\rm on}}{\gamma}[1-r(\infty)].

(S50)

In the final state, the inequality $0\leq 1-r(\infty)\leq\gamma/\beta_{\rm off}$ must hold because of the herd immunity condition. Therefore,

\displaystyle\frac{\partial r(\infty)}{\partial(\Delta r)}

\displaystyle\leq

\displaystyle 0,

(S51)

holds. As $\Delta r$ monotonically increases with respect to $t_{\rm off}$ for a fixed $t_{\rm on}$ , the final size is smaller if the intervention is longer. Equivalently, $r(\infty)$ monotonically increases with respect to $r(t_{\rm off})$ for a fixed $r(t_{\rm on})$ .

Note that there exists the upper bound in $r(t_{\rm off})$ . Let this upper bound be $\tilde{r}(r(t_{\rm on}))$ . It satisfies Eq. (7) in the main text:

\displaystyle\tilde{r}

\displaystyle=

\displaystyle 1-\exp\left[-\frac{\beta_{\rm off}-\beta_{\rm on}}{\gamma}r(t_{\rm on})\right]\exp\left[-\frac{\beta_{\rm on}}{\gamma}\tilde{r}\right],

(S52)

which is obtained for $\tilde{r}=r(\infty)$ in Eq. (6). This upper bound is observed in Figs. 1(F) and 2(F) in the main text.

Finally, let us discuss $r(t_{\rm on})$ and $r(t_{\rm off})$ which minimizes the final size with fixed $R_{0,\rm on}$ . The monotonicity of $r(\infty)$ with respect to $r(t_{\rm off})$ implies that the lower bound of the final size $r(\infty)$ with fixed $r(t_{\rm on})$ is obtained if $r(t_{\rm off})$ is replaced with $\tilde{r}(r(t_{\rm off}))$ . Therefore, we focus on the final size equation in which $r(t_{\rm off})$ is replaced with $\tilde{r}$ ,

\displaystyle r(\infty)

\displaystyle=

\displaystyle 1-\exp\bigg{[}-\frac{\beta_{\rm off}}{\gamma}r(\infty)\bigg{]}\exp\bigg{\{}\frac{\beta_{\rm off}-\beta_{\rm on}}{\gamma}\bigg{[}\tilde{r}(r(t_{\rm on}))-r(t_{\rm on})\bigg{]}\bigg{\}}.

(S54)

To minimize $r(\infty)$ in this equation, we should maximize $\tilde{r}(r(t_{\rm on}))-r(t_{\rm on})$ with respect to $r(t_{\rm on})$ . This condition is given by

	$\displaystyle\frac{\partial}{\partial r(t_{\rm on})}[\tilde{r}-r(t_{\rm on})]$	$\displaystyle=$	$\displaystyle\frac{(\beta_{\rm off}-\beta_{\rm on})(1-\tilde{r})}{\gamma-\beta_{\rm on}(1-\tilde{r})}-1$		(S55)
		$\displaystyle=$	$\displaystyle 0.$		(S56)

This can be satisfied for

	$\displaystyle\tilde{r}$	$\displaystyle=$	$\displaystyle 1-\frac{\gamma}{\beta_{\rm off}}=1-\frac{1}{R_{0,\rm off}},$		(S57)
	$\displaystyle r(t_{\rm on})$	$\displaystyle=$	$\displaystyle\frac{\gamma}{\beta_{\rm off}-\beta_{\rm on}}\left[\log\left(\frac{\beta_{\rm off}}{\gamma}\right)-\beta_{\rm on}\left(\frac{\beta_{\rm off}-\gamma}{\gamma\beta_{\rm off}}\right)\right].$		(S58)

By substituting these equations into Eq. (S54), we obtain the final size equation:

\displaystyle r(\infty)

\displaystyle=

\displaystyle 1-\frac{\gamma}{\beta_{\rm off}}\exp\bigg{[}\frac{\beta_{\rm off}}{\gamma}-1\bigg{]}\exp\bigg{[}-\frac{\beta_{\rm off}}{\gamma}r(\infty)\bigg{]}.

(S59)

Interestingly, this transcendental equation has the solution $r(\infty)=\tilde{r}=1-\frac{1}{R_{0,\rm off}}$ if $R_{0,\rm on}<R_{0,\rm off}\log\left(R_{0,\rm off}\right)/(R_{0,\rm off}-1)$ (Eq. (8) in the main text). In other words, one can achieve herd immunity with the theoretically minimal removed fraction by intervening once if the intervention reproduction number is sufficiently small.

S6 Appendix: Numerical results for the final size

Figures 1(B)–(E) and 2(B)–(E) in the main text present the numerical results of the final size $r(\infty)$ by varying the intervention onset time $t_{\rm on}$ with the constant basic reproduction number without the intervention $R_{0,{\rm off}}=2$ and intervention duration $\Delta t=60$ days. As suggested theoretically in S5 Appendix, the final size with the interventions (Figs. 1(B)–(E), 2(B)–(E) in the main text) is smaller than that without the intervention (Fig. 1(A) and 2(A) in the main text).

The final size is minimized in intermediate $t_{\rm on}$ in both cases (Figs. 1(D) and 2(D) in the main text). This result is intuitively understood as follows. The second wave occurs if the intervention onset is early (Figs. 1(B) and 2(B) in the main text). Conversely, the removed fraction is sufficiently enough before the intervention and the second wave is not observed if the onset is late (Figs. 1(E) and 2(E) in the main text).

For sufficiently large $\Delta t$ and small enough $R_{0,\rm on}$ , the final size reaches close to $1-\frac{1}{R_{0,\rm off}}=0.5$ , which is the lower bound of the final size necessary to achieve herd immunity, for $r(t_{\rm on})\approx 0.295$ (Fig. 2(D) in the main text), which corresponds to $t_{\rm on}=42.681$ days. It should be noted that a maximum of $\Delta r=r(t_{\rm off})-r(t_{\rm on})$ can be found for the intermediate $r(t_{\rm on})$ (Figs. S4 (left) and S5 (left)). Equation (6) in the main text (Eq. (6)) suggests that a larger $\Delta r$ yields a smaller $r(\infty)$ . Because $\tilde{r}$ is largest at intermediate $t_{\rm on}$ , the final size is smallest at intermediate $t_{\rm on}$ (see S5 Appendix).

As suggested by Eq. (6), the final size depends on $\Delta r$ but is independent of $r(t_{\rm on})$ (Figs. S4(A), S6(B), S5(A), and S7(B)).

Weak intervention (large $R_{0,\rm on}=1.4$ )

The final size can be minimized with $t_{\rm on}=0$ , the early implementation of the intervention, in this case (Fig. S4). The timing of the intervention that minimizes the final size is different from that which minimizes the maximum fraction of infected individuals. In general, the optimal timing of the intervention depends on the objective function to be minimized. The theoretical prediction in Eq. (6) (Eq. (6)), that the final size depends on $\Delta r$ , is numerically verified (Fig. S4(A)). For a constant intervention duration $\Delta t$ (Fig. S4(B)), the final size is minimized with intermediate $t_{\rm on}$ (Fig. 1(D) in the main text).

Strong intervention (small $R_{0,\rm on}=0.7$ )

The simulation results (Fig. 2(B)–(E) in the main text) suggest that the final size is minimized with an intermediate starting time for small $R_{0,\rm on}$ . Indeed, as shown in Fig. S5, the final size can be minimized at $t_{\rm on}\neq 0$ . The minimized final size is closer to $r(\infty)=1-\frac{1}{R_{0,\rm off}}=0.5$ , which is the minimum fraction of removed individuals necessary to achieve herd immunity, compared with the case of $R_{0,\rm on}=1.4$ (Fig. 2(D) in the main text).

S7 Figure: Contours of $i_{\rm max}$ and $r(\infty)$

As mentioned in the main text, the maximum fraction of infected individuals $i_{\rm max}$ depends linearly on $\Delta r$ if it appears after the intervention (case (i)). If the maximum appears during or at the onset of the intervention (cases (ii) and (iii)), $i_{\rm max}$ depends on $r(t_{\rm on})$ and does not depend on $\Delta r$ . Therefore, the contours in the $(r(t_{\rm on}),\Delta r)$ plane are the horizontal and vertical lines for the former and latter cases, respectively. To visualize them explicitly, Figs. S6(A) and S7(A) plots the five contours for $i_{\rm max}$ . As a theoretical prediction, they consist of horizontal and vertical lines.

Equation (6) theoretically predicts that the final size $r(\infty)$ with the intervention depends on $\Delta r$ only, which implies that the contours are horizontal lines (Figs. S6(B) and S7(B)).

S8 Appendix: Dependence of $\bar{i}_{\rm max}$ on $R_{0,\rm on}$

We want to minimize the maximum fraction of infected individuals $\bar{i}_{\rm max}$ for a fixed $R_{0,\rm on}$ by varying $t_{\rm on}$ and $t_{\rm off}$ . Note that $i_{\rm max}$ is decreasing function in $\Delta r$ (Eq. (S19)) for case (i), and $i_{\rm max}$ does not depend on $\Delta r$ in regions (ii) and (iii). Therefore, for a fixed $r(t_{\rm on})$ , $\bar{i}_{\rm max}$ is achieved either in region (ii) or (iii). Here, $\bar{i}_{\rm max}$ can be analytically or semi-analytically calculated by evaluating $i_{\rm max}$ on the boundaries between (i) and (ii) or between (i) and (iii). We show below that the dependence of $\bar{i}_{\rm max}$ on $R_{0,\rm on}$ changes, in which the boundary gives $\bar{i}_{\rm max}$ , leading to its non-trivial dependence on $R_{0,\rm on}$ .

Weak intervention ( $R_{0,\rm on}\geq R_{0,\rm on}^{\ast}\approx 1.23$ )

For large $R_{0,\rm on}$ , $\bar{i}_{\rm max}$ is achieved at the boundary between regions (i) and (ii), where the peaks of the infected individuals during and after the intervention are equal (Fig. 1(F) in the main text). To discuss the parameter regions where this boundary exists, it is necessary to take into account the upper bound of removed individuals $\tilde{r}(r(t_{\rm on}))$ . If $r(t_{\rm off})=r(t_{\rm on})+\Delta r$ at the boundary between regions (i) and (ii) given by Eq. (S41) is larger than $\tilde{r}$ , this boundary is not observed at $t_{\rm on}$ . In particular, it should be verified whether $\tilde{r}(0)$ is less than $\Delta r$ in Eq. (S41) for $r(t_{\rm on})=0$ . Let $R_{0,{\rm on}}^{\ast}$ be the value of $R_{0,{\rm on}}$ , below which the boundary between (i) and (ii) does not exist. For $R_{0,{\rm on}}\leq R_{0,{\rm on}}^{\ast}$ , the boundary crosses $\tilde{r}$ at a non-zero intervention onset time $r(t_{\rm on})$ , because the boundary is decreasing in $r(t_{\rm on})$ (Eq. (S41)) and $\tilde{r}$ is increasing in $r(t_{\rm on})$ for small $r(t_{\rm on})$ .

For $R_{0,\rm on}\geq R_{0,\rm on}^{\ast}$ , where there is a boundary between regions (i) and (ii) at $t_{\rm on}=0$ (Fig. 1(F) in the main text), $\Delta r$ can be the largest at this boundary at $t_{\rm on}=0$ . Equation (S17) suggests that a larger $\Delta r$ yields a smaller $i_{\rm max}$ , and $\bar{i}_{\rm max}$ is given by substituting $r(t_{\rm on})=0$ into Eq. (S30) as follows:

\displaystyle\bar{i}_{\rm max}

\displaystyle=

\displaystyle 1-\frac{1}{R_{0,\rm on}}\bigg{[}1+\log\left(R_{0,\rm on}\right)\bigg{]}.

(S60)

This equation is used as the theoretical curve for $R_{0,\rm on}\geq R_{0,\rm on}^{\ast}$ in Fig. 3 in the main text. This equation implies that $\bar{i}_{\rm max}$ is increasing in $R_{0,\rm on}$ for $R_{0,\rm on}\geq R_{0,\rm on}^{\ast}$ . See the next subsection for the value of $R_{0,\rm on}^{\ast}$ .

Medium-intervention ( $R_{0,\rm on}^{\ast}\geq R_{0,\rm on}\geq R_{0,\rm on}^{(1)}\approx 1.08$ )

If $R_{0,\rm on}$ is smaller than $R_{0,\rm on}^{\ast}$ , region (ii) does not give the global maximum of the fraction of infected individuals for $t_{\rm on}=0$ . As a result, a larger $\Delta r$ , which indicates a smaller $i_{\rm max}$ , is realized for $t_{\rm on}\neq 0$ . The upper bound of $r(t_{\rm off})$ at the boundary between regions (i) and (ii) equals $\tilde{r}(r(t_{\rm on}))$ . Substituting the condition for the boundary between (i) and (ii) (Eq. (S41)) in the limit $r(t_{\rm off})\rightarrow\tilde{r}$ into the condition for $\tilde{r}$ (Eq. (7)), we obtain

\displaystyle\tilde{r}

\displaystyle=

\displaystyle 1-e^{-A},

(S61)

where

\displaystyle A

\displaystyle:=

\displaystyle 1+\frac{1}{R_{0,\rm off}-R_{0,\rm on}}\big{[}R_{0,\rm off}\log\left(R_{0,\rm on}\right)-R_{0,\rm on}\log\left(R_{0,\rm off}\right)\big{]}.

(S62)

Substituting this into Eq. (S41) again, $r(t_{\rm on})$ at the boundary can be explicitly given as

\displaystyle r(t_{\rm on})

\displaystyle=

\displaystyle\frac{R_{0,\rm on}}{R_{0,\rm on}-R_{0,\rm off}}\bigg{[}(1-e^{-A})-\frac{1}{R_{0,\rm on}}A\bigg{]}.

(S63)

Substituting these into Eq. (S17) or (S18), we obtain

\displaystyle\bar{i}_{\rm max}

\displaystyle=

\displaystyle e^{-A}-\frac{1}{R_{0,\rm off}-R_{0,\rm on}}\log\left(\frac{R_{0,\rm off}}{R_{0,\rm on}}\right),

(S64)

which is the theoretical curve for $R_{0,\rm on}^{\ast}\geq R_{0,\rm on}\geq R_{0,\rm on}^{(1)}$ in Fig. 3 in the main text. Note that $\bar{i}_{\rm max}$ is a decreasing function in $R_{0,\rm on}$ . Therefore, $\bar{i}_{\rm max}$ is minimized at $R_{0,\rm on}=R_{0,\rm on}^{\ast}$ . The condition for $R_{0,\rm on}^{\ast}$ is obtained by equating Eqs. (S60) and (S64). It is difficult to calculate $R_{0,\rm on}^{\ast}$ explicitly, and it was numerically obtained as $R_{0,\rm on}^{\ast}\approx 1.23$ in the present parameter setting.

Strong intervention ( $R_{0,\rm on}\leq R_{0,\rm on}^{(1)}$ )

For small $R_{0,\rm on}$ , $\bar{i}_{\rm max}$ is determined differently. It is understood that a peak in region (ii) does not give $i_{\rm max}$ for $R_{0,\rm on}\leq R_{0,\rm on}^{(1)}$ , where $R_{0,\rm on}^{(1)}\approx 1.08$ .

In this case, the boundary between regions (i) and (iii) satisfying $r(t_{\rm off})=\tilde{r}(r(t_{\rm on}))$ gives the condition for $\bar{i}_{\rm max}$ (Fig. 2(F) in the main text). This condition cannot be solved explicitly with respect to $\bar{i}_{\rm max}$ . Therefore, we give the parametric equations for $\bar{i}_{\rm max}$ and $R_{0,\rm on}$ in terms of $r(t_{\rm on})$ . The parametric equation for $\bar{i}_{\rm max}$ is

\displaystyle\bar{i}_{\rm max}

\displaystyle=

\displaystyle 1-e^{-\frac{\beta_{\rm off}}{\gamma}r(t_{\rm on})}-r(t_{\rm on}).

(S65)

Next, we derive the parametric equation for $R_{0,\rm on}$ . To this end, let us rewrite (LABEL:eq:boundary_i_iii) as

\displaystyle\Delta r

\displaystyle=

\displaystyle\frac{B(r(t_{\rm on}))}{R_{0,\rm off}-R_{0,\rm on}},

(S66)

where

\displaystyle B(r(t_{\rm on}))

\displaystyle=

\displaystyle R_{0,\rm off}[e^{-\frac{\beta_{\rm off}}{\gamma}r(t_{\rm on})}+r(t_{\rm on})]-\bigg{[}1+\log\left(R_{0,\rm off}\right)\bigg{]}.

(S67)

Substituting this into Eq. (7), we obtain a transcendental equation for $\tilde{r}$ :

\displaystyle\tilde{r}

\displaystyle=

\displaystyle 1-\exp\left[-R_{0,\rm off}\tilde{r}\right]\exp\left[B(r(t_{\rm on}))\right],

(S68)

which can be solved as

\displaystyle\tilde{r}

\displaystyle=

\displaystyle 1+\frac{1}{R_{0,\rm off}}W_{-1}\bigg{(}-R_{0,\rm off}\exp\left[B(r(t_{\rm on}))-R_{0,\rm off}\right]\bigg{)},

(S69)

where $W_{-1}(\cdot)$ denotes the branch of the Lambert $W$ function which gives the real value other than the principal branch. Finally, we obtain the parametric equation for $R_{0,\rm on}$ as

	$\displaystyle R_{0,\rm on}=R_{0,\rm off}-\frac{B(r(t_{\rm on}))}{1+\frac{1}{R_{0,\rm off}}W_{-1}\left(-R_{0,\rm off}\exp\left[B(r(t_{\rm on}))-R_{0,\rm off}\right]\right)-r(t_{\rm on})}.$
			(S70)

Equations (S65) and (S70) are used to plot the theoretical curves for $R_{0,\rm on}\leq R_{0,\rm on}^{(1)}$ in Fig. 3 in the main text.

S9 Appendix: Optimization problems with respect to final size

Two scenarios for the optimization problems in the hybrid nonlinear dynamical systems are presented in the main text. Here, we discuss the problems using the final size as the objective functions.

Minimizing final size with a constraint in intervention duration

Targeting the final size to be minimized, we discuss two constraints based on Figs. S4 and S5. As in the case for the minimization problem of the maximum fraction of infected individuals, the final size is minimized at an intermediate $t_{\rm on}$ with the constraint of a constant $\Delta t$ (Figs. 1(D) and 2(D) in the main text, and Figs. S4 (right) and S5 (right)). Note that this $t_{\rm on}$ is different from the one which minimizes $i_{\rm max}$ . This result is intuitively understood in the right panels of these figures. The contour of the final size depends on $\Delta r$ . For a constant intervention duration, $\Delta r$ is largest with intermediate $r(t_{\rm on})$ , which maximizes the final size with this constraint.

Minimizing intervention duration with a constraint on the final size

It is possible to consider this rule to minimize the intervention duration along the contour of the final size in Figs. S4 and S5. This corresponds to intervening while keeping $\Delta r$ constant. There exists an intermediate $t_{\rm on}$ that minimizes the final size in this setting as well.

Analytical estimation of maximum fraction of infected individuals with one-shot non-pharmaceutical intervention in a hybrid epidemic model

Abstract

keywords:

Introduction

Materials and methods

Model

Non-pharmaceutical intervention

Analysis

Maximum fraction of infected individuals

Final size with intervention

Results

Weak intervention (large R0,on=1.4R_{0,\rm on}=1.4)

Strong intervention (small R0,on=0.7R_{0,\rm on}=0.7)

Maximum infected fraction i¯max\bar{i}_{\rm max} versus reproduction number under the intervention R0,onR_{0,\rm on}

Intervention strategies

Minimizing imaxi_{\rm max} with a constraint in the intervention duration

Minimizing intervention duration with a constraint in imaxi_{\rm max}

Conclusion

Competing interests

Author’s contributions

Acknowledgements

References

Additional Files

Additional file 1 — Supplementary texts

Supplementary materials: Analytical estimation of maximum fraction of infected individuals with one-shot non-pharmaceutical intervention in a hybrid epidemic model

S1 Method: Derivation of Eq. (4) [6, 7, 8, 9]

S2 Method: Analytical expressions for the maximum fraction of infected individuals

Derivation

(i) Maximum fraction after the intervention

(ii) Maximum fraction during the intervention

(iii) Maximum fraction at the onset of the intervention

(iv) Maximum fraction before the intervention

S3 Method: Derivation of boundaries between regions of different maxima

Boundaries between regions (i) and (ii)

Boundaries between (i) and (iii)

Boundaries between (ii) and (iii)

Boundaries between (iii) and (iv)

S4 Method: Derivation of the final size equation with the intervention

S5 Appendix: Features of the final size equation

S6 Appendix: Numerical results for the final size

Weak intervention (large R0,on=1.4R_{0,\rm on}=1.4)

Strong intervention (small R0,on=0.7R_{0,\rm on}=0.7)

S7 Figure: Contours of imaxi_{\rm max} and r​(∞)r(\infty)

S8 Appendix: Dependence of i¯max\bar{i}_{\rm max} on R0,onR_{0,\rm on}

Weak intervention (R0,on≥R0,on∗≈1.23R_{0,\rm on}\geq R_{0,\rm on}^{\ast}\approx 1.23)

Medium-intervention (R0,on∗≥R0,on≥R0,on(1)≈1.08R_{0,\rm on}^{\ast}\geq R_{0,\rm on}\geq R_{0,\rm on}^{(1)}\approx 1.08)

Strong intervention (R0,on≤R0,on(1)R_{0,\rm on}\leq R_{0,\rm on}^{(1)})

S9 Appendix: Optimization problems with respect to final size

Minimizing final size with a constraint in intervention duration

Minimizing intervention duration with a constraint on the final size

Weak intervention (large $R_{0,\rm on}=1.4$ )

Strong intervention (small $R_{0,\rm on}=0.7$ )

Maximum infected fraction $\bar{i}_{\rm max}$ versus reproduction number under the intervention $R_{0,\rm on}$

Minimizing $i_{\rm max}$ with a constraint in the intervention duration

Minimizing intervention duration with a constraint in $i_{\rm max}$

Weak intervention (large $R_{0,\rm on}=1.4$ )

Strong intervention (small $R_{0,\rm on}=0.7$ )

S7 Figure: Contours of $i_{\rm max}$ and $r(\infty)$

S8 Appendix: Dependence of $\bar{i}_{\rm max}$ on $R_{0,\rm on}$

Weak intervention ( $R_{0,\rm on}\geq R_{0,\rm on}^{\ast}\approx 1.23$ )

Medium-intervention ( $R_{0,\rm on}^{\ast}\geq R_{0,\rm on}\geq R_{0,\rm on}^{(1)}\approx 1.08$ )

Strong intervention ( $R_{0,\rm on}\leq R_{0,\rm on}^{(1)}$ )