Home | About | Journals | Submit | Contact Us | Français |

**|**Scientific Reports**|**PMC5587719

Formats

Article sections

Authors

Related links

Sci Rep. 2017; 7: 10663.

Published online 2017 September 6. doi: 10.1038/s41598-017-08828-8

PMCID: PMC5587719

U.S. Naval Research Laboratory, Code 6792, Plasma Physics Division, Nonlinear Systems Dynamics Section, Washington, DC 20375 United States

Jason Hindes, Email: lim.yvan.lrn@rtc.sednih.nosaj.

Received 2017 May 8; Accepted 2017 July 21.

Copyright © The Author(s) 2017

We propose an analytical technique to study large fluctuations and switching from internal noise in complex networks. Using order-disorder kinetics as a generic example, we construct and analyze the most probable, or optimal path of fluctuations from one ordered state to another in real and synthetic networks. The method allows us to compute the distribution of large fluctuations and the time scale associated with switching between ordered states for networks consistent with *mean*-*field* assumptions. In general, we quantify how network heterogeneity influences the scaling patterns and probabilities of fluctuations. For instance, we find that the probability of a large fluctuation near an order-disorder transition decreases exponentially with the *participation ratio* of a network’s principle eigenvector – measuring how many nodes effectively contribute to an ordered state. Finally, the proposed theory is used to answer how and where a network should be targeted in order to optimize the time needed to observe a switch.

Network science is highly interdisciplinary, and when combined with other fields such as statistical physics and nonlinear dynamics, provides a useful framework to address fundamental questions regarding complex systems^{1–3}. Network approaches have provided quantifiable results in diverse applications in nearly every field of science and engineering, from biological networks^{4} to climate networks^{5}, information networks^{6}, infrastructure networks^{7}, and social networks^{8}. Consequently, much progress has been made in understanding the role of topology in many collective processes in complex systems, including in adaptive and co-evolving networks, where the topological dynamics is itself a function of a network’s state^{9}. Popular examples where interaction structure is understood to strongly influence behavior are the spread of infectious diseases^{10}, the dynamics of neural systems^{11}, the synchronization of coupled oscillators^{12}, the patterns of voters^{13}, and the collective motion of networked swarms^{14}.

However, many theoretical results in network dynamics rely on deterministic, and “mean field” limits of some simple model^{2}. Though useful, such approaches typically ignore noise and dynamical fluctuations that are inherent in virtually all of the aforementioned examples. Although inherent noise may be considered small in large networks, the existence and observation of large fluctuations can result in drastic change in a network’s dynamics^{15}. Therefore, some recent efforts have been made to study the relationship between network dynamics and noise^{16–18}. It has been demonstrated that the interplay between complex topology and noise can alter well-known scaling laws and patterns for fluctuations^{19}, as well as provide new control mechanisms that take advantage of noise-induced phenomena including switching and extinction^{15, 20}.

An important class of statistical physics models for capturing many processes, and where noise is relevant, are spin systems, in which nodes in a network take on discrete states^{1, 2}. Typically, the probability (or probability of changing in time) of any configuration of states depends on the configuration’s energy^{21}. Such an approach has been useful for understanding opinion formation and dynamics^{22, 23}, rumor spreading^{24}, as well for machine learning on networks (e. g., Boltzmann machines)^{25} and network inference^{26}. In spin dynamics, often two limits are considered. The first entails zero temperature, where the energy of a network is minimized and generally does not fluctuate between configurations^{22, 27}. The second involves the average state of a network at finite temperatures and noise, where an order-disorder transition is observed and analyzed in the limit where the network size tends to infinity^{23}. However, all real networks are finite systems. Thus, fluctuations between distinct metastable configurations arise – effectively changing the collective order due to noise. Example systems are social networks, where fluctuations are known to be important^{28, 29}, and hence switching from one majority opinion to another is possible. Yet, many open questions remain about how switching occurs in complex networks as a result of random fluctuations.

On the other hand, the role of noise is reasonably well understood in simple well-mixed and spatially homogenous contexts^{30–33}. It has been demonstrated in many works that noise and collective dynamics can couple in such a way as to induce a large fluctuation – effectively driving a system to switch from one collective behavior to another. If the fluctuation is a *rare event*, then the process is captured by a most probable, or optimal path (OP) – where all others are exponentially less likely to occur^{34, 35}. In such cases, the OP is describable in an analytical-mechanics formalism familiar from classical physics^{30}, which has been used to elegantly describe a variety of rare phenomena including: fixation in evolutionary games^{36}, extinction of disease in homogenous populations^{37}, switching in self-regulating genes^{38}, large velocity fluctuations in propagating fronts^{39}, viral clearance^{40}, irreversible fluctuations in electronic circuits^{41}, and switching in quantum mechanical oscillators^{42}. Our strategy is to find the OP in general network configurations given a general spin (opinion) dynamics, and use it to understand the dynamical switching pathway between ordered (majority) states, the average time needed to see a switch, and the distribution of large fluctuations in both real and synthetic networks.

This report moves well beyond both the deterministic and small fluctuation limits of order-disorder dynamics in complex networks by addressing how large fluctuations occur. Our approach enables us to construct and analyze the approximate OP through a given network – reducing a very high dimensional stochastic and rare process to a single trajectory. Consequently, we are able to compute several topologically dependent quantities that have so far eluded analysis in network science: the distribution, shape, and time scale of large fluctuations. Beyond computation, we show that just above an order-disorder transition, the probability of a large fluctuation decreases exponentially with the participation ratio of a network^{43}, and hence is exponentially sensitive to topological heterogeneity– e.g., the ratio of second to fourth moments of a network’s degree distribution. Moreover, we find two quantitatively distinct scaling patterns for the distribution of large fluctuations, in which low eigenvector-centrality fluctuations predominate at high order (large majorities), and high centrality at low order (small majorities). Finally, we demonstrate with several examples on a Facebook network^{44} how the formalism is useful for designing controls that optimally leverage noise in order to minimize the time scale for switching – answering where and at what rate a network should be targeted in order to induce a switch.

