Adaptive Fitness Landscape for Replicator Systems: To Maximize or not to Maximize

TL;DR

This paper investigates when the adaptive landscape metaphor accurately describes evolutionary dynamics in replicator systems and explores the structure of mean fitness even when it is not maximized.

Contribution

It identifies conditions under which the adaptive landscape metaphor holds exactly and analyzes the implications of mean fitness structure on system dynamics.

Findings

Adaptive landscape accurately models some replicator dynamics.

Conditions for exact correspondence between dynamics and fitness maximization.

Insights into mean fitness landscape structure and its evolutionary implications.

Abstract

Sewall Wright's adaptive landscape metaphor penetrates a significant part of evolutionary thinking. Supplemented with Fisher's fundamental theorem of natural selection and Kimura's maximum principle, it provides a unifying and intuitive representation of the evolutionary process under the influence of natural selection as the hill climbing on the surface of mean population fitness. On the other hand, it is also well known that for many more or less realistic mathematical models this picture is a sever misrepresentation of what actually occurs. Therefore, we are faced with two questions. First, it is important to identify the cases in which adaptive landscape metaphor actually holds exactly in the models, that is, to identify the conditions under which system's dynamics coincides with the process of searching for a (local) fitness maximum. Second, even if the mean fitness is not maximized in the process of evolution, it is still important to understand the structure of the mean fitness manifold and see the implications of this structure on the system's dynamics. Using as a basic model the classical replicator equation, in this note we attempt to answer these two questions and illustrate our results with simple well studied systems.

Equations

$${\dot{u}}_i=u_i(m_i-\bar{m}(\mathbf{u})),i=1,\dots ,n,$$

$$\bar{m}(\mathbf{u})=\sum _{i=1}^n{m}_i{u}_i,$$

$$S_n=\{\mathbf{u}\in {\mathbf{R}}^n:u_i\ge 0,\sum _{i=1}^n{u}_i=1\}.$$

$$\dot{\bar{m}}(t)=\sum _{i=1}^n{m}_i{u}_i^2(t)-{\left(\sum _{i=1}^n{m}_i{u}_i(t)\right)}^2=\text{Var}(t)\ge 0,$$

$$\Sigma =\{(\bar{m},\mathbf{u}):\bar{m}(\mathbf{u})=\sum _{i=1}^n{m}_i{u}_i,\mathbf{u}\in S_n\},$$

$${\dot{u}}_i=u_i(a_i{u}_i-\bar{m}(\mathbf{u})),i=1,\dots ,n,$$

$$\bar{m}(\mathbf{u})=\sum _{i=1}^n{a}_i{u}_i^2.$$

$$\dot{\bar{m}}(t)\ge 0,$$

$$\dot{u}_{i}=u_{i}\Bigl{(}\sum_{j=1}^{n}a_{ij}u_{j}-\sum_{k,l=1}^{n}a_{kl}u_{k}u_{l}\Bigr{)}$$

$$\dot{u}_{i}=u_{i}\bigl{(}({\boldsymbol{Au}})_{i}-{\boldsymbol{u}}\cdot{\boldsymbol{Au}}\bigr{)},\quad i=1,\ldots,n,$$

$$\bar{m}(\mathbf{u})=\mathbf{u}\cdot \mathbf{Au}=\sum _{i,j=1}^n{a}_{ij}{u}_i{u}_j.$$

$$\Sigma=\Bigl{\{}(\bar{m},{\boldsymbol{u}})\colon\bar{m}({\boldsymbol{u}})=\sum_{i,j=1}^{n}a_{ij}u_{i}u_{j},\,{\boldsymbol{u}}\in S_{n}\Bigr{\}}.$$

$$\mathbf{A}=\frac{1}{2}(\mathbf{A}+{\mathbf{A}}^{\top })+\frac{1}{2}(\mathbf{A}-{\mathbf{A}}^{\top })=\mathbf{B}+\mathbf{C},$$

$$\bar{m}(\mathbf{u})=\mathbf{Au}\cdot \mathbf{u}=\mathbf{Bu}\cdot \mathbf{u}.$$

