
==== Front
Nat Commun
Nat Commun
Nature Communications
2041-1723
Nature Publishing Group UK London

38849343
49227
10.1038/s41467-024-49227-8
Article
Physical networks as network-of-networks
http://orcid.org/0009-0001-6878-4038
Pete Gábor gabor.pete@renyi.hu

12
Timár Ádám 13
Stefánsson Sigurdur Örn 3
http://orcid.org/0000-0002-6112-8273
Bonamassa Ivan 4
http://orcid.org/0009-0001-4680-3349
Pósfai Márton posfaim@ceu.edu

4
1 https://ror.org/03vw74f64 grid.423969.3 0000 0001 0669 0135 HUN-REN Alfréd Rényi Institute of Mathematics, Budapest, Hungary
2 https://ror.org/02w42ss30 grid.6759.d 0000 0001 2180 0451 Department of Stochastics, Institute of Mathematics, Budapest University of Technology and Economics, Műegyetem rkp. 3., H-1111 Budapest, Hungary
3 https://ror.org/01db6h964 grid.14013.37 0000 0004 0640 0021 University of Iceland, Reykjavík, Iceland
4 https://ror.org/02zx40v98 grid.5146.6 0000 0001 2149 6445 Department of Network and Data Science, Central European University, Vienna, Austria
7 6 2024
7 6 2024
2024
15 488214 6 2023
27 5 2024
© The Author(s) 2024, corrected publication 2024
2024
https://creativecommons.org/licenses/by/4.0/ 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 Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted 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 licence, visit http://creativecommons.org/licenses/by/4.0/.
Physical networks are made of nodes and links that are physical objects embedded in a geometric space. Understanding how the mutual volume exclusion between these elements affects the structure and function of physical networks calls for a suitable generalization of network theory. Here, we introduce a network-of-networks framework where we describe the shape of each extended physical node as a network embedded in space and these networks are bound together by physical links. Relying on this representation, we introduce a minimal model of network growth and we show for a general class of physical networks that volume exclusion induces heterogeneity in both node volume and degree, with the two becoming correlated. These emergent properties strongly affect the dynamics on physical networks: by calculating their Laplacian spectrum as a function of the coupling strength between the nodes we show that degree-volume correlations suppress the role of hubs as early spreaders in diffusive dynamics. We apply the network-of-networks framework to describe several real systems and find properties analog to the minimal model networks. The prevalence of these properties points towards general growth mechanisms that do not depend on the specifics of the systems.

In analysis of physically embedded complex networks, their nodes are usually considered as localized spheres connected by links, neglecting possible differences of nodes spatial shapes. The authors develop a representation of physical networks that captures arbitrary node shapes to characterize structural and dynamical network properties.

Subject terms

Complex networks
Applied mathematics
100010663 EC | EU Framework Programme for Research and Innovation H2020 | H2020 Priority Excellent Science | H2020 European Research Council (H2020 Excellent Science - European Research Council) 772466-NOISE 100010663 EC | EU Framework Programme for Research and Innovation H2020 | H2020 Priority Excellent Science | H2020 European Research Council (H2020 Excellent Science - European Research Council) 810115-DYNASNET 100010663 EC | EU Framework Programme for Research and Innovation H2020 | H2020 Priority Excellent Science | H2020 European Research Council (H2020 Excellent Science - European Research Council) 810115-DYNASNET Icelandic Research Fund grant No. 239736-051Icelandic Research Fund grant No. 239736-051100010663 EC | EU Framework Programme for Research and Innovation H2020 | H2020 Priority Excellent Science | H2020 European Research Council (H2020 Excellent Science - European Research Council) 810115-DYNASNET issue-copyright-statement© Springer Nature Limited 2024
==== Body
pmcIntroduction

The building blocks of physical networks are extended objects that do not intersect each other, resulting in non-trivial geometric layouts1, link entanglement2 and emergent correlations between physical and network structure3. Yet, these works model nodes as localized spheres connected by extended tube-like links, an assumption that does not necessarily reflect the structure of most real-world physical networks. In the connectome, for example, nodes represent neurons with non-trivial dendritic shapes, and links are point-like synapses4. A similar picture emerges for molecular networks such as the cytoskeleton, mitochondrial networks, or fiber materials, where nodes are extended molecular strands and bonds between them are localized5–7, as well as in the wood-wide-web, where extended tree roots and mycelia connect to form a complex underground network8,9. Therefore, the sphere-tube paradigm often falls short of describing physical networks, calling for a more general framework to cope with the complex shape of nodes and links.

In this work, we develop a network-of-networks representation of physical networks that is able to capture arbitrary node shapes10,11 and allows us to characterize both structural and dynamical properties of networks. Relying on the network-of-networks framework, we introduce a model that grows physical networks from fractal segments. Analytically solving the model, we show that physicality induces heterogeneity in both the physical and the network properties of the nodes and that the two become strongly correlated. These correlations also affect the dynamics on the networks: generalizing the combinatorial Laplacian to physical networks12–14, we show that fast dynamical modes associated to hubs (and corresponding to the tail of the Laplacian spectra) are suppressed by the emergent correlations between node volume and degree. The usefulness of the mathematical tools we develop in this paper goes beyond the specifics of the model, and we demonstrate this by applying our framework to several real physical networks, including a recently collected data set describing more than ~ 20, 000 neurons of the adult fruit fly’s brain15. In doing so, we identify positive node degree-volume correlations similar to our minimal growth model, and we show that these have an analog effect on the Laplacian spectrum of the connectome. The fact that degree-volume correlations emerge in a minimal model while also prevalent in real systems suggests a general mechanism behind such correlations that does not depend on the complex details of the growth of real networks.

Results

Network-of-networks representation

