Connecting complex networks to nonadditive entropies

Boltzmann-Gibbs statistical mechanics applies satisfactorily to a plethora of systems. It fails however for complex systems generically involving strong space-time entanglement. Its generalization based on nonadditive $q$-entropies adequately handles a wide class of such systems. We show here that scale-invariant networks belong to this class. We numerically study a $d$-dimensional geographically located network with weighted links and exhibit its 'energy' distribution per site at its quasi-stationary state. Our results strongly suggest a correspondence between the random geometric problem and a class of thermal problems within the generalised thermostatistics. The Boltzmann-Gibbs exponential factor is generically substituted by its $q$-generalisation, and is recovered in the $q=1$ limit when the nonlocal effects fade away. The present connection should cross-fertilise experiments in both research areas.


Introduction
Boltzmann-Gibbs (BG) statistical mechanics constitutes one of the pillars of contemporary theoretical physics. As such is has uncountable successes for a great variety of physical systems. However, when the system constituents have a generically strong space-time entanglement, this theory does not apply. Such is the case already pointed in 1902 by Gibbs himself, namely when the standard partition function diverges, e.g., gravitation. It is in this context that it was proposed in 1988 1 the generalisation -hereafter referred to as nonextensive statistical mechanics -of the BG theory based on nonadditive entropies, namely S q = k 1−∑ i p q i q−1 , which recovers S BG = −k ∑ i p i ln p i in the q → 1 limit. The composition of two probabilistically independent systems A and B yields straightforwardly S q (A + B)/k = [S q (A)/k] + [S q (B)/k] + (1 − q)[S q (A)/k][S q (B)/k]. As we see, the BG entropic additivity is recovered when q = 1. The fundamental advantage associated with q = 1 is that, for strongly correlated systems, it enables, as illustrated in 2 , the preservation of the extensivity of the thermodynamic entropy, mandated in all circumstances by the Legendre structure of classical thermodynamics.
In parallel with the above, the study of complex networks has been intensified around the world [3][4][5][6][7] . Networks can be found everywhere. Society is formed by humans linked through relationships. The Internet is a set of devices communicating with each other. The brain is formed by neurones communicating through synapses. All these completely different systems can be translated onto a simple set of nodes (or sites) and edges (or links) obeying some connection rule, and the tools of network science can be successfully used to study them. Typical applications of this area can be found in classical and quantum internet 8,9 , medicine 10, 11 , neuroscience 12 , and sociology 13,14 . It was thought, during more than a decade, that most of real networks were purely scale-free meaning that the distribution of the number of links in the network follows a power-law distribution. It was recently argued that most real networks are not pure scale-free 15 , paving the way for new possibilities to describe them.
During the initial years, network science and nonextensive statistical mechanics were seen as completely different areas. But meaningful connections started in 2005 [16][17][18][19][20][21][22] . It is nowadays known that the degree distribution of asymptotically scale-free networks at the thermodynamic limit is of the form P(k) ∝ e −k/κ q , where the q-exponential function is defined as . This form, more precisely p q (ε i ) = e −β q ε i q /Z q , optimizes the entropy S q under appropriate canonical constraints, ε i being the site energy and β q the inverse temperature; the BG weight is recovered at the q → 1 limit. This thermostatistical approach has been successfully applied in a wide diversity of areas, such as long-range-interacting Hamiltonian systems 23 , vortices in type II supercondutors 24 , cold atoms 25 , granular matter 26 , high-energy physics experiments Sample of a N = 100 network for (d, α A , α G , η, w 0 ) = (2, 1, 5, 1, 1). As can be seen, for this choice of parameters, hubs (highly connected nodes) naturally emerge in the network. Each link has a specific width w i j and the total energy ε i associated to the site i will be given by half of the sum over all link widths connected to the site i (see zoom of site i).
on Earth 27 and observations in the outer space 28,29 , civil engineering 30 , and for predicting COVID-19 peaks around the world 31,32 .
In this work, we introduce and study a geographically located d-dimensional network model. One of the main characteristics of this model is the possibility to control the long/short range nature of the interactions between the sites. BG statistics completely fails to describe systems that interact at long-range, and many theories have been proposed to approach this regime. The present model introduces a new property for this class of systems. In addition to the fact that Euclidean distances (d i j ) between the sites are relevant, it also takes into account the weights (w i j ) of the links (see Fig. 1) and associates them to the 'energy' (ε i ) of each site. Due to that new ingredient, we could compute the energy distribution of the ever growing network. This distribution turns out to have the functional form of the q-generalised BG distribution for nonextensive systems, based on nonadditive entropy. These numerical results strongly suggest a neat correspondence between the random network geometrical problem and a particular thermostatistical problem within the generalised theory.

