Evolutionary games on minimally structured populations

TL;DR

This paper introduces a minimal hierarchical population structure model to analyze evolutionary games, revealing effects of spatial structure independent of detailed topology and uncovering mechanisms that sustain cooperation and biodiversity.

Contribution

The model simplifies population structure to two parameters, enabling analysis of spatial effects without detailed topology, and uncovers new mechanisms for cooperation and biodiversity.

Findings

Existence of a transition favoring cooperation as benefit-to-cost ratio decreases.

Constant influx of defectors can sustain cooperation.

Reproduction of effects like biodiversity maintenance and oscillations seen in spatial models.

Abstract

Population structure induced by both spatial embedding and more general networks of interaction, such as model social networks, have been shown to have a fundamental effect on the dynamics and outcome of evolutionary games. These effects have, however, proved to be sensitive to the details of the underlying topology and dynamics. Here we introduce a minimal population structure that is described by two distinct hierarchical levels of interaction. We believe this model is able to identify effects of spatial structure that do not depend on the details of the topology. We derive the dynamics governing the evolution of a system starting from fundamental individual level stochastic processes through two successive meanfield approximations. In our model of population structure the topology of interactions is described by only two parameters: the effective population size at the local scale and the relative strength of local dynamics to global mixing. We demonstrate, for example, the existence of a continuous transition leading to the dominance of cooperation in populations with hierarchical levels of unstructured mixing as the benefit to cost ratio becomes smaller then the local population size. Applying our model of spatial structure to the repeated prisoner's dilemma we uncover a novel and counterintuitive mechanism by which the constant influx of defectors sustains cooperation. Further exploring the phase space of the repeated prisoner's dilemma and also of the "rock-paper-scissor" game we find indications of rich structure and are able to reproduce several effects observed in other models with explicit spatial embedding, such as the maintenance of biodiversity and the emergence of global oscillations.

Figures

a). More importantly, diversity was also maintained for finite local population sizes if global mixing was present. Simulations of the time evolution of $\rho({\bf x})$ also revealed a Hopf bifurcation leading to the oscillation of the global average as $\mu$ was increased above a critical value $\mu_{\rm c}$ depending on $N$ (Fig.4a). These results show that previous results obtained from simulations of populations constrained to different lattice topologies can be considered universal in the sense that not only lattices, but any population structure that can be approximated by two distinct i

d,e). Further, for any given value of $N$ and $\mu$ the global configuration is described by a unique $I$ due to the presence of diffusion. For appropriate values of the parameters the global average converges to a stationary value in phase (ii) similarly to case of explicit spatial embedding (Fig.4c).

Equations

$${\pi }_k={\pi }_{\mathrm{base}}+\sum _{j=1}^d\frac{{\alpha }_{kj}(n_j-{\delta }_{kj})}{N-1},$$

$$T_{ik}=\frac{n_i}{N}\frac{{\pi }_k{n}_k}{\bar{\pi }N},$$

$${\dot{x}}_k=a_k(\mathbf{x})+\sum _{j=1}^{d-1}{c}_{kj}(\mathbf{x}){\xi }_j(t),$$

$$a_k(\mathbf{x})=\sum _{j=1}^d(T_{jk}-T_{kj})=x_k\frac{{\pi }_k(\mathbf{x})-\bar{\pi }(\mathbf{x})}{\bar{\pi }(\mathbf{x})},$$

$${\hat{T}}_{ik}^l=\frac{n_i^l}{N}\left(\frac{{\pi }_k^l{n}_k^l+\mu ⟨{\sigma }_k{n}_k⟩}{{\sum }_{k=1}^d({\pi }_k^l{n}_k^l+\mu ⟨{\sigma }_k{n}_k⟩)}\right)=\frac{n_i^l}{N}\left(\frac{{\pi }_k^l{n}_k^l+\mu ⟨{\sigma }_k{n}_k⟩}{N({\bar{\pi }}^l+\mu \overline{⟨\sigma ⟩})}\right),$$

$${\bar{T}}_{ik}^l=\frac{n_i^l}{N}\left(\frac{{\bar{\pi }}^l}{{\bar{\pi }}^l+\mu }\frac{{\pi }_k^l{n}_k^l}{{\bar{\pi }}^{l}N}+\frac{\mu }{{\bar{\pi }}^l+\mu }\frac{⟨n_k⟩}{N}\right).$$

$${\dot{x}}_k^l={\hat{a}}_k({\mathbf{x}}^l,⟨\mathbf{x}⟩)+\sum _{j=1}^{d-1}{\hat{c}}_{kj}({\mathbf{x}}^l,⟨\mathbf{x}⟩){\xi }_j(t),$$

