Applied Mathematics, 2013, 4, 1682-1693 Published Online December 2013 (http://www.scirp.org/journal/am) http://dx.doi.org/10.4236/am.2013.412229
Deterministic and Stochastic Schistosomiasis Models with General Incidence Stanislas Ouaro, Ali Traoré Laboratoire d’Analyse Mathématique des Equations (LAME), Unité de Formation et de Recherche en Sciences Exactes et Appliquées, Département de Mathématiques, Ouagadougou, Burkina Faso Email:
[email protected],
[email protected],
[email protected] Received September 11, 2013; revised October 11, 2013; accepted October 18, 2013 Copyright © 2013 Stanislas Ouaro, Ali Traoré. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
ABSTRACT In this paper, deterministic and stochastic models for schistosomiasis involving four sub-populations are developed. Conditions are given under which system exhibits thresholds behavior. The disease-free equilibrium is globally asymptotically stable if R 0 1 and the unique endemic equilibrium is globally asymptotically stable when R 0 1 . The populations are computationally simulated under various conditions. Comparisons are made between the deterministic and the stochastic model. Keywords: Computational Simulation; General Incidence; Reproduction Number; Schistosomiasis Model; Stochastic Differential Equation
1. Introduction Schistosomiasis (Bilharzia or snail fever) remains a serious public health problem in the tropics and subtropics. There are about 240 million people worldwide infected and more than 700 million people live in endemic areas [1]. It probably originated around the Great Lakes of Central Africa and possibly spread to other parts of Africa, the west Indies and South America [2]. Human schistosomiasis is a parasitic disease caused by five species of the genus Schistosoma or flatworms: Schistosoma mansoni, Schistosoma intercalatum, Schistosoma japonicum, Schistosoma mekongi and Schistosoma haematobium. The three most widespread are Schistosoma japonicum, Schistosoma haematobium and Schistosoma mansoni. The large number of parasite eggs retained in the host tissues rather than shed in stools can cause the inflammatory reactions. Indeed, a study conducted in Sudan indicated that schistosomiasis depleted the blood hemoglobin levels to the extent that oxygen flowed to the muscles and the brain and impaired physical activity [3]. Several mathematical models for schistosomiasis disease have been done ([4-15] and the references therein). Stochastic differential equation (SDE) model is a natural generalization of ordinary differential equation (ODE) model. SDE models became increasingly more popular Open Access
in mathematical biology ([4,16,17] and the references therein). In [4], a SDE model for transmission of schistosomiasis was analyzed. That model assumes that human births and deaths are neglected. So, the computational work is involved in a computation of B that one requires other schemes in which we solve an initial value problem. In this paper, we derive an equivalent stochastic model for schistosomiasis infection which involves four subpopulations. This model allows the recruitments and the natural deaths of human. Additionally, the incidence forms are taken as general as possible. The paper is structured as follows: In Section 2, we present the schistosomiasis model while Section 3 provides the basic properties of the model. Section 4 provides an analysis of the full model. We construct a stochastic differential schistosomiasis model in Section 5. In Section 6, we derive an equivalent stochastic model. We describe a numerical method to solve the equivalent stochastic model in Section 7. Computational simulations are performed in Section 8 and finally, in the last section, we end with a conclusion.
2. Mathematical Model Here we consider a relatively isolated community where there are no immigration or emigration. Additionaly, we AM
S. OUARO, A. TRAORÉ
suppose that all newborns are susceptible and the infection does not imply death or isolation of the different hosts. Then, according to Figure 1, we obtain the following system of four differential equations dH s dt bh d h H s f Si , H s H i , dH i f S , H d H H , i s h i i dt dS s b d S g H , S , s s s i s dt dS i g H i , S s d s Si , dt
1683
is no infected in the human and snail populations, then the incidence functions are equal to zero. The incidence functions are also equal to zero when there is no susceptible in the human and snail populations.
3. Model Basic Properties
(1)
In this section, we study the basic properties of the solutions of system (1). Theorem 3.1 [18] Let L : n be a differentiable function and let X x be a vector field. Let a and let G x n : L x a be such that L x 0 for all x L1 a x n : L x a . If X x L x 0 for all x L1 a then, the set G is positively invariant. Lemma 3.2 The system Equation (1) preserve the positivity of solutions. proof. We start by proving that H s , H i , Ss , Si 4 ; H i 0 is positively invariant. For that, let x H s , H i , S s , Si and L x H i . L x is differentiable and L x 0, 1, 0, 0 0 for all x L1 0 x 4 : L x 0 . The vector field on the set H s , H i , Ss , Si 4 ; H i 0 is
where H s t is the size of the susceptible human population; H i t is the size of the infected human population; S s t is the size of the susceptible snails; Si t is the size of the infected snails; bh is the recruitment rate of human hosts; bs is the recruitment rate of snails; d h is the per capita natural death rate of human hosts; d s is the per capita natural death rate of snails; is the treatment rate of human hosts; f is the infection function of susceptible human due to infected snails; g is the infection function of susceptible snails due to infected human. We assume that the functions f and g satisfies the following assumptions. H1 f and g are nonnegative C1 functions on the nonnegative orthant. H2 For all H s , H i , S s , Si 4 , f 0, H s f Si , 0 0 and g 0, S s g H i , 0 0 . Remark 2.1 f and g are two incidence functions which explain the contact between two species. Therefore, f and g are nonnegative. Note also that when there
bh f Si , H s d h H s f Si , H s . X x bs d s S s d s Si
(2)
We obtain X x L x f Si , H s 0 . Then, using Theorem 3.1, we conclude that H s , H i , Ss , Si 4 ; H i 0 is positively invariant. Similarly we prove that H s , H i , S s , Si 4 ; H s 0 , H s , H i , Ss , Si 4 ; Ss 0 and H s , H i , Ss , Si 4 ; Si 0 are positively invariant Lemma 3.3 Each nonnegative solution of model system (1) is bounded.
Figure 1. Transfer diagram. Open Access
AM
S. OUARO, A. TRAORÉ
1684
Proof. Let H t H s t H i t and S t S s t Si t . Then, adding the first two equations and the second two Equations of (1) respectively, we get H bh d h H and S bs d s S . (3)
Proof. We use the method in [20] to compute the reproduction number R 0 . Let X X 1 , X 2 , where X 1 H s , S s is the healthy population and X 2 = H i , Si , the infected population. H f Si , H s d h H i H i X 1 i Si g H i , S s d s Si X1 , X 2 X1 , X 2 .
According to [19], it follows that
bh b H 0 h e dh t dh dh b b and S t s S 0 s e ds t , ds ds H t
where H 0 H s 0 H i 0 and S 0 S s 0 Si 0 . Thus, b b lim H t h and lim S t s t t dh ds Corollary 3.4 The set
(4) Then,
F
and V
(5)
K H s t , H i t , S s t , Si t 4 ; 0 H s t Hi t
bh b , 0 S s t Si t s dh ds
is invariant and attracting for system (1). Theorem 3.5 For every nonzero, nonnegative initial value, solutions of model system (1) exist for all time t 0. Proof. Local existence of solutions follows from standard arguments since the right-hand side of system (1) is locally Lipschitz. Global existence follows from the a priori bounds.
4. Reproduction Number and Equilibria Analysis In this section, we define the reproduction number and analyze the eventual equilibria. For that, let us denote f f f1 Si , H s , f 2 Si , H s , Si H s g1
g g H i , Ss , g2 H i , Ss . H i S s
System (1) has a disease-free equilibrium (DFE) given by b b EDFE H s0 , 0, S s0 , 0 h , 0, s , 0 . d d s h
Proposition 4.1 The reproduction number for system (1) is R0
Open Access
f1 0, H s0 g1 0, S s0 d s dh
.
f1 0, H s 0
0 X1 , 0 X 2 g1 0, S s d X1 , 0 h X 2 0
0 . d s
Hence,
V 1
1 d h 0
0 0 1 and FV g 0, S 1 s 1 ds d h
f1 0, H s ds . 0
The basic reproduction number is defined as the dominant eigenvalue of matrix FV 1 [20]. Therefore, R0
f1 0, H s0 g1 0, S s0 d s dh
The basic reproduction number evaluate the average number of new infections generated by a single infected individual in a completely susceptible population. Using Theorem 2 in [20], the following result is established. Lemma 4.2 If R 0 1 , then EDFE is locally asymptotically stable. Proof. Using the assumption H2, it follows that f 2 0, H s 0 and g 2 0, S s 0 for all H s and S s . Then, the linearization of system (1) at EDFE gives the following equation
d h d s d s d h
(6)
f1 0, H s0 g1 0, S s0 0.
We can see that all solution of Equation (6) have negative real part. Indeed, the Equation (6) has negative roots d h , d s and other roots are given by
d s d h
f1 0, H s0 g1 0, S s0 .
(7)
Now, we suppose that there exists at least one root 1 of the Equation (7) which have positive real part. According to Proposition 4.1, it follows that AM
S. OUARO, A. TRAORÉ
f1 0, H s0 g1 0, S s0 R 02 d s d h .
Adding the first two and the last two equations of system (9) respectively, we get
Since R 0 < 1 , we obtain
f1 0, H s0 g1 0, S s0 d s d h 1 d s 1 d h .
This show that the Equation (7) cannot have roots with positive real part. Hence, EDFE is locally asymptotically stable according to Theorem 2 in [20] Now, let us analyze the global behavior of the DFE. The global stability of the DFE will be studied using the basic reproduction number R 0 . We first make the following additional assumption. H3 For all H s , H i , S s , Si 4 , f Si , H s f1 0, H s0 Si and g H i , S s g1 0, S s0 H i . Theorem 4.3 The disease-free equilibrium is globally asymptotically stable in K whenever R 0 < 1 . Proof. The proof is based on using a comparison theorem [21]. Note that the equations of the infected components in system (1) can be expressed as
dH i dt H F V i , dSi Si dt
(8)
where F and V , are defined in the proof of Proposition 1. Using the fact that all the eigenvalues of the matrix F V have negative real parts, it follows that the linearized differential inequality system (8) is stable whenever R 0 < 1 . Indeed, we can see that all eigenvalues of the matrix F V satisfy Equation (7). By the same reasoning as in the proof of Lemma 4.2 we deduce that, the eigenvalues of the matrix F V have negative real parts since R 0 < 1 . Thus, H i t , Si t 0, 0 as t for the system (8). Consequently, by a standard comparison theorem [21], H i t , Si t 0, 0 as t and substituting H i Si 0 into system (1) gives H s H s0 and S s S s0 as t . Thus, H s t , H i t , Ss t , Si t H s0 , 0, Ss0 , 0 as t for R 0 < 1 . Therefore, EDFE is globally asymptotically stable if R 0 < 1 . Lemma 4.4 If R 0 1 then, EDFE is unstable and there exists a unique endemic equilibrium E H s , H i , S s , Si . Proof. Let E H s , H i , S s , Si . It follows that E is a positive equilibrium if and only if
bh d h H s f Si , H s H i 0, f Si , H s d h H i H i 0, bs d s S s g H i , S s 0, g H i , S s d s Si 0. Open Access
1685
bh d h H s d h H i 0 and bs d s S s d s Si 0. Then, Hs
bh d h H i b d s Si and S s s . dh ds
Let b dh H i h Si , H i f Si , h dh b d s Si dh H i , g H i , s d s Si . ds We can see that the function h is continuous and so system (9) hold whenever h Si , H i 0 . Therefore any b b zero of h in the set 0, s 0, h where H i and ds dh Si are positives corresponds to an endemic equilibrium. b b According to H2, h 0, 0 0 and h s , h 0 . d s dh Then, in order to have a solution of equation h 0 in bs bh 0, 0, , h must increase at 0. This leads to ds dh dh 0, 0 0, dX 2
where X 2 is defined in the proof of Proposition 4.1. We can see that
dh f1 0, H s0 d h , g1 0, S s0 d s . dX 2
It follows from inequality (10) that
f1 0, H s0 d h 0
and g1 0, S s0 d s 0
d 0 g 0, S d 0.
or f1 0, H and
1
Which implies that (9)
(10)
0 s
h
0 s
s
f1 0, H s0 g1 0, S s0 ds dh
1 ,
that is
R0 1 . To study the global behaviour of this endemic equilibrium, we second make this additional assumption. H4 For all H s , H i , S s , Si 4 , AM
S. OUARO, A. TRAORÉ
1686
f Si , H s Hs S 1 i f Si , H s Si Hs g Hi , Ss Ss H and 1 i. g Hi , Ss Hi Ss
and
(11)
Theorem 4.5 When 0 and R 0 1 , the endemic equilibrium E is globally asymptotically stable in K . Proof. If we consider the system (1) when R 0 1 , there exists a unique endemic equilibrium E . We now establish the global asymptotic stability of this endemic equilibrium. Evaluating both sides of (1) at E with 0 gives bh d h H s f Si , H s , d h H i f Si , H s , bs d s S s g H i , S s , d s Si g H i , S s .
(12)
h x x 1 ln x
(13)
H t Vhs t g H i , S s H s h s , Hs H t Vhi t g H i , S s H i h i , Hi (14) Ss t Vss t f Si , H s S s h , Ss S t Vsi t f Si , H s Si h i . Si We can see that h : * * has the strict global minimun h 1 0 . Thus, Vhs 0, Vhi 0, Vss 0, Vsi 0 with equality if and only if H s H s , H i H i , S s S s and Si Si . We will study the behaviour of the Lyapunov function V t Vhs Vhi Vss Vsi .
(15)
We can see that V t 0 with equality if and only if H s t Hs
Let
1,
Hi t Hi
1,
Ss t Ss
1 and
Si t Si
1.
The derivatives of Vhs , Vhi , Vss and Vsi will be calculated separately and then combined to get the dV desired quantity . dt
H dH H dVhs g H i , S s 1 s s g H i , S s 1 s bh d h H s f Si , H s . dt H s dt Hs
Using the first Equation of (12) to replace bh gives H dVhs g H i , S s 1 s d h H s H s f S i , H s f Si , H s dt Hs 2 H s H s g H , S f S , H 1 f Si , H s 1 H s d h g H i , S s i s i s f S , H H Hs i s s 2 H s H s g H , S f S , H 1 H s f Si , H s f Si , H s H s . d h g H i , S s i s i s H f S , H f S , H H Hs i s i s s s
By adding and substracting the quantity 1 ln
f Si , H s H s , we obtain f Si , H s H s
H s H s g H , S f S , H H s 1 ln H s dVhs d h g H i , S s i s i s H dt Hs Hs s f Si , H s f Si , H s f S i , H s H s f Si , H s H s 1 ln 1 ln f Si , H s f Si , H s f Si , H s H s f Si , H s H s 2 H s H s g H , S f S , H h H s h f Si , H s h f Si , H s H s . d h g H i , S s i s i s H f S , H f S , H H Hs i s i s s s 2
Open Access
(16)
AM
S. OUARO, A. TRAORÉ
Next, we calculate
1687
dVhi . dt
H dH H dVhi g H i , S s 1 i i g ( H i , S s ) 1 i f Si , H s d h H i H dt dt i Hi f Si , H s H H g H i , S s 1 i f Si , H s dh Hi i . Hi f Si , H s H i Using the second Equation of (12) to replace d h H i gives H dVhi f Si , H s g H i , S s 1 i dt Hi
f Si , H s H i f Si , H s H i
By adding and substracting the quantity 1 ln
H f Si , H s H i f Si , H s f Si , H s g H i , S s 1 i . H i f Si , H s H i f Si , H s
f Si , H s H i , we obtain f Si , H s H i
f Si , H s H i f Si , H s H i 1 ln f Si , H s H i f Si , H s H i H f S , H H f S , H f S , H f Si , H s i s i s i s i . f Si , H s g H i , S s h i h h 1 ln f Si , H s H H f , S H f S , H f S , H i s i i i s i s H dVhi H f Si , H s g H i , S s i 1 ln i dt H H i i
Now, we calculate
(17)
dVss . dt S dS S dVss f Si , H s 1 s s f Si , H s 1 s bs d s S s g H i , S s . dt S s dt Ss
Using the third Equation of (12) to replace bs gives S dVss f Si , H s 1 s d s S s S s g H i , S s g H i , S s dt Ss
d s f Si , H s d s f Si , H s
S
s
Ss
2
g H ,S S i s 1 s f Si , H s g H i , S s 1 g H i , Ss Ss
2
S g H i , Ss g H i , Ss Ss f Si , H s g H i , S s 1 s Ss g H i , Ss g H i , Ss Ss
Ss
S
s
Ss Ss
By adding and substracting the quantity 1 ln
g H i , Ss Ss g H i , Ss Ss
.
, we obtain
Ss Ss f S , H g H , S Ss 1 ln Ss dVss d s f Si , H s i s i s S dt Ss Ss s 2
g H ,S g H i , Ss g H i , Ss Ss g H i , Ss Ss i s 1 ln 1 ln g H i , Ss g H i , S s g H i , S s S s g H i , Ss Ss d s f Si , H s
Open Access
S
s
Ss Ss
2
S f Si , H s g H i , S s h s Ss
g H ,S g H ,S S i s i s s h h S g H , S g H , S i s i s s
(18) . AM
S. OUARO, A. TRAORÉ
1688
After that, we evaluate
dVsi . dt
S dVsi g H i , S s f Si , H s 1 i dt Si
S dS dVsi f Si , H s 1 i i dt Si dt S f Si , H s 1 i g H i , S s d s Si Si S f Si , H s 1 i Si
f Si , H s g H i , S s
g H i , S s Si g H i , S s Si
S g H ,S S g H ,S i s i s i . 1 i Si g H i , S s Si g H i , S s
g H i , Ss S d s Si i . g H i , Ss Si g H i , Ss
By adding and substracting the quantity g H i , S s Si 1 ln , we obtain g H i , S s Si
Using the last Equation of (12) to replace d s Si gives
S dVsi S f Si , H s g H i , S s i 1 ln i dt Si Si g H ,S S g H i , S s Si g H i , S s g H i , Ss i s i 1 ln 1 ln S g H i , S s Si g H i , S s g H i , S s g H i , Ss i S g H ,S S i s i f Si , H s g H i , S s h i h g H i , S s Si Si
(19)
g H , S i s h . g H i , Ss
Combining equations (16)-(19), we obtain 2 H s H s d f S , H Ss Ss g H , S f S , H h Si dV d h g H i , S s i s i s S s i s Hs Ss dt i 2
g H ,S S H i s s h i h H S i g H i , Ss s
f S , H H S i s i h h s f Si , H s H i S s
Since the function h is monotone on each side of 1
g H ,S S Hs i s i h h H S s g H i , Ss i
5. Stochastic Differential Equation Model To derive a stochastic model, we apply a similar procedure to that described in Allen [23]. Here, we Open Access
(20)
H h i . Hi
neglect the possibility of multiple events of order t . The possible changes in the populations over a short time t , concern individual births, deaths and transformation. These changes are produced in Table 1, together with their corresponding probability. Let’s denote this change T by H s , H i , S s , Si . 2 Neglecting terms of the order t , the mean of system (0.1) is given by 2
dV 0, (21) dt for all H s , H i , S s , Si K with equality only for H s H s , H i H i , S s S s and Si Si . Hence, the endemic equilibrium E is the only positively invariant set of the system (0.1) contained in H s , H i , S s , Si 4 ; H s H s , H i H i , S s S s , Si Si . Then, it follows that E is globally asymptotically stable on K [22].
.
and is minimized at 1, H4 implies
f S , H H g H ,S S S i s i s s s h i and h h f Si , H s H s S S i g H i , Ss s
Since h 0 , then
f S , H H i s s h H f S H , s i s
bh d h H s f Si , H s H i f Si , H s d h H i H i E Pii t t. bs d s S s g H i , S s i 1 g H i , S s d s Si (22) 10
AM
S. OUARO, A. TRAORÉ
6. Equivalent Stochastic Differential Equation Model
Table 1. Possible changes in the population. Change
1 1,0,0,0
Probability
2 1,0,0,0 3 1,1,0,0
6 0,0,1,0
P2 d h H s t
T
P3 f Si , H s t
T
4 0, 1,0,0 5 1, 1,0,0
T
P4 d h H i t
T
P5 H i t
7 0,0, 1,0
10 0,0,0,0
P7 d s S s t
T
P8 g H i , S s t
T
9 0,0, 0, 1
H s u1 u2 u3 u4 , H i u3 u4 u5 , S s u6 u7 u8 , Si u8 u9 ,
P6 bs t
T
8 0,0, 1,1
In this section, we develop a stochastic model to examine the changes occured on each vector individually. Unless we say otherwise, we adopt the vectors defined in the previous section but here the Poisson processes P are used to establish the different probabilities as mean times before occurence. Then, we have
P1 bh t
T
P9 d s Si t
T
u1 P bh t , u2 P d h H s t ,
P10 1 i 1Pi
Further, the covariance matrix of system (1) is given by
E T
B11 B PiiiT 12 0 i 1 0 10
B12 B22
0
0 0
B33 B34
0
0 0 t Bt , B34 B44 (23)
where B11 bh d h H s f Si , H s H i , B22 f Si , H s d h H i H i , B12 f Si , H s H i , B33 bs d s S s g Si , H s , B34 g Si , H s and B44 g Si , H s d s Si .
Note that it has been proved in [23] that the changes
are normally distributed. Then,
Y t t Y t Y t t Y t t Bt ,
where i N 0,1 for i 1, 2,, 4 . Furthermore, as t 0 , Y t converges strongly to the solution of the stochastic system dY t dW t , Y t B Y t dt dt
(24)
where Y t H s , H i , S s , Si and W t is the fourdimensional Wiener process [23]. The computational work of Equation (24) involves the calculation of the quantity B Y t at each time step [24,25]. This procedure is often cumbersome. In the next section, we derive an equivalent stochastic model which seem to be easier to implement. T
Open Access
(25)
where
9
T
1689
u3 P f Si , H s t , u4 P H i t , u5 P d h H i t , u6 P bs t , u7 P d s S s t , u8 P g H i , S s t , u9 P d s Si t .
We now normalize the Poisson process to get H s bh t bh t 1 d h H s t d h H s t 2 f Si , H s t f Si , H s t 3 H i t H i t 4 , H f S , H t f S , H t i s i s 3 i d H t d H t H t H t , h i 5 i i 4 h i S s bs t bs t 6 d s S s t d s S s t 7 g H i , S s t g H i , S s t 8 , Si g H i , S s t g H i , S s t 8 d s Si t d s Si t 9 , (26) where i N 0,1 for i 1, 2, ,9 . Then, as t 0 , the system (26) converges to the following Itô stochastic differential equation [24] dH s bh d h H s f Si , H s H i dt bh dW1 d h H s dW2 f Si , H s dW3 H i dW4 , dH f S , H d H H dt f S , H dW i s h i i i s 3 i H i dW4 d h H i dW5 , dS s bs d s S s g H i , S s dt bs dW6 d s S s dW7 g H i , S s dW8 , dSi g H i , S s d s Si dt g H i , S s dW8 d s Si dW9 . (27)
System (27) can be rewritten as follows. AM
S. OUARO, A. TRAORÉ
1690
dY t dW t Y t G , dt dt
G3
(28)
f Si , H s , G4 H i ,
where Y t and are the same as in system (24), W is the nine-dimensional Wiener process and G is defined by
G5 d h H i , G6 bs ,
0 0 0 0 0 G1 G2 G3 G4 0 0 G3 G4 G5 0 0 0 0 G , 0 0 0 0 0 G6 G7 G8 0 0 0 0 0 0 0 G8 G9 0 (29)
and G9 d s Si .
G7 d s S s , G8 g H i , S s
7. Computational Results In this section, computational results are given for the stochastic system (28). We use the Euler-Maruyama method with 1000 sample to solve the SDE model (28). Let h be a specified time step. The Euler-Maruyama numerical method for system (28) is given by:
where G1 bh , G2 d h H s ,
H sk 1 H sk h bh d h H s f Si , H s H i h H ik 1 H ik h f Si , H s d h H i H i h S
k 1 s
f Si , H s 3, k H i 4, k ,
bh1, k d h H s 2, k
f Si , H s 3, k H i 4, k d h H i 5, k ,
S h bs d s S s g H i , S s h k s
Sik 1 Sik h g H i , S s d s Si h
for k 0,1, 2, until the maximum time is reached. Here, i ,k for i 1, 2, ,9 and k 0,1, 2, are normally distributed numbers with zero mean and unit variance.
8. Computational Simulations In this section, computational simulations are given for the models described in the previous sections. Here, two different cases of computational simulations were studied: in the first case R 0 1 while in the second case R 0 1 . In the computations, the functions f and g are chosen as follows f Si , H s Si H s and g H i , S s H i S s (mass action), where is the rate infection of hosts by cercaria released by an infected snail and is the infection rate of snails by miracidia of the parasite eggs of an infected human. The parameters values used are given as: d h 0.014 (human beings average life span 70 years 1 d h ), bh 14 (the carrying capacities of human population is 1000), d s 0.2 (snails average life span is 2 - 8 years), bs 400 (the carrying capacities of snails population is 2000), 0.0004 [26], 0.000027 [26]. The treatment rate in the first case is 0.1 which give R 0 0.97 . In the second case the treatment rate is chosen as 0.05 to have R 0 1.29 . The initial values of the population sizes are taken as H s 0 600 , H i 0 400 , S s 0 1300 , Si 0 700 . The time step h was chosen as h 0.3 year and the final time was taken as 300 years. These figures are produced by Matlab. Open Access
bs 6, k d s S s 7, k g H i , S s 8, k ,
g H i , S s 8, k d s Si 9, k ,
Figure 2 illustrates the deterministic model (1) and the equivalent stochastic model (28) when R 0 1 . We can see that, in the Figure 2 the trajectory of de- terministic and stochastic graphs are approximately the same behaviour. Indeed, the parasite extinction is effective if R0 1 . Figure 3 depict the deterministic model (1) and the equivalent stochastic model (28) when R 0 1 . In Figure 3 we can observe that, the features of deterministic graph are similar to those of the stochastic graph, namely the parasite extinction is unlikely for this two models when R 0 1 . The Table 2 and Table 3 below contains the results that supports the comparison we made about the Table 2. Mean and standard deviation for the ODE epidemic model (1) and the SDE epidemic model (28) at t = 300 years where R 0 = 0.97 . Models
ODE (1)
SDE (28)
Variables X i
E X i
σ Xi
Hs
992.25
0
Hi
7.74
0
Ss
1968.63
0
Si
31.36
0
Hs
993.21
32.49
Hi
6.44
9.44
Ss
1985.71
719.28
Si
26.12
38.78
AM
S. OUARO, A. TRAORÉ
1691
Figure 2. Deterministic and Equivalent stochastic models for R 0 = 0.97 .
Figure 3. Deterministic and Equivalent stochastic models for R 0 = 1.29 .
deterministic model (1) and the equivalent stochastic model (28).
9. Conclusion An ordinary differential equation model and a corresponding equivalent stochastic differential equation Open Access
model for schistosomiasis were studied. Four subpopulations were modelled: susceptible human, infected human, susceptible snails and infected snails. A threshold number was exhibited. Computational simulations were presented. The behavior of the deterministic and equivalent stochastic models are approximately the same. For the deterministic and equivalent stochastic models AM
S. OUARO, A. TRAORÉ
1692
Table 3. Mean and standard deviation for the ODE epidemic model (1) and the SDE epidemic model (28) at t = 300 years where R 0 = 1.29 . Models
ODE (0.1)
SDE (0.28)
Variables X i
E X i
σ Xi
Hs
813.52
0
Hi
186.47
0
Ss
1456.71
0
Si
543.28
0
Hs
818.21
42.05
Hi
181.56
36.85
Ss
1474.13
451.02
Si
532.71
148.20
128. [11] G. Macdonald, “The Dynamics of Helminth Infections, with Special Reference to Schistosomiasis,” Transactions of the Royal Society of Tropical Medicine and Hygiene, Vol. 59, No. 5, 1965, pp. 489-506. http://dx.doi.org/10.1016/0035-9203(65)90152-5 [12] I. Nasell, “A Hybrid Model of Schistosomiasis with Snail Latency,” Theoretical Population Biology, Vol. 10, No. 1, 1976, pp. 47-69. http://dx.doi.org/10.1016/0040-5809(76)90005-8 [13] M. E. J. Woolhouse, “On the Application of Mathematical Models of Schistosome Transmission Dynamics. I. Natural Transmission,” Acta Trop. 49, 1241-1270, (1991). http://dx.doi.org/10.1016/0001-706X(91)90077-W [14] Woolhouse, M.E.J., On the application of mathematical models of schistosome transmission dynamics. II. Control. Acta Tropica, Vol. 50, 1992, pp. 189-204. http://dx.doi.org/10.1016/0001-706X(92)90076-A
the disease dies out for R 0 1 and the disease limits to an endemic equilibrium when R 0 1 .
[15] Y. Yang and D. Xiao, “A Mathematical Model with Delays for Schistosomiasis,” Chinese Annals of Mathematics, Vol. 31B, No. 4, 2010, pp. 433-446. http://dx.doi.org/10.1007/s11401-010-0596-1
REFERENCES
[16] L. J. S. Allen, “An Introduction to Stochastic Processes with Applications to Biology,” Prentice-Hall, Englewood Cliffs, 2003.
[1]
World Health Organization, 2012. http://www.who.int/schistosomiasis/en/
[2]
P. Jordan and G. Webbe, “Human Schistosomiasis,” Chales Thomas, Springfield, 1969.
[3]
P. L. Rosenfield, “The Management of Schistosomiasis,” Resources for the Future, Washington DC, 1979.
[4]
E. J. Allen and H. D. Victory, “Modelling and Simulation of a Schistosomiasis Infection with Biological Control,” Acta Tropica, Vol. 87, No. 2, 2003, pp. 251-267. http://dx.doi.org/10.1016/S0001-706X(03)00065-2
[5]
R. M. Anderson and R. M. May, “Prevalence of Schistosome Infections within Molluscan Populations: Observed Patterns and Theoretical Predictions,” Parasitology, Vol. 79, No. 1, 1979, pp. 63-94. http://dx.doi.org/10.1017/S0031182000051982
[6]
R. M. Anderson and R. M. May, “Helminth Infections of Humans: Mathematical Models, Population Dynamics, and Control,” Advances in Parasitology, Vol. 24, 1985, pp. 1-101. http://dx.doi.org/10.1016/S0065-308X(08)60561-8
[7]
J. E. Cohen, “Mathematical Models of Schistosomiasis,” Annual Review of Ecology and Systematics, Vol. 8, 1977, pp. 209-233. http://dx.doi.org/10.1146/annurev.es.08.110177.001233
[8]
Z. Feng, C. Li and F. A. Milner, “Schistosomiasis Models with Density Dependence and Age of Infection in Snail Dynamics,” Mathematical Biosciences, Vol. 177-178, 2002, pp. 271-286.
[9]
Z. Feng, C. Li and F. A. Milner, “Schistosomiasis Models with Two Migrating Human Groups,” Mathematical and Computer Modelling, Vol. 41, 2005, pp. 1213-1230.
[10] Z. Feng and F. A. Milner, “A New Mathematical Model of Schistosomiasis,” In: Innovation Applied Mathematics, Vanderbilt University Press, Nashville, 1998, pp. 117Open Access
[17] T. C. Gard, “Introduction to Stochastic Differential Equations,” Marcel Dekker, New York, 1987. [18] J. M. Bony, “Principe du Maximum, Inégalité de Harnack et Unicité du Problème de Cauchy Pour les Opérateurs Elliptiques Dégénérés,” Annales de L’institut Fourier (Grenoble), Vol. 19, 1969, pp. 277-304. [19] G. Birkhoff and G. C. Rota, “Ordinary Differential Equations,” Ginn, Boston, 1982. [20] P. Van den Driesche and J. Watmough, “Reproduction Numbers and Subthreshold Endemic Equilibria for the Compartmental Models of Disease Transmission,” Mathematical Biosciences, Vol. 180, No. 1-2, 2002, pp. 2948. http://dx.doi.org/10.1016/S0025-5564(02)00108-6 [21] V. Lakshmikantham, S. Leela, A. A. Martynyuk, “Stability Analysis of Nonlinear Systems,” Marcel Dekker, New York, 1989. [22] J. P. LaSalle, “The Stability of Dynamical Systems,” SIAM, Philadelphia, 1976. http://dx.doi.org/10.1137/1.9781611970432 [23] E. J. Allen, “Stochastic Differential Equations and Persistence Time of Two Interacting Populations,” Dynamics of Continuous, Discrete and Impulsive Systems, Vol. 5, 1999, pp. 271-281. [24] E. J. Allen, L. J. S. Allen, A. Arciniega and P. E. Greenwood, “Construction of an Equivalent Stochastic Differential Equation Models,” Stochastic Analysis and Applications, Vol. 26, No. 2, 2008, pp. 274-297. http://dx.doi.org/10.1080/07362990701857129 [25] E. J. Allen, J. Baglama and S. K. Boyd, “Numerical Approximation of the Product of the Square Root of a Matrix with a Vector,” Linear Algebra and Its Applications, Vol. 310, No. 1-3, 2000, pp. 167-181. http://dx.doi.org/10.1016/S0024-3795(00)00068-9 AM
S. OUARO, A. TRAORÉ [26] Z. Feng, A. Eppert, F. A. Milner and D. J. Minchella, “Estimation of Paramaters Governing the Transmission Dynamics of Schistosomes,” Applied Mathematics Let-
Open Access
1693
ters, Vol. 17, No. 10, 2004, pp. 1105-1112. http://dx.doi.org/10.1016/j.aml.2004.02.002
AM
