Abstract
Echinococcus multilocularis, a major cause of echinococcosis in human, is a parasitic sylvatic disease between two major hosts in a predator-prey relation. A new model for the transmission dynamics of Echinococcus multilocularis in the population of red foxes and voles with environment as a source of infection is formulated and rigorously analyzed. The model is used to access the impact of treatment on red foxes and environmental disinfection as control strategies on the disease dynamics. The control reproduction number is computed and is used to rigorously prove the local and global dynamics of models’ equilibria. Using available data on Echinococcus, elasticity indices and partial rank correlation coefficients of control reproduction number and cumulative new cases in red foxes and voles are computed. Parameters that have high influence locally and globally are identified. Numerical experiments indicate that administering disinfection of environment only induces more positive impact than applying treatment only on red foxes in controlling the infection. Generally, interventions towards treating red foxes and environmental disinfection could be sufficient in tackling transmission of disease in the populations.
1. Introduction
Echinococcus multilocularis (EM) is a parasitic taeniid tapeworm and one of the six species of the Echinococcus genus. The disease is sylvatic and has indirect life cycle between two major hosts in a predator-prey interaction [1–3]. Adult EM inhabits the small intestine of canines (such as red foxes) which are regarded as definitive hosts and produces eggs that are released to the environment, typically through faeces. After oral ingestion of eggs by rodents (such as voles), regarded as intermediate hosts, a larval stage (metacestode) develops in any of the internal organs (liver, kidney, heart, etc). The mature metacestodes are capable of producing numerous protoscoleces, each having the potential to develop into an adult EM when a definitive host preyed on an intermediate host, and the cycle continues [1, 2, 4, 5]. In some parts of the world, other wild canids (such as coyotes, raccoon dogs, and wolves) can also serve as definitive hosts while other animals (like sheep, domestic dogs, cats, and rats) can be regarded as intermediate hosts [1, 4, 6]. The parasite (EM) causes alveolar echinococcosis in humans regarded as accidental hosts and is characterized with a tumour-like, destructive growth with the potential of causing high fatality rate [1, 4, 7]. The disease was initially confined to a certain part of the globe; however, as researches indicated, it spread all over the globe especially in rural nomadic communities that are economically less privileged, geographically and/or behaviourally detached to certain extent from healthcare systems [8–12]. Furthermore, the disease has high prevalence in red foxes population (1%–76.7%) [3, 5, 13–16] and low prevalence in rodents (voles) (0.4%–30%) [3, 5, 13, 15]. Due to the prevalence of the disease in intermediate hosts and the difficulties involved in treating definitive hosts in a given community, control can take longer times and in many cases may last indefinitely.
Modelling approaches, and in particular mathematical modelling, can give an insight into the biology and epidemiology of diseases in terms of revealing facts on data gaps, understanding the interaction between organisms, predicting future, control, and quantifying certain unmeasurable quantities such as force of infection, basic reproduction number, etc [2]. Mathematical models of Echinococcus multilocularis have been developed in literature to attend to must of the above assertions. For instance, in [5], Ishikawa et al. proposed a mathematical model that described the transmission of Echinococcus multilocularis in both the definitive (foxes) and intermediate host (voles) populations in Hokkaido. They quantitatively studied seasonal transition in the prevalence of EM in foxes and the risk of infection with human alveolar echinococcosis. Roberts and Aubert [3] developed a simple mathematical model in an attempt to determine the likely effect of combining treatment for EM infection in red foxes and voles in France with the existing vaccination campaign against rabies. It was shown that if the prevalence of the EM in foxes is low in a particular region, the parasite can be eradicated or controlled. Wang et al. [17] proposed a model for the transmission of echinococcosis in dogs, livestock, and human populations to explore effective control and preventive measures in Xianjing. Basic reproduction number, on which the dynamics of model was completely determined, has been estimated, and sensitivity analysis was carried out based on data relevant to the study area. In [18], a model that took into account, the contribution of domestic and stray dogs on the transmission of the parasite in humans was proposed. The global dynamics revealed that, without disposing the stray dogs, the disease became endemic even if the domestic dogs are controlled.
It is a known fact that environment aids the transmission dynamics of Echinococcus multilocularis, as discussed in the above reviews. However, this important component is neglected in the modelling process of the disease. In this paper, we establish a noble mathematical model for the transmission dynamics of Echinococcus multilocularis within the population of foxes as definitive hosts and the voles as intermediate hosts with concentration of parasites in the environment as a source of infection for intermediate hosts. Considering the high prevalence of EM in red fox population (1–76.7%) compared to rodent population as reported in previous studies, we incorporate treatment as control strategy on the infected red foxes. As a result of the treatment, we assume that recovered red foxes will acquire long-time immunity and as such will not return to susceptible. Furthermore, disinfection or cleaning of environment to reduce the concentration of the disease is also incorporated as a second control strategy. We carry out rigorous analysis on the computation of the basic control reproduction number, a threshold quantity used, to determine the existence and stability dynamics of equilibria. Furthermore, using data available from literature, we conduct elasticity indices of parameters on the control reproduction number and global sensitivity analysis using partial rank correlation coefficients of control reproduction number and cumulative new infections on the two hosts populations. Numerical simulations are conducted to support analytic results and effects of control strategies on the model. This paper is organized as follows. The model formulation, equations, and flow diagram are presented in Section 2. Basic properties of the model on existence, uniqueness, positivity, and boundedness of solutions are discussed in Section 3. Furthermore, existence and global stabilities with systematic calculation of control reproduction number are presented in Section 4. Numerics which comprise of elasticity index, global sensitivity analysis, numerical simulations, and effects of control strategies are presented in Section 5. Finally, we present concluding remarks in Section 6.
2. Model Formulation
The total population of red foxes, which is assumed constant (birth and death rates, , are assumed equal) in the environment at time t, denoted by is divided into susceptible (), exposed (), infected (), and recovered () subpopulations so that
The susceptible population is increased by recruitment of foxes by birth or immigration at rate and is decreased when it preyed with searching efficiency s on an infected vole containing protoscoleces in hydatid cysts [5] with probability of becoming infectious. The exposed fox population is increasing by the number of susceptibles that preyed on infected voles and is decreasing by progression to infected population and natural deaths, at rates and , respectively. The infected fox population is increased by progression of exposed foxes and decreased as a result of treatment and natural deaths at rates and , respectively. The population of recovered is increased by the treated infected foxes and decreased by natural death at rates and , respectively. We assume here that treated red foxes have either permanent or long-lasting immunity to the parasite and hence will not return to susceptible population.
Similarly, the total population of voles, also assume constant (birth and death rates, , are assumed equal) in the environment at time t, denoted by is subdivided into susceptible (), exposed (), and infected () subpopulations so that
The population of susceptible voles is increased by birth at rate and decreased by infection from the concentration of parasites in the environment at the rate and natural deaths at rate , which is also applicable to all the subpopulations of voles. Furthermore, the concentration of Echinococcus in the environment is increased by shedding of the parasites by infected foxes at rate and decreased by disinfection or cleaning of environment at rate . Based on the above descriptions, the model can be described completely by the following system of ordinary differential equations, which follow from the schematic diagram shown in Figure 1:with initial conditions
3. Basic Properties of the Model
The model (3)–(10) monitors the dynamics of red foxes and voles populations, and all its associated parameters are assumed nonnegative. Hence, we now present the basic results for the properties of the model.
Theorem 1. The following region is positively invariant for the model (3): , where , and .
Proof. The detailed proof of Theorem 1 is presented in Appendix A.
Equations (A.14)–(A.18) establish the boundedness of total populations for red foxes, voles, and concentration of parasites, respectively, and by extension verifies the boundedness of subpopulations. Thus, region is positively invariant. Hence, in this region, the model (3)–(10) is considered to be mathematically and epidemiologically well posed, and therefore, the dynamics of the model can be studied in .
4. Existence and Stability of Equilibria
4.1. Disease-Free Equilibrium (DFE)
The DFE of the model is obtained by equating the right-hand sides of model equations (3)–(10) to zero as follows:
4.2. Calculation of Control Reproduction Number
The basic reproduction number in epidemic models is an important threshold value that quantifies the infection risk in order to effectively control the disease. Furthermore, it plays a vital role in stability analysis of equilibria of the models. It can be derived using the next-generation matrix approach [19]. However, when there is intervention, it is referred as the control reproduction number. For detailed computation of the control reproduction number, refer Appendix B. Therefore, the basic control reproduction number, denoted by , is given by the following equation:
It is worth stating that is the basic control reproduction number, which represents the number of secondary infection cases generated by introducing at least one infective agent into the population that is assumed wholly susceptible. This number is obtained from the contribution of average number of secondary infections through fox-to-environment-to-vole transmission , voles-to-fox transmission , and environment-to-voles transmission as a result of one infectious subject during its infectious period.
4.2.1. Stability of DFE
Theorem 2. The disease-free equilibrium is locally asymptotically stable when and unstable if .
Proof. The local asymptotic stability of DFE is established using Theorem 2 in [19].
The epidemiological implication of the result in Theorem 2 is that the Echinococcus m. can be eliminated from the populations when if the initial sizes of subpopulations are within the basin of attraction of the DFE. In order to guarantee the total elimination of the disease irrespective of the initial population started with in , it is necessary to prove the global asymptotic stability (in ) of the DFE, which is shown below. Here, we use the matrix-theoretic method as described in [20].
Theorem 3. The DFE of the model (3)–(10) given by (12) is globally asymptotically stable (GAS) in when . If , the DFE is unstable, the system is uniformly persistent, and there exists at least one endemic equilibrium in the interior of .
Proof. Let and matrices be given as above. The matrix is computed to beObserve that the matrix is reducible (second, fifth, and sixth columns are the nonzero columns); hence, Theorem 2.2 of [20] is not applicable, and instead we use the conditions of Theorem 2.1 in [20] to construct the Lyapunov function. Define and so that the dynamics of infected compartments can be expressed as follows:wheresince and inside .
We define the Lyapunov function as , where is the left eigenvector of the matrix with respect to the eigenvalue . Thus,Multiplying and equating (17), we haveClearly, it can be seen that and , and we can deduce thatIf we take , one solution of is as follows:The Lyapunov function L is now given asDifferentiating L along the solutions of infected compartments in (3)–(10) givesUsing (13) and (21), simplifying and rearranging, we haveFrom (23), and , and then and . Therefore, if , it implies that . Furthermore, implies either , or . Thus, the largest invariant set where is the singleton (). Therefore, by LaSalle’s invariance principle [21], is GAS in if . Furthermore, if in (23), the first term is positive, while the second and third terms will be zero in when ; therefore, , and hence, is unstable. Using the argument in Theorem 2.2 of Shuai and van den Driessche [20], it can be shown that the instability of when implies that the system is uniformly persistent in , thus implying the existence of at least one positive endemic equilibrium.
4.3. Existence and Global Stability of Endemic Equilibrium
The existence of endemic equilibrium follows from the argument in Theorem 3. In the presence of disease in the community, the endemic equilibrium is obtained by setting the right-hand sides of equations (3)–(10) to zero, and thus,where
Theorem 4. If , then the unique endemic equilibrium of model (3)–(10) is globally asymptotically stable (GAS) in .
Proof. To prove the uniqueness and global stability of , we apply the method of graph-theoretic as described in Section 3 of [20]. Detailed proof of the theorem is also given in Appendix C.
5. Numerics: Elasticity Indices, Numerical Simulations, and Control Strategies
In this section, we use the parameter values in Table 1 with the aim of illustrating the theoretical results and quantifying the control measures for Echinococcus multilocularis.
5.1. Elasticity Indices
As evident from the expression of basic control reproduction number, in (13), it is interesting to know qualitatively and estimate quantitatively how perturbations of associated parameters have influence on . In order to achieve this, we determine the normalized forward sensitivity index as introduced in Chitnis et al. [22], otherwise called elasticity indices [23] of parameters on . This quantity for with respect to a parameter is defined as follows:
In Table 2, we compute the elasticity indices of with respect to the parameters at static baseline values as indicated in Table 1 and arranged in the descending order of magnitudes. The computations indicate equal influence of six parameters associated with incidences for transmission of the parasites (), disinfection rate (), and concentration of Echinococcus in the environment (K). Most importantly, treatment of red foxes () has the second largest value followed by incubation rates of voles and foxes in that order.
5.2. Global Sensitivity Analysis
From our result in Section 5.1, it is obvious that the local sensitivity analysis on could not explicitly differentiate the most influential parameters and thus the need for global sensitivity analysis. We adapt the approach in [24] to analyze the global sensitivity of the parameters on both and cumulative new cases in rodents and red foxes, respectively. Using the method of partial rank correlation coefficients (PRCC), as described and implemented in [25], we carry out the global sensitivity analysis of 9 parameters on the control reproduction number and cumulative new infections in the populations of both red foxes and rodents. The main objective is to determine the most influential parameters for the purpose of control and extent of infectivity in the two populations. To compute the PRCC values, we used the MatLab R2017b with ranges of parameters in Table 2 divided into 1000 sample sizes, and the results are displayed in Figure 2. The parameter with the PRCC value far away from zero indicates the more influential parameter is on both and cumulative new cases.
(a)
(b)
(c)
In Figure 2(a), the global sensitivity of parameters on is depicted. It can be seen that the rates of cleaning/disinfecting the environment () and rate of treating red foxes () have the most global influence on , followed by rate of red foxes contribution of E. multilocularis to the environment (). The global sensitivity of parameters on the cumulative number of new cases for red foxes is also displayed in Figure 2(b) which indicates that the incubation rate in red foxes () has the highest global influence, followed by the rate of searching efficiency of red foxes (s) and probability that an infected vole preyed on infects a red fox () in that order. Lastly, in Figure 2(c), the global sensitivity of parameters on cumulative new infection cases in rodents indicates that the transmission rate from environment to rodents () is the most global influential parameter, followed by the incubation rate in rodents (). From the global sensitivity analysis, for control purposes, it can be suggested that more emphasis should be given to cleaning/disinfecting the environment, for example, by removing carcass and administering praziquantel to red foxes.
5.3. Numerical Simulations
Figure 3 depicts the global stability of disease-free equilibrium as proved in Theorem 3 with different initial conditions, where the numbers of foxes, voles, and concentration of E. multilocularis converges asymptotically to the equilibrium point using different initial conditions. The parameter values in Table 1 are used so that the control reproduction number . It can be seen that all disease compartments () converge asymptotically to zero while the noninfected compartments () converge to their respective total populations.
In Figure 4, the time evolution for number of red foxes, voles, and the concentration of Echinococcus multilocularis for model (3)–(10) is illustrated using the parameter values in Table 1, except for depicting GAS of endemic equilibrium as proved in Theorem 4. It can be observed that populations of foxes, voles, and concentration of parasite converge asymptotically to their respective endemic equilibrium points irrespective of the initial population started with. Here, the control reproduction number . Hence, the disease will persist in the community.
5.4. Effects of Control Strategies on
In this section, the effects of treating red foxes () and the cleaning or disinfecting the environment () are going to be explored. So here it will be interesting to see how values of infectious red foxes and/or the will change as these parameters are varied when other parameters are fixed at the baseline values.
The effect of treatment-only on red foxes is illustrated in Figure 5(a) using the baseline parameter values in Table 1 except for . The value of is varied from to ; as a result, the number of infectious red foxes is reduced from 126 to 19, respectively, as displayed by the graphs. Similarly, the effect of controlling the red foxes by cleaning and/or disinfecting environment-only is shown in Figure 5(b) with fixed parameter values except , and the value of varied from to . It can be seen that the cumulative number of infectious red foxes decrease from 134 with to 8, , respectively. It is evident that the implementation of either of the two control strategies may not be adequate in eradicating the parasite completely from the community. Therefore, when the control strategies are administered simultaneously, as depicted in Figure 5(c), the cumulative number of infectious red foxes decreases from 125 with to 0, with , respectively. Hence, combining the two control strategies is more effectively followed by environmental cleaning/disinfection and treatment of red foxes. The later results agree with our elasticity indices in Section 5.1 for the two parameters ( and ).
(a)
(b)
(c)
Given that and are the control parameters in the model, it is important to see how varies as the two parameters are varied with others fixed using contour plots. The objective is to estimate values of and that will ensure disease eradication (making ) as stated in Theorem 2. The results are displayed as contour curves of as a function of treatment on red foxes () and cleaning/disinfection of environment () at fixed baseline values in Figure 6(a). The least values of and that will ensure parasites eradication are estimated to be 0.025 and 0.01 so that or 0.005 and 0.05 with , respectively. Furthermore, to access the impact of combined treatment and reduced contribution of parasite to environment by red foxes, contour plots of as function of the control strategies with varying rate of contribution by red foxes to the environment () are displayed in Figures 6(b)–6(d). The figures show remarkable increase in the associated control reproduction number with increase in rate of contribution of parasites to environment by red foxes. In Figure 6(b), low control strategies are needed if the rate of contribution () is very small to ensure almost total eradication of the parasites, with range of and . In Figure 6(c), with high contribution rate (), the control strategies must also be high to lower the value of with . However, when the rate of contribution () is moderate, in Figure 6(b), the control strategies must be in reciprocal combinations (low treatment rate versus high disinfection rate and vice versa) putting the range of with .
(a)
(b)
(c)
(d)
6. Concluding Remarks
A new global deterministic model for the transmission of Echinococcus multilocularis in the population of red foxes and voles with environment as a source of infection is formulated and used to access the impact of control strategies on the disease dynamics. Moreover, sensitivity analysis is carried out to determine the parameters that have influence on the control reproduction number and cumulative new infectious cases of red foxes and rodents. We start by investigating the basic properties of the model to ascertain its worthiness mathematically and epidemiologically. The major findings of the study are outlined as follows:(1)The disease-free equilibrium of the model is obtained and used to systematically determine the basic control reproduction number (). Furthermore, using a matrix-theoretic method, the DFE is globally asymptotically stable whenever is less than unity. The implication of this result is that infection of the parasite can be control in the community if can be reduced and maintain below unity.(2)When the control reproduction number exceeds unity, a unique endemic equilibrium exists, and using a graph-theoretic method, it is shown to be globally asymptotically stable if is greater than unity. Hence, the infection will persist in the community when the control strategies fail to reduce below unity.(3)Elasticity indices of on key parameters are computed, and the results indicate equal influence of six parameters associated with incidence for transmission of the disease (), cleaning/disinfection (), and concentration of parasite in the environment (K), followed by rate of treatment on foxes () and incubation rates of voles () and foxes () in that order.(4)Having noticed that the local sensitivity analysis on could not differentiate explicitly the most influential parameter(s) of the model, a global sensitivity using PRCC is conducted. From the simulations, the two control parameters: rate of cleaning/disinfecting the environment () and rate of treating red foxes have the most global influence on , followed by rate of red foxes contribution of E. multilocularis to the environment (). The global sensitivity of parameters on the cumulative number of new cases for red foxes indicates that the incubation rate in red foxes () has the highest global influence, followed by rate of searching efficiency of red foxes (s) and probability that an infected vole preyed on infects a red fox () in that order. On the contrary, the global sensitivity of parameters on cumulative new infection cases in rodents shows that the transmission rate from environment to rodents () is the most global influential parameter, followed by incubation rate in rodents ().(5)Numerical simulations with baseline values and varying and , have indicated that administering disinfection of environment only induce more positive control impact (leaving only 8 infected red foxes) compared to applying treatment only on red foxes (leaving about 19 infected red foxes). However, administering the two control strategies induce the most positive impact by treating all the infected red foxes.(6)Contour curves are used to estimate the least values of two control parameters that will ensure disease eradication, i.e., the value of to be reduced and kept below unity with baseline values of other parameters fixed. These values of and are estimated at 0.025 and 0.01 so that or 0.005 and 0.05 with , respectively. Moreover, contour plots of as function of the control strategies with the varying rate of contribution by red foxes to the environment () are computed. The results show remarkable increase in the associated control reproduction number with increase in rate of contribution of parasites to environment by red foxes.
It is worth remarking that amongst the seven most influential parameters globally (), five are associated directly with the red foxes, and this justifies our choice of treatment on the population. This also explains the high prevalence of the disease on red foxes as reported in literature.
Appendix
A. Proof of Theorem 1
A.1. Existence and Uniqueness of Solutions
Let be a column vector in so that Define to be the vector valued function in , where are the right-hand sides of model (3)–(10). The model system with initial conditions can therefore be expressed as follows:
Therefore, using a standard theorem of the dynamical system [26], f is the Lipschitz continuous in . Hence, there exist a unique solution of (3)–(10) for all times .
A.2. Positivity of Solutions
Here, we prove that any solution with the nonnegative initial condition that start in will remain there (i.e., it is positive at all times ).
Let be any solution of model (3)–(10) with initial conditions (11).
Consider equation (1):and using the variable separable method, we have
From (A.3), it can be seen that
From (7),so that using the same method as above, we get
Hence,
For the positivity of other variables, we proceed by the method of contradiction as follows: Suppose by contradiction that the conclusion is not true, then there exists a time such that
Suppose , then from (4), we haveso thatwhich contradicts our earlier assumption. Hence, no such exists. Therefore, for all . Suppose also that , from (5), it implies
Hence,which leads to contradiction. Therefore, for all .
Similar approaches can be followed to show that , and for all . This completes the proof. Therefore, all solutions of the model (3)–(10) are nonnegative for all .
A.3. Boundedness of Solutions
Here, we prove that the solutions are bounded. For the red foxes population, we add equations (3)–(6), and thus,
Since , this implies thatand hence bounded.
From the voles population, also adding equations (7)–(9), we will get
Again, since , this givesand hence bounded above.
Lastly, for the concentration of Echinococcus, since , from (10), we have
Hence,
B. Calculation of Control Reproduction Number
Using the approach in [19], we define the matrix F as the generation of secondary infectious cases and V as the matrix of transition rates between the infected compartments. The infected compartments are listed in order for easy application of the method as . Hence, F and V are calculated as follows:where .
The inverse of V is obtained as follows:and hence,
Therefore, the basic control reproduction number, denoted by , is obtained by taking the spectral radius of the matrix, and thus,
C. Proof of Theorem 4
Proof. At equilibrium, equations of system (1) can be expressed as follows:LetDifferentiating for along solutions of system (3)–(10) and using the inequality for with equality only if , we haveFrom the derivatives of for , it can be seen that all and the digraph does not have any cycle, as shown in Figure 7. Using the application of Theorems 3.3 and 3.5 of [20], so that gives a Lyapunov function for system (3)–(10), where ’s can be determined as follows:
Since , .
Given , .
From , .
From , Therefore,From ,And when ,It can easily be verified that {} is the only invariant set in the interior of that can satisfy and thus the uniqueness of . Using this Lyapunov function in combination with LaSalle’s invariance principle [21] completes the proof that is unique and globally asymptotically stable in if .
Data Availability
The type of data used in our research are numerical and all the data used are obtained from published literature which are adequately cited therein. Where data are not available (not cited), we use reasonable estimate of the value concern.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
A. S. Hassan acknowledges with thanks the motivation from Research Support Directorate, University of South Africa, for the Research Grant Fellowship number 409000.