key: cord-0642434-fb2igfbk
authors: Bliman, Pierre-Alexandre; Duprez, Michel; Privat, Yannick; Vauchelet, Nicolas
title: Optimal immunity control by social distancing for the SIR epidemic model
date: 2020-06-10
journal: nan
DOI: nan
sha: 3b26f118f0379b1197e1db8a65a50018eb3f5aa9
doc_id: 642434
cord_uid: fb2igfbk

Until a vaccine or therapy is found against the SARS-CoV-2 coronavirus, reaching herd immunity appears to be the only mid-term option. However, if the number of infected individuals decreases and eventually fades only beyond this threshold, a significant proportion of susceptible may still be infected until the epidemic is over. A containment strategy is likely the best policy in the worst case where no vaccine or therapy is found. In order to keep the number of newly infected persons to a minimum, a possible strategy is to apply strict containment measures, so that the number of susceptible individuals remains close to herd immunity. Such an action is unrealistic since containment can only last for a finite amount of time and is never total. In this article, using a classical SIR model, we determine the (partial or total) containment strategy on a given finite time interval that maximizes the number of susceptible individuals over an infinite horizon, or equivalently that minimizes the total infection burden during the curse of the epidemic. The existence and uniqueness of the optimal strategy is proved and the latter is fully characterized. If applicable in practice, such a strategy would lead theoretically to an increase by 30% of the proportion of susceptible on an infinite horizon, for a containment level corresponding to the sanitary measures put in place in France from March to May 2020. We also analyze the minimum intervention time to reach a fixed distance from herd immunity, and show the relationship with the previous problem. Simulations are provided that illustrate and validate the theoretical results.

Important efforts have been devoted to the application of optimal control theory to the control of infectious diseases. Many works have considered the simple SIR model without demography that recently became even more famous in a broad community worldwide. Such a model readṡ

The state variables S, I, R correspond respectively to the proportions of susceptible, infected and removed individuals in the population. Notice that the sum of the derivatives of the three state variables is zero, so the sum of the variables remains equal to 1 if initially S 0 + I 0 + R 0 = 1. One may therefore describe the system solely with the equations (1a)-(1b). The total population is constant and no demographic effect (births, deaths) is modeled, as they are not relevant to the time scale to be taken into account in reacting to an outbreak. The positive parameter β in the infection term βSI accounts at the same time for the rate of encounters between the individuals and the probability of transmitting the infection during each of theses encounters. The positive parameter ν describes the rate at which the infected are removed (typically by recovery) by time unit. Cured people are supposed to have acquired permanent immunity. An important issue, abundantly discussed these days, is herd immunity. The latter occurs naturally when a large proportion of the population has become immune to the infection. Mathematically, it is defined as the value of S below which the number of infected decreases. For the SIR model (1), one hasİ = (βS − ν)I, and this number reads

(2)

While the number of infected decreases when S(·) S herd , epidemics continue to consume susceptible and to generate new infections once the immunity threshold has been crossed.

The basic reproduction number R 0 := β/ν = 1 S herd governs the behaviour of the system departing from its initial value. When R 0 < 1, no epidemic may occur, whereas when R 0 > 1 epidemic occurs if S 0 > 1 R0 . In such a case, I reaches a peak and then goes to zero due to the immunity effect previously presented. From now on, we assume S 0 > S herd and

Obviously, every solution of (1) is nonnegative and thus,Ṙ 0 Ṡ at any instant, so S may only decrease, while R may only increase. We infer that S remains bounded and therefore, the integral +∞ 0 βS(t)I(t) dt converges, and in any case the following limit exists: 

We remark that, after passing the collective immunity threshold, that is while t is large enough, one has S herd S(t) S ∞ .

In order to illustrate herd immunity, Table 1 displays, for several values of R 0 , the value of the herd immunity threshold and the number of susceptibles that remain after the fading out of the outbreak, in the case of an initially naive population (S 0 ≃ 1). The proportion of infections occurred after the overcome of the immunity is also shown.

Apart from medical treatment, there are generally speaking three main methods to control human diseases. Each of them, alone or in conjunction with the others, gave rise to applications of optimal control. We introduce a modified SIR system, where the control inputs u 1 (·), u 2 (·), u 3 (·) materialize these three methods.Ṡ = −u 3 βSI − u 1 , S(0) = S 0 , (5a)

R 0 1.5 2 2.5 2.9 3 3.5 S herd 0.67 0.50 0.40 0.34 0.33 0.29 S ∞ 0.42 0.20 0.11 0.067 0.059 0.034 S herd − S ∞ 1 − S ∞ 43% 37% 33% 30% 29% 27% Table 1 : Herd immunity level S herd and asymptotic susceptible proportion S ∞ for initial value S 0 ≃ 1 for several values of the basic reproduction number R 0 . The value of S herd comes from (2), and S ∞ may be deduced from the fact that S ∞ − 1 R0 ln S ∞ ≈ 1, see Lemma 2 below. The ratio

