The central role of peripheral nodes in directed network dynamics

Many social, technological, and biological systems with asymmetric interactions display a variety of collective phenomena, such as opinion formation and synchronization. This has motivated much research on the dynamical impact of local and mesoscopic structure in directed networks. However, the unique constraints imposed by the global organization of directed networks remain largely undiscussed. Here, we control the global organization of directed Erdős–Rényi networks, and study its impact on the emergence of synchronization and ferromagnetic ordering, using Kuramoto and Ising dynamics. In doing so, we demonstrate that source nodes – peripheral nodes without incoming links – can disrupt or entirely suppress the emergence of collective states in directed networks. This effect is imposed by the bow-tie organization of directed networks, where a large connected core does not uniquely ensure the emergence of collective states, as it does for undirected networks.


Results
Structure. While the global organization of directed networks (cf. Fig. 1) is characteristically feedforward (from sources to sinks), the overall connectivity of the core is feedback. On the one hand, CORE nodes can both dynamically adjust to and influence the state of their neighbors (through local pair-wise interactions), but on the other, they are subject to the influence of sources. Here, we study how the number of SOURCE and CORE nodes, and the inter-connectivity of the IN and CORE components affect the balance between the core's ability to support collective dynamics (through feedback) and the influence of sources. To this end, we consider a simple toy model based on directed Erdős-Rényi (ER) networks, where the structural features of interest are determined by the mean in-degree 〈 〉 q in 3 (i.e. the mean number of incoming links), starting with the emergence of the core in the network at 〈 〉 = q 1 in (percolation point). When determining how sources influence the collective dynamics of the core, we must consider how CORE nodes connect to each other, the number of SOURCE nodes relative to that of CORE nodes, and how the former connect to the latter. Connections from SOURCE to CORE nodes may be direct (through direct links) or indirect (through paths to other IN nodes that are themselves directly connected to CORE nodes). For any directed network where the IN component is entirely composed of SOURCE nodes, the dynamics of the SOURCE nodes are formally equivalent to external fields, which act directly on any number of CORE nodes. In the presence of indirect connections from SOURCE nodes to CORE nodes, the above-mentioned equivalence may no longer apply, depending both on the intra-connectivity of the IN component and the dynamical rules themselves. Here, we simply focus on the fraction of SOURCE nodes -N SOURCE /N -relative to the fraction of CORE nodes -N CORE /N CORE I N O U T SOURCE SINK Figure 1. Bow-tie architecture of a directed network, composed of a core and its periphery. The CORE is the largest strongly connected component (largest set of nodes reachable from each other through a sequence of directed links). The periphery comprises the IN and OUT components -the sets of nodes in sequences of directed links leading into and out of the CORE respectively 1-3,43 -and a hierarchy of tendrils (such as those connecting into and out of the OUT and IN components, respectively) and tubes (directly connecting the IN and OUT components) 36 . As indicated by the arrowheads, the overall connectivity of this architecture is feedforward, from SOURCE nodes (IN nodes without in-links) to SINK nodes (OUT nodes without out-links).
www.nature.com/scientificreports www.nature.com/scientificreports/ -and on the number of links between the IN and CORE components -L IN-CORE /L -in ER networks with N nodes and = 〈 〉 L N q in links. Above the percolation point -〈 〉 = q 1 in -the fraction of CORE nodes -N CORE /N -increases monotonically with 〈 〉 q in , as shown in Fig. 2(a), and eventually comprises the majority of the network's nodes. However, the fraction of SOURCE nodes -N SOURCE /N -remains of the same order of magnitude as the fraction of CORE nodes - in , as confirmed both through simulations and analytic calculations, and shown in Fig. 2(a). The influence of this finite fraction of SOURCE nodes is ultimately exerted through a fraction of IN-CORE links -L IN-CORE /L -which is also of the same order of magnitude as the fraction of CORE-CORE links -L CORE-CORE /L -i.e. non-vanishingly small, as shown in Fig. 2(b). Individual CORE nodes adjust their state both through links from other CORE nodes, and from IN nodes, so that their dynamics are partly determined by the balance between their external in-degree q in ext (the number of in-links from IN nodes) and their internal in-degree q in int (the number of in-links from other CORE nodes), where = + q q q in in in ext i nt . For a typical (average) CORE node, we may consider the ratio between the average external and internal in-degree as a measure of the node's susceptibility to the influence of the IN component. This is equivalent to analyzing the density of the inter-connectivity between IN and CORE components, relative to the density of the CORE's own connectivity. Symbolically, , and may therefore significantly impact the collective dynamics of the CORE given specific dynamical rules and parameters. In addition, the brief survey of directed networks presented in Table 1 below suggests that these structural features are equally significant (non-vanishingly small) in a broad range of directed networks.
Dynamics. Next, we investigated how the above-discussed structural features impact the emergence of synchronization and magnetization in the Kuramoto and Ising models, by controlling the corresponding dynamical parameters -the coupling strength K and the temperature T -along with the structural control parameter -the mean in-degree 〈 〉 q in .
Kuramoto Model. The Kuramoto model is a phenomenological model of synchronization between coupled phase oscillators 25,[33][34][35] . On a directed network with N nodes, the instantaneous change in the phase θ n of oscillator n is governed by:   where ω n is the oscillator's natural frequency, K is the coupling strength, and A mn is the network's adjacency matrix element, which is 1 if node m has an out-link to node n, and 0 otherwise. The state of the system is described by the complex order parameter ( ) is the average phase. For sufficiently large K, the system becomes partially synchronized: a small finite fraction of oscillators lock at frequency Ω, such that their individual phases θ = Ω + t c n n for any phase shift ψ π ψ π − < < + /2 c / 2 n . The phases of unlocked oscillators simply drift, i.e. they are time-wise uncorrelated with the average phase. The extent of synchronization (ordering) in the system is conventionally described by the amplitude of the complex order parameter, the order parameter By definition, the order parameter in equation (4) is an average over all N nodes in the network, so it is easily redefined over a subset X with N X nodes (e.g. the CORE), and interpreted as the extent of synchronization over the subset, denoted r X . In the interest of consistency, this notation will be used throughout the present work for all quantities of interest.
The large-time average of equation (4) is presented in Fig. 3  www.nature.com/scientificreports www.nature.com/scientificreports/ networks is presented in Fig. 3(b), where it is shown that the extent of synchronization is limited by the topology of the network when  K 1. The data in Fig. 3(a) show that synchronization emerges simultaneously in the network, the CORE, and the OUT component at 〈 〉 .  q 1 5 in . This suggests that synchronization is supported by the CORE, as expected from the latter's feedback connectivity, and is further substantiated by the increase in synchronization with 〈 〉 > . q 1 5 in , as the CORE spans an increasingly larger fraction of the network, and becomes more densely connected. In addition, Fig. 3(a) shows that the IN component is unsynchronized at all values of 〈 〉 q in , in the limit of large coupling (  K 1). In this limit, the extent of synchronization in the network is determined by 〈 〉 q in , as shown in Fig. 3(b). Thus, the overall feedforward connectivity of the network, the absence of synchronization in the IN component, and the increase in synchronization among CORE nodes with 〈 〉 > . q 1 5 in -as the IN component becomes less densely connected to the CORE -suggest that the IN component may disrupt synchronization in the CORE, and by extension, completely suppress synchronization in the CORE for ≤ 〈 〉 < .
, the network does not have a CORE, and synchronization between a finite fraction of the network's nodes is therefore impossible.
Ising Model. The Ising model describes ferromagnetism through the interaction of spins s n and s m at sites n and m. Here, we analyze the time evolution of the dynamics of a spin system in contact with a heat bath at temperature T 30 . Each spin s n experiences the effect of its neighboring spins as a local field  in is magnetized at any temperature < T T C , as shown in Fig. 4(a). In fact, T C increases monotonically with the CORE's link density, and is therefore characteristic of the internal connectivity of the CORE. The more densely connected the CORE, the larger the range of temperatures for which it can support magnetization, which is maximal in the limit when → T 0, independently of the CORE's connectivity, as shown in Fig. 4(b). In the presence of the IN component, the T C of the CORE is significantly different, and magnetization is disrupted/suppressed, as shown in Fig. 4(a,c), respectively. Although spins in the IN component are at the same temperature as the CORE, they can reduce the range of temperatures at which the CORE supports magnetization, or indeed suppress magnetization altogether, even in the limit where → T 0. The feedforward organization of the IN component, and the latter's impact on ordering among CORE nodes with Kuramoto and Ising dynamics hint at the role of SOURCE nodes in limiting the emergence of ordering in the CORE, as does the data in Fig. 4(b,d). Given the overall feedforward connectivity of the bow-tie architecture, we could expect the OUT component to be fully ordered in the absence of the IN component, driven by core dynamics. However, a comparison of Fig. 4(b,d) shows this is not the case. At 〈 〉 = .
q 1 1 in , magnetization is significantly lower in the OUT component than in the CORE. The size of the CORE decreases approaching the network's percolation point ( , while the number of tendrils increases 36 , including those that connect SOURCE nodes to OUT nodes (cf. Fig. 1), potentially explaining the above-mentioned discrepancy.
Pair-correlation functions: measuring the dynamical influence of SOURCE nodes. Whether direct or indirect, the influence of SOURCE nodes is experienced by individual CORE nodes as an external field through links from IN nodes. Any strategy to control both of these structural features is heavily dependent on the intra-connectivity of the IN component. For example, removing a SOURCE node may simply create one or more in its place, and removing any number of IN-CORE links does not control for the number of existing SOURCE nodes. Here, we indirectly control these features by simply removing fractions of IN nodes (f IN ) uniformly at random, and calculate the CORE's response to these structural changes. In the Kuramoto model, this response can be calculated in the rotating frame (at group velocity Ω), using the pair-correlation function for some observation time τ  1, for any time-dependent observable Q(t). In the Ising model, the pair-correlation function is similarly defined The dependence of the synchronization order parameter r CORE and C CORE on f IN is presented in Fig. 5(a,b) for unsynchronized networks, and in Fig. 5(c) for a partially synchronized network, i.e. for networks with 〈 〉 < . In networks with 〈 〉 < . q 1 5 in , synchronization appears only after the removal of a critical fraction f C IN of randomly selected IN nodes, as shown in Fig. 5(a,b). Once the CORE is partially synchronized, the removal of further IN nodes enhances synchronization, similarly to the removal of any fraction of IN nodes from networks with 〈 〉 > . q 1 5 in , where synchronization is linearly enhanced by the removal of randomly selected IN nodes, as shown in Fig. 5(c). The absence of any peaks in C CORE for partially synchronized networks and its decaying behavior with f IN further suggest that the observed enhancement in synchronization results from more CORE nodes being drawn into the existing synchronized group. In unsynchronized networks, the emergence of synchronization at f C IN is accompanied by a clear peak in C CORE , signaling the emergence of long-range order that is characteristic of a second order phase transition, in this case from an unsynchronized phase to a synchronized phase, and confirming that synchronization is supported by the CORE and disrupted by the IN component. For networks where 〈 〉 < .