$${\mathbf{U}}^{\top }\mathbf{B}\mathbf{U}=\mathbf{\Lambda }=\mathrm{diag}({\lambda }_1,\dots ,{\lambda }_n),$$

$$\bar{m}(\mathbf{w})=\sum _{i=1}^k{\lambda }_i^{+}{w}_i^2+\sum _{j=k+1}^n{\lambda }_j^{-}{w}_j^2,$$

$$W_n=\{\mathbf{w}\in {\mathbf{R}}^n:\mathbf{w}\cdot {\mathbf{U}}^{\top }1=1,\mathbf{Uw}\ge 0,1=(1,\dots ,1)\in {\mathbf{R}}^n\}.$$

$$\sum _{i=1}^k{\lambda }_i^{+}{w}_i^2⟶\max ,\mathbf{w}\in W_n.$$

$$\mathcal{L}(\mathbf{u})=\mathbf{Bu}\cdot \mathbf{u}-\mu (\mathbf{u}\cdot 1-1),\mu \in \mathbf{R}.$$

$$\mathbf{Bu}=\mu 1.$$

$$\mathbf{Bw}\cdot \mathbf{w}<0(>0)$$

$$\mathbf{w}\cdot 1=0$$

$${\mathbf{r}}_1=(1,-1,0,\dots ,0),{\mathbf{r}}_2=(1,0,-1,\dots ,0),\dots ,{\mathbf{r}}_{n-1}=(1,0,\dots ,-1).$$

$$\mathbf{Bw}\cdot \mathbf{w}=\sum _{k=1}^{n-1}{\beta }_{ij}{\xi }_i{\xi }_j,{\beta }_{ij}={\mathbf{Br}}_i\cdot {\mathbf{r}}_j.$$

$$\hat{\mathbf{u}}\cdot \mathbf{Au}>\mathbf{u}\cdot \mathbf{Au}$$

$$ϵ(\mathbf{w}\cdot \mathbf{A}\hat{\mathbf{u}})+{ϵ}^2(\mathbf{w}\cdot \mathbf{Aw})<0.$$

$$\mathbf{w}\cdot \mathbf{Aw}=\mathbf{w}\cdot \mathbf{Bw}<0,$$

$${\mathbf{c}}_i=\sum _{k=1}^n{\mu }_{ik}{\mathbf{b}}_k,{\mathbf{b}}_k=(b_{k1},\dots ,b_{kn}),{\mathbf{c}}_i=(c_{i1},\dots ,c_{in}),$$

$$\sum _{k=1}^n{\mu }_{ik}=0,i=1,\dots ,n.$$

$$\mathbf{A}\hat{\mathbf{u}}=\mathbf{B}\hat{\mathbf{u}}+\mathbf{C}\hat{\mathbf{u}}=\bar{m}(\hat{\mathbf{u}})1.$$

Full text

Keywords:

Adaptive fitness landscape, population fitness, Fisher’s theorem of natural selection, Kimura’s maximum principle

AMS Subject Classification:

92D10, 92D15

1 Introductory remarks and mathematical problem statement

The infamous Fisher’s theorem of natural selection asserts that “the rate of increase in [mean] fitness of any organism at any time is equal to its genetic variance in fitness at that time”[3]. This theorem, for which Fisher himself did not provide any mathematical justification, can be indeed proved within the simplest possible framework of natural selection acting on a haploid population with $n$ alleles in the absence of mutations and recombination, assuming that we deal with an infinite population. Let $u_{i}$ denote the frequency of the $i$-th allele and let $m_{i}$ be the corresponding Malthusian fitness, then a very simple replicator equation

[TABLE]

where

[TABLE]

describes the change with time of the vector of frequencies ${\boldsymbol{u}}(t)=(u_{1}(t),\ldots,u_{n}(t))\in S_{n}$ that belongs for any time moment $t$ to the simplex

[TABLE]

