A Schistosomiasis Model with Praziquantel Resistance

A compartmental model is established for schistosomiasis with praziquantel resistance. The model considers the impact of genetic resistance and drug treatment on the transmission of schistosomiasis. We calculate the basic reproductive number and discuss the existence and stability of disease-free equilibrium, boundary equilibrium, and coexistence equilibrium. Our analysis shows that regardless of whether drug treatment leads to the emergence of resistance, once the impact of genetic resistance is larger, the resistant strain will be dominant, which is detrimental to the control of schistosomiasis. In addition, once the proportion of human with drug-resistant strain produced by drug treatment is larger, the number of human and snails with resistant strain is larger. This is not a good result for drug treatment with praziquantel.


Introduction
Currently, treatment of human beings infected by schistosomiasis primarily focuses on chemotherapy with praziquantel (PZQ).PZQ appeared as a new schistosomicidal compound during the 1970s [1].In recent years, PZQ has become the drug of choice in most endemic areas because of its efficacy, its ease of administration, its tolerable side-effects, and its cost [1].Although the effectiveness of PZQ against schistosomiasis is well documented, the precise mode of action of the drug has not clearly defined [2].It is reported that the chemotherapy of many helminth infections is complicated by the occurrence of drug resistance and drug tolerance (a natural resistance) to certain anthelmintics [2].Not surprisingly, recent epidemiological evidence suggests the emergence of PZQ-resistant\tolerant schistosomes [1,3,4].Resistance is defined as a genetically transmitted loss of sensitivity in a parasite population that was previously sensitive to a given drug [2].Tolerance is an innate insusceptibility of a parasite to a drug, with the caveat that the parasite must not have been previously exposed to the drug [2].The first report of possible PZQ resistance came from an intensive focus in northern Senegal, where the drug had produced very low cure rates (18-39%) [5,6].And snails collected in the area carried schistosoma strains.When tested in the laboratory, those snails had a decreased susceptibility to PZQ [7,8].
Additional PZQ-resistant evidence was collected in Egypt [4,9].Preliminary studies have begun on these isolates to identify genetic, physiological, and morphological characteristics associated with PZQ resistance, and some of these may find use as markers for monitoring whether or not resistance is developing in endemic areas, where the drug is used [10].
Many papers have reported that drug treatment results in the emergence of schistosome resistance to PZQ [11][12][13].Drug treatment can remove drug-susceptible parasites in infected human beings, while resistant parasites survive.However, many investigations find that traits of PZQ resistance of Schistosoma mansoni are dominant inheritance [14][15][16].The resistant worms can reproduce and pass the resistant genes to the next generation.Furthermore, the resistance of Schistosoma mansoni to PZQ can be expressed in eggs, miracidia, cercariae, adults, and all stages of development [14][15][16].In other words, some definitive hosts carrying resistant schistosomes can infect snails and make those snails carry resistant schistosomes.On the contrary, some snails carrying resistant schistosomes can infect definitive hosts and also make those definitive hosts carry resistant schistosomes.It goes on, the control of schistosomiasis will face enormous difficulties.Therefore, it is necessary to study the impact of this schistosome genetic resistance on the transmission of schistosomiasis.
In previous schistosomiasis models, resistant problems have been studied by considering that the resistance of schistosomiasis is due to drug treatment [11][12][13].In [11,13], the authors proposed a multistrain schistosome model including sensitive and resistant parasite strains.Their goal was to infer the impact of drug treatment on the maintenance of schistosome genetic diversity.In their assumptions, the drugsensitive parasite strain had an additional per capita death rate, , due to treatment.For a parasite strain that had developed drug resistance with a resistance level  ( > 1), this treatment-related death rate was assumed to be reduced by the factor  to /.Their results implied that higher treatment rate could allow for coexistence between sensitive and resistant parasite strains.In [12], the authors formulated a deterministic model with multiple strains of schistosomes in order to explore the role of drug treatment in the maintenance of a polymorphism of parasite strains that differed in their resistance levels.And snails infected by parasite strains were divided into multistrain subclasses according to the different level   of parasite strains.Analysis of the model showed that the likelihood that resistant strains would increase in frequency depended on the interplay between their relative fitness, the cost of resistance, and the degree of selection pressure exerted by drug treatments.
Motivated by [11][12][13], we establish a new model considering hosts with sensitive and resistant strains in this paper.Our purpose is mainly to study the impact of drug treatment and genetic resistance on the transmission of schistosomiasis.
Our paper is organized as follows.In Section 2, we establish a mathematical model with praziquantel resistance and obtain basic reproductive number and existence of equilibria.And then the stability of disease-free equilibrium is obtained in Section 3. Section 4 devotes to stability analysis of boundary equilibria.In Section 5, stability analysis of endemic equilibrium is performed.