$${\hat{a}}_k({\mathbf{x}}^l,⟨\mathbf{x}⟩)=\frac{x_k^l({\pi }_k({\mathbf{x}}^l)-\bar{\pi }({\mathbf{x}}^l))+\mu (⟨x_k⟩-x_k^l)}{\bar{\pi }({\mathbf{x}}^l)+\mu },$$

$$\dot{\rho }(\mathbf{x})=-\nabla \left\{\hat{\mathbf{a}}(\mathbf{x},⟨\mathbf{x}⟩)\rho (\mathbf{x})\right\}+\frac{1}{2}{\nabla }^2\left\{\hat{\mathbf{b}}(\mathbf{x},⟨\mathbf{x}⟩)\rho (\mathbf{x})\right\},$$

$$\left(\begin{array}{cc}b-c& -c\\ b& 0\end{array}\right),$$

$$\left(\begin{array}{ccc}0& -ϵ& ϵ\\ ϵ& 0& -ϵ\\ -ϵ& ϵ& 0\end{array}\right),$$

$$\left(\begin{array}{ccc}Rm& Sm& Rm\\ Tm& Pm& T+P(m-1)\\ Rm-c& S+P(m-1)-c& Rm-c\end{array}\right),$$

$${\hat{T}}_{ik}^{ml}=\frac{n_i^{ml}}{N}\left(\frac{{\pi }_k^{ml}{n}_k^{ml}+{\mu }^{(1)}{⟨}^{(1)}{\sigma }_k^{m{l}^{\prime }}{⟩}_{l^{\prime }}+{\mu }^{(2)}⟨{⟨}^{(2)}{\sigma }_k^{m^{\prime }{l}^{\prime }}{⟩}_{l^{\prime }}{⟩}_{m^{\prime }}}{{\sum }_{k=1}^d({\pi }_k^{ml}{n}_k^{ml}+{\mu }^{(1)}{⟨}^{(1)}{\sigma }_k^{m{l}^{\prime }}{⟩}_{l^{\prime }}+{\mu }^{(2)}⟨{⟨}^{(2)}{\sigma }_k^{m^{\prime }{l}^{\prime }}{⟩}_{l^{\prime }}{⟩}_{m^{\prime }})}\right),$$

Full text

I Introduction

The dynamics of Darwinian evolution is intrinsically frequency dependent, the fitness of individuals is tightly coupled to the type and number of competitors. Evolutionary dynamics acts, however, on populations, not individuals and as a consequence depends on not only population composition, but also population size and structure. Evolutionary game theory came about as the result of the realization that frequency dependent fitness introduces strategic aspects to evolution Fisher ; JMS ; HofbauerBOOK . More recently the investigation of the evolutionary dynamics of structured populations, where individuals only compete with some subset of the population, e.g. their neighbors in space or more generally in some graph NowakNAT92 ; LiebermanNAT , has lead to the recognition that the success of different strategies can be greatly influenced by the topology of interactions within the population. Fundamental differences were found – compared to well-mixed populations, where individuals interact with randomly chosen partners – in models that describe the evolution of cooperation (variants of the prisoner’s dilemma game NowakNAT92 ; Axelrod ; AxelrodBOOK ; NowakPNAS94 ; HauertPRL02 ) or deal with the maintenance of biodiversity in the context of competitive cycles (variants of the rock-paper-scissors game HofbauerBOOK ; KerrNAT ; NowakNATnw ; KirkupNAT ; CzaranPNAS ; LenskiPNASComm ).

In order to investigate the coevolutionary dynamics of games on structured populations the full set of connections between a potentially very large number of individuals must be specified. This is only possible by reducing the number of degrees of freedom considered, either through postulating a highly symmetrical (such as lattices NowakNAT92 ; NowakPNAS94 ; SzaboPRE98 ; SzaboClique ; GTphy ; Hui07 ; Hui05 ) or fundamentally random connection structure (such as some random graph ensemble SantosPRL ; OhtsukiNAT06 ). The question of how one goes about the task of reducing the number of degrees of freedom – of choosing the relevant parameters to describe the population structure constrained to which individuals undergo evolution – is not trivial. Both the explicit spatial as well as the random graph ensemble approach have clear precedents in condensed matter physics and network theory, respectively. It is not, however, clear which – if either – approach best describes natural populations of cyclically competing species or societies composed of individuals playing the prisoner’s dilemma game.

