Modeling Cholera Epidemiology Using Stochastic Differential Equations

Abstract

In this study, we extend Codeço’s classical SI-B epidemic and endemic model from a deterministic framework into a stochastic framework. Then, we formulated it as a stochastic differential equation for the number of infectious individuals $It$ under the role of the aquatic environment. We also proved that this stochastic differential equation (SDE) exists and is unique. The reproduction number, $R_0$, was derived for the deterministic model, and qualitative features such as the positivity and invariant region of the solution, the two equilibrium points (disease-free and endemic equilibrium), and stabilities were studied to ensure the biological meaningfulness of the model. Numerical simulations were also carried out for the stochastic differential equation (SDE) model by utilizing the Euler-Maruyama numerical method. The method was used to simulate the sample path of the SI-B stochastic differential equation for the number of infectious individuals $It$, and the findings showed that the sample path or trajectory of the stochastic differential equation for the number of infectious individuals $It$ is continuous but not differentiable and that the SI-B stochastic differential equation model for the number of infectious individuals $It$ fluctuates inside the solution of the SI-B ordinary differential equation model. Another significant feature of our proposed SDE model is its simplicity.

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

The fight against cholera is far from over; it, therefore, becomes very reasonable to try and tackle the cholera infection also from theoretical and numerical points of view. It must be emphasized that looking at the high death rates of cholera, any research that aims at improving the success rates becomes crucial.

Cholera is an acute intestinal infection caused by the Vibrio cholerae bacterium being ingested in contaminated water or food which is characterized by extreme diarrhea and vomiting. Cholera has a brief incubation period, ranging from one to five days. The bacteria Vibrio cholerae produces a toxin known as enterotoxin that dehydrates and prevents the human body from absorbing liquids which can lead to death if treatment is delayed within two to three hours of infection [1]. According to Crooks and Atesmachew [2], getting portable drinking water and basic environmental hygiene is a major problem in Africa and Asia where cholera cases are on a continuous rise.

Cholera is one of the oldest diseases that continue to harm people of all ages, generating epidemics and pandemic outbreaks notwithstanding continued efforts to control its transmission. There exist a number of environmental factors that contribute to the spread of cholera infections. Because Vibrio cholerae can travel around in the water, any change in the hydrological cycle has the potential to modify the concentration of the pathogens in the water. Droughts and floods can boost or decrease the transmission process depending on the amount of rain and its seasonal nature [3].

The symptoms of cholera include painless watery stools, extreme vomiting, irregular heartbeat, and low blood pressure [4]. Most people do not fall sick with infection of cholera when exposed to Vibrio cholerae, but they can still infect other susceptible individuals through polluted water, as they shed the pathogen in their stool for 7 to 14 days. According to Akor [5], most people who are sometimes exposed have mild or asymptomatic symptoms, and at other times, symptoms are very severe. About one in every 20 infected individuals develops severe diarrhea with vomiting which could lead to dehydration.

The transmission of cholera can be direct or indirect [6]. The direct (human-to-human) transmission of cholera occurs when the infected person contact, engages in sexual activity with, or bites other infected individuals, whereas indirect (environment-to-human) transmission of cholera occurs when infected individuals consume Vibrio cholerae bacteria through contaminated waters and food [7].

As cholera epidemics become a global health burden in recent decades, researchers have paid more attention to cholera epidemiology. The interactions of Vibrio cholerae bacteria with its host and other pathogens in the environment have revealed that the dynamics of cholera are far more complicated than imagined previously [3]. Several mathematical models have been presented in the past by different authors to study the complicated epidemic and endemic nature of cholera. The study done by Capasso and Paveri-Fontana [8] in the Mediterranean formulated a deterministic mathematical model to investigate the outbreak of cholera that occurred in $1973$. Codeco extended and explicitly included the role of the aquatic environment in the dynamics of cholera in 2001, based on the work of Capasso and Paveri-Fontana [8]. The study done by Hntsa and Kahsay [9] focused on a mathematical model of the dynamical behavior of a fractional order model of cholera. Opoku and Afriyie [10] proposed and developed a dynamic mathematical model to analyze the role of the environment and control measures on cholera transmission dynamics involving the human population and the population of bacteria. Wang and Modnak [7] developed and analyzed an epidemiological model of cholera incorporating control measures, which is an extension of the formulated model by Mukandavire et al. [11] to better gain an understanding of the complex dynamics of cholera which included medicinal treatment, vaccination, and water sanitation effects. Hntsa and Kahsay [9] formulated a deterministic mathematical epidemic model to examine the influence of adequate cholera prevention strategies on the disease’s dynamics. According to Mark et al. [12], lytic bacteriophage specific for Vibrio cholerae may decrease the severity of epidemics of cholera by killing bacteria both in the infected persons and in the reservoir, based on environmental and epidemiological studies of cholera epidemics in Dhaka. Ochoche [13] proposed and formulated a mathematical model for transmission dynamics of cholera control with a water treatment control strategy. Fatima et al. [4] formulated and analyzed a mathematical epidemic model for cholera control in Nigeria that differed from earlier cholera models. In order to control cholera epidemics, this research model integrates treatment, water hygiene, and environmental cleanliness. Togbenon and Moyo [14] proposed and analyzed a mathematical model for cholera transfer as a control strategy with a quarantine class and vaccination parameter. Posny and Wang [15] proposed in a periodic environment, a deterministic compartmental model for the dynamics of cholera. Seasonal variation is incorporated into a broad formulation for the pathogen concentration and the incidence in their model. According to Fakai et al. [6], epidemics of cholera have been on the rise and more than 250,000 cholera cases are recorded worldwide annually. Flood, drought, and river height are all factors that influence the outbreak of cholera. Draught and flood, according to Codeco [3], are likely to have a complex impact on the dynamics of cholera. Flooding can wash infected feces and sewage into rivers, disrupting water delivery and worsening hygiene conditions. Nyabadza et al. [16] developed a deterministic mathematical model to investigate the mechanisms of cholera transmission in the face of scarce resources which includes nonlinear recovery rates. Lemos-paião et al. [17] researched an epidemic cholera model optimum treatment with control in Portugal. They proposed in the study a quarantine-treated cholera mathematics model which employs an ideal control problem and controls a fraction of infected individuals who will be treated in quarantine until their recovery will be completed, with a reduced number of infected persons. Peter et al. [18] developed a deterministic mathematical model to examine the degree of sensitivity of certain factors that aid in cholera transmission and management. Nyaberi and Malonza [19] researched the trends of transmission dynamics of cholera through health education and quarantine treatment as an epidemic control measure, modeling disease transmission using ordinary differential equations with an appropriate simple reproductive number, $R_0$, measured with the use of next-generation matrices. Abdulai [20] conducted research on the dynamics of the fractional order cholera transmission behavior model in Ghana with the Atanackovic and Stankovic numerical methods. Mwasa and Tchuenche [21] worked on a SIR model combining programs for public health awareness, vaccination, quarantine, and treatment as preventive measures to curb the illness.

Other researchers on cholera transmission dynamics include Crooks and Atesmachew and Lemaitre et al. [2, 22], to mention but a few.

According to Liang et al. [23], mathematical models as a tool in analyzing and predicting dynamical behavior in biological systems have been successfully used in the past decay. Therefore, in this study, a stochastic differential equation model would be used which provides an effective and efficient tool to unravel the role of the aquatic environment in the transmission dynamics and a better understanding of the spread of cholera infection even under uncertainties. In this study, we extend the deterministic model developed in [3] by converting it to a stochastic model. To understand the flow and prevent cholera infection in a better way, we first study the deterministic model by deriving the $R_0$, qualitative features such as positivity and invariant region of the solution, and stabilities of both the disease-free and endemic equilibriums to ensure the biological meaningfulness of the model. Next, we pass to the stochastic differential equation model and show the existence and uniqueness of the model. Finally, numerical simulations are carried out using the Euler-Maruyama scheme and analyzed with the aid of MATLAB, and the findings showed that the sample path of the stochastic differential equation for the number of infectious individuals is continuous but not differentiable.

The rest of the study is organized as follows: in Section 2, the Codeco [3] cholera model is described and reformulated. A qualitative analysis of the model is also discussed. In Section 3, a stochastic differential equation (SDEs) is formulated from transition probabilities. The existence and uniqueness of the stochastic differential equation model are also discussed. In Section 4, we use MATLAB software to investigate the numerical simulation results. Finally, in Section 5, we present our discussions and conclusions.

2. Model Formulation

The cholera model formulated by [3] is a system that includes both the environment and human population, with the environment-to-human transmission route represented by a logistic function. This model incorporated explicitly the environmental component, namely, the concentration of Vibrio cholerae bacteria in water $B$, into a regular SIR system to create an epidemiological model of a combined human-to-environment SI-B. In this study, we extend this model to a stochastic model by adding a “noise” to the deterministic model developed by [3]. Thus, we add a Brownian motion (Wiener’s process) and the intensity or impact of the stochastic environmental factors on the deterministic model.

