Comparison principles and applications to mathematical modelling of vegetal meta-communities

This article partakes of the PEGASE project the goal of which is a better understanding of the mechanisms explaining the behaviour of species living in a network of forest patches linked by ecological corridors (hedges for instance). Actually we plan to study the effect of the fragmentation of the habitat on biodiversity. A simple neutral model for the evolution of abundances in a vegetal metacommunity is introduced. Migration between the communities is explicitely modelized in a deterministic way, while the reproduction process is dealt with using Wright-Fisher models, independently within each community. The large population limit of the model is considered. The hydrodynamic limit of this split-step method is proved to be the solution of a partial differential equation with a deterministic part coming from the migration process and a diffusion part due to the Wright-Fisher process. Finally, the diversity of the metacommunity is adressed through one of its indicator, the mean extinction time of a species. At the limit, using classical comparison principles, the exchange process between the communities is proved to slow down extinction. This shows that the existence of corridors seems to be good for the biodiversity.


Introduction
This article partakes of a research program aimed at understanding the dynamics of a fragmented landscape composed of forest patches connected by hedges, which are ecological corridors. When dealing with the dynamics of a metacommunity at a landscape scale we have to take into account the local competition between species and the possible migration of species.
We are interested here in the mathematical modelling of two species, on two forest patches linked by some ecological corridor. We model the evolution by a splitting method, performing first the exchange process (see the definition of the corresponding Markov chain in the sequel) on a small time step, and then we perform independently on each station a birth/death process according to the Wright-Fisher model, and we reiterate.
Our first mathematical result is to compute the limit equation of this modelling when the time step goes to 0 and the size of the population diverges to ∞. This issue, the hydrodynamic limit, i.e. to pass from the mesoscopic scale to the macroscopic one received increasing interest in the last decades (see for instance in various contexts [2], [11], [17]). As our main results on extinction times do not require the convergence in law of the processes, instead of using a martingale problem ( [8]), we prove directly the convergence of operators towards a diffusion semi-group ( [10]). We find a deterministic diffusion-convection equation, where the drift comes from the exchange process, while the diffusion comes from the limit of the Wright-Fisher process. We point out here that the fact that the diffusion operator L d satisfies a non standard comparison principle (or a maximum principle) is instrumental: first the comparison principle ensures the uniqueness of the limit of the approximation process and then the definition of the Feller diffusion process. Then this comparison principle yields our second result that is concerned with the comparison of the extinction time of one species for a system with exchange and a system without exchanges. Assuming that the discrete extinction time converges, we prove that the limit is solution of the equation −L d τ = 1. Taking advantage once again of comparison principles, we prove that the exchange process slows down the extinction time of one species. Thus, the fragmentation of the habitat seems to be good to the biodiversity.
This article outlines as follows. In a second section we describe the modelling at mesoscopic scale. We couple a Wright-Fisher model for the evolution of the abundances together with an exchange process. The third section is devoted to the large population limit of the discrete process. In a fourth section we discuss the issues related to the extinction time; we compare the extinction time of one species with and without exchange process. In a final section we draw some conclusion and prospects for ecological issues, and we address the question of convergence in law for our model.
2. The mathematical model 2.1. Modelling the exchange between patches. Consider two patches that have respectively the capacity to host (N 1 , N 2 ) individuals, to be chosen into two different species α and β. Set (y n 1 , y n 2 ) for the numbers of individuals of type α, respectively in patch 1 and 2, at time nδt, i.e. after n iterations and δt is the time that will be defined below.
The exchange process is then simply modelled by where κ is the instantaneous speed of exchanges and d = N 2 N 1 represents the distortion between the patches (the ratio between the hosting capacities); we may assume without loss of generality that d ≤ 1. With this modelling, and assuming that κδt ≤ 1, it is easy to check that is mapped into itself, i.e stable, by the exchange process.
• The total population of individuals of type α, y n 1 + y n 2 , is conserved. • If we start with only individuals of species α (respectively β) then we remain with only individuals from α (respectively β); this reads (N 1 , N 2 ) → (N 1 , N 2 ) (respectively (0, 0) → (0, 0)). Set x = (x 1 = y 1 N 1 , x 2 = y 2 N 2 ) belonging to D = [0, 1] 2 for the population densities of a species α on two separate patches and x n = (x n 1 , x n 2 ) for these densities at time nδt. Then we have alternatively This reads also x n+1 = Ax n where A is a stochastic matrix. Consider now the piecewise constant càdlàg process with jumps X → AX at each time step δt. In other words, for any continuous function f defined on D = [0, 1] 2 then P ex δt (f )(x) = f (Ax), where P ex δt is the transition kernel of the exchange process.

