Invasion and adaptive evolution for individualbased spatially structured populations
Abstract
The interplay between space and evolution is an important issue in population dynamics, that is in particular crucial in the emergence of polymorphism and spatial patterns. Recently, biological studies suggest that invasion and evolution are closely related. Here we model the interplay between space and evolution starting with an individualbased approach and show the important role of parameter scalings on clustering and invasion. We consider a stochastic discrete model with birth, death, competition, mutation and spatial diffusion, where all the parameters may depend both on the position and on the trait of individuals. The spatial motion is driven by a reflected diffusion in a bounded domain. The interaction is modelled as a trait competition between individuals within a given spatial interaction range. First, we give an algorithmic construction of the process. Next, we obtain large population approximations, as weak solutions of nonlinear reactiondiffusion equations with Neumann’s boundary conditions. As the spatial interaction range is fixed, the nonlinearity is nonlocal. Then, we make the interaction range decrease to zero and prove the convergence to spatially localized nonlinear reactiondiffusion equations, with Neumann’s boundary conditions. Finally, simulations based on the microscopic individualbased model are given, illustrating the strong effects of the spatial interaction range on the emergence of spatial and phenotypic diversity (clustering and polymorphism) and on the interplay between invasion and evolution. The simulations focus on the qualitative differences between local and nonlocal interactions.
MSC 2000 subject classifications: primary 60J85, 60K35, 92D15; secondary 92D25, 35K60.
Key words and phrases. Spatially structured population, adaptive evolution, stochastic individualbased process, birthanddeath point process, reflected diffusion, mutation and selection, nonlinear reactiondiffusion equation, nonlocal and local interactions, clustering and polymorphism, invasion and evolution.
1 Introduction
The spatial aspect of population dynamics is an important ecological issue that has been extensively studied (Murray [31], Durrett and Levin [12], Tilman and Kareiva [44], McGlade [26], Dieckmann et al. [9]). It is in particular crucial in environmental problems, such as spatial invasions and epidemics (Mollison [30], Murray [31], Rand et al. [36], Tilman and Kareiva [44], Lewis and Pacala [25]), and clustering or agglomeration of the population, i.e. the organization as isolated patches (Hassel and May [19], Hassel and Pacala [20], Niwa [32], Flierl et al. [14], Young et al. [45]). The combination of space and phenotype is also known for a long time to have important effects on population dynamics (Mayr [27], Endler [13]). In particular, it can strongly favor the coexistence of several types of individuals and the emergence and stability of polymorphism (Durrett and Levin [11], Dieckmann and Doebeli [7]). More recently, several biological studies (Thomas et al. [43], Phillips et al. [33]) observed that classical models could underestimate the invasion speed and suggested that evolution and invasion are closely related. Namely, the evolution of morphology can have strong impact on the expansion of invading species, such as insect species ([43]) or cane toads ([33]). In this context, the study of spacerelated traits, such as dispersal speed (Pr vost [34], Desvillettes et al. [6]), or sensibility to heterogeneously distributed resources (Bolker and Pacala [3], Grant and Grant [18]), is fundamental.
In this paper, we propose and construct stochastic and deterministic population models describing the interplay between evolution and spatial structure. We show how helpful these models can reveal to understand and predict several specific behaviors concerning clustering and invasion.
We study the dynamics of a spatially structured asexual population with adaptive evolution, in which individuals can move, reproduce with possible phenotypic mutation, or die of natural death or from the competition between individuals. The individuals are characterized both by their position and by one or several phenotypical adaptive traits (such as body size, rate of food intake, age at maturity or dispersal speed). The interaction is modelled as a trait competition between individuals in some spatial range. Our approach is based on a stochastic microscopic description of these individuals’ mechanisms, involving both space and traits. This approach has already been developed in simpler ecological contexts. For populations undergoing dispersal, Bolker and Pacala [2, 3] and Dieckmann and Law [8], offered the first microscopic heuristics and simulations. Their individualbased model has been rigourously constructed in Fournier and Méléard [15]. If one thinks of the dispersion in the physical space as a mutation in a trait space, this model translates into an evolutionary model. The generalization to adaptive population with general mutation and competition phenomena is achieved by Champagnat, Ferrière and Méléard [5, 4]. In these papers, different large population deterministic or stochastic approximations have been obtained, depending on several scalings on the microscopic parameters.
The basic mechanisms of the population dynamics we consider combine spatial motion and evolutionary dynamics (Section 2). The birth, mutation and death parameters of each individual depend on its position and trait. An offspring, appearing at the same position as its progenitor, usually inherits the trait value of the latter, except when a mutation causes the offspring to take an instantaneous mutation step at birth to new trait values. As soon as it is alive, an individual moves in the spatial domain according to a reflected diffusion process. Moreover, each individual dies because of natural death or is eliminated in the competition (selecting the fittest traits) between individuals living in a given spatial range .
Section 3 starts with the algorithmic construction of a stochastic Markov point process whose generator captures the individual migration and ecological mechanisms in the population. Then the existence of this measurevalued process and its martingale properties are proved under some moment condition on the initial data.
Next (Section 4), we study approximations of this model based on largepopulation limits. We consider a large number of individuals at initial time and assume that a fixed amount of available resources has to be partitioned between individuals. When tends to infinity, the conveniently normalized point process converges to a deterministic finite measure, solution of a nonlinear nonlocal integrodifferential equation with Neumann’s boundary conditions, parameterized by the spatial range. The proof is based on the martingale properties of the process and on limit theorems for measurevalued jump processes. We moreover prove that for sufficiently smooth and nondegenerate diffusion coefficients, assuming that the initial condition has a density, the limiting measure has at each time a density with respect to the Lebesgue measure. That is due to the regularizing effect of the reflected diffusion process. The proof mainly uses analytic tools, and is based on the mild formulation of the limiting nonlinear equation.
In Section 5, we study the behavior of this density function as the interaction range tends to . We show its convergence to the solution of a spatially local nonlinear integrodifferential equation with Neumann’s boundary conditions. This equation has been introduced and studied in Prévost [34] in an analytic point of view (see also Desvillettes et al. [6]). In this spatially local case, numerical simulations by finite element methods are given and show the influence of diffusion and mutation parameters on the invasion of the domain by the population.
In Section 6, we give simulations of the microscopic
process illustrating the timedependent interplay between space and
adaptation. We address the effect of the population size, and the
crucial role of the interaction range with respect to spatial
organization (clustering) and polymorphism. We focus on the
qualitative differences between nonlocal and local interactions. In a
first example, we show that, when migrations and mutations are not too
strong, a large interaction range induces a spatial organization of
the population as a finite set of isolated clusters, as assumed in
classical metapopulation models ([11]). Such a spatial
organization is related to the ecological notion of “niches”
(different types of individuals settle different regions of space,
Roughgarden [39]). Conversely, for sufficiently small
interaction range, the clustering phenomenon is no more observed.
Next, we propose another example where a similar phase transition
occurs for spatial clustering and in which the critical interaction
range can be identified. In our last example, we investigate a model
describing the invasion of a species with evolving dispersal speed (as
in [6]). The diffusion coefficient and the trait are
assumed to be proportional and a triangular invasion pattern is
observed, indicating that the invasion front is
composed of faster individuals ([33]).
Notation
The individuals live in the closure of a bounded domain of
of class and their trait values belong
to a compact set of .
 For , we denote by the outward
