The effect of increasing returns to scale in public goods investment on threshold values of cooperation under social exclusion mechanism

Physica A 532 (2019) 121866

Contents lists available at ScienceDirect

Physica A journal homepage: www.elsevier.com/locate/physa

The effect of increasing returns to scale in public goods investment on threshold values of cooperation under social exclusion mechanism ∗

∗

Ji Quan a , , Junjun Zheng b , , Xianjia Wang b , Xiukang Yang a a b

School of Management, Wuhan University of Technology, Wuhan 430070, China School of Economics and Management, Wuhan University, Wuhan 430072, China

highlights • • • •

Increasing returns to scale in PGG under social exclusion mechanism is studied. Possible stochastic stable equilibriums (SSEs) of the model are investigated. The range of parameters for the system to choose different SSEs is studied. The cooperation behavior of the population is revealed in the stochastic system.

article

info

Article history: Received 21 January 2019 Received in revised form 18 April 2019 Available online 24 June 2019 Keywords: Increasing returns to scale Social exclusion Stochastic stable equilibrium The evolution of cooperation

a b s t r a c t This paper establishes a stochastic evolutionary public goods game (PGG) model with exclusion-type strategies in a well-mixed and finite size population. The effect of increasing-returns-to-scale in public goods investment is considered. We focus on the stochastic stable equilibrium (SSE) of the evolutionary system. By numerical experiments, we observe that there are two types of SSE in the system, namely, ‘‘All individuals choosing defection’’ (‘‘All D’’ state) and ‘‘the coexistence of cooperation and exclusion’’ (‘‘C+E’’ states); and three types of phases, namely, choosing the ‘‘All D’’ state with probability one (‘‘D’’ phase), choosing the ‘‘C+E’’ states with probability one (‘‘C+E’’ phase) and choosing the ‘‘All D’’ and ‘‘C+E’’ states with non-zero probabilities (‘‘C+E+D’’ phase). We study the combined effects of four parameters (probability of exclusion success, unit exclusion cost, increasing-returns-to-scale coefficient, investment amplification factor) on the phase selection of the system. We get the boundary curves or surfaces for each pair or set of combined parameters, thus conditions for the evolution of cooperative strategies depending on these parameters can be given. Notably, we observe that in a well-mixed population with exclusion-type strategies, the increasing-returns-to-scale coefficient cannot change the parameter boundary between the ‘‘C+E’’ phase and the ‘‘C+E+D’’ phase within our parameters. Corresponding simulation experiments are also carried out in a structured population, and it is found that the increasing-returns-toscale effect can greatly reduce critical values of the amplification factor under which cooperation can emerge, and the increase in the increasing-returns-to-scale coefficient can induce the system to transfer directly from the ‘‘D’’ phase to the ‘‘C+E’’ phase, which is different from the well-mixed situation. This research can help us better understand the emergence of cooperation in the PGG with increasing-returns-to-scale effect under social exclusion mechanism. © 2019 Elsevier B.V. All rights reserved.

∗ Corresponding authors. E-mail addresses: [email protected] (J. Quan), [email protected] (J. Zheng). https://doi.org/10.1016/j.physa.2019.121866 0378-4371/© 2019 Elsevier B.V. All rights reserved.

2

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

1. Introduction Cooperation phenomena have prevailed since the birth of human society or other earlier biological populations. However, how to explain the emergence and evolution of cooperation remains a puzzle to human beings [1,2]. This topic has been the focus of both social and natural sciences [3]. In game theory, we usually use the prisoner’s dilemma game (PDG) or the public goods game (PGG), respectively, to describe the dilemma of two-person or multi-person cooperation, which provides simple and profound models for the study of cooperation [4–6]. Evolutionary game theory, which maps the individual’s payoff in the games into fitness and is based on the idea of survival of the fittest, provides a dynamic analysis and theoretically feasible method for studying the evolution of cooperation [7]. We believe that the study of the evolution of cooperation in the framework of evolutionary games requires two types of mechanisms. One is individual’s decision-making mechanism in an interactive environment, and the other is the interaction mechanism between individuals. Evolutionary game theory assumes that individuals are bounded rational, and individuals adjust their strategies through feedback and learning. In fact, in evolutionary games, individuals make decisions according to certain rules and the individual’s decision-making mechanism is to define this rule. Making decisions according to rules is fundamentally different from traditional decision-making theory based on revenue maximization. At present, many kinds of decision rules based on individual learning have been proposed [8–24]. Individual’s decisionmaking rules in the games will affect the cooperation behavior of the population [25]. The population can show a higher level of cooperation under some decision-making mechanisms [18]. The interaction mechanism between individuals defines the interaction rules and payoff structure of individuals. The five main cooperation mechanisms that Nowak [26] summarized can be categorized into different interaction rules. As an example, the network reciprocity mechanism changes rules of interaction by introducing a network structure between individuals [27–32]. Under the network reciprocity mechanism, the macro-cooperative phenomenon caused by individual micro-behavior is a complex system problem, which can be analyzed by statistical physics methods [33–35]. In addition to theoretical models, some recent studies have further explored cooperative mechanisms through behavioral experiments [36–38]. We focus on the direct reciprocity mechanism in this paper. Specifically, we introduce exclusion-type strategies in the traditional PGG model. Exclusion can be seen as a special punishment imposed by excluders on defectors. The same as traditional punishing strategies, excluders who implement this punishing action also have to pay an extra cost. The difference is that the defectors who are punished will not pay a fixed fine, but be expelled from the group with a certain probability, and individuals being driven out cannot share the cooperative benefit of the group. One specific example is that for drivers who violate traffic rules can be deprived of their right to use public roads by revoking their driver’s license. Like other punishment [23,39–52] or reward [37,53–57] mechanisms, social exclusion is also a mechanism that has been proven to be effective in promoting cooperation in human public goods experiments [58–60]. On one hand, there is evidence from experimental studies that exclusion is widespread in some biological populations [61,62] and human social systems [63]. On the other hand, exclusion-type strategy itself is a second-order social dilemma. Who will implement this strategy, or how this strategy can be evolved in a self-organized system has always been the focus of academic attention [64–69]. In addition, in the standard PGG model, the total investment of all cooperative individuals is evenly allocated after amplification, and the amplification factor is independent of the number of cooperative individuals. However, in many real-life situations, the value of public goods increases non-linearly with the number of investors. Based on this observation, some typical nonlinear PGG models have been proposed [70,71]. We consider another type of nonlinear relationship between the return and the sum of contributions in this paper. We observe that in some types of multiperson collaborative activities, the greater the number of people working together, the faster the value of the products obtained from cooperation will increase. A typical example is cooperative hunting [72]. Suppose the hunters form a circle for hunting, and the length of the region that each hunter can control is l. When the number of cooperative hunters is n, the perimeter of the region that can be controlled is nl, then the radius is nl/2π , and the area is n2 l2 /4π . Obviously, the area of the controllable region is proportional to the square of the number of hunters. Assuming that the number of prey is proportional to the area of the region, then the return of cooperation will be proportional to the square of (and not linear to) the number of cooperators. In economics, the academic terminology ‘‘increasing returns to scale’’ is usually used to describe this situation [72,73]. Thus, if increasing-returns-to-scale effect in the investment is taken into account under social exclusion mechanism, further research is needed on how the population’s cooperative behavior and the system’s stable equilibrium will change. This is the main motivation for this study. Different from the Moran-process-based rule which uses the fixation probability to analyze the conditions of mutual intrusion between different strategies, or other methods based on embed-Markovchain which only focus on the homogeneous states in which all individuals choose the same strategy [65–68,74–76], we describe the state of the system as a continuous time quasi-birth-and-death process. Moreover, we calculated the limit distributions of all possible states, through which, we can analyze the evolutionary stable state of the system. This analysis method can obtain the stochastic stable equilibrium of the system [51,77].

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

3

2. Model

2.1. PGG with increasing returns to scale under social exclusion mechanism

We assume a finite population of size M, and the three strategy types, namely, cooperation, defection and exclusion are well mixed. Let X , Y and W (X + Y + W = M) denote the number of each type, respectively. The exclusion-type strategy invests in the public goods game, and also excludes defectors in his group. Exclusion behavior is costly, which will bring additional cost to the excluders, but it can prevent free-riders from sharing their investment income. Defectors who are excluded from the group cannot share the benefit of the public goods. In the model, we assume that an excluder cannot expel a defector successfully with certainty, but only with probability β (0 < β ≤ 1). Let cE denote unit cost for each excluder excluding each defector. Moreover, we consider a positive correlation between the number of investors and the amplification factor of the game in our model. Let α (α ≥ 1) denote the coefficient of increasing-returns-to-scale of the public goods. In the following, we analyze the expected payoff of each type strategy in the PGG with increasing-returns-to-scale under social exclusion mechanism in this finite size population. Each time, N individuals are sampled randomly from the population to participate in the PGG. Let l and j denote the number of excluders and defectors, respectively, in the sampled group. Let s denote the number of defectors who escape from all excluders in the game. For each defector, whether it is excluded from the group obeys a 0–1 distribution(with ) parameter p1 . Thus, when j l ̸ = 0, the probability that s defectors escape from all excluders in the group is p1 s (1 − p1 )j−s , where p1 = s (1 − β )l is the probability of a defector escaping from all excluders. Without loss of generality, suppose the unit r(N −j)α investment cost of cooperators equals one. In this situation, the payoff of cooperators is N −j+s − 1. Calculating expected payoff of cooperators with respect to random variable s, we obtain

∑j

s=0

()

j r(N −j)α p1 s (1 − p1 )j−s [ N −j+s − 1]. s

For an(existing probability ( ) ) ( cooperator, the )/ ( ) that there are l(l > 0) excluders and j defectors in the sampled group is W Y M −1−W −Y M −1 . Calculating expected payoff of cooperators again with respect to random l j N −1−l−j N −1 ( )( )( ) W Y M −1−W −Y () ∑N −1 ∑N −1−l l ∑j j N −1−l−j j r(N −j)α ( ) variables l(l > 0) and j, we obtain × s=0 p1 s (1 − p1 )j−s [ N −j+s − 1]. l=1 j=0 s M −1 N(− 1 )/( ) M −1−W M −1 The probability of no excluders in the group is p(l = 0) = . In this situation, the payoff of N −1 N −1 r(N −j)α

the probability that there are j defectors in the sampled group cooperators ( ) ( is N − 1.) For / (an existing cooperator, ) Y M −1−W −Y M −1−W is , and calculating expected payoff of cooperators with respect to random j N −1−j N −1 ( ) ( )( ) M −1−W Y M −1−W −Y ∑N − 1 j α N −1 N −1−j ( ) ( ) variable j and l = 0, we obtain [ r(NN−j) − 1]. j=0 M −1 M −1−W N −1 N −1 Thus considering the two situations of l ̸ = 0 and l = 0, the total expected payoff of the cooperation-type strategies is

( )( )( ,Y ,W ) π (X = C

N −1 N −1−l ∑ ∑

W l

Y j

( l=1

( +

j=0

M −1−W N −1

(

)

M −1 N −1

N −1 ∑

Y j

( j=0

)

)

×

j ( ) ∑ j

s

p1 s (1 − p1 )j−s [

s=0

( )(

)

M −1 N −1

M −1−W −Y N −1−l−j

M −1−W −Y N −1−j

M −1−W N −1

)

r(N − j)α N −j+s

− 1] .

) [

r(N − j)α N

− 1]

(1)

4

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

Similarly, we can get expected payoffs of the defection-type and exclusion-type strategies as follows.

( )( (X ,Y ,W )

πD

=

W l

N −1 N −1−l ∑ ∑

( l=1

( +

j=0

M −1−W N −1

(

M −1 N −1

( ,Y ,W ) π (X = E

)(

Y −1 j

N −1 N −1−l ∑ ∑

)

j=0

× p1 ×

M −1 N −1

N −1 ∑

)(

M −W −Y N −1−j

(

M −1−W N −1

j=0

)( )( Y j

j ( ) ∑ j

s

p1 s (1 − p1 )j−s [

r(N − j − 1)α

s=0

Y −1 j

)

W −1 l

)

)

(

( l=0

M −W −Y N −1−l−j

M −W −Y N −1−l−j

r(N − j − 1)α N

] .

(2)

− 1 − cE j].

(3)

) [

)

N −j+s

]

)

)

M −1 N −1

×

j ( ) ∑ j

s s=0

p2 s (1 − p2 )j−s [

r(N − j)α N −j+s

2.2. Stochastic evolutionary dynamics and stochastic stable equilibrium The number of each strategy-type will evolve with the amount of its payoff. In order to describe the evolutionary process, we introduce a multi-dimensional stochastic process z(t) = (z1 (t), z2 (t), M − z1 (t) − z2 (t)) to denote the system state, where the three components represent the number of cooperators, defectors and excluders, respectively, in the population at time t. For convenience, we abbreviate it as (z1 (t), z2 (t)). The state space of the system is S = {(i, j)|0 ≤ (M +1)(M +2) i + j ≤ M ; i, j ∈ N}, and there are |S | = elements in the state space. Each time step, one individual in the 2 population adjusts its strategy. The adjustment of the individual’s strategy leads to the change of the system state. The three assumptions about bounded rationality of individuals in the population, namely, inertia, myopic and mutation in literature [11] are adopted in our model. Owing to inertia, we can assume that more than two individuals adjusting their strategies at the same time is impossible. Myopic means that when an individual chooses a strategy, it only considers the current payoff, regardless of the payoff in the future. Mutation refers to the possibility that individuals may choose a non-optimal strategy with a small probability because of the complex decision-making environment and the limited nature of individual cognitive ability. Based on the above assumptions, when the system state is (z1 (t), z2 (t)) = (i, j) ∈ S, the transfer rate of strategy x towards strategy y can be described as ,j) (i,j) p(ix→ − πx(i,j) )+ , x, y ∈ {C , D, E }, x ̸= y, y = ε + κ · ( πy

(4) (i,j) x

(i,j) y ,

where f = max(f , 0), ε > 0 is a small positive number, κ > 0. As an example, when π >π the individual in (i, j) state has a more incentive to move from strategy y to strategy x, but because of mutation, the transfer rate of strategy x (i,j) to strategy y is px→y = ε . Thus, parameter ε can be seen as the noise intensity in the environment, and κ can be seen as the speed at which the individual responds to the environment and we fix κ = 1 in our discussion. (i,j,0) ,j,k) (i,0,k) does not make sense. In this When i, j and k takes zero, respectively, the corresponding π (0 πD and πE C situation, the payoffs of cooperation, defection and punishment type strategies are defined as the average payoff of the population. Let I = (i, j), I ′ = (i′ , j′ ). Due to homogeneity, let pI ,I ′ (t) denote the probability of the system that transfers from state I to state I ′ after time t. That is, +

pI ,I ′ (t) = p{z(s + t) = (i′ , j′ )|z(s) = (i, j)}, ∀s > 0.

(5)

According to the system evolutionary rules and the transition rate between different strategies, after an adequately (i,j) small time t, the probabilities of the system transfers from state I = (i, j) to states (i − 1, j) and (i − 1, j + 1) are pC →E t + o(t) (i,j) (i,j) and pC →D t + o(t), respectively; transfers from state I = (i, j) to states (i, j − 1) and (i + 1, j − 1) are pD→E t + o(t) and (i,j) (i,j) (i,j) pD→C t + o(t), respectively; transfers from state I = (i, j) to states (i + 1, j) and (i, j + 1) are pE →C t + o(t) and pE →D t + o(t), ∑ (i,j) respectively; the probability to keep the same state is 1 − x,y∈{C ,D,E } px→y t − o(t), where o(t) is a high order infinitesimal x̸ =y

of t when t is adequately small. As ε > 0, this process is ergodic. According to the properties of the stochastic process, when t → +∞, the limit of pI ,I ′ (t) exists and it is independent of the initial state I. Let lim pI ,I ′ (t) = vIε′

t →+∞

(6)

Then vIε′ is the limit distribution of the two-dimensional stochastic process reaching arbitrary state I ′ (I ′ ∈ S) when the system noise is ε . According to the limit distribution, it is possible to determine the stable states of the system under

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

5

Fig. 1. Limit distributions of the stable states with different α for fixed r = 2, cE = 0.1, β = 0.7. The system chooses the ‘‘All D’’ state with probability 1 when α = 1, which corresponds to the ‘‘D’’ phase; whereas chooses the ‘‘All D’’ state with probability 0.035 and chooses ‘‘C+E’’ states with probability 0.965 when α = 1.3, which corresponds to the ‘‘C+E+D’’ phase; and choose the ‘‘All D’’ state with probability 0.032 and choose ‘‘C+E’’ states with probability 0.968 when α = 1.5, which also corresponds to the ‘‘C+E+D’’ phase. The increase in α from 1 to 1.3 changes the phase of the system, but the continued increase in α from 1.3 to 1.5 does not change the phase of the system, whereas only changes the probability of the system selecting each state.

arbitrary noise intensity. Further, when the noise parameter is gradually reduced to zero, let lim vIε′ = vI ′

ε→0+

(7)

Based on vI ′ , we can determine the limit states of the system and their probability distributions. According to Young’s description in reference [78], state I ′ ∈ S is stochastically stable if and only if vI ′ > 0. 3. Results and discussion We focus on the stochastic stable states of the evolutionary system and the corresponding probability of the system to choose each stable state. We utilize the Gauss–Seidel iterative algorithm introduced in literature [79] to calculate the limit distribution of the stochastic process. By a large number of numerical experiments under arbitrary parameters, we find that only states of (0, M , 0), (i, 0, M − i)(0 ≤ i ≤ M) may be stochastically stable, where state (i, j, k) denote the number of cooperators, defectors and excluders, respectively, in the population. For instance, (0, M , 0) corresponds to the ‘‘All D’’ state and (i, 0, M − i) corresponds to the ‘‘C+E’’ states. That is to say, there are following three different phases in this system. (i) Selecting the ‘‘All D’’ state with probability 1 (denote as ‘‘D’’); (ii) Selecting ‘‘C+E’’ states with probability 1 (denote as ‘‘C+E’’); (iii) Selecting the ‘‘All D’’ and ‘‘C+E’’ states with non-zero probabilities (denote as ‘‘C+E+D’’). In the following, we fix M = 30, N = 5 to investigate the impacts of parameters (α, r , cE , β ) on the probabilities of the system to choose different stable states, and thus to identify the range of parameters that the system reaches to each phase. Due to i + j + k = M for each state (i, j, k), the limit distribution of the stable states can be represented by a three-dimensional simplex graph or called S3. Figs. 1 to 3 utilize S3 to show the corresponding limit distributions under fixed parameter-combinations of (r , cE , β ) with three different values of α . As shown in Fig. 1, for small r, small cE , and relatively large β , e.g. (r , cE , β ) = (2, 0.1, 0.7), the system will choose the ‘‘All D’’ state with probability 1 when α = 1, which corresponds to the ‘‘D’’ phase; whereas choose the ‘‘All D’’ and ‘‘C+E’’ states with probability 0.035 and 0.965, respectively, when α = 1.3, which corresponds to the ‘‘C+E+D’’ phase; and choose the ‘‘All D’’ and ‘‘C+E’’ states with probability 0.032 and 0.968, respectively, when α = 1.5, which also corresponds to the ‘‘C+E+D’’ phase. The increase in α from 1 to 1.3 changes the phase of the system, but the further increase in α from 1.3 to 1.5 does not change the phase of the system, and only changes the probability of the system selecting each state. The probabilities of the system selecting (M , 0, 0) (‘‘All C’’ state) and other states that near the ‘‘All C’’ state are significantly increased when α increased from 1.3 to 1.5. For large r, small cE , and relatively large β , e.g. (r , cE , β ) = (4, 0.1, 0.7) as in Fig. 2, the system will choose ‘‘C+E’’ states with probability 1 for any values of α ≥ 1, which corresponds to the ‘‘C+E’’ phase. Although the system phases are the same, the probability of the system selecting (M , 0, 0) (‘‘All C’’ state) and other states that near the ‘‘All C’’ state increases when α increases. For relatively large r, large cE , and small β , e.g. (r , cE , β ) = (3.5, 1, 0.2) as in Fig. 3, the results are somewhat different. As seen in the figure, the system is in the ‘‘D’’ phase when α = 1 and in the ‘‘C+E+D’’ phase when α = 1.3 or 1.5. The difference is that for the ‘‘C+E+D’’ phase, relatively larger probabilities in the ‘‘C+E’’ states are located on the ‘‘All C’’ side. Intuitively, we can judge that parameters α, r , β will have a positive impact on the cooperation, whereas parameter cE will have a negative impact. However, we need to know exactly the range and extent of the impacts of these parameters on the system to choose different phases. By numerical experiments under arbitrary parameter combinations of (α, r , cE , β ), we got the critical relationship of different parameter combinations for the system to choose different phases. Fig. 4 shows

6

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

Fig. 2. Limit distributions of the stable states with different α for fixed r = 4, cE = 0.1, β = 0.7. The system chooses ‘‘C+E’’ states with probability 1 for any values of α ≥ 1, which corresponds to the ‘‘C+E’’ phase. Although the system is in the same phase in this situation, the probabilities of the system to choose (M , 0, 0) (‘‘All C’’ state) and other states near the ‘‘All C’’ state increase when α increases.

Fig. 3. Limit distributions of the stable states with different α for fixed r = 3.5, cE = 1, β = 0.2. The system reaches the ‘‘D’’ phase when α = 1 and the ‘‘C+E+D’’ phase when α = 1.3 or 1.5. For the ‘‘C+E+D’’ phase, relatively larger probabilities in the ‘‘C+E’’ states are located on the ‘‘All C’’ side.

the results under four fixed values of α . In each sub-figure, there are two boundary surfaces corresponding to (cE , β, r) which divide the 3D space into three parts, and the three parts of the space correspond to the three phases of the system. As seen in the figure, the spatial region of the ‘‘D’’ phase shrinks as α increases, and the spatial region of the ‘‘C+E+D’’ phase expands accordingly, whereas the spatial region of the ‘‘C+E’’ phase is relatively stable as α changes. In order to observe the critical values for each parameter combination more clearly, we select one pair of parameters each time and fix the other parameter to project the three-dimensional image onto a two-dimensional plane for observation. Fig. 5 shows the two-dimensional parameter regions of (cE , r) for the system to choose different phases when β is fixed at 0.3, 0.6 and 0.9, respectively. From left to right in the sub-figures, as β increases, the spatial region of the ‘‘D’’ phase (blue area) shrinks and the spatial region of the ‘‘C+E’’ phase (red area) expands. From top to bottom in the sub-figures, as α increases, the spatial region of the ‘‘D’’ phase shrinks and the spatial region of the ‘‘C+E+D’’ phase (yellow area) expands significantly. Fig. 6 shows the two-dimensional parameter regions of (β, r) for the system to choose different phases when cE is fixed at 0.1, 0.3 and 0.5, respectively. The effect of cE on the boundary curve between different phases can be seen from left to right in the sub-figures, and the effect of α can be seen from top to the bottom. Fig. 7 shows the two-dimensional parameter regions of (cE , β ) for the system to choose different phases when r is fixed at 2.2, 3.2 and 4.2, respectively. From Figs. 5 to 7, we found an interesting conclusion that the value of parameter α only affects the boundary between ‘‘D’’ and the ‘‘C+E+D’’, and it does not affect the boundary between ‘‘C+E+D’’ and ‘‘C+E’’. For any combination of other parameters (β, cE , r), the increase in α can only change the system from the ‘‘D’’ phase to the ‘‘C+E+D’’ phase, but it cannot change the system from the ‘‘C+E+D’’ phase to the ‘‘C+E’’ phase. In other words, it is not possible to eliminate the state in which the system selects ‘‘All D’’ simply by increasing the value of α . To further verify this conclusion, we let α increase continuously from 1 to 1.5. Fig. 8 shows the two-dimensional parameter regions of (α, r) for the system to choose different phases when cE is fixed at 0.1, 0.3 and 0.5, respectively, and β is fixed at 0.3, 0.6 and 0.9, respectively. From the figure we can see that the boundary between the yellow area (C+D+E) and the red area (C+E) is a vertical line, which further verifies the above conclusion.

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

7

Fig. 4. Three-dimensional parameter regions (r , cE , β ) corresponding to the three different phases with four different values of α . The red surface is the boundary between the two phases of ‘‘C+E’’ and ‘‘C+E+D’’; and the blue surface is the boundary between the two phases of ‘‘C+E+D’’ and ‘‘D’’. The spatial region of the ‘‘D’’ phase shrinks as α increases; the spatial region of the ‘‘C+E+D’’ phase expands accordingly, whereas the spatial region of the ‘‘C+E’’ phase relatively stable as α changes.

Fig. 5. Two-dimensional parameter regions (r , cE ) corresponding to the three different phases with different combinations of α and β . The blue, yellow, and red areas indicate the ‘‘D’’, ‘‘C+E+D’’ and ‘‘C+E’’ phase, respectively. From left to right, as β increases from 0.3 to 0.9, the spatial region of the ‘‘D’’ phase shrinks and the spatial region of the ‘‘C+E’’ phase expands. From top to bottom, as α increases from 1 to 1.5, the spatial region of the ‘‘D’’ phase shrinks and the spatial region of the ‘‘C+E+D’’ phase expands significantly.

8

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

Fig. 6. Two-dimensional parameter regions (r , β ) corresponding to the three different phases with different combinations of α and cE . From left to right, as cE increases from 0.1 to 0.5, the spatial region of the ‘‘D’’ phase expands and the spatial region of the ‘‘C+E’’ phase shrinks. From top to bottom, as α increases from 1 to 1.5, the spatial region of the ‘‘D’’ phase shrinks and the spatial region of the ‘‘C+E+D’’ phase expands.

We suspect that this conclusion is related to the introduction of the exclusion-type strategy. To verify this conjecture, we did the same numerical calculations for a simplified model with only C-type and D-type strategies. We found that an increase in the increasing-returns-to-scale coefficient would change the system from the ‘‘D’’ phase to the ‘‘C’’ phase. Fig. 9 shows the parameter regions (α, r) corresponding to two different phases for the simplified model. It can be seen from the above conclusions that in a well-mixed population with exclusion-type strategies, the increasing-returns-to-scale effect will have a significant impact on the equilibrium selection of the system. However, within our parameters, the increasing-returns-to-scale coefficient cannot change the parameter boundary between the ‘‘C+E’’ phase and the ‘‘C+E+D’’ phase. In order to verify the impact of this parameter on the emergence of cooperation in structured populations, inspired by literature [80–82], we only performed simulation experiments under some parameters to compare. We utilize a 100 × 100 square lattice with periodic conditions. Initially, three types of strategies are randomly distributed at each node. At each time step, an individual i is randomly selected, and then one of its neighbors j is randomly selected. Individual j randomly changes its strategy by the probability determined by the Fermi function p(si → sj ) = 1+exp[(π1 −π )/τ ] , where si and πsj denote the strategy and payoff of individual i, respectively, and τ denotes sj

si

the noise intensity. We fix τ = 1 in our simulation. In each iteration, this process is repeated 10,000 times to ensure that each individual has an opportunity on average to adjust its strategy. Through 10 independent simulation experiments, we found that the introduction of the increasing-returns-to-scale parameter α can greatly reduce critical values of parameter r under which cooperation can emerge. Figs. 10 and 11 show that when r = 1.5, cE = 0.1, the parameter combination β = 0.3, α = 1.3 or β = 0.6, α = 1.1 will cause the system to reach the ‘‘C + E’’ phase, that is, all individuals choose a cooperative strategy. Moreover, in our simulation, the increase in α can lead the system to transfer directly from the ‘‘D’’ phase to the ‘‘C+E’’ phase, which is different from the well-mixed situation. Fig. 12 shows the dynamics of the frequencies of the three strategy types over time for the same combination of parameters. All simulation results have been validated on larger scale networks of 200 × 200 and 500 × 500 using the procedure proposed in literature [83]. 4. Conclusions In this paper, we studied a stochastic evolutionary PGG model with exclusion-type strategies in a well-mixed and finite size population. The effect of increasing returns to scale in public goods investment is considered. By numerical experiments, we demonstrated how parameters of probability of exclusion success (β ), unit exclusion cost (cE ), coefficient of increasing returns to scale (α ) and investment amplification factor (r) can have combined impacts on the evolution of

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

9

Fig. 7. Two-dimensional parameter regions (β, cE ) corresponding to the three different phases with different combinations of α and r. From left to right, as r increases from 2.2 to 4.2, the spatial region of the ‘‘D’’ phase shrinks and the spatial region of the ‘‘C+E’’ phase expands. From top to bottom, as α increases from 1 to 1.5, the spatial region of the ‘‘D’’ phase shrinks and the spatial region of the ‘‘C+E+D’’ phase expands. The changes exhibit nonlinear characteristics, and they are only significant in certain places.

Fig. 8. Two-dimensional parameter regions (α, r) corresponding to the three different phases with different combinations of β and cE . Notably, the boundary between the ‘‘C+D+E’’ phase and the ‘‘C+E’’ phase is a vertical line, which further verifies that the increase in α can only change the system from the ‘‘D’’ phase to the ‘‘C+E+D’’ phase, but it cannot change the system from the ‘‘C+E+D’’ phase to the ‘‘C+E’’ phase.

10

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

Fig. 9. Two-dimensional parameter regions (α, r) corresponding to different phases when there are only cooperation and defection two strategy types. The blue and red areas indicate the ‘‘D’’ and ‘‘C’’ phase, respectively. An increase in α can change the system from the ‘‘D’’ phase to the ‘‘C’’ phase in this simplified model.

Fig. 10. Spatiotemporal distribution of the three strategies in the PGG at four different Monte Carlo steps (MCS) for different α , when r = 1.5, cE = 0.1 and β = 0.3. The blue, yellow and red color indicate the D, E and C type strategies, respectively. (a) α = 1; (b) α = 1.1; (c) α = 1.3; (d) α = 1.5. When α = 1 or 1.1, the system reaches the ‘‘D’’ phase rapidly. When α = 1.3, the D-strategy disappears after relatively long-time iterations and the system reaches the ‘‘C+E’’ phase eventually. When α = 1.5, the system reaches the ‘‘C+E ’’ phase rapidly. The increase in α can lead the system to transfer directly from the ‘‘D’’ phase to the ‘‘C+E’’ phase.

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

11

Fig. 11. Spatiotemporal distribution of the three strategies in the PGG at four different Monte Carlo steps (MCS) for different α , when r = 1.5, cE = 0.1 and β = 0.6. (a) α = 1; (b) α = 1.1; (c) α = 1.3; (d) α = 1.5. When α = 1, the system reaches the ‘‘D’’ phase. When α = 1.1, 1.3 or 1.5, the system reaches the ‘‘C+E’’ phase, but has different iteration convergence times. The increase in α also can lead the system to transfer directly from the ‘‘D’’ phase to the ‘‘C+E’’ phase.

cooperation. The effects of these parameters on cooperation are different and non-linear. From the results we can draw the following conclusions. (i) The value of parameter α only affects the boundary between the ‘‘D’’ phase and the ‘‘C+E+D’’ phase, and does not affect the boundary between the ‘‘C+E+D’’ phase and the ‘‘C+E’’ phase. For any combination of other parameters (β, cE , r), the increase in α can only change the system from the ‘‘D’’ phase to the ‘‘C+E+D’’ phase, but cannot change the system from the ‘‘C+E+D’’ phase to the ‘‘C+E’’ phase. In other words, it is not possible to eliminate the state that the system selects ‘‘All D’’ simply by increasing the value of α . (ii) The value of parameter β has a large impact on which phase the system enters. The increase in β can not only change the system from the ‘‘D’’ phase to the ‘‘C+E+D’’ phase, but also change the system from the ‘‘C+E+D’’ phase to the ‘‘C+E’’ phase. For example, for parameter combination cE = 0.1, r = 3, α = 1, the system will enter the ‘‘D’’ phase when β = 0.3, whereas enter the ‘‘C+E+D’’ phase when β = 0.6, and enter the ‘‘C+E’’ phase when β = 0.9. When β is relatively small (β < 0.4), for any combination of other parameters (α, cE , r), the system cannot enter the ‘‘C+E’’ phase and it is not possible to eliminate the ‘‘All D’’ state by adjusting other parameters. (iii) The value of parameter cE has a great influence on which phase the system enters. The decrease in cE can not only change the system from the ‘‘D’’ phase to the ‘‘C+E+D’’ phase, but also change the system from the ‘‘C+E+D’’ phase to the ‘‘C+E’’ phase. For example, for parameter combination β = 0.7, r = 3.5, α = 1.1, the system will enter the ‘‘D’’ phase when cE = 0.5, whereas enter the ‘‘C+E+D’’ phase when cE = 0.3, and enter the ‘‘C+E’’ phase when cE = 0.1. (iv) The value of parameter r also has a great influence on which phase the system enters. The increase in r can not only change the system from the ‘‘D’’ phase to the ‘‘C+E+D’’ phase, but also change the system from the ‘‘C+E+D’’ phase to the ‘‘C+E’’ phase. When α is relatively large (α = 1.3), there is critical value r, when r ≥ r, the ¯ ¯ system will not enter the ‘‘D’’ phase for any combination of parameters (cE , β ). In some cases, there is a correlation between the parameters. For instance, we always have to increase the cost of exclusion in order to increase the probability of exclusion success, in which case, parameters cE and β are positively

12

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

Fig. 12. Frequencies of the three types of strategies over time for different combinations of α and β when r = 1.5, cE = 0.1. From left to right, α increases from 1 to 1.5, and from top to bottom, β increases from 0.3 to 0.6. The corresponding convergence time of the system under different parameter combinations can be seen in this figure.

correlated. The above results are analyzed with independent parameters, and the results are also applicable to parameterrelated situations. Our model can also be extended to structured populations. We carried out corresponding simulation experiments for some parameter combinations on a homogeneous square lattice, and found that the introduction of the increasing-returns-to-scale effect can greatly reduce critical values of the amplification factor under which cooperation can emerge. We also found that the increase in the increasing-returns-to-scale coefficient can lead the system to transfer directly from the ‘‘D’’ phase to the ‘‘C+E’’ phase, which is different from the well-mixed situation. This research can help us better understand the emergence of cooperation in the PGG with increasing-returns-to-scale effect under social exclusion mechanism. Acknowledgments This research was supported by the National Natural Science Foundation of China (No. 71871173, 71771181, 71871171, 71501149) and the Fundamental Research Funds for the Central Universities, China (WUT: 2019VI029). References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19]

D. Kennedy, C. Norman, What don’t we know? Science 309 (2005) 75. E. Pennisi, How did cooperative behavior evolve, Science 309 (2005) 93. V. Capraro, M. Perc, Grand challenges in social physics: In pursuit of moral behavior, Front. Phys. 6 (2018) 107. A.W. Tucker, A two-person dilemma, Read. Games Inf. (1950) 7–8. G. Hardin, The tragedy of the commons, Science 162 (1968) 1243–1248. H. Hamburger, N-person prisoner’s dilemma, J. Math. Sociol. 3 (1973) 27–48. J. Maynard Smith, Evolution and the Theory of Games, Cambridge University, Cambridge, 1982. M. Kandori, G.J. Mailath, R. Rob, Learning, mutation, and long run equilibria in games, Econometrica 61 (1993) 29–56. D. Fudenberg, D. Levine, Learning in games, Eur. Econ. Rev. 42 (1998) 631–639. R. Jimenez, H. Lugo, M.S. Miguel, Gradual learning and the evolution of cooperation in the spatial continuous prisoner’s dilemma, Eur. Phys. J. B 71 (2009) 273–280. M. Amir, S.K. Berninghaus, Another approach to mutation and learning in games, Games Econ. Behav. 14 (1996) 19–43. L. Messerschmidt, A.P. Engelbrecht, Learning to play games using a PSO-based competitive learning approach, IEEE Trans. Evol. Comput. 8 (2004) 280–288. X.J. Chen, F. Fu, L. Wang, Effects of learning activity on cooperation in evolutionary prisoner’s dilemma game, Internat. J. Modern Phys. C 19 (2008) 1377–1387. J. Quan, S. Yang, X.J. Wang, Learning algorithm and the cooperation behavior of continuous prisoner’s dilemma game on complex networks, J. Inf. Comput. Sci. 10 (2013) 3031–3041. A. Amado, W. Huang, P.R.A. Campos, F.F. Ferreira, Learning process in public goods games, Phys. A 430 (2015) 21–31. P. Li, H.B. Duan, Evolution of cooperation driven by incremental learning, Phys. A 419 (2015) 14–22. P.H. Liu, J. Liu, Contribution diversity and incremental learning promote cooperation in public goods games, Phys. A 486 (2017) 827–838. J. Quan, X. Yang, X. Wang, Spatial public goods game with continuous contributions based on particle swarm optimization learning and the evolution of cooperation, Phys. A 505 (2018) 973–983. A.J. Bladon, T. Galla, Learning dynamics in public goods games, Phys. Rev. E 84 (2011).

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

13

[20] M.N. Burton-Chellew, H.H. Nax, S.A. West, Payoff-based learning explains the decline in cooperation in public goods games, Proc. R. Soc. B 282 (2015). [21] S. Lamba, Social learning in cooperative dilemmas, Proc. R. Soc. B 281 (2014). [22] H.H. Nax, M. Perc, Directional learning and the provisioning of public goods, Sci. Rep. 5 (2015) 8010. [23] J. Quan, X. Yang, X. Wang, Continuous spatial public goods game with self and peer punishment based on particle swarm optimization, Phys. Lett. A 382 (2018) 1721–1730. [24] L. Shi, C. Shen, Y. Geng, C. Chu, H. Meng, M. Perc, S. Boccaletti, Z. Wang, Winner-weaken-loser-strengthen rule leads to optimally cooperative interdependent networks, Nonlinear Dynam. (2019) http://dx.doi.org/10.1007/s11071-11019-04772-11076. [25] W. Zhen, M. Jusup, S. Lei, J.H. Lee, Y. Iwasa, S. Boccaletti, Exploiting a cognitive bias promotes cooperation in social dilemma experiments, Nature Commun. 9 (2018) 2954. [26] M.A. Nowak, Five rules for the evolution of cooperation, Science 314 (2006) 1560–1563. [27] C. Xia, X. Li, W. Zhen, M. Perc, Doubly effects of information sharing on interdependent network reciprocity, New J. Phys. 20 (2018) 075005. [28] J.H. Qin, Y.M. Chen, W.M. Fu, Y. Kang, M. Perc, Neighborhood diversity promotes cooperation in social dilemmas, IEEE Access 6 (2018) 5003–5009. [29] Y. Li, J. Zhang, M. Perc, Effects of compassion on the evolution of cooperation in spatial social dilemmas, Appl. Math. Comput. 320 (2018) 437–443. [30] Z. Wang, L. Wang, A. Szolnoki, M. Perc, Evolutionary games on multilayer networks: a colloquium, Eur. Phys. J. B 88 (2015) 124. [31] C. Shen, C. Chu, L. Shi, M. Jusup, M. Perc, Z. Wang, Coevolutionary resolution of the public goods dilemma in interdependent structured populations, Epl 124 (2018) 48003. [32] J. Quan, Y. Zhou, M. Zhang, C. Tang, X. Wang, The impact of heterogeneous scale return coefficient between groups on the emergence of cooperation in spatial public goods game, J. Stat. Mech. Theory Exp. (2019) 043402. [33] Z. Wang, C.T. Bauch, S. Bhattacharyya, A. d’Onofrio, P. Manfredi, M. Perc, N. Perra, M. Salathé, D. Zhao, Statistical physics of vaccination, Phys. Rep. 664 (2016) 1–113. [34] M.R. D’Orsogna, M. Perc, Statistical physics of crime: A review, Phys. Life Rev. 12 (2015) 1–21. [35] M. Perc, J.J. Jordan, D.G. Rand, Z. Wang, S. Boccaletti, A. Szolnoki, Statistical physics of human cooperation, Phys. Rep. 687 (2017) 1–51. [36] Z. Wang, M. Jusup, R.W. Wang, L. Shi, Y. Iwasa, Y. Moreno, J. Kurths, Onymity promotes cooperation in social dilemma experiments, Sci. Adv. 3 (2017) e1601444. [37] C. Du, D. Jia, L. Jin, L. Shi, The impact of neutral reward on cooperation in public good game, Eur. Phys. J. B 91 (2018) 234. [38] X. Li, M. Jusup, Z. Wang, H. Li, L. Shi, B. Podobnik, H. Stanley, S. Havlin, S. Boccaletti, Punishment diminishes the benefits of network reciprocity in social dilemma experiments, Proc. Natl. Acad. Sci. USA 115 (2018) 30–35. [39] A. Szolnoki, M. Perc, Effectiveness of conditional punishment for the evolution of public cooperation, J. Theoret. Biol. 325 (2013) 34–41. [40] B. Rockenbach, M. Milinski, The efficient interaction of indirect reciprocity and costly punishment, Nature 444 (2006) 718–723. [41] X. Chen, A. Szolnoki, M. Perc, Probabilistic sharing solves the problem of costly punishment, New J. Phys. 16 (2014) 083016. [42] D. Helbing, A. Szolnoki, M. Perc, G. Szabó, Punish, But not too hard: how costly punishment spreads in the spatial public goods game, New J. Phys. 12 (2010) 083005. [43] M. Perc, A. Szolnoki, Self-organization of punishment in structured populations, New J. Phys. 14 (2012) 43013–43025. [44] D. Helbing, A. Szolnoki, M. Perc, G. Szabó, Defector-accelerated cooperativeness and punishment in public goods games with mutations, Phys. Rev. E 81 (2010) 057104. [45] A. Szolnoki, G. Szabo, L. Czako, Competition of individual and institutional punishments in spatial public goods games, Phys. Rev. E 84 (2011) 046106. [46] A. Szolnoki, G. Szabo, M. Perc, Phase diagrams for the spatial public goods game with pool punishment, Phys. Rev. E 83 (2011) 036101. [47] J.H. Fowler, Altruistic punishment and the origin of cooperation, Proc. Natl. Acad. Sci. USA 102 (2005) 7047–7049. [48] C. Hauert, A. Traulsen, H. Brandt, M.A. Nowak, K. Sigmund, Via freedom to coercion: The emergence of costly punishment, Science 316 (2007) 1905–1907. [49] X. Chen, T. Sasaki, M. Perc, Evolution of public cooperation in a monitored society with implicated punishment and within-group enforcement, Sci. Rep. 5 (2015) 17050. [50] M. Perc, A. Szolnoki, A double-edged sword: Benefits and pitfalls of heterogeneous punishment in evolutionary inspection games, Sci. Rep. 5 (2015) 11027. [51] J. Quan, W. Liu, Y. Chu, X. Wang, Stochastic evolutionary voluntary public goods game with punishment in a quasi-birth-and-death process, Sci. Rep. 7 (2017) 16110. [52] J. Quan, Y. Chu, W. Liu, X. Wang, X. Yang, Stochastic evolutionary public goods game with first and second order costly punishments in finite populations, Chin. Phys. B 27 (2018) 060203. [53] P.A.I. Forsyth, C. Hauert, Public goods games with reward in finite populations, J. Math. Biol. 63 (2011) 109–123. [54] T. Sasaki, T. Unemi, Replicator dynamics in public goods games with reward funds, J. Theoret. Biol. 287 (2011) 109–114. [55] K. Sigmund, C. Hauert, M.A. Nowak, Reward and punishment, Proc. Natl. Acad. Sci. USA 98 (2001) 10757–10762. [56] T. Sasaki, S. Uchida, X. Chen, Voluntary rewards mediate the evolution of pool punishment for maintaining public goods in large populations, Sci. Rep. 5 (2015) 8917. [57] X.W. Wang, S. Nie, L.L. Jiang, B.H. Wang, S.M. Chen, Role of delay-based reward in the spatial cooperation, Phys. A 465 (2017) 153–158. [58] D. Masclet, Ostracism in work teams: a public good experiment, Int. J. Manpow. 24 (2003) 867–887. [59] M. Cinyabuguma, T. Page, L. Putterman, Cooperation under the threat of expulsion in a public goods experiment, J. Public Econ. 89 (2005) 1421–1435. [60] F.P. Maier-Rigaud, P. Martinsson, G. Staffiero, Ostracism and the provision of a public good: experimental evidence, J. Econ. Behav. Organ. 73 (2010) 387–395. [61] A.P. Melis, B. Hare, M. Tomasello, Chimpanzees recruit the best collaborators, Science 311 (2006) 1297–1300. [62] R. Bshary, A.S. Grutter, Punishment and partner switching cause cooperative behaviour in a cleaning mutualism, Biol. Lett. 1 (2005) 396–399. [63] R. Kurzban, M.R. Leary, Evolutionary origins of stigmatization: the functions of social exclusion, Psychol. Bull. 127 (2001) 187–208. [64] T. Sasaki, S. Uchida, The evolution of cooperation by social exclusion, Proc. R. Soc. B 280 (2013) 20122498. [65] K. Li, R. Cong, T. Wu, L. Wang, Social exclusion in finite populations, Phys. Rev. E 91 (2015) 042810. [66] K. Li, R. Cong, L. Wang, Cooperation induced by random sequential exclusion, Epl 114 (2016) 58001. [67] L. Liu, X. Chen, A. Szolnoki, Competitions between prosocial exclusions and punishments in finite populations, Sci. Rep. 7 (2017) 46634. [68] X.K. Sui, B. Wu, L. Wang, Rationality alters the rank between peer punishment and social exclusion, Epl 121 (2018) 38003. [69] A. Szolnoki, X. Chen, Alliance formation with exclusion in the spatial public goods game, Phys. Rev. E 95 (2017) 052316. [70] A. Szolnoki, M. Perc, Impact of critical mass on the evolution of cooperation in spatial public goods games, Phys. Rev. E 81 (2010) 057101. [71] M. Archetti, I. Scheuring, Coexistence of cooperation and defection in public goods games, Evolution (N Y) 65 (2011) 1140–1148. [72] H. Ye, S. Chen, J. Luo, F. Tan, Y.M. Jia, Y.F. Chen, Increasing returns to scale: The solution to the second-order social dilemma, Sci. Rep. 6 (2016) 31927.

14

J. Quan, J. Zheng, X. Wang et al. / Physica A 532 (2019) 121866

[73] N.J. Haering, Interjurisdictional redistribution and public goods with increasing returns to scale, Public Choice 95 (1998) 321–329. [74] T. Wu, L. Wang, Adaptive play stabilizes cooperation in continuous public goods games, Phys. A 495 (2018) 427–435. [75] J. Wang, B. Wu, X.J. Chen, L. Wang, Evolutionary dynamics of public goods games with diverse contributions in finite populations, Phys. Rev. E 81 (2010) 056103. [76] J. Wang, X.J. Chen, L. Wang, Effects of migration on the evolutionary game dynamics in finite populations with community structures, Phys. A 389 (2010) 67–78. [77] J. Quan, W. Liu, Y. Chu, X. Wang, Stochastic dynamics and stable equilibrium of evolutionary optional public goods game in finite populations, Phys. A 502 (2018) 123–134. [78] P. Young, The evolution of conventions, Econometrica 61 (1993) 57–84. [79] W.J. Stewart, Introduction to the Numerical Solution of Markov Chains, Princeton University Press, Princeton, 1994. [80] M. Perc, J. Gómez-Gardeñes, A. Szolnoki, L.M. Floría, Y. Moreno, Evolutionary dynamics of group interactions on structured populations: a review, J. R. Soc. Interface 10 (2013) 20120997. [81] A. Szolnoki, X. Chen, Benefits of tolerance in public goods games, Phys. Rev. E 92 (2015) 042813. [82] X. Chen, A. Szolnoki, Individual wealth-based selection supports cooperation in spatial public goods games, Sci. Rep. 6 (2016) 32802. [83] M. Perc, Stability of subsystem solutions in agent-based models, Eur. J. Phys. 39 (2018) 014001.