Faster network disruption from layered oscillatory dynamics

Various systems in nature and engineered applications can be modeled as individual dynamical units, interacting with one another in a complex way such as neurons spiking together in the brain² or rotating masses of power generators evolving at the same frequency in large-scale transmission networks.³ Collective phenomena taking place in these systems such as synchronization are enabled by both the internal dynamics and the interaction within the individual units.⁴ Another remarkable feature of such systems is multistability. Thanks to nonlinearities in the interaction, numerous equilibria might exist, each of them having their own basin of attraction.^5,6 Due to external perturbations such as environmental noise or faults occurring at some components, the whole system might undergo transitions between equilibrium points. In many instances, this is not desirable as both the transient dynamics and the new equilibrium might threaten the correct operation of the system or even cause damages.⁷ An important task is to predict such a transition as they could disrupt the desired functioning state of the networked system. The latter question has been investigated mostly in coupled systems made of a single layer of interaction, usually subjected to spatially uncorrelated noise.^8–11 However, many applications require more involved coupling structures in order to correctly describe them.¹² A natural extension of the single layer framework is to include multiple layers of dynamical systems.^13,14 In such layered dynamics, noise acting on one sub-network propagates to other layers with system-specific correlations.¹ In other words, when going through a layer, the noise acquires some specific statistical properties, which is then injected into other layers. Including such a system-specific correlation in the noise might increase the level of fluctuations and, thus, the risk of transitions between equilibria. This is depicted in Fig. 1, where uncorrelated noise is acting on the first layer [blue in panel (a)], which is then transmitted to the second one [red in panel (a)]. In this setting, the second layer is subjected to correlated noise, which we describe below. Both layers are made of a single cycle network, which has multiple equilibria as illustrated in Fig. 1(b). Due to the noise injected in both of them, the two layers are pushed away from their initial equilibrium and, depending on the noise amplitude and direction, eventually leave their initial basin of attraction. Transitions between basins of attraction can be detected by changes in the winding numbers in each layer $q1$ and $q2$ [see Fig. 1(c)]. One clearly sees that, while the first layer does not operate any transition, the second one exits its initial basin of attraction multiple times. Similar transitions were observed for a single layer system^8,10 and related to the eigenvalues of the network Laplacian matrix and the distance between the initial fixed point and the closest saddle point with a single unstable direction.¹¹ Here, I consider layered systems where each layer has its own individual units, nonlinearly coupled on a complex network. The different layers interact together via a coupling function. It was recently shown that noise acting on a layer might be strongly amplified in the other layers depending on the network connectivities.¹ In addition to being amplified, the noise structure seems to align with the lowest-lying eigenmodes of the network Laplacian. Building on these results, I investigate the first escape time from the initial basin of attraction in one layer due to noise originating from connected layers. Remarkably, correlations in space and time of this type of noise seem more efficient in driving the system through a transition to another basin than uncorrelated noise.

This paper is organized as follows. In Sec. II, I introduce the two-layer system of nonlinearly coupled oscillators considered, recall previous results,¹ and discuss escape from the basin of attraction. Numerical results are presented in Sec. III and conclusions in Sec. IV.

In order to investigate the effect of noise in layered networks, I consider a simple two-layer system made of Kuramoto oscillators¹⁵ as the one shown in Fig. 1(a). Note that the results below may apply to a larger set of diffusively coupled systems. Their dynamics is governed by $2n$ coupled differential equations

where degrees of freedom and natural frequencies, respectively, in layers 1 and 2 are denoted as ${xk}$⁠, ${yk}$ and $ωk(1)$⁠, $ωk(2)$⁠. The undirected coupling network in the $l$th layer is given by the adjacency matrix elements $bij(l)≥0$⁠, and $fi$ is a coupling function between the two layers. The noise acting on the first layer is encoded in $ηi$ that is taken as white and uncorrelated in space, i.e., $⟨ηi(t)ηj(t′)⟩=δijη02δ(t−t′)$⁠. This model made of Kuramoto oscillators is a nonlinear version of the one used in Ref. 1. Thanks to the sine coupling and provided that $bij(l)$’s are sufficiently large compared to the distributions of $ωi(l)$’s, Eq. (1) may have multiple stable or unstable fixed points $({xk(0)},{yk(0)})$⁠, depending on the coupling networks. In both cases, they satisfy

