Absorbing phase transition in the coupled dynamics of node and link states in random networks

We present a stochastic dynamics model of coupled evolution for the binary states of nodes and links in a complex network. In the context of opinion formation node states represent two possible opinions and link states a positive or negative relation. Dynamics proceeds via node and link state update towards pairwise satisfactory relations in which nodes in the same state are connected by positive links or nodes in different states are connected by negative links. By a mean-field rate equations analysis and Monte Carlo simulations in random networks we find an absorbing phase transition from a dynamically active phase to an absorbing phase. The transition occurs for a critical value of the relative time scale for node and link state updates. In the absorbing phase the order parameter, measuring global order, approaches exponentially the final frozen configuration. Finite size effects are such that in the absorbing phase the final configuration is reached in a characteristic time that scales logarithmically with system size, while in the active phase, finite-size fluctuation take the system to a frozen configuration in a characteristic time that grows exponentially with system size. There is also a finite-size topological transition associated with group splitting in the network of these final frozen configurations.

node is considered in most studies of opinion formation, as for example in 27,28 . A step beyond these studies are those, either in opinion formation (friendly or unfriendly link) 15 or epidemics (active or inert link) 29,30 , in which both nodes and links have a state, but the state of the link is either fixed or, else, determined by the state of the connecting nodes, so that, still, there is no coupled dynamics of node and link states.
The dynamical interplay between the states of the nodes and the state of the links has been considered in the context of an Ising model 18 and a susceptible-infected model 17 . In both cases the dynamics to approach structural balance is described by energy minimization of an appropriate Hamiltonian in a fully connected network. A different approach is that of Carro et al. 16 in the context of language competition (node state: language preference, link state: language use) where a genuine non-equilibrium dynamics with no Hamiltonian minimization is implemented in a complex network. This leads to a wide range of asymptotic states, including long-lived dynamical states. Here we also consider a non-equilibrium dynamics in the context of opinion formation and we focus on binary or pairwise relations: Nodes can be in either of two states or opinions and links can be positive (friendly) or negative (unfriendly). Satisfactory pairwise relations are those of friendly links connecting nodes in the same state or unfriendly links connecting nodes in different states. The basic assumption is that unsatisfying pairs evolve to satisfying ones either by updating the state of the link or by updating the state of one of the nodes (see Fig. 1). We note that the node-link state dynamical coevolution considered here is essentially different of the one considered in coevolving voter models 31,32 . In this latter case a link has no state, but the network evolves dynamically with a link connecting nodes in different states being rewired with some probability to connect nodes in the same state. Although one can associate a binary state of the link as existing or non existing, the dynamics is local in our present model, with a change of the state of the link only between neighbors. However, there is a non local network dynamics of random rewiring in the coevolving voter model.
As a main result of our study, we find a phase transition from a dynamically active state to an absorbing frozen configuration in which all pairwise relations are satisfactory. This transition shares many features with the one found in the coevolving voter model 31 . Here, the absorbing phase occurs beyond a critical value p c of the parameter that measures the ratio of time scales for links and nodes updates. This implies that, despite a dynamical rule of local convergence, a global convergence to the absorbing state only occurs when the states of the links evolve fast enough in comparison with the evolution of the states of the nodes. The two phases separated by p c can be also characterized dynamically: From general random initial conditions and for > p p c the evolution of an order parameter that measures global ordering indicates an exponential approach to the absorbing state, while for < p p c the order parameter falls into a plateau value characterizing partial ordering of the system. While the above scenario is valid in the thermodynamic limit where the number of nodes tends to infinite, finite-size effects introduce new time scales: In the absorbing phase the frozen configuration is reached in a characteristic time that scales logarithmically with system size, while in the active phase finite-size fluctuations bring the system to a frozen configuration in a characteristic time that scales exponentially with system size. In addition, we find a transition associated with the network topology of the frozen configurations, whether they have been reached in the absorbing phase or as a result of a finite-size fluctuation in the active phase. This is a finite-size transition that disappears in the thermodynamic limit. It manifests itself via group splitting in the network and we refer to it as a finite-size topological transition. This paper is organized as follows. After introducing our dynamical model, we discuss a mean-field rate equation approach, suitable for homogeneous random networks 33 , that is based on the assumption that each node has exactly μ neighbors. These equations, that neglect fluctuations and are only valid in the thermodynamic limit, predict a continuous transition between absorbing and active phases and allow the calculation of p c as a function of the mean degree of the network μ. Next we report Monte Carlo simulations of the model on Erdős-Rényi networks that confirm the rate equation predictions, and describe the dynamical properties of the active and frozen phase, including finite-size effects. We finally summarize our results.

Results
Coupled evolution of node and links property in the imitating process. Let us consider a connected network defined as a set of nodes and links. The nodes represent individuals and the links, understood as undirected connections, indicate a relation between the nodes. Each node holds a binary state variable whose value represents one of two possible opinions. In the figures those two possible values are indicated by a dark (blue) or white color. The links between nodes represent one of the two possible types of relationship: friendly (attraction) and unfriendly (repulsion). In the figures they are indicated, respectively, by a continuous or a dashed line. According to the aforementioned interpretation, we consider that friendly links between pairs of nodes Figure 1. We present in this figure all six possible configurations of pairs constructed by two types of link states, friendly (solid line) and unfriendly (dashed line), and two types of node states, blue opinion and white opinion. Thus, the pairs a, c and e are unsatisfying while pairs b, d and f are satisfying. Moreover, we also depict the dynamical rules that turn unsatisfying pairs into satisfying ones through node or link updates. p is the probability of link update and − p 1 , the complementary probability, is the probability of node update.
holding the same opinion or unfriendly links between pairs of nodes holding different opinions are satisfying links (or form satisfying pairs), while friendly links between nodes holding different opinions or unfriendly links between nodes holding the same opinion are unsatisfying. All possible situations are displayed in Fig. 1: pairs a, c and e are unsatisfying, while pairs b, d and f are satisfying. Our basic assumption is that people in a community act in order to maximize their level of satisfaction, a desirable option from the psychological viewpoint. Hence, individuals belonging to an unsatisfying pair connection take action to turn it into satisfying. To convert unsatisfying pairs into satisfying ones, individuals can either change their opinions or alter the link state. Consider, for instance, the unsatisfying a pair in Fig. 1. It can become satisfying if one of the individuals changes its opinion state from white to blue or, alternatively, if they decide to change their link state from unfriendly to friendly. We consider that these two options happen with probability − p  1 and p, respectively. A similar scenario holds for the unsatisfying link c. In the case e we need to assign a probability − p (1 )/2 to the node holding the blue opinion changing it to white and a probability − p (1 )/2 to the node holding the white opinion changing it to blue. In the pairs a, c, it does not matter which one node changes opinion, as it always leads to an f pair.
In a Monte Carlo implementation of these dynamical rules, a link is chosen at random from all existing links. If the corresponding pair is satisfying, nothing happens; otherwise, if the pair is unsatisfying, it is converted into satisfying by applying the rules, with their respective probabilities, displayed in Fig. 1. A Monte Carlo step, as usual, is defined as a number of consecutive link selections equal to the total number of links existing in the system. Most of our results concern an uncorrelated Erdős-Rényi network with N nodes and average degree (number of links per node) μ such that the total number of links is . When all pairs are satisfying, no further evolution is possible and the system is dynamically frozen in an absorbing state. Note that an update that converts a pair from unsatisfying to satisfying by changing an individual node state might change the status of another pair, to which the node involved in the update also belongs to, from satisfying to unsatisfying. Therefore, we ask the question of under which conditions the dynamical rules defined above lead to an absorbing, all pairs being satisfying, global state.
In the frozen state, the densities ρ a , ρ c , ρ e are zero and hence in order to determine whether the frozen state has been reached in a particular realization, we focus on the time evolution of the link densities ρ t c d e f { , , , , , }. This evolution has been analyzed, either using the Monte Carlo procedure explained before or by a set of approximate rate equations, derived in the Method section.
In the Erdős-Rényi Network, where links amongst nodes are generated randomly with probability μ/N, it is possible to relate the different link densities with the number n(t) of nodes holding the white opinion. The relations are Rate equation, fixed points and critical line. As described in the Method section, we use a mean-field approximation to derive the rate equations of the aforementioned dynamics, based on the update rules sketched in Fig. 1. The main assumption in the derivation is that each node has exactly μ neighbors distributed randomly amongst all possible nodes. The mean-field treatment is exact in the thermodynamic limit (an infinite system size N) of an all-to-all network where all nodes are connected to each other, but it is only an approximation for other networks that we consider in our simulations. This is the case of an Erdős-Rényi network with a Poisson distribution, of mean degree μ, for the number of neighbors. Furthermore, as the rate equations are deterministic, they can not describe the finite-size fluctuations observed in the numerical simulations. The rate equations are six coupled differential equations for the time evolution of the six types of pair densities ρ ρ ρ ρ ρ ρ { , , , , , } a b c d e f and include a combination of linear and nonlinear terms of those variables, see Eq. (6). The linear terms reflect the direct change of densities in any type of update, either node or link, and the nonlinear terms are the indirect consequence of node update, by which the status of other links connected to the updated node will also change. The latter effect plays a crucial role in the evolution of the system. We first identify the fixed points of our 6th-order dynamical system, found by setting all time derivatives equal to zero. It turns out that there are two independent sets of fixed points that can be easily obtained by a standard software of mathematical analysis. In the first set all densities take a well defined, non-null, value www.nature.com/scientificreports www.nature.com/scientificreports/ while in the second set there are no unsatisfying pairs, but there is an arbitrariness in the values of the densities of the satisfying links (always verifying the normalization condition ρ ρ ρ ). In the first solution, Eq. (2), the densities of the pair connections reach asymptotically a non-null plateau depending on p and μ and independent of initial conditions. However, in the second set of solutions, Eq. (3), there are no unsatisfying pair connections, and the densities of the satisfying ones depend on the initial conditions.
The condition ρ ≥ 0 e st determines that the first solution is only relevant for A linear stability analysis 34 shows that the first solution, Eq. (2), is always stable (negative eigenvalues of the linearized equations) whenever it leads to non-negative densities, and that the set, Eq. (3), is marginally stable (eigenvalues equal to zero). Therefore for μ < p p ( ) c there are always active links in the asymptotic state and, consequently, continuous changes of the microscopic state. We say that the system stays in an active, or dynamical, phase. However, for μ > p p ( ) c , the densities of all active pairs tend to zero and the dynamics reaches a frozen, or absorbing, phase. As order parameter distinguishing one phase from another we choose the density ρ e st , which is zero in the absorbing phase and positive in the active phase. Note that, in turn, the condition μ ≥ p ( ) 0 c requires μ > 5/2. If μ ≤ 5/2 the first solution does not exist and the dynamics always leads to a frozen state. As it follows from Eq. (2), as p approaches the critical value p c from below, the order parameter tends to zero as with a critical exponent β = 1, a continuous phase transition. As this coincides with the mean-field critical exponent of directed percolation 35 , we speculate that the active-frozen transition that we just described can be categorized under the class of directed percolation, but much more extensive numerical simulations would be needed in order to be conclusive in this point. These conclusions are based on the study of the rate equations and are strictly valid only in the thermodynamic limit for the all-to-all network. In a finite system, an active state will display fluctuations of the order parameter around its mean value. Due to the stochastic nature of the dynamics, there will be always a fluctuation that takes an active state into a frozen one and, from there on, all microscopic dynamics stops. As discussed in detail in the next sections, the likeness of such fluctuation tends to zero with increasing system size and the average time to reach the frozen state diverges exponentially with system size. phase transition. We have carried out extensive Monte-Carlo simulations of the dynamical rules presented in Fig. 1 on an Erdős-Rényi network with mean degree μ and system size N. Once the network, nodes and links, has been constructed, we assign randomly a state (white/blue) to each node and then a state to each link (friendly/ unfriendly). The initial density of white nodes is x 0 and that of friendly links is  0 . Note that the all-to-all network corresponds to μ = − N 1. Some representative trajectories of the order parameter ρ t ( ) e can be seen in Fig. 2 for = N 400 for some system parameters leading to the active phase, panel (a), or to the frozen phase, panel (b). In all curves we have taken the same values of = .
x 0 5 0 and = .  0 5 0 but different realizations of the networks, initial conditions and the dynamics. Observe the dispersion in the different curves in the active phase and that the microscopic dynamics continues in this phase for the whole range of time displayed in the figure. In the frozen phase, there is no further dynamics when the density of unsatisfied links reaches zero. Note, however, that there is also a dispersion in the times it takes the different realizations of the dynamics to reach the frozen state. p 0 8 and μ = 12, leading to an active phase where the evolution continues for the whole range of time displayed. In panel (b) we take = .
p 0 8 and μ = 6, leading to a frozen phase with a zero density of unsatisfying links. All trajectories are generated using the Monte Carlo method described in the main text and, in both panels, start out from the same value of the density of white nodes, = .
x 0 5 0 , and the density of friendly links, = . In panels (a) and (c) we take μ = 399 which, for the Erdős-Rényi network, means that all nodes are connected to each other. In this case, as expected, the agreement between the simulations and the rate equations is very good. In panels (b) and (d) we take μ = 6, but still observe a good agreement between simulations and rate equations. For μ = 399, the chosen value of = .
p 0 9 is below the critical one μ < = . p p ( ) 0 996 c and, consequently, the dynamics leads to an active state. For μ = 6, on the contrary, the value of = .
p 0 9 is above the critical one μ > = .  Fig. 4( st st in the μ p ( , ) plane in a color code yielding plots which are indistinguishable from the corresponding Fig. 4(b,e), respectively, for ρ 〈 〉 f st in the frozen phase, A more detailed analysis of the dependence of x st and hence of ρ b st , ρ d st and ρ f st on the initial condition, in the frozen phase has been performed in Fig. 6(a,b). We set a value of x 0 (initial density of white nodes) and generate many initial conditions varying the value of the initial density of friendly links ∈  (0, 1) 0 . Each of these microscopic configurations evolves to for fixed μ = 6 have been plotted in Fig. 6(a,b), respectively, for different values of x 0 .
We have also carried out numerical simulations in regular lattices in one and two dimensions with nearest neighbors (results not shown). The qualitative phenomenology is the same than the one described previously: p 0 9 and start out from same value, = .
x 0 1 0 and = .   The existence of absorbing states and the ergodicity of the stochastic dynamical rules imply that for a finite system there is always a fluctuation that will take the system towards one of the absorbing states. Therefore, for a finite system, the ultimate fate is to end in the absorbing phase. The key point is to analyze the dependence of the time τ to reach the absorbing state on system size. τ is a random variable and we present in panels (a) and (b) of Fig. 8 its probability density function (pdf) τ f ( ) for two values of μ p ( , ) corresponding to the active phase and several values of the system size N. It appears from these figures that, at least for large τ, the pdf can be fitted by an exponential form τ ∼ τ τ − 〈 〉 f e ( ) / . The average value of the distribution τ 〈 〉 is plotted in panels (e) and (f) Fig. 8 as a function of the system size N, showing an exponential dependence τ 〈 〉 ∼ α e N . This exponential dependence indicates that the decay to the absorbing state becomes very rare for increasing system size. From the point of view of Statistical Mechanics, in the limit → ∞ N , the active state remains forever and represents a genuine macroscopic phase.  p 0 8 and initial condition = .
x 0 5 0 and = .  0 5 0 but different realizations of the network and the dynamics. In both cases, the averages use 100 runs that did not end in the frozen phase. In panel (b) we observe an exponential behavior ρ 〈 〉 ∝ τ − e e t/ 0 with a small dependence of τ 0 with system size. In panel (a), the fluctuations of ρ 〈 〉 e around the well defined value decrease with system size.
www.nature.com/scientificreports www.nature.com/scientificreports/ The active phase is remarkable in the sense that it indicates the failure of the local dynamical rule. While the evolution is dictated by a tendency to reduce the unsatisfying pairs, the final state is one of coexistence of all types of links. This is reminiscent of other dynamical models, the most notable being that of Axelrod 36,37 that predicts that local convergence can generate global polarization in an agent-based model of dissemination of culture.
From the topological point of view, the active phase is characterized by a continuously evolving, and apparently disordered, structure of nodes and links, see Fig. 9(a,b). The topology of the corresponding frozen phases that appear at a later time due a to a finite-size fluctuation, see Fig. 9(c,d), will be analyzed in another section.
Frozen phase. When the link-update probability is larger than the critical value μ > p p ( ) c the system falls into the absorbing or frozen phase. At variance with the system in the active phase, the densities of unsatisfied continuously decrease during the time evolution and never reach a plateau from which they eventually escape. Therefore, there is a continuous decay towards the absorbing phase, contrarily to the decay of the active phase that was produced by a rare fluctuation. This is evidenced in panel (b) of Fig. 7 where we plot the time evolution of ρ t ( ) e for different system sizes. We observe an exponential decay ρ ∼ τ − t e ( ) e t/ 0 with a very small dependence of τ 0 on the system size N and approaching a limiting value relatively close to the prediction of the rate equation. The time to reach the absorbing state, τ, is also a random variable characterized by a pdf τ f ( ). As shown in panels (c) and (d) Fig. 8, corresponding to two different points in the frozen phase, the tail of τ f ( ) can still be fitted by an exponential function τ − e t/ 2 , but this function now presents a well defined maximum located at τ 1 , the characteristic time for decay unto the frozen state. This characteristic time is plotted in panels (g) and (h) of Fig. 8 showing a logarithmic increase with system size τ ∼ N log 1 . It is possible to relate the exponential decay observed in ρ t ( ) e with this logarithmic dependence. The transition to the absorbing state will occur at the time τ 1 when the density of unsatisfying pairs falls below a value of order 1/N, i.e. ρ τ ∼ , a logarithmic dependence on system size, as observed. Snapshot of typical frozen configurations are displayed in Fig. 10(a,b). In the next section we discuss the possible topological structures of the frozen phases.
Finite-size topological transition. A fully satisfying, absorbing, configuration obtained in an Erdős-Rényi network displays a transition associated with some structure that can be described as "group splitting". By that, we mean that the nodes organize in several groups, defining a group as a set of nodes holding the same opinion and connected by friendly links among themselves and by unfriendly links to the members of other groups. This group splitting structure appears both when the absorbing configuration has been reached from a finite-size fluctuation of an active phase, μ < p p ( ) c , see Fig. 9(c,d), or when it corresponds to the frozen phase in parameter space, μ > p p ( ) c , see Fig. 10(a,b). The number of components is an important topological property of a network. We characterize the absorbing configurations by the number of groups N g . This number is a stochastic variable that depends on the initial con- www.nature.com/scientificreports www.nature.com/scientificreports/ dition, the system parameters and the particular realization of the dynamics. We denote by f N ( ) g #groups the probability distribution function of this variable as: and δ is the Dirac-delta function. The shape of these discrete distributions for small μ can be approximated by a continuous Gaussian shape, as shown in Fig. 11. What we have observed is that there is a transition in which f N ( ) g #groups changes from a unimodal distribution at = N 2 g (meaning that all realizations end up Figure 9. In panel (a) we display a snapshot of a configuration of the dynamics with parameters = .
p 0 3 and μ = 4 (active phase). In this configuration, all types of pair connections exist and evolve in time. This dynamics further evolves until a finite-size fluctuation takes it to the absorbing state displayed in panel (c), where there are no unsatisfying links. As it is evident from the figure, the system divides in more than two friendly groups with different opinions, a situation that corresponds to μ μ < ≈14 split for these parameters. In panel (b) we show an also active configuration for = .
p 0 7 and μ μ = > 16 split . This configuration evolves until a fluctuation, for large times, takes it to the absorbing one displayed in panel (d). In this absorbing configuration there are only two groups, a situation corresponding to μ μ > split . In all cases displayed, the number of nodes is = N 30. Figure 10. In panel (a) we display a final absorbing configuration of the dynamics in the case of parameters = . p 0 8 and μ = 4 (a point in parameter space corresponding to the frozen phase) that contains more than two groups, a situation corresponding to μ μ < split . Panel (b), with parameters = . p 0 95, μ = 16 (also in the frozen phase), shows a frozen configuration with two groups, as it corresponds to μ μ > split . In both cases, the number of nodes is = N 30. (2019) 9:9726 | https://doi.org/10.1038/s41598-019-45937-y www.nature.com/scientificreports www.nature.com/scientificreports/ in two groups) to a wide distribution in which different realizations reach a state in which two large groups coexist with a varying number of smaller groups. Examples of two groups-splitting can be seen in Figs 9(d) and 10(b), whereas more than two groups are displayed in Figs 9(c) and 10(a). The topological transition appears when crossing the line μ p ( ) split in parameter space μ p ( , ) such that for μ μ > split the system always fall into exactly two friendly group with different opinion, = N 2 g , but for μ μ < split , the system splits in more than two friendly groups, The exact location of the transition line μ p ( ) split depends on the initial condition and on system size N. The simulations indicate that for the initial condition = .
the value μ split is roughly independent on p, see panels (a) and (b) in Fig. 11. Furthermore the transition point μ split grows with system size, see panel (c) of Fig. 11, and we speculate that it tends to infinity with N, in such a way that in the thermodynamic limit, a typical absorbing configuration always contains more than two groups. For a different set of initial conditions, it is not true that the line μ p ( ) split is independent of p, but the same conclusion is reached about the disappearance of the two-groups phase in the large N limit.
The size G of the largest white group is also a random variable described by the corresponding pdf f G ( ) giant .
The mean value and variance, 〈 〉 G , σ G [ ] 2 of that distribution depend, besides p and μ, on the initial densities of white opinions x 0 and friendly links  0 , in a similar functional form that the final density 〈 〉 x st displayed in Fig. 6(c,d). For a given initial condition, 〈 〉 G and σ G [ ] 2 increase linearly with system size N (not shown). In Fig. 12 we show that f G ( ) giant can be well represented by a Gaussian distribution. Figure 11. We plot in different panels the probability density function of the number of groups for fixed p but different μ for the initial condition = . x The distributions with the same color in panel (a and b) correspond to the same μ but different p and are, roughly speaking, indistinguishable, implying that, for these particular initial conditions, the mean value and the variance of f groups # do not depend on p but on μ. In addition we can see that in panels (a and b) with = N 30 splitting occurs for μ μ <  14 split but in panel (c) with = N 250 it occurs for μ μ <  22 split , providing evidence that the critical value for splitting increases with system size. p 0 8 and μ = 6. The results for different system sizes have been rescaled by defining ξ = − 〈 〉 G G ( ) /σ G [ ], and can be clearly fitted by a Gaussian distribution (solid line). In all cases the initial condition is = . x www.nature.com/scientificreports www.nature.com/scientificreports/

summary and Discussion
We have introduced a model of opinion formation in the context of the study of coupled dynamics of node and link states in a complex network: We postulate that friendly/unfriendly links can affect the process of changing opinions so that friends like to have the same opinion and unfriendly relations are satisfied with different opinions. We have proposed a dynamical rule for the evolution of unsatisfied pairwise relations to satisfactory relations by either node or link updates. The relevant parameter of the problem p is the probability for link update instead of a node update in the local dynamics rule. By a mean-field rate equation analysis, corroborated by Monte Carlo simulations, we find an absorbing continuous phase transition from a frozen to a dynamically active phase occurring for a critical value of p. In spite of a dynamical rule of local convergence, global convergence to the satisfactory absorbing phase does not occur for slow link update. In the active phase, the densities of all possible pairwise relations fluctuate around well defined values that depend on μ and p. However, finite-size fluctuations take the system to a frozen configuration, but this occurs in a characteristic time that grows exponentially with system size. In the frozen phase, the system orders dynamically, with an order parameter decaying exponentially to zero. For a finite system, the characteristic time to reach the absorbing state in this phase grows logarithmically with system size. The final frozen configurations reached continuously in the absorbing phase or those reached by a finite-size fluctuation in the active phase show a group structure such that the links within a group are friendly and the links between groups are unfriendly. There is a finite-size topological transition between a two-group and a multi-group structure of those final frozen configurations. An interesting feature of these results is the polarized configuration of the final absorbing states in both the active and absorbing phases (Figs 9 and 10). In these configurations, two large and opponent groups with different opinion emerge out of a disordered configuration, illustrating the process of social polarization [21][22][23][24] in opinion formation. We note that this polarization is in contrast to the consensus configuration observed in the final absorbing state of the coevolving voter model 31 in its active phase. A challenging open question for future research is how these transitions are modified when the network evolves dynamically as in a coevolving voter model 31 : An unsatisfying pair could evolve to a satisfying pair by rewiring a link in the network with some probability.

Method
Rate equation of coupled evolution of node and link in imitating process in the mean-field approximation. To predict the behavior of the different densities of pairs as a function of time, we derive the rate equations of the dynamic sketched in Fig. 1 on a network in which each nodes has exactly μ links as