Mathematical Model
According to different level of parasite strains, we divide infected hosts into sensitive and resistant strains.Considering resistance and inheritance of resistance, new resistant strains are composed of two parts.We classify definitive and intermediate hosts as susceptible, sensitive and resistant in the following: (i)   (), the population of susceptible human; (ii)   (), the population of human infected with sensitive parasite strain; (iii)   (), the population of human infected with resistant parasite strain; (iv)   (), the population of susceptible snail host; (v)   (), the population of snail host carrying sensitive parasites; (vi)   (), the population of snail host carrying resistant parasites.
We follow some of the available models for schistosomiasis [17][18][19] and assume that the reproduction rate of hosts is constant, and we ignore the recovery class of host since the life span of infected is short in comparison to that of human [18,20].All parameters in the model are assumed to be nonnegative constants: (ix) , the rate of treatment for infected human; (x) , the proportion of human with drug-resistant strain produced by treatment; (xi) , represents the impact of inheritance and the cost of resistance on transmission rate, we assume that  ≤ 1 since the transmission rate is reduced due to resistance [12]; (xii)  1 , represents the impact of resistance on diseaseinduced death rate of human; (xiii)  2 , represents the impact of resistance on diseaseinduced death rate of snails.
Then, we have a model with the form: Using standard methods, it is easy to see that disease free equilibrium  0 = (Λ 1 / 1 , 0, 0, Λ 2 / 2 , 0, 0) always exists.Let According to the concept of next generation matrix [21] and the formula of the basic reproductive number for ODE compartmental models [22], if we let one can calculate that the eigenvalues of the next generation matrix  −1 are given by Then, it follows that the basic reproductive number for the system (1) is given by To obtain other equilibria, we let the right-hand side of (1) equal to zero and obtain If   = 0, we have   = 0 and a formula of   as follows: Note that there exists a boundary equilibrium with only resistant type, given as  0 = (  , 0,   ,   , 0,   ), where If drug treatment does not lead to drug resistance, that is,  = 0, then when  0 > 1, we can obtain the other boundary equilibrium with only sensitive type, given as  0 = (  ,   , 0,   ,   , 0), where Now, we study existence of coexistence equilibrium for the system (1).From ( 5)- (10), we obtain Following ( 6), (7), (9), and (10), we have In the case that  = 0, (15) leads to Note that  0 =  0 equals to  1 /  (14) represents the existence of coexistence equilibrium in the form of a line.
The following section shows that the basic reproductive number  0 provides a threshold condition for schistosoma extinction in (1).