We consider a system of *N* nodes interacting through a network. The network is represented by a real-valued matrix, A, where *A*
_{ij} gives the influence strength of node *i* on node *j*. At any instant, each node *i* is in one of two possible opinion (spin) states, characterized by ${s}_{i}\phantom{\rule{-.25em}{0ex}}\in \phantom{\rule{-.25em}{0ex}}\{\phantom{\rule{-.25em}{0ex}}-\phantom{\rule{-.25em}{0ex}}\mathrm{1,}\phantom{\rule{.25em}{0ex}}\mathrm{1\}}$
^{16, 21}. Nodes can change state by interacting with their nearest neighbors in the network, such that *s*
_{i} evolves *stochastically* in time with a probability that depends on *s*
_{i} and *s*
_{j} of neighbors, for *A*
_{ij}≠0. We study a simple interaction rule motivated from statistical physics, where each node has a tendency to align its opinion with neighboring opinions such that the energy, *E*_{i} ≡ − *s*_{i}∑_{j}*A*_{ij}*s*_{j}, is minimized^{2, 21, 45}. We choose the kinetics to be a continuous-time Glauber dynamics– a Markov process with a transition rate for each node:

$$\mathrm{Rate}\phantom{\rule{.25em}{0ex}}({s}_{i}\to -\phantom{\rule{-.25em}{0ex}}{s}_{i})=\frac{\alpha}{1+{e}^{-\phantom{\rule{-.25em}{0ex}}2\beta {E}_{i}}}+\phantom{\rule{-.25em}{0ex}}{f}_{i},$$

1

where *β* is an inverse temperature that measures the ratio of energy to thermal noise, *f*
_{i} is a local spontaneous flipping rate^{16}, and *α* is a rate constant that determines the units of time; we take *α*=1 without loss of generality. Qualitatively, minimizing energy pushes the network toward complete order (unanimous majority or consensus), while thermal noise and spontaneous flipping inject fluctuations that tend to break up order. Because of this interplay, a majority order emerges as long as *β* is above some critical value (see supplementary information section two, SI. ^{2}, for derivation)^{23, 46}. Below, our results are to be compared with the stochastic process defined by Eq. (^{1}).

In general, the network dynamics is governed by a master equation for the *N*-node probability distribution, *ρ*(**s**, *t*), which is high-dimensional and difficult to analyze in its entirety. Therefore, we seek a simplified description of the dynamics by first considering the opinion *density* at each node in a large ensemble. The ensemble consists of *C* identical networks with the same *A*, but *independent* realizations of the dynamics given by Eq. (^{1}). The opinion density in the ensemble is ${m}_{i}\equiv {\sum}_{c=1}^{C}{s}_{i,c}/C$, where *s*
_{i,c} is the state of node *i* in realization *c* of the stochastic process. Ultimately, we are interested in the limit *C*→∞, or continuous density. Our goal is to find an approximate master equation for the network ensemble distribution, *P*(**m**, *t*), that is a function of the densities alone, and extract a particular solution relevant to large fluctuations. As we will see, such a solution will correspond to the OP. To find it, we must consider the transition rates for *m*
_{i} and make a *mean-field* approximation:

$$\begin{array}{rcl}{R}_{i}^{\pm}(\mathbf{m})\equiv \mathrm{R}\mathrm{a}\mathrm{t}\mathrm{e}({m}_{i}\to {m}_{i}\pm 2/C)& =& \sum _{c}\frac{1}{2}(1\mp {s}_{i,c})[\frac{1}{1+{e}^{\mp 2\beta \sum _{j}{A}_{ij}{s}_{j,c}}}+{f}_{i}]\\ & & \approx \phantom{\rule{thinmathspace}{0ex}}\frac{C}{2}(1\mp {m}_{i})[\frac{1}{1+{e}^{\mp 2\beta \sum _{j}{A}_{ij}{m}_{j}}}+{f}_{i}].\end{array}$$

2

The mean-field approximation replaces *s*
_{i,c} by its ensemble average *m*
_{i} – effectively neglecting correlations between neighbors in the network. The result is a master equation that describes a simplified stochastic process in terms of the opinion density at each node,

$$\begin{array}{lll}\hfill \frac{\partial P}{\partial t}(\mathbf{m},t)& \hfill =\hfill & \sum _{i}[{R}_{i}^{+}\left(\mathbf{m}-\frac{2}{C}{\mathbf{1}}_{i}\right)P\left(\mathbf{m}-\frac{2}{C}{\mathbf{1}}_{i},t\right)-{R}_{i}^{+}(\mathbf{m})P(\mathbf{m},t)\hfill \\ \hfill & \hfill & +\phantom{\rule{.25em}{0ex}}{R}_{i}^{-}\left(\mathbf{m}+\frac{2}{C}{\mathbf{1}}_{i}\right)P\left(\mathbf{m}+\frac{2}{C}{\mathbf{1}}_{i},t\right)-{R}_{i}^{-}\left(\mathbf{m}\right)P(\mathbf{m},t)],\hfill \end{array}$$

3

where 1_{i} = 〈0_{1}, 0_{2}, *…*, 0_{i−1}, 1_{i}, 0_{i+1}, *…*〉. In Sec. 3 we analyze Eq. (^{3}) and compare to simulations of Eq. (^{1}) on several networks.

When *β* is above threshold, *P*(**m**, *t*) is peaked around one of two ordered equilibrium states, **m**(*t*)≈±**m***. Dynamically, a finite network fluctuates around equilibrium, **m***, after an initial transient, for a long time until a large fluctuation occurs, which carries the network to the opposite ordered state, −**m**
***, as shown in Fig. 1(a). Such order switches are rare events in large networks (*N* ≫ 1), and we expect them to be encoded in the tails of *P*(**m**, *t*). In particular, if **m** corresponds to a large deviation from **m***, we expect an exponential reduction in probability, as demonstrated in Fig. 1(b). Therefore, it is convenient to constrain our search for solutions of the master equation to the exponentially-distributed tail that is relevant for rare events, since Eq. (^{3}) contains too much information to be useful in practice.

