Mathematical modelling of the dynamics of typhoid fever and two modes of treatment in a Health District in Cameroon

1.   Introduction

Typhoid fever (TF) is one of the most common infectious diseases in South Asia and Sub-Saharan Africa, with over 9 million cases per year ^[1,2]. It is caused by the bacterium Salmonella typhi, which is usually spread through contaminated food or water. Persons with typhoid fever carry the bacteria in their bloodstream and intestinal tract. The development of the disease follows three phases: the incubation phase, the invasion phase (acute), and the chronic phase. The incubation period varies between 7 and 21 days ^[3]. The acute phase is characterized by the gradual onset of fever, headaches, dizziness, insomnia, epistaxis, and sometimes constipation ^[4]. In the chronic phase, the infected may have encephalitis or inflammation of the brain, dehydration, weakness, abdominal hemorrhage due to severe intestinal perforations, and other serious complications that may lead to death^[2].

In Sub-Saharan Africa, it is estimated that more than $80\%$ of the population use medicinal plants as their main source of treatment ^[5,6]. Antibiotics have been the primary therapy for typhoid fever and also the cornerstone of modern medicine. However, the increased resistance of some bacteria to antibiotics may be one of the reasons why patients seek traditional treatment. One principal reason typhoid fever patients may prefer traditional medicine over modern medicine is the nonexistence of modern healthcare facilities, especially in rural areas. In areas where these facilities exist, the antibiotics are expensive, as the patients are mainly low-income earners. In contrast, traditional medicine is relatively cheap and readily available. The effectiveness of traditional medicine is recognized in several countries and has made important contributions to modern medicine ^[5,7]. A report in ^[8] says, WHO has supported clinical trials, leading 14 countries to issue marketing authorization for 89 traditional medicine products that met international and national requirements for registration. It is important to note that since the 1980s, traditional medicine has been benefitting from such authorizations ^[9,10].

Sometimes, some of the patients taking antibiotics experience a relapse of typhoid fever after an initial recovery. Relapse is the return of a disease or the signs and symptoms of a disease after a period of improvement or treatment. This phenomenon may be due to incomplete treatment or from antibiotics resistance due to the abuse of prescribed antibiotics. The re-emergence of a disease may either be due to relapse or reinfection. Reinfection occurs when a patient after; treatment becomes infected again, meanwhile, relapse is the recurrence of the same infection. Relapse dynamics models account for the possibility that individuals who have been treated for typhoid fever may have a reactivation of the same Salmonella typhi bacteria.

Models incorporating treatment relapse are particularly important concerning TF. Many mathematical models have been developed to study the spread of typhoid fever ^{[12,13,14,15,16,17,18,19]}. The authors of ^[13], proposed a human-to-human (direct) and environment-to-human (indirect) transmission model, in which susceptible individuals get infected with Salmonella typhi at a rate proportional to the susceptible population and the environmental bacteria concentration at a constant rate. They proved that sanitation, vaccination, and treatment of symptomatic and asymptomatic infected individuals are efficient measures, to reduce the severity of the disease. Considering the latter control measures, ^[12] presented a model to investigate the outcome of the disease in Mbandjock, a town in Cameroon. Through calibration, they obtained a control reproduction number of $2.4750$, meaning that TF was endemic in that region. Similar to ^[13], the authors of ^[18] proposed a model with both direct and indirect transmissions, but also considered in addition the influence of vaccination. Through a sensitivity analysis, they found that the proportion of susceptible unvaccinated immigrants had a high impact on the infected populations. Thus, people entering the town were to be checked, to ensure they were vaccinated against TF. In ^[11], a non-autonomous model was proposed accounting for temporal variations and it was shown that the human-to-human infection rate had a significant impact on the reproduction number, while the environment-to-human infection rate and the bacteria shedding rate had an effect on long-cycle infections.

Mushayabasa ^[14], proposed an $S\, E\, I\, I_c\, Q\, R$ compartmental model, with the susceptible individuals $(S)$, exposed individuals $(E)$, symptomatic infected individuals $(I)$, chronic enteric carriers $(I_c)$, quarantined symptomatic persons and chronic enteric carriers $(Q)$, and recovered individuals $(R)$. He assumed that the population reduced due to natural death $(\mu)$ and the acquired infection rate $(\lambda)$ expressed as a function of $I$ and $I_c$. The disease-free equilibrium was shown to be globally asymptotically stable whenever the reproductive number was less than unity. The transmission of the disease was through human contact only and he concluded that in the event of an outbreak of typhoid in the community the disease can be effectively controlled if optimal intervention strategies are implemented. In ^[19], a mathematical model that explored the dynamics of typhoid fever transmission with particular focus on the impact of treatment relapse is presented. The study formulated a deterministic mathematical model to analyze both direct and indirect transmission modes of typhoid infection. The findings reveal that limited efficacy of antibiotics and relapse response significantly influence the spread of typhoid infection.

