The Dynamics of Vector-Host Feeding Contact Rate with Saturation : A Case of Malaria in Western Kenya

In this study, we develop an expression for a saturated mosquito feeding rate in an SIS malaria model to determine its effect on infection and transmission dynamics of malaria in the highlands of Western Kenya. The basic reproduction number is established as a sharp threshold that determines whether the disease dies out or persists in the population. Precisely, if 0  0 1   , the disease-free equilibrium is globally asymptotically stable and the disease always dies out and if , there exists a unique endemic equilibrium which is globally stable and the disease persists. The contribution of the saturated contact rate to the basic reproduction number and the level of the endemic equilibrium are also analyzed. 0 1  


Introduction
Malaria is an infectious disease caused by a parasite of the genus, Plasmodium.It is transmitted between human hosts by female anopheles mosquitoes as they seek blood meal for their eggs development.When a mosquito bites an infected person, a small amount of blood is taken in, which contains microscopic malaria parasites.When the mosquito takes its next blood meal, these parasites mix with the mosquito saliva and are injected into the person being bitten and the transmission process is perpetuated.
Malaria constitutes a big health problem especially within sub-saharan Africa and Asia.World Malaria Report 2011 estimates that it causes between 250 -260 million infections and more than a million deaths (mostly among children in Africa), annually.In Kenya, reports show that despite the many control strategies to eliminate malaria, it has re-emerged and increased in incidence.The disease continues to wreck havoc on millions especially from the poor countries [1].
Vector abundance in Western Kenya is driven by temperature variation, ecosystem characteristics and human activities.The population varies depending on the site, the season and the species of the vector.Some sites in Western Kenya have 12.7 fold indoor resting densities during the long rainy season (March-June) and 23.3 fold during the dry season (January-March) [2].This implies that the vector population is never constant as assumed in many models.On the other hand host population changes due to seasons and economic activities, natural deaths and death due to diseases like malaria and migration to urban centers and other regions for greener pastures.For our model to capture the reality of the epidemics in Western Kenya highlands, we assume that the host and mosquito populations change with time.
Most malaria models assume a constant human biting rate in their models, which means that hosts are freely available whenever a mosquito wants to bite, but in practice, this is more of a simplifying assumption.Research shows that for small host population, this rate is proportional to the host population size, and for large host population, it is constant [3,4].The feeding cycle of a mosquito involves, host-seeking, feeding, resting, site-seeking, oviposition and host seeking resumes [5,6].The probability of finding a host and successfully obtaining a blood meal depends on many factors.Among them is human avoidance and defensive behaviour [7,8].
If the mosquito survives this process, and has a successful blood meal, it rests, finds a larval habitat, oviposits and continues host seeking.Since this process drives malaria transmission, it is necessary to address the specific form of the mosquito-human contact process.Arditi [9], argues that the rates of (successful predation) contact between a predator and a prey is most properly a function of the ratios of their proportions.This would fit into malaria mosquitoes which inhabit homesteads and other areas where human hosts are available, like farms and urban areas [10].It is also clear that this contact rate does not increase without bound, as the predator-prey ratio increases, this is because once a mosquito is fed, it rests before ovipositing, to resume host seeking and biting again [11].
When the predator-prey ratio value is low, the contact rate will be limited by the predators ability to find the prey, on the other hand if the ratio is high, the contact rate is limited by the predators satiation (desired predation rate) [4,9].For malaria, the contact rate takes a similar course where the mosquito bites will increase as a function of host-vector ratio until the ratio reaches a critical level [12].
Saturation models are also not lacking in literature.A cholera model with saturation in the incidence was proposed by Capasso [13].They argue that when there is a real threat to infected people becoming cautious and taking preventive measures which control further infection.Heesterbeek [14] formulated a saturated individual contact rate in relationships such as courting and marriage, where they assume that the population mixes randomly.Zu and Ma [15] analyzed a SEIR epidemic model whose latent period is described by delay and included a saturated incidence rate.Zhang and Ma [16] studied a SEIR model with saturation in contact rates and did a thorough analysis of its global dynamics.
In 2010, Ming and Li [17] formulated vector borne disease model, where they argue that increasing the density of the susceptible hosts with respect to infected ones leads to Holling type II saturation on the force of infection of host.In this model the biting rate of vectors and both populations are assumed to be constant.Further the model neglects disease related deaths a very crucial factor in malaria infection.
A model by [18] on dengue with variable human population is formulated and analyzed for both local and global dynamics.Ngwa et al. [19] analyzed the stability of a malaria model, with disease deaths, recovery and variable host and vector populations.However they assumed that the biting rate of vectors is constant hence their infection term is the one described in [20].Realising then the need to predict the dynamics and transmission of malaria with great precision, we are motivated to engage in this study, as we pay particular attention to saturation in mosquito feeding habits and the varying host and vector populations.
The rest of the paper is subdivided as follows.Section 2 covers vector-host contact with saturation.In Section 3, the saturated contact process model is formulated.Section 4 is dedicated to the existence of equilibiria, while Section 5 studies the stability of the Disease Free and the Endemic equilibrium.Finally in Section 6 we give some results on numerical simulation.

Vector-Host Contact with Saturation
Vector-host contact results from the need for mosquitoes to obtain a blood meal for their eggs development.A given vector's biting rate is limited by both host population density and its own feeding frequency [21].
Therefore the per vector biting rate should increase as a function of the host-vector ratio until the ratio reaches a critical threshold, which we denote v , above which, biting rate saturates and the average vector can feed at its preferred rate v (contacts per vector per time).Below this threshold, the relative scarcity of hosts constrains the rate at which a vector can feed on the given type of hosts (it must seek other sources).

Q b
We assume that an average host can receive bites at a maximum rate h beyond which it successfully defends itself against the vector (including leaving the place altogether) [12].Then this threshold density ratio is given by b There are two ways of modeling the biting rate which increases for small host population and then approaches a maximum for large populations.The first is using a smooth verhulst-type function of the form here v is the preferred biting rate, and A, the population density ratio at which there is saturation, and it measures how soon the saturation occurs.Dietz [22] and Ming [17] have used an equivalent form of saturation [23].

b 50%
The other way of modeling the vector biting rate saturation is by using a continuous and piecewise function with a "switch point", v deleanating the boundary between the two ranges as used by kribs in [23,24].This function takes the form where v is the maximum preferred biting rate and v is the threshold density ratio.There are other Hollingtype responses suggested in [14].The difference between the two saturation models is the saturation sharpness, that is, how quickly the per vector biting rate levels off as the hosts become plentiful, that a vector can feed at its preferred rate, v .The saturation is gradual for the func- The sharpness is crucial for our malaria model because mosquito populations can change drastically in a short span of time [23].We also choose the predator-prey type or response since the biting in malaria transmission is vector initiated.The fact that the second model is also easy and captures more dynamics adds to the reason we apply it here in this malaria model.
From the reasons above, we shall assume a saturated contact process as used in [23,24] with the so-called Holling Type 1 form.Under this assumption, the per-vector contact rate can be described as a function of the hostvector density ratio v (many host per vector) the rate completely saturates at the maximum desired biting rate f z  b , while for (i.e.few hosts per vector) , and the rate rises linearly with the hostvector ratio.The later is our interest in this study since the saturated contact process   v f z b  has been used in the classical Ross model for malaria [20].Substituting which can be rewritten in the form   min , .
We note that for the current host density the maximum number of vectors that can effectively bite hosts at one given time is v v .Therefore the parameter v is an important determinant in our model.It determines which of the two population densities is driving the biting contact rate.

Q N Q
To examine the rate of appearance of new malaria infections from the rate of mosquito feeding contacts, we have to take into account the probability of infection resulting from an effective contact where one party (host or vector) is infected with malaria parasite and the other is not.Let be the probability that such a contact between an infected vector and an uninfected host results in infecting the host, as the proportion of blood meal contacts between infected hosts and uninfected vectors which result in an infected vector.Further if are the susceptible and infectious hosts respectively and are the susceptible and infectious vectors respectively, the new infections will be given as defined in [12].
• For the Hosts In many malaria models, saturation has been assumed to be constant [20,25,26].Here we consider the density dependent biting rate, and the populations ratio plays a vital role in the transmission.If the vectors-host ratio is low, the bites are few and hence the probability of transmission reduces too.As the ratio increases the infection will rise as a function of the ratio until it reaches the threshold and it becomes a constant.We wish to model this change in biting rate as the vector-host ratio changes and study its effect on the basic reproduction ratio.

The Model Equations
The model we derive here is mathematically equivalent to the classical Ross model [20].The saturation in contact processes will address how the infection rates depend on this ratio.We assume that mosquito has a variable population growth such that birth v v    .For the human population, we assume a density dependent mortality rate, such that the total population vary with time and is modified by a logistic equation that include dis-ease induced deaths.
A description of the variables and parameters used in the model follows in Tables 1 and 2 respectively.
The dynamics of our model will be governed by the following set of equations: The term h in the susceptible hosts compartment corresponds to a constant recruitment of susceptible hosts respectively.Infective period of mosquitoes ends with their death due to their relatively short life-cycle so we do not have recovery or immune term in the vector equations [27,28].
All the parameters in the model are non negative and the model equations are well posed.For initial values

 
, , , , , S I S I N N in , the solutions exist and remains in the region for all .
6   t 0  In the absence of disease the host population dynamics is given by h h h h .In this kind of demographic structure, the total human and mosquito population size For ease of studying the system, we let setting , and the equation now takes the form   , , , . .
which is defined in feasible region (i.e.where the model makes biological sense) , 0 , 0, 0   denotes the non-negative cone of including its lower dimensional faces.It is clear that 6


 is positively invariant with respect to (3).We denote the boundary and the interior of  by and    respectively.

Global Stability of the Disease Free Equilibrium
We recall the equations of the model and use the relation and , to study the system

A Compact Positively Invariant Set
Using Barrier theorems (e.g.[29,30]) we prove that the following set

 
, , , 0 , 0 is a positively invariant compact set for system (4).Moreover K is a global attractor on the nonnegative orthant .
4   Since the ODE is Lipschitz, it is sufficient to check that the vector field induced by the system is either tangent or entering K on the boundary K .
Clearly we have the following implications: 1) and 0 0 2) ; 0 0 we have ; 0 0 6) When and The preceding relations prove that all trajectories tends to K , which ends the proof of our claim.This also implies that all the trajectories are forward bounded.

Reduction of the Model
We remark that in the last equation only the state is appearing.Hence our system is a triangular system.Using Vidyasagar theorem on v N K ([31], Theorem A.4 given in annex) we can reduce the stability study to the stability of the equivalent system This system is considered now on A similar argument, as in the preceding section, shows that  is a global attractor on the nonnegative orthant 3   for system (5), as shown in Figure 2.

Basic Reproduction Ratio
For system (5) there exists a disease free equilibrium (DFE), which is (0,0, in .According to the technique by van den Driessche-Watmough [32], we define the transmission vector and the vector of other transfers defined on the "infected components" The Jacobian of these vectors computed at the DFE gives The spectral radius of I I N denote an equilibria.We have the following system of equations The following relations are satisfied : • From (6a) ; Replacing these values in (6a) gives


We see that if In conclusion, if 0 , we have an unique endemic equilibrium in the interior of .

Global Stability of the DFE
The DFE is on the boundary of .

Theorem 4.1
The disease-free equilibrium of system ( 5) is globally asymptotically stable in the nonnegative orthant if and only if and is unstable if .
the unstability of the DFE is a consequence of the theorem in [32,33].

 
To prove the global stability, we consider a Lasalle-Lyapunov function [34][35][36] on the positively invariant compact set.By this definition we mean that on V  , is continuous and nonnegative.The function V is not positive definite.We define We take advantage on the fact that the ODE (5 ) can be written , then, the derivative along the trajectories is given by For the inequality in the last row we assume on .
Finally, with , since Now the largest invariant set contained in the set E defined by is certainly contained in the set of points for which or .We have two situations: To be in an invariant set with has for consequences that and consequently , hence to be at the DFE; • To be in an invariant set with implies either , thus to be at the DFE, or , which is cannot be contained in an invariant set by the last equation system (5).
We have proved that the largest invariant set contained in and in the positively invariant set 0 V    is the singleton constituted by the DFE.By Lasalle's theorem [35,36] this proves that the DFE is globally asymptotically stable in , hence in the nonnegative orthant. 

Stability of the Endemic Equilibrium
In this section we will prove that, when 0 then the unique endemic equilibrium which is asymptotically sta-ble.
Proof.We consider the Jacobian of system ( 5) We can see that our system is not monotone and cannot be monotone for any orthant of , since we have the entry of J negative and the entry positive.
The Jacobian computed at the EE, using the relation of ( 6), can be expressed as We will examine the characteristic polynomial of

 
, , The Routh-Hurwitz criteria ensures that   , , We have now to compute the sum of the 2 a 2 2 principal minors.

h h h h h h h v v v h h h v v h h h h h v h h v v
The first three requirements of Routh-Hurwitz's criteria are satisfied.We have now to conduct a long calculation to prove 1 2 Let compute This proves the asymptotic stability of the endemic equilibrium.

Global Stability of the Endemic Equilibrium, when There Is No Disease Induced Mortality
We suppose in this section that 0 h   .In this case we can reduce the system to a two dimensional system, thanks to Vidyasagar's theorem.The original system is then, for the stability point of view, equivalent to In this case this system is a Ross system, and the global stability of the endemic equilibrium is well known, see for example [25,35,37].

Estimating v Q
The saturation threshold ratio v is an important parameter that determines which of the two populations, that is the host and the mosquito, is driving malaria in- fection.We will assume that the indoor resting density for mosquitoes is 2.5 for each household that has approximately 6.5 persons (Githeko, personal communicaz   1 .This implies from our model that there are few host for each vector and the relative scarcity of hosts constrains the rate at which an average vector can feed.The biting rate therefore will be driven by the host and not the vector population.For its survival it may be forced to seek other sources for its blood meal.Reduced host biting translates to lower malaria infection and a reduced basic reproduction number.Then 0 0.718 , the infection dies out, and the infectious population goes to zero as shown in Figures 3(a) and (c).

Using the values in
implies that v , and the vector feeds at its preferred rate, v .When v , h must be greater than v , so the biting rate will depend on the vector and not the host population.This would result to high vector host contact and high malaria transmission.In this case and malaria will persist in the population  In this study, we developed a vector f Dis ration model for the spread of malaria with disease induced deaths and varying human and host populations.
We have shown that the two populations drive the entire infection process through the threshold population density ratio v Q , which plays a vital role in the basic reproduction number.Our model captures the natural fluctuations known to occur in mosquito and host populations in malaria dynamics, and the effect of varying contacts between the vector and the host.The inclusion of diseaseinduced deaths was also of importance.Having in mind that majority of malaria deaths occur in children in Kenya, Mathematical analysis was done to establish that in the absence of the disease, a disease-free equilibrium will always exist if 0 1    .In the presence of the disease, that is when 0 1 , an endemic equilibrium is established with the infec opulations greater than zero.We observe from Figure 3 that a decrease in v Q increases   and vice versa.When the human popul-tious p 0 ation is very low, mosquitoes will turn to other bloodmeal source, and malaria transmission goes down.Then v Q controls the magnitude of malaria transmission.This implies that the best methods of controlling malaria in the highlands should target the adult mosquito, its biting habits and alternative sources of blood meals.Our results are consistent with results in literature that, 0   is a threshold that completely determines the global dynamics of disease transmission [38].

2 f
sudden for the model described by see Figure1.

Figure 1 .
Figure 1.A graphical comparison between f 1 and f 2 .


by natural birth.The transmission term πv h h h v I b S N  corresponds to frequency dependent infection of susceptible hosts by infectious mosquitoes, on infection they move to the infectious compartment.The infected hosts who recover h h I  become susceptible again as malaria has no permanent immunity.The last terms , deaths of the susceptible, infected
Then, from the last expression we obtain

Q 3 h
as shown in Figures3(b) and (d).If however the host irritation from vector bites is high, then we assume a

Figure 3 .= 2 .
Figure 3. (a) Variation of susceptible host and vector populations at the Disease Free equilibrium.Q v = 0.3 and 0  R = 0.718 < 1; (b) variation of susceptible host and vector populations at the Endemic equilibrium.Q v = 0.09 and 0  R = 2.495 > 1; (c) variation of Infected Host, Infected Vector, Total Host and Total Vector populations at the Disease Free equilibrium.Q v = 0.3 and 0  R = 0.718 < 1; (d) variation of Infected Host, Infected Vector, Total Host and Total Vector populations at the Endemic Equilibrium.Q v = 0.09 and 0  R = 2.495 > 1.

Table 1 . Variables used in the model related to infection contact process. Variable Definition
h S Susceptible Host Population Density v S Susceptible Vector Population Density h I Infectious Host Population Density v I Infectious Vector Population Density h N Total Host Population Density (Constant) v N Total Vector Population Density (Variable)

Table 2 . Parameters used in the model related to infection contact process.
h  Hosts Rate of Recovery v Q Vector-Host Ratio above Which Per-Vector Biting Saturates h b Host Irritability Biting Threshold v b Preferred (max.)Vector Feeding Rate , h 

Table 3
we simulate system 2, and show the dynamics of various populations in the figures below.