(**a**) Switching in a Facebook network^{44} between meta-stable ordered states. Average opinion, weighted by eigenvector centrality (*η*
_{i} for node *i*), is shown versus time. Parameters are *βλ*=1.37 and *f*
_{i}=0.02 **...**

These observations suggest extracting a solution of Eq. (^{3}) with exponential, or Wentzel–Kramers–Brillouin (*WKB*) form, *P*(**m**, *t*)=*ae*
^{−CS(m,t)}. The WKB solution for the ensemble distribution can be viewed as a product of independent and identical distributions for each realization in the ensemble. Thus, we can approximate the probability distribution for states in a single realization, *ρ*(**s**, *t*), by

$$\rho (\mathbf{s},t)\cong \rho (\mathbf{m},t)=b{e}^{-\phantom{\rule{-.10em}{0ex}}\mathbf{S}(\mathbf{m},t)}\mathrm{.}$$

4

Predictions from Eq. (^{4}) are in good agreement with simulations on a real Facebook network^{44} – shown in red in Fig. 1(b).

We can find the continuous-density solution for *P*(**m**, *t*) in the large-ensemble limit by substituting the WKB ansatz into Eq. (^{3}), Taylor expanding Eq. (^{3}) in powers of the small parameter 1/*C*, and neglecting terms of $\mathcal{O}\mathrm{(1/}C)$ or smaller (see SI. ^{1} for expansion). As is customary, this approximation converts the master equation into a familiar Hamilton-Jacobi equation (HJE) from analytical mechanics^{30, 33}:

$$\frac{\partial S}{\partial t}+H(\mathbf{m},\partial S/\partial \mathbf{m}\mathrm{)=0,}$$

5

where *S* and *H* are called the Action and Hamiltonian, respectively^{33}. Once *S* is found from Eq. (^{5}), the distribution of large fluctuations is determined. Following analytical mechanics, we define a momentum, *p*_{i} = ∂*S*/∂*m*_{i}, and conveniently express the Hamiltonian,

$$\begin{array}{l}\phantom{\rule{-.25em}{0ex}}\phantom{\rule{-.25em}{0ex}}H(\mathbf{m},\mathbf{p})=\sum _{i}\phantom{\rule{0.25em}{0ex}}[\frac{1}{2}\mathrm{(1}-{m}_{i}\left)\right({e}^{2{p}_{i}}-\mathrm{1)}\left(\frac{1}{1+{e}^{-2\beta \phantom{\rule{-.25em}{0ex}}\sum _{j}\phantom{\rule{-.25em}{0ex}}{A}_{ij}{m}_{j}}}+{f}_{i}\right)\\ \phantom{\rule{1em}{0ex}}\phantom{\rule{1em}{0ex}}\phantom{\rule{1em}{0ex}}\phantom{\rule{1em}{0ex}}+\phantom{\rule{.25em}{0ex}}\frac{1}{2}\mathrm{(1}+{m}_{i}\left)\right({e}^{-2{p}_{i}}-\mathrm{1)}\left(\frac{1}{1+{e}^{2\beta \sum _{j}\phantom{\rule{-.25em}{0ex}}{A}_{ij}{m}_{j}}}+{f}_{i}\right)]\mathrm{.}\end{array}$$

6

A crucial result of the WKB approximation, which makes the approach worthwhile, is that solutions of the HJE extremize *S*, when expressed as the integral^{30}

$$S(\mathbf{m},t)={\int}_{{t}_{0}}^{t}[\mathbf{p}\cdot \frac{d\mathbf{m}}{d{t}^{\mathrm{\prime}}}-H(\mathbf{m},\mathbf{p})]d{t}^{\mathrm{\prime}}={\int}_{\mathbf{m}(t={t}_{0})}^{\mathbf{m}}\phantom{\rule{thinmathspace}{0ex}}\mathbf{p}\cdot d\mathbf{m}-{\int}_{{t}_{0}}^{t}H(\mathbf{m},\mathbf{p})d{t}^{\mathrm{\prime}},$$

7

where **m**(*t*) and **p**(*t*) are determined from Hamilton’s equations of motion below, Eqs (^{8}–^{9}). Since *S* is minimized, the probability of a stochastic path associated with *S* is maximized. In summary, when the transition between ordered states in a network is exponentially rare, there is a *least-action* path in (**m**, **p**) phase-space that connects the two, which is a local maximum in probability, and therefore corresponds to the OP through a network. The Action along the OP gives *ρ*(**m**, *t*) from Eqs (^{4}) and (^{7}).

Finally, just as in analytical mechanics, a convenient approach for computing the OP is to solve Hamilton’s equations of motion for the system, $\partial H/\partial {p}_{i}={\dot{m}}_{i}$ and $\partial H/\partial {m}_{i}=-\phantom{\rule{-.25em}{0ex}}{\dot{p}}_{i}$:

$${\dot{m}}_{i}=\frac{\mathrm{(1}\phantom{\rule{-.25em}{0ex}}-{m}_{i}){e}^{2{p}_{i}}}{1\phantom{\rule{-.25em}{0ex}}+{e}^{-2\beta \sum _{j}{A}_{ij}{m}_{j}}}-\frac{\mathrm{(1}\phantom{\rule{-.25em}{0ex}}+{m}_{i}){e}^{-2{p}_{i}}}{1\phantom{\rule{-.25em}{0ex}}+{e}^{2\beta \sum _{j}{A}_{ij}{m}_{j}}}+{f}_{i}\mathrm{[(1}-{m}_{i}){e}^{2{p}_{i}}-\mathrm{(1}+{m}_{i}){e}^{-2{p}_{i}}],$$

8