The model
Our growing d-dimensional network starts with one site at the origin. We then stochastically locate a second site (and then a third, a fourth, and so on up to N) through a probability p(r) ∝ 1/r d+α G (α G > 0), where r ≥ 1 is the Euclidean distance from the newly arrived site to the center of mass of the pre-existing cluster; α G is the growth parameter and d = 1, 2, 3 is the dimensionality of the system (large α G yields geographically concentrated networks).
The site i = 1 is then linked to the site j = 2. We sample a random number w i j from a distribution P(w) that will give us the corresponding link weight. Each site will have a total energy ε i that will depend on how many links it has, noted k i , and the widths {w i j } of those links. At each time step, the site i only has access to its local energy ε i defined as: The value of ε i will directly affect the probability of the site i to acquire new links. Indeed, from this step on, the sites i = 3, 4, ...

2/8
will be linked to the previous ones with probability where d i j is the Euclidean distance between i and j, where j runs over all sites linked to the site i. The attachment parameter α A controls the importance of the distance in the preferential attachment rule (2). When α A 1 the sites tends to connect to close neighbours, whereas α A 0 tends to generate distant connections all over the network. Notice that, while the network size increases up to N nodes, the variables k i and ε i (number of links and total energy of the i-th node; i = 1, 2, 3 . . . , N) also increase in time (see Fig. 1 for a sample of the ever growing network).
If we consider the particular case , thus recovering the usual preferential attachment rule. Consequently, the present model recovers the one in 19-21 as a particular instance. Note that, if we additionally consider the particular case α A = 0, we recover the standard Barabási-Albert model with Π i ∝ k i 5, 6 . We are considering here the case where w is given by the following stretched-exponential distribution: which satisfies ∞ 0 dw P(w) = 1. As particular cases of Eq. (3) we have: η = 1, which corresponds to an exponential distribution, η = 2, which corresponds to a half-Gaussian distribution, and η → ∞, which corresponds to an uniform distribution within w ∈ [0, w 0 ].

Results
Our focus here is to analyse the energy distribution p(ε) of the N 1 network. We have in fact analyzed a large amount of typical cases in the space (d, α A , α G , w 0 , η), and have systematically found the same results for d = 1, 2, 3 within the intervals . Similarly to previous works 16,[19][20][21] , p(ε) does not depend on α G ; also, it does not depend independently on d and α A , but only, remarkably, on the ratio α A /d; p(ε) also depends on w 0 and η (see Fig. 2(a)-(d)). Because of these features, and without loss of generality, we have once for ever fixed α G = 1, and d = 2. The simulations were done for 10 3 realisations of size N = 10 5 , which was verified to be enough for observing the asymptotic distribution p(ε) with high precision.
We know that the signature of the Boltzmannian systems is the presence of exponentials and Gaussians distributions. Similarly, the nonextensive systems based on the entropy S q can be recognised by the emergence of q-exponentials and q-Gaussians distributions. We have here found that, independent of the choice of (d, α A , α G , w 0 , η), the 'energy' distribution p(ε) associated with the network is invariably well fitted by the following q-exponential: where p q (ε) represents the generalisation, within nonextensive statistical mechanics, of the BG energy weight with ε, β q and Z q playing respectively the roles of energy, inverse temperature and normalisation factor (see Fig. 3). Note that, when q → 1, we do recover the standard Boltzmann distribution since e −β 1 ε 1 ≡ e −β ε . This result exhibits an interesting emergence of correspondence between a random network geometric problem and a particular case within generalised thermostatistics. This fact definitively reminds the Kasteleyn-Fortuin theorem 33 , which establishes an important isomorphism between the bond-percolation geometric problem and the q Potts → 1 limit of the q Potts -state Potts ferromagnet.
In Fig. 4(a,d) we show the numerical results for the index q as function of α A /d. This result is consistent with 19,20 , where the behaviour of q characterises the existence of three regimes. As can be seen, q is constant and equal to 4/3 in the range 0 ≤ α A /d ≤ 1. This interval describes the regime of strong-long-range interactions characterised by the highest value of q. In the interval 1 < α A /d 5 we have the moderately-long-range interaction regime, characterised by q smaller than 4/3 but still greater than 1. In this regime q displays no abrupt transition from 4/3 to 1 but instead it decreases exponentially with α A /d through the function e 1−α A /d 19, 20 . This behaviour exhibits that the BG regime was not yet reached. In the last regime, α A /d 5, the Boltzmannian-like regime finally emerges and q = 1. Similar results for q were found in 17 for a gas-like network model. In the Fig. 4 (b, c,e) we show similar results for the parameter β q which equals β q 0 in the range 0 ≤ α A /d ≤ 1 and, then, exponentially decreases with α A /d; β q increases with η and decreases with w 0 . However, if we plot (β q − β q ∞ )/(β q 0 − β q ∞ ), all curves collapse as a function of α A /d. Moreover, we verify in Fig. 4(f) that the effective 'temperature' decreases with increasing q. There is no thermodynamical prescription which would impose that. In 30 , for instance, both possibilities are in fact observed.

