Density dependent replicator-mutator models in directed evolution

TL;DR

This paper studies a density-dependent replicator-mutator model in directed evolution, providing a complete mathematical analysis of its solutions, long-term behavior, and a novel socially modulated mutation model exhibiting concentration and diffusion phenomena.

Contribution

It introduces a comprehensive analysis of a density-dependent replicator-mutator model, including explicit solutions and long-term dynamics, and proposes a new socially modulated mutation model.

Findings

Explicit solution expressions derived

Long-term behavior characterized

Socially modulated mutation model shows concentration and diffusion phenomena

Abstract

We analyze a replicator-mutator model arising in the context of directed evolution [23], where the selection term is modulated over time by the mean-fitness. We combine a Cumulant Generating Function approach [13] and a spatio-temporal rescaling related to the Avron-Herbst formula [1] to give of a complete picture of the Cauchy problem. Besides its well-posedness, we provide an implicit/explicit expression of the solution, and analyze its large time behaviour. As a by product, we also solve a replicator-mutator model where the mutation coefficient is socially determined, in the sense that it is modulated by the mean-fitness. The latter model reveals concentration or anti diffusion/diffusion phenomena.

Figures

Equations

$$\frac{\partial u}{\partial t}(t,x)={\sigma }^2\frac{{\partial }^{2}u}{\partial x^2}(t,x)+u(t,x)\frac{x-\overline{x}(t)}{\overline{x}(t)},t>0,x\in \mathbb{R},$$

$$\overline{x}(t):={\int }_{\mathbb{R}}xu(t,x)dx,$$

$$\frac{\partial v}{\partial t}(t,x)={\sigma }^2\overline{x}(t)\frac{{\partial }^{2}v}{\partial x^2}(t,x)+v(t,x)(x-\overline{x}(t)),t>0,x\in \mathbb{R},$$

$${\int }_{\mathbb{R}}{u}_0(x)dx={\int }_{\mathbb{R}}{v}_0(x)dx=1,$$

$$M^{\prime }(t)=1-M(t),M(0)=1,$$

$$M^{\prime }(t)=\overline{v}(t)(1-M(t)),M(0)=1,$$

$$\frac{\partial u}{\partial t}=\underset{\text{mutations}}{\underbrace{{\sigma }^2\frac{{\partial }^{2}u}{\partial x^2}}}+\underset{\text{selection}}{\underbrace{u\left(\mathcal{W}(x)-{\int }_{\mathbb{R}}\mathcal{W}(y)u(t,y)dy\right)}}.$$

$$\frac{\partial u}{\partial t}=J\ast u-u+u\left(x-{\int }_{\mathbb{R}}yu(t,y)dy\right),$$

$$\frac{\partial u}{\partial t}(t,x)=(x-\overline{x}(t))u(t,x),$$

$$u(t+1,x)=\frac{xu(t,x)}{\overline{x}(t)}$$

$$u(t+1,x)-u(t,x)=\frac{x-\overline{x}(t)}{\overline{x}(t)}u(t,x).$$

$$\frac{\partial u}{\partial t}(t,x)=\frac{x-\overline{x}(t)}{\overline{x}(t)}u(t,x),$$

$$u(t+1,x)=\frac{[g(l)x+(1-g(l))\overline{x}(t)]u(t,x)}{\overline{x}(t)}$$

$$u(t+1,x)-u(t,x)=g(l)\frac{x-\overline{x}(t)}{\overline{x}(t)}u(t,x).$$

$${\int }_{\mathbb{R}}f(x)dx=1,$$

$$\underset{x\to \pm \infty }{\lim }f(x)e^{\alpha |x|}=0,\forall \alpha >0,$$

$$m_0:={\int }_{\mathbb{R}}xf(x)dx>0.$$

$$\overline{u}(t)=\sqrt{m_0^2+2{\sigma }^2{t}^2+2{\int }_0^t{C_0}^{\prime \prime }\left({\int }_0^s\frac{dy}{\overline{u}(y)}\right)ds}\text{lower}0.4pt\text{vrule}height=13.88889pt,depth=-11.11116pt$$

$$u(t,x)=w\left(t,x+2{\int }_0^t{\int }_0^s\frac{{\sigma }^2}{\overline{u}(\tau )}d\tau ds\right)\exp \left(-t+x{\int }_0^t\frac{ds}{\overline{u}(s)}+{\int }_0^t{\sigma }^2{\left({\int }_0^s\frac{d\tau }{\overline{u}(\tau )}\right)}^{2}ds\right),$$

$$V(t)$$