2.1. The Deterministic Model Equations

In this model, we consider a deterministic compartmental human population and a Vibrio cholera bacteria population. The total human population is divided into three subclasses which are the susceptible population $S$, the cholera infectious population $I$, and a Vibrio cholerae bacteria population $B$. Susceptibles are recruited with the rate of $n$ through birth, susceptible can get cholera with the rate of $aB/K+B$, where $a$ is the contact rate with untreated water and $B/K+B$ is the probability of an individual to catch cholera. Infected individuals also contribute to the enhancement of the Vibrio cholera bacteria population through excretion at the rate of $e$. The infected individual may recover at the rate of $r$. In the aquatic environment, the bacteria population also grows at a rate determined by environmental factors, and b is the size of bacteria in the aquatic reservoir. A susceptible individual dies naturally at a rate of $\mu$. The above description of the model is plotted in Figure 1, while the variables and parameters are defined in Tables 1 and 2, respectively.

[figure(s) omitted; refer to PDF]

Table 1

Definition of variables in Codeco’s model.

Variable	Definition
$$S$$	Number of susceptible individuals
$$I$$	Number of infected individuals
$$B$$	Toxigenic concentration of Vibrio cholerae in water (cell/ml).

Table 2

Description of parameters in Codeco’s model.

Parameter	Definition
$$H$$	Total human population
$$n$$	Birth and death rates of humans (${\text{day}}^{-1}$)
$$a$$	Exposure rate of individuals to contaminated water (${\text{day}}^{-1}$)
$$K$$	Vibrio cholerae concentration in water that yields 50% of catching cholera (cell/ml)
$$r$$	Recovery rate of individuals (${\text{day}}^{-1}$)
$$\text{nb}-\text{mb}$$	Difference between the growth rate and loss rate of Vibrio cholerae in the aquatic reservoir (${day}^{-1}$)
$$e$$	Infected individual’s contribution to the Vibrio cholerae bacteria population in the aquatic reservoir (cell/ml ${\text{day}}^{-1}{\text{individual}}^{-1}$)

To extend the deterministic model in [3] to a stochastic one, we reformulate the basic model as follows.

The above assumptions of the model lead to the system of ordinary differential equations shown below as in [3]. $$\begin{array}{}\text{(1)}& \begin{array}{c}\frac{dS}{\text{dt}}=nH-S-\frac{a\text{BS}}{K+B}\hfill \\ \frac{dI}{\text{dt}}=\frac{a\text{BS}}{K+B}-rI\hfill \\ \frac{dB}{\text{dt}}=B\text{nb}-\text{mb}+eI\hfill \end{array},\end{array}$$with initial conditions $S0=S_0,I0=I_0,B0=B_0$.

2.2. Basic Properties of the Deterministic Model

For our model to make sense, it is necessary to show at least that this $\text{SI}-B$ model has a solution, and also the solution will remain within $0,H$ whenever it starts from there.

2.2.1. Positivity of the Solutions

In order for our model to be realistic, solutions will have to be nonnegative at all times for all $t\ge 0.$ We show that every state variable in the system equations will remain nonnegative.

Theorem 1.

If our initial values of the parameters are $S0\ge 0,I0\ge 0,B0\ge 0$, then the solution set $St,It,Bt$ is nonnegative for all $t\ge 0.$

Proof.

Let $t_{\ast }=\sup t>0:St>0,It>0,Bt>0$

Since $St$, $It$, and $Bt$ are continuous, we deduce that $t_{\ast }>0.$ If $t_{\ast }=+\infty ,$ then positivity holds but if $0<t_{\ast }<+\infty ,$$St=0$ or $It=0$, or $Bt=0$. Now, consider the first equation of model (1), $$\begin{array}{}\text{(2)}& \frac{dSt}{\text{dt}}=nH-\frac{aBtSt}{K+Bt}-nSt.\end{array}$$

Rewritten as $$\begin{array}{}\text{(3)}& \frac{dSt}{\text{dt}}+n+\lambda tSt=nH,\end{array}$$where $\lambda t=\alpha Bt/K+Bt$.

So, by integrating Equation (3), we have $$\begin{array}{}\text{(4)}& \frac{d}{\text{dt}}St{e}^{\text{nt}+\int \lambda \tau d\tau }=nH{e}^{\text{nt}+\int \lambda \tau d\tau }.\end{array}$$

Thus, $$\begin{array}{}\text{(5)}& S{t}_{\ast }{e}^{\text{n}{\text{t}}_{\ast }+{\int }_0^{t_{\ast }}\lambda \tau d\tau }-S0={\int }_0^{t_{\ast }}{e}^{\int n+\lambda t\text{dt}}.nH\text{dt}.\end{array}$$

This implies $$\begin{array}{c}\text{(6)}& S{t}_{\ast }=S0e^{-\text{n}{\text{t}}_{\ast }+{\int }_0^{t_{\ast }}\lambda \tau d\tau }+e^{-\text{n}{\text{t}}_{\ast }+{\int }_0^{t_{\ast }}\lambda \tau d\tau }\times {\int }_0^{t_{\ast }}{e}^{\int n+\lambda t\text{dt}}.nH\text{dt},S{t}_{\ast }=T_1\text{S}0+T_1{\int }_0^{t_{\ast }}{e}^{\int n+\lambda t\text{dt}}.nH\text{dt}>0,\end{array}$$where $T_1=e^{-\text{n}{\text{t}}_{\ast }+{\int }_0^{t_{\ast }}\lambda \tau d\tau }>0,$$S0>0,$ and from the definition of $t_{\ast }$ above, we have $Bt>0$. Therefore, the solution $S{t}_{\ast }>0$ and, hence, $S{t}_{\ast }\ne 0.$

Consider the second equation of model (1), $$\begin{array}{}\text{(7)}& \frac{dIt}{\text{dt}}=\frac{aBtSt}{K+Bt}-rIt,\end{array}$$which can be written as $$\begin{array}{}\text{(8)}& \frac{dIt}{\text{dt}}+rIt=\lambda tSt,\end{array}$$where $\lambda t=\alpha Bt/K+Bt$.

Also, by integrating Equation (8), we have $$\begin{array}{}\text{(9)}& \frac{d}{\text{dt}}It{e}^{\int \text{rdt}}=\lambda tSt{e}^{\int \text{rdt}}.\end{array}$$

Thus, $$\begin{array}{}\text{(10)}& I{t}_{\ast }{e}^{{\text{rt}}_{\ast }}-I0={\int }_0^{t_{\ast }}{e}^{\int \text{rdt}}.\lambda tSt\text{dt}.\end{array}$$

This implies, $$\begin{array}{c}\text{(11)}& I{t}_{\ast }=I0e^{-\text{r}{\text{t}}_{\ast }}+e^{-\text{r}{\text{t}}_{\ast }}{\int }_0^{t_{\ast }}{e}^{\int \text{rdt}}.\lambda tSt\text{dt},I{t}_{\ast }=T_{1}I0+T_1{\int }_0^{t_{\ast }}{e}^{\int \text{rdt}}.\lambda tSt\text{dt}>0,\end{array}$$where $T_1=e^{-\text{r}{\text{t}}_{\ast }}>0,I0>0,$ and from the above, $St>0:$ then, the solution $I{t}_{\ast }>0$ and, hence, $I{t}_{\ast }\ne 0.$

Consider the third equation of model (1), $$\begin{array}{c}\text{(12)}& \frac{dBt}{\text{dt}}=eIt-\text{mb}-\text{nb}Bt,\frac{dBt}{\text{dt}}+\text{mb}-\text{nb}Bt=eIt,\\ \frac{d}{\text{dt}}Bt{e}^{\int \text{mb}-\text{nb}dt}=e^{\int \text{mb}-\text{nb}dt}.eIt.\end{array}$$

Thus, $$\begin{array}{}\text{(13)}& B{t}_{\ast }{e}^{\text{mb}-\text{nb}{t}_{\ast }}-B0={\int }_0^{t_{\ast }}{e}^{\int \text{mb}-\text{nb}dt}.eIt\text{dt}.\end{array}$$

This implies, $$\begin{array}{c}\text{(14)}& B{t}_{\ast }=B0e^{-\text{mb}-\text{nb}{t}_{\ast }}+e^{-\text{mb}-\text{nb}{t}_{\ast }}{\int }_0^{t_{\ast }}{e}^{\int \text{mb}-\text{nb}\text{dt}}.eIt\text{dt},B{t}_{\ast }=k_{1}B0+k_1{\int }_0^{t_{\ast }}{e}^{\int \text{mb}-\text{nb}\text{dt}}.eIt\text{dt}>0,\end{array}$$where $k_1=e^{-\text{mb}-\text{nb}{t}_{\ast }}>0$ for $\text{mb}-\text{nb}>0$, $B0>0$ and from the above, $It>0;$ then the solution $B{t}_{\ast }>0,$ and, hence, $B{t}_{\ast }\ne 0.$

