DYNAMICAL BEHAVIORS OF AN ECHINOCOCCOSIS EPIDEMIC MODEL WITH DISTRIBUTED DELAYS

. In this paper, a novel spreading dynamical model for Echinococcosis with distributed time delays is proposed. For the model, we ﬁrstly give the basic reproduction number R 0 and the existence of a unique endemic equilibrium when R 0 > 1. Furthermore, we analyze the dynamical behaviors of the model. The results show that the dynamical properties of the model is completely determined by R 0 . That is, if R 0 < 1, the disease-free equilibrium is globally asymptotically stable, and if R 0 > 1, the model is permanent and the endemic equilibrium is globally asymptotically stable. According to human Echinococcosis cases from January 2004 to December 2011 in Xinjiang, China, we estimate the parameters of the model and study the transmission trend of the disease in Xinjiang, China. The model provides an approximate estimate of the basic reproduction number R 0 = 1 . 23 in Xinjiang, China. From theoret- ic results, we further ﬁnd that Echinococcosis is endemic in Xinjiang, China. Finally, we perform some sensitivity analysis of several model parameters and give some useful measures on controlling the transmission of Echinococcosis .

1. Introduction. Echinococcosis, which is often referred to as a hydatid disease, is a parasitic disease that affects both humans and other mammals, such as sheep, dogs, rodents and horses [3]. There are three different forms of Echinococcosis found in humans, each of which is caused by the larval stages of different species of the tapeworm of genus Echinococcus. The first of the three and also the most common form found in humans is cystic Echinococcosis, which is caused by Echinococcus granulosus. The second is alveolar Echinococcosis, which is caused by Echinococcus multilocularis and the third is polycystic Echinococcosis, which is caused by Echinococcus vogeli and very rarely, Echinococcus oligarthus. Alveolar and polycystic Echinococcosis are rarely diagnosed in humans and are not as widespread as cystic Echinococcosis. Thus, we focus on cystic Echinococcosis in this paper.
Echinococcus granulosus is an extremely small tapeworm, only 4-6 mm in length in its adult stage. This stage occurs in carnivore "definitive hosts," generally dogs, wolves, coyotes, or wild dogs [17]. When mature, the adult tapeworm produces eggs that are passed in feces of the infected carnivore. The eggs contaminate foliage or vegetation which may be eaten by grazing animals, the "intermediate hosts," generally ruminants such as sheep, goats, or cattle, although an extremely wide range of potential intermediate hosts exists. After the intermediate host ingests tapeworm eggs, the eggs hatch into tiny embryos in the small intestine of that host. These small mobile forms penetrate the intestinal mucosa, enter the circulatory tract, and are then transported via the bloodstream to various filtering organs in that host. The embryo lodges in those organs and there develops into a larval cystic form termed the "hydatid cyst". Eventually a thick wall develops around the hydatid cyst, and numerous tiny tapeworm heads are produced via asexual division within "brood capsules" from the germinal lining of each hydatid cyst. If parts of the viscera containing hydatid cysts are eaten by suitable carnivore hosts, then the tiny protoscolices are liberated from the brood capsules, evaginate their head region with its suckers and hooks, and move actively to eventually attach to the intestinal mucosa. There they grow within approximately 40 days to the adult stage to the tapeworm, which is then capable of producing eggs. Humans act as intermediate hosts for Echinococcus granulosus, and are infected when they ingest tapeworm eggs from the definitive host. The eggs may be eaten in foods such as vegetables, fruits or herbs, or drunk in contaminated water. The life cycle of Echinococcus granulosus is shown in Figure 1.  [20].
Echinococcosis is widely distributed and highly endemic in Xinjiang, China. Between 2004 and 2011, a total of 6086 human Echinococcosis cases were reported in Xinjiang, China (see Figure 2 [22]).
The average infection rates in sheep, cattle and dogs were 54.20%, 15.21% and 24.52%, respectively [9].  samples of cattle and sheep in fix point slaughter houses (see Table 2 [11]). Meanwhile, Echinococcus granulosus infection for dogs was also examined in various parts of Xinjiang, China (see Table 2 [23]). Mathematical modeling is an important method of studying the spread of infectious disease qualitatively and quantitatively. Models can be used to understand how an infectious disease spreads in the real world, and how various complexities affect the dynamics. There have been some modeling studies on different aspects of Echinococcosis, and most of them are statistical models. A base-line survey was carried out on the transmission dynamics of Echinococcus granulosus, Taenia hydatigena and Taenia ovis in sheep in Uruguay by Cabrera et al. in [4]. The Echinococcus granulosus was only relative stable and the basic reproduction ratio was estimated about 1.2 in an endemic steady state. The effects of the different transmission levels of the parasites on potential control strategies were discussed. Torgerson et al. [16] examined the abundance and prevalence of infection of Echinococcus granulosus in cattle and sheep in Kazakhstan. Maximum likelihood techniques were used to define the parameters and their confidence limits in the model and the negative binomial distribution was used to define the error variance in the observed data. Lahmar et al. [8] examined the abundance and prevalence of infection of Echinococcus granulosus in camels in Tunisia.
Recently, dynamic models of Echinococcosis have been discussed in [19,21]. The work [19] proposed a deterministic model to study the transmission dynamics of Echinococcosis in Xinjiang, China. The results demonstrated that the dynamics of the model were completely determined by the basic reproductive number R 0 . The authors showed that Echinococcosis was endemic in Xinjiang, China with the current control measures. In [21] two mathematical models, the baseline model and the intervention model, were proposed to study the transmission dynamics of Echinococcosis. Applying these mathematical models to Qinghai Province, China, the authors showed that the infection of Echinococcosis was in an endemic state. In [19,21], the authors considered the density of Echinococcus eggs as a compartment in their models. However, it is hard to obtain the parameters of Echinococcus eggs in actual epidemiological surveys. To solve this problem, we take Echinococcus eggs in the environment as vectors and propose a novel spreading model for Echinococcosis with distributed time delays in this paper. We will evaluate the basic reproduction number R 0 and analyze the dynamical behaviors of the model. According to human Echinococcosis cases from January 2004 to December 2011 in Xinjiang, China, we will estimate the parameters of the model and study the transmission trend of the disease in Xinjiang, China. We will also carry out some sensitivity analysis of the basic reproduction number R 0 in terms of various model parameters. Finally, we shall suggest some strategies for control Echinococcosis infection in Xinjiang, China.
The article is organized as follows. The model is presented in Section 2. The basic properties on the positivity and boundedness of solutions, computation of the basic reproduction number and the existence of endemic equilibrium for the model are discussed in Section 3. In Section 4, we establish the global stability of the diseasefree equilibrium for the model. In Section 5, we will apply the theory of permanence for infinite-dimensional systems to obtain the permanence of the model. The global stability theorem of the endemic equilibrium is stated and proved in Section 6. In Section 7, we estimate the parameters of the model and study the transmission trend of the disease in Xinjiang, China. A brief discussion is given in Section 8. The dogs population is subdivided into two classes: the susceptible population and the infected population denoted by S D (t) and I D (t), respectively. For livestock population, we separate the total livestock population into two classes: susceptible and infectious denoted by S L (t) and I L (t), respectively. Humans are infected when they ingest tapeworm eggs from the definitive host. We separate the total population into three classes: susceptible, exposed and infectious denoted by S H (t), E H (t) and I H (t), respectively.
To model the Echinococcosis infection delay from infected livestock to susceptible dogs, we let τ be the random variable that describes the time between infection offal thrown in the environment and eaten by dogs with a probability distribution f 1 (τ ).
To model the Echinococcosis infection delay from infected dogs to susceptible livestock, we let τ be the random variable that is the time between eggs produced by infected dogs and ingestion by susceptible livestock with a probability distribution f 2 (τ ). h 2 describes average life expectancy for Echinococcus eggs. f 2 (τ ) is assumed to be non-negative and continuous on [0, h 2 ] satisfying h2 0 f 2 (τ )dτ = 1. When infected dogs produced Echinococcus eggs at time t − h 2 , the number of infected livestock at time t is given by Humans act as an accidental intermediate host for Echinococcosis. Thus, similarly to the infection of livestock, the number of infected human at time t is given by where f 3 (τ ) is a probability distribution. Because the range of movement for humans and livestock is different, we choose distribution function in different forms. Hence, distribution f 3 (τ ) is different from f 2 (τ ). Similarly, f 3 (τ ) is assumed to be nonnegative and continuous on [0, h 2 ] satisfying h2 0 f 3 (τ )dτ = 1. Based on the above discussions, a dynamical model of Echinococcosis with distributed time delays can be written as: (1) All parameters in model (1) are assumed to be positive constants. For the dog population, A 1 describes the annual recruitment rate; d 1 is the natural death rate; σ denotes the recovery rate of transition from infected to non-infected dogs, including natural recovery rate and recovery due to anthelmintic treatment; For the livestock population, A 2 is the annual recruitment rate; d 2 is the death rate; For human population, A 3 is the annual recruitment rate; d 3 is the natural death rate; µ is the disease-related death rate; 1/ω denotes the incubation period of infected individuals; γ denotes the recovery rate. Motivated by biological background of model (1), the initial conditions for model (1) are given as follows In this paper, for any integer n > 0 we denote by 3. Basic properties. For any ε > 0, we define the region Γ ε as follows.
Firstly, on the positivity and ultimate boundedness of solutions for model (1), we have the following result.
(ii) All solutions of model (1) with initial conditions (2) ultimately turn into the region Γ ε as t → ∞.
The proof of Lemma 3.1 is simple, we hence omit it here.
Remark 3.2. According to Lemma 3.1, all feasible solutions of model (1) enter or remain in the region Γ ε as t is large enough. Therefore, the dynamics of model (1) can be considered only in Γ ε .
Simple algebraic calculation shows that model (1) always has a disease-free equi- Following the idea in [24,18], we can introduce the basic reproduction number for model (1) as follows .

