Optimal Immunity Control and Final Size Minimization by Social Distancing for the SIR Epidemic Model

Abstract

The aim of this article is to understand how to apply partial or total containment to SIR epidemic model during a given finite time interval in order to minimize the epidemic final size, that is the cumulative number of cases infected during the complete course of an epidemic. The existence and uniqueness of an optimal strategy are proved for this infinite-horizon problem, and a full characterization of the solution is provided. The best policy consists in applying the maximal allowed social distancing effort until the end of the interval, starting at a date that is not always the closest date and may be found by a simple algorithm. Both theoretical results and numerical simulations demonstrate that it leads to a significant decrease in the epidemic final size. We show that in any case the optimal intervention has to begin before the number of susceptible cases has crossed the herd immunity level, and that its limit is always smaller than this threshold. This problem is also shown to be equivalent to the minimum containment time necessary to stop at a given distance after this threshold value.

Appendix A-Implementation Issues

We provide here an insight of the numerical methods used in Sect. 3. The codes are available on:

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

1.1 A.1. Solving Problem $\mathcal {P}_{\alpha ,T}$ by a Direct Approach

In order to check the consistency of the results of Theorem 2.3, we solved directly the optimal control problem $\mathcal {P}_{\alpha ,T}$. Our approach rests upon the use of gradient like algorithms, which necessitates the computation of the differential of $S_\infty$ in an admissible direction⁴h. According to the proof of Theorem 2.2 (see Sect. 4.3), this differential reads

$$\begin{aligned} DS_\infty (u)\cdot h=\int _0^T\left( \gamma -\beta S\left( p_I-p_S\right) \right) Ih \,dt, \end{aligned}$$

where $(p_S,p_I)$ denotes the adjoint state, solving the backward adjoint system (14a)–(14b). Thanks to this expression of $DS_\infty (u)\cdot h$, we deduce a simple projected gradient algorithm to solve numerically the optimal control problem $\mathcal {P}_{\alpha ,T}^{\varPhi }$, then $\mathcal {P}_{\alpha ,T}$. The algorithm is described in Algorithm 1. The projection operator $\mathbb {P}_{\mathcal {U}_{\alpha ,T}}$ is given by

$$\begin{aligned} \mathbb {P}_{\mathcal {U}_{\alpha ,T}}u(t)=\min \left\{ \max \{u(t),\alpha \},1\right\} , \quad \text {for a.e. }t\in [0,T]. \end{aligned}$$

1.2 A.2. Solving Problem $\mathcal {P}_{\alpha ,T}$ Thanks to the Theoretical Results

Taking advantage of the theoretical results given in Theorem 2.3, we then considered a simpler algorithm, based on the solution of Problem $\widetilde{\mathcal {P}}_{\alpha ,T}$ by bisection method. This method is described in Algorithm 2.

Fig. 8

Graph of the cost value j with respect to $T_0$ for $T=200$ (left) and $T=300$ (right)

All computations shown in Sect. 3 indicate, as expected, that the optimal trajectories computed by the two methods do coincide.

For illustrative purpose, we provide in Fig. 8 the graph of the cost $j(T_0)$ defined in (18) that corresponds to the one-dimensional optimization Problem $\widetilde{\mathcal {P}}_{\alpha ,T}$. As can be seen, the cost is not convex.