A structural transition in physical networks.

In many physical networks, including neurons in the brain1,2, three-dimensional integrated circuits3 and underground hyphal networks4, the nodes and links are physical objects that cannot intersect or overlap with each other. To take this into account, non-crossing conditions can be imposed to constrain the geometry of networks, which consequently affects how they form, evolve and function. However, these constraints are not included in the theoretical frameworks that are currently used to characterize real networks5-7. Most tools for laying out networks are variants of the force-directed layout algorithm8,9-which assumes dimensionless nodes and links-and are therefore unable to reveal the geometry of densely packed physical networks. Here we develop a modelling framework that accounts for the physical sizes of nodes and links, allowing us to explore how non-crossing conditions affect the geometry of a network. For small link thicknesses, we observe a weakly interacting regime in which link crossings are avoided via local link rearrangements, without altering the overall geometry of the layout compared to the force-directed layout. Once the link thickness exceeds a threshold, a strongly interacting regime emerges in which multiple geometric quantities, such as the total link length and the link curvature, scale with the link thickness. We show that the crossover between the two regimes is driven by the non-crossing condition, which allows us to derive the transition point analytically and show that networks with large numbers of nodes will ultimately exist in the strongly interacting regime. We also find that networks in the weakly interacting regime display a solid-like response to stress, whereas in the strongly interacting regime they behave in a gel-like fashion. Networks in the weakly interacting regime are amenable to 3D printing and so can be used to visualize network geometry, and the strongly interacting regime provides insights into the scaling of the sizes of densely packed mammalian brains.

To lay out physical networks, we must arrange the links and the nodes in such a way to avoid crossing each other, while minimizing the total link length. In other words, we must find the shortest path for each link, even when the straight path is obstructed by other nodes and links, a problem equivalent to stretching a rubber band between flexible obstacles (See Fig. 1, also SI 3.1 for proof [13]). To achieve this, we propose a model in which the forces governing the motion of the nodes and links is determined by the gradient of the total potential energy, where V is the total elastic potential of all links l = 1, ..., L. Each link is modeled as an elastic cylinder with radius r, experiencing both internal elastic forces and short-range external repulsive forces from other links and nodes. V captures the node-link interactions at link endpoints; the non-crossing condition is ensured by a short-range repulsive force in V (node-node interaction) and V (link-link interaction) modeled as short-range Gaussian potentials whose strength is set by A and A. In (1) s is the length parameter of link l and represents the position of a point along the center of the link at time t; is the position of node i; r is the range of the node-node repulsive force; k is the elastic constant of the links. The potential energy (1) is inspired by models employed in self-avoiding polymer chains [14] and manifold dynamics [15], but given the constraints induced by the network structure, (1) has different terms and displays behavior unique to networks. With V = 0 and replacing V with the elastic energy of a spring, (1) reduces to the potential energy of FDL with short-range node repulsion. The lowest energy solution of (1) can lead to sharp bending of some links, which we avoid using a Gay-Berne potential [16] employed in polymer physics (SI 4). Finally, we embed the network in a high viscosity medium, allowing it to relax to a low energy state without oscillations. Therefore, the node and link positions follow the first order gradient descent equations of motion, where λ and λ are the node and link friction constants (SI 3.6). We use FDL to set the initial node positions and explore two versions of the model with different constraints: (i) In the Elastic Link Model (ELI), which corresponds to the node positions are kept fixed and only the links are allowed to reorganize themselves; (ii) In the Fully Elastic Model (FUEL) we assume hence both nodes and links are free to move.

Figure 1

Modeling Framework for Avoiding Link Crossings:

(A) We model each link as a stretched, flexible rubber band, corresponding to many short springs connected to each other, pulled apart by elastic forces F. The links exert a repulsive force F on each other that falls sharply at radii larger than r. While in FDL the links cross each other (left figure), in ELI and FUEL such crossings are prohibited (right figure). (B) A small network with N = 6 nodes laid out with FDL (left), resulting in multiple link crossings (red links). The right plot shows how the network laid out by ELI that resolves the crossings. (C) Finding the final layout of a lattice with also showing the evolution of the total link length over time during the simulation. We started from a random layout and used simulated annealing to find the final layout. The thermal noise helps links pass through each other and resolve crossings.

