A discrete-time risk-structured model of cholera infections in Cameroon

Abstract

In a recent paper, Che et al. [5] used a continuous-time Ordinary Differential Equation (ODE) model with risk structure to study cholera infections in Cameroon. However, the population and the reported cholera cases in Cameroon are censored at discrete-time annual intervals. In this paper, unlike in [5], we introduce a discrete-time risk-structured cholera model with no spatial structure. We use our discrete-time demographic equation to ‘fit’ the annual population of Cameroon. Furthermore, we use our fitted discrete-time model to capture the annually reported cholera cases from 1987 to 2004 and to study the impact of vaccination, treatment and improved sanitation on the number of cholera infections from 2004 to 2019. Our discrete-time cholera model confirms the results of the ODE model in [5]. However, our discrete-time model predicts a decrease in the number of cholera cases in a shorter period of cholera intervention (2004–2019) as compared to the ODE model's period of intervention (2004–2022).

Keywords:

• Cholera
• risk-structure model
• cholera intervention
• discrete-time
• Cameroon

1. Introduction

According to the World Health Organization (WHO), cholera is an acute intestinal infection caused by ingestion of food or water contaminated with the bacterium Vibrio cholerae [26]. In many parts of the world, cholera remains a significant threat to public health [25]. It continues to devastate impoverished populations with limited access to medication, clean water and proper sanitation amenities [16]. Annually, about 1.3 to 4.0 million cholera cases and 21, 000 to 143, 000 cholera-induced deaths are reported worldwide [1,26]. In recent years, there have been several reports of major cholera outbreaks. For example, in Haiti from 2010 to 2012 there were more than 530,000 reported cases [16]. Major outbreaks were also reported in Sierra Leone (2012), Ghana (2011), Nigeria (2010), Vietnam (2009), Zimbabwe (2008), and India (2007). Cameroon, a cholera endemic African country experienced large cholera outbreaks in 1991, 1996, 1998, 2004, 2010 and 2011 [19]. In 2017, 1,227,391 cases and 5,654 deaths of cholera were reported from34 countries [24].

Transmission of cholera can be indirect, from water infected with the bacterium to susceptible humans, or direct, from infectious humans to susceptible humans [18,20,23]. Several continuous-time mathematical ODE models have been developed, in an effort to gain a deeper understanding of cholera transmission dynamics, see [23] and the citations therein. In a recent paper, Che et al. introduced a continuous-time ODE low–high risk-structured cholera model and used it to capture the annual reported cholera infections in Cameroon from 1987 to 2017 [5].

The cholera census data used in [5] were reported annually [27]. Thus, unlike the ODE cholera model in [5], to capture the annual census data in Cameroon, in this paper, we introduce a discrete-time low-high risk-structured cholera model of Cameroon. That is, as in the reported annual cholera data of Cameroon, we assume that the unit time interval is 1 year and we use the knowledge of the population sizes at time t to predict the population sizes at time $(t+1)$. We will use our discrete-time model's demographic equation to ‘fit’ the population of Cameroon, and then use the fitted discrete-time cholera model to capture reported annual cholera cases in Cameroon from 1987 to 2004. Furthermore, as in [5], we will use our model to study the impact of three cholera intervention strategies (vaccination, treatment, and sanitation) on the number of cholera infections in Cameroon from 2004. Also, we will use our model to predict future cholera cases in Cameroon.

Both the continuous-time and discrete-time cholera models of Cameroon are based on observations and assumptions of the cholera disease dynamics in Cameroon. However, unlike the continuous-time model, the discrete-time model is also based on the observation that the cholera cases in Cameroon are reported annually. Thus, the difference equations of our novel discrete-time model are a closer fit to the observed cholera cases in Cameroon than the ordinary differential equations of the continuous-time choleramodel.

Discrete-time disease epidemic models are ubiquitous in the literature [4,10,13,21,22,28]. For example, van den Driessche and Yakubu, in [22], introduced a discrete-time SIR-B cholera model with direct and indirect cholera transmission pathways. Our discrete-time risk-structured cholera model, unlike the model in [22], partitions the entire Cameroon population into two cholera risk groups; ‘low’ and ‘high’ cholera risk groups.

The rest of the paper is organized as follows: In Section 2, we identify the low and high cholera risk groups in Cameroon and then introduce our discrete-time risk-structured model with no disease intervention. We compute the discrete-time demographic threshold parameter, ${\mathcal{R}}_d$, and then use our model's discrete-time demographic equation to ‘fit’ the census data of the population of Cameroon. We compute the basic reproduction number, ${\mathcal{R}}_0$, and then use the fitted discrete-time risk-structured model to capture reported cholera cases in Cameroon from 1987 to 2004. Also, in Section 2, we use sensitivity analysis to study the impact of our model parameters on ${\mathcal{R}}_d$, ${\mathcal{R}}_0$, and on the total number of our model's predicted cholera cases. In Section 3, we study the impact of vaccination, treatment and sanitation as cholera intervention strategies on the number of cholera infections in Cameroon from 2004 to 2019. We compute the effective reproduction number, ${\mathcal{R}}_E$, and use sensitivity analysis to study the impact of our model parameters on ${\mathcal{R}}_E$ and on the total number of predicted cholera cases by our model with intervention. In Section 4, we compare the results of our discrete-time cholera model of Cameroon to that of the corresponding continuous-time ODE Cholera model in [5]. Our conclusions are in Section 5.

2. Model formulation

In [5], Che et al. partitioned the ten regions of Cameroon, shown in Figure 1, into ‘low’ and ‘high’ risk cholera regions, see Table 1. Unlike in [5], we use the two risk groups to construct a discrete-time structured cholera model. Also, in [5], Che et al. used their scaled cholera model to capture scaled values of the 1987–2017 Cameroon's yearly reported cholera cases from WHO [27]. The scaled values are obtained by dividing the yearly reported cholera cases from [27] by the Cameroon population census data [17], see Figure 2. As in [5], we will use our discrete-time model to capture the actual cholera cases. First, we introduce the discrete-time demographic equation and use it to capture the 1987–2017 Cameroon's population census data from the World Bank [17].

Figure 1. Map of Cameroon showing the ten regions and neighbouring countries [15].

Figure 2. Scaled values of the yearly reported cases of cholera in Cameroon from 1987 to 2017. The scaled values are obtained by dividing the yearly reported cholera cases from [27] by the Cameroon population census data [17].

Table 1. Low- and high-risk regions of Cameroon.

Group	Regions in Cameroon
Low risk	Adamawa, South, North West, East
High risk	Far North, North, Littoral, South West, Centre, West

2.1. Demographic equation and Cameroon population

To introduce the discrete-time demographic equation, we let $N(t)$ denote the total population of Cameroon, at time $t=0,1,2,\dots$. The population of Cameroon is modelled by the difference equation (1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) where the parameters β, natural birth rate, and δ, natural death rate, are determined from the Cameroon's population census data [17]. In Model (1(1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) ), time is in years. The first term denotes the recruitment (birth) function of individuals into the population per year and the second time is the proportion of those that survive per year. We take t = 0 to correspond to 1987, t = 1 is 1988, and so on. Consequently, the initial condition $N(0)=10,730,000$ is the population of Cameroon in 1987, rounded to the nearest ten thousand (see [17]). Model (1(1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) ) has a dimensionless demographic threshold parameter ${\mathcal{R}}_d=\frac{\beta }{\delta }.$The total population goes extinct when ${\mathcal{R}}_d<1$ while it persists when ${\mathcal{R}}_d>1$.

In [5], Che et al. estimated that the birth rate per year, $\beta =4.1232\times {10}^{-2}$, and the death rate per year, $\delta =1.3301\times {10}^{-2}$. Consequently, as in [5], ${\mathcal{R}}_d=3.0999$. Hence, $\beta >\delta$ and ${\mathcal{R}}_d>1$ imply that the population of Cameroon is experiencing exponential growth.

Let $\eta =\beta +1-\delta .$Then $\eta =1.0279>1$ is the constant per capital growth rate.

A scatter plot of the Cameroon population census data in [17] and a plot of the solution $N(t)$ of (1(1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) ) using the estimates of β and δ are shown in Figure 3. From Figure 3, the solution $N(t)$ of Model (1(1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) ) fits the Cameroon population census data with relative error $=8.2652\times {10}^{-3}$, where as in [5] $\text{Relative error}=\frac{\Vert YP(\beta ,\delta )-Y{P}^{\ast }{\Vert }_2}{\Vert Y{P}^{\ast }{\Vert }_2},$$YP(\beta ,\delta )$ is the vector of the number of yearly population given by Model (1(1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) ) and $Y{P}^{\ast }$ is the vector of the actual number of the yearly population [17]. The scalar value $\Vert YP(\beta ,\delta )-Y{P}^{\ast }{\Vert }_2$ is the square root of the sum of the squares of the differences between the model results and the actual population data at the corresponding times, and $\Vert Y{P}^{\ast }{\Vert }_2$ is the square root of the sum of the squares of the actual population data for each year.

Figure 3. Scatterplot of Cameroon census data for years 1987–2017 [17] versus ‘fitted’ Model (1(1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) ) solution. The model solution $N(t)$ is shown as asterisks, and the dots are the data. $N(t)$ fits the Cameroon population data with relative error $=8.2652\times {10}^{-3}$.

2.2. Low–high-risk structured cholera model

Based on the schematic diagram shown in Figure 4, we introduce a simple risk-stratified discrete-time model of cholera infections in Cameroon with no explicit spatial structure. The list of the model variables and their descriptions are in Table 2, while the parameters are in Table 3.

Figure 4. Flow diagram for the risk-structured model. Solid lines represent flow between compartments, the straight dashed line represents the infectious class shedding pathogen into the environment and the curved dashed line represents the source of infections resulting from susceptibles interacting with the pathogen.

Table 2. Model variables, descriptions and units.

Variable	Description	Unit
$S_l$	Susceptible individuals with low risk	individuals
$E_l$	Exposed individuals with low risk	individuals
$I_l$	Infectious individuals with low risk	individuals
$R_l$	Recovered individuals with low risk	individuals
$B_l$	Pathogen concentration in water environment in low-risk population	cells ml${}^{-1}$
$N_l$	Subpopulation with low risk	individuals
$S_h$	Susceptible individuals with high risk	individuals
$E_h$	Exposed individuals with high risk	individuals
$I_h$	Infectious individuals with high risk	individuals
$R_h$	Recovered individuals with high risk	individuals
$B_h$	Pathogen concentration in water environment in high risk population	cells ml${}^{-1}$
$N_h$	Sub-population with high risk	individuals
N	Total population	individuals

Table 3. Model parameters, descriptions and values.

Parameter	Description	Value	Source
β	natural birth rate	$4.1232\times {10}^{-2}$ year${}^{-1}$	[5]
δ	non-disease related death rate	$1.3301\times {10}^{-2}$ year${}^{-1}$	[5]
${\theta }_l$	fraction of $S_l$ that become exposed via direct contact	$1.9998\times {10}^{-3}$	m. est.
${\theta }_h$	fraction of $S_h$ that become exposed via direct contact	$1.8009\times {10}^{-3}$	m. est.
${\rho }_l$	transmission rate between $S_l$ and $I_l$	0.2630	[22]
${\rho }_h$	transmission rate between $S_h$ and $I_h$	0.3072	m. est.
${\beta }_l$	contact rate between $S_l$ and $B_l$	0.5000	m. est.
${\beta }_h$	contact rate between $S_h$ and $B_h$	0.7900	m. est.
${\tau }_l$	fraction of $R_l$ that can become susceptible	0.1400	m. est.
${\tau }_h$	fraction of $R_h$ that can become susceptible	0.1000	m. est.
${\sigma }_l$	fraction of $E_l$ that can progress to infectious class	0.4100	m. est.
${\sigma }_h$	fraction of $E_h$ that can progress to infectious class	0.0033	m. est.
${\gamma }_l$	fraction of $I_l$ that can progress to recovered class	0.7000	m. est.
${\gamma }_h$	fraction of $I_h$ that can progress to recovered class	0.8000	m. est.
${\xi }_l$	shedding rate of pathogens by $I_l$	122 cells ml${}^{-1}$ year${}^{-1}$ individuals${}^{-1}$	m. est.
${\xi }_h$	shedding rate of pathogens by $I_h$	120 cells ml${}^{-1}$ year${}^{-1}$ individuals${}^{-1}$	m. est.
${\mu }_l$	fraction of $I_l$ that can die due to cholera infection	$7.2672\times {10}^{-3}$	[5]
${\mu }_h$	fraction of $I_h$ that can die due to cholera infection	0.2020	[5]
λ	net decay rate of pathogen	0.2160	m. est.
${\alpha }_{lh}$	fraction of $S_l$ that transition to $S_h$	0.1800	m. est.
${\alpha }_{hl}$	fraction of $S_h$ that transition to $S_l$	$1.2000\times {10}^{-4}$	m. est.
${\nu }_{lh}$	fraction of $E_l$ that transition to $E_h$	0.5800	m. est.
${\nu }_{hl}$	fraction of $E_h$ that transition to $E_l$	0.5700	m. est.
${\omega }_{lh}$	fraction of $R_l$ that transition to $R_h$	0.8524	[5]
${\omega }_{hl}$	fraction of $R_h$ that transition to $R_l$	0.7752	[5]

M. est. is an abbreviation for multistart estimation.

At time $t=0,1,2,\dots ,$ we denote by $S_l(t),E_l(t)$, $I_l(t)$, and $R_l(t)$ the population sizes of the susceptible, exposed, infectious, and recovered low-risk individuals, respectively. The total population size of the low-risk subpopulation at time t is $N_l(t)=S_l(t)+E_l(t)+I_l(t)+R_l(t).$Also, we denote by $S_h(t),E_h(t)$, $I_h(t)$, and $R_h(t)$ the population sizes of the susceptible, exposed, infectious, and recovered high-risk individuals at time $t=0,1,2,\dots ,$ respectively. Similarly, the total population size of the high-risk subpopulation at time $t=0,1,2,\dots ,$ is $N_h(t)=S_h(t)+E_h(t)+I_h(t)+R_h(t).$The total population of individuals in the two risk groups at time $t=0,1,2,\dots ,$ is $N(t)=N_l(t)+N_h(t).$At time $t=0,1,2,\dots ,$ we denote by $B_l(t)$ and $B_h(t)$ the concentration of pathogens in the contaminated water in the low-risk and high-risk environments, respectively.

