Mathematical analysis of cholera typhoid co-infection transmission dynamics

1. Introduction

Cholera, an acute gastro-intestinal water-borne infection, is caused by the bacterium Vibrio Cholerae, V. cholerae O1 or O139. Some of the symptoms are vomiting and diarrhea. If treatment is delayed, it can lead to severe dehydration and death within a few hours. The disease has two modes of transmission: direct and indirect transmission. Direct transmission (human-human) is very uncommon, whilst indirect transmission (environment-human), which occurs through the ingestion of contaminated food or water [1], is more frequent. Known estimates of the incubation period for the cholera disease is 1.4 days [2]. On the other hand, the Salmonella Typhi bacteria is responsible for causing the life threatening typhoid fever disease. Cholera and typhoid fever have the same transmission modes. The recticuloendothelical system, the intestinal lymphoid, and the gall bladder are severely damaged by the typhoid fever disease. Once a susceptible individuals is infected by the disease, roughly 19 days are required for the disease to incubate within the host [3].

Mathematical models have been used for the past decades to give insights into the transmission dynamics of co-infections within the human population. Akinyi et al. [4], showed that whenever the basic reproduction number is lowered to below one, then the malaria and the pneumonia cases will be reduced in a model of malaria-pneumonia coinfection. Onyinge et al. [5] modeled the co-dynamics of pneumonia and HIV, and they showed that the model was mathematically and epidemiologically sound; Mushayabasa et al. [6] modeled malaria-typhoid co-infection and demonstrated that a typhoid outbreak will inevitably lead to a spike in the malaria cases.

A number of mathematical models on typhoid have been proposed by a number of researchers. Mushayabasa [7], modeled how vaccines can help mitigate the spread of typhoid in Ghana. Pitzer et al. [8], extended the work in Mushayabasa [7] by applying the model to South Asia. Khan et al. [9], studied the typhoid disease with a saturated incidence rate.

To the best of our knowledge, the co-dynamics of typhoid and cholera have not been investigated in the literature. A recent outbreak of these two infections in Zimbabwe prompted this theoretical inquiry into how these infections interact. Due to the complicated nature of the co-infection model, we begin our analysis by studying the underlying sub-models; namely, the cholera only and the typhoid only sub-models. For each of the models, a number of pertinent questions are investigated. The questions explored include: Which factors in the models are key to decreasing the prevalence of each disease and the co-infection? Within the population, are the infections in competition with each other, or are they symbiotic? The implications of the results to the public health are discussed.

The paper is arranged as follows; the development of the model and the properties of the basic reproduction number are established in Section 2. Section 3 contains the stability analysis of the model at the fixed points. Numerical simulations and parameter estimations are done in section 4. Section 5 concludes the articles.

2. Methodology

2.1. Model development

Our typhoid cholera co-infection model partitions the human population N(t), at time t, into a susceptible class S(t), a cholera infection class I_c(t), a typhoid infection class I_t(t), a coinfection class I_ct(t), a cholera recovery class R_c(t), a typhoid recovery class R_c(t), and a coinfection recovery class R_ct(t). Thus,

$\begin{array}{l}N(t)=S(t)+I_c(t)+R_c(t)+I_t(t)+R_t(t)+I_{ct}(t)+R_{ct}(t).\end{array}$

The bacterial concentration of Salmonella Typhi, B_t(t), and Vibrio Cholerae, B_c(t), in the environment are incorporated into the model as well. The formulation of this model is an extension to the work carried out by Matsebula et al. [10].

Since the incubation periods of the two infections are different, we assume that dually infected individuals can only transmit either cholera or typhoid but not both infections simultaneously. Transmission of cholera to susceptible individuals occurs in one of two routes—the direct transmission route (human-to-human) and the indirect transmission route (envirnment-to-human). The rates of the transmission routes, respectively, are given by

$\begin{array}{l}{\lambda }_{c_1}=\frac{{\beta }_{c_1}(I_c+{\eta }_c{I}_{ct})}{N},{\lambda }_{c_2}=\frac{{\beta }_{c_2}{B}_c}{B_c+{\kappa }_c}.\end{array}$

The parameter β_c1 denotes the person-to-person cholera transmission. The effective contact rate for cholera multiplied by the probability of cholera transmission per contact gives the person to person cholera transmission. The modification parameter η_c, accounts for the relative infectiousness of individuals in class I_c relative to individuals in class I_ct. We assume that η_c ∈ (0, 1). This assumption is motivated by the fewer numbers of co-infected individuals as compared to those infected with cholera only. The parameter β_c2 denotes the environment-to-humans per capita contact rate and the Vibrio Cholerae in the contaminated environment, whilst the parameter κ_c denotes the half saturation constant of the Vibrio Cholerae. The half saturation constant is the bacterial concentration that is required to support half of the maximum rate, β_c2.

Similarly, the transmission of typhoid to susceptible individuals occurs in one of two routes—the direct transmission route (human-to-human) and the indirect transmission route (envirnment-to-human). The rates of the transmission routes, respectively, are given by

$\begin{array}{l}{\lambda }_{t_1}=\frac{{\beta }_{t_1}(I_t+{\eta }_t{I}_{ct})}{N},{\lambda }_{t_2}=\frac{{\beta }_{t_2}{B}_t}{B_t+{\kappa }_t}.\end{array}$

The parameter β_t1 denotes the effective person-to-person typhoid transmission rate. The effective contact rate for typhoid multiplied by the probability of typhoid transmission per contact gives the person to person typhoid transmission. The modification parameter η_t, accounts for the relative infectiousness of individuals in class I_t relative to individuals in class I_ct. We also assume that η_t ∈ (0, 1) following the assumptions given in the cholera infection dynamics. The per capita contact rate between the susceptibles and Salmonela typhi is represented by β_t2, and the half saturation constant for λ_c2 is κ_t.

Transmission of cholera to typhoid infected individuals occurs in one of two routes—the direct transmission route (human-to-human) and the indirect transmission route (environment-to-human). The rates of the transmission routes, respectively, are given by

$\begin{array}{l}{\lambda }_{c_3}=\frac{{\beta }_{c_3}(I_c+{\eta }_c{I}_{ct})}{N},{\lambda }_{c_4}=\frac{{\beta }_{c_4}{B}_c}{B_c+{\kappa }_c}.\end{array}$

The parameter β_c3 denotes the effective person-to-person cholera transmission rate of individuals in class I_t. The parameter β_c4 denotes the environment-to-humans per capita contact rate for individuals in class I_t and the Vibrio Cholerae in the contaminated environment. Transmission of typhoid to cholera infected individuals occurs in one of two routes—the direct transmission route (human-to-human) and the indirect transmission route (environment-to-human). The rates of the transmission routes, respectively, are given by

$\begin{array}{l}{\lambda }_{t_3}=\frac{{\beta }_{t_3}(I_t+{\eta }_t{I}_{ct})}{N},{\lambda }_{t_4}=\frac{{\beta }_{t_4}{B}_t}{B_t+{\kappa }_t}.\end{array}$

The parameter β_t3 denotes the person-to-person typhoid transmission rate of individuals in class I_c. The parameter β_t4 denotes the environment-to-humans per capita contact rate for individuals in class I_c and the Salmonella Typhi in the contaminated environment.

Infected individuals in classes I_c, I_t and I_ct experience disease related death at rates given, respectively by δ_c, δ_t and δ_ct. Individuals in the infectious states I_c and I_t, respectively, excrete Vibrio Cholerae bacteria and Salmonella Typhi bacteria into the environment at rates α_c and α_t. Coinfected individuals shed Vibrio Cholerae and Salmonella Typhi into the environments at rates θ_c and θ_t, respectively. Infection is assume to confer temporary immunity. The cholera and typhoid immunity wanes at rates ρ_c, ρ_t and ρ_ct.