Based on the definition, $t_{\ast }$ is not finite which means $t_{\ast }=+\infty$; therefore, the solutions of the Codeco model system (1) are always positive.

2.2.2. Invariant Region for the Deterministic Model

In this subsection, we obtain the invariant region of the model Equation (1).

Theorem 2.

The solutions for the model system (1) are contained and remain in the region $\Omega =S,I,BS\ge 0,I\ge 0,0\le S+I\le \Delta /n,0\le B\le e\Delta /\text{nq}$, for all time $t\ge 0$.

Proof.

Consider the total human population, $Ht=St+It$, its time derivative satisfies, $$\begin{array}{c}\text{(15)}& \frac{dHt}{\text{dt}}=\frac{dSt}{\text{dt}}+\frac{dIt}{\text{dt}},\frac{dHt}{\text{dt}}=nHt-\frac{aBtSt}{K+Bt}-nSt+\frac{aBtSt}{K+Bt}-rIt.\end{array}$$

Let $nH$ the rate of recruitment into the susceptible class be $\Delta$ and $S=H-I.$$$\begin{array}{c}\text{(16)}& \frac{dHt}{\text{dt}}=\Delta -nHt-nIt-rIt,\frac{dHt}{\text{dt}}\le \Delta -nHt.\end{array}$$

By integrating, we obtain $$\begin{array}{}\text{(17)}& Ht\le \frac{\Delta }{n}+\text{k}{\text{e}}^{-nt},\end{array}$$where $k$ is a constant. Initially, at $t=0,$$H0-\Delta /n\le k.$

Therefore, $$\begin{array}{}\text{(18)}& Ht\le \frac{\Delta }{n}+H0-\frac{\Delta }{n}{e}^{-\text{nt}}.\end{array}$$

Thus,

$\underset{t⟶\infty }{\lim }Ht\le \Delta /n$, which implies that $0\le Ht\le \Delta /n$

Similarly, $$\begin{array}{}\text{(19)}& \frac{dBt}{\text{dt}}=eI-qB\le e\frac{\Delta }{n}-qB,\end{array}$$where $q=\text{mb}-\text{nb}.$

By integrating, we obtain $$\begin{array}{}\text{(20)}& Bt\le \frac{e\Delta }{\text{nq}}+\text{k}{\text{e}}^{-qt},\end{array}$$where $k$ is a constant. Initially, at $t=0,$$$\begin{array}{}\text{(21)}& B0-\frac{e\Delta }{\text{nq}}\le k.\end{array}$$

Therefore, $$\begin{array}{}\text{(22)}& Bt\le \frac{e\Delta }{\text{nq}}+B0-\frac{e\Delta }{\text{nq}}{e}^{-\text{qt}}.\end{array}$$

Thus,

$\underset{t⟶\infty }{\lim }Bt\le e\Delta /\text{nq},$ which implies that $0\le Bt\le e\Delta /\text{nq}$.

Therefore, the feasibility of the solution set of the system of Equation (1) is in the region $\Omega =S,I,BS\ge 0,I\ge 0,0\le S+I\le \Delta /n,0\le B\le e\Delta /\text{nq}$; hence, the model is well posed and biologically meaningful.

2.3. Existence and Stability of the Disease-Free Equilibrium

Here, the existence of the equilibrium state of the model is discussed.

We set $dS/\text{dt}=dI/\text{dt}=dB/\text{dt}=0$ at the equilibrium state.

Therefore, by setting the model equations of the system (1) to zero, we have $$\begin{array}{}\text{(23)}& nH-S-\frac{a\text{BS}}{K+B}=0,\text{(24)}& \frac{a\text{BS}}{K+B}-rI=0,\text{(25)}& \text{nb}-\text{mb}B+eI=0.\end{array}$$

There are no infections at the disease-free state; thus, $I=0$. Substituting this into Equation (25), we have $\text{nb}-\text{mb}B=0⟹B=0$ provided $\text{nb}-\text{mb}\ne 0.$ Therefore, putting $B=0$ into Equation (23), we have $nH-S=0⟹n=0$ and $S=H.$

Therefore, a disease-free equilibrium state exists and is given by $E_{0}H,0,0.$

Now, we analyze the disease-free equilibrium’s stability. Let us consider taking the Jacobian of the model Equations (23)–(25) to prove that they are locally asymptotically stable around the equilibrium point as in [19].

When it comes to disease spread, local asymptotic stability means that if there is a small change or perturbation on the system, the system will still return to the disease-free equilibrium.

2.4. Basic Reproductive Number for the Deterministic Model

The basic reproductive number $R_0$ is calculated from the disease compartments as follows:

$\mathcal{F}:$ rate of appearance of new infection$=aBS/K+B/0$

$V_i^{-}:$ transfer rate of disease out of the disease compartment $=\begin{array}{c}rI\\ 0\end{array}$

$V_i^{+}:$ transfer rate of infection into the disease compartment by other means $$\begin{array}{c}\text{(26)}& =\begin{array}{c}0\\ \text{nb}-\text{mb}B+eI\end{array},V_i=V_i^{-}-V_i^{+}=\begin{array}{c}rI\\ \text{mb}-\text{nb}B-eI\end{array}.\end{array}$$

$F=\partial F{E}_0/\partial x_j:$ the Jacobian of $\mathcal{F}$ with respect to disease compartments $x_j$ evaluated at disease-free equilibrium. The model’s disease-free equilibrium point is obtained by setting $$\begin{array}{c}\text{(27)}& \frac{dSt}{\text{dt}}=\frac{dIt}{\text{dt}}=\frac{dBt}{\text{dt}}=0,nHt-\frac{aBtSt}{K+Bt}-nSt=0,\\ \text{(28)}& \frac{aBtSt}{K+Bt}-rIt=0,\text{nb}-\text{mb}Bt+eIt=0.\end{array}$$

Now, from system above, we have, $$\begin{array}{c}\text{(29)}& S^{\ast }=\frac{nH}{a{B}^{\ast }/K+B^{\ast }+n},I^{\ast }=\frac{a{B}^{\ast }{S}^{\ast }}{rK+B^{\ast }},\\ B^{\ast }=\frac{e{I}^{\ast }}{\text{nb}-\text{mb}}.\end{array}$$

From the equations above,

$I^{\ast }=0$ and $\text{nb}-\text{mb}{B}^{\ast }=0⟹B^{\ast }=0$ provided $\text{nb}-\text{mb}\ne 0.$ Therefore, putting $B^{\ast }=0$ into $S^{\ast }=nH/a{B}^{\ast }/K+B^{\ast }+n$, we have $S^{\ast }=H$.

Therefore, $E_0=St,It,Bt=H,0,0$. $$\begin{array}{}\text{(30)}& F=\begin{array}{cc}\frac{\partial m_1}{\partial I}& \frac{\partial m_1}{\partial B}\\ \frac{\partial m_2}{\partial I}& \frac{\partial m_2}{\partial B}\end{array},\end{array}$$where $$\begin{array}{c}\text{(31)}& m_1=\frac{aBtSt}{K+Bt}-rIt.m_2=nb-mbBt+eIt\end{array}$$

Then, the Jacobian matrix is given by the following equation: $$\begin{array}{}\text{(32)}& F=\begin{array}{cc}0& \frac{aSK}{{K+B}^2}\\ 0& 0\end{array}.\end{array}$$

At disease-free equilibrium, $$\begin{array}{c}\text{(33)}& S=H,B=0,\\ F=\begin{array}{cc}0& \frac{aH}{K}\\ 0& 0\end{array}.\end{array}$$

$V=\partial V_i{E}_0/\partial x_j:$ the Jacobian of $V_i$ with respect to disease compartments $x_j$ evaluated at DFE, $$\begin{array}{c}\text{(34)}& V=\begin{array}{cc}r& 0\\ -e& \text{mb}-\text{nb}\end{array},V^{-1}=\frac{1}{r\text{mb}-\text{nb}}\begin{array}{cc}\text{mb}-\text{nb}& 0\\ e& r\end{array},\\ V^{-1}=\begin{array}{cc}\frac{1}{r}& 0\\ \frac{e}{r\text{mb}-\text{nb}}& \frac{1}{\text{mb}-\text{nb}}\end{array}.\end{array}$$

Therefore, $$\begin{array}{c}\text{(35)}& \text{F}{\text{V}}^{-1}=\begin{array}{cc}0& \frac{aH}{K}\\ 0& 0\end{array}\begin{array}{cc}\frac{1}{r}& 0\\ \frac{e}{r\text{mb}-\text{nb}}& \frac{1}{\text{mb}-\text{nb}}\end{array},\text{F}{\text{V}}^{-1}=\begin{array}{cc}\frac{eaH}{Kr\text{mb}-\text{nb}}& \frac{aH}{K\text{mb}-\text{nb}}\\ 0& 0\end{array}.\end{array}$$

