A Multi-Patch Model of Vibrio cholera Transmission with Food Insecurity and Environmental Coupling: From Theory to Application

1. Introduction

Cholera remains a persistent public health concern in many developing regions, especially where access to clean water, adequate nutrition, and basic sanitation is limited. In Sub-Saharan Africa, Vibrio cholerae continues to cause recurrent epidemics, leading to significant morbidity and mortality, particularly among vulnerable populations. Food insecurity and malnutrition, exacerbated by poverty and environmental stressors such as droughts or floods, may amplify the risk and severity of cholera outbreaks (FAO 2021; UNICEF 2022). Poor nutritional status weakens the immune system and may increase susceptibility to enteric infections, creating a vicious feedback loop between undernutrition and infectious diseases (Guerrant et al. 2013; Bhutta et al. 2008). Environmental factors also play a central role in cholera transmission. V. cholerae can survive and proliferate in aquatic reservoirs, with its persistence influenced by temperature, salinity, rainfall, and contamination from human activities (Bertuzzo et al. 2010; de Magny et al. 2008). In such contexts, the transmission of cholera is not only the result of direct human-to-human contact but also strongly coupled to environmental and social dynamics. Consequently, modeling approaches that incorporate both environmental transmission pathways and socio-economic vulnerabilities are essential to better understand disease persistence and control (Che et al. 2019). While numerous models of cholera transmission have been proposed (Codeço 2001; Hartley et al. 2006; Tien and Earn 2010b), few have integrated the spatial heterogeneity of both population vulnerability and environmental exposure across multiple communities or regions. Even fewer have addressed the compounded effects of food insecurity on cholera dynamics. Yet, these two factors malnutrition and environmental contamination can interact in ways that strongly affect both the spread and severity of the disease, particularly in decentralized or poorly connected health systems.

This work is motivated by the following question: how does food insecurity influence cholera transmission and persistence at the population level, particularly when environmental contamination and spatial connectivity are accounted for? To address this, we propose a novel multi-patch compartmental model that captures the interplay between cholera transmission, food availability, and spatial structure. The model explicitly considers two distinct classes of susceptible individuals those in food security and those in food insecurity each with different vulnerabilities to infection and mortality. It also includes compartments for acutely infected individuals, chronically infected carriers, and an environmental reservoir representing water contamination. Here, we consider a simple extension of the classical SIR model with water compartment W by adding a foodborne transmission to obtain the resulting “SCIWR-F”model, allowing for both person-person, person-water-person and person-food-person transmission. Our goal is to provide a theoretical and numerical framework that allows us to investigate the impact of food insecurity on the burden of cholera, identify threshold parameters (such as basic reproduction numbers, ${\mathcal{R}}_0$) associated with disease persistence or elimination, and explore scenarios where food-related stress may give rise to complex dynamics such as backward bifurcation or multiple endemic equilibria. Numerical simulations calibrated to World bank data both malnutrition and waterborne disease burden are high used to support the theoretical results and quantify the potential health losses due to the coupling of food stress and epidemic dynamics. In particular, we show that food insecurity, acting as a vulnerability amplifier, may increase both the force of infection and the cholera-induced mortality, thereby altering the stability landscape of the system. The presence of multiple patches further reveals the importance of spatial feedbacks, with poorly connected regions acting as sources or sinks of infection. This study highlights the critical need to address food security not only as a development goal but also as a core component of epidemic preparedness and resilience planning.

The paper is organized as follows. In Section 4, we present the compartmental model and its assumptions. A detailed analysis of the life cycle of V. cholerae and the structure of food dynamics are also given. Section 5 contains the mathematical analysis, including the computation of ${\mathcal{R}}_0$ and equilibrium properties. In Section 6, we present numerical simulations based on data from Douala and its surrounding areas. Conclusions and perspectives are given in Section 7.

2. Life cycle and biological background of Vibrio cholerae

Vibrio cholerae is a Gram-negative, facultative anaerobic bacterium that thrives in aquatic environments, especially in estuarine, brackish, and coastal waters (Lizárraga-Partida et al. 2009). Its life cycle is strongly influenced by ecological and environmental conditions such as temperature, salinity, and nutrient availability (Colwell and Huq 2004; Nelson et al. 2009). The bacterium alternates between two major phases: a free-living stage in the environment and a parasitic stage within human hosts. In its environmental phase, V. cholerae is capable of surviving in both planktonic and biofilm-associated forms. It adheres to biotic surfaces such as copepods, phytoplankton, and chitinous exoskeletons of aquatic invertebrates (Lipp et al. 2002; Huq et al. 1995), facilitating persistence in nutrient-limited waters. Environmental survival is further enhanced by the transition to a viable but non culturable state under unfavorable conditions such as low temperature or nutrient deprivation (Oliver 2005)). The bacterium can remain viable in this state for extended periods and regain infectivity when conditions become favorable.

Upon ingestion of contaminated water or food, V. cholerae enters the human gastrointestinal tract, where it must overcome gastric acidity and colonize the small intestine. This colonization involves chemotaxis, mucin penetration, and the expression of key virulence factors such as the toxin-coregulated pilus and cholera toxin (Kaper et al. 1995). Within the host, the bacterium multiplies rapidly, leading to massive fluid loss through diarrhea, which in turn contributes to environmental recontamination . Excreted bacteria are often in a hyperinfectious state for several hours after being shed, with significantly increased infectivity compared to environmental strains (Codeço 2001; Merrell et al. 2002). This hyperinfectious phase plays a critical role in epidemic amplification. Once reintroduced into the environment, the bacteria return to their aquatic phase, completing the cycle. The environmental-to-human-to-environment loop, modulated by seasonality, temperature, and human behavior, defines the full transmission cycle of V. cholerae. Understanding this cycle is essential to accurately model both the environmental persistence and the outbreak dynamics of cholera.

3. A revision of old models: SIWR-F

In a quest to improve existing models of cholera, we undertook a revision of these model with the aim of making them more realistic. Indeed, cholera transmission is not only waterborne through the ingestion of contaminated water or contact with bodily fluids from infected individuals but also foodborne, via the interaction with food/biomass, a route that current models have yet to capture. By introducing a term of the form ${\beta }_W\left(1+{\eta }_F\left(F\right)\right)·W$, we propose a novel coupling between contamination dynamics and transmission risk, bridging ingestion-based exposure and envronmental pathways. Moreover, if $F=0$, we obtain classicals force of infection such as (Tien and Earn 2010a,b). Then, we propose this model (1) with Flowchart given in Figure 1, variables and parameters in Table 1.

$$\begin{array}{cc}\hfill \dot{S}& =\mu N-b_W(1+{\eta }_F\left(F\right))WS-b_{I}SI-\mu S,\hfill \end{array}$$

(1a)

$$\begin{array}{cc}\hfill \dot{I}& =b_W(1+{\eta }_F\left(F\right))WS+b_{I}SI-(\gamma +\mu )I,\hfill \end{array}$$

(1b)

$$\begin{array}{cc}\hfill \dot{W}& =\alpha I-\xi W,\hfill \end{array}$$

(1c)

$$\begin{array}{cc}\hfill \dot{R}& =\gamma I-\mu R,\hfill \end{array}$$

(1d)

$$\begin{array}{cc}\hfill \dot{F}& =rF(1-{\displaystyle \frac{F}{K}})-aF,\mathrm{with}{\eta }_F\left(F\right)={\displaystyle \frac{F}{K}}.\hfill \end{array}$$

(1e)

By giving ${\displaystyle s=\frac{S}{N},i=\frac{I}{N},r=\frac{R}{N},w=\frac{\xi }{\alpha N}W}$ and $f={\displaystyle \frac{F}{K}}$, one has

$$\begin{array}{cc}\hfill \dot{s}& =\mu -{\beta }_W(1+f)ws-{\beta }_{I}si-\mu s,\hfill \end{array}$$

(2a)

$$\begin{array}{cc}\hfill \dot{i}& ={\beta }_W(1+f)ws+{\beta }_{I}si-(\gamma +\mu )i,\hfill \end{array}$$