In each risk group $k\in \{l,h\}$, we assume that a fraction ${\theta }_k\in (0,1)$ of susceptible individuals, become exposed via contact with infectious individuals (direct transmission) with probability ${\hat{\phi }}_{I_k}\left(\frac{I_k(t)}{N(t)}\right)=1-{\phi }_{I_k}\left(\frac{I_k(t)}{N(t)}\right),$and a fraction ${\hat{\theta }}_k=1-{\theta }_k\in (0,1)$ of susceptible individuals, become exposed via contact with contaminated water (indirect transmission) with probability ${\hat{\phi }}_{B_k}\left(\frac{B_k(t)}{N(t)}\right)=1-{\phi }_{B_k}\left(\frac{B_k(t)}{N(t)}\right).$The ‘escape from direct infection’ and ‘escape from indirect infection’ functions are, respectively, the nonlinear decreasing smooth concave-up functions ${\phi }_{I_k}:[0,\mathrm{\infty })\to [0,1]$and ${\phi }_{B_k}:[0,\mathrm{\infty })\to [0,1],$where ${\phi }_{I_k}(0)={\phi }_{B_k}(0)=1$. That is, for all $m\in \{I_k,B_k\}$, ${\phi }_m^{\prime }<0$ and ${\phi }_m^{″}\ge 0.$

If, for example, the interactions between susceptible individuals and infectious individuals or susceptible individuals and the contaminated water are modelled as Poisson processes, then (2) ${\phi }_{I_k}\left(\frac{I_k(t)}{N(t)}\right)=e^{-{\rho }_k\left(\frac{I_k(t)}{N(t)}\right)}$(2) and (3) ${\phi }_{B_k}\left(\frac{B_k(t)}{N(t)}\right)=e^{-{\beta }_k\left(\frac{B_k(t)}{N(t)}\right)},$(3) where the parameters ${\rho }_k$ and ${\beta }_k$ are the transmission rates between susceptible individuals and infectious individuals, and susceptible individuals and contaminated water, respectively (see Table 3).

For each unit time interval $[t,t+1]$, and in each risk group $k\in \{l,h\}$, the fraction of exposed individuals that can progress to the infectious class is ${\sigma }_k\in (0,1)$ and the fraction that remains in the exposed class is $(1-{\sigma }_k)$. Infectious individuals can contaminate the water by shedding the pathogen at a positive rate ${\xi }_k$ per unit time interval $[t,t+1]$. We assume that in Cameroon, the number of cholera infections from cholera infected corpses is very small and can be ignored. Cholera is a treatable disease. In Cameroon cholera is usually treated with oral rehydration salt solutions and antibiotics. In our model, for each unit time interval $[t,t+1]$, the fraction of infectious individuals that can progress to the recovered class is ${\gamma }_k\in (0,1)$ and the fraction that remains in the infectious class is $(1-{\gamma }_k)$. The fraction of infectious individuals that can die from cholera infection is ${\mu }_k\in (0,1)$. We assume that $\delta <1-{\mu }_k$, so that per unit time interval $[t,t+1]$, the fraction of infectious individuals that die from other causes does not exceed the fraction that do not die from cholera infection. Recovered individuals do not have permanent immunity, and so the fraction that can become susceptible to the infection is ${\tau }_k\in [0,1)$. In the pathogen class, a fraction $g\in [0,1)$ of cholera pathogen is added (growth) and a fraction $d\in (0,1)$ is removed from the water (decay).

Our model allows for the movement of individuals between the two risk groups. For each unit time interval $[t,t+1]$, the fraction of susceptible individuals that can move from the low-risk group to the high-risk group is ${\alpha }_{lh}\in (0,1)$, and the fraction of susceptible individuals that can move from the high-risk group to the low-risk group is ${\alpha }_{hl}\in (0,1)$. The fraction of exposed individuals that can move from the low-risk group to the high-risk group is ${\nu }_{lh}\in (0,1)$, and the fraction of exposed individuals that can move from the high-risk group to the low-risk group is ${\nu }_{hl}\in (0,1)$. We assume that ${\nu }_{lh}<1-{\sigma }_l$, so that the fraction of exposed individuals that can move from the low-risk group to the high-risk group does not exceed the fraction that does not become infectious in the low-risk group. Also, we assume that ${\nu }_{hl}<1-{\sigma }_h$, so that the fraction of exposed individuals that can move from the high-risk group to the low-risk group does not exceed the fraction that does not become infectious in the high-risk group. The fraction of recovered individuals that can move from the low-risk group to the high-risk group is ${\omega }_{lh}\in (0,1)$, and the fraction of exposed individuals that can move from the high-risk group to the low-risk group is ${\omega }_{hl}\in (0,1)$. We assume that ${\omega }_{lh}<1-{\tau }_l$, so that the fraction of recovered individuals that can move from the low-risk group to the high-risk group does not exceed the fraction that does not lose their immunity in the low-risk group. We assume that ${\omega }_{hl}<1-{\tau }_h$, so that the fraction of recovered individuals that can move from the high-risk group to the low-risk group does not exceed the fraction that does not lose their immunity in the high-risk group. However, we assume that infectious individuals are too sick to move and are confined to their respective risk groups. The simple mathematical model for the spread of the disease in the two risk groups with movements of non-infectious individuals between the groups is thus given by the following system of difference equations: (4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) with initial conditions: $S_l(0)>0,E_l(0)\ge 0,I_l(0)\ge 0,R_l(0)\ge 0,B_l(0)\ge 0,$and $S_h(0)>0,E_h(0)\ge 0,I_h(0)\ge 0,R_h(0)\ge 0,B_h(0)\ge 0.$In Figure 4, we have set $\begin{array}{rl}{a}_l& =(1-{\alpha }_{lh})(1-\delta ),c_l={\sigma }_l(1-\delta ),\\ q_l& ={\gamma }_l(1-\delta -{\mu }_l),p_l=(1-{\tau }_l-{\omega }_{lh})(1-\delta ),\\ f_l& ={\xi }_l(1-\delta -{\mu }_l),a_h=(1-{\alpha }_{hl})(1-\delta ),\\ c_h& ={\sigma }_h(1-\delta ),q_h={\gamma }_h(1-\delta -{\mu }_h),\\ p_h& =(1-{\tau }_h-{\omega }_{hl})(1-\delta ),f_h={\xi }_l(1-\delta -{\mu }_l),\\ a_{lh}& ={\alpha }_{lh}(1-\delta ),b_{lh}={\nu }_{lh}(1-\delta ),\end{array}$and $c_{lh}={\omega }_{lh}(1-\delta ),a_{hl}={\alpha }_{hl}(1-\delta ),b_{hl}={\nu }_{hl}(1-\delta ),c_{hl}={\omega }_{hl}(1-\delta ).$We note that all the terms in our structured Model (4(4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) ) are nonnegative. In fact, starting with initial conditions $S_l(0)>0,E_l(0)\ge 0,I_l(0)\ge 0,R_l(0)\ge 0,B_l(0)\ge 0,$and $S_h(0)>0,E_h(0)\ge 0,I_h(0)\ge 0,R_h(0)\ge 0,B_h(0)\ge 0,$we have $S_l(t)>0,E_l(t)\ge 0,I_l(t)\ge 0,R_l(t)\ge 0,B_l(t)\ge 0,$and $S_h(t)>0,E_h(t)\ge 0,I_h(t)\ge 0,R_h(t)\ge 0,B_h(t)\ge 0,$for all $t=0,1,2,\dots .$

For each risk group $k\in \{l,h\}$, in the absence of shedding of the Vibrio cholerae, the pathogen equations in Model (4(4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) ) reduce to the following difference equation $B_k(t+1)=(g+1-d)B_k(t).$Consequently, when d<g then the pathogen maintains itself in the environment via geometric growth. In the absence of shedding, to exclude pathogen population explosion in the Cameroon environment, we assume that d>g. That is, $\lambda =g+1-d<1$ is the net decay rate of the pathogen. Summing the equations for $S_l,S_h,I_h,I_l,E_h,E_l,R_h$ and $R_l$, we obtain (5) $N(t+1)=(\beta +(1-\delta ))N(t)-({\mu }_l{I}_l(t)+{\mu }_h{I}_h(t)).$(5) In the absence of cholera infections, $I_l(t)=I_h(t)=0$, and (5(5) $N(t+1)=(\beta +(1-\delta ))N(t)-({\mu }_l{I}_l(t)+{\mu }_h{I}_h(t)).$(5) ) reduces to (1(1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) ), the demographic equation of Cameroon. Let ${\kappa }_l=\frac{\lambda }{{\xi }_l}\text{and}{\kappa }_h=\frac{\lambda }{{\xi }_h}.$To consider proportions,we let $\begin{array}{rl}{s}_l(t)& =\frac{S_l(t)}{N(t)},e_l(t)=\frac{E_l(t)}{N(t)},i_l(t)=\frac{I_l(t)}{N(t)},r_l(t)=\frac{R_l(t)}{N(t)},b_l(t)=\frac{{\kappa }_l{B}_l(t)}{N(t)},\\ s_h(t)& =\frac{S_h(t)}{N(t)},e_h(t)=\frac{E_h(t)}{N(t)},i_h(t)=\frac{I_h(t)}{N(t)},r_h(t)=\frac{R_h(t)}{N(t)},b_h(t)=\frac{{\kappa }_h{B}_h(t)}{N(t)},\end{array}$and $n_l(t)=\frac{N_l(t)}{N(t)},n_h(t)=\frac{N_h(t)}{N(t)},$for $t=0,1,2,\dots .$

Also, we let $\varphi (t)=\frac{N(t)}{N(t+1)}=\frac{1}{\eta -({\mu }_l{i}_l(t)+{\mu }_h{i}_h(t))},$for $t=0,1,2,\dots .$ Model (4(4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) ) becomes (6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) where $\begin{array}{rl}{n}_l(t)& =s_l(t)+i_l(t)+e_l(t)+r_l(t),\\ n_h(t)& =s_h(t)+i_h(t)+e_h(t)+r_h(t),\end{array}$and $n_l(t)+n_h(t)=1,$for $t=0,1,2,\dots$.

We can rewrite Equation (5(5) $N(t+1)=(\beta +(1-\delta ))N(t)-({\mu }_l{I}_l(t)+{\mu }_h{I}_h(t)).$(5) ) as $\begin{array}{rl}N(t+1)& =(\beta -\delta )N(t)+S_l(t)+E_l(t)+R_l(t)+S_h(t)+E_h(t)+R_h(t)\\ & +(1-{\mu }_l)I_l(t)+(1-{\mu }_h)I_h(t).\end{array}$In Cameroon, $\beta >\delta$, ${\mu }_l<\delta +{\mu }_l<1$ and ${\mu }_h<\delta +{\mu }_h<1$. So, $N(t+1)>0.$ Moreover, in Cameroon, the population is growing exponentially. So, $N(t+1)>N(t)>0,N(t+1)>N_l(t)>0\text{and}N(t+1)>N_h(t)>0,$for $t=0,1,2,\dots$. Thus, $0<\varphi (t)<1$. With $s_l(t)>0,e_l(t)\ge 0,i_l(t)\ge 0,r_l(t)\ge 0,b_l(t)\ge 0,$and $s_h(t)>0,e_h(t)\ge 0,i_h(t)\ge 0,r_h(t)\ge 0,b_h(t)\ge 0,$we see that in Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ), $s_l(t+1)>0,e_l(t+1)\ge 0,i_l(t+1)\ge 0,r_l(t+1)\ge 0,b_l(t+1)\ge 0,$and $s_h(t+1)>0,e_h(t+1)\ge 0,i_h(t+1)\ge 0,r_h(t+1)\ge 0,b_h(t+1)\ge 0,$for $t=0,1,2,\dots$.

Since $i_l(t),i_h(t)\in [0,1]$ and $\varphi (t)<1$, for $t=0,1,2,\dots ,$ the equation for $b_l(t+1)$ and $b_h(t+1)$ in Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) implies $b_l(t+1)<\lambda (1-\delta -{\mu }_l)+\lambda b_l(t)$and $b_h(t+1)<\lambda (1-\delta -{\mu }_h)+\lambda b_h(t).$To obtain the feasible region of solutions to Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ), we consider the first-order linear difference equations $b_l(t+1)=\lambda (1-\delta -{\mu }_l)+\lambda b_l(t)$and $b_h(t+1)=\lambda (1-\delta -{\mu }_h)+\lambda b_h(t)$Recall that $0<\lambda <1$. Let $b_{\mathrm{\infty }}^l=\frac{\lambda (1-\delta -{\mu }_l)}{1-\lambda }$and $b_{\mathrm{\infty }}^h=\frac{\lambda (1-\delta -{\mu }_h)}{1-\lambda }.$For each risk group $k\in \{l,h\}$, we let ${\mathbf{x}}_k=(s_k,e_k,i_k,r_k,b_k)\in {\mathbf{R}}_{+}^5.$Assuming $b_k(t)\le b_{\mathrm{\infty }}^k$ for all $t=0,1,2,\dots$, the feasible region of solutions of Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) is the compact set $\mathrm{\Omega }=\{({\mathbf{x}}_l,{\mathbf{x}}_h)\in {\mathbf{R}}_{+}^{10}|0\le b_l\le b_{\mathrm{\infty }}^l\text{and}0\le b_h\le b_{\mathrm{\infty }}^h\}.$

