A Schistosomiasis Model with Diffusion Effects

In this paper, we propose a schistosomiasis model in which two human groups share the water contaminated by schistosomiasis and migrate each other. The dynamical behavior of the model is studied. By calculation, the threshold value is given, which determines whether the disease will be extinct or not. The existence and global stability of the parasite-free equilibrium and the locally stability of the endemic equilibrium are discussed. Numerical simulations indicate that the diffusion from the mild endemic village to severe endemic village is benefit to control schistosomiasis transmission; otherwise it is bad for the disease control.


Introduction
Schistosomiasis is frequently a serious health problem, which was first described by Theodor Bilharz in 1851, after whom the disease was initially named bilharzia [1].The WHO has recently identified schistosomiasis as the second most important human parasitic disease in the world, after malaria [2].The infection is endemic in approximately 70 countries with about 200 million people affected worldwide [3], and resulting in about 200,000 deaths annually [4].Despite major advances in its control that have lead to substantial decreases in morbidity and mortality, schistosomiasis continues to spread to new geographic areas [5].Although significant progress has been made in chemotherapy with safer and more effective drugs, these cannot prevent the high reinfection rates of schistosomes, and there have been dramatic recurrences in both its prevalence and associated morbidity [6].
During their complex developmental cycle, schistosomes alternate between a mammalian host and a snail host through the medium of fresh water.Mammals are infected by free-swimming larval forms of the parasite called cercariae.These larvae enter through the skin, and mature through different larval stages while circulating through the blood to the lungs before entering the hepatic portal system as mature males and females.They release thousands of eggs daily, which are discharged in the faeces after a damaging passage through the intestinal wall.Once into the fresh water, the eggs hatch and produce free-swimming miracidia, which infect amphibious snails from the genus Oncomelania.The miracidia reproduce asexually through sporocyst stages within these intermediate hosts, resulting in the production of many free-swimming cercariae [7]- [10].MacDonald (1965) was the first to use simple mathematical models to study the transmission dynamics of schistosomiasis [11].The earliest models of schistosomiasis described the population sizes of both humans and snails to be constant [11] [12].In [11] [13] [14], authors considered that models were based on describing the dynamics of transmission between man and snails.Previous several models focused on the interactions between one group of human hosts and schistosomes in a contaminated water resource(for example [15] [16]).However, in realistic situations, the contaminated water might be shared by several human groups.In [15], Feng et al. proposed a model that described the disease dynamics involved two migrated human groups.They also analyzed the mathematical properties of the systems.Meanwhile, they established models with multiple human groups and found some structurally similarities between the models involved two human groups and those involved n groups.
Incidence rate plays an important role in the modeling of epidemic dynamics.In many epidemic models, the bilinear incidence rate SI β and the standard incidence rate SI N β are frequently used.The saturated incidence rate ( ) β implicits the infection force of the schistosomiasis and ( ) describes the psychological effect or inhibition effect from the behavioral change of the susceptible individuals with the increase of the infective individuals.It seems more reasonable than the bilinear incidence rate SI β , and it is a good approximation if the number of available partners is large enough and everybody could not make more contacts than is practically feasible, and includes the behavioral change and crowding effect of the infective individuals and prevents the unboundedness of the contact rate [17]- [19].In this paper, we develop a new mathematical model with saturated incidence function and diffusion effect.In many literatures [20]- [22], the diffusion effect is studied.Numerical simulations demonstrate that the diffusion effect is an important parameters for epidemic transmission or species survival.
In order to keep the model manageable, Feng et al. assumed that the disease-induced death rate of snails 0 s d = in [10].Previous studies suggested that the disease-induced death rate of snails s d was an important parameter in the study of population dynamics [23].In this paper, we investigate firstly a schitosomiasis model with saturated incidence and diffusion effect, in which the disease-induced death rate of snails s d is taken into consideration.Further, by the spectral radius theory, we get the threshold value 0 R , below which the parasites die out, and above which the disease persists.When the threshold 0 1 R = , we consider that the model may produce a bifurcation.And we study that exchange of stability between disease-free and endemic equilibria at bifurcation point.
This paper is organized as follows.In Section 2, we introduce model formulation.In Section 3, we analyze equilibria states of model.The basic reproduction number of the model is determined and the stability of the equilibria is studied.Numerical simulations and control strategies are presented in Section 4. Finally, we summarize and discuss the results in Section 5.