Discussion
The results presented in the preceding sections show that the global organization of directed networks constrains the emergence of collective phenomena. Specifically, it was shown that the emergence of ordered states in the Kuramoto and Ising models is supported by the CORE and disrupted by the IN component. Here, we discuss how the uncorrelated dynamics of SOURCE nodes are responsible for this effect. By definition, SOURCE nodes do not have in-links, and are therefore unaffected by feedback (cf. Fig. 1), so that SOURCE nodes with randomly distributed initial states remain uncorrelated for all time, acting as a source of noise or fluctuations.
The exact output of the IN component is determined not only by the dynamics of SOURCE nodes, but also by its internal connectivity, which is beyond the scope of this discussion. Nonetheless, it is clear that the number of  www.nature.com/scientificreports www.nature.com/scientificreports/ direct links from SOURCE nodes to the CORE increases with 〈 〉 q in (cf. Fig. 2). Additionally, it has been recently shown that, in directed uncorrelated random complex networks, the number of nodes s in the finite in-component of any node scales as − − ⁎ s e s s 3/2 / , where s* is a characteristic parameter that depends on 〈 〉 q in 38 . This means that only a small number of IN nodes can be reached by traveling backwards from any other IN node, and at least one of these must be a SOURCE node. Above the network's percolation threshold, the number of SOURCE nodes is always of the same order as the total number of IN nodes. Given the effective absence of reciprocal links and structural correlations that characterize directed ER networks, the above information suggests that paths from individual SOURCE nodes to the CORE are largely disjoint, even when the links to the CORE are not direct.
Let us now consider the specific case of the Ising model. According to equation (5), each SOURCE node i experiences a null local field = h 0 i as a consequence of not having in-links. From equation (6), it then follows that SOURCE nodes are independently and identically found in spin up states ( = s 1 i ) and spin down states ( = − s 1 i ) with probability 1/2, acting as a source of fluctuations that will affect all downstream nodes, at any temperature T. In Fig. 5, the CORE is shown to fully magnetize in the absence of the IN component, when → T 0, so that all CORE nodes are in the same spin state. Upon restoring the IN component, and with it q in ext links to a given CORE node j, the latter will experience fluctuations in its local field h j . If the number of links from CORE nodes ≤ q q in int in ext , the fluctuations will inevitably lead to a configuration of spins where = h 0 j or < h 0 j . Both outcomes will frustrate the dynamics of CORE node j, which will change its spin state. While the exact probability of such an event is outside the scope of this discussion, in Fig. 6, we consider a simplified picture, where a single SOURCE node i is directly linked to a CORE node j, and = = q q 1 in int in ext , a frequent configuration in directed ER networks with ≤ 〈 〉 ≤ q 1 2 in . In the above simplified picture, the dynamics of node j are frustrated with probability 1/2, as depicted in Fig. 6(a,b). When node j is frustrated, it may then frustrate another CORE node k in a similar manner, as depicted in Fig. 6(c). Clearly, the fraction of CORE nodes which can be frustrated in this manner is dictated by the internal connectivity of the CORE, the number of links from IN nodes to individual CORE nodes, and the number of SOURCE nodes themselves. Given the similar dependence of magnetization on the above-mentioned structural features identified in the preceding section, frustration presents itself as the putative mechanism behind this dependence, driven by the intrinsic dynamics of SOURCE nodes. When a critical fraction of CORE nodes becomes frustrated, long-range order is broken, suppressing magnetization altogether. www.nature.com/scientificreports www.nature.com/scientificreports/ In the Kuramoto model, SOURCE nodes with randomly-distributed natural frequencies act as a source of incoherence: each SOURCE node i rotates steadily at its own natural frequency ω i , as a consequence of not having in-links, with a phase shift determined by its initial phase (cf. equation (2)). If the coupling strength K is sufficiently larger than the width of the frequency distribution ω g( ), any downstream node j with a single in-link from node i will be driven at frequency ω i . A simple stability analysis of equation (2) also shows that node j is driven at frequency ω i if the number of in-neighbors rotating at ω i exceeds the magnitude of the remaining in-neighbor's combined angular velocity. On the one hand, this is demonstrative of the difficulty in predicting the exact output of the IN component to the CORE, other than when the IN component is composed of linear disjoint chains or trees with distinct SOURCE nodes at their root. On the other hand, it reveals that a CORE node j can be prevented from joining a group of locked CORE nodes simply by receiving an appropriate number of links from IN nodes. The data presented in Fig. 5 show that the CORE fully synchronizes in the absence of the IN component, and by averaging equation (2) over the CORE it follows that its nodes are rotating at the mean natural frequency ω 〈 〉. So, upon restoring the IN component, a single SOURCE node i rotating at a frequency other than ω 〈 〉 can prevent a CORE node j from joining the synchronized group. Moreover, as ω i is narrowly distributed about ω 〈 〉, the desynchronizing influence directly experienced by different CORE nodes will also vary, which may compound the desynchronizing effect.

