# Extremal Coupled Map Lattices

\[

## Abstract

We propose a model for co-evolving ecosystems that takes into account two levels of description of an organism, for instance genotype and phenotype. Performance at the macroscopic level forces mutations at the microscopic level. These, in turn, affect the dynamics of the macroscopic variables. In some regions of parameter space, the system self-organises into a state with localised activity and power law distributions.

PACS: 05.45.+b, 64.60.Lx, 87.10.+e

\]

Complex extended systems showing a lack of scale in their features appear to be widespread in nature, being as diverse as earthquakes , creep phenomena , material fracturing , fluid displacement in porous media , interface growth , river networks  and biological evolution . At variance with equilibrium statistical mechanics, these systems do not need any fine tuning of a parameter to be in a critical state. To explain this behaviour, Bak, Tang and Wiesenfeld introduced the concept of self-organised criticality (SOC) through the simple sand-pile model . In recent years, several models with extremal dynamics have been shown to exhibit SOC in the presence of a sufficiently decorrelated signal .

Several years ago, Bak and Sneppen (BS) proposed a SOC model  for the co-evolution of species. In the BS model each species occupies a site on a lattice. Each site $`j`$ is assigned a fitness, namely a random number between $`0`$ and $`1`$. At each time step in the simulation the smallest fitness is found. Then the fitnesses of the minimum and of its nearest neighbours are updated according to the rule

$$f_{n+1}=F(f_n)$$
(1)   
that assigns a new fitness $`f_{n+1}`$ at time $`n+1`$ to the chosen lattice site. This corresponds to the extinction of the less fit and its impact on the ecosystem. Indeed, in the original BS model, the function $`F`$ is just a random function with a uniform distribution between 0 and 1. The system reaches a stationary critical state in which the distribution of fitnesses is zero below a certain threshold and uniform above it (the actual value of the threshold depends on the updating rule). It has been shown  that the exact nature of the updating rule is not relevant. Indeed, the use of a chaotic map instead of a random update, preserves the universality class (even if the final distribution may be altered).

As a result of the dynamical rules, the BS system exhibits sequences of causally connected evolutionary events called avalanches . The number of avalanches $`N`$ follows a power law distribution $`N(s)s^\tau `$ where $`s`$ is the size of the avalanche and $`\tau 1.07`$ . Other quantities also exhibit a power law behaviour with their own critical exponents. Prominent among them are the first and all return time distributions of activity (a site is defined as active when its fitness is the minimum one),

$$P_f(t)t^{\tau _f},P_a(t)t^{\tau _a},$$
(2)   
where $`\tau _f1.58`$ and $`\tau _a0.42`$ .

In the BS model, two basic ingredients are needed for SOC to occur :   
i) Order from extremal dynamics (minimum rule).   
ii) Disorder from the updating rule (stochastic or otherwise).

An oversimplification in the BS model is apparent: each species is described by a single variable. The minimum rule is applied to this variable, and so is the effect of mutations. In natural evolving systems, however, at least two—interacting—levels of organisation are present, and both play a role in the evolution. Mutations occur at a microscopic level, namely at the molecular level of the genotype. This affects a macroscopic level, defining the phenotype. Natural selection acts on the phenotype, allowing the survival of some and the extinction of others, according to the observed power law distributions.

In this paper we propose a new model, namely an Extremal Coupled Map Lattice (ECML), that takes into account this two level structure. The ecosystem is represented as an ensemble of $`N`$ species arranged on a one-dimensional lattice. Each species is described by means of macroscopic variable $`x_i`$ subject to a nonlinear dynamics and coupled to its nearest neighbours. In general, $`x_i`$ can be identified with a population, or some other function of the phenotype. The control parameter $`\lambda _i`$ of the nonlinear dynamics is identified with the microscopic level (genotype). To find the exact function connecting these two levels of description is beyond the scope of this simplified model. For this reason, we chose a logistic map for the independent evolution of $`x_i`$ with $`\lambda _i`$ acting as the nonlinear parameter in the map . With this in mind, the evolution of each species is given by:

$$x_{n+1}^i=(1ϵ)f(x_n^i)+\frac{ϵ}{2}\left[f(x_n^{i1})+f(x_n^{i+1})\right]$$
(3)   
where $`f(x)`$ is the logistic map $`f(x)=\lambda x(1x)`$ (in general any chaotic function should do). Each site has its own $`\lambda _i`$, extracted from a fixed distribution $`g(\lambda )`$ . The local coupling of strength $`ϵ_i`$ emulates the ecological interaction between neighbouring species. A general description should include inhomogeneities in $`ϵ`$. Here, for simplicity we limit our analysis to the homogeneous case, $`ϵ_iϵ`$.

We consider this CML as the substrate on which evolution takes place. The parameter $`\lambda _i`$, that determines the behaviour of $`x_i`$, is regarded as the microscopic level, the genotype, and is subject to mutation. We suppose that, through mutation, a species is able to alter this strategy to adapt to the environment defined by the collective behaviour. For this we propose an extremal mechanism, akin to the BS model . The species that, at each time step, has the minimum value of $`x`$, is considered the candidate to mutation. Its $`\lambda `$ is replaced by a new value drawn from the distribution $`g(\lambda )`$.

