# Adaptation and migration of a population between patches

###### Abstract

A Hamilton-Jacobi formulation has been established previously for phenotypically structured population models where the solution concentrates as Dirac masses in the limit of small diffusion. Is it possible to extend this approach to spatial models? Are the limiting solutions still in the form of sums of Dirac masses? Does the presence of several habitats lead to polymorphic situations?

We study the stationary solutions of a structured population model, while the population is structured by continuous phenotypical traits and discrete positions in space. The growth term varies from one habitable zone to another, for instance because of a change in the temperature. The individuals can migrate from one zone to another with a constant rate. The mathematical modeling of this problem, considering mutations between phenotypical traits and competitive interaction of individuals within each zone via a single resource, leads to a system of coupled parabolic integro-differential equations. We study the asymptotic behavior of the stationary solutions to this model in the limit of small mutations. The limit, which is a sum of Dirac masses, can be described with the help of an effective Hamiltonian. The presence of migration can modify the dominant traits and lead to polymorphic situations.

Key-Words: Structured populations, phenotypical and spatial structure, Hamilton-Jacobi equation, viscosity solutions, Dirac concentrations, stationary solutions

AMS Class. No: 35B25, 47G20, 49L25, 92D15

## 1 Introduction

Non-local Lotka-Volterra equations arise in models of adaptive evolution of phenotypically structured populations. These equations have the property that the solutions concentrate generally, in the limit of small diffusion, on several isolated points, corresponding to distinct traits. Can we generalize these models by adding a spatial structure? How do the dominant traits evolve if we introduce a new habitat? To understand the interaction of ecological and evolutionary processes in population dynamics, spatial structure of the communities and adaptation of species to the environmental changes, it is crucial to dispose mathematical models that describe them jointly. We refer to [19] and the references therein for general literature on the subject. In this manuscript we consider a model where several distinct favorable habitable zones are possible. Population dynamics models structured by spatial patches have been studied using both deterministic and probabilistic methods (see for instance [30, 1]). Our model, in the case of two patches, is indeed very close to the one studied in [30] where the authors use adaptive dynamics theory (adaptive dynamics is a theory, based on dynamical systems and their stability, to study population dynamics [11]). Here we model similar phenomena, by adding a spatial structure to an earlier known integro-differential model describing the darwinian evolution. Integro-differential models have the advantage that the mutations can be considered directly in the model without assuming a separation of time scales of evolution and ecology. The present work provides a general description of the asymptotic stationary solutions, in the general case where two or several patches are possible.

We study the asymptotic behavior of solutions of a system of coupled elliptic integro-differential equations with small diffusion terms. These solutions are the stationary solutions to a parabolic system describing the dynamics of a population density. The individuals are characterized by phenotypical traits, that we denote by . They can move between two or several patches, which are favorable habitable zones, with constant rates (that we denote by and in the case of two patches). The mathematical modeling is based on the darwinian evolution and takes into account mutations and competition between the traits. There is a large literature for mathematical modeling and analysis on the subject of adaptive evolution, we refer the interested reader to [16, 15, 11, 12, 22, 10]. Here, we represent the birth and death term by a net growth term that is different in each patch, for instance because of a change in the temperature, and depends on the integral parameter , which corresponds to the the pressure exerted by the whole population within patch on the resource. To model the mutations, we use Laplace terms with a small rate that is introduced to consider only rare mutations. We study the asymptotic behavior of stationary solutions as the mutation rate goes to . The asymptotic solutions are generally concentrated on one or several Dirac masses. We describe the position and the weight of these Dirac masses using a Hamilton-Jacobi approach.

The time-dependent model, in the case of two patches, is written as

(1) |

with

(2) |

Such models, without the structure in space, have been derived from stochastic individual based models in the limit of large populations (see [7, 6]). This manuscript follows earlier works on parabolic Lotka-Volterra type equations to study concentration effects in models of phenotypically structured populations, that are based on a Hamilton-Jacobi formulation (see [12, 25, 2, 20]). The novelty of our work is that we add a spatial structure to the model by considering a finite number of favorable habitable zones. We thus have a system instead of a single equation. A Hamilton-Jacobi approach in the case of systems has also been introduced in [5] for an age structured model. See also [4] for a study of stationary solutions of the latter system. The Hamilton-Jacobi approach can also be used in problems other than adaptive evolution to prove concentration phenomena. See for instance [27, 26, 24] where related methods have been used to study the motion of motor proteins.