2.2.
Wright-Fisher reproduction model. On each patch we now describe the death/birth process that is given by the Wright-Fisher model. The main assumption is that the death/birth process on one patch is independent of the other one. Consider then the first patch that may host N 1 individuals. The Markov chain is then defined by the transition matrix, written for Since the two Wright-Fisher processes are independent, the corresponding transition kernel reads for any function f defined on D = [0, 1] 2 . Notice that P wf δt is a two-variable version of the usual Bernstein polynomials. In the sequel, we will also use the notation B N (f ) and write for the sake of conciseness 2.3. The full disrete model. Starting from the state x = (x 1 , x 2 ), during a time step, we apply first the exchange process and then the Wright-Fisher reproduction process. In this way, the sequence of random variables x n is a Markov chain with state space 1} and the transition kernel reads as follows

From discrete model to continuous one
We consider the same scaling as for the Wright-Fisher usual model, that is N 1 δt = 1. We set N = N 1 in the sequel to simplify the notations. We may consider either the càdlàg process associated to the reproduction-exchange discrete process defined by x t = x n if nδt ≤ t < (n + 1)δt or the continuous piecewise linear function x t such that x t = x n for t = nδt. We consider an analogous interpolation in space in order to deal with function that are defined on [0, T ] × D where T > 0 is given.
For a given continuous function f that vanishes at (0, 0) and (1, 1), we now define the sequence of functions (5) u We may also use analogously u N (t, The fonctions u N and u N represent the average densities of the species at a macroscopic level. If X N is the Lagrangian representation of the densities, then u N represents the densities in Eulerian variables. and with initial data u(0, x) = f (x).

Remark 3.2.
We may have proved that the càdlàg process associated to the reproductionexchange process u N converges to a diffusion equation. We will discuss this in the sequel. Besides, we prove the convergence results for a sufficiently smooth f , and we will extend in the sequel the definition of a mild solution to the equation for functions f in the Banach space E = {f ∈ C(D); f (0, 0) = f (1, 1) = 0}. The theory for Markov diffusion process and the related PDE equations is well developed in the litterature (see [1], [9], [16] and the references therein). The particularity of our diffusion equation is that the boundary of the domain is only two points.

3.2.
Proof of Theorem 3.1. The proof of the theorem is divided into several lemmata. The first lemma describes in a way how the discrete process is close to a martingale.
Lemma 3.3. The conditional expectation of the discrete reproduction-exchange process is

Proof
Using the properties of the Bernstein polynomials, Then the proof of the lemma is completed, observing that A − Id = o(1).
The following lemma is useful to prove that x t and x t are close.
Lemma 3.4. There exists a constant C such that , then the following conditional expectation reads We expand the ℓ 2 norm in R 2 as We first have by linearity and by the Lemma 3.3 above that that completes the proof of the lemma.
The next statement is a consequence of the inequality |x t − x t | ≤ |x n − x n+1 | for t ∈ (nδt, (n + 1)δt) and of the previous lemma The processes x t and x t are asymptotically close, i.e. there exists a constant C such that As a consequence, when looking for the limit when N diverges towards +∞ of the process, we may either work with x t or x t . The next lemma is a compactness result on the bounded sequence u N defined in (5).
Lemma 3.6. There exists a constant C that depends on ||f || lip and on T such that for any, x, y in D and s, t in [0, T ], Remark 3.7. Since the constants C do not depend on N we can infer letting N → ∞ some extra regularity results for u, assuming that f is Lipschitz.

Proof
We begin with the first estimate. Introduce n such that nδt ≤ t < (n + 1)δt. Set y t for the process that starts from y = y 0 .
therefore, proving the first inequality for u N (which amounts to controlling |x n − y n |)) will imply the inequality for u N . Due to the properties of Bernstein's polynomials we have that where ω(f, 1 N ) is the modulus of continuity of f . Then, using that ||A − Id|| ≤ CN −1 , we infer that Iterating in time we have that, The other derivative is similar and then we infer from this computation that the first inequality in the statement of Lemma 3.6 is proved. We now proceed to the proof of the second one. Introduce the integers m, n such that mδt ≤ s < (m + 1)δt and nδt ≤ t < (n + 1)δt. Using that On the one hand, by Lemma 3.4 we have that the first term in the right hand side of (11) is bounded by above by C(m−n) N . On the other hand, using that the E(y j |x j ) = (A−Id)x j then Since ||Id − A|| ≤ CN −1 then the right hand side of (12) is also bounded by above by . This completes the proof of the lemma.
Thanks to Ascoli's theorem, up to a subsequence extraction, u N converges uniformly to a continuous function u(t, x). We now prove that u is solution of a diffusion equation whose infinitesimal generator is defined as the limit of N (P wf δt P ex δt − Id). Lemma 3.8. Consider f a function of class C 2 on D that vanish at (0, 0) and (1, 1). Then where L d is defined in Theorem 3.1.