The generation rate of Vibrio Cholerae is $g_c{B}_c\left(1-\frac{B_c}{k_c}\right)$, and its growth is enhanced by cholera infected individuals and the coinfected individuals that are shedding into the environment. The generation rate of Salmonella Typhi is $g_t{B}_t\left(1-\frac{B_t}{k_t}\right)$ and its growth is enhanced by typhoid infected individuals and the coinfected individuals that are shedding into the environment. We assume that the Vibrio Cholerae and the Salmonella Typhi bacteria in the environment are respectively removed by interventions such as improved sanitation and treatment of contaminated environments at rates μ_c and μ_t. The parameter Λ represents the recruitment into the susceptibles, while the parameter μ represents the natural death rate. It is assumes that individuals mix homogeneously and that they are indistinguishable in each of the classes. The model diagram is shown in Figure 1.

FIGURE 1

Figure 1. The cholera typhoid co-infection compartmental model. For the concise presentation of our model flow diagram, we make use of the following expressions: ${\chi }_1=g_c{B}_c\left(1-\frac{B_c}{k_c}\right)+{\alpha }_c{I}_c+{\theta }_c{I}_{ct}$, ${\chi }_2=g_t{B}_t\left(1-\frac{B_t}{k_t}\right)+{\alpha }_t{I}_t+{\theta }_t{I}_{ct}$, ${\lambda }_c={\lambda }_{c_1}+{\lambda }_{c_2}$, ${\lambda }_t={\lambda }_{t_1}+{\lambda }_{t_2}$, ${\lambda }_1={\lambda }_{t_3}+{\lambda }_{t_4}$, ${\lambda }_2={\lambda }_{c_3}+{\lambda }_{c_4}$.

The dynamical system associated with the schematic diagram in Figure 1 is;

$\begin{array}{l}\frac{dS}{dt}=\Lambda -({\lambda }_{c_1}+{\lambda }_{c_2}+{\lambda }_{t_1}+{\lambda }_{t_2})S-\mu S+{\rho }_c{R}_c+{\rho }_t{R}_t+{\rho }_{ct}{R}_{ct},\\ \frac{d{I}_c}{dt}=({\lambda }_{c_1}+{\lambda }_{c_2})S-({\lambda }_{t_3}+{\lambda }_{t_4})I_c-(\mu +{\delta }_c+{ϵ}_c)I_c,\\ \frac{d{I}_t}{dt}=({\lambda }_{t_1}+{\lambda }_{t_2})S-({\lambda }_{c_3}+{\lambda }_{c_4})I_t-(\mu +{\delta }_t+{ϵ}_t)I_t,\\ \frac{d{I}_{ct}}{dt}=({\lambda }_{t_3}+{\lambda }_{t_4})I_c+({\lambda }_{c_3}+{\lambda }_{c_4})I_t-(\mu +{\delta }_{ct}+{ϵ}_{ct})I_{ct},\\ \frac{d{R}_c}{dt}={ϵ}_c{I}_c-(\mu +{\rho }_c)R_c,\\ \frac{d{R}_t}{dt}={ϵ}_t{I}_t-(\mu +{\rho }_t)R_t,\\ \frac{d{R}_{ct}}{dt}={ϵ}_{ct}{I}_{ct}-(\mu +{\rho }_{ct})R_{ct},\\ \frac{d{B}_c}{dt}=g_c{B}_c\left(1-\frac{B_c}{k_c}\right)+{\alpha }_c{I}_c+{\theta }_c{I}_{ct}-{\mu }_c{B}_c,\\ \frac{d{B}_t}{dt}=g_t{B}_t\left(1-\frac{B_t}{k_t}\right)+{\alpha }_t{I}_t+{\theta }_t{I}_{ct}-{\mu }_t{B}_t,& (1)\end{array}$

with initial conditions

$\begin{array}{l}S(0)=S_0>0,B_c(0)=B_{c0}\ge 0,B_t(0)=B_{t0}\ge 0,\\ I_c(0)=I_{c0}\ge 0,I_t(0)=I_{t0}\ge 0,I_{ct}(0)=I_{ct0}\ge 0,\\ R_c(0)=R_{c0}\ge 0,R_t(0)=R_{t0}\ge 0,R_{ct}(0)=R_{ct0}\ge 0.\end{array}$

2.2. Cholera only model

We define the cholera only model as the model obtained from setting all the typhoid classes and its associated parameters to zero. We thus have the following

$\begin{array}{l}\frac{dS}{dt}=\Lambda -({\tilde{\lambda }}_{c_1}+{\lambda }_{c_2})S-\mu S+{\rho }_c{R}_c,\\ \frac{d{I}_c}{dt}=({\tilde{\lambda }}_{c_1}+{\lambda }_{c_2})S-q_c{I}_c,\\ \frac{d{R}_c}{dt}={ϵ}_c{I}_c-(\mu +{\rho }_c)R_c,\\ \frac{d{B}_c}{dt}=g_c{B}_c(1-\frac{B_c}{k_c})+{\alpha }_c{I}_c-{\mu }_c{B}_c,& (2)\end{array}$

where

$\begin{array}{l}{\tilde{\lambda }}_{c_1}=\frac{{\beta }_{c_1}{I}_c}{N_c},q_c=\mu +{\delta }_c+{ϵ}_c,N_c=S+I_c+R_c,\end{array}$

with initial conditions

$\begin{array}{l}S(0)=S_0>0,B_c(0)=B_{c0}\ge 0,I_c(0)=I_{c0}\ge 0,\\ R_c(0)=R_{c0}\ge 0.\end{array}$

2.2.1. Boundedness and non-negative trajectories

We argue that model (Equation 2) yields bounded-non negative-trajectories in this section provided the initial conditions are non-negative.

Theorem 1. All solutions of the cholera only sub-model (Equation 2) are non-negative if all the initial conditions are non-negative.

Proof. Define t₁ = sup{t > 0|S(τ₁) > 0, I_c(τ₁) ≥ 0, R_c(τ₁) ≥ 0, B_c(τ₁) ≥ 0, ∀τ₁ ∈ [0, t]}. It follows t₁ > 0 since

S₀ > 0, I_c0 ≥ 0, R_c0 ≥ 0, B_c0 ≥ 0.

Assume t₁ < ∞, then S(t₁) > 0, I_c(t₁) = 0, R_c(t₁) = 0, B_c(t₁) = 0. Applying variation of constants to

$\frac{dS}{dt}=\Lambda -({\tilde{\lambda }}_{c_1}+{\lambda }_{c_2})S-\mu S+{\rho }_c{R}_c,$

yields

$S(t_1)={\displaystyle {\int }_0^{t_1}f(r)\exp \left(-{\displaystyle {\int }_0^{t_1}P(x)dx}\right)}dr+S_0\exp \left(-{\displaystyle {\int }_0^{t_1}P(x)dx}\right),$

where $P(x)=({\tilde{\lambda }}_{c_1}+{\lambda }_{c_2}+\mu )$ and f(r) = Λ + ρ_cR_c. Clearly,

$S(t_1)>0$

Since f(r) > 0 and P(x) > 0 when x, r ∈ [0, t₁]. Similarly, I_c(t₁) > 0 and R_c(t₁) > 0. This produces a contradiction, hence t₁ = ∞.