Now, we calculate the eigenvalues of the matrix $\begin{array}{cc}eaH/Kr\text{mb}-\text{nb}& aH/K\text{mb}-\text{nb}\\ 0& 0\end{array}$ to determine the basic reproduction number, $R_0$ which is defined as the spectral radius or the dominant eigenvalue of the matrix as in [7]. We compute this by $\begin{array}{cc}eaH/Kr\text{mb}-\text{nb}& aH/K\text{mb}-\text{nb}\\ 0& 0\end{array}-\text{I}\lambda =0$, where $I$ is a $2\times 2$ identity matrix. Hence, $$\begin{array}{c}\text{(36)}& \begin{array}{cc}\frac{eaH}{Kr\text{mb}-\text{nb}}& \frac{aH}{K\text{mb}-\text{nb}}\\ 0& 0\end{array}-\begin{array}{cc}\lambda & 0\\ 0& \lambda \end{array}=0,\begin{array}{cc}\frac{eaH}{Kr\text{mb}-\text{nb}}-\lambda & \frac{aH}{K\text{mb}-\text{nb}}\\ 0& -\lambda \end{array}=0,\\ \frac{eaH}{Kr\text{mb}-\text{nb}}-\lambda \lambda =0.\end{array}$$

Either $eaH/Kr\text{mb}-\text{nb}-\lambda =0$ or $\lambda =0\Rightarrow {\lambda }_1=eaH/Kr\text{mb}-\text{nb}$ and ${\lambda }_2=0$

The eigenvalues of ${\text{FV}}^{-1}$ are $eaH/Kr\text{mb}-\text{nb},0.$ Clearly, ${\lambda }_1$ is the dominant eigenvalue and becomes the basic reproduction number $R_0$ of the model as in [7]. $$\begin{array}{}\text{(37)}& R_0=\frac{eaH}{Kr\text{mb}-\text{nb}}.\end{array}$$

With respect to the cholera disease, the basic reproduction number ($R_0$) describes the expected number of cholera infections generated by one cholera case in a susceptible population [7].

Theorem 3.

The disease-free equilibrium points of the system (1) is asymptotically stable if and only if $R_0<1$ [19].

Proof.

The Jacobian matrix is given by $$\begin{array}{}\text{(38)}& J{S}_0,I_0,B_0=\begin{array}{ccc}-n-\frac{aB}{K+B}& 0& -\frac{aSK}{{K+B}^2}\\ \frac{aB}{K+B}& -r& \frac{aSK}{{K+B}^2}\\ 0& e& \text{nb}-\text{mb}\end{array}.\end{array}$$

Now, at disease-free equilibrium point, $S_0,I_0,B_0=H,0,0$. $$\begin{array}{}\text{(39)}& J{S}_0,I_0,B_0=\begin{array}{ccc}-n& 0& -\frac{aH}{K}\\ 0& -r& \frac{aH}{K}\\ 0& e& \text{nb}-\text{mb}\end{array}.\end{array}$$

The equilibrium point is asymptotically stable if the following condition (Routh-Hurwitz) holds for polynomial $P$ and its determinant.

The Jacobian matrix has characteristic equation given as, $$\begin{array}{c}\text{(40)}& PJ{S}_0,I_0,B_0-\lambda I=\begin{array}{ccc}-n-\lambda & 0& -\frac{aH}{K}\\ 0& -r-\lambda & \frac{aH}{K}\\ 0& e& \text{nb}-\text{mb}-\lambda \end{array}=0,-n-\lambda \begin{array}{cc}-r-\lambda & \frac{aH}{K}\\ e& \text{nb}-\text{mb}-\lambda \end{array}-0\begin{array}{cc}0& \frac{aH}{K}\\ 0& \text{nb}-\text{mb}-\lambda \end{array}-\frac{aH}{K}\begin{array}{cc}0& -r-\lambda \\ 0& e\end{array}=0,\\ -n-\lambda -r-\lambda \text{nb}-\text{mb}-\lambda -\frac{aeH}{K}=0,\\ -n-\lambda {\lambda }^2+\text{mb}-\text{nb}+r\lambda +r\text{mb}-\text{nb}-\frac{\text{ae}H}{K}=0,\\ {\lambda }^3+n+\text{mb}-\text{nb}+r{\lambda }^2-n+r\text{nb}-\text{mb}-\text{rn}+\frac{\text{ae}H}{K}\lambda -\text{rn}\text{nb}-\text{mb}+\frac{\text{nae}H}{K}=0,\\ \text{(41)}& {\lambda }^3+n+\text{mb}-\text{nb}+r{\lambda }^2+n+r\text{mb}-\text{nb}+rn-\frac{\text{ae}H}{K}\lambda +\text{rn}\text{mb}-\text{nb}-\frac{naeH}{K}=0.\end{array}$$

It is observed that the Jacobian matrix characteristic equation has three roots.

We can rewrite the characteristic Equation (40) as $P\lambda ={\lambda }^3+A{\lambda }^2+B\lambda +C=0$, where $$\begin{array}{c}\text{(42)}& A=r+n+\text{mb}-\text{nb},B=\text{nr}+r+n\text{mb}-\text{nb}-\frac{\text{ae}H}{K}=n\text{mb}-\text{nb}+r+r\text{mb}-\text{nb}1-R_0,\\ C=\text{nr}\text{mb}-\text{nb}-\frac{\text{nae}H}{K}=\text{nr}\text{mb}-\text{nb}1-R_0,\\ \text{AB}-C=\text{mb}-\text{nb}+n+rn\text{mb}-\text{nb}+r+r\text{mb}-\text{nb}1-R_0.\end{array}$$

Therefore, the eigenvalues that correspond to the equilibrium $E_0$ are $$\begin{array}{c}\text{(43)}& {\lambda }_1=-n,{\lambda }_{2,3}=\frac{-\text{mb}-\text{nb}+r\pm \sqrt{{\text{mb}-\text{nb}+r}^2-4r\text{mb}-\text{nb}1-R_0}}{2}.\end{array}$$

We have $A>0$, $B>0$, $C>0$ and $AB-C>0$ and according to the Routh-Hurwitz criterion, all the characteristic equation’s roots have negative real part when $R_0<1$ and $\text{mb}-\text{nb}>0$. Hence, the disease-free equilibrium (DFE) point is stable given that $R_0<1$ [19].

2.5. Existence and Stability of the Endemic Equilibrium

$I\ne 0$ will be used to find the endemic state of the system of Equation (1). We now investigate if the endemic equilibrium condition of the system is stable. $$\begin{array}{}\text{(44)}& nH-S-\frac{a\text{BS}}{K+B}=0,\text{(45)}& \frac{a\text{BS}}{K+B}-rI=0,\text{(46)}& \text{nb}-\text{mb}B+eI=0,\end{array}$$

From Equation (46), $$\begin{array}{}\text{(47)}& B=\frac{eI}{\text{mb}-\text{nb}}.\end{array}$$

From Equation (45), $$\begin{array}{}\text{(48)}& S=\frac{K+BrI}{aB}.\end{array}$$

From Equation (44), $$\begin{array}{}\text{(49)}& nH-Sn+\frac{aB}{K+B}=0.\end{array}$$

Substitute Equation (48) into Equation (49). $$\begin{array}{c}\text{(50)}& nH-\frac{B+KrI}{aB}n+\frac{aB}{B+K}=0,nH-rI\frac{nB+K}{aB}+1=0.\end{array}$$

Substitute Equation (47) into Equation (50). $$\begin{array}{c}\text{(51)}& nH-rI\frac{K+eI/\text{mb}-\text{nb}}{aeI/\text{mb}-\text{nb}}n+1=0,\text{nae}H-\text{nr}K\text{mb}-\text{nb}-\text{er}In+a=0,\\ I=\frac{\text{nae}H-\text{nr}Kmb-nb}{\text{er}n+a}.\end{array}$$

Supposedly, for $I>0,$$$\begin{array}{c}\text{(52)}& \frac{\text{nae}H-\text{nr}K\text{mb}-\text{nb}}{\text{er}n+a}>0,\frac{\text{ae}H}{rK\text{mb}-\text{nb}}>1.\end{array}$$

Therefore, $R_0>1$, where $R_0$ is the reproduction number and is written as $R_0=\text{ae}H/rK\text{mb}-\text{nb}.$

We observed from the analysis that if $R_0>1$ and $\text{mb}-\text{nb}>0$, there exists a positive endemic equilibrium.

Theorem 4.

The endemic equilibrium state $E^{\ast }$ is locally asymptotically stable if $R_0>1.$

Proof.

