Minimizing the epidemic final size while containing the infected peak prevalence in SIR systems.

Mathematical models are critical to understand the spread of pathogens in a population and evaluate the effectiveness of non-pharmaceutical interventions (NPIs). A plethora of optimal strategies has been recently developed to minimize either the infected peak prevalence ( I P P ) or the epidemic final size ( E F S ). While most of them optimize a simple cost function along a fixed finite-time horizon, no consensus has been reached about how to simultaneously handle the I P P and the E F S , while minimizing the intervention's side effects. In this work, based on a new characterization of the dynamical behaviour of SIR-type models under control actions (including the stability of equilibrium sets in terms of herd immunity), we study how to minimize the E F S while keeping the I P P controlled at any time. A procedure is proposed to tailor NPIs by separating transient from stationary control objectives: the potential benefits of the strategy are illustrated by a detailed analysis and simulation results related to the COVID-19 pandemic.

Introduction

Infectious disease outbreaks are latent threats to humankind — killing annually millions worldwide (Abumalloh et al., 2021, Hernandez-Vargas et al., 2019). During an outbreak, public health agencies aim to limit the spread of pathogen by non-pharmaceutical interventions (NPIs), including the implementation of lockdowns of varying intensity and geographic scope (Ferguson et al., 2020). While an effective vaccine is developed to counter the pathogens and the new variants that could emerge during the pandemic, new waves of infections may still take place, saturating public health capacities (Contreras & Priesemann, 2021).

A crucial aspect for policymakers during pandemics is to set up and remove intervention measures while avoiding the collapse of public health capacities and the economy. From a theoretical control perspective, this problem falls into the classic framework of optimal control (Lewis, Vrabie, & Syrmos, 2012). The major aim of this control strategy is to obtain the best possible performance by using the least control actions, in a kind of cause–effect balance. However, in the case of SIR-type systems based on the seminal work (Kermack & McKendrick, 1927), such a cause–effect separation is not so clear, and evidence can be collected showing that the use of direct and simple control objectives – i.e., to minimize the number of infected individuals at every time or to minimize the final total number of infected – largely produces sub-optimal solutions (Abbasi, 2020, Alamo et al., 2021, Hale et al., 2021, Punzo, 2022, Rypdal et al., 2020).

Several approaches have been proposed to find an optimal-control-based intervention for SIR models. The two main metrics to measure the disease impact are Di Lauro, Kiss, and Miller (2021): the infected peak prevalence, (maximal fraction of infected individuals along time), which is closely related to the health systems capacity, and the epidemic final size, (total final fraction infected). A first result – when minimizing either or – is that optimal solutions can be obtained for the ubiquitous single-interval intervention (Sadeghi, Greene, & Sontag, 2021), i.e., a fixed reduction of the reproduction number (by reducing the infection rate) for a given period of time , with . In Morris, Rossine, Plotkin, and Levin (2021) and Sadeghi et al. (2021) rigorous analyses are made to show how to find the optimal single-interval control action that minimizes . The focus is put on the optimal intervention starting and finishing times and, as stated in Morris et al. (2021), even when theoretical optimal (and near-optimal) interventions are found, they are not resilient to errors in timing. The main problem of this strategy is, however, that it does not account for the other severity index, the . Thus, the total number of infected (and the total number of deceases) is far to be minimized.

Similar approaches (and results) are presented in Bliman and Duprez, 2021, Di Lauro et al., 2021 and Ketcheson (2021), but minimizing the . In these studies, it is found that an optimal intervention also exists, but in a rather unimplementable context. According to Bliman and Duprez (2021), the best policy consists in leaving the system in open-loop until the susceptible/non-infected fraction of individuals approaches the herd immunity threshold and then, at a particular time, implementing the hardest possible intervention. As it can be inferred, this strategy has two main drawbacks: any small error in the timing produces a performance drastically different from the optimal one and, more importantly, the is unacceptably large (since the system is left in open-loop for a long time before acting). The latter point is mentioned in Di Lauro et al. (2021), where it is said that ‘to minimize the total number of infected, the intervention should start close to the peak’. The former point, on the other hand, is demonstrated through simulations in Ketcheson (2021), where a slightly different intervention from the optimal produces severe sub-optimal results. Another (practical) optimal control approach that led to a more realistic scheduling can be found in Köhler et al. (2021), where a model predictive control (MPC) is proposed based on the SIDARTHE (Giordano et al., 2020), and the control objective consists in minimizing the current number of infected individuals (and fatalities) and the time of isolation.

In any case, the common factor in all the recent literature – at the best of the author’s knowledge – is that no conclusive results are shown concerning which is the best policy to simultaneously minimize the , the , and intervention’s severity. In this article, we show that the key point to achieve – or, at least, to arbitrarily approximate – such a goal is the way the optimization problem is posed, by properly separating transient and stationary regimes. Based on a new set-based dynamic analysis of the SIR-type models, a different perspective to formulate the optimal control problem is presented. Instead of considering the control objective of minimizing the or the , the susceptible are directly steered to the (open-loop) herd immunity, since this threshold represents the minimal at steady-state, for any finite-time intervention. Furthermore, since it is independent of the , the is maintained under an upper bound (computed to cope with the health system capacity), while the only quantity to be minimized is the strength and time of the NPI. As demonstrated by several simulation results (related to the spread of COVID-19 in France, during 2020, Bliman and Duprez (2021)), this strategy seems to be general enough to provide a confidence baseline to policymakers in the critical task of decision-making in a pandemic context.

Review of control SIR model

The SIR epidemic model (Kermack & McKendrick, 1927) describes the fractions of susceptible and infectious individuals in a population, at time . New infections occur proportional to at a transmission rate , and infectious individuals recover or die at a rate . NPIs reduce the effective transmission rate, , below its value in the absence of intervention (which is considered fixed). By rescaling the time by , the SIR model can be written in non-dimensional form as (Bertozzi, Franco, Mohler, Short, & Sledge, 2020): where denotes the time-varying reproduction number fulfilling , with , being the starting and ending intervention time ( is assumed to be finite since social intervention has always an end), and the minimal and maximal values for the reproduction number, respectively ( and correspond to non-intervention and maximal intervention, respectively; the case is not considered, since a perfect full lockdown is not possible).

Susceptible and infectious are positive and constrained to be in the set , for all . Particularly, denoting the epidemic outbreak time, it is assumed that , with ; i.e., the fraction of susceptible individuals is smaller than, but close to , and the fraction of infectious is close to zero at .

No-intervention dynamical analysis

Assume first that , for , which represents the no-intervention (or open-loop) scenario. The solution of (1) for – which was analytically determined in Harko, Lobo, and Mak (2014) – depends on and the initial conditions . is a decreasing function for all , while is decreasing for all if and, if , initially increases, then reaches a global maximum, and finally decreases to zero. In this latter case, the peak of , , is reached at , when , and depends on initial conditions and : Condition implies that , where is the ‘herd immunity’ (i.e., the value of under which cannot longer increase).

Define now and , which depend on , and . By taking for the solutions proposed in Harko et al. (2014), we obtain . Furthermore, following a similar procedure to Abuin, Anderson, Ferramosca, Hernandez-Vargas, and Gonzalez (2020) for in-host models, is given by where is the Lambert function (Pakes, 2015). The , which as the , is a function of and , is then given by

The following lemma states the maximum of over .

Consider system (1) with initial conditions , for some , and fixed. Then, for , with , the maximum of occurs at and is given by (particularly, by , if ).

See the Appendix.

Next, some properties of , for different values of its arguments, are given.

Consider system (1) with arbitrary initial conditions , for some and . Then: (i) and . (ii) For and fixed , , decreases with , and . This means that the closer is to from above, the closer will be to from below. (iii) For and fixed , , increases with , and . This means that the closer is to from below, the closer will be to , from below. (iv) For any fixed and , decreases with and . (v) . If and , , for any value of (note that for ).

The proof of Property 2.1 is omitted for brevity. Fig. 1 shows how behaves for different initial conditions.

Fig. 1