We aim to describe physical networks embedded in some substrate or medium. In its most general form, the substrate is represented by a graph S, and each physical node i is an extended object occupying a subgraph Vi⊂S. To capture volume exclusion, we do not allow nodes to overlap, i.e., Vi∩Vj=∅ for i ≠ j. Two nodes i and j may form a link (i, j) if they occupy adjacent sites in S. The physical layout P of the network is a network-of-networks, i.e., it is the union of physical nodes, where each node is a network itself, together with the bonds forming the connections between the nodes (Fig. 1a). The layout P is a physical realization of the combinatorial network G, where node i of G corresponds to the physical node Vi, and nodes i and j are connected if there is a bond between Vi and Vj in P (Fig. 1b). Though the substrate S can represent any available space, here we focus on substrates that are d dimensional cubic lattices with linear size L and periodic boundary conditions. Note that network representations of this kind are employed in the graph drawing literature with the focus on algorithms that embed a given combinatorial network into S16. Here, we are interested in physical networks P growing in S, the emergent relation between P and G, and its consequences on the dynamics on the network.Fig. 1 Network-of-networks representation.

a Each physical node Vi is a subgraph of the substrate S. Physical nodes cannot overlap, i.e., Vi∩Vj=∅ for i ≠ j. The physical layout P (dashed area) is a network-of-networks: it is the union of physical nodes Vi together with the bonds connecting them. b The combinatorial network G is a coarse-grained representation of P capturing the connections between the nodes without the physical structure.

Network growth

To study the effect of physicality on network evolution, we define a model of network growth relying on the network-of-networks representation. We start with an empty S and we place a single physical node V0 occupying a subset of the sites. We add the rest of the nodes iteratively: At time step t we add a new node Vt that is seeded at a random unoccupied site and grows until it hits an already existing node Vs and a link (t − s) is formed. The growth of node Vt is driven by some random or deterministic process; and we assume that the physical nodes are characterized by a fractal dimension df ∈ [1, d]17,18. We add N physical nodes or until all of S is occupied; in the latter case we call the physical network saturated.

Since the total volume of the network increases over time, later nodes hit the network at higher rates, and the typical size of nodes decreases. Hence, we expect that nodes added early have a higher degree than nodes added in the final stages of the network evolution, both because they are larger and they have more time to collect connections. This suggests that to analytically characterize the evolution of the physical network two ingredients have to be considered: (i) network growth, i.e., nodes are added sequentially to the system and (ii) externally limited node growth, i.e., the nodes grow until they hit the already existing network. We show that these two ingredients lead to the emergence of power law combinatorial networks with degree exponents γ≤3.

We start the analytical treatment of the model by estimating the probability pij that two randomly placed physical nodes Vi and Vj intersect. If the two boxes containing the physical nodes have side length li ≫ lj, respectively, and the larger node Vi intersects the box containing the smaller node Vj, then, by dimension count, the two nodes overlap with positive probability if df≥d/2. We can tile the lattice with (L/lj)d boxes with side length lj, and the number of such boxes intersected by Vi is ~(li/lj)df. Therefore the intersection probability is1 pij~lidfljd−dfLd~vivjd/df−1Ld,

where vi=∣Vi∣~lidf is the volume of node i. If, however, df < d/2 and lj≤li ≪ L, then the nodes avoid each other with high probability. In this case, the intersection probability will have the meanfield behavior, well-approximated by the probability of selecting the sites of Vi and Vj uniformly from S, i.e., pijMF~vivj/Ld, which is independent of df and agrees with Eq. (1) for df = d/2.

Using the same box-counting argument that led to Eq. (1), the probability that a node added at time t intersects any existing node s < t is approximately given by ∑s<tpst=vtd/df−1Vt−1/Ld, where Vt−1 is the total volume of nodes s < t. A key observation is that the size of node t increases until it hits the existing network, meaning that vt increases until ∑s<tpst ≈ 1, allowing us to estimate the volume of node t as2 vt≈Vt−1/Ld−dfd−df.

Equation (2) allows us to express the evolution of the expected total volume via the recursion Vt+1 = vt+1 + Vt with initial condition V0 = v0. Using a continuous time approximation, we obtain3 Vt≈Lddd−dft+cLdd−dfd~LdtLd1−df/d,

where c is a constant depending on v0. A natural choice for the latter is that the first node spans the entire available space, i.e., v0~Ldf, in which case c is independent of L. Equation (3) predicts that Nsat, the number of nodes when the network saturates, scales as Nsat ~ Ld, meaning that the average node volume 〈v〉 remains constant in the L → ∞ large system limit. Therefore, the physical layout P is optimal in the sense that no physical representation of a combinatorial network of Nsat nodes can fit into a smaller volume than ~ Nsat ~ Ld. It is noteworthy that the model achieves this bound despite the fact that the nodes grow randomly.

In light of Eq. (3), we can now calculate the expected degree of the physical nodes in the combinatorial network G. In the continuous time approximation, the volume of the newly added node vt is provided by the time derivative of Vt, i.e., vt~t/Ld−df/d. Hence, following Eq. (1), the expected degree of node t after the addition of N nodes is4 kt(N)=1+vtLd∑s=t+1Nvsd/df−1~vt⋅NLdddf,

where the proportionality is valid for t ≪ N. This means that the volume occupied by large nodes (i.e., nodes that were added early) in the physical layout is proportional to their degree in the combinatorial network.

We finally calculate the complementary cumulative degree distribution P(k)=1−1N∑t:kt≥kN1, finding that P(k) ~ k−(γ−1) with exponent5 γ=1+ddf.

For df≤d < 2df the degree exponent falls in the range 2≤γ < 3. In the mean-field regime, the degree exponent can be obtained by substituting d/df with 2, yielding γMF = 3. Note that the upper critical dimension of the physical network depends on the kinetic growth of the nodes. For example, growing nodes along a straight trajectory in a random direction generates nodes with df = 1; therefore, the networks fall in the mean-field regime df≤d/2 even for embedding dimension d = 2.

In the above model, each physical node grows starting from a random location following some growth process. We stress that our calculations hold for a general class of node growth algorithms, the crucial assumption being that if the boxes around two random walk pieces intersect, then with uniformly positive probability the trajectories also intersect, which implies a level of isotropy of node growth. As a counterexample, consider nodes that always grow in one direction along one of the axes. Such nodes will run parallel to each other, avoiding intersection, hence the resulting network will be a collection of disconnected chains. If, however, the nodes grow along straight lines but in random directions, thus restoring isotropy on average, then the resulting network has a power law degree distribution (SI Sec. S1.4).