We can summarise these simple dynamical rules by

1. Evaluate expression (3).
2. Find the site with the absolute minimum fitness on the lattice (this site will be called the “active” site).
3. Change the value of $`\lambda `$ of the active site.
4. Go to step 1.

In Fig. 1 we show the space-time picture of the active sites for four different sets of parameters. In all cases we take $`g(\lambda )`$ uniform within the range $`(\lambda _0,\lambda _0+\mathrm{\Delta }\lambda )`$ and zero outside, with $`\mathrm{\Delta }\lambda =0.1`$. This four cases can be considered paradigmatic of the behaviours observed in a wide region of parameter space. Fig. 1(a) and (b) correspond to a rather high value of the coupling, $`ϵ=0.5`$. Fig. 1(c) and (d) have, instead, $`ϵ=0.1`$. The value of $`\lambda _0`$ is $`3.9`$ in (b) and (d), $`3.7`$ in (a), and $`3.8`$ in (c). Depending on the value of the parameters, the behaviour of the successive minima can be classified in one of the following four categories.

* discontinuous lines at $`i=\text{constant}`$: The activity is concentrated in a few sites of the system, and in first approximation the same sites remain active as time goes by.
* uniform: The activity is spread all over the system without any apparent order.
* “worms”: The activity is mostly concentrated in a few sites of the system. The activity, however, wanders as time goes by following a random pattern.
* clusterized: The activity is spread all over the system, but in clusters. At variance with (b), at any given time, one sees regions in space where there is no activity.

Let us discuss in more detail Fig. 1(c). The worms are created during the transient. Once created, each worm wanders in space till it encounters a second worm. At that moment they annihilate each other. We have observed that the inactive regions correspond to a periodic pattern in space . The worms, in turn, correspond to defects in the periodicity. Eventually, the last two surviving worms merge, and from there on the activity remains concentrated on a single “fat” worm.

For each one of these cases, we have computed the first return time distribution. This corresponds to the distribution of times between consecutive activity in the same site. We observe that when the successive minima are distributed uniformly in space, the first return time distribution is close to an exponential (see Fig. 2(b)). As soon as the activity is not uniformly distributed, be it lines or irregular clusters, the first return time distribution exhibits a power law decay (asymptotically in case (d)).

To make this picture more quantitative, we present in Fig. 3 a phase diagram, in $`\lambda `$-$`ϵ`$ space ($`\mathrm{\Delta }\lambda `$ is held constant). Phase I is characterised by the presence of power law behaviour in the first return time distribution, $`P(\tau )\tau ^\alpha `$. The value of the exponent of this distribution is not constant within the region. Phase II, on the other hand, is characterised by a non power law behaviour (close to exponential) in the first return time distribution (see Fig. 2(b)). This is tantamount to saying that in phase I there is clusterization of the activity, that is absent in phase II. In Fig. 3 we indicate also a third phase (III), which corresponds to an asymptotic power law (see Fig. 2(d)), with an exponent that is different from that in region I. The border of this region (that we show only schematically) is not as well defined as the border between regions I and II.

This picture holds, qualitatively, even in the absence of the evolutionary dynamics. Indeed, if in every time step we track the position of the minimum, but do not change the value of its $`\lambda `$ (this corresponds, effectively, to setting $`\mathrm{\Delta }\lambda =0`$), we still obtain the four abovementioned regimes and the corresponding first return time distributions. As $`\mathrm{\Delta }\lambda `$ changes, the value of the exponent in the first return time distribution (in region I) changes as well. Moreover, the area of region II increases with increasing $`\mathrm{\Delta }\lambda `$.

The presence of evolutionary dynamics (i. e. $`\mathrm{\Delta }\lambda >0`$), has yet another effect. The distribution of $`\lambda `$ evolves in time, from an originally uniform to a stationary non-uniform one. Indeed, this stationary distribution is peaked close to $`\lambda _0`$ and monotonically decresases for larger values of $`\lambda `$: The system has “self-organized”. Since the extremal dynamics favors higher values of $`x`$, and the higher the value of $`\lambda `$ the more likely small values of $`x`$ are , then large values of $`\lambda `$ are more likely to be updated. The exact shape of this distribution depends on the values of the parameters .

Summarizing, the dynamics of extremal coupled maps exhibits a behavior usually associated with criticality:

* First return time distribution is a power law (even in the absence of extremal dynamics)
* there is self-organization in the $`\lambda `$ space
* the active regions are localized in space, very much like the avalanches in SOC models (BS).

It is worth emphasizing that this behavior (except the self organization in $`\lambda `$ space) remains even in the absence of extremal dynamics. In particular, the first return time distribution exhibits a power law behavior. This implies that power law distributions are compatible with non-extremal systems, as has been previously observed by Newman an coworkers , in the context of noise-driven systems.

Both SOC and CML models have been used separately in the past to describe evolving ecosystems. To our knowledge, the model presented here is the first synthesis of both ideas that includes the level structure (genotype/phenotype) of real organisms.