Considering $\bar{m}$ as a function of $t$ and calculating the derivative of this function along the orbits of our system, we get

[TABLE]

which can be taken literally as a mathematical form of Fisher’s theorem.

Moreover, system (1.1) can be also used to illustrate another important metaphor of theoretical population genetics, — namely, Wright’s adaptive landscape [11]. In words, Wright envisioned a surface of the mean fitness and argued that the mean fitness maxima represent the so-called adaptive picks, to which natural selection leads evolving populations. For system (1.1) we can give an illustration of this process by considering how the mean population fitness, which, according to the calculations above, is non-decreasing, changes with time. To be precise and introduce notations that we will use in the following, let $\gamma_{t}$ be an orbit of dynamical system (1.1) (see Fig. 1.1(a)), whose behavior is quite simple because, assuming there is unique maximum $m_{1}=\max\{m_{1},\ldots,m_{n}\}$, we find ${\boldsymbol{u}}(t)\to{\boldsymbol{e}}_{1}$, where ${\boldsymbol{e}}_{1}=(1,0,\ldots,0)$. The mean fitness surface is given by

[TABLE]

and to each $\gamma_{t}\in S_{n}$ there corresponds a curve $\Gamma_{t}$ on $\Sigma$, which in this particular case will be approaching the highest vertex on the hyperplane $\bar{m}({\boldsymbol{u}})=\sum_{i=1}^{n}m_{i}u_{i}$ that also must satisfies the constraint ${\boldsymbol{u}}\in S_{n}$ (see Fig. 1.1).

Of course, Sewall Wright had a more general picture in mind. For example he admitted the existence of several adaptive picks, such that the population evolves in the direction of a nearest one, which may be only local and not global maximum. We can illustrate this idea using yet another extremely simple mathematical model, namely, the model of an autocatalytic growth. Now

[TABLE]

where now to guarantee that the simplex $S_{n}$ is invariant, we have that the mean population fitness is defined as

[TABLE]

In this case, if we still stick to a haploid interpretation of our model, we now deal with a frequency dependent selection, albeit with its simplest case. Again, as in the previous example, it is straightforward to show that

[TABLE]

and all but one alleles disappear. Now, however, which allele will survive depends also on the initial conditions, as well as the parameters of the systems. Note that each vertex of the simplex is asymptotically stable and each has a basin of attraction of non-zero measure (see Fig. 1.2).

Taken together, Wright’s adaptive landscape and Fisher’s theorem of natural selection can, in some sense, to be taken as a mathematical manifestation of a “layman’s” understanding of how natural selection works. In short, the natural selection can be regarded as a process of fitness maximization (much more biological context can be found in [1]). On the other hand, it is well known in population genetics that the fitness is maximized only in very simple models, such as (1.1) or (1.2), see an extensive discussion in [1] and references therein. It is not our goal in this short note to discuss various implications of these two seemingly mutually inconsistent points of view at the evolutionary process. We agree, however, with the general premise of [1] that it is important to reconcile them and, in particular, find maximization principles that apply conditionally. In other words, it is important to understand under which conditions our theoretical mathematical models imply fitness maximization.

For example, one such well known condition for the replicator equation of the form

[TABLE]

is that the matrix ${\boldsymbol{A}}=[a_{ij}]$ that defines this equation is symmetric, ${\boldsymbol{A}}={\boldsymbol{A}}^{\top}$. In the language of the game theory we deal in this case with the so-called partnership games. It is known that in this case the mean fitness $\bar{m}({\boldsymbol{u}})=\sum_{i,j=1}^{n}a_{ij}u_{i}u_{j}$ is 1) non-decreasing, since it is a strict Lyapunov function, and, even more, 2) it is a potential for the replicator equation, which hence becomes in this case a gradient system. In the language of theoretical population genetics it is said that in this case the state of the gene pool evolves under the influence of natural selection that not only the fitness increases, but also at a maximal rate (this is often called Kimura’s maximum principle [8]).