Theorem 2. All solutions of the cholera only sub-model (Equation 2) are bounded within Ω whenever g_c ≥ μ_c.

Proof. The time derivative of the population for the cholera model (Equation 2) is bounded above by

$\begin{array}{l}\frac{d{N}_c}{dt}=\Lambda -\mu N_c-{\delta }_c{I}_c\le \Lambda -\mu N_c,\end{array}$

Upper bounds for the human population, N_c(t), are obtained by integrating the separable differential inequality as follows,

$\begin{array}{l}{N}_c\le \frac{\Lambda -Mexp(-\mu t)}{\mu }\le \frac{\Lambda }{\mu }.\end{array}$

By extension, Λ/μ is also the upper bound for each of the human classes. Whereas, owing to I_c ≤ N_c ≤ Λ/μ, an upper bound for the bacterial classes can be obtained as follows,

$\begin{array}{ll}\frac{d{B}_c}{dt}=g_c{B}_c(1-\frac{B_c}{k_c})+{\alpha }_c{I}_c-{\mu }_c{B}_c\le g_c{B}_c(1-\frac{B_c}{k_c})+{\alpha }_c\frac{\Lambda }{\mu }-{\mu }_c{B}_c.& (3)\end{array}$

From inequality (3), if

$\begin{array}{ll}{B}_c\ge {\alpha }_c\frac{\Lambda }{\mu },& (4)\end{array}$

Where ${\alpha }_c\frac{\Lambda }{\mu }$ is the maximum shedding rate from the cholera infected individuals, then

$\begin{array}{l}\frac{d{B}_c}{dt}\le (g_c-{\mu }_c)B_c-\frac{g_c}{k_c}{B}_c^2+B_c\\ =(g_c-{\mu }_c+1)B_c(1-\frac{g_c{B}_c}{k_c(g_c-{\mu }_c+1)}).& (5)\end{array}$

The constant

$\begin{array}{ll}\frac{k_c(g_c-{\mu }_c+1)}{g_c},& (6)\end{array}$

is the upper bound for the differential inequality (Equation 5) since (Equation 5) is the logistic growth model with carrying capacity (Equation 6). For some t ≥ 0, (α_c + θ_c)Λ/μ is an upper bound for B_c whenever (Equation 4) is false, whilst B_c is bounded above by Equation (6) for the rest of the time points in the domain of B_c if (Equation 4) is true. Thus, in both cases, $B_c\le max\left\{\frac{k_c(g_c-{\mu }_c+1)}{g_c},{\alpha }_c\frac{\Lambda }{\mu }\right\}$.

Within the feasible region,

$\begin{array}{l}{\Omega }_c=\{(S,I_c,R_c,B_c)|0\le N_c\le \frac{\Lambda }{\mu }],\\ B_c\in \left[0,\max \{\frac{k_c(g_c-{\mu }_c+1)}{g_c},{\alpha }_c\frac{\Lambda }{\mu }\}\right]\},\end{array}$

We have summarized the results on the boundedness and positivity of the solutions of the cholera only sub-model 2.

2.2.2. The stability of the disease free equilibrium and the reproduction number, ${ℛ}_C$

The disease free equilibria of system (Equation 2) is given by

$\begin{array}{l}{\text{x}}_0=(S,I_c,R_c,B_c)=\left(\frac{\Lambda }{\mu },0,0,0\right),\\ {\text{x}}_1=(S,I_c,R_c,B_c)=(c_1,0,0,c_2),\end{array}$

where $c_1=\frac{\Lambda (g_c-{\mu }_c)k_c+g_c{\kappa }_c}{(\mu +{\beta }_{c_2})(g_c-{\mu }_c)k_c+\mu {\kappa }_{c}gc}$ and $c_2=\frac{g_c-{\mu }_c}{g_c}{k}_c$. The equilibrium x₁ > 0 if g_c > μ_c. The Jacobian of dynamical system (Equation 2) is given by

$\begin{array}{l}J({\text{x}}_0)=\left(\begin{array}{cccc}-\mu & -{\beta }_{c_1}& {\rho }_c& -\frac{\Lambda {\beta }_{c_2}}{\mu {\kappa }_c}\\ 0& {\beta }_{c_1}-q_c& 0& \frac{\Lambda {\beta }_{c_2}}{\mu {\kappa }_c}\\ 0& ϵ& -(\mu +{\rho }_c)& 0\\ 0& {\alpha }_c& 0& g_c-{\mu }_c\end{array}\right).\end{array}$

The dynamical system (Equation 2) is locally asymptotically stable if all four of its eigenvalues have negative real parts. Two of the eigenvalues for the Jacobian, J, are λ₁ = −μ and λ₂ = −(μ + ρ_c). The other two eigenvalues for J are the eigenvalues from the sub-matrix

$\begin{array}{l}\bar{J}=\left(\begin{array}{cc}{\beta }_{c_1}-q_c& \frac{\Lambda {\beta }_{c_2}}{\mu {\kappa }_c}\\ {\alpha }_c& g_c-{\mu }_c\end{array}\right).\end{array}$

The characteristic equation for matrix $\bar{J}$ is ${\lambda }^2+{\bar{\nu }}_1\lambda +{\bar{\nu }}_2$, where

$\begin{array}{l}{\overline{\nu }}_1=-\left((g_c-{\mu }_c)+({\beta }_{c_1}-q_c)\right),\\ {\overline{\nu }}_2=({\beta }_{c_1}-q_c)(g_c-{\mu }_c)\left(1-\frac{{\alpha }_c{\beta }_{c_2}\Lambda }{({\beta }_{c_1}-q_c)(g_c-{\mu }_c){\kappa }_c\mu }\right)\\ =({\beta }_{c_1}-q_c)(g_c-{\mu }_c)\left(1-{ℛ}_C\right),\end{array}$

and

$\begin{array}{l}{ℛ}_C=\frac{{\alpha }_c{\beta }_{c_2}\Lambda }{{\kappa }_c\mu {\mu }_c{q}_c(1-{ℛ}_h)(1-{ℛ}_b)},{ℛ}_h=\frac{{\beta }_{c_1}}{q_c}{ℛ}_b=\frac{g_c}{{\mu }_c}.\end{array}$

The constants ${ℛ}_b$ and ${ℛ}_h$ are the bacterial regeneration threshold and the human-to-human sub reproduction number, respectively. The constant ${ℛ}_C$ is the so-called basic reproduction number for the system (Equation 2). Clearly, ${ℛ}_h,{ℛ}_b<1$ or ${ℛ}_h,{ℛ}_b>1$ if and only if ${ℛ}_C>0$.

It follows from the Routh Hurwitz criterion that the two eigenvalues of $\bar{J}$ have negative real parts if ${\bar{\nu }}_1,{\bar{\nu }}_2>0$. It is easy to see that ${\bar{\nu }}_1,{\bar{\nu }}_2>0$ if ${ℛ}_h<1$, ${ℛ}_h<1$ and ${ℛ}_C<1$. Hence, a positive basic reproduction number for system (Equation 2) that is less than unity implies that the system is locally asymptotically stable at the disease free equilibrium.

2.2.3. Endemic equilibrium

Setting the derivatives of the classes to zero gives the endemic equilibrium for the cholera only sub-model (Equation 2). Let $\lambda ={\tilde{\lambda }}_{c_1}+{\lambda }_{c_2}$.

$\begin{array}{ll}0=\Lambda -({\lambda }^{*}+\mu )S^{*}+{\rho }_c{R}_c^{*},& (7)\end{array}$