$$\overline{u}(t)\{\begin{array}{ll}\ge \sqrt{2{\sigma }^2}\text{lower}0.4pt\text{vrule}height=6.44444pt,depth=-5.15558ptt& \text{ in any case}\\ \sim \sqrt{2{\sigma }^2}\text{lower}0.4pt\text{vrule}height=6.44444pt,depth=-5.15558ptt& \text{ if }\sup \mathrm{supp}(u_0)<+\infty ,\end{array}$$

$${\int }_0^{+\infty }\frac{dy}{\overline{u}(y)}=+\infty .$$

$$V(t)\{\begin{array}{ll}\ge 2{\sigma }^{2}t& \text{ in any case}\\ \sim 2{\sigma }^{2}t& \text{ if }\sup \mathrm{supp}(u_0)<+\infty ,\end{array}$$

$${\psi }^{\prime \prime }(x)-x\psi (x)=0,x\in \mathbb{R}.$$

$$a(t):=\frac{a_0}{1+2a_0{\sigma }^{2}t},m(t):=\{\begin{array}{ll}+\sqrt{2{\sigma }^2{t}^2+\frac{2t}{a_0}+m_0^2}\text{lower}0.4pt\text{vrule}height=7.95523pt,depth=-6.36421pt& \text{ when }{m}_0>0\\ \pm \sqrt{2{\sigma }^2{t}^2+\frac{2t}{a_0}}\text{lower}0.4pt\text{vrule}height=7.95523pt,depth=-6.36421pt& \text{ when }{m}_0=0\\ -\sqrt{2{\sigma }^2{t}^2+\frac{2t}{a_0}+m_0^2}\text{lower}0.4pt\text{vrule}height=7.95523pt,depth=-6.36421pt& \text{ when }{m}_0<0.\end{array}$$

$$u(t,x):=\sqrt{\frac{a(t)}{2\pi }}\text{lower}0.4pt\text{vrule}height=13.22221pt,depth=-10.57782pt{e}^{-a(t)\frac{(x-m(t){)}^2}{2}}$$

$${\partial }_{t}u(t,x)$$

$${\partial }_{xx}u(t,x)$$

$$u(t,x)\frac{x-\overline{x}(t)}{\overline{x}(t)}$$

$$-\frac{a^{\prime }(t)}{2}={\sigma }^2{a}^2(t),$$

$$m^2(t)=2{\sigma }^2{t}^2+\frac{2t}{a_0}+m_0^2,$$

Full text

1. 1 Introduction

2. 2 Main results

3. 3 Special solutions

4. 3.1 Steady states

5. 3.2 Gaussian solutions

6. 4 The CGF approach

7. 5 The solution implicitly/explicitly

8. 6 When the mutation coefficient is socially determined

9. 6.1 Gaussian solutions: concentration vs. extinction

10. 6.2 General solutions

1. Introduction

In this paper, we analyze a mathematical model of a directed evolution process which is density-dependent. The model in consideration is derived in [23] (see below for details) and is given by the following nonlocal equation

[TABLE]

where the nonlocal term is given by

[TABLE]

and will also be denoted $\overline{u}(t)$, to remind its dependence on the solution itself. We will prove the well-posedness of the Cauchy problem associated with (1.1), obtain the solution through an implicit/explicit expression, and analyze its long time behaviour.

As a by product of our analysis of (1.1), we will collect results on the well-posedness and the long time behaviour of the solution to the Cauchy problem associated with equation

[TABLE]

where $\overline{x}(t):=\int_{\mathbb{R}}xv(t,x)\,dx$ will also be denoted $\overline{v}(t)$.

Throughout this work, we assume that the initial data $u_{0}$ and $v_{0}$, are non-negative and have unitary mass

[TABLE]

so that, formally, $\int_{\mathbb{R}}u(t,x)\,dx=\int_{\mathbb{R}}v(t,x)\,dx=1$ for later times. Indeed, if we formally integrate (1.1) over $x\in\mathbb{R}$, we see that the total mass $M(t):=\int_{\mathbb{R}}u(t,x)\,dx$ solves the Cauchy problem

[TABLE]

so that, by the Cauchy-Lipschitz theorem, $M(t)\equiv 1$ for all $t>0$. Hence, $u(t,\cdot)$ is a probability density on $\mathbb{R}$ and $\overline{x}(t)=\overline{u}(t)$ is its mean. As far as (1.3) is concerned, we reach

[TABLE]

so that, by Gronwall’s lemma, $M(t)=1$ as long as “ $\overline{v}(t)$ is meaningful ”. We refer to [1] for situations where, in some sense, $\overline{v}(t)$ blows up in finite time, thus leading to extinction of the solution, which contradicts the formal conservation of the mass.

Let us now comment on the emergence of equations (1.1) and (1.3), termed as replicator-mutator models, in evolutionary biology and their relevance in biotechnology.

Replicator-mutator models aim at describing Darwinian evolutionary processes, whose fundamental elements are mutations and selection. Originally introduced by Taylor and Jonker [22], the replicator dynamics was developed for symmetric games in order to describe the evolution (determined by the payoff matrix) of the frequencies of strategies in a population, see [16] for a complete derivation. Nevertheless, such an approach neglects the effect of mutations. As an attempt to fill this gap, replicator-mutator models have been developed, starting with the work of Kimura [17] and followed by models considering different types of mutations, both in the local [6], [1, 2], [4] and nonlocal [17], [12], [8, 9, 10], [13, 14] cases.

It is important to mention the diversity of research areas where this type of model is applied: population genetics  [15], game theory [7], auto-catalytic reaction networks [21] and language evolution [18]. As pointed out by Schuster and Sigmund [20], in the ordinary differential equation case, several evolutionary models in different biological fields lead independently to the same class of replicator dynamics, showing some sort of universal structure; this idea is also developed in [19], where authors show how apparently very different formulations of evolutionary dynamics are part of a single unified framework given by the replicator-mutator equation.

When mutations are modeled by the local diffusion operator, and under the constrain of mass, the replicator-mutator equation typically takes the form

[TABLE]

In the context of evolutionary biology, $u(t,x)$ represents the density of population, at time $t$, per unit of phenotypic trait on a one-dimensional phenotypic trait space. The function $\mathcal{W}(x)$ stands for the fitness of the phenotype $x$ and models the individual reproductive success. Hence the nonlocal term $\overline{u}(t)=\overline{\mathcal{W}}(t)=\int_{\mathbb{R}}\mathcal{W}(y)u(t,y)\,dy$ is the mean fitness at time $t$, and can be seen as a Lagrange multiplier for the mass to be conserved, thus yielding an equation on the frequencies.

Due to their nonlocal structure, the analysis of replicator-mutator equations often requires new approaches compared with local-pde (e.g. comparison principle does not hold), both from the analytic and numerical viewpoint. In the framework of (1.7), we mention the works [2], [4] for the cases of quadratic or confining fitness functions, but now stick to the case where the fitness linearly depends on the phenotypic trait, say $\mathcal{W}(x)=x$.

In this setting, the authors of [1] proved that the solution of the replicator mutator Cauchy problem (1.7) with linear fitness can be written explicitly in terms of the solution of the Heat equation. On the other hand, when a probability density (or kernel) $J$ models mutations, the equation is recast

[TABLE]

for which a Cumulant Generating Functions (CGF) approach has been developed in [13]: it turns out that the CGF satisfies a first order non-local partial differential equation that can be explicitly solved, thus giving access to many information such as mean fitness (or mean trait since $\mathcal{W}(x)=x$), variance, position of the leading edge. In the present paper, we shall combine these two techniques to first reach an implicit expression of the mean fitness $\overline{u}(t)$ of the solution to (1.1), and then to obtain a so called implicit/explicit expression of the solution $u(t,x)$ itself. Additionally, this procedure allow us to develop a simple numerical strategy for solving the Cauchy problem associated to (1.1).

In this work, and in contrast with (1.7), we consider the case when the fitness function is modulated by $\overline{u}(t)=\overline{x}(t)$, the mean-fitness at time $t$, as can be seen in (1.1). Let us now make a short description of the emergence of model (1.1) in [23] to which we refer for further details.

The original mutator model, yielding the aforementioned replicator-mutator model (1.7), is the continuous time model

[TABLE]

and considers the so-called Malthusian fitness, that is the rate of growth of a particular genotype.

On the other hand, in the case of non-overlapping generations, one may start from the discrete time model

[TABLE]

for the change in the so-called Darwinian fitness (success in propagating genes to the next generation) distribution. The Darwinian fitness being nonnegative, equation (1.9) is supplemented with $u(t=0,\cdot)$ supported in $[0,+\infty)$. Model (1.9) is immediately recast

[TABLE]

Formally, at least for small times and narrow distributions, the above is approached by the continuous in time selection model

[TABLE]

which becomes (1.1) after inserting the effect of mutations through a diffusion term. Notice however that, in this derivation of (1.1), the fitness of an individual with trait $x$ depends only on $x$ and thus selection is not density dependent.

In [23], the authors consider directed evolution, that is a laboratory technique that mimics natural evolution and can be used for example to acquire proteins with new or improved properties. More precisely, they consider a high-throughput compartmentalized directed evolution process. Genotypes inside the compartments have different phenotypes. They not only pool their activity but also share the total number of produced copies, which makes the selection density dependent. In this process, the analogous of (1.9) is given by

[TABLE]

where the constant $g(+\infty)=0<g(l)<1=g(0)$ depends on $l$, a Poisson parameter measuring the occupancy of compartments. The analogous of (1.10) is then

[TABLE]

Compared to (1.10), the presence of the coefficient $g(l)<1$ in (1.13) is the revelator of the density dependent selection. Now, if $g(l)$ is small (meaning $l$ large), the process is slowed down and the validity of (1.11) as a continuous in time approximation should hold in much more cases than previously, meaning for larger time periods but also for more various shapes of distribution.

Hence, the compartmentalization that takes place in directed evolution (but also in natural systems like viruses with multiple infections) and the associated frequency dependence are the key tools to reach the continuous time model (1.1), starting from a discrete time model written in terms of the Darwinian fitness, see [23].

Last, our second main focus, namely equation (1.3), corresponds to the replicator-mutator model (1.7), but with the additional effect that the mutations are frequency-dependent: the diffusion coefficient is modulated by the mean trait $\overline{u}(t)=\overline{x}(t)$ and is thus “socially determined”.

2. Main results

We here state our main results on (1.1), those on (1.3) will be gathered in the final section, where we transfer our developments on the solution $u$ of (1.1) to the solution $v$ of (1.3).

We first need to define an admissible set where to look after solutions. We denote by $\mathcal{A}$, the set of non-negative functions $f\in L^{1}(\mathbb{R})$ such that

[TABLE]

which decrease faster than any exponential function, that is

[TABLE]

and, last, such that

[TABLE]

Remark* 2.1**.*

Notice that the case $m_{0}<0$ follows by symmetry from the case $m_{0}>0$, whereas the case $m_{0}=0$ is singular, as clear from the equation and the Gaussian case studied in subsection 3.2. Notice also that assumption (2.1) could be relaxed by only assuming limit $x\to+\infty$ when $m_{0}>0$ (or $x\to-\infty$ when $m_{0}<0$), the relevant tail being the right one as already observed in related situations by [1] or [13]. This would require to adapt $(iii)$ in the below definition by adding some integrability properties as $x\to-\infty$. For simplicity, in this work, we stick to (2.1).

Definition 2.2** (Notion of solution).**

Let $u_{0}\in\mathcal{A}$ be given. We say that $u=u(t,x)$ is a global solution of (1.1) starting from $u_{0}$ if

1. $(i)$

$u\in\mathscr{C}^{\infty}((0,+\infty)\times\mathbb{R})$, 2. $(ii)$

for all $t\geq 0$, $u(t,\cdot)\in\mathcal{A}$, 3. $(iii)$

for all $t>0$, $\partial_{t}u(t,\cdot)$, $\partial_{x}u(t,\cdot)$ and $\partial_{xx}u(t,\cdot)$ decrease faster than any exponential function, in the sense of (2.1), 4. $(iv)$

$u$* solves (1.1) in the classical sense,* 5. $(v)$

$u(t,\cdot)\to u_{0}$* in $L^{1}(\mathbb{R})$, as $t\to 0$.*

Our main result consists in the well-posedness of the Cauchy problem (1.1) with, moreover, an implicit/explicit expression of the solution.

Theorem 2.3** (The solution of the Cauchy problem).**

Let $u_{0}\in\mathcal{A}$ be given. Then there is a unique global solution of (1.1) starting from $u_{0}$ (in the sense of Definition 2.2). Moreover, its mean $\overline{u}(t)$ is implicitly (and uniquely) determined by

[TABLE]

where $m_{0}=\int_{\mathbb{R}}xu_{0}(x)dx>0$ and $C_{0}(z):=\ln\left(\int_{\mathbb{R}}u_{0}(x)e^{zx}dx\right)$, $z\geq 0$, is the cumulant generating function (cgf) of $u_{0}$. Last, $u(t,x)$ is given by

[TABLE]

where $w=w(t,y)$ is the solution of the Heat equation $\partial_{t}w=\sigma^{2}\partial_{yy}w$ starting from $u_{0}$.

The proof is based in a combination of the approaches of [13] and [1]. In the course of the proof, we collect the expressions and some estimates on the mean $\bar{u}(t)$ and the variance

[TABLE]

of the solution at time $t$.

Corollary 2.4** (Long time behaviour).**

Let $u_{0}\in\mathcal{A}$ be given. Then the mean $\overline{u}(t)$ of the global solution given by Theorem 2.3 satisfies

[TABLE]

and, in any case,

[TABLE]

The variance $V(t)$ of the global solution satisfies

[TABLE]

where ${\rm supp}(u_{0})$ denotes the support of $u_{0}$.

The above shows that, as $t\to+\infty$, the solution moves to the right since $\overline{u}(t)\to+\infty$, and is flattening since $V(t)\to+\infty$. In particular, for one side compactly supported initial data, the propagation (asymptotically) occurs at constant speed $\textstyle\sqrt{2\sigma^{2}\,}$ , which is in contrast with the acceleration phenomenon proved in [1]. This is also true for Gaussian data as will be noticed in subsection 3.2. The role of (2.2) is to provide “a kind of upper bound” on $\overline{u}(t)$ when, for example, the initial tails are still lighter than any exponential but heavier than Gaussian, e.g. of the magnitude $e^{-x^{\alpha}}$ ($1<\alpha<2$) or even $e^{-x\ln x}$ as $x\to+\infty$.

The paper is organized as follows. In Section 3, we investigate special solutions, in particular Gaussian ones which provide a preliminary understanding of equation (1.1). In Section 4, we begin the proof of our main result, Theorem 2.3, by using the cgf approach introduced in [13]. The proof of Theorem 2.3 and its corollary are completed in Section 5 thanks to the methodology of [1]. As a by product of our analysis, we collect a numerical strategy for solving the Cauchy problem. Last, Section 6 is devoted to obtain the companion results on (1.3) by expressing $v(t,x)$ in terms of $u(t,x)$.

3. Special solutions

In this section, we investigate special solutions to (1.1), in particular non-negative, integrable steady states and Gaussian solutions. In the first case we prove the non existence of such steady states. This is due to the particular form of the fitness function and is an analogous result to that one obtained by Alfaro an Carles in [1]. The situation is different in the case of a confining fitness function [4, 2] where the existence and uniqueness of a non-negative stationary solution is ensured and corresponds to the ground-state of the Schrödinger Hamiltonian where the potential is the opposite of fitness function. In the second case, Gaussian solutions are computed explicitly by solving the differential equations describing the evolution of the mean and the inverse of the variance.

3.1. Steady states

Proposition 3.1** (Steady state).**

There is no nontrivial non-negative steady state $\phi=\phi(x)$ solving (1.1) and satisfying $\phi(\pm\infty)=0$.

Proof.

A steady state with $\overline{x}\neq 0$ must solve $\phi^{\prime\prime}(x)+\frac{1}{\sigma^{2}}\frac{x-\overline{x}}{\overline{x}}\phi(x)=0$ for all $x\in\mathbb{R}$. Letting $\psi(x):=\phi(\overline{x}-{(\sigma^{2}\overline{x})}^{1/3}\,x)$, this is recast

[TABLE]

Hence $\psi(x)$ is a linear combination of the Airy functions ${\rm Ai}(x)$ and ${\rm Bi}(x)$. From $\psi(+\infty)=0$, we deduce that $\psi(x)$ is a multiple of ${\rm Ai}(x)$. Hence either it is trivial, or it changes sign. ∎

3.2. Gaussian solutions

We investigate the propagation of a Gaussian initial data, which is relevant for biological applications. In the context of evolutionary genetics, families of Gaussian solutions for nonlinear and nonlocal equations can be found in [6], [1, 2]. In a different context involving logarithmic non-linearities, we also refer to [5], [11] for the Schrödinger equation and to [3] for the Heat equation.

Proposition 3.2** (Propagation of Gaussian initial data).**

For $a_{0}>0$ and $m_{0}\in\mathbb{R}$, let us define

[TABLE]

Then

[TABLE]

solves (1.1) in $(0,+\infty)\times\mathbb{R}$ and starts from u_{0}(x)=\mathchoice{{\hbox{\displaystyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=9.3322pt,depth=-7.46579pt}}}{{\hbox{\textstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=6.55832pt,depth=-5.24669pt}}}{{\hbox{\scriptstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=5.05275pt,depth=-4.04222pt}}}{{\hbox{\scriptscriptstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=5.05275pt,depth=-4.04222pt}}}e^{-a_{0}\frac{(x-m_{0})^{2}}{2}}.

Proof.

From the ansatz (3.2), we compute

[TABLE]

We plug the above into equation (1.1), and identify the $x^{2}$, the $x^{1}$ and the $x^{0}$ coefficients to obtain three differential equations. The first one is

[TABLE]

whose solution, starting from $a_{0}$, is given by the first part in (3.1). Using (3.3), we see that the two other equations reduce to $a(t)m(t)m^{\prime}(t)=1$, which is solved as

[TABLE]

and thus the second part in (3.1). ∎

Notice that functions $m(t)$ and $V(t)=\frac{1}{a(t)}$ are respectively the mean and the variance of the density $u(t,\cdot)$. Hence Proposition 3.2 shows that, starting from a Gaussian profile, the solution remains a Gaussian function, is asymptotically propagating at constant speed and flattening since m(t)\sim\pm\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=6.44444pt,depth=-5.15558pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=6.44444pt,depth=-5.15558pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.44165pt,depth=-2.75334pt}}}\,t, $V(t)\sim 2\sigma^{2}t$, as $t\to+\infty$, see Figure 1. Notice that the direction of propagation is determined by the initial value $m_{0}$: towards the right when $m_{0}>0$, the left when $m_{0}<0$. The case $m_{0}=0$ is singular because of the multiplicity of solutions to (3.1), two Gaussian solutions emerge, propagating to the left and to the right, see Figure 1.

