# Mathematical Modeling and Analysis of Transmission Dynamics and Control of Schistosomiasis

Authors: Ebrima Kanyi; Ayodeji Sunday Afolabi; Nelson Owuor Onyango
Journal: Journal of Applied Mathematics
 (2021)

---

## Abstract

This paper presents a mathematical model that describes the transmission dynamics of schistosomiasis for humans, snails, and the free living miracidia and cercariae. The model incorporates the treated compartment and a preventive factor due to water sanitation and hygiene (WASH) for the human subpopulation. A qualitative analysis was performed to examine the invariant regions, positivity of solutions, and disease equilibrium points together with their stabilities. The basic reproduction number,R0, is computed and used as a threshold value to determine the existence and stability of the equilibrium points. It is established that, under a specific condition, the disease-free equilibrium exists and there is a unique endemic equilibrium when R0>1. It is shown that the disease-free equilibrium point is both locally and globally asymptotically stable provided R0<1, and the unique endemic equilibrium point is locally asymptotically stable whenever R0>1 using the concept of the Center Manifold Theory. A numerical simulation carried out showed that at R0=1, the model exhibits a forward bifurcation which, thus, validates the analytic results. Numerical analyses of the control strategies were performed and discussed. Further, a sensitivity analysis of R0 was carried out to determine the contribution of the main parameters towards the die out of the disease. Finally, the effects that these parameters have on the infected humans were numerically examined, and the results indicated that combined application of treatment and WASH will be effective in eradicating schistosomiasis.


---

## Body

## 1. Introduction

Schistosomiasis is a disease that is common in the tropics, and it is caused by a group of parasitic worms known as schistosomes. These collections of parasitic worms are, in some instances, called blood flukes. Several species of schistosomes cause illnesses in human beings. The three major ones areSchistosoma haematobium, Schistosoma mansoni, and Schistosoma japonicum. These species of schistosomes can be found in parts of Africa, Asia, and South America, [1]. Schistosoma mansoni and Schistosoma japonicum mainly cause diseases in the liver and bowels but Schistosoma haematobium mostly affects the urinary and genital areas [2].The life cycles of schistosomes are very complex and involve fresh water and (mostly) two hosts: a definite human host (also cattle in the case ofSchistosoma japonicum) and an intermediate snail host. Schistosomes reproduce both sexually and asexually [2]. The sexual reproduction takes place in the definite hosts, and the asexual reproduction takes place in an intermediate snail host. The life cycle begins in the human portal system. Schistosome eggs are released to freshwater by infected humans through their feces or urine. These eggs hatch into larvae once the eggs enter freshwater. These larvae are known as miracidia. The miracidia, within a day, find and infect the intermediate host, the snails. After a period of development in the snail, the miracidia further mature and transform into cercariae through asexual reproduction and are released into the water in the presence of light. The cercariae (with a life expectancy of up to 2 days) will further swim and penetrate into the skin of a definite host (humans and also cattle in the case of Schistosoma Japonicum). The cercariae find their way to the liver through the lungs and develop into matured adult worms. Thereafter, they move into the blood vessels and the organs around it. These worms can live in these organs for many years. The adult schistosomes can lay up to 300 eggs within a definite host daily, and these eggs can penetrate into the bladder or intestines. Some of these eggs get released during urination or defecation and may find their way back into fresh waters, and their life cycles will begin afresh. However, other eggs will be trapped in some of these organs, and their activities cause damage to the organs. It takes between 4 and 6 weeks for an exposed person to become infectious [2].The disease is contracted by humans through contact with infested freshwater which is inhabited by infected snails; there is no direct human to human transmission. People with schistosomiasis can suffer from rash or itchy skin, cough, fever, blood in urine or stool, muscle ache, stunted growth, etc. Later in their lives, they can experience malfunctioning of some major organs like the liver, intestines, lungs, and bladder due to damage caused by the eggs trapped in their body. These eggs can be found in the spinal cord and even the brain in some rare cases. It may result in paralysis, inflammation of the spinal cord, or seizures. Children who are reinfected with this disease may experience malnutrition, anemia, and even difficulty of learning [1].There is still no commercially available vaccine for schistosomiasis. However, praziquantel is the WHO-approved medication for the treatment of the disease. It works by killing the adult worms, hence preventing the production of new eggs in the human host, but does not prevent reinfection as treated individuals can easily get reinfected [3]. Improved sanitation and snail control are also mechanisms used to put the disease under control. Snail control involves the use of molluscicides or biological control methods like use of competitor snails (see [4, 5]); however, their effectiveness varies and the disease continues to spread.Mathematical models can be used as guiding tools to effectively study the spread of diseases, and they have over the years helped both policymakers and public health officials in their decision-making processes with respect to key intervention programs. Its role in the study and understanding of the transmission dynamics as well as the effectiveness of the various control strategies of many infectious diseases has been considered to be very effective. Several mathematical models were designed over the years to study the transmission dynamics and control of the spread of schistosomiasis [4–11].A review of a number of existing mathematical models of schistosomiasis, their merits, and demerits was carried out by [12]. [13] suggested ways by which controls could be incorporated into mathematical models in order to eradicate or reduce the prevalence of schistosomiasis in an endemic area. [14] presented some of the measures adopted for the control schistosomiasis in the People’s Republic of China. [15] integrated a mathematical model for schistosomiasis to obtain a set of nonlinear integral equations. The Contraction Mapping Principle was then used to establish the existence of a unique solution of the schistosomiasis model.The transmission dynamics of schistosomiasis among humans, cattle, and snails was presented in a mathematical model by [16]. Sensitivity analysis of the model parameters on the basic reproduction number showed that cattle-snail transmission of S. japonicum played an important role in the endemicity of the disease in Hubei Province. [17] performed stability analysis on a schistosomiasis model, and the result indicated that the prevalence of schistosomiasis was affected by the latent period of infection and that effective public health education campaign on drug treatment of the disease could help to curtail it.Schistosomiasis is declared by the WHO as one of the most deadly of the neglected tropical diseases. Despite numerous intervention programs to curb the spread of the disease, it still remains endemic in many tropical and subtropical countries particularly in the poverty-stricken African countries [1]. Mathematical models have been used to model the control and spread of many diseases [18, 19]. Hence, there is the need to develop a mathematical model that will help to explain the spread and control of schistosomiasis. To the very best of our knowledge, the proposed model is the first to include the treated human class in a schistosomiasis model and has incorporated the exposed period for both human and snail subpopulations alongside a parameter for water, sanitation, and hygiene (WASH). In order to capture the complex nature of the disease dynamics and proffer effective strategies for its control, compartments that have to do with the incubation periods for the hosts that are infected but not releasing pathogen into the environment alongside the treated humans should be incorporated. Consequently, a deterministic model that divides the populations into susceptible humans and snails, infected humans and snails, exposed humans and snails, and treated humans, with the dynamics of free living miracidia and cercariae together with their interactions is proposed.


## 2. Model Formulation

Here, a deterministic model that describes the transmission dynamics of schistosomiasis is proposed. We considered the fact that once a pathogen is introduced within a population, it divides the population into fractions called compartments or components. The human subpopulation is divided into four compartments: susceptible humansSh, exposed humans Eh, infected humans Ih, and treated humans Th; the snail subpopulation into three compartments: susceptible snails Ss, exposed snails Es, and infected snails Is; and compartments for the free living miracidia Nm and free living cercariae Nc based on the life cycle of the parasite.The susceptible humans comprise those individuals capable of contracting the disease but have yet to do so. The exposed humans include those individuals who contract cercariae, undergoing and incubation period and have yet to begin producing (miracidia) eggs. The infected humans are those individuals that are symptomatic and have (miracidia) eggs appearing in their feces. The treated humans include individuals who have received treatment. Praziquantel works only by killing the adult worms; thus, only the infected humans are moved to the treated class, and the treated humans are assumed not to be infectious, i.e., they do not produce eggs for miracidia. The susceptible snails are those snails that have the possibility of getting exposed to miracidia but have yet to do so. The exposed snails are those snails that get infected by miracidia but have yet to shed cercariae. The infected snails are those snails that contract miracidia, survived the exposed period, and are shedding cercariae. Moreover, a parameter for prevention through water sanitation and hygiene (WASH) has been incorporated. Table1) gives a description of various parameters used in the model, and Figure 1) is the schematic diagram for the proposed model.Table 1
A table showing the description of the model parameters and parameter values.
ParameterDescriptionValueSourceαhRecruitment rate of humans254d−1EstimateαsRecruitment rate of snails3000d−1[27]βhRate of transmission of humans from susceptible to exposed0.09753Lcer−1d−1[28]βsRate of transmission of snails from susceptible to exposed0.615Lmir−1d−1[7, 29]ξhPreventive factor due to WASH0-1EstimateμhNatural death rate of humans0.00004379d−1EstimateμsNatural death rate of snails0.000569d−1[7]μmNatural death rate of miracidia0.9d−1[7]μcNatural death rate of cercariae0.004d−1[7]δhDeath rate of humans due to infection0.000274-0.000913d−1EstimateδsDeath rate of snails due to infection0.0004012d−1[7]σhRate of transmission of humans from exposure to infection0.0238-0.0286d−1EstimateσsRate of transmission of snails from exposure to infection0.0286-0.0357d−1EstimateγhTransmission rate of humans from infection to treatment0.03d−1[29]λ1Rate individuals produce miracidia6.96mirhost−1d−1[9]λ2Rate snails produce cercariae2.6cerhost−1d−1[7]ρTreatment efficacy0.8Assumedm0Miracidia saturation constant1×108[7]c0Cercariae saturation constant9×107[7]εLimitation of the growth velocity0.2[7]Figure 1
A compartmental chart for the mathematical model.The model is derived based on these assumptions:(1)
There is no mother-to-child transmission of the disease among humans, and treated humans have no immunity(2)
The rates of depletion of the parasites by humans or snails have no effect on the population dynamics of the parasites(3)
Contact with free living cercariae is the only mean through which susceptible humans become exposed to the disease(4)
There are no introduction of infectious (humans or snails) into the studied population(5)
Treated humans do not produce eggs for the production of miracidia(6)
Infected snails do not get treated and, hence, do not recoverwhereωh=1−ξhβh/c0+εNc and ωs=βs/m0+εNm.The model is represented by the following system of ODEs:(1)dShdt=αh+ρTh−1−ξhβhNcShc0+εNc−μhSh,dEhdt=1−ξhβhNcShc0+εNc−σhEh−μhEh,dIhdt=σhEh−γhIh−δhIh−μhIh,dThdt=γhIh−ρTh−1−ρδhTh−μhTh,dNmdt=λ1Ih−μmNm,dSsdt=αs−βsNmSsm0+εNm−μsSs,dEsdt=βsNmSsm0+εNm−σsEs−δsEs−μsEs,dIsdt=σsEs−δsIs−μsIs,dNcdt=λ2Is−μcNc.With initial conditionsSh0≥0,Eh0≥0,Ih0≥0,Th0≥,Nm0≥0,Ss0≥0,Es0≥0,Is0≥0,andNc0≥0.