The network defined by (1)–(3) has an uneven potential energy landscape [17] with a very large number of local minima, hence identifying the globally optimal configuration is NP hard (SI 3.7). We therefore use simulated annealing [18] to approach an energetically favorable local minimum (SI 3.7). The computational complexity of the model is discussed in SI 8.3. Figure 1C shows how FUEL finds the correct 3D configuration of a lattice, helped by the noise to tunnel through the finite potential walls and escape local minima.

As FDL ignores the physical dimensions of the nodes and links, it suffers from multiple link and node crossings (SI 2). The number of such conflicts increases linearly with r (Fig. 2A), a behavior analytically predicted by a geometric model (SI 2). To avoid these conflicts, we applied ELI and FUEL to several networks with different topologies (regular lattice, random, scale-free), spanning a range of sizes and link densities. We find that the obtained layouts undergo a geometrical layout transition as we increase the link thicknesses (Fig. 2E, F, G, H).

Figure 2

Crossover in Network Layouts:

(A) The number of link crossings in FDL grows linearly with r, saturating at very high r. A proper layout must resolve this increasing number of conflicts. (B) The average link length remains largely constant in the weakly interacting regime, but grows linearly (dashed gray line) in the strongly interacting regime. (C) The average link curvature rises slowly in the weakly interacting regime, then falling linearly in the strongly interacting regime. (D) The relaxation time grows significantly near the transition point. The transition occurs over a finite range of r/r (Shaded area A–D), regardless of the system size, offering evidence of a crossover. ELI and FUEL geometries: ELI (E, F, orange) and FUEL (G, H, blue) layout for a BA network [29] with N = 20, m = 2. When both ELI (E) and FUEL (G) layouts are similar to FDL. At larger r links bend to avoid each other. (F) In ELI at large r, links don’t fit inside the region containing the nodes and make outward arcs. (H) As nodes are free to move in FUEL, at large r the layout behaves more gently than ELI (link skeletons shown in white in (F) and (H)). (I) In the weakly interacting regime the links are thin and the layout radius is approximately the radius R of the bounding sphere surrounding N balls of radius r. (L) At larger r/r, thick links exclude each other and their volume dominates the volume of the layout. (J) The order parameter (scaling exponent of versus r/r for networks with different N and L and geometries (random, scale-free). (K) Rescaling the ratio of r/r using (4) collapses the transition point.

Weakly Interacting Regime: For small link thicknesses, r, the ELI and FUEL layouts are largely indistinguishable from the initial FDL layout. Indeed, at low r, the average link length is independent of r even as the link thickness increases by orders of magnitude (Fig. 2B). This is unexpected, given that there is a tenfold increase in the number of conflicts (potential link crossings) in this regime (Fig. 2A). The unchanged indicates that ELI and FUEL avoid the increasing number of conflicts by small local bending of the links, without the need to noticeably alter A similar behavior is seen for the average curvature of the links, finding that it shows modest changes throughout the weakly interacting regime (Fig. 2C), indicating that despite multiple local bending necessary to avoid conflicts the links remain largely straight. Note that the behavior of in this regime is model dependent, as the movement of nodes in FUEL offers avenues to avoid crossings, requiring less curving. Altogether, we find that in the weakly interacting regime local link rearrangements are sufficient to solve the multiple conflicts FDL suffers from.

Strongly Interacting Regime: Once r exceeds a critical value, we observe a dramatic change in the layout geometry (Fig. 2F, H). In ELI, with fixed node positions, the links must take long convoluted routes outside the network to reach their end-nodes, as they are unable to find sufficient empty space between the nodes. This change in the link structure is particularly visible on the skeleton of the layout (white links in Fig. 2F, H). In FUEL, with flexible node positions, the links reach their destination by pushing the nodes away from each other. These changes alter the behavior of which in the strongly interacting regime increases linearly with r, and induce rapid changes in the average curvature, at After the transition, the curvature decreases as 1/r. It is remarkable, however, that despite the different mechanisms the two models employ, the scaling of and in the strongly interacting regime in ELI and FUEL is universal, independent of the network topology. The observed linear increase in and the decrease in is consistent with isometric scaling, indicating that the layouts in the strongly interacting regime are structurally similar to each other if we rescale them by r (SI 5.1).

We can reveal the origin of the observed transition in the layout geometry by estimating the transition point When the links are much thinner than the node repulsion range r, the layout is dominated by the repulsive forces between the nodes, that together occupy the volume (SI 10). When the volume occupied by the links becomes comparable to V the layout must change to accommodate the links, inducing a transition from the weakly interacting regime to the strongly interacting regime. To calculate the transition point, the volume of all nodes and links need to be taken into account (Sec. SI 10), obtaining where and and the average is taken over the degree distribution. In scale-free and random networks, in the limit we obtain (SI 10), Using that in many real and model networks we have we obtain implying that in the limit i.e. the weakly interacting regime is absent in the thermodynamic limit. In other words, for large networks the crossings are so numerous that they cannot be ignored. Consequently, FDL and all currently used layout tools that ignore link crossings are expected to fail for large networks, as their layout is dominated by conflicts.

While networks with different N and L transition at different r/r ratios, if we scale r/r with the layout transition occurs near one for all networks. Using the scaling exponent of the average link length as the order parameter, the data collapse of Fig. 2K confirms the validity of (4). The fact that networks of rather different topologies (Scale-free, random, lattices and random geometric graphs (Sec. SI 11)) show similar dependence of the transition point on r, suggests that the transition of Fig. 2 independent of the network topology and the degree distribution. After rescaling of the transition point, networks of various sizes collapse on the same curve with the same width (Fig. 2K, see SI 10 and SI 11).

Finite-size scaling indicates that the layout transition occurs over a finite range of r/r, regardless of the network size (Sec. SI 11), indicating that we are observing a crossover [19, 20] from a mean-field behavior (exponent ) to a scaling behavior Indeed, the weakly interacting regime of both ELI and FUEL is well described by FDL, which ignores the interactions between the links and the possible link crossings are resolved as local perturbations. Yet, this weakly interacting regime disappears in the thermodynamic limit. For only the strongly interacting regime survives that, dominated by strong link-link interactions, displays universal scaling.

The observed crossover also alters the network’s physical properties. For example, a network’s response to external forces is captured by the Cauchy stress tensor [21] (SI 6), which depends on the physical and material properties of nodes and links. In the weakly interacting regime, the links are largely straight, hence the node terms V and V dominate the total stress. As the nodes are surrounded by a varying number of nodes, the stress does not spread uniformly in all directions, but has shear (off-diagonal) stress components, a common feature of solids. In the strongly interacting regime, the links fill up the space, hence the link contributions dominate resulting in a diagonal total stress tensor (SI 6). In other words, we predict that networks in the strongly interacting regime display a fluid or gel-like response to external stress. To test the validity of the predicted solid-gel transition, we compress the networks generated by FUEL in the y direction and measure the tensile forces (Fig. 3A, SI 6). We again, observe a crossover at predicted by (4) from a roughly constant stress in the weakly interacting regime to a monotonically increasing stress in the strongly interacting regime (Fig. 3B). Furthermore, as we rotate the network, we find that the total stress ratio displays large fluctuations in the weakly interacting regime, a common behavior in anisotropic solids. The fluctuations vanish at the transition point reaching the hydrostatic ratio (Fig. 3C), as expected for gels under pressure.

Figure 3

Stressed Networks and 3D Printing Networks:

(A) The tensile stress build-up in nodes and links as a result of compressing the network between two walls. Arrows indicate tensile stress (cyan indicates parallel to to the direction of compression, and green is perpendicular to it). The networks are colored based on the total amount of stress. In the weakly interacting regime (left) the stress is concentrated in the nodes and in the strongly interacting regime (right), almost all the stress is in the links. (B) The parallel stress of scale-free and random network topologies as a function of link thickness. Since the definition of x,y,z is frame-dependent, we measure the forces for 50 random network orientations. (C) The ratio of parallel and transverse tensile stress error bars correspond to one standard deviation around the mean, calculated over the 50 random orientations. In the weakly interacting regime the ratio depends on layout orientation, a solid-like feature. In the strongly interacting regime, the fluctuations of decay to zero. 3D Printing Networks: The Flavor Network [23] with N = 184, L = 716, represents a network of ingredients sharing flavor compounds. (D) An FDL layout applied to this network results in multiple conflicts (red). The inset highlights a densely connected region containing dairy products, so highly overlapping that it is impossible to discern the underlying network. (E) Laying out the flavor network using FUEL, the conflicts disappear, unveiling the inner structure of network communities. We 3D printed the resulting flavor network using a commercial 3D printer, allowing us to inspect the full internal structure of the network.

In summary, the layout geometry of physical networks is characterized by two distinct regimes: a weakly interacting regime, where the overlap between the nodes and links is avoided via local link rearrangements, and a strongly interacting regime, whose layout is shaped by excluded link volume interactions. Networks in the weakly interacting regime are solid-like, while those in the strongly interacting regime behave like a gel. The observed transition between the two regimes is unique to three dimensions. Indeed, as links are one dimensional objects, the non-crossing condition results in knot-like constraints in 3D, preventing the links from passing through each other. In D ≥ 4, knots of one dimensional objects can be untied [22], so the non-crossing conditions will not constraint the embedding. Therefore, D = 3 is the lowest dimension where links can avoid each other by bending, and also the highest dimension where they cannot pass by each other without breaking or tunneling.

Both phases have potentially important applications. The weakly interacting regime is the one traditionally used in network visualizations, and can help us 3D print networks, offering novel ways to explore the inner structure of a complex system. We illustrate this in Fig. 3D, where we applied FDL to the Flavor Network [23]. Given the high link density, 2D visualizations suffer from visual cluttering, making only a fraction of the links visible [23]. A 3D layout offers more clarity, but given the dense communities, FDL leads to excessive node and link overlaps (inset, Fig. 3D), obstructing the true wiring of the network. Applying FUEL, and choosing r so that the layout stays in the weakly interacting regime, the obtained layout geometry reveals the underlying communities present in the network. This layout is amenable to 3D printing, as shown in Fig. 3E, offering unprecedented opportunities to interact with the network and to directly inspect its inner structure. Note that the finding that for large N the weakly interacting regime disappears, indicates that the toolset introduced here, capable of resolving crossings, is inevitable as we aim to visualize large networks, when link crossings dominate the layout.

The strongly interacting regime has direct relevance to the brain, a three dimensional physical network where the close-packing of the axons is critical to their ability to form synapses [24, 25]. Isometric scaling laws observed in rodent brains [26] indicate that volume V and the area A of the white matter relate to each other as [26]. This implies that in these networks the average neuron length scales as with the axon thickness, a behavior accurately predicted by the strongly interacting regime (Fig. 2B). Indeed, if we adapt the conventional approach and describe anatomical regions as nodes, and axon bundles connecting these regions as links, the thickness of the axon bundles (r) is comparable to the size of the anatomical regions, supporting the prediction of the empirical scaling that these brain networks are in the strongly interacting regime. Thus, (1)–(3) offers an appropriate modeling framework to capture the layout geometry of dense neuronal networks, generating an economic spatial layout [27, 28] while respecting the non-crossing conditions each axon must obey [1].