4. The CGF approach

In this section, we assume that we are equipped with a solution of (1.1) starting from $u_{0}\in\mathcal{A}$. We define the Cumulative Generating Function (cgf)

[TABLE]

of such a solution, in the spirit of [13]. From Definition 2.2, $C(t,z)$ is well defined and smooth on $(0,+\infty)\times[0,+\infty)$. We shall derive an implicit expression for $C(t,z)$, and then for the value $\overline{u}(t)$ of the mean. This crucial step will enable us to complete the analysis in the next section, which is much in the spirit of [1].

The following computations are validated by the notion of solution adopted in Definition 2.2 and, possibly, the dominated convergence theorem. Observe that

[TABLE]

and that

[TABLE]

Notice in particular that the mean is reached through

[TABLE]

whereas the variance is reached through

[TABLE]

Hence, multiplying equation (1.1) by $e^{zx}$, integrating over $x\in\mathbb{R}$ and finally dividing by $\int_{\mathbb{R}}u(t,x)e^{zx}$, we obtain the following nonlocal first order pde

[TABLE]

where we have used integration by parts to get

[TABLE]

Furthermore, since $\int_{\mathbb{R}}u(t,x)dx=1$, the condition $C(t,0)=0$ must hold for any $t\geq 0$. As a result we are facing the problem