## 3. Model Analysis

### 3.1. Invariant Region and Positivity of Solutions

Since a population model always requires positive solutions, it is relevant that all solutions be nonnegative. Consider the model as described by the model system equation (1). The total human subpopulation at any time, t, is given by(2)Nht=Sht+Eht+Iht+Tht,dNhdt=αh−μhNh−δhIh+1−ρTh.In the absence of mortality due to infection, we have(3)dNhdt≤αh−μhNh.Thus,Nh≤αh/μh.Similarly, the total snail subpopulation is given by(4)Nst=Sst+Est+Ist,dNsdt=αs−μsNs−δsEs+Is.In the absence of mortality due to infection, we get(5)dNsdt≤αs−μsNs.So,Ns≤αs/μs.From equation (18) of the state system equation (1), we have(6)dNmdt=λ1Ih−μmNm,Nm≤αhλ1μhμm.From equation (25) of the system equation (1), we have(7)dNcdt=λ2Is−μcNc,Nc≤αsλ2μsμc.Thus,Nh, Ns, Nm, and Nc are bounded, and the region(8)Ω=Sh,Eh,Ih,Th,Nm,Ss,Es,Is,Nc∈ℝ+9:Nh≤αhμh,Nm≤αhλ1μhμm,Ns≤αsμs,Nc≤αsλ2μsμc,is positively invariant, and accordingly, the system equation (1) is epidemiologically meaningful and analytically well posed in the region Ω.Theorem 1.
Suppose that(9)Sh0,Eh0,Ih0,Th0,Nm0,Ss0,Es0,Is0,Nc0∈ℝ+9∪0,is the initial condition for the model system equation (1), then the set of solutions(10)Sht,Eht,Iht,Tht,Nmt,Sst,Est,Ist,Nct,is nonnegative ∀t≥0.Proof.
Given that the set of initial conditions is nonnegative, consider the model system equation (1). Then,(11)dShtdt=αh+ρTh−ωhNcSh−μhSh,dShtdt≥−ωhNc+μhSh.
Separating the variables and integrating, we have(12)Sht≥Sh0e∫−ωhNcdt−μht≥0,∀t≥0.
Similarly,(13)Eht≥Eh0e−σh+μht≥0,∀t≥0,Iht≥Ih0e−γh+δh+μht≥0,∀t≥0,Tht≥Th0e−ρ+1−ρδ+μht≥0,∀t≥0,Nmt≥Nm0e−μmt≥0,∀t≥0,Sst≥Ss0e∫−ωsNmdt−μst≥0,∀t≥0,Est≥Es0e−σs+δs+μst≥0,∀t≥0,Ist≥Is0e−δs+μst≥0,∀t≥0,Nct≥Nc0e−μct≥0,∀t≥0.


### 3.2. Equilibrium Points

At the equilibrium point, all states are considered to have a constant solution. So, it is obtained by settingdSh∗/dt=dEh∗/dt=dIh∗/dt=dTh∗/dt=dNm∗/dt=dSs∗/dt=dIs∗/dt=dNc∗/dt=0. It is given by(14)E=Sh∗,Eh∗,Ih∗,Th∗,Nm∗,Ss∗,Es∗,Is∗,Nc∗.Due to the complexity of the model system 1, all other state variables at the steady state are expressed in terms of the steady state ofIh which is denoted as Ih∗. Accordingly,(15)Sh∗=αhρ+1−ρδh+μh+γhρIh∗c0μcδs+μsσs+δs+μsβs+εμsλ1Ih∗+m0μsμm+αsβsσsλ1λ2εIh∗ρ+1−ρδh+μh1−ξhβh+εμhαsβsσsλ1λ2Ih∗+c0μhμcδs+μsσs+δs+μsβs+εμsλ1Ih∗+m0μsμm,Eh∗=γh+δh+μhIh∗σh,Th∗=γhIh∗ρ+1−ρδh+μh,Nm∗=λ1Ih∗μm,Ss∗=αsm0μm+ελ1Ih∗βs+εμsλ1Ih∗+m0μsμm,Es∗=αsβsλ1Ih∗σs+δs+μsβs+εμsλ1Ih∗+m0μsμm,Is∗=αsβsσsλ1Ih∗δs+μsσs+δs+μsβs+εμsλ1Ih∗+m0μsμm,Nc∗=αsβsσsλ1λ2Ih∗μcδs+μsσs+δs+μsβs+εμsλ1Ih∗+m0μsμm,and a combination of the above equation (15) and the second equation of the model system 1, the following relation was obtained for Ih∗:(16)Ih∗Ih∗−AR02−1B+C=0, where(17)A=m0c0μhμsμmμcσh+δhγh+δh+μhδs+μsσs+δs+μsρ+1−ρδh+μh,B=1−ξhβhαsβsσsλ1λ2sigmahγh1−ρδh+μh+ρ+1−rhoδh+μhσhδs+μs+μhγh+δh+μh,C=λ1μhσh+μhγh+δh+μhρ+1−ρδh+μhαsβsσhσsλ2ε+c0μcδs+μsbetas+εμsσs+δs+μs.


#### 3.2.1. Disease-Free Equilibrium Point

At the disease-free equilibrium, there are no infections or recovery, and thus, no miracidia or cercariae is produced. Accordingly, at this point,Ih∗ in the system equation (16) must be zero. Hence, the point is obtained as(18)E0=Sh0=αhμh,Eh0=0,Ih0=0,Th0=0,Nm0=0,Ss0=αsμs,Es0=0,Ih0=0,Nc0=0,where Ih∗=0.


#### 3.2.2. Endemic Equilibrium Point

At the endemic equilibrium point, all disease states in equation (14) are considered to be positive, and consequently, Ih∗ must be greater than zero for all the other states to be positive. Therefore, the unique positive endemic equilibrium point exists only when Ih∗ in the system equation (16) is positive, which is also in this case the same as R0>1. And thus, it is the point(19)E∗=Sh∗∗,Eh∗∗,Ih∗∗,Th∗∗,Nm∗∗,Ss∗∗,Es∗∗,Is∗∗,Nc∗∗,where(20)Sh∗∗=Sh∗,Eh∗∗=Eh∗,Ih∗∗=AR02−1B+C,Th∗∗=Th∗,Nm∗∗=Nm∗,Ss∗∗=Ss∗,Es∗∗=Es∗,Is∗∗=Is∗,Nc∗∗=Nc∗,whenever R0>0.


### 3.3. Stability Analysis

#### 3.3.1. Basic Reproduction Number,R0

The basic reproduction number is defined as the average number of secondary infections caused by the emergence of an infectious individual into a completely susceptible population [20]. It is denoted by R0. We considered the next generation matrix approach according to Chavez et al. [21] (for details, see Appendix A). The next generation matrix is denoted by G and is obtained as(21)G=0000λ2σs1−ξhβhαhc0μhμcδs+μsσs+δs+μsλ21−ξhβhαhc0μhμcδs+μs1−ξhβhαhc0μhμc000000000000000000000λ1σhβsαsm0μsμmσh+μhγh+δh+μhλ1βsαsm0μsμmγh+δh+μh0βsαsm0μsμm00000000000000000.Accordingly, the basic reproduction number,R0, which is the dominant eigenvalue of G, is obtained as(22)R0=αh1−ξhβhσhαsβsσsλ1λ2c0m0μhμsμmμcσh+μhγh+δh+μhδs+μsσs+δs+μs.For simplicity, sinceR0 denotes the possible number of uninfected a particular infected individual (in this case either human or snail) can infect, we can rewrite R0 as(23)R0=R0HS·R0SH, where R0HS=αh1−ξhβhσhλ1/m0μhμmσh+μhγh+δh+μh represents the number of snails an infected human is capable of infecting through the production of miracidia and R0SH=αsβsσsλ2/c0μsμcδs+μsσs+δs+μs denotes the number of humans an infected snail is capable of infecting through the production of cercariae.


#### 3.3.2. Stability of Disease-Free Equilibrium Point

