Assessing the Impact of Increasing Antimicrobial Resistance of Vibrio cholerae on the Future Trends of Cholera Epidemic

Cholera, an acute intestinal infection caused by the bacterium Vibrio cholerae, remains a major public health problem in many parts of Africa, Asia, and Latin America. A mathematical model is developed, to assess the impact of increasing antimicrobial resistance of Vibrio cholerae on the future trends of the cholera epidemic. Equilibrium states of the model are determined and their stabilities have been examined. The impacts of increasing antimicrobial resistance of Vibrio cholerae on the future trends of cholera epidemic have been investigated through the reproductive number. Numerical results are provided to support analytical findings.


Introduction
The global burden of cholera remains substantial. In 2005, 131,943 cases and 2,272 deaths were reported to the WHO (World Health Organization), and recently major, sustained epidemics have been reported in many parts of Africa, such as Zimbabwe and Nigeria to mention a few [1]. These statistics are gross underestimates, as many cholera-endemic countries do not report cholera to the WHO, including Zimbabwe, which has one of the highest rates of cholera in the world. More realistic estimates of the global burden of cholera mortality place the figure at 100,000-150,000 deaths per year [2]. Cholera is contracted by ingestion of food or water contaminated with the Gram-negative bacterium Vibrio cholerae. The bacteria pass through the human gastric acid barrier into the small intestine where they colonize, multiply, and begin to secrete cholera toxin. The primary treatment for cholera is rehydration with oral or intravenous fluids [3]. For severe cases, antimicrobial agents may reduce the volume and duration of diarrhea [3]. Tetracyclines (e.g., doxycycline), fluoroquinolones (e.g., ciprofloxacin), macrolides (e.g., erythromycin), and trimethoprim/sulfamethoxazole have commonly been used to treat cholera [3].
Antimicrobial drug resistance can undermine the success of antimicrobial therapy. V. cholerae becomes drug resistant by exporting drugs through efflux pumps, chromosomal mutations, or developing genetic resistance via the exchange of conjugative plasmids, conjugative transposons, and integrons [4]. Several reports have documented tetracycline-and fluoroquinolone-resistant V. cholerae, and multidrug resistance is increasing [4]. The aim of this study is to assess the impact of an increase on antimicrobial resistance of Vibrio cholerae on the future trends of the cholera epidemic, with the aid of a simple mathematical model. Mathematical models have become an important tool in describing the dynamics of the spread of an infectious disease and the effects intervention strategies have on its dynamics. Modeling of epidemiological phenomena has a very long history with the first model for small pox formulated by Daniel Bernoulli in 1760 [5].
A brief survey on previous studies provides the context of this paper. Various theoretical studies have been carried out on mathematical modeling of cholera transmission dynamics, focusing on a number of different issues, see [6][7][8][9][10][11], to mention a few. King et al. (2008) [11] proposed a two-path cholera model including a class for severe infections as well as a class for mild or inapparent infections. Log-likelihood estimates using 1900's mortality data in 26 Bengal districts were generated by King et al. [11] and indicate that the districts of Bogra and Calcutta are well described by the two-path model, yet differ greatly in several parameter estimates, most notably the estimated proportion of infections resulting in severe symptoms (approximately 0.24 and 0.02, resp.). In addition to emphasizing the need for two infectious classes, the estimates show that natural immunity to cholera may wane within a year, as opposed to the previously suspected period of 3-10 years. Neilan et al. (2010) [8] proposed a mathematical model for investigating the optimal intervention strategies for cholera. Their work managed to provide a framework for designing costeffective strategies for diseases with multiple intervention methods. More recently, Das and Mukherjee (2012) [6] proposed a mathematical model to investigate the role of lytic bacteriophage in the cyclic behaviour of cholera outbreaks. Results from their study suggest that, regardless whether phage exists in the environment or not, cholera persists in the environment. The paper is structured as follows. The cholera transmission model is formulated in Section 2. Analytical results for the model are presented in Section 3. The future trend and projection profiles of the proposed model are presented in Section 4. A short discussion rounds up the paper.

Model Formulation
Our objective is to formulate a model for cholera that includes drug resistant and drug sensitive bacteria. Thus, based on the individual's epidemiological status the total human population N is made up of the following subgroups: susceptible individuals S, individuals infected with cholera drug sensitive strain I 1 , and individuals infected with drug resistance bacteria I 2 . The total bacteria population B constitutes the drug sensitive bacteria B 1 and drug resistant bacteria B 2 . The susceptible human population (S) becomes infected with cholera through ingestion of bacteria contaminated water or food which contains either drug sensitive strain bacteria or drug resistant strain at rate λ 1 = ϑB 1 /(K 1 + B 1 ) and λ 2 = ϑB 2 /(K 2 + B 2 ), for drug sensitive and drug resistant, respectively. The parameter ϑ denotes the ingestion rate, K 1 , K 2 captures the bacteria concentrations measured with respect to their infection doses, for drug sensitive strain and drug resistant strain, respectively. Assuming homogeneous mixing and constant recruitment of the human population into the susceptible class at rate Λ, the model takes the form Natural mortality μ is assumed to occur in all human population classes at a constant rate; ω represents the net death rate of vibros; φ 1 , φ 2 defines the average contribution of each cholera infected individual to the aquatic population of V. cholerae, for those infected with drug sensitive strain and drug resistant strain, respectively. Treated individuals recover at rate r 1 (for those infected with drug sensitive strain) and r 2 (for those infected with drug resistant strain). We assume that a proportion 1 − p of individuals infected with drug sensitive strain on treatment will successfully recover from the disease, and the remainder p develop resistance due to chromosomal mutations or genetic resistance. Individuals infected with drug sensitive bacteria are assumed to be susceptible to infection by the drug resistance bacteria at rate λ 2 . Further, untreated infected individuals may succumb to cholera related mortality at rate δ 1 and δ 2 , for individuals infected with drug sensitive strain and drug resistant strain, respectively. The model flow diagram is depicted in Figure 1.
For system (1), the first octant in the state space is positively invariant and attracting; that is, solutions that start where all the variables are nonnegative remain there. Thus, system (1) will be analyzed in a suitable region which is positively invariant and attracting. Existence, uniqueness, and continuation results for system (1) holds in this region.

Disease-Free Equilibrium and Its
Stability. System (1) has a disease-free E 0 , given by Following Driessche and Watmough [12], the reproductive number for system (1), defined as the spectral radius of an irreducible or primitive nonnegative matrix is given by Biologically R 1 and R 2 measure the average number of secondary infections generated by a cholera sensitive strain and resistant strain bacteria, respectively, when introduced in a disease-free population. The reproductive number R 0 controls the number of equilibria of system (1). If R 0 ≤ 1, then the disease-free equilibrium E 0 is the only equilibrium in Φ.
If R 0 > 1, then there exists a unique endemic equilibrium E * . Using Theorem 2 in van den Driessche and Watmough (2002) [12], the following result is established.
Following Kamgang and Sallet (2008) [13], we write system (1) in the form on the positively invariant set Φ ⊂ R n1+n2 + . Here x 1 = S and x 2 = (I 1 , I 2 , B 1 , B 2 ). Here x 1 ∈ R 1 + denotes the uninfected components and x 2 ∈ R 4 + denoting the number of infected components, and x * 1 = E 0 . We have to prove that the following conditions hold.
(H 1 ) The system is defined on a positively invariant set Φ of the nonnegative orthant. The system is dissipative on Φ.
(H 2 ) The subsystem is Metzler (a Metzler matrix is a matrix with off-diagonal entries nonnegative [14]) and irreducible for any given x ∈ Φ.

(H 4 ) There exists an upper-bound matrix
, the points where the maximum is realized are contained in the disease-free submanifold).
We express the subsystem This is a linear system which is globally asymptotically stable at the equilibrium Λ/μ corresponding to E 0 , which satisfies conditions H 1 and H 2 . The matrix A 2 (x) is given by which is an irreducible matrix. The upper bound for x ∈ Φ is given by which is not attained on Φ and which is not the corresponding block in the Jacobian matrix of the system at E 0 . Thus we obtain only a sufficient condition. Condition H 5 requires that α(A) ≤ 0; thus

ISRN Biomathematics
The maximum of two quantities in (9) is the reproductive number for system (1). Thus, E 0 is globally asymptotically stable whenever R 0 ≤ 1. We summarize the result in Theorem 2.