$\begin{array}{ll}0={\lambda }^{*}{S}^{*}-q_c{I}_c^{*},& (8)\end{array}$

$\begin{array}{ll}0={ϵ}_c{I}_c^{*}-(\mu +{\rho }_c)R_c^{*},& (9)\end{array}$

$\begin{array}{ll}0=g_c{B}_c^{*}\left(1-\frac{B_c^{*}}{k_c}\right)+{\alpha }_c{I}_c^{*}-{\mu }_c{B}_c^{*}.& (10)\end{array}$

From Equation (8),

$S^{*}=\frac{q_c}{{\lambda }^{*}}{I}_c^{*}.$

From Equation (7),

$R_c^{*}=\frac{1}{{\rho }_c}\left[\frac{(\lambda *+\mu )q_c{I}_c^{*}}{{\lambda }^{*}}-\Lambda \right].$

Consider (Equation 9),

$ϵI_c^{*}-\frac{(\mu +\rho c)}{{\rho }_c}\left[\frac{(\lambda *+\mu )q_c{I}_c^{*}}{{\lambda }^{*}}-\Lambda \right]=0,$

therefore

$I_c^{*}=\frac{{\lambda }^{*}\Lambda (\mu +\rho c)}{q_c(\lambda *+\mu )(\mu +\rho c)-{ϵ}_c{\rho }_c{\lambda }^{*}}.$

Given that

${\lambda }^{*}=\frac{{\beta }_{c_1}{I}_c^{*}}{S+I_c+R_c}+\frac{{\beta }_{c_2}{B}_c}{B_c+{\kappa }_c}.$

Using (Equation 10), we have a quadratic equation in B_c of the form

${\bar{\nu }}_2{B}_c^2+{\nu }_1{B}_c+{\nu }_0=0,$

where

$\begin{array}{l}{\nu }_2=g_c\left[q_c(\lambda *+\mu )(\mu +\rho c)-{ϵ}_c{\rho }_c{\lambda }^{*}\right],\\ {\nu }_1=-{\mu }_c{\kappa }_c({ℛ}_b-1){\nu }_2,{\nu }_0=-{\lambda }^{*}\Lambda {\kappa }_c{\alpha }_c(\mu +\rho c),\end{array}$

with

${ℛ}_b=\frac{g_c}{{\mu }_c}.$

Clearly, ν₀ < 0, ν₁ < 0 if ${ℛ}_b>1$. Since

$\begin{array}{ll}{B}_c=\frac{-{\nu }_1\pm \sqrt{{\nu }_1^2-4{\nu }_2{\nu }_0}}{2{\nu }_2},& (11)\end{array}$

it follows that if ${\nu }_2<0,{ℛ}_b>1$, then it follows from Descartes' rule of signs that B_c has no positive roots, and if ${\nu }_2>0,{ℛ}_b>1$, then B_c has only one positive root. We shall call the positive root $B_c^{+}$.

Let

$\begin{array}{l}B=\frac{B_c^{+}}{B_c^{+}+{\kappa }_c}.\end{array}$

Then

$\begin{array}{l}{\lambda }_{c_2}^{*}=\frac{{\beta }_{c_2}{B}_c^{+}}{B_c^{+}+{\kappa }_c}={\beta }_{c_2}B.\end{array}$

We have an expression for λ^* such that

$\begin{array}{ll}{a}_2{\lambda }^{*2}+a_1{\lambda }^{*}+a_0=0& (12)\end{array}$

where,

$\begin{array}{l}{a}_2=\mu +{ϵ}_c+{\rho }_c>0,a_1=q_c(\mu +{\rho }_c)-({\beta }_{c_1}(\mu +{\rho }_c)\\ +B{\beta }_{c_2}(\mu +{ϵ}_c+{\rho }_c)),\\ a_0=-B{q}_c{\beta }_{c_2}(\mu +{\rho }_c)<0.\end{array}$

Since

$\begin{array}{l}\lambda =\frac{-a_1\pm \sqrt{a_1^2-4a_2{a}_0}}{2a_2},\end{array}$

it follows that if a₁ > 0, then the polynomial (Equation 12) has a positive root, and if a₁ < 0, then the polynomial (Equation 12) has a positive root. So the polynomial (Equation 12) will always have one positive root.

So system (Equation 2) has a unique endemic equilibrium if ${ℛ}_b>1.$

Remark 1. Due to the symmetric structure of the cholera only and typhoid only sub-models, the typhoid only sub-model has similar structural results to those obtained for the cholera only sub-model. To avoid repetition, we have not shown the analysis of the typhoid only sub-model.

2.3. Cholera-typhoid co-infection model

We study the full co-infection model (Equation 1) in this section.

2.3.1. Non-negative trajectories and boundedness

We prove in this subsection that model (Equation 1) has non-negative trajectories.

Theorem 3. All solutions of the co-infection model (Equation 1) are non-negative if all the initial conditions are non-negative.

Proof. Define t₁ = sup{t > 0|S(τ₁) > 0, I_c(τ₁) ≥ 0, I_t(τ₁) ≥ 0, I_ct(τ₁) ≥ 0, R_c(τ₁) ≥ 0, R_t(τ₁) ≥ 0, R_ct(τ₁) ≥ 0, B_c(τ₁) ≥ 0, B_t(τ₁) ≥ 0, ∀τ₁ ∈ [0, t]}. It follows that t₁ > 0 since S₀ > 0, I_c0 ≥ 0, I_t0 ≥ 0, I_ct0 ≥ 0, R_c0 ≥ 0, R_t0 ≥ 0, R_t0 ≥ 0, B_c0 ≥ 0, B_t0 ≥ 0.

Assume t₁ < ∞, then S(t₁) > 0, I_c(t₁) = 0, I_t(t₁) = 0, I_ct(t₁) = 0, R_c(t₁) = 0, R_t(t₁) = 0, R_ct(t₁) = 0, B_c(t₁) = 0, B_t(t₁) = 0. Applying variation of constants to

$\frac{dS}{dt}=\Lambda -({\lambda }_{c_1}+{\lambda }_{c_2}+{\lambda }_{t_1}+{\lambda }_{t_2})S-\mu S+{\rho }_c{R}_c+{\rho }_t{R}_t+{\rho }_{ct}{R}_{ct}$

yields

$S(t_1)={\displaystyle {\int }_0^{t_1}f(r)\exp }\left(-{\displaystyle {\int }_0^{t_1}P(x)dx}\right)dr+S_0\exp \left(-{\displaystyle {\int }_0^{t_1}P(x)dx}\right),$

Where P(x) = λ_c1 + λ_c2 + λ_t1 + λ_t2 + μ and f(r) = Λ + ρ_cR_c + ρ_tR_t + ρ_ctR_ct. Clearly,

$S(t_1)>0$

since f(r) > 0 and P(x) > 0 when x, r ∈ [0, t₁]. Similarly, I_c(t₁) > 0, I_t(t₁) > 0, I_ct(t₁) > 0, R_t(t₁) > 0 and R_c(t₁) > 0. This produces a contradiction, hence t₁ = ∞.

Theorem 4. All solutions of the co-infection model (Equation 1) are bounded within Ω whenever g_c ≥ μ_c and g_t ≥ μ_t.

Proof. Since δ_ct(I_c + I_ct) ≥ 0, it follows that the upper bound for the time derivative of the total human population, N(t), is

$\begin{array}{l}\frac{dN}{dt}=\Lambda -\mu N-{\delta }_{ct}(I_c+I_{ct})\le \Lambda -\mu N.\end{array}$