Despite the avalanche of studies on the dynamics of TF, there is little literature on the role of traditional medicine in the treatment of TF. There are some studies on typhoid fever-treating herbs ^[7], which demonstrated that extracts from these plants have an inhibiting impact on Salmonella typhi. Herbal medicine is a feature and at times synonymous with African traditional medicine. The role of traditional medicine in the treatment of COVID-19 has been discussed in ^[16]. The results prescribed the combination of modern medicine with traditional medicine instead of one medicine as a stand-alone treatment. The mathematical model in the study did not consider possible relapse. However, the COVID-19 virus and Salmonella typhi are diseases with different modes of transmission and hence their models are dissimilar.

To the best of our knowledge, no typhoid fever model that incorporates modern treatment and traditional treatment has been proposed in literature so far. This paper contributes to the existing literature by proposing a novel mathematical deterministic model that incorporates these modes of treatment and treatment relapse due to resistance to antibiotics. The structure of the model model described and formulated in Section 2, to an extent, is a modified version of Mushayabasa's in ^[14]. Essentially, we have replaced the quarantine compartment ($Q$) with two compartments for treatment modes (modern and traditional). In addition, we assume that the bacterial infection is indirect (environment-to-human). Together with real data, we study and make predictions about the future dynamics of the disease in the locality of Penka-Michel in Cameroon. We theoretically analyze our model and use the fminsearchbnd and ode45 functions in MATLAB for calibration and numerical simulations, respectively.

The rest of the paper is organised as follows. In Section 2, the main assumptions of the model are presented with a detailed description of the model. Section 3 presents the quantitative and qualitative analysis of the model. In Section 4, we calibrate and validate the model. Section 5 deals with the sensitivity analysis and numerical simulations. The conclusion of the paper is provided in Section 6.

2.   Model formulation

2.1.   Main assumptions and variables

In this study, the following assumptions are made for our model. These assumptions are crucial in determining the effectiveness and tractability of the model.

Assumption 1. Infected individuals can choose either modern or traditional treatment. However, individuals under modern treatment can switch to traditional therapy. This assumption is motivated by the resistance of Salmonella to some modern treatment protocols ^[20].

Assumption 2. Although human-to-human transmission of typhoid is possible, it is less probable ^[21]. Therefore, only indirect means of transmission are considered, i.e., transmission through contaminated food, poor sanitation and unsafe drinking water.

Assumption 3. Individuals in the chronic phase of the disease can get serious complications like internal bleeding and encephalitis (brain inflammation) ^[22]. We assume that this category of patients will choose modern treatment only.

Assumption 4. It is known that after recovering from typhoid fever, individuals do develop some level of immunity. Recovery cases develop antibodies and remain immune against the disease for at least three years ^[23]. Thus, for simplicity, we assume that the recovered individuals enjoy permanent immunity.

Assumption 5. We consider a possible relapse of typhoid fever for individuals receiving modern treatment only. In our case, a relapse occurs when patients taking antibiotics fall back to the acute infection stage, after a period of improvement.

Assumption 6. The shedding rate of individuals under treatment is negligible. That is, we have well-informed patients who are taking all the precautionary measures to avoid spreading the disease.

Assumption 7. Under treatment, the mortality rate due to typhoid is less than 1% ^[24]. We consider this value negligible and ignore the mortality rate for patients receiving any of the treatments.

Assumption 8. The individuals in each group have an equal natural death rate.

The epidemiological model we develop in this study has eight compartments for the transmission dynamics of the disease. We denote by $N(t)$ the human population at time $t$, which has been divided into mutually exclusive compartments as follows:

● $S: = S(t)$: Susceptible individuals.

● $E: = E(t)$: Exposed individuals. This refers to the individuals infected with Salmonella but in a latent period.

● $I_a: = I_a(t)$: Infected individuals in the acute phase of the disease.

● $I_c: = I_c(t)$: Infected individuals in the chronic phase of the infection. The symptoms of these individuals develop very slowly, so the disease persists for a long time in their bodies.

● $M_m: = M_m(t)$: Individuals under modern treatment.

● $M_t: = M_t(t)$: Individuals under traditional treatment.

● $R: = R(t)$: Recovered individuals.