$$\begin{array}{c}{\dot{p}}_{i}=\frac{\frac{1}{2}({e}^{2{p}_{i}}-1)}{1+{e}^{-2\beta \sum _{j}{A}_{ij}{m}_{j}}}-\frac{\frac{1}{2}({e}^{-2{p}_{i}}-1)}{1+{e}^{2\beta \sum _{j}{A}_{ij}{m}_{j}}}-\beta \sum _{j}\phantom{\rule{thickmathspace}{0ex}}{A}_{ji}[\frac{(1-{m}_{j})({e}^{2{p}_{j}}-1)-(1+{m}_{j})({e}^{-2{p}_{j}}-1)}{{({e}^{\beta \sum _{k}{A}_{jk}{m}_{k}}+{e}^{-\beta \sum _{k}{A}_{jk}{m}_{k}})}^{2}}]\\ \phantom{\rule{2em}{0ex}}+\phantom{\rule{thinmathspace}{0ex}}\frac{{f}_{i}}{2}[{e}^{2{p}_{i}}-{e}^{-2{p}_{i}}].\end{array}$$

9

It is important to notice that if one takes **p** ≡ 0 in Eqs (^{8}–^{9}) the “quenched mean field” equations are derived for a kinetic Ising model as a special case which ignores fluctuations^{46, 47}. Similar findings have been shown recently for epidemic dynamics^{48}. Therefore, the OP formalism naturally generalizes deterministic approaches for network dynamics to include large fluctuations. In practice, all that is needed to find the OP are appropriate boundary conditions for Eqs (^{8}–^{9}). These are derivable directly from the distribution as explained in Sec. 3.2.

By considering histograms of **m**(*t*) from time-series data, e.g., Fig. 1(a), we can see that the expected exponential distribution of large fluctuations appears. Moreover, simple inspection gives us the boundary conditions for solving Eqs (^{8}–^{9}). Since the distribution takes a maximum value at the equilibrium **m*** (satisfying $\dot{\mathbf{m}}=\dot{\mathbf{p}}=\mathbf{0}$), the initial boundary condition for a large fluctuation is **m**(*t*=0)=**m*** and $\partial S/\partial \mathbf{m}=\mathbf{0}=\mathbf{p}(t=\mathrm{0)}$. Because we are interested in large fluctuations that lead to a switch, the final boundary condition is similarly **m**(*t*→∞)=−**m**
*** and **p**(*t*→∞)=**0**
^{30, 33, 49}. Therefore, by solving Eqs (^{8}–^{9}) subject to *zero-momentum boundary conditions*, we determine the OP: a single trajectory that gives the probability of a large fluctuation to **m** within logarithmic accuracy, ln *ρ*(**m**) ≈ − *S*(**m**). We note that since the distribution is nearly constant in time, the network Action is time-independent ∂*S*/∂*t* = 0 = *H* ∀*t*. Therefore, *S*(**m**) is equal to the line integral of the momentum along the OP, from Eq. (^{7}).

Optimal paths can be computed numerically from Eqs (^{8} and ^{9}) with zero-momentum boundary conditions using quasi-Newton methods^{49} and dimension-reduction techniques based on a spectral decomposition of *A* (see SI. ^{6–7} for details on numerical approaches). In general, we find two distinct segments of the OP: an activation segment with **p**≤0, requiring noise to carry the network from **m**=**m*** to **m**=**0** (shown in Fig. 3), and a deterministic segment with **p**=**0** that leads from **m**=0 to **m**=−**m**
***. Since the deterministic segment has *p*=**0**, the probability distribution is predicted to be flat from **m**=0 to **m**=−**m***, as demonstrated in Fig. 1(b). In addition to the distribution, the shapes of large fluctuations are predicted and can be compared to simulations. For example, Fig. 2 shows OP projections and prehistory trajectory-density plots from many stochastic realizations of Eq. (^{1}) for two networks. We see that the OP for each network lies near the maximum of the corresponding heat maps. Importantly, since the OP predicts the most likely sequence of changes in opinion density, **m**-projections at various points along a path (such as in Fig. 2) give us insight into how and where noise acts in a network during a large fluctuation.

Momentum versus opinion density for several Facebook network eigenvector centralities. The centralities are (from red to black) *η*
_{i}=0.0058, 0.0384, 0.0582, 0.0685, 0.0748, 0.0802 and 0.0864. The shaded area illustrates the contribution **...**

Prehistory heat maps showing the density of the final *N* stochastic events in 400 stochastic simulations that resulted in switching for two networks. Trajectories are projected into two network bins, denoted by subscripts *b* and *b*′. (**a**) Network **...**

Beyond finding a prescription for computing the OP numerically, we are interested in understanding how the path structure is related to topological properties of a network. The OP dependencies can be found analytically in certain limiting cases of interest, such as near threshold. A general method used throughout Sec. 3.2 is to expand Eqs (^{8} and ^{9}) around the OP boundaries – giving the local functional forms of **m**, **p**, and *S*(**m**). We show in SI. ^{4} that when **f**=**0** and the distance to threshold is small, *δ* ≡ *β**λ* − 1 ≳ 0 (where *λ* is the largest eigenvalue of *A*), the OP depends on the principle right and left eigenvectors of *A*, ** η** and

$$\mathbf{m}(h)=\mathit{\eta}{\delta}^{\mathrm{1/2}}\mathit{h}\sqrt{\mathrm{3/}\sum _{j}{\zeta}_{j}{\eta}_{j}^{3}},$$

10

$$\mathbf{p}(h)=\mathit{\zeta}{\delta}^{\mathrm{3/2}}h(h-\mathrm{1)(}h+\mathrm{1)}\sqrt{\mathrm{3/}\sum _{j}{\zeta}_{j}{\eta}_{j}^{3}}\mathrm{.}$$

11

*h* is a unit-length parameter along the activation segment, $h\equiv {m}_{i}/{m}_{i}^{\u204e}\forall \phantom{\rule{-.25em}{0ex}}i$. In general, the principle right eigenvector satisfies ** η**=