Using separation of variables, we obtain the following upper bound for the human population,

$\begin{array}{l}N\le \frac{\Lambda -Mexp(-\mu t)}{\mu }\le \frac{\Lambda }{\mu }.\end{array}$

This upper bound for the population implies that each of the classes are also bounded above by the same constant Λ/μ. Since I_c, I_ct ≤ Λ/μ, it follows that the upper bound for the bacterial concentration of Vibros Cholerae is bounded above by

$\begin{array}{l}\frac{d{B}_c}{dt}=g_c{B}_c\left(1-\frac{B_c}{k_c}\right)+{\alpha }_c{I}_c+{\theta }_c{I}_{ct}-{\mu }_c{B}_c\\ \le g_c{B}_c\left(1-\frac{B_c}{k_c}\right)+({\alpha }_c+{\theta }_c)\frac{\Lambda }{\mu }-{\mu }_c{B}_c.& (13)\end{array}$

From inequality (Equation 13), if

$\begin{array}{lll}\hfill B_c\ge ({\alpha }_c+{\theta }_c)\frac{\Lambda }{\mu },& & (14)\end{array}$

then

$\begin{array}{ll}\hfill \frac{d{B}_c}{dt}\le (g_c-{\mu }_c)B_c-\frac{g_c}{k_c}{B}_c^2+B_c& \\ \hfill =(g_c-{\mu }_c+1)B_c\left(1-\frac{g_c{B}_c}{k_c(g_c-{\mu }_c+1)}\right).& (15)\end{array}$

The constant

$\begin{array}{ll}\frac{k_c(g_c-{\mu }_c+1)}{g_c},& (16)\end{array}$

is the upper bound for the differential inequality (Equation 15) since (Equation 15) is the logistic growth model with carrying capacity (Equation 16). For some t ≥ 0, (α_c + θ_c)Λ/μ is an upper bound for B_c whenever (Equation 14) is false, whilst B_c is bounded above by Equation (16) for the rest of the time points in the domain of B_c if (Equation 14) is true. The constant $({\alpha }_c+{\theta }_c)\frac{\Lambda }{\mu }$ is the maximum shedding rate from the cholera infected individuals and dually infected individuals. In both cases, $B_c\le max\left\{\frac{k_c(g_c-{\mu }_c+1)}{g_c},({\alpha }_c+{\theta }_c)\frac{\Lambda }{\mu }\right\}.$

$\begin{array}{ll}J({\text{X}}_0)=\left(\begin{array}{ccccccccc}-\mu & -{\beta }_{c_1}& -{\beta }_{t_1}& -({\beta }_{c_1}{\eta }_c+{\beta }_{t_1}{\eta }_t)& {\rho }_c& {\rho }_t& {\rho }_{ct}& -\frac{\Lambda {\beta }_{c_2}}{\mu {\kappa }_c}& -\frac{\Lambda {\beta }_{t_2}}{\mu {\kappa }_t}\\ 0& {\beta }_{c_1}-q_c& 0& {\eta }_c{\beta }_{c_1}& 0& 0& 0& \frac{\Lambda {\beta }_{c_2}}{\mu {\kappa }_c}& 0\\ 0& 0& {\beta }_{t_1}-q_t& {\eta }_t{\beta }_{t_1}& 0& 0& 0& 0& \frac{\Lambda {\beta }_{t_2}}{\mu {\kappa }_t}\\ 0& 0& 0& -(\mu +{\delta }_{ct}+{ϵ}_{ct})& 0& 0& 0& 0& 0\\ 0& {ϵ}_c& 0& 0& -(\mu +{\rho }_c)& 0& 0& 0& 0\\ 0& 0& {ϵ}_t& 0& 0& -(\mu +{\rho }_t)& 0& 0& 0\\ 0& 0& 0& {ϵ}_{ct}& 0& 0& -(\mu +{\rho }_{ct})& 0& 0\\ 0& {\alpha }_c& 0& {\theta }_c& 0& 0& 0& g_c-{\mu }_c& 0\\ 0& 0& {\alpha }_t& {\theta }_t& 0& 0& 0& 0& g_t-{\mu }_t\end{array}\right),& (17)\end{array}$

Within the feasible region,

$\begin{array}{l}{\Omega }_{co}=\{(S,I_c,I_t,I_{ct},R_c,R_t,R_{ct},B_c,B_t)\in {ℝ}_{+}^9|0\le N\le \frac{\Lambda }{\mu },\\ B_c\in [0,\max \{\frac{k_c(g_c-{\mu }_c+1)}{g_c},({\alpha }_c+{\theta }_c)\frac{\Lambda }{\mu }\}],\\ B_t\in [0,\max \{\frac{k_t(g_t-{\mu }_t+1)}{g_t},({\alpha }_t+{\theta }_t)\frac{\Lambda }{\mu }\}]\},\end{array}$

We have summarized the results on the boundedness and positivity of the solutions to the co-infection model (Equation 1).

2.3.2. Stability analysis of the disease free equilibrium and reproduction number, ${ℛ}_0$

We find the conditions required for the disease free equilibrium for dynamical system (Equation 1) to be locally asymptotically stable in this section. The disease free equilibria of dynamical system (Equation 1) are

$\begin{array}{l}{\text{X}}_0=(S,I_c,I_t,I_{ct},R_c,R_t,R_{ct},B_c,B_t)=\left(\frac{\Lambda }{\mu },0,0,0,0,0,0,0,0\right),\\ {\text{X}}_1=(S,I_c,I_t,I_{ct},R_c,R_t,R_{ct},B_c,B_t)=(c_3,0,0,0,0,0,0,c_4,c_5),\\ {\text{X}}_2=(S,I_c,I_t,I_{ct},R_c,R_t,R_{ct},B_c,B_t)=(c_6,0,0,0,0,0,0,c_4,0),\\ {\text{X}}_3=(S,I_c,I_t,I_{ct},R_c,R_t,R_{ct},B_c,B_t)=(c_7,0,0,0,0,0,0,0,c_5),\end{array}$

where

$\begin{array}{l}{c}_3=\frac{\Lambda ((g_c-{\mu }_c)k_c+g_c{\kappa }_c)((g_t-{\mu }_t)k_t+g_t{\kappa }_t)}{\mu ((g_c-{\mu }_c)k_c+g_c{\kappa }_c)((g_t-{\mu }_t)k_t+g_t{\kappa }_t)+{\beta }_{c_2}{g}_t(g_c-{\mu }_c)k_c+{\beta }_{t_2}{g}_c(g_t-{\mu }_t)k_t},\end{array}$

$\begin{array}{l}{c}_4=\frac{g_c-{\mu }_c}{g_c}{k}_c,c_5=\frac{g_t-{\mu }_t}{g_t}{k}_t,\\ c_6=\frac{\Lambda ((g_c-{\mu }_c)k_c+g_c{\kappa }_c)(g_t{\kappa }_t)}{\mu ((g_c-{\mu }_c)k_c+g_c{\kappa }_c)(g_t{\kappa }_t)+{\beta }_{c_2}{g}_t(g_c-{\mu }_c)k_c},\\ c_7=\frac{\Lambda (g_c{\kappa }_c)((g_t-{\mu }_t)k_t+g_t{\kappa }_t)}{\mu (g_c{\kappa }_c)((g_t-{\mu }_t)k_t+g_t{\kappa }_t)+{\beta }_{t_2}{g}_c(g_t-{\mu }_t)k_t}.\end{array}$

It is observed that the disease free equilibria, X₁, X₂, X₃, are always unstable due to the condition, g_c ≥ μ_c, that is requires for their existence.