As an example let us consider colicin producing bacteria, that play the so called ”rock-paper-scissors” (RPS) game (for details see below). This system has recently been the subject of two experimental studies aimed at demonstrating the role of structured populations in the maintenance of diversity. In the first study KerrNAT ; NowakNATnw bacteria were cultured in vitro in Petri dishes, effectively restricting competition between bacteria to neighbors on the (2D) Petri dish surface (Fig.1 top left), while in the second experiment KirkupNAT in vivo bacterial colonies were established in co-caged mice and their development was subsequently followed. In the case of the first experiment the analogy with explicit 2D spatial embedding (present by construction) is clear (Fig.1 bottom left). The population structure of the second experiment is, however, clearly different. The bacteria in individual mice can be readily considered as locally well-mixed populations, the coevolutionary dynamics of which reduces in the standard meanfield limit to a system of non-linear differential equations (the adjusted replicator equations TraulsenPRL ). As the experiments show, however, migration of bacteria between mice may also occur – resulting in the observed cyclic presence of the three strains in individuals. There are two distinct scales of mixing present in the system. Bacteria within each mice compete with each other forming local populations – an unstructured neighborhood composed of individual bacteria, while also being exposed to migrants from mice with whom they share the cage, together forming a global population – an unstructured neighborhood composed of individual local populations (Fig.1 top and bottom right). This setup is referred to in the ecology literature – albeit in significantly different contexts – as a ”structured metapopulation” MetapopBOOK ; Hanski where structured here refers to the detailed consideration of the population dynamics of the individual populations (often called ”patches”) comprising the metapopulation and is also related to the finite island models of population genetics Pannell .

The above example of co-caged mice is not unique, we may readily think of other ecological or sociological examples where an approximation with hierarchical scales of mixing with no internal structure can be relevant (such as human societies with two distinct scales of mixing present, the first within individual nations the between them at an international level). We have, also, recently used a similar approach to construct a model of genetic exchange among bacteria of the same species (the bacterial equivalent of sex) with which we were able to take into account the effects of spatial and temporal fluctuations in a manner that can explain the benefit of such genetic exchange at the level of the individual SzollosiGEN .

In this paper we construct a hierarchical meanfield theory where the two distinct (i.e. local and global) scales of mixing are each taken into account in terms of two separate meanfield approximations and fluctuations resulting from finite population size on the local scale of mixing are also considered. We subsequently explore the similarities and differences between this and other models of structured populations in the case of the ”rock-paper-scissors” and prisoner’s dilemma games. Through these examples we suggest that our approach allows the separation of the effects of structured populations on coevolutionary dynamics into effects which are highly sensitive to and dependent on the details of the topology and those which only require the minimal structure present in our approximation and can consequently (in terms of sensitivity to the details of the topology) be considered more robust.

II Hierarchical Meanfield Theory for Two Distinct Scales

Let us consider an evolutionary game between $d$ types (strategies) described by the $d\times d$ payoff matrix $A$ with elements $\alpha_{kj}$. Assuming finite and constant population size, natural selection can be described at the level of the individual by the so called the Moran process Moran , during which at each time step an individual is selected randomly from the population to be replaced (death) by the offspring of an individual that is chosen proportional to its fitness to reproduce (birth). This models a population in equilibrium, where the time scale of the population dynamics is set by the rate at which ”vacancies” become available in the population. The fitness of each individual depends on the payoff received from playing the game described by $A$ with competitors (an individual of type $k$ receiving a payoff $\alpha_{kj}$ when playing with an individual of type $j$). In well-mixed populations, individuals can be considered to come into contact (compete) with equal probability with any member of the population excluding themselves – this allows one to calculate the fitness of an individual of type $k$ in a meanfield manner, yielding

[TABLE]

where $n_{k}$ is the number of individuals of type $k$ in the population, $\sum_{k=1}^{d}n_{k}=N$ is the size of the population, $\pi_{\rm base}$ is some baseline fitness and the Kronecker delta symbol $\delta_{kj}$ is equal to unity if $k=j$ and is zero otherwise. From this we may calculate the transition probabilities of our stochastic process, i.e., the probability of an individual of type $i$ being replaced by an offspring of an individual of type $k$ is given by

[TABLE]

where $\bar{\pi}=\sum_{k=1}^{d}\pi_{k}n_{k}/N$. The state of any population is completely described by the frequency of the different strategies $x_{k}=n_{k}/N$. Due to the normalization $\sum_{k=1}^{N}x_{k}=1$, the values of $x_{k}$ are restricted to the unit simplex $S_{d}$ HofbauerBOOK . For $d=2$ this is the interval $[0,1]$, $S_{3}$ is the triangle with vertices $\{(1,0,0),(0,1,0),(0,0,1)\}$ while $S_{4}$ is a tetrahedron etc.

As Traulsen *et al. *have recently shown TraulsenPRL ; TraulsenPRE for sufficiently large, but finite populations the above stochastic process can be well approximated by a set of stochastic differential equations combining deterministic dynamics and diffusion (population drift) referred to as Langevin dynamics:

[TABLE]

where the effective deterministic terms $a_{k}({\bf x})$ are given by

[TABLE]

