The Use of Non-Standard Finite Difference Schemes to Solve the DAMP and SIT Models

Joshua A Mwasunda, Eunice W Mureithi, Nyimvua Shaban

Journal of Mathematical Sciences and Applications

The Use of Non-Standard Finite Difference Schemes to Solve the DAMP and SIT Models

Joshua A Mwasunda1, Eunice W Mureithi2, Nyimvua Shaban2,

1Department of Mathematics, Mkwawa University College of Education, Tanzania

2Department of Mathematics, University of Dar es Salaam, Tanzania

Abstract

Sterile insect technique (SIT) is a method of biological control that uses sterile male insects to reduce the reproductive rate of a species of target insect. The method relies on the release of sterile or treated males in order to reduce the native population of insects. We propose the model that governs the dynamics of the anopheles mosquito population, and then modify to incorporate the sterile insect technique as an intervention to curtail the reproduction of mosquitoes. The nonstandard finite difference numerical schemes and simulations for these models are provided. The results indicate that sterile technique with frequent and high rate of release can be an alternative to chemical control tools in the fight against malaria.

Cite this article:

  • Joshua A Mwasunda, Eunice W Mureithi, Nyimvua Shaban. The Use of Non-Standard Finite Difference Schemes to Solve the DAMP and SIT Models. Journal of Mathematical Sciences and Applications. Vol. 3, No. 2, 2015, pp 25-32. http://pubs.sciepub.com/jmsa/3/2/2
  • Mwasunda, Joshua A, Eunice W Mureithi, and Nyimvua Shaban. "The Use of Non-Standard Finite Difference Schemes to Solve the DAMP and SIT Models." Journal of Mathematical Sciences and Applications 3.2 (2015): 25-32.
  • Mwasunda, J. A. , Mureithi, E. W. , & Shaban, N. (2015). The Use of Non-Standard Finite Difference Schemes to Solve the DAMP and SIT Models. Journal of Mathematical Sciences and Applications, 3(2), 25-32.
  • Mwasunda, Joshua A, Eunice W Mureithi, and Nyimvua Shaban. "The Use of Non-Standard Finite Difference Schemes to Solve the DAMP and SIT Models." Journal of Mathematical Sciences and Applications 3, no. 2 (2015): 25-32.

Import into BibTeX Import into EndNote Import into RefMan Import into RefWorks

At a glance: Figures

1. Introduction

Malaria is a vector-borne disease caused by the Plasmodium parasite and transmitted between humans through the bite of the female Anopheles mosquito. Four species of protozoan parasite of the plasmodium genus - P. falciparum, P. vivax, P. ovale, and P. malariae - cause malaria in humans. Though malaria caused by P.vivax is the most common, but it is, however, malaria caused by P. falciparum that is most dangerous and sufficient to cause death [1]. Most people who die from malaria are African children below the age of 5 years and pregnant women. The dreaded disease is recognised as difficult to eradicate and its control is possible only with coordinated efforts of the general public, healthcare personnel and government and international agencies [2].

Chemicals have been and are still extensively used all over the world to control wild insect populations. However, there are significant challenges, including the increase insect resistance to insecticides, environmental contamination, effects on non-target organisms and the selection of resistance hampering its effectiveness [3]. As viable alternative, non-polluting methods also known as biological control tools are more explored, with a special focus on the ecology and behaviour of the involved species. One of the most promising methods is the sterile insect technology (SIT). For example, we refer the reader to [4] and [5] for details.

This technology (SIT) can also be used to control anopheles mosquito, the vector responsible for malaria transmission [6]. In this case, male mosquitoes are bred and then exposed to enough gamma radiation to render them sterile. The sterile males are then released repeatedly into the environment in large numbers in order to mate with the native female anopheles mosquitoes preventing production of offspring (e.g. see [7, 8]).

This paper contributes to this work by modifying the models by [6, 9] and [10] to include the release of sterile male mosquitoes. Non-standard finite difference schemes are developed to simulate the models. The use of these schemes is motivated by the fact that standard schemes may give spurious and unstable numerical solutions which depend strongly on time step-size. Also, standard schemes may give results which do not preserve positivity property. This paper is organised as follows: In Sections 2 and 3 DAMP and SIT models are formulated and discussed. Non-standard finite difference schemes for DAMP and SIT are developed in Sction 6. Numerical simulations are implemented in Section 5 and Section 6 gives concluding remarks.

2. DAMP Model Formulation

The model divides population at any time t is into four compartments namely; immature population, I; the non-laying female mosquitoes, Y; fertilized egg laying mosquitoes, F and male mosquitoes, M.

We consider a female mosquito to be in the Y compartment starting from its emergence from pupae until when her gonotrophic cycle has begun, (that is the time of mating and taking the first blood meal) which takes typically 3-4 days.

A female mosquito needs to mate successfully only once and get a blood meal before it starts laying eggs during the gonotrophic cycle. The compartments for the dynamics of anopheles mosquitoes are shown in Figure 1.

Figure 1. A Compartmental Model for Wild Mosquito Population

The population in immature stage grows to adult mosquitoes at a per capita rate; a fraction f of such emerging population represents females and the remaining fraction , males. The per capita mortality rates for the immature stage, the young females not yet laying eggs, the mating fertilized females, and the wild male mosquitoes are denoted by, , and respectively.

The net oviposition rate per female mosquito is proportional to their density, though it is also regulated by the carrying capacity effect depending on the occupation of the available breeding sites. In this model the per capita oviposition rate is assumed to be , where K is the carrying capacity related to the amount of available nutrients and space, and is the average amount of eggs laid per fertilized female per day. Thus, the rate (per day) of laying eggs in the breeding sites is .

The male mosquito can mate practically throughout all its life. Since the female needs one successful mating there is an overabundance of males. Therefore, in general, it is reasonable to assume that the waiting time for mating does not depend on the number of males M in the sense that if M is increased further this rate remains the same. For the model this means that the transfer rate from compartment Y to compartment F is independent of M.

The description above leads to the following system of ordinary differential equations for the DAMP:

(1)

The DAMP model (1) has two equilibrium points, the trivial equilibriumwhere the wild mosquitoes are absent and the non-trivial

withwhere the wild mosquitoes survive and R is the basic offspring number defined by .

The trivial equilibrium point is locally asymptotically stable (LAS) if and unstable otherwise while the non-trivial equilibrium point is locally asymptotically stable whenever and unstable otherwise.

3. The SIT Model Formulation

In this section we extend the basic model for the dynamics of anopheles mosquito population by introducing compartments of treated males and of females that would be laying sterile (not hatching). This formulation yields a SIT model, in which we study the effect of sterile male release for control of the wild mosquito population. The death rate, of treated male mosquitoes depends on the procedure used during sterilization and we denote the death rate of females that would not be hatching eggs by.

Flows from Y to F and from Y to U compartments depend mainly on the number of encounters of females with native and sterile males, and on the corresponding mating rates. Under the assumption that the mosquitoes in the compartments M and are equally likely to mate, a mating female mosquito has probability to be with wild mosquito and probability to be with a sterile mosquito. Hence the transfer rate from the compartment Y splits into transfer rate of to compartment F of fertile females that would be laying eggs that hatch into young mosquitoes and a transfer rate of to compartment U of infertile females that would be laying sterile (not hatching) eggs.

Note that the total mating rate remains unchanged by the introduction of the sterile mosquitoes. This indeed has to be the case since, as mentioned earlier, an increase of the number of males does not change the mating rate.

Similar to the immature stage for natural population of mosquitoes, the growth rate of sterile male mosquitoes is assumed to be regulated by a carrying capacity effect. This feature takes into account the limiting capacity (of laboratories, for instance) to produce and release sterilized male mosquitoes, which allows the release rate of sterile male insects to be described by where C is the maximum capacity related to the sterile insect production. Note that is the per capita release rate where represents the intrinsic release rate.

Based on the model descriptions (see Figure 2) and assumptions, we have the following equations:

(2)
(3)
(4)
(5)
(6)
(7)

Note that equation (6) can be considered as a logistic release of sterile male mosquitoes. Equation (7) for mating unfertilized females can be decoupled from the system (2) – (7) and equation (6) can be easily solved, which has the exact solution given by

where B is a constant of integration which can be found from the initial condition

Equation (6) for the sterile male mosquitoes has two equilibrium points: and. Note that equation (6) is biologically feasible if and only if. The value gives rise to two equilibrium points. One is the trivial equilibrium where both natural and sterile male mosquitoes are absent and the other is

With where only the natural mosquitoes survive. When, the SIT model reduces to the DAMP model and therefore the stability nature of and is the same as and respectively.

For , we have three equilibrium points. The trivial equilibrium with and the non-trial equilibriums

and

with

and the parameter R is the basic offspring number where and are the conditions for biological existence of the non-trivial equilibriums.

The remaining decoupled mating unfertilized females is given by

implying that with being substituted by and to produce the equilibrium points and respectively, and .

The parameter W measures the ratio between the number of sterile males and the number of natural male mosquitoes in equilibrium. If W is sufficiently high, the next generation of wild mosquitoes would be much more lower than the actual one since a considerable proportion of eggs would not hatch. If sterile male mosquitoes are released for a long period of time, this pattern would drive the natural mosquito population to zero.

The stability nature of the non-trivial equilibriums indicates that is always unstable andlocally asymptotically stable.

4. The NSFD for DAMP and SIT Models

A continuous dynamical system given by system of ordinary differential equations can be solved by using standard numerical methods such as Runge–Kutta methods [11]. However, these methods may give spurious solutions and numerical instabilities that depend strongly on the time step-size. They may also give results that do not preserve the positivity property. New procedures known as Non-standard Finite Difference (NSFD) schemes were developed to obtain schemes that are dynamically consistent and their solutions preserve the physical properties of the approximated differential system for arbitrary time step-sizes [12]. Such properties include conservation law, positivity, monotonicity, replication of the fixed points and their stability.

A finite difference scheme is called nonstandard finite difference method if at least one of the following conditions is met:

(i) In the discrete derivative, the traditional denominator is replaced by a nonnegative function such that

(ii) Nonlinear and negative linear terms that occur in the differential equation are approximated in a non-local way. For example, may be approximated by or by, by and by, where is any given parameter and x and y represent a given phenomenon (see [12, 14] and [15]).

In this section, we develop NSFD numerical schemes for DAMP and SIT systems (1) and (2) - (7). The new NSFD methods preserve both the positivity of the solutions and the stability of the equilibriums of the corresponding DAMP and SIT systems. In addition, the designed numerical approximations allow us to solve the discrete systems explicitly, which increases the efficiency of the methods.

Assume that all the death rates are equal denoted byand let the total number of natural population be given by. The DAMP model becomes:

The discrete model is given as:

Since we do not have an exact scheme we need to solve for the total population V as well. Note that if, we have an exact scheme with

Expressing in Gauss-Seidel-like structure, we have:

(8)

Proposition 1: The NSFD scheme preserves the positivity property.

Proof: Consider the Gauss-Seidel-like system (8) above. Using the initial conditions , , , it can be seen clearly from (10) that and . If we substitute and into equation for results into and substituting and into equation for F gives.

Assuming that,,.Then using the same argument, we have. Again this will result into.Therefore the positivity property holds for all

Proposition 2: The discrete DAMP model has the same equilibrium points as those of the continuous system.

We can use system (8) to find the fixed points of the system by lettingand similarly with Y, F and M.

(9)

Implying that

(10)

Thus, the NSFD scheme reduces to DAMP model after setting the derivatives equals to zero and therefore it will have the same fixed points and as shown in the continuous case (since). The equilibrium solution for V can be obtained by solving.

For the trivial equilibrium point,, which implies that as well. For the non-trivial equilibrium point , we need to solve (for I and F different from zero), implying that .

4.1. The Stability Analysis of the Discrete Model

The stability properties for the discrete model can be studied by analysing the Jacobian matrix of the Gauss-Seidel-like structure and applying the Jury stability criterion [16]. The Jury stability criterion is a method of determining the stability of a linear discrete time system by analysis of the coefficients of its characteristic polynomial. It is the discrete time analogue of the Routh-Hurwitz stability criterion. The Jury stability criterion requires that the system poles are located inside the unit circle centred at the origin, while the Routh-Hurwitz stability criterion requires that the poles are in the left half of the complex plane.

The discrete approximation (NSFD) for SIT model is given by:

Expressing discrete SIT model in Gauss-Seidel-like structure, we have

Summing up the right equations Y, F, M, and U gives the total number of wild flying mosquitoes . Now, we find the denominator function for treated male mosquitoes using the following procedures:

First, the fixed points of equation (6) need to be obtained by solving equation Thus the fixed points are and. Then the roots of equation (6) need to be obtained using the formula .

But .

Thus and

, and the denominator function is defined by . Thus the denominator function for treated male mosquitoes is given by where.

For other equations we use the denominator function

5. Numerical Simulations

In this section, simulations of the discrete (NSFD) DAMP and SIT models are carried out using Matlab software. The main aim is to verify some of the analytical results on the stability of the systems. Beside verification of our analytical findings, these numerical solutions are very important from practical point of view. The basic model for the dynamics of anopheles mosquito population is simulated in the absence of release of treated male mosquitoes and then this model is simulated after release of treated male mosquitoes, and then the effects of varying the release rate () are observed.

The vital parameters are,,,,and f as indicated in the table below. The values of the parameters are kept fixed at K=120 and C=100 where K and C represents the maximum carrying capacity for immature stage and sterile insect production respectively. The figures are plotted using the parameter values in table 5 and the initial conditions (estimated initial average values of the population) are given for each graph. The initial time is taken to be and the step size

Table 1. Parameter Values Used in Simulations

5.1. Simulation without SIT Control

The simulation of DAMP model has been done to find out the dynamics of the anopheles mosquito population before release of treated male mosquitoes.

The results show that, the immature population increases/decreases with time until it reaches its equilibrium level depending on the initial condition as shown by Figure 3(a).

The same trend as for immature population occurs for young female population, fertile females and wild males as indicated by Figures 3(b), 3(c) and 3(d) respectively. This shows that as long as there is no release of treated male mosquitoes to control the population, there will exist a wild population of natural mosquitoes since the average number of secondary female mosquitoes produced by a single female mosquito R > 1 (R = 14.7465). This result supports the theorem on local stability of nontrivial equilibrium point.

5.2. Simulations with SIT Control

When sterile males are released, all the subpopulations are affected. We simulate the SIT discrete model using initial conditions: and. Figure 4(a) shows the number of released sterile male mosquitoes with time. This subpopulation increases with time until when it reaches its equilibrium. It can be observed in Figure 4(b) that immature population initially increases to a maximum and then starts to decrease with increase in time until it goes to extinction as a result of the released treated male mosquitoes.

Figure 4(a). The graph representing exact solution for (Treated males)

Figure 4(c), shows the sum of the young females, fertile females, wild males and unfertile females that are grouped together as wild flying mosquito population. The graph indicates that after the release of treated males, the wild mosquito population initially grows to a maximum and then starts to decrease with time to until it goes to extinction.

Figure 4(c). The graph for (Wild flying population) after release of

The reason for the trends in Figures 4(b) and 4(c) is that, when sterile male mosquitoes are released into the environment they will mate with the native female mosquitoes. A native female that mates with a sterile male will lay eggs, but the eggs will not hatch. Hence the immature population and wild mosquito population will decrease with time eventually reaching extinction.

5.3. Variation of Population for Different Values of Release Rate

The simulation of SIT model when there are changes to the release rate value has been done to find out the dynamics of the population. We simulate this model using initial condition: and two different release rates; and by keeping all other parameter values fixed and observe the nature of resulting graphs. Figure 5 shows the variation of immature population and wild flying mosquito population for different values of release rate (α).

The graph indicates that the control of the wild-type mosquito population is highly dependent on the rate () at which the sterile males are released, with only high release rates giving sufficient control for a short period of time. This is because the sterile males become larger in size as compared to wild males such that the probability of sterile males to mate with fertile female mosquitoes is greater than the probability of wild males to mate with fertile female mosquitoes. In turn, most of the fertile female mosquitoes will lay eggs that cannot hatch to become larvae and pupae.

Figure 5. SIT Model for Different Values of Release Rate

6. Conclusion

In this study, we have presented both the DAMP and SIT models. The main objective of the study was to assess the effect of SIT for the control of anopheles mosquito using NSFD schemes. Numerical results show that before release of sterile male mosquitoes the wild mosquito population was approaching the non-trivial equilibrium. Soon after release, there was a major change in trend whereby the immature and wild flying population started to decrease to zero. With increase in release of sterile males, it is was shown that we can control the wild mosquito population for a short period of time due to the fact that the probability of sterile males to mate with wild females will be higher compared to that between wild males and females. Therefore the eggs that will be laid by unfertilized female mosquitoes will not hatch. In general the model predicts that when (where is the death rate for treated males) extinction of the wild mosquitoes depends on its initial population size and the rate of release of sterile males, with only high release rate giving sufficient control. Based on the model of this study, it is proposed that future work should consider:

1. Finding better methods of sterilizing and releasing sterile males that do not affect their behavior as well as fitness to mate with wild females.

2. Comparing the SIT approach with standard chemical vector control for anopheles mosquito population.

Acknowledgements

J.A. Mwasunda gratefully acknowledges for the study leave and financial support from MUCE.

References

[1]  World Health Organization (WHO): Media Centre Malaria Factsheet, No.94.October, 2011 http://www.who.int/ mediacentre / factsheets/fs094/en/ [Accessed October 22, 2011].
In article      
 
[2]  World Health Organization (WHO), Media centre – Malaria Factsheet, No.94. April 2010. http: // www.who.int/mediacentre/ factsheets/ fs094/en/[Accessed August 28, 2011].
In article      
 
[3]  Dorta,D.M., Vasuki, V and Rajavel, A, Evaluation of organophosphorus and synthetic pyrethroid insecticides against six vectormosquitoes species. Revistade SaúdePública 1993, 27: 391-7.
In article      PubMed
 
[4]  Knipling, E.F, Sterile insect technique as a screwworm control measure: the concept and its development, in: O.H.Graham (Ed.), Symposium on Eradication of the Screwworm from the United States and Mexico, 62, Misc. Publ.Entomol. Soc. America, College Park, MD, p. 4, 1985.
In article      
 
[5]  Bartlett,A.C and Staten, R.T, The sterile release method and other genetic control strategies, in: E.B. Radcliffe, W.D.Hutchison (Eds.), Radcliffe’s IPM World Textbook, University of Minesota, St. Paul, MN, Available at: http://ipmword.umn.edu, 1996.
In article      
 
[6]  Anguelov, R.,Dumont, Y andLubuma,J, Mathematical modelling of sterile insect technology for control of anopheles mosquito, Computers and Mathematics with Applications.
In article      
 
[7]  Steinau, R, Tips from a real pest Control expert-Sterile Insect Technique: www.asktheexterminator.com/Do_It_ Yourself_Pest_Control/Sterile_Insect_ Technique.shtml, 2007.
In article      
 
[8]  Mwasunda, JA, Modelling the Effect of Sterile Insect Technology for control of Anopheles mosquito population in Tanzania, MSc. Dissertation , University of Dar es Salaam, 2012.
In article      
 
[9]  Esteva, Land Yang, H.M, Mathematical model to assess the control of Aedesaegypti mosquitoes by the sterile insect technique, Mathematical Biosciences, 2005, 198, 132-147.
In article      View Article  PubMed
 
[10]  Esteva, L and Yang, H. M, Control of Dengue Dengue Vector by the Sterile Insect Technique Considering Logistic Recruitment, TEMA Tend. Mat. Apl. Comput., 7, No. 2, 259-268, 2006.
In article      View Article
 
[11]  Dobromir, T. D and Hristo, V. K, Nonstand-ard finite-difference methods for predator–prey models with general functional response, J. Mathematics and Computers in Simulation 2007, 78 1-11.
In article      
 
[12]  Mickens, R.E Advances in the applications of nonstandard finite difference schemes.World Scientific, Singapore, 1994.
In article      
 
[13]  Dumont, Y and Lubuma, J Non-standard finite difference methods for vibro-impact problems, Proc.R. Soc. London, 461A, 2005, 1927-1950.
In article      View Article
 
[14]  Mickens, R..E, Advances in the applications of nonstandard finite difference schemes. World Scientific, Singapore. 1994.
In article      
 
[15]  Mickens, R, Non-standard Finite Difference models of Differential Equation. World Scientific, Singapore 2005.
In article      
 
[16]  Kasim, M. A, Stability of Real-Time Systems, Computer Engineering Department, Philadelphia University, 2011.
In article      
 
  • CiteULikeCiteULike
  • MendeleyMendeley
  • StumbleUponStumbleUpon
  • Add to DeliciousDelicious
  • FacebookFacebook
  • TwitterTwitter
  • LinkedInLinkedIn