In general, since each node’s contribution to the Action is equal to the line integral of the momentum, we expect the contribution to increase with increasing eigenvector centrality. This pattern can be seen in Fig. 3, where the area under the darker curves (corresponding to higher centrality) contains the area under the lighter curves. Near threshold, we can calculate the line integral explicitly by substituting Eqs (^{10} and ^{11}) into Eq. (^{7}): $S(\mathbf{m})={\int}_{\mathbf{m}\u204e}^{\mathbf{m}}\mathbf{p}\cdot \mathbf{d}\mathbf{m}\mathbf{\prime}\approx {\sum}_{i\mathrm{=1}}^{i=N}\phantom{\rule{-.25em}{0ex}}{\int}_{1}^{h}{p}_{i}(h\prime )[d{m}_{i}/dh\prime ]dh\prime $ or

$$S(\mathbf{m}(h))=\frac{3{\delta}^{2}}{4\sum _{j}{\zeta}_{j}{\eta}_{j}^{3}}{\mathrm{(1}-{h}^{2})}^{2}\mathrm{.}$$

12

The expression is interesting, since for a symmetric network we find that *S* is a function of the fourth moment of the eigenvector-centrality distribution, or inverse participation ratio^{43}. To understand what this implies we first consider the case of a homogeneous complete graph, where ${\eta}_{i}\phantom{\rule{-.25em}{0ex}}=\phantom{\rule{-.25em}{0ex}}{\zeta}_{i}\phantom{\rule{-.25em}{0ex}}=\phantom{\rule{-.25em}{0ex}}\mathrm{1/}\sqrt{N}\mathrm{.}$ In this case the probability of a large fluctuation to zero consensus, and therefore the probability of switching, is to logarithmic accuracy $\mathrm{ln}\phantom{\rule{.10em}{0ex}}\rho (\mathbf{0})=-\phantom{\rule{-.25em}{0ex}}3N{\delta}^{2}\mathrm{/4}$ (for *N**δ*^{2} ≫ 1). On the other hand, for a random network without degree correlations the standard annealed-network approximation gives ${\eta}_{i}={\zeta}_{i}={k}_{i}/\sqrt{N\langle {k}^{2}\rangle},$ given a degree *k*
_{i} for node *i* and a network average of *k*
^{2}, *k*
^{2}. Hence, $\mathrm{ln}\phantom{\rule{.10em}{0ex}}\rho (\mathbf{0})=-\phantom{\rule{-.25em}{0ex}}3N{\delta}^{2}{\langle {k}^{2}\rangle}^{2}\mathrm{/4}\langle {k}^{4}\rangle $. If the the network is composed of nodes with only degree *k*, the switching probability is predicted to equal the complete graph’s value. However for a heterogeneous network, e.g., when the distribution of nodes with degree *k* scales like *k*
^{−γ}, the additional topological factor, 〈*k*^{2}〉^{2}/〈*k*^{4}〉, can be significantly smaller than 1 depending on the cutoff of the distribution. Generally then, given the same distance to threshold and system size the switching probability is predicted to increase exponentially with increasing topological heterogeneity^{15}. This will be particularly the case if *γ*≤5, which can be contrasted to other kinetic systems such as epidemics^{48}.

Near threshold a network is highly disordered, and we may wonder how large fluctuations behave as the order increases with *β* far above threshold. In fact, the scaling patterns for the very largest fluctuations to small **m**≈**0** maintain a similar form. In SI. ^{5} we show that *m*
_{i} and *p*
_{i} remain proportional to centrality in the tail of the distribution by expanding Eqs (^{8} and ^{9}) around the origin (assuming *f*
_{i}=*f* and *A* is symmetric). Hence, the relative probabilities for observing small ordering are simple functions of a network’s average-squared opinion density, $\langle {m}^{2}\rangle =(\mathbf{m}\cdot \mathbf{m})/N\gtrsim \mathrm{0:}$

$$\rho (\mathbf{m})/\rho (\mathbf{m}{}^{\phantom{\rule{-0.15em}{0ex}}\mathbf{\prime}})=\mathrm{exp}\phantom{\rule{thinmathspace}{0ex}}[\frac{N}{2}(\frac{\beta \lambda -1-2f}{1+2f})(\u27e8{m}^{2}\u27e9-\u27e8{{m}^{\mathrm{\prime}}}^{2}\u27e9)].$$

13

The generic scaling with *m*
^{2} occurs because *p*
_{i} tends to the same slope with respect to *m*
_{i} along the OP for all nodes when **m**≈**0**; the scaling is demonstrated in Fig. 3 with a black-dashed line. Notably Eq. (^{13}) implies that the tail of *ρ*(**m**) is an exponential with a rate that increases *linearly* with the network size and distance to threshold, *βλ*−1−2*f*, but decreases with *f*. The latter effect entails an additional broadening of the distribution in the presence of spontaneous flipping that is independent of network topology.

In contrast, the scaling at high order is significantly different, with fluctuations that are very sensitive to node centrality. In particular, when a network is near consensus, **m** ≈ **m**^{⁎} ≲ 1, lower-centrality nodes have fluctuations that are *exponentially* larger than higher-centrality nodes along the OP, $[{m}_{i}-{m}_{i}^{\u204e}]/\phantom{\rule{.25em}{0ex}}[{m}_{j}-{m}_{j}^{\u204e}]\sim [{\eta}_{i}/{\eta}_{j}]{e}^{2\beta \lambda {\sum}_{l}\phantom{\rule{-.25em}{0ex}}{\eta}_{l}{m}_{l}^{\u204e}[{\eta}_{j}-{\eta}_{i}]},$ and are distributed according to a simple Gaussian:

$$\rho (\mathbf{m})=b\phantom{\rule{0.15em}{0ex}}\mathrm{exp}\phantom{\rule{0.15em}{0ex}}[\sum _{i}-\frac{1}{8}{({m}_{i}-{m}_{i}^{\ast})}^{2}{e}^{2\beta \lambda {\eta}_{i}{\sum}_{l}{\eta}_{l}{m}_{l}^{\ast}}].$$

14