We are interested in the equilibria of 1 limited to a bounded domain, that are given by solutions of the following system

(3) |

where is a ball of radius with center in and is the unit normal vector, at the point , to the boundary of . The Neumann boundary condition is a way to express that mutants cannot be born in .

To formulate our results we introduce the assumptions we will be using throughout the paper. We assume that, there exist positive constants and such that

(4) |

Moreover there exist positive constants , , and such that, for all and ,

(5) |

(6) |

(7) |

We use the Hopf-Cole transformation

(8) |

and replace the latter in the system satisfied by to obtain

(9) |

We prove the following

###### Theorem 1.1

The function is indeed an effective Hamiltonian that contains information from the two patches and helps us in Theorem 1.2 to describe the support of the weak limits of as . We can interpret as the fitness of the system in the limit of (see [23] for the definition of fitness).

The difficulty here is to find appropriate regularity estimates on , that we obtain using the Harnack inequality [3] and the Bernstein method [9]. To prove convergence to the Hamilton-Jacobi equation, we are inspired from the method of perturbed test functions in homogenization [14].

The above information on the limit of allows us to describe the limit of the densities as vanishes. We prove

###### Theorem 1.2

Assume 4–7. Consider a subsequence such that and converge uniformly to and goes to , as , with solution of 10. Let , for , converge weakly in the sense of measures to along this subsequence. We have

(12) |

with

(13) |

Moreover, we have

(14) |

in the sense of distributions. The above condition is coupled by

(15) |

Theorem 1.2 provides us with a set of algebraic constraints on the limit, which allows us to describe the latter. In particular, if the support of , for , is a set of distinct points: , 14 implies that

(16) |

with

(17) |

Furthermore, the weights satisfy the normalization condition

(18) |

Condition 17 means that the vector is the eigenvector corresponding to the largest eigenvalue of the matrix at the point , which is . Thereby 17 implies once again that .

We point out that since , for , is such that the fitness vanishes on the support of and is negative outside the support, we can interpret as evolutionary stable distribution of the model. In adaptive dynamics, evolutionary stable distribution (ESD) corresponds to a distribution that remains stable after introduction of small mutants (see [21, 13, 17] for a more detailed definition). See also [10, 28] for related works on stability and convergence to ESD for trait-structured integro-differential models.

The set of assumptions in Theorem 1.2 allows us to describe the asymptotics of the stationary solutions, in the limit of rare or small mutations. In Section 5 we provide some examples where based on this information we can describe the asymptotics. In particular, we notice that the introduction of a new environment can lead to dimorphic situations. We refer to [8] for a related work using the Hamilton-Jacobi approach, where polymorphic situations can also appear in a model with multiple resources.

The paper is organized as follows. In Section 2 we prove some bounds on and some regularity properties on that allow us to pass to the limit as and derive the Hamilton-Jacobi equation with constraint. Theorem 1.1 is proved in Section 3. Using the results obtained on the asymptotic behavior of we prove Theorem 1.2 in Section 4. In Section 5 we provide some examples where the information given by Theorem 1.1 and Theorem 1.2 allows us to describe the limit. The asymptotic behavior of the stationary solutions in a more general framework, where more than two habitable zones are considered, is given in Section 6. Finally in Section 7 we present some numerical simulations for the time-dependent problem and compare them with the behavior of stationary solutions.

## 2 Regularity results

###### Lemma 2.1

###### Remark 1

###### Proof 1

We prove the result by contradiction. We suppose that (the case with , and the inequalities from below can be treated following similar arguments). We multiply the first equation in 3 by , integrate, and use 4 to obtain

Using now 5, 6 and the fact that we deduce that, for small enough,

and thus

Now we multiply the equations in 3 by , integrate and add them and use 4 to obtain

From 5, 6 and the above bounds on and it follows that

which is not possible if is small enough. We conclude that .

###### Theorem 2.2

###### Proof 2

(i) We define

From 3 we have

(22) |

Moreover, from 5, 6 and 19 we have, for ,

Therefore the coefficients of the linear elliptic system 22 are bounded uniformly in . It follows from the classical Harnack inequality ([3], Theorem 8.2) that there exists a constant such that for all such that and for ,

Rewriting the latter in terms of and and replacing by we obtain

(ii) To prove the Lipschitz bounds, we use the Bernstein method (see [9]). We assume that

(23) |