Endemic Equilibrium and Its Stability Analysis.
There are three possible endemic equilibrium states for system (1), namely, (a) if R 1 > 1 and R 2 < 1, then there exists a drug sensitive strain induced endemic equilibrium point; (b) if R 1 < 1 and R 2 > 1, then there exists a drug resistant strain induced endemic equilibrium point; (c) if R 1 > 1 and R 2 > 1, then there coexists a drug sensitive and drug resistant strain induced endemic equilibrium point.
If either case (a) or case (b) exists, then the global stability of resultant system follows from [15]. Now we investigated the local stability for case (c), with the aid of the centre manifold theory. When case (c) exists, then the endemic equilibrium for system (1) in terms of the forces of infection is given by Before starting our main results, we give the following lemma which will be useful in the subsequent section.
Theorem 3 (see [16]). Consider the following general system of ordinary differential equations with a parameter φ, Without loss of generality, it is assumed that 0 is an equilibrium for System (11) for all values of the parameter φ, that is, f (0, φ) = 0 for all φ and assumed 0)) is the linearization of system (11) around the equilibrium 0 with φ evaluated at 0. Zero is a simple eigenvalue of A and other eigenvalues of A have negative real parts; (A2) matrix A has a right eigenvector w and a left eigenvector v corresponding to the zero eigenvalue.
Let f k be the Kth component of f and The local dynamics of (11) around 0 are totally governed by a and b.
To apply this method, the following change of variables is required. Let S = x 1 , I 1 = x 2 , I 2 = x 3 , B 1 = x 4 , and B 2 = x 5 , further, by using vector notation Linearizing system (13) about E 0 gives from (8) it can be shown that Consider a case R 2 = 1 (i.e., R 1 < R 2 = 1). Suppose, further, that ϑ = ϑ * is chosen as a bifurcation parameter. Solving for ϑ = ϑ * , we have Thus, the linearized system of the transformed equation (13) with ϑ = ϑ * chosen as the bifurcation parameter has a simple zero eigenvalue. Hence, it can be shown that the Jacobian of (13) at ϑ = ϑ * has right and left eigenvectors (associated with the zero eigenvalue) given below.

Eigenvectors of J.
It can be shown that the Jacobian J of system (13) at ϑ = ϑ * has a right eigenvector (corresponding to the zero eigenvalue) given by w = (w 1 , w 2 , w 3 , w 4 , w 5 ) T , where Further, J has a left eigenvector (associated with the zero eigenvalue) given by

Computations of the Bifurcation Coefficients a and b.
For the sign of a, it can be shown that the associated nonvanishing partial derivatives of F are (20) It follows from (19) that For the sign of b, it can be shown that the associated nonvanishing partial derivatives of F are It follows from (22) Based on the results in (21), (23), and Theorem 3, we establish the following results.
Theorem 4. The endemic equilibrium E * is locally asymptotically stable for R 0 > 1, but close to 1.

Impact of Drug Resistance on the Dynamics of Cholera.
The reproduction number of an epidemic (the mean number of secondary cases infected by a single infectious case) is a key parameter for the analysis of infectious diseases because it summarizes the potential transmissibility of the disease and indicates whether an epidemic is under control. Now, we investigate the impact of increasing R 2 on R 1 , through partial differentiation, as follows: Results in (24) suggest that an increase on antimicrobial resistance of Vibrio cholerae may lead to an increase in cholera prevalence in the community.

Population Level Effects
In order to support analytical findings in this study, we now simulate system (1) using the parameters in Table 1. The average life span 1/μ was set to be 50 years.
Numerical results in Figure 2 suggests that whenever cholera sensitive strain and resistant strain coexist in the community, then in the long run (averagely 25 days or less according to Figures 2(a) and 2(c)) cumulative resistant strains will outnumber cumulative sensitive strains. In   Table 1, and the following assumed initial conditions S = 20000, I 1 = 500, I 2 = 500, and B 1 = B 2 = 1000, 0000, 000. Figure 2(b) we note that if we have more drug resistant cases than drug sensitive cases on the outbreak, then they will be a pronounced increase of cumulative drug resistant strain cases on a period, averagely of 10 days or less. Here (Figure 3) we note that an increase of p will lead to a decrease on the susceptible population and cumulative infectious drug sensitive cases, while cumulative infectious drug resistant cases will rise. These results further highlight the impact drug resistance has on the long-term dynamics of cholera.
Simulations in Figure 4 illustrate the effect of varying bacteria ingestion rate ϑ, on cumulative cholera cases. Figures  4(a) and 4(b) depict the associated infection levels when cholera drug sensitive and drug resistant coexist. Figure 4(c) illustrates cumulative infections in the absence of drug resistant strain. For ϑ > 0.5 we note that infection levels have a sharp increase and they attain their maximum in a period of 20 days or less, despite whether there is coexistence of drug sensitive and drug resistant or not, then decreases sharply, until they attain stability. Further analysis of Figure 4 reveals that for low bacteria ingestion rate (ϑ < 0.5) the cumulative infections increase gradually and attain the maximum after a period of 20 days, and then decrease gradually until stability is attained.