To confirm the stability of the endemic equilibrium state, we used the system of Equation (1), so substituting $B^{\ast }$ for $B$, $S^{\ast }$ for $S$ and taking $E^{\ast }=S^{\ast },I^{\ast },B^{\ast }$ as the endemic equilibrium, the Jacobian matrix of the system Equation (1) becomes, $$\begin{array}{}\text{(53)}& J{E}^{\ast }=\begin{array}{ccc}-n-W& 0& -U\\ W& -r& U\\ 0& e& -V\end{array},\end{array}$$where $$\begin{array}{c}\text{(54)}& W=\frac{{aB}^{\ast }}{K+B^{\ast }}=\frac{na{R}_0-1}{a+n{R}_0},U=\frac{a{S}^{\ast }K}{{K+B^{\ast }}^2}=\frac{ra+n\text{mb}-\text{nb}}{ea+n{R}_0},\\ V=\text{mb}-\text{nb}.\end{array}$$

The characteristic polynomial of the matrix $J{E}^{\ast }$ is $$\begin{array}{}\text{(55)}& \det J{E}^{\ast }-\lambda I=a_0{\lambda }^3+a_1{\lambda }^2+a_2\lambda +a_3,\end{array}$$where $$\begin{array}{c}\text{(56)}& a_0=1,a_1=n+r+W+V,\\ a_2=\text{nr}+nV+rW+\text{VW}-eU,\\ a_3=\text{nr}V+r\text{VW}-\text{ne}U.\end{array}$$

The Routh-Hurwitz criterion [19] requires $a_0>0,$$a_1>0,$$a_2>0,$$a_3>0$, and $a_1{a}_2-a_0{a}_3>0$ as the necessary conditions for the endemic equilibrium state to be locally asymptotically stable, i.e, all roots of the characteristic polynomial above have negative real parts. Clearly $a_0>0,$$a_1>0,$ since $n,$$a>0$ and $R_0>1$. Similarly, $a_2,a_3>0,$ since $n,$$e,$$r,$$a,$$\text{mb}-\text{nb}>0$ and $R_0>1$ if and only if $\text{nr}+nV+rW+\text{VW}>eU$ and $nrV+r\text{VW}>\text{ne}U.$

Now, it remains to show that $a_1{a}_2-a_0{a}_3>0.$ Thus, $$\begin{array}{}\text{(57)}& n+r+W+V\text{nr}+nV+rW+\text{WV}-eU-\text{nr}V+r\text{WV}-\text{ne}U=n^{2}r+n^{2}V+nrW+n\text{WV}+\text{n}{\text{r}}^2+r^{2}W+nrW+n\text{WV}+r{W}^2+W^{2}V+\text{nr}V+n{V}^2+r\text{WV}+\text{W}{\text{V}}^2-\text{re}U-eWU-eVU>0.\end{array}$$

Therefore, when $R_0>1$, the model system (1) has a unique endemic equilibrium $E^{\ast }$, and it is stable if $n^{2}r+n^{2}V+nrW+n\text{WV}+n{r}^2+r^{2}W+\text{nr}W+n\text{WV}+r{W}^2+W^{2}V+\text{nr}V+n{V}^2+r\text{WV}+\text{W}{\text{V}}^2>reU+e\text{WU}+e\text{VU}$, and $\text{mb}-\text{nb}>0$.

3. Transition Probabilities

Allen [24] developed the Ito stochastic differential equations from transition probabilities approach, which is based on the diffusion process. The stochastic differential equation to be formed is in the form, $$\begin{array}{}\text{(58)}& dXt=fXt,tdt+gXt,tdWt,\end{array}$$where $fXt,t$ is the drift (deterministic part of the model), $gXt,t$ is the diffusion part given as $gXt,t=V^{1/2}$, $V$ is the covariance to order $\Delta t$, $V\Delta t$ is the approximate covariance matrix, and $Wt$ is the vector of independent Wiener’s process.

From the deterministic equations above, we formulate the stochastic differential equation as $$\begin{array}{}\text{(59)}& dXt=fXt,tdt+gXt,tdWt,\end{array}$$where $fXt,t={\rm E}\Delta X/\Delta t$, and $gXt,t=\sqrt{V}=\sqrt{{\rm E}\Delta X{\Delta X}^T/\Delta t}$ [24].

Throughout this study, we assume that time is a continuous variable and the state variables $St$ and $It$ are continuous random variables. Now, let $\Delta S=St+\Delta t-St$ and $\Delta I=It+\Delta t-It.$ In addition, we also assume that the change (transition) of random variables $St$ and $It$ is approximately normally distributed, $\Delta St~N\mu s\Delta t,{\sigma }^{2}s\Delta t$ and $\Delta It~N\mu I\Delta t,{\sigma }^{2}I\Delta t$ for small time intervals $\Delta t$.

Let $Xt={X_1,X_2}^T$, where $X_1$ and $X_2$ correspond to the number of individuals $St,It\in 0,H$, respectively, and $Bt$ is the concentration of Vibrio cholerae bacteria in the aquatic environment. But $Bt$ was not considered as it is not compartment occupancy as the susceptible and the infected; as a result, its transition is not taken into consideration in formulating the transitions.

The following are needed in forming a stochastic differential equation model: the expectation ${\rm E}\Delta X$ and the covariance ${\rm E}\Delta X{\Delta X}^T$ which is a matrix. To compute the expectation and the covariance matrix, the possible changes or the transitions along with their associated probabilities are first computed as shown in Table 3 below. The transition probabilities are formulated from the deterministic model.

Table 3

Transition probabilities.

Possible changes	Probabilities	Description
$$\Delta X_1={10}^T$$	$$P_1=nH\Delta t$$	Birth of a susceptible
$$\Delta X_2={-11}^T$$	$$P_2=\frac{a\text{BS}}{K+B}\Delta t$$	Susceptible becomes infected
$$\Delta X_3={-10}^T$$	$$P_3=nS\Delta t$$	Susceptible dies natural death
$$\Delta X_4={0-1}^T$$	$$P_4=rI\Delta t$$	Infected recovers

In Table 3, the expectation is computed as follows: $$\begin{array}{}\text{(60)}& {\rm E}\Delta X=\sum _{i=1}^4{P}_i\Delta X_i=P_1\Delta X_1+P_2\Delta X_2+P_3\Delta X_3+P_4\Delta X_4,\end{array}$$substituting the values of $P_i,\Delta X_i$ and ${X_1,X_2}^T={St,It}^T$ into Equation (60), $$\begin{array}{}\text{(61)}& {\rm E}\Delta X=P_1\begin{array}{c}1\\ 0\end{array}+P_2\begin{array}{c}-1\\ 1\end{array}+P_3\begin{array}{c}-1\\ 0\end{array}+P_4\begin{array}{c}0\\ -1\end{array}=\begin{array}{c}{P}_1\\ 0\end{array}+\begin{array}{c}-P_2\\ P_2\end{array}+\begin{array}{c}-P_3\\ 0\end{array}+\begin{array}{c}0\\ -P_4\end{array}=\begin{array}{c}{P}_1-P_2-P_3\\ P_2-P_4\end{array}=\begin{array}{c}nH\Delta t-\frac{aBS}{K+B}\Delta t-nS\Delta t\\ \frac{aBS}{K+B}\Delta t-rI\Delta t\end{array}\therefore {\rm E}\Delta X=\begin{array}{c}nH-S-\frac{aBS}{K+B}\\ \frac{aBS}{K+B}-rI\end{array}\Delta t.\end{array}$$

And the covariance matrix is also computed as follows: $$\begin{array}{}\text{(62)}& {\rm E}\Delta X{\Delta X}^T=\sum _{i=1}^4{P}_i\Delta X_i{{\Delta X}_i}^T=P_1\Delta X_1{\Delta X_1}^T+P_2\Delta X_2{\Delta X_2}^T+P_3\Delta X_3{\Delta X_3}^T+P_4\Delta X_4{\Delta X_4}^T,\end{array}$$substituting the values of $P_i,\Delta X_i$ and ${X_1,X_2}^T={St,It}^T$ into Equation (62), $$\begin{array}{}\text{(63)}& {\rm E}\Delta X{\Delta X}^T=P_1\begin{array}{c}1\\ 0\end{array}1\hspace{1em}0+P_2\begin{array}{c}-1\\ 1\end{array}-1\hspace{1em}1+P_3\begin{array}{c}-1\\ 0\end{array}-1\hspace{1em}0+P_4\begin{array}{c}0\\ -1\end{array}0\hspace{1em}-1=P_1\begin{array}{cc}1& 0\\ 0& 0\end{array}+P_2\begin{array}{cc}1& -1\\ -1& 1\end{array}+P_3\begin{array}{cc}1& 0\\ 0& 0\end{array}+P_4\begin{array}{cc}0& 0\\ 0& 1\end{array}=\begin{array}{cc}{P}_1+P_2+P_3& -P_2\\ -P_2& P_2+P_4\end{array}=\begin{array}{cc}nH\Delta t+\frac{a\text{BS}}{K+B}\Delta t+nS\Delta t& -\frac{a\text{BS}}{K+B}\Delta t\\ -\frac{a\text{BS}}{K+B}\Delta t& \frac{a\text{BS}}{K+B}\Delta t+rI\Delta t\end{array}\therefore {\rm E}\Delta X{\Delta X}^T=\begin{array}{cc}nH+S+\frac{a\text{BS}}{K+B}& -\frac{a\text{BS}}{K+B}\\ -\frac{a\text{BS}}{K+B}& \frac{a\text{BS}}{K+B}+rI\end{array}\Delta t.\end{array}$$