● $B: = B(t)$: The concentration of Salmonella bacteria in the environment.

$N(t) = S(t)+ E(t)+ I_a(t)+ I_c(t)+ M_m(t)+ M_t(t)+ R(t),$ is the total human population at time $t$.

2.2.   The model equations

In this section, we derive the necessary equations for the construction of our mathematical model while also defining the inherent variables.

The susceptible individuals are recruited at a constant number $\pi$, through birth or immigration. The susceptible population diminishes either due to natural death at the rate $\mu$, or due to Salmonella bacteria infection. We assume that the infection pressure $\lambda$, the rate at which susceptible individuals acquire the infection after their exposure to the disease, is given by the Holling Type Ⅱ functional response (see ^[21]). That is

where $\beta$ is the ingestion rate (that is, the product of the contact rate between individuals and the environment, and the percentage of Salmonella successfully ingested). The quantity $k_b$ is the concentration of bacteria in what is being consumed which gives 50% chance of getting infected. The fraction $B/(B+k_b)$, measures the probability of individuals getting infected through contaminated food or water. It is clear that $\lambda$ is an increasing function of the concentration of Salmonella (B) and when $B$ is large enough ($B\gg k_b$), $\lambda$ saturates at the constant value $\beta$. If $B\ll k_b, \, \lambda$ will grow linearly with $B$. Thus, the dynamic equation of the susceptible compartment is depicted by:

Individuals in the exposed compartment exit this class either through death due to natural causes at the rate $\mu$, or enter the acute infection compartment at the rate $\alpha$. This leads to the dynamic equation of compartment $E$ as:

The acute infected compartment $(I_a)$ is increased by movements of individuals from the exposed group at the rate $\alpha$, or from the modern treatment class as a result of treatment failure (antibiotics resistance), at the rate $\theta_m$. Individuals leaving compartment $(I_a)$ enters either compartment $(M_m)$ at the rate $\gamma_{am},$ the chronic infection stage $(I_c)$ at the rate $\kappa$, the traditional medicine treatment class $(M_t)$ at the rate $\gamma_{at},$ or these infected individuals die due to the disease at the rate $\delta_a$ or die due to natural causes at the rate $\mu$. The dynamic equation of $I_a$ is given by:

The infected individuals in the chronic stage compartment are replenished by the $\kappa I_a$ individuals who left compartment $I_a$. The number of individuals in $I_c$ is reduced by choosing the modern treatment, at the rate $\omega_{cm},$ by natural mortality or by disease-induced mortality, at the rate $\delta_c.$ In compartment $I_c$, the dynamic equation is:

The acute infected individuals move into the traditional treatment compartment at the rate $\gamma_{at}$, while those receiving modern treatment switch to this compartment at the rate $\psi_{mt}$. The population of $M_t$ decreases either by recovery at a rate $\sigma_t,$ or by natural mortality at a rate $\mu$. The dynamic equation of $M_t$ is thus given as:

Compartment $M_m$ is supplied by compartment $I_a$ and $I_c$ at the rates $\gamma_{am}$ and $\omega_{cm}$, respectively. This population $M_m$ decreases in four different ways: Switching to traditional treatment at a rate $\psi_{mt}$, having a relapse at the rate $\theta_m$, recovering at a rate $\sigma_m$, and death due to natural mortality at the rate $\mu$. This leads to the following dynamic equation of $M_m$:

Patients in compartments $M_m$ and $M_t$ die due natural causes at the rate $\mu$ and recover at the rates $\sigma_m$ and $\sigma_t$, respectively. This population decreases by natural mortality. The dynamic equation of $R$ is:

We assume that Salmonella bacteria have a natural growth rate $r$ in the environment. Infected individuals in the compartments $I_a$ and $I_c$ shed bacteria in the environment at the rates $\eta_a$ and $\eta_c$, respectively. The bacterial decay rate in the environment is $\mu_b$. Thus, the dynamic equation of the concentration of bacteria is given as:

Table 1 presents the descriptions of all the parameters of interest in our model. The Flow diagram of the model is given in Figure 1. Putting all the equations together give System (2.6) of non-linear differential equations, which is the main model of this paper.

subject to the initial conditions

For mathematical convenience, we set

Throughout this paper, we assume $\mu_b > r$ in order for $k_6$ to be positive. The biological significance of this assumption is that without shedding of Salmonella typhi in the environment, the concentration of the bacteria will decrease exponentially ^[25].

3.   Theoretical analysis of the model

3.1.   Well-posedness of the model

