• Open Access

Binary-State Dynamics on Complex Networks: Pair Approximation and Beyond

James P. Gleeson*

  • MACSI, Department of Mathematics and Statistics, University of Limerick, Limerick, Ireland

  • *james.gleeson@ul.ie

Phys. Rev. X 3, 021004 – Published 29 April, 2013

DOI: https://doi.org/10.1103/PhysRevX.3.021004

Abstract

A wide class of binary-state dynamics on networks—including, for example, the voter model, the Bass diffusion model, and threshold models—can be described in terms of transition rates (spin-flip probabilities) that depend on the number of nearest neighbors in each of the two possible states. High-accuracy approximations for the emergent dynamics of such models on uncorrelated, infinite networks are given by recently developed compartmental models or approximate master equations (AMEs). Pair approximations (PAs) and mean-field theories can be systematically derived from the AME. We show that PA and AME solutions can coincide under certain circumstances, and numerical simulations confirm that PA is highly accurate in these cases. For monotone dynamics (where transitions out of one nodal state are impossible, e.g., susceptible-infected disease spread or Bass diffusion), PA and the AME give identical results for the fraction of nodes in the infected (active) state for all time, provided that the rate of infection depends linearly on the number of infected neighbors. In the more general nonmonotone case, we derive a condition—that proves to be equivalent to a detailed balance condition on the dynamics—for PA and AME solutions to coincide in the limit . This equivalence permits bifurcation analysis, yielding explicit expressions for the critical (ferromagnetic or paramagnetic transition) point of such dynamics, that is closely analogous to the critical temperature of the Ising spin model. Finally, the AME for threshold models of propagation is shown to reduce to just two differential equations and to give excellent agreement with numerical simulations. As part of this work, the Octave or Matlab code for implementing and solving the differential-equation systems is made available for download.

Popular Summary

Article Text

I. INTRODUCTION

Binary-state dynamics on complex networks are frequently used as simple models of social interaction . Each node in a network is considered to be in one of two possible states (e.g., susceptible or infected, inactive or active) at each moment of time. Switching of nodes to the opposite state occurs stochastically, with probabilities that depend on the state of the updating node and on the states of its neighbors. Examples include models for competition between two opinions in a network-based population, such as the voter model and the majority-vote model . Models for rumor propagation, adoption of new technologies, or the spreading of behavior are often based on generalizations of disease-spread or individual-threshold dynamics. Recently, data from online social networks have been mined to help increase the understanding of how individuals are influenced by their network neighbors ; further insight is provided by experiments where the network topology and social interactions are carefully controlled . Binary-state dynamics provide the simplest test for newly developed models ; interest is often focused on the critical parameters where the behavior of the model changes qualitatively, as at the paramagnetic or ferromagnetic phase transition that occurs in the Ising spin model at the critical temperature .

All these models, and many others, can be considered as special cases of a general formulation for stochastic binary-state dynamics. Letting be the degree (number of network neighbors) of a node in an unweighted, undirected, uncorrelated, static network, we denote by the number of neighbors of that are in state 1, i.e., infected or active or spin-up, depending on the model context. If node itself is in state 0 (i.e., susceptible or inactive or spin-down), then the rate is defined by letting be the probability that will switch from state 0 to state 1 in an infinitesimally small time interval . On the other hand, if node is in state 1, then is the probability that node will switch to state 0 within the time interval . We call the functions and the infection rate and recovery rate, respectively; collectively, they are also called the transition rates or spin-flip rates of the dynamics. These rates depend (only) on the degree of node and on the number of neighbors of that are in state 1; we also assume that the rates and do not vary in time. In Sec.  (see Table ), we show that many stochastic binary-state social-interaction models can be expressed in terms of suitable and functions . Deterministic threshold models —where nodes change state (become active) only when the number of active neighbors exceeds a predefined threshold level—may also be put into this form (see Sec. ).

TABLE I.

Infection and recovery rates for some examples of binary-state dynamics on networks: is the node’s degree, is its number of infected neighbors, and is the mean degree of the network. The parameters of the various models are described in Sec. . Symmetric models are those with rates that obey condition ; equilibrium models obey condition .

Process or model Infection rate Recovery rate Symmetric model Equilibrium model
SIS No No
SI 0 No No
Bass 0 No No
Kirman No No
Voter Yes No
Link-update voter Yes No
Nonlinear voter on four-regular random graph Yes No
Language model Yes, if No
Majority vote Yes No
Ising Glauber Yes Yes
Ising Metropolis Yes Yes
Threshold 0 No No

Approximations for macroscopic quantities of interest—such as the expected fraction of active nodes in the network at a given time—are needed in order to understand how the combination of network topology (e.g., the degree distribution ) and the microscopic dynamics (defined by the rates and ) affects the emergent dynamics at the macroscopic level. In the limit of infinite network size, systems of deterministic differential equations may be derived to describe such macroscopic quantities at various levels of approximation. The most common analytical approach for dynamics on complex networks is mean-field (MF) theory. Mean-field theories are typically derived under a number of assumptions, the most important of which for the current discussion is the assumed lack of dynamical correlations . Under this assumption, the neighbors of a node are considered to be active (i.e., in state 1) with a probability that is determined by the overall fraction of active nodes (of same degree) across the network as a whole. In particular, the probability of a neighbor of being active is assumed to be independent of the state of node : Thus, there is assumed to be no correlation between the state of and the state of ’s neighbors.

Examples of MF theories for binary-state dynamics are found in Refs. , and the validity of the assumptions of mean-field theories on networks is considered in detail in Refs. . Mean-field theories, although relatively simple to derive and solve, are often quite inaccurate, especially on sparse networks (with low mean degree ) and near to critical points of the dynamics. Pair-approximation (PA) theories improve upon the accuracy of MF theories—at the cost of extra complexity in the resulting system of differential equations—by including dynamical correlations at a pairwise level . In addition to the fraction of active nodes in each degree class (as in the MF), PA theories account for the fraction of edges in the network that link active nodes to active nodes (called edges here), inactive nodes to inactive nodes ( edges), and inactive to active nodes ( edges). However, dynamical correlations beyond nearest neighbors are not captured by PA theories. If node is inactive, for example, the probability that each of its neighbors is in the active state is approximated in PA theory by the relative preponderance of edges over edges in the whole network. Moreover, each neighbor of is considered to be independent of every other neighbor of ; thus, if is the number of ’s neighbors that are in the active state, is assumed to have a binomial distribution.

Higher-order accuracy, beyond the PA level, has been obtained recently using compartmental models or approximate master equations (AMEs) . These approaches consider the state of nodes and their immediate neighbors, generating large systems of differential equations that better capture the (nonbinomial) distributions of the number of active neighbors of updating nodes. The increased complexity of the differential-equation systems gives improved accuracy over PA and MF, particularly near critical points of the dynamics .

In this paper, we investigate the possibility of, for certain classes of dynamics, reducing the dimensionality of the full AME system to obtain a simpler set of differential equations, but without any loss of accuracy. Ideally, the reduced-dimension system would permit analytical solutions, but, even if closed-form solutions are not possible, the smaller system is less computationally expensive to solve. We restrict our attention to static networks but note that AME approaches have also been successfully used for models where the network structure coevolves with the dynamics . We also assume here that the networks are generated by the configuration model , and so they have zero clustering (transitivity) in the infinite-size limit. These simplifications enable us to focus on the types of dynamics for which the AME can be exactly reduced to a simpler (e.g., pair-approximation) system; we anticipate that complicating factors such as nonzero clustering can be studied in subsequent work, for example, by employing the models of Refs. .

The remainder of the paper is organized as follows. In Sec. , we give examples of binary-state dynamics, and in Sec.  we briefly review the approximate-master-equation approach of Ref. , showing how the equations of PA and MF theories can be derived as systematic approximations of the full AME. In Sec. , the AME and PA solutions are shown to give identical results for certain macroscopic quantities within the class of dynamics we call monotone. In Sec. , a more general class of dynamics (corresponding to equilibrium spin systems that obey a detailed balance) is shown to have AME and PA solutions that are equal in the steady-state limit but not at finite . Focusing on equilibrium models with up-down symmetry, this result enables us to derive—in Sec. —an analytical result for the critical point of such dynamics, i.e., the bifurcation point marking the phase transition from paramagnetic (disordered) to ferromagnetic (ordered) behavior, as in the Ising model on complex networks . Finally, in Sec. , we show that (monotone) threshold models may be accurately approximated using an extended version of the AME, and that the large system of equations may be reduced to just two differential equations to determine the expected fraction of active nodes. Details of some mathematical derivations are contained in the Appendixes.

II. EXAMPLES OF BINARY-STATE DYNAMICS

Here, we briefly examine some of the binary-state models that can be described using transition rates and ; see Table  and the schematic Fig. . Note that in many of the models described here, a single randomly chosen node is updated at each time step, for which the time increment is , where is the number of nodes in the network (and we take the limit). In such a model, is the probability that a spin-down node, when selected for updating, flips to spin-up, while is the corresponding probability of a node flipping from spin-up to spin-down .

FIG. 1

Schematic of the infection rate and recovery rate as defined in the text. Here—and in Figs.  and —white nodes are susceptible or inactive or spin-down and black nodes are infected or active or spin-up. Examples of the dependence of the rates on the number of infected neighbors include (a) the susceptible-infected-susceptible disease-spread model, (b) Glauber dynamics for the Ising model, and (c) a stylized “friend-saturation” model. Case (c) is motivated by the data analysis of social contagion for Digg stories in Fig. 4 of Ref. , where the probability of, for example, becoming infected (voting on a Digg story) initially increases with the number of infected neighbors but then saturates and eventually decreases, giving a function that is roughly symmetric about .

The susceptible-infected-susceptible (SIS) disease-spread model assumes that infected individuals transmit disease to each of their network neighbors at a rate . Thus, if a susceptible node has infected neighbors, the probability of it remaining susceptible for a time interval is . Therefore, its probability of infection is , and in the limit this probability is , giving the rate for SIS dynamics in Table . The recovery rate in the SIS model is a constant , and in the case of , recovery is impossible; this special case is called the SI model.

The Bass model for diffusion of technological innovations is an extension of the SI model and may be similarly specialized to a network context. In the Bass model, nodes move from the susceptible (nonadopted) state to the infected (adopted) state with the rate . The parameter represents individual action, independent of the influence of neighbors, while the parameter quantifies the “herding behavior,” whereby individuals copy the actions of their neighbors: In this case, they do so by adopting the innovation when of their neighbors have already done so. Transitions from the adopted state to the nonadopted state are not permitted, so . Kirman’s ant colony model for choice between stock market trading strategies has similar herding effects but with switching permitted between both states (strategies). A node with of its neighbors infected has neighbors in the susceptible state, so the rate for transitions of such a node from infected to susceptible is . Note that mean-field theory [see Eq.  below] for the complete-graph case —where all nodes are connected to all other nodes—gives the standard population-level version of the Bass and Kirman models .

In the voter model , a node is chosen at random in each time step ( ) and it adopts the same state of one of its (randomly chosen) neighbors. If of the node’s neighbors are in the active state, the node thus becomes active with probability and becomes inactive with probability , these probabilities being the fractions of its neighbors in the respective states. Many variants of the voter model have been studied . In the link-update voter model, for example, an edge (rather than a node) is chosen at random in each time step, and one of the end nodes of the chosen edge then copies the state of the other. Nonlinear voter models and other variants have been examined by a number of authors . In Ref. , for example, the current state of the node is incorporated into the probabilities for switching state, lending an effective inertia to the dynamics: See the transition rates for this case (on four-regular random graphs) in Table . The language model of Ref.  has been examined in detail in Ref. ; in this model, the two states represent the primary language choice of a person (node), and the probability of switching states is proportional to the fraction of speakers in the locality, raised to the power , multiplied by the status parameter or of the respective language. ( for the symmetric case of equal-status languages; then gives the voter model.)

Many models of opinion dynamics are based on the classical Ising model of magnetic spins. Here, each node is in either the spin-up or spin-down state, and transitions occur according to dynamical rules that minimize the Hamiltonian of ferromagnetic interactions . Letting represent the temperature of the heat bath and the ferromagnetic-interaction parameter, the transition rates for Glauber and Metropolis dynamics are shown in Table . The majority-vote model is a nonequilibrium spin model, with spins tending to align with the local neighborhood majority but with a noise parameter giving the probability of misalignment.