represents the proportion of susceptible that occur after passing the collective immunity threshold. The column in bold corresponds to R 0 found in [SKL + 20] for the SARS-CoV-2 in France before the lockdown of March-May 2020.

A first class of interventions consists in vaccination or immunization [Aba74, MW74, DB80, Gre88, Beh00, GS09, HD11, HZ13, AER + 14, YWZZ15, LT15, BBSG17, BBDMG19, Shi19]. It consists in transferring individuals from the S compartment to the R one: this is the effect of the additive input u 1 in (5). The members of the latter may not be anymore "recovered" stricto sensu, but the key point in the model is that they are excluded from the infective process. A second class of measures corresponds to screening and quarantining of infected [Aba73, Wic75, Wic79, Beh00, AI12, ZWW13, BBSG17, BBDMG19]. It may be modelled by transfer of individuals from the I compartment to the R one, as done by the input u 2 in (5). Last, it is possible to reduce transmission through health promotion campaigns or lockdown policies [Beh00, BBSG17, MRPL20] , as done by the multiplicative input u 3 in (5). Of course these methods may be employed jointly [HD11] .

Other modelling frameworks have also been considered. More involved models called SEIR and SIRS have been analyzed in [Beh00, GS09] , and SIR model structured by the infection age in [AI12] . Constraints on the number of infected persons that the public health system can accommodate or on available resources, particularly in terms of vaccination, were studied [ZWW13, YWZZ15, BBSG17, BBDMG19]. Economical considerations may be aggregated to the epidemiological model [KS20, PS20, AAL20]. Ad-hoc models for tackling emergence of resistance to drugs issues have been introduced [JDM13, JDHWM13], as well as framework allowing to study revaccination policies [KGZ15, KZ16] . Optimization of vaccination campaigns for vector-borne diseases have also been considered [Shi19] . The references cited above are limited to deterministic differential models, but discrete-time models and stochastic models have also been used e.g. in [Aba73, Aba74, Gre88] .

The costs considered in the literature are usually integral costs combining an "outbreak size" (the integral of the number of infected, or of the newly infected term, or the largest number of infected) and the input variable on a given finite time horizon. Some minimal time control problems have also been considered [AER + 14, BBSG17, BBDMG19]. Few results consider infinite horizon [Beh00] . Qualitatively, optimal solutions attached to the vaccination or isolation protocols are in general bang-bang 1 , with an intervention from the very beginning. Bounds on the number of switching times between the two modes (typically zero or one) are sometimes provided. By contrast, protocols of the third type (health promotion campaigns or lockdown policies) usually provide bang-bang optimal solutions with transmission reduction beginning after a certain time.

The present article is dedicated to the optimal control issue of obtaining, by enforcing social distancing, the largest value for S ∞ , the limit number of susceptible individuals at infinity. This means that we act on System (5) through the function u 3 (·) and set u 1 (·) = u 2 (·) = 0. This aim is also equivalent to minimize the total infection burden minimize the final number of removed ν T 0 I(t)dt, see (5a)-(5b). We thus seek here to determine how close to the herd immunity threshold it is possible to stop the spread of the disease, in the case where no vaccine or treatment is found to modify its evolution.

In the case of interventions on a infinite horizon, a possible action is to let the susceptible proportion reach the collective immunity level (u 3 = 1), and then impose total lockdown (u 3 = 0). This situation is illustrated in Figure 1 (right) for parameters of the Covid-19 in France in 2020 borrowed from [SKL + 20] and given in Table 2 below. We will focus here on interventions on a given finite horizon through possibly partial lockdown, i.e. with u 3 bounded from below. This setting is a prototypal problem, voluntarily unrealistic since for instance, no constraint on the peak value of the number of infected is considered. In a nutshell, the aim is to determine what best result may be obtained in terms of reduction of the total cumulative number of infected individuals, by applying lockdown of given maximal intensity and duration, in the worst conditions where no medical solution is discovered to stop earlier the epidemic spread. This models public interventions on the transmission rate by measures like social distanciation, restraining order, lockdown and so on, imposed on finite time horizon. Table  2 and u 1 (·) = u 2 (·) = 0. Left: no action (u 3 (·) = 1). Right: bang-bang action with one switch at the epidemic peak (u 3 (·) = 1 on (0, 62) and u 3 (·) = 0 otherwise).

For better readability, all results are exposed in Section 2. It is shown in Section 2.1 that, at least for strong enough lockdown measures and long enough intervention time T , it is possible to stop the epidemics arbitrarily close past the herd immunity. In Section 2.2, we provide and analyze the optimal control law that leads asymptotically, through an intervention of duration T , to the largest number of susceptible individuals. Last, we show in Section 2.3 that this optimal strategy coincides with a time minimal policy. For the sake of readability, all the proofs are postponed to Section 3. These results are numerically illustrated in Section 4.