Stability Analysis of the Disease Free Equilibrium
In this section, we will analyze stability of the disease free equilibrium of the model ( 1).The stability of the disease free equilibrium determines whether schistosomiasis will be permanent in an uninfected population.The following result shows that schistosome will go extinct if  0 < 1.
Theorem 2. The disease free equilibrium  0 of the system (1) is locally asymptotically stable if  0 < 1 and unstable if  0 > 1.
Proof.The Jacobian matrix for the system (1) is given by Then, the eigenvalues of  0 are − 1 , − 2 and roots of the following equations: Note that  0 < 1 equals to  0 < 1 and  0 < 1, which leads to Hence, if  0 < 1, then all the roots of (20) have negative real parts.Hence, using the Routh-Hurwitz criterion, we can obtain that the disease free equilibrium  0 of the system (1) is locally asymptotically stable if  0 < 1 and unstable if  0 > 1.
Now, we turn to the study of the global stability of the disease free equilibrium of the model (1) by using Metzler matrix theory and the technique of Kamgang and Sallet [23].
Consider systems of the following form: where + , and  and  are  1 .We denote by  = ( 1 ,  2 ) the state of the system and ( * 1 , 0) is a disease free equilibrium on a positively invariant set Ω ⊂ . Now rewrite (21) as For the system (22), we make the following assumptions.
(h 1 ) The system is defined on the positively invariant set Ω of the nonnegative orthant.The system is dissipative on Ω.
) is globally asymptotically stable at the equilibrium  * 1 on the canonical projection of Ω on  (h 4 ) There exists an maximum matrix  2 , then for any  ∈ Ω such that  2 =  2 (),  ∈ (h 5 ) ( 2 ) ≤ 0, that is, the greatest real part of eigenvalues of  2 is nonnegative.
For convenience, we state two lemmas due to Kamgang and Sallet [23].Next, we discuss the global stability of the disease free equilibrium  0 of the system (1) using the above two Lemmas.

Lemma 3. If the above hypotheses, h
From the system (1), we know This proves that the set is a compact positively invariant absorbing set contained in the nonnegative orthant.Thus, the system (1) is dissipative on Ω because the trajectories of (1) are forward bounded.Now, we will study the system (1) on Ω.
We set for system (1)  1 = (  ,   ),  2 = (  ,   ,   ,   ), and As in [23], we express the subsystem as This is a linear system, and its unique equilibrium ( 1 / 1 ,  2 / 2 ) (corresponding to the disease free equilibrium of ( 1)) is globally asymptotically stable, hence the assumptions (h 1 ) and (h 2 ) are satisfied.The matrix  2 () is given by As required by hypothesis h 3 , for any  ∈ Ω, the matrix  2 () is irreducible.Now, let us check (h 4 ).There is a maximum which is uniquely realized in Ω if   = Λ 1 / 1 and   = Λ 2 / 2 , which corresponds to the disease free equilibrium.This maximum matrix is then  2 , the subblock of the Jacobian matrix at the disease free equilibrium, corresponding to the matrix  2 ().The matrix  2 is given by Therefore, we are in the situation of Lemma 4, where the maximum is attained at the disease free equilibrium.The hypothesis (h 5 ) requires that ( 2 ) ≤ 0. Writing  2 as a block matrix  2 = ( Since  is already a Metzler stable matrix, the condition ( 2 ) ≤ 0 is equivalent to the condition ( −  −1 ) ≤ 0 [23], where ) . ( Then, Hence, the condition ( −  −1 ) ≤ 0 is equivalent to  0 ≤ 1.We have seen that the hypotheses (h 1 ), (h 2 ), (h 3 ), (h 4 ), and (h 5 ) are satisfied.Then, by Lemma 4, we have the following result.
Theorem 5.The disease free equilibrium  0 of the system (1) is globally asymptotically stable if  0 ≤ 1.

Stability Analysis of the Boundary Equilibria
In this section, we turn to study stability of the two boundary equilibria.From Lemma 1, we know that if  0 > 1, there exists a boundary equilibrium with only resistant type  0 .Through calculations, we can obtain the characteristic equation as following: Hence, the eigenvalue of  0 are roots of the following equations: Here, From (12), we can obtain It is easy to see that  2 0 / 2 0 = ( 2  4 / 2 )/ 1  3 .Then,  0 <  0 equals to  1  3 >  1  2     , which implies that the roots of (31) have negative real parts if Here,  3 > 0 and Hence, we can obtain It follows from Routh-Hurwitz criterion that all roots of (32) have negative real parts if  0 > 1. Summering above analyses, we have the following result.
Now, we study the global stability of the boundary equilibrium  0 .Consider the Lyapunov function where . The Lyapunov derivative is This implies that the sensitive type dies out if  0 ≤ 1.Then the largest compact invariant set of the system (1) in the set Using the LaSalle-Lyapunov theorem, we know that all trajectories in Ω eventually tend to Ω 1 as  → ∞.Then, we only need to study the dynamical behavior of (1) in Ω 1 .At this time, (1) reduces to the following system To show that all trajectories of (39) in the interior of Ω 1 approach the point (  ,   ,   ,   ) corresponding to the boundary equilibrium  0 , consider the Lyapunov function where positive constants  1 and  2 are defined in the following.It is easy to see that  2 ≥ 0 for (  ,   ,   ,   ) ∈ Ω 1 , and  2 = 0 ⇔ (  ,   ,   ,   ) = (  ,   ,   ,   ).Hence, the function  2 is positive definite with respect to the point (  ,   ,   ,   ).
Computing the derivative of  2 along solutions of system (39), we have Substituting Λ 1 = It is easy to see that where  is an arbitrary positive number.Substituting   =   ,   =   ,   =   , and   =   into the first equation of system (39), we obtain Now, we turn to the other boundary equilibrium.From Lemma 1, we know that if  = 0 and  0 > 1, the boundary equilibrium with only sensitive type  0 exists.Through calculations, we can obtain the characteristic equation as following: Hence, the eigenvalue of  0 are roots of the following equations: Here, From ( 13), we can obtain that  1  2     =  2  4 .Then, it is easy to see that the roots of (48) have negative real parts if  0 >  0 .Similarly to the case of  0 , using Routh-Hurtwitz criterion, we can obtain that all roots of (49) have negative real parts if  0 > 1. Summering above analysis, we have the following result.Theorem 8.When  = 0, the boundary equilibrium  0 of the system (1) is locally asymptotically stable if  0 > 1 and  0 >  0 .Now, we study the global stability of the boundary equilibrium  0 .Consider the Lyapunov function where . The Lyapunov derivative is This implies that the resistant type dies out if  0 ≤ 1.Then, the largest compact invariant set of the system (1) in the set Using the LaSalle-Lyapunov theorem, we know that all trajectories in Ω eventually tend to Ω 2 as  → ∞.Then, we only need to study the dynamical behavior of (1) in Ω 2 .At this time, (1) reduces to the following system: To show that all trajectories of (53) in the interior of Ω 2 approach the point (  ,   ,   ,   ) corresponding to the boundary equilibrium  0 , consider the Lyapunov function where positive constants  1 and  2 are defined in the following.
Computing the derivative of  2 along solutions of system (39), we have Substituting Note that for the limiting system Hence, for all (  ,   ,   ,   ) ∈ Ω 2 , It is easy to see that Therefore, the only compact invariant subset of the set where  2 / = 0 is the point (  ,   ,   ,   ), corresponding to the boundary equilibrium  0 .By LaSalle's Invariance Principle,  0 is globally asymptotically stable if  0 > 1 and  0 ≤ 1.
Summering above analysis, we have the following result.

Stability Analysis of the Coexistence Equilibrium
In this section, we turn to study the local stability of the coexistence equilibrium   in the limiting system of (1) by using Krasnoselskii sublinearity trick [24], as in [25,26].In detail, if   = () is a system of differential equations and  * is an equilibrium point, then to prove the local asymptotical stability of  * is to prove that the linearized equation   =   ( * ) has no solutions of the form with  0 ∈   ,  ∈  and, Re  ≥ 0. This implies that the eigenvalues of the characteristic polynomial associated with the linearized equations have negative real part, that is, Re  < 0.Then, the coexistence equilibrium   is locally asymptotically stable.
Considering the limiting system where  1 = Λ 1 / 1 and  2 = Λ 2 / 2 .In this way, let  0 = ( 1 ,  2 ,  3 ,  4 ),   ∈ .Substituting a solution of the form (60) into the linearized system (61) of the coexistence equilibrium   , we obtain the following linear equations: which is equivalent to the system Moving all the negative terms to the left-hand side, after some manipulations we obtain the system where ) . (65) Note that the matrix  has nonnegative entries, and   = (  ,   ,   ,   ) satisfies To show that Re  < 0, we distinguish two cases:  = 0 and  ̸ = 0.In the first case, (62) is a homogeneous linear system.Through calculations, we have the determinant of ( 62) is where Then, we can obtain Note that the coexistence equilibrium   exists under the condition which contradicts the minimality of .Hence, Re  < 0. Summering above analysis, we have the following result.
The existence and stability of equilibria can be summered in Table 1.