As a final example, we consider threshold models , which are used to model the diffusion of fads, collective action, or adoption of innovations . Each node has a (frozen) threshold level, which may depend on the degree of the node. In each time step, a fraction of the nodes in the network is updated. When updated, a node compares the number of its active neighbors to its threshold and activates (with probability 1) if is greater than or equal to its threshold. Similarly to the Bass model of innovation diffusion, no recovery is usually permitted in this model (but see Ref.  for a recent extension that includes recovery). Note that a coordination game, modeling the diffusion of new behavior through a network, can also be written as a threshold model of this type . Threshold models, and the extension of the AME required to describe them, are considered below in Sec. .

III. APPROXIMATE MASTER EQUATIONS, PAIR APPROXIMATIONS, AND MEAN-FIELD THEORY

For completeness, we here briefly recapitulate the AMEs introduced in Ref. . These equations were derived by generalizing the approach used in Refs.  for SIS dynamics; see Appendix  for details. Let [respectively, ] be the fraction of -degree nodes that is susceptible (respectively, infected) at time and have infected neighbors. Then, the fraction of -degree nodes that is infected at time is given by 𝜌𝑘(𝑡) =𝑘𝑚=0𝑖𝑘,𝑚 =1 𝑘𝑚=0𝑠𝑘,𝑚, and the fraction of infected nodes in the whole network is found by summing over all 𝑘 classes: 𝜌(𝑡) =𝜌𝑘(𝑡) 𝑘𝑃𝑘𝜌𝑘(𝑡). (Recall that 𝑃𝑘 is the degree distribution, i.e., the probability that a randomly chosen node has 𝑘 neighbors.)

The approximate master equations for the evolution of 𝑠𝑘,𝑚(𝑡) and 𝑖𝑘,𝑚(𝑡) are

𝑑𝑑𝑡𝑠𝑘,𝑚=𝐹𝑘,𝑚𝑠𝑘,𝑚+𝑅𝑘,𝑚𝑖𝑘,𝑚𝛽𝑠(𝑘𝑚)𝑠𝑘,𝑚+𝛽𝑠(𝑘𝑚+1)𝑠𝑘,𝑚1𝛾𝑠𝑚𝑠𝑘,𝑚+𝛾𝑠(𝑚+1)𝑠𝑘,𝑚+1,
(1)
𝑑𝑑𝑡𝑖𝑘,𝑚=𝑅𝑘,𝑚𝑖𝑘,𝑚+𝐹𝑘,𝑚𝑠𝑘,𝑚𝛽𝑖(𝑘𝑚)𝑖𝑘,𝑚+𝛽𝑖(𝑘𝑚+1)𝑖𝑘,𝑚1𝛾𝑖𝑚𝑖𝑘,𝑚+𝛾𝑖(𝑚+1)𝑖𝑘,𝑚+1
(2)

for each 𝑚 in the range 0,,𝑘 and for each 𝑘 class in the network, with 𝛽𝑠 =𝑘𝑚=0(𝑘 𝑚)𝐹𝑘,𝑚𝑠𝑘,𝑚/ 𝑘𝑚=0(𝑘 𝑚)𝑠𝑘,𝑚, 𝛾𝑠 =𝑘𝑚=0(𝑘 𝑚)𝑅𝑘,𝑚𝑖𝑘,𝑚/ 𝑘𝑚=0(𝑘 𝑚)𝑖𝑘,𝑚, 𝛽𝑖 =𝑘𝑚=0𝑚𝐹𝑘,𝑚𝑠𝑘,𝑚/𝑘𝑚=0𝑚𝑠𝑘,𝑚, and 𝛾𝑖 =𝑘𝑚=0𝑚𝑅𝑘,𝑚𝑖𝑘,𝑚/𝑘𝑚=0𝑚𝑖𝑘,𝑚. Equations  and , with the time-dependent rates 𝛽𝑠, 𝛾𝑠, 𝛽𝑖, and 𝛾𝑖 (defined as nonlinear functions of 𝑠𝑘,𝑚 and 𝑖𝑘,𝑚), form a closed system of deterministic equations that can be solved numerically using standard methods; Octave or Matlab scripts are available from the author’s Web site . Assuming that a randomly chosen fraction 𝜌(0) of nodes is initially infected, the initial conditions are 𝑠𝑘,𝑚(0) =[1 𝜌(0)]𝐵𝑘,𝑚[𝜌(0)] and 𝑖𝑘,𝑚(0) =𝜌(0)𝐵𝑘,𝑚[𝜌(0)], where 𝐵𝑘,𝑚(𝑞) denotes the binomial factor

𝐵𝑘,𝑚(𝑞)=(𝑘𝑚)𝑞𝑚(1𝑞)𝑘𝑚.
(3)

For dynamics on a general network, with nonempty degree classes from 𝑘 =0 up to a cutoff 𝑘max, the number of differential equations in the system and is (𝑘max +1)(𝑘max +2), and so grows with the square of the largest degree. Some approximation is therefore necessary if it is desirable to reduce the AME to a lower-dimensional system. One possibility is to consider the parameters 𝑝𝑘(𝑡) [respectively, 𝑞𝑘(𝑡)], defined as the probability that a randomly chosen neighbor of a susceptible (respectively, infected) 𝑘-degree node is infected at time 𝑡. Noting that 𝑝𝑘(𝑡) can be expressed in terms of 𝑠𝑘,𝑚 as 𝑘𝑚=0𝑚𝑠𝑘,𝑚/𝑘𝑚=0𝑘𝑠𝑘,𝑚, an evolution equation for 𝑝𝑘 may be derived by multiplying Eq.  by 𝑚 and summing over 𝑚. The right-hand side of the resulting equation contains higher moments of 𝑠𝑘,𝑚, so a closure approximation is needed to proceed. If we make the ansatz that 𝑠𝑘,𝑚 and 𝑖𝑘,𝑚 are proportional to binomial distributions 𝑠𝑘,𝑚 (1 𝜌𝑘)𝐵𝑘,𝑚(𝑝𝑘) and 𝑖𝑘,𝑚 𝜌𝑘𝐵𝑘,𝑚(𝑞𝑘), we obtain the PA that consists of the 3𝑘max +1 differential equations

𝑑𝑑𝑡𝜌𝑘=𝜌𝑘𝑘𝑚=0𝑅𝑘,𝑚𝐵𝑘,𝑚(𝑞𝑘)+(1𝜌𝑘)𝑘𝑚=0𝐹𝑘,𝑚𝐵𝑘,𝑚(𝑝𝑘),𝑑𝑑𝑡𝑝𝑘=𝑘𝑚=0[𝑝𝑘𝑚𝑘][𝐹𝑘,𝑚𝐵𝑘,𝑚(𝑝𝑘)𝜌𝑘1𝜌𝑘𝑅𝑘,𝑚𝐵𝑘,𝑚(𝑞𝑘)]+𝛽𝑠(1𝑝𝑘)𝛾𝑠𝑝𝑘,𝑑𝑑𝑡𝑞𝑘=𝑘𝑚=0[𝑞𝑘𝑚𝑘][𝑅𝑘,𝑚𝐵𝑘,𝑚(𝑞𝑘)1𝜌𝑘𝜌𝑘𝐹𝑘,𝑚𝐵𝑘,𝑚(𝑝𝑘)]+𝛽𝑖(1𝑞𝑘)𝛾𝑖𝑞𝑘
(4)

for each 𝑘 class. The rates here are given by inserting the binomial ansatz into the general formulas, so that 𝛽𝑠, for example, is (1 𝜌𝑘)𝑚(𝑘 𝑚)𝐹𝑘,𝑚𝐵𝑘,𝑚(𝑝𝑘)/ (1 𝜌𝑘)𝑘(1 𝑝𝑘); the initial conditions are 𝜌𝑘(0) =𝑝𝑘(0) =𝑞𝑘(0) =𝜌(0).

A cruder, MF approximation results from replacing both 𝑝𝑘 and 𝑞𝑘 with 𝜔: 𝑠𝑘,𝑚 (1 𝜌𝑘)𝐵𝑘,𝑚(𝜔) and 𝑖𝑘,𝑚 𝜌𝑘𝐵𝑘,𝑚(𝜔), where 𝜔 =𝑘𝑧𝜌𝑘 is the probability that one end of a randomly chosen edge is infected, and 𝑧 =𝑘 is the mean degree of the network. Using this ansatz in the master equations yields a closed system of 𝑘max +1 differential equations for the fraction 𝜌𝑘 of infected 𝑘-degree nodes:

𝑑𝑑𝑡𝜌𝑘=𝜌𝑘𝑘𝑚=0𝑅𝑘,𝑚𝐵𝑘,𝑚(𝜔)+(1𝜌𝑘)𝑘𝑚=0𝐹𝑘,𝑚𝐵𝑘,𝑚(𝜔),
(5)

with 𝜌𝑘(0) =𝜌(0).

In Ref. , we showed that the AME system and generally gives improved accuracy over the approximations and . Moreover, the general Eqs.  and for pair approximation and mean-field theory reduce to previously studied equations when the specific rates 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚 are selected from Table . In this paper, we concentrate on finding dynamics for which the AME system can be reduced exactly to lower-dimensional systems, for instance, to the PA equations.

IV. MONOTONE DYNAMICS

We first consider the case where 𝑅𝑘,𝑚 =0 for all 𝑘 and 𝑚, with 𝐹𝑘,𝑚 nonzero for at least some arguments. The zero-recovery rate means that it is impossible for nodes to switch from the infected state to the susceptible state, and so we call this type of dynamics monotone . In this case, the 𝑠𝑘,𝑚 equations of the AME system are decoupled from the 𝑖𝑘,𝑚 equations, with Eq.  reducing to

𝑑𝑑𝑡𝑠𝑘,𝑚=𝐹𝑘,𝑚𝑠𝑘,𝑚𝛽𝑠(𝑘𝑚)𝑠𝑘,𝑚+𝛽𝑠(𝑘𝑚+1)𝑠𝑘,𝑚1,
(6)

and the fraction of infected nodes is given by 𝜌(𝑡) =1 𝑘𝑃𝑘𝑘𝑚=0𝑠𝑘,𝑚(𝑡). Thus, several important macroscopic quantities (although not all) can be found by solving Eq.  for 𝑠𝑘,𝑚(𝑡), without considering the 𝑖𝑘,𝑚(𝑡) variables.

Now, consider whether or not the pair-approximation ansatz 𝑠𝑘,𝑚(𝑡) =[1 𝜌𝑘(𝑡)]𝐵𝑘,𝑚[𝑝𝑘(𝑡)]—that was used to derive Eqs. —could be an exact solution of Eq. . Inserting the PA ansatz into Eq.  and dividing by (1 𝜌𝑘)𝐵𝑘,𝑚(𝑝𝑘) gives the condition

11𝜌𝑘𝑑𝜌𝑘𝑑𝑡+(𝑚𝑝𝑘𝑘𝑚1𝑝𝑘)𝑑𝑝𝑘𝑑𝑡=𝐹𝑘,𝑚𝛽𝑠(𝑘𝑚)+𝛽𝑠1𝑝𝑘𝑝𝑘𝑚,
(7)

where the identities

𝑑𝑑𝑝𝐵𝑘,𝑚(𝑝)=(𝑚𝑝𝑘𝑚1𝑝)𝐵𝑘,𝑚(𝑝)
(8)

and

(𝑘𝑚+1)𝐵𝑘,𝑚(𝑝)=1𝑝𝑝𝑚𝐵𝑘,𝑚(𝑝)
(9)

have been, respectively, utilized on the left-hand side and the right-hand side of Eq. .

Equation  can be viewed as a condition on the forms of the infection rate 𝐹𝑘,𝑚 for which the PA ansatz on 𝑠𝑘,𝑚 is an exact solution of the corresponding approximate master equation. Note that all terms in Eq.  are linear in 𝑚, so 𝐹𝑘,𝑚 is necessarily of the form

𝐹𝑘,𝑚=𝑐𝑘+𝑑𝑘𝑚
(10)

for (possibly 𝑘-dependent) constants 𝑐𝑘 and 𝑑𝑘. Using this form in confirms that the solutions of the AME for 𝑠𝑘,𝑚 are given by the PA ansatz, with 𝜌𝑘 and 𝑝𝑘 being the solutions of the reduced system