Numerical simulations

To test our analytical predictions, we numerically generate physical networks where nodes grow according to random walk trajectories. Specifically, we generate nodes using loop-erased random walks (LERWs), i.e., a trajectory that evolves as a simple random walk in which any loop is erased as soon as it is formed19–22. Here, we focus on the LERW, as it represents a tractable model ofnon-self-intersecting random trajectories with well-understood non-trivial critical properties. Its critical properties are studied both in the mathematics and physics literature23–26; for example, their fractal dimension in d = 2 is df = 5/422, in d = 3 it is df ≃ 1.6236(4)27,28, while its upper critical dimension is du = 4 where df = 2 with a logarithmic correction29. (See “METHODS” for further details.) We remark that our predictions are not specific to LERWs, in Sec. S1 of the Supplementary Information, we study various alternative growth processes.

Knowing the fractal dimensions of LERWs allows us to directly verify the predictions of Eqs. (3)–(5):

Volume evolution

Equation (3) predicts that the total volume of the physical network evolves as Vt~t1−df/d. Figure 2a shows the excellent agreement between the theoretical predictions and numerical simulations. Note that, as expected, in the mean-field regime d≥4 the network volume follows the classic diffusion growth Vt ~ t1/2. Figure 2b further corroborates the predicted scaling of the number of nodes at saturation, i.e., Nsat ~ Ld.Fig. 2 Evolution of LERW physical networks.

a The temporal evolution of the total volume Vt of physical networks, dashed lines represent the theoretical prediction Eq. (3). Networks built from LERWs in dimensions d≥4 fall into the mean-field regime. b The number of physical nodes in the saturated networks is proportional to the volume of the substrate ∣S∣ irrespective of df and d. c Node degree is proportional to node volume independently from df and d. d The complementary cumulative degree distribution function (CCDF) of the physical networks. Dashed lines indicate the predicted degree exponent γ = 1 + d/df≤3, while the dotted line corresponds to the mean-field behaviour, ﻿γMF = 3. Plots (a), (c) and (d) represent measurements of single networks with initial condition v0=Ldf and ∣S∣=106. In (b) markers are an average of 10 independent networks, error bars represent the standard error of the mean. Lines corresponding to different slopes are shifted to increase readability.

Degree-volume correlations

A second prediction is the emergence of degree-volume correlations, capturing the interplay between the physical layout P and the combinatorial network G. In particular, Eq. (4) predicts a linear proportionality between the node volume vi and degree ki, and we again find excellent agreement with simulations for all the tested dimensions (Fig. 2c).

Power law emergence

As a final test, we verify the emergence of power law scaling in the degree distribution of the combinatorial networks G. As anticipated in Eq. (5), the degree exponent depends on both the dimensionality of the embedding substrate, d, and the fractal dimension of the nodes, df. Figure 2d shows that numerical simulations confirm the predicted degree exponent γ = 1 + d/df for dimensions d < 4, while the mean-field exponent γMF = 3 is found for d≥4. In traditional models of combinatorial networks, heterogeneity typically arises from preferential attachment or some other optimization process. Our model is based on random growth without any explicit preference to create highly connected nodes; therefore, the uniform attachment tree may be considered as the non-physical counterpart of our model. Uniform attachment yields exponential degree distribution, hence the power law distribution observed here is a direct consequence of volume exclusion, which, together with the dynamic network growth, induces effective preferential attachment.

Physical network Laplacian

The layout P is a physical realization of the combinatorial network G. Traditional studies of dynamics on physical networks ignore the layout P and focus only on the role of G, thus prompting the question: does modeling dynamics on G accurately capture dynamics on physical networks? To explore this, we study the spectral properties of P and show that physical nodes emerge as functional units through timescale separation, yet even in this limit the structure of P continues to affect the dynamics. We focus on the Laplacian spectrum12 since it influences the behavior of several dynamical processes on networks30 including diffusion31,32, synchronization33 and it underlies the definition of several information-theoretic tools to analyze the multiscale functioning of networks10,14,34–37.

We study the problem by invoking, once again, the network-of-networks representation. In our setup, we assign a weight to each connection in P such that links within physical nodes have weight 1 and links connecting two physical nodes have weight w, capturing that in real physical networks bonds between nodes are often qualitatively different than those within nodes. The weighted Laplacian matrix of P occupying V sites is then QP=DP−AP, where AP is the V × V weighted adjacency matrix and DP is a diagonal matrix such that [DP]ss=∑u[AP]su is the sum of the weights of the links adjacent to site s in P. If we now set w = 0, the network-of-networks falls apart and each physical node becomes a separate connected component, resulting in a block-diagonal Laplacian QP(0)=diag(QV1,QV2,…,QVN), where QVi is the Laplacian of the physical node Vi. The Laplacian QP(0) has N zero eigenvalues corresponding to the N blocks (i.e., the physical nodes), hence we can assign an eigenvector ui(w = 0) to the i-th node such that [ui(0)]s=1/vi if site s is within node i, otherwise [ui(0)]s=0, where vi=∣Vi∣ is the volume of node i. Since linear combinations of these vectors are also eigenvectors, we can write the zero eigenvectors of QP as u(0)=Mu~, where M is the N × V membership matrix such that [M]si=1/vi if site s is part of node i, otherwise [M]si = 0, and u~ is any normalized N dimensional vector.

We can gain insights about the spectral properties of QP by working in the weak coupling regime w ≪ 1 and relying on perturbation theory. Following a treatment similar to the one adopted to study diffusion in multilayer networks38–40, we consider w a small perturbation and write QP(w)=QP(0)+wQP′, where QP′ is the Laplacian matrix of the subnetwork of P formed by the bonds between physical nodes. The characteristic equation, up to first order in w, becomes then6 QP(0)+wQP′u(0)+wu′≈λ(0)+wλ′u(0)+wu′.