(2b)

$$\begin{array}{cc}\hfill \dot{w}& =\xi (i-w),\hfill \end{array}$$

(2c)

$$\begin{array}{cc}\hfill \dot{r}& =\gamma i-\mu r,\hfill \end{array}$$

(2d)

$$\begin{array}{cc}\hfill \dot{f}& =\tilde{r}f(1-f)-af.\hfill \end{array}$$

(2e)

Assertion 3.1. 

In this study, we suppose that food is always present in the environment, in other words, $\forall t\ge ,F\left(t\right)>0$. This condition (r > a) ensures that the foodborne transmission component of cholera remains active at all times, which is the central focus of this work. Moreover, the case F = 0 reduces to the classical SIWR form already studied by (Tien and Earn 2010b).

Remark 3.1. 

Since the population is constant in this model, adding explicitly a term for food consumption by humans is of little interest.

3.1. Basic Properties of the SIWR-F Model

Before proceeding to qualitative analysis, we first verify that the SIWR-F model is mathematically well-posed and biologically meaningful. Since the right-hand sides of the system (2) are composed of multivariate polynomials and hence ${\mathcal{C}}^{\infty }$, then the system is locally Lipschitz. Therefore, by the Cauchy-Lipschitz theorem, there exists a unique local solution for any initial condition in ${\mathbb{R}}_{+}^5$.

Moreover, we need to ensure that the positive orthant is invariant under the flow of the system (2). Let $Z\left(t\right)=\left(s\right(t),i(t),w(t),r(t),f(t\left)\right)$ denote the solution vector of the system (2).

Lemma 3.1. 

Let $Z\left(0\right)\in {\mathbb{R}}_{+}^5$ with $s\left(0\right)>0$ and $i\left(0\right),w\left(0\right),r\left(0\right),f\left(0\right)\ge 0$. Then, the solution $Z\left(t\right)$ satisfies $Z\left(t\right)\in {\mathbb{R}}_{+}^5$ for all $t\ge 0$, and $s\left(t\right)>0$ for all $t\ge 0$.

Proof. 

The right-hand side of each equation in system (2) is such that if any variable reaches zero, its derivative becomes non-negative (except possibly for i and w, which depend on infections, but are fed by $s>0$ and the initial infection).

Take for instance $s\left(t\right)$: if $s\left(t_1\right)=0$ at some $t_1>0$, then from Equation (2a),

$$\dot{s}\left(t_1\right)=\mu >0,$$

contradicting the assumption that $s\left(t\right)$ reaches zero. The same reasoning applies to $f\left(t\right)$ and $w\left(t\right)$. Therefore, all state variables remain nonnegative and $s\left(t\right)$ remains strictly positive. □

We now show that the solution is uniformly bounded for all $t\ge 0$.

Lemma 3.2. 

The solution $Z\left(t\right)$ of system (2) is bounded $\forall t\ge 0$.

Proof. 

Since the total population is constant and the food dynamics are logistic, we can use standard comparison arguments. In particular, $s\left(t\right)+i\left(t\right)+r\left(t\right)\le 1$, and $f\left(t\right)\le 1$ for the normalized model. The bacteria compartment $w\left(t\right)$ satisfies the same reasoning. □

These properties confirm that the SIWR-F model is well-posed and suitable for epidemiological analysis.

Definition 3.1. 

The basic reproduction rate is usually the average number of newly infected that a vector can produce in a population made up entirely of susceptible individuals during it period of infection without any control.

Using technique provide by (van den Driessche and Watmough 2002), we write the next generation matrix at the disease free equilibrium as $F{V}^{-1}$, where the $ij$ entry of the matrix F is the rate at which infected individuals in compartment j produce new infections in compartment i, and the $jk$ entry of $V^{-1}$ is the average duration of stay in compartment j starting from k.

From (1), we have

$$F=\left(\begin{array}{cc}{\beta }_I& {\beta }_W(1+f^{*})\\ 0& 0\end{array}\right),\mathrm{with}{f}^{*}=1-{\displaystyle \frac{a}{r}}\mathrm{and}{V}^{-1}=\left(\begin{array}{cc}{\displaystyle \frac{1}{\gamma +\mu }}& 0\\ {\displaystyle \frac{1}{\gamma +\mu }}& {\displaystyle \frac{1}{\xi }}\end{array}\right)$$

(3)

Then, the basic reproductive number is given by:

$$\begin{array}{cc}\hfill {\mathcal{R}}_0& =\rho \left(F{V}^{-1}\right)\hfill \\ \hfill & ={\displaystyle \frac{{\beta }_I+{\beta }_W(1+f^{*})}{\gamma +\mu }}.\hfill \end{array}$$

(4)

To ensure that the effective control of the disease is not dependent on the initial size of V. cholerae concentration, a global stability result must be established for the disease free equilibrium point $x^0=(1,0,0,0)$.

Rappel 3.1. 

It is easy show that when $i\left(t\right)\to 0,w\left(t\right)\to 0.$

Lemma 3.3. 

If, ${\mathcal{R}}_0\le 1$, system (1) has a disease free equilibrium who is globally asymptotically stable and unstable when ${\mathcal{R}}_0>1$.

Proof. 

The local asymptotic stability follows clearly using (van den Driessche and Watmough 2002) (Theorem 2). Consider the following Lyapunov function using (Korobeinikov and Wake 2002) reasoning.

$$V(s,i)=s-1-ln\left(s\right)+i,s,i\ge 0.$$

(5)

Since $s\le 1,f\le f^{*}$ ,because the system (2) is in proportions. Yields,

$$\begin{array}{cc}\hfill {\displaystyle \frac{dV}{dt}}{|}_{DFE}& =\dot{s}\left(1-{\displaystyle \frac{1}{s}}\right)+\dot{i},\hfill \\ \hfill & =(\mu -{\beta }_W(1+f)ws-{\beta }_{I}si-\mu s)\left(1-{\displaystyle \frac{1}{s}}\right)+i(s({\beta }_W(1+f)+{\beta }_I)-(\gamma +\mu )),\hfill \\ \hfill & \le (\mu -\mu s)\left(1-{\displaystyle \frac{1}{s}}\right)+i\left({\displaystyle \frac{{\beta }_W(1+f^{*})-{\beta }_I}{\gamma +\mu }}-1\right),\hfill \\ \hfill & \le -{\displaystyle \frac{\mu }{s}}{(s-1)}^2+i(\gamma +\mu )({\mathcal{R}}_0-1),\hfill \\ \hfill & \le 0.\hfill \end{array}$$

(6)

Since, ${\displaystyle \frac{dV}{dt}}{|}_{DFE}=0$ only at the disease free equilibrium, thus the largest invariant set

$\mathsf{\Gamma }=\left\{(s,i)\in [0,1]\times [0,1]:{\displaystyle \frac{dV}{dt}}{|}_{DFE}=0\right\}$ is the singleton $(1,0)$. By (LaSalle 1976) invariant principle, one has that the disease free equilibrium is globally asymptotically stable in $\mathsf{\Gamma }$. This achieves the proof. □

3.2. Stability of the endemic equilibrium

Herein, we study the stability of an endemic equilibrium whenever it exists denote $x^{*}=(s^{*},i^{*},w^{*},f^{*})$ of system (2), where

$$(s^{*},i^{*},w^{*},f^{*})=\left({\displaystyle \frac{1}{{\mathcal{R}}_0}},{\displaystyle \frac{\mu }{\gamma +\mu }}(1-s^{*}),i^{*},1-{\displaystyle \frac{a}{r}}\right).$$

(7)

and ${\mathcal{R}}_0$ define in Eq.4.

We compute the Jacobian matrix (2) evaluated at $x^{*}$:

$${\displaystyle {J|}_{x^{*}}=\left(\begin{array}{cccc}-\frac{\mu }{s^{*}}& -{\beta }_{I}s& -{\beta }_W{s}^{*}(1+f^{*})& -{\beta }_W{s}^{*}{w}^{*}\\ \mu \left(-1+\frac{1}{s^{*}}\right)& -{\beta }_W(1+f^{*})s^{*}& {\beta }_W{s}^{*}(1+f^{*})& {\beta }_W{s}^{*}{w}^{*}\\ 0& \xi & -\xi & 0\\ 0& 0& 0& r(1-2f^{*})-a\end{array}\right).}$$

(8)

because, at equilibrium ($\dot{s}=0$),

$$\mu =s^{*}({\beta }_W(1+f^{*})w^{*}+{\beta }_I{i}^{*}+\mu )$$

This implies that

$${\displaystyle \frac{\mu }{s^{*}}}={\beta }_W(1+f^{*})w^{*}+{\beta }_I{i}^{*}+\mu$$

(9)

Moreover, for the case ($\dot{i}=0$), we have

$${\beta }_W(1+f^{*})s^{*}{w}^{*}+{\beta }_I{s}^{*}{i}^{*}-(\gamma +\mu )i^{*}=0,$$

Using the fact that at equilibrium, $w^{*}=i^{*}$, one has

$$({\beta }_W(1+f^{*})s^{*}+{\beta }_I{s}^{*})i^{*}=(\gamma +\mu )i^{*},$$

Then

$$\gamma +\mu ={\beta }_W(1+f^{*})s^{*}+{\beta }_I{s}^{*},$$

(10)

which implies

$${\beta }_I{s}^{*}-(\gamma +\mu )={\beta }_I{s}^{*}-{\beta }_W(1+f^{*})s^{*}-{\beta }_I{s}^{*}=-{\beta }_W(1+f^{*})s^{*}.$$

(11)

The characteristic polynomial is ${\lambda }^4+a_1{\lambda }^3+a_2{\lambda }^2+a_3\lambda +a_4,$ where

$$\begin{array}{cc}\hfill & a_1={\beta }_W(1+f^{*})s^{*}+{\displaystyle \frac{\mu }{s^{*}}}+\xi -(-r+a),\hfill \\ \hfill & a_2={\beta }_W(1+f^{*})(\mu -(-r+a)s^{*})+\mu {\beta }_I(1-s^{*})+{\displaystyle \frac{\mu }{s^{*}}}(\xi -(-r+a))-(-r+a)\xi ,\hfill \\ \hfill & a_3=\mu \left[(1-s^{*})\left({\beta }_W(1+f^{*})\xi +{\beta }_I(\xi -(-r+a))\right)-(-r+a)\left({\beta }_W(1+f^{*})+{\displaystyle \frac{\xi }{s^{*}}}\right)\right],\hfill \\ \hfill & a_4=\xi \mu (r-a)(1-s^{*})({\beta }_I+{\beta }_W(1+f^{*})).\hfill \end{array}$$

(12)

Since ${\mathcal{R}}_0>1$ and $s^{*}={\displaystyle \frac{1}{{\mathcal{R}}_0}}$, we have $(1-s^{*})>0$, so $a_i>0,\forall i=1,\dots ,4$. The Routh-Hurwitz criteria give that the endemic equilibrium is stable if $a_1,a_4>0$ and $a_1{a}_2-a_3>0,a_1{a}_2{a}_3-{\left(a_1\right)}^2{a}_4-{\left(a_3\right)}^2>0$. Clearly, we prove that $a_1{a}_2-a_3>0$. Indeed:

$$\begin{array}{cc}\hfill a_1{a}_2-a_3& =a_2\left({\beta }_W(1+f^{*})s^{*}+{\displaystyle \frac{\mu }{s^{*}}}\right)+a_2(\xi -(-r+a))-a_3,\hfill \\ \hfill & =a_2\left({\beta }_W(1+f^{*})s^{*}+{\displaystyle \frac{\mu }{s^{*}}}\right)+{\beta }_W\mu \xi (1+f^{*})s^{*}+{(-r+a)}^2\left({\beta }_W(1+f^{*})s^{*}+{\displaystyle \frac{\mu }{s^{*}}}+\xi \right)\hfill \\ \hfill & +\xi (r-a)\left({\beta }_W(1+f^{*})s^{*}+{\displaystyle \frac{\mu }{s^{*}}}+\xi \right)+{\displaystyle \frac{\mu {\xi }^2}{s^{*}}},\hfill \\ \hfill & >0.\hfill \end{array}$$

(13)

Thus, local stability of the endemic equilibrium is determined by the sign of the last condition