normal to the boundary at point .
 For a sufficiently smooth function and , we denote by the
scalar product .
 We denote by the space of measurable functions
of class in and bounded in satisfying
for all
and by the subspace of functions which are
moreover continuous in .
 For each , the norm on
is denoted by .
 We denote by the set of finite
measures on , endowed by the weak
topology, and by the subset of composed of all finite point measures, that is
where denotes the Dirac measure at . (If
, one obtains by extension the null measure). For any and for any measurable function
on , we write indifferently or . If , then
.
 We denote by
the Skorohod space of left limited and right continuous functions from
to
, endowed with the Skorohod topology.
 The constant will be a constant which can change from line to
line.
Remark 1.1
Let us remark that the space of functions with a vanishing normal derivative is dense, for the uniform norm, in . Indeed, let us consider the Cauchy problem for the parabolic differential equation with the boundary condition . Since is of class , we may apply SatoUeno [40] Theorem 2.1. There exists a smooth fundamental solution to this problem and each is the uniform limit of the sequence of functions with vanishing normal derivative, as tends to .
We easily extend this result and show that the space is dense in the space of continuous functions on .
2 The model
Let us now describe the evolutionary process we are interested in. The population will be described at any time by a finite point measure . Each individual, characterized by its position and trait , may move, give birth or die, as described below.