$c_{kj}({\bf x})$ are effective diffusion terms, that can also be expressed in terms of the transition probabilities as described in TraulsenPRE , and $\xi_{j}$ are delta correlated $\langle\xi_{k}(t)\xi_{j}(t^{\prime})\rangle=\delta_{kj}\delta(t-t^{\prime})$ Gaussian white noise terms. As $N\to\infty$ the diffusion term tends to zero as $1/\sqrt{N}$ and we are left with the modified replicator equation.

In the context of our hierarchical mixing model the topology of connections can be described by two parameters, the populations size at the local scale of mixing $N$, and a second parameter $\mu$, which tunes the strength of global mixing relative to the local dynamics. We take into account the second (global) scale of mixing – mixing among local populations – by introducing a modified version of the Moran process. In the modified process a random individual is replaced at each time step either with the offspring of an individual from the same population (local reproduction) or with an individual from the global population (global mixing). This is equivalent to considering the global population to be well-mixed at the scale of local populations.

Let us consider a global population that is composed of $M$ local populations of size $N$. In each local population vacancies become available that local reproduction and global mixing compete to fill. In any local population $l$ the probability of an individual of some type $k$ filling a new vacancy due to local reproduction must be proportional to the number of individuals of type $k$ multiplied by their fitness i.e. $\pi_{k}^{l}n_{k}^{l}$, where we consider $\pi_{k}^{l}$ to be determined only by interactions with individuals in the same local population according to equation (1). To describe the tendency of individuals of some type $k$ in local population $l$ to contribute to global mixing we introduce the parameters $\sigma_{k}^{l}$. The choice of appropriate $\sigma_{k}^{l}$ depends on the details of the global mixing mechanism, for systems where only the offspring of individuals mix globally it is proportional to the fitness of a given type, while for mechanisms such as physical mixing, by e.g. wind or ocean currents, it may be identical for each type. Irrespective of the details, however, the probability of an individual of some type $k$ filling in a new vacancy due to global mixing should be proportional to the global average of the number of individuals of type $k$ multiplied by their mixing tendency, which we denoted as $\langle\sigma_{k}n_{k}\rangle=\sum_{l=1}^{M}\sigma_{k}^{l}n_{k}^{l}/M$, and the strength of global mixing $\mu$. These consideration lead to the new transition probabilities:

[TABLE]

where $\bar{\pi}^{l}=\sum_{k=1}^{d}\pi^{l}_{k}n^{l}_{k}/N$ and $\bar{\langle\sigma\rangle}=\sum_{k=1}^{d}\langle\sigma_{k}n_{k}\rangle/N$.

We have found that the results presented below are qualitatively the same for both the fitness dependent choice of $\sigma_{k}^{l}=\pi_{k}^{l}$ and the fitness independent choice of $\sigma_{k}^{l}=1$. Therefore, in the following we restrict ourselves to the somewhat simpler fitness independent choice of $\sigma_{k}^{l}=1$, which can be considered to correspond to some form of physical mixing mechanism. The transition probabilities (5) then reduce to:

[TABLE]

We can see that after a vacancy appears either local reproduction occurs, with probability $\bar{\pi}^{l}/(\bar{\pi}^{l}+\mu)$, or global mixing, with probability $\mu/(\bar{\pi}^{l}+\mu)$. From (6) we may derive the Langevin equation describing the coevolutionary dynamics of population $l$ from the

[TABLE]

with the modified deterministic terms given by

[TABLE]

where the vector $\langle{\bf x}\rangle=\sum_{l=1}^{M}{\bf x}^{l}/M$ with components $\langle x_{k}\rangle=\sum_{l=1}^{M}x^{l}_{k}/M$ describes the frequencies of the individual types in the global population and the diffusion terms $\hat{c}({x}^{l},\langle{\bf x}\rangle)$ can be expressed in terms of the modified transition probabilities $\hat{T}^{l}_{ik}$ as above.

Equations (7) describe the coevolutionary dynamics of the global population through the coupled evolution of the $\{{\bf x}^{1},\dots,{\bf x}^{M}\}$ local populations. In the limit of a large number of local populations ($M\to\infty$) the distribution of the local populations over the space of population states (the simplex $S_{d}$) is described by a density function $\rho({\bf x})$ that is normalized over $S_{d}$, i.e., $\int_{S_{d}}\rho({\bf x})=1$. The time evolution of $\rho(\bf x)$ follows a $d-1$ dimensional advection-diffusion equation – the Fokker-Planck equation corresponding to eq. (7):

[TABLE]

with the global averages $\langle x_{k}\rangle=\int_{S_{d}}x_{k}\rho({\bf x})$ coupled back in a self-consistent manner into the deterministic terms $\hat{a}_{k}({\bf x},\langle{\bf x}\rangle)$ and the diffusion matrix $\hat{b}_{kj}({\bf x},\langle{\bf x}\rangle)=\sum_{i=1}^{d-1}\hat{c}_{ki}({\bf x},\langle{\bf x}\rangle)\hat{c}_{ij}({\bf x},\langle{\bf x}\rangle)$. For large local populations ($N\to\infty$) the diffusion term vanishes as $1/N$.