Theorem 2.
The disease-free equilibrium point,E0, is locally asymptotically stable provided R0<1 and unstable otherwise.Proof.
To verify the validity of the above proposition, see Theorem 2 in Van den Driessche and Watmough [22].
The method according to Chavez et al. [21] is used to determine the global stability of the disease-free equilibrium point within Ω∈ℝ+9. Hence, we denote the system equation (1) by(24)dXdt=FX,Y,dYdt=GX,Y,GX,0=0,where X=Sh,Ss comprises of the uninfected subpopulations and Y=Eh,Ih,Th,Nm,Es,Is,Nc comprises of the infected and infectious subpopulations. E0 is guaranteed to be globally asymptotically stable (GAS) if R0<1 and the following two conditions (C1) and (C2) hold:
(C1): FordX/dt=FX,0,X∗=E0 is GAS
(C2):GX,Y=AY−G∗X,Y,G∗X,Y≥0forX,Y∈Ω
whereA=DYGX∗,0 and Ω∈ℝ+9 is the invariant region of the model.Theorem 3.
The disease-free equilibrium point,E0=αh/μh,0,0,0,0,αs/μs,0,0,0∈Ω, is GAS for the system equation (1) provided that R0<1 and conditions (C1) and (C2) both hold.Proof.
From the system of equation (1), we have(25)FX,0=αh−μhShαs−μsSs.
First, suppose that the time,t⟶∞, we need to show that X⟶E0.
But we know that forFX,0,(26)dShdt=αh−μhSh,dShdt+μhSh=αh,Sh=αhμh.
Also,(27)dSsdt=αs−μsSs,dSsdt+μhSs=αs,Ss=αsμs.
Thus, all points with respect to this condition converge at(28)E0=αhμh,0,0,0,0,αsμs,0,0,0∈Ω.
Now, from the system equation (1), we know that(29)dYdt=GX,Y,where(30)GX,Y=ωhShNc−σh+μhEhσhEh−γh+δh+μhIhγhIh−ρ+1−ρδh+μhThλ1Ih−μmNmωsSsNm−σs+δs+μsEsσsEs−δs+μsIsλ2Is−μcNc,and GX,Y=AY−G∗X,Y, where(31)A=−σh+μh00000αh1−ξhβhc0μhσh−γh+δh+μh000000γh−ρ+1−ρδh+μh00000λ10−μm000000αsβsm0μs−σs+δs+μs000000σs−δs+μs000000λ2−μc,G∗X,Y=αhc0μh−Shc0+εNc1−ξhβhNc000αsm0μs−Ssm0+εNmβsNm00.
Clearly,G∗X,Y≥0,∀X,Y∈Ω, since(32)Nh≤αhμh,Ns≤αsμs⇒αhc0μh−Shc0+εNc≥0,αsm0μs−Ssm0+εNm≥0.
Moreover,A is an M-matrix, since the nondiagonal entries are nonnegative [21]. Therefore, the DFE, denoted E0, is globally asymptotically stable provided that R0<1.


#### 3.3.3. Stability of the Endemic Equilibrium Point

We consider the concept ofCenter Manifold Theory as described in Theorem 4.1 by Chavez and Song [23]. The theorem helps to determine the course of the bifurcation at a critical value of the parameter [24]. The theorem is reproduced in Appendix B.Theorem 4.
The endemic equilibrium point, based on the Theorem 4.1 by [23] is locally asymptotically stable for R0>1 (but near 1).Proof.
Now, consider the model system equation (1). Using Theorem 4.1 by Chavez and Song [23], let Sh=x1, Eh=x2, Ih=x3, Th=x4, Nm=x5, Ss=x6, Es=x7, Is=x8, and Nc=x9. Further, denote X=x1,x2,x3,x4,x5,x6,x7,x8,x9 and F=f1,f2,f3,f4,f5,f6,f7,f8,f9. Then, the model system equation (1) can be expressed in the form dX/dt=F, where(33)dx1dt=αh+ρx4−ωhx9x1−μhx1=f1,dx2dt=ωhx9x1−σhx2−μhx2=f2,dx3dt=σhx2−γhx3−δhx3−μhx3=f3,dx4dt=γhx3−ρx4−1−ρδh−μhx4=f4,dx5dt=λ1x3−μmx5=f5,dx6dt=αs−ωsx5x6−μsx6=f6,dx7dt=ωsx5x6−σsx7−δsx7−μsx7=f7,dx8dt=σsx7−δsx8−μsx8=f8,dx9dt=λ2x8−μcx9=f9.
The role of the parameterϕ in Theorem 4.1 by Chavez and Song [23] is played by βh with the critical value obtained from R0=1 [7, 24]. Hence, βh=βh∗, regarded as the bifurcation point, can be obtained by evaluating it at R0=1. So(34)βh∗=c0m0μhμsμmμcσh+μhγh+δh+μhδs+μsσs+δs+μsαh1−ξhαsβsσhσsλ1λ2.
Let the Jacobian of the transformed system equation (33) at E0 when βh=βh∗ is denoted by A. Then, it can be shown that 0 is a simple eigenvalue of the matrix A, where(35)A=−μh00ρ0000−αh1−ξhβh∗c0μh0−σh+μh000000αh1−ξhβh∗c0μh0σh−γh+δh+μh00000000γh−ρ+1−ρδh+μh0000000λ10−μm00000000−αsβsm0μs−μs0000000αsβsm0μs0−σsδs+μs00000000σs−δs+μs00000000λ2−μc,and if the vector w=w1,w2,w3,w4,w5,w6,w7,w8,w9T is a right eigenvector corresponding to the eigenvalue, 0 of A, then the components of w can be obtained as(36)w1=−αh1−ξhβh∗σhγh1−ρδh+μh+ρ+1−ρδh+μhσhδh+μh+μhγh+δh+μhc0μh2σh+μhγh+δh+muhρ+1−ρδh+μhw9,w2=αh1−ξhβh∗c0μhσh+μhw9,w3=αh1−ξhβh∗σhc0μhσh+μhγh+δh+μhw9,w4=αh1−ξhβh∗σhγhc0μhσh+μhρ+1−ρδh+μhγh+δh+μhw9,w5=αh1−ξhβh∗σhλ1c0μhμmσh+μhγh+δh+μhw9,w6=−αh1−ξhβh∗αsβsσhλ1m0c0μhμs2μmσh+μhγh+δh+μhw9,w7=αh1−ξhβh∗αsβsσhλ1m0c0μhμsμmσh+μhγh+δh+μhσs+δs+μsw9,w8=αh1−ξhβh∗αsβsσhσsλ1m0c0μhμsμmσh+μhγh+δh+μhdeltas+μsσs+δs+μsw9,w9=w9>0,
whereβh∗ is defined in equation (34).
Similarly, letv=v1,v2,v3,v4,v5,v6,v7,v8,v9 be a left eigenvector corresponding to the zero eigenvalue of A. Then, the components of v are(37)v1=v4=v6=0,v2=c0μhμcαh1−ξhβh∗v9,v3=c0μhμcσh+μhαh1−ξhβh∗v9,v5=c0μhμcσh+μhγh+δh+μhαh1−ξhβh∗λ1v9,v7=m0c0μhμsμmμcσh+μhγh+δh+μhαh1−ξhβh∗αsβsλ1v9,v8=m0c0μhμsμmμcσh+μhγh+δh+μhσs+δs+μsαh1−ξhβh∗αsβsσsλ1v9,v9=v9>0.
Now, to get the values ofa and b as in Theorem 4.1 by Chavez and Song [23], we have to get the second partial derivatives in the system equation (33). Clearly, we do not need the partial derivatives of f1,f4, and f6, since v1=v4=v6=0. Thus, the second (nonzero) partial derivatives from the linearized system equation (33) for a are(38)∂2f2∂x1∂x90,0=1−ξhβh∗c0,∂2f7∂x5∂x60,0=βsm0,and for b is(39)∂2f2∂x9∂βh∗0,0=1−ξhx10c0=αh1−ξhc0μh.
Therefore,(40)a=∑k,i,j=1nvkwiwj∂2fk∂xi∂xj0,0,a=v2w1w9∂2f2∂x1∂x9+v7w5w6∂2f7∂x5∂x6,a=−1−ξhβh∗μcσhγh1−ρδh+μh+ρ+1−ρδh+μhσhδh+μhγh+δh+μhc0μhσh+μhγh+δh+μhρ+1−ρδh+μhv9w92−αh1−ξhβh∗σh2λ1μcc0μhμmσh+μhγh+δh+μhρ+1−ρδh+μh2v9w92.
Soa<0.(41)b=∑k,i=1nvkwi∂2fk∂xi∂ϕ0,0,b=v2w9∂2f2∂x9∂βh∗0,0,b=μcβh∗v9w9,
Sob=μc/βh∗v9w9.
Clearly,a<0 and b>0, and with the help of condition (iv) in Theorem 4.1 by Chavez and Song [23], it is confirmed that the system follows a forward bifurcation at R0=1(see Figure 2), and hence, the endemic equilibrium is locally asymptotically stable for R0>1, but sufficiently closed to 1.Figure 2
A forward bifurcation showing an exchange of stability of the equilibrium points atR0=1.Generally, also based on the principles of the next generation matrix according to Van den Driessche and Watmough [22], the disease-free equilibrium point is stable when R0<1 and unstable when R0>1. To validate these analytic results, a numerical simulation was carried out. Figure 2 shows a graph of R0 against the steady state of Ih, denoted by Ih∗. The figure shows that at R0=1, the model system equation (1) exhibits a forward bifurcation. Forward bifurcations have been observed in many other epidemic models [7, 23, 25]. It shows an exchange of stability between the disease-free equilibrium point and the disease endemic equilibrium point. To this effect, it can be concluded that the model system equation (1) has a disease-free equilibrium point that is globally asymptotically stable when R0<1, and when R0>1, it posses an endemic equilibrium point that is locally asymptotically stable. This is actually a validation of the results from the analytic processes which also confirms that the disease-free equilibrium point is both locally and globally asymptotically stable when R0<1, and when R0>1, a locally asymptotically stable unique endemic equilibrium point exists.Figure 3
Evolution of each subpopulation with time.
(a)
Evolution of the human subpopulation with time(b)
Evolution of the snail subpopulation with time(c)
Evolution of the free living miracidia with time(d)
Evolution of the free living cercariae with timeFigure 4
Effects of WASH on the dynamics of humans and the pathogens.
(a)
Effects of varying WASH on susceptible humans(b)
Effects of varying WASH on the infected humans(c)
Effects of varying WASH on miracidia(d)
Effects of varying WASH on cercariaeFigure 5
Effects of varying treatments on infected humans, miracidia, and cercariae.
(a)
Effects of varying treatments on the infected humans(b)
Effects of varying treatments on miracidia(c)
Effects of treatments on cercariaeFigure 6
Effects of varying sensitive parameters on infected humans.
(a)
Effects of varyingαh on infected humans(b)
Effects of varyingαs on infected humans(c)
Effects of varyingβh on infected humans(d)
Effects of varyingβs on infected humans(e)
Effects of varyingσh on infected humans(f)
Effects of varyingλ1 on infected humans(g)
Effects of varyingλ2 on infected humans