[TABLE]

where $C_{0}(z)=\ln\left(\int_{\mathbb{R}}u_{0}(x)e^{zx}dx\right)$ is the cgf of the initial data $u_{0}$. Notice that, from a straightforward computation,

[TABLE]

Hence, from Cauchy Schwarz inequality, we see that ${C_{0}}^{\prime\prime}(z)\geq 0$ for all $z\geq 0$, and even ${C_{0}}^{\prime\prime}(z)>0$ for all $z\geq 0$ if $u_{0}\not\equiv 0$. This convexity property of cumulant generating functions is well-known in probability, and will be used in the following.

Fix $t>0$ and $z>0$. For $s$ such that $-t\leq s$ and $\int_{0}^{s}\frac{d\tau}{\overline{u}(t+\tau)}\leq z$, set

[TABLE]

which we differentiate to get

[TABLE]

As a result,

[TABLE]

From (4.4), we deduce that

[TABLE]

by Fubini-Tonelli theorem. We thus conclude that

[TABLE]

This is an implicit expression for $C(t,z)$ since it still involves $\overline{u}(t)=\partial_{z}C(t,0)$. Differentiating with respect to $z$ and evaluating at $z=0$, we reach another implicit formula for the mean value of the solution

[TABLE]

whereas differentiating twice with respect to $z$ and evaluating at $z=0$, we reach the variance of the solution