$$\begin{array}{cc}\hfill C& =a_1(a_2{a}_3-a_1{a}_4)-a_3^2,\hfill \\ \hfill & =a_1((1-s^{*})\xi \mu (r-a){\displaystyle \frac{\left({\beta }_W(1+f^{*})s^{*2}+\mu +s^{*}(\xi -(-r+a))\right)\left({\beta }_I+{\beta }_W(1+f^{*})\right)}{s^{*}}}+\hfill \\ \hfill & ({\beta }_W(1+f^{*})\mu -{\beta }_W(1+f^{*})s^{*}(-r+a)+{\beta }_I\mu (1-s^{*})-\xi (-r+a)+{\displaystyle \frac{\mu (\xi -(-r+a\left)\right)}{s^{*}}}{a}_3)-a_3^2,\hfill \\ \hfill & >0.\hfill \end{array}$$

(14)

Thus, the endemic equilibrium in system (2) is locally stable whenever it exists and we prove that it persists inside the population when it appears using (Thieme 1993) in Lemma.

The left panel in Figure 2 displays a forward bifurcation: as the basic reproduction number ${\mathcal{R}}_0$ crosses the threshold of 1, the disease-free equilibrium loses stability and an endemic equilibrium emerges. This confirms the analytical result (Lemma 3.3) that is a critical threshold for cholera persistence. The right panel presents the phase portrait of the system, showing the trajectories in the plane. It demonstrates that trajectories starting from different initial conditions converge to the endemic equilibrium (if ${\mathcal{R}}_0>1$). The vector field highlights the direction of movement, reinforcing the local and global stability properties derived in the present section.

3.3. Comparison of SIR, SIWR and SIWR-F dynamics

Figure 3. Comparison between SIR, SIWR and SIWR-F under cholera reported cases in Cameroon during 1976-2000 according to World Bank and WHO Data.

Figure 3. Comparison between SIR, SIWR and SIWR-F under cholera reported cases in Cameroon during 1976-2000 according to World Bank and WHO Data.

Figure 4. Simulation of system (2) when ${\mathcal{R}}_0>1$. We present dynamics of susceptible and infected during an epidemic.

Figure 4. Simulation of system (2) when ${\mathcal{R}}_0>1$. We present dynamics of susceptible and infected during an epidemic.

4. An extension : Multi-patch framework

However, this new modifications inside system (1) to correct the previous models (Anderson and May 1997, Eisenberg et al 2013, (Bertuzzo et al. 2010), (Tien and Earn 2010a), (Che et al. 2019)... It do not taking into account some others realities like migration, presence of asymptomatic carriers, non-linearity of diseases transmissions and vulnerability of susceptible due to poverty, under-nourrished or others. So, the goal of this section is to solve them.

4.1. Model Formulation

We need the following hypotheses:

Assertion 4.1. 

• The foodborne transmission channel mediates environmental contamination by V. cholerae, thus capturing both indirect environmental and ingestion-based exposure routes. Then, $\lambda (.)={\beta }_B·{\displaystyle \frac{F_i^c}{{\Delta }_1+F_i^c}}+{\beta }_H·{\displaystyle \frac{I_i+\epsilon C_i}{N_i^{q_p}}},$ where $\lambda (.)$ is the force of infection with ${\eta }_F={\beta }_F·{\displaystyle \frac{F_i}{{\Delta }_2+F_i}}$ and $q_p\in \{0,1\}$ i.e inside the principal node, the transmission look the standard incidence function ($q_p=1$) and the mass action law for their satellites ($q_p=0$).
• ${\varphi }_i=\left\{\begin{array}{cc}max(0,\delta (K-e{a}_i{F}_i))\hfill & i=1\hfill \\ max(0,\delta (e{a}_i{F}_i-K))\hfill & i=2\hfill \end{array}\right.$. For the transition to vulnerable individuals, it is assumed that vulnerability appears when the amount of nutriments absorbed (converted into biomass) is below than a threshold K (2000 calories.day⁻¹), i.e $e{a}_i{F}_i<K$, or even $K-e{a}_i{F}_i>0$. Individuals become vulnerable at rate ${\varphi }_i$, with δ the apparition rate of vulnerability.
• $d(.)=d+{\displaystyle \frac{\theta }{e{a}_i{S}_{2,i}+1}}$ where d is the mortality rate due to the disease and η the maximal mortality rate possible due to the severe food insecurity (if any food is absorb).

Here, we develop a mathematical model that describes the multi-patch dynamics of cholera transmission in a population subject to nutritional vulnerability and contaminated food. The model is formulated using a compartmental framework, and explicitly considers both spatial heterogeneity (through patch structure) and multiple transmission routes (direct and indirect). Individuals move between compartments according to their epidemiological status and are exposed to infection either via direct contact with infectious individuals or through ingestion of contaminated food.

The model rests on the following assumptions:

(i)

Individuals who are food insecure or nutritionally vulnerable are more likely to become infected and suffer severe outcomes if exposed to cholera.

(ii)

Individuals within a given patch are homogeneously mixed, but inter-patch coupling exists via human mobility and environmental contamination.

(iii)

Food contamination arises from bacterial load in the environment and is modulated by local food availability and hygiene conditions.

(iv)

Since the study focused on a relatively short period (January-April 2022 epidemic), long-term effects such as loss of immunity or natural mortality were neglected.

(v)

We assume that food and water ingestion occur jointly during meals, as is common in many societies worldwide. Therefore, waterborne and foodborne exposures are combined into a single ingestion-based transmission route. This simplifies the model compared to classical cholera frameworks where water is treated separately.

(vi)

Food contamination is assumed to result primarily from environmental exposure to V. cholerae, reflecting the dominant route of contamination observed in cholera-endemic settings. Direct human contamination (e.g., via food handling) is not explicitly modeled, as its contribution is generally secondary compared to water-related pathways.

Each spatial unit (or patch) contains human and environmental compartments. At time , individuals in patch are divided into the following disjoint classes: the Susceptible individuals with adequate nutrition ($S_{1,i}\left(t\right)$), the Susceptible individuals under nutritional stress or vulnerable ($S_{2,i}\left(t\right)$), Symptomatic infectious individuals ($I_i\left(t\right)$) , Asymptomatic carriers ($C_i\left(t\right)$) and the Recovered individuals with temporary immunity ($R_i\left(t\right)$). Environmental compartments include: Concentration of Vibrio cholerae in the local environment $\left(B_i\left(t\right)\right)$, Safe food availability ($F_i^s\left(t\right)$) and Contaminated food ($F_i^c\left(t\right)$). Individuals can move from $S_1$ to $S_2$ at rate ${\varphi }_1$ when their nutritional intake falls below a threshold K and return to $S_1$ at rate ${\varphi }_2$.

The total population in patch at time is:

$$N_i\left(t\right)=S_{1,i}\left(t\right)+S_{2,i}\left(t\right)+I_i\left(t\right)+C_i\left(t\right)+R_i\left(t\right).$$

(15)

The dynamics of cholera in each patch are governed by both local processes (e.g., infection, recovery, bacterial shedding) and inter-patch interactions (e.g., movement of individuals, environmental exchange). The infection force combines two mechanisms:

1.

Direct human-to-human transmission, proportional to the prevalence of infectious individuals.

2.

Food-borne transmission through ingestion of contaminated food, governed by a saturation function of disease.

Safe food availability is regenerated logistically and consumed by the population. The contamination of food depends on environmental bacterial concentration and decreases with better hygiene at rate $\omega$. The concentration of bacteria increases with shedding from both symptomatic and asymptomatic individuals and decays naturally or via dilution/migration who is temperature-dependent at rate ${\mu }_B\left(T\right)$. Contaminated food is produced at rate proportional to the product of bacteria and available food , and decays at a rate depending on sanitation and population size. Migration between patches is modeled via matrices that govern movement of susceptible, infected, or bacteria, allowing the disease to propagate spatially. The diagram of Figure 5 shows the overall flow of infection, food contamination, and vulnerability dynamics. Consider the following Flowchart Figure 5 with:

From the flowchart of diagram in Figure 5, the dynamics of the interaction between foods, the concentration of V. cholerae in the environment and human disease transmission is given by the following system of ordinary differential equations describing the temporal evolution of each class across all patches:

$$\begin{array}{cc}\hfill \dot{S_{1,i}}& =(1-\pi )\mathsf{\Lambda }+{\xi }_1{R}_i+{\varphi }_2{S}_{2,i}-\zeta \lambda (.)S_{1,i}-({\varphi }_1+\mu )S_{1,i}+\sum _{j=1}^n(m_{ji}^S{S}_{1,j}-m_{ij}^S{S}_{1,i}),\hfill \end{array}$$

(16a)

$$\begin{array}{cc}\hfill \dot{S_{2,i}}& =\pi \mathsf{\Lambda }+{\xi }_2{R}_i+{\varphi }_1{S}_{1,i}-\lambda (.)S_{2,i}-({\varphi }_2+\mu )S_{2,i}+\sum _{j=1}^n(m_{ji}^S{S}_{2,j}-m_{ij}^S{S}_{2,i}),\hfill \end{array}$$

(16b)

$$\begin{array}{cc}\hfill \dot{C_i}& =(1-p)\lambda (.)(S_{2,i}+\zeta S_{1,i})-({\gamma }_C+\mu )C_i+\sum _{j=1}^n(m_{ji}^C{C}_j-m_{ij}^C{C}_i),\hfill \end{array}$$

(16c)

$$\begin{array}{cc}\hfill \dot{I_i}& =p\lambda (.)(S_{2,i}+\zeta S_{1,i})-({\gamma }_I+d(.)+\mu )I_i+\sum _{j=1}^n(m_{ji}^I{I}_j-m_{ij}^I{I}_i),\hfill \end{array}$$

(16d)

$$\begin{array}{cc}\hfill \dot{R_i}& ={\gamma }_C{C}_i+{\gamma }_I{I}_i-(\mu +{\xi }_1+{\xi }_2)R_i+\sum _{j=1}^n(m_{ji}^R{R}_j-m_{ij}^R{R}_i),\hfill \end{array}$$

(16e)

$$\begin{array}{cc}\hfill {\dot{F_i}}^s& =r_i{F}_i\left(1-{\displaystyle \frac{F_i^s}{K_1}}\right)-(a_i+\omega )F_i^s-\rho B_i{F}_i^s,\hfill \end{array}$$

(16f)

$$\begin{array}{cc}\hfill \dot{F_i^c}& =\rho B_i{F}_i^s-(\omega +a_i)F_i^c,\hfill \end{array}$$

(16g)

$$\begin{array}{cc}\hfill \dot{B_i}& ={\alpha }_C{C}_i+{\alpha }_I{I}_i-{\mu }_{B_i}\left(T\right)B_i+\sum _{j=1}^n(m_{ji}^B{B}_j-m_{ij}^B{B}_i).\hfill \end{array}$$

(16h)

With initial conditions $S_1\left(0\right)>0,S_2\left(0\right)>0,C\left(0\right)>0,I\left(0\right)>0,B\left(0\right)>0,F^s\left(0\right)>0$ and $F^c\left(0\right)>0$. Table 2 and Table 3 present the variables and parameters of system (16), respectively. With $S_i=S_{1,i}+S_{2,i}$, $F_i=F_i^s+F_i^c$ and $N_i=S_{1,i}+S_{2,i}+C_i+I_i+R_i$.

5. Mathematical analysis

5.1. Basic properties of (16)

First of all, we need to establish that system (16) is well-posed to ensure that the model makes sense biologically.

Proposition 5.1. 

The nonnegative orthant ${\mathbb{R}}_{+}^{8n}$ is positively invariant under the flow of system (16) if initial conditions satisfy $S_{1,i}\left(0\right),S_{2,i}\left(0\right),C_i\left(0\right),I_i\left(0\right),F_i^s\left(0\right),F_i^c\left(0\right),B_i\left(0\right)>0$ for all i and all other variables are nonnegative, then $X\left(t\right)\ge 0$ for all $t\ge 0$ and for all variables $X\in \{S_{1,i},S_{2,i},C_i,I_i,R_i,F_i^s,F_i^c,B_i\}$.

Proof. 

Assume the initial conditions are nonnegative with $S_{1,i}\left(0\right),S_{2,i}\left(0\right),C_i\left(0\right),I_i\left(0\right),F_i^s\left(0\right),$

$F_i^c\left(0\right),B_i\left(0\right)>0$ for all $i=1,\dots ,n$. Let $t_1>0$ be the minimal time at which some component $X_k\left(t_1\right)=0$ for $X\in \{S_1,S_2,I,C,R,F^s,F^c,B\}$, and $X_k\left(t\right)>0$ for all $t<t_1$. We investigate which component can be the first to exit the nonnegative orthant.

Suppose that for some i, $S_{1,i}\left(t_1\right)=0$. From the corresponding equation,

$${\dot{S}}_{1,i}=(1-\pi )(\mathsf{\Lambda }+\xi R_i)+{\varphi }_2{S}_{2,i}+-\zeta \lambda (.)S_{1,i}-({\varphi }_1+\mu )S_{1,i}+\sum _{j=1}^n(m_{ji}^S{S}_{1,j}-m_{ij}^S{S}_{1,i}),$$

(17)

we observe that at $t_1$, all negative terms vanish due to $S_{1,i}\left(t_1\right)=0$, and the remaining terms ${\varphi }_2{S}_{2,i}\left(t_1\right)$ and ${\sum }_j{m}_{ji}^S{S}_{1,j}\left(t_1\right)$ are nonnegative and strictly positive since $S_{2,i}\left(t_1\right)>0$ and $S_{1,j}\left(t_1\right)>0$ for $j\ne i$ (as $S_{1,i}$ is the first to reach zero). Therefore,

$${\dot{S}}_{1,i}\left(t_1\right)>0,$$

And it comes that

$${\int }_0^{t_1}{\dot{S}}_{1,i}\left(t\right)dt>0,$$

which implies

$$\begin{array}{cc}\hfill & S_{1,i}\left(t_1\right)-S_{1,i}\left(0\right)>0⟹S_{1,i}\left(t_1\right)>S_{1,i}\left(0\right)>0.\hfill \end{array}$$

Finally, one obtains $S_{1,i}\left(t_1\right)>0$, that contradict the above hypothesis. The same reasoning applies for $S_{2,i}$ using its equation.

Now suppose $I_i\left(t_1\right)=0$ for some i, with $I_i\left(t\right)>0$ for $t<t_1$. Then from

$${\dot{I}}_i=p\lambda (S_{2,i}+\zeta S_{1,i})-({\gamma }_I+d(.)+\mu )I_i+\sum _j(m_{ji}^I{I}_j-m_{ij}^I{I}_i),$$

(18)

we have ${\dot{I}}_i\left(t_1\right)=p\lambda (S_{2,i}+\zeta S_{1,i})+{\sum }_j{m}_{ji}^I{I}_j>0$, since $S_{1,i},S_{2,i},I_j>0$ for all j and $t<t_1$. This again contradicts the assumption that $I_i\left(t_1\right)=0$ and decreasing.

Similar reasoning applies to the asymptomatic class $C_i$, where the derivative at $t_1$ takes the form:

$${\dot{C}}_i=(1-p)\lambda (S_{2,i}+\zeta S_{1,i})+\sum _j{m}_{ji}^C{C}_j>0.$$

(19)

For the recovered class $R_i$, we have

$${\dot{R}}_i={\gamma }_C{C}_i+{\gamma }_I{I}_i-(\mu +\xi )R_i+\sum _j(m_{ji}^R{R}_j-m_{ij}^R{R}_i).$$

(20)

At $t_1$, if $R_i\left(t_1\right)=0$, using similar process the remaining terms are nonnegative and strictly positive, leading to ${\dot{R}}_i\left(t_1\right)>0$.

Next, consider the food compartments. If $F_i^s\left(t_1\right)=0$, from

$${\dot{F}}_i^s=r_i{F}_i\left(1-{\displaystyle \frac{F_i^s}{K_1}}\right)-(a_i+\omega )F_i^s-\rho B_i{F}_i^s,$$

(21)

we have ${\dot{F}}_i^s\left(t_1\right)=r_i{F}_i>0$ , since $F_i^c<K_1$ and $F_i=F_i^s+F_i^c>0$.

Similarly, if $F_i^c\left(t_1\right)=0$,

$${\dot{F}}_i^c=\rho B_i{F}_i^s>0,$$

since $B_i,F_i^s>0$.

For the bacterial concentration $B_i$, the equation

$${\dot{B}}_i={\alpha }_C{C}_i+{\alpha }_I{I}_i-{\mu }_B\left(T\right)B_i+\sum _j(m_{ji}^B{B}_j-m_{ij}^B{B}_i)$$

(22)

gives ${\dot{B}}_i\left(t_1\right)>0$ because $C_i,I_i>0$ and $B_i=0$ eliminates the negative term.

Hence, no variable can be the first to exit the nonnegative orthant. Therefore, all components remain nonnegative for all $t\ge 0$, and $S_{1,i}\left(t\right),S_{2,i}\left(t\right),F_i^s\left(t\right)$ remain strictly positive. □

Lemma 5.1. 

The total human population $T_H\left(t\right)={\sum }_{i=1}^n{N}_i\left(t\right)$ and the total bacterial concentration $T_B\left(t\right)={\sum }_{i=1}^n{B}_i\left(t\right)$ are bounded for all $t\ge 0$.

Proof. 

Consider first the evolution of the total human population:

$$T_H\left(t\right)=\sum _{i=1}^n{N}_i\left(t\right)=\sum _{i=1}^n\left(S_{1,i}\left(t\right)+S_{2,i}\left(t\right)+C_i\left(t\right)+I_i\left(t\right)+R_i\left(t\right)\right).$$

(23)

By summing the differential equations for $S_{1,i},S_{2,i},C_i,I_i$ and $R_i$ over all $i=1,\dots ,n$, we obtain:

$${\displaystyle \frac{d{T}_H}{dt}}=\sum _{i=1}^n\left({\dot{S}}_{1,i}+{\dot{S}}_{2,i}+{\dot{C}}_i+{\dot{I}}_i+{\dot{R}}_i\right).$$

(24)

Observe that all inter-patch movement terms cancel out in the sum. For example, migration of $S_{1,i}$ from patch i to patch j appears with opposite signs in the equations for $S_{1,i}$ and $S_{1,j}$:

$$\sum _{i=1}^n\sum _{j=1}^n(m_{ji}^X{X}_j-m_{ij}^X{X}_i)=0.$$

(25)

Now, from the system (16), suppose

$${\mu }^H=\underset{i=1,\dots ,n}{min}{\mu }_i,$$

Thus,

$${\displaystyle \frac{d{T}_H}{dt}}\le {\mathsf{\Lambda }}^{*}-{\mu }^H{T}_H\left(t\right).$$

(26)

This implies that the total human population is non-increasing over time. And using Proposition 5.1, hence for all $t\ge 0$:

$$0\le T_H\left(t\right)\le max\left\{T_H\left(0\right),{\displaystyle \frac{{\mathsf{\Lambda }}^{*}}{{\mu }^H}}\right\}.$$

(27)

with

$${\mathsf{\Lambda }}^{*}=\sum _{i=1}^n\mathsf{\Lambda }.$$

Now consider the total bacterial concentration:

$$T_B\left(t\right)=\sum _{i=1}^n{B}_i\left(t\right).$$

(28)

From the bacterial dynamics:

$${\dot{B}}_i={\alpha }_C{C}_i+{\alpha }_I{I}_i-{\mu }_{B_i}\left(T\right)B_i+\sum _{j=1}^n(m_{ji}^B{B}_j-m_{ij}^B{B}_i),$$

(29)

the migration terms again cancel upon summing over i.

Thus:

$${\displaystyle \frac{d{T}_B}{dt}}=\sum _{i=1}^n({\alpha }_C{C}_i+{\alpha }_I{I}_i-{\mu }_{B_i}\left(T\right)B_i).$$

(30)

Let ${\alpha }_{max}=max({\alpha }_C,{\alpha }_I)$ (intuitively ${\alpha }_{max}={\alpha }_I$), $C_{max}\left(t\right)={max}_i{C}_i\left(t\right)$, and similarly for $I_i$ (because $C_i,I_i$ are bounded). Then:

$${\displaystyle \frac{d{T}_B}{dt}}\le {\alpha }_{max}\sum _{i=1}^n(C_i+I_i)-{\mu }_{min}{T}_B,$$

(31)

where ${\mu }_{min}={min}_{i,t}{\mu }_{B_i}\left(T\left(t\right)\right)>0$ (by hypothesis that $T\left(t\right)<T_{max}$ and temperature is bounded in real situations).

Let $K={\alpha }_{max}·{max}_t{T}_H\left(t\right)$. Then:

$${\displaystyle \frac{d{T}_B}{dt}}\le K-{\mu }_{min}{T}_B\left(t\right).$$

(32)

This is a linear differential inequality of the form:

$${\displaystyle \frac{d{T}_B}{dt}}+{\mu }_{min}{T}_B\le K.$$

(33)

Integrating, we find:

$$T_B\left(t\right)\le e^{-{\mu }_{min}t}\left(T_B\left(0\right)-{\displaystyle \frac{K}{{\mu }_{min}}}\right)+{\displaystyle \frac{K}{{\mu }_{min}}}.$$

(34)

Therefore, $T_B\left(t\right)$ is bounded for all $t\ge 0$, and in particular:

$$0\le T_B\left(t\right)\le max\left(T_B\left(0\right),{\displaystyle \frac{K}{{\mu }_{min}}}\right).$$

(35)

Combining both bounds, all human compartments and the bacteria concentration remain bounded over time. Since the other compartments are involved in equations with logistic or saturating growth, bounded by carrying capacities $K_1$ and food production rates $r_i$, they are also bounded.

Thus, all variables of system (16) are uniformly bounded, and the solution remains in a compact, positively invariant subset of ${\mathbb{R}}_{+}^{8n}$. □

Property 5.1. 

Let the initial conditions satisfy $X\left(0\right)\in \mathsf{\Omega }$, where

$$\mathsf{\Omega }=\left\{X\in {\mathbb{R}}_{+}^{8n}|\begin{array}{c}0\le T_H\left(t\right)\le M_H,\hfill \\ 0\le T_B\left(t\right)\le M_B,\hfill \\ 0\le F_i^s\left(t\right)+F_i^c\left(t\right)\le K_1,\forall i=1,\dots ,n\hfill \end{array}\right\}.$$

(36)

with,

$$M_H=max\left\{T_H\left(0\right),{\displaystyle \frac{{\mathsf{\Lambda }}^{*}}{{\mu }^H}}\right\}\mathrm{and}{M}_B=max\left(T_B\left(0\right),{\displaystyle \frac{K}{{\mu }_{min}}}\right).$$

Then, the region Ω is positively invariant under the flow of system (16). That is, any solution $X\left(t\right)$ with $X\left(0\right)\in \mathsf{\Omega }$ remains in Ω for all $t\ge 0$.

Proof. 

From Proposition 5.1, we know that all components of the solution remain non-negative, i.e., $X\left(t\right)\in {\mathbb{R}}_{+}^{8n}$ for all $t\ge 0$. Furthermore, from Property 5.1, the total human population satisfies $T_H\left(t\right)={\sum }_{i=1}^n{N}_i\left(t\right)\le M_H$, and the total bacterial concentration satisfies $T_B\left(t\right)={\sum }_{i=1}^n{B}_i\left(t\right)\le M_B$ for all $t\ge 0$.

Moreover, since the food compartment $F_i\left(t\right)=F_i^s\left(t\right)+F_i^c\left(t\right)$ evolves according to a logistic growth law with upper bound $K_1$, it follows that $F_i\left(t\right)\le K_1$ for each $i=1,\dots ,n$.

Hence, all conditions defining the set $\mathsf{\Omega }$ are preserved for all $t\ge 0$, which proves that $\mathsf{\Omega }$ is positively invariant under the flow of system (16). □

Theorem 5.1 

(Global existence and boundedness of solutions). Let the initial conditions $X\left(0\right)\in \mathsf{\Omega }\subset {\mathbb{R}}_{+}^{8n}$, and suppose that the right-hand side of system (16) is locally Lipschitz-continuous in X. Then, the system (16) admits a unique global solution $X\left(t\right)$ for all $t\ge 0$, which remains in the positively invariant set Ω.

Proof. 

The system (16) is a system of ordinary differential equations with locally Lipschitz-continuous right-hand side in ${\mathbb{R}}_{+}^{8n}$. Thus, by the Cauchy-Lipschitz theorem, a unique local solution $X\left(t\right)$ exists for some maximal time interval $[0,t_{max})$.

From Proposition 5.1, the solution remains non-negative for all $t\in [0,t_{max})$. Lemma 5.1 ensures that the total human population and bacterial concentration remain bounded, and Property 5.1 confirms that the compact positively invariant set $\mathsf{\Omega }$ is preserved under the flow of the system.

Therefore, the solution remains bounded and cannot blow up in finite time, implying that $t_{max}=+\infty$. Consequently, the local solution extends to a global solution. Since the solution stays in $\mathsf{\Omega }$ for all $t\ge 0$, the proof is complete. □

Remark 5.1. 

• We need movements within the same patch ($m_{ii}=0,\forall i$).
• In this paper, authors assimilate vulnerable individuals at individuals in food insecurity.
• Often, when we have patch model in case without movement ($m_{ji}^S=...=m_{ji}^B=0$), using (Arino and Driessche 2004) the basic reproduction number may be given by:

  $${\mathcal{R}}_0^{\mathrm{global}}=\underset{i=1,2,3,4}{max}{\mathcal{R}}_{0,i}.$$

  (37)

5.2. Disease-Free Equilibrium

At the disease-free equilibrium (DFE), the reduced system becomes:

$$\left\{\begin{array}{cc}{\dot{S}}_{1,i}& =(1-\pi )\mathsf{\Lambda }+{\xi }_1{R}_i+{\varphi }_2{S}_{2,i}-({\varphi }_1+\mu )S_{1,i}+\sum _{j=1}^n(m_{ji}^S{S}_{1,j}-m_{ij}^S{S}_{1,i}),\hfill \\ {\dot{S}}_{2,i}& =\pi \mathsf{\Lambda }+{\xi }_2{R}_i+{\varphi }_1{S}_{1,i}-({\varphi }_2+\mu )S_{2,i}+\sum _j(m_{ji}^S{S}_{2,j}-m_{ij}^S{S}_{2,i}),\hfill \\ {\dot{R}}_i& =-(\mu +{\xi }_1+{\xi }_2)R_i+\sum _j(m_{ji}^R{R}_j-m_{ij}^R{R}_i),\hfill \\ {\dot{F}}_i^s& =r_i{F}_i^s\left(1-{\displaystyle \frac{F_i^s}{K_1}}\right)-(a_i+\omega )F_i^s.\hfill \end{array}\right.$$

Let us define the following vectors:

$$S_1={(S_{1,1},\dots ,S_{1,n})}^{\top },S_2={(S_{2,1},\dots ,S_{2,n})}^{\top },R={(R_1,\dots ,R_n)}^{\top }$$

Let ${\mathcal{M}}^S$ and ${\mathcal{M}}^R$ be the movement matrices for susceptible and recovered individuals, respectively, and consider all the other parameters in matrix form such as

$${\varphi }_1=\mathrm{diag}\left({\varphi }_1^{\left(i\right)}\right),{\varphi }_2=\mathrm{diag}\left({\varphi }_2^{\left(i\right)}\right)$$

Property 5.2. 

${({\varphi }_1+{\varphi }_2+{\mu }^H-{\mathcal{M}}^S)}^{-1}$ is nonsingular and ${({\varphi }_1+{\varphi }_2+{\mu }^H-{\mathcal{M}}^S)}^{-1}>>0$.

Proof. 

The proof is guaranteed by (Arino et al. 2019), Proposition 3. □

The linearized system at the DFE becomes:

$$\left\{\begin{array}{c}{\displaystyle \frac{d}{dt}}\left(\begin{array}{c}{S}_1\\ S_2\end{array}\right)=\left(\begin{array}{c}(1-\pi )\mathsf{\Lambda }\\ \pi \mathsf{\Lambda }\end{array}\right)+\left(\begin{array}{cc}{\mathcal{M}}^S-({\varphi }_1+{\mu }^H)& {\varphi }_2\\ {\varphi }_1& -({\varphi }_2+{\mu }^H)+{\mathcal{M}}^S\end{array}\right)\left(\begin{array}{c}{S}_1\\ S_2\end{array}\right)=0\hfill \\ {\displaystyle \frac{dR}{dt}}=-\mu R+{\mathcal{M}}^{R}R\Rightarrow R^{*}\in ker\left({\mathcal{M}}^R\right)(\sin \mathrm{ce}{\mathcal{M}}^R\mathrm{is}\mathrm{irreducible})\hfill \end{array}\right.$$

Then, define:

$$A=\left(\begin{array}{cc}{\mathcal{M}}^S-({\varphi }_1+{\mu }^H)& {\varphi }_2\\ {\varphi }_1& -({\varphi }_2+{\mu }^H)+{\mathcal{M}}^S\end{array}\right),\mathrm{and}b=\left(\begin{array}{c}{b}_1\\ b_2\end{array}\right).$$

with:

$$b_1=\mathsf{\Lambda }(1-\pi ),b_2=\pi \mathsf{\Lambda }.$$

One has

$$S_1^{*}={\left(({\mathcal{M}}^S-{\mu }^H)({\mathcal{M}}^S-{\mu }^H-{\varphi }_1-{\varphi }_2)\right)}^{-1}((1-\pi )(-{\mathcal{M}}^S+{\mu }^H)+{\varphi }_2)\mathsf{\Lambda },$$

$$S_2^{*}={\left(({\mathcal{M}}^S-{\mu }^H)({\mathcal{M}}^S-{\mu }^H-{\varphi }_1-{\varphi }_2)\right)}^{-1}(\pi ({\mu }^H-{\mathcal{M}}^S)+{\varphi }_1)\mathsf{\Lambda }.$$

(38)

Then, the equilibrium is given by:

$${\mathcal{E}}_0^{(*)}=\left(S_1^{*},S_2^{*},\dots ,F_s^{*},0,0\right),$$

(39)

This expression defines the disease-free equilibrium distribution of susceptible individuals, accounting for nutritional transitions and mobility.

Consequence 5.1. 

Using the above result in Eq.39 without migration there exists two (02) diseases-free equilibria ${\mathcal{E}}_0^{\left(1\right)},{\mathcal{E}}_0^{\left(2\right)}\in {\mathbb{R}}_{+}^{8n}$ for system (16).

1.

A state of equilibrium with susceptible individuals but without food.

$${\mathcal{E}}_0^{\left(1\right)}=\left(S_1^0,S_2^0,\dots ,0\right),$$

2.

A state of equilibrium where there are individuals and food.

$${\mathcal{E}}_0^{\left(2\right)}=(S_1^0,S_2^0,0,0,0,F_s^0,0,0).$$

where $S_1^0={\displaystyle \frac{\mathsf{\Lambda }({\varphi }_2+\mu (1-\pi ))}{\mu (\mu +{\varphi }_1+{\varphi }_2)}},S_2^0={\displaystyle \frac{\mathsf{\Lambda }({\varphi }_1+\mu \pi )}{\mu (\mu +{\varphi }_1+{\varphi }_2)}}$ and $F_s^0=K\left(1-{\displaystyle \frac{a+\omega }{r_i}}\right)$.

5.3. Basic reproduction number

Let us consider the system (16) without inter-patch movements (i.e $m_{ji}^X=0$, for all variables X and all $i,j$). We find disease-free equilibria, and they occur when : $I=C=B=F^c=0$.

To compute this basic reproduction number, we use the method of Van den (van den Driessche and Watmough 2002). To do so, system (16) can be written in the following form:

$$\dot{X}=\mathcal{F}\left(X\right)-\mathcal{V}\left(X\right),$$

(40)

where $X=(C,I,F^s,B)$,

Then inside a satellite, we have:

$$\begin{array}{cccc}\hfill \mathcal{F}\left(X\right)& =\left(\begin{array}{c}(1-p)\left({\beta }_H·(I+\epsilon C)N\right)(S_2+\zeta S_1)\\ \\ p\left({\beta }_H·(I_i+\epsilon C)\right)(S_2+\zeta S_1)\\ 0\\ 0\end{array}\right)and\hfill & \hfill \mathcal{V}\left(X\right)& =\left(\begin{array}{c}({\gamma }_C+\mu )C\\ \\ ({\gamma }_I+d+\mu )I\\ -\rho B{F}^s+(a+\omega )F^c\\ -({\alpha }_{C}C+{\alpha }_{I}I)+{\mu }_B\end{array}\right).\hfill \end{array}$$

(41)

We point that $\mathcal{F}\left(X\right)$ is the speed of appearance of newly infected humans. These are the newly infected obtained by transmissions of all kinds. and $\mathcal{V}\left(X\right)$ is the rate of onset of new cases for reasons other than disease.

The Jacobian matrix of $\mathcal{F}\left(X\right)$ and $\mathcal{V}\left(X\right)$ at the disease free equilibrium point ${\mathcal{E}}_0^{\left(4\right)}$ are:

$${\displaystyle \frac{d\mathcal{F}}{dy}}=\left(\begin{array}{cccc}(1-p){\beta }_H\epsilon (S_2^{*}+\zeta S_1^{*})& (1-p){\beta }_H(S_2^{*}+\zeta S_1^{*})& (1-p){\beta }_B(S_2^{*}+\zeta S_1^{*})& 0\\ \\ p{\beta }_H\epsilon (S_2^{*}+\zeta S_1^{*})& p{\beta }_H(S_2^{*}+\zeta S_1^{*})& p{\beta }_B(S_2^{*}+\zeta S_1^{*})& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\end{array}\right),$$

(42)

and

$${{\displaystyle \frac{d\mathcal{V}}{dy}}{|}_{X={\mathcal{E}}_0^{\left(2\right)}}}^{-1}=\left(\begin{array}{cccc}{\displaystyle \frac{1}{{\gamma }_C+\mu }}& 0& 0& 0\\ 0& {\displaystyle \frac{1}{d+{\gamma }_I+\mu }}& 0& 0\\ {\displaystyle \frac{F_s^0{\alpha }_C\rho }{{\gamma }_C{\mu }_B\left(a+\omega \right)}}& {\displaystyle \frac{F_s^0{\alpha }_I\rho }{{\mu }_B\left(ad+a{\gamma }_I+d\omega +{\gamma }_I\omega \right)}}& {\displaystyle \frac{1}{a+\omega }}& {\displaystyle \frac{F_s^0\rho }{{\mu }_B\left(a+\omega \right)}}\\ {\displaystyle \frac{{\alpha }_C}{{\gamma }_C{\mu }_B}}& {\displaystyle \frac{{\alpha }_I}{{\mu }_B\left(d+{\gamma }_I\right)}}& 0& {\displaystyle \frac{1}{{\mu }_B}}\end{array}\right),$$

(43)

Then,

$${\displaystyle \mathcal{F}{\mathcal{V}}^{-1}=j\left(\begin{array}{cccc}\frac{m\left(F_s^0{\alpha }_C{\beta }_B\rho +{\beta }_Hϵ{\mu }_B\left(a+\omega \right)\right)}{{\gamma }_C{\mu }_B\left(a+\omega \right)}& \frac{m\left(F_s^0{\alpha }_I{\beta }_B\rho +a{\beta }_H{\mu }_B+{\beta }_H{\mu }_B\omega \right)}{{\mu }_B\left(ad+a{\gamma }_I+d\omega +{\gamma }_I\omega \right)}& \frac{{\beta }_{B}m}{a+\omega }& \frac{F_s^0{\beta }_{B}m\rho }{{\mu }_B\left(a+\omega \right)}\\ \frac{p\left(F_s^0{\alpha }_C{\beta }_B\rho +{\beta }_Hϵ{\mu }_B\left(a+\omega \right)\right)}{{\gamma }_C{\mu }_B\left(a+\omega \right)}& \frac{p\left(F_s^0{\alpha }_I{\beta }_B\rho +a{\beta }_H{\mu }_B+{\beta }_H{\mu }_B\omega \right)}{{\mu }_B\left(ad+a{\gamma }_I+d\omega +{\gamma }_I\omega \right)}& \frac{{\beta }_{B}p}{a+\omega }& \frac{F_s^0{\beta }_{B}p\rho }{{\mu }_B\left(a+\omega \right)}\\ 0& 0& 0& 0\\ 0& 0& 0& 0\end{array}\right).}$$

(44)

with

$$j=S_2^{*}+\zeta S_1^{*},\mathrm{and}m=1-p.$$

One has,

$$Sp\left(\mathcal{F}{\mathcal{V}}^{-1}\right)=\left\{0,R_0\right\}.$$

Using Theorem 2, (van den Driessche and Watmough 2002)

$${\mathcal{R}}_0=\rho \left(\mathcal{F}{\mathcal{V}}^{-1}\right)$$

Yields

$${\mathcal{R}}_0=(S_2^0+\zeta S_1^0)({\mathcal{R}}_0^H+{\mathcal{R}}_0^E).$$

(45)

$${\mathcal{R}}_0^H={\beta }_H\left({\displaystyle \frac{(1-p)\epsilon }{{\gamma }_C+\mu }}+{\displaystyle \frac{p}{{\gamma }_I+d(.)+\mu }}\right)\mathrm{and}{\mathcal{R}}_0^E={\beta }_B{\displaystyle \frac{F_s^0\rho }{{\mu }_B(a+\omega )}}\left({\displaystyle \frac{(1-p){\alpha }_C}{{\gamma }_C}}+{\displaystyle \frac{p{\alpha }_I}{{\gamma }_I+d}}\right).$$

(46)

Where ${\mathcal{R}}_0^H$ represents the transmission of the disease via human interactions.

Therefore, if a patch i, is isolated from the others i.e $m_{ji}^X=0$, for all variables X and all $i,j$, then the basic reproduction number in each patch is given by:

$${\mathcal{R}}_0=\left\{\begin{array}{cc}{\mathcal{R}}_0^{\left(1\right)}={\displaystyle \frac{S_2^0+\zeta S_1^0}{N_c}}({\mathcal{R}}_0^H+{\displaystyle \frac{{\mathcal{R}}_0^E}{\chi }})\hfill & ,\mathrm{for}\mathrm{the}\mathrm{principal}\mathrm{node}\hfill \\ {\mathcal{R}}_0^{\left(i\right)}=(S_2^0+\zeta S_1^0)({\mathcal{R}}_0^H+{\mathcal{R}}_0^E)\hfill & ,\mathrm{for}\mathrm{their}\mathrm{satellites}\hfill \end{array}\right.$$

(47)

with $N_c$ the total population inside the central node.

Theorem 5.2. 

Suppose that $m_{ij}=m_{ji}=m$ for all i,j=1,...,n. Then

$$\underset{i=1,\dots ,n}{min}{\mathcal{R}}_0^i\le {\mathcal{R}}_0\le \underset{i=1,\dots ,n}{max}{\mathcal{R}}_0^i$$

(48)

Proof. 

See (Arino and Driessche 2004). □

6. Model Application

This section presents numerical simulations that reflect local transmission patterns, spatial connectivity, and food-related vulnerability.

6.1. Numerical simulations

The initial conditions reflect a scenario where most individuals are susceptible, with a few infectious cases in the central node and no contamination in peripheral patches. However, to assess the importance of dynamically modeling nutritional vulnerability, we simulated the model under varying transition rate ${\varphi }_1$ between well-nourished and vulnerable individuals. Figure 6 shows that epidemic size and mortality are highly sensitive to these parameters. In particular, increasing the rate ${\varphi }_1$ at which individuals become vulnerable significantly amplifies both infection and death peaks. This justifies the inclusion of two susceptible compartments ($S_1$, $S_2$), as a single-compartment model with a static risk modifier would fail to capture these emergent dynamics.

Moreover, it confirms that nutritional transitions drive epidemic severity and that one class susceptible model would not capture these nonlinear effects.

Figure 7 illustrates periodic dynamics in the population, which may indicate the presence of a limit cycle in the model.

6.2. Case study of Douala and its Surroundings Areas

The modeling processes start inside two regions (SU, SW) and finish at (Littoral and South-west regions) cover an area of more than 100000 km². However, we focus the study on the city of Douala and its neighbors, which recorded one of the highest cholera case-fatality rates in Cameroon during the 2022 outbreak ($\approx 3\%$), and lies adjacent to the South-West region, the most affected area that epidemic. In a context marked by sociopolitical instability and high population displacement from conflict-affected “anglophone” regions, Douala serves as the main entry point for migrants from the South-West into “francophone” urban zones. This demographic pressure, combined with preexisting challenges related to water, sanitation, and food safety, creates conditions conducive to cholera transmission and persistence.

To capture both local dynamics and spatial interactions, we consider a multi-patch framework composed of four interconnected zones:

• a central node representing Douala,
• and three surrounding patches: Bonaberi, Bomono, and Yassa connected to Douala through human movement and food exchange.
• ${\mu }_B\left(T\right)=\widetilde{{\mu }_B}\left(1-\kappa {\displaystyle \frac{T-\overline{T}}{T_{max}-\overline{T}}}\right)$, where the temperature is in degree Celsius. $\kappa$ represents the dependency on temperature, ($\kappa =3$ according to (Bertuzzo et al. 2010)), $T_{max}$ and $\overline{T}$ correspond respectively to the maximum and mean temperature of Douala city over the 20 years.

This structure reflects the urban-peripheral gradient and allows us to account for spatial heterogeneity in exposure, vulnerability, and mobility, essential to understand the spread and control of cholera in the Douala metropolitan area.

7. Discussions and Perspectives

In this study, we developed a novel mathematical model to understand cholera transmission in the context of environmental contamination and food insecurity. The model progressively incorporates the key pathways of cholera spread from basic SIR dynamics, to waterborne transmission (SIWR), and finally to food/biomass mediated transmission via the proposed SIWR-F (2) framework.

Moreover, the multi-patch extension (16) captures the geographic and demographic heterogeneity of cholera dynamics, emphasizing the role of migration, patch connectivity, and spatially varying food availability. Simulations based on the Douala metropolitan area illustrate how vulnerable populations in connected patches can act as persistent sources or sinks for the disease.

In conclusion, this work provides a comprehensive framework to analyze cholera outbreaks in resource-limited settings, showing that food insecurity is not just a background condition but a core driver of epidemic dynamics. The model can inform targeted interventions that integrate food, water, and mobility data paying the way toward more resilient health systems in vulnerable regions.

As possible extensions of this work, one can assume without loss of generality to extend this model (16) over several years (given that cholera is a disease in several countries of sub-Saharan Africa) by change with parameters time dependent.