## 3.1. Invariant Region and Positivity of Solutions

Since a population model always requires positive solutions, it is relevant that all solutions be nonnegative. Consider the model as described by the model system equation (1). The total human subpopulation at any time, t, is given by(2)Nht=Sht+Eht+Iht+Tht,dNhdt=αh−μhNh−δhIh+1−ρTh.In the absence of mortality due to infection, we have(3)dNhdt≤αh−μhNh.Thus,Nh≤αh/μh.Similarly, the total snail subpopulation is given by(4)Nst=Sst+Est+Ist,dNsdt=αs−μsNs−δsEs+Is.In the absence of mortality due to infection, we get(5)dNsdt≤αs−μsNs.So,Ns≤αs/μs.From equation (18) of the state system equation (1), we have(6)dNmdt=λ1Ih−μmNm,Nm≤αhλ1μhμm.From equation (25) of the system equation (1), we have(7)dNcdt=λ2Is−μcNc,Nc≤αsλ2μsμc.Thus,Nh, Ns, Nm, and Nc are bounded, and the region(8)Ω=Sh,Eh,Ih,Th,Nm,Ss,Es,Is,Nc∈ℝ+9:Nh≤αhμh,Nm≤αhλ1μhμm,Ns≤αsμs,Nc≤αsλ2μsμc,is positively invariant, and accordingly, the system equation (1) is epidemiologically meaningful and analytically well posed in the region Ω.Theorem 1.
Suppose that(9)Sh0,Eh0,Ih0,Th0,Nm0,Ss0,Es0,Is0,Nc0∈ℝ+9∪0,is the initial condition for the model system equation (1), then the set of solutions(10)Sht,Eht,Iht,Tht,Nmt,Sst,Est,Ist,Nct,is nonnegative ∀t≥0.Proof.
Given that the set of initial conditions is nonnegative, consider the model system equation (1). Then,(11)dShtdt=αh+ρTh−ωhNcSh−μhSh,dShtdt≥−ωhNc+μhSh.
Separating the variables and integrating, we have(12)Sht≥Sh0e∫−ωhNcdt−μht≥0,∀t≥0.
Similarly,(13)Eht≥Eh0e−σh+μht≥0,∀t≥0,Iht≥Ih0e−γh+δh+μht≥0,∀t≥0,Tht≥Th0e−ρ+1−ρδ+μht≥0,∀t≥0,Nmt≥Nm0e−μmt≥0,∀t≥0,Sst≥Ss0e∫−ωsNmdt−μst≥0,∀t≥0,Est≥Es0e−σs+δs+μst≥0,∀t≥0,Ist≥Is0e−δs+μst≥0,∀t≥0,Nct≥Nc0e−μct≥0,∀t≥0.


## 3.2. Equilibrium Points

At the equilibrium point, all states are considered to have a constant solution. So, it is obtained by settingdSh∗/dt=dEh∗/dt=dIh∗/dt=dTh∗/dt=dNm∗/dt=dSs∗/dt=dIs∗/dt=dNc∗/dt=0. It is given by(14)E=Sh∗,Eh∗,Ih∗,Th∗,Nm∗,Ss∗,Es∗,Is∗,Nc∗.Due to the complexity of the model system 1, all other state variables at the steady state are expressed in terms of the steady state ofIh which is denoted as Ih∗. Accordingly,(15)Sh∗=αhρ+1−ρδh+μh+γhρIh∗c0μcδs+μsσs+δs+μsβs+εμsλ1Ih∗+m0μsμm+αsβsσsλ1λ2εIh∗ρ+1−ρδh+μh1−ξhβh+εμhαsβsσsλ1λ2Ih∗+c0μhμcδs+μsσs+δs+μsβs+εμsλ1Ih∗+m0μsμm,Eh∗=γh+δh+μhIh∗σh,Th∗=γhIh∗ρ+1−ρδh+μh,Nm∗=λ1Ih∗μm,Ss∗=αsm0μm+ελ1Ih∗βs+εμsλ1Ih∗+m0μsμm,Es∗=αsβsλ1Ih∗σs+δs+μsβs+εμsλ1Ih∗+m0μsμm,Is∗=αsβsσsλ1Ih∗δs+μsσs+δs+μsβs+εμsλ1Ih∗+m0μsμm,Nc∗=αsβsσsλ1λ2Ih∗μcδs+μsσs+δs+μsβs+εμsλ1Ih∗+m0μsμm,and a combination of the above equation (15) and the second equation of the model system 1, the following relation was obtained for Ih∗:(16)Ih∗Ih∗−AR02−1B+C=0, where(17)A=m0c0μhμsμmμcσh+δhγh+δh+μhδs+μsσs+δs+μsρ+1−ρδh+μh,B=1−ξhβhαsβsσsλ1λ2sigmahγh1−ρδh+μh+ρ+1−rhoδh+μhσhδs+μs+μhγh+δh+μh,C=λ1μhσh+μhγh+δh+μhρ+1−ρδh+μhαsβsσhσsλ2ε+c0μcδs+μsbetas+εμsσs+δs+μs.


### 3.2.1. Disease-Free Equilibrium Point

At the disease-free equilibrium, there are no infections or recovery, and thus, no miracidia or cercariae is produced. Accordingly, at this point,Ih∗ in the system equation (16) must be zero. Hence, the point is obtained as(18)E0=Sh0=αhμh,Eh0=0,Ih0=0,Th0=0,Nm0=0,Ss0=αsμs,Es0=0,Ih0=0,Nc0=0,where Ih∗=0.


### 3.2.2. Endemic Equilibrium Point

At the endemic equilibrium point, all disease states in equation (14) are considered to be positive, and consequently, Ih∗ must be greater than zero for all the other states to be positive. Therefore, the unique positive endemic equilibrium point exists only when Ih∗ in the system equation (16) is positive, which is also in this case the same as R0>1. And thus, it is the point(19)E∗=Sh∗∗,Eh∗∗,Ih∗∗,Th∗∗,Nm∗∗,Ss∗∗,Es∗∗,Is∗∗,Nc∗∗,where(20)Sh∗∗=Sh∗,Eh∗∗=Eh∗,Ih∗∗=AR02−1B+C,Th∗∗=Th∗,Nm∗∗=Nm∗,Ss∗∗=Ss∗,Es∗∗=Es∗,Is∗∗=Is∗,Nc∗∗=Nc∗,whenever R0>0.


## 3.2.1. Disease-Free Equilibrium Point

At the disease-free equilibrium, there are no infections or recovery, and thus, no miracidia or cercariae is produced. Accordingly, at this point,Ih∗ in the system equation (16) must be zero. Hence, the point is obtained as(18)E0=Sh0=αhμh,Eh0=0,Ih0=0,Th0=0,Nm0=0,Ss0=αsμs,Es0=0,Ih0=0,Nc0=0,where Ih∗=0.


## 3.2.2. Endemic Equilibrium Point

At the endemic equilibrium point, all disease states in equation (14) are considered to be positive, and consequently, Ih∗ must be greater than zero for all the other states to be positive. Therefore, the unique positive endemic equilibrium point exists only when Ih∗ in the system equation (16) is positive, which is also in this case the same as R0>1. And thus, it is the point(19)E∗=Sh∗∗,Eh∗∗,Ih∗∗,Th∗∗,Nm∗∗,Ss∗∗,Es∗∗,Is∗∗,Nc∗∗,where(20)Sh∗∗=Sh∗,Eh∗∗=Eh∗,Ih∗∗=AR02−1B+C,Th∗∗=Th∗,Nm∗∗=Nm∗,Ss∗∗=Ss∗,Es∗∗=Es∗,Is∗∗=Is∗,Nc∗∗=Nc∗,whenever R0>0.


## 3.3. Stability Analysis

### 3.3.1. Basic Reproduction Number,R0