Sensitivity
Analysis of R 0 . Sensitivity analysis of model parameters is very important to design and control strategies as well as a direction to future research. There are many methods [17] available for conducting sensitivity analysis such as differential analysis, response surface methodology, the Fourier amplitude sensitivity test (FAST) and other variance decomposition procedures, fast probability integration, and sampling-based procedures. In this section, sensitivity analysis based on Latin hypercube sampling (LHS) [17,18] has been performed using relevant parameters' values in Table 1. Sensitivity analysis assesses the amount and type of change inherent in the model as captured by the terms that define the reproductive number (R 0 ). If R 0 is very sensitive to a particular parameter, then a perturbation of the conditions that connect the dynamics to such a parameter may prove useful in identifying policies or intervention strategies that reduce epidemic prevalence. In this section the 8 ISRN Biomathematics  partial rank correlation coefficients (PRCCs) were calculated to estimate the correlation between values of R 0 and the model parameters across 1000 random draws from the empirical distribution of R 0 and its associated parameters. Results in Figure 5 suggest that reduction of the cholera bacteria through death has a significant impact on reducing the magnitude of the reproductive number. Although, there is a need to treat infected individuals, results on Figure 5 suggest that implementing intervention strategies targeting at reducing the bacteria population may have more impact on controlling cumulative cholera cases and related mortality since half saturation constant and human recruitment rate are the two parameters, which have the greatest impact on decreasing and increasing the magnitude of R 2 , respectively.  Figure 5: Partial rank correlation coefficients showing the effects of parameter variation on, R 2 using ranges in Table 1. Parameters with positive PRCCs will increase R 2 when they are increased, whereas parameters with negative PRCCs will decrease R 2 when they are increased.     Table 1. Parameters with positive PRCCs will increase R 1 when they are increased, whereas parameters with negative PRCCs will decrease R 1 when they are increased.
We now examine the dependence of the R 2 using the Latin hypercube sampling technique. Results in Figure 6 support the earlier findings in Figure 5, that high values of K 2 and Λ lead to low and high values, respectively, of R 2 .
Overall, results in Figure 7 are in agreement to the results in Figure 5; however it is worth noting that the inclusion of a proportion of individuals who develop drug resistance has an impact on increasing the magnitude of the reproductive number, though not highly significant compared to the other model parameters.

Discussion
As the causative agent of cholera, the bacterium Vibrio cholerae represents an enormous public health burden, especially in the developing countries around the world. A deterministic model for assessing the impact of increasing antimicrobial resistance of Vibrio cholerae on the long-term transmission dynamics of cholera has been proposed and analyzed. The reproductive number for the model has been computed. Qualitative analysis of the model reveals that the model has a globally asymptotically stable infection free equilibrium state whenever the associated reproductive number is less than unity. With the aid of centre manifold theory it has been established that for the reproductive number greater but close to unity the model has a locally asymptotically stable endemic equilibrium point. Comprehensive analytical and numerical results from the study suggest that an increase on cumulative antimicrobial resistance of Vibrio cholerae has a negative impact on controlling cholera prevalence within the community, more specifically, if the outbreak is constituted of higher drug resistant cases than drug sensitive cases.

Concluding Remarks
In conclusion, our model suggests that the level of cholera infection depends highly on the rate of human recruitment and bacteria ingestion rate. We note that for bacteria ingestion rate above 0.5 cumulative cholera cases have a sharp increase within the period of 20 days, while for bacteria have ingestion less than 0.5 cumulative cholera cases a gradual increase and attain maximum after 20 days, and they also decrease steadily until stability is attained. We also observed that if a cholera outbreak have higher drug resistant cases than drug sensitive cases on the outbreak, then there will be a pronounced increase of cumulative drug resistant strain cases on a period, averagely of 10 days or less.