Ising model Kuramoto model
Regardless of the exact output of the IN component, it is clear that SOURCE nodes can directly or indirectly impact the collective dynamics of the CORE. SOURCE nodes are a characteristic feature of the global organization of directed networks, and in this work it was shown that, at least for directed ER networks, when the dynamics of SOURCE nodes are uncorrelated, magnetization and synchronization can be suppressed. This represents an additional constraint when compared to undirected networks, where these ordered states emerge at the network's percolation threshold, as summarized in Fig. 7.
Given the evident importance of SOURCE nodes in determining the state of the CORE, these findings reinforce their importance for the controllability and resilience of systems abstracted onto directed networks. For example, a synchronized brain network may be desynchronized through the random destruction of synapses, e.g. in Alzheimer's disease, as the random destruction of links alters the global organization of the network 38 . Similarly, a targeted attack on low in-degree centrality nodes may allow an attacker to exert considerable influence on a network, by controlling the dynamics of newly-created SOURCE nodes. In general, local pairwise dynamics on bow-tie architectures, which may include intrinsic (node-specific), external (field), and stochastic contributions 23 , are mapped to CORE dynamics under the direct or indirect action of the effective field created by SOURCE nodes.
Despite the broad class of dynamics that are susceptible to the influence of SOURCE nodes, the body of work demonstrating bow-tie organization in real systems, and the applications of dynamical models under the action of external fields, any generalizations regarding the dynamics of real systems must be made with care. Firstly, the existence of heterogeneities in the CORE of real networks, such as hubs (scale-free networks) and reciprocal connections, may modify the nature of critical transitions and the values of the critical parameters themselves 32 . However, for particular phase transitions, e.g. percolation, one can also note that structural features such as degree-degree correlations and clustering can change the critical point and critical exponents, but not the nature of the transition 39,40 . Secondly, there are dynamics, e.g. social dynamics, which are more aptly described by agent-based models 19,41,42 . Thirdly, there are further relevant theoretical aspects to consider, such as boundary conditions, the nature of the interactions themselves, and the existence of external driving fields, which are common in real-world systems. In some sense, sources and sinks are a boundary in bow-tie architectures, and , and is now equally likely to be in a spin up or down state, cf. equation (6). Box (b) represents the latter outcome, and box (c) a further step with a similar outcome, where node j frustrates another CORE node k.
the boundary conditions are therefore automatically defined by the dynamics of sources and sinks (regardless of whether these are internally or externally driven). Changing these boundary conditions may produce significantly different collective behavior. For example, if all SOURCE oscillators are assigned the same natural frequency, they may drive synchronization in the CORE at the same frequency, or even cause the emergence of other macroscopic dynamics, analogously to what happens in the Kuramoto model driven by an external field 26,27 . Likewise, different types of interaction can produce different behavior on the same network structure e.g. the antiferromagnetic Ising model displays spin-glass like behavior caused by long loops (cycles) 32 .
In conclusion, we considered a toy model based on directed Erdős-Rényi networks, where the global organization of the network is determined by the mean in-degree 〈 〉 q in . Considering both Ising and Kuramoto dynamics, we found that the global organization of directed networks constrains the emergence of magnetized and synchronized phases. Unlike in undirected Erdős-Rényi networks, where these collective phenomena emerge at the network's percolation point, in directed Erdős-Rényi networks, SOURCE nodes with uncorrelated dynamics act to disrupt the collective dynamics of the CORE, delaying the onset of ordered states above the network's percolation point, where 〈 〉 = q 1 in . Magnetization was confirmed to emerge for  〈 〉 . q 1 9 in , and synchronization found to emerge for  〈 〉 . q 1 5 in . The clear impact of SOURCE nodes, a topological feature only found in directed networks, highlights the need to consider the global organization of directed networks when considering their dynamics, robustness, and controllability. Based on the discussion in the preceding paragraph, we believe that source nodes will produce similar effects in strongly heterogeneous, degree-degree correlated, and clustered networks with a bow-tie architecture, but further investigations are required.