Intuitively, when a network is near consensus, it is very improbable for high centrality nodes to change state. This pattern can be seen in Fig. 3 where the slope of *p*
_{i} with respect to *m*
_{i} along the OP is much steeper for higher centrality nodes when **m**≈**m***– corresponding to a more rapid decrease in probability from the maximum value. The quantitative scaling for high order is depicted in Fig. 3 with a blue-dashed line, and reflects the pattern that the standard deviation of each node’s distribution is exponentially decreasing with its centrality near consensus. Equation (^{14}) results from an expansion of Eqs (^{8} and ^{9}) near consensus, given in SI. ^{5}, and assumes *f*≈0 and *A*≈*λ*
*ηη*^{T}.

The final observable we consider in this report is the average switching time, *T*, since it quantifies the expected time scale over which the very largest fluctuations occur. Generally *T* takes the form:

$$\langle T\rangle =\tau (\beta ,f,A){e}^{S\mathrm{(0)}},$$

15

from the assumption that switching between ordered states has a rate, or inverse time, proportional to the probability, *ρ*(0)
^{33, 34}. For sufficiently large *S*, the exponential contribution dominates, and therefore ln*T*~*S*(**0**), as demonstrated in Fig. 4(a) for several networks. Given the dependence on *ρ*(0), our analysis indicates that switching times have exponential sensitivity to the topological properties which determine the Action, such as network heterogeneity.

Since we are able to predict large fluctuations with the OP theory, we may be interested in adding controls to the network dynamics in order to optimize an observable– e.g., minimize *T*. In this report, our approach will involve minimizing *S*(**0**) (a deterministic quantity derived from theory). However, it is important to note that the formalism ultimately entails using internal fluctuations in the network to do work, which would not be possible in the absence of noise^{15}. In this section we study the example of targeting a subset, *F* of nodes with spontaneous flipping, for which *f*
_{i}=*f*≠0 given *i**F*, while all other nodes have rate equal to zero. We consider several cases using the Facebook network^{44} as an illustration, and answer the following questions: what nodes should be targeted in order to minimize *T* (given a fixed number for targeting, |*F*|, and a fixed flipping rate, *f*), and whether targeting a larger subset of nodes at a lower rate, or a smaller set with a higher rate tends to minimize *T*.

Following the analysis above, we differentiate among nodes by eigenvector centrality, *η*
_{i} for node *i*. In order to first compare equal numbers of nodes for targeting, we have arranged nodes in a list according to increasing *η*
_{i}, and binned the list into roughly equally sized bins. In Fig. 4(b) ln *T* is shown as a function of the average *η*
_{i} in *F*, *η**F*, for several *f*. The size of *F* is fixed at 32 nodes. The *F* with highest *η*_{F} in Fig. 4(b) represents nodes near the maximum value of *η**F* in the network: i.e., the 32 nodes with highest centrality. The second highest *η**F* represents nodes with the next highest centrality, and so on (see SI. ^{7–8} for more details). Since the times decrease with *η*_{F}, we can see that it is optimal to target nodes with higher centrality. The finding makes intuitive sense, since high *η*
_{i} implies high importance for a given node^{3}, and therefore one might expect increased effect from flipping high-centrality nodes. Nevertheless, the scale of difference is important. Because of the exponential form of *T*, targeting nodes with the highest centrality can reduce ln*T* by 25% for the parameters shown, even though <1% of the network is controlled. Predictions are in good agreement with simulations, Fig. 4(b).