The nontrivial point here to see is that Kimura’s maximum principle becomes a precise theorem only if one replaces the usual Euclidian metric, generated by the standard dot product, with another Riemannian metric generated by the so-called Shahshahani inner product [9, 5]. This maximum principle can be extended for some particular cases of frequency depended selection [10], but in general we know that if ${\boldsymbol{A}}\neq{\boldsymbol{A}}^{\top}$ for the replicator equation, the mean fitness is not always maximized (actually, it can be shown that in this case the so-called KL-divergence between the initial and current distributions is minimized, see [7] for information–theoretic approach to the replicator equation), and are left with an open question what are the necessary or/and sufficient conditions on the matrix ${\boldsymbol{A}}$ of the replicator equation, under which the dynamics of its solutions coincides with mean fitness maximization. In this note we incept such a study, derive several such conditions and illustrate them with examples of well studied systems.

2 Replicator Equation

Consider the now classical replicator equation (see, e.g., [5, 6] and references therein)

[TABLE]

where ${\boldsymbol{u}}\in S_{n}$, matrix ${\boldsymbol{A}}=[a_{ij}]$ is given and the mean fitness is

[TABLE]

Here and below the dot denotes the standard dot product in $\mathbf{R}^{n}$: ${\boldsymbol{x}}\cdot{\boldsymbol{y}}=\sum_{i=1}^{n}x_{i}y_{i}$. Together with solutions to (2.1) we also consider the surface of the mean fitness:

[TABLE]

The orbits $\gamma_{t}$ of (2.1) correspond to the curves $\Gamma_{t}$ with directions on $\Sigma$, see the introductory section for simple examples (Figs. 1.1 and 1.2).

We represent matrix ${\boldsymbol{A}}$ as a sum of two matrices:

[TABLE]

that is, as a sum of symmetric (${\boldsymbol{B}}$) and skew-symmetric (${\boldsymbol{C}}$) matrices. From (2.2) it follows that the mean fitness depends only on the symmetric part of ${\boldsymbol{A}}$, which is a simple yet key observation in our analysis:

[TABLE]

This means that there exists an orthogonal transformation ${\boldsymbol{U}}$ such that

[TABLE]

where $\lambda_{i}$ are real eigenvalues of ${\boldsymbol{B}}$.

This orthogonal transformation ${\boldsymbol{u}}={\boldsymbol{U}}{\boldsymbol{w}}$ brings the quadratic form (2.3) to the canonical form

[TABLE]

where we denote $\lambda_{i}^{+}$ and $\lambda_{j}^{-}$ the positive and negative eigenvalues of ${\boldsymbol{B}}$ respectively. The same transformation maps the simplex $S_{n}$ into the convex set

[TABLE]

This means that the maximum value of the fitness surface $\Sigma$ can be obtained as a result of solving the problem of convex programming

[TABLE]

Solution to this problem always exists but can be non-unique. An analogous statement holds for the minimal value on the surface $\Sigma$.

Yes another approach to find the extrema on $\Sigma$ is to use the classical method of Lagrange multipliers. We write down Lagrange’s function

[TABLE]

The necessary condition for a local extremum takes the form

[TABLE]

The equality ${\boldsymbol{u}}\cdot{\boldsymbol{1}}=1$ implies that $\mu={\boldsymbol{B^{-1}{\boldsymbol{1}}\cdot{\boldsymbol{1}}}}$. Vector ${\boldsymbol{u}}\in S_{n}$ is a local maximum (minimum) if the condition of being negative definite (respectively, positive definite) for the quadratic form

[TABLE]

is satisfied for any ${\boldsymbol{w}}\in\mathbf{R}^{n}$ for which

[TABLE]

holds.

Consider a basis in the orthogonal complement to ${\boldsymbol{1}}$:

[TABLE]

Then ${\boldsymbol{w}}=\sum_{k=1}^{n-1}\xi_{k}{\boldsymbol{r}}_{k}$, and

[TABLE]

Therefore, checking the condition of sign defiteness of the form (2.5) boils down to checking the sign defitness of the form (2.7), which is generated by the matrix with elements $\beta_{ij}={\boldsymbol{Br}}_{i}\cdot{\boldsymbol{r}}_{j}$.