2.3. Structured model: DFE and ${\mathcal{R}}_0$

Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) has a trivial disease-free equilibrium point (DFE), $p_{0-}=(0,0,0,0,0,0,0,0,0,0)$and a non-trivial DFE, $\begin{array}{rl}{p}_{0+}& =(s_l^0,s_h^0,e_l^0,e_h^0,i_l^0,i_h^0,r_l^0,r_h^0,b_l^0,b_h^0)\\ & =(\frac{{\alpha }_{hl}}{{\alpha }_{lh}+{\alpha }_{hl}},\frac{{\alpha }_{lh}}{{\alpha }_{lh}+{\alpha }_{hl}},0,0,0,0,0,0,0,0).\end{array}$For Cameroon, ${\mathcal{R}}_d>1$ and the trivial equilibrium, $p_{0-}$, is unstable. Now, for $p_{0+}$ we use the next generation matrix (NGM) method to compute ${\mathcal{R}}_0$, the basic reproduction number of Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) [2,22]. Assuming that $\lambda (1-\delta -{\mu }_l)i_l(t)\varphi (t)$ and $\lambda (1-\delta -{\mu }_h)i_h(t)\varphi (t)$ are not new infections, we obtain that F, the matrix of new infections, is $\begin{array}{r}F=\left(\begin{array}{cccccc}0& F_{12}& F_{13}& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& F_{45}& F_{46}\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\end{array}\right),\end{array}$where $\begin{array}{rl}{F}_{12}& =\frac{-(1-{\alpha }_{lh})(1-\delta ){\theta }_l{\phi }_{I_l}^{\prime }(0)s_l^0}{\eta }\\ F_{13}& =\frac{-(1-{\alpha }_{lh})(1-\delta ){\hat{\theta }}_l{\phi }_{B_l}^{\prime }(0)s_l^0}{{\kappa }_l\eta }\\ F_{45}& =\frac{-(1-{\alpha }_{hl})(1-\delta ){\theta }_h{\phi }_{I_h}^{\prime }(0)s_h^0}{\eta }\\ F_{46}& =\frac{-(1-{\alpha }_{hl})(1-\delta ){\hat{\theta }}_h{\phi }_{B_h}^{\prime }(0)s_h^0}{{\kappa }_h\eta },\end{array}$and T, the matrix of transition terms, is $\begin{array}{rl}T& =(\begin{array}{ccc}{\displaystyle \frac{(1-{\sigma }_l-{\nu }_{lh})(1-\delta )}{\eta }}& 0& 0\\ {\displaystyle \frac{{\sigma }_l(1-\delta )}{\eta }}& {\displaystyle \frac{(1-{\gamma }_l)(1-\delta -{\mu }_l)}{\eta }}& 0\\ 0& {\displaystyle \frac{\lambda (1-\delta -{\mu }_l)}{\eta }}& {\displaystyle \frac{\lambda }{\eta }}\\ {\displaystyle \frac{{\nu }_{lh}(1-\delta )}{\eta }}& 0& 0\\ 0& 0& 0\\ 0& 0& 0\end{array}\\ & \begin{array}{ccc}{\displaystyle \frac{{\nu }_{hl}(1-\delta )}{\eta }}& 0& 0\\ 0& 0& 0\\ 0& 0& 0\\ {\displaystyle \frac{(1-{\sigma }_h-{\nu }_{hl})(1-\delta )}{\eta }}& 0& 0\\ {\displaystyle \frac{{\sigma }_h(1-\delta )}{\eta }}& {\displaystyle \frac{(1-{\gamma }_h)(1-\delta -{\mu }_h)}{\eta }}& 0\\ 0& {\displaystyle \frac{\lambda (1-\delta -{\mu }_h)}{\eta }}& {\displaystyle \frac{\lambda }{\eta }}\end{array}).\end{array}$Thus, ${\mathcal{R}}_0$, the spectral radius of $F(Id-T{)}^{-1}$, is (7) ${\mathcal{R}}_0=\frac{Q_{11}+Q_{44}+\sqrt{(Q_{11}-Q_{44}{)}^2+4Q_{14}{Q}_{41}}}{2},$(7) where $\begin{array}{rl}{Q}_{11}& =(F_{12}+F_{13}\frac{\lambda (1-\delta -{\mu }_l)}{\eta -\lambda })\left(\frac{-{\sigma }_l\eta (1-\delta )[\eta -(1-{\sigma }_h-{\nu }_{hl})(1-\delta )]}{D[\eta -(1-{\gamma }_l)(1-\delta -{\mu }_l)]}\right),\\ Q_{14}& =(F_{12}+F_{13}\frac{\lambda (1-\delta -{\mu }_l)}{\eta -\lambda })\left(\frac{-{\sigma }_l{\nu }_{hl}\eta (1-\delta {)}^2}{D[\eta -(1-{\gamma }_l)(1-\delta -{\mu }_l)]}\right),\\ Q_{41}& =(F_{45}+F_{46}\frac{\lambda (1-\delta -{\mu }_h)}{\eta -\lambda })\left(\frac{-{\sigma }_h{\nu }_{lh}\eta (1-\delta {)}^2}{D[\eta -(1-{\gamma }_h)(1-\delta -{\mu }_h)]}\right),\\ Q_{44}& =(F_{45}+F_{46}\frac{\lambda (1-\delta -{\mu }_h)}{\eta -\lambda })\left(\frac{-{\sigma }_h\eta (1-\delta )[\eta -(1-{\sigma }_l-{\nu }_{lh})(1-\delta )]}{D[\eta -(1-{\gamma }_h)(1-\delta -{\mu }_h)]}\right),\\ D& ={\nu }_{hl}{\nu }_{lh}(1-\delta {)}^2-[\eta -(1-{\sigma }_l-{\nu }_{lh})(1-\delta )][\eta -(1-{\sigma }_h-{\nu }_{hl})(1-\delta )].\end{array}$In Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) ${\mathcal{R}}_0<1$ implies the DFE, $p_{0+}$, is locally asymptotically stable and cholera is eliminated in Cameroon. However, ${\mathcal{R}}_0>1$ implies $p_{0+}$ is unstable and cholera invades the Cameroon population.

In Cameroon, reported cholera interventions began in 2004 [9]. To compute ${\mathcal{R}}_0$ for Cameroon, in the next subsection, we use the scaled values of pre-intervention Cameroon's yearly reported cholera cases from 1987–2004 (see Figure 2) to estimate the remaining Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) parameter values. We will use the fitted Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) to illustrate that in Cameroon, ${\mathcal{R}}_0>1$ and cholera is an endemic infectious disease.

2.4. Structured model: parameter estimation

As in Section 2.1, we take $\beta =4.1232\times {10}^{-2}$ year${}^{-1}$ and $\delta =1.3301\times {10}^{-2}$ year${}^{-1}$. For the transmission rate between susceptible and infectious individuals in the low-risk group, ${\rho }_l=0.2630$ [22]. For the fraction of infectious individuals in the low-risk group who can die due to cholera infection, we multiply the cholera-related death rate per day in [5] by 365 days to obtain ${\mu }_l=(1.9910\times {10}^{-5})\times 365$ per year, that is ${\mu }_l=7.2672\times {10}^{-3}$. For the fraction of infectious individuals in the high-risk group who can die due to cholera infection, we multiply the cholera-related death rate per day in [5] by 365 days to obtain ${\omega }_{lh}=(5.5329\times {10}^{-4})\times 365$ per year, that ${\mu }_h=0.2020$. For the fraction of recovered individuals that can transition from the low-risk group to the high-risk group, we multiply the transition rate per day in [5] by 365 days to obtain ${\omega }_{lh}=(2.3353\times {10}^{-3})\times 365$ per year, that is ${\omega }_{lh}=0.8524$. For the fraction of recovered individuals that can transition from the high-risk group to the low-risk group, we multiply the transition rate per day in [5] by 365 days to obtain ${\omega }_{hl}=(2.1238\times {10}^{-3})\times 365$ per year, that is ${\omega }_{hl}=0.7752$.

As in [5], we estimate the remaining model parameters using the multistart algorithm with the fmincon function in the MATLAB optimization toolbox. The fmincon function takes an iterative approach to solving an optimization problem that includes our parameter values as unknowns with the scalar value to be minimized being the relative error, (8) $\text{Relative error}=\frac{\Vert YI(m)-Y{I}^{\ast }{\Vert }_2}{\Vert Y{I}^{\ast }{\Vert }_2},$(8) where $YI(m)$ is the vector of the number of yearly new infections given by the structured Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) with parameter vector m, and $Y{I}^{\ast }$ is the corresponding vector of the scaled values of the yearly reported cholera cases [27]. The scaled values are obtained by dividing the yearly cholera cases by the yearly population data [17].

There are more job, educational and trade opportunities in the high risk regions than the low-risk regions. Consequently, more people tend to move from the low-risk regions to the high-risk regions in search of those opportunities. To capture this movement pattern in the structured Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ), as in [5], we choose ${\alpha }_{hl}<{\alpha }_{lh},{\nu }_{hl}<{\nu }_{lh}$ and ${\omega }_{hl}<{\omega }_{lh}$, (see Table 3).

For the initial conditions, we take $N(0)=10,730,000$. As in [5], we assume that the high-risk population is about three times that of the low-risk population. Consequently, for the low-risk population, we let $N_l(0)=2,682,500$, $S_l(0)=2,682,476$, $E_l(0)=0$, $I_l(0)=24$ and $R_l(0)=0$. Recall that $s_l(0)=\frac{S_l(0)}{N(0)},e_l(0)=\frac{E_l(0)}{N(0)},i_l(0)=\frac{I_l(0)}{N(0)},r_l(0)=\frac{R_l(0)}{N(0)},n_l(0)=\frac{N_l(0)}{N(0)}.$For the high-risk population, we let $N_h(0)=8,047,500$, $S_h(0)=8,047,430$, $E_h(0)=0$, $I_h(0)=70$ and $R_h(0)=0$. Similarly, $s_h(0)=\frac{S_h(0)}{N(0)},e_h(0)=\frac{E_h(0)}{N(0)},i_h(0)=\frac{I_h(0)}{N(0)},r_h(0)=\frac{R_h(0)}{N(0)},n_h(0)=\frac{N_h(0)}{N(0)}.$As noted in [5], $75\mathrm{\%}$ percent of the total population of Cameroon in 1987 is high risk. To ensure that a higher fraction of individuals in the $S_h$ class are infected than the individuals in the $S_l$ class, as in [5], we choose ${\beta }_l$, ${\beta }_h$, ${\rho }_l$ and ${\rho }_h$ so that ${\beta }_l<{\beta }_h$ and ${\rho }_l<{\rho }_h$ (see Table 3).

We do not have data for the 1987 pathogen concentrations in the low and high-risk regions, $b_l(0)$ and $b_h(0)$, respectively. As in the continuous-time ODE model, we take $b_l(0)=1.3840\times {10}^{-6},b_h(0)=1.6915\times {10}^{-6}.$Also, we let ${\phi }_{I_l},{\phi }_{I_h},{\phi }_{B_l}$ and ${\phi }_{B_h}$ be defined as in (2(2) ${\phi }_{I_k}\left(\frac{I_k(t)}{N(t)}\right)=e^{-{\rho }_k\left(\frac{I_k(t)}{N(t)}\right)}$(2) ) and (3(3) ${\phi }_{B_k}\left(\frac{B_k(t)}{N(t)}\right)=e^{-{\beta }_k\left(\frac{B_k(t)}{N(t)}\right)},$(3) ). Using pre-intervention observations for Cameroon population, with $b_l(0)=1.3840\times {10}^{-6}$, $b_h(0)=1.6915\times {10}^{-6}$, $\beta =4.1232\times {10}^{-2}$ year${}^{-1}$ and $\delta =1.3301\times {10}^{-2}$ year${}^{-1}$, we estimate the remaining parameter values for ${\theta }_l,{\theta }_h,{\rho }_h,{\beta }_l,{\beta }_h,{\sigma }_l,{\sigma }_h,{\xi }_l,{\xi }_h,\lambda ,{\tau }_l,{\omega }_{lh},{\omega }_{hl}{\tau }_h,{\nu }_{lh}\text{and}{\nu }_{hl}$(see Table 3). Using these parameter values and Equation (7(7) ${\mathcal{R}}_0=\frac{Q_{11}+Q_{44}+\sqrt{(Q_{11}-Q_{44}{)}^2+4Q_{14}{Q}_{41}}}{2},$(7) ), we estimate that ${\mathcal{R}}_0=1.1832>1$. In [5], Che et al. used their ODE model to obtain ${\mathcal{R}}_0=1.1803>1$ for cholera disease in Cameroon. Thus, as in [5], our discrete-time risk-structured cholera model correctly predicts cholera endemicity in Cameroon.

In [5], we obtained that with a relative error of 0.8637, the continuous-time ODE structured model captures the known cholera endemicity in Cameroon, while predicting a decreasing trend from 1987 to 1994 and an increasing trend from 1995 to 2004 in the pre-intervention reported number of cholera cases in Cameroon from 1987 to 2004. Figure 5 shows the resulting Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) output values of the yearly number of new cholera infections predicted by the scaled structured model in comparison with the scaled values of the actual yearly cholera reported cases in Cameroon from 1987–2004, see Figure 2. From Figure 5, we see that unlike the continuous-time ODE model in [5], with a relative error of 0.8441, our discrete-time risk-structured model captures the known cholera endemicity in Cameroon, and predicts an increasing trend from 1987 to 1989, 1990 to 1992, 1993 to 1995, 1997 to 1998 and 2003 to 2004, and a decreasing trend from 1989 to 1990, 1992 to 1993, 1995 to 1997 and 1998 to 2003.

Figure 5. Scatterplot of the yearly cholera disease cases in Cameroon for 1987-2004 from [27] plotted against simulation output using estimated parameters. With a relative error of 0.8441, the structured model captures the known cholera endemicity in Cameroon, and predicts an increasing trend from 1987 to 1989, 1990 to 1992, 1993 to 1995, 1997 to 1998 and 2003 to 2004, and a decreasing trend from 1989 to 1990, 1992 to 1993, 1995 to 1997 and 1998 to 2003.

2.5. Sensitivity analysis