Model (2.6) is about a human population and a bacterial population. It is important to check if this model has a solution and if the variables remain positive for a positive initial condition. This is ensured by the following result.

Theorem 3.1. Model (2.6) is a dynamic system on the biologically feasible compact domain:

Proof. According to Cauchy-Lipschitz theorem, Model (2.6) has a unique local solution as its right-hand side is locally Lipschitz.

To show that $S(t) \geq 0$ for all $t \geq 0$, we rewrite the first equation of ($2.6$) as follows:

Integrating the above equation from $0$ to $t$ gives

This leads to

Thus, if $S(0) > 0$, then $S(t) > 0$, for all $t\geq0$.

To show that $E(t), \; I_a(t), \; I_c(t), \; M_m(t), \; M_t(t), \; B(t)$ and $R(t)$ are positive when $E(0)\geq 0, \; I_a(0) \geq 0, \; I_c(0)\geq0, \; M_m(0)\geq 0, \; M_t(0)\geq 0, \; R(0)\geq 0, \; B(0)\geq 0,$ we consider the following sub-system

This can be rewritten as

where

Clearly, $A$ is a Metzler matrix. Hence System (3.2) is positively invariant in $\mathbb{R}_+^7.$

Furthermore, by adding the first seven equations of System (2.6), we have:

The application of Gronwall's inequality gives

Thus,

Finally, under the hypothesis that $N(0)\leq\pi/\mu$ and knowing that $I_a\leq N$, $I_c\leq N$, $M_m\leq N$, and $M_t\leq N$, the application of Gronwall's inequality once more leads to

Using the fact that the solutions of (2.6) are bounded on $\mathbb{R}_8^+,$ we then conclude that with a non-negative initial condition, the solution of System (2.6) remains non-negative and exists globally over time.

3.2.   Disease-free equilibrium and its stability

The unique DFE of System (2.6) is given by:

We use the next-generation matrix approach ^[26,27], to compute $\mathcal{R}_0$. According to ^[27], the vector of the new infections $\mathcal{F}$ and that of the remaining transfer terms $\mathcal{V}$ are respectively given by:

The next-generation matrix $FV^{-1}$, where $F$ and $V$ are the Jacobian matrices of $\mathcal{F}$ and $\mathcal{V}$ at the DFE, respectively, is given by:

where the expressions for $A_{12}, \; A_{13}, \; A_{14}$ and $A_{16}$ are given in Appendix A.

Hence

The relevance of the reproduction number $\mathcal{R}_0$ is established in the following result ^[27].

Lemma 3.1. The DFE $P_0$ of System (2.6) is locally asymptotically stable (LAS) if $\mathcal{R}_0 < 1$ and it is unstable if $\mathcal{R}_0 > 1$.

According to Lemma 3.1, the disease dies out when if $\mathcal{R}_0 < 1$ provided that the initial population sizes are in the basin of attraction of the DFE. For the global control of the disease, the global asymptotic stability needs to be proven.

Theorem 3.2. Assume $\mathcal{R}_0 < 1$, then the $DFE$ point $P_0$ of Model (2.6) is globally asymptotically stable (GAS) in $\varOmega$.

Proof. Let us consider the candidate Lyapunov function

where

Differentiating $V$ on both sides gives

That is,

The resulting equation is due to the fact that

and

Following (3.6), Equation (3.9) becomes

Since $S\leq N\leq\dfrac{\pi}{\mu} = S_0,$ and $\dfrac{1}{B+k_b} \leq \dfrac{1}{k_b}$, we have:

Thus, whenever $\mathcal{R}_0 < 1$, we have $dV/dt\leq 0$, and $V$ is a strict Lyapunov function in $\varOmega$. Hence, the disease-free equilibrium $P_0$ is globally asymptotically stable.

Figure 2 illustrates the GAS of the DFE for Model (2.6) when $\mathcal{R}_0 < 1$.

3.3.   Interior equilibrium of the model and its stability

Theorem 3.3. Model (2.6) has a unique positive equilibrium if and only if $\mathcal{R}_0 > 1$.

Proof. Let us denote by $\mathcal{\varepsilon}^{\ast} = (S^*, E^*, I_a^*, I_c^*, M_m^*, M_t^*, R^*, B^*)$ a nontrivial equilibrium for Model (2.6). Setting the right-hand side of Model (2.6) to zero gives:

with

Substituting the expressions in (3.12) in Eq (3.13) gives

That is

Hence, when $\mathcal{R}_0\leq 1$, the unique equilibrium is the DFE, while when $\mathcal{R}_0 > 1$, the DFE $P_0$ coexists with the unique endemic equilibrium $\mathcal{\varepsilon}^*$. □