[TABLE]

where the last estimate follows from [13, Lemma 4.5].

Next, differentiating (4.5) with respect to $t$, we get

[TABLE]

so that $\frac{d}{dt}\overline{u}^{2}(t)=2V(t)$ and thus

[TABLE]

Also, integrating $\frac{d}{dt}\overline{u}^{2}(t)=2V(t)$, we collect the very useful expression

[TABLE]

Remark* 4.1**.*

Notice that if u_{0}(x)=\mathchoice{{\hbox{\displaystyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=9.3322pt,depth=-7.46579pt}}}{{\hbox{\textstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=6.55832pt,depth=-5.24669pt}}}{{\hbox{\scriptstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=5.05275pt,depth=-4.04222pt}}}{{\hbox{\scriptscriptstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=5.05275pt,depth=-4.04222pt}}}e^{-a_{0}\frac{(x-m_{0})^{2}}{2}} is a Gaussian initial data as in subsection 3.2, we know that its cgf is

[TABLE]

Hence \overline{u}(t)=\mathchoice{{\hbox{\displaystyle\sqrt{m_{0}^{2}+2\sigma^{2}t^{2}+\frac{2}{a_{0}}t,}}\lower 0.4pt\hbox{\vrule height=9.49944pt,depth=-7.59958pt}}}{{\hbox{\textstyle\sqrt{m_{0}^{2}+2\sigma^{2}t^{2}+\frac{2}{a_{0}}t,}}\lower 0.4pt\hbox{\vrule height=7.95523pt,depth=-6.36421pt}}}{{\hbox{\scriptstyle\sqrt{m_{0}^{2}+2\sigma^{2}t^{2}+\frac{2}{a_{0}}t,}}\lower 0.4pt\hbox{\vrule height=5.59444pt,depth=-4.47557pt}}}{{\hbox{\scriptscriptstyle\sqrt{m_{0}^{2}+2\sigma^{2}t^{2}+\frac{2}{a_{0}}t,}}\lower 0.4pt\hbox{\vrule height=4.94304pt,depth=-3.95445pt}}} and $V(t)=\frac{1}{a_{0}}+2\sigma^{2}t$, which agrees with the values of $m(t)$ and $\frac{1}{a(t)}$ in Proposition 3.2.

