Diffusively-Coupled Rock-Paper-Scissors Game with Mutation in Scale-Free Hierarchical Networks

We present a metapopulation dynamic model for the diﬀusively-coupled rock-paper-scissors (RPS) game with mutation in scale-free hierarchical networks. We investigate how the RPS game changes by mutation in scale-free networks. Only the mutation from rock to scissors ( R -to- S ) occurs with rate μ . In the network, a node represents a patch where the RPS game is performed. RPS individuals migrate among nodes by diﬀusion. The dynamics are represented by the reaction-diﬀusion equations with the recursion formula. We study where and how species coexist or go extinct in the scale-free network. We numerically obtained the solutions for the metapopulation dynamics and derived the transition points. The results show that, with increasing mutation rate μ , the extinction of P species occurs and then the extinction of R species occurs, and ﬁnally only S species survives. Thus, the ﬁrst and second dynamical phase transitions occur in the scale-free hierarchical network. We also show that the scaling law holds for the population dynamics which suggests that the transition points approach zero in the limit of inﬁnite size.


Introduction
Recently, rock-paper-scissors (RPS) games have been extensively investigated in various fields, such as physics and ecology. RPS games are known as one of the cyclic dominance models  and capture the central features of cyclic systems, which are abundant in nature [26][27][28]. Spatial and temporal behaviors of RPS games are key to understanding not only biodiversity but also a variety of properties in cyclic systems. Stochastic cellular automaton models for the RPS games on lattices reveal that three species sometimes can coexist but sometimes cannot. e spatial RPS games have been extended to small-world networks [25]. Masuda and Konno have studied the RPS games without mutation on scale-free networks [29]. It has been shown that the network's heterogeneity can help to maintain the stable coexistence of species in a cyclic competition.
RPS games in networks have usually been studied by an agent-based model, namely, a node denotes an agent and a link means the interaction between agents. In previous papers, we have applied a metapopulation model to the RPS games where each node denotes a patch which contains many agents. e metapopulation model consists of spatially separated patches of species. Agents can migrate between patches. Small sets of patches are represented by a graph (small network). en, the metapopulation dynamics have been represented by reaction-diffusion equations on the graph. It has been shown that the dynamics are stabilized by the graph's heterogeneity [3].
Mutation is a basic mechanism in evolution. For example, in microorganisms, particular species can change phenotypes by mutation [30,31]. us, mutation affects the stability and population dynamics of competing species [1]. RPS games with mutation have been considered in many studies [1,[32][33][34][35][36]. Traditionally, mutation was introduced as a form that a particular species changes to another competing species with specific rates [1,32,33]. On the contrary, various mutations among three species where mutation occurs in certain targeted species have been considered [34,35]. ey showed such mutation can promote biodiversity by exhibiting multistability of certain survival states of subpopulation groups. Park considered a different form of mutation in RPS games where a new species emerges from a native species and showed the complex coexistence of the mutant and native species [36]. Nagatani et al. have studied the effect of mutation on the RPS games on small graphs [4]. It has been shown that the dynamical phase transitions induced by mutation can be divided into two stages with increasing mutation rate: the first transition from the coexistence of RPS species to the coexistence of R and S species, and the second transition from the coexistence of R and S species to only S species.
Until now, the dynamical stability for the coexistence of three species has been studied in networks. However, the dynamical phase transitions between coexistence and extinction have not been investigated for the RPS games with mutation in large-scale networks. It is poorly understood how the metapopulation dynamics and phase transitions change in a large-scale network. It is also not known how the phase transition points change in the limit of infinity system size. How does the network's structure affect the dynamical phase transition in the diffusively-coupled RPS game with mutation? In the scale-free network, can three species coexist or not? It is also unknown whether or not the scaling law holds. Scaling properties are important in the sense that they characterize the behavior of physical quantities near continuous phase transitions. It is believed, though not proven, that they are universal. We investigated the scaling properties since the properties are unknown in the dynamical phase transitions of RPS games with mutation.
In this paper, we study the metapopulation dynamics and dynamical phase transitions of the diffusively-coupled RPS game with mutation in a hierarchical network with scale-free and self-similar properties where the dynamics are represented by reaction-diffusion equations. e hierarchical structure of a network is described by the simple recursion formula for the reaction-diffusion equations. We carry out numerical calculations and obtain the numerical solutions at a steady state. We derive the dependence of RPS densities on the mutation rate. We show where and how the extinction of species occur in the network. We clarify the dependence of the transition points (thresholds) on the highest degree.

Mean-Field and Metapopulation Models
In the metapopulation dynamics of diffusively-coupled RPS games with mutation in a network, a node represents a patch (subpopulation) in which three species (rock, paper, and scissors) live and a link corresponds to a migration path on which three species move by diffusion. e three species are mixed well within a patch. e three species have a cyclic relation: species R eats S, S eats P, and P eats R. A mutation is introduced as R ⟶ S, where R individuals change to S at rate μ. Traditionally, there are two standard approaches in the dynamics of RPS games, either the Lotka-Volterra limit or the May-Leonard limit. Phase transitions depend on whether the total density is conserved or not. We use the mean-field dynamics based on the Lotka-Volterra limit because it is a conserved system. e mean-field dynamics (MFD) hold for each subpopulation. e MFD for the RPS games with mutation are given by where ρ R (t), ρ S (t), ρ P (t) are the densities of species R, S, and P at time t, respectively. For simplicity, we set all three victory rates to one. Without mutation (μ � 0), the RPS dynamics of equations (1)-(3) show oscillatory behavior (neutrally stable). However, the dynamics change to stable focuses by adding mutation (μ > 0). In the equilibrium state, RPS densities show the following behavior [4]. With increasing rate μ, paper density decreases and becomes zero at μ � 0.5. Rock density equals scissors density until μ � 0.5. For 0.5 < μ < 1, rock density decreases linearly with increasing μ and becomes zero at μ � 1, while scissors density increases linearly with increasing μ for 0.5 < μ < 1 and becomes one at μ � 1. With increasing mutation rate μ, furthermore, only scissors survive for μ > 1. us, the first phase transition between the coexistence of three species (R, P, and S species) and the coexistence of two species (R and P species) occur at μ � 0.5. e second phase transition between the coexistence of two species (R and P species) and a single species (S species) occurs at μ � 1.
Here, we consider the dynamics of three species in a network. e number of subpopulations is represented by N. e structure of the metapopulation is expressed by a network. Subpopulation (patch) i is expressed as node i in the network. e connection between two subpopulations is expressed by a link (edge). Many individuals walk randomly among nodes through links.
us, for the diffusively-coupled RPS game in a network, the dynamics of individuals are simply represented by reaction-diffusion equations.
In the metapopulation model of the RPS game with mutation, densities ρ R,i (t), ρ S,i (t), ρ P,i (t) of RSP species in node (subpopulation) i are described by 2 Complexity where D is the diffusion constant of individuals [3,4]. In equation (4), the first and second terms represent the inflow and outflow of R at node i.
to the increase of R at node i by migrating into the node while − Dρ R,i (t) corresponds to the decrease of R at node i by migrating out from the node. ere are two kinds of random walks: node-centric and edge-centric random walks [37]. On homogeneous lattices, the node-centric random walks agree with the edge-centric random walks. However, in heterogeneous networks, the node-centric random walks do not agree with the edge-centric random walks. In our model, a node corresponds to a patch. As described in the first and second terms of equations (4)- (6), we apply the node-centric random walks to the diffusion process in networks. e third and fourth terms represent the change of R species by the RPS game in subpopulation i. e fifth term represents the decrease of R species by mutation. Equations (5) and (6) are defined in the same way.

Metapopulation in Hierarchical Network
We extend the above metapopulation dynamics of the diffusively-coupled RPS game in a network to the case of a hierarchical network with scale-free and self-similar properties. e first, second, and third generations of the hierarchical network are shown respectively in Figures 1(a)-1(c). e n th generation is obtained through repetition. We study the effect of mutation rate μ on the metapopulation dynamics in this hierarchical network.
We derive the density equations of RPS species in the hierarchical network by applying equations (4)-(6) to metapopulations for the second generation of hierarchical network in Figure 1(b). When the initial values of densities are the same on all nodes, the time-varying densities behave symmetrically. e densities at the nodes on the upper-left side of Figure 1(b) are consistent with those on the upperright, bottom-left, and bottom-right sides. e densities of RSP species at node i for the second generation are represented by ρ 2,R,i (t), ρ 2,S,i (t), and ρ 2,P,i (t). en, the density equations at nodes 1 and 2 are given by dρ 2,P,2 (t) dt � 1 2 ρ 2,P,1 (t) − ρ 2,P,2 (t) + ρ 2,P,2 (t)ρ 2,R,2 (t) − ρ 2,S,2 (t) 2,P,2 (t).

Complexity
Here, coefficients C 1 (n, i), C 2 (n, i) (n > 3) are given by Recursion formula (22) is obtained for the coefficients of density equations (13)- (21). One can calculate the time evolution of RPS densities on all nodes in the hierarchical network by using equations (13)-(21) with (22).

Numerical Result
We carry out numerical calculations for the RPS game with mutation in the hierarchical network. We consider a static network for each generation. At each generation, the network structure is fixed during the numerical simulations. We calculate the densities of RPS individuals on all nodes by using density equations (7)- (22). For each simulation, the maximum iteration for each difference equation is 10 8 where the interval is Δt � 0.0001. In the diffusively-coupled RPS game with mutation, the sum of RPS individuals is conserved over the network. We set the total density as one, ρ 0 � 1.
First, we calculate the densities of RPS individuals for the second generation in Figure 1(b) by solving equations (7)- (12) numerically. e densities of RPS individuals at a steady state are obtained after a sufficient length of time. We omit the subscript representing the generation number n and add the subscript e indicating the steady state. e steady-state densities on node i are indicated by ρ R,i,e , ρ S,i,e , and ρ P,i,e . Figure 2(a) shows the plots of steady-state densities ρ R,i,e , ρ S,i,e , and ρ P,i,e (i � 1, 2) of RSP individuals against mutation rate μ of the second generation (n � 2) in Figure 1(b). Black, red, and green circles indicate the steadystate densities of RSP individuals, respectively. Small and large circles indicate the steady-state densities on nodes 1 and 2, respectively. P densities ρ P,1,e and ρ P,2,e at nodes 1 and 2 decrease with increasing mutation rate. R and S densities at nodes 1 and 2 increase with the mutation rate. RSP densities at node 1 are higher than those at node 2. When the mutation rate reaches a first critical value (transition point or threshold), P individuals disappear at nodes 1 and 2 simultaneously. Above the first transition point, P densities at nodes 1 and 2 are zero. e first transition point (threshold) 0.190 is lower than the mean-field result 0.5. By increasing the mutation rate, furthermore, S densities at nodes 1 and 2 increase and R densities at nodes 1 and 2 decrease. When the mutation rate reaches a second transition point, R individuals disappear at nodes 1 and 2 simultaneously. Above the second transition point, R densities at nodes 1 and 2 are zero. As a result, only S individuals survive. e second transition point (threshold) 0.383 is lower than the meanfield result 1.0. By increasing the mutation rate, furthermore, S densities at nodes 1 and 2 remain constant. S density at node 1 is higher than that at node 2. RPS individuals gather at node 1 at a higher degree than node 2.
Next, we calculate the densities of RPS individuals on the fourth generation by solving equations (13)-(22) numerically. e densities of RPS individuals at a steady state are obtained after a sufficient length of time. e steady-state densities are indicated by ρ R,i,e , ρ S,i,e , and ρ P,i,e (i � 1 − 8). Figure 2(b) shows the plots of steady-state densities ρ R,i,e , ρ S,i,e , and ρ P,i,e (i � 1, 2, 5) of RSP individuals against mutation rate μ of the fourth generation (n � 4). Black, red, and green circles indicate the steady-state densities of RSP individuals, respectively. Small and large circles indicate the steady-state densities on nodes 1 and 2, respectively. Crosses indicate the steady-state densities of RSP individuals on node 5. P densities ρ P,i,e (i � 1 − 8) at all nodes decrease with increasing mutation rate. R and S densities at all nodes increase with the mutation rate. RSP densities at node 1 are highest compared to those at other nodes. Each density of RSP individuals decreases with the decreasing node's degree. When the mutation rate reaches the first transition point, P individuals disappear at all nodes simultaneously. Above the first transition point, P densities at all nodes are zero. e first transition point (threshold) 0.020 is significantly lower than 0.190 of the second generation. By increasing the mutation rate, furthermore, S densities at all nodes increase and R densities at all nodes decrease. When the mutation rate reaches the second transition point, R individuals disappear at all nodes simultaneously. Above the second transition point, R densities at all nodes are zero. As a result, only S individuals survive. e second transition point (threshold) 0.040 is significantly lower than 0.383 of the second generation. By increasing the mutation rate, furthermore, S densities at all nodes remain constant. S density at node 1 is higher than those at other nodes. S density decreases with the decreasing node's degree. RPS individuals gather at node 1 at the highest degree.
In the limit of μ � 0, no mutation occurs. Steady-state RSP densities are the same. We study the dependence of steady-state densities at μ � 0 on the node's degree. Figure 3(a) shows the log-log plot of RPS densities against node's degree k at mutation rate μ � 0 for generation number n � 6 where the highest degree is 2 6 . Black, red, and green circles indicate RSP densities in the limit of μ � 0, respectively. e RSP densities are on the linear line with slope 1. e following scaling is obtained for RPS games without mutation in the hierarchical network: Complexity e RSP densities are proportional to the degree of nodes.
Also, only scissors individuals survive at a sufficiently high value of mutation rate. We study the dependence of steady-state S densities at μ � 0.5 on the node's degree. Figure 3(b) shows the log-log plot of S densities against node's degree k at mutation rate μ � 0.5 for generation number n � 6. Red circles indicate S densities at μ � 0.5. e S densities are on the linear line with slope 1. e following scaling is obtained for RPS games with a high value of mutation rate in the hierarchical network: e S densities are proportional to the degree of nodes. e first and second transition points decrease with increasing generation number n. We study the dependence of first and second transition points on highest degree k max (� 2 n ) at node 1. Figure 4 shows the log-log plot of first and second transition points μ c1 and μ c2 against highest degree k max . Black and red circles indicate transition points μ c1 and μ c2 , respectively. e transition points are on the linear lines. e slope of the first transition point agrees with that of the second transition point. e following scaling is obtained: e scaling exponent of the first and second transition points is lower than − 1. us, the first and second transition points decrease with the increasing highest degree and have scaling form (25). Because those points have negative scaling exponents, they approach zero when the generation number (maximal degree) approaches infinity.  Figure 2: (a) Plots of steady-state densities ρ R,i,e , ρ S,i,e , and ρ P,i,e (i � 1, 2) of RSP individuals against mutation rate μ for a diffusively-coupled RPS game with mutation for the second generation (n � 2) in Figure 1(a). Black, red, and green circles indicate the steady-state densities of RSP individuals, respectively. Small and large circles indicate the steady-state densities on nodes 1 and 2, respectively. (b) Plots of steadystate densities ρ R,i,e , ρ S,i,e , and ρ P,i,e (i � 1, 2, 5) of RSP individuals against mutation rate μ for the diffusively-coupled RPS game with mutation on the fourth generation (n � 4). Black, red, and green circles indicate the steady-state densities of RSP individuals, respectively. Small and large circles indicate the steady-state densities on nodes 1 and 2, respectively. Crosses indicate the steady-state densities of RSP individuals on node 5.  Figure 3: (a) Log-log plot of RPS densities against node's degree k at mutation rate μ � 0 for generation number n � 6 where highest degree is 2 6 . Black, red, and green circles indicate RSP densities in the limit of μ � 0, respectively. e RSP densities are on the linear line with slope 1. (b) Log-log plot of S densities against node's degree k at mutation rate μ � 0.5 for generation number n � 6. Red circles indicate S densities at μ � 0.5. e S densities are on the linear line with slope 1.

Complexity
By diffusion, densities on hubs increase while those of the other nodes decrease. e transition point increases as the density increases. In contrast, it decreases as the density decreases. Since the number of nonhub nodes is much larger than that of hub nodes, the transition point decreases as a whole. Finally, the density per node approaches zero with k max going to infinity. As a result, the transition points approach zero. In the limit of infinite number of generations, the mutation reaction becomes dominant and only S species survives.

Summary
e RPS games without mutation have been investigated by reaction-diffusion equations on one-dimensional and twodimensional lattices. In the previous paper, we have studied the effect of both graph's structures and mutation on metapopulation dynamics of the RPS games on various graphs [4]. By introducing R-to-S mutation into the metapopulation model of the RPS games, the first and second dynamical phase transitions of P and R extinctions occur. However, the diffusively-coupled RPS games with mutation have not been investigated in scale-free networks. Here, we studied how the scale-free hierarchical network affects the dynamical phase transitions of the three species. us, we presented reaction-diffusion equations for the metapopulation RPS game with mutation in the hierarchical network. e reaction-diffusion equations are described by means of recursion formula. By numerical calculations, we found that the first dynamical phase transition of P extinction occurs at all nodes simultaneously by increasing the mutation rate and the second dynamical phase transition of R extinction also occurs at all nodes simultaneously. We found that RPS densities show a scaling behavior which follows equation (23) at μ � 0 and S density shows a scaling behavior which follows equation (24) after the second phase transition. When the generation number approaches infinity, the first and second transition points approach zero. If R-to-S mutation occurs at a finite rate, only S species survive in an infinite size. Our study contributes to the further understanding of the coexistence of cyclic biological systems.

Data Availability
No data were used in this study.

Conflicts of Interest
e authors declare that they have no conflicts of interest. Complexity 7