In this section, we use sensitivity analysis to study the impact of each model parameter on the threshold parameters, ${\mathcal{R}}_d$ and ${\mathcal{R}}_0$ and the number of model predicted cholera cases. For a variable v and a parameter p, the analytic sensitivity index to the parameter p is defined as ${\mathrm{{\rm Y}}}_p^v=\frac{\mathrm{\partial }v}{\mathrm{\partial }p}\times \frac{p}{v}.$This index gives the proportional rate of change of v as p changes [6]. A positive index will mean that an increase in the variable p will lead to a proportional increase in v, and a negative index will mean an increase in p will lead to a proportional decrease in v. We have ${\mathrm{{\rm Y}}}_{\beta }^{{\mathcal{R}}_d}=\frac{\mathrm{\partial }{\mathcal{R}}_d}{\mathrm{\partial }\beta }\times \frac{\beta }{{\mathcal{R}}_d}=1$and ${\mathrm{{\rm Y}}}_{\delta }^{{\mathcal{R}}_d}=\frac{\mathrm{\partial }{\mathcal{R}}_d}{\mathrm{\partial }\delta }\times \frac{\delta }{{\mathcal{R}}_d}=-1.$So, as in [5], the sensitivity indices of ${\mathcal{R}}_d$, the demographic threshold of Cameroon, do not depend on either β or δ, the parameters of (1(1) $N(t+1)=\beta N(t)+(1-\delta )N(t)=(\beta +(1-\delta ))N(t),$(1) ).

2.5.1. Sensitivity indices of ${\mathcal{R}}_0$

Next, we compute the sensitivity indices for ${\mathcal{R}}_0$ with respect to each of the 25 structured Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) parameters of Table 3. Most of the sensitivity indices are complicated algebraic expressions with no obvious structure. As in [6], we evaluate the sensitivity indices at the baseline parameter values given in Table 3. Table 4 summarizes our results with the parameters ordered from the most sensitive to the least sensitive. Unlike the ODE model in [5], where the most sensitive parameter is ${\gamma }_h$, the recovery rate of infectious individuals in the high-risk group, in our discrete-time model, ${\mathrm{{\rm Y}}}_{{\beta }_h}^{{\mathcal{R}}_0}=+9.7958\times {10}^{-1}$, and the most sensitive parameter is ${\beta }_h$, the contact rate between the susceptible individuals and the cholera pathogen in the high-risk group (see Table 4). Thus, decreasing (or increasing) the most sensitive parameter ${\beta }_h$, by $10\mathrm{\%}$ decreases (or increases) ${\mathcal{R}}_0$ by $9.7958\mathrm{\%}$. Also, from Table 4, ${\mathrm{{\rm Y}}}_{{\sigma }_h}^{{\mathcal{R}}_0}=+9.6839\times {10}^{-1}$. So, another important parameter in our discrete-time model is ${\sigma }_h$, the fraction of exposed individuals in the high-risk group that can progress to the infectious class. Thus, increasing (or decreasing) ${\sigma }_h$ by $10\mathrm{\%}$ increases (or decreases) ${\mathcal{R}}_0$ by $9.6839\mathrm{\%}$. The least sensitive parameters in our discrete-time model are ${\tau }_l$ and ${\tau }_h$, the fraction of recovered individuals in the high-risk group who become susceptible, and $w_{lh}$ and $w_{hl}$, the fraction of recovered individuals who transition from the low-risk group to the high-risk group and from the high-risk group to the low-risk groups, respectively. Because ${\mathrm{{\rm Y}}}_{{\tau }_l}^{{\mathcal{R}}_0}={\mathrm{{\rm Y}}}_{{\tau }_h}^{{\mathcal{R}}_0}={\mathrm{{\rm Y}}}_{w_{lh}}^{{\mathcal{R}}_0}={\mathrm{{\rm Y}}}_{w_{hl}}^{{\mathcal{R}}_0}=0,$decreasing (or increasing) ${\tau }_l$, ${\tau }_h$, $w_{lh}$ or $w_{hl}$ by $10\mathrm{\%}$ does not affect ${\mathcal{R}}_0$.

Table 4. Sensitivity indices of ${\mathcal{R}}_0$ with respect to the 25 structured Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) parameters, evaluated at the baseline parameter values given in Table 3.

Parameter	Sensitivity Index
${\beta }_h$	$+9.7958\times {10}^{-1}$
${\sigma }_h$	$+9.6839\times {10}^{-1}$
${\nu }_{hl}$	$-8.4628\times {10}^{-1}$
${\gamma }_h$	$-7.0608\times {10}^{-1}$
${\nu }_{lh}$	$+4.5360\times {10}^{-1}$
${\sigma }_l$	$-4.2758\times {10}^{-1}$
β	$-2.8572\times {10}^{-1}$
λ	$+2.6603\times {10}^{-1}$
${\xi }_l$	$+2.6603\times {10}^{-1}$
${\xi }_h$	$+2.6603\times {10}^{-1}$
${\mu }_h$	$-2.9798\times {10}^{-1}$
${\alpha }_{lh}$	$-2.4234\times {10}^{-2}$
${\alpha }_{hl}$	$+2.1187\times {10}^{-2}$
${\beta }_l$	$+2.0418\times {10}^{-2}$
${\gamma }_l$	$-1.9069\times {10}^{-2}$
${\theta }_h$	$-1.7614\times {10}^{-3}$
δ	$-1.1009\times {10}^{-3}$
${\mu }_l$	$-1.7313\times {10}^{-4}$
${\theta }_l$	$-4.0767\times {10}^{-5}$
${\rho }_h$	$+5.9253\times {10}^{-6}$
${\rho }_l$	$+1.4623\times {10}^{-7}$
${\tau }_l$	0
${\tau }_h$	0
${\omega }_{lh}$	0
${\omega }_{hl}$	0

The parameters are ordered from the most sensitive to the least sensitive. The most sensitive parameter is ${\beta }_h$, the contact rate between susceptible individuals and the cholera pathogen in the high-risk group, while the least sensitive parameters are ${\tau }_l$ and ${\tau }_h$, the fraction of recovered individuals in the high risk group who become susceptible, and $w_{lh}$ and $w_{hl}$, the fraction of recovered individuals who transition from the low-risk group to the high-risk group and from the high-risk group to the low-risk groups, respectively.

2.5.2. Global sensitivity analysis

Next, as in the ODE model in [5], we perform a global sensitivity analysis for our model using Latin Hypercube Sampling (LHS) to sample the parameter space and Partial Rank Correlation Coefficients (PRCC) to evaluate the sensitivity of the outcome variable, the total number of model predicted cholera cases, to uncertainty in the input variables [3,7,12]. For parameter intervals, as in [5], we pick values $50\mathrm{\%}$ above and below the values in Table 3. Uniform probability distributions were used for each parameter interval. From the recommendation in [14], we take $N>4M/3$ draws of the LHS design, where M, the number of input parameters is 25 and N, the number of LHS draws is 50.

Table 5 gives the PRCCs and the p-values for each model parameter while Figure 6 shows the sensitivity of the number of cholera cases to changes in parameters in Table 3 as computed by the Latin Hypercube Sampling-Partial Rank Correlation Coefficient (LHS-PRCC) index. Based on the p-values, we observe that PRCC values for ${\beta }_h$, 0.5440, for ${\sigma }_h$, 0.7091 and for ${\xi }_h$, 0.7054 are statistically significant. As compared to the ODE model in [5], where $w_h$, the scaled contact rate between susceptible individuals and the cholera pathogen in the high-risk group shows a significant positive correlation with the total number of cholera cases and ${\gamma }_h$, the recovery rate of infectious individuals in the high-risk group shows a significant negative correlation with the total number of cholera cases, in our discrete-time model, ${\beta }_h$, the contact rate between susceptible individuals and the cholera pathogen in the high risk, ${\sigma }_h$, the fraction of exposed individuals that can progress to the infectious class and ${\xi }_h$, the shedding rate of pathogens by infectious individuals in the high risk group show a significant positive correlation with the total number of cholera cases. This is in agreement with our discrete-time model predictions. In the next section, we explore the effects of the three cholera intervention strategies, vaccination, treatment and improved sanitation on the predictions of the fitted structured discrete-time model.

Figure 6. Sensitivity of the number of cholera cases to changes in parameters in Table 3 as computed by the Latin Hypercube Sampling-Partial Rank Correlation Coefficient (LHS-PRCC) index.

Table 5. Model parameters, corresponding PRCC and corresponding p-values resulting from the sensitivity analysis.

Parameter	PRCC	P-value
β	$--0.0994$	0.6291
δ	0.0373	0.8566
${\theta }_l$	$--0.2322$	0.2538
${\theta }_h$	0.1575	0.4422
${\rho }_l$	0.2770	0.1706
${\rho }_h$	0.1665	0.4161
${\beta }_l$	0.0068	0.9738
${\beta }_h$	0.5440	$\mathbf{0.0041}$
${\tau }_l$	$--0.0087$	0.9662
${\tau }_h$	0.3196	0.1115
${\sigma }_l$	$--0.4137$	0.0356
${\sigma }_h$	0.7091	$\mathbf{0.0001}$
${\gamma }_l$	0.0195	0.9248
${\gamma }_h$	$--0.3400$	0.0892
${\mu }_l$	$--0.2410$	0.2356
${\mu }_h$	$--0.4392$	0.0248
λ	0.4398	0.0246
${\xi }_l$	$--0.3397$	0.0895
${\xi }_h$	0.7054	$\mathbf{0.0001}$
${\alpha }_{lh}$	$--0.0741$	0.7192
${\alpha }_{hl}$	0.2122	0.2980
${\nu }_{lh}$	0.1954	0.3389
${\nu }_{hl}$	$--0.2057$	0.3135
${\omega }_{lh}$	0.0366	0.8592
${\omega }_{hl}$	0.3252	0.1051

Significant p-values $(p<0.01)$ are bold-faced.

3. Discrete-time risk-structured model with vaccination, treatment and sanitation

As noted in Che et al., cholera intervention strategies such as vaccination, treatment and improved water sanitation have been used to reduce the number of cholera infections. In Cameroon, treatment and water sanitation were first reported in 2004, while vaccination was first reported in 2015 [8,9,11]. In the next section, we introduce these three cholera intervention strategies in Model (4(4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) ).

3.1. Formulation of structured model with vaccination, treatment and sanitation

In Model (4(4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) ), we make the following assumptions in each unit time interval $[t,t+1]$.

1. The fractions of high- and low-risk susceptible populations that are vaccinated are $v_h\in [0,1)$ and $v_l\in [0,1)$, respectively.

2. The fractions of susceptible individuals that are vaccinated do not exceed the fractions that do not transition to another risk group. That is, $v_h<1-{\alpha }_{hl}$ and $v_l<1-{\alpha }_{lh}$.

3. The fractions of the high- and low-risk infectious populations that receive therapeutic treatment are $m_h\in [0,1)$ and $m_l\in [0,1)$, respectively.

4. The fractions of infectious individuals that receive therapeutic treatment do not exceed the fractions that do not recover from the cholera infection. That is, $m_l<1-{\gamma }_l$ and $m_h<1-{\gamma }_h$.

5. The fractions of vibrios in the high and low-risk groups that die due to improved water sanitation are ${\zeta }_h\in [0,1)$ and ${\zeta }_l\in [0,1)$, respectively.

6. The fractions of vibrios that die due to improved water sanitation do not exceed the fraction that is not cleared. That is, ${\zeta }_l<1-d$ and ${\zeta }_h<1-d$.

Our assumptions and notation lead to the following structured model with vaccination, treatment and sanitation: (9) $\{\begin{array}{ll}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\left(\frac{B_l(t)}{N(t)}\right)\\ & +{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t){\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right)\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t),\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& =({\gamma }_l+m_l)(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)\\ & +{\omega }_{hl}(1-\delta )R_h(t)+v_l(1-\delta )S_l(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d-{\zeta }_l)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\frac{I_h(t)}{N(t)}\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)\\ & +{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right)\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t),\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& =({\gamma }_h+m_h)(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)\\ & +{\omega }_{lh}(1-\delta )R_l(t)+v_h(1-\delta )S_h(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d-{\zeta }_h)B_h(t)\end{array}$(9) with initial conditions: $S_l(0)>0,E_l(0)\ge 0,I_l(0)\ge 0,R_l(0)\ge 0,B_l(0)\ge 0,$and $S_h(0)>0,E_h(0)\ge 0,I_h(0)\ge 0,R_h(0)\ge 0,B_h(0)\ge 0.$When there are no cholera interventions and $v_h=v_l=m_h=m_l={\zeta }_h={\zeta }_l=0$, then Model (9(9) $\{\begin{array}{ll}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\left(\frac{B_l(t)}{N(t)}\right)\\ & +{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t){\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right)\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t),\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& =({\gamma }_l+m_l)(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)\\ & +{\omega }_{hl}(1-\delta )R_h(t)+v_l(1-\delta )S_l(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d-{\zeta }_l)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\frac{I_h(t)}{N(t)}\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)\\ & +{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right)\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t),\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& =({\gamma }_h+m_h)(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)\\ & +{\omega }_{lh}(1-\delta )R_l(t)+v_h(1-\delta )S_h(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d-{\zeta }_h)B_h(t)\end{array}$(9) ) reduces to Model (4(4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) ). Note that Model (4(4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) ) and Model (9(9) $\{\begin{array}{ll}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\left(\frac{B_l(t)}{N(t)}\right)\\ & +{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t){\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right)\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t),\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& =({\gamma }_l+m_l)(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)\\ & +{\omega }_{hl}(1-\delta )R_h(t)+v_l(1-\delta )S_l(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d-{\zeta }_l)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\frac{I_h(t)}{N(t)}\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)\\ & +{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right)\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t),\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& =({\gamma }_h+m_h)(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)\\ & +{\omega }_{lh}(1-\delta )R_l(t)+v_h(1-\delta )S_h(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d-{\zeta }_h)B_h(t)\end{array}$(9) ) share the same equation for $N(t+1)$, Equation (5(5) $N(t+1)=(\beta +(1-\delta ))N(t)-({\mu }_l{I}_l(t)+{\mu }_h{I}_h(t)).$(5) ). As in the scaling of Model (4(4) $\{\begin{array}{rl}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh})(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+(1-{\alpha }_{lh})(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\\ & \times \left(\frac{B_l(t)}{N(t)}\right)+{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh})(1-\delta )S_l(t)({\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t)\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& ={\gamma }_l(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)+{\omega }_{hl}(1-\delta )R_h(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl})(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\left(\frac{I_h(t)}{N(t)}\right)+(1-{\alpha }_{hl})(1-\delta )\\ & \times S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)+{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl})(1-\delta )S_h(t)({\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)+{\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t)\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& ={\gamma }_h(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)+{\omega }_{lh}(1-\delta )R_l(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d)B_h(t)\end{array}$(4) ), we let $\begin{array}{rl}{\kappa }_l& =\frac{g+1-d}{{\xi }_l}=\frac{\lambda }{{\xi }_l},{\kappa }_h=\frac{g+1-d}{{\xi }_h}=\frac{\lambda }{{\xi }_l},\\ s_l(t)& =\frac{S_l(t)}{N(t)},e_l(t)=\frac{E_l(t)}{N(t)},i_l(t)=\frac{I_l(t)}{N(t)},r_l(t)=\frac{R_l(t)}{N(t)},b_l(t)=\frac{{\kappa }_l{B}_l(t)}{N(t)},\\ s_h(t)& =\frac{S_h(t)}{N(t)},e_h(t)=\frac{E_h(t)}{N(t)},i_h(t)=\frac{I_h(t)}{N(t)},r_h(t)=\frac{R_h(t)}{N(t)},b_h(t)=\frac{{\kappa }_h{B}_h(t)}{N(t)}\end{array}$and $n_l(t)=\frac{N_l(t)}{N(t)},n_h(t)=\frac{N_h(t)}{N(t)},$for $t=0,1,2,\dots .$