Proof
Due to Taylor formula Using that the linear operator P wf δt is positive and bounded by 1 we then have The well-known properties of Bernstein polynomials (see [6]) entail that uniformly in x On the other hand, the operator P wf δt is the tensor product of two one-dimensional Bernstein operators. Then by Voronovskaya-type theorem (see [6]), for f ( 2d f x 2 x 2 . By density of the linear combinations of tensor products f 1 (x 1 )f 2 (x 2 ) this result extend to general f as (15) lim Denoting ∆ d the diffusion operator defined by the right hand side of (15), the Kolmogorov limit equation of our coupled Markov process is with initial data u(0, x) = f (x). Let us observe that u, the limit of E(f (x t )|x 0 = x), vanishes at two points (0, 0) and (1, 1) in the boundary ∂D.
We now complete the proof of the Theorem. Considering f such that the convergence in Lemma 3.8 holds. Then, for n ≤ tN < n + 1, N (s, .)))(x)ds.
Using the uniform convergence of u N , Lemma 3.8 and a recurrence on n we may prove that at the limit where we have omitted the variable x for the sake of convenience. We now state a result that ensures the uniqueness of a solution to the diffusion equation (18). Such a solution is a solution to the diffusion equation in a weak PDE sense. Introduce Remark 3.9. We precise here the regularity of the functions f in D(L d ). Since L d is a strictly elliptic operator on any compact subset of the interior of D then f is C 2 (D) ∩ C(D) (see [12]). The regularity of f up to the boundary is a more delicate issue (see [14], [15]). Besides, to determine exactly what is the domain of L d is a difficult issue. For PDEs the unbounded operator is also determinated by its boundary conditions. Here we have boundary conditions of Ventsel'-Vishik type, that are integro-differential equations on each side of the square linking the trace of the function f and its normal derivative. This is beyond the scope of this article.
• Parabolic version: Consider a func- We postpone the proof of this theorem until the end of this section. We point out that a comparison principle for L d is not standard since it requires only information on two points {(0, 0), (1, 1)} in ∂D and not on the whole boundary. Theorem 3.10 implies uniqueness of the limit solution. Therefore the whole sequence u N converge and the semigroup is well defined. Actually, setting S(t)f = u(t) we then have defined for smooth f the solution to a Feller semigroup (see [1]) as follows (1) S(0) = Id.
The second property comes from uniqueness, the last one passing to the limit in The third one is then simple. The third property allows us to extend the definition of S(t) to functions in E by a classical density argument. Then we have a Feller semigroup in E that satisfies the assumptions of the Hille-Yosida theorem (see [5]).

