CUBO A Mathematical Journal Vol.20, No¯ 02, (53–66). June 2018 http: // dx. doi. org/ 10. 4067/ S0719-06462018000200053 Optimal control of a SIR epidemic model with general incidence function and a time delays Moussa Barro, Aboudramane Guiro and Dramane Ouedraogo Département de Mathématiques, UFR des Sciences et Techniques, Université Nazi Boni, Laboratoire d’Analyse Mathématique et d’Informatique (LAMI) B.P. 1091 Bobo-Dioulasso, Burkina Faso. mousbarro@yahoo.fr, abouguiro@yahoo.fr, dramaneouedraogo268@yahoo.ca. ABSTRACT In this paper, we introduce an optimal control for a SIR model governed by an ODE system with time delay. We extend the stability studies of model (2.2) in section 2, by incorporating suitable controls. We consider two control strategies in the optimal control model, namely: the vaccination and treatment strategies. The model has a time delays that represent the incubation period. We derive the first-order necessary conditions for the optimal control and perform numerical simulations to show the effec- tiveness as well as the applicability of the model for different values of the time delays. These numerical simulations show that the model is sensitive to the delays representing the incubation period. RESUMEN En este art́ıculo, introducimos un control óptimo para un modelo SIR gobernado por un sistema de EDOs con retardo temporal. Extendemos los estudios de estabilidad del modelo (2) en la sección 2, incorporando controles apropiados. Consideramos dos estrategias de control en el modelo de control óptimo, llámense: las estrategias de vacunación y tratamiento. El modelo tiene un retardo en el tiempo que representa el peŕıodo de incubación. Derivamos las condiciones necesarias de primer orden para http://dx.doi.org/10.4067/S0719-06462018000200053 54 Moussa Barro, Aboudramane Guiro and Dramane Ouedraogo CUBO 20, 2 (2018) el control óptimo y realizamos simulaciones numéricas para mostrar la efectividad y también la aplicabilidad del modelo para diferentes valores de los retardos temporales. Estas simulaciones numéricas muestran que el modelo es sensible a los retardos que representan el peŕıodo de incubación. Keywords and Phrases: SIR, general incidence, delays, optimal control, epidemic models, Hamiltonian. 2010 AMS Mathematics Subject Classification: 34K19, 34K20, 49K25, 49K30, 65N06, 90C90. CUBO 20, 2 (2018) Optimal control of a SIR epidemic model with general incidence . . . 55 1 Introduction Mathematical modeling of population are often used to describe the dynamics of epidemic diseases. This is a fast growing research area and has been paying important roles in discovering relations be- tween species and their interactions. There have been many variations such as classical epidemiolog- ical models [11]. These models are based on the standard Susceptible-Infectious-Susceptible (SIS), Susceptible-Infectious-Recovered (SIR) and Susceptible-Exposed-Infectious-Recovered (SEIR) mod- els, which are determined according to the difference on the method of transmission, nature of the disease, those with short/long incubation period, killer/curable diseases, etc, and the response of the individuals to it, for instance, gaining transient/permanent immunity, dying from the disease, etc.[6, 16]. The main purpose of formulating a such epidemiological model is to understand the long-term behavior of the epidemic disease and to determine the possible strategies to control it. Differential equations, whether there are ordinary, delay, partial or stochastic are one of the main mathematical tools being used to formulate many epidemiological models. The focus in such epi- demiological models has been on the general incidence at which people move from the class of susceptible individuals to the class of infective individuals.these general incidence have been mod- eled mostly by using bilinear and Holling type of functional responses [10, 12]. On the other hand, optimal control has extensively been used a strategy to control the epidemic outbreaks [8]. The main idea behind using the optimal control in epidemics is to search for, among the available strategies, the most effective strategy that reduces the infection rate to a minimum level while optimizing the cost of deploying a therapy or preventive vaccine that is used for con- trolling the disease progression [18]. In terms of epidemic diseases, such strategies can include therapies, vaccines, isolation and educational campaigns [3, 5]. Mathematical models have become important tools in analyzing the spread and control of infec- tious diseases. The model formulation process clarifies assumptions, variables parameters. There have been many studies that have mathematically analyzed infectious diseases [4, 7, 15]. Recently, many control optimal models pertaining to epidemic disease to epidemic diseases have appeared in the literature. They include, but not limited to, delayed SIRS epidemic model [13], delayed SIR model [1], tuberculosis model [17], HIV model [9] and dengue fever [2]. In this paper, we consider an optimal control problem governed by a system of delay differential equations with general incidence function and time delays. The governing state equations of the optimal control are described in a SIR framework with a general incidence function and a time delays representing the incubation period. Then we derive first-order necessary conditions for ex- istence of the optimal control and develop a numerical method to solve them. The rest of this paper is organized as follows. In section 2, we give the statement of the optimal control problem. We derive the necessary conditions for existence of the optimal control in section 3. In section 4, we describe the numerical method and present the resulting numerical simulations. Finally, we discuss these results in section 5 along with some concluding remarks. 56 Moussa Barro, Aboudramane Guiro and Dramane Ouedraogo CUBO 20, 2 (2018) 2 Statement of the optimal control problem Compute the optimal pair of vaccination and treatment strategies (u1, u2) that would maximize the recovered population and minimize both the infected and susceptible population, and at the same time minimize the costs of applying the vaccination and treatment strategies. So we consider the optimal control problem of the form (see Eihab B. M. et al): min (u1,u2)∈U J(u1(t), u2(t)) =    S(T) + I(T) − R(T) + ∫T 0 (c1u 2 1(t) + c2u 2 2(t) + S(t) + I(t) − R(t))dt, (2.1) subject to the quation    Ṡ = B − µ1S − f(S, Iτ) − u1S, İ = f(S, Iτ) − (µ2 + γ)I − u2I, Ṙ = γI − µ3R. (2.2) The two functions u1(t) and u2(t) represent vaccination and treatment strategies. These control functions are assumed to be L∞(0, T) functions belonging to a set of admissible controls U = {(u1, u2) ∈ (L ∞ (0, T))2 : u1min ≤ u1(t) ≤ u1max, u2min ≤ u2(t) ≤ u2max}, where 0 ≤ u1min < u1max ≤ 1 and 0 ≤ u2min < u2max ≤ 1. The two constants c1 and c2 are weighted cost associated with the use of the controls u1(t) and u2(t), respectively. The state equations are formulated from an SIR model with general incidence model, where S(t), I(t), R(t) are the numbers of susceptible, infected and recovered individuals at time t, respectively. The parameters B is the recruitment rate, the death rates for the classes are µ1, µ2 and µ3, respectively. The average time spent in class I before recovery is 1/γ. For biological reasons, we assume that µ1 ≤ µ2 + γ; that is, removal of infectives is at least as fast as removal of susceptibles. The time delays τ represents the incubation period. That is to say, only susceptible individuals who got infected a time t − τ are able to communicate the disease at time t. As general as possible, the incidence function f must satisfy technical conditions. Thus, we assume that H1 f is non-negative C1 functions on the non-negative quadrant, H2 for all (S, I) ∈ R2+, f(S, 0) = f(0, I) = 0. Let us denote by f1 and f2 the partial derivatives of f with respect to the first and to the second variable The differential equation model described by (2.2) without controls (u1 = u2 = 0) has two equilibrium points: a disease-free equilibrium E0 given by E0 = ( B µ1 , 0, 0 ) CUBO 20, 2 (2018) Optimal control of a SIR epidemic model with general incidence . . . 57 and an endemic equilibrium E∗ = (S∗, I∗, R∗) where, S∗ = B − (µ2 + γ)I ∗ µ1 I∗ = I∗ R∗ = γ µ3 I∗ The basic reproduction number of (2.2) without controls is given by R0 = f2(S 0, 0) µ2 + γ It was proven that if R0 < 1, then the disease-free equilibrium is asymptotically stable and if R0 > 1 then it is unstable. On the other-hand, when the controls is not null (u1 6= 0 and or u2 6= 0), we have the SIR model (2.2). The disease-free equilibrium for system (2.2) is given by Ec0 = ( B µ1 + u1 , 0, 0 ) (2.3) whereas the endemic equilibrium E∗c is given by E∗c = ( B − (µ2 + γ + u2)I ∗ µ1 + u1 , I∗, γ µ3 I∗ ) The basic reproduction number Rc of system (2.2) is given by Rc = f2(S 0, 0) µ2 + γ + u2 and it is clear that when u1 → 0 and u2 → 0 then Rc → R0 3 Existence and characterization of the optimal control In this section, we discuss the existence of the optimal control and then construct the Hamiltonian of the optimal control problem to derive the first order necessary conditions for the optimal control. 3.1 Existence of optimal control To show the existence of the optimal control for the problem under consideration, we notice that the set of admissible controls U is, by definition, closed and bounded. It is also convex because [u1min, u1max] × [u2min, u2max] is convex in R 2. It is obvious that there is an admissible pair ((u1(t), u2(t))) for the problem. Hence, the existence of the optimal control comes as a direct result from the Filippove-Cesari theorem [14]. We therefore, have the following result: 58 Moussa Barro, Aboudramane Guiro and Dramane Ouedraogo CUBO 20, 2 (2018) Theorem 3.1. Consider the optimal control problem (2.1) subject to (2.2). Then there exists an optimal pair of controls (u∗ 1 , u∗ 2 ) and a corresponding optimal states (S∗, I∗, R∗) that minimizes the objective function J(u1, u2) over set of admissible controls U. Proof. To prove the existence of an optimal control pair, it is important to verify the following assertion. (1) The set of controls and corresponding state variables is nonempty. (2) The admissible set U is convex and closed. (3) The right-hand side of the state system (2.2) is bounded by a linear function in the state and control variables. (4) The integrand of the objective functional LS,I,R(u1, u2) is convex on the set U. The hessian matrix of LS,I,R(u1, u2) on U is done by : M = ( 2c1 0 0 2c2 ) , Sp(M) = {2c1, 2c2} ⊂ R ∗ +, then, LS,I,R(u1, u2) is strictly convex in U. (5) There exists constants ω1 > 0, ω2 and ρ > 1 such that the integrand LS,I,R(u1, u2) of the objective functional satisfies LS,I,R(u1, u2) ≥ ω1|(u1, u2)| ρ − ω2. LS,I,R(u1, u2) = c1u 2 1(t) + c2u 2 2(t) + S(t) + I(t) − R(t) ≥ min(c1, c2)(u 2 1(t) + u 2 2(t)) − R(t) R(t) is bounded because N = S + I + R i.e ∃α, β, α < R(t) < β, ∀t Let ω1 = min(c1, c2) and ω2 = β. We have, LS,I,R(u1, u2) ≥ ω1‖(u1; u2)‖ 2 − ω2. CUBO 20, 2 (2018) Optimal control of a SIR epidemic model with general incidence . . . 59 3.2 Characterization of optimal control In this subsection, we derive the first order necessary conditions for the existence of optimal con- trol, by constructing the Hamiltonian H and applying the Pontryagin’s maximum principle. To simplify the notations, we write x(t) = [S(t), I(t), R(t)] T , u(t) = [u1(t), u2(t)] T and λ(t) = [λ1(t), λ2(t), λ3(t)]. We denote by g(u(t), x(t)) the integrand part of the objective function (2.1).With these notations and terminologies, the Hamiltonian is given by H = H(u(t), x(t), λ(t)) H = g(u(t), x(t)) + λT (t).ẋ(t) = c1u 2 1 + c2u 2 2 + S + I − R + λ1 ( B − µ1S − f(S, Iτ) − u1S ) (3.1) +λ2 ( f(S, Iτ) − (µ2 + γ)I − u2I ) + λ3 ( γI − µ3R ) . Let χ[a,b](t) be the characteristic function defined by χ[a,b](t) =    1, if t ∈ [a, b], 0, otherwise. (3.2) Let u∗ = [u∗1, u ∗ 2] T be the optimal control and x∗(t) = [S∗(t), I∗(t), R∗(t)]T be the corresponding optimal trajectory. Then there exists λ(t) ∈ R3 such that the first order necessary conditions for the existence of optimal control are given by the equations ∂H ∂u (t) = 0, (3.3) dx dt (t) = ∂H ∂λ , (3.4) dλ dt (t) = − ∂H ∂x . (3.5) The optimality conditions: [ ∂H ∂u1 (t) ] u(t)=u∗(t) = 0, (3.6) [ ∂H ∂u2 (t) ] u(t)=u∗(t) = 0. (3.7) Simplifying (3.5) and (3.6), we obtain 2c1u ∗ 1 − Sλ1 = 0, (3.8) 2c2u ∗ 2 − λ2I = 0. (3.9) 60 Moussa Barro, Aboudramane Guiro and Dramane Ouedraogo CUBO 20, 2 (2018) Further simplification of (3.8) and (3.9) yields u∗1(t) = min { u1max; max { 0; S(t)λ1(t) 2c1 }} (3.10) and u∗2(t) = min { u1max; max { 0; I(t)λ2(t) 2c2 }} . (3.11) The state equations: given by the forms (2.2) The co-state equations: dλ1 dt (t) = − ∂H ∂S (t), dλ2 dt (t) = − [ ∂H ∂I (t) + χ[0,T−τ] ∂H ∂Iτ (t + τ) ] , dλ3 dt (t) = − ∂H ∂R (t), which when simplified, lead to dλ1 dt = −1 + (λ1(t) − λ2(t))f1(S, I) + (µ1 + u1)λ1(t), dλ2 dt = −1 + (λ1(t + τ) − λ2(t + τ))χ[0,T−τ](t)f2(S, I) + (µ2 + γ + u2)λ2(t) − γλ3(t), dλ3 dt = 1 + µ3λ3(t). The transversality conditions: λ1(T) = 1, λ2(T) = 1, λ3(T) = −1. Remark 3.2. It is noting that 1. The Hamiltonian function H is strongly convex in the control variables. 2. The right-hand sides of the state and co-state equations are Lipschitz continuous. 3. The set of the admissible controls U is convex 4 Numerical simulations In this section, we apply the above optimal control theory with consideration of its applicability. we discuss the discretization of the optimal control problem described and present the numerical CUBO 20, 2 (2018) Optimal control of a SIR epidemic model with general incidence . . . 61 results obtained through our simulations. The algorithm describing the approximation method to obtain the optimal control is the following algorithm inspired from [13]. The Algorithm used here is a numerical variation of forward Euler method with a step size h. We explicitly write the forward Euler method for the state and the adjoint. step1: for i = −m, ..., 0, do : Si = S0; Ii = I0; Ri = R0; u i 1 = 0; ui 2 = 0 end for for i = n, ..., n + m λi1 = 1; λ i 2 = 1; λ i 3 = −1 end for step2 :for i = 0, ..., n − 1 Si+1 = Si + h(B − µ1Si − βSiIi − u i 1Si) Ii+1 = Ii + h(βSiIi − (µ2 + γ)Ii − u i 2 Ii) Ri+1 = Ri + h(γIi − µ3Ri) λn−i−1 1 = λn−i 1 − h(−1 + (λn−i 1 − λn−i 2 )βIi+1 + (µ1 + u i 1 )λn−i 1 λn−i−1 2 = λn−i 2 − h(−1 + (λn+m−i 1 − λn+m−i 2 )χ[0,T−τ](tn−i)βSi+1 + (µ2 + γ + u i 2 )λn−i 2 − γλn−i 3 ) λn−i−1 3 = λn−i 3 − h(1 + µ3λ n−i 3 ) ui+1 1 = Si+1λ n−i 1 /2c1 ui+1 2 = Ii+1λ n−i 2 /2c2 end for step3 :for i = 1, ..., n, write S∗(ti) = Si, I ∗(ti) = Ii R ∗(ti) = Ri u ∗ 1 (ti) = u i 1 and u∗ 2 (ti) = u i 2 . 62 Moussa Barro, Aboudramane Guiro and Dramane Ouedraogo CUBO 20, 2 (2018) Comments Fig 1. represent the different dynamics of the susceptible population for different aspect of control. The Orange color represent the population when there is treatment but not vaccination (u1 = 0 and u2 6= 0). The blue curve represent the population when there are vaccination and treatment (u1 6= 0 and u2 6= 0). The green curve show the evolution of the susceptible population when there is just treatment but not vaccination (u1 = 0 and u2 6= 0). This show that, without vaccination so many people are exposed to disease. Fig 1. Dynamic of susceptible population with different aspect of control. CUBO 20, 2 (2018) Optimal control of a SIR epidemic model with general incidence . . . 63 Fig 2. represent the different dynamics of the infected population for different aspect of con- trol. The Orange color present the evolution of the infected population when there is treatment but not vaccination (u1 = 0 and u2 6= 0). The blue curve represent the population when there are vaccination and treatment (u1 6= 0 and u2 6= 0). The green curve show the evolution of the infected population when there is just treatment but not vaccination (u1 = 0 and u2 6= 0). Naturally, as many people are exposed to the disease without vaccination, we see the growth of the infected population. Fig 2. Dynamic of infected population with different aspect of control. 64 Moussa Barro, Aboudramane Guiro and Dramane Ouedraogo CUBO 20, 2 (2018) Fig 3. represent the different dynamics of the infected population for different aspect of con- trol. The Orange color present the evolution of the recovered population when there is treatment but not vaccination (u1 = 0 and u2 6= 0). The blue curve represent the population when there are vaccination and treatment (u1 6= 0 and u2 6= 0). The green curve show the evolution of the recov- ered population when there is just treatment but not vaccination (u1 = 0 and u2 6= 0). As many people are exposed to the disease without vaccination, indeed we see the growth of the recovered population. Fig 3. Dynamic of recovered population with different aspect of control. In conclusion, observing the figures, we can deduce that the strategy leading to the vaccination alone (u1 6= 0 and u2 = 0) should be preferable to the joint use of vaccination (u1 6= 0) and treatment (u2 6= 0). The optimal control strategy here shows that prevention is more effective for the eradication of the disease. CUBO 20, 2 (2018) Optimal control of a SIR epidemic model with general incidence . . . 65 5 Conclusion In this paper, we considered an optimal control problem for a SIR model with time delay (repre- senting the incubation period ) and general incidence function. The main idea developed here is the optimal control in epidemics in order to search among the available strategies, the most effi- science one that reduce the infection rate to a minimum level while optimizing the cost deploying a therapy and preventive vaccine that is used to control the disease progression. The two control functions u1(t) and u2(t), which represent the vaccination and the treatment strategies are subject to time delays before being effective. Then we formulated the objective function of the optimal control problem. We discussed the existence of the optimal control and then derived the first order necessary conditions for the optimal control through constructing the Hamiltonian and using the Pontryagin’s maximum principle to achieve our aim. Finally, to end our study, we do a numerical simulation to corroborate the theoretical results obtained. Acknowledgments The authors want to thank the anonymous referee for his valuable comments on the paper. Competing interests The authors declare that they have no competing interests. Author’s contribution Aboudramane Guiro provide the subject, wrote the introduction and the conclusion and verified some calculation. Moussa BARRO conceived the study and computed the equilibria and their local stabilities. Dramane OUEDRAOGO wrote mathematical formula, he bring up the control strategy and did all the calculus with the second author. All the authors read and approved the final manuscript. References [1] A. Abta, H. Laarabi, H. T. Alaoui , The hopf bifurcation analysis and optimal Control of a delayed SIR epidemic model, Int. J. Anal. (2014), 1-10. [2] D. Aldila, T. Gotz, E. Soewono, An optimal Control problem arising from a dengue disease transmission model, Math. Biosci. 242 (1)(2013), 9-16. [3] F. G. Ball, E. S. Knock, P.D. O’Neil Control of emerging infectious diseases using responsive imperfect vaccination and isolation , Math. Biosci. 216 (1) (2008), 100-113. [4] N. Becker, The use of epidemic models, Biometrics 35(1978) 295-305. 66 Moussa Barro, Aboudramane Guiro and Dramane Ouedraogo CUBO 20, 2 (2018) [5] C. Castilho , Optimal Control of an epidemic through educational campaigns , Electron. J. Differ. Equ. 2006 (2006), 1-11. [6] C. Chiyaka, W. Garira, S. Dube, Transmission model of endemic human malaria in a partially immune population. Math. Comput. Model. 46 (2007), 806-822. [7] K. Dietz, The first epidemic model: a historical note on P.D. En’ko. Aust. J. Stat. 30A (1988), 56-65. [8] H. Gaff, E.Schaefer, Optimal control applied to vaccination and treatment strategies for various epidemiological models, Math. Biosci. Eng. 6 (3) (2009), 469-492. [9] K. Hattaf, N. Yousfi, Optimal Control of a delayed HIV infection model with immune response using an efficient numerical method , Int. Sch. Res. Netw. (2012) (2012), 1-7. [10] H. W. Hethcote, P. van den Driessche; Some epidemiological models with nonlinear incidence. J.Math. Biol. 29 (1991), 271-287. [11] H. W. Hethcote, The mathematics of infectious, SIAM Rev.42 (2000), 599-653. [12] A. Kaddar; On the dynamics of a delayed SIR epidemic model with a modified saturated incidence rate, Electron. J. Differ. Equ. 13 (2009), 1-7. [13] H. Laarabi, A. Abta, K. Hattaf , Optimal Control of a delayed SIRS epidemic model with vaccination and treatment , Acta Biotheor. 63 (15) (2015), 87-97. [14] S. Nababan ; A filippov-type lemma for functions involving delays and its application to time- delayed optimal control problems,optim. Theory appl. 27 3 (1979), 357-376. [15] P. Ogren, C. F. Martin, Vaccination strategies for epidemics in highly mobile populations. Appl. Math. Comput. 127 (2002), 261-276. [16] S. Ruan, D. Xiao, J. C. Beier; On the delayed ross-macdonald model for malaria transmission, Bull. Math. Biol. 70 (2008), 1007-1025. [17] C. J. Silva, D. F. Torres, Optimal Control strategies for tuberculosis treatment: a case study in angola , Numer. Algebra Control Optim. 2 (3) (2012), 601-617. [18] G. Zaman, Y.H. Kang, J.H. Jung, Optimal treatment of an SIR epidemic model with time delay , Biosystems 98 (1) (2009), 43-50. Introduction Statement of the optimal control problem Existence and characterization of the optimal control Existence of optimal control Characterization of optimal control Numerical simulations Conclusion