In the following, I consider a simple coupling function $fi({xk},{yk})=d(xi−n−1∑jxj)$ where degrees of freedom in the first layer tune the natural frequencies in the second one, and $d$ serves as a tuning parameter. The dynamics of small deviations close to a stable fixed point is given by the Taylor expansion of Eq. (1) to the first order in $[δxi(t),δyi(t)]=[xi(t)−xi(0),yi(t)−yi(0)]$⁠, and read

which are diffusive linear systems. Indeed, one can define the weighted Laplacian matrices,

which depend on both the initial coupling networks in each layer and the considered equilibrium. Using these Laplacian matrices, Eq. (3) is conveniently rewritten in a vector form as

where $δx¯i=δxi−n−1∑jδxj$⁠. From the Laplacian properties of $Lij(1)({xi(0)})$ and $Lij(2)({yi(0)})$⁠, their eigenvalues satisfy: $0=λ1(k)<λ2(k)<⋯<λn(k)$ for $k=1,2$ with the first eigenvector being constant, i.e., $u1(k)=(1,…,1)/n$  for $k=1,2$⁠. The latter modes correspond to a symmetry of the system as any shift along it does not modify the dynamics or the fixed points. The second smallest eigenvalue $λ2(k)$ is usually called algebraic connectivity of the network, with its corresponding eigenvector $u2(k)$ called Fiedler mode. The larger the algebraic connectivity, the better connected is the network, in the sense that even the slowest modes are quickly damped. The solution of Eq. (4) is obtained by expanding the deviations over the eigenvectors as $δxi=∑αcα(1)uα,i(1)$⁠, $δyi=∑αcα(2)uα,i(2)$⁠.

It has been shown recently that, for linear systems governed by the dynamics of Eq. (4), spatially uncorrelated and white in time noise acting on the first layer $xi$’s might lead to amplified fluctuations in the second layer $yi$’s.¹ Indeed, the variance of the degrees of freedom in the first layer is generated by $ηi$’s, which are i.i.d. and white in time, and yields¹ (see  Appendix A for details)

where $⟨⋯⟩$ refers to the average. In contrast, the variance in the second layer is generated by the $xi$’s, which play the role of additive noise that is both correlated in space and time as shown by their two-point, time shifted correlator,¹ 

Using the latter correlator together with Eq. (4) allows us to calculate the variance in the second layer, which yields

For linearly coupled oscillators, the variance is exactly described by Eqs. (5) and (7). Depending on the overlap between the slowest eigenvectors of $Lij(1)({xi(0)})$ and $Lij(2)({yi(0)})$⁠, fluctuations in the second layer might be strongly amplified. In the particular case where the networks are the same in both layers, i.e., $Lij(1)({xi(0)})=Lij(2)({yi(0)})$⁠, Eq. (7) simplifies to

Due to the third power of $λα$ at the denominator, the variance is dominated by the slowest modes. So, not only the overall amplitude of the fluctuations might be large but also their shape follows those of the slowest modes. This is shown in Fig. 2. For both $xi$’s and $yi$’s, an important contribution to the variance comes from $u2,i2$ [see Eqs. (5) and (8)]. However, for the second one, this contribution dominates even more due to the denominator in Eq. (8). One clearly sees that nodes corresponding to the largest $u2,i2$'s (lighter colors) have the largest standard deviations in the second layer. Therefore, the second layer is perturbed away from its initial equilibrium mostly along the directions of the slowest eigenmodes.

While the linear system on the first line of Eq. (4) has a single fixed point (neglecting the symmetry along $u1(1)$⁠), it is not the case of nonlinearly coupled systems, which may have multiple stable fixed points each of them with their own basin of attraction. If the excursion from the initial equilibrium becomes too important, the system undergoes a transition to another basin of attraction corresponding to another fixed point or reaches an unstable region of the phase space. Obviously, large perturbations leading to important excursions result in short first escape times. However, not only the amplitude but also the direction of the perturbation matter. Indeed, basins of attraction are not isotropic,^5,16,17 and excursion along specific directions may lead to faster transitions. This is already pointed at by Fig. 1(c), where the amplitudes of the noise acting on layers 1 and 2 are the same. However, the second layer (in red) operates more transitions than the first one (in blue), which stays in its initial basin of attraction.