Remark 3.3.
To better understand the basic reproduction number, we rewrite it as The biological meaning of R 0 can be interpreted as follows, near the disease free equilibrium, the total number of livestock reaches at a stable state A 2 /d 2 and every livestock will be susceptible. Then susceptible livestock is infected by infectious dogs over its expected infectious period 1/(d 1 + σ) with probability β 2 . Furthermore, the total number of dogs reaches stable state A 1 /d 1 and every dog will be susceptible, which is infected by infectious livestock over its expected infectious period 1/d 2 with probability β 1 . The square root arises from the two 'generations' required for an infected dogs or livestock to 'reproduce' itself.
It can be proved that if R 0 > 1 then model (1) has a unique endemic equilibrium Therefore, we finally have the following theorem. 4. Global stability of E 0 . Regarding the stability of the disease-free equilibrium E 0 , we have the following result. (1) is globally asymptotically stable in Γ ε .
Proof. Notice that the last three equations are independent of the first four equations, we start by considering the first four equations, Since R 0 < 1, it follows that Hence, there is a small enough constant ε > 0 such that Consider the following Lyapunov functional calculating the derivative of V (t) along the solution of system (3), we havė Solving second equation, we have Hence, I D (t) ≡ 0. Otherwise, I D (t) will be unbounded on R. Since (3) is globally asymptotically stable. Next, we consider the last three equations Since I D (t) → 0 as t → ∞, we obtain that E H (t) → 0 as t → ∞ by the second equation of system (4). Furthermore, from the third equation of system (4) we also have I H (t) → 0 as t → ∞. Lastly, from the first equation of system (4) we finally have that S H (t) → A 3 /d 3 as t → ∞. Therefore, (A 3 /d 3 , 0, 0) is globally attractive with respect to system (4). Thus, according to the theory of asymptotic autonomous system [15], we finally obtain that the disease-free equilibrium E 0 is globally asymptotically stable for model (1) when R 0 < 1. This completes the proof of Theorem 4.1.

