![MathJax Logo](/templates/jsp/_style2/_tandf/pb2/images/math-jax.gif)
ABSTRACT
We analyse the asymptotic behaviour of integro-differential equations modelling N populations in interaction, all structured by different traits. Interactions are modelled by non-local terms involving linear combinations of the total number of individuals in each population. These models have already been shown to be suitable for the modelling of drug resistance in cancer, and they generalize the usual Lotka-Volterra ordinary differential equations. Our aim is to give conditions under which there is persistence of all species. Through the analysis of a Lyapunov function, our first main result gives a simple and general condition on the matrix of interactions, together with a convergence rate. The second main result establishes another type of condition in the specific case of mutualistic interactions. When either of these conditions is met, we describe which traits are asymptotically selected.
1. Introduction
1.1. Biological motivations
We are interested in the evolution of N populations of individuals, each of which is structured by a continuous phenotypic trait, also called trait. In each species the phenotype models some continuous biological characteristics (such as the size of the individual, the concentration of a protein inside it, etc). We shall consider both interactions inside a given population and between the populations and we assume that mutations can be neglected. Mathematical modelling and analysis of such ecological scenarios is one purpose of the field of adaptive dynamics, a branch of mathematical biology which aims at describing evolution among a population of individuals, see [Citation15,Citation29,Citation32] for an introduction to deterministic models.
The basis for our model stems from the logistic ODE where r is the net growth rate, dN the logistic death rate due to competition for nutrients and for space by direct or indirect inhibition of proliferation between individuals. Its natural extension to a density
of individuals of phenotype x (say in
) is a non-local logistic model
(1)
(1) with
the total number of individuals.
Following [Citation37], one might consider that these two terms combine both Darwinism selection through the intrinsic growth rate r, and Lamarckism induction through the logistic death rate since it depends on the environment. Note that these models can be derived from stochastic models at the individual level [Citation7,Citation16,Citation21], and more generally measure-valued functions n can be considered [Citation6,Citation22]. The asymptotic behaviour of the previous model (Equation1(1)
(1) ) and variants is analysed in [Citation20,Citation26,Citation32], and one important property among others is that solution typically tend to concentrate on a few phenotypes, a convergence to Dirac masses in mathematical terms. These models are thus successful at representing the survival of only a few phenotypes, which we will refer to as selected.
To account for more complex interactions, one may want to consider a more general non-local logistic term than , in the form
. The behaviour of such equations strongly depends on how localized the kernel is, and therefore so do the mathematical techniques to analyse them. Indeed, with an added diffusion term, a special case of this situation is the non-local Fisher-KPP equation. When the kernel is localized (small as soon as
is large), then the solutions typically remain bounded independently of the mutation rate [Citation23]: selection is no longer a feature. This property highlights how differently the solutions behave depending on the kernel, and that some choices are not appropriate for the ecological situation we are concerned with.
The non-locality actually implies that interaction of an individual of phenotype x with other individuals of phenotype y has the same strength
regardless of y, because individuals do not only necessarily compete with those which have close phenotypes. As such, our choice can serve as a case study to understand the effect of a blind competition across individuals, essentially mediated by the total density.
In [Citation20,Citation26], the biological motivation to use this type of models comes from drug resistance in cancer: the phenotype represents the level of resistance to a given drug and the authors argue that this might be a better approach than a discrete description of the phenotype taking only a finite number of values. Indeed, it can be correlated to some continuous biological characteristics, such as the intracellular concentration of a detoxication molecule, the activity of detoxifying enzymes in metabolizing the administered drug, or drug efflux transporters eliminating the drug [Citation9].
This model is further developed in [Citation26] to incorporate the healthy cell population. Neglecting mutations, it becomes a system of two integro-differential equations of the form
with
,
,
,
the total number of individuals in the cancer cell and healthy cell populations. The interspecific competition (between the two populations) is taken to be competitive, that is
,
, but below the intraspecific competition because each cell population is considered to belong to a different ecological niche:
(2)
(2) In the context of a system arises the central question of persistence (whether asymptotically both species remain), complementing that of identifying which phenotypes are selected. With the additional difficulty of control terms to represent chemotherapeutic drugs, the asymptotic properties of this model are elucidated in [Citation34], and assumption (Equation2
(2)
(2) ) happens to be crucial.
In the aforementioned work, carrying out the asymptotic analysis is instrumental for the application, which is the optimal control problem of minimizing the number of cancer cells through chemotherapy. Simulations suggest that the optimal strategy is to first let the cancer cell density concentrate on a sensitive phenotype, before using the maximum tolerated doses. The convergence to Dirac masses is at the cornerstone of the theoretical analysis of why this strategy happens to be optimal. It is also in sharp contradiction with the classical approach in the clinic, which relies on constant infusion of high doses.
These integro-differential models have therefore already proved their efficiency at helping understanding phenotypic heterogeneity in cancer. The mathematical results available for N=1 and N=2 for competitive interactions naturally call for generalizations on systems of interacting species with such non-local logical terms based on the total number of individuals. For instance, to study resistance in cancer, one may think also of different cancer subpopulations interacting with healthy cells and between them, each one of them being endowed with a specific drug resistance phenotype in a tumour 'bet hedging' strategy [Citation4]. These generalizations, in turn, might help both unravel general principles about the underlying ecological processes, and develop new mathematical techniques to analyse them.
1.2. The model
We consider N populations structured by respective phenotypes , where
is some compact subset of
, with
, for
. Although they model distinct quantities, we abusively denote all variables x to improve readability.
The model writes
(3)
(3) where, for
,
and
are functions in
,
is the total number of individuals in species i, and
are fixed (interaction) coefficients.
The initial conditions are
(4)
(4) where each initial density
is taken to be a non-negative function in
. From now on, we will call these equations integro-differential Lotka-Volterra equations.
The matrix , called matrix of interactions, describes the interactions between the populations: if
, the species j acts positively on the species i, and negatively if
. We will also say that the interaction between species i and j for
is:
mutualistic if
and
,
competitive if
and
,
predator-prey like if
.
Finally, we will say that the equations are competitive (resp., mutualistic) if (resp.,
) for all
.
Another interpretation of the equations is to see them as coupled logistic equations of the form
(5)
(5) In other words, the species i reacts to its environment through the non-local variable
defined for
by
(6)
(6) The terms
and
respectively stand for the net proliferation rate and logistic death rate of individuals in species i, of phenotype x.
We will also use the notation , with which the equations rewrite:
(7)
(7) These models generalize Lotka-Volterra ordinary differential equation (ODE) models [Citation1]: if the functions
,
are all constant (say equal to some
, and
), then after integration with respect to
, the equations boil down to
(8)
(8) which we will from now on refer to as classical N-dimensional Lotka-Volterra equations. Thus, another advantage of a logistic term directly defined by ρ is that it makes our model more tractable with respect to the corresponding already well understood ODE models. Conversely, the integro-differential model can be seen as a perturbation of the ODE one where the individuals among a given population are allowed to have different proliferation and death rates depending on their phenotype x.
Our goal is to understand the asymptotic behaviour of the solutions to these equations, both in terms of convergence at the level of the total number of individuals , and in terms of concentration at the level of the densities
. The first problem is usual in population dynamics while the second is specific to adaptive dynamics and consists of determining which traits asymptotically survive, taking over the whole population. These are then called Evolutionary Attractors, and the fact that it is the generic situation has been coined exclusion principle. Mathematically, this corresponds to a given density
converging to a sum of Dirac masses. For one Dirac mass only, concentration writes, for some
:
(9)
(9) as
, in the weak sense of measures.
The more precise aim of this paper is to study the global asymptotic stability (GAS) of what we will call coexistence steady-states, namely of possible with positive components (all species asymptotically survive) such that ρ converges to
, because we will see how it determines on which phenotypes the densities concentrate. When it is possible, we will investigate the speed at which convergence and concentration occur. An interesting question within the scope of this paper is also to see if a result of that type is sharp, i.e. to compare the assumptions needed to obtain global asymptotic stability in our generalized setting to those known for classical Lotka-Volterra equations.
At this stage, we did not make any restrictive assumptions on the matrix A. However, it will be clear from the results recalled below in the ODE case and the ones presented in Section 2, that answers to the previous questions are available when interspecific interactions are low compared to the intraspecific ones (reminiscent of (Equation2(2)
(2) )). Thus, we are covering the ecological scenario of each species i having its own niche, but inside which competition (if
) is blind because of the term
.
Notations. In what follows, will stand for the positive orthant in
, the set of vectors whose components are all positive, and we will write x>y when
. We will also use the usual ordering on the set of symmetric matrices: for A a real symmetric matrices, we denote
(resp., A>0) when A is positive semidefinite (resp., positive definite). Finally,
will denote the set of Radon measures supported in X.
1.3. State of the art
Classical Lotka-Volterra equations. The ODE system (Equation8(8)
(8) ) has been extensively studied, dating back to the pioneering works of Lotka and Volterra for two populations of preys and predators [Citation28,Citation36]. Since then, many contributions to the analysis of steady-states and their stability have been made, and we refer to [Citation31] for an introduction and to [Citation1] for a review.
Regarding the global asymptotic stability of coexistence steady-states, a very well-known result due to Goh [Citation19] states a simple and very general condition on the matrix which ensures convergence to the (unique) coexistence steady-state:
Theorem 1.1
[Citation19]
Assume that the equation
where
and
are the vectors
and
has a solution
in
. If there exists a diagonal matrix D>0 such that
then
is GAS in
and hence is the unique coexistence steady-state
for system (Equation8
(8)
(8) ).
A result also worth stating is that the mere existence of a unique coexistence steady-state is not enough for it to be GAS. Other steady-states on the boundary of can attract trajectories even in dimension N=2. Another possibility is the occurrence of chaotic behaviour even in low dimension as evidenced in [Citation35] for N=3. Finally, we mention the more recent work [Citation13], where the authors tackle the problem of GAS for some type of Lotka-Volterra ODEs with mutations. In particular, they obtain GAS of the coexistence steady-state in the case where the logistic variables
,
all coincide, that is, when they are equal to some variable
. In such a case, it is proved that the convergence to the equilibrium is exponential. The result of GAS is also extended to perturbations of this specific case.
Integro-differential Lotka-Volterra equations. The first question for such equations is the existence of a solution for all positive times. This obviously does not hold true in full generality since the ODE is a particular case. Let us first state an existence and uniqueness theorem.
Theorem 1.2
Assume that for a given
there exists
such that we have an a priori upper bound
for the functions
whenever they are defined. Then the Cauchy problem (Equation3
(3)
(3) ) and (Equation4
(4)
(4) ) has a unique solution
in
.
The proof is a straightforward generalization of that given in [Citation32, Theorem 2.4] for N=1, relying on the Banach Fixed Point Theorem.
In the case of a single equation, the asymptotic behaviour is well understood. For N=1, assuming to avoid blow-up, the equation is simply
where, without loss of generality, we have set
. The first result is that
converges.
Theorem 1.3
Assume some regularity on
and
Then, for any positive continuous initial condition
the function
is well defined on
and converges to
as
.
This, in turn, completely determines where concentrates.
Corollary 1.1
Under the previous hypotheses, viewed as a Radon measure on
concentrates on the set
as
. If this set is reduced to some
, we obtain in particular
as
in
equipped with its usual weak star topology.
A proof of this result can be found in [Citation34] in the special case , and it relies on proving that
is a bounded variation (BV) function on
. Let us stress that when the set on which
concentrates is not reduced to a singleton, the steady-state (at the level of
) is not unique. For example, if the set is made of two points, the repartition of the limiting density on these two points depends on the initial condition, see for example [Citation12]. This is why for this equation and the general equations considered here, there is no hope in proving general GAS results directly at the level of the densities
.
For a general logistic term and a single equation, the asymptotic behaviour is also analysed in detail in both [Citation14] and [Citation24]. In the latter, under some suitable assumptions on the kernel K, a Lyapunov functional is used to prove that some measure is GAS, in a very specific sense depending on K. Similar results can be found in [Citation8], where their counterpart for competitive classical Lotka-Volterra equations are also discussed.
In the case of integro-differential systems, however, much less is known about the asymptotic behaviour. A Lyapunov functional inspired by [Citation24] has been used successfully in [Citation34] to prove GAS for a competitive system of two populations which writes exactly as (Equation3(3)
(3) ). We also mention [Citation5] where an integro-differential system of two populations is analysed, and whose form does not fit in our framework. A particular triangular coupling structure allows the authors to perform an asymptotic analysis.
The paper is organized as follows. In Section 2, we explain how coexistence steady-states can be identified, allowing us to state rigorously what we mean by GAS for system (Equation3(3)
(3) ). We explain why, under the hypothesis of GAS, only some phenotypes are generically selected, and how to compute them. Then, we present the two main results about GAS for such equations. Section 3 is devoted to the proof of the first result, which applies for any type of interactions and relies on analysing a suitably designed Lyapunov functional. In the specific case of mutualistic interactions, our second main result gives alternative conditions sufficient for GAS. It is presented in Section 4. In Section 5, we conclude with several comments and open questions.
2. Possible coexistence steady-states and main results
For the rest of the article, we will work with the following assumptions:
(10)
(10) This will simplify statements, but we will be more specific below as to which data our results generalize.
2.1. Analysis of coexistence steady-states
In the context of this system of integro-differential equations, the expression ‘GAS in ’ must be defined. By that, we mean that there exists
such that, whatever the positive continuous initial conditions
are,
converges to
for all i.
First, let us explain how to compute the possible steady-states at the level of ρ, i.e. possible limits for positive continuous initial conditions. We will work with the following topological assumption on the sets
:
(11)
(11) where
stands for the Lebesgue measure on
and
for the open ball of centre x and radius η.
Assume that each converges to some
, in which case the exponential behaviour of
is asymptotically governed by
, the sign of which we can analyse as follows.
If this quantity is positive for some
, let us prove that
blows up in its neighbourhood, leading to the explosion of
. If
, there exists
such that by continuity
on
, and then
whether
or also if
thanks to (Equation11
(11)
(11) ). For
small enough and t large enough (say
) such that
, we can write:
with the right-hand side blowing up as
, which cannot be compatible with the convergence of
.
If
is negative globally on
, this clearly implies the extinction of species i, which is also incompatible with the convergence of
to a positive limit.
Remark 2.1
It is possible to relax the regularity on both the sets and the data
and
by working only with points which are both Lebesgue points of
and of Lebesgue density 1 for
, see [Citation17]. If the functions
are in
, one can indeed check that
on
.
The previous results motivate the following definition:
With this definition, a steady-state
exists if and only if the following assumption holds:
(12)
(12) which we assume from now on.
The previous computations also show that vanishes where
, which implies the following result:
Proposition 2.1
Assume that assumption (Equation12(12)
(12) ) holds, and that ρ converges to the coexistence steady-state
. Then,
viewed as a Radon measure, concentrates on the set
as
for all
.
If, for some i, is reduced to some
we obtain in particular
as
in
.
Densities of total mass
and concentrated on
are called Evolutionary Stable Distributions (ESD) in [Citation24].
Remark 2.2
In the hypothesis of global existence and convergence of ρ towards , the previous reasoning actually shows that the concentration is ensured as soon as
is bounded by below by a positive constant on a neighbourhood of one of the points of
. For more general hypotheses ensuring concentration, we refer to [Citation24].
Remark 2.3
If all the sets are reduces to some singletons
, then the dynamics of ρ are asymptotically governed by classical Lotka-Volterra equations concentrated in
, namely
as t goes to
. For a precise statement, see [Citation34].
2.2. Analysis of coexistence steady-states
Our first approach to prove GAS is to further generalize the approach of [Citation34] in dimension N. The main idea is to mix a Lyapunov functional which is inspired by the one designed in [Citation24] and the Lyapunov functional used for classical Lotka-Volterra equations [Citation19], which is the key tool to obtain Theorem 1.1. With some mild regularity assumptions on the data, we obtain the following result:
Theorem 2.1
Assume (Equation12(12)
(12) ) and that there exists a diagonal matrix D>0 such that DA is symmetric and DA<0. Then the solution to the Cauchy problem (Equation3
(3)
(3) ) and (Equation4
(4)
(4) ) is globally defined. Furthermore, the solution
to
is GAS
and hence, unique
.
We emphasize that there is no assumption on the type of interactions, i.e. on the sign of the coefficients of A. However, a consequence of this result is that A must be such that for all i, j. For this result to apply, interactions must therefore be stronger inside species than between them.
We also remark that our hypothesis is exactly the one exhibited in [Citation8] for competitive classical Lotka-Volterra equations. The analysis of the Lyapunov functional allows to determine a speed at which convergence to and concentration on a given set of phenotypes occur. In dimension 2, we also analyse more deeply the link between this condition and the one for classical Lotka-Volterra equations, which in most interesting cases happen to be equivalent.
Our second main result focuses on the special case of mutualistic interactions, and an informal statement of the theorem is the following.
Theorem 2.2
Assume (Equation12(12)
(12) ), that for
and that for some explicit
the matrix
defined by
and
for
is Hurwitz. Then the solution to the Cauchy problem (Equation3
(3)
(3) ) and (Equation4
(4)
(4) ) is globally defined. Furthermore, the solution
to
is GAS.
Again, this applies to the case of interspecific interactions being higher that intraspecific ones, because a Hurwitz matrix is a matrix such that all its eigenvalues have negative real part and it has to do with diagonally dominant matrices (see Section 4).
Because of the hypothesis of mutualism, the system is cooperative, and sub and supersolution techniques can succeed. More precisely, it is possible to prove that all functions are BV on
and this implies their convergence.
3. General interactions
3.1. Proof of the main theorem
In this section, we need slightly more regularity for the data, namely that the functions are Lipschitz continuous:
(13)
(13) We can now restate the first theorem:
Theorem 3.1
Assume (Equation12(12)
(12) ) and (Equation13
(13)
(13) ). Assume that there exists a diagonal matrix D>0 such that DA is symmetric and DA<0. Then the solution to the Cauchy problem (Equation3
(3)
(3) ) and (Equation4
(4)
(4) ) is globally defined.
Furthermore, the solution to
is GAS with
(14)
(14) Concentration of a given
occurs at speed
in the following sense:
(15)
(15) In particular, if
is reduced to a singleton
then
(16)
(16)
Proof.
First step: definition of the Lyapunov functional.
From (Equation12(12)
(12) ), Evolutionary Stable Densities exist and we can pick one of them: for
, we choose any measure
in
satisfying
, which is furthermore concentrated on
, i.e.
(17)
(17) We abusively write integration of functions g against measures μ as
. We also set
and define N functions
by
In what follows, we consider the following Lyapunov functional:
where the positive constants
are to be chosen later. The diagonal matrix of diagonal entries
is denoted by D.
Second step: computation and sign of the derivative.
For any i, we compute
The definition of
implies that the first term simplifies as follows
For the second term, the choice (Equation17
(17)
(17) ) leads to
The functions
are all non-positive by definition of
.
Defining the vector , we arrive at:
Thus, we end up with the expression
(18)
(18) Since the antisymmetric part of DA does not play any role, this can also be expressed
By assumption,
, from which we deduce
. Furthermore, the convergence of the term
to 0 is equivalent to that of ρ to
. However, we do not have the usual property
for Lyapunov functions, so that we cannot yet conclude.
Third step: estimates on .
Let
We are going to show that G is non-decreasing.
We denote by the canonical scalar product of two vectors u and v in
.
For
, the derivative of
is given by
leading to the bound
Put together, these estimates yield:
The last expression is equal to 0 if DA is symmetric, in which case G is non-decreasing as claimed.
The assumptions that DA is symmetric and that are equivalent to the assumption made for the theorem: DA is supposed to be a symmetric negative definite matrix.
As a consequence of the monotonicity of G, we get for all t. The left-hand side is the square of some norm on
, which means that ρ has to be bounded: these a priori bounds ensure the global definition of the solution to (Equation3
(3)
(3) ) and (Equation4
(4)
(4) ).
Fourth step: a lower estimate for V .
To estimate V from below, we need a uniform (in x) upper bound on the densities . Because of the regularity assumption (Equation13
(13)
(13) ), there exists C>0 such that:
The constant C can be chosen to be independent of t since the functions
are bounded.
This implies for a given i
Computing the integral, we write, thanks to the boundedness of
and
and (C has changed and is independent of t and x): for t large enough,
. The bound on V follows immediately:
(19)
(19)
Fifth step: convergence.
G bounds from above:
. Thus
thanks to the third step.
We can now write :
, consequently, each non-positive term it is composed of also vanishes like
as
.
In other words, and each
onverge to 0 as as well
. This is nothing but the two first statements (Equation14
(14)
(14) ) and (Equation15
(15)
(15) ).
For the last statement, we fix i and . We denote
, which is non-negative on
, and by assumption vanishes at
only. We choose a>0 small enough such that
on
. This enables us to write
Remark 3.1
The first interesting fact is that the convergence rate of G to 0, in , is almost optimal in many cases. Indeed, if the sets
are reduced to singletons, there cannot exist any
such that this sum vanishes like
. This comes from the fact that if it were true,
would be integrable on the half-line, which would imply the convergence of V . This is not possible since each
has to go to
as t goes to
.
This might seem contradictory with the exponential convergence rates obtained in [Citation13] for some classical Lotka-Volterra equations, but the Lyapunov functional gives us information on the speed of both phenomena in the sense defined above (through the function G) and it does not say whether one of the two is faster.
3.2. Sharpness in dimension 2
It is clear that if we can find D>0 diagonal such that DA is symmetric and DA<0, then . The condition that DA should be symmetric is constraining, especially if
in which case it imposes some polynomial equalities on the coefficients of the matrix A. In dimension 2, however, the result is sharp in various contexts, as stated in the following proposition.
Proposition 3.1
Assume N=2,
and
. Then the following conditions are equivalent.
there exists D>0 diagonal such that DA is symmetric and
there exists D>0 diagonal such that
.
Proof.
(i) implies (ii) as noticed before.
Now, let us assume (ii) and compute , which has positive determinant, i.e.
. If
,
, a contradiction.
Now, if (iii) holds, we take and
for which
is clearly symmetric negative definite, whence (i).
4. Cooperative case
4.1. Some facts about Hurwitz matrices
We now focus on the cooperative case, i.e. on the case where all off-diagonal elements of A are non-negative. We will also assume that the diagonal elements are negative, since otherwise blow-up clearly occurs: there is intra-spectific competition inside any given species. We will say that such a matrix is cooperative.
In this case, we can hope for stronger results at the level of the integro-differential system because we can use sub and super-solution techniques. For our purpose, the following result on ODEs is sufficient.
Lemma 4.1
For T>0 possibly
let
be a continous function on
locally Lipschitz in
uniformly in
. Denoting
further assume that for all
is non-decreasing with
for all
.
Assume that we have a solution z on of the following Cauchy problem:
(20)
(20) where
and a function y subsolution to the previous Cauchy problem, i.e.
(21)
(21) Then
on
.
When the matrix A is cooperative, it is possible to give an equivalent condition to the one required in Theorem 1.1 for GAS in classical Lotka-Volterra equations. Let us explain how, starting with the three following lemmas, the two first of which can be found in [Citation1].
Lemma 4.2
Let A be a cooperative matrix. Then it is Hurwitz if and only if it is negatively diagonally dominant, i.e. if and only if there exists a vector v>0 such that for
.
This first lemma will be useful in its own right in this section. A consequence is that
Lemma 4.3
If A is cooperative and
has a unique solution in
if and only if A is Hurwitz.
Finally, it comes from the theory of M-matrices (see [Citation33] for a review) that
Lemma 4.4
Let A be cooperative. Then A is Hurwitz if and only if there exists diagonal such that
.
An important consequence of these three lemmas is the following rephrasing of Theorem 1.1 for classical Lotka-Volterra equations.
Proposition 4.1
Assume that A is cooperative, r>0 and that the equations (Equation8(8)
(8) ) have a unique steady-state in
. Then the equations are globally defined and this steady-state is GAS.
Thus, in the context of cooperation between the species, the requirement that A is Hurwitz is somehow optimal to obtain a GAS coexistence steady-state, since it is already required to have its mere existence, a fact mentioned in [Citation18]. We will mainly work with this characterization (rather than the equivalent one given by Lemma 4.4 which we used for a general interaction matrix A) because the next results will lead us to modify the matrix A: analysing whether it is still Hurwitz or not is easier than checking this equivalent condition.
4.2. A priori bounds
For the remaining part of this section, we make the assumption that is positive on
for
, and we define the lower and upper bounds
,
.
Theorem 4.1
Assume that the matrix defined by
and
is Hurwitz. Then the solutions to (Equation3
(3)
(3) ) are globally defined and bounded.
Proof.
First remark that since is Hurwitz, then so is A from Lemma 4.2.
We integrate the equations with respect to x and bound them (through )
Since the diagonal elements are negative, the off-diagonal non-negative, we obtain
Thus, the vector
is a subsolution of the previous system which is nothing but classical Lotka-Volterra equations with interaction matrix
. Thanks to (4.1), the solutions to this system are bounded. Thus, so are those of the integro-differential one.
Remark 4.1
Note that the assumption that is Hurwitz reduces to A being Hurwitz in the case of constant coefficients. Thus, this result for boundedness is sharp, since it is exactly the one required for convergence to the coexistence steady-state when the equations at hand are classical Lotka-Volterra equations.
Using Theorem 4.1, we can thus define as the GAS steady-state for the system obtained in the previous proof, that is to the equations
In other words,
where
is the vector
. This means that we can write
(22)
(22) In a similar fashion to the previous proposition, bounding
away from 0:
leading to
(23)
(23) where
with B, a Hurwitz matrix defined by
,
for
and
.
4.3. GAS in the mutualistic case
We can now state the main result:
Theorem 4.2
Assume for all
and that the matrix
defined by
and
for
is Hurwitz.
Then A and B are also Hurwitz. Furthermore,
lies in
and it is GAS.
Proof.
The fact that , A and B are also Hurwitz is a direct consequence of Lemma 4.2.
being Hurwitz, the solutions are globally defined with ρ bounded thanks to Theorem 4.1.
A being Hurwitz, it is invertible and is in
thanks to Lemma 4.3.
Let us now prove that it is GAS. The idea is to prove that each is BV on
. Identifying the limit is straightforward, thanks to the reasoning made in Section 2.
For any i, we define and write
for readability. Since
, we obtain
Let
. The idea is that
is ‘mostly’ increasing, so we are interested in the negative part of
, denoted by
. For this quantity we have the (a.e.) bound
On the one hand,
On the other hand, for
,
Combining these two, we get
We fix
small enough and t large enough (say
) such that
is Hurwitz (where J is the matrix composed of ones only) and such that, for each
,
. The first requirement is easily derived from Lemma 4.2 since
is clearly cooperative and negatively diagonally dominant for
small enough. The second requirement comes from the lower and upper bounds for the functions
as stated in (Equation22
(22)
(22) ) and (Equation23
(23)
(23) ).
The resulting inequality is
so that
is a subsolution of the system with same initial conditions at
, given by
The solutions to this system go exponentially to 0 since
is Hurwitz.
For any i, we have thus proved that goes to 0 exponentially. Together with the fact that
is bounded from above, we conclude that it is BV on
. Indeed it holds true that a function
which is both bounded from above and such that
is BV on
, see [Citation32, Lemma 6.7]. Therefore, ρ converges (to
).
5. Conclusion
We have analysed the global asymptotic stability properties for integro-differential systems of N species structured by traits x belonging to different trait spaces . The coupling comes from a non-local logistic term, which is a linear combination of the total number of individuals
in each species. Theses systems generalize the usual Lotka-Volterra ODEs for which many stability analyses are available in the literature. Our main focus has been on the asymptotic behaviour of the functions
as
, especially towards equilibrium
with positive components, i.e. of persistence of all species. In Section 2, we explained how identifying it essentially determines the asymptotic behaviour of the underlying density
, namely the phenotypes on which the measures
concentrate in large time.
In Section 3, an adequate Lyapunov functional allowed us to state a general result relying solely on an assumption on the matrix A, regardless of the type of interactions. For N=2, this is essentially a sharp result, but becomes more restrictive for . This tool also provided us with convergence rates to equilibrium. In Section 4, we presented another strategy based on a BV bound which yielded a second result of global asymptotic stability, this time for mutualistic equations.
The result of Theorem 4.2 is partly less general than the one given in Theorem 3.1 because it requires a sign on the coefficients of A. However, the set of matrices which satisfy the hypothesis given in the last theorem is an open subset of the set of real matrices in any dimension. This in sharp contrast with the hypothesis of Theorem 3.1, which, as already mentioned, imposes some polynomial equalities on the coefficients of A as soon as
. In other words, for a small perturbation of a cooperative matrix for which GAS holds, GAS still holds. In particular, if one has weakly (but mutualistically) coupled equations, GAS holds, whereas Theorem 3.1 does not cover the case of any weakly coupled equations for general interactions, except for N=2.
In both cases, the assumptions fall within the class of matrices which cannot have off-diagonal coefficients which are too high compared to the diagonal ones. The present results thus apply to cases where interactions among individuals of a same species are not only blind because of the term , but also stronger than the interactions between species. In other words, each one of them has its own ecological niche inside which interactions are independent of how given phenotypes x and y are away from another.
Let us remark that the BV method would apply to more general functions , as long as they are increasing in the variables
,
. However, the Lyapunov functional used in Theorem 3.1 seems to be dependent on the linear coupling chosen here and it is an open problem to generalize our results for other settings. Another open question is about finding whether there are matrices A for which the underlying classical Lotka-Volterra equations converge to the coexistence steady state (for example such that there exists D>0 with
), but for which there is no GAS for the integro-differential system. Numerically at least, we could not build any such case.
We finally mention a natural extension of these integro-differential systems, which has drawn much attention in adaptive dynamics: that is when one adds a mutation term , and a typical single equation then writes
(24)
(24) This is usually done either through a second order elliptic operator like
with Neumann boundary conditions, or using an integral operator allowing for long-range mutations like
Finally, let us mention that in some more recent works an advection term is also considered [Citation10,Citation11]. The idea is to model stress-induced adaptation: individuals actively adapt to their environment and this can be thought of as an appropriate modelling of Lamarckism induction.
Most of the studies on model (Equation24(24)
(24) ) have aimed at understanding how small mutations affect the dynamics of the surviving phenotype (when one expects a single Dirac mass) with a vanishing mutation term
, after a proper rescaling of time [Citation2,Citation27,Citation30]. Without smallness assumptions on the mutations, non-trivial steady-states have been investigated in detail in [Citation3,Citation12,Citation25], and results of GAS have essentially been obtained for
. It would be interesting to investigate the asymptotic stability properties of steady states for more general kernels. To the best of our knowledge, general asymptotic results are indeed still elusive, even in the case explored in our paper, namely when
for d non constant. It is not clear whether the techniques developed in the present paper can be extended to that setting nor to systems of this form.
Disclosure statement
No potential conflict of interest was reported by the authors.
ORCID
Camille Pouchol http://orcid.org/0000-0002-1152-9896
References
- S. Baigent, Lotka-Volterra Dynamics: an introduction, preprint (2010).
- G. Barles, S. Mirrahimi, and B. Perthame, Concentration in Lotka-Volterra parabolic or integral equations: A general convergence result, Methods Appl. Anal. 16(3) (2009), pp. 321–340.
- O. Bonnefon, J. Coville, and G. Legendre, Concentration phenomenon in some non-local equation, preprint (2015). Available at arXiv:1510.01971.
- B. Brutovsky and D. Horvath, Structure of intratumor heterogeneity: Is cancer hedging its bets? (2013). arxiv, page 1307.0607.
- J.-E. Busse, P. Gwiazda, A. Marciniak-Czochra, Mass concentration in a nonlocal model of clonal selection. J. Math. Biol. 73 (2016), pp. 1–33. doi: 10.1007/s00285-016-0979-3
- J.A. Cañizo, J.A. Carrillo, and S. Cuadrado, Measure solutions for some models in population dynamics, Acta Appl. Math. 123(1) (2013), pp. 141–156. doi: 10.1007/s10440-012-9758-3
- N. Champagnat, R. Ferrière, and S. Méléard, From individual stochastic processes to macroscopic models in adaptive evolution, Stoch. Models 24(S1) (2008), pp. 2–44. doi: 10.1080/15326340802437710
- N. Champagnat, P.-E. Jabin, and G. Raoul, Convergence to equilibrium in competitive Lotka-Volterra and chemostat systems, Comptes Rendus Math. 348(23) (2010), pp. 1267–1272. doi: 10.1016/j.crma.2010.11.001
- R.H. Chisholm, T. Lorenzi, and J. Clairambault, Cell population heterogeneity and evolution towards drug resistance in cancer: Biological and mathematical assessment, theoretical treatment optimization, Biochimic. Biophys. Acta (BBA) - General Subjects 1860(11) (Nov 2016), pp. 2627–2645. doi: 10.1016/j.bbagen.2016.06.009
- R.H. Chisholm, T. Lorenzi, and A. Lorz, Effects of an advection term in nonlocal lotka-volterra equations, Commun. Math. Sci. 14 (2016), pp. 1181–8. doi: 10.4310/CMS.2016.v14.n4.a16
- R.H. Chisholm, T. Lorenzi, A. Lorz, A.K. Larsen, L.N. de Almeida, A. Escargueil, and J. Clairambault, Emergence of drug tolerance in cancer cell populations: An evolutionary outcome of selection, nongenetic instability, and stress-induced adaptation, Cancer Res. 75(6) (2015), pp. 930–939. doi: 10.1158/0008-5472.CAN-14-2103
- J. Coville, Convergence to equilibrium for positive solutions of some mutation-selection model, preprint (2013). Available at arXiv:1308.6471.
- J. Coville and F. Fabre, Convergence to the equilibrium in a Lotka-Volterra ODE competition system with mutations, preprint (2013). Available at arXiv:1301.6237.
- L. Desvillettes, P.E. Jabin, S. Mischler, and G. Raoul, On selection dynamics for continuous structured populations, Commun. Math. Sci. 6(3) (2008), pp. 729–747. doi: 10.4310/CMS.2008.v6.n3.a10
- O. Diekmann, A beginner's guide to adaptive dynamics, Banach Center Publications 63 (2004), pp. 47–86.
- O. Diekmann, M. Gyllenberg, and J. Metz, Steady-state analysis of structured population models, Theor. Popul. Biol. 63(4) (2003), pp. 309–338. doi: 10.1016/S0040-5809(02)00058-8
- L.C. Evans, R.F. Gariepy, Measure Theory and Fine Properties of Functions, CRC press, Boca Raton, 2015.
- B. Goh, Stability in models of mutualism. Am. Nat. 113 (1979), pp. 261–275. doi: 10.1086/283384
- B.S. Goh, Global stability in many-species systems. Am. Nat. 111 (1977), pp. 135–143. doi: 10.1086/283144
- J. Greene, O. Lavi, M.M. Gottesman, and D. Levy, The impact of cell density and mutations in a model of multidrug resistance in solid tumors, Bull. Math. Biol. 76(3) (2014), pp. 627–653. doi: 10.1007/s11538-014-9936-8
- J.M. Greene, D. Levy, K.L. Fung, P.S. Souza, M.M. Gottesman, and O. Lavi, Modeling intrinsic heterogeneity and growth of cancer cells, J. Theor. Biol. 367 (2015), pp. 262–277. doi: 10.1016/j.jtbi.2014.11.017
- M. Gyllenberg and G. Meszéna, On the impossibility of coexistence of infinitely many strategies, J. Math. Biol. 50(2) (2005), pp. 133–160. doi: 10.1007/s00285-004-0283-5
- F. Hamel and L. Ryzhik, On the nonlocal fisher–kpp equation: Steady states, spreading speed and global bounds, Nonlinearity 27(11) (2014), pp. 2735. doi: 10.1088/0951-7715/27/11/2735
- P.-E. Jabin and G. Raoul, On selection dynamics for competitive interactions, J. Math. Biol. 63(3) (2011), pp. 493–517. doi: 10.1007/s00285-010-0370-8
- H. Leman, S. Méléard, and S. Mirrahimi, Influence of a spatial structure on the long time behavior of a competitive Lotka-Volterra type system, Discrete Continuous Dyn. Syst. Ser. B. J. Bridging Math. Sci. 20(2) (2015), pp. 469–493. doi: 10.3934/dcdsb.2015.20.469
- A. Lorz, T. Lorenzi, M.E. Hochberg, J. Clairambault, and B. Perthame, Populational adaptive evolution, chemotherapeutic resistance and multiple anti-cancer therapies, ESAIM: Math. Model. Numer. Anal. 47(2) (2013), pp. 377–399. doi: 10.1051/m2an/2012031
- A. Lorz, S. Mirrahimi, and B. Perthame, Dirac mass dynamics in multidimensional nonlocal parabolic equations, Commun. Part. Diff. Equ. 36(6) (2011), pp. 1071–1098. doi: 10.1080/03605302.2010.538784
- A.J. Lotka, Elements of physical biology, (1924). Reprinted 1956 as elements of mathematical biology.
- J.A. Metz, O. Diekmann, The Dynamics of Physiologically Structured Populations, Vol. 68, Springer, Berlin, 2014.
- S. Mirrahimi and B. Perthame, Asymptotic analysis of a selection model with space, J. Math. Pures Appl. 104(6) (2015), pp. 1108–1118. doi: 10.1016/j.matpur.2015.07.006
- J.D. Murray, Mathematical Biology I: An Introduction, Vol. 17 of Interdisciplinary Applied Mathematics, 2002.
- B. Perthame, Transport Equations in Biology, Springer Science & Business Media, Basel, 2006.
- R.J. Plemmons, M-matrix characterizations. I-nonsingular M-matrices, Linear Algebra Appl. 18(2) (1977), pp. 175–188. doi: 10.1016/0024-3795(77)90073-8
- C. Pouchol, J. Clairambault, A. Lorz, E. Trélat, Asymptotic analysis and optimal control of an integro-differential system modelling healthy and cancer cells exposed to chemotherapy. J. Math. Pures Appl. 116 (2017), pp. 268–308. doi: 10.1016/j.matpur.2017.10.007
- R. Vance, Predation and resource partitioning in one predator – two prey model communities. Am. Nat. 112 (1978), pp. 797–813. doi: 10.1086/283324
- V. Volterra and M. Brelot, Leçons sur la théorie mathématique de la lutte pour la vie, Vol. 1, Gauthier-Villars, Paris, 1931.
- J.X. Zhou, A.O. Pisco, H. Qian, and S. Huang, Nonequilibrium population dynamics of phenotype conversion of cancer cells, PloS one 9(12) (2014), pp. e110714. doi: 10.1371/journal.pone.0110714