Optimal control of an epidemic model with a saturated incidence rate ∗

Abstract. In this study we consider a mathematical model of an SIR epidemic model with a saturated incidence rate. We used the optimal vaccination strategies to minimize the susceptible and infected individuals and to maximize the number of recovered individuals. We work in the nonlinear optimal control framework. The existence result was discussed. A characterization of the optimal control via adjoint variables was established. We obtained an optimality system that we sought to solve numerically by a competitive Gauss–Seidel like implicit difference method.


Introduction
The dynamics of the infectious diseases is an important research area in mathematical epidemiology.Understanding the transmission characteristics of infectious diseases in communities, regions and countries can lead to better approaches for decreasing the transmission of these diseases [1][2][3].
Epidemic models have been studied by many authors.Most of them are interested in the formulation of the incidence rate, i.e. the infection rate of susceptible individuals through their contacts with infected individuals (see, for example, [2,[4][5][6][7] and references therein).In order to model this disease transmission process, several authors employ the following incidence functions.The first one is the bilinear incidence rate βSI, where S and I are respectively the number of the susceptible and infected individuals in the population, and β is a positive constant.The second one is the saturated incidence rate of the form βSI/(1 + α 1 S), where α 1 is a positive constant.The effect of the saturation factor (refer to α 1 ) stems from epidemic control (taking appropriate preventive measures).The third one is the saturated incidence rate of the form βSI/(1 + α 2 I), where α 2 is a positive constant.In this incidence rate the number of the effective contacts between infected and susceptible individuals may saturate at high infecting levels due to the crowding of the infected individuals or due to the protection measures by the susceptible individuals.
In this study we consider a mathematical model of an SIR epidemic model with a modified saturated incidence rate: where S is the number of the susceptible individuals, I is the number of infected individuals, R is the number of the recovered individuals, A is the recruitment rate of the population, µ is the natural death of the population, α is the death rate due to disease, β is the transmission rate, α 1 and α 2 are the parameter that measure the inhibitory effect, and γ is the recovery rate of the infective individuals.

The optimal vaccination
Optimal control techniques are of great use in developing optimal strategies to control various kinds of diseases.To solve the challenges of obtaining an optimal vaccination strategy, we use optimal control theory [11].Our goal, then, is to reduce the numbers of susceptible and infected individuals and increase the number of recovered individuals.
For the optimal control problem, we consider the control variable u(t) ∈ U ad to be the percentage of susceptible individuals being vaccinated per unit of time.Here indicates an admissible control set.In this optimal problem, we assume a restriction on the control variable such as 0 ≤ u(t) ≤ u max [12], because vaccination of all of the entire susceptible individuals at one time is impossible.The physical meaning of the control variable in this problem is that if the low levels of the numbers of infected and susceptible individuals is reached we get more recovered individuals.In case of no vaccination, the number of infected and susceptible individuals increases while the numbers of recovered individuals decrease.The prefect time for vaccinations can bring the number of infected individuals down to a small level, the number of susceptible individuals begins to be built again, and more individuals recover from the infection.Now, we consider an optimal control problem to minimize the objective (cost) functional given by subject to • The first two terms in the functional objective represent benefit of S(t) and I(t) populations that we wish to reduce.A 1 and A 2 are positive constants to keep a balance in the size of S(t) and I(t), respectively.
• Naturally, each control incurs in some costs.Unfortunately, we do not have good data on the costs associated with vaccination control.Hence, we focus on the use of "relative" cost for the control.We use in the second term in the functional objective (as it is customary), the quadratic term (1/2)τ u 2 , where τ is a positive weight parameter which is associated with the control u(t) and the square of the control variable reflects the severity of the side effects of the vaccination [13].
There are other types of functions in literature (see, for example, [14]).As every vaccination requires a certain effort, the more natural is to consider quadratic function because it can make the analogy with the energy expended.
The objective of our work is to minimize the infected and susceptible individuals and, therefore, to maximize the total number of recovered individuals (vaccinated individuals are added directly to the recovered compartment), by using possible minimal control variables u(t).Susceptible individuals induce a percentage of vaccination control u(t) to protect against possible infection at a per unit time.Note that such a vaccine is protective against the primary infection.Therefore, we expect that the vaccine will be effective in the case of re-infection.So, recovered individuals, who are susceptible again, should provide the same percentage of vaccination control u(t) to protect from possible reinfection.

Existence of an optimal control
For existence, we consider a control system (2) with initial conditions.Then, we rewrite our system (2) in the following form: where φ = [S(t) and φ t denote derivative of φ with respect to time t.Equation ( 3) is a non-linear system with a bounded coefficient.
We set where Thus, it follows that the function D is uniformly Lipschitz continuous.From the definition of the control u(t) and the restriction on S(t), I(t) and R(t) ≥ 0, we see that a solution of the system (3) exists [15].
Let us go back to the optimal control problem (1)- (2).In order to find an optimal solution, first we find the Lagrangian and Hamiltonian for the optimal control problem (1)-(2).In fact, the Lagrangian of the optimal problem is given by L(S, I, u) = A 1 S(t) + A 2 I(t) + (1/2)τ u 2 (t).
Theorem 1.There exists an optimal control u * (t) such that subject to the control system (2) with initial conditions.
Proof.To prove the existence of an optimal control we use the result in [16].Note that the control and the state variables are nonnegative values.In this minimizing problem, the necessary convexity of the objective functional in u(t) is satisfied.The control space is also convex and closed by definition.The optimal system is bounded which determines the compactness needed for the existence of the optimal control.In addition, the integrand in the functional (1), A 1 S(t) + A 2 I(t) + (1/2)τ u 2 (t) is convex on the control u(t).Also, we can easily see that, there exist a constant ρ > 1, positive numbers ω 1 and ω 2 such that J(u(t)) ≥ ω 2 + ω 1 (|u| 2 ) ρ/2 .We conclude that there exists an optimal control.

Characterization of the optimal control
In the previous section we showed the existence of an optimal control, which minimize the functional (1) subject to the system (2).In order to derive the necessary conditions for this optimal control, we apply Pontryagin's maximum principle to the Hamiltonian H (4).If (x * (t), u * (t)) is an optimal solution of an optimal control problem, then there exists a non-trivial vector function λ(t) = (λ 1 (t), λ 2 (t), . . ., λ n (t)) satisfying the following equalities: which gives after derivation Now, we apply the necessary conditions to the Hamiltonian H (4).
Theorem 2. Let S * (t), I * (t) and R * (t) be optimal state solutions with associated optimal control variable u * (t) for the optimal control problem (1) and (2).Then, there exist adjoint variables λ 1 , λ 2 and λ 3 that satisfy where with transversality conditions Furthermore, the optimal control u * (t) is given by Proof.We use the Hamiltonian (4) in order to determine the adjoint equations and the transversality conditions.By putting S(t) = S * (t), I(t) = I * (t), R(t) = R * (t) and differentiating the Hamiltonian with respect to S, I and R, we obtain and by using the optimality conditions we find Nonlinear Anal.Model.Control, 2012, Vol.17, No. 4, 448-459 which gives Using the property of the control space, we obtain So the optimal control is characterized as The optimality system consists of the state system coupled with the adjoint system with the initial and transversality conditions together with the characterization of the optimal control.Utilizing the characterization of the optimal control, we have the following optimality system.

The improved GSS1 method
The resolution of the optimality system is created improving the Gauss-Seidel-like implicit finite-difference method developed by Gumel et al. [17], presented in [18] and www.mii.lt/NAdenoted GSS1 method.It consists on discretizing the interval [t 0 , t end ] the points t k = kl + t 0 (k = 0, 1, . . ., n), where l is the time step.Next, we define the state and adjoint variables S(t), I(t), R(t), λ 1 (t), λ 2 (t), λ 3 (t) and the control u(t) in terms of nodal points 3 and u k with S 0 , I 0 , R 0 , λ 0 1 , λ 0 2 , λ 0 3 and u 0 as the state and adjoint variables and the controls at initial time t 0 .S n , I n , R n , λ n 1 , λ n 2 , λ n 3 and u n as the state and adjoint variables and the controls at final time t end .
As it is well known, the approximation of the time derivative by its first-order forwarddifference is given, for the first state variable S(t), by We use GSS1 to adapt it to our case as following: we visit the variables one by one by blocking all other value to the most recently calculated where By applying an analogous technology, we approximate the time derivative of the adjoint variables by their first-order backward-difference and we use the appropriated scheme as follows: Hence, we can establish an algorithm to solve the optimality system and then to compute the optimal control.

Numerical results
Here we consider a general SIR epidemic model and all the parameter values are chosen hypothetically due to the unavailability of real world data.
We obtain Figs.1-4.The control individuals are marked by blue lines while the individual without control is marked by dashed black lines.In Fig. 1, we see that the population of the susceptible individuals sharply decreases in first week after that, it begins go to the stable state.
Figure 2 shows a significant difference in the number of infected individuals with and without control from the very beginning days of vaccination (in the first three days), after that it begins go to the stable state.
In Fig. 3, we see that the population of recovered individuals with vaccination increases rapidly from the third day and it goes to its stable state, we also see that the number of recovered individuals without control is very small.As expected, the population of the susceptible group decreases and consequently, the recovered group increases including the vaccinated susceptible group.

Conclusion
In this paper, we do not consider any special disease but our aim is to set up an optimal control problem relative to epidemic model with a modified saturated incidence rate, so it is to minimize the infected and susceptible populations and to maximize recovered populations.A comparison between optimal control and no control is presented.It is easy to see that the optimal control is much more effective for reducing the number of infected individuals.In order to illustrate the overall picture of the epidemic, the numbers of infected, susceptible and recovered individuals under the optimal control and no control are shown in figures.