In particular, condition (2.5) implies the following proposition, which also connects the notion of an evolutionary stable state with the fact that an equilibrium of (2.1) coincides with a local strict maximum of the fitness surface.

Proposition 2.1**.**

If an equilibrium ${\boldsymbol{\hat{u}}}\in S_{n}$ of the replicator equation (2.1) is evolutionary stable then it is a strict local maximum of the fitness surface $\Sigma$.

Proof.

${\boldsymbol{\hat{u}}}$ is evolutionary stable means (e.g., [5])

[TABLE]

for all ${\boldsymbol{u}}$ in some neighborhood of ${\boldsymbol{\hat{u}}}$. Take ${\boldsymbol{u}}={\boldsymbol{\hat{u}}}+\epsilon{\boldsymbol{w}},\,\epsilon>0$. Since we require ${\boldsymbol{u}}\in S_{n}$ then ${\boldsymbol{w}}\cdot{\boldsymbol{1}}=0$. The condition for the evolutionary stability implies

[TABLE]

Since ${\boldsymbol{A\hat{u}}}=\bar{m}({\boldsymbol{\hat{u}}}){\boldsymbol{1}}$, then the first term is zero and hence

[TABLE]

which concludes the proof. ∎

In other words, the existence of a strict maximum of the surface $\Sigma$ is a necessary condition for an equilibrium to be evolutionary stable.

Now let $\gamma_{t}$ be an orbit of (2.1) that approaches the equilibrium ${\boldsymbol{\hat{u}}}\in\operatorname{int}S_{n}=\{{\boldsymbol{u}}\in S_{n}\colon u_{i}>0,\,i=1,\ldots,n\}$. Consider the conditions under which this would correspond to the case when the curve $\Gamma_{t}$ on the surface $\Sigma$ approaches an extremal point.

Theorem 2.2**.**

An isolated equilibrium ${\boldsymbol{\hat{u}}}\in\operatorname{int}S_{n}$ of (2.1) coincides with an extremal point on the surface $\Gamma$ if and only if the vectors that compose the rows of the skew symmetric matrix ${\boldsymbol{C}}=[c_{ij}]$ can be represented as linear combinations of the vectors that compose the rows of ${\boldsymbol{B}}=[b_{ij}]$:

[TABLE]

and, additionally,

[TABLE]

Proof.

For the equilibrium point we have

[TABLE]

The solution to (2.4) is also a solution to (2.10) if and only if

[TABLE]

Indeed, the sufficiency of (2.11) is clear. Let us prove its necessity. We multiply the equalities (2.4) and $\eqref{eq2:10}$ by ${\boldsymbol{\hat{u}}}$ and use the fact ${\boldsymbol{C\hat{u}}}\cdot{\boldsymbol{\hat{u}}}=0$ to obtain ${\boldsymbol{B\hat{u}}}\cdot{\boldsymbol{\hat{u}}}=\bar{m}({\boldsymbol{\hat{u}}})=\mu$ and therefore (2.11) holds. Let us show now that in order for the solutions of

[TABLE]

be in the kernel of ${\boldsymbol{C}}$, it is necessary and enough for the conditions (2.8) and (2.9) to hold. Assume that (2.8) and (2.9) hold. That is

[TABLE]

Hence

[TABLE]

In the opposite direction, since the rows of ${\boldsymbol{B}}$ linearly independent, we can always write (2.8) with some $\mu_{ij}$. Let ${\boldsymbol{\hat{u}}}$ solve (2.10), therefore (2.13) holds, since ${\boldsymbol{C\hat{u}}}=0$, therefore, we must have ${\boldsymbol{\mathcal{M}}}{\boldsymbol{1}}=0$, that is, the equalities (2.9). ∎

Corollary 2.3**.**

Assume that conditions of Theorem 2.2 hold and ${\boldsymbol{\hat{u}}}\in\operatorname{int}S_{n}$ be a strict local maximum of the mean fitness. Then ${\boldsymbol{\hat{u}}}$ is an evolutionary stable state.