5. The solution implicitly/explicitly

In this section, we mainly complete the proof of Theorem 2.3 and of Corollary 2.4. As a by product, we present some numerical strategies for solving (1.1).

First, assume that we are equipped with a solution $u(t,x)$ of (1.1) starting from $u_{0}\in\mathcal{A}$. In particular we know from Section 4 that $\overline{u}(t)$ is given by (4.8). Following [1], we write

[TABLE]

This clearly defines a unique $w(t,y)$, $t\geq 0$, $y\in\mathbb{R}$, which has to satisfy

[TABLE]

We refer to [1] for more details on such a change of unknown function, which enables to reduce some replicator-mutator equations to the Heat equation.

Now, our main task is to show that the problem (4.8), is globally well-posed.

Lemma 5.1**.**

For any given $T>0$, there is a unique solution $\overline{u}\in\mathscr{C}^{0}([0,T])$ to

[TABLE]

Proof.

Recall that $m_{0}>0$. We define the following subspace of $\mathscr{C}^{0}([0,T])$, equipped with the $L^{\infty}$ norm,

[TABLE]

which is a Banach space. Now, we define the operator $\mathcal{L}:h\in X\mapsto q\in X$, and $q$ is given by

[TABLE]

The proof consists in mimicking that of the Banach fixed point theorem. In order to prove existence, we introduce the sequence $(q_{n})\in X^{\mathbb{N}}$ defined inductively by

[TABLE]

We have

[TABLE]

where $k:=\frac{1}{m_{0}^{3}}{\lVert{{C_{0}}^{(3)}}\rVert}_{L^{\infty}(0,T)}$. Then we straightforwardly prove by induction that, for any $n\geq 0$, any $0\leq s\leq T$,

[TABLE]

In particular, the series $\sum{{\left\lVert q_{n+1}-q_{n}\right\rVert}}_{L^{\infty}(0,T)}$ converges and therefore $(q_{n})$ converges uniformly in $\mathscr{C}^{0}(0,T)$ to some $q$ which is a fixed point of $\mathcal{L}$.

As far as uniqueness is concerned, if $q$ and $\tilde{q}$ are two fixed points of $\mathcal{L}$ then the argument to reach (5.5) also yields $\|q-\tilde{q}\|_{L^{\infty}(0,T)}\leq\|q-\tilde{q}\|_{L^{\infty}(0,T)}\frac{(kT^{2})^{n}}{n!}$. Letting $n\to+\infty$ enforces $q\equiv\tilde{q}$. ∎

Completion of the proof of Theorem 2.3 and Corollary 2.4.

Hence, any solution in the sense of Definition 2.2 has to be given by (5.1), where $\overline{u}(t)$ is given by Lemma 5.1. Notice that, from (5.3), $\overline{u}$ has to be smooth.

Conversely, it is now a matter of straightforward computations — based on (5.1),(5.2) and (5.3)— to check that this does provide the solution. Notice that, in particular, the initial data is understood in the sense of Definition 2.2 $(v)$ since, in the Cauchy problem for the heat equation (5.2), the initial data is understood in the sense $w(t,\cdot)\to u_{0}$ in $L^{1}(\mathbb{R})$ as $t\to 0$. This proves Theorem 2.3.

Last, the conclusions of Corollary 2.4 have been collected in the course of Section 4, except estimate (2.2) which will be obtained in the proof of Lemma 6.3, where its relevance will become much clearer. ∎

Numerical implications. Three major difficulties are encountered when setting up a numerical strategy for problem (1.1). The first one is the nonlocal nature of the equation. At this stage, two natural approaches can be considered: $i)$ use finite differences and approximate the nonlocal term by a Riemann sum; $ii)$ apply a splitting method, computing alternatively the solution of the non-local term and the resulting local reaction-diffusion equation.