In order to detect escapes, one can use the winding number defined on a cycle $c$ as

where $|.|[−π,π)$ brings the argument into the interval $[−π,π)$ and $zi$ is the degree of freedom of the $i$th oscillator in the cycle. The winding number is, therefore, an integer given by the number of times that the degrees of freedom wind around the origin along a cycle $c$ [see Fig. 1(b)]. On a planar meshed network, one can further define a winding vector $q$ whose components are given by the winding number on each cycle composing the network. Each stable fixed point of Eq. (1) is unambiguously identified by a pair of winding vectors $(q1,q2)$⁠, where $ql$ is the collection of winding numbers on each cycle of the $l$th layer. Therefore, by recording $(q1,q2)$⁠, one is able to detect basin’s escapes independently in each layer. This technique is useful in cases where one is able to measure the flows on each edge and, therefore, determine the winding numbers on each cycle. Another way of measuring first escape times is simply by letting the system relax at every time step and checking that the stable fixed point is the same as the original. This is for sure more precise in determining the first escape time; however, it is much more costly in terms of numerical simulations.

In this section, I compare the first escape time in both layers varying the settings. I first consider a system made of two times the same network. To disentangle the effects of amplitude and shape of the noise coming from the first layer, I investigate the situations where (i) the noise is not rescaled, i.e., large fluctuation might appear in the second layer; (ii) rescaled noise in the second layer so that noise amplitudes in both layers are the same on average, but their shapes are different. The latter is done by tuning the $d$ parameter in the inter-layer coupling function. More details are given below. Then, I briefly consider the case of different networks. The method to detect escape based on winding numbers is used for cycle networks, and a combination or the more computationally expensive one and the winding numbers for more complex networks (see Sec. II C). For simplicity, I consider identical oscillators, i.e., $ωi(l)=0$⁠, $∀i,l$⁠.

When the two networks are the same, the noise injected in the second layer has statistical properties that are directly connected to its structure [see Eq. (6)]. One, therefore, expects an important difference in the first escape times of the two layers. More precisely, for networks with low algebraic connectivity, fluctuations are amplified in the second layer, and, therefore, one expects shorter first escape times for the second layer compared to the first one. For the network with high algebraic connectivity, one expects the opposite, as fluctuations are reduced in the second layer.

Here, we set $d=1$ in the coupling function, meaning that the noise in the second layer has both a specific modal structure and is potentially strongly amplified or reduced. In this case, the two networks are the same and, therefore, share the same eigenvectors and eigenvalues. One expects enhances fluctuations in the second layer if the algebraic connectivity is less than one, i.e., $λ2≪1$⁠.¹ This is, in particular, the case for large cycle networks whose algebraic connectivity scales with $n$ as $λ2=2−2cos⁡(2π/n)$⁠. Here, I focus on such networks with amplification, as those with $λ2≫1$ have reduced fluctuations in the second layer and, thus, are less likely to leave the initial basin of attraction. Figure 3 gives the first escape times for both layers in a system made of two cycles. One clearly sees that the first escape time is shorter by several orders of magnitudes in the second layer compared to the first one for noise amplitude $η0<2.5$⁠. Then, for larger $η0$⁠, the first layer exits its initial basin of attraction faster than the second one. In this regime, the dynamics in the first layer is essentially dominated by the noise with only little influence coming from the coupling network. The noise injected in the second layer may, therefore, not have any amplification or specific structure in this case as the linearized calculation presented here does not apply anymore.

In order to differentiate the effect of noise amplitude and shape, one needs to rescale the noise coming from the first layer so that amplitudes are the same in both layers. The noise coming from the first layer and to which is subjected the second layer is given by $ξ=dδx¯$⁠. The noise variance in the first layer is $η02$⁠. Choosing $d2=2/∑αλα(1)−1$ produces an input noise in the second layer that has on average a variance of $η02$⁠, i.e.,