The Jacobian of the full system is

The dynamical system (Equation 1) is locally asymptotically stable if all nine of its eigenvalues have negative real parts. Five of the eigenvalues for the Jacobian, J, are λ₁ = −μ, λ₂ = −(μ + ρ_c), λ₃ = −(μ+ρ_t), λ₄ = −(μ+ρ_ct) and λ₅ = −(μ+δ_ct+ϵ_ct). The other four eigenvalues for J are the eigenvalues from the sub-matrix

$\begin{array}{l}\bar{J}=\left(\begin{array}{cccc}{\beta }_{c_1}-q_c& 0& \frac{\Lambda {\beta }_{c_2}}{\mu {\kappa }_c}& 0\\ 0& {\beta }_{t_1}-q_t& 0& \frac{\Lambda {\beta }_{t_2}}{\mu {\kappa }_t}\\ {\alpha }_c& 0& g_c-{\mu }_c& 0\\ 0& {\alpha }_t& 0& g_t-{\mu }_t\end{array}\right).\end{array}$

The characteristic equation for matrix $\bar{J}$ is $({\lambda }^2+{\nu }_1\lambda +{\nu }_2)({\lambda }^2+{\nu }_3\lambda +{\nu }_4)$, where

$\begin{array}{ll}{\nu }_1& =-((g_c-{\mu }_c)+({\beta }_{c_1}-q_c)),\\ {\nu }_2& =({\beta }_{c_1}-q_c)(g_c-{\mu }_c)(1-{ℛ}_C),\\ {\nu }_3& =-((g_t-{\mu }_t)+({\beta }_{t_1}-q_t)),\\ {\nu }_4& =({\beta }_{t_1}-q_t)(g_t-{\mu }_t)(1-{ℛ}_T),\end{array}$

and

$\begin{array}{l}{ℛ}_C=\frac{{\alpha }_c{\beta }_{c_2}\Lambda }{{\kappa }_c\mu {\mu }_c{q}_c(1-{ℛ}_h^c)(1-{ℛ}_b^c)},\\ {ℛ}_T=\frac{{\alpha }_t{\beta }_{t_2}\Lambda }{{\kappa }_t\mu {\mu }_t{q}_t(1-{ℛ}_h^t)(1-{ℛ}_b^t)},\\ {ℛ}_h^c=\frac{{\beta }_{c_1}}{q_c}{ℛ}_b^c=\frac{g_c}{{\mu }_c},{ℛ}_h^t=\frac{{\beta }_{t_1}}{q_t}{ℛ}_b^t=\frac{g_t}{{\mu }_t}.\end{array}$

The constants ${ℛ}_b^c$ and ${ℛ}_h^c$ are the bacterial regeneration threshold and the human-to-human sub reproduction number, respectively, for the cholera only sub-model. The constants ${ℛ}_b^t$ and ${ℛ}_h^t$ are the bacterial regeneration threshold and the human-to-human sub reproduction number, respectively, for the typhoid only sub-model. The constants ${ℛ}_C$ and ${ℛ}_T$ are the so-called basic reproduction numbers for the cholera only sub-model and the typhoid only sub-model, respectively. Clearly, ${ℛ}_h^c,{ℛ}_b^c<1$ or ${ℛ}_h^c,{ℛ}_b^c>1$ if and only if ${ℛ}_C>0$; Similarly, ${ℛ}_h^t,{ℛ}_b^t<1$ or ${ℛ}_h^t,{ℛ}_b^t>1$ if and only if ${ℛ}_T>0$.

We note that

$\begin{array}{l}{\nu }_2>({\beta }_{c_1}-q_c)(g_c-{\mu }_c)(1-max\left\{{ℛ}_C,{ℛ}_T\right\}),\\ {\nu }_4>({\beta }_{t_1}-q_t)(g_t-{\mu }_t)(1-max\left\{{ℛ}_C,{ℛ}_T\right\}).\end{array}$

Thus

${ℛ}_0=max\left\{{ℛ}_C,{ℛ}_T\right\}.$

The constant ${ℛ}_0$ is the basic reproduction number for the systems (Equation 1). It follows from the Routh Hurwitz criterion that the four eigenvalues of $\bar{J}$ have negative real parts if ν₁, ν₂, ν₃, ν₄ > 0. It is easy to see that ν₁, ν₂, ν₃, ν₄ > 0 if ${ℛ}_h^t<1,{ℛ}_b^t<1$, ${ℛ}_h^t<1,{ℛ}_b^t<1$ and ${ℛ}_0<1$. Hence, a positive basic reproduction number for system (Equation 1) that is less than unity implies that the system is locally asymptotically stable at the disease free equilibrium.

2.3.3. Impact analysis

In this section, we show how cholera affects typhoid, and through symmetry, we show how typhoid affects cholera.

The reproduction numbers for cholera and typhoid are

$\begin{array}{l}{ℛ}_C=\frac{{\alpha }_c{\beta }_{c_2}\Lambda }{{\kappa }_c\mu {\mu }_c{q}_c(1-{ℛ}_h^c)(1-{ℛ}_b^c)},\end{array}$

$\begin{array}{ll}{ℛ}_T=\frac{{\alpha }_t{\beta }_{t_2}\Lambda }{{\kappa }_t\mu {\mu }_t{q}_t(1-{ℛ}_h^t)(1-{ℛ}_b^t)},& (18)\end{array}$

respectively. These two reproduction numbers are dependent on each other. The constant, Λ/μ, allows for the expression of one reproduction number in terms of the other. From the second equation above, Equation (18), isolating, Λ/μ, yields

$\begin{array}{ll}{ℛ}_C={ℛ}_T\frac{{\alpha }_c{\beta }_{c_2}{\kappa }_t{\mu }_t{q}_t(1-{ℛ}_h^t)(1-{ℛ}_b^t)}{{\alpha }_t{\beta }_{t_2}{\kappa }_c{\mu }_c{q}_c(1-{ℛ}_h^c)(1-{ℛ}_b^c)}.& (19)\end{array}$

Differentiating ${ℛ}_C$ with respect to ${ℛ}_T$ gives

$\begin{array}{ll}\frac{\partial {ℛ}_C}{\partial {ℛ}_T}=\frac{{\alpha }_c{\beta }_{c_2}{\kappa }_t{\mu }_t{q}_t(1-{ℛ}_h^t)(1-{ℛ}_b^t)}{{\alpha }_t{\beta }_{t_2}{\kappa }_c{\mu }_c{q}_c(1-{ℛ}_h^c)(1-{ℛ}_b^c)}.& (20)\end{array}$

We conclude that an increase in cholera cases may be associated with an increase in typhoid cases, and an increase in typhoid cases may be associated with an increase in cholera cases. This conclusion is subject to the following conditions: firstly, the the bacterial regeneration threshold for both cholera and typhoid must be less than unity; secondly, the human-to-human sub reproduction number for both cholera and typhoid must also be less than unity. This result proves the symbiotic nature of the relationship between the typhoid disease and the Cholera disease.

3. Numerical simulations

In this section, we give a brief outline of the numerical results obtained in the investigation. Table 1 shows the parameters of the cholera typhoid co-infection model (Equation 1). The basic reproduction number, ${ℛ}_0$, obtained from the Table 1 is 1.4. The initial conditions used to produce the figures in this section were: S(0) = 99980, I_c(0) = 20, I_t(0) = 20, I_ct(0) = 20, R_c(0) = 0, R_t(0) = 0, R_ct(0) = 0, B_c(0) = 40000, B_t(0) = 40000. Note that all figures in this section are presented in the logarithmic scale since the range of some of the plots spanned several orders of magnitude.

