Global Dynamics Analysis of a Cholera Transmission Model with General Incidence and Multiple Modes of Infection

Daijun Li, Lianwen Wang^*

School of Mathematics and Statistics, Hubei Minzu University, Enshi, China.

DOI: 10.4236/jamp.2023.1111236                

Abstract

This paper develops an SIBR cholera transmission model with general incidence rate. Necessary and sufficient conditions for local and global asymptotic stability of the equilibria are established by Routh Hurwitz criterium, Lyapunov function, and the second additive composite matrix theorem. What is more, exploiting the DED is cover simulation tool, the parameter values of the model are estimated with the 1998-2021 cholera case data in China. Finally, we perform sensitivity analysis for the basic reproduction number to seek for effective interventions for cholera control.

Keywords

Cholera Transmission Model, Nonlinear Incidence Rate, Second Additive Composite Matrix, Lyapunov Function, Global Stability

1. Introduction

Cholera, a grave waterborne ailment caused by Vibrio cholerae, exhibits a remarkable ability to persist in certain aquatic environments for durations spanning three months to two years. Its clinical manifestations are characterized by intense diarrhea and vomiting, with severe cases leading to potentially fatal dehydration due to significant loss of bodily fluids and electrolytes. Categorized as a class A infectious disease in China, cholera is distinguished by its sudden onset and rapid transmission dynamics. The primary mode of transmission involves the interaction between humans and their environment, particularly through the ingestion of food or water contaminated by the Vibrio cholerae bacteria [1] [2]. Additionally, a secondary transmission route occurs through direct human-to-human interactions, specifically close contact with infected individuals [3] [4]. Globally, cholera imposes a substantial burden, resulting in 3 - 5 million new infections annually and causing 28,800 - 130,000 deaths each year. Despite a typically low mortality rate of under five percent, this rate can skyrocket to fifty percent in regions where access to treatment is limited. Children, especially in Africa and Southeast Asia, bear the brunt of this disease’s impact. Cholera continues to pose a significant public health threat in developing nations, attracting sustained attention and research interest from scholars [5]-[10].

Significantly, Wang et al. [8] introduced a distinctive SIBR cholera transmission model, encompassing both direct human-to-human and indirect environment-to-human transmission pathways:

$\{\begin{array}{l}{S}^{\prime }=\mu N-{\beta }_1\left(I\right)SI-\frac{{\beta }_2\left(I\right)SB}{B+K}-\mu S+\sigma R,\\ I^{\prime }={\beta }_1\left(I\right)SI+\frac{{\beta }_2\left(I\right)SB}{B+K}-\left(\gamma +\mu \right)I,\\ B^{\prime }={\beta }_3\left(I\right)SI-\delta B,\\ R^{\prime }=\gamma I-\left(\mu +\sigma \right)R.\end{array}$ (1.1)

In this framework, the entire human population, denoted by the constant size N, is categorized into distinct compartments: the susceptible (S), the infectious (I) and the recovered (R) [11] [12] [13]. Additionally, there is an auxiliary compartment B to quantify the concentration of Vibrio in contaminated water. The transmission rates and bacterial shedding rates are contingent upon the number of infectious individuals. Within this context, the parameters are defined as follows: $\mu$ represents the natural mortality rate, $\delta$ signifies the net mortality rate of bacteria, $\gamma$ denotes the recovery rate, $\sigma$ represents the rate of host immune loss, $g\left(I\right)$ represents the infection rate of patients, influenced by the virus’s spread function $f\left(I,B\right)$ . Furthermore, $\eta \left(I\right)$ signifies the host shedding rate. It stands to be reasonable that all the parameters keep nonnegative.

Furthermore, this model accounts for the influence of human behavior arising from health education, improved hygiene, and sanitation practices [14] [15]. Importantly, given that immunity gradually wanes, individuals who have recovered from cholera remain susceptible to reinfection after a certain period. Studies findings indicate that diminished exposure due to altered human behavior can lead to a reduction in the scale of both epidemics and endemic diseases. Building upon the concepts used in [8], we delve into a more comprehensive model. Specifically, we explore the properties of a generalized SIBR cholera transmission model with a broadened incidence function. The proposed model is expressed as follows:

$\{\begin{array}{l}{S}^{\prime }=\mu N-Sg\left(I\right)-Sf\left(I,B\right)-\mu S+\sigma R,\\ I^{\prime }=Sg\left(I\right)+Sf\left(I,B\right)-\left(\gamma +\mu \right)I,\\ B^{\prime }=\eta \left(I\right)-\delta B,\\ R^{\prime }=\gamma I-\left(\mu +\sigma \right)R.\end{array}$ (1.2)

Assume that the total population is constant $N=S+I+R$ . The model is based on the standard SIR (susceptibility-infection-recovery) compartment structure and has an additional compartment B indicating the concentration of Vibrio cholera in contaminated water. Based on the above assumptions, model (1.2) can be established [16] [17] [18].

The structure of this article unfolds as follows: Section 1 is the proof of the positivity and boundedness of solutions for the model (1.2). Section 2 is dedicated to establishing the existence and stability of equilibria, encompassing both the disease-free equilibrium and the positive equilibrium points. Some numerical simulations and sensitivity analyses are performed in Section 3 for the parameter $R_0$ , corroborating the theoretical analysis mentioned earlier. Finally, the article concludes in Section 4 with a concise discussion summarizing our findings.

2. Positive and Boundedness of Solutions

2.1. Positivity of Solutions

Theorem 2.1. Under nonnegative initial conditions, for $t>0$ , the solution $\left(S\left(t\right),I\left(t\right),B\left(t\right),R\left(t\right)\right)$ of model (1.2) is nonnegative.

Proof. Let $t=Sup\left\{t>0|S>0,I>0,B>0,R>0\right\}$ . Now, from the first equation of model (1.2), we obtain

$\frac{dS\left(t\right)}{dt}\ge \mu N-\left(g\left(I\right)+f\left(I,B\right)+\mu \right)S$ .

From the above equation, we can reduce that

$\begin{array}{l}\begin{array}{l}\frac{d}{dt}\left[S\left(t\right)\exp \left\{{\displaystyle {\int }_0^{t}g\left(I\right)d\varsigma }+{\displaystyle {\int }_0^{t}f\left(I,B\right)d\varsigma }+\mu t\right\}\right]\hfill \end{array}\\ \ge \mu N\exp \left\{{\displaystyle {\int }_0^{t}g\left(I\right)d\varsigma }+{\displaystyle {\int }_0^{t}f\left(I,B\right)d\varsigma }+\mu t\right\}.\end{array}$

Further,

$\begin{array}{l}S\left(t_1\right)\exp \left\{{\displaystyle {\int }_0^{t_1}g\left(I\right)d\varsigma }+{\displaystyle {\int }_0^{t_1}f\left(I,B\right)d\varsigma +\mu t_1}\right\}-S\left(0\right)\\ \ge {\displaystyle {\int }_0^{t_1}\mu N\exp \left\{{\displaystyle {\int }_0^{\varphi }g\left(I\right)d\varphi }+{\displaystyle {\int }_0^{\varphi }f\left(I,B\right)d\varphi +\mu \varphi }\right\}d\varphi }.\end{array}$

Therefore,

$\begin{array}{l}S\left(t_1\right)\ge S\left(0\right)\exp \left\{-{\displaystyle {\int }_0^{t_1}g\left(I\right)d\varsigma }-{\displaystyle {\int }_0^{t_1}f\left(I,B\right)d\varsigma -\mu t_1}\right\}\hfill \\ +\exp \left\{-{\displaystyle {\int }_0^{t_1}g\left(I\right)d\varsigma }-{\displaystyle {\int }_0^{t_1}f\left(I,B\right)d\varsigma +\mu t_1}\right\}\hfill \\ \times {\displaystyle {\int }_0^{t_1}\mu N\exp \left\{{\displaystyle {\int }_0^{\varphi }g\left(I\right)d\varphi }+{\displaystyle {\int }_0^{\varphi }f\left(I,B\right)d\varphi +\mu \varphi }\right\}d\varphi }\hfill \\ >0.\hfill \end{array}$

Similarly, we can obtain the bounds for the other components of the solution.

2.2. Boundedness of Solutions

Theorem 2.2. All solutions of model (1.2) are bounded.

Proof. The model (1.2) consists of two populations, namely human and pathogen. Therefore, we will break the model (1.2) into two parts, of which one involves the human population $\left(S,I,R\right)$ and the other the pathogen population B. According to model (1.2) we obtain

$\frac{d\left(S+I+R\right)}{dt}d\left(S+I+R\right)/dt=\mu \left(N-S-I-R\right).$

Further, from the first equation of the model (1.2), we have

$S^{\prime }\le N\left(\mu +\sigma \right)-S\left(\mu +\sigma \right).$

Hence, we conclude that $S<N$ . Now, from the last equation of the model, we deduce that $R^{\prime }\le rN-\left(r+\mu +\sigma \right)R$ . We can obtain $R\le rN/\left(r+\mu +\sigma \right)$ . According to the third equation of the model (1.2) and assumptions in reference [8], we can arrive at

$B^{\prime }={\alpha }_3-b_3{m}_3\left(I\right)-\delta B\le {\alpha }_3-\delta B.$

Therefore, $B\le {\alpha }_3/\delta$ . From the above discussion, it is clear from the above discussion that all solutions are bounded. Next, we obtain the feasible region for the human population as

${\Omega }_H=\left\{\left(S,I,R\right)|S+I+R=N,0\le S\le N,0\le R\le rN/\left(r+\mu +\sigma \right)\right\}.$

And the feasible region for pathogen population is

${\Omega }_B=\left\{B|0\le B\le {\alpha }_3/\delta \right\}.$

Define $\Omega ={\Omega }_H\times {\Omega }_B$ . Now, $\Omega$ is a positively invariant region for the model (1.2). Moreover, the model (1.2) is mathematically and epidemiologically well-posed with the method utilized in [19].

3. Main Content

3.1. The Existence of the Equilibria

The existence of equilibria is discussed below, we define

$\begin{array}{l}g\left(0\right)=0,f\left(0,B\right)\ge 0,f\left(I,0\right)=0,\eta \left(0\right)=0,{f^{\prime }}_I\left(0,0\right)=0,{f^{\prime }}_I\left(I,B\right)<f\left(I,B\right)/I,\\ {f^{\prime }}_I\left(I,B\right)<0,0\le {f^{\prime }}_B\left(I,B\right)\le f\left(I,B\right)/B,0\le {\eta }^{\prime }\left(I\right)\le \eta \left(I\right)/I,0\le g^{\prime }\left(I\right)\le g\left(I\right)/I.\end{array}$

Theorem 3.1 When $R_0>1$ , then model (1.2) has two equilibria $E^0$ and $E^{*}$ , when $R_0\le 1$ , model (1.2) has a unique equilibrium $E^0$ .

Proof. Due to $N=S+I+R$ , $R=\gamma I/\left(\mu +\sigma \right)S=N-I-R=N-\left(\mu +\sigma +\gamma \right)I/\left(\mu +\sigma \right)=\varphi \left(I\right)$ , and

$S=\frac{\left(\mu +\gamma \right)I}{g\left(I\right)+f(I,\eta \left(I\right)/\delta )}=\frac{\left(\mu +\gamma \right)}{h\left(I\right)}\triangleq \psi \left(I\right),h\left(I\right)=\left(g\left(I\right)+f(I,\eta \left(I\right)/\delta )\right)/I.$

Because of $S=\varphi \left(I\right)>0$ , then $0<I<\left(\mu +\sigma \right)N/\left(\mu +\sigma +\gamma \right)=I_{\max }$ . Thus, considering $S=\varphi \left(I\right)$ and $S=\psi \left(I\right)$ at the intersection where $[0,I_{\max })$ Thus, we can obtain

$\begin{array}{c}{h}^{\prime }\left(I\right)=\frac{1}{I^2}\left[g^{\prime }\left(I\right)I-g\left(I\right)+{f^{\prime }}_I\left(I,B\right)I+\left({f^{\prime }}_B\left(I,B\right){\eta }^{\prime }\left(I\right)\right)I/\delta -f\left(I,B\right)\right]\\ \le \frac{1}{I^2}\left[g^{\prime }\left(I\right)I-g\left(I\right)+{f^{\prime }}_B\left(I,B\right)\eta \left(I\right)/\delta -f\left(I,B\right)\right]\\ =\frac{1}{I^2}[g^{\prime }\left(I\right)I-g\left(I\right)+{f^{\prime }}_B\left(I,B\right)B-f\left(I,B\right)]\\ \le 0.\end{array}$

This indicates that ${\psi }^{\prime }\left(I\right)\ge 0$ , On the other hand ${\varphi }^{\prime }\left(I\right)\le 0,\varphi \left(0\right)=N,\varphi \left(I_{\max }\right)=0$ and $\psi \left(I_{\max }\right)>0$ , then

$\psi \left(0^{+}\right)=\underset{I\to 0^{+}}{\lim }\frac{\mu +\gamma }{h\left(I\right)}=\frac{\mu +\gamma }{g^{\prime }\left(0\right)+{f^{\prime }}_I\left(0,0\right)+{f^{\prime }}_B\left(0,0\right){\eta }^{\prime }\left(0\right)/\delta }=\frac{N}{T_1}.$

Therefore, when $T_1>1$ , $\varphi \left(0\right)>\psi \left(0^{+}\right)$ is the only one node, denoted as $I^{\ast }$ . When $T_1\le 1$ , $\varphi \left(0\right)\le \psi \left(0^{+}\right)$ , there is no node. With the utilization of next generation matrix method mentioned in [20] [21], the matrix sum F and V can be written as

$F=\left[\begin{array}{cc}N{g}^{\prime }\left(0\right)& N{f^{\prime }}_B\left(0,0\right)\\ {\eta }^{\prime }\left(0\right)& 0\end{array}\right],V=\left[\begin{array}{cc}\mu +\gamma & 0\\ 0& \delta \end{array}\right].$

Therefore, the basic reproduction number $R_0$ of the model can be obtained as follows

$R_0=\rho \left(M\right)=\frac{1}{2}\left[\frac{N{g}^{\prime }\left(0\right)}{\mu +\gamma }+\sqrt{{\left(\frac{N{g}^{\prime }\left(0\right)}{\mu +\gamma }\right)}^2+4\frac{N{f^{\prime }}_B\left(0,0\right)}{\delta }\frac{{\eta }^{\prime }\left(0\right)}{\mu +\gamma }}\right].$

At the same time, if disease control targets at a particular host type, a useful threshold is called the reproduction number T. The reproduction number defines the expected number of secondary infections due to a typical primary case in a fully susceptible population [22] [23]. It is an extension of the basic reproduction number $R_0$ . According to literature [23], it is concluded that $R_0<1\left(\ge 1\right)\iff T_1<1\left(\ge 1\right)$ . In the following analysis, we will use both and realize that the two are equivalent in characterizing the disease threshold dynamics.

3.2. Stability of Disease-Free Equilibrium E⁰

Theorem 3.2.1. When $R_0<1$ , the disease-free equilibrium point of model (1.2) is locally asymptotically stable. When $R_0>1$ , the disease-free equilibrium point of model (1.2) is unstable.

Proof. According to the model (1.2), we take $R=N-S-I$ . Thus, it can be seen that

$\{\begin{array}{l}{S}^{\prime }=\left(\mu +\sigma \right)N-Sg\left(I\right)-Sf\left(I,B\right)-S\left(\mu +\sigma \right)-\sigma I,\hfill \\ I^{\prime }=Sg\left(I\right)+Sf\left(I,B\right)-\left(\gamma +\mu \right)I,\hfill \\ B^{\prime }=\eta \left(I\right)-\delta B.\hfill \end{array}$ (3.1)

Then, its Jacobi matrix of (3.1) at $E^0$ is

$J=\left[\begin{array}{ccc}\text{}\text{}\text{}-g\left(I\right)-f\left(I,B\right)-\left(\mu +\sigma \right)& \text{}\text{}-S{g}^{\prime }\left(I\right)-S{f^{\prime }}_I\left(I,B\right)-\sigma & \text{}-S{f^{\prime }}_B\left(I,B\right)\text{}\text{}\text{}\\ g\left(I\right)+f\left(I,B\right)& \text{}\text{}S{g}^{\prime }\left(I\right)+S{f^{\prime }}_I\left(I,B\right)-\left(\gamma +\mu \right)& S{f^{\prime }}_B\left(I,B\right)\text{}\text{}\text{}\\ 0& {\eta }^{\prime }\left(I\right)& -\delta \end{array}\right].$

The characteristic equation is $\left(\lambda +\mu +\sigma \right)\left({\lambda }^2+a_1\lambda +a_2\right)=0$ .

As for $a_1=\gamma +\mu +\delta -N{g}^{\prime }\left(0\right)$ , $a_2=\delta \left(\gamma +\mu -N{g}^{\prime }\left(0\right)\right)-N{f^{\prime }}_B\left(0,0\right){\eta }^{\prime }\left(0\right)$ . It is easy to know that one of the eigenvalues is $\lambda =-\left(\mu +\sigma \right)$ , rest up to $\left({\lambda }^2+a_1\lambda +a_2\right)=0$ , take advantage of $T_1=N{g}^{\prime }\left(0\right)/\left(\mu +\gamma \right)+N{f^{\prime }}_B\left(0,0\right){\eta }^{\prime }\left(0\right)/\left(\mu +\gamma \right)\delta$ , $T_1<1$ , $a_1>0$ , and $a_2>0$ .

By the Hurwitz criterion, when $T_1<1$ , all the eigenvalues are negative, the disease-free equilibrium is locally asymptotically stable. When $T_1>1$ , the characteristic equation consists of one positive root and two negative roots, so the disease-free equilibrium is unstable, completing the proof.

Theorem 3.2.2. When $R_0\le 1$ , model (1.2) has a globally asymptotically stable disease-free equilibrium point.

Proof. Now we use the next generation matrix method to prove this theorem. Establish

$F_1=\left[\begin{array}{cc}N{g}^{\prime }\left(0\right)& N{f^{\prime }}_B\left(0,0\right)\\ 0& 0\end{array}\right],V_1=\left[\begin{array}{cc}\mu +\gamma & 0\\ -{\eta }^{\prime }\left(0\right)& \delta \end{array}\right].$

Because of $I\ge 0$ , $g\left(I\right)\ge 0$ . If and only if $I=0$ , $g^{″}\left(I\right)\le 0$ . So $g^{\prime }\left(I\right)I\le g\left(I\right)\le g^{\prime }\left(0\right)I$ . When $y={\left(I,B\right)}^T$ , model (3.1) satisfies

$\frac{dy}{dt}\le dy\left(F_1-V_1\right)y.$

Le $w=\left(N{g}^{\prime }\left(0\right),N{f^{\prime }}_B(0,0)\right)$ , due to $T_1=\rho \left(F_1{V}_1{}^{-1}\right)=\rho \left(V_1{}^{-1}{F}_1\right)$ we can prove $w{V}_1{}^{-1}{F}_1=T_{1}w$ . From [24] we define a Lyapunov function as follows $L=w{V}_1{}^{-1}y$ , then the derivative of L can be written as

$L^{\prime }=w{V}_1{}^{-1}\frac{dy}{dt}\le w{V}_1{}^{-1}\left(F_1-V_1\right)y=\left(T_1-1\right)wy.$

When $T_1\le 1$ , the disease-free equilibrium is globally asymptotically stable.

3.3. Stability of Positive Equilibrium E^*

Consider the differential equation $\stackrel{\dot{}}{x}=f(x),x\in D$ . Let $x\left(t,x_0\right)$ be the solution to this equation with the initial value $x\left(0,x_0\right)=x_0$ satisfying two hypotheses listed as follows [25].

(H₁) There exists a compact attractive subset $K\subset D$ ,

(H₂) Model (3.1) has a unique equilibrium $x\in D$ .

Lemma 3.3.1. [26] If $tr\left(M\left(E^{*}\right)\right)<0$ , $det\left(M\left(E^{*}\right)\right)<0$ , and $det\left(M^{\left[2\right]}\left(E^{*}\right)\right)<0$ , then all eigenvalues are negative real numbers.

Theorem 3.3.1 When $R_0>1$ , the positive equilibrium point $E^{*}$ is locally asymptotically stable.

Proof. The Jacobi matrix of model (3.1) at $E^{*}$ is

$M\left(E^{*}\right)=\left[\begin{array}{ccc}-g\left(I^{*}\right)-f\left(I^{*},B^{*}\right)-\left(\mu +\sigma \right)& -S{g}^{\prime }\left(I^{*}\right)-S{f^{\prime }}_I\left(I^{*},B^{*}\right)-\sigma & -S{f^{\prime }}_B\left(I^{*},B^{*}\right)\text{}\text{}\\ g\left(I^{*}\right)+f\left(I^{*},B^{*}\right)& S{g}^{\prime }\left(I^{*}\right)+S{f^{\prime }}_I\left(I^{*},B^{*}\right)-\left(\gamma +\mu \right)\text{}& S{f^{\prime }}_B\left(I^{*},B^{*}\right)\\ 0& {\eta }^{\prime }\left(I^{*}\right)& -\delta \end{array}\right].$

Therefore $tr\left(M\left(E^{*}\right)\right)=S{g}^{\prime }\left(I^{*}\right)+S{f^{\prime }}_I\left(I^{*},B^{*}\right)-g\left(I^{*}\right)-f\left(I^{*},B^{*}\right)-\left(\gamma +\mu +\delta \right)-\left(\gamma +\mu \right)$ , take advantage of $\left(\mu +\gamma \right)=\left(Sg\left(I\right)+Sf\left(I,B\right)\right)/I$ , it is easy to derive

$\begin{array}{c}tr\left(M\left(E^{*}\right)\right)=S{g}^{\prime }\left(I^{*}\right)-Sg\left(I^{*}\right)/I^{*}+S{f^{\prime }}_I\left(I^{*},B^{*}\right)-Sf\left(I^{*},B^{*}\right)/I^{*}\\ -g\left(I^{*}\right)-f\left(I^{*},B^{*}\right)-\left(\mu +\sigma +\delta \right)\\ <0.\end{array}$

The determinant of $M\left(E^{*}\right)$ is

$\begin{array}{c}det\left(M\left(E^{*}\right)\right)=\delta S/I^{*}\left(\mu +\sigma \right)\left(B{f^{\prime }}_B\left(I^{*},B^{*}\right)-f\left(I^{*},B^{*}\right)\right)\\ +\delta \left(\mu +\sigma \right)\left(S{g}^{\prime }\left(I^{*}\right)-Sg\left(I^{*}\right)/I^{*}+S{f^{\prime }}_I\left(I^{*},B^{*}\right)\right)\\ -\delta (g\left(I^{*}\right)+f\left(I^{*},B^{*}\right))\left(S{g}^{\prime }\left(I^{*}\right)+S{f^{\prime }}_I\left(I^{*},B^{*}\right)+\sigma \right)\\ <0.\end{array}$

We can write

$M^{\left[2\right]}\left(E^{*}\right)=\left[\begin{array}{ccc}{M}_{11}\left(E^{*}\right)& S{f^{\prime }}_B\left(I^{*},B^{*}\right)& S{f^{\prime }}_B\left(I^{*},B^{*}\right)\\ {\eta }^{\prime }\left(I^{*}\right)& -g\left(I^{*}\right)-f\left(I^{*},B^{*}\right)-\left(\mu +\sigma +\delta \right)& -S{g}^{\prime }\left(I^{*}\right)-S{f^{\prime }}_I\left(I^{*},B^{*}\right)-\sigma \\ 0& \text{}\text{}\text{}g\left(I^{*}\right)+f\left(I^{*},B^{*}\right)& S{g}^{\prime }\left(I^{*}\right)+S{f^{\prime }}_I\left(I^{*},B^{*}\right)-\left(\gamma +\mu +\delta \right)\text{}\text{}\text{}\end{array}\right].$

Among $M_{11}\left(E^{*}\right)=S{g}^{\prime }\left(I^{*}\right)+S{f^{\prime }}_I\left(I^{*},B^{*}\right)-g\left(I^{*}\right)-f(I^{*},B^{*})-\left(\sigma +\gamma +2\mu \right)$ . Calculate its determinant as

$\begin{array}{l}\det \left(M^{\left[2\right]}\left(E^{*}\right)\right)<\left(S{g}^{\prime }\left(I^{*}\right)+S{f^{\prime }}_I\left(I^{*},B^{*}\right)-\left(\gamma +\mu +\delta \right)-g\left(I^{*}\right)-f\left(I^{*},B^{*}\right)\right)\\ \left(\delta Sf\left(I^{*},B^{*}\right)/I^{*}-{\eta }^{\prime }\left(I^{*}\right)S{f^{\prime }}_B\left(I^{*},B^{*}\right)\right)\\ <0.\end{array}$

It’s made use of ${\eta }^{\prime }\left(I\right)S{f^{\prime }}_B\left(I,B\right)\le \eta \left(I\right)S{f^{\prime }}_B\left(I,B\right)/I\le \delta BS{f^{\prime }}_B\left(I,B\right)/I\le \delta Sf\left(I,B\right)/I$ . In conclusion, when $R_0>1$ , $E^{*}$ is locally asymptotic stability, the theorem is proven.

Lemma 3.3.2. [19] If the region D is simply connected and conditions (H₁) and (H₂) hold. When $q<0$ , the only internal equilibrium solution $E^{*}$ of the model $\stackrel{\dot{}}{x}=f(x)$ is globally asymptotically stable.

Since (H₁) is equivalent to the consistent persistence of model (3.1) [26] and bounded in the feasible domain $\stackrel{°}{\Gamma }$ . Then the consistent persistence of model (3.1) is equivalent to $E_0$ being unstable [27]. According to Theorem 2.3.1, when $R_{\text{0}}>\text{1}$ , $E_0$ is unstable. So, the following lemma holds.

Lemma 3.3.3. When $R_0>1$ , model (3.1) is consistent and persistent.

Theorem 3.3.2. When $R_0>0$ , the positive equilibrium point $E^{*}$ is globally asymptotically stable.

Proof. First, based on the $M(E^{*})$ obtained above, the Lyapunov function is established as follows.

$V_2\left(x,u\right)=\left\{\max \left|X\right|,I\left(\left|Y\right|+\left|Z\right|\right)/B\right\}$ .

According to model (3.1), it can be concluded that

$\{\begin{array}{l}{X}^{\prime }=M_{11}X+S{f^{\prime }}_B\left(I,B\right)\left(Y+Z\right),\hfill \\ Y^{\prime }={\eta }^{\prime }\left(I\right)X-\left(g\left(I\right)+f\left(I,B\right)+\left(\mu +\sigma +\delta \right)\right)Y\text{}-\text{}\left(S{g}^{\prime }\left(I\right)+S{f^{\prime }}_I\left(I,B\right)+\sigma \right)Z,\hfill \\ Z^{\prime }=\left(g\left(I\right)+f\left(I,B\right)\right)Y-\left(-S{g}^{\prime }\left(I\right)-S{f^{\prime }}_I\left(I,B\right)+\left(\gamma +\mu +\delta \right)\right)Z.\hfill \end{array}$ (3.2)

Then

$\begin{array}{l}{D}_{+}\left|X\right|\le M_{11}\left|X\right|+\frac{B}{I}S{f^{\prime }}_B\left(I,B\right)\left(\frac{I}{B}\left(\left|Y\right|+\left|Z\right|\right)\right),\hfill \\ D_{+}\left|Y\right|\le {\eta }^{\prime }\left(I\right)\left|X\right|-\left(g\left(I\right)+f\left(I,B\right)+\left(\mu +\sigma +\delta \right)\right)\left|Y\right|-\left(S{g}^{\prime }\left(I\right)+S{f^{\prime }}_I\left(I,B\right)+\sigma \right)\left|Z\right|,\hfill \\ D_{+}\left|Z\right|\le \left(g\left(I\right)+f\left(I,B\right)\right)\left|Y\right|-\left(-S{g}^{\prime }\left(I\right)-S{f^{\prime }}_I\left(I,B\right)+\left(\gamma +\mu +\delta \right)\right)\left|Z\right|.\hfill \end{array}$

The derivative of $V_2$ along the positive solution of model (3.2) can be simplified as the following differential inequality

$\begin{array}{l}{D}_{+}\left[\frac{I}{B}\left(\left|Y\right|+\left|Z\right|\right)\right]=\frac{I}{B}\left(\frac{I^{\prime }}{I}-\frac{B^{\prime }}{B}\right)\left[\left|Y\right|+\left|Z\right|\right]+\frac{I}{B}{D}_{+}\left[\left|Y\right|+\left|Z\right|\right]\\ \le \frac{I}{B}\left\{{\eta }^{\prime }\left(I\right)\left|X\right|-\left(\mu +\sigma +\delta \right)\left|Y\right|-\left(\gamma +\mu +\sigma +\delta \right)\left|Z\right|\right\}\\ +\frac{I}{B}\left(\frac{I^{\prime }}{I}-\frac{B^{\prime }}{B}\right)\left[\left|Y\right|+\left|Z\right|\right]\\ \le \frac{\eta \left(I\right)}{B}\left|X\right|+\frac{I}{B}\left(\frac{I^{\prime }}{I}-\frac{B^{\prime }}{B}-\gamma -\mu -\sigma -\delta \right)\left[\left|Y\right|+\left|Z\right|\right].\end{array}$

Obtained from model (3.1) that $\frac{I^{\prime }}{I}=\frac{Sg\left(I\right)}{I}+\frac{Sf\left(I,B\right)}{I}-\left(\gamma +\mu \right)$ and $\frac{B^{\prime }}{B}=\frac{\eta \left(I\right)}{B}-\delta$ , then

$\begin{array}{l}{g}_1\left(t\right)=S{g}^{\prime }\left(I\right)+S{f^{\prime }}_I\left(I,B\right)-g\left(I\right)-f\left(I,B\right)-\left(\sigma +\gamma +2\mu \right)+\frac{B}{I}S{f^{\prime }}_B\left(I,B\right)\\ =\frac{I^{\prime }}{I}+\left(\mu +\sigma \right)-g\left(I\right)-f\left(I,B\right)+S{g}^{\prime }\left(I\right)-\frac{Sg\left(I\right)}{I}+\frac{B}{I}S{f^{\prime }}_B\left(I,B\right)\\ -\frac{Sf\left(I,B\right)}{I}+S{f^{\prime }}_I\left(I,B\right)\\ \text{ <}\frac{I^{\prime }}{I}-\left(\mu +\sigma \right)\end{array}$ (3.3)

$\begin{array}{l}{g}_2\left(t\right)=\frac{\eta \left(I\right)}{B}+\left(\frac{I^{\prime }}{I}-\frac{B^{\prime }}{B}-\gamma -\mu -\sigma -\delta \right)\\ =\frac{B^{\prime }}{B}+\delta +\frac{I^{\prime }}{I}-\frac{B^{\prime }}{B}-\mu -\sigma -\delta \\ =\frac{I^{\prime }}{I}-\left(\mu +\sigma \right)\end{array}$ (3.4)

By using inequalities (3.3) and (3.4) we have

$D_{+}V(t)<\sup \left\{g_1\left(t\right),g_2\left(t\right)\right\}\le \frac{I^{\prime }}{I}-\left(\mu +\sigma \right)={\lambda }^{\prime }\left(t\right)-\left(\mu +\sigma \right).$

For satisfying that $\left(S\left(0\right),I\left(0\right),B\left(0\right)\right)\in K$ (K is the inner compact attractive set $\stackrel{°}{\Gamma }$ ) is the solution of the model (3.1) $X\left(t\right)=\left(S\left(t\right),I\left(t\right),B\left(t\right)\right)$ , there must be

$\frac{1}{t}{\displaystyle {\int }_0^t{D}_{+}V(t)ds}\le \frac{1}{t}{\displaystyle {\int }_0^t\frac{I^{\prime }}{I}-\left(\mu +\sigma \right)ds}=\frac{1}{t}\ln \frac{I^{\prime }\left(t\right)}{I\left(0\right)}-\left(\mu +\sigma \right).$

Thereby $q\le -\left(\mu +\sigma \right)<-\left(\mu +\sigma \right)/2<0$ . Therefore, when $R_0>1$ , $E^{*}$ is globally asymptotically stable. The proof is now complete.

4. Numerical Simulation

In this section, the DED is cover simulation tool [28] [29] was used to numerically simulate 24 case data from 1998 to 2021, sensitivity analysis for $R_0$ is conducted to reveal the influence degree on model outcomes. According to the assumption in reference [8], we obtain the basic regeneration number

$R_0=\frac{1}{2}\left[\frac{{\alpha }_{1}N}{\mu +\gamma }+\sqrt{{\left(\frac{{\alpha }_{1}N}{\mu +\gamma }\right)}^2+4\frac{{\alpha }_2{\alpha }_{3}N}{\delta \left(\mu +\gamma \right)K}}\right].$

where ${\alpha }_1$ , ${\alpha }_2$ , and ${\alpha }_3$ are the direct propagation rate, indirect propagation rate, and dropout rate, respectively. The biological meaning and the standard deviation of each parameters in model (1.2) are listed in Table 1.

Figure 1 illustrates the trend of cholera cases from 1998 to 2021 in China, indicating a decline in the transmission of this infectious disease. According to the simulation results, the transmission of cholera will eventually be gradually reduced and ultimately controlled.

In Figure 2, sensitivity analysis was applied. The basic reproduction number $R_0$ of model (1.2) may determine the transmissibility, severity, and outcome of the pandemic. In order to seek for effective disease control measures, we therefore shall be concerned with the effects of input parameters ${\alpha }_1$ , ${\alpha }_2$ , ${\alpha }_3$ , $\gamma$ , on $R_0$ . The results show that ${\alpha }_1$ , ${\alpha }_2$ , ${\alpha }_3$ and $\gamma$ are highly correlated with cholera infected persons. In order of relevance ${\alpha }_1$ , $\gamma$ , ${\alpha }_2$ , ${\alpha }_3$ . According to

Table 1. Parameters value.

Figure 1. Comparison of the reported cholera case data in China and the simulated solution $I(t)$ of model (1.2).

Figure 2. The PRCC values affecting the key parameters of $R_0$ are obtained.

numerical simulation, ${\alpha }_1$ , ${\alpha }_2$ , ${\alpha }_3$ are positively proportional to I, $\gamma$ is negatively proportional to I. Undoubtedly, reducing the coefficient of disease transmission ${\alpha }_1$ and ${\alpha }_2$ , such as epidemic prevention propaganda, isolation, sterilization, and wearing masks can effectively control the spread of cholera. On the other hand, shortening the disease course of disease $\gamma$ can reduce the number of infected individuals. Therefore, it is possible for policy-makers to use multiple control measures jointly during the influenza pandemic.

5. Conclusion

Our study delves deep into the intricacies of the SIBR cholera transmission model, incorporating multiple modes of infection and a generalized incidence function. Initially, we derive the expression for the basic reproduction number. Subsequently, employing the Routh-Hurwitz condition and constructing the Lyapunov function, we establish a pivotal insight: when $R_0\le 1$ , the disease-free equilibrium point is globally asymptotically stable. This implies that, absent any interventions, infectious diseases will eventually fade away. However, if $R_0>1$ , we ascertain that the endemic equilibrium point becomes globally asymptotically stable, indicating the disease’s perpetual presence. Notably, we address and surmount the constraints posed by the existing literature [8], demonstrating that infectious diseases persist under specific conditions denoted by $S{\left({\beta }_1\left(I\right)I\right)}^{\prime }\le \left(\gamma -\sigma \right)/2$ . Furthermore, through meticulous numerical simulations, we elucidate the implications of the generalized incidence transmission model on the proliferation and containment of infectious diseases. These simulations form a theoretical foundation, enabling the evaluation of the efficacy of disease control measures. In essence, our analysis furnishes invaluable insights into the dynamics of the cholera transmission model, underscoring the criticality of implementing robust control strategies to thwart outbreaks and curtail the disease’s spread.

Acknowledgements

This work is supported partially by National Natural Science Foundation of China (Nos.12001178).

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Tien, J.H. and Earn, D.J.D. (2010) Multiple Transmission Pathways and Disease Dynamics in a Waterborne Pathogen Model. Bulletin of Mathematical Biology, 72, 1506-1533.
[2]	Wang, Y. and Cao, J. (2015) Global Stability of General Cholera Models with Nonlinear Incidence and Removal Rates. Journal of the Franklin Institute, 352, 2464-2485.
[3]	Eisenberg, M.C., Shuai, Z., Tien, J.H., et al. (2013) A Cholera Model in a Patchy Environment with Water and Human Movement. Mathematical Biosciences, 246, 105-112.
[4]	Hartley, D.M., Morris Jr, J.G. and Smith, D.L. (2006) Hyperinfectivity: A Critical Element in the Ability of V. cholerae to Cause Epidemics? PLoS Medicine, 3, e7.
[5]	Shuai, Z. and Van den Driessche, P. (2011) Global Dynamics of Cholera Models with Differential Infectivity. Mathematical Biosciences, 234, 118-126.
[6]	Lipp, E.K., Huq, A. and Colwell, R.R. (2002) Effects of Global Climate on Infectious Disease: The Cholera Model. Clinical Microbiology Reviews, 15, 757-770.
[7]	Andrews, J.R. and Basu, S. (2011) Transmission Dynamics and Control of Cholera in Haiti: An Epidemic Model. The Lancet, 377, 1248-1255.
[8]	Wang, X., Gao, D. and Wang, J. (2015) Influence of Human Behavior on Cholera Dynamics. Mathematical Biosciences, 267, 41-52.
[9]	Duan, L. and Xu, Z. (2020) A Note on the Dynamics Analysis of a Diffusive Cholera Epidemic Model with Nonlinear Incidence Rate. Applied Mathematics Letters, 106, 106356.
[10]	Liao, S., Yang, W. and Fang, F. (2021) Traveling Waves for a Cholera Vaccination Model with Nonlocal Dispersal. Mathematical Methods in the Applied Sciences, 44, 5150-5171.
[11]	Posny, D., Wang, J., Mukandavire, Z., et al. (2015) Analyzing Transmission Dynamics of Cholera with Public Health Interventions. Mathematical Biosciences, 264, 38-53.
[12]	Tian, X., Xu, R. and Lin, J. (2019) Mathematical Analysis of a Cholera Infection Model with Vaccination Strategy. Applied Mathematics and Computation, 361, 517-535.
[13]	Cai, L.M., Modnak, C. and Wang, J. (2017) An Age-Structured Model for Cholera Control with Vaccination. Ap-plied Mathematics and Computation, 299, 127-140.
[14]	Lin, J., Xu, R. and Tian, X. (2018) Global Dynamics of an Age-Structured Cholera Model with Both Human-to-Human and Environment-to-Human Transmissions and Saturation Incidence. Applied Mathematical Modelling, 63, 688-708.
[15]	Shuai, Z., Tien, J.H. and Van den Driessche, P. (2012) Cholera Models with Hyperinfectivity and Temporary Immunity. Bulletin of Mathematical Biology, 74, 2423-2445.
[16]	Wang, Y. and Wei, J. (2013) Global Dynamics of a Cholera Model with Time Delay. International Journal of Biomathematics, 6, 1250070.
[17]	Wang, X., Chen, Y. and Song, X. (2019) Global Dynamics of a Cholera Model with Age Structures and Multiple Transmission Modes. International Journal of Biomathematics, 12, 1950051.
[18]	Sun, G.Q., Xie, J.H., Huang, S.H., et al. (2017) Transmission Dynamics of Cholera: Mathematical Modeling and Control Strategies. Communications in Nonlinear Science and Numerical Simulation, 45, 235-244.
[19]	Lin, J., Xu, R. and Tian, X. (2019) Transmission Dynamics of Cholera with Hyperinfectious and Hypoinfectious Vibrios: Mathematical Modelling and Control Strategies. Mathematical Biosciences and Engineering, 16, 4339-4358.
[20]	Diekmann, O., Heesterbeek, J.A.P. and Metz, J.A.J. (1990) On the Definition and the Computation of the Basic Reproduction Ratio R0 in Models for Infectious Diseases in Heterogeneous Populations. Journal of Mathematical Biology, 28, 365-382.
[21]	Van den Driessche, P. and Watmough, J. (2002) Reproduction Numbers and Sub-Threshold Endemic Equilibria for Compartmental Models of Disease Transmission. Mathematical Biosciences, 180, 29-48.
[22]	Heesterbeek, J.A.P. and Roberts, M.G. (2007) The Type-Reproduction Number T in Models for Infectious Disease Control. Mathematical Biosciences, 206, 3-10.
[23]	Roberts, M.G. and Heesterbeek, J.A.P. (2003) A New Method for Estimating the Effort Required to Control an Infectious Disease. Proceedings of the Royal Society of London. Series B: Biological Sciences, 270, 1359-1364.
[24]	Shuai, Z., Heesterbeek, J.A.P. and van Den Driessche, P. (2013) Extending the Type Reproduction Number to Infectious Disease Control Targeting Contacts between Types. Journal of Mathematical Biology, 67, 1067-1082.
[25]	Cai, L.M. and Li, X.Z. (2009) Analysis of a SEIV Epidemic Model with a Nonlinear Incidence Rate. Applied Mathematical Modelling, 33, 2919-2926.
[26]	Waltman, P. (2006) A Brief Survey of Persistence in Dynamical Systems. Delay Differential Equations and Dynamical Systems: Proceedings of a Conference in Honor of Kenneth Cooke Held in Claremont, California, Jan. 13-16, 1990. Berlin, Heidelberg: Springer Berlin Heidelberg.
[27]	Li, M.Y., Graef, J.R., Wang, L., et al. (1999) Global Dynamics of a SEIR Model with Varying Total Population Size. Mathematical Biosciences, 160, 191-213.
[28]	Hethcote, H.W. (2000) The Mathematics of Infectious Diseases. SIAM Review, 42, 599-653.
[29]	Wu, H., Miao, O., Warnes, G.R., et al. (2008) Dediscover: A Computation and Simulation Tool for HIV Viral Fitness Research. 2008 International Conference on BioMedical Engineering and Informatics, 1, 687-694.