that is the maximum is achieved at a point and for (the case where the maximum is achieved for can be treated by similar arguments). From the Neumann boundary condition in 9 we know that is an interior point of . We define and notice that

We now differentiate the first equation in 9 with respect to and multiply it by to obtain

From 23 we have

and thus

Moreover from 23 we have and . It follows that

Using again 9 we obtain

From 5, 6 and 19 we find that is uniformly bounded for . We conclude that is uniformly Lipschitz for .

## 3 Convergence to the Hamilton-Jacobi equation

In this section we prove Theorem 1.1.

###### Proof 3

Convergence to the Hamilton-Jacobi equation: From (ii) and (iii) in Theorem 2.2 we have that for , the families are uniformly bounded and Lipschitz. Therefore, from the Arzela-Ascoli Theorem we obtain that, along subsequences, and converge locally uniformly to some continuous functions , with . Moreover, from (i) in Theorem 2.2 we deduce that . Here we consider a subsequence of that converges to .

Let , be the largest eigenvalue of the matrix

and be the corresponding eigenvector. Since the non-diagonal terms in are strictly positive, using the Perron-Frobinius Theorem, we know that such eigenvalue exist and that and are strictly positive. We write

We prove that is a viscosity solution of

To this aim, suppose that has a maximum in . Then, we consider a sequence , such that as , and

is attained at the point and for (The case with can be treated similarly). In this case, we have in particular that

Using the latter and the viscosity criterion for the first equation in 9 we obtain that

(24) |

We notice that, by definition of and , we have

From the latter and by letting in 24 we deduce that

and thus is a subsolution of 10 in the viscosity sense. The supersolution criterion can be proved in a similar way.

## 4 Asymptotic behavior of stationary solutions

In this section we prove Theorem 1.2.

###### Proof 4

Support of : From 19, we deduce that, along subsequences and for , converges weakly to a measure . The fact that is a consequence of the Hopf-Cole transformation 8. To prove 12 it is enough to prove . To this aim following the idea in [25] we first prove that, for , are uniformly semi-convex. Recall that the smooth function is semiconvex with constant , if we have

Let

(25) |

The case where the minimum is achieved for can be treated similarly. We differentiate twice the first equation in 9 with respect to and obtain

From 25 we obtain that , and . Using 7 It follows that

Since , we have . We deduce that

and thus

This proves that , for are semiconvex functions with constant . By passing to the limit in we obtain that is also semiconvex with the same constant.

## 5 Examples of application

In 16–18 we give a description of , assuming that the support of , for , is a set of distinct points, i.e. is a sum of Dirac masses and does not have a continuous distribution. This is what we expect naturally in the models based on darwinian evolution. More precisely, from Volterra-Gause s competitive exclusion principle (see [18, 29]) it is known in theoretical biology that in a model with limiting factors (as nutrients or geographic parameters) at most distinct species can generally survive. Here we have two limiting factors, represented by and , that correspond to the environmental pressures in the two patches. We thus expect to observe only monomorphic or dimorphic situations. This is also the case in the numerical simulations represented in Section 7.

From 12 we know that the support of is included in the set of maximum points of , , with limits of . If now is such that, for fixed , the corresponding set consists of isolated points, it follows that the supports of and consist also of isolated points. We give an example below where has clearly this property.

###### Example 5.1

(monomorphism towards dimorphism) Consider a case with the following values for the parameters of the system

(26) |

with

Then the supports of and consist at most of two single points.

We first notice that in the case where there is no migration between patches (), from the results in [20], we know that in patch , the population concentrates in large time on the maximum points of with the limit of . Since is a quadratic function in , it has a unique maximum and thus is a single Dirac mass on this maximum point. However, allowing migration by taking positive values for and the population can become dimorphic. In Section 7 we give a numerical example where a dimorphic situation appears (see Figure 2). This is in accordance with the competitive exclusion principle since we have introduced a new limiting factor, which is the choice of habitable zones.

Next, we prove the result:

###### Proof 5 (Proof of Example 5.1. )

From 12 we have that the stationary solutions concentrate asymptotically on the maximum points of defined as below

with

(27) |

Since , we deduce that

(28) |

and

(29) |

For fixed , is a polynomial of order . Therefore it has at most two maximum points. It follows that consists of one or two distinct points.

###### Example 5.2

(An asymmetric case) We assume that the parameters are such that the support of , for , consists of isolated points, and we have