TABLE 1

Table 1. Parameter values used for numerical simulation.

Coupled with the parameters from Table 1, the sensitivity indices of the variables above are shown on Figure 2. Latin Hypercube sampling was utilized to generate the plot above (Figure 2). This method returns the correlation between the state variable I_ct and each of the model parameters, and it also returns the ranks of all these correlations (PRCC). The simulation was carried out over 1,000 runs. A parameter with a negative PRCC value means that parameter is negatively correlated with I_ct, whilst a parameter with a positive PRCC value represents a positive correlation between that parameter and I_ct. Relative to the current model parameters, we note that the coinfection class is most sensitive to changes to the person-to-person typhoid transmission rate, β_t1, and the correlation is positive between this parameter and the state variable. The typhoid induced death rate is that second most sensitive parameter to the coinfection class, and it is negatively correlated to the coinfection class. Due to the large number of parameters in model 1, we have opted to split the PRCC values into 4 equal sets, see Figure 2.

FIGURE 2

Figure 2. The correlation between the co-infected class and each of the model's parameter are shown in this bar graph (PRCC). (A) Shows the PRCC values for {Λ, β_c1, β_t1, β_c2, β_t2, κ_c, κ_t, μ}. (B) Shows the PRCC values for {ρ_c, ρ_t, ρ_ct, δ_c, δ_t, ϵ_c, ϵ_t, ϵ_ct}. (C) Shows the PRCC values for {g_c, g_t, k_c, k_t, α_c, α_t, μ_c, μ_t}. (D) Shows the PRCC values for {β_t3, β_t4, β_c3, β_c4, θ_c, θ_t, η_c, η_t}.

The contour map of ${ℛ}_0$ as a function of the typhoid recovery rate, ϵ_t, and the cholera recovery rate, ϵ_c is shown in Figure 3. Using the parameters from Table 1, the base case as well as the contour levels are also shown in Figure 3. The basic reproduction number ${ℛ}_0$ attains its global minimum if both the typhoid and cholera recovery rate are maximized. It is be observed that, locally, a reduction in the reproduction number, ${ℛ}_0$, —moving the base case to a lower contour level—is only achieved by increasing the cholera recovery rate. Since the reproduction number, ${ℛ}_0$, is the maximum of the reproduction numbers of the individual diseases, it follows that a reduction in the reproduction number, ${ℛ}_0$, means a reduction in the reproduction numbers of each of the diseases. Hence, locally, an increase in the cholera recovery rate will not only reduce the cholera reproduction number, ${ℛ}_c$, but it has the added benefit of indirectly reducing the reproduction number for typhoid, ${ℛ}_t$, as well. It is also observed that increasing the typhoid recovery rate exclusively will have no immediate benefits locally. This finding is consistent with the previous findings of an optimal treatment plan being centered around the recovery rate of cholera.

FIGURE 3

Figure 3. The contour map of the basic reproduction number, ${ℛ}_0$, as a function of the typhoid recovery rate, ϵ_t, and the cholera recovery rate, ϵ_c.

We show the trajectories of the three infectious classes of model (Equation 1). An initial surge in infections followed closely by an immediate recovery is shown in Figure 4. The phenomenon of waning immunity results in the smaller second wave of infections. The co-infected class is the only exception to this observation. We see the co-infected class reach a local minimum before the first surge in cholera only or typhoid only infections is reached. A possible reason for this is that, unlike the cholera and typhoid classes, the co-infected class does not recruit directly from the susceptible class. This is due to the fact that the cholera disease has a shorter incubation period than the typhoid disease. The incubation periods are 1.4 days for cholera [2] and 19 days for typhoid [3]. What is then observed in the co-infected class is a case of people leaving the class either through death or recovery coupled with the delayed recruitment into the class. All the diseases reach stability after the second waves of infection.

FIGURE 4

Figure 4. The trajectories of the infectious classes.

In order to understand how the diseases interact with each other, we vary the different recovery rates and observe how the prevalence of each of the infections change. In Figures 5A,B show the impact of varying the recovery rate of the co-infected on the cholera and typhoid prevalence, whilst (Figures 5C,D) show the impact of varying the recovery rates of cholera and typhoid on the prevalence of the co-infected individuals. Figure 5C shows a significant reduction in the co-infected class' prevalence when the cholera recovery rate is increased, whilst plots (Figure 5D) shows that this reduction is negligible when the typhoid recovery rate was increased. Figures 5A,B show that an increase in the co-infected class' recovery rate reduces the typhoid prevalence more than the cholera prevalence. The net effect is that an increased cholera recovery rate may be associated with a decreased prevalence of the co-infected individuals and a higher co-infected recovery rate. This in turn, produces a reduced typhoid prevalence. Given the current model parameters, this finding suggests that an optimal treatment plan for the two infections should primarily focus on increasing the cholera recovery rate as opposed to the typhoid recovery rate. This also underscores the point made earlier about the symbiotic nature of the two diseases.

FIGURE 5

Figure 5. Plots (A,B) show the cholera and typhoid prevalence, respectively, as the co-infection recovery rate, ϵ_ct, runs through {0.1, 0.6, 1.1, 1.6}. Plots (C,D) show the prevalence of the co-infected as the cholera and typhoid recovery rates are varied through the sets {0.07, 0.075, 0.08, 0.085} and {0.1, 0.2, 0.3, 0.4}, respectively.

4. Discussion and conclusion

In this article, we formulated and analyzed a theoretical model for the transmission dynamics of a cholera typhoid co-infection model. Through numerical simulations, we were able to verify a number of the results obtained analytically.

The birth and death rates of the bacteria are central to proving the boundedness and positivity of all three models—cholera only sub-model, typhoid only sub-model, and the full cholera typhoid co-infection model. For the cholera-only model, if the birth rate of the Vibrio Cholerae bacteria exceeds its death rate, then the cholera only model has non-negative and bounded trajectories. For the typhoid only model, if the birth rate of the Salmonella Typhi bacteria exceeds its death rate, then the typhoid only model has non-negative and bounded trajectories. For the full cholera typhoid co-infection model, if the birth rates of the Vibrio Cholerae bacteria and the Salmonella Typhi bacteria exceed their death rates, simultaneously, then the cholera-typhoid co-infection model has non-negative and bounded trajectories.

In analyzing the equilibria of the co-infection models, several key sights were discovered. We showed the existence of the disease free equilibria, by finding them, for all three models. Sufficient conditions for the existence of the endemic equilibria for the cholera only sub-model and the typhoid only sub-model were documented. We showed that if the reproduction number is less than one for the all the models, then the disease free equilibria are locally asymptotically stable, otherwise they are unstable. Global stability could not be guaranteed, both at the disease free equilibria and the endemic equilibria, in any of the models. Sensitivity analysis revealed the parameters in the model were at the heart of the spread of the cholera typhoid co-infection. The prevalence of cholera is decreased whenever η_t, β_t3, β_t4 are increased and/or β_c1, β_c2, ϵ_ct, η_c and θ_c are decreased. The prevalence of typhoid is decreased whenever η_c, β_c3, β_c4 are increased and/or β_t1, β_t2, ϵ_ct, η_t and θ_t are decreased.

From the impact analysis section, we found that an increase in cholera cases may be associated with an increased risk of typhoid and that an increase in typhoid cases may be associated with an increased risk of cholera. This result proves the symbiotic nature of the relationship between the typhoid disease and the cholera disease.