Discussion
All these results strongly suggest that the 'energy' distribution p(ε) of the network is given by the very same expression which q-generalises the Boltzmann-Gibbs weight when it is the nonadditive entropy S q which is optimised. Naturally, since , all curves collapse and exponentially decrease with α A /d > 1 (black straight line). Inset: β q 0 and β q ∞ exponentially vary with η; β q ∞ was estimated by fixing α A /d = 10. In (d) and (e) we present log-linear representations of the same data as in (a) and (c) respectively, exhibiting the exponential dependence of both q and β q on α A /d, when α A /d ≥ 1. ( f ) By eliminating the variable α A /d, we show 1/(w 0 β q ) as functions of q for the same set of data shown in the previous plots; q is related with 1/(w 0 β q ) through the equation 1/w 0 [β q ∞ + 9(β q 0 − β q ∞ )(q − 1) 2 ] that is valid for all values of (w 0 , η). For very large values of α A /d and the extreme regions η → 0 and η → ∞, the numerical precision needed to attain the stationary-state distribution is too high for our present computational effort and further analysis would be needed, which is out of the present scope. the present study is numerical, we can not exclude very minor corrections due to high-order metric-topological terms. The BG limit is rapidly reached when α A /d 5. Not less important, q and β q depend on α A and d only through the ratio α A /d; also, interestingly enough, none of them depends on α G . The fact that q depends only on α A /d means that this ratio uniquely determines the entropic nonadditivity universality class. The quantity β q also depends on η and w 0 . Consistently, w 0 , which characterises the width of the stretched-exponential distribution P(w), plays here the same role as T in usual thermal BG problems. This seemingly is the first time that, in a complex network, we identify a parameter which plays the role of an external parameter that we may fix at will, similarly to the manner in which we fix, in BG statistical mechanics, the temperature at which the thermally equilibrated system is placed. In all previous connections with random networks [16][17][18][19][20][21][22] , β q (sometimes noted 1/κ q ) was univocally related to q. A single value for β q for a given value of q is analogous to traditional critical points in BG statistics. In our present case, we have, for a fixed value of q, the freedom of making β q to vary, like T in BG thermal statistics.
The present results strongly support the conjecture of existence of a neat correspondence between geometrical random network (asymptotically) scale-invariant problems and the present specific class of many-body models within nonextensive statistical mechanics, constructed upon nonadditive entropies. This is analogous to the Kasteleyn-Fortuin theorem for the q Potts -state Potts model, whose q Potts → 1 limit rigorously corresponds to the bond percolation problem 33 , and also to the de Gennes celebrated isomorphism for the n-vector ferromagnetic model, whose n → 0 limit precisely corresponds to the self-avoiding random walk 34 , which constitutes a pillar in polymer physics. It is possible to think of a variety of applications of connections of the present kind, for example the maintaining budgets to be distributed among cities connected within a large regional network of roads. Each city could, for instance, receive a support proportional to the sum of the widths of the roads arriving to it.