We study the local stability of $\mathcal{\varepsilon}^{\ast}$ in Theorem 3.4.

Theorem 3.4. System (2.6) has a trans-critical bifurcation at $\mathcal{R}_0 = 1$, which is the bifurcation parameter. Moreover, the unique endemic equilibrium $\mathcal{\varepsilon}^{\ast}$ is LAS when $\mathcal{R}_0 > 1$.

Proof. Let us consider the case where $\mathcal{R}_0 = 1$, and choose $\beta = \beta^{\ast}$ as a bifurcation parameter.

The Jacobian matrix of System (2.6) at the DFE is given as,

System (2.6), with $\beta = \beta^*$ has a nonhyperbolic equilibrium point. The other eigenvalues have negative real parts. Therefore, the center manifold theory ^[28] can be applied to analyze the dynamics of System (2.6) near the bifurcation parameter $\beta^*$. The components of a right eigenvector $w = (w_1, w_2, w_3, w_4, w_5, w_6, w_7, w_8)^T$ of $J_{\beta^*}$ and a non-negative left-eigenvector $v = (v_1, v_2, v_3, v_4, v_5, v_6, v_7, v_8)^T$ of $J_{\beta^*}$ associated with the zero eigenvalue are, respectively, given as

and

The coefficients $a$ and $b$ as defined in ^[28] are:

and

Since $a < 0$ and $b > 0$, the model has a trans-critical bifurcation at $\mathcal{R}_0 = 1$ and the endemic equilibrium is LAS for $\mathcal{R}_0 > 1$, but close to $1$ ^[28].

Figure 3 illustrates the forward bifurcation for Model (2.6).

Theorem 3.4 proves the persistence of the disease when $\mathcal{R}_0 > 1$, for values close to $1$. To determine the asymptotic behavior of the model for higher values of $\mathcal{R}_0$, the global asymptotic stability of $\mathcal{\varepsilon}^{\ast}$ has to be proven. This is stated in Theorem 3.5.

Theorem 3.5. The endemic equilibrium $\varepsilon^{\ast}$ for Model $2.6$ is GAS provided that $\mathcal{R}_0 > 1$.

Proof. See Appendix B.

Figure 4 illustrates the persistence of typhoid and the global asymptotic stability of $\mathcal{\varepsilon}^{\ast}$ when $\mathcal{R}_0 > 1$ for trajectory plot when using the values of Table $2$ and $\mathcal{R}_0 > 1.$ From this figure, we can observe that the infected individuals and bacteria are always present in the population. This means that the trajectories converge to the endemic equilibrium point. Thus, whenever $\mathcal{R}_0 > 1,$ the disease persists in the host population as established in Theorem 3.3.

4.   Model calibration and sensitivity analysis

4.1.   Model calibration

In this section, we fitted our model to the weekly cumulative reported cases of TF in the Penka-Michel health district in Cameroon from 30 December 2019 to $1$ January 2023 (157 weeks). The number of infected cases during the first week was $48$. We used the following initial conditions for the human population: $S(0)$ = 120,000, $E(0)$ = 0, $I_a(0)$ = 10, and $I_c(0)$ = 20, $M_m(0)$ = 10, $M_t(0)$ = 8, $R(0) = 0$. The value $B(0)$ = 7371 was obtained by calibration as a parameter.

According to ^[29], the life expectancy in Cameroon is 54.4 years. Therefore, the estimated value for $\mu$ is $0.000353$ per week. For simplicity, an estimate for $\pi$ is $\pi = S_0\times \mu = 42.36$. The data were fitted using the nonlinear least squares algorithm implemented by the fminsearcbnd function in MATLAB. Plots of the model's fit and calibration are shown in Figure 5.

Figure 5a presents model fit to the cumulative cases of diagnosed typhoid fever, while Figure 5b shows the model's validation and prediction on future cases based on our model. The fitted parameters as well as those obtained from literature are presented in Table 2. The results in Figure 5a, show that our model is a very good fit for the typhoid disease dynamics in Penka-Michel. To validate our model, in Figure 5a, we plot the curve from 2 January 2023 to 5 November 2023 (from Week 158 to Week 201). The model gives an estimate of 12,950 cases of typhoid at the end of Week 201. We observe that this number is quite close to the 12,873 reported cases at that same period. Hence, the model can be used for predictions on the disease trend in the district. Extrapolation of the trend curve up to 6 November 2026 (Week 358), predicts a total number of 21,270 cases of typhoid before the end of 2026. This knowledge, together with a basic reproduction number $\mathcal{R}_0 = 1.2058 > 1$, calls for the urgent need to put control measures in place to overcome the present trend of the disease. However, a knowledge on the most influential parameters in the model is crucial for an effective control strategy.

4.2.   Sensitivity analysis

Sensitivity analysis has almost becoming an integral part of mathematical modeling. In infectious disease modeling, it is the technique commonly used to determine which parameters have a significant effect on the spread of an infectious disease. We have shown that $\mathcal{R}_0$ verifies the sharp threshold property, and hence, the control of the disease lies in the control of $\mathcal{R}_0$. We performed global sensitivity analysis on $\mathcal{R}_0$ to identify its most influential parameters. The sensitivity analysis actually assesses the weight and type of change on the basic reproduction number when the model parameters vary across the entire range of values as shown in Table 2. We used the Latin hypercube sampling (LHS) technique to run 2500 simulations and we compute thr partial rank correlation coefficients (PRCC) between $\mathcal{R}_0$ and each parameter of Model (2.6). Usually, parameters with high PRCC absolute values (> 0.5 or < -0.5) and with small $p$-values (<0.05) are considered most influential on $\mathcal{R}_0$ ^[32].

The PRCC computed values are displayed in Figure 6. This figure shows that $\mu_b, \, \gamma_{am}$ and $\beta$ are the most influential parameters in either decreasing or increasing the value of $\mathcal{R}_0$. The analysis shows that shedding rate of Salmonella typhi in the environment increases the burden of the disease. The modern medicine parameter is seen to be crucial in decreasing the value of $\mathcal{R}_0$. Other parameters that significantly reduce $\mathcal{R}_0$ whenever they increase are $\mu$ and $\delta_c$. Finally, though $\gamma_{at}$ and $\psi_{mt}$ are not among the most influential parameters, the results point out that taking either of the treatments or switching from modern treatment to traditional therapy will decrease the disease burden. This highlights the important role traditional medicine should play in the control of typhoid fever.

Using PRCC and scatterplots together provides a robust way to understand and visualize the relationships (monotonicity) between the model's parameters and outputs. A monotonic relationship means that the output consistently increases or decreases as the parameter changes. PRCC measures the strength and direction of the monotonic relationship between a model parameter and the model output while controlling for other parameters. Figure 7 depicts scatterplots for each parameter against the model output $\mathcal{R}_0$. The plots show the monotonicity between each parameter of our model and $\mathcal{R}_0$.

5.   Numerical simulations

In this section, we provide some simulations to illustrate how a change in the parameters influences the model's dynamics. Most of the parameters used are presented in Table 2. In Figure 8, we explore the impact of relapse over time by considering varying values of $\theta_m$: $\theta_m = 0.2$; $\theta_m = 0.4$, and $\theta_m = 0.6$. It is seen that as the value of $\theta_m$ increases, typhoid cases also increase. This is an indication that resistance to modern treatment leads to a relapse into the disease stage.

One way of avoiding the relapse of TF is to switch from modern treatment to traditional therapy. Figure 9, gives a picture of the dynamics as infected individuals receiving modern treatment switch to traditional treatment at the rates $\psi_{mt} = 0.2$; $\psi_{mt} = 0.4$, and $\psi_{mt} = 0.6$. The output shows that over time, there is a decrease in the number of cases as the rate of movement from the modern treatment class to the traditional medicine compartment increases. This is a pointer to the positive impact traditional medicine has on lessening the burden of the typhoid fever disease.

Though it is recommended to switch from modern treatment to traditional treatment, it is imperative to at least receive some treatment, be it modern or traditional. Figure 10 presents different scenarios that emphasize the importance of taking treatment against typhoid fever. The different subfigures (10a–10f) affirm that receiving any of the two modes of treatment will cause a decrease in the number of cases.

In Figure 11, we investigate the impact of Salmonella bacteria shedding in the environment. It is seen that when the concentration of the bacteria increases in the environment, the number of cases increases as well. This calls for the necessity of adopting measures to prevent the spread of the disease, such as good sanitation habits, eating safe food and drinking potable water.

The key parameters in our model by virtue of their influence and importance are: $\gamma_{am}, \; \gamma_{at}, \; \mu_b$ and $\omega_{cm}$. Studying the effect of these parameters on the basic reproduction number $(\mathcal{R}_0)$ is crucial in order to understand the disease dynamics. In our context, this will help in either choosing a particular mode of treatment or adopting a synergistic therapy approach in the treatment of typhoid fever. Figure 12 presents the effect of each of these parameters on $\mathcal{R}_0$, when all other parameter values are as shown in Table 2. The parameters $\gamma_{am}$ and $\gamma_{at}$ measure the rate at which individuals in the acute stage of disease choose the modern medicine and traditional medicine, respectively. Meanwhile $\omega_{cm}$ measures the rate at which patients in the chronic stage of typhoid choose modern treatment, and $\mu_b$ is the decay rate of the bacteria in the environment. Figure 12a shows that $\gamma_{am} = 2.07$ (30% relative error) is the minimum value for $\mathcal{R}_0 < 1$. With this value, $\mathcal{R}_0 = 0.9967 < 1$, which means the disease will be under control. Figure 12b gives a minimum value of $\gamma_{at} = 2.49$ (36.38% relative error). The corresponding value of the basic reproduction number is $\mathcal{R}_0 = 0.9904 < 1$. Figure 12c focuses on the effect of concentration of bacterium Salmonella typhi in the environment. The value of $\mu_b$ should be $\mu_{b} = 0.44$ (45.5% relative error), and this gives $\mathcal{R}_0 = 0.9914 < 1$. Interestingly, unlike the other three, the output in Figure 12d, indicates that there is no value of $\omega_{cm}$ for which $\mathcal{R}_0 < 1$. This tells us that with the choice of modern medicine alone, the disease will prevail. It is therefore obvious that symptoms of the disease will not disappear for chronically sick patients receiving modern treatment only. There is a need to explore other means through which the disease can be controlled.

Contour plots make it easy for us to simultaneously compare how different conditions on the influential and important parameters affect $\mathcal{R}_0$. We use contour plots to display the change of $\mathcal{R}_0$ in the 2-dimensional space parameter $(\omega_{cm}, \; \gamma_{at})$, $(\omega_{cm}, \; \mu_b)$, $(\gamma_{am}$, and $\gamma_{at})$, respectively. These plots help to in determining joint parameters values that produce a value of $\mathcal{R}_0 < 1$. The contour plots are given in Figure 13. Figure 13a depicts the correlation between the rate at which patients in the chronic stage of typhoid choose modern treatment and the rate at which infected persons in the acute phase choose traditional medicine and their effect on the reproduction number. It is seen that $\mathcal{R}_0 < 1$, if $(\omega_{cm}\geq 1.0, \gamma_{at}\geq 1.1)$. This is an indication that neither modern medicine alone nor traditional medicine alone can eradicate the disease. Figure 13b shows the correlation between the rate at which patients in the chronic stage of typhoid choose modern treatment and the decay rate of the bacteria in the environment and their effect on $\mathcal{R}_0$. The output demonstrates that to maintain $\mathcal{R}_0 < 1$, it is important to maintain the shedding rate of Salmonella typhi bacteria within the interval $[0.17, 0.32]$ and ensure that the rate at which patients choose modern treatment should be above 1.3767. Similarly, Figure 13c reveals that, it is important that the rate at which infected individuals choose modern medicine is above 2.8371, while the rate at which infected persons opt for traditional medicine must be at least 2.1786, in order for $\mathcal{R}_0$ to be less than unity. In summary, the plots emphasize the necessity to associate traditional medicine with modern medicine in the treatment of typhoid fever.

6.   Conclusions

Typhoid fever has remained a major public health concern in Cameroon, especially in communities where there are inadequate preventive measures against the spread of the disease such as, good sanitation habits, eating safe food, and drinking potable water. Resistance of the disease to antibiotics is worrisome. Mathematical models that will improve our understanding suggest that new combination therapies and good practices that culminate in tackling and reducing the disease burden are important. In this paper, we have proposed a typhoid fever mathematical model that takes into account the evolution of the disease from the incubation phase to the invasion (acute) and chronic phases. The key novelty of the proposed model is the inclusion of modern and traditional treatment components. This is done to assess the impacts of both modes of treatment on the Salmonella typhi bacteria. Theoretical analysis of the model was performed. We have proved that the basic reproduction number $\mathcal{R}_0$ is a sharp threshold that ensures the global asymptotic stability of the DFE when its value is less than one and endemicity of the disease otherwise. Thus, for efficient control of the disease, one has to reduce $\mathcal{R}_0$ to below one.

