Evolutionary Games of Multiplayer Cooperation on Graphs.

There has been much interest in studying evolutionary games in structured populations, often modeled as graphs. However, most analytical results so far have only been obtained for two-player or linear games, while the study of more complex multiplayer games has been usually tackled by computer simulations. Here we investigate evolutionary multiplayer games on graphs updated with a Moran death-Birth process. For cycles, we obtain an exact analytical condition for cooperation to be favored by natural selection, given in terms of the payoffs of the game and a set of structure coefficients. For regular graphs of degree three and larger, we estimate this condition using a combination of pair approximation and diffusion approximation. For a large class of cooperation games, our approximations suggest that graph-structured populations are stronger promoters of cooperation than populations lacking spatial structure. Computer simulations validate our analytical approximations for random regular graphs and cycles, but show systematic differences for graphs with many loops such as lattices. In particular, our simulation results show that these kinds of graphs can even lead to more stringent conditions for the evolution of cooperation than well-mixed populations. Overall, we provide evidence suggesting that the complexity arising from many-player interactions and spatial structure can be captured by pair approximation in the case of random graphs, but that it need to be handled with care for graphs with high clustering.

Introduction

Graphs are a natural starting point to assess the role of population structure in the evolution of cooperation. Vertices of the graph represent individuals, while links (edges) define interaction and dispersal neighborhoods. Classical models of population structure, such as island models [1, 2] and lattices [3, 4], often developed before the current interest in complex networks [5, 6], can all be understood as particular instances of graphs [7, 8]. More recently, the popularity of network theory has fueled a renewed interest in evolutionary dynamics on graphs, especially in the context of social behaviors such as cooperation and altruism [7-21].

When selection is weak on two competing strategies, such that fitness differences represent only a small perturbation of a neutral evolutionary process, a surprisingly simple condition for one strategy to dominate the other, known as the “sigma rule”, holds for a large variety of graphs and other models of spatially structured populations [22]. Such a condition depends not only on the payoffs of the game describing the social interactions, but also on a number of “structure coefficients”. These coefficients are functions of demographic parameters of the spatial model and of its associated update protocol, but are independent of the payoffs. In the case of two-player games, the sigma rule depends on a single structure coefficient σ. The larger this σ, the greater the ability of spatial structure to promote the evolution of cooperation or to choose efficient equilibria in coordination games [22]. Partly for this reason, the calculation of structure coefficients for different models of population structure has attracted significant interest during the last years [8, 21–27].

Despite the theoretical and empirical importance of two-player games, many social interactions involve the collective action of more than two individuals. Examples range from bacteria producing extracellular compounds [28-31] to human social dilemmas [32-36]. In these situations, the evolution of cooperation is better modeled as a multiplayer game where individuals obtain their payoffs from interactions with more than two players [37-43]. An example of such multiplayer games is the volunteer’s dilemma, where individuals in a group must decide whether to volunteer (at a personal cost) or to ignore, knowing that volunteering from at least one individual is required for a public good to be provided [44-46]. Importantly, such a multiplayer interaction cannot be represented as a collection of pairwise games, because changes in payoff are nonlinear in the number of co-players choosing a particular action.

Multiplayer games such as the volunteer’s dilemma can also be embedded in graphs, assuming, for instance, that nodes represent both individuals playing games and games played by individuals [47-49]. Most previous studies on the effects of graph structure on multiplayer game dynamics have relied on computer simulations [49]. However, similar to the two-player case, some analytical progress can be made if selection is assumed to be weak. In the multiplayer case, the sigma rule depends no longer on one, but on up to d − 1 structure coefficients, where d is the number of players [50]. Although exact formulas for structure coefficients of multiplayer games can be obtained for relatively simple models such as cycles [51], analysis has proved elusive in more complex population structures, including regular graphs of arbitrary degree. Indeed, extending analytical results on evolutionary two-player games on graphs to more general multiplayer games is an open problem in evolutionary graph theory [52].

Here, we contribute to this body of work by deriving approximate analytical expressions for the structure coefficients of regular graphs updated with a Moran death-Birth model, and hence for the condition of one strategy to dominate another according to the sigma rule. The expressions we find for the structure coefficients suggest that regular graphs updated with a Moran death-Birth model lead to less stringent conditions for the evolution of cooperation than those characteristic of well-mixed populations. Computer simulations suggest that our approximations are good for random regular graphs, but that they systematically overestimate the condition for the evolution of cooperation in graphs with more loops and higher clustering such as rings and lattices. In these cases, cooperation can be no longer promoted, but even be hindered, with respect to the baseline case of a population lacking spatial structure.

Methods

We consider stochastic evolutionary dynamics on a graph-structured population of size N. Each individual is located at the vertex of a regular graph of degree k. Individuals obtain a payoff by interacting with their k neighbors in a single d-person symmetric game (i.e., d = k+1). If j co-players play A, a focal A-player obtains a whereas a focal B-player obtains b, as indicated in Table 1.

Table 1

Payoffs to A-players and B-players.

Opposing A-players	0	1	…	j	…	d − 1
payoff to A	a0	a1	…	aj	…	ad−1
payoff to B	b0	b1	…	bj	…	bd−1

We model the stochastic evolutionary dynamics as a Markov process on a finite space state. More specifically, we consider a Moran death-Birth process [12, 14, 53] according to which, at each time step: (i) a random individual is chosen to die, and (ii) its neighbors compete to place a copy of themselves in the new empty site with probability proportional to 1 − w + w × payoff, where the parameter w measures the intensity of selection. Without mutation, such a Markov process has two absorbing states: one where all vertices are occupied by A-players and one where all vertices are occupied by B-players. Let us denote by ρ the fixation probability of a single A-player in a population of B-players, and by ρ the fixation probability of a single B-player in a population of A-players. We take the comparison of fixation probabilities, i.e. as a measure of evolutionary success [54] and say that A is favored over B if condition (1) holds.

Under weak selection (i.e., w ≪ 1) the condition for A to be favored over B holds if the sigma rule for multiplayer games [50] is satisfied, i.e., if where σ0, …, σ are the d structure coefficients (constants that depend on the population structure and on the update dynamics), and are differences between payoffs, which we will refer to in the following as the “gains from flipping”. The gains from flipping capture the change in payoff experienced by a focal individual playing B in a group where j co-players play A when all players simultaneously switch strategies (so that A-players become B-players and B-players become A-players). It turns out that the payoffs of the game only enter into condition (1) via the gains from flipping Eq (3), as the structure coefficients are themselves independent of a and b.

Structure coefficients are uniquely determined up to a constant factor. Setting one of these coefficients to one thus gives a single structure coefficient for d = 2 [22]. For d > 2, and in the usual case where structure coefficients are nonnegative, we can impose without affecting the selection condition (2). For our purposes, this normalization turns out to be more useful than setting one coefficient to one, as it allows us to rewrite the sigma rule Eq (2) as where f(j) ≡ f, and J is the random variable with probability distribution prescribed by the “normalized structure coefficients” . In light of condition (4), the sigma rule can be interpreted as stating that strategy A is favored over B if the expected gains from flipping are greater than zero when the number of co-players J is distributed according to the normalized structure coefficients. From this perspective, different models of population structure lead to different normalized structured coefficients and hence to different expected gains from flipping, which in turn imply different conditions for strategy A to be favored over B in a given multiplayer game [51]. For instance, a well-mixed population with random group formation updated with either a Moran or a Wright-Fisher process leads to normalized structure coefficients given by [39, 40]:

A normalized sigma rule such as the one given by Eq (4) holds for many spatial models and associated updating protocols [50, 51]. Here, we focus on the case of regular graphs updated with a Moran death-Birth process. We provide exact expressions for the case of cycles for which k = 2. For k ≥ 3, we bypass the difficulties of an exact calculation by using a combination of pair approximation [55, 56] and diffusion approximation [14]. Our approach implicitly assumes that graphs are equivalent to Bethe lattices (or Cayley trees) with a very large number of vertices (N ≫ k). In addition, weak selection intensities (wk ≪ 1) are also required for an implicit argument of separation of timescales to hold. In order to assess the validity of our analytical approximations, we implemented a computational model of a Moran death-Birth process in three different types of regular graphs (rings, random graphs, and lattices) with different degrees and estimated numerically the fixation probabilities ρ and ρ as the proportion of realizations where the mutant succeeded in invading the wild-type.

Results

Exact structure coefficients and sigma rule for cycles

Going beyond the complete graph representing a well-mixed population, the simplest case of a regular graph is the cycle, for which k = 2 (and consequently d = 3). In this case, we find the following exact expressions for the structure coefficients (S1 Text, Section 1):

For large N, the structure coefficients reduce to , and the sigma rule Eq (4) simplifies to

This is also the condition for the boundary between a cluster of A-players and a cluster of B-players to move in favor of A-players for weak selection [57] (Fig 1). Condition (7) implies that A can be favored over B even if A is strictly dominated by B (i.e., a < b for all j) as long as the payoff for mutual cooperation a2 is large enough so that a2 > b0+(b1 − a1); a necessary condition for this inequality to hold is that A strictly Pareto dominates B (i.e., a2 > b0). Such a result is impossible in well-mixed populations, where the structure coefficients Eq (5) prevent strictly dominated strategies from being favored by selection. Condition (7) provides a simple example of how spatial structure can affect evolutionary game dynamics and ultimately favor the evolution of cooperation and altruism.

Fig 1

Payoffs at the boundary of two clusters in the cycle.

Under weak selection, the cluster of A-players expands if the sigma rule a1+a2 > b1+b0 holds. As a player is never paired with two players of the opposite strategy, neither a0 nor b2 enter into this expression. This provides an intuition behind our analytical results in the simple case when the graph is a cycle.

Approximate structure coefficients and sigma rule for regular graphs of degree k ≥ 3

For regular graphs of degree k ≥ 3, we find that the structure coefficients can be approximated by (S1 Text, Section 2) where and

These expressions are nontrivial functions of the degree of the graph k and thus difficult to interpret. For instance, for k = 3, we obtain .

Promotion of multiplayer cooperation

The previous results hold for any symmetric multiplayer game with two strategies. To investigate the evolution of multiplayer cooperation, let us label strategy A as “cooperate”, strategy B as “defect”, and assume that, irrespective of the focal player’s strategy, the payoff of a focal player increases with the number of co-players playing A, i.e.,

This restriction on the payoffs is characteristic of “cooperation games” [51] in which playing A is beneficial to the group but might be costly to the individual. Well-known multiplayer games belonging to this large class of games include different instances of volunteer’s dilemmas [44, 46], snowdrift games [58], stag hunts [59], and many other instances of public, club, and charity goods games [43].

We are interested in establishing whether graph-structured populations systematically lead to structure coefficients that make it easier to satisfy the normalized sigma rule Eq (4) than well-mixed populations (the baseline case scenario of a population with no spatial structure) for any cooperation game satisfying condition (11). In other words, we ask whether a graph is a stronger promoter of cooperation than a well-mixed population. Technically, this is equivalent to asking whether the set of games for which cooperation is favored under a graph contains the set of games for which cooperation is favored under a well-mixed population, i.e., whether a graph is greater than a well-mixed population in the “containment order” [51]. A simple sufficient condition for this is that the difference in normalized structure coefficients, ςG − ςW, has exactly one sign change from − to + [51]. This can be verified for any N > 3 in the case of cycles (k = 2) by inspection of eqs (5) and (6). For large regular graphs of degree k ≥ 3 and hence multiplayer games with d ≥ 4 players, we checked the condition numerically by comparing eqs (5) and (8) for k = 3, …, 100. We find that ςG − ςW always has a single sign change from − to + and hence that, in the limit of validity of our approximations, regular graphs promote more cooperation than well-mixed populations for all games fulfilling Eq (11) (Fig 2). In the following, we explore in more detail the sigma rule for particular examples of multiplayer games.

Fig 2

Structure coefficients and containment order of cooperation.

Approximated (normalized) structure coefficients ς for large regular graphs of degree k = 5 updated with a Moran death-Birth process () and large well-mixed populations where groups of d = 6 players are randomly matched to play a game (). Since ςG − ςW has one sign crossing from − to +, the graph is greater in the containment order than the well-mixed population (denoted by ). Consequently, if the sigma rule holds for a well-mixed population with coefficients ςW, then it also holds for a graph-structured population with coefficients ςG, for any cooperation game.

Examples

Collections of two-player games

As a consistency check, let us consider the case where individuals play two-player games with their k neighbors and collect the payoffs of the different interactions. The two-player game is given by the payoff matrix

The payoffs for the resulting multiplayer game, which are just the sum of payoffs of the pairwise games, are then given by a = jα + (k − j)β and b = jγ + (k − j)δ. The sigma rule Eq (4) can hence be written as

We can show that (S1 Text, Section 2.9) so that condition (13) is equivalent to i.e., the sigma rule previously established for pairwise games in regular graphs [cf. Eq (24) in the Supplementary Material of Ref. [14]]. For a pairwise donation game (for which , , , δ = 0, where and are respectively the benefit and cost of donation) this reduces to the well-known rule [7, 14, 16].

Linear games

Suppose now that a and b are both linear functions of j. We can thus write , for some parameters , , and . When , such a game can be interpreted in terms of a social dilemma as follows. Cooperators each pay a cost in order to provide a benefit to each of their co-players; defectors receive the benefits but pay no cost. In addition to the benefit , cooperators also get an additional bonus per other cooperator in the group. This bonus can be positive or negative.

For such linear games, and by making use of Eq (14), the sigma condition simplifies to . When there is no bonus () the game is an additive prisoner’s dilemma [60] and we recover the condition . In the limit of large k, the sigma condition becomes .

Volunteer’s dilemma

As an example of a nonlinear multiplayer game satisfying condition (11), consider the volunteer’s dilemma [44, 45]. In such a game, one cooperator can produce a public good of value at a personal cost ; defectors pay no cost and provide no benefit. Payoffs are then given by for all j, b0 = 0, and for j > 0. The sigma rule, Eq (4), for the volunteer’s dilemma hence reduces to

For the cycle, we thus find which in the limit of large N reduces to . For large regular graphs of degree k ≥ 3, our approximations lead to

These conditions contrast with that for a large well mixed population, which is given by .

Suppose now that the cost of producing the public good is shared among cooperators [46]. Payoffs are then given by , b0 = 0 and for j > 0. In this case the sigma rule simplifies to

This leads to in the case of a finite cycle of size N and for a large cycle. Contrastingly, in a well-mixed population,

Computer simulations

To assess the validity of our approximations, we compare our analytical results with explicit simulations of evolutionary dynamics on graphs (Fig 3, N = 100; S1 Fig, N = 500). We implemented three different kinds of regular graphs: (i) random regular graphs, (ii) rings (generalized cycles in which each node is connected to k/2 nodes to the left and k/2 nodes to the right), and (iii) lattices (a square lattice with von Neumann neighborhood with k = 4, a hexagonal lattice with k = 6, and a square lattice with Moore neighborhood and k = 8). Analytical predictions are in good agreement with simulation results in the case of cycles (i.e., rings with k = 2, for which our expressions are exact) and for all random regular graphs that we explored. Contrastingly, for rings with k ≥ 4 and lattices, our approximations tend to underestimate the critical benefit-to-cost ratio beyond which the fixation probability of cooperators is greater than that of defectors. In other words, our analytical results seem to provide necessary but not sufficient conditions for cooperation to be favored. Such discrepancies stem from the fact that our analysis assumes graphs with no loops such as Cayley trees; the error induced by our approximations is more evident when looking at the actual fixation probabilities (S2 Fig, N = 100, S3 Fig, N = 500) and not just at their difference. As all graphs with k > 2 we considered do contain loops, such mismatch is expected—in particular for rings and lattices, which are characterized by high clustering. Perhaps more importantly, our simulations suggest that the critical benefit-to-cost ratio for the volunteer’s dilemma without cost sharing in rings and lattices with k ≥ 6 is greater than the corresponding values for random graphs and well-mixed populations. This illustrates a case in which a graph-structured population updated with a death-Birth process leads to less favorable conditions for the evolution of cooperation than a well-mixed population.

Fig 3

Simulations of evolutionary game dynamics on graphs (difference in fixation probabilities, N = 100).

The first row shows the type of (regular) graph for the particular case of k = 4, i.e., each node has exactly four neighbors. The second and third rows show simulation results for the volunteer’s dilemma without cost-sharing and with cost-sharing, respectively. Simulation data in the first column correspond to random regular graphs, in the second column to rings, and in the third column to lattices. The fixation probability of cooperators, ρ (defectors, ρ) was calculated as the fraction of runs where a single cooperator (defector) reached fixation out of 107 runs. Symbols show the difference between such fixation probabilities, as a function of the benefit-to-cost ratio , for different types and degrees of the graph. Lines indicate analytical predictions for the difference in fixation probabilities (left hand side of Eq (4) with normalized sigmas given by Eqs (6) or (8)). Dashed vertical lines the critical benefit-to-cost ratios above which we have ρ > ρ for well-mixed populations (right hand side of Eqs (16) or (19) with normalized sigmas given by Eq (5)). Parameters: population size N = 100, intensity of selection w = 0.01, payoff cost .

Discussion

We studied evolutionary multiplayer game dynamics on graphs, focusing on the case of a Moran death-Birth process on regular structures. First, we used a combination of pair approximation and diffusion approximation to provide analytical formulas for the structure coefficients of a regular graph, which together with the payoffs from the game determine when a strategy is more abundant than another in the limits of weak selection and weak mutation. Such a condition is valid for any symmetric multiplayer game, including the volunteer’s dilemma [44-46] and other multiplayer social dilemmas discussed in the recent literature [38, 41, 58, 59, 61]. The condition can be used to determine the specific conditions (in terms of the degree of the graph and the parameters of the game, such as payoff costs and benefits) under which cooperation will thrive. The structure coefficients also provide a way of comparing the graph with other population structures, such as the well-mixed population. In particular, and to the extent that our approximations are valid, graphs updated with a death-Birth process are more conducive to the evolution of cooperation than well-mixed populations for a large class of games (see condition (11)).

Second, we used numerical simulations to estimate the fixation probabilities and the difference in fixation probabilities of different strategies for particular examples of games (volunteer’s dilemma with and without cost sharing) and graphs (random regular graphs, rings, and lattices). Although simulations agree very well with the analytical approximations in the case of random regular graphs, discrepancies are evident in the case of rings and lattices, which are characterized by higher clustering and for which pair approximation is not sufficiently accurate. In these cases, the analytical approximations systematically overestimate the ability of a graph to promote the evolution of cooperation. Importantly, in the case of the volunteer’s dilemma without cost sharing and for rings or lattices of relatively large degree, the critical benefit-to-cost ratio above which cooperation is favored is greater, not smaller, than the corresponding value for a well-mixed population. Even though detrimental effects of spatial structure on cooperation have been previously noted in similar studies [62], our results are counterintuitive given the updating protocol and the intensity of selection we explored. Indeed, a death-Birth Moran process under weak selection would always favor cooperation (with respect to a well-mixed population of the same size) for any linear cooperation game, including any collection of two-player cooperation games. Our simulations show that this might not be the case when social dilemmas are instead modelled as nonlinear games such as the volunteer’s dilemma.

We used pair approximation and diffusion approximation to find approximate values for the structure coefficients, but other approaches can be used to obtain better estimates of them. In particular, coalescent theory [63] allows us to write the sigma rule in terms of selection coefficients (dependent on the payoffs of the game and the demographic parameters of the model) and expected coalescence times under neutrality [64, 65]; however, such expected coalescence times can be difficult to obtain exactly. Alternatively, for small graphs, the sigma rule and hence the structure coefficients can be explicitly calculated from the transition matrix of the evolutionary process (cf. Appendix C of Ref. [26]). Finally, we note that even in cases for which the structure coefficients are difficult to obtain by purely analytical means, they can be estimated numerically, either indirectly (by estimating the expected times to coalescence) or directly (by computing and comparing fixation probabilities).

For simplicity, we assumed that a focal player obtains its payoff from a single multiplayer game with its k immediate neighbors. Such assumption allowed us to consider multiplayer interactions on graphs in a straightforward way. However, this is in contrast with a common assumption of many studies of multiplayer spatial and network games in which a focal player’s total payoff is the sum of payoffs obtained in k+1 different games, one “centered” on the focal player itself and the other k centered on its neighbors [47-49]. As a result, focal players interact not only with first-order but also with second-order neighbors, which would lead to more intricate structure coefficients. For example, in this case the structure coefficients of a cycle are given by [51, 66]

These values are different from those we calculated under the assumption that individuals play a single game with first-order neighbors, given by Eq (6). For N > 4, the structure coefficients fulfill , meaning that our assumption of payoffs from a single game leads to less restrictive conditions for cooperation to be favored by selection. This observation is in line with previous results for pairwise games on graphs suggesting that the condition for the evolution of cooperation is optimized when interaction and replacement neighborhoods coincide [67], which corresponds to our assumption of individuals playing a single game. Future work should consider the calculation of structure coefficients for the cases where the payoff to a player also depends on games centered on neighbors and how the condition for the promotion of cooperation differs from the one resulting from our simplifying assumption.

We modelled social interactions as multiplayer matrix games with two discrete strategies (A and B) and obtained our results by assuming that selection is weak (w is small). Alternatively, one could model the same multiplayer game but assume instead that players can choose between two similar mixed strategies z and z + δ, where z and z + δ refer to the probability of playing A for each strategy, and δ is small [43, 68, 69]. In such a “δ-weak selection” scenario, and for any number of players, only a single structure coefficient is needed to identify conditions under which a higher probability of playing A is favored by natural selection. For transitive graphs of size N and degree k, this structure coefficient is given by [7, 25]

Exchanging the structure coefficient σ for the “scaled relatedness coefficient” κ of inclusive fitness theory via the identity κ = (σ − 1)/(σ+1) [65], we obtain [16]

With such a value, recent results on multiplayer discrete games in structured populations under δ-weak selection [43] can be readily applied to show that, for all cooperation games as we defined them and for a death-Birth protocol, A is favored over B more easily for a graph-structured population than for a well-mixed population, as long as N > k+1. Such prediction qualitatively coincides with the one obtained from our analytical approximations, but does not capture our numerical results for the volunteer’s dilemma in rings and lattices.

To sum up, we have shown that even for multiplayer games on graphs, which are routinely analyzed by simulation only, some analytical insight can be generated. However, fully accounting for the complexity of evolutionary multiplayer games in graphs with high clustering remains a challenging open problem.

Supplementary Methods.

Calculations for the exact structure coefficients of cycles, the approximate structure coefficients of graphs of degree k ≥ 3, and a short description of the computational model used for our simulations.

(PDF)

Click here for additional data file.

Simulations of evolutionary game dynamics on graphs (difference in fixation probabilities, N = 500).

Same as in Fig 3, but for a population size N = 500.

(PDF)

Click here for additional data file.

Simulations of evolutionary game dynamics on graphs (fixation probabilities, N = 100).

Open symbols show the fixation probability of cooperators (ρ) and filled symbols the fixation probability of defectors (ρ) as a function of the benefit-to-cost ratio , for different types and degrees of the graph. Lines indicate analytical predictions for the fixation probabilities. Parameters: population size N = 100, intensity of selection w = 0.01, payoff cost .

(PDF)

Click here for additional data file.

Simulations of evolutionary game dynamics on graphs (fixation probabilities, N = 500).

Same as in S2 Fig, but for a population size N = 500.

(PDF)

Click here for additional data file.