Model Formulation
In [16] where N, P, S, I, C denote the numbers of human hosts living in village, adult parasites that are hosted by human hosts in village, uninfected snails, infected snails and free-living cercaria, respectively.τ is infection-age, and ( ) , x t τ is the infection-age density of snails at time t.k is the clumping parameter which determines the degree of over-dispersion in the negative binomial distribution.The following parameters is used in system (1), all of them positive, h Λ is the recruitment rate of human hosts; s Λ is the recruitment rate of snails; h µ is the per capita natural death rate of human hosts; p µ is the per capita death rate of adult parasites; α is the disease-induced death rate of humans per parasite; σ is the effective treatment rate of human hosts; s µ is the per capita natural death rate of snails; s d is the disease-induced death rate of snails; ξ is the per capita (successful) rate of infection of snails by miracidia produced by one pair of adult parasites; β is the per capita (successful) rate of infection of humans by one cercaria; ( ) r τ is the releasing rate of cercariae, when the infection age is τ .
In [15], Feng el at.considered two neighboring villages sharing the same contaminated water resource and migrated between these two villages, and proposed the following model which based on the system (1).

(
) ( ) ( ) where , is the recruitment rate of human hosts of village i and , i j m is the immigration rate of human hosts from village i to village j, ( ) , 1, 2 i j = . For system (2), Feng el at.made the following assumptions: 1) the snails do not move; 2) the parasites are overdispersed; 3) they have negative binomial distributions among human hosts with clumping parameters i k ; 4) the releasing rate of cercariae is infection-age independent, i.e., ( ) r τ γ ≡ .Thus, ( ) ( ) In system (2), authors introduced the bilinear incidence rate NI β .Whereas the number of uninfected snails is limited within a certain time which contacted by the adult parasites.So the saturated incidence may be more suitable for the realistic situation.The following new model with the saturated incidence function is derived: where  is limitation of the growth velocity of infection of snails.In a contaminated water resource, many people are infected, which develops into chronic disease if not treated.Current control programs primarily focus on chemotherapy with Praziquantel, it is a new drug that is very effective, they can almost kill the adult parasites which reside within the patient.Thus, the disease-induced death rate of human hosts α is very small.For analysing the properties of the model, we let 0 α = .Then the first two equations become: The equilibrium points are obtained by setting the right-hand side of system (4) to zero, we solve the following system of equations: The unique solution of system ( 5) is ( ) , N N , which is globally asymptotically stable, where , with , 1, 2 i j = and i j ≠ .Therefore, we have the following four-dimensional limit system of system (3) which summarizes the above result.
( ) The existence and the uniqueness of solutions of system ( 6) can be proved by using standard methods (see, for example, [24]).

Equilibrium States
In this section, the equilibrium states of system (6) are discussed.The system (6) admits two steady states.We establish sufficient condition for the globally asymptotic stable of infection-free solution and for the permanence of the system (6).

Boundedness
The model (6) describes the dynamics of adult parasites and snail.It is important to prove that these populations are positive and bounded for 0 t > with any positive initial data.So we have the following results.
, , , P t P t S t I t is any solution of system (6), and ( ) ( ) From the first equation of system (6), we have It follows from the first and second equations of ( 6) and ( 7) that  is a ultimately upper bound of 1 P and 2 P , respectively.The proof is completed.The equilibrium states of the basic model are obtained by setting the right-hand side of system (6) to zero.The system (6) has two steady states of the disease-free equilibrium 0 E and the endemic equilibrium * E .