The basic reproduction number is defined as the average number of secondary infections caused by the emergence of an infectious individual into a completely susceptible population [20]. It is denoted by R0. We considered the next generation matrix approach according to Chavez et al. [21] (for details, see Appendix A). The next generation matrix is denoted by G and is obtained as(21)G=0000λ2σs1−ξhβhαhc0μhμcδs+μsσs+δs+μsλ21−ξhβhαhc0μhμcδs+μs1−ξhβhαhc0μhμc000000000000000000000λ1σhβsαsm0μsμmσh+μhγh+δh+μhλ1βsαsm0μsμmγh+δh+μh0βsαsm0μsμm00000000000000000.Accordingly, the basic reproduction number,R0, which is the dominant eigenvalue of G, is obtained as(22)R0=αh1−ξhβhσhαsβsσsλ1λ2c0m0μhμsμmμcσh+μhγh+δh+μhδs+μsσs+δs+μs.For simplicity, sinceR0 denotes the possible number of uninfected a particular infected individual (in this case either human or snail) can infect, we can rewrite R0 as(23)R0=R0HS·R0SH, where R0HS=αh1−ξhβhσhλ1/m0μhμmσh+μhγh+δh+μh represents the number of snails an infected human is capable of infecting through the production of miracidia and R0SH=αsβsσsλ2/c0μsμcδs+μsσs+δs+μs denotes the number of humans an infected snail is capable of infecting through the production of cercariae.


### 3.3.2. Stability of Disease-Free Equilibrium Point

Theorem 2.
The disease-free equilibrium point,E0, is locally asymptotically stable provided R0<1 and unstable otherwise.Proof.
To verify the validity of the above proposition, see Theorem 2 in Van den Driessche and Watmough [22].
The method according to Chavez et al. [21] is used to determine the global stability of the disease-free equilibrium point within Ω∈ℝ+9. Hence, we denote the system equation (1) by(24)dXdt=FX,Y,dYdt=GX,Y,GX,0=0,where X=Sh,Ss comprises of the uninfected subpopulations and Y=Eh,Ih,Th,Nm,Es,Is,Nc comprises of the infected and infectious subpopulations. E0 is guaranteed to be globally asymptotically stable (GAS) if R0<1 and the following two conditions (C1) and (C2) hold:
(C1): FordX/dt=FX,0,X∗=E0 is GAS
(C2):GX,Y=AY−G∗X,Y,G∗X,Y≥0forX,Y∈Ω
whereA=DYGX∗,0 and Ω∈ℝ+9 is the invariant region of the model.Theorem 3.
The disease-free equilibrium point,E0=αh/μh,0,0,0,0,αs/μs,0,0,0∈Ω, is GAS for the system equation (1) provided that R0<1 and conditions (C1) and (C2) both hold.Proof.
From the system of equation (1), we have(25)FX,0=αh−μhShαs−μsSs.
First, suppose that the time,t⟶∞, we need to show that X⟶E0.
But we know that forFX,0,(26)dShdt=αh−μhSh,dShdt+μhSh=αh,Sh=αhμh.
Also,(27)dSsdt=αs−μsSs,dSsdt+μhSs=αs,Ss=αsμs.
Thus, all points with respect to this condition converge at(28)E0=αhμh,0,0,0,0,αsμs,0,0,0∈Ω.
Now, from the system equation (1), we know that(29)dYdt=GX,Y,where(30)GX,Y=ωhShNc−σh+μhEhσhEh−γh+δh+μhIhγhIh−ρ+1−ρδh+μhThλ1Ih−μmNmωsSsNm−σs+δs+μsEsσsEs−δs+μsIsλ2Is−μcNc,and GX,Y=AY−G∗X,Y, where(31)A=−σh+μh00000αh1−ξhβhc0μhσh−γh+δh+μh000000γh−ρ+1−ρδh+μh00000λ10−μm000000αsβsm0μs−σs+δs+μs000000σs−δs+μs000000λ2−μc,G∗X,Y=αhc0μh−Shc0+εNc1−ξhβhNc000αsm0μs−Ssm0+εNmβsNm00.
Clearly,G∗X,Y≥0,∀X,Y∈Ω, since(32)Nh≤αhμh,Ns≤αsμs⇒αhc0μh−Shc0+εNc≥0,αsm0μs−Ssm0+εNm≥0.
Moreover,A is an M-matrix, since the nondiagonal entries are nonnegative [21]. Therefore, the DFE, denoted E0, is globally asymptotically stable provided that R0<1.


### 3.3.3. Stability of the Endemic Equilibrium Point

We consider the concept ofCenter Manifold Theory as described in Theorem 4.1 by Chavez and Song [23]. The theorem helps to determine the course of the bifurcation at a critical value of the parameter [24]. The theorem is reproduced in Appendix B.Theorem 4.
The endemic equilibrium point, based on the Theorem 4.1 by [23] is locally asymptotically stable for R0>1 (but near 1).Proof.
Now, consider the model system equation (1). Using Theorem 4.1 by Chavez and Song [23], let Sh=x1, Eh=x2, Ih=x3, Th=x4, Nm=x5, Ss=x6, Es=x7, Is=x8, and Nc=x9. Further, denote X=x1,x2,x3,x4,x5,x6,x7,x8,x9 and F=f1,f2,f3,f4,f5,f6,f7,f8,f9. Then, the model system equation (1) can be expressed in the form dX/dt=F, where(33)dx1dt=αh+ρx4−ωhx9x1−μhx1=f1,dx2dt=ωhx9x1−σhx2−μhx2=f2,dx3dt=σhx2−γhx3−δhx3−μhx3=f3,dx4dt=γhx3−ρx4−1−ρδh−μhx4=f4,dx5dt=λ1x3−μmx5=f5,dx6dt=αs−ωsx5x6−μsx6=f6,dx7dt=ωsx5x6−σsx7−δsx7−μsx7=f7,dx8dt=σsx7−δsx8−μsx8=f8,dx9dt=λ2x8−μcx9=f9.
The role of the parameterϕ in Theorem 4.1 by Chavez and Song [23] is played by βh with the critical value obtained from R0=1 [7, 24]. Hence, βh=βh∗, regarded as the bifurcation point, can be obtained by evaluating it at R0=1. So(34)βh∗=c0m0μhμsμmμcσh+μhγh+δh+μhδs+μsσs+δs+μsαh1−ξhαsβsσhσsλ1λ2.
Let the Jacobian of the transformed system equation (33) at E0 when βh=βh∗ is denoted by A. Then, it can be shown that 0 is a simple eigenvalue of the matrix A, where(35)A=−μh00ρ0000−αh1−ξhβh∗c0μh0−σh+μh000000αh1−ξhβh∗c0μh0σh−γh+δh+μh00000000γh−ρ+1−ρδh+μh0000000λ10−μm00000000−αsβsm0μs−μs0000000αsβsm0μs0−σsδs+μs00000000σs−δs+μs00000000λ2−μc,and if the vector w=w1,w2,w3,w4,w5,w6,w7,w8,w9T is a right eigenvector corresponding to the eigenvalue, 0 of A, then the components of w can be obtained as(36)w1=−αh1−ξhβh∗σhγh1−ρδh+μh+ρ+1−ρδh+μhσhδh+μh+μhγh+δh+μhc0μh2σh+μhγh+δh+muhρ+1−ρδh+μhw9,w2=αh1−ξhβh∗c0μhσh+μhw9,w3=αh1−ξhβh∗σhc0μhσh+μhγh+δh+μhw9,w4=αh1−ξhβh∗σhγhc0μhσh+μhρ+1−ρδh+μhγh+δh+μhw9,w5=αh1−ξhβh∗σhλ1c0μhμmσh+μhγh+δh+μhw9,w6=−αh1−ξhβh∗αsβsσhλ1m0c0μhμs2μmσh+μhγh+δh+μhw9,w7=αh1−ξhβh∗αsβsσhλ1m0c0μhμsμmσh+μhγh+δh+μhσs+δs+μsw9,w8=αh1−ξhβh∗αsβsσhσsλ1m0c0μhμsμmσh+μhγh+δh+μhdeltas+μsσs+δs+μsw9,w9=w9>0,
whereβh∗ is defined in equation (34).
Similarly, letv=v1,v2,v3,v4,v5,v6,v7,v8,v9 be a left eigenvector corresponding to the zero eigenvalue of A. Then, the components of v are(37)v1=v4=v6=0,v2=c0μhμcαh1−ξhβh∗v9,v3=c0μhμcσh+μhαh1−ξhβh∗v9,v5=c0μhμcσh+μhγh+δh+μhαh1−ξhβh∗λ1v9,v7=m0c0μhμsμmμcσh+μhγh+δh+μhαh1−ξhβh∗αsβsλ1v9,v8=m0c0μhμsμmμcσh+μhγh+δh+μhσs+δs+μsαh1−ξhβh∗αsβsσsλ1v9,v9=v9>0.
Now, to get the values ofa and b as in Theorem 4.1 by Chavez and Song [23], we have to get the second partial derivatives in the system equation (33). Clearly, we do not need the partial derivatives of f1,f4, and f6, since v1=v4=v6=0. Thus, the second (nonzero) partial derivatives from the linearized system equation (33) for a are(38)∂2f2∂x1∂x90,0=1−ξhβh∗c0,∂2f7∂x5∂x60,0=βsm0,and for b is(39)∂2f2∂x9∂βh∗0,0=1−ξhx10c0=αh1−ξhc0μh.
Therefore,(40)a=∑k,i,j=1nvkwiwj∂2fk∂xi∂xj0,0,a=v2w1w9∂2f2∂x1∂x9+v7w5w6∂2f7∂x5∂x6,a=−1−ξhβh∗μcσhγh1−ρδh+μh+ρ+1−ρδh+μhσhδh+μhγh+δh+μhc0μhσh+μhγh+δh+μhρ+1−ρδh+μhv9w92−αh1−ξhβh∗σh2λ1μcc0μhμmσh+μhγh+δh+μhρ+1−ρδh+μh2v9w92.
Soa<0.(41)b=∑k,i=1nvkwi∂2fk∂xi∂ϕ0,0,b=v2w9∂2f2∂x9∂βh∗0,0,b=μcβh∗v9w9,
Sob=μc/βh∗v9w9.
Clearly,a<0 and b>0, and with the help of condition (iv) in Theorem 4.1 by Chavez and Song [23], it is confirmed that the system follows a forward bifurcation at R0=1(see Figure 2), and hence, the endemic equilibrium is locally asymptotically stable for R0>1, but sufficiently closed to 1.Figure 2
A forward bifurcation showing an exchange of stability of the equilibrium points atR0=1.Generally, also based on the principles of the next generation matrix according to Van den Driessche and Watmough [22], the disease-free equilibrium point is stable when R0<1 and unstable when R0>1. To validate these analytic results, a numerical simulation was carried out. Figure 2 shows a graph of R0 against the steady state of Ih, denoted by Ih∗. The figure shows that at R0=1, the model system equation (1) exhibits a forward bifurcation. Forward bifurcations have been observed in many other epidemic models [7, 23, 25]. It shows an exchange of stability between the disease-free equilibrium point and the disease endemic equilibrium point. To this effect, it can be concluded that the model system equation (1) has a disease-free equilibrium point that is globally asymptotically stable when R0<1, and when R0>1, it posses an endemic equilibrium point that is locally asymptotically stable. This is actually a validation of the results from the analytic processes which also confirms that the disease-free equilibrium point is both locally and globally asymptotically stable when R0<1, and when R0>1, a locally asymptotically stable unique endemic equilibrium point exists.Figure 3
Evolution of each subpopulation with time.
(a)
Evolution of the human subpopulation with time(b)
Evolution of the snail subpopulation with time(c)
Evolution of the free living miracidia with time(d)
Evolution of the free living cercariae with timeFigure 4
Effects of WASH on the dynamics of humans and the pathogens.
(a)
Effects of varying WASH on susceptible humans(b)
Effects of varying WASH on the infected humans(c)
Effects of varying WASH on miracidia(d)
Effects of varying WASH on cercariaeFigure 5
Effects of varying treatments on infected humans, miracidia, and cercariae.
(a)
Effects of varying treatments on the infected humans(b)
Effects of varying treatments on miracidia(c)
Effects of treatments on cercariaeFigure 6
Effects of varying sensitive parameters on infected humans.
(a)
Effects of varyingαh on infected humans(b)
Effects of varyingαs on infected humans(c)
Effects of varyingβh on infected humans(d)
Effects of varyingβs on infected humans(e)
Effects of varyingσh on infected humans(f)
Effects of varyingλ1 on infected humans(g)
Effects of varyingλ2 on infected humans


## 3.3.1. Basic Reproduction Number,R0

The basic reproduction number is defined as the average number of secondary infections caused by the emergence of an infectious individual into a completely susceptible population [20]. It is denoted by R0. We considered the next generation matrix approach according to Chavez et al. [21] (for details, see Appendix A). The next generation matrix is denoted by G and is obtained as(21)G=0000λ2σs1−ξhβhαhc0μhμcδs+μsσs+δs+μsλ21−ξhβhαhc0μhμcδs+μs1−ξhβhαhc0μhμc000000000000000000000λ1σhβsαsm0μsμmσh+μhγh+δh+μhλ1βsαsm0μsμmγh+δh+μh0βsαsm0μsμm00000000000000000.Accordingly, the basic reproduction number,R0, which is the dominant eigenvalue of G, is obtained as(22)R0=αh1−ξhβhσhαsβsσsλ1λ2c0m0μhμsμmμcσh+μhγh+δh+μhδs+μsσs+δs+μs.For simplicity, sinceR0 denotes the possible number of uninfected a particular infected individual (in this case either human or snail) can infect, we can rewrite R0 as(23)R0=R0HS·R0SH, where R0HS=αh1−ξhβhσhλ1/m0μhμmσh+μhγh+δh+μh represents the number of snails an infected human is capable of infecting through the production of miracidia and R0SH=αsβsσsλ2/c0μsμcδs+μsσs+δs+μs denotes the number of humans an infected snail is capable of infecting through the production of cercariae.


## 3.3.2. Stability of Disease-Free Equilibrium Point

Theorem 2.
The disease-free equilibrium point,E0, is locally asymptotically stable provided R0<1 and unstable otherwise.Proof.
To verify the validity of the above proposition, see Theorem 2 in Van den Driessche and Watmough [22].
The method according to Chavez et al. [21] is used to determine the global stability of the disease-free equilibrium point within Ω∈ℝ+9. Hence, we denote the system equation (1) by(24)dXdt=FX,Y,dYdt=GX,Y,GX,0=0,where X=Sh,Ss comprises of the uninfected subpopulations and Y=Eh,Ih,Th,Nm,Es,Is,Nc comprises of the infected and infectious subpopulations. E0 is guaranteed to be globally asymptotically stable (GAS) if R0<1 and the following two conditions (C1) and (C2) hold:
(C1): FordX/dt=FX,0,X∗=E0 is GAS
(C2):GX,Y=AY−G∗X,Y,G∗X,Y≥0forX,Y∈Ω
whereA=DYGX∗,0 and Ω∈ℝ+9 is the invariant region of the model.Theorem 3.
The disease-free equilibrium point,E0=αh/μh,0,0,0,0,αs/μs,0,0,0∈Ω, is GAS for the system equation (1) provided that R0<1 and conditions (C1) and (C2) both hold.Proof.
From the system of equation (1), we have(25)FX,0=αh−μhShαs−μsSs.
First, suppose that the time,t⟶∞, we need to show that X⟶E0.
But we know that forFX,0,(26)dShdt=αh−μhSh,dShdt+μhSh=αh,Sh=αhμh.
Also,(27)dSsdt=αs−μsSs,dSsdt+μhSs=αs,Ss=αsμs.
Thus, all points with respect to this condition converge at(28)E0=αhμh,0,0,0,0,αsμs,0,0,0∈Ω.
Now, from the system equation (1), we know that(29)dYdt=GX,Y,where(30)GX,Y=ωhShNc−σh+μhEhσhEh−γh+δh+μhIhγhIh−ρ+1−ρδh+μhThλ1Ih−μmNmωsSsNm−σs+δs+μsEsσsEs−δs+μsIsλ2Is−μcNc,and GX,Y=AY−G∗X,Y, where(31)A=−σh+μh00000αh1−ξhβhc0μhσh−γh+δh+μh000000γh−ρ+1−ρδh+μh00000λ10−μm000000αsβsm0μs−σs+δs+μs000000σs−δs+μs000000λ2−μc,G∗X,Y=αhc0μh−Shc0+εNc1−ξhβhNc000αsm0μs−Ssm0+εNmβsNm00.
Clearly,G∗X,Y≥0,∀X,Y∈Ω, since(32)Nh≤αhμh,Ns≤αsμs⇒αhc0μh−Shc0+εNc≥0,αsm0μs−Ssm0+εNm≥0.
Moreover,A is an M-matrix, since the nondiagonal entries are nonnegative [21]. Therefore, the DFE, denoted E0, is globally asymptotically stable provided that R0<1.


## 3.3.3. Stability of the Endemic Equilibrium Point

We consider the concept ofCenter Manifold Theory as described in Theorem 4.1 by Chavez and Song [23]. The theorem helps to determine the course of the bifurcation at a critical value of the parameter [24]. The theorem is reproduced in Appendix B.Theorem 4.
The endemic equilibrium point, based on the Theorem 4.1 by [23] is locally asymptotically stable for R0>1 (but near 1).Proof.
Now, consider the model system equation (1). Using Theorem 4.1 by Chavez and Song [23], let Sh=x1, Eh=x2, Ih=x3, Th=x4, Nm=x5, Ss=x6, Es=x7, Is=x8, and Nc=x9. Further, denote X=x1,x2,x3,x4,x5,x6,x7,x8,x9 and F=f1,f2,f3,f4,f5,f6,f7,f8,f9. Then, the model system equation (1) can be expressed in the form dX/dt=F, where(33)dx1dt=αh+ρx4−ωhx9x1−μhx1=f1,dx2dt=ωhx9x1−σhx2−μhx2=f2,dx3dt=σhx2−γhx3−δhx3−μhx3=f3,dx4dt=γhx3−ρx4−1−ρδh−μhx4=f4,dx5dt=λ1x3−μmx5=f5,dx6dt=αs−ωsx5x6−μsx6=f6,dx7dt=ωsx5x6−σsx7−δsx7−μsx7=f7,dx8dt=σsx7−δsx8−μsx8=f8,dx9dt=λ2x8−μcx9=f9.
The role of the parameterϕ in Theorem 4.1 by Chavez and Song [23] is played by βh with the critical value obtained from R0=1 [7, 24]. Hence, βh=βh∗, regarded as the bifurcation point, can be obtained by evaluating it at R0=1. So(34)βh∗=c0m0μhμsμmμcσh+μhγh+δh+μhδs+μsσs+δs+μsαh1−ξhαsβsσhσsλ1λ2.
Let the Jacobian of the transformed system equation (33) at E0 when βh=βh∗ is denoted by A. Then, it can be shown that 0 is a simple eigenvalue of the matrix A, where(35)A=−μh00ρ0000−αh1−ξhβh∗c0μh0−σh+μh000000αh1−ξhβh∗c0μh0σh−γh+δh+μh00000000γh−ρ+1−ρδh+μh0000000λ10−μm00000000−αsβsm0μs−μs0000000αsβsm0μs0−σsδs+μs00000000σs−δs+μs00000000λ2−μc,and if the vector w=w1,w2,w3,w4,w5,w6,w7,w8,w9T is a right eigenvector corresponding to the eigenvalue, 0 of A, then the components of w can be obtained as(36)w1=−αh1−ξhβh∗σhγh1−ρδh+μh+ρ+1−ρδh+μhσhδh+μh+μhγh+δh+μhc0μh2σh+μhγh+δh+muhρ+1−ρδh+μhw9,w2=αh1−ξhβh∗c0μhσh+μhw9,w3=αh1−ξhβh∗σhc0μhσh+μhγh+δh+μhw9,w4=αh1−ξhβh∗σhγhc0μhσh+μhρ+1−ρδh+μhγh+δh+μhw9,w5=αh1−ξhβh∗σhλ1c0μhμmσh+μhγh+δh+μhw9,w6=−αh1−ξhβh∗αsβsσhλ1m0c0μhμs2μmσh+μhγh+δh+μhw9,w7=αh1−ξhβh∗αsβsσhλ1m0c0μhμsμmσh+μhγh+δh+μhσs+δs+μsw9,w8=αh1−ξhβh∗αsβsσhσsλ1m0c0μhμsμmσh+μhγh+δh+μhdeltas+μsσs+δs+μsw9,w9=w9>0,
whereβh∗ is defined in equation (34).
Similarly, letv=v1,v2,v3,v4,v5,v6,v7,v8,v9 be a left eigenvector corresponding to the zero eigenvalue of A. Then, the components of v are(37)v1=v4=v6=0,v2=c0μhμcαh1−ξhβh∗v9,v3=c0μhμcσh+μhαh1−ξhβh∗v9,v5=c0μhμcσh+μhγh+δh+μhαh1−ξhβh∗λ1v9,v7=m0c0μhμsμmμcσh+μhγh+δh+μhαh1−ξhβh∗αsβsλ1v9,v8=m0c0μhμsμmμcσh+μhγh+δh+μhσs+δs+μsαh1−ξhβh∗αsβsσsλ1v9,v9=v9>0.
Now, to get the values ofa and b as in Theorem 4.1 by Chavez and Song [23], we have to get the second partial derivatives in the system equation (33). Clearly, we do not need the partial derivatives of f1,f4, and f6, since v1=v4=v6=0. Thus, the second (nonzero) partial derivatives from the linearized system equation (33) for a are(38)∂2f2∂x1∂x90,0=1−ξhβh∗c0,∂2f7∂x5∂x60,0=βsm0,and for b is(39)∂2f2∂x9∂βh∗0,0=1−ξhx10c0=αh1−ξhc0μh.
Therefore,(40)a=∑k,i,j=1nvkwiwj∂2fk∂xi∂xj0,0,a=v2w1w9∂2f2∂x1∂x9+v7w5w6∂2f7∂x5∂x6,a=−1−ξhβh∗μcσhγh1−ρδh+μh+ρ+1−ρδh+μhσhδh+μhγh+δh+μhc0μhσh+μhγh+δh+μhρ+1−ρδh+μhv9w92−αh1−ξhβh∗σh2λ1μcc0μhμmσh+μhγh+δh+μhρ+1−ρδh+μh2v9w92.
Soa<0.(41)b=∑k,i=1nvkwi∂2fk∂xi∂ϕ0,0,b=v2w9∂2f2∂x9∂βh∗0,0,b=μcβh∗v9w9,
Sob=μc/βh∗v9w9.
Clearly,a<0 and b>0, and with the help of condition (iv) in Theorem 4.1 by Chavez and Song [23], it is confirmed that the system follows a forward bifurcation at R0=1(see Figure 2), and hence, the endemic equilibrium is locally asymptotically stable for R0>1, but sufficiently closed to 1.Figure 2
A forward bifurcation showing an exchange of stability of the equilibrium points atR0=1.Generally, also based on the principles of the next generation matrix according to Van den Driessche and Watmough [22], the disease-free equilibrium point is stable when R0<1 and unstable when R0>1. To validate these analytic results, a numerical simulation was carried out. Figure 2 shows a graph of R0 against the steady state of Ih, denoted by Ih∗. The figure shows that at R0=1, the model system equation (1) exhibits a forward bifurcation. Forward bifurcations have been observed in many other epidemic models [7, 23, 25]. It shows an exchange of stability between the disease-free equilibrium point and the disease endemic equilibrium point. To this effect, it can be concluded that the model system equation (1) has a disease-free equilibrium point that is globally asymptotically stable when R0<1, and when R0>1, it posses an endemic equilibrium point that is locally asymptotically stable. This is actually a validation of the results from the analytic processes which also confirms that the disease-free equilibrium point is both locally and globally asymptotically stable when R0<1, and when R0>1, a locally asymptotically stable unique endemic equilibrium point exists.Figure 3
Evolution of each subpopulation with time.
(a)
Evolution of the human subpopulation with time(b)
Evolution of the snail subpopulation with time(c)
Evolution of the free living miracidia with time(d)
Evolution of the free living cercariae with timeFigure 4
Effects of WASH on the dynamics of humans and the pathogens.
(a)
Effects of varying WASH on susceptible humans(b)
Effects of varying WASH on the infected humans(c)
Effects of varying WASH on miracidia(d)
Effects of varying WASH on cercariaeFigure 5
Effects of varying treatments on infected humans, miracidia, and cercariae.
(a)
Effects of varying treatments on the infected humans(b)
Effects of varying treatments on miracidia(c)
Effects of treatments on cercariaeFigure 6
Effects of varying sensitive parameters on infected humans.
(a)
Effects of varyingαh on infected humans(b)
Effects of varyingαs on infected humans(c)
Effects of varyingβh on infected humans(d)
Effects of varyingβs on infected humans(e)
Effects of varyingσh on infected humans(f)
Effects of varyingλ1 on infected humans(g)
Effects of varyingλ2 on infected humans


## 4. Numerical Results

Table1 gives the parameter values used in the numerical simulation of the model equation (1). Many of the parameters are obtained from existing literature and are rightly cited while some are taken from real life data. Other parameters are estimated based on convenience using what is generally known about the dynamics of the studied population and the disease. In particular, the model is designed to study the dynamics of schistosomiasis within the Gambian population. In 2020, The Gambia has a population of 2,486,945, with a birth rate of 37.272 per 1,000 population and life expectancy of up to 62.57 years. Thus, αh=37.272×2,486,945/1,000×365≈254 and μh=1/62.57×365≈0.00004379. The human exposed period for schistosomiasis, i.e., the time period of humans from having contact with free living cercariae to the period eggs begin to appear in the human feces ranges from 5 to 6 weeks [2]. Hence, σh is estimated to lie between 1/6×7≈0.0238 and 1/5×7≈0.0286. The incubation period for snails ranges from 4 weeks to 5 weeks [2]. Thus, σs lies between 1/5×7≈0.0286 and 1/4×7≈0.0357. The parasite can live within a definite human host for a period ranging from 3 years to 10 years [26]. So, δh is varied between values corresponding to 3 and 10 years. Consequently, δh lie between 1/10×365≈0.000274 and 1/3×365≈0.000913.Usingode45 in MATLAB, a numerical simulation of the model equation (1) was carried out with the help of the parameter values presented in Table 1, and the results illustrated graphically show the behavior of the various subpopulations: the human subpopulation (Figure 3(a)), the snail subpopulation (Figure 3(b)), and the subpopulations of the free living miracidia (Figure 3(c)) and cercariae (Figure 3(d)), against time (days).It can be observed from Figure3(a) that the susceptible humans have risen within the first few days but gradually fall afterwards as they continuously come in contact with free living cercariae and later became stable. We experienced a steady and then a quick rise in the exposed and infected human compartments due to the introduction of the disease, both of which later had a slow rise as the disease continue to spread. Equally, the treated human compartment shows a steady rise within the first few days and, afterwards, gains stability. The susceptible snail and the exposed snail compartments, from Figure 3(b), had a little rise within the first few days but become stable shortly and then fall gradually afterwards due to more of the susceptible snails getting in contact with free living miracidia and the exposed snails shedding cercariae. However, there is a high rise in the number of infected snails whose population became stable later on. The infected snails rise above both the susceptible and exposed snails causing the number of susceptible and exposed snails to fall gradually as a result of the fact that infected snails do not recover. Hence, it demonstrates the dominance of the cercariae shedding snails over all other snail compartments. The free living miracidia, from Figure 3(c), also rises within the first few days but reaches a peak and slowly changed course but continues to rise slowly as more infected individuals continue to release eggs for the production of miracidia. The dominance of the cercariae shedding snails explains why more and more cercariae are produced as seen in Figure 3(d). The cercariae subpopulation rises very quickly and gain stability after a number of days, and this can be interpreted to be the result of the dominance of the infected snails over all other snail compartments.


### 4.1. Effects of WASH

Here, we examined the effectiveness of the application of WASH on the dynamics of the human subpopulation and the dynamics of the free living miracidia and cercariae. It can be observed from the subplots in Figure4 that WASH as a preventive mechanism in controlling the spread of schistosomiasis is effective if nearly the entire human subpopulation practices it. The results indicate that if more and more individuals within the studied population are encouraged to practice WASH, it will mitigate the number of people contracting the disease, thereby lowering the number of infected individuals and, consequently, reducing the production rate of miracidia. However, it may effectively slow down the rate of production of cercariae, but it will not stop the number of cercariae produced from reaching the point of endemicity as infected snails will continue to produce free living cercariae. Thus, practicing WASH alongside water chlorination (which is considered to be effective in killing cercariae) should be considered.


### 4.2. Effects of Treatment

Here, we look into how increase in the number of treated infected humans can affect the number of infected humans and also the production rates of the pathogens. It can be observed from Figure5(a) that increasing the rate of treatment of humans will be efficient in reducing the number of infected humans. Further, treatment can help in reducing the production rate of miracidia as seen in Figure 5(b). However, it does not seem to have an impact in reducing the production rate of cercariae, and thus, by not combining it with preventive strategies like WASH, individuals will continue to get infected or reinfected (upon treatment) (see Figure 5(c)). This is an indication that the effectiveness of treatment is a short-term measure. As a result, the susceptible and treated individuals should be educated on WASH procedures to effectively curtail the spread of the disease.


## 4.1. Effects of WASH

Here, we examined the effectiveness of the application of WASH on the dynamics of the human subpopulation and the dynamics of the free living miracidia and cercariae. It can be observed from the subplots in Figure4 that WASH as a preventive mechanism in controlling the spread of schistosomiasis is effective if nearly the entire human subpopulation practices it. The results indicate that if more and more individuals within the studied population are encouraged to practice WASH, it will mitigate the number of people contracting the disease, thereby lowering the number of infected individuals and, consequently, reducing the production rate of miracidia. However, it may effectively slow down the rate of production of cercariae, but it will not stop the number of cercariae produced from reaching the point of endemicity as infected snails will continue to produce free living cercariae. Thus, practicing WASH alongside water chlorination (which is considered to be effective in killing cercariae) should be considered.


## 4.2. Effects of Treatment

Here, we look into how increase in the number of treated infected humans can affect the number of infected humans and also the production rates of the pathogens. It can be observed from Figure5(a) that increasing the rate of treatment of humans will be efficient in reducing the number of infected humans. Further, treatment can help in reducing the production rate of miracidia as seen in Figure 5(b). However, it does not seem to have an impact in reducing the production rate of cercariae, and thus, by not combining it with preventive strategies like WASH, individuals will continue to get infected or reinfected (upon treatment) (see Figure 5(c)). This is an indication that the effectiveness of treatment is a short-term measure. As a result, the susceptible and treated individuals should be educated on WASH procedures to effectively curtail the spread of the disease.


## 5. Sensitivity Analysis

Sensitivity indices permit us to estimate the perturbations in a variable with respect to changes in a parameter. The normalized forward sensitivity index of a variable in relation to a parameter is the proportion of the perturbations in the variable with respect to changes in the given parameter. Suppose that a variable is partially differentiable with respect to the parameter, then the sensitivity index can be defined using partial derivatives.Definition 5.
The normalized forward sensitivity index ofR0, which depends differentiably on a parameter ε, is defined by(42)SεR0=∂R0∂ε×εR0.
Table2 contains all model parameters that are partially differentiable with respect to R0, their values, and sensitivity indices.Table 2
A table showing the sensitivity index of each model parameter onR0.
ParameterParameter valueSensitivity indexαh254+0.5000αs3000+0.5000βh0.09753+0.5000βs0.615+0.5000σh0.0236+0.00092σs0.0286+0.0164λ16.96+0.5000λ22.6+0.5000γh0.03-0.4948μh0.00004379-0.5016μs0.000569-0.8029μm0.9-0.5000μc0.004-0.5000δh0.000274-0.0045δs0.0004012-0.2135


### 5.1. Interpretation of the Sensitivity Indices onR0

The sensitivity indices for the model parameters to the basic reproduction number,R0, are presented in Table 2. The most sensitive parameters include the recruitment rates of humans and snails and the production rate of free living miracidia and cercariae as produced by the infected humans and infected snails, respectively. In general, the results indicate that the parameters αh, αs, βh, βs, σh, σs, λ1, and λ2 have positive indices and, therefore, are each directly proportional to the value of R0. Thus, increasing the value of any of these parameters will lead to the disease remaining endemic within the studied population and vice versa. However, the parameters γh, δh, δs, μh, μs, μm, and μc have each a negative sensitivity index, which means that each of them varies inversely as the value of R0. Therefore, increasing the value of any of these parameters with negative index, while holding all other parameters constant, will reduce the value of R0 and, hence, contribute to the elimination of the disease and vice versa. For example, increasing the recruitment rate of either humans or snails by 10% will lead to a 5% increase on R0, and increasing the rate of treatment by 10% will result to a reduction of R0 by 4.948%.


### 5.2. Effects of Parameters on Infected Humans

Figure6 shows the contributions of the parameters αh, αs, βh, βs, σh, λ1, and λ2 on the infected human subpopulation.Figures6(a) and 6(b) indicate that reducing the recruitment rates of either humans or (most importantly) snails will significantly reduce the number of humans getting infected. The effects of reducing the forces of infection as seen in Figures 6(c) and 6(d) on the infected humans are an indication that vaccination or provision of quality water and improved sanitation system are highly essential in containing the spread of the disease. Figure 6(e) shows that reducing the rate exposed individuals become infected will significantly reduce the number of individuals that become infected eventually. This can be done through provision of early diagnostic materials that will be able to detect the disease in an exposed individual before the production of eggs and also invention of drugs that can treat the disease at the very onset of infection. Figures 6(f) and 6(g) show that reducing the production rate of either miracidia or cercariae will drastically reduce the number of infected individuals. Miracidia production can be reduced by treatment of infected individuals (see Figure 5(b)) or practice of WASH (see Figure 4(c)) or possibly through the introduction of vaccines for humans. Although miracidia can be reduced by treatment of infected individuals, the high possibility of them getting reinfected will make it very difficult to control it unless the individuals are vaccinated or sensitized on the key preventive mechanisms. Cercariae production, on the other hand, can be reduced through chlorination of fresh waters or snail control using chemical molluscicides.


## 5.1. Interpretation of the Sensitivity Indices onR0

The sensitivity indices for the model parameters to the basic reproduction number,R0, are presented in Table 2. The most sensitive parameters include the recruitment rates of humans and snails and the production rate of free living miracidia and cercariae as produced by the infected humans and infected snails, respectively. In general, the results indicate that the parameters αh, αs, βh, βs, σh, σs, λ1, and λ2 have positive indices and, therefore, are each directly proportional to the value of R0. Thus, increasing the value of any of these parameters will lead to the disease remaining endemic within the studied population and vice versa. However, the parameters γh, δh, δs, μh, μs, μm, and μc have each a negative sensitivity index, which means that each of them varies inversely as the value of R0. Therefore, increasing the value of any of these parameters with negative index, while holding all other parameters constant, will reduce the value of R0 and, hence, contribute to the elimination of the disease and vice versa. For example, increasing the recruitment rate of either humans or snails by 10% will lead to a 5% increase on R0, and increasing the rate of treatment by 10% will result to a reduction of R0 by 4.948%.


## 5.2. Effects of Parameters on Infected Humans

Figure6 shows the contributions of the parameters αh, αs, βh, βs, σh, λ1, and λ2 on the infected human subpopulation.Figures6(a) and 6(b) indicate that reducing the recruitment rates of either humans or (most importantly) snails will significantly reduce the number of humans getting infected. The effects of reducing the forces of infection as seen in Figures 6(c) and 6(d) on the infected humans are an indication that vaccination or provision of quality water and improved sanitation system are highly essential in containing the spread of the disease. Figure 6(e) shows that reducing the rate exposed individuals become infected will significantly reduce the number of individuals that become infected eventually. This can be done through provision of early diagnostic materials that will be able to detect the disease in an exposed individual before the production of eggs and also invention of drugs that can treat the disease at the very onset of infection. Figures 6(f) and 6(g) show that reducing the production rate of either miracidia or cercariae will drastically reduce the number of infected individuals. Miracidia production can be reduced by treatment of infected individuals (see Figure 5(b)) or practice of WASH (see Figure 4(c)) or possibly through the introduction of vaccines for humans. Although miracidia can be reduced by treatment of infected individuals, the high possibility of them getting reinfected will make it very difficult to control it unless the individuals are vaccinated or sensitized on the key preventive mechanisms. Cercariae production, on the other hand, can be reduced through chlorination of fresh waters or snail control using chemical molluscicides.


## 6. Conclusion

A mathematical model that explains the transmission dynamics of the disease involving humans, snails, and the dynamics of the free living miracidia and the free living cercariae has been presented. Attention is paid to controls within the human subpopulation by incorporating the treated human compartment and preventive mechanism through water sanitation and hygiene (WASH). The equilibrium points were obtained and analyzed for both local and global stabilities. Considering the basic reproduction number,R0, as the threshold and given an R0<1, it is established that the disease-free equilibrium point is both locally and globally asymptotically stable. The existence of a unique endemic equilibrium point when R0>1 was established and using the Center Manifold Theory, it is confirmed that the endemic equilibrium is locally asymptotically stable when R0>1 but close to 1. To validate these analytic results, a numerical simulation was carried out. The results showed that the model system follows a forward bifurcation at R0=1, and hence, the disease-free equilibrium is globally asymptotically stable with R0<1, and the endemic equilibrium is locally asymptotically stable when R0>1.The effectiveness of control by WASH and treatment of infected individuals was analyzed numerically. The results indicated that a combination of both controls is required for the proper eradication of the disease. However, WASH is a more effective control mechanism especially if it is exercised by the entire population.Using the normalized sensitivity index ofR0, we examined the contribution of the model parameters on the basic reproduction number, R0. The results demonstrate that the production rates of the free living miracidia and the free living cercariae are very effective contributors to R0 and hence very essential in the spread and control of the disease. Thus, control mechanisms that can significantly mitigate the production of either pathogen will most definitely reduce the disease endemicity. It was shown that this could be achieved through the application of the combined control measures, i.e., treatment of infected individuals and WASH. Hence, infected individuals should be treated, and at the same time, communities should be educated regarding best practices on schistosomiasis prevention. Adequate provision of clean water and proper sanitation facilities will curtail the spread of the disease.



---

*Source: 6653796-2021-05-19.xml*