For the sake of clarity, the demonstrations of all the results of this section are postponed to Section 3.

According to the introduction in the previous section, we consider in the sequel the following "SIR type" systemṠ (t) = −u(t)βS(t)I(t), t 0 I(t) = u(t)βS(t)I(t) − νI(t), t 0 (6) complemented with nonnegative initial data S(0) = S 0 , I(0) = I 0 such that S 0 + I 0 1. For given T > 0 and α ∈ [0, 1), the input control u will be assumed to belong to the so-called admissible set U α,T defined by

The constant T characterizes the duration of the intervention, and α its intensity (typically the strength of a lockdown procedure).

2.1 Toward an optimal control problem: reachable asymptotic immunity levels

The following result assesses the question of stopping the evolution exactly at, or arbitrarily close to, the herd immunity S herd defined by (2).

Proposition 1. Let α ∈ [0, 1) and T > 0. Assume that S 0 > S herd and consider

(i) If α α, then, for all ε ∈ (0, S herd ), there exist T > 0 and a control u ∈ U α,T such that the solution S to (6) associated to u satisfies

(ii) If α > α, then for all u ∈ U α,T the solution S to (6) associated to u satisfies

where S α is the solution to (6) associated to u ≡ α.

Moreover, the map α → lim t→+∞ S α (t) is strictly decreasing.

(iii) There is no T ∈ (0, +∞) and u ∈ U α,T such that the solution S to (6) associated to u satisfies lim t→+∞ S(t) = S herd .

The proof of Proposition 1 is given in Section 3.2.

Remark 1. The values given in Table 2 , which correspond to the sanitary measures put in place in France during the March to May 2020 containment period, satisfy S 0 > S herd and α α. According to Table 1 , there thus exists a containment strategy that increases the limit proportion of susceptible by 30%.

Conclusion of this section. To tackle the issue of reaching the closest state to the threshold of herd immunity with the help of control functions u(·), it seems relevant to maximize the limit of S(t) when t → +∞ for a given duration T of the intervention. This leads to consider the optimal control problem sup

with (S, I) the solution to (6) associated to u. We now make an observation, which will be crucial in the analysis (see details in Section 3.3). It turns out that the quantity S + I − 1 R0 ln(S) is constant on any time interval on which u(·) = 1 (see Lemma 2). Therefore, using the fact that lim t→+∞ I(t) = 0 and the monotonicity of

where (S, I) is the solution to (6) associated to u.

This section is devoted to the analysis results of the optimal control problem (P α,T ). The first result of this section reduces the study of the complete problem to that of a one dimensional optimization problem, whose unknown denoted T 0 stands for a switching time. For simplicity, given α ∈ [0, 1], T > 0 and T 0 ∈ [0, T ], we define the function u T0 ∈ U α,T by

Also, we denote (S T0 , I T0 ) the solution of (6) with u = u T0 .

The proof of Theorem 1 is given in Section 3.3.

Figure 2: Opt. control At this step, the result above proves that the optimal control belongs to a family of bang-bang functions parametrized by the switching time T 0 . This reduces the study to a rather simple one-dimensional optimization problem, at least from the numerical point of view.

We deduce that the switching time associated to the optimal solution of Problem (P α,T ) solves the 1D optimization problem sup

where j(T 0 ) := S ∞ (u T0 ), and u T0 is defined by (8). The next result characterizes the optimal T 0 for ( P α,T ) whenever α ∈ [0, 1]. To state it, one needs to introduce the function ψ given by

where (S T0 , I T0 ) denotes the solution to (6) with u = u T0 .

Theorem 2. Let T > 0, α ∈ [0, 1) and T * 0 the unique optimal solution to Problem ( P α,T ) (hence u T * 0 is the optimal solution to Problem (P α,T ) with the notation (8)). One has the following characterization:

In the particular case α = 0, one has T * 0 > 0 if, and only if T > 1 ν ln S0 S0−S herd , and in that case, T * 0 is the unique solution to the equation

The reduction to the optimal control ( P α,T ) and this characterization of the optimal T * 0 constitute useful tools to refine the numerical search of optimal controls. Ongoing work. We conclude this part by underlying that the optimal control problem we have solved is a prototypal one. In particular, it does not appear completely relevant since it does not take into account a cost related to the occupancy of hospital beds and the fact that this number of beds is limited. This leads, as a first approximation, to add to the problem of optimal control a constraint on the total number of infected persons. This will be the main topic of future work.

We have shown (see Theorem 1) that, provided that α is sufficiently small, there exists for every ε > 0 a time of control T such that, by acting on [0, T ], we have lim t→+∞ S(t) S herd − ε, i.e. at a distance ε of the herd immunity threshold. Then, we may wonder what is the minimal time of action such that the final value of susceptible is at a given distance of this threshold. This amounts to solve the following optimal control problem.

