Relay synchronization in multiplex networks

Relay synchronization in multiplex networks www.nature.com/scientificreports OPEN Relay synchronization in multiplex networks 1,2 1,2 3 3 4 I. Leyva , I. Sendiña-Nadal , R. Sevilla-Escoboza , V. P. Vera-Avila , P. Chholak & S. Boccaletti Received: 16 February 2018 Relay (or remote) synchronization between two not directly connected oscillators in a network is an Accepted: 21 May 2018 important feature allowing distant coordination. In this work, we report a systematic study of this Published: xx xx xxxx phenomenon in multiplex networks, where inter-layer synchronization occurs between distant layers mediated by a relay layer that acts as a transmitter. We show that this transmission can be extended to higher order relay configurations, provided symmetry conditions are preserved. By first order perturbative analysis, we identify the dynamical and topological dependencies of relay synchronization in a multiplex. We find that the relay synchronization threshold is considerably reduced in a multiplex configuration, and that such synchronous state is mostly supported by the lower degree nodes of the outer layers, while hubs can be de-multiplexed without affecting overall coherence. Finally, we experimentally validated the analytical and numerical findings by means of a multiplex of three layers of electronic circuits. Synchronization is one of the most important collective phenomena in nature. It can be observed in natural, 1–3 social and technological systems, and it became one of the most active research topics in network science . e Th huge amount of new data collected in the last years has permitted a higher resolution network representation of real systems. In particular, the inclusion of new features shaped multi-layer representations, i.e. approaches in which the network units are arranged in several layers, each one accounting for a different kind of interactions 4–7 among the nodes . Multi-layer structures determine scenarios where novel forms of synchronization are rele- 8,9 vant. Despite an analytical approach has been tackled in just a few particular cases , several synchronization sce- narios have been already addressed, as unidirectional coordination between layers , explosive synchronization 11,12 emerging from the interactions between dynamical processes in multiplex networks , complete synchroniza- 13,14 15–21 22 23,24 tion , cluster synchronization , intra-layer or inter-layer synchronization. Very recently, relay (RS) and remote synchronization (two very well known phenomena in chains, or small motifs, of coupled oscillators) have captured the attention of researchers. This form of synchronization is observed when two units of a network (identical or slightly different) synchronize despite not being directly linked, and due instead to the intermediation of a relay mismatched unit. The phenomenon has been experimentally detected in 25 26,27 lasers and circuits . In general, the relay units exhibit generalized or delay synchronization with the units they actually pace to synchrony . RS is of outstanding relevance in the brain: the thalamus acts as a relay between distant cortical areas through 29–32 the thalamo-cortical pathways, playing the role of a coordination hub that maintains the information flow . Complex structures and neuronal dynamics are implicated in this process involving not only simple, but higher 30,31 order relay paths, that transfer the information through multiple-step relay chains . Recently, remote syn- chronization has been addressed in the context of complex networks , revealing the extremely important role of 34–36 network structural and dynamical symmetries in the appearance of distant synchronization , as it was already 37,38 suggested by the observation of zero-lag delays between mirror areas of the brain . Nevertheless, the interplay between symmetry, dynamics and multi-layer structure remains still mostly unexplored. In this work, we perform a systematic study of inter-layer relay synchronization in a multiplex network, where distant layers synchronize their dynamics while their intra-layer motion remains unsynchronized. We con- sider generic high-order structures where multi-site relay pathways are verified. The dynamical and topological 1 2 Complex Systems Group & GISC, Universidad Rey Juan Carlos, Móstoles, Madrid, 28933, Spain. Center for Biomedical Technology, Universidad Politécnica de Madrid, Pozuelo de Alarcón, Madrid, 28223, Spain. Centro Universitario de los Lagos, Universidad de Guadalajara, Jalisco, 47460, Mexico. Department of Mechanical Engineering, Indian Institute of Technology Bombay, Powai, Mumbai, 400076, India. CNR-Institute of complex systems, Via Madonna del Piano 10, Sesto Fiorentino, 50019, Italy. Correspondence and requests for materials should be addressed to I.L. (email: inmaculada.leyva@gmail.com) SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 1 www.nature.com/scientificreports/ Figure 1. Schematic representation of a multiplex of 2M + 1 layers (here M = 2) labeled as k = −M, …, −1, 0, 1, …, M where each pair of layers k and −k (painted with the same color) are networks of identical oscillators with the same topology  and intra-layer coupling σ and whose dynamical state is described by the variable U and −k U , respectively. The multiplex is symmetric with respect to the layer k = 0 and the nodes are coupled to their replicas in the rest of layers with an inter-layer coupling strength λ. dependencies of the phenomena are studied, using perturbation stability analysis. The robustness of the relay synchronization against de-multiplexing the layers is reported, revealing the key role of low degree nodes in maintaining the layers coordination. Finally, the findings are experimentally validated in a multiplex network of electronic circuits. Results Model. We start by considering 2M + 1 layers (or networks), arranged as shown in Fig. 1. Each layer k, with k =− MM ,, …… 0, , , is formed by N dynamical systems (each of which being m-dimensional), whose states are kk k k km represented by the column vectors U =… {, uu ,, u }, with u ∈  , i = 1, …, N, and whose intra-layer inter- 12 N i k k actions are encoded through the Laplacian matrices  = {}  . The layer stack is symmetric with respect to k = 0 ij k −k in such a way that Laplacians  and  have the same structure. The dynamical systems are also paired: nodes +k −k belonging to layers U and U are identical to each other, and different (in some parameter) from the rest of the layers. Consequently, layer k = 0 has no counterpart, and acts as a relay between all layers situated above and below it. Layers are coupled in a multiplex cong fi uration, and the dynamical evolution of the system is described by the following equations: qk =+1 ≤M k kk qk UF =− () UG σλ ()  ⊗+ UH ()  ⊗− () UU k k N qk =−1 ≥−M (1) k k k kT mm where the functions FU () =… [( fu ), fu (),, fu ()] (with f :  → representing the vectorial func- k k 12 k kN k tions evolving each dynamical unit), are identical for the same |k|. G, H are the m × m matrices representing respectively the linear intra- (G) and inter- (H) layer coupling schemes.  is the N × N identity matrix, σ is the N k intra-layer coupling strength within layers k and −k, and λ is the inter-layer coupling strength. +k −k Due to the reflection symmetry of the system under study (i.e. as long as the U and U layers are identi- cal for all k), a synchronous inter-layer evolution (with layers evolving in a pairwise synchronized fashion, i.e. + −k k k′ where U = U ) at all k without necessarily implying U = U for k ≠ k′) is a solution of Eq. (1), independently of intra-layer synchronization (i.e. independently on whether the state of the systems within each layers are SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 2 www.nature.com/scientificreports/ k +k −k synchronized). Let therefore δU (t) = U (t) − U (t), with k = 1, …, M be the vector describing the difference between the dynamics of the paired layers. kk k k Considering the smallness of and expanding around the inter-layer solution up to δδ Uu =… {, δδ uu ,, } 12 N first order, one obtains a set of M linearized vector equations for the perturbations δU : k k ∼∼ k k k δσ Uf =− [( JJ UG )( ⊗− )( λδ 2) − J HU ()]δU k k kM qk =+1 +δ λ J HU () U qk =−1qk ≠ (2) where J denotes the Jacobian operator, δ is the Kronecker delta accounting for the boundary condition at k = M kM ±M (as the stack end layers U are only connected to the previous neighbor layer). The vector U = {} u describes k −k 0 the dynamical state of any of the k = 0, …, M layers at the synchronous state U = U ≠ U and, therefore, the whole dynamics is reduced to the dynamics of M + 1 layers. Such evolution at the node level is given by: qk =+1 ≤M k k k q k uf  =− () ug σλ  () uh  +− [(uh )(u )] ∑∑ i ki k ij j i i j qk =−10 ≥ (3) where i = 1 …, N is the node index, k, q = 0, 1, …, M, and g(u) = Gu and h(u) = Hu are the projections of the inter- and intra- later coupling operators to the node level. Notice that, since each paired layers k and −k is kk − inter-layer synchronized () UUU == , each layer acts therefore as a relay to the rest of the stack. The prob- lem consists now in solving MmN linear equation (2), together with solving in parallel the (M + 1)mN nonlinear equation (3) for . Although the total number of equations to compute is still of the same order as in Eq. (1), that u is (2M + 1)mN, the fact that MmN of those equations are linear results in a much faster computation of the dynamics. Then, calculating the maximum Lyapunov exponent (MLE) transverse to the manifold U as a function of the system parameters actually gives the necessary conditions for the stability of the synchronous solution: whenever MLE0 < , perturbations transverse to the manifold will die out, and the multi-relay synchronous solu- tion will be stable. In order to monitor the synchronization error between layers, we define the inter-layer synchronization error as, q k E =− lim uu () tt () dt, qk , ∑ i i T →∞T 0 (4) i=1 where ||·|| stands for the Euclidean norm and q, k are the layers’ indexes, such that E denotes the inter-layer −k,k synchronization error of mirror layers. Without lack of generality, In our numerical simulations we consider two types of topologies where layers are either (i) Erdös-Renyi (ER) or (ii) scale-free (SF, generated by means of the Barabási-Albert’s algorithm ), in all cases with N = 500. We classify the layer stacks regarding the topology sequence of each layer. For instance, a triplex where the three layers have ER topology will be denoted as EEE, and a system where two identical SF layers are mediated by a center ER will be denoted as SES. The nodes are chaotic Rössler oscillators , defined by the m = 3 state vector u = (x, y, z) whose autonomous evolution is given by f () uf == () u [, −− yz xa +. yz ,0 2( +− x 9)] and the heterogeneity between the layers is introduced kk − k through the parameter a , such that each layer develops a different chaotic dynamics. In our case study, the intra- T T and inter- layer coupling functions are set to be g() uG == u (0,0,) z and hu () == Hu (0,, y 0) respec- tively. These coupling schemes ensure that intra-layer synchronization is prevented when layers are isolated and not multiplexed (class I layers, according to the standard master stability function (MSF) classification established in ref. ) whereas multiplexed nodes along the layers can synchronize for a coupling strength λ above a given threshold (class II MSF) . Layers with identical topology. With the aim of determining whether relay synchronization can be achieved in a multiplex cong fi uration let us first consider the multiplex structure den fi ed by Eq. ( 1) for the case of 01 −1 three identical SF layers  ==  and where the parameters a = a = 0.2 for the outer layers and a = 0.3 1 −1 0 for the relay units of the central layer, although different selections of these parameters and topologies produce a similar behavior. Results are collected in Fig. 2, where the synchronization error between the outer layers E is plotted versus −1,1 the inter-layer coupling λ for several values of the intra-layer couplings σ and σ in the outer and relay layers 1 0 respectively, with σ = σ . In all cases, there is a critical coupling λ above which complete synchronization 1 0 between layers k = 1 and k = −1 occurs, that is, E = 0 is achieved, while the relay layer (k = 0) remains unsyn- −1,1 chronized to any of the two outer layers (k = 1, −1) as shown in the inset where E > 0 for any parameter choice. 0,1 In addition, the calculation of the corresponding MLE given by Eq. (2) (bottom panel of Fig. 2) confirms −1 1 * that the relay synchronous solution U = U reaches stability (MLE < 0) at the same critical λ where the error between the outer layers is zero, as indicated by the vertical lines. Therefore, one can conclude that inter-layer MLE is a useful tool for reducing the system’s dimensionality and use it for evaluation of the critical inter-layer coupling λ from now on. In order to better understand the different roles played by external and relay layers, we show in Fig.  3 the criti- cal inter-layer coupling value in the parameter region (σ , σ ), that is, when the intra-layer coupling σ is different 0 1 k SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 3 www.nature.com/scientificreports/ Figure 2. Relay synchronization in a triplex (M = 1) with identical SF layers (SSS configuration). (Main panel) Synchronization error between the outer layers (k = −1, and k = 1) E (see Eq. 4) as a function of the inter- −1,1 layer coupling λ for three different values of the intra-layer coupling σ = σ (see legend). The inset shows the 0 1 corresponding synchronization errors between the relay and one of the outer layers. (Bottom panel) Maximum 1 −1 Lyapunov exponent (MLE) of the relay synchronization manifold U = U as a function of λ for the same cases as in the main panel. Vertical lines mark the points where the MLE becomes negative. All points are averages of 10 network realizations with N = 500 and ⟨k⟩ = 4. See the main text for the relay and outer layer Rössler oscillators specifications. Figure 3. Relay synchronization in a triplex network with identical SF layers as a function of the intra-layer couplings for the relay (σ ) and outer (σ ) layers. (Left) Color map of the inter-layer coupling threshold λ for the 0 1 relay state (E = 0 and E ≠ 0) in the σ − σ parameter space. (Right) Inter-layer coupling threshold λ for the −1,1 0,1 0 1 relay state as a function of the coupling strength in the relay layer σ for a fixed value of σ = 1 (red dashed line in 0 1 left panel) and as a function of the coupling strength in the outer layers σ for a fixed value of σ = 1 (black dashed 1 0 line in left panel). Each point is an average of 10 SF network realizations with N = 500 and ⟨k⟩ = 8. for the relay and outer layers. It can be seen that the system’s ability to synchronize is practically unaltered with σ , while increasing σ makes the value of λ to drop drastically. This therefore reveals that multiplex relay synchroni - zation is much more sensitive to changes ae ff cting the mirror layers than to those arising in the transmission layer. Our results can be generalized to any number of layers. As an example, we report also the case M = 2, which corresponds to two outer layers above (k = 1, 2) and below (k = −1, −2) the relay layer (k = 0). We choose a = a = 0.2 and a = a = 0.3, and a = 0.25 for the central layer. The results stand for any other parameter −1 1 −2 2 0 choice. In Fig. 4 we plot the inter-layer synchronization errors E (void markers) and E (full markers), vs. −1,1 −2,2 the inter-layer coupling λ for several values of the intra-layer coupling σ. As in the triplex case, the critical λ at which complete inter-layer synchronization is achieved depends on σ, but it is the same for both pairs of layers, as E and E drop to zero simultaneously. In the inset we plot the inter-layer synchronization errors between −1,1 −2,2 the non-paired layers, E , E to check that they remain mutually incoherent. Therefore, we have verified that 0,1 1,2 relay synchronization also occurs in cascade for arbitrarily high-order multiplex systems, provided a structural and dynamical symmetry is conserved. Layers with non-identical topology. So far, we have dealt with multiplexes of pairwise identical layers. However, this condition is too strong a limitation to hope that it would capture and properly represent the case of many real systems. The next step needed for generalization is studying then the relay synchronization scenario in SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 4 www.nature.com/scientificreports/ Figure 4. Relay synchronization in a pentaplex (M = 2) with identical N = 500 ER layers (EEEEE configuration). The synchronization error between the two pair of outer layers E (empty symbols) and E −1,1 −2,2 (full symbols) is shown as a function of λ for three different values of the intra-layer coupling σ, being σ = σ , ∀k. The inset shows the synchronization errors between each one of the outer layers and the relay layer. The results are averaged over 10 different network realizations and initial conditions. Figure 5. Relay synchronization in a triplex with different layers. Inter-layer coupling threshold λ vs the intra- layer couplings σ = σ for (a) a mixed ER-SF-ER (ESE) and identical (EEE) configurations and (b) a mixed SF- 0 1 ER-SF (SES) and identical (SSS) configurations. the case in which the topology of the relay layer differs from that of the outer layers. In Fig.  5 we have reported the critical inter-layer coupling λ in two heterogeneous triplex cases: (a) a pair of Erdös-Rényi layers mediated by a scale-free relay layer (ESE situation) and (b) the opposite case where SF layers are connected through a ER layer (SES). Each case is compared with the topologically homogeneous EEE and SSS structures, respectively. For the sake of simplification and of better assessment of the role of the topology, we keep σ = σ . 0 1 Figure 5(a) shows that, for a large range of intra-layer couplings, the mediation of a SF relay facilitates the synchronization between the paired layers, since λ in the ESE case (void blue circles) is smaller than the one corresponding to the homogeneous case (EEE, full blue circles). On the contrary, a relay ER layer intermediat- ing between two outer SF layers (Fig. 5(b)) does not determine a significant difference as long as the intra-layer coupling strength is low, but increases the threshold λ for higher σ, as compared to the homogeneous SSS case. Robustness. In the previous Sections we have addressed the dependence of relay synchronization in a multi- plex on the dynamical and structural layer heterogeneity, and proved that the phenomenon still holds even when the intermediate layer has a completely different structure and dynamics than the mirrored ones. The present section is devoted instead to assess the robustness of relay synchronization against structural changes by means of a de-multiplexing process of the layers, that is, against performing a progressively shutting down of the inter-layer links such that a fraction of nodes in each layer is not linked to their counterparts in the other layers. In addition, we also investigated several cases to test the robustness against mismatches in the oscillators parameters to closely resemble real experimental conditions. Structural robustness. To closely check this process, we initially consider a 3-layer multiplex with identical topol- ogy (EEE or SSS). We choose the inter- and intra-layer couplings to guarantee a relay synchronous state with the layers fully multiplexed. Then, we proceed to disconnect one by one the inter-layer links according to the nodes degree ranking, both in the ascending and the descending order, and re-evaluate in every step the state of the relay synchronization by measuring the E error. Figure 6(a) reports the averaged evolution of E as a function −1,1 −1,1 of the number of multiplexed nodes, aer h ft aving performed the whole de-multiplexing process for 10 different network realizations. It can be seen that, starting from a situation with E = 0, the EEE multiplex configuration −1,1 (blue void markers) loses the synchronization immediately with just a few of inter-layer links being removed. On SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 5 www.nature.com/scientificreports/ Figure 6. Structural robustness of the network relay synchronization for identical layers. (a) Synchronization error between the outer layers E vs the decreasing number of connected relay lines for identical ER (blue −1,1 empty symbols) and SF (black solid symbols) layers. Relay lines are disconnected following a descending (circle symbols) or ascending (square symbols) node degree ranking of the outer layers (seed legend in (b)). Parameter values are N = 500, ⟨k⟩ = 8, λ = 0.23 and σ = σ = 0.8. (b) Number of multiplexed relay lines needed 0 1 to support a relay network state as a function of the intra-layer coupling strength σ = σ while keeping constant 0 1 λ = 0.23. The different curves are explained in the legend. All the results are averaged over 10 different network realizations. the other hand, relay synchronization is resilient in SSS triplex configurations also when more than 30% of the nodes are not multiplexed. A more detailed view can be obtained from Fig. 6(b), where the number of multiplexed nodes needed to support the relay synchronization is represented as a function of the intra-layer coupling σ = σ . As expected, 0 1 when the coupling is weak, all the N nodes need to be linked to preserve relay synchronization. However, as the interaction within the layers increases, the intra-layer connectivity helps to maintain a synchronous state despite an increasing number of nodes are being de-multiplexed without damaging the coherence between the outer layers. In Fig. 6(b), we can see that for both the EEE (blue void markers) and the SSS (black full markers) tri- plex configurations, removing the links between layers connecting nodes with higher degree (descending degree ranking, circle markers) is much more robust than following an ascending degree ranking (square markers). This is indeed a very interesting result: relay synchronization in a multiplex network is supported by the low degree nodes, while the hubs can be safely disconnected without perturbing the transmission. This is notably evidenced in the SSS case (black full squares) where after having removed the 40% of the inter-layer links connecting the highest degree nodes, the relay synchronization is still supported by the multiplex structure connected through the lowest degree nodes. Once we have singled out that, from our trial rankings, the descending degree ranking is the most conven- ient way to de-multiplex part of the network without loosing coherence, we proceed our study by evaluating the impact of having a relay layer with different topology from the outer layers, as we did in the previous Section 4. In this scenario, we have two possible descending degree rankings, the one dictated by the relay layer and the one dictated by the outer layers. The results are summarized in Fig.  7 where we plot, as in Fig. 6(b), the number of nodes that need to be linked to maintain synchronization as a function of σ = σ . For the sake of comparison, 0 1 we added the curves for the homogeneous EEE and SSS (full markers) multiplex configurations, together with the data for the mixed ESE and SES (void markers) layers. Notice that the chosen inter-layer coupling λ = 0.23 is well above threshold for all the cases, as it can be derived from Fig. 5. All the reported evidence indicates that the introduction of a relay layer with a topology different from that of the outer layers has little influence on the min- imum number needed to support the relay synchronization, as long as the first removed inter-layer connections correspond to the hubs in the outer layers (blue and back curves). Curiously, the alternative of using the relay layer topology to rank the degree of the nodes, destroys the coherence between the outer layers as soon as a tiny fraction of links is removed (red curves). er Th efore, the relay synchronization in a multiplex is very unstable if just a few links connecting nodes which are hubs in the relay layer are removed. Notice that this unlinking criterion is equivalent to randomly disconnect the multiplex. Therefore, the robustness of the relay synchrony relies mainly in the low degree nodes of the external layers. The relevance of the low degree nodes in controlling the dynamics 43,44 of complex networks has been pointed out in other contexts . Dynamical robustness. In order to explore a more realistic situation of a multiplex composed of non identical oscillators, we have taken further the generalization by introducing some heterogeneity in the node’s parameter a , that is, a = a . The node values of the outer layers’s were randomly picked from the interval a = [0.20, 0.21]. k k k,i ±1,i We prepared three different setups depending on how the parameter setting is introduced in the outer layers: i) a = a but aa ≠ , that is, oscillators are not identical within each external layer but the two external 1,i −1,i 1,ij 1, layers are symmetric; ii) a = a but a ≠ a , that is, layers are composed of identical oscillators but there is a 1,i 1,j 1,i −1,i mismatch between oscillators in layer k = 1 and their corresponding replicas in layer k = −1; and iii) SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 6 www.nature.com/scientificreports/ Figure 7. Structural robustness of the network relay synchronization for non identical layers. Number of multiplexed relay lines needed to support a relay network state as a function of the intra-layer coupling strength σ = σ while keeping constant the inter-layer coupling strength λ = 0.23 for mixed ER (ESE, empty circles) and 0 1 mixed SF (SES, empty squares) layer configurations. The relay lines are disconnected following a descending order of the outer node degrees and for comparison the corresponding values for identical layers -see Fig. 6(b) are plotted in solid symbols. The red solid (ESE) and empty (SES) triangles show the behavior when the relay lines are disconnected following the degree ranking of the relay layer. Figure 8. Dynamical robustness of the network relay synchronization as a function of the heterogeneity in the a parameter. Synchronization error between the outer layers E vs the inter-layer coupling λ for the three k −1,1 cases reported in the main text with a =. [0 20,0.21]: i) mismatched layers with oscillators within the outer ±1,i layers being identical but mismatched with their replicas (green squares); ii) symmetrical outer layers with nonidentical oscillators within the layers (full cyan circles); and iii) not symmetrical layers with nonidentical oscillators (purple diamonds). Inset is a zoom of the main plot. Simulations have been conducted with SF layers of identical topology with N = 500 and ⟨k⟩ = 8. a ≠ a ≠ a , that is, combining cases i) and ii), a fully random case, not preserving the symmetry neither the 1,i 1,j −1,i parameter homogeneity in each layer. Results are summarized in Fig. 8, where the synchronization error between the outer layers is plotted vs the inter-layer coupling λ. For comparison, the perfectly homogenous case with all layers having identical oscilla- tors (red circles) is also included. In the inset, a zoom of the main figure, it is clear that breaking the symmetry between the outer layers slightly deteriorates the synchronization, as it occurs in cases i) (green squares) and iii) (purple diamonds). However, introducing some heterogeneity in the a parameter (full cyan circles) in the layer’s oscillators does not modify the synchronization threshold at all, as long as each oscillator in one outer layer is identical to its replica in the mirrored layer. We will actually show that our experimental results can be framed within this latter case. Experimental validation Finally, we report experimental evidence of relay synchronization in a multiplex of nonlinear electronic circuits, with the setup sketched in Fig. 9 (left). The array is made of 21 Rössler-like circuits arranged in three layers of 7 nodes, with the relay layer having different topology as the outer layers. Each layer has two different electronic couplers, one for the coupling among nodes in the same layer (σ ) and the second for the interaction of each node with its replica in the other layers (λ ). The chaotic dynamics of the circuits is well approximated by the three variables (v , v , v ) obeying : 1 2 3 SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 7 www.nature.com/scientificreports/ Figure 9. (Left) Schematic representation of the experimental arrangement of three layers of electronic circuits. e b Th idirectional coupling is adjusted by means of three strips of digital potentiometers X9C103 (XDCP), the resistance is controlled through digital pulses sent by a DAQ (NI USB 6363). (Right) Graph structure used for the upper and lower layers (top) and for the relay layer (bottom). C = 1 nF C = 1 nF C = 1 nF σe,λe = [0−0.6] 1 2 3 R = 2MΩ R2 = 200kΩ R3 = 10kΩ R4 = 100kΩ R = 50kΩ R6 = 5MΩ R7 = 100kΩ R8 = 10kΩ R = 10kΩ R10 = 100kΩ R11 = 100kΩ R12 = 150kΩ R = 68kΩ R14 = 10kΩ R15 = 75kΩ R16 = 120kΩ 0 1 Rk =Ω 50 Rk =Ω 35 I = 0.7 V = 15 d ee c c Table 1. Parameter values of the chaotic dynamics of one Rössler like circuit as described in Eq. (5).   1 R R k k 1 k 1 k   v =− v ++ v v  1i 1i 2i 3i RC  R R   11 2 4 1 R 1 k k k −− σ av v e ∑ ij () 11 j i RC R 11 15 j=1     RR RR 1  k  k k 68 68    v =− −+ v 1 − v 2i 1i 2i  k  RC  RR RR  62  97  c 7     q=1   1   q k − λ vv − e ∑ 22 i i RC  R 62  16 q=−1    k k k 10  v =− −+ Gv v  3i () 13 i i RC  R    10 3 11 (5) where G is a nonlinear gain function given by: 1i   R R   14  14 0, if vI ≤ 1 +  + V   1id ee    R R   13 13 Gv () = 1i      R R R R R R  12 12 12 12 14 14      vV −− ee I + ,1 if vI > +  + V   1id 1id ee       R R  R R  R R      14 13 13 14 13 13 where the parameter values are gathered in Table 1. The reader is referred to ref. for a detailed description of the 9,23,24,46 experimental implementation of the Rössler-like circuit in the networks, and refs for previous realizations in different network configurations. Both the intra-layer σ and the inter-layer λ are set by means of the digital e e potentiometers X9C103, that working as voltage divisor for the maximum resistance (10k Ω), σ and λ is set to e e zero, this potentiometers are controlled through the digital ports (P0.0, P0.1, P0.2, P0.3) of a DAQ card. First that all we send all the coupling value to zero, aer 500 ft ms takes the sample of the time series of each networks, all the variables v of each oscillator enter to the DAQ card through the analogue ports (AI0, AI1, …, AI20) and saved in 2i SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 8 www.nature.com/scientificreports/ Figure 10. Experimental results of relay synchronization in a triplex network with non-identical layers, as a function of the intra-layer (σ ) and inter-layer (λ ) couplings. (Top) Colormap of the inter-layer synchronization e e errors between the outer layers E (left) and between one outer layer and the relay layer E (right) in the −1,1 0,1 σ − λ parameter space. The white contour line in each panel indicates the isoline for E and E respectively e e −1,1 0,1 equal to 0.12, error value taken as a reference. (Bottom) Inter-layer E , E synchronization errors as a −1,1 0,1 function of (left) λ for fixed σ = 0.5 (vertical continuous lines in the above panels) and (right) σ for fixed e e e λ = 0.5 (horizontal dashed lines). the PC for further analysis. Next, the coupling between the inter-layer (λ ) increases one step (0.01), digital pulses are sent to the potentiometers corresponding to that coupling and decreases the resistance in 100Ω each time it passes through this state, until the maximum value of λ is reached (Ω in potentiometers). When the entire run is finished, σ is increased by one step, and another cycle of λ is initiated. The whole procedure is repeated until each coupling reached its maximum value. The experiment is controlled from a PC with the LabVIEW software. The experimental results are summarized in Fig.  10. The top panels represent the averaged experimental inter-layer synchronization error for the outer layers E (left) and between the relay and one of the outer layers −1,1 E (right), for all the experimental range of intra-layer σ = [0, 0.6] and inter-layers λ = [0, 0.6] couplings. Even 0,1 e e though the system is unavoidably ae ff cted by noise and parameter mismatch in the electronic components, for high enough λ the value of E is well below E and therefore the inter-layer relay synchronization is verified e −1,1 0,1 in our experimental setup. Superimposed to the colormaps, we also have drawn the isoline for E = 0.12 in both panels (white lines), showing that the threshold λ value for which E and E are below the value of the isoline e −1,1 0,1 is always smaller in the E case. The fact that the perfect synchronization between the two outer layers is never −1,1 achieved agrees with our numerical predictions reported in the dynamical robustness section and cleary visual- ized in Fig. 8. For a clearer view, in the bottom left panel we have just plotted E and E as a function of λ for a fixed −1,1 0,1 e intra-layer coupling σ = 0.5 (corresponding to the blue and black dashed lines in the respective colormap panels in the upper part of Fig. 10), showing that E monotonically goes to zero and is always below E . −1,1 0,1 Finally, in the bottom-right panel of Fig. 10 we plot both errors, E and E , as a function of σ for a fixed −1,1 0,1 e value of the intra-layer coupling λ = 0.5 (vertical cuts in red and magenta in the colormap plots). That is done in order to show the effect of increasing the interaction in the intra-layer connectivity. Similarly to what observed in Fig. 5, promoting the topological difference between layers as σ increases rises the synchronization threshold. SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 9 www.nature.com/scientificreports/ Discussion The synchronous behavior of groups of units in a complex system is often a signature of a common functional involvement. Theoretically, this has been studied in the framework of modular networks , where nodes densely 48–55 connected among them in the same mesoscale structure usually share similar functions , or in the context of cluster synchronization associated to the existence of network symmetries and where modularity is not rele- 35,56 vant . However, most of these studies disregard the possible different nature of the involved nodes and links, a common feature in real systems. For instance, in the brain coherence is observed involving die ff rent cell types and electrical/chemical synapses. In addition, cluster and modular synchronization requires the direct connection between the synchro- nizing nodes. However, long distance coherence between complex mirrored structures mediated through non-synchronous differentiated ones plays a key role in the functioning of several real-world systems. Zero-lag 37,38 synchronization has been indeed observed between distant areas of the cortex , with the thalamus acting as a coordinator, and the transcendental role of symmetry in its dynamics has been lately pointed out in other contexts 57,58 like in the evolution of complex developmental systems . In this work, we have overcame these limitations by using a multilayer description in the study of distant synchronization in heterogeneous ensembles. Within this framework we have accounted for, besides considering different coupling functions between the dynamical units, the impact of having topologically different layers and heterogeneity in the node dynamics. We have implemented the concept of relay synchronization to the case of a multiplex network, showing that the intermediation of a relay layer can lead to inter-layer synchronization of a set of paired layers, both topologically and dynamically different from the transmitter. The phenomenon can be extended to indefinitely higher order relay configurations, provided a mirror symmetry is preserved in the multiplex. The coherent state is very robust to changes in the dynamics, topology, and even to strong multiplex disconnection. In this latter scenario, we show that the lower degree nodes in the synchronized outer layers are responsible for resilience of the synchronous state, while hubs can be safely de-mutiplexed. Finally, we experi- mentally validated our results in a multiplex network of three layers of electronic oscillators. Our results provide a new path for starting the study of the role of symmetries in setting long distance coherence in real systems. References 1. Barahona, M. & Pecora, L. M. Synchronization in small-world systems. Physical Review Letters 89, 054101 (2002). 2. Boccaletti, S., Latora, V., Moreno, Y., Chavez, M. & Hwang, D.-U. Complex networks: Structure and dynamics. Phys. Rep. 424, 175–308 (2006). 3. Arenas, A., Daz-Guilera, A., Kurths, J., Moreno, Y. & Zhou, C. Synchronization in complex networks. Phys. Rep. 469, 93–153 (2008). 4. Sorrentino, F. & Ott, E. Network synchronization of groups. Physical Review E 76, 056114 (2007). 5. Irving, D. & Sorrentino, F. Synchronization of dynamical hypernetworks: Dimensionality reduction through simultaneous block- diagonalization of matrices. Phys. Rev. E 86, 056102 (2012). 6. De Domenico, M. et al. Mathematical formulation of multilayer networks. Phys. Rev. X 3, 041022 (2013). 7. Boccaletti, S. et al. The structure and dynamics of multilayer networks. Phys. Rep. 544, 1–122 (2014). 8. Sorrentino, F. Synchronization of hypernetworks of coupled dynamical systems. New J. Phys. 14, 33035 (2012). 9. Aguirre, J., Sevilla-Escoboza, R., Gutiérrez, R., Papo, D. & Buldú, J. M. Synchronization of interconnected networks: The role of connector nodes. Phys. Rev. Lett. 112, 248701 (2014). 10. Gutiérrez, R., Sendiña-Nadal, I., Zanin, M., Papo, D. & Boccaletti, S. Targeting the dynamics of complex networks. Sci. Rep. 2, 396 (2012). 11. Zhang, X., Boccaletti, S., Guan, S. & Liu, Z. Explosive Synchronization in Adaptive and Multilayer Networks. Phys. Rev. Lett. 114, 038701 (2015). 12. Nicosia, V., Skardal, P. S., Arenas, A. & Latora, V. Collective phenomena emerging from the interactions between dynamical processes in multiplex networks. Phys. Rev. Lett. 118, 138302 (2017). 13. del Genio, C. I., Gómez-Gardeñes, J., Bonamassa, I. & Boccaletti, S. Synchronization in networks with multiple interaction layers. Science Advances 2 (2016). 14. Tang, L., Wu, X., Lü, J., Lu, J.-a. & D’Souza, R. M. Master stability functions for multiplex networks. arXiv preprint arXiv:1611.09110 (2016). 15. Huang, L., Park, K., Lai, Y.-C., Yang, L. & Yang, K. Abnormal synchronization in complex clustered networks. Phys. Rev. Lett. 97, 164101, https://doi.org/10.1103/PhysRevLett.97.164101 (2006). 16. Wang, X., Huang, L., Lai, Y.-C. & Lai, C. H. Optimization of synchronization in gradient clustered networks. Phys. Rev. E 76, 056113, https://doi.org/10.1103/PhysRevE.76.056113 (2007). 17. Guan, S., Wang, X., Lai, Y.-C. & Lai, C.-H. Transition to global synchronization in clustered networks. Phys. Rev. E 77, 046211, https://doi.org/10.1103/PhysRevE.77.046211 (2008). 18. Ma, X., Huang, L., Lai, Y.-C., Wang, Y. & Zheng, Z. Synchronization-based scalability of complex clustered networks. Chaos: An Interdisciplinary Journal of Nonlinear Science 18, 043109, https://doi.org/10.1063/1.3005782 (2008). 19. Louzada, V. H. P., Araújo, Na. M., Andrade, J. S. & Herrmann, H. J. Breathing synchronization in interconnected networks. Sci. Rep. 3, 3289 (2013). 20. Jalan, S. & Singh, A. Cluster synchronization in multiplex networks. EPL (Europhysics Letters) 113, 30002 (2016). 21. Singh, A., Jalan, S. & Boccaletti, S. Interplay of delay and multiplexing: Impact on cluster synchronization. Chaos 27, 043103 (2017). 22. Gambuzza, L. V., Frasca, M. & Gómez-Gardeñes, J. Intra-layer synchronization in multiplex networks. EPL 110, 20010 (2015). 23. Sevilla-Escoboza, R. et al. Inter-layer synchronization in multiplex networks of identical layers. Chaos 26, 065304 (2016). 24. Leyva, I. et al. Inter-layer synchronization in non-identical multi-layer networks. Phys. Rep. 7, 45475 (2017). 25. Fischer, I. et al. Zero-lag long-range synchronization via dynamical relaying. Phys. Rev. Lett. 97, 123902 (2006). 26. Bergner, A. et al. Remote synchronization in star networks. Phys. Rev. E 85, 026208 (2012). 27. Banerjee, R. et al. Enhancing synchrony in chaotic oscillators by dynamic relaying. Phys. Rev. E 85, 027201 (2012). 28. Gutiérrez, R. et al. Generalized synchronization in relay systems with instantaneous coupling. Phys. Rev. E 88, 052908 (2013). 29. Guillery, R. & Sherman, S. M. Thalamic Relay Functions and Their Role in Corticocortical Communication: Generalizations from the Visual System. Neuron 33, 163–175 (2002). 30. Sherman, S. M. The thalamus is more than just a relay. Current Opinion in Neurobiology 17, 417–422 (2007). 31. Mitchell, A. S. et al. Advances in Understanding Mechanisms of a Th lamic Relays in Cognition and Behavior. Journal of Neuroscience 34 (2014). 32. Vlasov, V. & Bifone, A. Hub-driven remote synchronization in brain networks. Phys. Rep. 7, 10403 (2017). SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 10 www.nature.com/scientificreports/ 33. Gambuzza, L. V., Frasca, M., Fortuna, L. & Boccaletti, S. Inhomogeneity induces relay synchronization in complex networks. Phys. Rev. E 93, 042203 (2016). 34. Nicosia, V., Valencia, M., Chavez, M., Daz-Guilera, A. & Latora, V. Remote Synchronization Reveals Network Symmetries and Functional Modules. Physical Review Letters 110, 174102 (2013). 35. Pecora, L. M., Sorrentino, F., Hagerstrom, A. M., Murphy, T. E. & Roy, R. Cluster synchronization and isolated desynchronization in complex networks with symmetries. Nat. Commun. 5, 4079 (2014). 36. Zhang, L., Motter, A. E. & Nishikawa, T. Incoherence-Mediated Remote Synchronization. Phys. Rev. Lett. 118, 174102 (2017). 37. Konigqt, P. & Singer, W. Visuomotor integration is associated with zero time-lag synchronization among cortical areas. Nature 385, 157 (1997). 38. Soteropoulos, D. S. & Baker, S. N. Cortico-cerebellar coherence during a precision grip task in the monkey. Journal of Neurophysiology 95, 1194–1206 (2006). 39. Erdös, P. & Rényi, A. On random graphs I. Publ. Math. Debrecen 6, 290–297 (1959). 40. Barabási, A.-L. & Albert, R. Emergence of scaling in random networks. Science 286, 509–512 (1999). 41. Rössler, O. An equation for continuous chaos. Phys. Lett. 57, 397–398 (1976). 42. Pecora, L. M. & Carroll, T. L. Master stability functions for synchronized coupled systems. Physical Review Letters 80, 2109 (1998). 43. Liu, Y.-Y., Slotine, J.-J. & Barabási, A.-L. Controllability of complex networks. Nature 473, 167–173 (2011). 44. Skardal, P. S. & Arenas, A. Control of coupled oscillator networks with application to microgrid technologies. Science Advances 1, e1500339 (2015). 45. Sevilla-Escoboza, R. & Buldú, J. M. Synchronization of networks of chaotic oscillators: Structural and dynamical datasets. Data in Brief 7, 1185–1189 (2016). 46. Sevilla-Escoboza, R. et al. Enhancing the stability of the synchronization of multivariable coupled oscillators. Phys. Rev. E 92, 032804 (2015). 47. Fortunato, S. Community detection in graphs. Phys. Rep. 486, 75–174, https://doi.org/10.1016/j.physrep.2009.11.002 (2010). 48. Arenas, A., Daz-Guilera, A. & Pérez-Vicente, C. J. Synchronization reveals topological scales in complex networks. Phys. Rev. Lett. 96, 114102, https://doi.org/10.1103/PhysRevLett.96.114102 (2006). 49. Boccaletti, S., Ivanchenko, M., Latora, V., Pluchino, A. & Rapisarda, A. Detecting complex network modularity by dynamical clustering. Phys. Rev. E 75, 045102, https://doi.org/10.1103/PhysRevE.75.045102 (2007). 50. Gómez-Gardeñes, J., Moreno, Y. & Arenas, A. Paths to synchronization on complex networks. Phys. Rev. Lett. 98, 034101, https:// doi.org/10.1103/PhysRevLett.98.034101 (2007). 51. Li, D. et al. Synchronization interfaces and overlapping communities in complex networks. Phys. Rev. Lett. 101, 168701, https://doi. org/10.1103/PhysRevLett.101.168701 (2008). 52. Almendral, J. A., Criado, R., Leyva, I., Buldú, J. M. & Sendiña-Nadal, I. Introduction to focus issue: Mesoscales in complex networks. Chaos: An Interdisciplinary Journal of Nonlinear Science 21, 016101, https://doi.org/10.1063/1.3570920 (2011). 53. Gutiérrez, R. et al. Emerging meso- and macroscales from synchronization of adaptive networks. Phys. Rev. Lett. 107, 234103, https://doi.org/10.1103/PhysRevLett.107.234103 (2011). 54. Rad, A. A. et al. Topological measure locating the effective crossover between segregation and integration in a modular network. Phys. Rev. Lett. 108, 228701, https://doi.org/10.1103/PhysRevLett.108.228701 (2012). 55. Prignano, L. & Daz-Guilera, A. Extracting topological features from dynamical measures in networks of kuramoto oscillators. Phys. Rev. E 85, 036112, https://doi.org/10.1103/PhysRevE.85.036112 (2012). 56. Sorrentino, F., Pecora, L. M., Hagerstrom, A. M., Murphy, T. E. & Roy, R. Complete characterization of the stability of cluster synchronization in complex dynamical networks. Science Advances 2, e1501737 (2016). 57. Ruvinsky, I. & Gibson-Brown, J. Genetic and developmental bases of serial homology in vertebrate limb evolution. Development 127, 5233–5244 (2000). 58. Smith, R. S. et al. A plausible model of phyllotaxis. Proc. Natl. Acad. Sci. USA 103, 1301–1306, https://doi.org/10.1073/ pnas.0510457103 (2006). Acknowledgements Work partly supported by the Spanish Ministry of Economy under projects FIS2013-41057-P and FIS2017- 84151-P. Authors acknowledge the computational resources and assistance provided by CRESCO, the supercomputing center of ENEA in Portici, Italy. R.S.E. acknowledges support from Secretaría de Educación Pública, PRODEP, grant number UDG-PTC-1289-DSA/103.5/16/10313 and SEP-CONACYT/CB-2016-01, grant number 285909. Author Contributions I.L., I.S.N., P.C. and S.B. conceived the study and devised the model. R.S. and V.V.A. carried out the experiments. I.L., I.S.N. and P.C. carried out the numerical simulations. All Authors wrote the Manuscript. Additional Information Competing Interests: The authors declare no competing interests. Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Cre- ative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not per- mitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/. © The Author(s) 2018 SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 11 http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Scientific Reports Springer Journals
Free
11 pages
Loading next page...
 