The above advection-diffusion equation (9) presents an intuitive picture of the coevolutionary dynamics of the population at a global scale. We can see that local populations each attempt to follow the trajectories corresponding to the deterministic replicator dynamics, while under the influence of two additional opposing forces: (i) global mixing, which attempts to synchronize local dynamics and (ii) diffusion resulting from finite population size effects, which attempts to smear them out over the simplex. The strength of these forces are tuned by two parameters $\mu$ and $N$, respectively.

If, further, the effects of synchronization are irrelevant, as for example in the case of populations where selection is externally driven by independent environmental fluctuations, we may replace the global population average with the time average of any single population. This is the approach we used in our study of genetic mixing in bacteria SzollosiGEN .

During our numerical investigations we found solving the advection-diffusion equation (9) numerically challenging, particularly in the $N\to\infty$ limit. We resorted instead to solving the coupled Langevin equations (7) for large $M=10^{4}-10^{5}$ to simulate the time evolution of $\rho({\bf x})$.

III Cooperation in populations with hierarchical levels of mixing

The evolution of cooperation is a fundamental problem in biology, as natural selection under most conditions favors individuals who defect. Despite of this, cooperation is widespread in nature. A cooperator is an individual who pays a cost $c$ to provide another individual with some benefit $b$. A defector pays no cost and does not distribute any benefits. This implies the payoff matrix

[TABLE]

where $b$ is the benefit derived from playing with a cooperator while $c$ is the cost for cooperation. From the perspective of evolutionary game theory, which equates payoff with fitness, the apparent dominance of defection is simply the expression of the fact that natural selection a priori selects for fitness of individuals and not the fitness of groups.

Defection dominates cooperation in any well-mixed population HofbauerBOOK . Population structure induced by spatial structure NowakNAT92 ; GTphy and more general networks of interactions SantosPRL ; OhtsukiNAT06 ; Csermely ) has, however, been found to facilitate the emergence and maintenance of cooperation. The mechanism responsible, termed spatial, or more generally, network reciprocityNowakSCI06 depends strongly on the details of local topology. In particular, it seems that lattice like connectivity structures where three-site clique percolation occurs SzaboClique and more general interaction graphs where the degree of nodes $k$ does not exceed the ratio of benefit to cost (i.e. $k<b/c$) OhtsukiNAT06 are required for cooperation to be favored.

Examining the effects of hierarchical mixing on the evolutionary dynamics of cooperation we found that a sharp, but continuous transition leads to the dominance of cooperation as the benefit to cost ratio becomes smaller then the local population size, i.e. $b/c<N$. If the benefit to cost ratio is larger then the local population size the global population is dominated by defectors. The mechanism leading to the dominance of cooperation arises due to the competition between local reproduction and global mixing. In local populations with lower average fitness – larger number of defectors – the influx of individuals from the global scale will be larger than in local populations with higher average fitness (cf. eq. (6) where the relative strength of the two terms on the left hand side depends on the sum of the average fitness of population $l$ and $\mu$). The crucial ingredient for cooperation to be successful is population drift introduced by finite local population size. It is biased influx coupled with drift that can result in cooperation being favored in the global population (Fig 2.).

IV The RPS Game

To explore the effects of hierarchical mixing in the context of games with three strategies we first turn to the case of the so called ”rock-paper-scissors” (RPS) game. In the original popular version of the game two players are afforded the chance to simultaneously display either rock (fist), paper (flat hand) or scissors (two fingers). If player one displays a flat hand while player two displays a fist, player one wins as paper wraps rock. Similarly scissors cut paper, and rocks smashes scissors. Several examples of this game have been found in nature (e.g. among lizards SinervoNAT ), but it is bacteria that have received the most experimental and theoretical attention.

In ecology the often high diversity among microbial organisms in seemingly uniform environments, referred to as the ”paradox of the plankton”, has been difficult to understand. Several models based on spatially explicit game theoretical models have been proposed to explain this diversity CzaranPNAS ; LenskiPNASComm ; KerrNAT ; NowakNATnw . These models are all variants of the RPS game played by colicin producing bacteria. Colicins are antibiotics produced by some strains of Echerichia coli. In experiments (see Fig.1) typically three strains are used: colicin producing (C), sensitive (S) and resistant (R). The coevolutionary dynamics of the three strains can be cast in terms of an RPS game, C strains kill S strains, but are outcompeted, by R strains, because toxin production involves the suicide of bacteria. The cycle is closed by S strains that outcompete R strains, because resistance requires mutant versions of certain membrane protein, which are less efficient than the wild type KerrNAT . Despite the cyclic dynamics colicin-producing strains cannot coexist with sensitive or resistant strains in a well-mixed culture, yet all three phenotypes are recovered in natural populations. Local dispersal (modeled as explicit spatial embedding) has widely been credited with promoting the maintenance of diversity in this system KerrNAT ; NowakNATnw ; CzaranPNAS ; LenskiPNASComm .

In its most symmetric form the RPS game is described by the payoff matrix

[TABLE]

and some $\pi_{\rm base}>\epsilon$. The dynamics of this game in an infinitely large well mixed population consists of neutral orbits along which the product $x_{\rm R}x_{\rm P}x_{\rm S}$ is conserved. For any finite $N$, however, fluctuations lead to the inevitable extinction of all but one of the strategies ReichenbachPRE . Spatial population structure can avert this reduction in diversity CzaranPNAS ; KerrNAT through the emergence of a stable fixed point at the barycentre of the simplex . The effect of the gradual randomization of different lattice topologies (where a small number of edges are randomly rewired) on the dynamics of the game has also been investigated. A Hopf bifurcation leading to global oscillations was observed SzolnokiPRE ; SzaboJPA as the fraction of rewired links was increased above some critical value.

Examining the dynamics of the symmetric RPS game in terms of our hierarchical meanfield approximation we observed that an internal fixed point emerged for $N\to\infty$ (Fig.3a). More importantly, diversity was also maintained for finite local population sizes if global mixing was present. Simulations of the time evolution of $\rho({\bf x})$ also revealed a Hopf bifurcation leading to the oscillation of the global average as $\mu$ was increased above a critical value $\mu_{\rm c}$ depending on $N$ (Fig.4a). These results show that previous results obtained from simulations of populations constrained to different lattice topologies can be considered universal in the sense that not only lattices, but any population structure that can be approximated by two distinct internally unstructured scales of mixing are sufficient for their existence. In the context of the ”paradox of the plankton” these results imply that aside of local dispersal (modeled as explicit spatial embedding) a minimal metapopulation structure (with local competition and global migration) can also facilitate the maintenance of diversity in cyclic competition systems.

V The Repeated Prisoner’s Dilemma Game

In the general formulation of the prisoner’s dilemma (PD) game, two players have the choice to cooperate or to defect. Both obtain some payoff $R$ for mutual cooperation and some lower payoff $P$ for mutual defection. If only one of the players defects, while the other cooperates, the defector receives the highest payoff $T$ and the cooperator receives the lowest payoff $S$. That is $T>R>P>S$ and defection dominates cooperation in any well-mixed population. New strategies become possible, however if the game is repeated, and players are allowed to chose whether to defect or cooperate based on the previous actions of the opponent. In the following we consider, similar to refs. NowakNAT04 and Imhof that recently examined the role of finite population size and mutation and finite population size, respectively in terms of the repeated PD game with three strategies: always defect (ALLD), always cooperate (ALLC), and tit-for-tat (TFT). TFT cooperates in the first move and then does whatever the opponent did in the previous move. TFT has been a world champion in the repeated prisoner’s dilemma ever since Axelrod conducted his celebrated computer tournaments AxelrodBOOK , although it does have weaknesses and may be defeated by other more complex strategies Molander .

Previous results indicate that if only the two pure strategies are present (players who either always defect or ones who always cooperate) explicit spatial embedding NowakNAT92 and some sufficiently sparse interaction graphs OhtsukiNAT06 ; Taylor07 allow cooperation to survive and the behavior of populations is highly sensitive to the underlying topology of the embedding SzaboClique . We have found that introducing global mixing into the PD game with only the two pure strategies present also allows cooperation to survive. The mechanism responsible for favoring cooperation in this case, however, depends on the details of the competition between local reproduction and global mixing. For more than two strategies these details are much less relevant and do not qualitatively influence the dynamics. We will, therefore, consider the delicate issues concerning the PD game with only the two pure strategies in a separate publication, and concentrate here on the repeated PD game with three strategies.

To investigate the effect of global mixing on the repeated PD game with three possible strategies: ALLD, ALLC and TFT following Imhof *et al. * Imhof we considered the payoff matrix:

[TABLE]

where the strategies are considered in the order ALLC, ALLD, TFT, $m$ corresponds to the number of rounds played and $c$ to the complexity cost associated with conditional strategies (TFT). The dynamics of this game has a single unstable internal fixed point and the state where each member of the population plays ALLD is the only nontrivial stable equilibrium (Fig.3b).

Introducing global mixing, between local well-mixed populations, however, causes new stationary states to emerge . Three phases can be identified: (i) ALLD wins (ii) large fraction of local populations in the ALLD corner maintains local cycles of cooperation defection and reciprocity through providing an influx of defectors that prevent TFT players from being outcompeted by ALLC playing individuals (iii) a self maintaining internal globally oscillating cycle emerges. The simplest scenario of two ($M=2$) deterministic ($N\to\infty$) local populations coupled by global mixing ($\mu>0$) already leads to the emergence of phase (ii) as demonstrated in Fig.3b while phase (iii) only emerges for larger $M$. For larger $M$ simulations show that in the limit of large local populations all global configurations with less than some maximum ratio of the populations $I$ on the internal cycle are stable in phase (ii). A transition from phase (ii) to (i) happens as $\mu$ is decreased below a critical value $\mu^{\rm ii\to i}_{\rm c}$ and $I$ approaches zero as $I=(1-\mu^{\rm ii\to i}_{\rm c}/\mu)$ (data not shown). This can be understood if we considered that near the transition point a critical proportion $C=\mu(1-I)$ of ALLD individuals needs to arrive to stabilize local cycles of cooperation defection and reciprocity. At the critical point $I=0$ and $\mu=\mu^{\rm ii\to i}_{\rm c}$ which implies $C=\mu^{\rm ii\to i}_{\rm c}$ giving $I=(1-\mu^{\rm ii\to i}_{\rm c}/\mu)$

Exploring the $N-\mu$ phase space (Fig.5) we see that the transition from phase (i) to (ii) becomes discontinuous for finite $N$ (Fig.4d,e). Further, for any given value of $N$ and $\mu$ the global configuration is described by a unique $I$ due to the presence of diffusion. For appropriate values of the parameters the global average converges to a stationary value in phase (ii) similarly to case of explicit spatial embedding (Fig.4c).

For very small ($\mu\ll\pi_{k}$ for all $k$) and very large ($\mu\ll\pi_{k}$) values of $\mu$ the global dynamics can be reduced to that of some well-mixed population where only ALLD persists (Fig 5.). For small $N$ we again have an effective well-mixed population – the only limit were defectors do not dominate is $N\to\infty$. In comparison with previous results of Imhof et al. we can see that evolutionary cycles of cooperation defection and reciprocity can be maintained not only by mutation, but also by population structures with hierarchical levels of mixing.

VI Discussion

While it is, of course, clear that the reduction of any realistic population structure to a manageable construction is always an approximation, it has not been clearly established what the relevant degrees of freedom are in terms of evolutionary dynamics. Meanfield approximations are a classic method of statistical and condensed matter physics and are routinely used to circumvent intractable combinatorial problems which arise in many-body systems. Cluster-meanfield approximations of sufficient precision GTphy ; Hui07 have been developed that adequately describe the evolutionary dynamics of explicitly structured populations through systematically approximating the combinatorial complexity of the entire topology with that of small motif of appropriate symmetry. The effects of more minimal effective topologies have, however, not been investigated previously. In the above we have shown that straightforward hierarchical application of the meanfield approximation (the assumption of a well-mixed system) surprisingly unveils a new level of complexity.

In the broader context of ecological and population genetics research on structured populations our model can be described as a metapopulation model. The term ’metapopulation’ is, however, often used for any spatially structured population Hanski , and models thereof. More restrictive definitions of the term are often implied in the context of ecology and population genetics literature.

The foundations of the classic metapopulation concept where laid down by Levin’s vision of a ”metapopulation” as a population of ephemeral local populations prone to extinction. A classic metapopulation persists, like an ordinary population of mortal individuals, in a balance between ’deaths’ (local extinctions) and ’births’ (establishment of new populations at unoccupied sites) Hanski . This classic framework is most wide spread in the ecology literature, a less often employed extension is the concept of a structured metapopulation where the state of the individual populations is considered in more detail, this is more similar to our concept of hierarchical mixing, but differs in considering the possibility of local extinctions.

The effects of finite population size and migration, which our model considers, has been of more central concern in the population genetics literature. The analog of Levin’s classic metapopulation concept is often referred to as the ’finite-island’ model Pannell the effective population genetic parameters describing which, have been explored in detailwhitlock . The study of the population genetics of spatially subdivided populations in fact predates Levin, Wright having emphasised the capacity of drift in small populations to bring about genetic differentiation in the face of selection and/or migration several decades priorPannell .