5.
Permanence. In this section, we will apply the theory of persistence for infinitedimensional dynamical system from [7] to obtain the permanence of model (1). Let X be a complete metric space with metric d. Suppose that T (t) : X → X, t ≥ 0, is a C 0 -semigroup on X. For any x ∈ X, the positive orbit γ + (x) through x is defined as γ + (x) = {T (t)x : t ≥ 0}, and its ω-limit set is ω(x) = {y ∈ X: there is a time sequence t k → ∞ such that lim k→∞ T (t k )x = y}. Let X 0 and X 0 be two nonempty subsets such that X = X 0 ∪ X 0 and X 0 ∩ X 0 = ∅. Further, we assume that C 0 -semigroup T on X satisfies We denote T 0 (t) = T (t)| X 0 and T 0 (t) = T (t)| X0 . We further assume that there is a global attractor A 0 for T 0 in X 0 , and let Furthermore, the definitions of the compactness, point dissipativity and the uniform persistence for C 0 -semigroup T on X, and the definitions of the isolated invariance and acyclic covering for a subset A ⊂ X, can be found in Kuang [7]. Now, we state the following lemma, which can be found in ( [7], Theorem 2.4).
(iii)Ã 0 is isolated and has an acyclic covering In the following, as an application of this lemma, we establish the result on the permanence of model (1).
Denote by ω(φ) the ω-limit set of the solution of model (3) starting at t = 0 with initial value φ ∈ X. Let Restricting system (3) to X 0 , we obtain It is easy to verify that system (6) has a unique equilibrium (A 1 /d 1 , A 2 /d 2 ), which is globally asymptotically stable. Therefore, we have Ω = {E 0 }, and E 0 is a covering of Ω, which is isolated (since E 0 is the unique equilibrium) and is acyclic (since there exists no solution in X 0 which links E 0 to itself). It remains to show that where W s (E 0 ) denotes the stable manifold of E 0 . Suppose that (7) For this > 0, by (8), there exists a t 1 > 0 such that Hence, by system (3), we have for all t ≥ t 1 . Consider the following auxiliary system Let (U (t), V (t)) be the solution of system (9) with initial conditions U (t 1 + θ) = I D (t 1 + θ) and V (t 1 + θ) = I L (t 1 + θ) for all θ ∈ [−τ, 0], by the comparison theorem of functional differential equations( [12], Chapter 5, Theorem 1.1) we have From (8), we further have Consider the following auxiliary function where Then, from (10) we have lim t→∞ W (t) = 0.
Calculating the derivative of W (t) along the solution of system (9), we havė Therefore, W (t) increases as t increases. This shows that W (t) does not tend to zero as t → ∞, which leads a contradiction. Therefore, (7) holds and system (3) satisfies all conditions of Lemma 2. Therefore, system (3) is permanent. Then there exists a constant η > 0 such that lim inf t→∞ I D (t) > η. For any positive number ε < η, there is a T > 0 such that I D (t) > η − ε, for all t ≥ T.
Next, we prove that system (4) is also permanent. From system (4), we have Using a similar argument as above, we can obtain that This implies that model (1) is permanent.
6. Global stability of E * . In this section, we will obtain the stability of endemic equilibrium of model (1). We will use the theory of cooperative systems to prove the global stability of system. Consider the functional differential equations where f : U → R n is a continuously differentiable map on the open subset U of C r . Denote x(t, φ) for the unique solution of (11) satisfying x 0 = φ and let [0, σ φ ] be its maximal interval of existence. Now we have the following lemma, which can be found in ( [26], Theorem 3.2).
Proof. We still start by considering system (3). It is easily proved that S D (t) + I D (t) → A 1 /d 1 and S L (t) + I L (t) → A 2 /d 2 as t → ∞. Therefore, in system (4) we can represent S D (t) and S L (t) by A 1 /d 1 − I D (t) and A 2 /d 2 − S L (t), respectively, then system (3) reduces to the following system with two equations: By Lemma 3.1, the dynamics of system (12) can be focused on the following restricted region We will use Lemma 6.1 to prove the global stability of system (12). We only need to verify all assumptions in Lemma 6.1 for system (12). Consider the following equation For any φ, ψ ∈ C([−h, 0], R 2 + ), we have Hence, for any ψ(0) ∈ (0, which imply that f is a cooperative map. For any φ ∈ C([−h, 0], R 2 + ) with φ(0) = 0, we have It is easy to see from the expression of f that f maps bounded sets in C([−h, 0], R 2 + ) into bounded sets in R 2 + . For any α ∈ (0, 1) and φ ∈ C([−h, 0], R 2 + ) we have Using the same argument, we can show that f 2 (αφ) ≥ αf 2 (φ). Therefore, the functional f is sublinear. For x ∈ R 2 + with x 0, we have and This shows that F is strictly sublinear. Since f (0) = 0 and for any φ ∈ C([−h, 0], R 2 + ), it is easy to see that df (0) satisfies assumption (R) in [26]. By computing DF (x), we have Clearly, DF (x) is irreducible for x ∈ Ω. Since When R 0 > 1, we have a 0 < 0. Hence, s(DF (0)) := max{Re λ : det(λI −DF (0)) = 0} > 0. By Lemma 6.1, one can conclude that the equilibrium (I * D , I * L ) of system (12) is globally asymptotically stable. According to the theory of asymptotic autonomous systems [15], we further obtain that endemic equilibrium (S * D , I * D , S * L , I * L ) is globally attractive for system (3).
Next, when t → ∞, the limiting system of equation (4) is Let we can transform equation (13) into the following system The Jacobian matrix of equation (14) at (0, 0, 0) is The characteristic equation is It is clear that a 0 , a 1 , a 2 > 0, and by simple computation, we have Hence, by Routh-Hurwitz criteria, all roots of φ(λ) = 0 have negative real parts. Furthermore, by the stability theory of linear differential equations, it follows that (0, 0, 0) is globally attractive for equation (14). Therefore, (S * H , E * H , I * H ) is also globally attractive for equation (13). According to theory of asymptotic autonomous systems [15], we finally obtain that endemic equilibrium E * (S * D , I * D , S * L , I * L , S * H , E * H , I * H ) of model (1) is globally asymptotically stable. 7. Numerical simulations. In this section, we shall estimate the parameters of model (1) and study the transmission trend of the disease in Xinjiang, China. From Xinjiang Center for Disease Control and Prevention (Xinjiang CDC) [22], we have obtained the monthly numbers of reported Echinococcosis cases from January 2004 to December 2011 (see Figure 2).
As an illustration, we take the parameters in Table 3 and give a few necessary comments: (a) The number of dogs in Xinjiang, China was estimated to be 2 millions in 2004 [25]. The average life expectancy of dogs was 12.5 years [1]. Thus, d 1 = 1/(12.5 × 12) = 0.0067. We estimated the capacity of dogs in Xinjiang, China to be 2 millions. Thus, A 1 ≈ 2000000 × 0.0067 = 1.34 × 10 4 . The adult Echinococcus  [13]. We select 6 months. So, the recovery rate σ = 1/6. (b) According to China Statistical Yearbook [10], the recruitment rate A 2 for livestock population was estimated to be 8.7×10 6 and the mortality rate of livestock d 2 was estimated to be 0.0275.
(c) The average survival time of larval cysts into the infection offal thrown in the environment is assumed to be 7-10 days [6]. We select h 1 = 1/3. Life expectancy is assumed to be 35 days for parasite eggs [17]. So h 2 = 35/30 = 1.17.
(d) The average life expectancy of people in Xinjiang Uygur Autonomous Region of China was 71.12 years in 2005 [14]. We take it as the current average life expectancy. Thus, d 3 = 1/(71.12 × 12) = 0.0012. According to the data of the Sixth Census of China [2], human annual birth population in Xinjiang, China was 2.4 × 10 5 , that is A 3 = (2.4 × 10 5 )/12 = 2 × 10 4 . The incubation period of Echinococcosis can be months to years or even decades. We select 14 years [13].
(e) The parameters β 1 , β 2 and β 3 are obtained by fitting the model to data. We first apply the Runge-Kutta method to obtain the numerical expression formula of model (1) with parameters β 1 , β 2 and β 3 , then by the least-square estimation the estimation values of β 1 , β 2 and β 3 are given through MATLAB implementation.
The numerical simulation of the model on the number of human Echinococcosis cases is shown in Figure 3. We can predict the general tendency of the epidemic in a long term according to the current situation, which is presented in Figure 4. It shows that the number of human Echinococcosis cases will increase steadily in the next 6 or 7 years, then reach the peak (about 120) in 2017 and begin a slow decline, and finally become stable.
Moreover, with these parameter values, we can roughly estimate the basic reproduction number R 0 = 1.23 for model (1) under the current circumstances in Xinjiang, China. This indicates that human Echinococcosis is endemic in Xinjiang, China: it stabilizes and is approaching its equilibrium. Therefore, if no further effective prevention and control measures are taken, the disease will not vanish. From Figure 5, we can see that when R 0 < 1 the number of Echinococcosis humans I H (t) tends to 0. On the contrary, when R 0 > 1, I H (t) tends to a stable equilibrium. We can see that R 0 is really the threshold for the establishment of the disease and the number of infections increases with the increase of R 0 .   Next, to find better control strategies for Echinococcosis infection, we perform some sensitivity analysis to determine the influence of parameters A 1 , d 1 , σ and β 1 on R 0 . From Figure 6 we can see that R 0 decreases if A 1 decreases, or d 1 increases, or σ increases, or β 1 decreases. From Figure 6(a), it is obvious that when A 1 is less than 6500, R 0 can be less than 1. However, the recruitment rate of dogs can achieve 13400 in Xinjiang, China. This indicates that human Echinococcosis in Xinjiang, China cannot be eradicated if recruitment rate of dogs cannot be controlled under 6500. From Figure 6(d), we can see that when β 1 is less than 3.5 × 10 −10 , R 0 can be less than 1. In order to decrease the parameter β 1 , we have the following measures. Dogs should be barred from slaughter houses and should not be fed uncooked offal. Infection carcasses and offal should be burned or buried. The parameter σ can be increased through increasing the frequency of dog anthelmintic.
8. Discussion. In this paper, in order to understand how Echinococcosis spreads in the real world and forecast the disease trends, we propose a novel spreading model for Echinococcosis with distributed time delays. We calculated the basic   Table 3 do not change, R 0 = 0.9932 and 1.2321, respectively.  Table 3 do not change.
reproduction number R 0 . The mathematical results show that the dynamics of the model is completely determined by R 0 . If R 0 < 1, the disease-free equilibrium is globally asymptotically stable, and if R 0 > 1, the model is permanence and endemic equilibrium is globally asymptotically stable.
With suitable parameter values, we estimated the basic reproduction number to be R 0 = 1.23 in Xinjiang, China for model (1). Moreover, sensitivity analysis demonstrates that Echinococcosis can be controlled with several strategies which can change four parameters A 1 , d 1 , β 1 and σ(see Figure 6). From Figure 3 we can see that our model matches the trend of monthly data for human Echinococcosis cases well. However, model (1) cannot simulate the oscillations of the data. Therefore, extension of our framework to include the stochastic disturbance is an interesting future work. We hope this work would lead to a better understanding of transmission of Echinococcosis.