Also, we let $\varphi (t)=\frac{N(t)}{N(t+1)}=\frac{1}{\eta -({\mu }_l{i}_l(t)+{\mu }_h{i}_h(t))},$for $t=0,1,2,\dots .$ In proportions, Model (9(9) $\{\begin{array}{ll}{S}_l(t+1)& =\beta N_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t){\theta }_l{\phi }_{I_l}\left(\frac{I_l(t)}{N(t)}\right)\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t)\widehat{{\theta }_l}{\phi }_{B_l}\left(\frac{B_l(t)}{N(t)}\right)\\ & +{\alpha }_{hl}(1-\delta )S_h(t)+{\tau }_l(1-\delta )R_l(t)\\ E_l(t+1)& =(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t){\theta }_l{\hat{\phi }}_{I_l}\left(\frac{I_l(t)}{N(t)}\right)\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )S_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}\left(\frac{B_l(t)}{N(t)}\right)\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )E_l(t)+{\nu }_{hl}(1-\delta )E_h(t),\\ I_l(t+1)& ={\sigma }_l(1-\delta )E_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)I_l(t)\\ R_l(t+1)& =({\gamma }_l+m_l)(1-\delta -{\mu }_l)I_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})R_l(t)\\ & +{\omega }_{hl}(1-\delta )R_h(t)+v_l(1-\delta )S_l(t)\\ B_l(t+1)& ={\xi }_l(1-\delta -{\mu }_l)I_l(t)+g{B}_l(t)+(1-d-{\zeta }_l)B_l(t)\\ S_h(t+1)& =\beta N_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\theta }_h{\phi }_{I_h}\frac{I_h(t)}{N(t)}\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\hat{\theta }}_h{\phi }_{B_h}\left(\frac{B_h(t)}{N(t)}\right)\\ & +{\alpha }_{lh}(1-\delta )S_l(t)+{\tau }_h(1-\delta )R_h(t)\\ E_h(t+1)& =(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\theta }_h{\hat{\phi }}_{I_h}\left(\frac{I_h(t)}{N_h(t)}\right)\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )S_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}\left(\frac{B_h(t)}{N_h(t)}\right)\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )E_h(t)+{\nu }_{lh}(1-\delta )E_l(t),\\ I_h(t+1)& ={\sigma }_h(1-\delta )E_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)I_h(t)\\ R_h(t+1)& =({\gamma }_h+m_h)(1-\delta -{\mu }_h)I(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )R_h(t)\\ & +{\omega }_{lh}(1-\delta )R_l(t)+v_h(1-\delta )S_h(t)\\ B_h(t+1)& ={\xi }_h(1-\delta -{\mu }_h)I_h(t)+g{B}_h(t)+(1-d-{\zeta }_h)B_h(t)\end{array}$(9) ) becomes (10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) where $\begin{array}{rl}{n}_l(t)& =s_l(t)+i_l(t)+e_l(t)+r_l(t),\\ n_h(t)& =s_h(t)+i_h(t)+e_h(t)+r_h(t),\end{array}$and $n_l(t)+n_h(t)=1,$for $t=0,1,2,\dots .$

In the next section, we compute ${\mathcal{R}}_E$, the effective reproduction number for Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ).

3.2. Discrete-time risk-structured model with vaccination, treatment and sanitation: DFE and ${\mathcal{R}}_E$

Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) has a trivial DFE $p_{0-}^E=(0,0,0,0,0,0,0,0,0,0),$and a non-trivial DFE $\begin{array}{rl}{p}_{0+}^E& =(s_l^c,s_h^c,e_l^c,e_h^c,i_l^c,i_h^c,r_l^c,r_h^c,b_l^c,b_h^c)\\ & =(s_l^c,s_h^c,0,0,0,0,r_l^c,r_h^c,0,0),\end{array}$where $\begin{array}{rl}{s}_l^c& =\frac{{\alpha }_{hl}(1-\delta )[\beta +{\tau }_h(1-\delta )][\beta +(1-\delta )({\tau }_l+{\omega }_{lh})]}{\mathrm{\Delta }}\\ & +\frac{[\beta +{\tau }_l(1-\delta )](1-\delta {)}^2{\omega }_{hl}({\alpha }_{hl}+v_h)}{\mathrm{\Delta }},\\ s_h^c& =\frac{{\alpha }_{lh}(1-\delta )[\beta +{\tau }_l(1-\delta )][\beta +(1-\delta )({\tau }_h+{\omega }_{hl})]}{\mathrm{\Delta }}\\ & +\frac{[\beta +{\tau }_h(1-\delta )](1-\delta {)}^2{\omega }_{lh}({\alpha }_{lh}+v_l)}{\mathrm{\Delta }},\\ r_l^c& =\frac{(1-\delta {)}^3{v}_h{\omega }_{hl}({\alpha }_{lh}+v_l)+v_l{\alpha }_{hl}(1-\delta {)}^2(\beta +(1-\delta )({\tau }_h+{\omega }_{hl}))}{\mathrm{\Delta }},\\ r_h^c& =\frac{(1-\delta {)}^2{v}_h{\alpha }_{lh}(\beta +(1-\delta )({\tau }_l+{\omega }_{lh}))+(1-\delta {)}^3{v}_l{\omega }_{lh}({\alpha }_{hl}+v_h)}{\mathrm{\Delta }},\end{array}$and $\mathrm{\Delta }={\mathrm{\Delta }}_1+{\mathrm{\Delta }}_1,$where $\begin{array}{rl}{\mathrm{\Delta }}_1& =[\beta +{\tau }_l(1-\delta )][v_h{\omega }_{hl}(1-\delta {)}^2+({\alpha }_{hl}+{\alpha }_{lh})(1-\delta )(\beta +(1-\delta )({\tau }_h+{\omega }_{hl}))]\\ & +{\omega }_{lh}(1-\delta )[({\alpha }_{lh}+{\alpha }_{hl}+v_l)(1-\delta )(\beta +{\tau }_h(1-\delta ))+(1-\delta {)}^2{v}_l({\alpha }_{hl}+v_h)]\end{array}$and $\begin{array}{rl}{\mathrm{\Delta }}_2& =(1-\delta {)}^3({\alpha }_{lh}+v_l)v_h{\omega }_{hl}+(\beta +(1-\delta )({\tau }_h+{\omega }_{hl}))v_l{\alpha }_{hl}(1-\delta {)}^2\\ & +(\beta +(1-\delta )({\tau }_l+{\omega }_{hl}))v_h{\alpha }_{lh}(1-\delta {)}^2.\end{array}$For Cameroon, ${\mathcal{R}}_d>1$ and the trivial DFE, $p_{0-}^E$, is unstable. For the non-trivial DFE, $p_{0+}^E$, as in the previous section, assuming that $\lambda (1-\delta -{\mu }_l)i_l(t)\varphi (t)$ and $\lambda (1-\delta -{\mu }_h)i_h(t)\varphi (t)$ are not new infections, we obtain that F, the matrix of new infections, is $\begin{array}{r}F=\left(\begin{array}{cccccc}0& F_{12}& F_{13}& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& F_{45}& F_{46}\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\end{array}\right),\end{array}$where $\begin{array}{rl}{F}_{12}& =\frac{-(1-{\alpha }_{lh}-v_l)(1-\delta ){\theta }_l{\phi }_{I_l}^{\prime }(0)s_l^c}{\eta }\\ F_{13}& =\frac{-(1-{\alpha }_{lh}-v_l)(1-\delta ){\hat{\theta }}_l{\phi }_{B_l}^{\prime }(0)s_l^c}{{\kappa }_l\eta }\\ F_{45}& =\frac{-(1-{\alpha }_{hl}-v_h)(1-\delta ){\theta }_h{\phi }_{I_h}^{\prime }(0)s_h^c}{\eta }\\ F_{46}& =\frac{-(1-{\alpha }_{hl}-v_h)(1-\delta ){\hat{\theta }}_h{\phi }_{B_h}^{\prime }(0)s_h^c}{{\kappa }_h\eta },\end{array}$and T, the matrix of transition terms, is $\begin{array}{r}T=\left(\begin{array}{cccccc}{T}_{11}& 0& 0& {\displaystyle \frac{{\nu }_{hl}(1-\delta )}{\eta }}& 0& 0\\ {\displaystyle \frac{{\sigma }_l(1-\delta )}{\eta }}& T_{22}& 0& 0& 0& 0\\ 0& {\displaystyle \frac{\lambda (1-\delta -{\mu }_l)}{\eta }}& {\displaystyle \frac{\lambda -{\zeta }_l}{\eta }}& 0& 0& 0\\ {\displaystyle \frac{{\nu }_{lh}(1-\delta )}{\eta }}& 0& 0& T_{44}& 0& 0\\ 0& 0& 0& {\displaystyle \frac{{\sigma }_h(1-\delta )}{\eta }}& T_{55}& 0\\ 0& 0& 0& 0& {\displaystyle \frac{\lambda (1-\delta -{\mu }_h)}{\eta }}& {\displaystyle \frac{\lambda -{\zeta }_h}{\eta }}\end{array}\right),\end{array}$where $\begin{array}{rl}{T}_{11}& =\frac{(1-{\sigma }_l-{\nu }_{lh})(1-\delta )}{\eta }\\ T_{22}& =\frac{(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)}{\eta }\\ T_{44}& =\frac{(1-{\sigma }_h-{\nu }_{hl})(1-\delta )}{\eta }\\ T_{55}& =\frac{(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)}{\eta },\end{array}$Thus, ${\mathcal{R}}_E$, the spectral radius of $F(Id-T{)}^{-1}$, is (11) ${\mathcal{R}}_E=\frac{Q_{11}+Q_{44}+\sqrt{(Q_{11}-Q_{44}{)}^2+4Q_{14}{Q}_{41}}}{2},$(11) where $\begin{array}{rl}{Q}_{11}& =(F_{12}+F_{13}\frac{\lambda (1-\delta -{\mu }_l)}{\eta -\lambda +{\zeta }_l})\left(\frac{-{\sigma }_l\eta (1-\delta )[\eta -(1-{\sigma }_h-{\nu }_{hl})(1-\delta )]}{D[\eta -(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)]}\right),\\ Q_{14}& =(F_{12}+F_{13}\frac{\lambda (1-\delta -{\mu }_l)}{\eta -\lambda +{\zeta }_l})\left(\frac{-{\sigma }_l{\nu }_{hl}\eta (1-\delta {)}^2}{D[\eta -(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)]}\right),\\ Q_{41}& =(F_{45}+F_{46}\frac{\lambda (1-\delta -{\mu }_h)}{\eta -\lambda +{\zeta }_h})\left(\frac{-{\sigma }_h{\nu }_{lh}\eta (1-\delta {)}^2}{D[\eta -(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)]}\right),\\ Q_{44}& =(F_{45}+F_{46}\frac{\lambda (1-\delta -{\mu }_h)}{\eta -\lambda +{\zeta }_h})\left(\frac{-{\sigma }_h\eta (1-\delta )[\eta -(1-{\sigma }_l-{\nu }_{lh})(1-\delta )]}{D[\eta -(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)]}\right),\\ D& ={\nu }_{hl}{\nu }_{lh}(1-\delta {)}^2-[\eta -(1-{\sigma }_l-{\nu }_{lh})(1-\delta )][\eta -(1-{\sigma }_h-{\nu }_{hl})(1-\delta )].\end{array}$${\mathcal{R}}_E<1$ implies the DFE, $p_{0+}^E$, is locally asymptotically stable and our intervention strategies have successfully eliminated cholera in Cameroon. However, ${\mathcal{R}}_E>1$ implies the DFE, $p_{0+}^E$, is unstable, cholera invades and our intervention strategies have not succeeded in eliminating the disease in Cameroon population.

3.3. Numerical exploration