Our hierarchical mixing model treats the coevolutionary dynamics of evolutionary games on structured populations in a manner similar to the most simple population genetic models of spatially subdivided populations, focusing on the parallel effects of selection, drift and migration. It goes beyond these models both in considering the effects of frequency dependent selection (and the strategic aspects of the evolutionary dynamics this implies) and in using a self-consistent approach to describe the global state of the subdivided population. Also, in order to maintain a connection with previous work on the effects of spatial structure on evolutionary games, which rely on Nowak’s concept of spatial games NowakNAT92 , with individuals restricted to interact, and hence compete, only with neighbours as defined by some topology of interaction, we develop our model from the level of the individual by introducing a modified version of the Moran process – and not by extending the Wright-Fisher process (which considers discrete generations and binomial sampling to account for finite population size). The effective population structure described by our hierarchical mixing model can be thought of as a population of individuals, interactions among which are specified by the edges of a hierarchically organized random graph. The fundamental difference in our picture is that the edges of this graph of interactions are not considered to be fixed, but are instead in a constant state of change, being present with a different probability between pairs of individuals who share the same local population and between pairs of individuals who do not (Fig.1.). We consider annealed randomness, which in contrast to the usual quenched picture of fixed edges is insensitive to the details of topology. Our approach we believe best facilitates the exploration of the effects of changing the relative strengths of drift and migration in the context of evolutionary games on structured populations.

Examining the effects of hierarchical mixing in the context of the evolution of robustness we demonstrated that biased influx coupled with drift can result in cooperation being favored, provided the ratio of benefit to cost exceeds the local population size. This result bears striking resemblance to that of Ohtsuki et al. OhtsukiNAT06 , who were able to calculate the fixation probability of a randomly placed mutant for any two-person, two-strategy game on a regular graph and found that cooperation is favored provided the ratio of benefit to cost exceeds the degree of the graph. Our results demonstrate that this rule extends to the minimal spatial structure induced by hierarchical levels of mixing.

Applying our model of spatial structure to the repeated prisoners dilemma revealed that a constant influx of defectors can help to stabilize cycles of cooperation, defection, and reciprocity through preventing the emergence of an intermittent period of ALLC domination in the population, which would present a situation that ”leaves the door wide open” to domination by defectors. While previous work has been done on the effects of ”forcing” cooperation SzaboForcing the idea that an influx of defectors can in fact stabilize the role of reciprocity in promoting cooperation has not been proposed previously. It seems highly unlikely that this mechanism can be explained in terms of kin or multilevel (group) selection, the similarities between which in structured populations have recently been the subject of intensive debate (see e.g. Killingback and Grafen or TraulsenML and Lehmann ). Kin selection can operate whenever interactions occurring among individuals who share a more recent common ancestor than individuals sampled randomly from the whole population Lehmann are relevant. In our case it is the interaction between defectors, arriving from the global scale, and TFT players present at the local scale that is important, and not the interaction between individuals in the local population, who may be thought of as sharing a recent common ancestor due to local dispersal. Also, while the concept of multilevel selection presents a promising framework for the study of evolution of cooperation, it must nonetheless be possible to derive it from ”first principles” – just as kin selection can be cast as an emergent effect of local dispersal.

While there has been considerable work on studying the evolutionary games on graphs and highly symmetric spatial structures very little attention has been paid to the effects of more minimal effective population structures, despite their widespread application in ecology and population genetics, fields from which evolutionary game theory was born and must ultimately reconnect with. We believe that the minimal population structure that such a hierarchical meanfield theory describes is potentially more relevant in a wide range of natural systems, than more subtle setups with a delicate dependence on the details and symmetries of the topology. We showed through two examples that such structure is sufficient for the emergence of some phenomena previously only observed for explicit spatial embedding, demonstrating the potential of our model to identify robust effects of population structure on the dynamics of evolutionary games that do not depend on the details of the underlying topology. The practical advantage of our approach, lies in its ability to readily determine whether or not some feature of a structured population depends on the topological details of local interactions.

Recent simulation result concerning the dynamics of public goods games on different population structures HauertSCI02 ; HauertPRL02 and experiments where global mixing in an RPS like bacteria-phage system lead to the emergence of a ”Tragedy of the commons” scenario KerrTCNAT should all be amicable to analysis in terms of our method.

VII Acknowledgments

This work was partially supported by the Hungarian Scientific Research Fund under grant No: OTKA 60665.

VIII Appendix

Our approach readily generalizes for an arbitrary number of hierarchical mixing levels. For three levels of mixing we may consider the global population to be comprised of $\mathcal{M}$ subpopulations each of which is in turn subdivided into $M$ local populations. With $m\in\{1,\cdots,\mathcal{M}\}$ running over subpopulations and $l\in\{1,\cdots,M\}$ over local populations the transition probabilities can be written as:

[TABLE]

where primed indices indicate the scale of mixing over which the average is taken, $\mu^{(1)}$ describes the strength of mixing, and the ${}^{(1)}\sigma_{k}^{ml}$ the tendencies of mixing among local populations within a subpopulation, while $\mu^{(2)}$ describes the strength of mixing, and the ${}^{(2)}\sigma_{k}^{ml}$ the tendencies of mixing among subpopulations in the global population.