The migration is described as a diffusion process normally reflected at the boundary of the domain . Biologists usually assume that the random behavior is isotropic, so the diffusion matrix is chosen with the form (Id is the identity matrix on ) and the nonnegative coefficient (depending on the position and the trait value ), is the diffusion coefficient. We moreover model the environment heterogeneity (resources, topography, external effects,…) by a drift term driven by a vector .

Births and mutations. We consider a population with asexual reproduction. An individual with position and trait can give birth either to a clonal child at rate , or to a mutant with trait at rate , both at position .

The death rate of an individual depends on its position and trait and on the spatial and phenotypic competition with the individuals located around . Let us call the range of this spatial interaction.
For a population , the death rate is given by
The function is assumed to be measurable on .
This interaction assumes that spatial and phenotypic interactions
are independent, which is realistic in many biological situations.
One could of course consider a more complicated interaction. Since
our ultimate goal is to make the spatial
interaction range go to zero, we have chosen this particular form.
Hypotheses (H):
1) The coefficients and
depend Lipschitz continuously on the position and measurably on the
trait, and there exist constants and such that
for all
(2.1) 
2) It is natural from a biological point of view to assume that all birth rates are bounded. There exists such that
(2.2) 
The kernel is nonnegative and symmetric in for each and
(2.3) 
3) There exists a positive constant such that
(2.4) 
4) For each , the spatial kernel is nonnegative and bounded and for each ,
5) The competition kernel is
nonnegative and bounded on .
Let us remark that if is proportional to , then (H4) means that the interaction is proportional to
the surface in around . This is a natural biological
assumption, especially if lies on the boundary of . We
will later assume that the measure weakly converges
to the Dirac measure as tends to .
Hypotheses (H) will be assumed in all the sequel. They imply in particular that for each and each ,
(2.5) 
which yields
(2.6) 
and the total jump rate for a population is bounded by
(2.7) 
We are interested in the evolution of the stochastic point process , taking its values in and describing the evolution of the population at each time . We define
standing for the number of living individuals at time , describing their positions (in ) and their trait values (in ).
The dynamics of the population can be roughly summarized as follows. The initial population is characterized by a measure at time , and any individual located at with trait at time has four independent exponential clocks: a “clonal reproduction” clock with parameter , a “mutant reproduction” clock with parameter , and a “mortality” clock with parameter . If the “mortality” clock of an individual rings, then this individual disappears; if the “clonal reproduction” clock of an individual rings, then it produces at the same location an individual with the same trait as itself; if the “mutant reproduction” clock of an individual rings, then it produces at the same location an individual with characteristics .
The living individuals evolve in the domain, according to diffusion processes with diffusion coefficient and drift , normally reflected at the boundary of .
The measurevalued process is a Markov process whose infinitesimal generator captures this dynamics. This generator is the sum of a jump part corresponding to the phenotypic evolution and of a diffusion part . The generator is defined for bounded and measurable functions from into and for by
(2.8) 
A standard class of cylindrical functions generating the set of bounded and measurable functions from into is the class of functions
(2.9) 
for bounded and measurable functions and .
For such functions , with and , the diffusive part of the generator can easily be deduced from Itô’s formula. Its form is similar to the one obtained in the whole space for branching diffusing processes (cf. RoellyRouault [38]) and is given by
(2.10) 
Hence,
(2.11) 
3 Construction of the particle system and martingale properties
In this section, we construct a Markov process on the path space with infinitesimal generator . Then we prove some martingale properties satisfied by this process, which are the key point to obtain large population approximations.
Let us firstly present an iterative construction of the process, which gives an effective simulation algorithm, if combined with a diffusion simulation step such as an Euler scheme for reflected diffusions (see Lépingle [24], Gobet [17] and Section 6).
The initial number of individuals is equal to some natural integer and the vector of random variables denotes the position and trait values of these individuals. More generally, we denote by the number of individuals at time and by the vector of their positions and traits. Let us introduce the following sequences of independent random variables, independent of .
 are dimensional Brownian motions,
 are uniform random variables on ,
 take values in with law ,
 are exponential random variables with law . (The constant is defined in (2.6)).
The system is obtained inductively for as described below. We set and . Assume that are given. If , then for all . If not, let

. Notice that represents the time between possible jumps for individuals and that gives an upperbound on the total jump rate for a population with individuals, as seen in (2.7).

On the timeinterval , the number of particles remains equal to , their trait values are equal to and their positions evolve according to the following stochastic differential equation with normal reflection : ,
(3.1) 
At time , one chooses at random an individual uniformly among the individuals living during the timeinterval . Its position and trait are .

If dies and . , then the individual

If , then the individual gives birth to an offspring with characteristics and .

If , then the individual gives birth to a mutant offspring with trait at the position , and .

If , nothing happens and .

The total number of individuals at time is equal to
, and .
This stochastic individualbased process can be rigorously expressed as solution of a stochastic differential equation driven by dimensional Brownian motions and the valued multivariate point process
associated with the birth, mutation and death of individuals. We will prove its existence on , deduced from moment properties, and develop some martingale properties that we will use below.
Let us consider . For each function , we define the process as solution of the stochastic differential equation
(3.2) 
where , and have been defined previously.
By Remark 1.1, the knowledge of for is enough to characterized the finite measurevalued process .
We introduce the canonical filtration
where is the Borel field on .
Lemma 3.1
The measure
is the (predictable) compensator of the multivariate point process .
Proof. For , a regular version of the conditional law of with respect to is given by the measure
The conclusion is thus a consequence of [21] Theorem 1.33
p.136.
Using Lemma 3.1 and Itô’s formula, one can immediately show that any solution of (3.2), such that , is a Markov process with infinitesimal generator defined by (2.11). Moreover, we also deduce the following existence, moment and martingale properties.
Proposition 3.2
1) Assume Hypotheses (H) and that .
Then for each and the process defined by (3.2) is well defined on .
2) If furthermore for some , , then for each
Proof. We firstly prove 2). For each integer , define . A simple computation using (3.2), and dropping the nonpositive death terms, gives
Gronwall’s lemma implies that for any , there exists a constant independent of , such that . One easily deduces that tends a.s. to infinity when tends to infinity and next, Fatou’s lemma yields .
Point 1) is a consequence of point 2). Indeed, one builds the
solution step by step. One only has to check
that the sequence of jump instants goes to infinity a.s. as tends to infinity. But this follows from .
The following martingale properties are the key point to study large population approximations.
Theorem 3.3
Assume Hypotheses (H) and that for some , .
1) Then, for and such that for all , , the process
is a càdlàg martingale starting from . It is in particular true for .
2) The process defined for by
(3.3) 
is a càdlàg martingale starting from with predictable quadratic variation
(3.4) 
Proof. Point 1) is immediate. For point 2), we first assume that . Applying point 1) with (or (3.2) and Lemma 3.1) leads to . Then one applies 1) again with , and thus
(3.5) 
is a càdlàg martingale. In another hand, Itô’s formula allows us to compute from (3.3): the process
(3.6) 
is a càdlàg martingale. Comparing (3.5) and (3.6)
leads to (3.4). The extension to the case where
is straightforward, noticing
that .
4 Large population approximation for a fixed interaction range
We are now interested in deterministic approximations of the population point process when the size of the population increases. We assume in this section that the interaction range is fixed.
Let us consider a sequence of initial measures belonging to . For each , we keep all parameters unchanged, except the competition kernel. We assume that for each ,
(4.1) 
This assumption has a biological interpretation. In a case of fixed amount of available global resources, a large system of individuals may only exist if the biomass of each interacting individual scales as , which implies that the interaction effect between two individuals scales as as well. The parameter can also be interpreted as scaling the resources available, so that the renormalization of reflects the decrease of competition for resources.
We assume that the sequence converges, as tends to infinity. The size of the population is then of order and will stay at this order (or at a smaller order) during finite timeintervals, since birth rates are bounded. Hence, our aim is to study the asymptotic behavior, as tends to infinity, of the càdlàg process
(4.2) 
taking values in .
The process is a Markov process with generator . An easy computation, for and , gives that
(4.3) 
and (using (4.1))
(4.4) 
We deduce from Theorem 3.3 the following martingale properties.
Lemma 4.1
Let be fixed and assume that for some , . For all function , the process