Lymphatic Filariasis Transmission and Control: A Mathematical Modelling Approach

Lymphatic filariasis has an effect on almost 120 million individuals all over the world. The disease may cause a chronic morbidity if the persons who are infected are left untreated. It is endemic in many parts of tropical countries. To prevent worldwide parasite transmission, the World Health Organization initiated the Global Programme to Eliminate Lymphatic Filariasis (GPELF) by eliminating filarial parasites from their human hosts (Molyneux & Zagaria, 2002). Various GPELF implementations are done in many participating countries. In 2004 alone there were more than thirty countries have started elimination program and this number is still rising. Various degrees of success have emerged as a result of the implementation of this program. Although it was reported that in some places the program has interrupted the transmission, in many other places the program could not stop the transmission of the disease (WHO, 2005). It has been argued that strategic choices and operational or biological factors contribute to the success or failure of the program. In general, it is difficult to evaluate the success or the failure of a health program, especially in the beginning of the program. A mathematical model provides useful tools for planning and evaluation of control program in disease elimination (Goodman, 1994). In our earlier work (Supriatna et al., 2009) we develop a mathematical model for the transmission of Lymphatic Filariasis disease in Jati Sampurna, Indonesia. In Indonesia, the disease is already alarming. For example, the incidence of filariasis in Jati Sampurna (a district in the West Java province) is more than 1%. Within less than five years since the date of the publication confirming that Jati Sampurna is an endemic area, almost all regions nearby Jati Sampurna, and other relatively far distance areas are affected by the disease, and some of them are also categorized as endemic areas. Other cases of filarial prevalence are reported outside Java island, such as in Alor islands (the province of Nusa Tenggara Timur). On Alor islands, both B. timori and W. bancrofti are circulated, with a prevalence of up to 20% (Supali et al., 2002). Indonesia joined the GPELF since 2001 and implemented administration of a single dose regimen of diethylcarbamazine (DEC) and albendazole in endemic areas (Krentel et al., 2006). Our previous model tries to capture the effectiveness of this scenario in the attempt of controlling the spread of the disease, inspired by the transmission of the disease in Jati Sampurna. The model assumes that acute infected humans are infectious and treatment is given to a certain number of acute infected humans found from screening process. The screening is


Introduction
Lymphatic filariasis has an effect on almost 120 million individuals all over the world.The disease may cause a chronic morbidity if the persons who are infected are left untreated.It is endemic in many parts of tropical countries.To prevent worldwide parasite transmission, the World Health Organization initiated the Global Programme to Eliminate Lymphatic Filariasis (GPELF) by eliminating filarial parasites from their human hosts (Molyneux & Zagaria, 2002).Various GPELF implementations are done in many participating countries.In 2004 alone there were more than thirty countries have started elimination program and this number is still rising.Various degrees of success have emerged as a result of the implementation of this program.Although it was reported that in some places the program has interrupted the transmission, in many other places the program could not stop the transmission of the disease (WHO, 2005).It has been argued that strategic choices and operational or biological factors contribute to the success or failure of the program.In general, it is difficult to evaluate the success or the failure of a health program, especially in the beginning of the program.A mathematical model provides useful tools for planning and evaluation of control program in disease elimination (Goodman, 1994).In our earlier work (Supriatna et al., 2009) we develop a mathematical model for the transmission of Lymphatic Filariasis disease in Jati Sampurna, Indonesia.In Indonesia, the disease is already alarming.For example, the incidence of filariasis in Jati Sampurna (a district in the West Java province) is more than 1%.Within less than five years since the date of the publication confirming that Jati Sampurna is an endemic area, almost all regions nearby Jati Sampurna, and other relatively far distance areas are affected by the disease, and some of them are also categorized as endemic areas.Other cases of filarial prevalence are reported outside Java island, such as in Alor islands (the province of Nusa Tenggara Timur).On Alor islands, both B. timori and W. bancrofti are circulated, with a prevalence of up to 20% (Supali et al., 2002).Indonesia joined the GPELF since 2001 and implemented administration of a single dose regimen of diethylcarbamazine (DEC) and albendazole in endemic areas (Krentel et al., 2006).Our previous model tries to capture the effectiveness of this scenario in the attempt of controlling the spread of the disease, inspired by the transmission of the disease in Jati Sampurna.The model assumes that acute infected humans are infectious and treatment is given to a certain number of acute infected humans found from screening process.The screening is done every time a new chronic reported.The treated acute individuals are assumed to be remains susceptible to the disease.The model is analyzed and it is found a condition for the existence and stability of the endemic equilibrium.A well known rule of thumb in epidemiological model, that is, the endemic equilibrium exists and stable if the basic reproduction number is greater than one, is established.Moreover, it is also shown that if the level of screening is sufficiently large, current medical treatment strategy will be able to reduce the long-term level of incidences.However, in practice it is not realistic for the following reasons.One important concept in mathematical epidemiology regarding transmission of a disease is the basic reproduction number.It measures the number of new infections caused by an infective during the life time of the infective.Although our previous model is able to gain some insights on how the provision of a medical treatment can reduce the level of disease incidence, however it is worth to note that the basic reproduction number does not depend on the level of the treatment.It means that the treatment, no matter how large it is, will not be able to annihilate the endemicity of the disease.This is some what surprising and unexpected, because normally, in many epidemiological models, any medical treatment should reduce the basic reproduction number.Our earlier work shows that the medical treatment given in the model scenario cannot eliminate the disease, in terms of reducing the basic reproduction number.Our previous model has also ignored an important factor in the transmission stage, namely the time delay.The model has assumed that once an individual infected, he/she become infectious without any delay.Nonetheless, the reproduction number can be reduced by giving additional treatments, such as reducing the biting rate and mosquito's density.This suggests that there should be a combination of treatment to eliminate the disease.In this chapter we review our earlier model of the filariasis transmission and a new model based on the earlier work is developed and analysed.The chapter gives a step by step improvement of our previous model.We do not carry out a heavy mathematical analysis instead some simulations of the models are presented.Finally, some interpretations are derived from the results.

Mathematical model with no time delay in infection period
To formulate the model we use the assumptions that initially the human population is virgin, i.e. there is no infection, and the total population of human is constant.We assume that there is an invasion by few infective individuals of either human or mosquitoes.There is only one species of worm and one species of mosquito, and there is no vertical transmission of the disease, either in human or mosquitoes populations.The human population is divided into three subpopulations, susceptible H S , infected-carrier A and infected-chronic K , with the total number of the population given by H N .We assume that once a human individual is infected then without any delay the individual becomes infectious.However, we strictly assume that transmission to the mosquitoes is only from the acute population.All chronic individuals are isolated perfectly.This strict assumption will be relaxed in some simulation later on. and successful rate of filarial transmission from human to susceptible mosquitoes V p .If the medical treatment is quantified by n number of people screened by the health authority, for every single chronic found, with the successful probability of the treatment 0 p , then the governing differential equations describing the mathematical model of the disease transmission are given by the following equations: We can evaluate the effectiveness of the medical treatment n in managing the disease within the presumed policy, by inspecting its appearance in the endemic equilibrium and in the basic reproduction number.From the model, by assuming the host and vector populations are constant, so that , we found the endemic and non-endemic equilibria of the model related to the basic reproduction number 2 0 () We also establish a theorem saying that "if 0 1 R  then the endemic equilibrium of the system is locally asymptotically stable, otherwise it is unstable".The details of the derivation can be seen in Supriatna et al. (2009).In terms of controlling the disease it means that we should keep the basic reproduction number as low as possible so that it is lower than the unity by adjusting the level of the treatment n.The basic reproduction number is obtained using the next generation matrix (see Diekmann & Heesterbeek, 2000).It is worth to note that the basic reproduction number does not depend on the level of screening n, and hence, current presumed method of treatment does not annihilate the endemicity of the disease.This is partially because of the re-susceptibility of the treated population.However, our earlier work show that it indeed reduces the number of the acute population in the longterm as shown in the following section.

Numerical examples for the model with no delay time in infection period
To facilitate some interpretation regarding the results in our previous work, we present numerical examples using the parameters shown in Table 1.The simulation uses Powersim Constructor Ver.2.5d with the program listing equivalent to basic model of equations ( 1) to ( 5) is provided in the Appendices.Powersim code for other models in the preceding section can be easily modified from this basic model.We give two examples: the first example assumes that a virgin population is invaded by acute infected human (via human immigration) and the second example assume that a virgin population is invaded by infected mosquitoes (e.g. a container un-intentionally transporting infected mosquitoes from an endemic area).
Table 1.The main values of parameters used in the numerical examples Figure 1 depicts the following scenario.Suppose that a population is initially virgin and stays at its equilibrium.We assume that it is then invaded by 10 acute infected human individual, with all the mosquitoes are also virgin.Using the parameter values given in Table 1, we obtain the value of the basic reproduction number is 3.02, which means that the disease will increase if there is no intervention.Figure 1 shows the dynamics when there is no treatment (n=0).The effect of the values of the parameters on the basic reproduction number is clear from equation ( 6).However its effect on the dynamics and the endemic equilibrium is not so obvious.Figure 2 shows the same dynamics as in Figure 1, with an addition that in the 25 th year after the invasion of infective individuals there is a medical treatment with n=200.Figure 3 shows the same dynamics as in Figure 1, but here the treatment is carried out as early as the 5 th year after the invasion with only 100 screening (n=100).These figures reveal that an early average treatment is better than a late huge treatment.
The scenario in Figures 1 to 3 assumes that the medical treatment given to the infected persons does not affect the transmission parameters given in Table 1 other than the screening parameter n.The screening parameter n does not appear in the basic reproduction number formula (6).Hence, this treatment does not affect the endemic status of the disease.
In reality, there are some treatments that could alter the values of the disease transmission parameters.For example, if we assume that some portion of the population is treated by giving them some insect repellent, then the biting rate b could be altered.Let us assume that an effective insect repellent could decrease the biting rate to 50% of its current level.Figure 4 shows the dynamic when there exist this effective insect repellent, and used from the 5 th year in the absence of the medical treatment (n=0) and Figure 5 shows the same scenario as in the previous figure but in the presence of the medical treatment with n=100 given by the same time as the insect repellent provision.Compared to the case when there is no insect repellent (Figure 1), the introduction of the insect repellent is significantly reduces the level of the disease outbreak (Figure 4) and in the same time reduces the endemic level of the disease (changing the value of the basic reproduction number from 3.02 to 1.51).Meanwhile, if we also apply the medical treatment with only average treatment (n=100), then the level of the outbreak is relatively the same, but apparently with a shorter period of the outbreak (Figure 5).   .The dynamics of infected population when there is an effective insect repellent which changes the biting rate to its 50% of the current level with no medical treatment in the 5 th year after the disease invasion (n=0).

Infected Population Sizes
Other scenarios could also be considered.Some are already known to be ineffective if only applied solely, such as fogging (Soewono & Supriatna, 2002) and other still unexplored, such as newly developed method for shortening mosquitoes life expectancy (Turley et al., 2009).Supposed that with some ways we can reduce the mosquito life expectancy down to 50 % of the existing level (from 30 days as in Table 1 to 15 days).Figure 6 shows its .The dynamics of infected population when there is an effective insect repellent which changes the biting rate to its 50% of the existing level combining with average medical treatment in the 5 th year after the disease invasion (n=100).

Infected Population Sizes
dynamics which is the same as the dynamics in Figure 4.This is not surprising considering the form of the basic reproduction number (equation ( 6)), in which the decrease of biting rate acts the same as the decrease of the mosquitoes life expectancy (equivalently the increase of the mosquitoes mortality rate V  ).If we decrease both values, i.e. the values of the biting rate and the life expectancy, then their effect in reducing the basic reproduction number doubled, such as shown by Figure 7, resulting in the value of the basic reproduction number to be less than one (only 0.755), which means the disappearance of the disease is guaranteed.Even in the absence of medical treatment, Figure 8 shows that if we do this strategy before one year has elapsed then the disease does not have any chance to grow.This suggests that preventive action is better than curative action.
In the previous example we assume that invasion is done by infected human.Next in the f o l l o w i n g e x a m p l e w e a s s u m e t h a t i n v a s i o n i s d o n e b y i n f e c t e d m o s q u i t o e s f r o m a n endemic area.Considering the short distance of the mosquito flight, we can assume that this invasion happens un-deliberately, for example via container and other transportation modes.However, considering the stability theorem of the endemic equilibrium point in our previous work (Supriatna et al., 2009), we expect that the long term behaviour of the disease transmission dynamics would be the same as in the first example.In other words, there is an independence of initial values, such as illustrated by Figure 9, in which we assume that there are 100 infected mosquitoes invades the virgin population as described in the first example (Figure 1).Fig. 6.The dynamics of infected population when there is an intervention which changes the mosquitoes life expectancy to its 50% of the current level with no medical treatment in the 5 th year after the disease invasion (n=0). .The dynamics of infected population when there is an intervention which changes both the mosquitoes life expectancy and the biting rate to their 50% level with no medical treatment in the 5 th year after the disease invasion (n=0). .The dynamics of infected population when there is an intervention which changes both the mosquitoes life expectancy and the biting rate to their 50% level done before one year after the disease invasion has elapsed (n=0).

Mathematical model with delay time in infection period
In this section, we take a delay time into account and re-analysed the resulting model by introducing a new compartment to mimic the presence of the delay time.This is done by adding a sub acute or minor acute compartment m A into the previous system.The sub acute population has a lower force of infection than the acute population considering their worm burden status, and it might be not infectious yet.This is reflected by a lower successful rate of filarial transmission 2 V p from the sub acute to susceptible mosquito population, compared to the successful rate of filarial transmission from the acute population 1 V p .In this case, we can consider the sub acute compartment consist of exposed or latent individuals.Individuals stay in sub acute compartment with the sojourn time 1/ before they leave to the acute compartment.The system of equations takes form as the following, 0 ()

Numerical examples for the model with delay time in infection period
As in the previous section, we provide a simulation for the model of equations ( 7) to (11) to gain some insights.The parameters are the same as before unless it is stated explicitly.Compared to Figure 1, in which there are 10 acute infected human initially, Figure 10 shows that the present of time delay, by assuming that the sojourn time in the sub acute compartment is 5 years (hence  is 1/5) with the probability of transmission to the mosquitoes is only 10% of the probability of the acute compartment (hence 2 V p is 0.01), has an effect on significantly delaying the accumulation of the chronic and reducing the number of acute human population.However, the total infectious ( A + m A ) in Figure 10 is slightly greater than the total infectious ( A ) in Figure 1.We can also simulate if in fact we were unable to perfectly isolate the chronic individuals, hence there is a transmission from a portion of them to the mosquitoes.We would expect the transmission rate from the chronic is far greater than the one from the acute population, say the transmission is more certain considering the worm burden carried by them.One of .The dynamics of infected population when there is no medical treatment after the invasion of human comprising of 10 acute individuals.Here we assume that there is no sub acute individual, initially.Fig. 11.The dynamics of infected population as with all parameters as in Figure 10, with an addition that infection also occurs from the chronic by assuming there is no perfect isolation.Fig. 13.The dynamics of chronic population when there is an early medical treatment with n=100 together with the various reduction of the biting rate up to a certain level.

Numbers of Chronic
the realisations is shown in Figure 11.The figure reveals that the peak of the outbreak is higher and reached earlier compared to that in Figure 10.Note that in the early years, there is an iceberg phenomenon, in which the number of chronic is far less than the number of acute.This indicates that early treatment is better than late treatment.Suppose that we administer a medical treatment as in the previous section, measured by the number of screening n. Figure 12 shows various regimes of treatment done continuously since the beginning of the course of the epidemic.Figure 13 shows that a low level of medical treatment combined with the high reduction of biting rate (e.g. up to 15% of the original biting rate) performs better than that resulting from high level of medical treatment with no reduction in biting rate.Suppose now that we have another scenario of treatment, that is giving prophylaxis to all the populations (set n=0, simply to evaluate the effectiveness of this prophylaxis).The prophylaxis works by inhibiting the growth of the worms inside human, say by delaying the recruitment into the acute population A from the sub acute population m A .Technically this is done by varying the values of the transition rate  (or equivalently the sub acute sojourn time 1/ ) in the model.Figure 14 shows the effect of delay for various sojourn time due to the effect of the prophylaxis application.It seems that all the graphs increase exponentially (upper figure), but in fact at the end they end up to their stable equilibrium (lower figure) with different speed and different peak.This indicates that controlling the density of worm inside the body of infective human is effective in reducing the number of filarial infection.

Numbers of Chronic
The model assumes that the delivery of prophylaxis has a result in a constant effect over time, which doesn't reflect the reality.To increase the realism, we should consider the decrease of prophylaxis effectiveness by modifying or refining the model.Nevertheless, we still can apply the current model by only believing the short-term prediction given by the model, say only in one to two years prediction and use it as guidance in a periodic delivery of a mass drug administration program.
The introduction of a single exposed compartment is not without a problem.Getz and Lloyd-Smith (2006) showed that a single exposed compartment will produce an exponentially-distributed sojourn time in the exposed stage.Referring to our delay model (equations ( 8) and ( 9)), this distribution has mean at 1/ while its modus is at 0 , which is a poor match to the real distribution of latent periods.Plant and Wilson (1986) pointed out that the drawback can be resolved by introducing a distributed delay or staging delay time approach comprising of k classes of sub acute or exposed individual.This approach gives a gamma-distributed total time of individuals staying in the exposed class with mean 1 / and variance 2 1/( ) k .Note that a fixed time delay 1 / is obtained whenever the number of delay stages k approaches the infinity.In this part we use this approach (see also Getz and Lloyd-Smith (2006)) to our delay model by introducing multiple exposed compartments which is more appropriate to the disease like filariasis which has more than one different exposed stages.The general model is the same as equations ( 7) to (11) except that equations ( 8) and ( 9) are replaced by The system is much more complex since it consists of 15 differential equations compared to just 6 differential equations in the previous model.However, numerical example in Figure 15 shows that for 10 k  (and also for any 1 k  ), the simpler model of equations ( 7) to (11), qualitatively, is a good approximation of the more realistic model of the same equations but with equations ( 8) and ( 9) are replaced by equations ( 12) to ( 14).Initially, the prediction of simpler model ( 1K in Figure 15) slightly overestimates, but then after a certain years it begins to underestimate, the "true" numbers of chronic individuals ( 10 K in Figure 15).However in the long-term both model produce the same equilibrium point (not shown in the Figure ).7) to ( 11) and the more realistic staging delay time model of the same equations but with equations ( 8) and ( 9) are replaced by equations ( 12) to (14).

Conclusion
In this chapter we review a mathematical model of filarial transmission in human and in mosquitoes.Some simulations are carried out to obtain some insights regarding the transmission and possible actions to control the transmission.Some refinement of the model could be done in many directions to increase the realism of the model and to obtain a more accurate prediction.New directions may include the evolutionary, sosio-economics, and climatology aspects of the disease (Levin, 2002).
In the evolutionary issues of epidemiology, some agents of diseases may develop resistance to certain drug.It is worth to explore how this affects the transmission of the diseases.In many situations, especially in developing countries, there always competing interests related to limited resources and budget.There are many other important diseases, other than filariasis, needs for attention.Choosing the right priorities are among the concerns of health managers and authorities.In the absence of sufficient health budget it is important to address questions like the long term consequences when the treatment is terminated, either purportedly, e.g. because the budget is re-allocated to a higher priority health problem (to other endemic places of the same disease or to other disease problems) or inadvertently (due to the decreasing compliance of the program implementation).This is an example of sosio-economics issues in epidemiology (Supali et al., in prep.).Climate change also regarded as a factor contributes to current emerging and re-emerging infectious diseases.For example, since the global temperature is rising then suitable habitat for mosquitoes becomes wider.It is reported that many parts in the globe of previously free from mosquito is now invaded by incoming mosquitoes.To obtain a better prediction of global filarial transmission, this climatology aspect also should be considered.We believe that there are many other venues are possible for future research in mathematical aspect of filariasis transmission.

Acknowledgment
Part of the research is funded by the Indonesian Government through the scheme of Penelitian Hibah Kompetensi 2012.

Fig. 1 .
Fig.1.The dynamics of infected population when there is no medical treatment after the invasion of 10 infected human.

Fig. 2 .Fig. 3 .
Fig.2.The dynamics of infected population when there is a medical treatment in the 25 th year with n=200.
Fig.4.The dynamics of infected population when there is an effective insect repellent which changes the biting rate to its 50% of the current level with no medical treatment in the 5 th year after the disease invasion (n=0).
Fig.5.The dynamics of infected population when there is an effective insect repellent which changes the biting rate to its 50% of the existing level combining with average medical treatment in the 5 th year after the disease invasion (n=100).
Fig. 7.The dynamics of infected population when there is an intervention which changes both the mosquitoes life expectancy and the biting rate to their 50% level with no medical treatment in the 5 th year after the disease invasion (n=0).
Fig.8.The dynamics of infected population when there is an intervention which changes both the mosquitoes life expectancy and the biting rate to their 50% level done before one year after the disease invasion has elapsed (n=0).

Fig. 9 .
Fig. 9.The dynamics of infected population when there is no medical treatment after the invasion of 100 infected mosquitoes into a totally virgin population.The figure is similar to Figure 1 in which invasion is done by 10 infected human.
Fig.10.The dynamics of infected population when there is no medical treatment after the invasion of human comprising of 10 acute individuals.Here we assume that there is no sub acute individual, initially.

Fig. 12 .
Fig. 12.The dynamics of chronic population when there is an early medical treatment with various values of n.

Fig. 14 .
Fig.14.The dynamics of chronic population when prophylaxis is given to the whole population with various effects to the sub acute sojourn time (equivalent to the reciprocal value of the recruitment rate from sub acute to acute population).

Fig. 15 .
Fig. 15.The dynamics of chronic population predicted by the simple model of equations (7) to (11) and the more realistic staging delay time model of the same equations but with equations (8) and (9) are replaced by equations (12) to (14).