We calibrated the model using weekly reported cumulative cases of typhoid fever in the Penka-Michel health district in Cameroon from 30 November 2019 to 1 January 2023. The model provides a good fit for data, and this means the predictions will be more reliable and accurate over time. The model was validated using weekly cumulative cases from 2 January 2023 to 5 November 2023. It is common knowledge that a well-calibrated model will have an $\mathcal{R}_0$ value that closely matches the observed transmission dynamics of a disease. We obtained $\mathcal{R}_0 = 1.2058$, which is an indication that the disease would persist if the necessary control measures are not put in place. In addition, the model predicts that by the end of 2026, more than 21,270 reported cases will be registered. Numerical results suggest that a reduction in the number of cases will be achieved if patients are advised to at least receive some treatment (modern or traditional) and eventually to switch from modern to traditional treatment if symptoms persist. Conversely, we have found that relapses due to antibiotics resistance will increase the level of the disease.

The results from the proposed model imply that integrating traditional medicine into conventional medicine will be effective in the treatment of typhoid fever and prevent antibiotics resistance. Traditional medicine is a useful adjuvant therapy for patients who choose modern medicine. Moreover, it is common practice in Cameroon that many typhoid patients, after receiving antibiotics, continue their treatment at home with traditional medicine. Hence, adopting a synergistic therapy approach in the treatment of typhoid fever will significantly mitigate typhoid disease cases.

In future research, the proposed model can be extended to include the optimal control of vaccination with relapse and reinfection.

Use of AI tools declaration

The authors declare they have not used artificial intelligence (AI) tools in the creation of this article.

Acknowledgments

The third author acknowledges Prof. Jean Lubuma for his support through the National Research Foundation (NRF) under the Competitive Programme for Rated Researchers (CPRR: Grant no. 138013) and the Science Faculty Start-up Funds for Research and the Postdoctoral Programme at the University of Witwatersrand, Johannesburg. The authors will like to thank the following two staff members (Mr. Gaston Keunang and Dr. Yvonne Wamba Mbaduet) of the Health District in Penka-Michel, Cameroon, for generously sharing with us the data used in this article. Finally, the authors thank the anonymous reviewers and the editor for their insightful remarks which helped to improve the quality of the manuscript.

Conflict of interest

The authors declare no conflict of interest.

Appendix

Appendix A.   Next–generation matrix coefficients

The coefficients of the next-generation matrix are given as follows:

Appendix B.   Proof of Theorem 3.5

Let us consider the Volterra candidate Lyapunov function,

where $a_i, \; i = 1, \cdots, 5$ are positive numbers to be determined. We use the following function $H(c) = 1 -c +\ln(c)$, which is negative when $c > 0$, and is equal to zero for $c = 1$. The time derivative of $Z$ along the trajectories of System (2.6) yields

The following relations can be shown from (2.6):

Let us denote $\lambda S = f(B, S)$, i.e., $\lambda^{\ast}S^{\ast} = f(B^{\ast}, S^{\ast})$. To prove the global asymptotic stability of the endemic equilibrium, we use the following inequality.

As

one has,

Substituting Eq (B.3) into Eq (B.2) leads to (B.14).

According to (B.6), we have

Now using Eqs (B.7) to (B.13), we have

Let $l(x, x^{\ast}) : = \dfrac{x}{x^{\ast}}-\ln\dfrac{x}{x^{\ast}}\geq 1 > 0$ for $x > 0,$ and $x^{\ast} > 0$. Therefore,

Using the fact that $-\left\lbrace a_4\gamma_{at}I_a^{\ast}+a_4\psi_{mt}M_m^{\ast} \right\rbrace l(M_t, M_t^{\ast}) = -k_5M_t^{\ast}l(M_t, M_t^{\ast}) < 0$, we get

We choose $a_i$ such that the expressions in the square brackets vanish. That is, $a_i$ are solutions of the system

Then by (B.18) and (B.19), we have:

Equation (B.22) leads to

Consequently,

Let us consider (B.20) and (B.21). Here $a_2$ and $a_5$ verify the system

Using Eq (B.25), we have

For $a_4 = 1,$ we get:

These values of $a_i, \; i = 1, \cdots, 5$ imply that $\dfrac{dZ}{dt}\leq 0.$

Furthermore, the equality $\dfrac{dZ}{dt} = 0$ holds only for

Thus, $\{\mathcal{\varepsilon}^{\ast}\}$ is the largest positive invariant set which is contained in the set

Hence, it follows from LaSalle's invariance principle ^[33] that any solution of Eq (2.6) with an initial condition in $\varOmega$ converges to the endemic equilibrium point $\varepsilon^{\ast},$ as $t \longrightarrow \infty.$ Therefore, the positive equilibrium $\varepsilon^{\ast}$ is globally asymptotically stable if $\mathcal{R}_0 > 1.$