Proof.

We have for any ${\boldsymbol{u}}\in S_{n}$ that

[TABLE]

On the other hand

[TABLE]

By assumption

[TABLE]

therefore

[TABLE]

∎

Corollary 2.4**.**

Assume that conditions of Theorem 2.2 hold. Then the graph of the curve $\Gamma_{t}$ on the surface $\Sigma$ has neither local minima nor local maxima for ${\boldsymbol{u}}(t)\neq{\boldsymbol{\hat{u}}},\,{\boldsymbol{u}}(t)\in\operatorname{int}S_{n}$.

Proof.

Taking the derivative of the mean fitness along the orbits of (2.1), we get

[TABLE]

Looking for a contradiction, assume that ${\boldsymbol{u}}^{\ast}\in\operatorname{int}S_{n}$ be a local extremum of $\bar{m}$ for some $t=t_{\ast}$ and therefore the condition (2.4) holds with ${\boldsymbol{u}}={\boldsymbol{u}}^{\ast}$. Consider the vector

[TABLE]

Since the simplex is convex, ${\boldsymbol{u}}_{\lambda}$ belongs to $\operatorname{int}S_{n}$ for any admissible $\lambda$. On the other hand, from (2.4) it follows

[TABLE]

and conditions (2.8) and (2.9) imply

[TABLE]

Similarly, for ${\boldsymbol{\hat{u}}}$ we have

[TABLE]

which means that

[TABLE]

and hence $\bar{m}$ is constant along the interval corresponding to $0\leq\lambda\leq 1$, which implies that $\mu=\bar{m}$ and hence ${\boldsymbol{u}}^{\ast}={\boldsymbol{\hat{u}}}$. ∎

3 Examples

In this section we consider two simple examples that illustrate the theoretical results we obtained in the previous section.

Example 3.1**.**

Consider a hypercycle equation [2] in $S_{3}$ with the matrix

[TABLE]

The coordinates of the internal equilibrium are given by

[TABLE]

Simple calculations yield

[TABLE]

The extremal point of the fitness surface is given by

[TABLE]

where

[TABLE]

We can use equation (2.7) to check that this point is actually a local maximum. We obtain the quadratic form

[TABLE]

which is negative definite if

[TABLE]

In particular, the last condition holds if, for example, $k_{1}=0.25,k_{2}=0.3,k_{3}=0.35$. In this case

[TABLE]

whereas the fitness surface maximum is at the point ${\boldsymbol{u}}_{m}=(0.2692,0.3846,0.3462)$, moreover

[TABLE]

which implies that the curve on the surface $\Sigma$ can both increase and decrease while the system moving along trajectory $\gamma_{t}$. Let us show in this example that if the conditions (2.8), (2.9) hold then ${\boldsymbol{u}}$ and ${\boldsymbol{\hat{u}}}$ coincide.

We have

[TABLE]

Conditions (2.8), (2.9) satisfied if and only if $k_{1}=k_{2}=k_{3}=k$. In this case, putting $k=1$, the eigenvalues of ${\boldsymbol{B}}$ are $\lambda_{1}=1,\lambda_{2}=\lambda_{3}=1/2$. The orthogonal transformation ${\boldsymbol{u}}={\boldsymbol{Uw}}$ that brings the mean fitness to the canonical form $w_{1}^{2}-\frac{1}{2}(w_{2}^{2}+w_{3}^{2})$ takes the form

[TABLE]

The condition ${\boldsymbol{u}}\cdot{\boldsymbol{1}}={\boldsymbol{w}}\cdot{\boldsymbol{U}}^{\top}{\boldsymbol{1}}=1$ implies the constraint $w_{1}=\sqrt{3}/3$ and simplex $S_{3}$ is being transformed into the convex set

[TABLE]

That it, in new coordinates, the fitness surface if given by the two-sheeted hyperboloid

[TABLE]