Perturbations around λ(0) = 0 lead to N eigenvalues that are O(w), while the rest of the eigenvalues are constant in leading order (Fig. 3a). This means that on the 1/w timescale, diffusion-like dynamics on the physical network are captured by the N slow eigenmodes. We obtain these from Eq. (6) (see “METHODS”), yielding7 V−1/2QGV−1/2u~=λ′u~,

where QG is the N × N Laplacian matrix of the combinatorial network G and V is an N × N diagonal matrix such that its diagonal elements are [V]ii = vi. We call the volume-normalized Laplacian the physical network Laplacian Qphys=V−1/2QGV−1/2.Fig. 3 Laplacian of LERW physical networks.

a For decreasing weight w the eigenvalues of QP separate into two groups: eigenvalues corresponding to the zero eigenmodes of QP(w=0) decay as ~ w (blue), while the rest converge to a constant value (teal). (b–d) Comparing the spectrum of the physical Laplacian Qphys=V−1/2QGV−1/2, the randomized physical Laplacian Qphysrand=Vrand−1/2QGVrand−1/2, and the Laplacian of the combinatorial network QG. b, c The heterogeneous node volume distribution and the correlation between node degree and volume significantly reduce the largest eigenvalues of the spectra. d Heterogeneous node volumes alone explain the reduction of the algebraic connectivity λ2. Eigenvalues are calculated for d = 2 and L = 10 in (a) and L = 100 in (b). In (c) and (d) markers are an average of 10 independent networks, error bars represent the standard error of the mean.

Equation (7) is a key relation for understanding the dynamics on physical networks since it allows to characterize the dynamics on P on the timescale 1/w in a coarse-grained way: after integrating out the fast modes corresponding to eigenvalues λ(w) ≫ w, the state of each physical node Vi is given by a single variable, while the coupling between the nodes is provided by the combinatorial network G. However, the combinatorial Laplacian QG is not sufficient to capture the dynamics, and we must normalize QG by the volume of the nodes, as shown in Eq. (7). This means that physical networks with the same combinatorial network but different layout can have drastically different dynamical properties. For example, if nodes have approximately the same size, i.e., vi ≈ V/N, then the physical layout only affects the overall timescale, otherwise the Laplacian spectrum is determined by QG. If, however, node sizes are heterogeneously distributed, normalizing by volume will also have a heterogeneous effect on the eigenvalues.

Application to the physical network growth model

We showed above that physical networks generated by our network growth model are characterized by heterogeneous node-volume distribution and proportionality between the degree and the volume of nodes (Fig. 2). To probe the effect of this emergent correlation, we shuffle the volume of the nodes of a LERW physical network to remove the correlation between network and physical structure. We then compare the spectrum of the volume-normalized Laplacian Qphys=V−1/2QGV−1/2 to its randomized version Qphysrand=Vrand−1/2QGVrand−1/2 and to the Laplacian spectrum of the combinatorial network QG. Figure 3b shows that the spectrum of QG has a heavy tail characterized by the same γ exponent of Eq. (5), as expected for combinatorial networks with power law degree distributions12. Adding heterogeneous but uncorrelated node sizes does not influence the tail while taking into account the degree-volume correlation of nodes removes the heavy tail and leads to a rapidly decaying spectrum. In power law networks, the eigenvector corresponding to the largest eigenvalue λN of QG is typically concentrated on the node with the largest degree41,42. In our model, the largest degree node also has the largest volume; therefore normalizing by node volume V−1/2QGV−1/2 significantly lowers λN. Since node sizes are heterogeneously distributed, with high probability, we associate volume ~ 1 to the highest degree node after randomization. Hence, the eigenvalue λN of QG is largely unaffected by the randomized normalization (Fig. 3c). At the other end of the spectrum, controlling the long-time mixing of the dynamics, the eigenvector associated with the algebraic connectivity λ2 typically spans the entire network. Figure 3d shows that taking node volumes into account slows the dynamics down; however, degree-volume correlations do not significantly affect λ2.

Note that positive degree-volume correlations, responsible for the suppression of the tail of the Laplacian spectrum, naturally arise in minimum-volume physical realizations of combinatorial networks. Any combinatorial network G has many possible physical realizations P, a minimum volume realization is a P that minimizes the total volume of the network. Consider node i∈G with degree ki; in any possible P, the physical realization of node i must have volume at least proportional to ki, otherwise it is unable to support ki connections. Therefore, we expect positive degree-volume correlations in minimum-volume physical layouts. This means that any physical network generation process that minimizes total volume – either explicitly or as an emergent property, like in our model – is characterized by positive degree-volume correlations and hence that the spectrum of Qphys is similarly affected by physicality as in our model.

Real physical networks

We identified the degree-volume correlations and the profile of the Laplacian spectrum as important features of physical networks that can emerge even in the simplest models. To measure these properties, we do not need a detailed description of the layout of a physical system – we only need the combinatorial network and a list of the node volumes, allowing us to describe very large and complex physical networks. As a case study, we investigate a recently published data set providing the three-dimensional layout of more than 20,000 neurons of the brain of an adult fruit fly and the location of more than 13 million synapses connecting them (Fig. 4a)15. Although our simple growth model does not attempt to capture the myriad of complex mechanisms shaping brain development, we find that the fruit fly brain is characterized by similar emergent properties as the model networks. Figure 4b shows, for example, that the multiplicity-weighted node degree, i.e., the number of synapses a neuron has, can be approximated by a power law γff ≈ 2.3, albeit with an exponential cutoff43,44. We also find a strong positive correlation between the weighted degree and the volume of the nodes (Fig. 4c).Fig. 4 Fruit fly brain network.