The second difficulty is the unboundedness of the domain. To manage it, one usually imposes artificial boundary conditions and solves numerically the equation on a sufficiently large domain, approximating the true solution by the latter.

The third complexity comes from the propagation of the solution, at least linear in view of (4.7), making it difficult to track it over time.

Previously, we gave an implicit/explicit construction of the solution $u=u(t,x)$ of the Cauchy problem associated with equation (1.1), which actually provides a strategy for solving numerically the problem.

The idea is to initially find an approximation of the nonlocal term $\overline{u}(t)$ on a time interval $[0,T],\ T>0$. This step is performed using the fixed-point iteration sequence (5.4), for which the estimate (5.5) provides a convergence rate. This requires in particular to compute (analytically or numerically) the cumulant generating function $C_{0}(z)$ of the initial datum $u_{0}(x)$.

The next step consists in calculating the solution $w=w(t,x)$ of the heat equation with initial datum $u_{0}(x)$ and evaluating it in the values indicated by formula (5.1). Obviously, in the case where an explicit solution $w=w(t,x)$ of the heat equation can be found, artificial boundary conditions for solving the heat equation are not needed, leading to a better numerical solution. On the other hand, in the case where no closed form for $w=w(t,x)$ is available, it is important to previously solve numerically the heat equation in a sufficiently large spatial domain. For instance, if we want to solve numerically (1.1) in $[x_{-},x_{+}]\times[0,T]\subset\mathbb{R}\times\mathbb{R}$, assuming $m_{0}>0$ (the case $m_{0}<0$ following by symmetry), then it is clear from (5.1) that the function $w$ needs to be evaluated in space up to $x=x_{+}+2\int_{0}^{t}\int_{0}^{s}\frac{\sigma^{2}}{\overline{u}(\tau)}d\tau ds$ and since $\overline{u}(t)\geq m_{0}>0$, then we should solve the heat equation at most in $[x_{-},x_{+}+\tfrac{\sigma^{2}T^{2}}{m_{0}}]\times[0,T]$.

In figure 2(b), we have represented the numerical solution obtained by this method for the initial condition $u_{0}(x)=\mathds{1}_{[\tfrac{1}{2},\tfrac{3}{2}]}(x)$, so that $m_{0}=1$ and $C_{0}(z)=\ln\left(\frac{2e^{z}\sinh(z/2)}{z}\right)$, and $\sigma^{2}=1$; in this case, $w$ is known to be given by

[TABLE]

6. When the mutation coefficient is socially determined

This section is devoted to the analysis of equation (1.3) supplemented with a non-negative initial data $v_{0}\in L^{1}(\mathbb{R})$ satisfying $\int_{\mathbb{R}}v_{0}(x)\,dx=1$, and decreasing faster than any exponential data in the sense of (2.1). We denote $m_{0}:=\int_{\mathbb{R}}xv_{0}(x)\,dx$. Interestingly, one can avoid the assumption $m_{0}>0$ in subsection 6.1 and obtain anti-diffusing/diffusing solutions that are mathematically interesting.

6.1. Gaussian solutions: concentration vs. extinction

Proposition 6.1** (Propagation of Gaussian initial data).**

For $a_{0}>0$ and $m_{0}\in\mathbb{R}$. Then there is a unique Gaussian

[TABLE]

which solves, at least locally in time, (1.3) and starts from u_{0}(x)=\mathchoice{{\hbox{\displaystyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=9.3322pt,depth=-7.46579pt}}}{{\hbox{\textstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=6.55832pt,depth=-5.24669pt}}}{{\hbox{\scriptstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=5.05275pt,depth=-4.04222pt}}}{{\hbox{\scriptscriptstyle\sqrt{\frac{a_{0}}{2\pi},}}\lower 0.4pt\hbox{\vrule height=5.05275pt,depth=-4.04222pt}}}e^{-a_{0}\frac{(x-m_{0})^{2}}{2}}. Its behaviour depends strongly on the parameters $a_{0}>0$, $m_{0}\in\mathbb{R}$ and $\sigma>0$.

1. $(i)$

Assume m_{0}<-\frac{1}{a_{0}\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}}. Then $a(t)$ blows up in finite time: there is $0<T^{\star}<+\infty$ such that

[TABLE] 2. $(ii)$

Assume m_{0}=-\frac{1}{a_{0}\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}}. Then $a(t)$ and $m(t)$ are global and

[TABLE] 3. $(iii)$

Assume m_{0}>-\frac{1}{a_{0}\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}}. Then $a(t)$ and $m(t)$ are global and

[TABLE]