Discussion
In this paper we established a new schistosomiasis model.In contrast to previous schistosomiasis models with drug resistance, the model established in this study consider many aspects.First, snail is considered as a variable in the model since the resistance of schistosoma to PZQ can be expressed in snails [14][15][16].Second, previous models considered that  resistance was caused by drug treatment, while a large number of the literature show that some cases are due to inheritance [1,3,4].Therefore, in this study, we consider the reasons for resistance are drug treatment and genetic.Last, in previous models a resistance level () was used to discount the treatment rate.But the resistance level could not be measured, and its value was assumed.This model consider what percentage of infected human after treatment will recover and what percentage will emerge drug resistance.In medicine, from the occurrence of cases the value of this ratio can be identified.For example, this ratio was given in [27,28] (1 −  = 0.28 ∼ 0.609).Therefore, it is easy to operate.In addition, we separate the discussion of the case that treatment will cause resistance ( ̸ = 0) and treatment will not cause resistance ( = 0).The reproductive number  0 and  0 of the sensitive and resistant strains are given, respectively.It is easy to see that  0 is a decreasing function of , and  0 is an increasing function of .If the basic reproductive number  0 = max{ 0 ,  0 } of the model (1) is less than 1, one can prove the stability of the disease free equilibrium.This means that the spread of schistosomiasis can be effectively controlled.When the basic reproductive number  0 is greater than 1, we first consider the case that  = 0.If  0 =  0 > 1, the two strains can coexist.There is a line of coexistence equilibria in this case.The infected human will evolve to one of them with higher reproductive number (see Table 1).
When drug treatment can not cause resistance ( = 0), that is, the new resistant strain is due to the inheritance of resistance, the sensitive strain will dominate if the treatment rate is smaller and the impact of the inheritance of resistance is smaller such that  0 > 1 ≥  0 .This result accords to the results of previous models.On the other hand, if the treatment rate and the impact of the inheritance of resistance are both larger such that  0 > 1 ≥  0 , the resistant strain will dominate.This shows that although the assumption is that drug treatment does not result in the emergence of drug resistance, once the treatment rate is greater than a value, and the impact of genetic resistance is larger, there will still be the emergence of resistant strain, and the resistant strain is dominant.This further shows that genetic resistance has a great impact on the system.
When drug treatment can cause resistance ( ̸ = 0), we can show that the sensitive strain either does not appear or coexist with the resistant strain under certain condition.If  0 > 1, there is only resistant strain.It can be seen that, regardless of whether drug treatment leads to the emergence of resistance, once the impact of genetic resistance is larger, resistant strain will be dominant, which is detrimental to the control of schistosomiasis.
Finally, from the formula of the coexistence equilibrium   , it is easy to see that the value of the resistant strain is increased with the value of .This means once the proportion of human with drug-resistant strain produced by drug treatment is larger, the number of human and snails with resistant strain is larger.This is not a good result for drug treatment with praziquantel.
Hence, for poor treatment, there are two possible reasons: drug therapy and genetic.An important priority in developing new control strategies is to search new drug targets, in combination with selection of genetic methods such as that viable vaccine candidates.And there is already a need for alternative drugs to treat PZQ-resistant schistosomiasis, such as already exists in northern Senegal [10].

Table 1 :
Existence and stability of equilibria.The sign "∃ !" means the existence and uniqueness.LAS means locally asymptotically stable, and GAS means globally asymptotically stable.