/lp/springer_journal/relay-synchronization-in-multiplex-networks-9HVaebE2mn
Publisher
Nature Publishing Group UK
Copyright
Copyright © 2018 by The Author(s)
Subject
Science, Humanities and Social Sciences, multidisciplinary; Science, Humanities and Social Sciences, multidisciplinary; Science, multidisciplinary
eISSN
2045-2322
D.O.I.
10.1038/s41598-018-26945-w
Publisher site
See Article on Publisher Site

Abstract

www.nature.com/scientificreports OPEN Relay synchronization in multiplex networks 1,2 1,2 3 3 4 I. Leyva , I. Sendiña-Nadal , R. Sevilla-Escoboza , V. P. Vera-Avila , P. Chholak & S. Boccaletti Received: 16 February 2018 Relay (or remote) synchronization between two not directly connected oscillators in a network is an Accepted: 21 May 2018 important feature allowing distant coordination. In this work, we report a systematic study of this Published: xx xx xxxx phenomenon in multiplex networks, where inter-layer synchronization occurs between distant layers mediated by a relay layer that acts as a transmitter. We show that this transmission can be extended to higher order relay configurations, provided symmetry conditions are preserved. By first order perturbative analysis, we identify the dynamical and topological dependencies of relay synchronization in a multiplex. We find that the relay synchronization threshold is considerably reduced in a multiplex configuration, and that such synchronous state is mostly supported by the lower degree nodes of the outer layers, while hubs can be de-multiplexed without affecting overall coherence. Finally, we experimentally validated the analytical and numerical findings by means of a multiplex of three layers of electronic circuits. Synchronization is one of the most important collective phenomena in nature. It can be observed in natural, 1–3 social and technological systems, and it became one of the most active research topics in network science . e Th huge amount of new data collected in the last years has permitted a higher resolution network representation of real systems. In particular, the inclusion of new features shaped multi-layer representations, i.e. approaches in which the network units are arranged in several layers, each one accounting for a different kind of interactions 4–7 among the nodes . Multi-layer structures determine scenarios where novel forms of synchronization are rele- 8,9 vant. Despite an analytical approach has been tackled in just a few particular cases , several synchronization sce- narios have been already addressed, as unidirectional coordination between layers , explosive synchronization 11,12 emerging from the interactions between dynamical processes in multiplex networks , complete synchroniza- 13,14 15–21 22 23,24 tion , cluster synchronization , intra-layer or inter-layer synchronization. Very recently, relay (RS) and remote synchronization (two very well known phenomena in chains, or small motifs, of coupled oscillators) have captured the attention of researchers. This form of synchronization is observed when two units of a network (identical or slightly different) synchronize despite not being directly linked, and due instead to the intermediation of a relay mismatched unit. The phenomenon has been experimentally detected in 25 26,27 lasers and circuits . In general, the relay units exhibit generalized or delay synchronization with the units they actually pace to synchrony . RS is of outstanding relevance in the brain: the thalamus acts as a relay between distant cortical areas through 29–32 the thalamo-cortical pathways, playing the role of a coordination hub that maintains the information flow . Complex structures and neuronal dynamics are implicated in this process involving not only simple, but higher 30,31 order relay paths, that transfer the information through multiple-step relay chains . Recently, remote syn- chronization has been addressed in the context of complex networks , revealing the extremely important role of 34–36 network structural and dynamical symmetries in the appearance of distant synchronization , as it was already 37,38 suggested by the observation of zero-lag delays between mirror areas of the brain . Nevertheless, the interplay between symmetry, dynamics and multi-layer structure remains still mostly unexplored. In this work, we perform a systematic study of inter-layer relay synchronization in a multiplex network, where distant layers synchronize their dynamics while their intra-layer motion remains unsynchronized. We con- sider generic high-order structures where multi-site relay pathways are verified. The dynamical and topological 1 2 Complex Systems Group & GISC, Universidad Rey Juan Carlos, Móstoles, Madrid, 28933, Spain. Center for Biomedical Technology, Universidad Politécnica de Madrid, Pozuelo de Alarcón, Madrid, 28223, Spain. Centro Universitario de los Lagos, Universidad de Guadalajara, Jalisco, 47460, Mexico. Department of Mechanical Engineering, Indian Institute of Technology Bombay, Powai, Mumbai, 400076, India. CNR-Institute of complex systems, Via Madonna del Piano 10, Sesto Fiorentino, 50019, Italy. Correspondence and requests for materials should be addressed to I.L. (email: inmaculada.leyva@gmail.com) SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 1 www.nature.com/scientificreports/ Figure 1. Schematic representation of a multiplex of 2M + 1 layers (here M = 2) labeled as k = −M, …, −1, 0, 1, …, M where each pair of layers k and −k (painted with the same color) are networks of identical oscillators with the same topology  and intra-layer coupling σ and whose dynamical state is described by the variable U and −k U , respectively. The multiplex is symmetric with respect to the layer k = 0 and the nodes are coupled to their replicas in the rest of layers with an inter-layer coupling strength λ. dependencies of the phenomena are studied, using perturbation stability analysis. The robustness of the relay synchronization against de-multiplexing the layers is reported, revealing the key role of low degree nodes in maintaining the layers coordination. Finally, the findings are experimentally validated in a multiplex network of electronic circuits. Results Model. We start by considering 2M + 1 layers (or networks), arranged as shown in Fig. 1. Each layer k, with k =− MM ,, …… 0, , , is formed by N dynamical systems (each of which being m-dimensional), whose states are kk k k km represented by the column vectors U =… {, uu ,, u }, with u ∈  , i = 1, …, N, and whose intra-layer inter- 12 N i k k actions are encoded through the Laplacian matrices  = {}  . The layer stack is symmetric with respect to k = 0 ij k −k in such a way that Laplacians  and  have the same structure. The dynamical systems are also paired: nodes +k −k belonging to layers U and U are identical to each other, and different (in some parameter) from the rest of the layers. Consequently, layer k = 0 has no counterpart, and acts as a relay between all layers situated above and below it. Layers are coupled in a multiplex cong fi uration, and the dynamical evolution of the system is described by the following equations: qk =+1 ≤M k kk qk UF =− () UG σλ ()  ⊗+ UH ()  ⊗− () UU k k N qk =−1 ≥−M (1) k k k kT mm where the functions FU () =… [( fu ), fu (),, fu ()] (with f :  → representing the vectorial func- k k 12 k kN k tions evolving each dynamical unit), are identical for the same |k|. G, H are the m × m matrices representing respectively the linear intra- (G) and inter- (H) layer coupling schemes.  is the N × N identity matrix, σ is the N k intra-layer coupling strength within layers k and −k, and λ is the inter-layer coupling strength. +k −k Due to the reflection symmetry of the system under study (i.e. as long as the U and U layers are identi- cal for all k), a synchronous inter-layer evolution (with layers evolving in a pairwise synchronized fashion, i.e. + −k k k′ where U = U ) at all k without necessarily implying U = U for k ≠ k′) is a solution of Eq. (1), independently of intra-layer synchronization (i.e. independently on whether the state of the systems within each layers are SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 2 www.nature.com/scientificreports/ k +k −k synchronized). Let therefore δU (t) = U (t) − U (t), with k = 1, …, M be the vector describing the difference between the dynamics of the paired layers. kk k k Considering the smallness of and expanding around the inter-layer solution up to δδ Uu =… {, δδ uu ,, } 12 N first order, one obtains a set of M linearized vector equations for the perturbations δU : k k ∼∼ k k k δσ Uf =− [( JJ UG )( ⊗− )( λδ 2) − J HU ()]δU k k kM qk =+1 +δ λ J HU () U qk =−1qk ≠ (2) where J denotes the Jacobian operator, δ is the Kronecker delta accounting for the boundary condition at k = M kM ±M (as the stack end layers U are only connected to the previous neighbor layer). The vector U = {} u describes k −k 0 the dynamical state of any of the k = 0, …, M layers at the synchronous state U = U ≠ U and, therefore, the whole dynamics is reduced to the dynamics of M + 1 layers. Such evolution at the node level is given by: qk =+1 ≤M k k k q k uf  =− () ug σλ  () uh  +− [(uh )(u )] ∑∑ i ki k ij j i i j qk =−10 ≥ (3) where i = 1 …, N is the node index, k, q = 0, 1, …, M, and g(u) = Gu and h(u) = Hu are the projections of the inter- and intra- later coupling operators to the node level. Notice that, since each paired layers k and −k is kk − inter-layer synchronized () UUU == , each layer acts therefore as a relay to the rest of the stack. The prob- lem consists now in solving MmN linear equation (2), together with solving in parallel the (M + 1)mN nonlinear equation (3) for . Although the total number of equations to compute is still of the same order as in Eq. (1), that u is (2M + 1)mN, the fact that MmN of those equations are linear results in a much faster computation of the dynamics. Then, calculating the maximum Lyapunov exponent (MLE) transverse to the manifold U as a function of the system parameters actually gives the necessary conditions for the stability of the synchronous solution: whenever MLE0 < , perturbations transverse to the manifold will die out, and the multi-relay synchronous solu- tion will be stable. In order to monitor the synchronization error between layers, we define the inter-layer synchronization error as, q k E =− lim uu () tt () dt, qk , ∑ i i T →∞T 0 (4) i=1 where ||·|| stands for the Euclidean norm and q, k are the layers’ indexes, such that E denotes the inter-layer −k,k synchronization error of mirror layers. Without lack of generality, In our numerical simulations we consider two types of topologies where layers are either (i) Erdös-Renyi (ER) or (ii) scale-free (SF, generated by means of the Barabási-Albert’s algorithm ), in all cases with N = 500. We classify the layer stacks regarding the topology sequence of each layer. For instance, a triplex where the three layers have ER topology will be denoted as EEE, and a system where two identical SF layers are mediated by a center ER will be denoted as SES. The nodes are chaotic Rössler oscillators , defined by the m = 3 state vector u = (x, y, z) whose autonomous evolution is given by f () uf == () u [, −− yz xa +. yz ,0 2( +− x 9)] and the heterogeneity between the layers is introduced kk − k through the parameter a , such that each layer develops a different chaotic dynamics. In our case study, the intra- T T and inter- layer coupling functions are set to be g() uG == u (0,0,) z and hu () == Hu (0,, y 0) respec- tively. These coupling schemes ensure that intra-layer synchronization is prevented when layers are isolated and not multiplexed (class I layers, according to the standard master stability function (MSF) classification established in ref. ) whereas multiplexed nodes along the layers can synchronize for a coupling strength λ above a given threshold (class II MSF) . Layers with identical topology. With the aim of determining whether relay synchronization can be achieved in a multiplex cong fi uration let us first consider the multiplex structure den fi ed by Eq. ( 1) for the case of 01 −1 three identical SF layers  ==  and where the parameters a = a = 0.2 for the outer layers and a = 0.3 1 −1 0 for the relay units of the central layer, although different selections of these parameters and topologies produce a similar behavior. Results are collected in Fig. 2, where the synchronization error between the outer layers E is plotted versus −1,1 the inter-layer coupling λ for several values of the intra-layer couplings σ and σ in the outer and relay layers 1 0 respectively, with σ = σ . In all cases, there is a critical coupling λ above which complete synchronization 1 0 between layers k = 1 and k = −1 occurs, that is, E = 0 is achieved, while the relay layer (k = 0) remains unsyn- −1,1 chronized to any of the two outer layers (k = 1, −1) as shown in the inset where E > 0 for any parameter choice. 0,1 In addition, the calculation of the corresponding MLE given by Eq. (2) (bottom panel of Fig. 2) confirms −1 1 * that the relay synchronous solution U = U reaches stability (MLE < 0) at the same critical λ where the error between the outer layers is zero, as indicated by the vertical lines. Therefore, one can conclude that inter-layer MLE is a useful tool for reducing the system’s dimensionality and use it for evaluation of the critical inter-layer coupling λ from now on. In order to better understand the different roles played by external and relay layers, we show in Fig.  3 the criti- cal inter-layer coupling value in the parameter region (σ , σ ), that is, when the intra-layer coupling σ is different 0 1 k SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 3 www.nature.com/scientificreports/ Figure 2. Relay synchronization in a triplex (M = 1) with identical SF layers (SSS configuration). (Main panel) Synchronization error between the outer layers (k = −1, and k = 1) E (see Eq. 4) as a function of the inter- −1,1 layer coupling λ for three different values of the intra-layer coupling σ = σ (see legend). The inset shows the 0 1 corresponding synchronization errors between the relay and one of the outer layers. (Bottom panel) Maximum 1 −1 Lyapunov exponent (MLE) of the relay synchronization manifold U = U as a function of λ for the same cases as in the main panel. Vertical lines mark the points where the MLE becomes negative. All points are averages of 10 network realizations with N = 500 and ⟨k⟩ = 4. See the main text for the relay and outer layer Rössler oscillators specifications. Figure 3. Relay synchronization in a triplex network with identical SF layers as a function of the intra-layer couplings for the relay (σ ) and outer (σ ) layers. (Left) Color map of the inter-layer coupling threshold λ for the 0 1 relay state (E = 0 and E ≠ 0) in the σ − σ parameter space. (Right) Inter-layer coupling threshold λ for the −1,1 0,1 0 1 relay state as a function of the coupling strength in the relay layer σ for a fixed value of σ = 1 (red dashed line in 0 1 left panel) and as a function of the coupling strength in the outer layers σ for a fixed value of σ = 1 (black dashed 1 0 line in left panel). Each point is an average of 10 SF network realizations with N = 500 and ⟨k⟩ = 8. for the relay and outer layers. It can be seen that the system’s ability to synchronize is practically unaltered with σ , while increasing σ makes the value of λ to drop drastically. This therefore reveals that multiplex relay synchroni - zation is much more sensitive to changes ae ff cting the mirror layers than to those arising in the transmission layer. Our results can be generalized to any number of layers. As an example, we report also the case M = 2, which corresponds to two outer layers above (k = 1, 2) and below (k = −1, −2) the relay layer (k = 0). We choose a = a = 0.2 and a = a = 0.3, and a = 0.25 for the central layer. The results stand for any other parameter −1 1 −2 2 0 choice. In Fig. 4 we plot the inter-layer synchronization errors E (void markers) and E (full markers), vs. −1,1 −2,2 the inter-layer coupling λ for several values of the intra-layer coupling σ. As in the triplex case, the critical λ at which complete inter-layer synchronization is achieved depends on σ, but it is the same for both pairs of layers, as E and E drop to zero simultaneously. In the inset we plot the inter-layer synchronization errors between −1,1 −2,2 the non-paired layers, E , E to check that they remain mutually incoherent. Therefore, we have verified that 0,1 1,2 relay synchronization also occurs in cascade for arbitrarily high-order multiplex systems, provided a structural and dynamical symmetry is conserved. Layers with non-identical topology. So far, we have dealt with multiplexes of pairwise identical layers. However, this condition is too strong a limitation to hope that it would capture and properly represent the case of many real systems. The next step needed for generalization is studying then the relay synchronization scenario in SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 4 www.nature.com/scientificreports/ Figure 4. Relay synchronization in a pentaplex (M = 2) with identical N = 500 ER layers (EEEEE configuration). The synchronization error between the two pair of outer layers E (empty symbols) and E −1,1 −2,2 (full symbols) is shown as a function of λ for three different values of the intra-layer coupling σ, being σ = σ , ∀k. The inset shows the synchronization errors between each one of the outer layers and the relay layer. The results are averaged over 10 different network realizations and initial conditions. Figure 5. Relay synchronization in a triplex with different layers. Inter-layer coupling threshold λ vs the intra- layer couplings σ = σ for (a) a mixed ER-SF-ER (ESE) and identical (EEE) configurations and (b) a mixed SF- 0 1 ER-SF (SES) and identical (SSS) configurations. the case in which the topology of the relay layer differs from that of the outer layers. In Fig.  5 we have reported the critical inter-layer coupling λ in two heterogeneous triplex cases: (a) a pair of Erdös-Rényi layers mediated by a scale-free relay layer (ESE situation) and (b) the opposite case where SF layers are connected through a ER layer (SES). Each case is compared with the topologically homogeneous EEE and SSS structures, respectively. For the sake of simplification and of better assessment of the role of the topology, we keep σ = σ . 0 1 Figure 5(a) shows that, for a large range of intra-layer couplings, the mediation of a SF relay facilitates the synchronization between the paired layers, since λ in the ESE case (void blue circles) is smaller than the one corresponding to the homogeneous case (EEE, full blue circles). On the contrary, a relay ER layer intermediat- ing between two outer SF layers (Fig. 5(b)) does not determine a significant difference as long as the intra-layer coupling strength is low, but increases the threshold λ for higher σ, as compared to the homogeneous SSS case. Robustness. In the previous Sections we have addressed the dependence of relay synchronization in a multi- plex on the dynamical and structural layer heterogeneity, and proved that the phenomenon still holds even when the intermediate layer has a completely different structure and dynamics than the mirrored ones. The present section is devoted instead to assess the robustness of relay synchronization against structural changes by means of a de-multiplexing process of the layers, that is, against performing a progressively shutting down of the inter-layer links such that a fraction of nodes in each layer is not linked to their counterparts in the other layers. In addition, we also investigated several cases to test the robustness against mismatches in the oscillators parameters to closely resemble real experimental conditions. Structural robustness. To closely check this process, we initially consider a 3-layer multiplex with identical topol- ogy (EEE or SSS). We choose the inter- and intra-layer couplings to guarantee a relay synchronous state with the layers fully multiplexed. Then, we proceed to disconnect one by one the inter-layer links according to the nodes degree ranking, both in the ascending and the descending order, and re-evaluate in every step the state of the relay synchronization by measuring the E error. Figure 6(a) reports the averaged evolution of E as a function −1,1 −1,1 of the number of multiplexed nodes, aer h ft aving performed the whole de-multiplexing process for 10 different network realizations. It can be seen that, starting from a situation with E = 0, the EEE multiplex configuration −1,1 (blue void markers) loses the synchronization immediately with just a few of inter-layer links being removed. On SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 5 www.nature.com/scientificreports/ Figure 6. Structural robustness of the network relay synchronization for identical layers. (a) Synchronization error between the outer layers E vs the decreasing number of connected relay lines for identical ER (blue −1,1 empty symbols) and SF (black solid symbols) layers. Relay lines are disconnected following a descending (circle symbols) or ascending (square symbols) node degree ranking of the outer layers (seed legend in (b)). Parameter values are N = 500, ⟨k⟩ = 8, λ = 0.23 and σ = σ = 0.8. (b) Number of multiplexed relay lines needed 0 1 to support a relay network state as a function of the intra-layer coupling strength σ = σ while keeping constant 0 1 λ = 0.23. The different curves are explained in the legend. All the results are averaged over 10 different network realizations. the other hand, relay synchronization is resilient in SSS triplex configurations also when more than 30% of the nodes are not multiplexed. A more detailed view can be obtained from Fig. 6(b), where the number of multiplexed nodes needed to support the relay synchronization is represented as a function of the intra-layer coupling σ = σ . As expected, 0 1 when the coupling is weak, all the N nodes need to be linked to preserve relay synchronization. However, as the interaction within the layers increases, the intra-layer connectivity helps to maintain a synchronous state despite an increasing number of nodes are being de-multiplexed without damaging the coherence between the outer layers. In Fig. 6(b), we can see that for both the EEE (blue void markers) and the SSS (black full markers) tri- plex configurations, removing the links between layers connecting nodes with higher degree (descending degree ranking, circle markers) is much more robust than following an ascending degree ranking (square markers). This is indeed a very interesting result: relay synchronization in a multiplex network is supported by the low degree nodes, while the hubs can be safely disconnected without perturbing the transmission. This is notably evidenced in the SSS case (black full squares) where after having removed the 40% of the inter-layer links connecting the highest degree nodes, the relay synchronization is still supported by the multiplex structure connected through the lowest degree nodes. Once we have singled out that, from our trial rankings, the descending degree ranking is the most conven- ient way to de-multiplex part of the network without loosing coherence, we proceed our study by evaluating the impact of having a relay layer with different topology from the outer layers, as we did in the previous Section 4. In this scenario, we have two possible descending degree rankings, the one dictated by the relay layer and the one dictated by the outer layers. The results are summarized in Fig.  7 where we plot, as in Fig. 6(b), the number of nodes that need to be linked to maintain synchronization as a function of σ = σ . For the sake of comparison, 0 1 we added the curves for the homogeneous EEE and SSS (full markers) multiplex configurations, together with the data for the mixed ESE and SES (void markers) layers. Notice that the chosen inter-layer coupling λ = 0.23 is well above threshold for all the cases, as it can be derived from Fig. 5. All the reported evidence indicates that the introduction of a relay layer with a topology different from that of the outer layers has little influence on the min- imum number needed to support the relay synchronization, as long as the first removed inter-layer connections correspond to the hubs in the outer layers (blue and back curves). Curiously, the alternative of using the relay layer topology to rank the degree of the nodes, destroys the coherence between the outer layers as soon as a tiny fraction of links is removed (red curves). er Th efore, the relay synchronization in a multiplex is very unstable if just a few links connecting nodes which are hubs in the relay layer are removed. Notice that this unlinking criterion is equivalent to randomly disconnect the multiplex. Therefore, the robustness of the relay synchrony relies mainly in the low degree nodes of the external layers. The relevance of the low degree nodes in controlling the dynamics 43,44 of complex networks has been pointed out in other contexts . Dynamical robustness. In order to explore a more realistic situation of a multiplex composed of non identical oscillators, we have taken further the generalization by introducing some heterogeneity in the node’s parameter a , that is, a = a . The node values of the outer layers’s were randomly picked from the interval a = [0.20, 0.21]. k k k,i ±1,i We prepared three different setups depending on how the parameter setting is introduced in the outer layers: i) a = a but aa ≠ , that is, oscillators are not identical within each external layer but the two external 1,i −1,i 1,ij 1, layers are symmetric; ii) a = a but a ≠ a , that is, layers are composed of identical oscillators but there is a 1,i 1,j 1,i −1,i mismatch between oscillators in layer k = 1 and their corresponding replicas in layer k = −1; and iii) SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 6 www.nature.com/scientificreports/ Figure 7. Structural robustness of the network relay synchronization for non identical layers. Number of multiplexed relay lines needed to support a relay network state as a function of the intra-layer coupling strength σ = σ while keeping constant the inter-layer coupling strength λ = 0.23 for mixed ER (ESE, empty circles) and 0 1 mixed SF (SES, empty squares) layer configurations. The relay lines are disconnected following a descending order of the outer node degrees and for comparison the corresponding values for identical layers -see Fig. 6(b) are plotted in solid symbols. The red solid (ESE) and empty (SES) triangles show the behavior when the relay lines are disconnected following the degree ranking of the relay layer. Figure 8. Dynamical robustness of the network relay synchronization as a function of the heterogeneity in the a parameter. Synchronization error between the outer layers E vs the inter-layer coupling λ for the three k −1,1 cases reported in the main text with a =. [0 20,0.21]: i) mismatched layers with oscillators within the outer ±1,i layers being identical but mismatched with their replicas (green squares); ii) symmetrical outer layers with nonidentical oscillators within the layers (full cyan circles); and iii) not symmetrical layers with nonidentical oscillators (purple diamonds). Inset is a zoom of the main plot. Simulations have been conducted with SF layers of identical topology with N = 500 and ⟨k⟩ = 8. a ≠ a ≠ a , that is, combining cases i) and ii), a fully random case, not preserving the symmetry neither the 1,i 1,j −1,i parameter homogeneity in each layer. Results are summarized in Fig. 8, where the synchronization error between the outer layers is plotted vs the inter-layer coupling λ. For comparison, the perfectly homogenous case with all layers having identical oscilla- tors (red circles) is also included. In the inset, a zoom of the main figure, it is clear that breaking the symmetry between the outer layers slightly deteriorates the synchronization, as it occurs in cases i) (green squares) and iii) (purple diamonds). However, introducing some heterogeneity in the a parameter (full cyan circles) in the layer’s oscillators does not modify the synchronization threshold at all, as long as each oscillator in one outer layer is identical to its replica in the mirrored layer. We will actually show that our experimental results can be framed within this latter case. Experimental validation Finally, we report experimental evidence of relay synchronization in a multiplex of nonlinear electronic circuits, with the setup sketched in Fig. 9 (left). The array is made of 21 Rössler-like circuits arranged in three layers of 7 nodes, with the relay layer having different topology as the outer layers. Each layer has two different electronic couplers, one for the coupling among nodes in the same layer (σ ) and the second for the interaction of each node with its replica in the other layers (λ ). The chaotic dynamics of the circuits is well approximated by the three variables (v , v , v ) obeying : 1 2 3 SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 7 www.nature.com/scientificreports/ Figure 9. (Left) Schematic representation of the experimental arrangement of three layers of electronic circuits. e b Th idirectional coupling is adjusted by means of three strips of digital potentiometers X9C103 (XDCP), the resistance is controlled through digital pulses sent by a DAQ (NI USB 6363). (Right) Graph structure used for the upper and lower layers (top) and for the relay layer (bottom). C = 1 nF C = 1 nF C = 1 nF σe,λe = [0−0.6] 1 2 3 R = 2MΩ R2 = 200kΩ R3 = 10kΩ R4 = 100kΩ R = 50kΩ R6 = 5MΩ R7 = 100kΩ R8 = 10kΩ R = 10kΩ R10 = 100kΩ R11 = 100kΩ R12 = 150kΩ R = 68kΩ R14 = 10kΩ R15 = 75kΩ R16 = 120kΩ 0 1 Rk =Ω 50 Rk =Ω 35 I = 0.7 V = 15 d ee c c Table 1. Parameter values of the chaotic dynamics of one Rössler like circuit as described in Eq. (5).   1 R R k k 1 k 1 k   v =− v ++ v v  1i 1i 2i 3i RC  R R   11 2 4 1 R 1 k k k −− σ av v e ∑ ij () 11 j i RC R 11 15 j=1     RR RR 1  k  k k 68 68    v =− −+ v 1 − v 2i 1i 2i  k  RC  RR RR  62  97  c 7     q=1   1   q k − λ vv − e ∑ 22 i i RC  R 62  16 q=−1    k k k 10  v =− −+ Gv v  3i () 13 i i RC  R    10 3 11 (5) where G is a nonlinear gain function given by: 1i   R R   14  14 0, if vI ≤ 1 +  + V   1id ee    R R   13 13 Gv () = 1i      R R R R R R  12 12 12 12 14 14      vV −− ee I + ,1 if vI > +  + V   1id 1id ee       R R  R R  R R      14 13 13 14 13 13 where the parameter values are gathered in Table 1. The reader is referred to ref. for a detailed description of the 9,23,24,46 experimental implementation of the Rössler-like circuit in the networks, and refs for previous realizations in different network configurations. Both the intra-layer σ and the inter-layer λ are set by means of the digital e e potentiometers X9C103, that working as voltage divisor for the maximum resistance (10k Ω), σ and λ is set to e e zero, this potentiometers are controlled through the digital ports (P0.0, P0.1, P0.2, P0.3) of a DAQ card. First that all we send all the coupling value to zero, aer 500 ft ms takes the sample of the time series of each networks, all the variables v of each oscillator enter to the DAQ card through the analogue ports (AI0, AI1, …, AI20) and saved in 2i SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 8 www.nature.com/scientificreports/ Figure 10. Experimental results of relay synchronization in a triplex network with non-identical layers, as a function of the intra-layer (σ ) and inter-layer (λ ) couplings. (Top) Colormap of the inter-layer synchronization e e errors between the outer layers E (left) and between one outer layer and the relay layer E (right) in the −1,1 0,1 σ − λ parameter space. The white contour line in each panel indicates the isoline for E and E respectively e e −1,1 0,1 equal to 0.12, error value taken as a reference. (Bottom) Inter-layer E , E synchronization errors as a −1,1 0,1 function of (left) λ for fixed σ = 0.5 (vertical continuous lines in the above panels) and (right) σ for fixed e e e λ = 0.5 (horizontal dashed lines). the PC for further analysis. Next, the coupling between the inter-layer (λ ) increases one step (0.01), digital pulses are sent to the potentiometers corresponding to that coupling and decreases the resistance in 100Ω each time it passes through this state, until the maximum value of λ is reached (Ω in potentiometers). When the entire run is finished, σ is increased by one step, and another cycle of λ is initiated. The whole procedure is repeated until each coupling reached its maximum value. The experiment is controlled from a PC with the LabVIEW software. The experimental results are summarized in Fig.  10. The top panels represent the averaged experimental inter-layer synchronization error for the outer layers E (left) and between the relay and one of the outer layers −1,1 E (right), for all the experimental range of intra-layer σ = [0, 0.6] and inter-layers λ = [0, 0.6] couplings. Even 0,1 e e though the system is unavoidably ae ff cted by noise and parameter mismatch in the electronic components, for high enough λ the value of E is well below E and therefore the inter-layer relay synchronization is verified e −1,1 0,1 in our experimental setup. Superimposed to the colormaps, we also have drawn the isoline for E = 0.12 in both panels (white lines), showing that the threshold λ value for which E and E are below the value of the isoline e −1,1 0,1 is always smaller in the E case. The fact that the perfect synchronization between the two outer layers is never −1,1 achieved agrees with our numerical predictions reported in the dynamical robustness section and cleary visual- ized in Fig. 8. For a clearer view, in the bottom left panel we have just plotted E and E as a function of λ for a fixed −1,1 0,1 e intra-layer coupling σ = 0.5 (corresponding to the blue and black dashed lines in the respective colormap panels in the upper part of Fig. 10), showing that E monotonically goes to zero and is always below E . −1,1 0,1 Finally, in the bottom-right panel of Fig. 10 we plot both errors, E and E , as a function of σ for a fixed −1,1 0,1 e value of the intra-layer coupling λ = 0.5 (vertical cuts in red and magenta in the colormap plots). That is done in order to show the effect of increasing the interaction in the intra-layer connectivity. Similarly to what observed in Fig. 5, promoting the topological difference between layers as σ increases rises the synchronization threshold. SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 9 www.nature.com/scientificreports/ Discussion The synchronous behavior of groups of units in a complex system is often a signature of a common functional involvement. Theoretically, this has been studied in the framework of modular networks , where nodes densely 48–55 connected among them in the same mesoscale structure usually share similar functions , or in the context of cluster synchronization associated to the existence of network symmetries and where modularity is not rele- 35,56 vant . However, most of these studies disregard the possible different nature of the involved nodes and links, a common feature in real systems. For instance, in the brain coherence is observed involving die ff rent cell types and electrical/chemical synapses. In addition, cluster and modular synchronization requires the direct connection between the synchro- nizing nodes. However, long distance coherence between complex mirrored structures mediated through non-synchronous differentiated ones plays a key role in the functioning of several real-world systems. Zero-lag 37,38 synchronization has been indeed observed between distant areas of the cortex , with the thalamus acting as a coordinator, and the transcendental role of symmetry in its dynamics has been lately pointed out in other contexts 57,58 like in the evolution of complex developmental systems . In this work, we have overcame these limitations by using a multilayer description in the study of distant synchronization in heterogeneous ensembles. Within this framework we have accounted for, besides considering different coupling functions between the dynamical units, the impact of having topologically different layers and heterogeneity in the node dynamics. We have implemented the concept of relay synchronization to the case of a multiplex network, showing that the intermediation of a relay layer can lead to inter-layer synchronization of a set of paired layers, both topologically and dynamically different from the transmitter. The phenomenon can be extended to indefinitely higher order relay configurations, provided a mirror symmetry is preserved in the multiplex. The coherent state is very robust to changes in the dynamics, topology, and even to strong multiplex disconnection. In this latter scenario, we show that the lower degree nodes in the synchronized outer layers are responsible for resilience of the synchronous state, while hubs can be safely de-mutiplexed. Finally, we experi- mentally validated our results in a multiplex network of three layers of electronic oscillators. Our results provide a new path for starting the study of the role of symmetries in setting long distance coherence in real systems. References 1. Barahona, M. & Pecora, L. M. Synchronization in small-world systems. Physical Review Letters 89, 054101 (2002). 2. Boccaletti, S., Latora, V., Moreno, Y., Chavez, M. & Hwang, D.-U. Complex networks: Structure and dynamics. Phys. Rep. 424, 175–308 (2006). 3. Arenas, A., Daz-Guilera, A., Kurths, J., Moreno, Y. & Zhou, C. Synchronization in complex networks. Phys. Rep. 469, 93–153 (2008). 4. Sorrentino, F. & Ott, E. Network synchronization of groups. Physical Review E 76, 056114 (2007). 5. Irving, D. & Sorrentino, F. Synchronization of dynamical hypernetworks: Dimensionality reduction through simultaneous block- diagonalization of matrices. Phys. Rev. E 86, 056102 (2012). 6. De Domenico, M. et al. Mathematical formulation of multilayer networks. Phys. Rev. X 3, 041022 (2013). 7. Boccaletti, S. et al. The structure and dynamics of multilayer networks. Phys. Rep. 544, 1–122 (2014). 8. Sorrentino, F. Synchronization of hypernetworks of coupled dynamical systems. New J. Phys. 14, 33035 (2012). 9. Aguirre, J., Sevilla-Escoboza, R., Gutiérrez, R., Papo, D. & Buldú, J. M. Synchronization of interconnected networks: The role of connector nodes. Phys. Rev. Lett. 112, 248701 (2014). 10. Gutiérrez, R., Sendiña-Nadal, I., Zanin, M., Papo, D. & Boccaletti, S. Targeting the dynamics of complex networks. Sci. Rep. 2, 396 (2012). 11. Zhang, X., Boccaletti, S., Guan, S. & Liu, Z. Explosive Synchronization in Adaptive and Multilayer Networks. Phys. Rev. Lett. 114, 038701 (2015). 12. Nicosia, V., Skardal, P. S., Arenas, A. & Latora, V. Collective phenomena emerging from the interactions between dynamical processes in multiplex networks. Phys. Rev. Lett. 118, 138302 (2017). 13. del Genio, C. I., Gómez-Gardeñes, J., Bonamassa, I. & Boccaletti, S. Synchronization in networks with multiple interaction layers. Science Advances 2 (2016). 14. Tang, L., Wu, X., Lü, J., Lu, J.-a. & D’Souza, R. M. Master stability functions for multiplex networks. arXiv preprint arXiv:1611.09110 (2016). 15. Huang, L., Park, K., Lai, Y.-C., Yang, L. & Yang, K. Abnormal synchronization in complex clustered networks. Phys. Rev. Lett. 97, 164101, https://doi.org/10.1103/PhysRevLett.97.164101 (2006). 16. Wang, X., Huang, L., Lai, Y.-C. & Lai, C. H. Optimization of synchronization in gradient clustered networks. Phys. Rev. E 76, 056113, https://doi.org/10.1103/PhysRevE.76.056113 (2007). 17. Guan, S., Wang, X., Lai, Y.-C. & Lai, C.-H. Transition to global synchronization in clustered networks. Phys. Rev. E 77, 046211, https://doi.org/10.1103/PhysRevE.77.046211 (2008). 18. Ma, X., Huang, L., Lai, Y.-C., Wang, Y. & Zheng, Z. Synchronization-based scalability of complex clustered networks. Chaos: An Interdisciplinary Journal of Nonlinear Science 18, 043109, https://doi.org/10.1063/1.3005782 (2008). 19. Louzada, V. H. P., Araújo, Na. M., Andrade, J. S. & Herrmann, H. J. Breathing synchronization in interconnected networks. Sci. Rep. 3, 3289 (2013). 20. Jalan, S. & Singh, A. Cluster synchronization in multiplex networks. EPL (Europhysics Letters) 113, 30002 (2016). 21. Singh, A., Jalan, S. & Boccaletti, S. Interplay of delay and multiplexing: Impact on cluster synchronization. Chaos 27, 043103 (2017). 22. Gambuzza, L. V., Frasca, M. & Gómez-Gardeñes, J. Intra-layer synchronization in multiplex networks. EPL 110, 20010 (2015). 23. Sevilla-Escoboza, R. et al. Inter-layer synchronization in multiplex networks of identical layers. Chaos 26, 065304 (2016). 24. Leyva, I. et al. Inter-layer synchronization in non-identical multi-layer networks. Phys. Rep. 7, 45475 (2017). 25. Fischer, I. et al. Zero-lag long-range synchronization via dynamical relaying. Phys. Rev. Lett. 97, 123902 (2006). 26. Bergner, A. et al. Remote synchronization in star networks. Phys. Rev. E 85, 026208 (2012). 27. Banerjee, R. et al. Enhancing synchrony in chaotic oscillators by dynamic relaying. Phys. Rev. E 85, 027201 (2012). 28. Gutiérrez, R. et al. Generalized synchronization in relay systems with instantaneous coupling. Phys. Rev. E 88, 052908 (2013). 29. Guillery, R. & Sherman, S. M. Thalamic Relay Functions and Their Role in Corticocortical Communication: Generalizations from the Visual System. Neuron 33, 163–175 (2002). 30. Sherman, S. M. The thalamus is more than just a relay. Current Opinion in Neurobiology 17, 417–422 (2007). 31. Mitchell, A. S. et al. Advances in Understanding Mechanisms of a Th lamic Relays in Cognition and Behavior. Journal of Neuroscience 34 (2014). 32. Vlasov, V. & Bifone, A. Hub-driven remote synchronization in brain networks. Phys. Rep. 7, 10403 (2017). SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 10 www.nature.com/scientificreports/ 33. Gambuzza, L. V., Frasca, M., Fortuna, L. & Boccaletti, S. Inhomogeneity induces relay synchronization in complex networks. Phys. Rev. E 93, 042203 (2016). 34. Nicosia, V., Valencia, M., Chavez, M., Daz-Guilera, A. & Latora, V. Remote Synchronization Reveals Network Symmetries and Functional Modules. Physical Review Letters 110, 174102 (2013). 35. Pecora, L. M., Sorrentino, F., Hagerstrom, A. M., Murphy, T. E. & Roy, R. Cluster synchronization and isolated desynchronization in complex networks with symmetries. Nat. Commun. 5, 4079 (2014). 36. Zhang, L., Motter, A. E. & Nishikawa, T. Incoherence-Mediated Remote Synchronization. Phys. Rev. Lett. 118, 174102 (2017). 37. Konigqt, P. & Singer, W. Visuomotor integration is associated with zero time-lag synchronization among cortical areas. Nature 385, 157 (1997). 38. Soteropoulos, D. S. & Baker, S. N. Cortico-cerebellar coherence during a precision grip task in the monkey. Journal of Neurophysiology 95, 1194–1206 (2006). 39. Erdös, P. & Rényi, A. On random graphs I. Publ. Math. Debrecen 6, 290–297 (1959). 40. Barabási, A.-L. & Albert, R. Emergence of scaling in random networks. Science 286, 509–512 (1999). 41. Rössler, O. An equation for continuous chaos. Phys. Lett. 57, 397–398 (1976). 42. Pecora, L. M. & Carroll, T. L. Master stability functions for synchronized coupled systems. Physical Review Letters 80, 2109 (1998). 43. Liu, Y.-Y., Slotine, J.-J. & Barabási, A.-L. Controllability of complex networks. Nature 473, 167–173 (2011). 44. Skardal, P. S. & Arenas, A. Control of coupled oscillator networks with application to microgrid technologies. Science Advances 1, e1500339 (2015). 45. Sevilla-Escoboza, R. & Buldú, J. M. Synchronization of networks of chaotic oscillators: Structural and dynamical datasets. Data in Brief 7, 1185–1189 (2016). 46. Sevilla-Escoboza, R. et al. Enhancing the stability of the synchronization of multivariable coupled oscillators. Phys. Rev. E 92, 032804 (2015). 47. Fortunato, S. Community detection in graphs. Phys. Rep. 486, 75–174, https://doi.org/10.1016/j.physrep.2009.11.002 (2010). 48. Arenas, A., Daz-Guilera, A. & Pérez-Vicente, C. J. Synchronization reveals topological scales in complex networks. Phys. Rev. Lett. 96, 114102, https://doi.org/10.1103/PhysRevLett.96.114102 (2006). 49. Boccaletti, S., Ivanchenko, M., Latora, V., Pluchino, A. & Rapisarda, A. Detecting complex network modularity by dynamical clustering. Phys. Rev. E 75, 045102, https://doi.org/10.1103/PhysRevE.75.045102 (2007). 50. Gómez-Gardeñes, J., Moreno, Y. & Arenas, A. Paths to synchronization on complex networks. Phys. Rev. Lett. 98, 034101, https:// doi.org/10.1103/PhysRevLett.98.034101 (2007). 51. Li, D. et al. Synchronization interfaces and overlapping communities in complex networks. Phys. Rev. Lett. 101, 168701, https://doi. org/10.1103/PhysRevLett.101.168701 (2008). 52. Almendral, J. A., Criado, R., Leyva, I., Buldú, J. M. & Sendiña-Nadal, I. Introduction to focus issue: Mesoscales in complex networks. Chaos: An Interdisciplinary Journal of Nonlinear Science 21, 016101, https://doi.org/10.1063/1.3570920 (2011). 53. Gutiérrez, R. et al. Emerging meso- and macroscales from synchronization of adaptive networks. Phys. Rev. Lett. 107, 234103, https://doi.org/10.1103/PhysRevLett.107.234103 (2011). 54. Rad, A. A. et al. Topological measure locating the effective crossover between segregation and integration in a modular network. Phys. Rev. Lett. 108, 228701, https://doi.org/10.1103/PhysRevLett.108.228701 (2012). 55. Prignano, L. & Daz-Guilera, A. Extracting topological features from dynamical measures in networks of kuramoto oscillators. Phys. Rev. E 85, 036112, https://doi.org/10.1103/PhysRevE.85.036112 (2012). 56. Sorrentino, F., Pecora, L. M., Hagerstrom, A. M., Murphy, T. E. & Roy, R. Complete characterization of the stability of cluster synchronization in complex dynamical networks. Science Advances 2, e1501737 (2016). 57. Ruvinsky, I. & Gibson-Brown, J. Genetic and developmental bases of serial homology in vertebrate limb evolution. Development 127, 5233–5244 (2000). 58. Smith, R. S. et al. A plausible model of phyllotaxis. Proc. Natl. Acad. Sci. USA 103, 1301–1306, https://doi.org/10.1073/ pnas.0510457103 (2006). Acknowledgements Work partly supported by the Spanish Ministry of Economy under projects FIS2013-41057-P and FIS2017- 84151-P. Authors acknowledge the computational resources and assistance provided by CRESCO, the supercomputing center of ENEA in Portici, Italy. R.S.E. acknowledges support from Secretaría de Educación Pública, PRODEP, grant number UDG-PTC-1289-DSA/103.5/16/10313 and SEP-CONACYT/CB-2016-01, grant number 285909. Author Contributions I.L., I.S.N., P.C. and S.B. conceived the study and devised the model. R.S. and V.V.A. carried out the experiments. I.L., I.S.N. and P.C. carried out the numerical simulations. All Authors wrote the Manuscript. Additional Information Competing Interests: The authors declare no competing interests. Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Cre- ative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not per- mitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/. © The Author(s) 2018 SCIeNtIfIC REpo R tS | (2018) 8:8629 | DOI:10.1038/s41598-018-26945-w 11

Journal

Scientific ReportsSpringer Journals

Published: Jun 5, 2018

References

You’re reading a free preview. Subscribe to read the entire article.


DeepDyve is your
personal research library

It’s your single place to instantly
discover and read the research
that matters to you.

Enjoy affordable access to
over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month

Explore the DeepDyve Library

Search

Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly

Organize

Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.

Access

Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals.

Your journals are on DeepDyve

Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more.

All the latest content is available, no embargo periods.

See the journals in your area

DeepDyve

Freelancer

DeepDyve

Pro

Price

FREE

$49/month
$360/year

Save searches from
Google Scholar,
PubMed

Create lists to
organize your research

Export lists, citations

Read DeepDyve articles

Abstract access only

Unlimited access to over
18 million full-text articles

Print

20 pages / month

PDF Discount

20% off