The findings in this investigation come with some limitations. The most glaring of all is the lack of data to fit the model to. Our model also fails to take into account the highly seasonal nature of each of the diseases. For the two infections, fear has a significant impact on the transmission dynamics. Future work should also be able to account for the effects of fear in the transmission dynamics of both infections. Notwithstanding these limitations, we believe that the findings of this investigation can still be useful to policy makers in containing an outbreak of these two diseases.

Data availability statement

The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding author.

Author contributions

LM prepared the manuscript and contributed to the analysis, whilst FN formulated the model, proofread the manuscript, and also contributed to the analysis. Both authors contributed to the article and approved the submitted version.

Acknowledgments

The authors would like to acknowledge the financial support received from the University of Johannesburg during the production of this manuscript.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

1. Brachman PS, Evans AS. Bacterial Infections of Humans: Epidemiology Control. 3rd ed. Springer (1998). Available online at: http://gen.lib.rus.ec/book/index.php?md5=ae64651f2efa7ade60c537325be29cab

Google Scholar

2. Azman AS, Rudolph KE, Cummings DA, Lessler J. The incubation period of cholera: a systematic review. J Infect. (2013) 66:432–8. doi: 10.1016/j.jinf.2012.11.013

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Taylor A, Santiago A, Gonzalez-Cortes A, Gangarosa EJ. Outbreak of typhoid fever in Trinidad in 1971 traced to a commercial ice cream product. Am J Epidemiol. (1974) 100:150–7. doi: 10.1093/oxfordjournals.aje.a112017

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Akinyi OC, Mugisha J, Manyonge A, Ouma C, Maseno K. Modelling the impact of misdiagnosis and treatment on the dynamics of malaria concurrent and co-infection with pneumonia. Appl Math Sci. (2013) 7:6275–96. doi: 10.12988/ams.2013.39521

CrossRef Full Text | Google Scholar

5. Onyinge DO, Ongati NO, Odundo F. Mathematical model for co-infection of HIV/AIDS and pneumonia with treatment. Int J Scientific Eng Appl Sci. (2016) 2:106–11. Available online at: http://ir.jooust.ac.ke:8080/xmlui/handle/123456789/2643

Google Scholar

6. Mushayabasa S, Bhunu CP, Mhlanga NA. Modeling the transmission dynamics of typhoid in malaria endemic settings. Appl Appl Math. (2014) 9:121–40. Available online at: https://digitalcommons.pvamu.edu/aam/vol9/iss1/9

Google Scholar

7. Mushayabasa S, Bhunu CP, Ngarakana-Gwasira ET. Assessing the impact of drug resistance on the transmission dynamics of typhoid fever. Comput Biol J. (2013) 2013:1–13. doi: 10.1155/2013/303645

CrossRef Full Text | Google Scholar

8. Pitzer VE, Bowles CC, Baker S, Kang G, Balaji V, Farrar JJ, et al. Predicting the impact of vaccination on the transmission dynamics of typhoid in South Asia: A mathematical modeling study. PLoS Negl Trop Dis. (2014) 8:e2642. doi: 10.1371/journal.pntd.0002642

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Khan MA, Parvez M, Islam S, Khan I, Shafie S, Gul T. Mathematical analysis of typhoid model with saturated incidence rate. Adv Stud Biol. (2015) 7:65–78. doi: 10.12988/asb.2015.41059

CrossRef Full Text | Google Scholar

10. Matsebula L, Nyabadza F, Mushanyu J. Mathematical analysis of typhoid fever transmission dynamics with seasonality and fear. Commun Math Biol Neurosci. (2021) 2021:36. doi: 10.28919/cmbn/5590

CrossRef Full Text | Google Scholar

11. Codeço CT. Endemic and epidemic dynamics of cholera: the role of the aquatic reservoir. BMC Infect Dis. (2001) 1:1. doi: 10.1186/1471-2334-1-1

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Bertuzzo E, Casagrandi R, Gatto M, Rodriguez-Iturbe I, Rinaldo A. On spatially explicit models of cholera epidemics. J R Soc Interface. (2010) 7:321–33. doi: 10.1098/rsif.2009.0204

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Mari L, Bertuzzo E, Righetto L, Casagrandi R, Gatto M, Rodriguez-Iturbe I, Rinaldo A. Modelling cholera epidemics: the role of waterways, human mobility and sanitation. J R Soc Interface. (2012) 9:376–88. doi: 10.1098/rsif.2011.0304

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Miller Neilan RL, Schaefer E, Gaff H, Fister KR, Lenhart S. Modeling optimal intervention strategies for cholera. Bull Math Biol. (2010) 72:2004–18. doi: 10.1007/s11538-010-9521-8

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Mushayabasa S. Impact of vaccines on controlling typhoid fever in Kassena-Nankana district of upper east region of Ghana: insights from a mathematical model. J Modern Math Stat. (2011) 5:54–9. doi: 10.3923/jmmstat.2011.54.59

CrossRef Full Text

16. Sepulveda J, Gomez-Dantes H, Bronfman M. Cholera in the Americas: an overview. Infection. (1992) 20:243–8. doi: 10.1007/BF01710787

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Shuai Z, Tien JH, van den Driessche P. Cholera models with hyperinfectivity and temporary immunity. Bull Math Biol. (2012) 74:2423–45. doi: 10.1007/s11538-012-9759-4

PubMed Abstract | CrossRef Full Text | Google Scholar

18. King AA, Ionides EL, Pascual M, Bouma MJ. Inapparent infections and cholera dynamics. Nature. (2008) 454:877–80. doi: 10.1038/nature07084

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Sanches RP, Ferreira CP, Kraenkel RA. The role of immunity and seasonality in cholera epidemics. Bull Math Biol. (2011) 73:2916–31. doi: 10.1007/s11538-011-9652-6

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Okosun K, Makinde OD. Modelling the impact of drug resistance in malaria transmission and its optimal control analysis. Int J Phys Sci. (2011) 6:6479–87. doi: 10.5897/IJPS10.542

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Mutua JM, Wang FB, Vaidya NK. Modeling malaria and typhoid fever co-infection dynamics. Math Biosci. (2015) 264:128–44. doi: 10.1016/j.mbs.2015.03.014

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Hartley DM, Morris Jr JG, Smith DL. Hyperinfectivity: a critical element in the ability of V. cholerae to cause epidemics? PLoS Med. (2006) 3:e7. doi: 10.1371/journal.pmed.0030007

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Munro PM, Colwell RR. Fate of Vibrio cholerae O1 in seawater microcosms. Water Res. (1996) 30:47–50. doi: 10.1016/0043-1354(95)00137-A

CrossRef Full Text | Google Scholar

24. Blayneh K, Cao Y, Kwon HD. Optimal control of vector-borne diseases: treatment and prevention. Discrete Continuous Dyn Syst B. (2009) 11:587. doi: 10.3934/dcdsb.2009.11.587

CrossRef Full Text | Google Scholar

25. Hendrix TR. The pathophysiology of cholera. Bull N Y Acad Med. (1971) 47:1169–80.

Google Scholar

26. Adetunde IA. Mathematical models for the dynamics of typhoid fever in Kassena-Nankana district of upper east region of Ghana. J Modern Math Stat. (2008) 2:45–9. Available online at: https://medwelljournals.com/abstract/?doi=jmmstat.2008.45.49

27. Mushayabasa S. Modeling the impact of optimal screening on typhoid dynamics. Int J Dyn Control. (2016) 4:330–8. doi: 10.1007/s40435-014-0123-4

CrossRef Full Text | Google Scholar