Function is bounded from above by (, light red plane). Furthermore, reaches its maximum, given by , at , . (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Equilibrium characterization and stability

The equilibrium set of system (1), with and initial conditions , is given by . Next, a key theorem concerning the asymptotic stability of a subset of is introduced.

Asymptotic Stability of

Consider system (1) with and constrained by . Then, the set with the herd immunity, is the unique asymptotically stable (AS) set of system (1) , with a domain of attraction given by .

See the Appendix.

A corollary of Theorem 1, concerning the properties of , is presented next.

Consider system (1) with arbitrary initial conditions , for some . Then: (i) Set is a subset of (for , ), and its size depends on , but not on the initial conditions. (ii) Subsets of are stable but not attractive (i.e., is AS as a whole, but no subset of it is AS). This is particularly true for the state . (iii) If , . Then, and the so called healthy equilibrium with , lies in , and so it is stable, but not attractive (any small makes the system to converge to ) with . (iv) If , set can be divided into two sets, , where is is AS, for any value of . However, if , it contains an unstable set, .

It is noteworthy that the key dynamic behaviours described in this section also hold for complex and realistic SIR-type models, as stated in Sadeghi et al. (2021).

Epidemic control

NPIs (such as social distancing, isolation measures, mask wearing, etc.) are the typical measures that policymakers implement to control epidemics when vaccination (effectiveness and distribution) is not enough, to lessens the disease transmission rate or, directly, parameter in system (1). Assuming now that , is no longer given by Eq. (2). However, as the final intervention time is finite, Eq. (4) still describes . The following lemma gives an upper bound for , for any .

Consider system (1) with initial conditions , , such that , and . Then, (i) the system converges to with , being the herd immunity corresponding to no intervention and, (ii), the only way to achieve is with an producing and , i.e., the system achieving a quasi steady-state (QSS), at . The cases in which does not approach are: (a) if and , then a and, finally, the system will converge to (the greater is with respect to , the smaller will be ), (b) if and , then will be close to (the smaller is with respect to , the smaller will be ), and (c) if does not approach zero (i.e., no QSS conditions is reached at ), then no matter which value takes, will be smaller than (the farther is from from above or from below, the smaller will be ).

See the Appendix.

Lemma 3.1 is a simple but strong result concerning any kind of NPIs, interrupted at a finite time. On one hand, it establishes that the minimal is completely determined by the epidemic itself (the value of ) and, provided that no immunization (by vaccination and/or development of the individual’s immune system) is considered, it cannot be modified by NPIs. On the other hand, Lemma 3.1 establishes that must be reached as a QSS condition (i.e., with ), since otherwise will decrease after the NPI is interrupted at , thus producing .

To show the latter properties, consider system (1) with ( days−1 and days−1), , and (COVID-19 circulation in France, March–May, 2020 Bliman & Duprez, 2021). Fig. 2 shows the phase portraits of four different controls . First, long term strategies (solid lines) are implemented (one strong and the other soft) at some time , being and , up to a large time , such that the system reaches a QSS. In the first case (blue line) and and given that is significantly greater than and is still positive, a second wave occurs that steers to a small value (). In the second case (red line), and and since is already under (it cannot grow for ), so is again significantly small (). In a second stage, short terms strategies (dashed lines) are considered (one strong and the other soft). In the first case (blue dashed line), and , which leads to , with significantly small. In the second case (red dashed line), and , which leads, again, to , showing that any control interrupted before the QSS produces a small .

Fig. 2

Long term (solid line) and short term (dashed line) controls. Initial state: empty circle, : filled circles, Lyapunov function level curves (Appendix): green lines. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Redefining the control objective

Now, the idea is to find an that arbitrarily reduces the while maintaining . To have a first insight into the answer, consider the integral of equation ((1).b), , where is a constant determined by the initial values . Then, given that , with , it follows that . Now, taking the limits for , and recalling that for , it follows that . This equality means that, even when varies over time, only determines the area under the curve of , , but not its peak . In other words, it is possible to minimize the and also keep the under a maximal value imposed by the health system capacity, as long as the condition is respected. This separation of steady-state (minimize the ) and transient (keep the low) objectives can be then exploited to avoid unnecessary competitions between them, which would produces suboptimal solutions, as done in Bliman and Duprez, 2021, Di Lauro et al., 2021, Ketcheson, 2021 and Köhler et al., 2021, Morris et al., 2021. Given that the side effects of the NPIs should be minimized as well, we propose two different control objectives, a primary and a secondary one.

Epidemiological and Social/Economic Control Objectives

Consider system (1) with initial conditions , , such that , and . Consider also that a maximal value for , , is established according to the health system capacity. Then, the epidemiological control objective (ECO) consists in steering to , as , while maintaining , for all . Furthermore, the social/economic control objective (SECO) consists in minimizing the Social Distancing Index , provided that the ECO is achieved.

Single-interval intervention

We want to find first the simplest that fulfils the ECO of Definition 3.1. A single-interval intervention, , is defined as , for , and , for , where is a fixed value of intervention. Denote by the value of producing , when applied at time , for a large enough . By making , we obtain , which is a decreasing function of . Now, denote by the value of guaranteeing that , for all . By making , we (implicitly) obtain , which is also a decreasing function of . Finally, by merging the latter conditions, it is possible to define a (so-called) goldilocks intervention:

Goldilocks Single-interval Intervention

The goldilocks single-interval intervention is defined by a starting time, , fulfilling condition , and the fixed reproduction number value, .

The goldilocks single-interval intervention allows us to establish the following theorem:

Consider system (1) with initial conditions , , such that , and . Consider a given . Then, if for and there exists a goldilocks single-interval intervention, it is the only one that arbitrarily approaches the ECO, as .

See the Appendix.

We now resume the example of the previous section to evaluate the goldilocks single-interval intervention. Considering and , we obtain days and . Fig. 3 (left), shows (upper plot, blue solid line), (upper plot, red solid line) and (lower plot, solid line) for a period of time of days. Fig. 3 (right) shows the corresponding phase portrait (blue solid line), and the level curves of Lyapunov functions, (green curves) and (red curves), corresponding to and , respectively. The goldilocks strategy picks the red curve ‘guiding’ the system exactly to , while for . The performance indexes of the strategy are: , , .

Fig. 3

and time evolution (left), and phase portrait, (right). Goldilocks single-interval (solid line) and ‘wait, maintain, suspend’ interventions (dashed line). System with (black dashed line) and (black dotted–dashed line). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

The existence of the goldilocks single-interval intervention depends on and , and is not analysed here, for the sake of brevity. In any case, goldilocks interventions should be understood just as a first-step approach, since they can hardly be applied in realistic cases.

‘Wait, maintain, suspend’ strategy

Another strategy that accounts for the ECO of Definition 3.1 and avoids the problem of the existence of a solution for any is the ‘wait, maintain, suspend’ (the intervention) strategy (Morris et al., 2021): , for , , for , and for , where , is the time at which a threshold condition (specified later on) is reached, and is the fixed intervention that, if started at (and applied for a large enough), produces .

Time is considered now as the time at which the open-loop system reaches . At this time, the control action is applied to system (1), making constant for the period . As a result, decreases linearly for (since ). Now, if time is not large enough, may increase for , reaching a peak that overpasses , violating this way the control objective . On the other hand, if is too large, may decrease under , violating the control objective . The next theorem establishes that does exists, such that the ‘wait, maintain, suspend’ strategy fulfils the ECO.

Consider system (1) with initial conditions , , such that , and . Consider also a given . Then, there exist some such that the ‘wait, maintain, suspend’ strategy produces , for , that arbitrarily approach the ECO, as .

See the Appendix.

Fig. 3 (left), shows (upper plot, blue dashed line), (upper plot, red dashed line) and (lower plot, dashed line) for a period of time of days, of the ‘wait, maintain, suspend’ intervention. Fig. 3 (right) shows the corresponding phase portrait (blue dashed line). The times are given by , days, and and the ECO is reached, since and for . The performance indexes of this strategy are: , , .

As before, this intervention strategy is rather unrealistic, since the control action varies continuously in the interval and the severity of the intervention is extremely high.

Optimal control strategy

Now, by taking advantage of the analysis of the previous strategies the following optimal control problem, , is proposed, which accounts also for the SECO: where is a large enough (possibly infinite) time that covers the whole dynamic of the epidemic and is the set defined in Section 2. Conditions forces to be smaller than the externally imposed maximum at every time , while and , with being an arbitrary small detectable value, forces the system to reach a QSS at , with . The key point of Problem is that the ECO is imposed by constraints while the SECO is achieved by optimality, so the competition between them is avoided. The next theorem establishes that Problem is well-posed and achieves the ECO and the SECO.

Consider system (1) with initial conditions , , such that , and . Consider also a given . Then, the solution of Problem , denoted as , produces for that arbitrarily approaches the ECO, as , and minimizes the SDI (SECO).

See the Appendix.

Fig. 4 (left) shows (upper plot, blue line), (upper plot, red line) and (lower plot, blue line), of the optimal intervention , with days. A significant improvement is obtained in terms of the SECO: the SDI drop from in the previous strategies to , while the ECO is (practically) reached. As a particularity, separates the epidemiological objectives over time: first, it handles the (from to days) and, once cannot further increase, it tries to reach , at steady-state, to minimize . The performance indexes are: , , .

Fig. 4

and time evolution (left), and phase portrait (right). Optimal control . System with and . (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Although not simulated here, it can be shown that any other optimal control problem (i.e., the one minimizing a cost or ) will necessarily produce results far from the optimal one (for any combination of the weights and ), since in these cases the SECO and the ECO compete between them.

It is noteworthy that the optimal solution is an important theoretical result to settle down the need of more complex control strategies to account for real epidemics. The main practical limitations are: (i) Problem is in open loop (no state update is considered, as it is done in an MPC); (ii) the model is still a simplistic one; (iii) varies continuously over time while real NPIs are prescribed by taking a few possible values, each one applied for bounded periods; (iv) control actions affect in rather simplistic (linear) forms. Overall, is intended to elucidate some fundamental aspects of the optimization of SIR-type models, and not to formalize an applicable strategy.

Conclusions and future work

In this work, a new set-based equilibrium and stability characterization of SIR-type models was made. Based on this characterization, the and the were analysed, and concrete epidemiological objectives involving both indexes were proposed. It is shown that simple non-pharmaceutical intervention strategies exist that (theoretically) accomplish such objectives. If social/economic side effects of the NPIs are also considered – which leads to a non-trivial optimal control problem – a solution can be found that, in addition to achieving the epidemiological objective, minimizes the side effects. Future work will include the consideration of more complex and realistic models and the design of a proper feedback controller able to account for uncertainties/disturbances over short, updated, time horizons. MPC appears to be the right framework to account for such a challenge.