In Equation (61), $$\begin{array}{}\text{(64)}& \frac{E\Delta X}{\Delta \text{t}}=\begin{array}{c}nH-S-\frac{a\text{BS}}{K+B}\\ \frac{a\text{BS}}{K+B}-rI\end{array}=fXt,t,\end{array}$$

and in Equation (63), $$\begin{array}{}\text{(65)}& \frac{E\Delta X{\Delta X}^T}{\Delta t}=\begin{array}{cc}nH+S+\frac{aBS}{K+B}& -\frac{a\text{BS}}{K+B}\\ -\frac{a\text{BS}}{K+B}& \frac{a\text{BS}}{K+B}+rI\end{array}=VXt,t.\end{array}$$

Now, we compute $V^{1/2}$.

According to Allen [24], in a 2-dimensional system $g=V^{1/2}$ can be computed exactly, and it is computed as follows: $$\begin{array}{}\text{(66)}& g=V^{1/2}=\frac{1}{\beta }\begin{array}{cc}\delta +\theta & \rho \\ \rho & \omega +\theta \end{array},\end{array}$$where $\theta =\sqrt{\delta \omega -{\rho }^2}$ and $\beta =\sqrt{\delta +\omega +2\theta }$ with $$\begin{array}{c}\text{(67)}& \delta =nH+S+\frac{a\text{BS}}{K+B},\rho =-\frac{\text{aBS}}{\text{K}+\text{B}},\\ \omega =\frac{a\text{BS}}{K+B}+rI.\end{array}$$

Therefore, $\theta =\sqrt{nH+S+a\text{BS}/K+Ba\text{BS}/K+B+rI-{a\text{BS}/K+B}^2}.$$$\begin{array}{}\text{(68)}& \theta =\sqrt{nH+S\frac{a\text{BS}}{K+B}+rI\frac{a\text{BS}}{K+B}+nH+SrI},\end{array}$$and $\beta =\sqrt{nH+S+2a\text{BS}/K+B+rI+2\theta }$.

Therefore, the stochastic differential equation for the dynamics of the cholera infection is given as follows: $dXt=fXt,tdt+gXt,tdWt$ with the initial condition $X0=X_0$ and $Wt={W_{1}t,W_{2}t}^T,W_{1}t,\text{and}{W}_{2}t$ represent independent Wiener’s process. For each compartment, the following differential equations are obtained: $$\begin{array}{}\text{(69)}& \begin{array}{c}dSt=nH-S-\frac{a\text{BS}}{K+B}dt+\frac{\delta +\theta }{\beta }d{W}_{1}t+\frac{\rho }{\beta }d{W}_{2}t\hfill \\ dIt=\frac{a\text{BS}}{K+B}-rIdt+\frac{\rho }{\beta }d{W}_{1}t+\frac{\omega +\theta }{\beta }d{W}_{2}t\hfill \end{array}.\end{array}$$

The stochastic differential equations above are known as stochastic differential equations SI-B model.

Now, we model the populations each in system (69) to a single dimension Brownian motion (Wiener’s process). Thus, we write equations of system (69) in their simplified form.

Let us integrate the first equation of system (69), $$\begin{array}{c}\text{(70)}& {\int }_0^{t}dSs={\int }_0^{t}nH-S-\frac{a\text{BS}}{K+B}\text{ds}+{\int }_0^t\frac{\delta s+\theta s}{\beta s}d{W}_{1}s+{\int }_0^t\frac{\rho s}{\beta s}d{W}_{2}s,St={\int }_0^{t}nH-S-\frac{a\text{BS}}{K+B}\text{ds}+{\int }_0^t\frac{\delta s+\theta s}{\beta s}d{W}_{1}s+{\int }_0^t\frac{\rho s}{\beta s}d{W}_{2}s.\end{array}$$

Now, we define $$\begin{array}{}\text{(71)}& Mt={\int }_0^t\frac{\delta s+\theta s}{\beta s}d{W}_{1}s+{\int }_0^t\frac{\rho s}{\beta s}d{W}_{2}s.\end{array}$$

According to Greenhalgh et al. [25], the above is a martingale in terms of filtration and can be written in its quadratic variation as follows: $$\begin{array}{}\text{(72)}& Mt={\int }_0^t\frac{{\delta s+\theta s}^2}{\beta s^2}\text{ds}+{\int }_0^t\frac{\rho s^2}{\beta s^2}ds={\int }_0^t\frac{\rho s^2+\delta s^2+2\delta s\theta s+\theta s^2}{\beta s^2}\text{ds}={\int }_0^t\frac{\begin{array}{c}{-a\text{BS}/K+B}^2+{nH+S+a\text{BS}/K+B}^2+2nH+S+a\text{BS}/K+B\theta +\\ {\sqrt{nH+Sa\text{BS}/K+B+rIa\text{BS}/K+B+nH+SrI}}^2\end{array}}{{\sqrt{nH+S+2a\text{BS}/K+B+rI+2\theta }}^2}\text{ds}={\int }_0^t\frac{nH+S+a\text{BS}/K+BnH+S+2a\text{BS}/K+B+rI+2\theta }{nH+S+2a\text{BS}/K+B+rI+2\theta }\text{ds}={\int }_0^{t}nH+S+\frac{a\text{BS}}{K+B}\text{ds}.\end{array}$$

Martingale representation theorem allows the above equation to be written as Ito integral [25] in terms of Brownian motion as $$\begin{array}{}\text{(73)}& Mt={\int }_0^t\sqrt{nH+S+\frac{a\text{BS}}{K+B}}\text{dw}s.\end{array}$$

As a result, $$\begin{array}{}\text{(74)}& dSt=nH-St-\frac{aBtSt}{K+Bt}\text{dt}+\sqrt{nH+St+\frac{aBtSt}{K+Bt}}\text{dw}t.\end{array}$$

Similarly, the same process could be applied to the second equation of system (69) to get the following: $$\begin{array}{c}\text{(75)}& dIt=\frac{aBtSt}{K+Bt}-rIt\text{dt}+\sqrt{\frac{aBtSt}{K+Bt}+rIt}d{W}_{3}t,dIt=\frac{aBtH-It}{K+Bt}-rIt\text{dt}+\sqrt{\frac{aBtH-It}{K+Bt}+rIt}d{W}_{3}t.\end{array}$$

Hence, the system of stochastic differential Equations (74) and (75) describe how the susceptible and the infected population change with respect to time for $Ht>0.$ However, this system can be more simply described by the stochastic differential equatio. $$\begin{array}{}\text{(76)}& dIt=\frac{aBtH-It}{K+Bt}-rItdt+\sqrt{\frac{aBtH-It}{K+Bt}+rIt}d{W}_{3}t,\end{array}$$where $St+It=Ht>0$.

3.1. Existence and Uniqueness for the Stochastic Differential Equations

In this section, in order for the stochastic differential equation model (74) and (75) to make sense, we need to show at least that this model does not only have a unique solution but also exist.

Assume that the coefficients in the system of stochastic differential equation, $$\begin{array}{}\text{(77)}& dXt=f_{i}Xt,tdt+\sum _{i=1}^n\sum _{j=1}^m{g}_{ij}Xt,tdWt,\end{array}$$where $Xt={X_{1}t,X_{2}t}^T,Wt={W_{1}t,W_{2}t}^T,f_{i}Xt,t$ is a 2-dimensional vector with entries $f_{i}x,t$ and $g_{\text{ij}}Xt,t$ is a $2\times 2$ matrix with entries $g_{\text{ij}}x,t$ satisfy the following Lipschitz and growth conditions in the equations below for some constant $k<\infty ,$ and for all $t\in \mathfrak{R}$ and $x,y\in {\mathfrak{R}}^2$ as in Ogunlade et al. [26]. $$\begin{array}{c}\text{(78)}& \begin{array}{c}fx,t-fy,t\le kx-y\hfill \\ gx,t-gy,t\le kx-y\hfill \\ \begin{array}{c}fx,t\le kx\hfill \\ gx,t\le kx\hfill \end{array}\hfill \end{array},g=\sqrt{\sum _{i=1}^n\sum _{j=1}^m{g}_{ij}{x}^2},\\ f=\sqrt{\sum _{i=1}^n{f}_i{x}^2,}\end{array}$$