Pre-intervention Model (6(6) $\{\begin{array}{rl}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh})(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh})(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh})(1-\delta )s_l(t)({\theta }_l{\hat{\phi }}_{I_l}(i_l(t))+\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t)))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+\lambda b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl})(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl})(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl})(1-\delta )s_h(t)({\theta }_h{\hat{\phi }}_{I_h}(i_h(t))+{\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t)))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+\lambda b_h(t)]\varphi (t),\end{array}$(6) ) and Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) share the common pre-intervention parameter values of Table 3. As in [5], we use Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) to assess the impact of the seven cholera intervention strategies, described in Table 7, on the number of cholera infections in Cameroon. As in Section 2.4, using scaled values of the actual yearly reported cholera cases (as shown on Figure 2) and the multistart algorithm with the fmincon function, we estimated that $\begin{array}{rl}& v_l=0.6710,v_h=0.6891,m_l=0.0080,m_h=0.0850,\\ & {\zeta }_l=0.0030\text{and}{\zeta }_h=0.0070\end{array}$(see Table 6).

Table 6. Intervention parameters, their descriptions and values.

Parameter	Description	Value	Source
$v_l$	fraction of low-risk susceptible individuals that are vaccinated	0.6710	multistart est.
$v_h$	fraction of high-risk susceptible individuals that are vaccinated	0.6891	multistart est.
$m_l$	fraction of low-risk infectious individuals that are treated	0.0080	multistart est.
$m_h$	fraction of high-risk infectious individuals that are treated	0.0850	multistart est.
${\zeta }_l$	fraction of pathogen in low-risk group that dies due to sanitation	0.0030	multistart est.
${\zeta }_h$	fraction of pathogen in high-risk group that dies due to sanitation	0.0070	multistart est.

Multistart est. is the abbreviation for multistart estimation.

Table 7. Intervention strategies and their descriptions [5].

Intervention strategy	Description
1	Vaccination only
2	Treatment only
3	Sanitation only
4	1 and 2
5	2 and 3
6	1 and 3
7	1, 2 and 3

From Table 7, intervention strategies 1, 2 and 3 are single policies of cholera vaccination only, treatment only or improved sanitation only, respectively. However, intervention strategies 4, 5 and 6 are double policies of 1 and 2 only, or 2 and 3 only, or 1 and 3 only, respectively, while intervention strategy 7 is a combination of all the three single policies.

Intervention Strategy 1: The fraction of the low-risk susceptible population that is vaccinated is less than the fraction of high-risk susceptible population that is vaccinated.

Cholera vaccination in Cameroon was first reported in 2015. To illustrate the implementation of strategy 1 only from 2015 to 2019, we let $m_l=m_h={\zeta }_l={\zeta }_h=0,v_l=0.6710\text{and}{v}_h=0.6891,$in Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ). Figure 7(a) shows that Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) captures the decreasing effect of the vaccination only intervention strategy. Furthermore, Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) predicted that in 2019, with vaccination only, the percentage decrease in the number of cholera cases is $97.2621\mathrm{\%}$. With vaccination only, the effective reproduction number ${\mathcal{R}}_E=1.0028$.

Figure 7. Single Intervention Strategies. (a) Impact of vaccination. ${\mathcal{R}}_E=1.0028$, and the number of cholera infections is decreased by $97.2621\mathrm{\%}$. (b) Impact of treatment. ${\mathcal{R}}_E=1.1005$, and the number of cholera infections is decreased by $90.7952\mathrm{\%}$. (c) Impact of sanitation. ${\mathcal{R}}_E=1.1732$, and the number of cholera infections is decreased by $88.1657\mathrm{\%}$.

Intervention Strategy 2: The fraction of the low-risk infectious population that is treated is less than the fraction of high-risk infectious population that is treated.

Cholera treatment in Cameroon was first reported in 2004. To illustrate the implementation of strategy 2 only from 2004 to 2019, we let ${\zeta }_l={\zeta }_h=v_l=v_h=0,m_l=0.0080\text{and}{m}_h=0.0850,$in Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ). Figure 7(b) shows that Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) captures the decreasing effect of the treatment only intervention strategy. Furthermore, Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) predicted that in 2019, with treatment only, the percentage decrease in the number of cholera cases is $90.7952\mathrm{\%}$. With treatment only, the effective reproduction number ${\mathcal{R}}_E=1.1005$.

Intervention Strategy 3: The fraction of the pathogen in low risk group that dies due to improved water sanitation is less than the fraction of pathogen in high-risk group that dies due to improved water sanitation.

Improved water sanitation in Cameroon was first reported in 2004. To illustrate the implementation of strategy 3 only from 2004 to 2019, we let $m_l=m_h=v_l=v_h=0,{\zeta }_l=0.0030\text{and}{\zeta }_h=0.0070,$in Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ). Figure 7(c) shows that Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) captures the decreasing effect of the improved water sanitation only strategy. Furthermore, Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) predicted that in 2019, with improved water sanitation only, the percentage decrease in the number of cholera cases is $88.1657\mathrm{\%}$. With improvement in water sanitation only, the effective reproduction number ${\mathcal{R}}_E=1.1732$.

In Cameroon, with our choice of parameters, none of the single intervention policies reduces ${\mathcal{R}}_E$ to a value below 1. However, as in [5], vaccination only policy decreases the number of cholera cases most and has the lowest single policy ${\mathcal{R}}_E$ value. Next, we consider dual policies, and the combined policy of all three strategies.

Intervention Strategy 4: The fraction of susceptible individuals that are vaccinated and the fraction of infectious individuals that are treated in the low-risk group are less than the fraction of susceptible individuals that are vaccinated and the fraction of infectious individuals that are treated in the high-risk group, respectively.

To illustrate the implementation of strategy 4 from 2015 to 2019, we let $v_l=0.6710,v_h=0.6891,m_l=0.0080,m_h=0.0850$and ${\zeta }_l={\zeta }_h=0,$in Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ). Figure 8(a) shows that Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) captures the decreasing effect of the vaccination and treatment strategy. Furthermore, Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) predicted that in 2019, with vaccination and treatment, the percentage decrease in the number of cholera cases is $99.7119\mathrm{\%}$. With vaccination and treatment, the effective reproduction number ${\mathcal{R}}_E=0.9909$.

Figure 8. Double Intervention Strategies. (a) Impact of vaccination and treatment. ${\mathcal{R}}_E=0.9909$, and the number of cholera infections is decreased by $99.7119\mathrm{\%}$. (b) Impact of treatment and sanitation. ${\mathcal{R}}_E=1.0912$, and the number of cholera infections is decreased by $98.7505\mathrm{\%}$. (c) Impact of vaccination and sanitation. ${\mathcal{R}}_E=0.9989$, and the number of cholera infections is decreased by $97.5075\mathrm{\%}$.

Intervention Strategy 5: The fraction of infectious individuals that are treated and the fraction of pathogen that die due to sanitation in the low-risk group are less than the fraction of infectious individuals that are treated and the fraction of pathogen that die due to sanitation in the high-risk group, respectively.

To illustrate the implementation of strategy 5 from 2004 to 2019, we let $\begin{array}{rl}& v_l=v_h=0,m_l=0.0080,m_h=0.0850,\\ & {\zeta }_l=0.0030\text{and}{\zeta }_h=0.0070,\end{array}$in Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ). Figure 8(b) shows that Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) captures the decreasing effect of the treatment and sanitation strategy. Furthermore, Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) predicted that in 2019, with treatment and improved water sanitation, the percentage decrease in the number of cholera cases is $98.7505\mathrm{\%}$. With treatment and sanitation, the effective reproduction number ${\mathcal{R}}_E=1.0912$.

Intervention Strategy 6: The fraction of susceptible individuals that are vaccinated and the fraction of pathogen that die due to improved water sanitation in the low-risk group are less than the fraction of susceptible individuals that are vaccinated and the fraction of pathogen that die due to improved water sanitation in the high-risk group,respectively.

To illustrate the implementation of strategy 6 from 2015 to 2019, we let $\begin{array}{rl}& m_l=0=m_h,v_l=0.6710,v_h=0.6891,\\ & {\zeta }_l=0.0030\text{and}{\zeta }_h=0.0070,\end{array}$in Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ). Figure 8(c) shows that Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) captures the decreasing effect of the vaccination and sanitation strategy. Furthermore, Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) predicts that by 2019, with vaccination and improved water sanitation, the percentage decrease in the number of cholera cases is $97.5075\mathrm{\%}$. With vaccination and sanitation, the effective reproduction number ${\mathcal{R}}_E=0.9989$. As in [5], under the double policy, a vaccination and treatment policy can lead to a highest percentage decrease in the number of cholera cases and the smallest ${\mathcal{R}}_E$.

Intervention Strategy 7: The fraction of susceptible individuals that are vaccinated, the fraction of infectious individuals that are treated and the fraction of pathogen that dies due to improved water sanitation in the low-risk group are less than the fraction of susceptible individuals that are vaccinated, the fraction of infectious individuals that are treated and the fraction of pathogen that dies due to improved water sanitation in the high risk group, respectively.

To implement strategy 7 from 2015 to 2019, we let $\begin{array}{rl}& v_l=0.6710,v_h=0.6891,m_l=0.0080m_h=0.0850,\\ & {\zeta }_l=0.0030\text{and}{\zeta }_h=0.0070,\end{array}$in Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ). Figure 9 shows that Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) captures the decreasing effect of the combined strategy. Furthermore, Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) predicted that in 2019, with the combined strategy, the percentage decrease in the number of cholera cases is $99.9565\mathrm{\%}$. With this strategy, the effective reproduction number ${\mathcal{R}}_E=0.9872$. Table 8 summarizes estimates of the reproduction numbers and percentage decreases in cholera cases for the 7 different intervention strategies over the periods of intervention.

Figure 9. Impact of vaccination, treatment and sanitation. ${\mathcal{R}}_E=0.9872$, and the number of cholera infections is decreased by $99.9565\mathrm{\%}$.

Table 8. Reproduction numbers and percentage decreases in cholera cases for the seven different intervention strategies over the periods of intervention.

Intervention Strategy	Reproduction number	Decrease in cholera cases
Strategy 1 (2015–2019)	1.0028	$97.2621\mathrm{\%}$
Strategy 2 (2004–2019)	1.1005	$90.7952\mathrm{\%}$
Strategy 3 (2004–2019)	1.1732	$88.1657\mathrm{\%}$
Strategies 1 and 2 (2015–2019)	0.9909	$99.7119\mathrm{\%}$
Strategies 2 and 3 (2004–2019)	1.0912	$98.7505\mathrm{\%}$
Strategies 1 and 3 (2015–2019)	0.9989	$97.5075\mathrm{\%}$
Strategies 1, 2 and 3 (2015–2019)	0.9872	$99.9565\mathrm{\%}$

In Table 9, we summarize the corresponding ODE cholera model results from [5]. From Tables 8 and 9, we see that the ODE and discrete-time cholera intervention strategies results are consistent. However, the discrete-time model predicts a decrease in the number of cholera cases in a shorter period of cholera intervention (2015–2019, see Table 8) as compared to the ODE model's period of intervention (2015–2022, see Table 9).

Table 9. Reproduction numbers and percentage decreases in cholera cases for the seven different intervention strategies over the periods of intervention [5].

Intervention strategy	Reproduction number	Decrease in cholera cases
Strategy 1 (2015–2022)	1.0018	$98.6166\mathrm{\%}$
Strategy 2 (2004–2022)	1.1761	$90.5225\mathrm{\%}$
Strategy 3 (2004–2022)	1.1767	$86.5359\mathrm{\%}$
Strategies 1 and 2 (2015–2022)	0.9982	$99.6735\mathrm{\%}$
Strategies 2 and 3 (2004–2022)	1.1725	$98.6191\mathrm{\%}$
Strategies 1 and 3 (2015–2022)	0.9987	$98.7498\mathrm{\%}$
Strategies 1, 2 and 3 (2015–2022)	0.9952	$99.7280\mathrm{\%}$

3.4. Sensitivity indices of ${\mathcal{R}}_E$

As before, we compute the sensitivity indices of ${\mathcal{R}}_E$ with respect to each of the 31 Model (10(10) $\{\begin{array}{ll}{s}_l(t+1)& =[\beta n_l(t)+(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\phi }_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\phi }_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +{\alpha }_{hl}(1-\delta )s_h(t)+{\tau }_l(1-\delta )r_l(t)]\varphi (t)\\ e_l(t+1)& =[(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t){\theta }_l{\hat{\phi }}_{I_l}(i_l(t))\\ & +(1-{\alpha }_{lh}-v_l)(1-\delta )s_l(t)\widehat{{\theta }_l}{\hat{\phi }}_{B_l}(\frac{1}{{\kappa }_l}{b}_l(t))\\ & +(1-{\sigma }_l-{\nu }_{lh})(1-\delta )e_l(t)+{\nu }_{hl}(1-\delta )e_h(t)]\varphi (t)\\ i_l(t+1)& =[{\sigma }_l(1-\delta )e_l(t)+(1-{\gamma }_l-m_l)(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ r_l(t+1)& =[{\gamma }_l(1-\delta -{\mu }_l)i_l(t)+(1-\delta )(1-{\tau }_l-{\omega }_{lh})r_l(t)\\ & +{\omega }_{hl}(1-\delta )r_h(t)+v_l(1-\delta )s_l(t)\\ & +m_l(1-\delta -{\mu }_l)i_l(t)]\varphi (t)\\ b_l(t+1)& =[\lambda (1-\delta -{\mu }_l)i_l(t)+(\lambda -{\zeta }_l)b_l(t)]\varphi (t)\\ s_h(t+1)& =[\beta n_h(t)+(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\phi }_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\phi }_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +{\alpha }_{lh}(1-\delta )s_l(t)+{\tau }_h(1-\delta )r_h(t)]\varphi (t)\\ e_h(t+1)& =[(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\theta }_h{\hat{\phi }}_{I_h}(i_h(t))\\ & +(1-{\alpha }_{hl}-v_h)(1-\delta )s_h(t){\hat{\theta }}_h{\hat{\phi }}_{B_h}(\frac{1}{{\kappa }_h}{b}_h(t))\\ & +(1-{\sigma }_h-{\nu }_{hl})(1-\delta )e_h(t)+{\nu }_{lh}(1-\delta )e_l(t)]\varphi (t)\\ i_h(t+1)& =[{\sigma }_h(1-\delta )e_h(t)+(1-{\gamma }_h-m_h)(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ r_h(t+1)& =[{\gamma }_h(1-\delta -{\mu }_h)i_h(t)+(1-{\tau }_h-{\omega }_{hl})(1-\delta )r_h(t)\\ & +{\omega }_{lh}(1-\delta )r_l(t)+v_h(1-\delta )s_h(t)\\ & +m_h(1-\delta -{\mu }_h)i_h(t)]\varphi (t)\\ b_h(t+1)& =[\lambda (1-\delta -{\mu }_h)i_h(t)+(\lambda -{\zeta }_h)b_h(t)]\varphi (t),\end{array}$(10) ) parameters, using the baseline parameter values given in Table 3 and that of the estimated intervention parameters in Table 6 (see Table 10). From these results, unlike the ODE model in [5], where the most sensitive parameter with the intervention model is ${\gamma }_h$, the recovery rate of infectious individuals in the high-risk group, the most sensitive parameter with our discrete-time intervention model is $v_l$, the fraction of low-risk susceptible individuals that are vaccinated. From Table 10, ${\mathrm{{\rm Y}}}_{v_l}^{{\mathcal{R}}_E}=-5.0919$. Thus, decreasing (or increasing) $v_l$ by $10\mathrm{\%}$ increases (or decreases) ${\mathcal{R}}_E$ by $50.9190\mathrm{\%}$. Another important parameter is ${\alpha }_{lh}$, the fraction of susceptible individuals that transition from the low-risk group to the high-risk group. From Table 10, ${\mathrm{{\rm Y}}}_{{\alpha }_{lh}}^{{\mathcal{R}}_E}=-1.4049$. Thus, decreasing (or increasing) ${\alpha }_{lh}$ by $10\mathrm{\%}$ increases (or decreases) ${\mathcal{R}}_E$ by $14.0490\mathrm{\%}$. Unlike the ODE model in [5], where the least sensitive parameter is ${\zeta }_l$, the death rate of the cholera pathogen due to improved water sanitation in the low-risk group, the least sensitive parameter in our model with intervention is ${\rho }_h$, the transmission rate between susceptible individuals and infectious individuals in the high-risk group. From Table 10, ${\mathrm{{\rm Y}}}_{{\rho }_h}^{{\mathcal{R}}_E}=+1.2192\times {10}^{-7}$. Thus, decreasing (or increasing) ${\rho }_h$ by $10\mathrm{\%}$ decreases (or increases) ${\mathcal{R}}_E$ by $1.2192\times {10}^{-6}\mathrm{\%}$. Among the six intervention parameters, unlike in [5], where the most sensitive parameter is $v_h$, the vaccination rate of susceptible individuals in the high-risk group, and the least sensitive parameter is ${\zeta }_l$, the death rate of the cholera pathogen due to improved water sanitation in the low-risk group, the most sensitive intervention parameter in our discrete-time model is $v_l$, the fraction of susceptible individuals in the low-risk group that are vaccinated, and the least sensitive parameter is ${\zeta }_h$, the fraction of the pathogen in the high-risk group that dies due to improved water sanitation.

Table 10. Sensitivity indices of ${\mathcal{R}}_E$ with respect to the 31 parameters for the structured model with intervention, evaluated at the baseline parameter values in Table 3 and that of the estimated intervention parameters in Table 6.

Parameter	Sensitivity index
$v_l$	$-5.0919$
${\alpha }_{lh}$	$-1.4049$
${\tilde{\beta }}_l$	$+9.8001\times {10}^{-1}$
${\xi }_l$	$+9.8001\times {10}^{-1}$
${\gamma }_l$	$-9.0561\times {10}^{-1}$
${\tau }_l$	$+6.6261\times {10}^{-1}$
${\omega }_{hl}$	$+5.05536\times {10}^{-1}$
${\omega }_{lh}$	$-4.8373\times {10}^{-1}$
λ	$+2.6503\times {10}^{-1}$
${\sigma }_l$	$+1.5086\times {10}^{-1}$
β	$-9.0633\times {10}^{-2}$
${\nu }_{lh}$	$-6.6352\times {10}^{-2}$
${\nu }_{hl}$	$+6.2599\times {10}^{-2}$
${\tau }_h$	$-5.1289\times {10}^{-2}$
$v_h$	$+4.5563\times {10}^{-2}$
${\tilde{\beta }}_h$	$+1.9983\times {10}^{-2}$
${\xi }_h$	$+1.9983\times {10}^{-2}$
$m_l$	$-1.0350\times {10}^{-2}$
${\sigma }_h$	$+1.3939\times {10}^{-2}$
${\gamma }_h$	$-1.3379\times {10}^{-2}$
${\mu }_l$	$-1.0074\times {10}^{-2}$
${\mu }_h$	$-5.6375\times {10}^{-3}$
${\zeta }_l$	$-3.6077\times {10}^{-3}$
δ	$+2.0881\times {10}^{-3}$
${\theta }_l$	$-1.9567\times {10}^{-3}$
$m_h$	$-1.4215\times {10}^{-3}$
${\alpha }_{hl}$	$+1.9379\times {10}^{-4}$
${\zeta }_h$	$-1.7081\times {10}^{-4}$
${\theta }_h$	$-3.5930\times {10}^{-5}$
${\rho }_l$	$+7.0446\times {10}^{-6}$
${\rho }_h$	$+1.2192\times {10}^{-7}$

The parameters are ordered from the most sensitive to the least sensitive. The most sensitive parameter with the intervention model is $v_l$, the fraction of susceptible individuals in the low-risk group that are vaccinated, while the least sensitive parameter is ${\rho }_h$, the transmission rate between susceptible individuals and infectious individuals in the high-risk group.

3.5. Global sensitivity analysis

As before, we perform a global sensitivity analysis for our discrete-time model with an intervention using LHS and PRCC. Now, our number of input parameters is 31 and the number of LHS draws is 80. Table 11 gives the PRCCs and the p-values for each model parameter, while Figure 10 shows the sensitivity of the number of model predicted cholera cases to changes in parameters of Tables 11 and 6 as computed by the LHS-PRCC index. Based on the p-values, we observe that PRCC values for ${\beta }_l$ and ${\gamma }_h$ are statistically significant. In the ODE model in [5], the parameters $w_l$ (scaled contact rate between susceptible individuals and cholera pathogen in the high-risk group), ${\tau }_h$ (rate at which recovered individuals become susceptible individuals in the high-risk group) and ${\omega }_{lh}$ (transition rate of recovered individuals from low-risk group to high-risk group) show a significant positive correlation with the total number of model predicted cholera cases. Meanwhile, ${\gamma }_l$ and ${\gamma }_h$ (the rates at which infectious individuals become recovered individuals in the low and high-risk groups, respectively), λ (the net decay rate of the cholera pathogen) and ${\omega }_{hl}$ (transition rate of recovered individuals from high-risk group to low risk group) show a significant negative correlation with the total number of the ODE model predicted cholera cases. In our discrete-time model, ${\beta }_l$, the contact rate between susceptible individuals and the cholera pathogen in the low-risk group shows a significant positive correlation with the total number of our discrete-time model predicted cholera cases. On the other hand, λ, the net decay rate of pathogen shows a significant negative correlation with the total number of our discrete-time model predicted cholera cases. This is in agreement with our discrete-time model predictions. That is, increasing ${\beta }_l$ increases the total number of model predicted cholera cases while increasing λ decreases the total number of model predicted cholera cases.

Figure 10. Sensitivity of the number of cholera cases to changes in parameters in Tables 3 and 6 as computed by the Latin Hypercube Sampling-Partial Rank Correlation Coefficient (LHS-PRCC) index.

Table 11. Model parameters, corresponding PRCC and corresponding p-values resulting from the sensitivity analysis.

Parameter	PRCC	P-value
β	$--0.1002$	0.4885
δ	$--0.0199$	0.8909
${\theta }_l$	0.0445	0.7588
${\theta }_h$	$--0.1166$	0.4202
${\rho }_l$	0.0780	0.5903
${\rho }_h$	0.1002	0.4887
${\beta }_l$	0.4589	$\mathbf{0.0008}$
${\beta }_h$	0.0131	0.9278
${\tau }_l$	$--0.1545$	0.2840
${\tau }_h$	$--0.1802$	0.2105
${\sigma }_l$	0.2087	0.1458
${\sigma }_h$	0.1175	0.4166
${\gamma }_l$	0.2956	0.0372
${\gamma }_h$	0.1837	2017
${\mu }_l$	$--0.0816$	0.5731
${\mu }_h$	$--0.2007$	0.1621
λ	$--0.4813$	$\mathbf{0.0004}$
${\xi }_l$	0.2217	0.1218
${\xi }_h$	0.3155	0.0257
${\alpha }_{lh}$	$--0.2486$	0.0818
${\alpha }_{hl}$	$--0.2706$	0.0573
${\nu }_{lh}$	$--0.2831$	0.0464
${\nu }_{hl}$	0.2706	0.0574
${\omega }_{lh}$	0.1761	0.2213
${\omega }_{hl}$	0.2172	0.1297
$v_l$	$--0.0640$	0.6586
$v_h$	0.1406	0.3301
$m_l$	$--0.1576$	0.2745
$m_h$	$--0.1258$	0.3842
${\zeta }_l$	$--0.2302$	0.1078
${\zeta }_h$	$--0.0431$	0.7664

Significant p-values $(p<0.01)$ are bold-faced.

4. Discrete-time model versus the continuous-time ODE model

Both our discrete-time and continuous-time models are Cameroon data-driven parametrized models. However, unlike the continuous-time model, the discrete-time model is based on the discrete-time censoring of the population and cholera data of Cameroon. From Figure 5, we see that unlike the continuous-time ODE model in [5] that captures the known cholera endemicity in Cameroon with a relative error of 0.8637, with a smaller relative error of 0.8441, our discrete-time risk-structured model also captures the known cholera endemicity in Cameroon. Moreover, unlike the continuous-time ODE model in [5] that predicted a decreasing trend from 1987 to 1994 and an increasing trend from 1995 to 2004 in the pre-intervention reported number of cholera cases in Cameroon from 1987 to 2004, our discrete-time model predicts the observed oscillatory nature in the reported cholera cases in Cameroon.

From our sensitivity analysis, we obtain that, as in the ODE model in [5], in our discrete-time model ${\mathcal{R}}_d$ is not sensitive to changes to the parameters of the demographic equation. For ${\mathcal{R}}_0$, unlike the ODE model in [5], where the most sensitive parameter is ${\gamma }_h$, the recovery rate of infectious individuals in the high-risk group, in our discrete-time model, the most sensitive parameter is ${\beta }_h$, the contact rate between the susceptible individuals and the cholera pathogen in the high-risk group. In [5], the least sensitive parameters are ${\tau }_l$ and ${\tau }_h$, the parameters associated with loss of immunity, and ${\omega }_{lh}$ and ${\omega }_{hl}$, transition rates of recovered individuals between the two risk groups. In our discrete-time model, the least sensitive parameters are ${\tau }_l$ and ${\tau }_h$, the fraction of recovered individuals in the high-risk group who become susceptible, and $w_{lh}$ and $w_{hl}$, the fraction of recovered individuals who transition from the low-risk group to the high-risk group and from the high-risk group to the low-risk groups, respectively. Another important parameter that has a high impact on ${\mathcal{R}}_0$ in our discrete-time model is ${\sigma }_h$, the fraction of exposed individuals in the high-risk group that can progress to the infectious class. In the ODE model in [5], $w_h$, the scaled contact rate between susceptible individuals and the cholera pathogen in the high-risk group shows a significant positive correlation with the total number of cholera cases and ${\gamma }_h$, the recovery rate of infectious individuals in the high-risk group shows a significant negative correlation with the total number of cholera cases. In our discrete-time model, ${\beta }_h$, the contact rate between susceptible individuals and the cholera pathogen in the high-risk, ${\sigma }_h$, the fraction of exposed individuals that can progress to the infectious class and ${\xi }_h$, the shedding rate of pathogens by infectious individuals in the high-risk group show a significant positive correlation with the total number of cholera cases.

Under the single intervention policies of vaccination only or treatment only or improved sanitation only, in Table 8, we obtained that, as in [5], a vaccination only policy has the smallest value of the effective reproduction number ${\mathcal{R}}_E$ and the highest percentage decrease in the number of cholera cases. As in [5], under the double intervention policies of vaccination and treatment only or treatment and improved sanitation only or vaccination and improved sanitation only, in Table 8 we obtained that the double policy of vaccination and treatment has the smallest value for ${\mathcal{R}}_E$ and the highest percentage decrease in the number of cholera cases. As in [5], each of the three strategies of either vaccination and treatment, or vaccination and improved sanitation, or the combined strategy of vaccination, treatment and improved sanitation is capable of eliminating cholera in Cameroon with the combined strategy having the lowest value for the effective reproduction number, ${\mathcal{R}}_E$, and the highest percentage decrease in the number of cholera cases. As in [5], the combined strategy of vaccination, treatment and improved sanitation lead to the smallest number of cholera infections in Cameroon. From Tables 8 and 9, our discrete-time model results confirm that of the ODE model. However, since the discrete-time model is built on the discrete-time census data it predicts a decrease in the number of cholera cases in a shorter time interval of cholera intervention (2004–2019, see Table 8) as compared to the ODE model's time interval of intervention (2004–2022, see Table 9).

Again from our sensitivity analysis, we obtain that for ${\mathcal{R}}_E$, unlike the ODE model in [5], where the most sensitive parameter with the intervention model is ${\gamma }_h$, the recovery rate of infectious individuals in the high-risk group and the least sensitive parameter is ${\zeta }_l$, the death rate of the cholera pathogen due to improved water sanitation in the low-risk group, in our discrete-time model with intervention, the most sensitive parameter is $v_l$, the fraction of susceptible individuals in the low-risk group that is vaccinated, while the least sensitive parameter is ${\rho }_h$, the transmission rate between susceptible individuals and infectious individuals in the high-risk group. Another important parameter is ${\alpha }_{lh}$, the fraction of susceptible individuals that transition from the low-risk group to the high-risk group. Among the six intervention parameters, unlike in [5], where the most sensitive parameter is $v_h$, the vaccination rate of susceptible individuals in the high-risk group, and the least sensitive parameter is ${\zeta }_l$, the death rate of the cholera pathogen due to improved water sanitation in the low-risk group, the most sensitive parameter in our discrete-time model is $v_l$, the fraction of susceptible individuals in the low-risk group that is vaccinated, and the least sensitive parameter is ${\zeta }_h$, the fraction of the pathogen in the high-risk group that dies due to improved water sanitation.

In the ODE model in [5], the parameters $w_l$ (scaled contact rate between susceptible individuals and cholera pathogen in the low-risk group), ${\tau }_h$ (the rate at which recovered individuals become susceptible individuals in the high-risk group) and ${\omega }_{lh}$ (transition rate of recovered individuals from low-risk group to high-risk group) show a significant positive correlation with the total number of model predicted cholera cases. Meanwhile, ${\gamma }_l$ and ${\gamma }_h$ (the rates at which infectious individuals become recovered individuals in the low and high-risk groups, respectively), λ (the net decay rate of the cholera pathogen) and ${\omega }_{hl}$ (transition rate of recovered individuals from high-risk group to low-risk group) show a significant negative correlation with the total number of model predicted cholera cases. In our discrete-time model, ${\beta }_l$, the contact rate between susceptible individuals and the cholera pathogen in the low-risk group shows a significant positive correlation with the total number of our discrete-time model predicted cholera cases. On the other hand, λ, the net decay rate of the pathogen shows a significant negative correlation with the total number of our discrete-time model predicted cholera cases. These global sensitivity results are in agreement with our discrete-time model predictions.

In Table 12, we summarize the results for our discrete-time model and the continuous-time ODE models. Our results show that both models lead to different results even when they have almost the same ${\mathcal{R}}_0$. Both models are based on observations and assumptions of the cholera disease dynamics in Cameroon. However, unlike the continuous-time model, the discrete-time model is also based on the observation that the Cameroon data are reported annually (at discrete-time intervals). Thus, the difference equations of the discrete-time model are a closer fit to the observed cholera cases in Cameroon than the ordinary differential equations of the continuous-time cholera model.

Table 12. Results for discrete-time model and continuous-time ODE model.

	Discrete-Time Model	Continuous-Time ODE Model
${\mathcal{R}}_0$	1.1832	1.1803
Relative error without intervention	0.8441	0.8637
Trend in number of cholera cases	Oscillatory	Decreasing and Increasing
Sensitivity indices of ${\mathcal{R}}_0$	Most sensitive parameter is ${\beta }_h$	Most sensitive parameter is ${\gamma }_h$
Sensitivity indices of ${\mathcal{R}}_0$	Least sensitive parameters are ${\tau }_l,{\tau }_h$ and $w_{lh},w_{hl}$	Least sensitive parameters are ${\tau }_l,{\tau }_h$ and ${\omega }_{lh},{\omega }_{hl}$
Sensitivity of parameters on number cholera cases without intervention	Most significant parameters are ${\beta }_h,{\sigma }_h,{\xi }_h$	Most significant parameters are $w_h,{\gamma }_h$
Single intervention policy	Vaccination only smallest ${\mathcal{R}}_E$ with highest percent decrease from 2015 to 2019	Vaccination only smallest ${\mathcal{R}}_E$ with highest percent decrease from 2015–2022
Double intervention policy	Vaccination and Treatment smallest ${\mathcal{R}}_E$ with highest percent decrease from 2015 to 2019	Vaccination and Treatment smallest ${\mathcal{R}}_E$ with highest percent decrease from 2015 to 2022
Combined intervention policy	shorter period from 2004 to 2019 with $99.9565\mathrm{\%}$ decrease	longer period from 2004 to 2022 with $99.7280\mathrm{\%}$ decrease
Sensitivity indices of ${\mathcal{R}}_E$	Most sensitive parameter is $v_l$	Most sensitive parameter is ${\gamma }_h$
Sensitivity indices of ${\mathcal{R}}_E$	Least sensitive parameter is ${\rho }_h$	Least sensitive parameter is ${\zeta }_l$
Sensitivity of parameters on number cholera cases with intervention	Most significant parameters are ${\beta }_l,\lambda$	Most significant parameters are $w_l,{\tau }_h$ and ${\omega }_{lh},{\gamma }_l,{\gamma }_h,\lambda ,{\omega }_{hl}$

5. Conclusion

Usually, populations are censored at discrete-time intervals. To capture this discrete censoring of populations in a disease model, we introduce a simple low-high-risk discrete-time structured parametric cholera model of Cameroon with no spatial structure but with direct (human-to-human) and indirect (contaminated water-to-human) infection pathways. We use the discrete-time demographic equation (disease-free equation) to fit the total annual population of Cameroon and then use the fitted discrete-time risk-structured model to capture the reported cholera cases in Cameroon from 1987–2004. In Figure 3, we obtain that, as in the ODE model in [5], our demographic equation captures the observed exponential population growth of Cameroon. The demographic threshold for the fitted model, ${\mathcal{R}}_d=3.0999$. The basic reproduction number for our discrete-time model is ${\mathcal{R}}_0=1.1832$, and that of the ODE cholera model in [5] with a similar risk-structure is 1.1803. Thus, the two models have approximately the same value for ${\mathcal{R}}_0$ and both predict cholera endemicity in Cameroon. From Figure 5, unlike the ODE model that predicted a decreasing trend from 1987 to 1994 and an increasing trend from 1995 to 2004 in the pre-intervention reported number of cholera cases in Cameroon from 1987 to 2004, our fitted discrete-time risk-structured cholera model predicts an increasing trend from 1987 to 1989, 1990 to 1992, 1993 to 1995, 1997 to 1998 and 2003 to 2004, and a decreasing trend from 1989 to 1990, 1992 to 1993, 1995 to 1997 and 1998 to 2003, in the pre-intervention reported number of cholera cases in Cameroon from 1987 to 2004.

As in [5], we perform sensitivity analysis to determine the impact of each model parameter on the threshold parameters ${\mathcal{R}}_d$ and ${\mathcal{R}}_0$ and on the total number of our model's predicted cholera cases. As in the ODE model in [5], we obtained that ${\mathcal{R}}_d$ is not sensitive to changes to the parameters of the demographic equation. For ${\mathcal{R}}_0$, unlike the ODE model in [5], where the most sensitive parameter is ${\gamma }_h$, the recovery rate of infectious individuals in the high-risk group, in our discrete-time model, the most sensitive parameter is ${\beta }_h$, the contact rate between the susceptible individuals and the cholera pathogen in the high-risk group. In [5], the least sensitive parameters are ${\tau }_l$ and ${\tau }_h$, the parameters associated with loss of immunity, and ${\omega }_{lh}$ and ${\omega }_{hl}$, transition rates of recovered individuals between the two risk groups. In our discrete-time model, the least sensitive parameters are ${\tau }_l$ and ${\tau }_h$, the fraction of recovered individuals in the high-risk group who become susceptible, and $w_{lh}$ and $w_{hl}$, the fraction of recovered individuals who transition from the low-risk group to the high risk group and from the high-risk group to the low-risk groups, respectively. Another important parameter that has a high impact on ${\mathcal{R}}_0$ in our discrete-time model is ${\sigma }_h$, the fraction of exposed individuals in the high-risk group that can progress to the infectious class. In the ODE model in [5], $w_h$, the scaled contact rate between susceptible individuals and the cholera pathogen in the high-risk group shows a significant positive correlation with the total number of cholera cases and ${\gamma }_h$, the recovery rate of infectious individuals in the high risk group shows a significant negative correlation with the total number of cholera cases. In our discrete-time model, ${\beta }_h$, the contact rate between susceptible individuals and the cholera pathogen in the high risk, ${\sigma }_h$, the fraction of exposed individuals that can progress to the infectious class and ${\xi }_h$, the shedding rate of pathogens by infectious individuals in the high risk group show a significant positive correlation with the total number of choleracases.

As in [5], using our fitted low-high discrete-time risk-structured cholera model, we study the impact of the following seven cholera intervention strategies on the number of reported cholera cases in Cameroon:

• Vaccination only

• Treatment only

• Improved Sanitation only

• Vaccination and Treatment

• Treatment and Improved Sanitation

• Vaccination and Improved Sanitation

• Vaccination, Treatment and Improved Sanitation

From Table 8, we see that ${\mathcal{R}}_E<{\mathcal{R}}_0$. So, as in [5], adopting any of the intervention strategies leads to a smaller reproduction number. Under the single intervention policies of vaccination only or treatment only or improved sanitation only, in Table 8, we obtained that, as in [5], a vaccination only policy has the smallest value of the effective reproduction number ${\mathcal{R}}_E$ and the highest percentage decrease in the number of cholera cases. As in [5], under the double intervention policies of vaccination and treatment only or treatment and improved sanitation only or vaccination and improved sanitation only, in Table 8 we obtained that the double policy of vaccination and treatment has the smallest value for ${\mathcal{R}}_E$ and the highest percentage decrease in the number of cholera cases. In Figures 7 and 8, we illustrate that implementation of any of the single or double policies lead to a decreasing trend in the number of cholera infections in Cameroon. The dual strategies of either vaccination and treatment or vaccination and improved sanitation or the combined strategy of vaccination, treatment and improved sanitation reduce the basic reproduction number of Cameroon from 1.1832 to 0.9909, 1.1832 to 0.9989 and 1.1832 to 0.9872, respectively, and the number of cholera cases by $99.7119\mathrm{\%}$, $97.5075\mathrm{\%}$ and $99.9565\mathrm{\%}$, respectively. Thus, as in [5], each of these three strategies is capable of eliminating cholera in Cameroon with the combined strategy having the lowest value for the effective reproduction number, ${\mathcal{R}}_E$, and the highest percentage decrease in the number of cholera cases. Overall, the combined strategy of vaccination, treatment and improved sanitation lead to the smallest number of cholera infections in Cameroon (Figure 9).

From Tables 8 and 9, our discrete-time model results confirm that of the ODE model. However, since the discrete-time model is built on the discrete-time census data it predicts a decrease in the number of cholera cases in a shorter time interval of cholera intervention (2004–2019, see Table 8) as compared to the ODE model's time interval of intervention (2004–2022, see Table 9).

We perform sensitivity analysis to determine the impact of each model parameter on the threshold ${\mathcal{R}}_E$ and on the total number of predicted cholera cases by our model with intervention. For ${\mathcal{R}}_E$, unlike the ODE model in [5], where the most sensitive parameter with the intervention model is ${\gamma }_h$, the recovery rate of infectious individuals in the high-risk group and the least sensitive parameter is ${\zeta }_l$, the death rate of the cholera pathogen due to improved water sanitation in the low-risk group, in our model the most sensitive parameter in our model is $v_l$, the fraction of susceptible individuals in the low-risk group that is vaccinated, while the least sensitive parameter is ${\rho }_h$, the transmission rate between susceptible individuals and infectious individuals in the high-risk group. Another important parameter is ${\alpha }_{lh}$, the fraction of susceptible individuals that transition from the low-risk group to the high-risk group. Among the six intervention parameters, unlike in [5], where the most sensitive parameter is $v_h$, the vaccination rate of susceptible individuals in the high-risk group, and the least sensitive parameter is ${\zeta }_l$, the death rate of the cholera pathogen due to improved water sanitation in the low-risk group, the most sensitive parameter in our discrete-time model is $v_l$, the fraction of susceptible individuals in the low-risk group that are vaccinated, and the least sensitive parameter is ${\zeta }_h$, the fraction of the pathogen in the high-risk group that dies due to improved water sanitation.

In the ODE model in [5], the parameters $w_l$ (scaled contact rate between susceptible individuals and cholera pathogen in the low-risk group), ${\tau }_h$ (rate at which recovered individuals become susceptible individuals in the high-risk group) and ${\omega }_{lh}$ (transition rate of recovered individuals from low-risk group to high-risk group) show a significant positive correlation with the total number of model predicted cholera cases. Meanwhile, ${\gamma }_l$ and ${\gamma }_h$ (the rates at which infectious individuals become recovered individuals in the low and high-risk groups, respectively), λ (the net decay rate of the cholera pathogen) and ${\omega }_{hl}$ (transition rate of recovered individuals from high-risk group to low-risk group) show a significant negative correlation with the total number of model predicted cholera cases. In our discrete-time model, ${\beta }_l$, the contact rate between susceptible individuals and the cholera pathogen in the low-risk group shows a significant positive correlation with the total number of our discrete-time model predicted cholera cases. On the other hand, λ, the net decay rate of the pathogen shows a significant negative correlation with the total number of our discrete-time model predicted cholera cases. These global sensitivity results are in agreement with our discrete-time model predictions.

Our results show that using ODE and difference equations to model populations that are censored at discrete-time intervals can lead to different results even when the different models have almost the same ${\mathcal{R}}_0$. Both the continuous-time and discrete-time cholera models of Cameroon are based on observations and assumptions of the cholera disease dynamics in Cameroon. However, unlike the continuous-time model, the discrete-time model is also based on the observation that the cholera cases in Cameroon are reported annually. Thus, the difference equations of our novel discrete-time model are a closer fit to the observed cholera cases in Cameroon than the ordinary differential equations of the continuous-time cholera model. As in [5], our discrete-time models are not spatially structured. When data become available and our model is fully parameterized, a study of a discrete-time spatially structured cholera model with cost analysis of intervention strategies in Cameroon will be an interesting extension of this work.

Acknowledgments

We also appreciate travel awards to support travel expenses from the Society of Industrial and Applied Mathematics, and Biology and Medicine through Mathematics. We thank the anonymous reviewers for a thorough reading and helpful suggestions that improved our exposition.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

We would like to acknowledge the support of the Graduate School of College of Arts and Science, Howard University through their award of Just-Julian Graduate Research Assistantship. This work was also partially supported by the National Science Foundation [grant number CCF-1522054].