which is a conus with the vertex at the point $w_{2}=w_{3}=0$. In this particular case the equilibrium point coincides with the maximum of the fitness surface. Any orbit of the system corresponds to the curve $\Gamma_{t}$ on the surface $\Sigma$, which climbs the conus for $t\to\infty$ (see also Fig. 3.1 for a graphical illustration).

Example 3.2**.**

It is quite interesting to consider examples of fitness surfaces for which there are several local extrema and moreover periodic solutions exist.

Consider a replicator equation in $S_{4}$ with the matrix (we take this example from [4])

[TABLE]

It is known that for $|\mu|<1$ the system is permanent, that is, its solutions are separated from the boundaries of the simplex, and for $\mu=0$ Poincaré–Andronov–Hopf bifurcation happens with appearance of a stable limit cycle for $\mu>0$. If $|\mu|>1$ then the system stops being permanent and its omega limit set is the heteroclinic cycle that contains the vertices of $S_{4}$.

In this case we have

[TABLE]

and the eigenvalues of ${\boldsymbol{B}}$ are $\lambda_{1}=1-\mu,\lambda_{2,3}=\mu,\lambda_{4}=-(1+\mu)$. The orthogonal transformation is explicitly given by

[TABLE]

which yields

[TABLE]

From the condition ${\boldsymbol{w}}\cdot{\boldsymbol{U}}^{\top}{\boldsymbol{1}}=1$ we find that $w_{1}=1/2$. Simplex $S_{4}$ is mapped to the convex set

[TABLE]

If $|\mu|<1$ then the maximum value of the fitness surface is reached at $w_{4}=0$. The set $W_{4}$ is this case represents a square in the plane $(w_{2},w_{3})$ with the vertices at $p_{i}=\{\pm\sqrt{2}/4,\pm\sqrt{2}/4\}$, and at all fours vertices of this square the fitness surface attains local maxima

[TABLE]

The internal equilibrium ${\boldsymbol{\hat{u}}}=(1/4,1/4,1/4,1/4)$ in the new variables is ${\boldsymbol{\hat{w}}}=(1/2,0,0,0)$, moreover

[TABLE]

and we know that this equilibrium for $-1<\mu<0$ is asymptotically stable. In this case the fitness surface at this point has the global maximum, since $\bar{m}({\boldsymbol{\hat{w}}})>\bar{m}(p_{i})$. In general we have a system of five local summits, among which the central one has the highest fitness, and while the orbits $\gamma_{t}$ approach the equilibrium, the curves $\Gamma_{t}$ on $\Sigma$ can both increase and decrease while approaching ${\boldsymbol{\hat{w}}}$.

If $0<\mu<1$ then the curves $\Gamma_{t}$ will be approaching a closed curve that surrounds ${\boldsymbol{\hat{w}}}$, moreover, the maxima at the points $p_{i},i=1,\ldots,4$ will be outside of this curve. On simplex $S_{4}$ the points $p_{i}$ correspond to the points $P_{1}=(0,0,1/2,1/2),P_{2}=(1/2,0,0,1/2),P_{3}=(1/2,1/2,0,0),P_{4}=(0,1/2,1/2,0)$. If $\mu>1$ then $P_{i}$ still give the maximal values of $\bar{m}({\boldsymbol{u}})$, however, $\bar{m}({\boldsymbol{\hat{w}}})<0$ and now it is a global minimum of the fitness surface. This case corresponds to the heteroclinic cycle.

Finally, if $\mu<-1$ then ${\boldsymbol{\hat{u}}}$ is unstable and the maximum value of the form $\bar{m}({\boldsymbol{w}})$ is reached when $w_{2}=w_{3}=0$. It is reached at the points $Q_{1,2}=(1/2,0,0,\pm 1/2)$, and

[TABLE]

and hence $\bar{m}(Q_{1,2})>\bar{m}({\boldsymbol{\hat{w}}})$. In this case we also have a heteroclinic cycle.

A schematic illustration of the system’s dynamics and the corresponding structure of the mean fitness surface is given in Fig. 3.2.