Mathematical analysis of a lymphatic filariasis model with quarantine and treatment

Background Lymphatic filariasis is a globally neglected tropical parasitic disease which affects individuals of all ages and leads to an altered lymphatic system and abnormal enlargement of body parts. Methods A mathematical model of lymphatic filariaris with intervention strategies is developed and analyzed. Control of infections is analyzed within the model through medical treatment of infected-acute individuals and quarantine of infected-chronic individuals. Results We derive the effective reproduction number, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\mathcal {R}_{0},$\end{document}R0, and its interpretation/investigation suggests that treatment contributes to a reduction in lymphatic filariasis cases faster than quarantine. However, this reduction is greater when the two intervention approaches are applied concurrently. Conclusions Numerical simulations are carried out to monitor the dynamics of the filariasis model sub-populations for various parameter values of the associated reproduction threshold. Lastly, sensitivity analysis on key parameters that drive the disease dynamics is performed in order to identify their relative importance on the disease transmission.


Background
Lymphatic Filariasis commonly known as elephantiasis is a globally neglected tropical parasitic disease caused by a thread-like worms of the Filarioidea type (Wuchereria bancrofti, Brugia malayi and Brugia timori) [1]. The most common of these, Wuchereria bancrofti, is a round worm that mainly infects the lymphatic system. Lymphatic filariasis involves asymptomatic, acute and chronic conditions with the majority of infections being asymptomatic [2]. Nearly 1.4 billion people in 73 countries worldwide are threatened by the disease of which over 120 million individuals are currently infected [3]. The round worm (nematode) is spread from person to person via a mosquito vector and infected individuals can suffer from chronic conditions such as lymphedema, elephantiasis and, in men, swelling of the scrotum called hydrocele  [1,2]. A description of the microfilariae life cycle is depicted in Fig. 1.
Lymphatic filariasis is still a major public health problem in Africa, South America and Asia despite existing knowledge of the disease pathology and global treatment campaign [3] with drugs such as Diethylcarbamazine plus Albendazole and Ivermectin plus Albendazole that kill the microfilariae and some of the adult worms. There is not enough evidence on effectiveness of the drug Albendazole, alone or in combination, for killing or interrupting transmission of threadlike worms that cause lymphatic filariasis [4]. Some studies have shown that treatment with Doxycycline could completely kill microfilariae [2].
Eradication of this disease has been a great challenge [3,5]. Thus, investigating the impact of combined intervention strategies of treatment and quarantine of chronically infected persons is viable. Chronically infected individuals may not transmit infection when quarantined from the rest of the population [6,7]. Insecticide treatedbed nets (ITNs) and sleeping in indoor residual sprayed houses (IRS) could reduce contact between humans (especially microfilariae carriers) and mosquito vectors [8].  the disease, despite uncertainty about the microfilarial prevalence threshold level below which transmission cannot be sustained even in the absence of any treatments [3,9].
Compartmental mathematical models of lymphatic filariasis abound in the literature [6,[9][10][11]. Two general simulation models of lymphatic filariasis transmission and control used to support decision-making arethe population-based deterministic model (EPIFIL) [12,13] and -the individual-based stochastic model (LYMFASIM) [14]. However these models have some limitations as they do not account for intervention measures such as quarantine. While EPIFIL uses a constant forceof-infection and accounts for the impact of age structure of the human community [13], LYMFASIM accounts for the role of the immune system in regulating parasite numbers [15]. Luz et al., [16] noted that mathematical modeling of transmission dynamics and cost-effectiveness of neglected diseases can help to maximize the utility of the limited available resources. Bhunu and Mushayabasa [10] considered treatment as the only intervention strategy in their model, while Ottesen et al. [8] presented strategies and tools to control transmission and morbidity of lymphatic filariasis. Although various transmission and control mathematical models of lymphatic filariasis abound in the literature, our proposed model is seemingly new as it includes latent stage, treatment and quarantine of chronically infected persons [8,17]. The latent stage is included in the model because of different developmental stages the worm undergoes in human and mosquito populations. The proposed compartmental model is not exhaustive, and here are some limitations: no density dependent and species-specific parasite prevalence [18], additional mortality experienced by infected mosquitoes as a result of carrying filarial infection [19]. Pichon [20] noted that mosquito density-dependent mortality may be associated with increased infection intensity within the mosquito and mass drug administration may lead to an increase in survival of the mosquito population and hence to an increase in transmission in the long-term [20].
In the following sections, we formulate and analyse a deterministic model with two key control measures: quarantine and treatment. Key parameters that influence transmission are identified via sensitivity analysis of the model. Finally, some parameter values are assumed within realistic ranges to support the analytical results, but with one caveat that the model outcomes are not compared with real data.

Methods: model formulation and description
Human and mosquito populations are divided based on their lymphatic filariasis status. Human sub-populations are susceptible humans S h (t), latent stage (not showing signs of lymphatic filariasis) E h (t), infected-acute stage I ha (t) and infected-chronic stage I hc (t), with the total human population given by The model is formulated with the assumption that no infection exists at the initial stage, and there is no vertical transmission in both human and mosquito populations [17,21]. In addition, the model considers one species of worm and one species of mosquito. We also assume that the transmission to mosquito population is from infected-acute and infected-chronic individuals despite the quarantine of some infected-chronic individuals.
The mosquito population is divided into three subgroups: susceptible S v (t), exposed E v (t) and infected I v (t), with the total mosquito population given by The recruitment rate of human population is h , while v is the recruitment rate of the mosquito population. The natural death rates of human and mosquito populations are μ h and μ v respectively. These death rates are proportional to the number of each individual or mosquito class. The biting rate of the mosquitoes to humans is β. The microfilariae which are found in lymphatic vessels and lymphatic nodes infect susceptible mosquitoes when a mosquito bites infected-acute and infected-chronic individuals at a rate where ϑ v is the success rate of microfilariae transmission from human to susceptible mosquitoes and θ ∈ (0, 1) accounts for the reduced number of adult microfilariae in humans due to treatment and quarantine of the infectedchronic individuals. The vector will ingest microfilarial differently when it bites humans in the I ha and I hc stages. The rate of release ϑ v θ of microfilarial by I hc (t) into the vector is different from the rate of release of microfilarial by I ha (t). Therefore the microfilariae ingested by vectors during a blood meal depend on the density of microfilariae in humans [11]. Thereafter, susceptible mosquitoes enter the exposed class E v (t). During this stage, the microfilariae develop into infective filariform larvae to become infectious, and hence these mosquitoes move into the infected class I v (t) at rate α v . The larvae infect the susceptible human host during a subsequent blood meal by the infected mosquitoes at a rate where ϑ h is the success rate of transmission of infective filariform larvae from infected mosquitoes I v (t) biting susceptible individuals during a blood meal. Individuals during the exposed stage have infective filariform larvae which migrate to lymphatic vessels and lymphatic nodes, and develop into adult worms. The latent individuals progress to infected-acute individuals at a rate α h when the microfilariae develop into adults which remain in the lymphatic vessels and lymphatic nodes. Furthermore, the infected-acute individuals who progress to chronic condition are quarantined at symptomatic rate κ to join the infected-chronic class. Individuals in the infected-acute class, I ha , are screened by health personnel at the rate n and treated at the rate ϕ. Treated infected-acute individuals join the susceptible class due to temporary immunity at a rate π = ϕn. Figure 2 provides a graphical interpretation of the lymphatic compartmental model (3). Based on our model description, assumptions, definitions of the state variables and parameters in Table 1, the proposed SEIS lymphatic filariasis model satisfies the following system of nonlinear ordinary differential equations:

Invariant region
Both the model state variables and parameters are assumed non-negative for all time t ≥ 0. Let be any solution of the system with non-negative initial conditions. Applying Similarly, it can be shown that the feasible solutions on the mosquito population given by Eq.
Therefore, from (4) and (5), the possible solutions of model (3) will enter the the positively invariant region

Positivity of the state variables
Since is a positively invariant set under the flow induced by model (3), we now show that every solution with initial condition in R 7 remains in that region for t > 0.
The first equation of model (3) gives Intergrating with respect to t gives Therefore, Hence S h is always positive for t > 0. The second equation of model (3) gives Similarly it can be shown that

Existence and stability of steady-state solutions
The disease-free equilibrium (DFE) of the lymphatic filariasis model (3) denoted by E 0 is given by The effective reproduction number is obtained by using the next generation matrix [23]. Let The effective reproduction number is the spectral radius ρ(FV −1 ) and the resulting expression is given by which is the number of secondary lymphatic filariasis infections caused by one infectious individual/mosquito during the infectious period in a completely susceptible population. The effective reproduction number is not only important for describing how fast the disease could spread, but can also provide information for controlling and preventing the spread of the disease [24].

Local stability of the disease-free equilibrium
Local stability of the DFE can be established from Theorem 2 in [23].

Lemma 2
The DFE for the lymphatic filariasis model (3) is locally asymptotically stable if R 0 < 1 and unstable when R 0 > 1.

Global stability of the disease-free equilibrium
The system of Eq. (3) is broken into subsystems such that X 1 = (S h , S v ) which denotes the number of susceptible individuals and susceptible mosquitoes, and Y 1 = (E h , I ha , I hc , E v , I v ) which denotes the number of exposed and infected individuals and mosquitoes. Hence the model system (3) now reduces to This could further be simplified by identifying X 1 with (X 1 , 0) and Hence we obtain the reduced system dX 1 dt = F(X 1 , 0) as Therefore asymptotically stable equilibrium for the reduced system dX 1 dt = F(X 1 , 0). This is verified by integrating the first equation of the reduced system (7) with respect to t.
Similarly integrating the second equation of the reduced system (7) with respect to t, Further G(X 1 , Y 1 ) satisfies the two conditions given as assumptions H3 and H4 in [25] namely: G(X 1 , 0) = 0 and G(X 1 , Here we assume a steady-state value of the total human Thus the DFE is globally asymptotically stable sinceG(X 1 , Y 1 ) is non-negative. The global stability of E 0 excludes any possibility of the phenomenon of backward bifurcation, that is the co-existence of a stable disease-free equilibrium with a stable endemic equilibrium [26].

Existence of endemic equilibria
To determine the existence of an equilibrium for which filariasis is endemic in the population defined by (3) is solved in terms of the force of infection at steady-state (λ * h ), given by Solving the system at an arbitrary equilibrium, we have . Substituting I * v and the above solutions into the expression for λ * h , we have This can be written as where The root λ * h = 0 corresponds to the DFE and its stability has already been established in Lemma 2. It is clear that A > 0, and B > 0 if and only if R 0 > 1. Thus the linear system would have a unique positive solution given by λ * h = B/A. The components of the endemic equilibrium, E * 0 , are then determined by substituting λ * h = B/A. For R 0 < 1, B < 1. Thus, the force of infection (λ * h ) at steady-state is negative (which is biologically meaningless). Hence, the model system has no positive equilibria in this case i.e when λ * h < 0. This result is summarized in the following Lemma.

Lemma 3
The lymphatic filariasis model has a unique endemic equilibrium whenever R 0 > 1, and no endemic equilibrium otherwise.

Sensitivity analysis and model simulations
The following local sensitivity analysis is closely related to that in [27]. Expressions for the sensitivity indices of the endemic equilibrium are complex, and since our focus is on disease transmission and not prevalence, we neither derive expressions nor numerically calculate sensitivity indices of the endemic equilibrium.

Sensitive indices of R 0
The sensitivity indices allow us to measure the relative change in a state variable when a parameter changes [27]. The normalized forward sensitivity index [27] of a variable, ψ, that depends differentiably on a parameter, p, is defined as: Next, we evaluate the sensitivity indices at the parameter values given in Table 2. The resulting sensitivity indices are shown in Table 3.
The most sensitive parameter to R 0 is the mosquitoes natural death rate, μ v (ϒ R 0 μ v = −1.360). This is followed by the mosquito per capita biting rate, β (ϒ R 0 β = 1). Reducing this parameter would have a huge effect on filarisis transmission regardless of other parameter values. We have that ϒ R 0 β = 1, then decreasing (or increasing) β by 10% decreases (or increases) R 0 by 10%.
The parameters are ordered from most to least sensitive for each R 0 Other key parameters include the success rate of microfilariae transmission from humans to susceptible mosquitoes, ϑ v , as well as the success rate of transmission of infective larvae from infected biting mosquitoes to susceptible individuals during a blood meal, ϑ h . The sensitivity indices of these two parameters are equal and independent of other system parameters. With ϒ R 0 α v = 8.19 x 10 −5 , the progression rate of mosquitoes from exposed to infectious state, α v , is the least sensitive parameter. For ethical reasons, one should only attempt to decrease human death rate and for this reason, we do not calculate the sensitivity index related to individuals natural death.
For the mosquito recruitment rate v , the filariasis reproduction number, R 0 , increases as v increases. If the mosquito recruitment increases, so does the mosquito death rate because the environment can only support a certain number of mosquitoes. Further, when the mosquito recruitment rate v is equal to the death rate μ v , the mosquito population is at equilibrium. If 1/ v is the life span of the mosquitoes, then increasing v reduces their life span. Reducing the life span of the vector population reduces R 0 , as more infected mosquitoes die before they become infectious.

Analysis of R 0
The objective here is to determine, using the threshold quantity R 0 , whether or not activities provided to quarantined infected-chronic individuals (modelled by κ) and treatment (modelled by ϕ) of infected-acute individuals can lead to the elimination of lymphatic filariasis in the community. It is evident from (6) that Thus, a sufficient effective quarantine programme (morbidity management and disability prevention activities) that focuses on quarantining infected individuals in the I ha stage (at a high rate, κ → 1) can lead to effective disease control if it results in the right hand side of (12) being less than unity. Likewise, for an effective treatment program, the right hand side of (13) should be less than unity. The profiles of R 0 as a function of the quarantine rate, κ, and treatment rate, ϕ, are depicted in Fig. 3a and  3b.
For the set of parameter values used, the strategy that focuses on treating the infected-acute individuals alone can dramatically reduce R 0 from around R 0 = 0.264 to R 0 = 0.077. The quarantine strategy increases R 0 = 9.795 x 10 −3 to R 0 = 0.329. Thus, lymphatic filariasis in the community could be reduced more slowly in the latter case, but will be eliminated faster in the former case. Figure 4a shows that the combined strategy of activities provided to quarantined infected-chronic individuals with an effective treatment of infected-acute individuals reduces R 0 to values far below unity than when each strategy is applied singly. Figure 4b shows that with an effective treatment strategy, increasing the quarantine rate does not necessary reduce the burden of filariasis in the community.

R 0 as a function κ and ϕ
The lymphatic filariasis burden in the community is evaluated by computing the partial derivatives of R 0 with respect to the quarantine and treatment parameters (κ and ϕ respectively). This gives Consider the case when ϕ = 0 (there is no treatment but only quarantine). It follows that

Lemma 4
The targeted quarantine strategy of infectedacute individuals will have a positive impact if θ < I , no impact if θ = I , and will have detrimental impact in θ > I .
Thus, we have Note that if conditions (16) and (17) are invalid, then application of the activities provided to quarantined infected-chronic individuals and treatment strategies would increase the burden of filariasis in the community (since it increases R 0 ). Treatment would increase the disease burden if it fails to reduce the infectiousness of those treated below a certain threshold (θ > T if treatment of infected-acute individuals is targeted or θ < I if quarantine is targeted).

Model simulations
Numerical simulations of the model system (3) are carried out using Wolfram Mathematica 9.0 to illustrate some of the analytical results. Parameter values used for the model simulations are provided in Table 2, some of these were obtained from the literature [27][28][29][30] while others were assumed (within realistic range) for the purpose of simulations. The dynamics of the human and mosquito populations when both treatment and quarantine are employed, are depicted in Fig. 5a and 5b, respectively. The effects of increasing the infected-acute individuals quarantine rate as well as scaling up the treatment of acute-infected individuals on the dynamics of the whole population are explored. Two scenarios are investigated: (1) we assume that a population is invaded by infected mosquitoes, and (2) we assume that a population is invaded by acuteinfected humans.

Case 1.
Assume that a population is initially disease-free and stays at equilibrium. This population is then infiltrated by 10 infected individuals in the latent state. At this point we assume that all the mosquitoes are also disease-free. Using the parameters in Table 2, Fig. 6a and 6b show the dynamics of the human and vector populations when there is no treatment or quarantine (selective treatment) (ϕ = κ = 0). The basic reproduction number is 0.482, which means that the disease will eventually  Table 2 for a the human, b the vector sub-populations over time  Fig. 7a and 7b show the same dynamics as Fig. 5a and 5b when there is treatment of infected-acute individuals but no quarantine. Figure 8a and 8b show the effects of using quarantine of infected-chronic individuals as the only control measure. A strategy that uses both treatment and quarantine (selective treatment) is better than using only treatment or quarantine).
The cases in Figs. 6a, 7 and 8b assume that the treatment given to acute-infected individuals does not affect the transmission parameters other than the screening parameter n and treatment parameter ϕ. The basic reproduction number depends on these two parameters and so they do affect initial disease transmission and consequently the endemic status of the disease. If patients are given medication that reduces the microfilariae in the lymphatic vessels and nodes, then the success rate of microfilariae transmission from humans to susceptible mosquitoes will decrease. Effective treatment is expected to decrease the number of microfilariae, and thus its transmission success rate from humans to susceptible mosquitoes could greatly be be altered. Assuming that individuals are further protected by some insect repellent, then the mosquito biting rate β could also be impacted. Consider a medication that could decrease the success rate of microfilariae transmission from humans to susceptible mosquitoes by 50% of its current level (ϑ v = 0.05). Figure 9a and 9b show the dynamics when this effective treatment exists.
Assume that insect repellent could decrease the insect biting rate β by 50% of its current level. Figure 10a and 10b show the dynamics when such a repellent is available and appropriately used. Figure 10a shows the case when the repellent is used in the absence of medical treatment (ϕ = 0) and Fig. 10b shows the same scenario but in the presence of medical treatment (ϕ = 0.3) administered at the same time as the provided insect repellent. Compared to the case when the repellent is used alone (Fig. 10a), coupling treatment with the use of a repellent spray significantly reduces the disease outbreak (Fig. 10b) and at Other intervention strategies such as using insecticide treated bed nets, indoor residual spraying (which both can be built into the parameters ϕ and β) and shortening the mosquitoes life span (increasing μ v ) can also be considered. Consider a situation where the life span of the mosquitoes is reduced by 25% of the existing level with both treatment and effective repellent. The population dynamics of humans after reducing the mosquito life span is the same as the dynamics in Fig. 10a and 10b. However, the dynamics of the mosquito population will differ significantly for the two cases. If we decrease both values of the mosquito biting rate and the life expectancy, then there is a 50% reduction of the reproduction number (decrease from R 0 = 0.105 to R 0 = 0.052). For these parameter values, the eradication of the disease is guaranteed. In the absence of medical treatment, the basic reproduction number is slightly greater (R 0 = 0.072) but still less than unity. This implies that if the strategy of shortening Fig. 9 Effects of treatment. The dynamics of the a human, b vector sub-populations when there is medical treatment (prevention chemotheraphy) and the treatment reduces ϑ by 50% the mosquito lifespan before a certain period of time has elapsed is applied, then the lymphatic filariasis disease is potentially bound to die out.

Case 2.
In the previous section we assumed that a disease-free population is invaded by acute-infected humans. We now consider the case where a virgin population is invaded by 70 exposed and 30 infected mosquitoes from an endemic area. Figure 11a and 11b show the dynamics of both the human and vector populations after the invasion of the infected mosquitoes with no medication or quarantine. Employing the same intervention strategies as before, treating acute-infected and isolating some, Fig. 12a and 12b depict the human and mosquito populations dynamics. Any of these strategies significantly reduces the basic reproduction number.

Discussion and conclusions
A mathematical model of the transmission dynamics of lymphatic filariasis incorporating both the human and Fig. 10 Other intervention strategies. The dynamics of the human sub-population when the a repellent is used in the absence of treatment, b repellent is used with treatment mosquito vector was formulated and stability of equilibria and sensitivity analysis were investigated. Numerical simulations were provided to support the theoretical results. Control of infections was analyzed through two intervention strategies, namely medical treatment (prevention chemotherapy) and quarantine (selective treatment for morbidity control). The model system was globally stable and thus the phenomenon of backward bifurcation was never observed [26].
By evaluating the sensitivity indices of the reproduction numbers, we were able to identify parameters for which the model system was most sensitive. We found that the mosquito death rate was the most sensitive parameter. By also analyzing the basic reproduction number, it was shown that combined intervention strategies could lead to lymphatic filariasis elimination in the community.
The proposed model is not exhaustive and can be refined and/or extended in various ways. For instance, the emergence of drug-resistant strains of pathogens is an increasing threat to eradication of infectious diseases.  Aggressive treatment might lead to drug resistance and it is worth exploring how this could affect the transmission dynamics of the disease. Also, patients' compliance could be incorporated into the model system by assuming that only a small portion of individuals in the treatment class adhere to complete treatment, while a small proportion that do not adhere move quickly to the drug resistant class. Model extension could also address climate change since it is considered as a contributor to re-emergence of vector-borne diseases. Heavy rains and global temperature rising provide a conducive habitat for mosquitoes. Future studies could include these external factors and also consider co-infections of individuals with two types of worms.
For mathematical tractability we made several assumptions. Therefore our results are based on the formulation of the model. However, However the research undertaken enables us to gain valuable insights into lymphatic filariasis and the effectiveness of intervention strategies being implemented.