so that both layers are subjected to noises of similar variance on average but yet have different shapes and time correlations [see Eq. (6)]. This is illustrated in Fig. 4 for a system made of two cycles, where the noise in the second layer is rescaled. Interestingly, one observes that the first escape time in the second layer is still shorter than in the first one for $η0<1.7$⁠. This means that the structure of noise as discussed in Sec. II B favors faster escapes from the initial basin of attraction compared to uncorrelated noise.

For two different networks, noise transmission essentially depends on the overlap between $uα,l(1)$ and $uβ,l(2)$ for the slowest modes. If the overlap is low, one expects that the noise is transmitted mostly without amplification and, therefore, the first escape times should be similar in both layers. In the other situation where the overlap is significant, one expects different first escape times in the two layers.

Similar effects as in Sec. III A 1 can be observed when networks are different in the two layers. Here, I consider two networks generated following Watts–Strogatz¹⁸ procedure with first and second nearest neighbor couplings and low rewiring probabilities. Doing so, the algebraic connectivities of both networks are smaller than one, i.e., $λ2(1)=0.033$ $λ2(2)=0.084$⁠, and the overlap between their slowest eigenvectors is sufficiently different from zero, i.e., $∑iu2,i(1)u2,i(2)=0.4$⁠. Therefore, fluctuations are expected to be amplified in the second layer. Indeed, Fig. 5 shows the first escape time in both layers, and the second layer leaves the initial basin of attraction many orders of magnitudes faster than the first one for low values of $η0$⁠.

Following the same steps as in Sec. III A 2 to emphasize the effect of the shape of the noise, I rescale the noise such that Eq. (10) is satisfied. This allows us to have injected noises of similar variance on average in each layer. Remarkably, Fig. 6 shows that the first escape time in the second layer is still shorter than in the first layer, despite the rescaling of the noise for $η0<2.6$⁠. The shape that the noise acquires when passing through the first layer seems to be more efficient in driving the second layer outside its initial basin of attraction, as in the previous case of two times the same network.

Noisy layered systems can exhibit amplified fluctuation patterns depending on their connectivity. Here, I showed that noise originally injected in one layer may induce faster basin escape in connected layers. This is both due to the amplification of the noise amplitude and the system-specific correlations that the noise acquires while going through the first layer. Indeed, from Eq. (6), one sees that the noise in the second layer is correlated in both space and time with clear dependence on the network structure. For networks with low algebraic connectivity, I numerically showed that the first escape time is shorter in the two cases where (i) fluctuations are amplified in the second layer and (ii) noise in the second layer is rescaled in order to have the same variance in both layers. While point (i) is rather intuitive, i.e., larger fluctuations lead to shorter first escape times, point (ii) is more involved. Indeed, this indicates that noise with spatial and temporal correlations given by Eq. (6) selects directions that enable faster exits from the initial basin of attraction.

Even though I consider Kuramoto oscillators, the results presented here may apply to other nonlinear coupling functions as long as the system has an equilibrium point. Further studies should consider time-correlated noise injected in the first layer and the effect of inertia.

This work was supported by U.S. DOE/OE as part of the DOE Advanced Sensor and Data Analytics Program and by the Laboratory Directed Research and Development program of Los Alamos National Laboratory under Project Nos. 20220797PRD2 and 20220774ER.

The authors have no conflicts to disclose.

Melvyn Tyloo: Conceptualization (lead); Data curation (lead); Investigation (lead); Methodology (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (lead).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

In order to calculate $⟨δxi2⟩$⁠, one needs the solution to the first line of Eq. (4). The latter is obtained by expanding the $δxi$'s over the eigenvectors of $L(1)$⁠, i.e., $δxi=∑αcα(1)uα,i(1)$⁠, which yields the differential equations,

The solution of the latter equations is given by

where the initial condition is $cα(1)(t=0)=0$$∀α$⁠. Then, the long time limit of the variance in the first layer is given by

where in the last equality, I used $⟨ηi(t)ηj(t′)⟩=δijη02δ(t−t′)$ and the orthogonality between the eigenvectors. Performing the integral and taking the limit of a long time gives the results of Eq. (5).