where C_{m}:=\frac{2a_{0}m_{0}\mathchoice{{\hbox{\displaystyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.277pt,depth=-3.42162pt}}}{{\hbox{\textstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.277pt,depth=-3.42162pt}}}{{\hbox{\scriptstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.01193pt,depth=-2.40956pt}}}{{\hbox{\scriptscriptstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}+\mathchoice{{\hbox{\displaystyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=2.25555pt,depth=-1.80446pt}}}}{4a_{0}\mathchoice{{\hbox{\displaystyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.277pt,depth=-3.42162pt}}}{{\hbox{\textstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.277pt,depth=-3.42162pt}}}{{\hbox{\scriptstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.01193pt,depth=-2.40956pt}}}{{\hbox{\scriptscriptstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}}>0 and C_{V}:=\frac{2a_{0}m_{0}\mathchoice{{\hbox{\displaystyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.277pt,depth=-3.42162pt}}}{{\hbox{\textstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.277pt,depth=-3.42162pt}}}{{\hbox{\scriptstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.01193pt,depth=-2.40956pt}}}{{\hbox{\scriptscriptstyle\sqrt{\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}+\mathchoice{{\hbox{\displaystyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=2.25555pt,depth=-1.80446pt}}}}{2\mathchoice{{\hbox{\displaystyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2,}}\lower 0.4pt\hbox{\vrule height=2.25555pt,depth=-1.80446pt}}}a_{0}}>0.

Remark* 6.2**.*

In Figure 3, we have subdivided the phase plane $(m_{0},a_{0})$ into regions corresponding to the cases $(i),(ii)$ and $(iii)$.

Case $(i)$ corresponds to a concentration phenomena in finite time: the Gaussian solution converges, in finite time, to a Dirac mass centered at $x=m(T^{\star})<0$. See Figure 4.

Case $(ii)$ corresponds to a concentration phenomena in infinite time: the Gaussian solution converges, at large times, to a Dirac mass centered at $x=0$. See Figure 5.

In case $(iii)$, if $m_{0}<0$ then, in contrast with cases $(i)$ and $(ii)$, the mean of the solution manages to “cross” the value zero. In case $(iii)$ convergence to a accelerating (m(t)\sim C_{m}e^{\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}t}) and flattening (V(t)\sim C_{V}e^{\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}t}) profile always occurs at large times. Moreover, the variance of the solution is decreasing before the mean reaches zero, and then starts to increase after the mean crosses zero, which reveals an anti-diffusion/diffusion phenomenon, see Figure 6. On the other hand, if $m_{0}\geq 0$, the solution does not anti-diffuse and only diffuses while the mean tends to infinity, see Figure 7.

Proof.

As in the proof of Proposition 3.2, we plug the ansatz (6.1) into (1.3) and arrive at

[TABLE]

so that $m^{\prime\prime}(t)=2\sigma^{2}m(t)$, $m(0)=m_{0}$, $m^{\prime}(0)=\frac{1}{a_{0}}$ which is globally solved as

[TABLE]

From the first equation in (6.2) we deduce

[TABLE]

as long as blow-up does not occur. We now easily see that in case $(i)$ blow up occurs at time

[TABLE]

In case $(ii)$ blow up does not occur, m(t)=m_{0}e^{-\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}t} and a(t)=a_{0}e^{\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}t}. In case $(iii)$ blow up does not occur and conclusion follows from (6.3) and (6.4). ∎

6.2. General solutions

Here we consider a non-negative initial data $v_{0}\in L^{1}(\mathbb{R})$ satisfying $\int_{\mathbb{R}}v_{0}(x)\,dx=1$, (2.1) and $m_{0}:=\int_{\mathbb{R}}xv_{0}(x)\,dx>0$, and investigate possible solutions $v=v(t,x)$ of (1.3) starting from $v_{0}$. From Theorem 2.3, we are equipped with the unique solution $u=u(t,x)$ of equation (1.1) starting from $v_{0}$. Roughly speaking, the understanding of (1.3) now reduces to that of the o.d.e. Cauchy problem

[TABLE]

Indeed, defining

[TABLE]

we have $\overline{v}(t)=\overline{u}(\varphi(t))$ and one checks that, since $u$ solves (1.1) (and starts from $v_{0}$), $v=v(t,x)$ solves (1.3) (and starts from $v_{0}$) as long as the solution $\varphi(t)$ of (6.5) exists.

Lemma 6.3**.**

The solution $\varphi(t)$ of the o.d.e. Cauchy problem (6.5) is global.

Proof.

If $\varphi$ blows up in finite time, say at time $T>0$, then we obtain from the o.d.e. that

[TABLE]

Hence we deduce from (4.5) that

[TABLE]

But, from (4.8), we know that \frac{y}{\overline{u}(y)}\leq\frac{1}{\mathchoice{{\hbox{\displaystyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\textstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=4.51111pt,depth=-3.6089pt}}}{{\hbox{\scriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=3.15778pt,depth=-2.52625pt}}}{{\hbox{\scriptscriptstyle\sqrt{2\sigma^{2},}}\lower 0.4pt\hbox{\vrule height=2.40916pt,depth=-1.92735pt}}}} and thus

[TABLE]

which contradicts $\frac{1}{\overline{u}}\in L^{1}(0,+\infty)$. This concludes the proof of the lemma and, clearly, of estimate (2.2). ∎

Since $\overline{u}$ is smooth, so is $\varphi$. Also, $\varphi^{\prime}(t)\geq m_{0}>0$ and, in particular, $\varphi(t)$ tends to $+\infty$ as $t\to+\infty$. Hence $\varphi:(0,+\infty)\rightarrow(0,+\infty)$ is a smooth diffeomorphism. In other words, (6.6) shows that there is a one-to-one correspondence between solutions to (1.3) and that of (1.1). We omit the full details but this clearly concludes the analysis of the Cauchy problem associated with (1.3), whose surprising qualitative behaviours have already been underlined in subsection 6.1.

Acknowledgements. The authors are very grateful to Anton Zadorin for sharing his expertise on the models appearing in [23], and for a precise relecture of the preprint, in particular the introduction. M. Alfaro is supported by the ANR i-site muse, project michel 170544IA (n*∘* ANR idex-0006). M. Veruete is supported by the National Council for Science and Technology of Mexico.