𝑑𝑑𝑡𝜌𝑘=(1𝜌𝑘)(𝑐𝑘+𝑘𝑝𝑘𝑑𝑘),𝑑𝑑𝑡𝑝𝑘=𝑑𝑘𝑝𝑘(1𝑝𝑘)+‾‾‾𝛽𝑠(1𝑝𝑘),
(11)

where ‾‾‾𝛽𝑠 ={𝑘𝑃𝑘(1 𝜌𝑘)𝑘(1 𝑝𝑘)[𝑐𝑘 +(𝑘 1)𝑝𝑘𝑑𝑘]}/ [𝑘𝑃𝑘(1 𝜌𝑘)𝑘(1 𝑝𝑘)]. For the special case of 𝑧-regular random graphs or Bethe lattices (i.e., 𝑃𝑘 =𝛿𝑘,𝑧 for integer 𝑧), system can be solved analytically to give the explicit solution for the infected fraction as

𝜌(𝑡)=1[1𝜌(0)]𝑒(𝑐+𝑧𝑑)𝑡[1[1𝜌(0)]𝑑(𝑧2)𝑐+𝑑(𝑧2)×(1𝑒[𝑐+𝑑(𝑧2)]𝑡)]𝑧/(𝑧2).
(12)

In Fig. , we show results for SI dynamics (see Table ), which are monotone, and of form with 𝑐𝑘 =0 and 𝑑𝑘 =𝜆 for all 𝑘. Note the exact match between the AME and PA solutions for 𝜌(𝑡) and for the fraction of 𝑆-𝐼 edges [see Eqs.  and ] in the network and the excellent agreement with numerical simulation on networks of size 𝑁 =105 .

FIG. 2

SI dynamics (with transmission rate 𝜆 =1) on a truncated scale-free network, with degree distribution given by 𝑃𝑘 𝑘2.5 for degrees in the range 3 𝑘 20 and with 𝑃𝑘 =0 otherwise. The initial infected fraction is 𝜌(0) =102. Here, and in all subsequent figures, the symbols show the results of numerical simulations on networks of size 𝑁 =105, using the time step 𝑑𝑡 =104. The results are averages over 24 realizations; the error bars indicate the mean ±1 standard deviation (but note that the error bars are smaller than the symbols in this example).

The parameters 𝑑𝑘 in may be negative (provided that all 𝐹𝑘,𝑚 values are non-negative). As an example, consider the Bass diffusion model on a 𝑧-regular random graph with 𝑐 =1 and 𝑑 =1/𝑧. This rate might model, for example, the adoption by indie music fans of a new band, where the attractiveness of the band decreases as the number of neighbors 𝑚 who have already “jumped on the bandwagon” increases; see Ref.  for another approach to this aspect of social contagion modeling, which they call “limited imitation contagion.” The expected fraction of fans 𝜌(𝑡) is given explicitly by Eq. ; note that for these parameters the entire network does not become infected; see Fig. .

FIG. 3

Bass diffusion dynamics on a four-regular random graph (or Bethe lattice, 𝑃𝑘 =𝛿𝑘,4), with parameters 𝑐 =1 and 𝑑 =1/4 and initial condition 𝜌(0) =102. The AME and PA solutions for the fraction of adopted nodes are both given by the explicit formula .

It is interesting to note that, while, in the monotone case , the AME solutions for 𝑠𝑘,𝑚(𝑡) are exactly reproduced by the PA equations, the corresponding 𝑖𝑘,𝑚(𝑡) variables are not necessarily equal to their respective PA ansatz 𝜌𝑘(𝑡)𝐵𝑘,𝑚[𝑞𝑘(𝑡)]. Focusing on the SI model, the effect of the nonbinomial 𝑚 distribution in 𝑖𝑘,𝑚 is not visible in Fig. , as the quantities shown there can be written in terms of 𝑠𝑘,𝑚 only; indeed, it is necessary to examine connected triples of nodes to demonstrate this effect . Figure  shows the fraction of connected triples (defined by choosing a node at random, then randomly choosing two of its neighbors) that are of type 𝑆-𝐼-𝑆, the chosen (middle) node being infected while both chosen neighbors are susceptible. In the AME formulation, this fraction is given by 2𝑘𝑃𝑘𝑘𝑚=0(𝑘 𝑚)(𝑘 𝑚 1)𝑖𝑘,𝑚/𝑘(𝑘 1), whereas the corresponding PA version is 2𝑘𝑃𝑘𝜌𝑘𝑘(𝑘 1)(1𝑞𝑘)2/𝑘(𝑘 1). The differences seen in Fig. —and the close fit of AME results to numerical simulations—indicate that the match between AME and PA results in Figs.  and is due to the binomial 𝑚 distribution 𝑠𝑘,𝑚 for susceptible nodes, but the corresponding distribution 𝑖𝑘,𝑚 for infected nodes is not binomial in 𝑚.

FIG. 4

Fraction of triples that are of the 𝑆-𝐼-𝑆 type, for SI dynamics on a three-regular random graph.

V. GENERAL DYNAMICS AND EQUILIBRIUM SPIN MODELS

We consider in this section the more general case of nonmonotone dynamics, where none of the rates 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚 are identically zero. Inserting the PA ansatz for 𝑠𝑘,𝑚 and 𝑖𝑘,𝑚 into the AME Eqs.  and , we find, after dividing by (1 𝜌𝑘)𝐵𝑘,𝑚(𝑝𝑘), the condition

𝐹𝑘,𝑚+𝑅𝑘,𝑚𝜌𝑘1𝜌𝑘𝐵𝑘,𝑚(𝑞𝑘)𝐵𝑘,𝑚(𝑝𝑘)=𝑐(1)𝑘+𝑐(2)𝑘𝑚,
(13)

where 𝑐(1)𝑘 and 𝑐(2)𝑘 represent combinations of terms that depend on 𝑡 and 𝑘 but are independent of 𝑚. Similarly, using the PA ansatz in the AME equation for 𝑖𝑘,𝑚 yields

𝐹𝑘,𝑚1𝜌𝑘𝜌𝑘𝐵𝑘,𝑚(𝑝𝑘)𝐵𝑘,𝑚(𝑞𝑘)𝑅𝑘,𝑚=𝑐(3)𝑘+𝑐(4)𝑘𝑚,
(14)

where, like Eq. , the left-hand side is an exponential function of 𝑚 [because of the presence of the 𝐵𝑘,𝑚 functions; see Eq. ], while the right-hand side is linear in 𝑚. Bearing in mind that the transition rates 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚 are time independent, it is not generally possible to simultaneously solve Eqs.  and to obtain constant transition rates. We conclude that in this general case, AME solutions are not exactly equal, for all times, to the corresponding PA solutions.

However, an important special case occurs if we restrict our attention to the steady-state limit 𝑡 . The AME solutions and PA solutions can be identical in the limit 𝑡 (despite being different at finite 𝑡) if Eqs.  and are simultaneously satisfied in the steady state. The 𝑚 dependence of the left- and right-hand sides then requires that

𝐹𝑘,𝑚𝑅𝑘,𝑚=‾‾‾𝜌𝑘1‾‾‾𝜌𝑘𝐵𝑘,𝑚(‾‾‾𝑞𝑘)𝐵𝑘,𝑚(‾‾‾𝑝𝑘),
(15)

and 𝑐(1)𝑘 =𝑐(2)𝑘 =𝑐(3)𝑘 =𝑐(4)𝑘 =0, where the overbar denotes the steady-state limit of the corresponding variable, e.g., ‾‾‾𝑝𝑘 =lim𝑡𝑝𝑘(𝑡).

The conditions 𝑐(1)𝑘 =𝑐(2)𝑘 =0 can be shown to imply that ‾‾‾𝑝𝑘/(1 ‾‾‾𝑝𝑘) =‾‾‾𝛽𝑠/‾‾𝛾𝑠, and the 𝑘 independence of the right-hand side implies that ‾‾‾𝑝𝑘 must, in fact, be independent of 𝑘:

‾‾‾𝑝𝑘=‾‾‾𝑝for all𝑘.
(16)

Similarly, an analysis of the conditions 𝑐(3)𝑘 =𝑐(4)𝑘 =0 shows that the steady-state values of 𝑞𝑘 must also be independent of 𝑘; i.e., ‾‾‾𝑞𝑘 =‾‾‾𝑞 for all 𝑘.

Replacing ‾‾‾𝑝𝑘 and ‾‾‾𝑞𝑘 in Eq.  with these simplifications, we can rewrite the condition on the rates as

𝐹𝑘,𝑚𝑅𝑘,𝑚=𝑏𝑘𝑎𝑚,
(17)

where

𝑏𝑘=‾‾‾𝜌𝑘1‾‾‾𝜌𝑘(1‾‾‾𝑞1‾‾‾𝑝)𝑘
(18)

and

𝑎=‾‾‾𝑞(1‾‾‾𝑝)‾‾‾𝑝(1‾‾‾𝑞).
(19)

It is demonstrated in Appendix  that Eq.  is precisely the condition for microscopic reversibility (detailed balance) of the stochastic dynamics, i.e., for the dynamics to correspond to an equilibrium spin-flip model. Moreover, we show in Appendix  that the steady solution of the AME for transition rates obeying condition can be fully specified in a simple manner, by the equation

‾‾‾𝜌𝑘=𝑏𝑘𝑏𝑘+(1‾‾‾𝑝+‾‾‾𝑝𝑎)𝑘,
(20)

where ‾‾‾𝑝 is a solution of the algebraic equation

‾‾‾𝑝(1‾‾‾𝑝+‾‾‾𝑝𝑎)1‾‾‾𝑝2+‾‾‾𝑝2𝑎=𝑘𝑘𝑧𝑃𝑘𝑏𝑘𝑏𝑘+(1‾‾‾𝑝+‾‾‾𝑝𝑎)𝑘.
(21)

As an example of dynamics obeying condition , we introduce a modified model of the susceptible-infected-susceptible type, with

𝐹𝑘,𝑚=𝜆𝑚and𝑅𝑘,𝑚=𝜇
(22)

for positive constants 𝜆 and 𝜇. Like the standard SIS disease-spread model (see Table ), recovery of infected nodes occurs at a constant rate 𝜇. However, in contrast to the linear-in- 𝑚 dependence of the SIS infection rate 𝐹𝑘,𝑚, here the probability of infection is assumed to grow exponentially with 𝑚. While admittedly artificial, such dynamics might find an application in fitting data from social contagion experiments on networks, where the dependence of the rates of 𝑚 is not yet clear . Figure  shows that the AME solutions and PA solutions indeed agree as 𝑡 , despite being different for finite 𝑡.

FIG. 5

Modified SIS-type model described by the rates on a three-regular random graph. The parameters are 𝜆 =6 and 𝜇 =3; the initial state has 𝜌(0) =102.

Other important examples of equilibrium models obeying are given by the Glauber and Metropolis dynamics for the Ising spin model (see Table  and Fig. ), which both have 𝑎 =𝑒4𝐽/𝑇 and 𝑏𝑘 =𝑒(2𝐽/𝑇)𝑘 in Eq. . These models are examples of spin models with up-down symmetry and are considered in greater detail below in Sec. .

FIG. 6

Glauber dynamics for the Ising spin model dynamics on a Poisson (Erdös-Rényi) random graph with mean degree 𝑧 =7. The interaction parameter 𝐽 is set to 1, the temperature 𝑇 is 2/log(2.5) 2.18, and the initial fraction of spin-ups is 0.51.

VI. MODELS WITH UP-DOWN SYMMETRY

Models with up-down symmetry have dynamics that are invariant to the swap of state labels (susceptible to infected, and vice versa) for all nodes. In Refs. , this symmetry is called “ 𝑍2 symmetry.” It is characteristic of the voter model and other opinion dynamics models; see the fourth column of Table . In terms of the transition rates, the symmetry condition implies that

𝑅𝑘,𝑚=𝐹𝑘,𝑘𝑚for𝑚=0,,𝑘and for all𝑘.
(23)

Note that for such dynamics, the AME system and is invariant under the change of variables 𝑠𝑘,𝑚 𝑖𝑘,𝑘𝑚 and 𝑖𝑘,𝑚 𝑠𝑘,𝑘𝑚. Since the expected fraction of degree- 𝑘 nodes can be written as 𝜌𝑘 =𝑘𝑚=0𝑖𝑘,𝑚 =1 𝑘𝑚=0𝑠𝑘,𝑚, this symmetry condition implies that a solution exists of the AME with 𝜌𝑘(𝑡) =1/2 for all 𝑡. However, this solution may not be the only one possible: Depending on the initial condition 𝜌(0) and on the parameter regime, other solutions of the AME may also be found.