The Disease-Free Equilibrium
At the disease-free state, there is no adult parasitrs and infected snails and hence no infection in the host and the intermediate host.Thus, the system (6) has a disease-free equilibrium ( ) ( ) , , , 0, 0, , 0 , E P P S I S = = where 0  .
s s S µ = Λ In many epidemic models, the basic reproductive number 0 R is a key parameter.It refers to the expected number of secondary infections during the entire period of infectiousness in a completely susceptible population [25].Following the idea in [26], we give the basic reproductive number for system (6).Rewrite system (6) as following form: , , Y P P I = . S denotes the number of uninfected snails, while components of Y represent the number of adult parasites that are hosted by human hosts in Village ( ) , and infected snails, respectively.Following the symbol in [26], we compute matrixes A, M and D as ( ) ( ) where A F V = − .Obviously, 0 F ≥ and 0 V > is a diagonal matrix.The basic reproductive number is the spectral radius (dominant eigenvalue) of the matrix 1 FV − , that is, ( ) where ( ) .
We know that 1 R presents the schistosomiasis transmission coefficient in village 1, and 2 R represents the schistosomiasis transmission coefficient in village 2.
From above discussion, we have following result.Theorem 3. The disease-free equilibrium point 0 E is locally asymptotically stable if 0 1 R < and unstable if 0 1 R > .
Next, we give two conditions which guarantee the global asymptotic stability of the disease-free state.
(H1) For ( ) , , For system (6), we have , 1 and A is given in (8).It is clear that ( ) It is easy to see that the conditions (H1) and (H2) hold.According to the result of literature [26], we have the following result.
Theorem 4. The disease-free equilibrium ( ) is globally asymptotically stable provided that 0 1 R < and the assumptions (H1) and (H2) are satisfied.

The Endemic Equilibrium
First, we show the existence of the unique endemic equilibrium ( ) , , , E P P S I = when 0 1 R > .Ex- pressing in terms of * I , we can derive from system (6) as follows.
By solving (10) for * I we get one of the solutions as * 0 I = which corresponds to the disease-free equilibrium.For 0 1 R > implies that 0 B > .Since 0 C > , then the endemic equilibrium exists.The results of the existence of the endemic equilibrium of system ( 6) can be summarized in the following lemma.
Lemma 5.The system (6) always has a disease-free equilibrium and a unique endemic equilibrium when Center Manifold Theory [19] has been used to determine the local stability of a nonhyperbolic equilibrium, we now employ the Center Manifold Theory to establish the local asymptotic stability of the endemic equilibrium.In order to apply the Center Manifold Theory, we make the following change of variables.Let = .Now we use the vector notation ( ) , , X x x x x = .Then the system ( 6) is written in following form ( ) Evaluating the Jacobian matrix of system (11) at the disease-free equilibrium, it can be shown that the reproduction number is β as the bifurcation parameter.Considering the case 0 1 R = and solving for 1 β , we get ( ) .
We notice that the linearized system (11) of the transformed equation with , , , u u u u u = , where ( ) The left eigenvector of ( ) We now use the following lemma whose proof is found in [27].
Lemma 6.Consider the following general system of ordinary differential equations with a parameter φ , where 0 is an equilibrium of the system, that is ( ) 0, 0 f φ = for all φ and assume A1: is the linearization of system (12) around the equilibrium 0 with φ evaluated at 0. Zero is a simple eigenvalue of A and other eigenvalues of A have negative real parts; A2: Matrix A has a right eigenvector u and a left eigenvector v corresponding to the zero eigenvalue.Let k f be the kth component of f and , , 1 , 1 0, 0 , 0, 0 , The local dynamics of ( 12) around 0 are totally governed by a and b. 1) 0, 0 a b > > .when 0 φ < with 1 φ  , 0 is locally asymptotically stable, and there exists a positive unstable equilibrium; when 0 1 φ <  , 0 is unstable and there exists a negative and locally asymptotically stable equilibrium; 2) 0, 0 a b < < .when 0 φ < with 1 φ  , 0 is unstable; when 0 1 φ <  , 0 is locally asymptotically stable, and there exists a positive unstable equilibrium; 3) 0, 0 a b > < .when 0 φ < with 1 φ  , 0 is unstable, and there exists a locally asymptotically stable negative equilibrium; when 0 1 φ <  , 0 is stable and a positive unstable equilibrium appears; 4) 0, 0 a b < > .When φ changes from negative to positive, 0 changes its stability from stable to unstable.
Correspondingly a negative unstable equilibrium becomes positive and locally asymptotically stable.We now compute a and b, for system (11), the associated non-zero partial derivatives of  at the disease free equilibrium ( ) Substituting the above expressions into (13), we get