a The neurons have complex three-dimensional shapes. Two intertwined neurons (teal and yellow) are connected by synapses (red circles). b The tail of the weighted degree distribution can be approximated with the power law with γff ≈ 2.3 with an exponential cutoff. c Similar to the network growth model, the volume of the nodes v is strongly correlated with their degree k. The markers indicate binned degree-volume averages, where the data points are binned based on node volume. The shading represents a kernel density estimate of the joint v-k distribution, the dashed line indicates the least squares fit of the power law scaling. d The effect of the node degree-volume correlation on the Laplacian spectrum in the brain network is analog to the effect of correlations in the network growth model. The spectra of Qphys and Qphysrand are shifted to the right to allow direct comparison of the tails of the distributions. e, f A visualization of the leading eigenvectors u~N of the combinatorial and physical Laplacians. The color of each node i corresponds to ∣u~N(i)∣, the weight of the leading eigenvector at node i, and the size of the nodes is a linear function of their degree. (e) The eigenvector u~NG is concentrated on the node with the highest degree iG. f Due to degree-volume correlations, the u~Nphys is concentrated on node iphys, which has the maximum degree-volume ratio in the network. g The weight of the leading eigenvectors ∣u~NG(i)∣ and ∣u~Nphys(i)∣ decays exponentially as a function of the geodesic distance from iG and iphys, respectively. Error bars indicate the standard error of the mean and are typically smaller than the marker size. Node locations in (e) and (f) are generated based on G, and do not correspond to the actual physical locations.

To compare the spectrum of the combinatorial Laplacian QG and the physical network Laplacian Qphys=V−1/2QGV−1/2, we measure volume in units such that the mean node volume is unity, i.e., 〈v〉 = 1 (see “METHODS” for further details). Calculating the leading eigenvalues of QG and Qphys, we find that λNG/λNphys≈32.7, indicating that degree-volume correlations greatly suppress the modes of the dynamics that spread the fastest, similarly to model networks. This is further supported by Fig. 4d, showing again that physicality suppresses the tail of the spectrum.

To further probe the role of degree-volume correlations, we calculate the leading eigenvectors u~N of QG and Qphys. Figure 4e, g show that, as expected for heterogeneous combinatorial networks, u~NG is concentrated on the largest hub iG in the network, and the weight of the eigenvector decays exponentially as the geodesic distance from iG in G. This means that, without taking physicality into account, the largest degree node is also the earliest spreader of diffusive dynamics. For the physical Laplacian Qphys we find a different picture: u~Nphys is again concentrated on a single node iphys; this node, however, is not the largest hub. The leading eigenvector instead is centered on a node that balances high degree and low volume: node iphys is the 159th largest degree node and is at the top 15 percentile of the volume distribution. In fact, node iphys is the node that maximizes the degree-volume ratio, i.e., iphys = argmaxiki/vi. This means that degree-volume correlations not only slow down spreading dynamics in physical networks, but also change the identity of the early spreaders.

Here we chose to focus on the fruit fly brain network as it represents one of the largest and most detailed maps of physical networks available; however, our framework is not specific to neural networks. In the Sec. S2 of the Supplementary Information, we analyze four additional real systems: a network describing the cavities of a porous material, a neural network of a nematode, a river network, and a vascular network. In each case, we find positive degree-volume correlations and that these correlations suppress the tail of the Laplacian spectra. The fact that the physical and network properties of nodes become intertwined in such a diverse set of real networks, and also in the simplest models, indicates a general mechanism behind the emergence of degree-volume correlations that do not depend on the details of the individual networks.

Discussion

Physical networks are complex networks that have a complex three-dimensional layout. The network-of-networks framework naturally lends itself to representing these systems: representing nodes as physically embedded networks allows us to capture arbitrary node shapes and complex wiring. Here, we relied on the network-of-networks framework to characterize both model and real physical networks. We identified correlations between node degree and volume as a prevalent feature of physical networks: We analytically showed that degree-volume correlations emerge in a minimal network growth model, in fact, we provided arguments that such correlations naturally arise through any growth process that minimizes network volume. We also showed that positive degree-volume correlations are generally present in real systems. These correlations have important consequences on dynamics unfolding on physical networks: the tail of the physical Laplacian spectrum is suppressed by the large volume of hubs. More broadly, these results vividly demonstrate that traditional methods of network science focusing on combinatorial networks cannot fully describe physical networks and that their three-dimensional layout must be accounted for.

Our work opens new avenues for physical network research in several ways. First, by establishing the connection between physical networks and network-of-networks, we allow future work to leverage the rich literature of multi-layer networks to characterize physical systems10,11. For example, multi-layer centrality measures can be used to quantify the importance of physical nodes45–47. Second, previous work on physical networks relies on methods that require a full description of their spatial layout and, therefore, are often limited to systems of a few hundred nodes1–3. In contrast, the quantities we studied can be measured relying on the combinatorial network and a list of node volumes, allowing the characterization of large-scale physical networks without the need of the full three-dimensional layout. For example, we can tune the volume of the nodes to systematically study how physical layout affects the Laplacian spectrum. Finally, the simple growth model and its analytical description can serve as the starting point for the exploration of additional growth mechanisms that characterize neural networks and other physical networks. For example, future work may study branched nodes, long-range interactions that guide the growth of physical nodes, or the expansion of available space by modeling the evolution of the underlying substrate.

Methods

Loop-erased random walks

In our network growth model, we can generate physical nodes with any stochastic or deterministic process that produces a growing fractal embedded in Zd. Standard self-avoiding walks are traditionally used to model polymers obeying volume exclusion and, therefore, represent a natural choice to model node growth17. However, the naïve kinetic version of the self-avoiding walk traps itself in two and three dimensions at finite length21, making it a poor candidate for constructing large physical networks. Instead, we focus on loop-erased random walks (LERW): a LERW evolves as a simple random walk, except when it intersects itself, we delete the loop that it created and continue the walk20. This guarantees that the final physical node does not intersect itself and that the walk never gets trapped. Alternatively, the LERW can be defined as a special case of Laplacian-random walks, where transition probabilities are defined by a harmonic function48,49. This alternative construction does not require deleting loops, hence is more realistic as a growth model. The LERW has attractive mathematical properties making it amenable to analytical treatment. For example, Wilson’s algorithm uses iterative LERWs to construct a uniform spanning tree (UST) of any graph24. In fact, the physical network our algorithm constructs is a UST of the S substrate together with a partition identifying the nodes. Future work may exploit this connection between USTs and LERW physical networks, together with known results in dimensions d = 2 and d > 422,50, to rigorously prove some of the results presented here.

Perturbation of the physical Laplacian

To obtain the slow eigenmodes, we match the first-order terms of Eq. (6) and substitute u(0)=Mu~, so that8 QP(0)u′+QP′Mu~=λ′Mu~.

Multiplying from the left by the transpose of the membership matrix M we get9 MTQP(0)u′+MTQP′Mu~=λ′MTMu~.

The ith row of MT is the trivial eigenvector ui(w = 0) corresponding to physical node i; therefore MTQP(0) is all zeros and MTM is the N × N identity matrix, leading to Eq. (7) in the text.

The fruit fly connectome

We study the Hemibrain data set, which describes a portion of the central brain of the fruit fly, Drosophila melanogaster15. The physical layout of the connectome is provided by the detailed three-dimensional shape of each neuron and the location of the synapses between them. The corresponding combinatorial network contains 21,662 nodes representing neurons and 13,603,750 links representing synapses. Synaptic partners are connected through approximately 5 synapses on average, and the maximum number of synapses between two neurons is 6039. In our calculations, we treat the combinatorial network as a weighted and undirected network, where the weight of the link (i, j) is equal to the number of synapses between neurons i and j. Note that we only require the combinatorial network and the volume of each node for our calculations; therefore, the detailed physical layout of the connectome is, in fact, not needed.

Note that the Hemibrain data set covers a large portion of, but not the entire, fruit fly brain. Since degree and volume are local properties of the nodes, we expect that the results presented here would not change significantly if the entire connectome were to be considered.

Degree distribution

We find that the weighted degree distribution has a heavy tail, which can be approximated by a power law with γff ≈ 2.3 for degrees ≥1058 with an exponential cutoff; the power law fit, however, cannot be distinguished from a lognormal fit on the same range43,44.

Laplacian spectrum

Comparing the spectrum of the combinatorial Laplacian QG and the volume-normalized Laplacian Qphys=V−1/2QGV−1/2 carries a level of ambiguity: QG does not depend on the node volumes, while changing the unit of volume multiplies the spectrum of Qphys by a constant. To meaningfully compare the two spectra, (i) we think of QG as a physical Laplacian where all nodes have unit volume, and (ii) we set the mean node volume in Qphys to unity, i.e., 〈v〉 = 1. With this choice of units, any difference in the eigenvalues is due to the heterogeneous distribution of node volumes in the physical network and not to a global shift caused by the choice of units.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

Supplementary information

Supplementary Information

Reporting summary

Peer Review File

Supplementary information

The online version contains supplementary material available at 10.1038/s41467-024-49227-8.

Acknowledgements

I.B., M.P., and Á.T. were funded by ERC grant No. 810115-DYNASNET. ÁT and SÖS acknowledge partial support from the Icelandic Research Fund, grant No. 239736-051. GP was funded by the ERC Consolidator Grant No. 772466-NOISE.

Author contributions

M.P. developed and performed the numerical simulations. M.P. and I.B. performed the data analysis. G.P., Á.T., S.Ö.S., I.B., and M.P. contributed to the analytical results and the conceptual design of the study. M.P. was the lead writer of the manuscript.

Peer review

Peer review information

Nature Communications thanks Jianxi Gao, and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. A peer review file is available.

Funding

Open access funding provided by HUN-REN Alfréd Rényi Institute of Mathematics.

Data availability

Data to reproduce the figures is available at https://github.com/posfaim/physnets_as_net-o-nets.

Code availability

Code to generate random networks and reproduce the figures is available at https://github.com/posfaim/physnets_as_net-o-nets51.

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.
==== Refs
References

1. Dehmamy N Milanlouei S Barabási A-L A structural transition in physical networks Nature 2018 563 676 680 10.1038/s41586-018-0726-6 30487615
Dehmamy, N., Milanlouei, S. & Barabási, A.-L. A structural transition in physical networks. Nature 563, 676–680 (2018).30487615
2. Liu Y Dehmamy N Barabási A-L Isotopy and energy of physical networks Nat. Phys. 2021 17 216 222 10.1038/s41567-020-1029-z
Liu, Y., Dehmamy, N. & Barabási, A.-L. Isotopy and energy of physical networks. Nat. Phys. 17, 216–222 (2021).
3. Pósfai M Impact of physicality on network structure Nat. Phys. 2024 20 142 149 10.1038/s41567-023-02267-1
Pósfai, M. et al. Impact of physicality on network structure. Nat. Phys. 20, 142–149 (2024).
4. Bullmore E Sporns O Complex brain networks: graph theoretical analysis of structural and functional systems Nat. Rev. Neurosci. 2009 10 186 198 10.1038/nrn2575 19190637
Bullmore, E. & Sporns, O. Complex brain networks: graph theoretical analysis of structural and functional systems. Nat. Rev. Neurosci. 10, 186–198 (2009).19190637
5. Viana MP Mitochondrial fission and fusion dynamics generate efficient, robust, and evenly distributed network topologies in budding yeast cells Cell Syst. 2020 10 287 297 10.1016/j.cels.2020.02.002 32105618
Viana, M. P. et al. Mitochondrial fission and fusion dynamics generate efficient, robust, and evenly distributed network topologies in budding yeast cells. Cell Syst. 10, 287–297 (2020).32105618
6. Fletcher DA Mullins RD Cell mechanics and the cytoskeleton Nature 2010 463 485 492 10.1038/nature08908 20110992
Fletcher, D. A. & Mullins, R. D. Cell mechanics and the cytoskeleton. Nature 463, 485–492 (2010).20110992
7. Picu, C. R. Network Materials: Structure and Properties. (Cambridge University Press, Cambridge, 2022).
8. Simard SW Net transfer of carbon between ectomycorrhizal tree species in the field Nature 1997 388 579 582 10.1038/41557
Simard, S. W. et al. Net transfer of carbon between ectomycorrhizal tree species in the field. Nature 388, 579–582 (1997).
9. Steidinger BS Climatic controls of decomposition drive the global biogeography of forest-tree symbioses Nature 2019 569 404 408 10.1038/s41586-019-1128-0 31092941
Steidinger, B. S. et al. Climatic controls of decomposition drive the global biogeography of forest-tree symbioses. Nature 569, 404–408 (2019).31092941
10. De Domenico M Mathematical formulation of multilayer networks Phys. Rev. X 2013 3 041022
De Domenico, M. et al. Mathematical formulation of multilayer networks. Phys. Rev. X 3, 041022 (2013).
11. Bianconi, G. Multilayer Networks: Structure and Function. (Oxford University Press, Oxford, 2018).
12. Van Mieghem, P. Graph Spectra for Complex Networks. (Cambridge University Press, Cambridge, 2010).
13. Castellano C Pastor-Satorras R Relating topological determinants of complex networks to their spectral properties: structural and dynamical effects Phys. Rev. X 2017 7 041024
Castellano, C. & Pastor-Satorras, R. Relating topological determinants of complex networks to their spectral properties: structural and dynamical effects. Phys. Rev. X 7, 041024 (2017).
14. Villegas, P., Gili, T., Caldarelli, G., and Gabrielli, A. Laplacian renormalization group for heterogeneous networks. Nat. Phys. 19 , 445–450, (2023).
15. Scheffer LK A connectome and analysis of the adult drosophila central brain Elife 2020 9 e57443 10.7554/eLife.57443 32880371
Scheffer, L. K. et al. A connectome and analysis of the adult drosophila central brain. Elife 9, e57443 (2020).32880371
16. Tamassia, R. editor. Handbook of Graph Drawing and Visualization. (CRC Press, Boca Raton, Fl, (2013).
17. Vicsek, T. Fractal Growth Phenomena. 2nd Edn (World Scientific, Singapore,1992).
18. Bunde, A. and Havlin, S. (eds). Fractals and Disordered Systems. (Springer Berlin, Heidelberg (2012).
19. de Gennes P-G Exponents for the excluded volume problem as derived by the wilson method Phys. Lett. A 1972 38 339 340 10.1016/0375-9601(72)90149-1
de Gennes, P.-G. Exponents for the excluded volume problem as derived by the wilson method. Phys. Lett. A 38, 339–340 (1972).
20. Lawler GF A self-avoiding random walk Duke Math. J. 1980 47 655 693 10.1215/S0012-7094-80-04741-9
Lawler, G. F. A self-avoiding random walk. Duke Math. J. 47, 655–693 (1980).
21. Pietronero L Survival probability for kinetic self-avoiding walks Phys. Rev. Lett. 1985 55 2025 10.1103/PhysRevLett.55.2025 10031990
Pietronero, L. Survival probability for kinetic self-avoiding walks. Phys. Rev. Lett. 55, 2025 (1985).10031990
22. Schramm O Scaling limits of loop-erased random walks and uniform spanning trees Isr. J. Math. 2000 118 221 288 10.1007/BF02803524
Schramm, O. Scaling limits of loop-erased random walks and uniform spanning trees. Isr. J. Math. 118, 221–288 (2000).
23. Niemeyer L Pietronero L Wiesmann HJ Fractal dimension of dielectric breakdown Phys. Rev. Lett. 1984 52 1033 10.1103/PhysRevLett.52.1033
Niemeyer, L., Pietronero, L. & Wiesmann, H. J. Fractal dimension of dielectric breakdown. Phys. Rev. Lett. 52, 1033 (1984).
24. Wilson, D. B. Generating random spanning trees more quickly than the cover time. In Proceed of the 28th ACM Theory of computing, 296–303. 10.1145/237814.237880 (1996).
25. Lawler GF Schramm O Werner W Conformal invariance of planar loop-erased random walks and uniform spanning trees Ann. Probab. 2004 32 939 995 10.1214/aop/1079021469
Lawler, G. F., Schramm, O. & Werner, W. Conformal invariance of planar loop-erased random walks and uniform spanning trees. Ann. Probab. 32, 939–995 (2004).
26. Wiese K Fedorenko AA Field theories for loop-erased random walks Nucl. Phys. B 2019 946 114696 10.1016/j.nuclphysb.2019.114696
Wiese, K. & Fedorenko, A. A. Field theories for loop-erased random walks. Nucl. Phys. B 946, 114696 (2019).
27. Agrawal H Dhar D Distribution of sizes of erased loops of loop-erased random walks in two and three dimensions Phys. Rev. E 2001 63 056115 10.1103/PhysRevE.63.056115
Agrawal, H. & Dhar, D. Distribution of sizes of erased loops of loop-erased random walks in two and three dimensions. Phys. Rev. E 63, 056115 (2001).
28. Grassberger P Scaling of loop-erased walks in 2 to 4 dimensions J. Stat. Phys. 2009 136 399 404 10.1007/s10955-009-9787-0
Grassberger, P. Scaling of loop-erased walks in 2 to 4 dimensions. J. Stat. Phys. 136, 399–404 (2009).
29. Wilson DB Dimension of the loop-erased random walk in three dimensions Phys. Rev. E 2010 82 062102 10.1103/PhysRevE.82.062102
Wilson, D. B. Dimension of the loop-erased random walk in three dimensions. Phys. Rev. E 82, 062102 (2010).
30. Barrat, A., Barthelemy, M., and Vespignani, A. Dynamical Processes on Complex Networks. (Cambridge University Press, Cambridge, 2008).
31. Masuda N Porter MA Lambiotte R Random walks and diffusion on networks Phys. Rep. 2017 716 1 58 10.1016/j.physrep.2017.07.007
Masuda, N., Porter, M. A. & Lambiotte, R. Random walks and diffusion on networks. Phys. Rep. 716, 1–58 (2017).
32. De Domenico M Biamonte J Spectral entropies as information-theoretic tools for complex network comparison Phys. Rev. X 2016 6 041062
De Domenico, M. & Biamonte, J. Spectral entropies as information-theoretic tools for complex network comparison. Phys. Rev. X 6, 041062 (2016).
33. Arenas A Díaz-Guilera A Kurths J Moreno Y Zhou C Synchronization in complex networks Phys. Rep. 2008 469 93 153 10.1016/j.physrep.2008.09.002
Arenas, A., Díaz-Guilera, A., Kurths, J., Moreno, Y. & Zhou, C. Synchronization in complex networks. Phys. Rep. 469, 93–153 (2008).
34. Boguna M Network geometry Nat. Rev. Phys. 2021 3 114 135 10.1038/s42254-020-00264-4
Boguna, M. et al. Network geometry. Nat. Rev. Phys. 3, 114–135 (2021).
35. Ghavasieh A Stella M Biamonte J De Domenico M Unraveling the effects of multiscale network entanglement on empirical systems Commun. Phys. 2021 4 129 10.1038/s42005-021-00633-0
Ghavasieh, A., Stella, M., Biamonte, J. & De Domenico, M. Unraveling the effects of multiscale network entanglement on empirical systems. Commun. Phys. 4, 129 (2021).
36. Villegas P Gabrielli A Santucci F Caldarelli G Gili T Laplacian paths in complex networks: Information core emerges from entropic transitions Phys. Rev. Res. 2022 4 033196 10.1103/PhysRevResearch.4.033196
Villegas, P., Gabrielli, A., Santucci, F., Caldarelli, G. & Gili, T. Laplacian paths in complex networks: Information core emerges from entropic transitions. Phys. Rev. Res. 4, 033196 (2022).
37. Ghavasieh A De Domenico M Generalized network density matrices for analysis of multiscale functional diversity Phys. Rev. E 2023 107 044304 10.1103/PhysRevE.107.044304 37198772
Ghavasieh, A. & De Domenico, M. Generalized network density matrices for analysis of multiscale functional diversity. Phys. Rev. E 107, 044304 (2023).37198772
38. Gomez S Diffusion dynamics on multiplex networks Phys. Rev. Lett. 2013 110 028701 10.1103/PhysRevLett.110.028701 23383947
Gomez, S. et al. Diffusion dynamics on multiplex networks. Phys. Rev. Lett. 110, 028701 (2013).23383947
39. Sole-Ribalta A Spectral properties of the laplacian of multiplex networks Phys. Rev. E 2013 88 032807 10.1103/PhysRevE.88.032807
Sole-Ribalta, A. et al. Spectral properties of the laplacian of multiplex networks. Phys. Rev. E 88, 032807 (2013).
40. Radicchi F Arenas A Abrupt transition in the structural formation of interconnected networks Nat. Phys. 2013 9 717 720 10.1038/nphys2761
Radicchi, F. & Arenas, A. Abrupt transition in the structural formation of interconnected networks. Nat. Phys. 9, 717–720 (2013).
41. Pastor-Satorras R Castellano C Distinct types of eigenvector localization in networks Sci. Rep. 2016 6 18847 10.1038/srep18847 26754565
Pastor-Satorras, R. & Castellano, C. Distinct types of eigenvector localization in networks. Sci. Rep. 6, 18847 (2016).26754565
42. Hata S Nakao H Localization of laplacian eigenvectors on random networks Sci. Rep. 2017 7 1 11 28127051
Hata, S. & Nakao, H. Localization of laplacian eigenvectors on random networks. Sci. Rep. 7, 1–11 (2017).28127051
43. Clauset A Shalizi CosmaRohilla Newman MarkEJ Power-law distributions in empirical data SIAM Rev. 2009 51 661 703 10.1137/070710111
Clauset, A., Shalizi, CosmaRohilla & Newman, MarkE. J. Power-law distributions in empirical data. SIAM Rev. 51, 661–703 (2009).
44. Alstott J Bullmore ED Plenz D powerlaw: a python package for analysis of heavy-tailed distributions PloS one 2014 9 e85777 10.1371/journal.pone.0085777 24489671
Alstott, J., Bullmore, E. D. & Plenz, D. powerlaw: a python package for analysis of heavy-tailed distributions. PloS one 9, e85777 (2014).24489671
45. Halu A Mondragón RaúlJ Panzarasa P Bianconi G Multiplex pagerank PloS one 2013 8 e78293 10.1371/journal.pone.0078293 24205186
Halu, A., Mondragón, RaúlJ., Panzarasa, P. & Bianconi, G. Multiplex pagerank. PloS one 8, e78293 (2013).24205186
46. Solé-Ribalta, A., De Domenico, M., Gómez, S., and Arenas, A. Centrality rankings in multiplex networks. In Proceedings of the 2014 ACM conference on Web science, 149–155. 10.1145/2615569.2615687 (2014).
47. Iacovacci J Rahmede C Arenas A Bianconi G Functional multiplex pagerank Europhys. Lett. 2016 116 28004 10.1209/0295-5075/116/28004
Iacovacci, J., Rahmede, C., Arenas, A. & Bianconi, G. Functional multiplex pagerank. Europhys. Lett. 116, 28004 (2016).
48. Lyklema JW Evertsz C Pietronero L The Laplacian random walk EPL (Europhys. Lett.) 1986 2 77 10.1209/0295-5075/2/2/001
Lyklema, J. W., Evertsz, C. & Pietronero, L. The Laplacian random walk. EPL (Europhys. Lett.) 2, 77 (1986).
49. Lawler GF Loop-erased self-avoiding random walk and the Laplacian random walk J. Phys. A: Math. Gen. 1987 20 4565 10.1088/0305-4470/20/13/056
Lawler, G. F. Loop-erased self-avoiding random walk and the Laplacian random walk. J. Phys. A: Math. Gen. 20, 4565 (1987).
50. Bhupatiraju, S., Hanson, J., and Járai, A. A. Inequalities for critical exponents in d-dimensional sandpiles. 10.48550/arXiv.1602.06475 (2017).
51. Pósfai, Márton posfaim/physnets_as_net-o-nets: physical networks as network-of-networks. 10.5281/zenodo.11140782, (2024).