Focusing first on equilibrium spin models with up-down symmetry, which obey condition in addition to , we investigate the stability of the symmetric solution with ‾‾‾𝜌𝑘 =1/2. First, note that there is only a single parameter in these models: Putting 𝑚 =𝑘/2 for even 𝑘, or, for odd 𝑘, 𝑚 =(𝑘 1)/2 and then 𝑚 =(𝑘 +1)/2 into Eq.  and imposing condition immediately yields the necessary relation

𝑏𝑘=𝑎𝑘/2
(24)

between the parameters of equilibrium spin models. Using the steady-state solution of Sec. , it is possible to show that a critical value 𝑎𝑐 of the parameter 𝑎 exists: In the language of dynamical systems, this value is a (pitchfork) bifurcation point . For parameter values 𝑎 with 𝑎 <𝑎𝑐, the symmetric solution 𝜌 =1/2 is stable, meaning that if 𝜌(0) is close to 1/2, the steady-state solution will be ‾‾‾𝜌 =1/2. In spin models (where the magnetization can be written as 𝑀 =|2𝜌 1|), this regime is the paramagnetic (disordered) phase; for opinion models in this regime, the two opinions coexist equally on the network. However, if the parameter 𝑎 exceeds the critical value 𝑎𝑐, then the symmetric solution 𝜌 =1/2 is unstable, and two other stable solutions, symmetric about 𝜌 =1/2, exist. This regime is the ferromagnetic (ordered) phase; for opinion models, one of the two opinions dominates the other. The critical value 𝑎𝑐 gives the phase transition point, and the behavior of ‾‾‾𝜌 near 𝑎𝑐 can be determined from Eq.  (see Appendix  for details) in a very similar fashion to the analysis of the Ising model in Refs. ; see also Ref. . The results of such an analysis (see Appendix ) may be summarized as follows. If the degree distribution 𝑃𝑘 possesses a finite fourth moment 𝑘4 =𝑘𝑘4𝑃𝑘, then the phase transition is of the mean-field type, with a critical parameter

𝑎𝑐=(𝑘2𝑘22𝑘)2
(25)

and with ‾‾‾𝜌 1/2 ±(𝑎𝑎𝑐)1/2 as 𝑎 𝑎𝑐 from above. Following Refs.  (see Appendix ), if the network has a scale-free degree distribution 𝑃𝑘 𝑘𝛾 as 𝑘 , then, for exponents 𝛾 in the range 2 <𝛾 <3, the critical point is 𝑎𝑐 =1, with ‾‾‾𝜌 1/2 ±(𝑎𝑎𝑐)1/(3𝛾), while, for exponents with 3 <𝛾 <5, we have 𝑎𝑐 given again by Eq. , but with near-criticality scaling of ‾‾‾𝜌 1/2 ±(𝑎𝑎𝑐)1/(𝛾3) as 𝑎 𝑎𝑐. The case 𝛾 =3 shows an infinite order transition at 𝑎 =1, as is discussed in Ref. .

As is mentioned at the end of Sec. , the Ising model (Glauber or Metropolis dynamics) is of type , with temperature 𝑇 related to the parameter 𝑎 via 𝑇 =4𝐽/ln𝑎, so 𝑎 =1 corresponds to infinite temperature. A social-influence model of this type might be given, for example, by transition rates with the properties

𝑅𝑘,𝑚=𝐹𝑘,𝑘𝑚=𝐹𝑘,𝑚for𝑚=0,,𝑘and all𝑘.
(26)

Here, all 𝐹𝑘,𝑚 values for 𝑚 =0 to 𝑘2 are free parameters; for example, the rates 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚 could be given by Fig. , which is motivated by the data analysis in Fig. 4 of Ref. . In any model satisfying , the parameter 𝑎 of is equal to 1, and by the results above the model is in the paramagnetic (disordered) phase for all network topologies. However, Eq.  shows that 𝑎 =1 is near the critical point of the system if the degrees in the network are very heterogeneous (so that 𝑘2 𝑘), and indeed the model is poised precisely at criticality ( 𝑎𝑐 =1) on scale-free networks with infinite-variance degree distributions.

For symmetric models that do not obey the equilibrium condition , the PA and AME solutions are different, even as 𝑡 ; see Fig.  for an example using majority-vote dynamics. The solutions of the PA equation are therefore of limited usefulness, but it is nevertheless interesting that for 𝑧-regular graphs an expression for the critical point can be explicitly obtained. In Appendix , we show that the PA critical point of symmetric models on such networks occurs precisely when the steady-state PA parameter 𝑝 has the value

‾‾‾𝑝𝑐=𝑧22𝑧2,
(27)

and the location of the critical point in parameter space is given by the solutions 𝐹𝑧,𝑚 (for 𝑚 =0,,𝑧) of the implicit equation

𝑧𝑚=0(12𝑚𝑧)𝐹𝑧,𝑚𝐵𝑧,𝑚(𝑧22𝑧2)=0.
(28)

Using the infection rate 𝐹𝑘,𝑚 for the majority-vote model (see Table ), for example, in Eq. , gives an explicit expression for the PA critical noise parameter 𝑄:

𝑄𝑐=[1𝑧12𝑚=0(12𝑚𝑧)𝐵𝑧,𝑚(𝑧22𝑧2)𝑧𝑚=𝑧+12(12𝑚𝑧)𝐵𝑧,𝑚(𝑧22𝑧2)]1.
(29)

Equation  can similarly be used to obtain analytical expressions for the PA critical points in other models on 𝑧-regular random graphs, e.g., Fig. 12 of Ref. , Fig. 1 of Ref. , or Fig. 1 of Ref. .

FIG. 7

Majority-vote model on a three-regular random graph, with 𝜌(0) =0.55 and noise parameter 𝑄 =0.07. This nonequilibrium model does not obey condition , and the AME and PA solutions are not equal, even as 𝑡 , in contrast to Figs.  and .

VII. THRESHOLD DYNAMICS

Threshold models are often used to model the propagation of fads or collective action through a population . Each node has a (frozen) threshold level; these thresholds may be chosen at random (e.g., from a Gaussian distribution) or assigned in some other way (perhaps depending on the degree of the node). In the asynchronous-updating version of these models, a fraction 𝑑𝑡 of nodes is randomly chosen for updating in each time step . When chosen for updating, an inactive node becomes active (with probability 1) if 𝑚, the number of its neighbors that are active, exceeds the node’s threshold . Once activated, nodes cannot subsequently become inactive, so the dynamics are monotone; cf. Sec. . Unlike in Eq.  of Sec. , however, the transition rate is not linear in 𝑚; in fact, it is given by

