## Abstract

It has recently been demonstrated that ecological feedback mechanisms can facilitate the emergence and maintenance of cooperation in public goods interactions: the replicator dynamics of defectors and cooperators can result, for example, in the ecological coexistence of cooperators and defectors. Here we show that these results change dramatically if cooperation strategy is not fixed but instead is a continuously varying trait under natural selection. For low values of the factor with which the value of resources is multiplied before they are shared among all participants, evolution will always favour lower cooperation strategies until the population falls below an Allee threshold and goes extinct, thus evolutionary suicide occurs. For higher values of the factor, there exists a unique evolutionarily singular strategy, which is convergence stable. Because the fitness function is linear with respect to the strategy of the mutant, this singular strategy is neutral against mutant invasions. This neutrality disappears if a nonlinear functional response in receiving benefits is assumed. For strictly concave functional responses, singular strategies become uninvadable. Evolutionary branching, which could result in the evolutionary emergence of cooperators and defectors, can occur only with locally convex functional responses, but we illustrate that it can also result in coevolutionary extinction.

## 1. Introduction

Public goods games for groups of interacting individuals have widely been used to investigate the possibility of emergence and maintenance of cooperation. In these games, individual participants can typically choose how much of their money, or resources in general, they invest into a common pool. These resources are then multiplied with the factor *ρ* > 1 and shared equally among all participants. On the population level, the best outcome is achieved when all participants invest all of their resources into the common good. This situation may, however, be vulnerable to the invasion of defectors. Let *s* denote the investment of a focal individual. The payoff that this focal individual receives from itself in a game with *N* participants is equal to (*ρ*/*N* − 1)*s*. It is thus well known that in a well-mixed model, cooperation is favoured only if *ρ* > *N* (e.g. Hauert *et al.* 2006).

Several mechanisms promoting cooperation in case *ρ* < *N* have been proposed. If the participants of the game are genetic relatives, helping others becomes more beneficial owing to kin selection (Hamilton 1964). Even among non-related individuals, cooperation can be maintained owing to direct reciprocity if the participants are likely to meet again. Various forms of reputation and social norms (Ohtsuki & Iwasa 2006) can result in indirect reciprocity and promote cooperation. Various spatial settings without complete mixing may result in the assortment of cooperators, which can also promote cooperation.

Hauert *et al.* (2006) suggested that ecological feedback mechanisms can promote cooperation. In their set-up, the groups of maximum *N* individual participants are chosen randomly in a well-mixed population, but because the probability of finding a participant is usually less than one, the effective group size *S* is typically less than *N*. If the population size is low such that the effective population size is small enough, cooperation will be favoured, also resulting in the increase in the population size, until the population size becomes so large that increasing cooperation is no longer beneficial. Hauert *et al.* (2006) investigated a natural extension of replicator dynamics (Hofbauer & Sigmund 1998) of defectors and cooperators with a fixed cooperative investment. Depending on parameter values and initial conditions, replicator dynamics may result in the extinction of the whole population, in the coexistence of defectors and cooperators either in a stable equilibrium or in a periodic orbit, or in the presence of cooperators without defectors (see also Hauert *et al.* 2008). They concluded that ecological feedback can maintain cooperation in public goods games.

In this article, we relax the assumption of fixed amounts of cooperative investment, and instead treat them as continuously varying traits under natural selection. Modelling traits as continuous variables has a long history in evolutionary biology. A real biological example of such a continuous public good is bacteria producing siderophores, which turn iron into bioavailable forms (Griffin *et al.* 2004; Kuemmerli *et al.* 2009). A mutation has a quantitative effect on the evolving trait. If a mutant performs better than other present traits, it will become more abundant, and possibly oust the other traits. Adaptive dynamics (Metz *et al.* 1996; Geritz *et al.* 1997, 1998) provides tools for the investigation of long-term evolution of continuous strategies in realistic ecological models.

In addition to analysing continuously varying traits, we make the generalization that payoffs affect the birth rate through a functional response. In §2, we present the generalized model, together with some concepts and methods of adaptive dynamics. In §3, we investigate first the monomorphic population dynamics of the model, and then monomorphic evolutionary dynamics. As some situations cause disruptive selection, we also investigate polymorphic adaptive dynamics. The analytical proofs of various results are given in the electronic supplementary material.

## 2. Model and methods

### (a) Ecological public goods game

Assume that there are *k* different strategies present in the population, and let *x*_{i} be the density of individuals with cooperation strategy *s*_{i}, *i* = 1, 2, … *k*. We assume that individuals are randomly distributed in a habitat, and the population is well mixed. Furthermore, let *z* denote the proportion of empty space, thus
2.1
Following Hauert *et al.* (2006), we assume that individuals are randomly chosen to form a group of at most *N* individuals. This assumption is motivated by events where cooperation is possible occurring locally, which thus involve only a small proportion of the whole habitat. For example, warning signals about an arriving predator are observed (and also useful) only locally. The chance that an individual finds itself in a group of size *S* ≤ *N* is given by the chance to find *S* − 1 interaction partners, where *S* − 1 is binomially distributed, *S* ∼ Bin(*N* − 1, 1 − *z*)
2.2
Provided that *S* − 1 interaction partners have been found, the numbers of partners *X*_{i} with strategy *s*_{i} are multinomially distributed, i.e. for we have
2.3
Equation (2.3) means also that the random variables *X*_{i} are binomially distributed: *X*_{i} ∼ Bin(*S* − 1, *x*_{i}/(1 − *z*)).

Following Hauert *et al.* (2006), we assume that the total amount of invested resources are multiplied with the factor *ρ* and then shared equally among all partners, and this has a positive effect on the birth rate. The cooperation strategy of an individual has also a direct negative linear effect on the birth rate of itself only. Death rate *δ* is assumed to be constant. Generalizing the model of Hauert *et al.* (2006), we assume that the positive effect to the birth rate occurs through a functional response, described by the function *h*(*y*). We assume that it is a positive, continuous and strictly increasing function with *h*(0) = 0. Some biological situations, such as group foraging, could also result in stepwise benefits (Bach *et al.* 2006). By setting *h*(*y*) = *y*, we obtain the original model by Hauert *et al.* (2006). With the assumptions above, for given numbers *α* = (*α*_{1}, … , *α*_{k}) of different types of partners **s** = (*s*_{1}, … , *s*_{k}), the birth rate of an individual with cooperation strategy *s*_{mut} is equal to
2.4

The average birth rate with *S* − 1 interaction partners is obtained by summing over all possible combinations of partners
2.5
where **x** = (*x*_{1}, … , *x*_{k}). In the case of a linear functional response *h*(*y*) = *y*, equation (2.5) can be written as
2.6
because the birth rate function (2.4) is then linear with respect to *α*_{i} and *s*_{i}. The final average payoff is obtained, when equation (2.5) is summed over the distribution of the number of interaction partners.
2.7
The threshold of requiring two individuals for public goods production (Hauert *et al.* 2006) is the mechanism behind the Allee effect in monomorphic population dynamics (figure 1). Replicator dynamics for other threshold values has been studied by Pacheco *et al.* (2009). Like Hauert *et al.* (2006), we assume that reproduction can only occur in empty space. Therefore, the *per capita* population growth rate is
2.8

### (b) Adaptive dynamics

Invasion fitness (Metz *et al.* 1992) is a central concept in adaptive dynamics. It is defined as the long-term exponential growth rate *r*(*s*_{mut}, *E*_{res}) of a rare mutant *s*_{mut} in an environment *E*_{res} set by the resident. If the invasion fitness of a mutant is positive, it is able to grow in population size. Therefore, it may invade and become a new resident itself. Random small mutations are assumed to happen rarely, so that the population dynamics has settled to an attractor before the next mutation happens. These mutation–invasion events result in the change of the strategy of the individuals constituting the population. These events define a trait-substitution sequence. Each element of such a sequence is either a strategy or a set of several strategies, replacing the phenotype(s) that were previously present. Although each trait-substitution sequence is a realization of a stochastic process, there are analytical methods of adaptive dynamics that tell what will happen to all possible trait-substitution sequences in the long run.

If the invasion fitness of any mutant is non-positive in the environment set by a specific resident, then this resident strategy is uninvadable, and called an evolutionarily stable strategy (Maynard Smith 1976). In this case, the resident strategy is necessarily a (local) fitness maximum and thus the selection gradient, i.e. the derivative of invasion fitness with respect to the strategy of the mutant is zero at such points,
2.9
and the second derivative (∂^{2}/∂*s*_{mut}^{2})*r*(*s*_{mut}, *E*_{res})|_{smut = sres} < 0 is negative. Strategies for which the selection gradient is zero are called evolutionarily singular strategies (Geritz *et al.* 1998). Such a singular strategy *s** is convergence stable or an evolutionary attractor if the repeated invasion of nearby mutant strategies into resident strategies will lead to the convergence of resident strategies towards *s** (Christiansen 1991). If an evolutionary attractor is also evolutionarily stable, it is called a continuously stable strategy (Eshel 1983), and it is a feasible final outcome of an evolutionary process. In case the second derivative (∂^{2}/∂*s*_{mut}^{2})*r*(*s*_{mut}, *E*_{res})|_{smut}_{= sres} > 0 is positive for a monomorphically attracting singular strategy, this singular strategy is not unbeatable, and evolutionary branching (Geritz *et al.* 1998) may occur because of disruptive selection. This phenomenon can result in the evolutionary emergence of cooperators and defectors in an initially monomorphic population (Doebeli *et al.* 2004; Brännström & Dieckmann 2005).

We can illustrate the monomorphic adaptive dynamics by drawing pairwise invasibility plots (PIPs) (Matsuda 1985; Van Tienderen & De Jong 1986), where the sign of the invasion fitness is displayed in dependence on resident and mutant strategies (figures 2 and 3*a*). Because the resident population is at a population-dynamical equilibrium, a mutant population with the same strategy as the resident will have a neutral growth rate. Therefore, there is a zero contour at the diagonal *s*_{mut} = *s*_{res} in PIPs. Evolutionarily singular strategies lie at those points where other zero contours cross the diagonal.

## 3. Results

Hauert *et al.* (2006) assumed that there are two strategies present, a defector with *s* = 0 and a cooperator with fixed *s* > 0, and studied their population dynamics and possibilities of coexistence. In this article, we assume instead that the cooperation strategy is a continuously varying trait under natural selection. We will next use adaptive dynamics (Metz *et al.* 1996; Geritz *et al.* 1997, 1998) to study the evolutionary consequences of this assumption.

### (a) Monomorphic population dynamics

Let us begin with a simplifying assumption that there is only one resident phenotype (*k* = 1, *x*_{1} = 1 − *z*). First we need to find the population-dynamical attractors of the resident population. Consider a focal individual with strategy *s*_{mut} = *s*_{res}. Since everybody in the population has the same strategy *s*_{res}, the amount of received benefits (2.4)
3.1
does not depend on the number of players in the game, provided that the game is played. The game is not played if no partners are found for the focal individual, which occurs with probability *z*^{N−1}. In this model, this assumption is the mechanism behind the Allee effect (Allee *et al.* 1949), which means that individuals benefit from the presence of conspecifics, and these benefits are lost as population densities decline. Based on the reasoning above, the *per capita* growth rate (equation 2.8) of a monomorphic resident is
3.2

For a given *s*_{res}, the equilibrium population sizes *x* = 1 − *z* are obtained by finding solutions to *r*(*s*_{res}, *N*, *z*) = 0, see figure 1.

Since *r*(*s*_{res}, *N*, 1) = −*δ* < 0, the extinction equilibrium *z* = 1, *x* = 1 − *z* = 0 is always stable. Using the fact that (1 − *z*^{N−1})*z* reaches its maximum at *z** = *N*^{1/(1−N)}, we can show (see the electronic supplementary material, proposition A.1) that there exist two positive equilibria if and only if
3.3

Strategies, for which the left-hand side of the inequality (3.3) is equal to zero, are boundaries of viability, denoted as *s*_{bif}. Note that the population is not viable for *s*_{res} = 0. The larger equilibrium is stable, and the smaller one is unstable. In the case of a linear functional response *h*(*y*) = *y*, condition (3.3) reduces to
3.4
and the value of the stable positive equilibrium is an increasing function of the resident strategy *s*_{res} (figure 1*a*), see also proposition A.2 in the electronic supplementary material. In the case of Holling type II functional response *h*(*y*) = *y*/(1 + *αy*), the population is viable only with intermediate cooperation strategies *s*_{res}. The population size at the stable equilibrium reaches its maximum with (see figure 1*b,c*).

### (b) Monomorphic adaptive dynamics

Assume now that the resident population with strategy *s*_{res} is at the stable positive equilibrium, where *z* is the proportion of empty space. Then investigate how an individual with strategy *s*_{mut} will perform in the environment set by the resident. The invasion fitness of the mutant (equation 2.8) now becomes
3.5

The fitness gradient *D*(*s*_{res}) = (∂/∂s_{mut})*r*(*s*_{mut}, *s*_{res}, *N*, *z*)|_{smut = sres} gives the direction into which strategies are expected to change with mutations of small size. Although singular strategies *s** can in principle be solved from *D*(*s**) = 0, we do not obtain a simple explicit expression for them. Some general results can, however, be obtained. We can determine whether a singular strategy *s** is locally uninvadable from the sign of the second derivative
3.6
Since *p*_{S}(*N*, *z*) > 0, the sign of equation (3.6) is determined by the sign of the second derivative of the functional response function *h*″(*ρs**). The singular strategy *s** is thus locally uninvadable (evolutionarily stable) if the functional response function *h* is locally concave, *h*″(*ρs**) < 0. In the case of a linear functional response *h*(*y*) = *y*, the second derivative *h*″(*y*) ≡ 0, and invasion fitness (equation 3.5) is a linear function of the strategy of the mutant. In case *h* is locally convex, *h*″(*ρs**) > 0, the singular strategy is not locally uninvadable, and if the singular strategy is also convergence stable, evolutionary branching will occur.

As expected based on the results above, the PIP in figure 2*a* illustrating a case with linear functional response has a vertical contour line, and the convergence stable singular strategy *s** ≈ 0.980 is neutral against mutant invasions. Because for Holling type II functional response *h*″(*y*) < 0, any singular strategy is uninvadable. The PIP in figure 2*b* has an expected shape: there is one convergence stable singular strategy *s** ≈ 0.595, which is uninvadable, and thus an evolutionary endpoint. Also note that these results are robust against large mutational steps.

Surprisingly, these two shapes are not the only possible shapes of a PIP in the model with a linear or concave functional response. As illustrated in figure 3*a*, it is possible that the fitness gradient is negative for all viable strategies, and there is no singular strategy at all. In this case, whatever the strategy of the resident is, a mutant with a lower cooperation strategy has positive fitness, and can at least initially increase in population size. Therefore, the strategy of the evolving population will decrease owing to natural selection (figure 3*b*) until it is near the extinction boundary. But even then, a mutant with a lower strategy value has positive fitness. This mutant is, however, not viable alone, and the invasion by this ‘kamikaze’ mutant will take the population away from the positive equilibrium to the extinction equilibrium (figure 3*c*). In other words, evolutionary suicide (Ferrière 2000) happens.

When does evolutionary suicide occur in the investigated model? We have already noticed that at the extinction boundary *s*_{bif}, we have *z** = *N*^{1/(1−N)}. By substituting this into the fitness gradient
3.7
we obtain that *D*(*s*_{bif}) < 0 and thus evolutionary suicide occurs if and only if
3.8
provided that the initial strategy of the population is near the lower boundary of viability *s*_{bif}. Note that *N*^{N/(N−1)} − *N* > 2 for *N* ≥ 3 and *N*^{N/(N−1)} − *N* = 2 for *N* = 2. For the linear functional response *h*(*y*) = *y*, condition (3.8) becomes *ρ* < *N*^{N/(N−1)} − *N*. In other words, evolutionary suicide occurs for values of *ρ* that are lower than the boundary given by equation (3.8), but still large enough that the population is viable for some strategies (figure 4).

For some types of the functional response function, it is possible, but also quite tedious, to prove that the local evolutionary properties of the lower boundary of viability also determine global evolutionary dynamics in the model. When equation (3.8) holds and thus the lower boundary of viability is locally evolutionarily attracting, it is also globally evolutionarily attracting. Therefore, there are no singular strategies and evolutionary suicide happens with any initial condition. Furthermore, when equation (3.8) does not hold, and thus the lower boundary of viability is evolutionarily repelling, there exists a unique singular strategy, which is globally evolutionarily attracting. Thus, we obtain analytically a full classification of monomorphic adaptive dynamics. The covered types include Holling type I and II, but not Holling type III functional response (Holling 1959).

### Theorem 3.1.

*Assume that the functional response function is either linear (h(y) = y) or concave (h″(y) ≤ 0). Furthermore, assume that with the given parameters, the population is viable for some strategies s*

_{res}

*and thus there necessarily exists a lower boundary of viability s*

_{bif}> 0.

*Then, we have the following classification of evolutionary dynamics:*

(*1*) *If equation (3.8) holds and thus D(s*

_{bif}) < 0,

*the fitness gradient is negative for all viable s*

_{res}.

*Therefore, evolutionary suicide occurs from any initial condition*(

*figure*3).

(*2*) *Otherwise, there exists a unique globally attracting singular strategy. It is uninvadable if the functional response function h is strictly concave* (*figure* 2*b*). *If it is linear h (y) = y, the singular strategy is neutral against mutant invasions* (

*figure*2

*a*).

*Proof.* See the electronic supplementary material. □

How much are we expected to observe cooperation in the investigated model? To answer this question, we have plotted in figure 5 the numerical values of the evolutionarily attracting singular strategies for a range of parameters. We observe that increasing the maximal group size will decrease the singular strategy or even result in evolutionary suicide. Our intuitive explanation for this is that for a given strategy, increasing the maximal group size will increase the population size *x* at equilibrium, as can be seen from equation (3.2). In the new situation, benefits of cooperation are shared among more individuals, which is expected to select for decreased cooperation.

How about the effect of the multiplication factor *ρ*? Intuitively, one could say that increasing *ρ* should always select for increased cooperation, because it directly increases the benefits of cooperation. However, we can see from equation (3.2) that increasing *ρ* will also increase the population size *x* at equilibrium, which is expected to have a negative effect on the benefits of cooperation. Depending on the relative strenghts of these two effects, cooperation may either increase or decrease. As we observe from figure 5, this is indeed the case. For lowest values of *ρ*, either the evolving population is not viable or evolutionary suicide occurs. When the values of *ρ* increase, a singular strategy appears and it first decreases, reaches a minimum and starts to increase again. The final increase is rather steep in case of linear trade-off, but the Holling type II functional response diminishes the highest benefits, and the observed increase in the singular strategy is essentially smaller. Intuitively, one could also say that increasing *α* increases the costs of cooperation, and thus cooperation should decrease. Our numerical explorations suggest that mostly this is the case, but not for all parameter values. Since the population size at equilibrium also decreases, the total effect may increase cooperation (not illustrated).

### (c) Polymorphic adaptive dynamics

So far, we have investigated only monomorphic adaptive dynamics of cooperation, which was relatively easy to investigate because all attractors of a monomorphic resident population in the investigated model are equilibria. This is not the case with a polymorphic resident: Hauert *et al.* (2006) observed that a dimorphic population of cooperators and defectors can coexist ecologically in an equilibrium or in a limit cycle. Therefore, we need to extend fitness calculation to cover limit cycles as well (Metz *et al.* 1992). Let *x*_{i} (*t*) denote the density of individuals with cooperation strategy *s*_{i} at time *t*. Now the polymorphic fitness of a mutant with strategy *s*_{mut} is
3.9

What is then the expected evolutionary outcome? Let us first look at the case when the functional response function is strictly concave, *h*″(*y*) < 0 for all *y* ≥ 0. In this case, the singular strategy, when it exists, is evolutionarily stable. In such a case, there is no need to investigate polymorphic adaptive dynamics in detail because the general theory of adaptive dynamics tells us that selection would cause a polymorphic population to become monomorphic at least if the initial strategies are near the singular strategy. Furthermore, by differentiating equation (3.9) twice, we actually observe that under the present assumptions the polymorphic fitness function is concave with respect to the strategy of the mutant. Therefore, the invasion fitness (equation 3.9) can be zero for at most two strategy values. This means that at most, two strategies can coexist ecologically, and only mutants with intermediate strategies can invade. Therefore, any polymorphic population will become monomorphic.

In the case of a linear functional response *h*(*y*) = *y*, by using equations (2.6) and (2.7), polymorphic invasion fitness (equation 3.9) can be written as
3.10
In other words, not only the monomorphic invasion fitness (equation 3.5), but also the polymorphic invasion fitness (equation 3.9) is a linear function of the strategy of the mutant *s*_{mut} when *h*(*y*) = *y*. This fact holds with any number *k* of resident strategies and any type of resident attractor. Note that necessarily *r*(*s*_{i}, **s**, **x**, *N*, *z*) = 0 for all *i* = 1, 2, … , *k*. Therefore, a direct consequence of the linearity of the fitness function is that *r*(*s*_{mut}, **s**, **x**, *N*, *z*) = 0 for all *s*_{mut} when *k* ≥ 2. In other words, any coalition of several coexisting strategies is neutral against mutant invasions. Therefore, although a monomorphic population can become polymorphic when *ρ* is large enough, there is no disruptive selection to cause evolutionary branching. Note also that any small change in the functional response function making it either locally convex or concave will annihilate the linearity of the fitness function. Therefore, the linear case is structurally unstable.

Let us next investigate the evolutionary dynamics of cooperation in case the functional response function is locally convex, such as the Holling type III function *h*(*y*) = (*y*^{2}/(1 + *αy*^{2})). It is convex for and concave otherwise. Note that it is not realistic for a functional response function to be convex for all *y*. In our numerical investigations, we observed that the effect of the multiplication factor *ρ* (figure 6*a*) is rather similar to the other cases investigated above: for small values of *ρ*, evolutionary suicide occurs, and when *ρ*, is increased, a convergence stable singular strategy appears. This singular strategy is, however, not evolutionarily stable, and thus evolutionary branching is expected to happen (figure 6*b*). Only when *ρ* is further increased, the singular strategy becomes evolutionarily stable.

We can use the polymorphic fitness function (equation 3.9) to investigate what is expected to happen as a result of evolutionary branching. Especially, we can draw the set of protected and unprotected dimorphisms together with the direction of the dimorphic selection gradient (∂/∂*s*_{mut})*r*(*s*_{mut}, **s**, **x**, *N*, *z*)|_{smut = si}. Note that since the numbering of strategies is arbitrary, figure 6*d* is symmetric across the diagonal. Here, an ecologically coexisting strategy pair is called a protected dimorphism if the strategies are mutually invadable (dark grey in figure 6*d*), although the coexistence is not globally protected (extinction equilibrium is always stable). Light grey corresponds to unprotected coexisting strategy pairs for which the less cooperating strategy is not viable alone, and thus its invadability is meaningless.

From figure 6*d* we observe that the strategy of the less cooperating resident is expected to decrease, and the strategy of the more cooperating resident is expected to increase, until the dimorphic strategy pair reaches the boundary of dimorphic coexistence at the left side of the upper left quadrant (or the mirror image). Coexistence mostly occurs in a stable equilibrium, but near this boundary a Hopf bifurcation occurs, and the two resident strategies coexist in a limit cycle. Even then, selection pressure is divergent, and the strategies will evolve to the extinction boundary, where the limit cycle disappears, and the whole population goes extinct. In other words, we are expected to observe an event of dimorphic evolutionary suicide (figure 6*c,e*). For some other parameter values, evolutionary branching may result in evolutionarily stable coexistence of a cooperator and a defector (figure 6*f*, uninvadability not illustrated).

Doebeli *et al.* (2004) introduced the term ‘tragedy of the commune’ for a paradox of cooperation in which ‘selection does not always generate the evolution of uniform and intermediate investment levels but may instead lead to an asymmetric stable state, in which some individuals make high levels of cooperative investment and others invest little or nothing.’ Our results show that such low investment of a part of the evolving population may even result in the extinction of the whole population, which would not occur if the population had remained monomorphic.

## 4. Conclusion

In this article, we have investigated the evolution of cooperation in an ecological model, which is an extended version of the model presented by Hauert *et al.* (2006). They investigated a natural extension of the replicator equation dynamics of defectors and cooperators with a fixed cooperation strategy. In their conclusion, they wrote ‘One extension that seems particularly worthwhile is the continuously varying cooperative strategies, rather than just the all-or-nothing strategies cooperate and defect, and how ecological feedback mechanisms affect the evolutionary dynamics of continuous cooperative investments.’ This is precisely what we have done here, with the additional extension that receiving benefits may also have a nonlinear increasing effect of birth rates (functional response) instead of a linear effect.

We indeed observed interesting new insights and scenarios as Hauert *et al.* (2006) expected, to the extent that they change the interpretation of the main results by Hauert *et al.* (2006): although eco-evolutionary feedback can in principle stabilize cooperation at intermediate frequencies of cooperators and defectors, we observed that in the original model with linear functional response, the ecological coexistence of cooperators and defectors is only neutrally stable against invasions by mutants. Furthermore, the coexistence of cooperators and defectors will not emerge starting from an initially monomorphic population because monomorphic singular strategies are neutrally evolutionarily stable. A convergence stable singular strategy should be invadable in order to observe evolutionary branching.

An even more striking observation is that an initially viable population may evolve towards lower and lower investment strategies until the population reaches a bifurcation point where a stable population-dynamical equilibrium collides with an unstable one. Even then, a mutant with a lower cooperation strategy has positive fitness and can invade. However, this invasion pushes the population below the Allee threshold (Allee *et al.* 1949) to extinction.

Already Hardin (1968) in his article about tragedy of the commons argued that selfish behaviour may lead to the ruin of the common good. The observed evolutionary suicide (Ferrière 2000) in this mathematical model is an extreme example of this phenomenon. Not only does cooperation steadily decrease in the population, but reaching a population-dynamical bifurcation point results in the extinction of the whole population. Evolutionary suicide has been observed in a wide variety of models (Matsuda & Abrams 1994; Gyllenberg & Parvinen 2001; Gyllenberg *et al.* 2002; Webb 2003; Parvinen 2005, 2007; Rankin & López-Sepulcre 2005); see Parvinen (2005) for a review of the phenomenon and presentation of the relevant theory. We observed that evolutionary suicide will occur if the multiplication factor *ρ* is low enough. Note that the assumption that at least two participants are required for the cooperative game to be played is the mechanism behind the Allee effect in this model. Without this, or some other assumption resulting in a discontinuous transition to extinction, evolutionary suicide cannot occur. It has namely been shown by Gyllenberg *et al.* (2002) that if the resident attractor gradually approaches the extinction equilibrium, when the resident strategy approaches the boundary of viability, then evolutionary suicide is not possible. But note that it is not a sufficient condition, as can be seen from figure 2*a,b*.

Analysis of the extension with a nonlinear functional response revealed a structural instability of the linear case. For a strictly concave functional response, any monomorphic singular strategy is properly uninvadable. Furthermore, at most two strategies can coexist ecologically, and only mutants with intermediate strategies can invade. Therefore, any polymorphic population will become monomorphic.

Although a small perturbation can turn a linear function to a convex one, it is not realistic for a functional response function to be globally convex. It can be locally convex, such as Holling type III functional response. If *s** is a singular strategy, and the functional response is convex at *ρs**, thus *h*″(*ρs**) > 0, the singular strategy is invadable by nearby mutants. If the singular strategy is a (monomorphic) evolutionary attractor, it is an evolutionary branching point. Therefore, the strategy of a monomorphic population will initially approach *s**, but eventually a successful mutant will not oust the resident, and the evolving population becomes dimorphic. The evolving population will thus consist of two groups, and the strategies of these groups will at least initially evolve further away from each other. Such evolutionary branching provides a mechanism for the evolutionary emergence of cooperators and defectors (Doebeli *et al.* 2004; Brännström & Dieckmann 2005). However, we illustrated that also such dimorphic population can experience evolutionary suicide. The existence of an evolutionary branching point does not guarantee the evolutionarily stable coexistence of cooperators and defectors.

We conclude that the analysis of ecological models of cooperation is not complete without an evolutionary analysis. Furthermore, evolutionarily stable coexistence of cooperators and defectors require some other mechanisms than those presented by Hauert *et al.* (2006). Ecological dynamics as such is not sufficient, but it may be possible in an ecological model with a mechanism resulting, for example, in direct or indirect reciprocity, assortment of cooperators or kin selection. Also the risk of the population evolving to extinction should be taken into account in empirical studies.

## Acknowledgements

The author wishes to thank anonymous referees for constructive comments and the Academy of Finland for financial support (grant number 128323).

## Footnotes

- Received January 28, 2010.
- Accepted March 16, 2010.

- © 2010 The Royal Society