3.3.
Proof of the comparison principle. We begin with the comparison principle for the parabolic operator. We use that C 2 (D) is dense in D(L d ), i.e. that any function u in D(L d ) can be approximated in E by smooth functions u k up to the boundary, and such that Lu k converges uniformly on any compact subset ofD. We then prove the comparison principle for smooth functions and we conclude by density.
Consider u as in the statement of the Theorem for a C 2 initial data f . Consider ε small enough.
We then have v t (t 0 , Let us observe that if ε is chosen small enough Second case: x 0 belongs to ∂D but the four corners. We may assume that x 0 = (0, x 2 ) the other cases being similar. We have that We then have as in (20) We now conclude. since v is nonnegative we have Letting ε goes to 0 completes the proof. Let us prove the elliptic counterpart of the result for a smooth function u (we also proceed by density). Set as above v(x) = u(x) + εψ(x) + ε 2 θ(x). Introduce x 0 where v achieves its minimum, i.e. v(x 0 ) = min D v(x).. First if x 0 belongs to the interior of D, then L d v(x 0 ) > 0 and we have a contradiction. We disprove the case where x 0 belongs to the boundary but {(0, 0), (1, 1)} exactly as in the evolution equation case. Assume first that x 0 belongs to ∂D but the four corners; for instance x 0 = (0, x 2 ).
Then L d (ψ + εθ)(x 0 ) =< 0 gives the contradiction. Assume then that x 0 = (0, 1).  1) are absorbing states. These two absorbing states correspond to the extinction of a species. Let us introduce the hitting time Θ N that is the random time when the Markov chain reaches the absorbing states, i.e. the extinction time. Since the restriction of the chain to the non absorbing states is irreducible and since there is at least one positive transition probability from the non absorbing states to the absorbing states then this hitting time is almost surely finite. This result is standard for Markov chains with finite state space (see [4], [7] and the references therein).
Let U be the complement of the trapping states (0, 0) and (1, 1). Consider the vector T N defined as the conditional expectation (T N ) j N ∈U = E j N (Θ N ) of this hitting time and denote byP N orP the restriction of the transition matrix to U .
Then for x ∈ U , denoting P x the conditional probability, we have using Markov property and time translation invariance This is equivalent to We are now interested in the limit of T N when N diverges towards ∞. Let us recall that for the one dimensional Wright-Fisher process the expectation of the hitting time starting from x converges towards the entropy H(x) (see [16]) defined by The entropy is a solution to the equation − x(1−x) 2 H xx = 1 that vanishes at the boundary. The proof, that can be found in Section 10 of [9], uses probability tools like the convergence in distribution of the processes and the associated stochastic differential equation. We believe that the same kind of tools would give the convergence of τ N in dimension two but this is beyond the scope of this article. Besides, for the sake of completeness we provide a proof for the convergence in distribution of our processes in Section 5.2 below.
Set now τ N for the polynomial of degree N in x 1 and x 2 that interpolates T N at the points of the grid. We have Theorem 4.1 (Extinction time). When N diverges to +∞ the sequence τ N converges towards τ that is solution to the elliptic equation −L d τ = 1.
Assuming the convergence of τ N , the proof of the theorem is straightforward by passing to the limit in (23) using Lemma 3.8.

Remark 4.2.
We expect the function τ to be smooth up to the boundary but at the two points (0, 0) and (1, 1). We admit here this result. This allows us to use the previous comparison result.

4.2.
Exchanges slow down extinction. Consider now a single patch whose hosting capacity is N 1 + N 2 = (d + 1)N for N = 1 δt . The limit equation for the classical Wright-Fisher related process is Then the corresponding extinction time for the Wright-Fisher process without ex- is the corresponding averaged starting density (see [16]) and where H is the entropy defined above (24). We shall prove in the sequel

Proof
We point out that to check that −L d satisfies the comparison result is not obvious (see Theorem 3.10). We first observe that the entropy (24) vanishes at the boundary points {(0, 0), (1, 1)}. Setting τ (x 1 , x 2 ) = g(z), we have We then have Observing that by a mere computation we have that τ is a subsolution to the equation.

4.3.
More comparison results. We address here the issue of the convergence of the limit extinction time τ = τ d,κ defined in Section 4 when κ or d converges towards 0. This extinction time depends on the starting point x.

Proof
Consider here the function . This function vanishes at x = (0, 0) and x = (1, 1) and satisfies Then V is a subsolution to the equation −L d τ = 1 and by the comparison principle V ≤ τ d,κ everywhere. Letting κ → 0 completes the proof of the Proposition.
Proposition 4.5. Assume κ be fixed. Then that is the extinction time for one patch.

Proof
We begin with .
Let us observe that due to (28) The strategy is to seek a supersolution X to the equation −L dX = 1 d that is bounded when d converges to 0. We first have, using the entropy function H 2 (x 1 , Therefore, since we have we obtain, for d small enough to have (1 + 2κ)d < 1, ).
Using the estimate we have that if d is small enough depending on κ then −L d (H 2 + D) ≥ 1 2d . Using the comparison principle we then have that and we conclude by letting d converge to 0 since τ converges towards H(x 1 ).

Miscellaneous results and comments
5.1. Discussion and prospects for ecological issues. To begin with, we have introduced a split-step model that balances between the local reproduction of species and the exchange process between patches. This split-step model at a mesoscopic scale converges towards a diffusion model whose drifts terms come from the exchanges. This has been also observed for instance in [19].
Here we deal with a neutral metacommunity model with no exchange with an external pool. Hence the dynamics converge to a fixation on a single species for large times. The average time to extinction of species is therefore an indicator of biodiversity. Here for our simple neutral model, Theorem 4.3 provides a strong reckon that the exchange process is good for the biodiversity. In some sense, the presence of two patches allows each species to establish itself during a larger time lapse.
In a forthcoming work we plan to numerically study a similar model but with more than two patches and several species. We plan also to calibrate this model with data measured in the south part of Hauts-de-France. The main interest is to assess the role of ecological corridors to maintain biodiversity in an area. The question of the benefit of maintaining hedges arises when the agricultural world works for their removal to enlarge the cultivable plots. This is one of the issue addresses by the Green and Blue Frame in Hauts-de-France.

Convergence un distribution.
We address here the convergence in law/in distribution of the infinite dimensional processes related to the x t N . This is related to the convergence of the process towards the solution of a stochastic differential equations; we will not develop this here. Following [18] or [13], it is sufficient to check the tightness of the process and the convergence of the finite m-dimensional law.
Dealing with x t N instead of x t N , the second point is easy. Indeed, Theorem 3.1 implies the convergence of the m-dimensional law for m = 1. We can extend the result for arbitrary m by induction using the Markov property. For the tightness, we use the so-called Kolmogorov criterion that is valid for continuous in time processes (see [18] chapter 2 and [13] chapter 14); this criterion reads in our case This is a consequence of the following discrete estimate, since x t N is piecewise linear with respect to t, Proposition 5.1. There exists a constant C such that for any m < n (37) E(|x n − x m | 4 ) ≤ C |n − m| 2 N 2 . Proof First step: using that x n is close to a true martingale.
Let us set A = Id − κ N M = Id − B. Introduce z 0 = x 0 and z n = x n + B k<n x k . Then since E(x n+1 |x n ) = x n − Bx n , we have that z n is a martingale. Moreover we have the estimate, for 0 ≤ m < n Second step: computing the fourth moment.
To begin with we observe that, due to (38) Therefore we just have to prove that (37) is valid with z n replacing x n . We introduce the increment y j = z j+1 − z j . We then expand as follows, setting |.| and (., .) respectively for the euclidian norm and the scalar product in R 2 .
Since y l is independent of the past, if for instance l > max(i, j, k) then E((y i , y j )(y k , y l )) = 0. Therefore, (40) reads also Third step: handling D 4 and D 3 .
The key estimate reads as follows Fourth step: handling D 1 and D 2 .
Using the conditional expectation we have Due to (42) and Cauchy-Schwarz inequality E(|y k | 2 |x k )) = O(N −1 ) and it follows We now handle D 2 exactly as we did for D 1 . This completes the proof.