𝐹𝐤,𝑚={0if𝑚<𝑀𝐤1if𝑚𝑀𝐤
(30)

to reflect the deterministic activation of a node (once it is chosen for updating) when 𝑚 exceeds the threshold level 𝑀𝐤. We have introduced here the vector 𝐤 to encode two properties that define a class of node: their degree 𝑘 (a scalar) and their type 𝑟, which together determine the threshold 𝑀𝐤 for such nodes. The types are assumed to be from a discrete set of possibilities: All nodes of type 𝑟 =1, for example, might have the same threshold 𝑀1, with all nodes of type 𝑟 =2 having a common threshold 𝑀2, with 𝑀2 𝑀1. In this way, the set of all nodes may be partitioned into disjoint sets that are labeled by their degree and their type; in mathematical notation, we combine these labels into a two-vector, defining 𝐤 ={𝑘,𝑟} for the 𝐤 class of nodes. All nodes in the same 𝐤 class have the same degree and the same type, and therefore all share the same threshold 𝑀𝐤. We generalize the degree distribution 𝑃𝑘 to the distribution 𝑃𝐤, which gives the probability that a randomly chosen node has vector 𝐤 (i.e., has degree 𝑘 and type 𝑟). For example, if the thresholds of the nodes are randomly chosen, independent of their degrees, then the 𝑃𝐤 distribution can be written as 𝑃𝐤 =𝑃𝑘𝑃𝑟, where 𝑃𝑘 is the degree distribution and 𝑃𝑟 is the probability that a node is of type 𝑟. By taking the discrete set of types to be sufficiently large, it is possible to approximate a continuous distribution of types or thresholds with a desired level of accuracy. With this extended notation, the AME approach can be generalized in an obvious manner, essentially replacing the scalar degree 𝑘 by the vector 𝐤 as appropriate in Eq. , to yield

𝑑𝑑𝑡𝑠𝐤,𝑚=𝐹𝐤,𝑚𝑠𝐤,𝑚𝛽𝑠(𝑘𝑚)𝑠𝐤,𝑚+𝛽𝑠(𝑘𝑚+1)𝑠𝐤,𝑚1for𝑚=0,,𝑘,
(31)

with the rate 𝛽𝑠 given by 𝛽𝑠 = [𝐤𝑃𝐤𝑘𝑚=0(𝑘 𝑚)𝐹𝐤,𝑚𝑠𝐤,𝑚]/[𝐤𝑃𝐤𝑘𝑚=0(𝑘 𝑚)𝑠𝐤,𝑚]. Note that the sums here are over all 𝐤 classes, i.e., over all degrees 𝑘 and all types 𝑟: 𝐤 𝑘𝑟.

In Ref.  (see also Ref. ), it is argued that for no-recovery threshold models of the type and described above, the fraction 𝜌(𝑡) of active nodes at times 𝑡 can be found by solving just two differential equations:

𝑑𝑑𝑡𝜌=(𝜙)𝜌,𝑑𝑑𝑡𝜙=𝑔(𝜙)𝜙,with𝜙(0)=𝜌(0)=𝐤𝑃𝐤𝜌𝐤(0),
(32)

where

(𝜙)=𝐤𝑃𝐤[𝜌𝐤(0)+[1𝜌𝐤(0)]𝑚𝑀𝐤𝐵𝑘,𝑚(𝜙)]
(33)

and

𝑔(𝜙)=𝐤𝑘𝑧𝑃𝐤[𝜌𝐤(0)+[1𝜌𝐤(0)]𝑚𝑀𝐤𝐵𝑘1,𝑚(𝜙)].
(34)

Here, 𝜌𝐤(0) is the fraction of nodes with vector 𝐤 that are activated (infected) at time 𝑡 =0; as in Ref. , we generalize the usual infected fraction 𝜌(0)—that is used elsewhere in this paper—to allow for possible dependence on the degree or type of the nodes chosen to “seed” the contagion.

In Appendix , we demonstrate that Eqs.  and reduce to through an exact solution of given by

𝑠𝐤,𝑚(𝑡)=[1𝜌𝐤(0)]𝐵𝑘,𝑚(𝜙)for𝑚<𝑀𝐤.
(35)

The distribution of 𝑠𝐤,𝑚 for 𝑚 𝑀𝐤 is, in general, not of the binomial form , but it can nevertheless be given explicitly, as is detailed in Appendix . The nonbinomial form of 𝑠𝐤,𝑚 means that the reduced-dimension system is not precisely of the PA type , but it enables an efficient and very accurate solution of threshold-dynamics models; see the example in Fig. .

FIG. 8

Threshold model on a five-regular random graph, with all nodes having the same threshold 𝑀 =2; the initial condition is 𝜌(0) =0.05. The AME solution is identical to the solution of the two-dimensional system .

VIII. CONCLUSIONS

We have pointed out that many stochastic binary-state dynamical systems on networks can be described using the transition rates 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚 that were introduced in Sec. . The main results of this paper are the identification of dynamics for which the full approximate-master-equation system and can be reduced, without loss of accuracy, to a lower-dimensional system of equations, occasionally even yielding closed-form solutions [such as Eqs.  and ]. We showed in Sec.  that the pair-approximation system exactly matches the AME results for nodes in one state (e.g., the susceptible state) if the dynamics are monotone (i.e., the recovery rate 𝑅𝑘,𝑚 is identically zero) and the infection rate is linear in the number of infected neighbors ( 𝐹𝑘,𝑚 =𝑐𝑘 +𝑑𝑘𝑚). The classical SI model and the Bass diffusion model are of this type—see Figs.  and —and for Bethe lattices ( 𝑧-regular random graphs) Eq.  explicitly gives the expected fraction of infected nodes. Interestingly, the PA solution is not exactly equal to the AME for nodes in the infected state; see Fig. .

In Sec. , we showed that if both sets of transition rates 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚 are nonzero, then the AME solutions cannot, in general, be given exactly by the PA ansatz for all time 𝑡. However, in the special case defined by Eq. , which corresponds to equilibrium stochastic dynamics, i.e., those dynamics that obey a detailed balance (microscopic reversibility) condition, the AME solutions reduce to the PA solution in the limit 𝑡 ; see Figs.  and . This property allows us to perform bifurcation analysis for systems with up-down ( 𝑍2) symmetry [i.e., for dynamics obeying condition Eq. ], giving the explicit expression for the critical point where the steady-state limit of 𝜌 changes from the disordered (paramagnetic) state ‾‾‾𝜌 =1/2 to an ordered (ferromagnetic) state with ‾‾‾𝜌 1/2. The Glauber and Metropolis dynamics for the Ising spin model are important examples that obey both conditions and , and our results reproduce the critical Ising temperature that was found using very different methods (replica trick, recursion approach) in Refs. . The analysis of Refs.  for the critical behavior on scale-free networks can also be applied here and gives results described in Sec.  and Appendix . We highlight the fact that the critical point 𝑎 =1 for networks with infinite-variance degree distributions is attainable in plausible social-influence models of this type, as is the case 𝑎 <1, although for the Ising model 𝑎 =1 corresponds to infinite temperature and the regime 𝑎 <1 is unphysical. We also give an explicit expression [Eq. ] for the critical point (pitchfork bifurcation point) predicted by pair approximation of models with up-down symmetry—including nonequilibrium dynamics that obey but not necessarily —on 𝑧-regular random graphs, with the caveat that this PA critical point is not identical to the AME critical point except for the class of equilibrium models that obey as well as . Finally, in Sec. , we showed that the AME approach can be extended to include threshold models of adoption or fad diffusion and that the AME system can be reduced to a system of only two differential equations; see Eqs.  and Fig. .

The exact agreement of AME and PA solutions—whether for all time as in Sec.  or in the steady state as in Sec. —implies that higher-order correlations (beyond nearest neighbor) are correctly captured by PA in the cases we have identified. Indeed, the agreement of the theory curves with numerical simulation results in all these cases (for all time in Figs.  and , and as 𝑡 in Figs.  and ) is excellent. We interpret this agreement as indicating that for the classes of dynamics where PA is equal to AME, the results of PA are essentially exact. In contrast, in cases where PA and AME are not equal, as in Fig. , for example, it is necessary to solve the full AME system to obtain high-accuracy approximations.

The present work is focused only on infinite, uncorrelated networks with negligible levels of clustering (transitivity). Generalizing the AME approach to clustered network models and/or to networks with degree-degree correlations remains a challenge, as the added complexities will lead to even larger differential-equation systems than needed for configuration-model networks. Nevertheless, we hope the insights gained here may assist in generating and analyzing pair approximations for dynamics on such networks. Another direction for further work is to consider the effects of finite 𝑁 upon the various approximations used here. Although the local (node-level) dynamics are stochastic, the differential equations derived for the emergent dynamics are deterministic because we assume the 𝑁 limit. On smaller networks, stochastic effects will be important even at the global level, and different approaches—such as branching processes —will be required to describe the variability of results across realizations.

This work was partly funded by Science Foundation Ireland (Grants No. 11/PI/1026 and No. 09/SRC/E1780) and by the European Commission through FET-Proactive Project PLEXMATH (FP7-ICT-2011-8; Grant No. 317614). I acknowledge the SFI/HEA Irish Centre for High-End Computing (ICHEC) for the provision of computational facilities, and Sean Lyons, Davide Cellai, Sergey Melnik, Adam Hackett, and Peter Fennell for helpful discussions.

Appendix A: DERIVATION OF MASTER EQUATIONS

We consider binary-state dynamics on static, undirected, connected networks in the limit of infinite network size (i.e., 𝑁 , where 𝑁 is the number of nodes in the network). For convenience, we call the two possible node states susceptible and infected, as is common in disease-spread models. However, this approach also applies to other binary-state dynamics, such as spin systems , where each node may be in the +1 (spin-up) or the 1 (spin-down) state. The networks have degree distribution 𝑃𝑘 and are generated by the configuration model . Dynamics are stochastic and are defined by infection and recovery probabilities 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚, which depend on the degree 𝑘 of a node and on the current number 𝑚 of infected neighbors of the node; see Sec. .

We now proceed to derive the approximate master equations for dynamics of this type, closely following the approach used in Refs.  for SIS dynamics. Let 𝑆𝑘,𝑚 (respectively, 𝐼𝑘,𝑚) be the set of nodes that are susceptible (respectively, infected), have degree 𝑘, and have 𝑚 infected neighbors. To quantify the size of these sets, define 𝑠𝑘,𝑚(𝑡) [respectively, 𝑖𝑘,𝑚(𝑡)] as the fraction of 𝑘-degree nodes that are susceptible (respectively, infected) at time 𝑡 and have 𝑚 infected neighbors. Then, the fraction 𝜌𝑘(𝑡) of 𝑘-degree nodes that are infected at time 𝑡 is given by

𝜌𝑘(𝑡)=𝑘𝑚=0𝑖𝑘,𝑚=1𝑘𝑚=0𝑠𝑘,𝑚,
(A1)

and the fraction of infected nodes in the whole network is found by summing over all 𝑘 classes:

𝜌(𝑡)=𝑘𝑃𝑘𝜌𝑘(𝑡).
(A2)

If a randomly chosen fraction 𝜌(0) of nodes is initially infected, then the initial conditions for 𝑠𝑘,𝑚 and 𝑖𝑘,𝑚 are easily seen to be

𝑠𝑘,𝑚(0)=[1𝜌(0)]𝐵𝑘,𝑚[𝜌(0)],𝑖𝑘,𝑚(0)=𝜌(0)𝐵𝑘,𝑚[𝜌(0)],
(A3)

where 𝐵𝑘,𝑚(𝑞) is the binomial distribution introduced in Eq. . Note that we can also calculate the number of edges of various types using this formalism. For example, the number of edges in the network that join a susceptible node to an infected node (we call these 𝑆-𝐼 edges for short) can be expressed in two equivalent ways:

𝑁𝑘𝑃𝑘𝑘𝑚=0𝑚𝑠𝑘,𝑚or𝑁𝑘𝑃𝑘𝑘𝑚=0(𝑘𝑚)𝑖𝑘,𝑚.
(A4)

The first of these expressions, for example, follows from noting that in a sufficiently large network there are 𝑁𝑃𝑘 nodes of degree 𝑘, of which a fraction 𝑠𝑘,𝑚 is susceptible and has 𝑚 infected neighbors. Each such node contributes 𝑚 edges to the total number of 𝑆-𝐼 edges. Similar expressions may also be given for the number of 𝑆-𝑆 and 𝐼-𝐼 edges in the network. We note that the equivalence of the two expressions in is preserved by the evolution equations that are described below.

Next, we examine how the size of the 𝑆𝑘,𝑚 set changes in time. We write the general expression

𝑠𝑘,𝑚(𝑡+𝑑𝑡)=𝑠𝑘,𝑚(𝑡)𝑊(𝑆𝑘,𝑚𝐼𝑘,𝑚)𝑠𝑘,𝑚𝑑𝑡+𝑊(𝐼𝑘,𝑚𝑆𝑘,𝑚)𝑖𝑘,𝑚𝑑𝑡𝑊(𝑆𝑘,𝑚𝑆𝑘,𝑚+1)𝑠𝑘,𝑚𝑑𝑡+𝑊(𝑆𝑘,𝑚1𝑆𝑘,𝑚)𝑠𝑘,𝑚1𝑑𝑡𝑊(𝑆𝑘,𝑚𝑆𝑘,𝑚1)𝑠𝑘,𝑚𝑑𝑡+𝑊(𝑆𝑘,𝑚+1𝑆𝑘,𝑚)𝑠𝑘,𝑚+1𝑑𝑡
(A5)

to reflect all the transitions whose rate is linear in 𝑑𝑡 (all other state transitions are negligible in the 𝑑𝑡 0 limit); see Fig. . Here, 𝑊(𝑆𝑘,𝑚 𝐼𝑘,𝑚)𝑑𝑡, for example, is the probability that a node in the 𝑆𝑘,𝑚 set at time 𝑡 moves to the 𝐼𝑘,𝑚 set by time 𝑡 +𝑑𝑡. It is clear from the definitions above that

𝑊(𝑆𝑘,𝑚𝐼𝑘,𝑚)=𝐹𝑘,𝑚and𝑊(𝐼𝑘,𝑚𝑆𝑘,𝑚)=𝑅𝑘,𝑚.
(A6)

A node moves from the 𝑆𝑘,𝑚1 set to the 𝑆𝑘,𝑚 set if it remains susceptible, while one of its susceptible neighbors becomes infected. Note that this event means that an 𝑆-𝑆 edge changes to an 𝑆-𝐼 edge. If we suppose that 𝑆-𝑆 edges change to 𝑆-𝐼 edges at a (time-dependent) rate 𝛽𝑠, we can write

𝑊(𝑆𝑘,𝑚𝑆𝑘,𝑚+1)=𝛽𝑠(𝑘𝑚)and𝑊(𝑆𝑘,𝑚1𝑆𝑘,𝑚)=𝛽𝑠(𝑘𝑚+1),
(A7)

since nodes in the 𝑆𝑘,𝑚 set have 𝑘 𝑚 susceptible neighbors, while those in the 𝑆𝑘,𝑚1 set have 𝑘 𝑚 +1 susceptible neighbors. To calculate 𝛽𝑠, we count the number of 𝑆-𝑆 edges in the network at time 𝑡 and then count the number of edges that switch from being 𝑆-𝑆 edges to 𝑆-𝐼 edges in the time interval 𝑑𝑡; the probability 𝛽𝑠𝑑𝑡 is given by taking the ratio of the latter number to the former, i.e.,

𝛽𝑠𝑑𝑡=𝑘𝑃𝑘𝑘𝑚=0(𝑘𝑚)𝐹𝑘,𝑚𝑠𝑘,𝑚𝑑𝑡𝑘𝑃𝑘𝑘𝑚=0(𝑘𝑚)𝑠𝑘,𝑚.
(A8)

A similar approximation is used to define 𝛾𝑠, the rate at which 𝑆-𝐼 edges change to 𝑆-𝑆 edges because of the recovery of an infected node,

𝛾𝑠=𝑘𝑃𝑘𝑘𝑚=0(𝑘𝑚)𝑅𝑘,𝑚𝑖𝑘,𝑚𝑘𝑃𝑘𝑘𝑚=0(𝑘𝑚)𝑖𝑘,𝑚,
(A9)

and we then write

𝑊(𝑆𝑘,𝑚𝑆𝑘,𝑚1)=𝛾𝑠𝑚and𝑊(𝑆𝑘,𝑚+1𝑆𝑘,𝑚)=𝛾𝑠(𝑚+1).
(A10)

Taking the limit 𝑑𝑡 0 of Eq.  gives the master equation for the evolution of 𝑠𝑘,𝑚(𝑡) (see Fig. ):

𝑑𝑑𝑡𝑠𝑘,𝑚=𝐹𝑘,𝑚𝑠𝑘,𝑚+𝑅𝑘,𝑚𝑖𝑘,𝑚𝛽𝑠(𝑘𝑚)𝑠𝑘,𝑚+𝛽𝑠(𝑘𝑚+1)𝑠𝑘,𝑚1𝛾𝑠𝑚𝑠𝑘,𝑚+𝛾𝑠(𝑚+1)𝑠𝑘,𝑚+1,
(A11)

where 𝑚 is in the range 0,,𝑘 for each 𝑘 class in the network (and adopting the convention 𝑠𝑘,1 𝑠𝑘,𝑘+1 0). Applying identical arguments, mutatis mutandis, to the set 𝐼𝑘,𝑚, we derive the corresponding system of equations for 𝑖𝑘,𝑚(𝑡):

𝑑𝑑𝑡𝑖𝑘,𝑚=𝑅𝑘,𝑚𝑖𝑘,𝑚+𝐹𝑘,𝑚𝑠𝑘,𝑚𝛽𝑖(𝑘𝑚)𝑖𝑘,𝑚+𝛽𝑖(𝑘𝑚+1)𝑖𝑘,𝑚1𝛾𝑖𝑚𝑖𝑘,𝑚+𝛾𝑖(𝑚+1)𝑖𝑘,𝑚+1
(A12)

for 𝑚 =0,,𝑘 and for each 𝑘 class in the network, with time-dependent rates 𝛽𝑖 and 𝛾𝑖 defined through 𝑠𝑘,𝑚 and 𝑖𝑘,𝑚 as

𝛽𝑖=𝑘𝑃𝑘𝑘𝑚=0𝑚𝐹𝑘,𝑚𝑠𝑘,𝑚𝑘𝑃𝑘𝑘𝑚=0𝑚𝑠𝑘,𝑚and𝛾𝑖=𝑘𝑃𝑘𝑘𝑚=0𝑚𝑅𝑘,𝑚𝑖𝑘,𝑚𝑘𝑃𝑘𝑘𝑚=0𝑚𝑖𝑘,𝑚.
(A13)
FIG. 9

Schematic of transitions to or from the 𝑆𝑘,𝑚 and 𝐼𝑘,𝑚 sets, as described in Eqs. . For each set, the central (ego) node is shown along with some of its neighbors: White nodes are susceptible or inactive or spin-down, and black nodes are infected or active or spin-up. See also Fig. 1 of Ref. .

The approximate master equations and , with the time-dependent rates 𝛽𝑠, 𝛾𝑠, 𝛽𝑖, and 𝛾𝑖 (defined as nonlinear functions of 𝑠𝑘,𝑚 and 𝑖𝑘,𝑚), form a closed system of deterministic equations that, along with initial conditions , can be solved numerically using standard methods . Note that the evolution equations are completely prescribed by the functions 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚, and so this method can be applied to any stochastic dynamical process defined by transition rates 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚. For the SIS model, Eqs.  and were derived in Ref.  (see also Refs. ), with additional terms to study the adaptive rewiring of the network.

Appendix B: MICROSCOPIC REVERSIBILITY IN EQUILIBRIUM MODELS

We consider transition rates 𝐹𝑘,𝑚 and 𝑅𝑘,𝑚 for which a detailed balance holds, i.e., for which the dynamics exhibit microscopic reversibility . Figure  shows a pair of connected nodes: Node 1—on the left—has degree 𝑘1, and node 2—on the right—has degree 𝑘2. We assume that 𝑚1 neighbors of node 1, other than node 2, are in the active (infected) state; similarly, 𝑚2 is the number of infected neighbors, other than node 1, of node 2.

FIG. 10

Rates for state transitions of a connected pair of nodes, assuming that the states of all their other neighbors remain unchanged. Node 1 (left-hand node) has degree 𝑘1 and has 𝑚1 infected neighbors that are not shown. Node 2 (right-hand node) has degree 𝑘2 and has 𝑚2 neighbors—other than node 1—that are infected.

Starting from the state at the top of Fig. , where both node 1 and node 2 are in the susceptible state, we consider possible cycles of state transitions for the node pair, assuming that 𝑚1 and 𝑚2 remain unchanged. The transition from the {𝑆,𝑆} state (top of the figure) to the {𝑆,𝐼} state (right of the figure), for example, occurs at the rate 𝐹𝑘2,𝑚2, since node 2 becomes infected while it has 𝑚2 infected neighbors. The transition from {𝑆,𝐼} (right of the figure) to {𝐼,𝐼} (bottom of the figure) occurs at the rate 𝐹𝑘1,𝑚1+1, since node 1 is required to become infected at a time when it has 𝑚1 +1 infected neighbors. The other rates are derived similarly.

For microscopic reversibility, it is necessary that the product of the rates around a closed cycle is independent of the direction of rotation around the cycle . From Fig. , this condition means we must have

𝐹𝑘2,𝑚2𝐹𝑘1,𝑚1+1𝑅𝑘2,𝑚2+1𝑅𝑘1,𝑚1=𝐹𝑘1,𝑚1𝐹𝑘2,𝑚2+1𝑅𝑘1,𝑚1+1𝑅𝑘2,𝑚2.
(B1)

Rearranging gives the condition

𝐹𝑘1,𝑚1+1𝐹𝑘1,𝑚1𝑅𝑘1,𝑚1𝑅𝑘1,𝑚1+1=𝐹𝑘2,𝑚2+1𝐹𝑘2,𝑚2𝑅𝑘2,𝑚2𝑅𝑘2,𝑚2+1.
(B2)

Since 𝑘1 and 𝑘2 (and 𝑚1 and 𝑚2) are independent, this condition requires that

𝐹𝑘,𝑚+1𝐹𝑘,𝑚𝑅𝑘,𝑚𝑅𝑘,𝑚+1=𝑎for all𝑘and all𝑚𝑘,
(B3)

where 𝑎 is a constant. Defining 𝑢𝑚 =𝐹𝑘,𝑚/𝑅𝑘,𝑚, this equation can be written as the first-order difference equation

𝑢𝑚+1=𝑎𝑢𝑚,
(B4)

with the solution 𝑢𝑚 =𝑎𝑚𝑢0. In terms of the rates, this solution means that detailed balance requires that

𝐹𝑘,𝑚𝑅𝑘,𝑚=𝑎𝑚𝐹𝑘,0𝑅𝑘,0.
(B5)

Identifying 𝑏𝑘 as 𝐹𝑘,0/𝑅𝑘,0 yields Eq. .

Appendix C: DERIVATION OF EQUATION 

Here, we manipulate Eqs.  and for equilibrium models to produce the solution and the algebraic equation  for ‾‾‾𝑝. Starting with Eq. , we can solve for ‾‾‾𝑞, yielding

‾‾‾𝑞=‾‾‾𝑝𝑎1‾‾‾𝑝+‾‾‾𝑝𝑎.
(C1)

Solving for ‾‾‾𝜌𝑘 yields

‾‾‾𝜌𝑘=𝑏𝑘𝑏𝑘+(1‾‾‾𝑞1‾‾‾𝑝)𝑘.
(C2)

Inserting expression into gives the desired Eq. .

Next, we consider the identity , for which the PA ansatz gives

𝑘𝑃𝑘(1‾‾‾𝜌𝑘)𝑘‾‾‾𝑝𝑘=𝑘𝑃𝑘‾‾‾𝜌𝑘𝑘(1‾‾‾𝑞𝑘).
(C3)

Since ‾‾‾𝑝𝑘 and ‾‾‾𝑞𝑘 are independent of 𝑘 for the case considered here—see Eq. —we obtain the simpler expression

‾‾‾𝑝(1‾‾‾𝜔)=‾‾‾𝜔(1‾‾‾𝑞),
(C4)

where ‾‾‾𝜔 =𝑘𝑧‾‾‾𝜌𝑘. Solving Eq.  for ‾‾‾𝜔 and using Eq.  gives

‾‾‾𝜔=‾‾‾𝑝(1‾‾‾𝑝+‾‾‾𝑝𝑎)1‾‾‾𝑝2+‾‾‾𝑝2𝑎,
(C5)

and equating this expression to 𝑘𝑘𝑧𝑃𝑘‾‾‾𝜌𝑘, with Eq. , gives Eq. .

Appendix D: CRITICAL POINT OF EQUILIBRIUM MODELS WITH UP-DOWN SYMMETRY

Equilibrium spin models with up-down symmetry are discussed in Sec. , and it is shown there that they constitute a special case of the class of models defined by relation . The up-down symmetry imposes condition on the parameters in the solution , so that ‾‾‾𝜌𝑘 can be expressed in the form

‾‾‾𝜌𝑘=11+(𝑎1‾‾‾𝑝+‾‾‾𝑝𝑎)𝑘.
(D1)

Recall from Sec.  that a symmetric solution of the AME with 𝜌𝑘 =1/2 exists for all 𝑡; comparison with shows that this solution corresponds to the case with ‾‾‾𝑝 = 1/(1 +𝑎). Inserting this solution into expression for ‾‾‾𝜔—the probability that one end of a randomly chosen edge is infected—gives ‾‾‾𝜔 =1/2 for the symmetric solution.

Next, we investigate the possibility of other solutions lying near the symmetric solution; in the language of dynamical systems, we construct the normal form of the bifurcation . We choose ‾‾‾𝜔 as the order parameter and let ‾‾‾𝜔 =12 +𝜖, for small 𝜖, to explore the neighborhood of the symmetric solution. Equation  can then be inverted to give ‾‾‾𝑝 in terms of 𝜖, and inserting the resulting expression into Eq.  leads to a self-consistent equation for 𝜖:

𝜖+12=𝑘𝑘𝑧𝑃𝑘[1+(2𝜖+𝑎+4(1𝑎)𝜖2(1+2𝜖)𝑎)𝑘]1.
(D2)

After some rearrangement, this equation can be written as

𝜖=12𝑘𝑘𝑧𝑃𝑘tanh[𝑘2ln(2𝜖+𝑎+4(1𝑎)𝜖2(1+2𝜖)𝑎)],
(D3)

and, in the case 𝑎 =𝑒4𝐽/𝑇, this expression is very similar to Eq. (13) of Ref.  and Eq. (28) of Ref. , which were derived—using very different methods (the recursion method and the replica trick)—for the Ising model.

It follows that the structure of the solutions of Eq.  near the critical point (bifurcation point) can be analyzed—for the class of models obeying both and —using the same methods for as used in Refs. . If 𝑃𝑘 has a finite fourth moment, for example, the right-hand side of can be expanded as a Taylor series in small 𝜖, giving

𝜖𝑐1𝜖+𝑐3𝜖3+𝒪(𝜖5)as𝜖0,
(D4)

where 𝑐1 =𝑘22𝑧(1 1𝑎), and 𝑐3 is a coefficient involving moments of 𝑃𝑘 up to the fourth. Equation  possesses the symmetric solution 𝜖 =0 for all parameter values, but solutions with nonzero 𝜖 can also be found, if leading-order terms balance such that

𝑐11+𝑐3𝜖2=0.
(D5)

The 𝜖-independent coefficient 𝑐1 1 in this equation is negative for small values of the parameter 𝑎, but it is positive when 𝑎 exceeds the value 𝑎𝑐 given by Eq. , while 𝑐3 is negative at 𝑎 =𝑎𝑐. Thus, Eq.  has real-valued solutions for 𝜖, provided that 𝑎 𝑎𝑐, and so 𝑎𝑐 marks the pitchfork bifurcation point (the critical point).

The case 𝑃𝑘 𝑘𝛾 can be examined in the same fashion as in Refs. . The small- 𝜖 expansion of has leading-order terms of the form

𝜖𝑐1𝜖+𝑐3𝜖3+𝑐0𝜖𝛾2as𝜖0.
(D6)

If 𝛾 is in the range 2 <𝛾 <3, the 𝜖𝛾2 term dominates both the linear and cubic terms, and the critical point is determined by the vanishing of the coefficient 𝑐0: This critical value occurs at 𝑎 =𝑎𝑐 =1. For exponents in the range 3 <𝛾 <5, the linear term dominates the 𝜖𝛾2 term, and the critical point is again but with a different scaling near criticality. For 𝛾 >5, we recover the case .

Appendix E: PA CRITICAL POINT OF SYMMETRIC MODELS ON 𝑧-REGULAR RANDOM GRAPHS

For models with the symmetry but not necessarily possessing the equilibrium property , we focus here on the steady state of the PA equations  on 𝑧-regular graphs. In particular, we suppose that the dynamics (through the rates 𝐹𝑘,𝑚) depend on a parameter 𝑄, and we derive the condition determining the critical value (the bifurcation point) of this parameter.

We begin by noting that the property means that ‾‾‾𝜌 =1/2 is always a steady-state solution of . Using this value and property on the right-hand side of the first of the PA equations  gives the steady-state relation

0=12𝑧𝑚=0𝐹𝑧,𝑧𝑚𝐵𝑧,𝑚(‾‾‾𝑞)+12𝑧𝑚=0𝐹𝑧,𝑚𝐵𝑧,𝑚(‾‾‾𝑝)=12𝑧𝑚=0𝐹𝑧,𝑚[𝐵𝑧,𝑚(1‾‾‾𝑞)𝐵𝑧,𝑚(‾‾‾𝑝)],
(E1)

which is satisfied—as are the other steady PA equations—if ‾‾‾𝑞 =1 ‾‾‾𝑝 in the symmetric state. Using this condition to replace ‾‾‾𝑞 in the second of the PA equations gives, after some manipulation, the implicit relation for the steady-state value of ‾‾‾𝑝:

𝑧𝑚=0(12𝑚𝑧)𝐹𝑧,𝑚𝐵𝑧,𝑚(‾‾‾𝑝)=0.
(E2)

Next, we consider the possibility of steady states that are distinct from the symmetric state solution discussed above; we introduce the symbol 𝜎 as a convenient shorthand for the symmetric state: 𝜎 ={‾‾‾𝜌 =1/2, ‾‾‾𝑝 solves , ‾‾‾𝑞 =1 ‾‾‾𝑝}. For a general—possibly nonsymmetric—steady state, the first of the PA equations is

0=‾‾‾𝜌𝑧𝑚=0𝐹𝑧,𝑚𝐵𝑧,𝑚(1‾‾‾𝑞)+(1‾‾‾𝜌)𝑧𝑚=0𝐹𝑧,𝑚𝐵𝑧,𝑚(‾‾‾𝑝).
(E3)

We now treat ‾‾‾𝜌, ‾‾‾𝑝, and ‾‾‾𝑞 as implicit functions of 𝑄, the parameter defining the rates 𝐹𝑧,𝑚. Differentiating both sides of Eq.  with respect to 𝑄, evaluating at the symmetric state 𝜎, using the relations ‾‾‾𝑞|𝜎 =1 ‾‾‾𝑝|𝜎, and, from ,

𝑄‾‾‾𝑞| | | |𝜎=𝑄(1‾‾‾𝜌‾‾‾𝜌‾‾‾𝑝)| | | |𝜎=4‾‾‾𝜌𝑄| | | |𝜎‾‾‾𝑝|𝜎‾‾‾𝑝𝑄| | | |𝜎,
(E4)

leads eventually to the relation

0=2‾‾‾𝜌𝑄[1+𝑧22‾‾‾𝑝11‾‾‾𝑝]𝑧𝑚=0𝐹𝑧,𝑚𝐵𝑧,𝑚(‾‾‾𝑝),
(E5)

where all terms are evaluated at the symmetric state 𝜎. For this equation to be true, we must either have ‾‾‾𝜌/𝑄 =0 or the term in square brackets must vanish. (Note that the third term is a sum of positive terms, so it cannot be zero.) The symmetry-breaking bifurcation points correspond to the vanishing of the square bracketed term, and this condition gives the critical value of ‾‾‾𝑝:

‾‾‾𝑝𝑐=𝑧22𝑧2.
(E6)

Finally, inserting this value into the general steady-state solution gives the criticality condition on the rate 𝐹𝑧,𝑚 ( =𝑅𝑧,𝑧𝑚) of the PA dynamics.

Appendix F: DERIVATION OF REDUCED-DIMENSION EQUATIONS FOR THRESHOLD MODELS

Our goal here is to demonstrate that Eqs.  and reduce to Eqs.  through an exact solution of given by .

We proceed to insert the ansatz into the AME for 𝑚 <𝑀𝐤 (for which 𝑚 values we have 𝐹𝑘,𝑚 =0). The left-hand side of then gives, after some rearrangement,

˙𝑠𝐤,𝑚=[1𝜌𝐤(0)][𝑚𝜙𝑘𝑚1𝜙]𝐵𝑘,𝑚(𝜙)˙𝜙,
(F1)

where dots denote time derivatives. Using the identity

𝐵𝑘,𝑚1(𝜙)=1𝜙𝜙𝑚𝑘𝑚+1𝐵𝑘,𝑚(𝜙)
(F2)

on the right-hand side of Eq.  yields the condition

˙𝜙=𝛽𝑠(1𝜙)
(F3)

on the function 𝜙(𝑡) for the ansatz to be an exact solution of . From the initial condition on 𝑠𝐤,𝑚, we also obtain

𝜙(0)=𝜌(0)=𝐤𝑘𝑧𝑃𝐤𝜌𝐤(0).
(F4)

Comparing and , we see that it remains for us to show that

𝛽𝑠=𝑔(𝜙)𝜙1𝜙.
(F5)

We first prove a preliminary and rather general result, as follows. Multiplying the AME by (𝑘 𝑚)𝑃𝐤 and summing over 𝑚 =0,,𝑘 and over the 𝐤 classes gives

𝑑𝑑𝑡𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚=𝐤𝑃𝐤𝑚𝐹𝐤,𝑚(𝑘𝑚)𝑠𝐤,𝑚𝛽𝑠𝐤𝑃𝐤𝑚[(𝑘𝑚)2𝑠𝐤,𝑚(𝑘𝑚)(𝑘𝑚+1)𝑠𝐤,𝑚1].
(F6)

Note that the second term on the right-hand side is a telescoping series that reduces to

𝛽𝑠𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚
(F7)

and that the definition of 𝛽𝑠 enables us to express the first term on the right-hand side as

𝛽𝑠𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚.
(F8)

Therefore, can be solved for 𝛽𝑠 to yield

𝛽𝑠=12𝑑𝑑𝑡𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚=12𝑑𝑑𝑡ln(𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚).
(F9)

Rewriting as 𝛽𝑠 =𝑑𝑑𝑡ln(1 𝜙) and comparing with gives the result that 𝐤𝑃𝐤𝑚(𝑘 𝑚)𝑠𝐤,𝑚 and (1𝜙)2 are equal, up to a multiplicative constant. Using the initial conditions on 𝑠𝐤,𝑚 and 𝜙 determines the constant, giving

𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚=𝑧(1𝜙)2.
(F10)

We now use this result and the definition of 𝛽𝑠 to prove . For the infection probabilities , 𝛽𝑠 is given by

𝛽𝑠=𝐤𝑃𝐤𝑚𝑀𝐤(𝑘𝑚)𝑠𝐤,𝑚𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚=𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚𝐤𝑃𝐤𝑚<𝑀𝐤(𝑘𝑚)𝑠𝐤,𝑚𝐤𝑃𝐤𝑚(𝑘𝑚)𝑠𝐤,𝑚=𝑧(1𝜙)2𝐤𝑃𝐤𝑚<𝑀𝐤(𝑘𝑚)[1𝜌𝐤(0)]𝐵𝑘,𝑚𝑧(1𝜙)2,
(F11)

where the final line uses twice and the ansatz . Dividing the numerator and denominator by 𝑧(1 𝜙) and using the identities (𝑘 𝑚)𝐵𝑘,𝑚(𝜙) =𝑘(1 𝜙)𝐵𝑘1,𝑚(𝜙) and 𝑚<𝑀𝐤𝐵𝑘1,𝑚(𝜙) =1 𝑚𝑀𝐤𝐵𝑘1,𝑚(𝜙) leads to

𝛽𝑠=𝑔(𝜙)𝜙1𝜙,
(F12)

as claimed.

Finally, let us show that the fraction 𝜌(𝑡) of active nodes given in the AME as 𝜌(𝑡) =1 𝐤𝑃𝐤𝑚𝑠𝐤,𝑚 obeys Eqs.  and . Multiplying Eq.  by 𝑃𝐤 and summing over 𝑚 and over 𝐤 classes gives

𝐤𝑃𝐤𝑚˙𝑠𝐤,𝑚=𝐤𝑃𝐤𝑚𝐹𝐤,𝑚𝑠𝐤,𝑚𝛽𝑠𝐤𝑃𝐤𝑚[(𝑘𝑚)𝑠𝐤,𝑚(𝑘𝑚+1)𝑠𝐤,𝑚1],
(F13)

where the second term on the right-hand side is easily seen to telescope to zero. Thus, we have

˙𝜌(𝑡)=𝐤𝑃𝐤𝑚𝐹𝐤,𝑚𝑠𝐤,𝑚=𝐤𝑃𝐤𝑚𝑠𝐤,𝑚𝐤𝑃𝐤𝑚<𝑀𝐤𝑠𝐤,𝑚=1𝜌𝐤𝑃𝐤[1𝜌𝐤(0)][1𝑚𝑀𝐤𝐵𝑘,𝑚(𝜙)],
(F14)

using . It is straightforward to verify that this equation reduces to Eqs.  and .

For completeness, note that the distribution of 𝑠𝐤,𝑚 for 𝑚 𝑀𝐤 is, in general, not of the binomial form , which applies only to 𝑚 values below the threshold 𝑀𝐤. To obtain the values of 𝑠𝐤,𝑚 for these 𝑚 values, note that Eq.  has a solution giving 𝑠𝐤,𝑚 explicitly in terms of 𝑠𝐤,𝑚1 by using the integrating factor

exp(𝑡+(𝑘𝑚)𝛽𝑠𝑑𝑡)=𝑒𝑡(1𝜙)(𝑘𝑚)
(F15)

and using for 𝛽𝑠. Then, we have

𝑠𝐤,𝑚(𝑡)=𝑠𝐤,𝑚(0)𝑒𝑡(1𝜙(𝑡)1𝜙(0))𝑘𝑚+(𝑘𝑚+1)×𝑒𝑡(1𝜙)𝑘𝑚𝑡0𝑒𝜏[1𝜙(𝜏)](𝑘𝑚)×𝛽𝑠(𝜏)𝑠𝐤,𝑚1(𝜏)𝑑𝜏,
(F16)

which can be solved recursively for 𝑚 =𝑀𝐤 to 𝑚 =𝑘.

References (86)

  1. C. Castellano, S. Fortunato, and V. Loreto, Statistical Physics of Social Dynamics, Rev. Mod. Phys. 81, 591 (2009).
  2. C. Castellano, Social Influence and the Dynamics of Opinions: The Approach of Statistical Physics, Managerial and decision economics : MDE 33, 311 (2012).
  3. A. Barrat, M. Barthélemy, and A. Vespignani, Dynamical Processes on Complex Networks (Cambridge University Press, Cambridge, England, 2008).
  4. M. E. J. Newman, Networks: An Introduction (Oxford University Press, Oxford, England, 2010).
  5. T. M. Liggett, Interacting Particle Systems (Springer, New York, 1985).
  6. V. Sood and S. Redner, Voter Model on Heterogeneous Graphs, Phys. Rev. Lett. 94, 178701 (2005).
  7. M. J. de Oliveira, Isotropic Majority-Vote Model on a Square Lattice, J. Stat. Phys. 66, 273 (1992).
  8. L. F. C. Pereira and F. G. B. Moreira, Majority-Vote Model on Random Graphs, Phys. Rev. E 71, 016123 (2005).
  9. R. Pastor-Satorras and A. Vespignani, Epidemic Spreading in Scale-Free Networks, Phys. Rev. Lett. 86, 3200 (2001).
  10. A. L. Hill, D. G. Rand, M. A. Nowak, and N. A. Christakis, Infectious Disease Modeling of Social Contagion in Networks, PLoS Comput. Biol. 6, e1000968 (2010).
  11. H. P. Young, Innovation Diffusion in Heterogeneous Populations: Contagion, Social Influence, and Social Learning, Am. Econ. Rev. 99, 1899 (2009).
  12. M. Granovetter, Threshold Models of Collective Behavior, Am. J. Sociology 83, 1420 (1978).
  13. D. J. Watts, A Simple Model of Global Cascades on Random Networks, Proc. Natl. Acad. Sci. U.S.A. 99, 5766 (2002).
  14. D. Centola, V. M. Eguíluz, and M. W. Macy, Cascade Dynamics of Complex Propagation, Physica (Amsterdam) 374A, 449 (2007).
  15. P. S. Dodds, K. D. Harris, and C. M. Danforth, Limited Imitation Contagion on Random Networks: Chaos, Universality, and Unpredictability, Phys. Rev. Lett. 110, 158701 (2013).
  16. D. M. Romero, B. Meeder, and J. Kleinberg, Proceedings of the 20th International Conference on the World Wide Web: Differences in the Mechanics of Information Diffusion Across Topics (ACM, New York, 2011), p. 695.
  17. G. Ver Steeg, R. Ghosh, and K. Lerman, Proceedings of the Fifth International AAAI Conference on Weblogs and Social Media (AAAI Press, Menlo Park, CA, 2011).
  18. D. Easley and J. Kleinberg, Networks, Crowds, and Markets (Cambridge University Press, New York, 2010).
  19. S. González-Bailón, J. Borge-Holthoefer, A. Rivero, and Y. Moreno, The Dynamics of Protest Recruitment through an Online Network, Sci. Rep. 1, 197 (2011).
  20. E. Bakshy, J. M. Hofman, W. A. Mason, and D. J. Watts, Proceedings of the 4th ACM International Conference on Web Search and Data Mining (ACM, New York, 2011), p. 65.
  21. N. O. Hodas and K. Lerman, How Visibility and Divided Attention Constrain Social Contagion, arXiv:1205.2736.
  22. D. Centola, The Spread of Behavior in an Online Social Network Experiment, Science 329, 1194 (2010).
  23. J. Shao, S. Havlin, and H. E. Stanley, Dynamic Opinion Model and Invasion Percolation, Phys. Rev. Lett. 103, 018701 (2009).
  24. F. J. Pérez-Reche, J. J. Ludlam, S. N. Taraskin, and C. A. Gilligan, Synergy in Spreading Processes: From Exploitative to Explorative Foraging Strategies, Phys. Rev. Lett. 106, 218701 (2011).
  25. S. N. Dorogovtsev, A. V. Goltsev, and J. F. F. Mendes, Ising Model on Networks with an Arbitrary Distribution of Connections, Phys. Rev. E 66, 016104 (2002).
  26. M. Leone, A. Vázquez, A. Vespignani, and R. Zecchina, Ferromagnetic Ordering in Graphs with Arbitrary Degree Distribution, Eur. Phys. J. B 28, 191 (2002).
  27. S. Galam, Y. Gefen, and Y. Shapir, Sociophysics: A New Approach of Sociological Collective Behaviour, J. Math Sociology 9, 1 (1982).
  28. F. Schweitzer and L. Behera, Nonlinear Voter Models: The Transition from Invasion to Coexistence, Eur. Phys. J. B 67, 301 (2009).
  29. C. Castellano and R. Pastor-Satorras, Universal and Nonuniversal Features of the Generalized Voter Class for Ordering Dynamics in Two Dimensions, Phys. Rev. E 86, 051123 (2012).
  30. F. Vazquez and C. López, Systems with Two Symmetric Absorbing States: Relating the Microscopic Dynamics with the Macroscopic Behavior, Phys. Rev. E 78, 061127 (2008).
  31. Dynamical correlations should not be confused with degree-degree (structural) correlations, which are a property of the network rather than of the dynamics. In this paper, we assume zero structural correlations in the networks (uncorrelated configuration model [4, 32, 33]).

  32. E. A. Bender and E. R. Canfield, The Asymptotic Number of Labeled Graphs with Given Degree Sequences, J. Comb. Theory Ser. A 24, 296 (1978).
  33. B. Bollobás, A Probabilistic Proof of an Asymptotic Formula for the Number of Labelled Regular Graphs, Eur. J. Combinatorics 1, 311 (1980).
  34. C. Castellano and R. Pastor-Satorras, Zero Temperature Glauber Dynamics on Complex Networks, J. Stat. Mech. 05 (2006) P05001.
  35. J. P. Gleeson, S. Melnik, J. A. Ward, M. A. Porter, and P. J. Mucha, Accuracy of Mean-Field Theory for Dynamics on Real-World Networks, Phys. Rev. E 85, 026106 (2012).
  36. S. Melnik, A. Hackett, M. A. Porter, P. J. Mucha, and J. P. Gleeson, The Unreasonable Effectiveness of Tree-Based Theory for Networks with Clustering, Phys. Rev. E 83, 036112 (2011).
  37. K. T. D. Eames and M. J. Keeling, Modeling Dynamic and Network Heterogeneities in the Spread of Sexually Transmitted Diseases, Proc. Natl. Acad. Sci. U.S.A. 99, 13330 (2002).
  38. S. A. Levin and R. Durrett, From Individuals to Epidemics, Phil. Trans. R. Soc. B 351, 1615 (1996).
  39. M. J. de Oliveira, J. F. F. Mendes, and M. A. Santos, Nonequilibrium Spin Models with Ising Universal Behaviour, J. Phys. A 26, 2317 (1993).
  40. M. Taylor, P. L. Simon, D. M. Green, T. House, and I. Z. Kiss, From Markovian to Pairwise Epidemic Models and the Performance of Moment Closure Approximations, J. Math. Biol. 64, 1021 (2012).
  41. F. Vazquez and V. M. Eguíluz, Analytical Solution of the Voter Model on Uncorrelated Networks, New J. Phys. 10, 063011 (2008).
  42. F. Vazquez, X. Castelló, and M. San Miguel, Agent Based Models of Language Competition: Macroscopic Descriptions and Order-Disorder Transitions, J. Stat. Mech. 04 (2010) P04007.
  43. R. Dickman, Kinetic Phase Transitions in a Surface-Reaction Model: Mean-Field Theory, Phys. Rev. A 34, 4246 (1986).
  44. V. Marceau, P.-A. Noël, L. Hébert-Dufresne, A. Allard, and L. J. Dubé, Adaptive Networks: Coevolution of Disease and Topology, Phys. Rev. E 82, 036116 (2010).
  45. J. Lindquist, J. Ma, P. van den Driessche, and F. H. Willeboordse, Effective Degree Network Disease Models, J. Math. Biol. 62, 143 (2011).
  46. J. P. Gleeson, High-Accuracy Approximation of Binary-State Dynamics on Networks, Phys. Rev. Lett. 107, 068701 (2011).
  47. T. Petermann and P. De Los Rios, Cluster Approximations for Epidemic Processes: A Systematic Description of Correlations Beyond the Pair Level, J. Theor. Biol. 229, 1 (2004).
  48. R. Durrett, J. P. Gleeson, A. L. Lloyd, P. J. Mucha, F. Shi, D. Sivakoff, J. E. S. Socolar, and C. Varghese, Graph Fission in an Evolving Voter Model, Proc. Natl. Acad. Sci. U.S.A. 109, 3682 (2012).
  49. M. Taylor, T. J. Taylor, and I. Z. Kiss, Epidemic Threshold and Control in a Dynamic Network, Phys. Rev. E 85, 016103 (2012).
  50. M. E. J. Newman, Random Graphs with Clustering, Phys. Rev. Lett. 103, 058701 (2009); J. C. Miller, Percolation and Epidemics in Random Clustered Networks, Phys. Rev. E 80, 020901(R) (2009).
  51. J. P. Gleeson, Bond Percolation on a Class of Clustered Random Networks, Phys. Rev. E 80, 036107 (2009).
  52. A. Hackett, S. Melnik, and J. P. Gleeson, Cascades on a Class of Clustered Random Networks, Phys. Rev. E 83, 056107 (2011).
  53. L. Hébert-Dufresne, P.-A. Noël, V. Marceau, A. Allard, and L. J. Dubé, Propagation Dynamics on Networks Featuring Complex Topologies, Phys. Rev. E 82, 036115 (2010).
  54. The case of symmetric models, where the dynamics are unchanged by simultaneous flipping of all nodes’ spins, is considered in detail in Sec. VI.

  55. N. T. J. Bailey, The Mathematical Theory of Infectious Diseases (Griffin, London, 1975).
  56. R. M. Anderson and R. M. May, Infectious Diseases in Humans (Oxford University Press, Oxford, England, 1992).
  57. T. E. Harris, Contact Processes on a Lattice, Ann. Probab. 2, 969 (1974).
  58. F. M. Bass, A New Product Growth Model for Consumer Durables, Management Science 15, 215 (1969).
  59. P. S. Dodds and D. J. Watts, Universal Behavior in a Generalized Model of Contagion, Phys. Rev. Lett. 92, 218701 (2004).
  60. D. J. Watts and P. S. Dodds, Influentials, Networks, and Public Opinion Formation, J. Consumer Res. 34, 441 (2007).
  61. A. Kirman, Ants, Rationality, and Recruitment, The Quarterly Journal of Economics 108, 137 (1993).
  62. S. Alfarano, T. Lux, and F. Wagner, Estimation of Agent-Based Models: The Case of an Asymmetric Herding Model, Computational Economics 26, 19 (2005).
  63. V. Gontis, A. Kononovicius, and S. Reimann, The Class of Nonlinear Stochastic Models as a Background for the Bursty Behavior in Financial Markets, Adv. Compl. Syst. 15, 1250071 (2012).
  64. J. Molofsky, R. Durrett, J. Dushoff, D. Griffeath, and S. Levin, Local Frequency Dependence and Global Coexistence, Theor. Popul. Biol. 55, 270 (1999).
  65. D. M. Abrams and S. H. Strogatz, Modelling the Dynamics of Language Death, Nature (London) 424, 900 (2003).
  66. P. L. Krapivsky, S. Redner, and E. Ben-Naim, A Kinetic View of Statistical Physics (Cambridge University Press, New York, 2010).
  67. R. J. Glauber, Time-Dependent Statistics of the Ising Model, J. Math. Phys. (N.Y.) 4, 294 (1963).
  68. N. Metropolis, A. W. Rosenbluth, M. N. Rosenbluth, A. H. Teller, and E. Teller, Equation of State Calculations by Fast Computing Machines, J. Chem. Phys. 21, 1087 (1953).
  69. J. P. Gleeson and D. J. Cahalane, Seed Size Strongly Affects Cascades on Random Networks, Phys. Rev. E 75, 056103 (2007).
  70. S. Morris, Contagion, Rev. Econ. Stud. 67, 57 (2000).
  71. A. Montanari and A. Saberi, The Spread of Innovations in Social Networks, Proc. Natl. Acad. Sci. U.S.A. 107, 20196 (2010).
  72. Octave or Matlab files for implementing and solving the AME, PA, and MF systems in Eqs. (1)–(5) are available for download from the author’s Web site, www.ul.ie/gleesonj. The documentation includes commands to reproduce the results shown in the figures of this paper and in Ref. [46].
  73. The case where , with some nonzero, is obviously equivalent, up to the swapping of state labels, to the case considered here, and so is not considered separately.

  74. The networks used in numerical simulations are sufficiently large enough to ensure that there are multiple nodes of every degree class considered in the AME: This condition avoids finite-size effects because of having small numbers of nodes described by a population-level quantity such as .

  75. S. H. Strogatz, Nonlinear Dynamics and Chaos (Perseus, Cambridge, MA, 2000).
  76. S. N. Dorogovtsev, A. V. Goltsev, and J. F. F. Mendes, Critical Phenomena in Complex Networks, Rev. Mod. Phys. 80, 1275 (2008).
  77. J.-M. Drouffe and C. Godrèche, Phase Ordering and Persistence in a Class of Stochastic Processes Interpolating between the Ising and Voter Models, J. Phys. A 32, 249 (1999).
  78. A. Galstyan and P. Cohen, Cascading Dynamics in Modular Networks, Phys. Rev. E 75, 036109 (2007).
  79. The extreme case is called synchronous updating [15, 81, 86]; however, our interest is in the asynchronous-updating case, with infinitesimally small, which gives continuous-time dynamics in the limit .

  80. Chapter 19 of Ref. [18] shows how a coordination game played on a network—as a model of direct-benefit effects from the diffusion of new behavior [70]—can also be expressed as a threshold model of the type considered here.

  81. J. P. Gleeson, Cascades on Correlated and Modular Random Networks, Phys. Rev. E 77, 046117 (2008).
  82. A. W. Hackett, Ph.D. thesis, University of Limerick, 2011.
  83. P. A. Noël, A. Allard, L. Hébert-Dufresne, V. Marceau, and L. J. Dubé, Propagation on Networks: An Exact Alternative Perspective, Phys. Rev. E 85, 031118 (2012).
  84. The main approximation here is to assume that the edge-state transition rate is the same for all edges in the network, regardless of their local neighborhood—the same assumption is made for the other rates , , and . See also the explanation in Ref. [44] for the SIS case.

  85. P.-A. Noël, B. Davoudi, R. C. Brunham, L. J. Dubé, and B. Pourbohloul, Time Evolution of Epidemic Disease on Finite and Infinite Networks, Phys. Rev. E 79, 026101 (2009).
  86. S. Gómez, A. Arenas, J. Borge-Holthoefer, S. Meloni, and Y. Moreno, Discrete-Time Markov Chain Approach to Contact-Based Disease Spreading in Complex Networks, Europhys. Lett. 89, 38009 (2010).