Methods
The ER toy-model is based on ensembles of 50 directed ER networks with mean in-degree 〈 〉 q in . Each network is generated from its undirected counterpart by first assigning a direction to all links ( → i j, where < i j) and then reversing it with probability 1/2. The undirected ER networks were generated by creating links between all possible labeled pairs of = N 10 5 nodes with probability 〈 〉 − q N 2 /( 1) in . Note that in this kind of ER network, two nodes are connected by a single directed link, and reciprocal links are absent by construction. Directed ER networks can also be built by forming any possible directed link with probability ~1/N. However, this will lead to the same network structure as above in the limit where  N 1, since the probability of forming reciprocal links becomes vanishingly small (~1/N 2 ). For an example of how bow-tie architectures can also be built with finite fractions of both reciprocal and single unidrectional links see 43 .
All structural and dynamical quantities presented are averaged over the above-mentioned ensemble. The time evolution of the Kuramoto and Ising models layered on each network were simulated from random initial conditions: the phase of each node θ n was drawn uniformly at random between −π and π (Kuramoto model), and the spin s n of each node was drawn binomially at random with probability = p 1/2 (Ising model). All time-averaged quantities were calculated over a time window in the steady state. The natural frequencies ω in the Kuramoto model were drawn normally at random, with mean ω 〈 〉 = 0 and standard deviation σ = 1. The results of analytical calculations presented in Fig. 2 were obtained following standard methods (see 38 and the references therein for further details.) In particular, for any directed uncorrelated random network with N nodes, and an in-degree (q in ) out-degree (q out ) distribution P q q ( , ) in out , the fraction of SOURCE nodes N SOURCE /N can be calculated from the probability of randomly selecting a SOURCE node. For any node with out-degree q out , this is equal to the probability of selecting a node with = q 0 in , given by = P q q ( 0, ) in out , and that at least one of its out-links leads to the CORE, given by − y 1 c q out , where y c is the probability that an out-link leads to a finite component. Taking into account all possible degrees, where y c can be determined self-consistently 38 . The average number of IN-CORE links can also be calculated using this formalism, considering a node randomly selected with probability P q q ( , ) in out , and the probability − that it receives m in-links from nodes in finite components, given its in-degree, the probability that an in-link comes from a finite component x c 38 , and the probability that at least one out-link leads to a CORE node − y 1 c q out . To ensure that such a node belongs to the CORE, it is sufficient to ensure its m in-links from the IN component account at most for q i − 1 of its total in-links, i.e. that at least one of its in-links is from another CORE node. The average number of links 〈 〉 m may then be explicitly calculated, and the number of IN-CORE links = 〈 〉 -L Nm IN CORE . Summing over m (with some in-between simplifications), where x c , like y c , is also determined self-consistently.