On the other hand, with a different control scheme targeting higher *η*
_{i} alone may be sub-optimal. Another approach is to start with the 32 highest-centrality nodes (i.e., the control with the smallest ln*T* in Fig. 4(b), and increase/decrease the size of *F* by adding/removing nodes with lower *η*
_{i}. An example is shown in Fig. 4(c) (blue circles) where |*F*| nodes with the highest *η*
_{i} are targeted. In order to keep the total rate of flipping constant, *f*|*F*| – a proxy for the amount of work done by the controller, *f* must vary accordingly. We can see that for the Facebook network it is more optimal to target a larger set of nodes at a slower rate, than fewer nodes at a higher rate, even though a larger set implies decreasing *η*_{F}. Similar to Fig. 4(b), a reduction by a factor of 25% in ln*T* is predicted and observed, despite targeting at most <2% of the network. Finally, it is interesting that minimizing *S*(**0**), with the controls considered, is strongly correlated with reducing the amount of order in the network. To demonstrate, we redo the second control by picking an *f* so that *m*
^{*2} is constant as we vary *F*. The result is a slowly varying *T*, illustrated in Fig. 4(c) (green diamonds).

There is much interest in understanding the relationship between dynamics in complex systems and their underlying topology. However, most theoretical results pertain to deterministic limits, where noise is ignored. Therefore, analytical and computational tools are needed to understand how noise and dynamics interact, especially when the interplay causes a large qualitative change in a network’s collective dynamics. In this report, we have developed an approach based on WKB techniques applied to finite networks, which allowed us to analyze large fluctuations in order-disorder transitions driven by internal noise. The work went well beyond steady state and threshold analysis – concerning itself with a *global* dynamical object, the optimal path of large fluctuations through a network. By computing the optimal path we were able to predict the probabilities and time scales associated with large fluctuations, including the largest fluctuations that entailed a noise-induced switch between deterministically stable states. Simulations in both real and synthetic networks showed good agreement with theory.

The optimal path approach applied to networks has several advantages, which were detailed in this report. First, it provides a way to predict large fluctuations that naturally reproduces popular mean-field results as a special *zero-momentum limit*. Because of this structure, we expect more accurate network-approximation techniques, such as those that include dynamical correlations, to be generalized in a similar way – and for many other dynamical processes in networks. For instance, WKB techniques can be used to analyze large fluctuations in nonlinear systems with time delays^{50}, colored noise^{51}, and memory effects^{52}. Second, it allows one to analytically quantify the scaling patterns of large fluctuations on topology, including the exponential sensitivity of fluctuation-probabilities to topological heterogeneity and the multi-step structure of large fluctuations from highly ordered states through heterogeneous networks. Third, optimal control of noisy network dynamics is reduced to deterministic control of a mechanical analog. By minimizing the network *Action* derived from the optimal path theory, one can construct the optimal way to leverage internal noise so as to maximize the probability of network switching. This was demonstrated in a Facebook network where minimum-Action controls correctly predicted large exponential reductions in the average switching time, by targeting an optimal subset of the network that represented less than two percent of the total. Fourth, we expect that our approach will be useful for current avenues of research in network science and many new applications, such as network inference from data in the presence of large fluctuations.

J.H. is a National Research Council postdoctoral fellows. I.B.S. was supported by the U.S. Naval Research Laboratory funding (N0001414WX00023) and office of Naval Research (N0001416WX00657) and (N0001416WX01643). We are very grateful to L.B. Shaw and L. Mier-y-Teran-Romero for useful discussions.

Author Contributions

Both J.H. and I.B.S. contributed to the methods and analysis, and wrote the manuscript.

The authors declare that they have no competing interests.

**Electronic supplementary material**

**Supplementary information** accompanies this paper at doi:10.1038/s41598-017-08828-8

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

1. Dorogovtsev SN, Goltsev AV, Mendes JFF. Critical phenomena in complex networks. Rev. Mod. Phys. 2008;80 doi: 10.1103/RevModPhys.80.1275. [Cross Ref]

2. Barrat, A., Barthélemy, M. & Vespignani, A. *Dynamical Processes on Complex Networks* (Cambridge University Press, 2008).

3. Newman, M. E. J. *Networks*: *An Introduction* (Oxford University Press, 2010).

4. Proulx SR, Promislow DEL, Phillips PC. Network thinking in ecology and evolution. Trends. Ecol. Evolut. 2005;20 doi: 10.1016/j.tree.2005.04.004. [PubMed] [Cross Ref]

5. Berezin Y, Gozolchiani A, Guez O, Halvin S. Stability of climate networks with time. Sci. Rep. 2012;2 doi: 10.1038/srep00666. [PMC free article] [PubMed] [Cross Ref]

6. Timár G, Goltsev AV, Dorogovtsev SN, Mendes JFF. Mapping the structure of directed networks: Beyond the bow-tie diagram. Phys. Rev. Letts. 2017;118 doi: 10.1103/PhysRevLett.118.078301. [PubMed] [Cross Ref]

7. Bagler G. Analysis of the airport network of india as a complex weighted network. Physica A. 2008;387 doi: 10.1016/j.physa.2008.01.077. [Cross Ref]

8. Weng L, Menczer F, Ahn Y-Y. Virality prediction and community structure in social networks. Sci. Rep. 2013;3 doi: 10.1038/srep02522. [PMC free article] [PubMed] [Cross Ref]

9. Gross, T. & Sayama, H. *Adaptive networks* (Springer, 2009).

10. Pastor-Satorras R, Castellano C, Van Mieghem P, Vespignani A. Epidemic processes in complex networks. Rev. Mod. Phys. 2015;87 doi: 10.1103/RevModPhys.87.925. [Cross Ref]

11. Yan Y, et al. Nonequilibrium landscape theory of neural networks. Proc. Natl. Acad. Sci. USA. 2013;110 doi: 10.1073/pnas.1310692110. [PubMed] [Cross Ref]

12. Pecora LM, Sorrentino F, Hagerstrom AM, Murphy TE, Roy R. Cluster synchronization and isolated desynchronization in complex networks with symmetries. Nat. Commun. 2014;5 doi: 10.1038/ncomms5079. [PubMed] [Cross Ref]

13. Fernández-Gracia J, Suchecki K, Ramasco JJ, Miguel MS, Eguíluz VM. Is the voter model a model for voters? Phys. Rev. Letts. 2014;112 doi: 10.1103/PhysRevLett.112.158701. [PubMed] [Cross Ref]

14. Hindes J, Szwaykowska K, Schwartz IB. Hybrid dynamics in delay-coupled swarms with mothership networks. Phys. Rev. E. 2016;94 doi: 10.1103/PhysRevE.94.032306. [PubMed] [Cross Ref]

15. Hindes J, Schwartz IB. Epidemic extinction and control in heterogeneous networks. Phys. Rev. Letts. 2016;117 doi: 10.1103/PhysRevLett.117.028302. [PubMed] [Cross Ref]

16. Carro A, Toral R. & Miguel, M. S. The noisy voter model on complex networks. Sci Rep. 2016;6 doi: 10.1038/srep24775. [PMC free article] [PubMed] [Cross Ref]

17. Ching ESC, Tam HC. Reconstructing links in directed networks from noisy dynamics. Phys. Rev. E. 2017;95 doi: 10.1103/PhysRevE.95.010301. [PubMed] [Cross Ref]

18. Böttcher L, Lukovi M, Nagler J, Havlin S, Herrmann HJ. Failure and recovery in dynamical networks. Sci. Rep. 2017;7 doi: 10.1038/srep41729. [PMC free article] [PubMed] [Cross Ref]

19. Assaf M, Mobilia M. Metastability and anomalous fixation in evolutionary games on scale-free networks. Phys. Rev. Letts. 2012;109 doi: 10.1103/PhysRevLett.109.188701. [PubMed] [Cross Ref]

20. Wells DK, Kath WL, Motter AE. Control of stochastic and induced switching in biophysical networks. Phys. Rev. X. 2015;5 [PMC free article] [PubMed]

21. Castellano C, Fortunato S, Loreto V. Statistical physics of social dynamics. Rev. Mod. Phys. 2009;81 doi: 10.1103/RevModPhys.81.591. [Cross Ref]

22. Castellano C, Pastor-Satorras R. Zero temperature glauber dynamics on complex networks. J. Stat. Mech. Theor. Exp. 2006;5

23. Dorogovtsev SN, Goltsev AV, Mendes JFF. Ising model on networks with an arbitrary distribution of connections. Phys. Rev. E. 2002;66 doi: 10.1103/PhysRevE.66.016104. [PubMed] [Cross Ref]

24. Ostilli M, et al. Statistical mechanics of rumour spreading in network communities. Procedia Comput. Sci. 2010;1 doi: 10.1016/j.procs.2010.04.262. [Cross Ref]

25. Tanaka T. Mean-field theory of boltzmann machine learning. Phys. Rev. E. 1998;58 doi: 10.1103/PhysRevE.58.2302. [Cross Ref]

26. Zeng H-L, Aurell E, Alava M, Mahmoudi H. Network inference using asynchronously updated kinetic ising model. Phys. Rev. E. 2010;83 doi: 10.1103/PhysRevE.83.041135. [PubMed] [Cross Ref]

27. Schneider-Mizell CM, Sander LM. A generalized voter model on complex networks. J. Stat. Phys. 2009;136 doi: 10.1007/s10955-009-9757-6. [Cross Ref]

28. Xiong F, Lui Y. Opinion formation on social media: An empirical apprach. Chaos. 2014;24 doi: 10.1063/1.4866011. [PubMed] [Cross Ref]

29. Acemoglu D, Como G, Fagnani F, Ozdaglar A. Opinion fluctuations and disagreement in social networks. Math. Oper. Res. 2013;38 doi: 10.1287/moor.1120.0570. [Cross Ref]

30. Dykman MI, Mori E, Ross J, Hunt PM. Large fluctuations and optimal paths in chemical kinetics. J. Chem. Phys. 1994;100 doi: 10.1063/1.467139. [Cross Ref]

31. Ovaskainen O, Meerson B. Stochastic models of population extinction. Trends Ecol. Evol. 2010;25 doi: 10.1016/j.tree.2010.07.009. [PubMed] [Cross Ref]

32. Kamenev A, Meerson B, Shklovskii B. How colored environmental noise affects population extinction. Phys. Rev. Letts. 2008;101 doi: 10.1103/PhysRevLett.101.268103. [PubMed] [Cross Ref]

33. Assaf M, Meerson B. Wkb theory of large deviations in stochastic populations. J. Phys. A: Math. Theor. 2017;50 doi: 10.1088/1751-8121/aa669a. [Cross Ref]

34. Assaf M, Meerson B. Extinction of metastable stochastic populations. Phys. Rev. E. 2010;81 doi: 10.1103/PhysRevE.81.021116. [PubMed] [Cross Ref]

35. Lindley BS, Shaw LB, Schwartz LB. Rare event extinction on stochastic networks. Europhys. Lett. 2014;108 doi: 10.1209/0295-5075/108/58008. [Cross Ref]

36. Mobilia M, Assaf M. Fixation in evolutionary games under non-vanishing selection. Europhys. Lett. 2010;91 doi: 10.1209/0295-5075/91/10002. [Cross Ref]

37. Schwartz IB, Billings L, Dykman M, Landsman AS. Predicting extinction rates in stochastic epidemic models. J. Stat. Mech.: Theory Exp. 2010;9

38. Assaf M, Roberts E, Luthey-Schulten Z. Determining the stability of genetic switches: explicitly accounting for mrna noise. Phys Rev Lett. 2011;106 doi: 10.1103/PhysRevLett.106.248102. [PubMed] [Cross Ref]

39. Meerson B, Sasorov PV. Emergence of fluctuating traveling front solutions in macroscopic theory of noisy invasion fronts. Phys. Rev. E. 2011;84 doi: 10.1103/PhysRevE.84.030101. [PubMed] [Cross Ref]

40. Chaudhury S, Perelson AS, Sinitstyn NA. Spontaneous clearance of viral infections by mesoscopic fluctuations. PLOS ONE. 2012;7 doi: 10.1371/journal.pone.0038549. [PMC free article] [PubMed] [Cross Ref]

41. Luchinsky DG, McClintock PVE. Irreversibility of classical fluctuations studied in analougue electronic circuits. Nature. 1997;389 doi: 10.1038/38963. [Cross Ref]

42. Lin ZR, Nakamura Y, Dykman MI. Critical fluctuations and the rates of interstate switching near the excitation threshold of a quantum parametric oscillator. Phys. Rev. E. 2015;92 doi: 10.1103/PhysRevE.92.022105. [PubMed] [Cross Ref]

43. Pastor-Satorras R, Castellano C. Distinct types of eigenvector localization in networks. Sci. Rep. 2016;6 doi: 10.1038/srep18847. [PMC free article] [PubMed] [Cross Ref]

44. Mcauley, J. & Leskovec, J. Learning to discover social circles in ego networks. *NIPS* (2012).

45. Ben-Naim, E., Krapivsky, P. L. & Redner, S. *A Kinetic View of Statistical Physics* (Cambridge University Press, 2010).

46. Lynn, C. W. & Lee, D. D. Maximizing influence in an ising network: A mean-field optimal solution. *NIPS* (2013).

47. Goltsev AV, Dorogovtsev SN, Oliveiram JG, Mendes JFF. Localization and spreading of diseases in complex networks. Phys. Rev. Letts. 2012;109 doi: 10.1103/PhysRevLett.109.128702. [PubMed] [Cross Ref]

48. Hindes J, Schwartz IB. Epidemic extinction paths in complex networks. Phys. Rev. E. 2017;95 doi: 10.1103/PhysRevE.95.052317. [PubMed] [Cross Ref]

49. Lindley BS, Schwartz IB. An iterative action minimizing method for computing optimal paths in stochastic dynamical systems. Physica D. 2013;255 doi: 10.1016/j.physd.2013.04.001. [Cross Ref]

50. Schwartz IB, Billings L, Carr TW, Dykman MI. Noise-induced switching and extinction in systems with delay. Phys. Rev. E. 2015;91 doi: 10.1103/PhysRevE.91.012139. [PubMed] [Cross Ref]

Articles from Scientific Reports are provided here courtesy of **Nature Publishing Group**

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's national science library in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |