Global dynamics of an age-structured cholera model with multiple transmissions, saturation incidence and imperfect vaccination

ABSTRACT In this paper, an age-structured cholera model with multiple transmissions, saturation incidence and imperfect vaccination is proposed. In the model, we consider both the infection age of infected individuals and the biological age of Vibrio cholerae in the aquatic environment. Asymptotic smoothness is verified as a necessary argument. By analysing the characteristic equations, the local stability of disease-free and endemic steady states is established. By using Lyapunov functionals and LaSalle's invariance principle, it is proved that the global dynamics of the model can be completely determined by basic reproduction number. The study of optimal control helps us seek cost-effective solutions of time-dependent vaccination strategy against cholera outbreaks. Numerical simulations are carried out to illustrate the corresponding theoretical results.


Introduction
Cholera is an acute disease caused by Vibrio cholerae O-group 1 or O-group 139, which can give rise to acute diarrhoea and vomit. The World Health Organization (WHO) estimates that there are 1.3-4 million cholera cases per year with about 21,000-143,000 deaths all over the world. Beginning in April 2017, a major cholera epidemic occurred in Yemen, with about 500,000 reported cases and 2000 deaths [29]. In order to better understand the transmission dynamics of cholera and provide some valuable insights on the prevention and control, some cholera models have been proposed (see, for example, [15,19,24,25,30]). In [25], Tien and Earn pointed out that new infections arise both through exposure to contaminated water (indirect transmission), as well as by human-to-human transmission pathway (direct transmission).
In [10], Hartley et al. found that short-lived, hyperinfectious state of vibrios decay in a matter of hours into a state of lower infectiousness and incorporated this hyperinfectious state into a cholera model to provide a much better fit with the observed epidemic pattern of cholera. Besides, Neilan et al. [18] formulated a mathematical model to include two classes of bacterial concentrations, one is hyperinfectious and another is less-infectious. Furthermore, in [17], Mukandavire et al. considered that human-to-human transmission is assumed to be a very fast transmission process with a lower infectious dose as a result of immediate water or food contamination by hyperinfectious vibrios from freshly passed human stool. Besides, Mukandavire et al. pointed out that if vibrios have been in the environment for a sufficiently long period (anywhere from 5 to 18 h), they are no longer hyperinfectious [17], namely, in a hypoinfectious state.
However, it was assumed in the models mentioned above that direct and indirect transmission rates, recovery rate of infected individuals, the contribution rate of each infected individual to the concentration of V. cholerae and net death rate of V. cholerae are invariant with time. In [1], Brauer et al. considered the following age-structured cholera model: with the boundary conditions where corresponding variables and parameters in system (1) are described in Table 1. System (1) was also investigated by Yang and Qiu [26] and Wang and Zhang [27]. It was assumed in [1] that incidence rates are bilinear, which regards the infection rate per density of infected individuals or per concentration of V. cholerae as a constant. Actually, incidence rate is influenced by the inhibition effect from behavioural change of susceptible individuals and the crowding effect of vibrios. In [4], Capasso and Serio introduced a saturation incidence rate βI/(1 + αI), where βI measures the infection force of the disease and 1/(1 + αI) measures the inhibition effect and the crowding effect. There have been several works on cholera models with saturation incidence in the literature (see, for example, [5,20]). As for the infection rate for per concentration of V. cholerae, experimental studies [21] indicated that the probability of infection depends on the concentration of vibrios in the contaminated water. In [6], Codeço introduced a new form βP/(k + P) to measure the effect of saturation, where β is the contact rate with contaminated water and P/(k + P) is the probability of individuals to develop cholera. This saturation incidence form has been used in some cholera models (see, for example, [19,24,30,31]).
To prevent and control cholera, some preventive strategies, such as vaccination, have been widely used in several cholera outbreaks and gained great efficiencies. The availability of effective oral cholera vaccines have renewed interests in the use of vaccines for cholera control with overwhelming evidence supporting their effectiveness in endemic settings [22]. In 2010, WHO recommended the use of vaccines in cholera endemic settings and preemptively during outbreaks. In [19], Posny et al. presented a cholera epidemiological model with three types of control measures: vaccination, treatment and sanitation, which is described by the following system: where S(t), V(t), I(t) and R(t) denote the densities of susceptible individuals, vaccinated individuals, infected individuals and recovered individuals, respectively. While, P(t) denotes the concentration of V. cholerae. N is the total population and all individuals are born and die naturally at rate μ. Susceptible or vaccinated individuals acquire cholera infection either by ingesting environmental vibrios from contaminated aquatic reservoirs or through human-to-human transmission, at rates λ e = (1 − ρ)β e P/(k + P) and λ h = (1 − ρ)β h I, respectively. Here, ρ = εp is the sanitation-induced preventability to cholera infection which is a product of the sanitation efficacy ε and compliance p. Susceptible individuals are vaccinated at a rate φ, while, the vaccine has a degree of losing protection efficacy that is denoted by σ . Infected individuals are treated at a rate τ , and some recover naturally at a rate γ . Infected individuals contribute to V. cholerae in the aquatic environment at a rate ξ , and vibrios have a net death rate δ. In addition, water sanitation leads to the death of vibrios at a rate ν. Motivated by the works of Posny et al. [19], Brauer et al. [1] and Capasso and Serio [4], in the present paper, we are concerned with both human-to-human and environment-tohuman transmissions, saturation incidence and imperfect vaccination on the transmission dynamics of cholera. To this end, we consider the following differential equations: with the boundary conditions and the initial condition where L 1 + (0, ∞) is the space of functions on (0, ∞) that are positive and Lebesgue integrable. In system (3), the infection rate of infected individuals is given by saturation infection rate β 1 (a)i(a, t)/[1 + αi(a, t)], where α is the saturation infection rate coefficient. Besides, the infection rate for per concentration of V. cholerae takes the form where k is the concentration of V. cholerae in contaminated water that yields 50% chance of catching cholera. Corresponding flowchart of cholera transmission in system (3) is depicted in Figure 1. Variables and parameters in system (3) are described in Table 1. Thereinto, A, μ, α and k are nonnegative and bounded, while, φ and σ are nonnegative and less than unity. Denote function space X = R + × R + × L 1 By the standard theory of age-structured model [11,28], it can be verified that system (3) with boundary conditions (4) and initial condition (5) has a unique nonnegative solution. Thus, system (3) generates a continuous semi-flow : R + × X → X , which takes the form (t, This paper is organized as follows. In Section 2, some preliminaries are given for later analysis. In Section 3, we investigate the asymptotic smoothness of the semi-flow { (t)} t≥0 . Next, we study the existence of disease-free and endemic steady states and calculate basic reproduction number in Section 4. In Section 5, the local asymptotic stability of each of steady states is established by analysing the distribution of roots of characteristic equations. In Section 6, we discuss the global asymptotic stability of each of steady states by using suitable Lyapunov functionals and LaSalle's invariance principle. We carry out a study of optimal control to seek cost-effective solutions of control strategies for cholera in Section 7.
In Section 8, we present numerical simulations to illustrate theoretical results and obtain the optimal control solution by Forward-Backward Sweep Method. The paper ends with a conclusion in Section 9.

Preliminaries
Before analysing the global dynamics of system (3), we make the following reasonable assumptions based on biological significance, which hold throughout the paper.
is point dissipative and attracts all points in X .
Proof: From (6), we obtain that It follows from (7) that d dt Noting that ρ 1 (0) = 1 and dρ 1 (a)/da = −θ(a)ρ 1 (a), we have d dt Similarly, we obtain d dt From (9) and the first two equations of (3), we get It follows from the variation of constants formula that Hence, for any solution of (3) satisfying x 0 ∈ , the following inequality holds From (10) and (12), we have Similarly, by the variation of constants formula, we obtain that Adding (11) and (13) yields which implies that (t, x 0 ) ∈ holds for ∀t ≥ 0, x 0 ∈ . Moreover, it follows from (14) that lim sup t→∞ t (x 0 ) X ≤ A/μ 0 for any x 0 ∈ X . Therefore, is point dissipative and attracts all points in X . This completes the proof.
From the proof of Proposition 2.1, it is not difficult to verify the following result that will be useful in next section. Corollary 2.1: If x 0 ∈ X and x 0 X ≤ K for some constant K ≥ A/μ 0 , then the following statements hold for t ≥ 0: for ∀t 1 , t 2 ∈ , in which denotes an interval defined on R + .

Asymptotic smoothness
Noting that system (3) is an infinite dimensional dynamical system, the asymptotic smoothness of solutions is the prerequisite before discussing the stability of each of steady states. Hence, in this section, we investigate the asymptotic smoothness of the semi-flow { (t)} t≥0 generated by system (3).

Definition 3.1: ([23]) For any nonempty and closed bounded set
The following lemmas are useful in proving the asymptotic smoothness of the semi-flow

following two conditions hold:
(I) There exists a continuous function u : From Proposition 3.13 in [28], we can obtain the following theorem by applying Lemmas 3.1 and 3.2.

Theorem 3.1: The semi-flow { (t)} t≥0 generated by system (3) is asymptotically smooth.
Proof: (7) and (8), we have This completes the proof of (I) in Lemma 3.1. From (7), we have It follows from (1) and (2) in Corollary 2.1 that, for any solution of system (3) with x 0 ∈ , Note that It follows from (1) and (2) in Corollary 2.1 that and According to (3) in Corollary 2.1, we have Finally, it follows from (17) It is not difficult to find that the right side of the above inequality converges uniformly to 0 as h → 0. Hence, for t ≥ 0, ϕ(t, ·) maps any bounded sets of into sets with compact closure in X . Similarly, we can obtain that p 2 (b, t) satisfies similar conditions in Lemma 3.2. Hence, the semi-flow { (t)} t≥0 is asymptotically smooth.

The existence of steady states
It is easy to see that system (3) always has a disease-free steady state E 0 (S 0 , V 0 , 0, 0, 0), where If system (3) has an endemic steady state it must satisfy the following equations: It follows from (21) and (22) that According to (23), solving (21e) yields Substituting (22) and (24) into (21a) and (21b), we have where Noting that 0 < S * + V * < A/μ, namely, A − μS * − μV * = 0, Equation (25) becomes Obviously, g belongs to a two-variable function. In order to establish the existence of the endemic steady state, we have to prove that there is a unique positive solution to the following system of two equations: where N is the same function defined in Equation (25) and This can be seen that h(0) = 0 and h(x) is a monotone increasing function because of It is sufficient to prove that h(x) = 1 has a unique positive solution if h(A/(μ + φ)) > 1.
We therefore obtain the basic reproduction number as follows which represents the average number of new infections generated by a single newly infectious individual during the full infectious period. Thus, if R 0 > 1, system (3) has a unique endemic steady state E * (S * , V * , i * (·), R * , p * (·)).

Theorem 5.2:
If R 0 > 1, the endemic steady state E * is locally asymptotically stable.

Global asymptotic stability
In this section, we investigate the global stability of each of steady states to system (3) by using suitable Lyapunov functionals and LaSalle's invariance principle. Define a function g(x) = x − 1 − ln x, which is always positive except for x = 1 where g(x) = 0. Theorem 6.1: If R 0 < 1, the disease-free steady state E 0 of system (3) is globally asymptotically stable in X .

Proof: Let (S(t), V(t), i(a, t), R(t), p(b, t)) be any positive solution of system (3) with boundary conditions (4) and initial condition (5). Define
where ε 1 (a) and η 1 (b) will be determined later. Calculating the derivative of V 1 (t) along positive solutions of system (3) yields Noting Using integration by parts, we obtain that Similarly, we have which have the following properties From (43)-(46), it follows thaṫ  p(b, t)) :V 1 (t) = 0}. Noting that if R 0 < 1, E 0 is locally asymptotically stable, thus we obtain the global asymptotic stability of E 0 directly from LaSalle's invariance principle.
Before investigating the global stability of the endemic steady state E * to system (3), defineā and Y = R + × R + ×Ỹ.

Optimal control strategy
From Section 6, we know that if R 0 > 1, a cholera outbreak will take place and the disease will persist. Therefore, it is critically important to investigate how to effectively control cholera and how to achieve the disease control with typically limited resources.
Control strategies, such as quarantine, vaccination, treatment, and sanitation, can realize the control of cholera to varying degree and at different cost. Thereinto, vaccination strategy is the most efficient control strategy to prevent, control and eradicate cholera. Next, we will find a cost-effective vaccination strategy to control cholera epidemic. We rewrite system (3) as where the constant is introduced to balance the different units between time and age. We also rewrite the vaccination rate as u(t), a function depending on time. The initial and boundary conditions for system (61) are given in (4) and (5). We consider this system on a time interval [0, t max ]. The control set is defined as = {u(t) ∈ L ∞ (0, t max )|0 ≤ u(t) ≤ u max }, where u max denotes the upper bound due to the practical limitation on the vaccination that can be implemented within a given population and a given time period. We aim to minimize the number of infected individuals and corresponding cost of the strategy during the course of an epidemic. Define the objective functional as where a 1 , a 2 and a 3 are the weight constants of infected individuals, vaccination strategy and side effects of vaccination strategy, respectively. The square of the control variable shows the severity of the side effects. The minimization process is subject to the state equations in (61), where S,V,i,p are the state variables, while, u is the control variable.
Following [2], [7] and [8], we construct the optimal control model through the combination of the state equations, the adjoint equations, and the optimality condition. Denote a sensitivity function F = (S, V, i, p) and a solution map u → F = F(u). Based on results in [2], the map is differentiable and the sensitivity function F is defined by the Gateaux derivative: with boundary conditions and initial condition Denote Then, the adjoint system associated with control u(t) and state variables S(t), V(t), i(a, t), Following [8] and [3], the optimality condition is found as Numerical techniques for optimal control problems can often be classified as either direct or indirect. In terms of disease control, indirect methods, such as Forward-Backward Sweep Method, approximate solutions to optimal control problems by numerically solving the boundary value problem for the differential-algebraic system generated by the Maximum Principle [13].

Numerical simulations
In this section, we want to illustrate the theoretical results for system (3) by numerical simulations. Furthermore, by Forward-Backward Sweep Method, we obtain the optimal control strategy and show the graph trajectories of infected individuals with optimal control and without optimal control.
Usually, the course of cholera is 3-7 days. After 6 days of quarantine since the symptoms disappeared, the faecal vibrios of infected individuals are negative on three consecutive checks, then one can think that the disease is cured. Hence, we set the length of infection age as 12.5 days. In the first 3 days, infected individuals become infectious gradually; in the middle 3 ∼ 7 days, infected individuals are in a hyperinfectious state; in the last days, infected individuals are cured or dead due to the disease, thus the transmission coefficient  [10] decreases to 0. The specific function is shown as follows: As for the contribution rate of each infected person to the concentration of V. cholerae, the change law is associated with β 1 (a). Since the susceptible individuals newly ingest the vibrios, they cannot produce vibrios immediately; 3 days after the onset of cholera, the contribution rate keeps in a high level; in the last days, because of the treatment or death, the contribution rate become smaller. Therefore, we set ξ(a) as ξ m a, 7d ≤ a < 12.5d.
Besides, V. cholerae can survive in river water, well water, or sea water for 1-3 weeks, while in fresh fish or shellfish for 1-2 weeks. Similarly, we set the length of biological age as 12.5 days. In the first 5 h, vibrios from freshly passed human stool are hyperinfectious; if vibrios have been in the environment from 5 to 18 h, they are no longer hyperinfectious [17], namely, in a hypoinfectious state; in the last hours, vibrios become non-infectious gradually. Corresponding function is shown as follows: From the practical situation of cholera epidemic, we consider if the infected individuals receive treatment too late, the recovery rate reduces sharply. Thus, γ (a) is set as Note that θ(a) contains the rates of recovery, disease-induced death, birth and natural death, while, the disease-induced death often happens in the middle course of cholera. Hence, θ(a) is chosen as The net death rate of V. cholerae increases after 5 hours due to their weak resistance in the environment, then δ p (b) follows that First, select Case 1 in Table 2 as the parameter values of system (3). Through direct calculation, we find that basic reproduction number R 0 is near 0.9877 and less than unity. From Theorem 5.1, we obtain that the disease-free steady state is locally asymptotically stable.
In Figure 2, we observe that susceptible individuals S(t) and vaccinated individuals V(t) converge to S 0 = 117.58 and V 0 = 897.10, respectively, while, infected individuals i(a, t) and V. cholerae p(b, t) converge to 0. Then, select Case 2 in Table 2 as the parameter values of system (3). We find that basic reproduction number R 0 is near 2.8722 and greater than unity. The endemic steady state is approximately calculated as S * = 2676.25, V * = 4984.97, i * (0) = 183.58, p * (0) = 7869.49. From Theorem 5.2, we obtain that the endemic steady state is locally asymptotically stable. From Figure 3, we observe that infected individuals i(a, t) and V. cholerae p(b, t) converge to the endemic steady state.

Optimal control solution
The numerical results associated with optimal control are obtained based on Forward-Backward Sweep Method [12], which has been developed for age-structured cholera models by [3]. A rough outline of the algorithm is given below. First, break the time interval [t 0 , t 1 ] into pieces with N+1 points. Here, x = (x 1 , . . . , x N+1 ) and λ = (λ 1 , . . . , λ N+1 ) are the vector approximations for the states in system (61) and the adjoint in system (65), respectively.
(1) Make an initial guess for u over the interval;  When all steps are complete, we obtain the optimal control strategy, which can be seen in Figure 4 (a). On account of medical technology and cost, vaccination control strategy has its limitation, thus we set u max = 80%. We observe that vaccination strategy could be reduced 140 days later from the beginning of the cholera outbreaks, which saves much cost of vaccination. In Figure 4 (b), we compare the graph trajectories of infected individuals with optimal vaccination strategy and original vaccination strategy (with the same parameter values to Case 1 in Table 2). It is clear that infected individuals have been reduced due to optimal vaccination control strategy, where the cost of optimal and original vaccination strategies are close. From Figure 4, we suggest that taking vaccination strategy at the very beginning of cholera outbreaks can reduce the number of infected individuals remarkably, which are also cost-effective optimal strategy.   Figure 4 is the graph trajectory of vaccination control strategy; Figure 4 is the graph trajectories of i(a, t) with optimal vaccination strategy and original vaccination strategy.

Conclusion
In this paper, we have considered a cholera model including both human-to-human and environment-to-human transmissions, saturation incidence and imperfect vaccination. By a complete mathematical analysis, the threshold dynamics of the model was established and it can be fully determined by basic reproduction number. If R 0 < 1, the disease-free steady state E 0 is locally and globally asymptotically stable; if R 0 > 1, the endemic steady state E * is locally and globally asymptotically stable. The study of optimal control helps us seek cost-effective solutions of time-dependent vaccination control strategy against cholera outbreaks. Numerical simulations vividly illustrate our main results of stability analysis for system (3). Furthermore, we obtain the optimal solution by Forward-Backward Sweep Method.
At the beginning of cholera epidemic, hyperinfectious vibrios freshly-shed from infected individuals play an important role on cholera transmission, due to that they are likely to come into contact with other individuals [10]. Therefore, the strategy of vaccination at the very beginning of the onset can effectively control the cholera epidemic. Besides, other control strategies, such as quarantine and sanitation, should be implemented as a supplement, which can better prevent and control cholera. It is worth mentioning that, after vaccination, the vaccine efficacy is waning as time goes by. Hence, some vaccinated individuals become susceptible again, which we leave for further investigation.

Disclosure statement
No potential conflict of interest was reported by the authors.