Minimal time problem: for ε > 0, determine the minimal time of action T * > 0 such that the optimal final number of susceptible individuals satisfies

We recall that S * ∞,α,T is defined in Theorem 1. The following result answers to this question by noting that solving this problem is equivalent to solve Problem (P α,T ).

Theorem 3. Assume α α (defined in (7)) and let ε > 0.

Let T * ε > 0 be the solution to the minimal time problem above and denote u * ε the corresponding control function. Then, u * ε is the unique solution of Problem (P α,T ) determined in Theorems 1 and 2 associated to T = T * .

Conversely, let T > 0 and S * ∞,α,T the maximum of Problem (P α,T ). Then, T is the minimal time of intervention such that S ∞ (u) S * ∞,α,T for some u ∈ U α,T .

3 Proofs of the main results

Before proving the main results, we will provide some useful elementary properties of the state variables (S, I) solving (6) whose role in the sequel will be central. To this aim, it is convenient to introduce the function Φ γ defined for any γ > 0

Let us start with a preliminary result regarding the value of Φ γ along the trajectories of System (6).

along any trajectory of system (6). In particular, if u is constant on a non-empty, possibly unbounded, interval, then the function Φ R0u (S(·), I(·)) is constant on this interval along any trajectory of System (6).

Proof. The proof of (12) follows from straightforward computations. Indeed along every trajectory of (6), one has

The second part of the statement is an obvious byproduct of this property, by setting γ = R 0 u.

Lemma 1 allows to characterize the value of the limit of S at infinity, as now stated.

Lemma 2. Let u ∈ U α,T . For any trajectory of (6), the limit S ∞ (u) of S(t) at infinity exists and is the unique solution in the interval [0, 1/R 0 ] of the equation

where Φ R0 is given by (11).

Proof. Any input control u from U α,T is equal to 1 on [T, +∞). Hence, applying Lemma 1 with u = 1 on this interval yields (13), by continuity of Φ R0 and because of (3). Moreover, Eq. (13) has exactly two roots. Indeed, this follows by observing that the mapping S → Φ R0 (S, 0) is first decreasing and then increasing on (0, +∞), with infinite limit at 0 + and +∞, and minimal value at S = 1/R 0 = S herd , equal to 1 R0 (1 + ln R 0 ). We conclude by noting that the limit S ∞ (u) cannot be larger than S herd : otherwise there would exist ε > 0 such thatİ > ε > 0 andṠ < −εβS for T large enough, so that S would tend to zero at infinity, yielding a contradiction. It follows that the value of S ∞ (u) is thus the smallest root of (13).

Remark 2 (On the control in infinite time). Notice that, in the quite unrealistic situation where one is able to act on System (6) up to an infinite horizon of time, then the optimal strategy to maximize S ∞ (u) is to consider the constant control function u α (·) = α on [0, +∞). Indeed, according to Lemma 1,

for all t ∈ [0, T ]. Hence, Φ αR0 (S 0 , I 0 ) Φ αR0 (S ∞ (u), 0), with equality if, and only if u = u α = α a.e. on R + . Since S ∞ (u) is maximal whenever Φ αR0 (S ∞ (u), 0) is minimal, we get that the optimal strategy in that case corresponds to the choice u = u α = α a.e. on R + . Moreover, it is notable that this maximal value S ∞ (u α ) is computed by solving the nonlinear equation Φ αR0 (S ∞ , 0) = Φ αR0 (S 0 , I 0 ). Therefore, an easy application of the implicit functions theorem yields that the mapping α → S ∞ (u α ) is decreasing.

Let us start with (iii). According to (6), one hasİ −νI, thus I(t) I 0 e −νt > 0 for all t > 0. Then,Ṡ < 0 and S is decreasing. Assume by contradiction that, for all t > 0, we have S herd < S(t).

Since I satisfies (6), it follows thatİ > 0 on (T, +∞) and then,Ṡ −βI(T )S on (T, +∞). We thus infer that S(t) e −βI(T )(t−T ) S(T ) for t > T , and thus S(t) → 0 as t → +∞, which is a contradiction.

Let us now show (i). For all α ∈ [0, 1), we denote by u α (·) the control equal to α on (0, ∞). Thanks to the same argument by contradiction as above, one has S ∞ (u α ) 1/(αR 0 ). Let us denote by (S uα , I uα ) the solution to System (6) associated to u α . Lemma 1 shows that the function t → Φ αR0 (S uα (t), I uα (t)) is conserved, and we infer that S ∞ (u α ) solves the equation

Using the expression of α,