Then, for each $x_0\in {\mathfrak{R}}^2$, the system of stochastic differential Equation (77) has a unique solution in which $X0=x_0.$

Now, consider Equations (69) and (74), $$\begin{array}{c}\text{(79)}& dSt=nH-St-\frac{aBtSt}{K+Bt}dt+\sqrt{nH+St+\frac{aBtSt}{K+Bt}}\text{dw}t,dIt=\frac{aBtH-It}{K+Bt}-rIt\text{dt}+\sqrt{\frac{aBtH-It}{K+Bt}+rIt}d{W}_{3}t,\end{array}$$with the following: $$\begin{array}{}\text{(80)}& X0=S0,I0,B0,\text{(81)}& f_1=nH-St-\frac{aBtSt}{K+Bt},\text{(82)}& f_2=\frac{aBtH-It}{K+Bt}-rIt.\end{array}$$

Then, a constant $M>0$ exists such that $$\begin{array}{c}\text{(83)}& \frac{\partial f_1}{\partial S}=-n-\frac{aB}{K+B}\le M,\frac{\partial f_1}{\partial I}=0,\\ \frac{\partial f_2}{\partial S}=0,\\ \frac{\partial f_2}{\partial I}=-\frac{aB}{K+B}-r\le M.\end{array}$$

The diffusion matrix’s elements are continuously differentiable.

Also, for Equations (74) and (75) as a system, $$\begin{array}{}\text{(84)}& f=\sqrt{\sum _{i=1}^2{f}_i{x}^2}\text{and}g=\sqrt{\sum _{i=1}^2\sum _{j=1}^2{g}_{ij}{x}^2},\end{array}$$where $f=\sqrt{{nH-S-a\text{BS}/K+B}^2+{a\text{BS}/K+B-rI}^2},$and $$\begin{array}{}\text{(85)}& g=\sqrt{nH+S+2\frac{a\text{BS}}{K+B}+rI}.\end{array}$$

Both $f_i$ and $g_{\text{ij}}$ are continuously differentiable at $S0,I0,B0$ and, hence, satisfy the Lipschitz condition by the mean value theorem for calculus. The norms are bounded because they exist. As a result, the drift and the diffusion matrices are bounded and, hence, satisfy the conditions for uniqueness and existence of the solution criteria [26].

4. Numerical Results and Discussion

In this section, we carry out numerical simulations of the stochastic differential equation model and ordinary differential equation model. The aim of these simulations is to demonstrate numerically the stochastic fluctuations of the infective of the Codeco cholera model for a given initial conditions and a population size.

To achieve this, the Euler-Maruyama scheme was implemented in MATLAB to integrate the model and the individual sample path behavior of the stochastic differential equations (SDEs) models compared to their deterministic solution. In our simulations, one infective is considered and introduced into the population. The Euler-Maruyama scheme is one of the numerical schemes for determining sample paths of stochastic differential equations [27]. We approximate our stochastic differential equation model in the Euler-Maruyama scheme as follows: $St+\text{dt}=St+fSt,tdt+gSt,tWt+\text{dt}-Wt$ and $It+\text{dt}=It+fIt,t\text{dt}+gIt,tWt+dt-Wt$, where $\text{dt}$is the time step.

Throughout the study, we use $\text{time}=100$ days. Thus, the time axis is the number of time steps. The actual total time of $100\text{dt}=1$ and the time step $\text{dt}=0.01.$

The values of our model parameters are based on published epidemiological data shown in Table 4.

Table 4

Parameter values for model simulation.

Parameter	Values	Unit
$$n$$	0.0001	$${\text{day}}^{-1}$$
$$a$$	0.5	$${\text{day}}^{-1}$$
$$K$$	10^6	$$\text{cell}/\text{ml}$$
$$r$$	0.2	$${\text{day}}^{-1}$$
$$\text{mb}-\text{nb}$$	0.33	$${\text{day}}^{-1}$$
$$e$$	10	$$\text{cell}/\text{ml}\text{da}{\text{y}}^{-1}{\text{person}}^{-1}$$

The sample paths of the SI-B stochastic differential, Equations (69) and (74), are shown below in Figures 2–7.

[figure(s) omitted; refer to PDF]

In Figure 2, we observed that at finite time, the susceptible population drops gradually to zero, and the entire individual in the population gets infected.

In Figure 3(b), we tried to demonstrate the impact of the Vibrio cholerae bacterium concentration in the water supply on the number of infected individuals. With Vibrio cholerae bacterium concentration of $B0=100\text{cell}/\text{ml}$, $I0=1$ and a step size of $dt=0.01$, we observed that the $SI-B$ Stochastic Differential Equations (SDE) model sample path of $It$ is continuous but not differentiable. The nowhere differentiability indicates a Wiener process property. The results from the stochastic model also indicate that there will not be much increase in cholera disease (few infections), and within a few days, it dies out of the community since the Vibrio cholerae bacterium concentration in the water supply is small. According to Doldersum [28], cholera is a dose-dependent disease that requires ${10}^4$ cells for an infection. We also, in Figure 3, display the numerical results to compare both the deterministic model and the stochastic one when the basic reproduction number of the deterministic model, $R_0=0.7576$. We see that the infectious population in the stochastic case goes extinct earlier than that of the deterministic one. Furthermore, the stochastic behavior of the curves shows real life as compared to the deterministic one. Thus, stochastic models are better due to their incorporation of white noise or stochastic environmental factors than their deterministic counterparts.

In Figure 4, we observe that when the basic reproduction number is greater than one, the infection persists in both the deterministic case and the stochastic one. We also see that the infectious population in the stochastic case fluctuates randomly which shows a real life behavior, while in the deterministic case, the random fluctuations were not observed. In Figure 4, we can further say that the stochastic solutions are more realistic than deterministic ones.

In Figure 5, we tried to demonstrate the impact of the Vibrio cholerae bacterium concentration in the water supply on the number of infected individuals by simulating the stochastic differential equation model and the corresponding deterministic model. With Vibrio cholerae bacterium concentration of $B0=100cell/ml$, $I0=1$ and a step size of $dt=0.01$, we observed that the $\text{SI}-B$ stochastic differential equations (SDE) model sample path of $It$ fluctuates in the solution of the $\text{SI}-B$ ordinary differential equations (ODE) model. The results from both the stochastic model and the deterministic model indicate that there will not be much increase in cholera disease since the Vibrio cholerae bacterium concentration in the water supply is small. According to [28], cholera is a dose-dependent disease, which requires ${10}^4$ cells for an infection.

In Figure 6, we tried to demonstrate the impact of the Vibrio cholerae bacterium concentration in the aquatic reservoir on the number of infected populations. The numerical results were obtained by varying the value of the Vibrio cholerae bacterium concentration, $Bt$ while keeping other parameters constant. The result from the stochastic model maintained their perturbing property due to the randomness behavior. However, the overall outcome is that the number of infected individuals increases significantly with increasing value of Vibrio cholerae bacterium concentration in the water supply. We also observed that the $\text{SI}-B$ stochastic differential equation model’s sample path of $It$ is continuous but not differentiable. The nowhere differentiability indicates a Wiener process property.

Therefore, we can infer that, when the Vibrio cholerae bacterium concentration in the water supply is increasing, and other parameters are kept constant, the cholera disease transmission expands in the community.

In Figure 7, we tried to demonstrate the impact of the Vibrio cholerae bacterium concentration, $Bt$, in the water supply on the number of infected individuals. The numerical result was obtained by varying the value of $Bt$ while keeping other parameters constant. In both the stochastic differential equation model and the deterministic model (black dashed curve), when the value of $Bt$ increased from $100\text{cell}/\text{ml}$ to $1000\text{cell}/\text{ml}$, there is a significant and regular increase in the number of infected individuals. We also observed that the $\text{SI}-B$ stochastic differential equation (SDE) model’s sample path of $It$ fluctuates in the solution of the $\text{SI}-B$ ordinary differential equation (ODE) model. Therefore, increase in bacteria level in the aquatic environment leads to an increase in infections.

5. Discussions and Conclusions

In this work, Codeco’s work on modeling cholera outbreak and endemic under the influence of the aquatic environment is reviewed and extended to stochastic model using transition probabilities. A stochastic differential equation model is designed from the deterministic model and both investigated for the dynamics of cholera transmission. The stochastic model is a 2-dimensional diffusion process of the susceptible and the infected classes. Our focus is on the interaction of the pathogens from the environment to human and the shedding of bacteria from the infected individuals into the environment.

For the deterministic model, a basic reproductive number $R_0$ was obtained which predicts whether the disease is eradicated or remain in the given population. We also, via mathematical analysis, shown that the model solutions are all nonnegative and bounded in a given region and, thus, exist uniquely. The disease-free equilibrium and the endemic equilibrium are, respectively, derived as $R_0<1$ and $R_0>1$, and the stability of both is examined. It is found that the disease-free equilibrium state exists and is locally asymptotically stable; hence, when $R_0<1$, there will be no cholera outbreak in the community. It was also found that there is an asymptotically stable positive endemic equilibrium locally; so, if $R_0>1$, there will be an outbreak of disease in the community. At zero infections and zero toxigenic bacteria, the disease-free state stability is achieved.