For the sign of b, it is associated with the following non-vanishing partial derivatives of  , Thus, 0 a < and 0 b > .According to Lemma 6, item (iv), we can yield the following result which only holds for 0 1 R > , but close to 1. Theorem 7. The unique endemic equilibrium * E is locally asymptotically stable for 0 R near 1.In summary, model ( 6) has a disease-free equilibrium which is globally asymptotically stable when 0 1 R < , and a unique endemic equilibrium point when 0 1 R > .The unique endemic equilibrium is locally asymptotically stable at least near 0 1 R = .We use numerical simulations to show the existence and stability of endemic equilibrium.

Numerical Simulations and Control Strategies
In this section, in order to understand our results more intuitively, some numerical simulations of system (6) that support and extend the conclusions of previous sections are carried out.We use year as unit of time, and choose the parameters 1 0.0005 β = In Figure 1, we show the relationship between the threshold 0 R and adult parasites 1 2 P P + for the mathematical model (6).It is easy to see that 0 1 R = is a bifurcation point, and the adult parasites 1 2 P P + are stable eventually, when the threshold 0 R increases.Otherwise, the adult parasites 1 2 P P + are extinct.It implies that the threshold 0 R is greater than unit, the schistosomiasis will be endemic.Figure 1 and Figure 2 show that if the threshold 0 R is less than unit, the schistosomiasis will be extinct.To see the relative effect of migration in each village, we plot the curved surface of the relationship between 1,2 m , 2,1 m and 0 R .From Figure 3 and Figure 4, we can observe that 0 R decreases dramatically when 1,2 m increases and 2,1 m is fixed with a small number, and 0 R increases sharply when 2,1 m increases and 1,2 m is fixed with a small number.This implies that the migration from severe endemic village to mild endemic village is bad for disease control.
In Figure 5, we consider the infection rates 1 β , and 2 β as the control factors.We plot the curved surface of the threshold 0 R as a function of β and 2 β .We observe that the threshold 0 R decreases dramatically when 2 β and 1 β decrease.It means that decreasing infection rates is helpful to prevent schistosomiasis transmission.

Conclusion and Discussion
As a kind of the tropical diseases, schistosomiasis continues to be a significant public health threat in the world.Following the pioneering work of Feng et al. [16] on modeling schistosomiasis, we establish and analyzed a schistosomiasis model with diffusion effect and saturated incidence function, in which two groups of human share the water contaminated by schistosomiasis and migrate each other.we derived the basic reproduction number 0 R and proved that the disease-free equilibrium is globally asymptotically stable when 0 1 R < , and the unique endemic equilibrium is locally asymptotically stable for 0 R is larger than 1 and near 1.Our results indicate that the diffusion rates and the infection rates play an important role in the determination of the permanence and extinction of schistosomiasis.The diffusion from the mild endemic village to severe endemic village is benefit to control schistosomiasis transmission.
In realistic situations, there might be several human groups sharing the contaminated water resource.Only considering the model with two human groups is insufficient, we expect a similar to work in higher-dimensional systems with n human groups and migration.It can be guessed that the model with n human groups has similar mathematical properties to two human groups.

=
, has a simple zero eigenvalue.Hence, Center Manifold Theory can be used to analyze the dynamics of (13Theorem 4.1 in Castillo-Chavez and Song[27], it can be shown that the Jacobian matrix at

Figure 1 .
Figure 1.The relationship between the threshold 0 R and 1 2 P P + for system (6).

Figure 3 .
Figure 3.It shows sensitive figure that the relationship between the threshold 0 R and migrated rate 1,2 m , 2,1 m .

Figure 5 .
Figure 5.It shows sensitive figure that the relationship between the threshold 0 R and migrated rate 1 β , 2 β .
, Feng et al. proposed a schistosomiasis model with age dependence: