AVIAN INFLUENZA DYNAMICS IN WILD BIRDS WITH BIRD MOBILITY AND SPATIAL HETEROGENEOUS ENVIRONMENT

In this paper, we propose a mathematical model to describe the avian influenza dynamics in wild birds with bird mobility and heterogeneous environment incorporated. In addition to establishing the basic properties of solutions to the model, we also prove the threshold dynamics which can be expressed either by the basic reproductive number or by the principal eigenvalue of the linearization at the disease free equilibrium. When the environment factor in the model becomes a constant (homogeneous environment), we are able to find explicit formulas for the basic reproductive number and the principal eigenvalue. We also perform numerical simulation to explore the impact of the heterogeneous environment on the disease dynamics. Our analytical and numerical results reveal that the avian influenza dynamics in wild birds is highly affected by both bird mobility and environmental heterogeneity.


1.
Introduction. Aquatic birds such as Anseriformes (ducks, geese and swans) and Charadriiformes (gulls, terns and waders) are the major reservoir of all influenza A viruses, including the highly pathogenic H5N1 AI virus transmitted to humans [9,17,22,36]. Understanding of the ecology of avian influenza (AI) virus and its dynamics in wild birds is useful in predicting influenza dynamics in human population and devising control strategies.
In the context of AI dynamics in wild birds, in addition to transmission of AI from bird to bird, another highly efficient route of transmission is through the excretion of AI virus by infected birds, followed by ingestion of virus in the drinking water of uninfected birds [4,12,13,37]. Laboratory experiments have shown that the persistence of AI virus in water depends on environmental factors such as temperature, pH and salinity [4,5,6,28,37]. Numerous mathematical modeling studies have already highlighted the importance of such environmental factors in AI dynamics among wild birds [3,10,23,24,29,32]. In particular, [32] has shown that even small differences in the environmental condition between two locations (for example, differences in temperature of less than 8 o C) can produce significantly different AI dynamics. In addition, mobility is a common strategy for birds to occupy seasonal habitats [2]. Because of their mobility, birds come across various environmental conditions. Therefore, it is quite obvious to raise the question of how such spatial heterogeneity affects AI dynamics in mobile aquatic wild birds.
In this study, we propose a dynamic model for the transmission of AI among wild birds that incorporates both mobility of the birds and spatial heterogeneity in the environmental condition. We introduce into the model a spatial diffusion of the bird population and use the experimentally-determined dependence of virus persistence in water upon temperature, pH and salinity. We focus on how spatial diffusion and environmental heterogeneity affect the basic reproductive number and threshold dynamics of the system.
One of the main technical difficulties in our analysis is the lack of compactness of solution maps of the model system (see also [15]). This is because one of the equations, which describes the dynamics of virus particles in water, loses its diffusion term as the diffusion of virus particles is negligible compared to the birds' mobility. To overcome this difficulty, we first prove that the solution maps associated with a linearized system around the disease-free equilibrium are κ-condensing, where κ is the Kuratowski measure of non-compactness (see, e.g., [7]). By a generalized Krein-Rutman Theorem, we can show that the principal eigenvalue of the associated eigenvalue problems exists. Next, we prove the solution maps associated with our model system are κ-contracting. Thus, we conclude that the solution maps admit a connected global attractor by appealing to some existing results in [21].
The basic reproductive number for an infectious disease, conventionally denoted by R 0 , is an important index in epidemiology which predicts whether an infectious disease will die out or persist in the host population. For models given by infinite dimensional systems, such as our model system with spatial structure, identifying R 0 is usually not trivial. For our model, by making use of the abstract results on this topic in [31], we are able to find the so called next generation operator L, and thereby, identify R 0 as the spectral radius of L. We show that R 0 plays a threshold role in the sense that when R 0 < 1, the disease dies out from the birds population while when R 0 > 1, the disease remains persistent. We also carry out model simulations to observe how AI dynamics vary with diffusion rate and spatial heterogeneity of the environmental conditions. The rest of the paper is organized as follows. The model is formulated in Section 2. The model analysis and simulation results are presented in Sections 3 and 4, respectively. Finally, we state the conclusions of the paper in Section 5.

Model.
Suppose Ω ⊂ R n is a bounded domain which is the habitat of a host (for the AI virus) bird species. We divide the total bird population into susceptible (S), infected (I) and recovered (R) groups. We further consider AI virus concentration in water as V . We introduce the birds' mobility by the spatial diffusion terms in the model. Following [32], we incorporate environmental effects via the spatiallyvarying viral decay rate c(x), which depends on environmental conditions such as temperature, pH and salinity. The model we consider is as follows: ∂R(x, t) ∂t with (x, t) ∈ Ω × (0, ∞). Here, D is the diffusion coefficient and ∆ is the Laplace operator. As included in the model, the susceptible birds get infected by direct bird-to-bird transmission at rate β I SI, and by indirect fecal-oral transmission at rate β V SV . The parameters λ, d, γ and η represent the rate of recruitment of susceptible birds, the rate of natural death, the rate of recovery from infection and the rate of immunity loss, respectively. Since AI virus is generally non-pathogenic in wild birds [22], we have ignored disease caused deaths in our model. Infected birds shed virus particles in their feces at rate α.
Here, we use the homogeneous Neumann boundary condition and initial conditions where ∂ ∂ν denotes the differentiation along the outward normal ν to ∂Ω. Here we assume that the habitat dependent parameter c(x) is strictly positive and continuous onΩ.
Furthermore, the semiflow Ψ t : X + → X + is point dissipative and the positive orbits of bounded subsets of X + for Ψ t are bounded.

3.2.
The basic reproductive number. The basic reproductive number, which is defined as the average number of secondary infections generated by a single infected individual introduced into a completely susceptible population, is one of the important quantities in epidemiology. In this subsection, we will identify the basic reproductive number for the model system (7)- (9).
Obviously E 0 = (w * , 0, 0, 0) is the disease free equilibrium of (7) where w * = λ d . Note that by Lemma 3.2, w * is the unique positive equilibrium of (16) and is globally asymptotically stable in C(Ω, R). Linearizing system (7)-(9) at E 0 , we get the following system for the infection related variables u 2 and u 4 : We first consider the following generalized version of the system (20): (21) leads to the following associated eigenvalue problem: It is easy to see that the system (21) is co-operative, but its solution map is not compact since the second equation in (21) has no diffusion term. The following lemma deals with the existence of the principal eigenvalue of (22).
has a principal eigenvalue, denoted by µ(H) which is associated with a strongly positive eigenfunction.
By the same argument as that in the proof of Lemma 3.5 in the next subsection, we can show that for each t > 0, Π t is an κ-contraction on Y in the sense that for any bounded set B in Y, wherec := min x∈Ω c(x) > 0 and κ is the Kuratowski measure of non-compactness as defined in (28). From the discussions above, it is easy to see that the solution map Π t generated by (21) is κ-condensing in the sense that Note that (21) is a cooperative system. By the generalized Krein-Rutman Theorem (see, e.g., [16,Lemma 2.2]) and [11,Chapter II.14], the equation (22) has a principal eigenvalue, denoted by µ(H), with an associated eigenvector ψ * = (ψ * 2 , ψ * 4 ) 0.
With the above preparation, we can now employ the ideas and theory in [8,18,31,33,35] to the linearized system (20) to define the basic reproductive number for the system (7)- (9). Assume that population is near the disease free equilibrium (w * , 0, 0, 0). Let ϕ := (ϕ 2 , ϕ 4 ) be the spatial distribution of (u 2 , u 4 ) and S(t)ϕ be the solution semiflow generated by the following linear system: Let T 2 (t) and T 4 (t) be the semigroup defined in (11) and (13), respectively. From the first two equations of (23), it follows that u 2 (·, t, ϕ) = T 2 (t)ϕ 2 and That is, It then follows that S(t) is a positive C 0 -semigroup on C(Ω, R 2 ) and S(t)ϕ represents the spatial distribution of u 2 and u 4 at time t > 0.
Let C be the positive linear operator on C(Ω, R 2 ) defined by

AVIAN INFLUENZA DYNAMICS IN SPATIAL HETEROGENEOUS ENVIRONMENT 2835
where Then, at time t > 0 and location x, there will be C 2 (S(t)ϕ)(x) individuals added per unit time into the u 2 compartment, and hence C(S(t)ϕ)(x) accounts for the infection force at time t and location x. Thus, the spatial distribution of total new infected individuals caused by the initial infective distribution ϕ = (ϕ 2 , ϕ 4 ) is Along the line of [29] where no spatial factor is considered, we can define the next generation operator L by By [31], the basic reproductive number for system (7)- (9) is given by the spectral radius of L, that is, Also, by the general results in [31] and the same arguments as in [35, Lemma 2.2], we have the following conclusion on R 0 and µ(w * ). Proposition 1. R 0 − 1 and µ(w * ) have the same sign.

Threshold dynamics.
In this subsection, we show that R 0 is, in fact, a threshold index for disease persistence. Since the last equation in (7) has no diffusion term, its solution map Ψ t is not compact. In order to overcome this problem, we introduce the Kuratowski measure of non-compactness (see [7]), κ, which is defined by κ(B) := inf{r : B has a finite cover of diameter < r}, for any bounded set B. We set κ(B) = ∞ whenever B is unbounded. It is easy to see that B is precompact (i.e.,B is compact) if and only if κ(B) = 0.
Further, Ψ t is κ-contracting on X + in the sense that Proof. It is easy to see that u 4 (·, t, φ) satisfies the following equations:

NAVEEN K. VAIDYA, FENG-BIN WANG AND XINGFU ZOU
Then Motivated by the discussion above, we define the following operators: and for any bounded set B ⊂ X + . It is easy to see that Q(t) : X + → X + is compact for each t > 0 and hence κ(Q(t)B) = 0, ∀ t ≥ 0. It is easy to see that there exists a real numberc := min x∈Ω c(x) > 0 such that c(·) ≥c and it then follows that and hence L(t) ≤ e −ct . Consequently, Thus, Ψ t is κ-contraction of order e −ct on X + . This implies Ψ t is κ-contracting on X + .
Theorem 3.6. Ψ t admits a connected global attractor on X + .
Proof. By Lemma 3.3, it follows that Ψ t is point dissipative on X + and that the positive orbits of bounded subsets of X + for Ψ t are bounded. Furthermore, Ψ t is κ-contracting on X + by Lemma 3.5. By [21, Theorem 2.6], Ψ t has a global attractor that attracts each bounded set in X + .
The following results will play an important role in establishing the persistence of (7)- (9).

AVIAN INFLUENZA DYNAMICS IN SPATIAL HETEROGENEOUS ENVIRONMENT 2837
Proof. It is easy to see that u 2 and u 3 satisfy the following inequality: From the last equation of (7), we get that This implies that part (ii) is valid.
From (15), (16), Lemma 3.1 and Lemma 3.2, it follows that there exists a t 1 > 0 such that u 2 (x, t) ≤ 2w * , ∀ t ≥ t 1 . From the last equation of (7), we get wherec := min x∈Ω c(x). Thus, there exists t 2 ≥ t 1 such that u 4 (x, t) ≤ 4 αw * c , t ≥ t 2 . The first equation of (7) gives ∂u1(x,t) ∂t which completes the proof. Now we prove the main result of this section, which shows that R 0 is a threshold index for disease persistence.

3.4.
Homogeneous environment: c(x) = c. In this subsection, we discuss the special case when c(x) = c, a positive constant. The space independence will allow us to obtain more explicit results on extinction and persistence of the disease (see also [38,Lemma 3.1]).

AVIAN INFLUENZA DYNAMICS IN SPATIAL HETEROGENEOUS ENVIRONMENT 2841
shows that µ * is the maximum root of equation (38), and hence where Next, we show that the basic reproductive number R 0 is indeed identical to R * 0 in this case. To this end, we use the idea in [35, Theorem 2.1] to consider the following operator obtained by perturbing L: where > 0 is a constant and Then L is a strongly positive linear operator. Moreover, by Ω Γ(x, y, t)dy = 1, ∀ x ∈ Ω, t > 0, it is easy to see that Since each element of matrix J is positive, its spectral radius r(J ) is an eigenvalue corresponding to which there is a positive eigenvector in R 2 . Straightforward calculation shows that If we can prove that r(L ) = r(J ), for > 0, (40) then by letting → 0 + , we would obtain Unfortunately, we are unable to prove (40). The main difficulty is that, like L, L is not compact; otherwise the uniqueness of the positive eigenvector associated with r(L ) would immediately confirm (40). Note that the formulas (39) and (41) coincide with the conclusion of Proposition 1. This makes us conjecture that (40) holds.
Below we provide an alternative way to confirm (41). Assume that initially there is no virus in the host population and a single infective individual is brought into the host population with the probability of landing at the location x being ϕ 2 (x) ≥ 0 (hence 0 ≤ ϕ 2 (x) ≤ 1 and Ω ϕ 2 (x) dx = 1). This corresponds to ϕ 4 = 0. By (25), the spatial distribution of total new infected individuals caused by such a single infective individual is given by Summing up the above over the spatial domain, we then obtain the total new infected individuals caused by the single infective individual as Making use of Ω Γ(x, y, t)dx = Ω Γ(x, y, t)dy = 1 and Ω ϕ 2 (x) dx = 1, as well as changing the orders of integrals in the above, we can show that the above quantity involving multiple intergrals actually is equal to On the other hand, by the biological definition of the basic reproductive number, (43) is nothing but R 0 . This confirms (41).  Table 1). It is quite difficult to obtain the actual function c depending on the spatial variable x. In fact, it varies from location to location and also from time to time. For illustration purposes, we consider a linearly decay function of temperature, i.e. T (x), to represent the spatial variation of the environment: T (0) = T 0 and T (1) = T 1 stand for temperatures of a warmer place and a cooler place, respectively, while the temperature in between is given by T (x) = (T 1 − T 0 )x + T 0 . Then as obtained in [32], the decay rate of viral particles, c(x) is given by the following relation: where a and b are constants (See Table 1). By taking different values of D and different temperature profiles, T , we study how the diffusion and the spatial heterogeneity of the environmental condition impact AI dynamics. Numerical solutions are obtained by using the method of lines to derive a system of ODEs, which are solved using MATLAB software.
As a case study, we consider spatial temperature variation in Canada [1]. As given in [1], the mean annual temperature in Canada varies from 13 o C at some places along the southern border to -18 o C in the north. Therefore we take T 0 = 13 o C and T 1 = −18 o C for our base case computation, while the base case value of the scaled diffusion coefficient is fixed at D = 5 × 10 −4 . We simulate the model for a six month time-frame, and analyze the spatial distribution of avian influenza prevalence (%) at the end of six months. Initially, birds are assumed to live locally in the middle of the domain (Fig. 1). We then observe how the spatial diffusion and the spatial environmental heterogeneity will drive the dynamics of avian influenza among these birds. We emphasize that this case study is not intended as a model of AI in Canada, but simply illustrates the pattern and magnitude of the effects we are studying.  Figure 1. Distribution of AI prevalence (%) among wild birds in a six month period for a model without spatial diffusion and spatial environmental heterogeneity (case 1), with spatial environmental heterogeneity only (case 2), with spatial diffusion only (case 3), and with both spatial diffusion and spatial environmental heterogeneity (case 4). The dotted line indicates the initial distribution.
We present spatial distributions of AI prevalence (%) predicted by the model for different cases (Fig. 1): (i) neither spatial diffusion nor spatial environmental heterogeneity, (ii) with spatial environmental heterogeneity only, (iii) with spatial diffusion only, and (iv) with both spatial diffusion and spatial environmental heterogeneity. We can clearly see a significant effect of both diffusion and environmental heterogeneity. As expected, without spatial diffusion, AI accumulates locally. While the AI prevalence without diffusion shown in Fig. 1 is at T = 5 o C, the level of AI prevalence is set by the value of T considered: the lower the temperature the higher the prevalence level, as the infectious viruses persist longer in cold temperatures. When diffusion is introduced without any environmental heterogeneity, AI spreads out symmetrically on both sides from the initial location. AI becomes asymmetrically distributed across the spatial domain if spatial environmental heterogeneity is taken into account in addition to diffusion.
We now fix the temperature distribution as T 0 = 13 o C, ∆T = T 1 − T 0 = −31 o C and observe how the spatial distribution of AI prevalence changes with the diffusion   (Fig. 2). A higher diffusion coefficient has tendency to increase AI prevalence towards the boundary while decreasing the prevalence in the middle as seen in Fig. 2. However, due to the spatial environmental heterogeneity, the effect of diffusion is more pronounced at the cooler boundary than the warmer boundary and the prevalence remains higher at cooler places, again due to longer persistence of viruses at low temperatures.
In Fig. 3, we show a spatial distribution of AI prevalence for different temperature profiles while fixing the diffusion constant at D = 5 × 10 −4 . Here, in each case we force the temperature in the middle of the domain to be always the same (5 o C), and then from the middle of the domain the temperature linearly increases to the x = 0 boundary and linearly decreases to the x = 1 boundary so that ∆T = 5, 15, 25 and 35 o C. AI prevalence is clearly affected by spatial environmental heterogeneity with a positive correlation between heterogeneity in environment and heterogeneity in the AI distribution. The results shown in Fig. 3 reveal that this effect is more sensitive at cooler places than warmer places.

5.
Conclusion. In this paper, we presented a transmission dynamic model of avian influenza among wild birds. The novelty of the model is that it includes both spatial diffusion of birds and spatial heterogeneity of the environment, which are critical in understanding AI dynamics and devising control strategies. The environmental heterogeneity was introduced into the model based on the experimentally-observed  dependence of virus persistence on the environmental factors such as temperature, pH and salinity. Mathematical analysis of the model allowed us to achieve a formula for the basic reproductive number and a threshold condition for the disease to die out. We found that the reproductive number is independent of the diffusion coefficient in the absence of environmental heterogeneity. However, diffusion comes into play to define the reproductive number due to spatial heterogeneity in the environmental condition.
We performed model simulations for various diffusion constants and environmental conditions. Our results show that the dynamics of AI prevalence among wild birds is highly affected by both bird diffusion and environmental heterogeneity. While diffusion has a tendency to spread AI across a larger space, the environmental heterogeneity brings an asymmetrical nature to the AI distribution. In our model, we have introduced only spatial heterogeneity. However, the environmental condition, for example temperature, varies widely even within a short period. Therefore, further extension of our work would be to analyze the effects of spatiotemporal variations of the environmental conditions on AI dynamics in wild birds. In c(x) -0.12 [32] b In c(x) 5.10 [32]