For the transition probabilities stochastic differential equations (SDEs), we determined the existence and uniqueness of the solutions using mean-value theorem of calculus.

The Euler-Maruyama numerical method is used to simulate the sample trajectories of the stochastic differential equation model via numerical simulations. The findings show that the sample paths of the stochastic differential equation model fluctuate in the solution of the corresponding deterministic model and are continuous but not differentiable which is a Wiener process property. Also, it is observed from the graphs that cholera outbreak is independent of the number of infected individual but on ingestion and discharge of the bacteria into the aquatic environment and the infectious population decreases, while the toxigenic Vibrio cholerae bacteria concentration in water remains low any time the infected individual’s contribution to the aquatic reservoir is small. This keeps the reproductive number, $R_0$, less than a unit and greater than a unit any time the concentration of toxigenic Vibrio cholerae in water supply and each infected individual’s contribution to the aquatic reservoir increases in the population. We also find out that the Vibrio cholerae bacteria concentration in the water depends mostly on the rate at which people are exposed to contaminated water supply and on each infected individual’s contribution to the aquatic reservoir.

Therefore, cholera transmission dynamics may also be studied applying stochastic differential equation (SDEs) models which allows for the inclusion of randomness. According to Keeling et al. [29], real world problems such as diseases experience stochasticity in terms of opportunities for transmission.

Acknowledgments

The authors funded the study themselves.

References

[1] A. McElroy, P. K. Townsend, Medical Anthropology in Ecological Perspective, 2009.

[2] A. T. Crooks, B. H. Atesmachew, "An agent-based modeling approach applied to the spread of cholera," Environmental Modelling and Software, vol. 62, pp. 164-177, DOI: 10.1016/j.envsoft.2014.08.027, 2014.

[3] C. T. Codeço, "Endemic and epidemic dynamics of cholera: the role of the aquatic reservoir," BMC Infectious Diseases, vol. 1 no. 1,DOI: 10.1186/1471-2334-1-1, 2001.

[4] S. Fatima, I. Krishnarajah, M. Z. A. M. Jaffar, M. B. Adam, "A mathematical model for the control of cholera in Nigeria," Research Journal of Environmental and Earth Sciences, vol. 6 no. 6, pp. 321-325, DOI: 10.19026/rjees.6.5776, 2014.

[5] A. Akor, Handbook on Management of Tropical Disease in Developing Countries. A Guide to the Diagnosis and Treatment of Common Diseases in the Tropics, 2007.

[6] S. A. Fakai, M. O. Ibrahim, A. M. Siddiqui, "A deterministic mathematical model on cholera dynamics and some control strategies," International Journal of Scientific Engineering and Technology, vol. 8 no. 3, pp. 1115-1118, 2014.

[7] J. Wang, C. Modnak, "Modeling cholera dynamics with controls," Canadian Applied Mathematics Quarterly, vol. 19 no. 3, pp. 255-273, 2011.

[8] V. Capasso, S. Paveri-Fontana, "A mathematical model for the 1973 cholera epidemic in the European Mediterranean region," Revue d'epidemiologie et de sante publique, vol. 27 no. 2, pp. 121-132, 1979.

[9] K. H. Hntsa, B. N. Kahsay, "Analysis of cholera epidemic controlling using mathematical modeling," International Journal of Mathematics and Mathematical Sciences, vol. 2020,DOI: 10.1155/2020/7369204, 2020.

[10] N. K. O. Opoku, C. Afriyie, "The role of control measures and the environment in the transmission dynamics of cholera," Abstract and Applied Analysis, vol. 2020,DOI: 10.1155/2020/2485979, 2020.

[11] Z. Mukandavire, S. Liao, J. Wang, H. Gaff, D. L. Smith, J. G. Morris, "Estimating the reproductive numbers for the 2008-2009 cholera outbreaks in Zimbabwe," Proceedings of the National Academy of Sciences of the United States of America, vol. 108 no. 21, pp. 8767-8772, DOI: 10.1073/pnas.1019712108, 2011.

[12] J. A. Mark, S. M. Faruque, J. J. Mekalanos, B. R. Levin, "Modeling the role of bacteriophage in the control of cholera outbreaks," Proceedings of the National Academy of Sciences of the United States of America, vol. 103 no. 12, pp. 4652-4657, DOI: 10.1073/pnas.0600166103, 2006.

[13] J. M. Ochoche, "A mathematical model for the transmission dynamics of cholera with control strategy. International," Journal of Science and Technology, vol. 2 no. 11, pp. 797-803, 2013.

[14] N. J. E. H. A. Togbenon, E. Moyo, "Modeling and analysis of cholera dynamics with vaccination," vol. 7 no. 1,DOI: 10.12691/ajams-7-1-1, 2019.

[15] D. Posny, J. Wang, "Modelling cholera in periodic environments," Journal of Biological Dynamics, vol. 8 no. 1,DOI: 10.1080/17513758.2014.896482, 2014.

[16] F. Nyabadza, J. M. Aduamah, J. Mushanyu, "Modelling cholera transmission dynamics in the presence of limited resources," BMC Research Notes, vol. 12 no. 1,DOI: 10.1186/s13104-019-4504-9, 2019.

[17] A. P. Lemos-paião, C. J. Silva, D. F. M. Torres, "An epidemic model for cholera with optimal control treatment," Journal of Computational and Applied Mathematics, vol. 318, pp. 168-180, DOI: 10.1016/j.cam.2016.11.002, 2017.

[18] O. Peter, A. Ayoade, A. Abioye, A. Vivtor, C. Akpan, "Sensitivity analysis of the parameters of a cholera model," Journal of Applied Sciences and Environmental Management, vol. 22 no. 4, pp. 477-481, DOI: 10.4314/jasem.v22i4.6, 2018.

[19] H. O. Nyaberi, D. M. Malonza, "Mathematical model of cholera transmission with education campaign and treatment through quarantine," Journal of Advances in Mathematics and Computer Science, vol. 32 no. 3,DOI: 10.9734/JAMCS/2019/v32i330145, 2019.

[20] T. Abdulai, Using Atanackovic and Stankovic Numerical Method to Investigate Fractional Order Cholera Model, 2015.

[21] A. Mwasa, J. M. Tchuenche, "Mathematical analysis of a cholera model with public health interventions," Biosystems, vol. 105 no. 3, pp. 190-200, DOI: 10.1016/j.biosystems.2011.04.001, 2011.

[22] J. Lemaitre, D. Pasetto, J. Perez-saez, C. Sciarra, J. F. Wamala, "Rainfall as a driver of epidemic cholera: comparative model assessments of the effect of intra-seasonal precipitation events," Acta Tropica, vol. 190, pp. 235-243, DOI: 10.1016/j.actatropica.2018.11.013, 2019.

[23] Y. Liang, D. Greenhalgh, X. Mao, "A stochastic differential equation model for the spread of HIV amongst people who inject drugs," Computational and Mathematical Methods in Medicine, vol. 2016,DOI: 10.1155/2016/6757928, 2016.

[24] E. Allen, Modeling with Itô Stochastic Differential Equations, 2007.

[25] D. Greenhalgh, Y. Liang, X. Mao, "SDE SIS epidemic model with demographic stochasticity and varying population size," Applied Mathematics and Computation, vol. 276, pp. 218-238, DOI: 10.1016/j.amc.2015.11.094, 2016.

[26] T. O. Ogunlade, O. M. Ogunmiloro, S. N. Ogunyebi, G. E. Fatoyinbo, J. O. Okoro, O. R. Akindutire, O. Y. Halid, A. O. Olubiyi, "On the effect of vaccination, screening and treatment in controlling typhoid fever spread dynamics: deterministic and stochastic applications," Mathematics and Statistics, vol. 8 no. 6, pp. 621-630, DOI: 10.13189/ms.2020.080601, 2020.

[27] K. Burrage, P. M. Burrage, T. Tian, "Numerical methods for strong solutions of stochastic differential equations: an overview," Proceedings of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences, vol. 460, pp. 373-402, DOI: 10.1098/rspa.2003.1247, 2004.

[28] T. Doldersum, The Role of Water in Cholera Diffusion: Improvements of a Cholera Diffusion Model for Kumasi, Ghana, [M.S. Thesis], 2013.

[29] M. J. Keeling, P. Rohani, B. Pourbohloul, "Modeling infectious diseases in humans and animals:modeling infectious diseases in humans and animals," Clinical Infectious Diseases, vol. 47 no. 6, pp. 864-865, DOI: 10.1086/591197, 2008.