It follows that for η > 0 small enough, there exists T > 0 such that for each t > T ,İ u α (t) (βα(1 + η) ν β − ν)I uα (t) < 0. By using a Gronwall lemma, one infers that I uα (t) → 0 as t → +∞. Then, for all k ∈ N * , there exists T k > 0, such that |S uα (T k ) − S herd | 1/k and I uα (T k ) 1/k. Consider u k := α1 (0,T k ) + 1 (T k ,∞) and let us denote by (S u k , I u k ) the solution to System (6) associated to u k . By continuity of Φ R0 ,

Let us finally prove (ii). We show that α → S ∞ (u α ) is strictly decreasing. Let α 1 , α 2 ∈ [0, 1) such that α 1 < α 2 and (S uα 1 , I uα 1 ) and (S uα 2 , I uα 2 ) the solutions to System (6) associated to u α1 and u α2 , respectively. Using Lemma 1, t → Φ α1R0 (S uα 2 (t), I uα 2 (t)) is strictly increasing, hence

Thanks to the equations satisfied by (S uα 1 , I uα 1 ) and (S uα 2 , I uα 2 ), one has S ∞ (u α1 ), S ∞ (u α2 ) S herd /α 1 . Since x → Φ α1R0 (x, 0) in strictly decreasing on (0, S herd /α 1 ), we deduce that S ∞ (u α2 ) < S ∞ (u α1 ). This concludes the proof since S ∞ (u α ) = S herd .

Solving (P α,T ) involves the resolution of an ODE system on an infinite horizon, and it is quite convenient to consider an equivalent version of this problem involving an ODE system on a bounded horizon. A key point for this is that, according to Lemma 2, S ∞ solves Eq. (13). Furthermore, since the mapping [0, 1/R 0 ] ∋ S → Φ R0 (S, 0) is decreasing, maximizing S ∞ is equivalent to minimize Φ R0 (S ∞ , 0). Combining all these observations yields that the optimal control problem is equivalent to the following version we will investigate hereafter:

and (S, I) solves the controlled system (6) associated to the control function u(·).

Proof of Theorem 1. For better readability, the proof of Theorem 1 is decomposed into several steps.

Step 1: existence of an optimal control We will prove the existence of an optimal control for the equivalent problem (P Φ α,T ). Let (u n ) n∈N be a maximizing sequence for Problem (P Φ α,T ). Since (u n ) n∈N is uniformly bounded, we may extract a subsequence still denoted (u n ) n∈N with a slight abuse of notation, converging towards u * for the weak-star topology of L ∞ (0, T ). It is moreover standard that U α,T is closed for this topology and therefore, u * belongs to U α,T . For n ∈ N, let us denote (S n , I n ) the solution to the SIR model (6) associated to u = u n . A straightforward application of the Cauchy-Lipschitz theorem yields that (Ṡ n ,İ n ) n∈N is uniformly bounded. By applying Ascoli's theorem, we may extract a subsequence still denoted (S n , I n ) n that converges towards (S * , I * ) in C 0 ([0, T ]). As usually, we consider an equivalent formulation of System (6), where (S n , I n ) can be seen as the unique fixed point of an integral operator. We then pass to the limit and show that (S * , I * ) solves the same equation where u has been replaced by u * . By combining all these facts with the continuity of Φ R0 , we then infer that (J Φ (u n )) n∈N converges up to a subsequence to J Φ (u * ), which gives the existence.

Step 2: optimality conditions and bang-bang property We will again establish these properties for the equivalent problem (P Φ α,T ). In what follows, for the sake of simplicity, we will consider and denote by u a solution to Problem (P Φ α,T ) and by (S, I), the associated pair solving System (6). Observe first that integrating (12) in Lemma 1, one has J Φ (u) = Φ R0 (S(T ), I(T )) = Φ R0 (S 0 , I 0 ) + ν It is standard to write the first order optimality conditions for such kind of optimal control problem. To this aim, we use the so-called Pontryagin maximum principle (see e.g. [LM67] ) and introduce the Hamiltonian H defined on R 4 by 

where w denotes the Lipschitz-continuous switching function given by w = S(p S − p I ). By using (6), one computeṡ

We will now prove that the optimal control can be written as u T0 defined in (8), for some T 0 ∈ [0, T ). From (14b), one has p S (T ) = p I (T ) = 0, and thus w(T ) = 0 > − 1 R0 . According to (15), this implies u(·) = α on a certain maximal interval [T 0 , T ], for some T 0 ∈ (0, T ), by continuity of w. By inserting the relation w − 1 R0 holding on (T 0 , T ) in the equation (14a) satisfied by p I , we deduce thatṗ I ν(p I − 1) on (T 0 , T ). Since p I (T ) = 0, the Gronwall lemma yields 

Since ε ∈ (0, νS(T 0 )(1 − p I (T 0 )), this implies that w is strictly increasing in V T0 . Therefore, there exists a maximal open interval (T 1 , T 0 ) with T 1 ∈ [0, T 0 ) on which w < − 1 R0 , and therefore on which u = 1. As a consequence, the left-derivative of w at T 0 exists and readṡ

We will in fact show that T 1 = 0. In other words, the control u can be written as (8).

To this aim, let us assume by contradiction that there exists T 1 ∈ (0, T 0 ) such that w(T 1 ) = − 1 R0 = w(T 0 ) and w < − 1 R0 on (T 1 , T 0 ). Observing that w is differentiable on (T 1 , T 0 ) and using Rolle's theorem yields the existence of τ ∈ (T 1 , T 0 ) such thatẇ(τ ) = 0.

Note that, according to (14a), one hasẇ = −Sṗ I , a.e. in (0, T ). Since S(τ ) > 0, one has alsȯ p I (τ ) = 0. Using the fact that u = 1 on (T 1 , T 0 ), this means that the point (p I (τ ), w(τ )) is a steady-state of the systemẇ = −βSw − νSp I ,ṗ I = βw + νp I .

According to the Cauchy-Lipschitz theorem, we infer that (p I , w) is constant on [τ, T 0 ] and therefore,ṗ I (T − 0 ) =ẇ(T − 0 ) = 0, which is in contradiction with (17). As a conclusion, T 1 = 0, and u can be written as (8).

Step 3: monotonicity of S * ∞,α,T Let T > 0 and 0 α α < 1. It is straightforward that Uα ,T ⊂ U α,T , and then S * ∞,α,T S * ∞,α,T . It follows that the map [0, 1) ∋ α → max u∈Uα,T S ∞ (u) is nonincreasing.

Let us show that the map T → max u∈Uα,T S ∞ (u) is nondecreasing. Let 0 < T T , 0 α < 1, and denote by u * ∈ U α,T the control realizing the maximum S * ∞,α,T . Since u * = 1 in (T,T ), one has u * ∈ U α,T . Thus S * ∞,α,T S ∞ (u * ) = S * ∞,α,T .

Step 4: uniqueness of the optimal control The demonstration of the uniqueness of the optimal control is achieved in the proof of Theorem 2 below, by demonstrating the uniqueness of the optimal switching point T 0 . Except this point, the demonstration of Theorem 1 is now complete.

Let us decompose this proof into several steps. We assume that α > 0, the case α = 0 will be considered in the last step. In the whole proof, we will deal with control functions u T0 as defined in formula (8).

Step 1: necessary first order optimality conditions on T 0 Let u = u T0 be an optimal control for problem (P α,T ), with T 0 be the associated optimal switching time. Let us introduce the criterion j Φ given by

where (S T0 , I T0 ) is the solution corresponding to the control u T0 , as previously defined. For the sake of simplicity, we omit these subscripts in the sequel. Again, it is equivalent to maximize j and to minimize j Φ . By using Lemma 2, one has

where c 0 = I 0 + S 0 − ν β ln S 0 . According to (18), we infer that S solves the systeṁ

with the initial data S(0) = S 0 . Using (18), one gets

so that the cost function reads

The next lemma allows to compute the derivative of j Φ with respect to T 0 . For the sake of clarity, its proof is postponed to the end of this section. Thanks to this result, we may compute

By noting that

we have for the function ψ defined in (9):

We deduce that j ′ Φ (T 0 ) = 0 is equivalent to ψ(T 0 ) = 0.

Step 2: Zeros of j ′ Φ and uniqueness of the optimal switching time According to (6), one has for any t ∈ (T 0 , T ), I(t) = I(T 0 ) exp t T0 (αβS(s) − ν) ds . Then, using the expression of ψ given in (20), it follows that

Introducing ϕ : [0, T ] → R defined by ϕ(s) = exp T s (αβS(s) − ν) ds , the last expression writes simply

Differentiating this identity with respect to T 0 yields

As a consequence of Lemma 3, both terms in the previous formula are negative, and ψ ′ (T 0 ) < 0. The function ψ is thus decreasing on (0, T 0 ). Moreover, ψ(T ) = −1 < 0. Therefore, if ψ(0) < 0, then (21), or equivalently j ′ Φ (T 0 ) = 0, has no solution, and thus T * 0 = 0. Conversely, if ψ(0) 0, then (21) admits a unique solution T * 0 which is the unique critical point of j Φ , and of j. In particular, in the case ψ(0) = 0, one has T * 0 = 0. We also deduce that the function j Φ is nonincreasing on (0, T * 0 ) and increasing on (T * 0 , T ), and so is j.

Step 3: Upper bound on T 0 Let us now prove that S T (T * 0 ) S herd , where, in accordance with the notation defined after (8), (S T , I T ) is the solution to System (6) for u(·) = 1. Assume by contradiction that S T (T * 0 ) < S herd , i.e. T * 0 > (S T ) −1 (S herd ). Since I is non-increasing on ((S T ) −1 (S herd ), ∞), we have I T ((S T ) −1 (S herd )) > I T (T * 0 ). Due to the fact (obtained from Lemma 1) that

Since S T (T * 0 ) ∈ (0, S herd ) and x → x − 1 R0 ln(x) is decreasing on (0, S herd ), the inequality (22) implies that S T (T * 0 ) > S herd , which leads to a contradiction. We have thus demonstrated that S T (T * 0 ) S herd .

Step 4: The case α = 0 Let us finally deal with the case "α = 0". Using the fact that S is constant on (T 0 , T ), we deduce that p S (t) = 0 and p I (t) = 1 − e ν(T0−T ) for all t ∈ (T 0 , T ). A commutation occurs at T 0 if, and only if,

The function S is nonincreasing, thus there exists T 0 > 0 satisfying this relation only if S 0 > S herd 1−e −νT which is equivalent to T > 1 ν ln S0 S0−S herd . If this is the case, then, since t → S(t) is nonincreasing and T 0 → S herd 1−e ν(T 0 −T ) is increasing, there exists a unique T 0 satisfying the relation (23). We also remark that (23) is equivalent to (21).

To achieve the proof of Theorem 2, it now remains to prove Lemma 3.

Proof of Lemma 3. Using the notation S T0 previously defined, one has (see (19b)) on (T 0 , T )

and, at T 0 , S T0 (T 0 ) is defined thanks to (19a) by

By differentiating (24) with respect to T 0 , one infers

which is the second identity in Lemma 3. Furthermore, using (19b), one has

Differentiating this relation with respect to T 0 yields for all t ∈ (T 0 , T ]

Let us simplify this latter identity. Observe first that, because of (18), one has for all t ∈ (T 0 , T ]

By using at the same time the change of variable v = S(t) and the identity c 0 + ν β 1 − 1 α ln S(T 0 )− S(t) + ν αβ ln S(t) = I(t), holding true for any t ∈ (T 0 , T ], we infer that

Combining all these facts leads to, for all t ∈ (T 0 , T ],

Therefore, we arrive at, for all t ∈ (T 0 , T ],

which is the first identity of the statement. This achieves the proof of Lemma 3, and consequently of Theorem 2.

Let T * be the minimal time associated to the Minimal time problem. Let u * α,T * be the unique solution of Problem (P α,T ) associated to T = T * . Assume by contradiction that u * α,T * does not solve the minimal time problem, i.e. S ∞ (u * α,T * ) < S herd − ε. Then for each u ∈ U α,T * one has

which is in contradiction with the fact that T * solves the minimal time problem. Conversely, let T > 0 and u * realizing the maximum of (P α,T ), i.e. S ∞ (u * ) = S * ∞,α,T . Let T * be the minimal time of intervention such that S ∞ (u) S * ∞,α,T for some u ∈ U α,T * . Since

Hence, by definition of S * ∞,α,T , one has v = u * . But v does not have the form of the minimal solution in Theorem 1, unless T * = T .

This section is devoted to numerical illustrations around the optimal control problem (P α,T ). The codes are available on:

https://github.com/michelduprez/optimal-immunity-control.git

We use the parameter values given in Table 2 , coming from [SKL + 20] and corresponding to the lockdown conditions in force in France between March 17th and May 11th 2020. We suppose that, on the total population of 6.7e7 persons in France, there are no removed individuals and 1000 infected individuals at the initial time, i.e. R 0 = 0 and I 0 = 1e3/6.7e7.

Solutions to ODEs are computed with the help of a Runge-Kutta fourth-order method and the research of the optimal controls is conducted thanks to the Algorithms described below.

Our approach rests upon the use of gradient like algorithms needing the computation of the differential of S ∞ in an admissible direction 3 h. According to the proof of Theorem 1 (see below Initial proportion of susceptible 1 − I 0 I 0 Initial proportion of infected 1e3/N R 0 Initial proportion of removed 0 Section 3.3), this differential reads

where (p S , p I ) denotes the adjoint state, solving the backward adjoint system (14a)-(14b). Thanks to this expression of DS ∞ (u) · h, we deduce a simple projected gradient algorithm to solve numerically the optimal control problem (P Φ α,T ), then (P α,T ). The algorithm is described in Algorithm 1. The projection operator P Uα,T is given by

Algorithm 1 Solving Problem (P α,T ) by projected gradient descent

Compute (S k , I k ) solution to the primal system (6) 3:

Compute (p S,k , p I,k ) solution to the dual system (14a)-(14b)

Compute u k+1 := P Uα,T (u k − ρ k DS ∞ (u k )) where ρ k is the step, chosen variable such that S ∞ (u k+1 ) − S ∞ (u k ) < 0. 5: end while 6: return u k+1

We recall that it is equivalent to maximize j and to minimize j Φ . According to Theorem 1, solutions u * and T * 0 to (P α,T ) and P α,T , respectively, satisfy u * = u T * 0 , moreover j Φ is decreasing on (0, T * 0 ) and increasing on (T * 0 , T ). Taking advantage of this property, we also test a second algorithm based on the bisection method. This method is described in Algorithm 2.

Let us compare the numerical solutions obtained with these approaches for each problem (P α,T ) and ( P α,T ). Optimal trajectories denoted respectively (S u * , I u * , R u * ) and (S T * 0 , I T * 0 , R T * 0 ) associated to the optimal control problems (P α,T ) and ( P α,T ) are provided on Figures 3 and 4 for the parameter choices α ∈ {0, α lock } and T = 100. To capture the full behavior of the trajectories, the computation window is [0, 200] . One recovers all the theoretical results: the solutions to (P α,T ) and ( P α,T ) are as expected and do coincide.

By using Lemma 2, it is easy to determine numerically the optimal value S * ∞,α,T by solving the equation

Algorithm 2 Solving Problem ( P α,T ) (then (P α,T )) by bisection method Require: k ∈ N * 1: Initialization: T 0,min,0 = 0 and T 0,max,0 = T 2: for i = 1, . . . , k do T 0,min,i+1 = T 0,min,i and T 0,max = T 0,right 10:

end if 11: end for 12: return T 0,k+1 := (T 0,max,k+1 + T 0,min,k+1 )/2, u k+1 := u T 0,k+1

This allows to investigate numerically on Fig. 5 and 6 the dependency of S * ∞,α,T and T * 0 with respect to the parameters T and α. On Fig. 5 , for T = 400, we observe numerically that the lower bound α given in Proposition 1 and below which S can get as close as we want to S herd over an infinite horizon, is optimal (α ≈ 0.56). In particular, for lockdown conditions similar to the ones in effect in France between March and May 2020 (α ≈ 0.231), it appears that it is possible to come as close as we want to the optimal bound S herd of inequality (4). Interestingly, we observe on Fig.  6 -left that when α is not small enough we have T 0 = 0 for T large enough.

We also mention that Fig 5-left gives the solution to the minimal time problem. Indeed, from Theorem 3, a control u * is optimal for (P α,T ) iff it is optimal for the minimal time problem. Then, given ε > 0 and α α, the minimal time of action such that the final value of susceptible is at a distance ε of S herd is obtained by computing the intersection of the curves in Fig 5- left with the horizontal line S herd − ε. As expected, when α is too large, i.e. when the lockdown is insufficient, the solution stays far away from S herd (see Fig 5-right) . Figure 7 represents the function j Φ for different values of α and T . We recover the fact that j Φ is decreasing then increasing stated in Theorem 2. The fact that j Φ is non-convex highlights the difficulty of the considered optimal control problem. [Gre88] David Greenhalgh, Some results on optimal control applied to epidemics, Mathematical Biosciences 88 (1988), no. 2, 125-158.

[GS09] Holly Gaff and Elsa Schaefer, Optimal control applied to vaccination and treatment strategies for various epidemiological models, Mathematical Biosciences & Engineering 6 (2009), no. 3, 469.

A simple planning problem for covid-19 lockdown

Andris Abakuks, An optimal isolation policy for an epidemic

Optimal immunisation policies for epidemics

Free terminal time optimal control problem of an SIR epidemic model with vaccination

Optimal screening in structured SIR epidemics

Optimal control of epidemic size and duration with limited resources

Time-optimal control strategies in SIR epidemic models

Optimal control of deterministic epidemics

A synthesis problem for the optimal control of epidemics

Optimal control of epidemics with limited resources

Optimal vaccination strategies for an influenza epidemic model

Optimal treatment profile during an influenza epidemic

Optimality of a time-dependent treatment profile during an epidemic

Optimal regime of double vaccination

Optimal control of an epidemic through social distancing

Optimization of immunocorrection of collective immunity

Foundations of optimal control theory

Globally optimal vaccination policies in the SIR model: smoothness of the value function and uniqueness of the optimal strategies

Optimal, near-optimal, and robust epidemic control

On the optimal control of a deterministic epidemic

The optimal COVID-19 quarantine and testing policies

Optimal dengue vaccination strategies of seropositive individuals

Estimating the burden of SARS-CoV-2 in France

Optimal isolation policies for deterministic and stochastic epidemics

Optimal immunization rules for an epidemic with recovery

Optimal vaccination policy and cost analysis for epidemic control in resource-limited settings

Optimal isolation strategies of emerging infectious diseases with limited resources

The authors express their sincere acknowledgment to working groups Maths4Covid19 and Op-timCovid19 for fruitful discussions and in particular their colleagues Luis Almeida (CNRS UMR 7598, LJLL, France), Emmanuel Franck (INRIA Grand-Est and IRMA Strasbourg, France), Sidi-Mahmoud Kaber (Sorbonne Université, LJLL, France), Grégoire Nadin (CNRS UMR 7598, LJLL, France), Benoît Perthame (Sorbonne Université, LJLL, France).