Fluctuations of motifs and non self-averaging in complex networks. A self- vs non-self-averaging phase transition scenario
eepl draft
Fluctuations of motifs and non self-averaging in complex networksA self- vs non-self-averaging phase transition scenario
M. Ostilli Cooperative Association for Internet Data Analysis, San Diego Supercomputer Center, UCSD, San Diego, CA
PACS – Structures and organization in complex systems
PACS – Networks and genealogical trees
PACS – Fluctuation phenomena, random processes, noise, and Brownian motion
Abstract – Complex networks have been mostly characterized from the point of view of thedegree distribution of their nodes and a few other motifs (or modules), with a special attentionto triangles and cliques. The most exotic phenomena have been observed when the exponent γ ofthe associated power law degree-distribution is sufficiently small. In particular, a zero percolationthreshold takes place for γ <
3, and an anomalous critical behavior sets in for γ <
5. In thisLetter we prove that in sparse scale-free networks characterized by a cut-off scaling with thesistem size N , relative fluctuations are actually never negligible: given a motif Γ, we analyzethe relative fluctuations R Γ of the associated density of Γ, and we show that there exists aninterval in γ , [ γ , γ ], where R Γ does not go to zero in the thermodynamic limit, where γ ≈ k min and γ ≈ k max , k min and k max being the smallest and the largest degree of Γ, respectively.Remarkably, in ( γ , γ ) R Γ diverges, implying the instability of Γ to small perturbations. [email protected] Introduction. –
In the last decade, several complexnetworks models have been proposed to explain and repro-duce the widespread presence of real-world networks [1–4].Observed real-world networks are the output of certainrandom processes. Therefore, a complex network modelshould reproduce not only the same observed averages butalso the same observed sample to sample fluctuations, ifavailable. If data about fluctuations are not available, themodel should remain maximally random around the ob-served averages [5] or make use of a minimal number ofassumptions (a null model approach).In all branches of physics, fluctuations have played acrucial role in the understanding of the underlying phe-nomena. It is not excessive to say that any observable inphysics does not have any objective meaning without theevaluation of its fluctuations. Observables in networks arenot an exception. Particularly important are the motifs,also known as modules, patterns, or communities (Fig.1), i.e. , sub-graphs on which the functionality of the net-work largely depends [6–9]. In this Letter we provide a Here the term “community” is meant in a broad sense. In generala community is a motif, but not vice-versa [6]. first systematic analysis of the fluctuations of the densityof motifs in an analytically treatable model of networks,the hidden variable model [10–14], and we show that thesefluctuations, in certain regions of the model parameters,are far from being negligible. In fact, we discover the ex-
Fig. 1: Examples of motifs. The labels in Γ may specify thenumber of links of the motif; for motifs made of k fully con-nected nodes, k -cliques, we also use the symbol Γ kc ; d in Γ d stands for the presence of a diagonal; Γ kc × stands for two k -cliques sharing a common node. istence of a phase transition scenario with regions whererelative fluctuations are negligible, separated from regionswhere the relative fluctuations diverge. The practical con-p-1 a r X i v : . [ c ond - m a t . d i s - nn ] F e b . Ostilli 1sequences of this analysis are dramatic. A large classof real networks can be mapped toward suitable hiddenvariable models. For each real network realization sucha mapping requires measuring the hidden variables { h i } associated to each node. The sequence { h i } , however, inthe process that produces the real network, can vary, anddifferent real network realizations in general are character-ized by different sequences { h i } . It is this variability of thesequence { h i } that generates extremely large fluctuations.In particular, two large network realizations will presenttwo totally different community structures. More in gen-eral, large fluctuations translate in a effective instabilityof the substructures composing the network to small per-turbations. For the same reason, in simulating complexnetworks, if we are in the high fluctuation region, in orderto have a fair evaluation of the average of the density of amotif we need to make use of a very demanding statistics.Complex networks studies have been mostly focused onthe analysis of local averages and variances, but without asystematic study of fluctuations. Analysis of correlationsin complex networks have been in fact mainly confined tothe degree of adjacent nodes [3,4,11], or to the presence ofloops and cliques [15–17] which are a manifestation of cor-relations. Particular attention has been paid to the “con-figuration model” [1, 3], i.e. , the “uncorrelated” networkmodel” originating from all possible graphs constrainedto satisfy a given degree distribution exactly (hard ver-sion) [1, 3, 18, 19], or in average (soft version) [10–12, 20].It is well known that, in the thermodynamic limit, in theconfiguration model the presence of loops of finite lengthis negligible, provided the exponent γ of the associatedpower law degree-distribution P ( k ) ∼ k − γ is sufficientlylarge [3, 15]. In turn, this has led perhaps to the erro-neous conclusion that, in any synthetic or real network,fluctuations are always negligible for sufficiently large γ .In fact, exotic phenomena are believed to occur only when γ <
5, where an anomalous critical behavior sets in, es-pecially when γ <
3, where a zero percolation thresholdtakes place due to a divergent second moment of P ( k )[21]. In particular, fluctuations of motifs are assumed tobe negligible for γ enough large. In this Letter, by usingthe framework of hidden-variable models, we prove thatin sparse networks ( γ >
2) characterized by a cut-off scal-ing with the system size N , fluctuations are actually nevernegligible: given a motif Γ, we analyze the relative fluc-tuations R Γ of the density of Γ, and we show that thereexists an interval ( γ , γ ) where R Γ diverges in the ther-modynamic limit, where γ ≈ k min and γ ≈ k max , k min and k max being the smallest and the largest degree of Γ.As a consequence, in ( γ , γ ), measuring the density ofΓ in simulations is a hard problem, and Γ is unstable tosmall perturbations, a fact that in turn provides a key to The expression “uncorrelated network model” might be mislead-ing, but we keep it using for historical reasons. A network is saidto be “uncorrelated” if the degrees k and k (cid:48) at the ends of a linkare independent random variables. However, a lack of degree-degreecorrelations does not imply the absence of other correlations. understand the stability/instability of communities [6]. The hidden-variable scheme. The choice of thecut-off. –
Given N nodes, hidden variable models aredefined in the following way: i) to each node we associatea hidden variable h drawn from a given probability den-sity function (PDF) ρ ( h ); ii) between any pair of nodes,we assign, or not assign, a link, according to a given prob-ability p ( h, h (cid:48) ), where h and h (cid:48) are the hidden variablesassociated to the two nodes. The probability p ( h, h (cid:48) ) canbe any function of the h ’s, the only requirement being that0 ≤ p ( h, h (cid:48) ) ≤
1. It has been shown that, when p ( h, h (cid:48) )has the following form (or similar generalizations) p ( h, h (cid:48) ) = (cid:18) k s hh (cid:48) (cid:19) − , k s = (cid:112) N ¯ k, (1)and ¯ k is the wanted average degree, for large N , the actualdegree k of the nodes of the network realized with theabove scheme are distributed according to ρ with actualaverage degree equal to ¯ k . In particular, if we choose thefollowing PDF having support in [ h min , h max ] ρ ( h ) = a h − γ , h max ≥ h ≥ h min > , (2)with γ >
2, the degree-distribution of the resulting net-work will be a power law with exponent γ and, for N sufficiently large, the normalization constant a and theso called structural cut-off k s are a = ( γ − / ( h − γ min ),and k s = (cid:112) N h min ( γ − / ( γ − h max (cid:28) k s ,correlations of the generated network are negligible, and p ( h, h (cid:48) ) (cid:39) hh (cid:48) /k s , while for h max (cid:29) k s correlations can beimportant. The choice of the cut-off h max is in principlearbitrary. However, most of real-world networks show thatthe maximal degree scales according to the so called nat-ural cut-off: k max ∼ h nc = N / ( γ − . As a consequence,in several models of complex networks it was assumed thechoice h max = h nc , justified as empirical. We think how-ever that such an approach is wrong: the fact that in mostof the real-world networks k max ∼ h nc is due to a proba-bilistic effect, is not due to a rigid upper bound k ≤ h nc .In fact, by using order-statistics one finds that by drawing N degree values from a power law with exponent γ , thehighest degree in average scales just as (cid:104) k max (cid:105) ∼ N / ( γ − [22, 23]. More precisely, it is possible to prove that thePDF for the rescaled random variable k max /N / ( γ − isalso a power law with exponent γ [24]. Power law distri-butions always lead to important fluctuations. It is thenclear that empirical observations of k max must be takenwith care: k max is not a self-averaging variable and sam-ples in which k max (cid:29) N / ( γ − , even if extremely rare,do exist and, as we shall see, have dramatic effects onthe fluctuations of motifs. We stress that here we fol-low a null-model approach. We do not claim that the de-gree of all real networks must have a cut-off scaling with N ; there might be of course many other possible scalingswhose value depend on the details of the system (relatedto physical, biological, or economical constraints). How-ever, if the only information that we have from a givenp-2 self- vs non-self-averaging phase transition scenario in complex networksreal network of size N is that i) the degree obeys a power-law distribution with exponent γ , and ii) highest degreesscale in average as N / ( γ − , forcing the model to have aspecific cut-off other than N would introduce a bias. Infact, order statistics tells us that a lower cut-off wouldproduce (cid:104) k max (cid:105) (cid:28) N / ( γ − . This observation leads us tochoose h max ∼ N for the hidden-variable scheme (1)-(2).More precisely, if we consider as target degree distribu-tion a power-law P ( k ) ∝ k − γ with finite support k ≤ N ,in order to reproduce its characteristics from the hidden-variable model, we need to use a cut-off h max = O ( N λ )with λ ≥
1. In fact, any other choice implies a difference inthe scaling of the moments (cid:104) k n (cid:105) between the hidden vari-able model (1)-(2) and the target distribution P ( k ). Fig.2 shows this for the second moment. Similar plots hold forhigher moments. In conclusion, the minimal cut-off of themodel (1)-(2) able to reproduce the correct scaling of allthe moments of the target degree distribution P ( k ) is just h max = O ( N λ ) with λ = 1. In this paper we set therefore h max = N . We stress that with this choice highest degreeswill be still order N / ( γ − , but only on average. N < k > λ=2/3λ=0.8λ=1λ=1.2 From Target P(k)
Fig. 2: Behavior of (cid:10) k (cid:11) vs the system size N for γ = 2 . P ( k ) = ak − γ with k ≤ N , whereas theother plots correspond to the hidden variable model (1)-(2)with different choices of the cut-off h max = N λ : λ = 1 / ( γ −
1) =2 / λ = 0 . λ = 1,and λ = 1 .
2. For higher values of λ , the plots saturate toa curve that, on this scale, is indistinguishable from the case λ = 1 .
2. The plots of the hidden-variable model have beencalculated by numerical evaluation of the involved integrals: (cid:10) k (cid:11) (cid:39) N (cid:82) h max h min dhdh (cid:48) dh (cid:48)(cid:48) ρ ( h ) ρ ( h (cid:48) ) ρ ( h (cid:48)(cid:48) ) p ( h, h (cid:48) ) p ( h (cid:48) , h (cid:48)(cid:48) ) (seebelow for a more detailed analysis of these techniques). Fluctuations of Motifs. –
Given the parameters N , h min , ¯ k , and γ , the above hidden-variables scheme pro-duces an ensemble of networks which, in terms of a fewcharacteristics, like statistics of the degree and motifs, arein part representative of many real-world networks withthose given parameters. In the following we will indicatethe ensemble averages with the bracket symbol (cid:104)·(cid:105) . Theaverages are built by following the above steps (i) and (ii) of the hidden-variables scheme. Notice that each time wegenerate a network realization, we need to draw N hiddenvariables from the PDF ρ ( h ), and N ( N − / p ( h, h (cid:48) ). In terms of the adjacency matrix a i,i , tak-ing value 0 or 1 for the presence or not of a link betweennodes i and j , steps (i) and (ii) give (cid:104) a i,j (cid:105) = (cid:90) dh i dh j ρ ( h i ) ρ ( h j ) p ( h i , h j ) . (3)We will indicate by n Γ the density of the motif Γ ina network realization. As is known [11, 25], for γ > C = (cid:104) n Γ3 (cid:105) / (cid:104) n Γ2 (cid:105) .For example, for γ (cid:29) (cid:104) n Γ2 (cid:105) = O (1), while (cid:104) n Γ3 (cid:105) = O (1 /N ). More in general, the more the motif isclustered, the smaller is its density. Yet, for finite N , andfor any motif Γ, clustered or not, by tuning the param-eters h min , ¯ k and γ , one can set, within some freedom,a desired value of (cid:104) n Γ (cid:105) . However, as we shall see, thesample-to-sample fluctuations of n Γ can be unexpectedlylarge. Fluctuations of n Γ must be compared with the cor-responding average of n Γ , therefore we are going to analyzethe following standard ratio R Γ = (cid:10) n Γ (cid:11) − (cid:104) n Γ (cid:105) (cid:104) n Γ (cid:105) . (4)In general R Γ will strongly depends on Γ, N and γ . Whenlim N →∞ R Γ = 0 the network is said to be self-averagingwith respect to the motif density n Γ . In practical terms,when this occurs, even one single sample is enough toget by simulations an accurate estimation of the average (cid:104) n Γ (cid:105) , provided N is large enough. The behavior of R Γ with respect to the network size N is therefore of crucialimportance: if the network is not self-averaging withrespect to some motif Γ, the number of samples necessaryto get a good estimation of (cid:104) n Γ (cid:105) in simulations will haveto grow with N or, from another perspective, it is hardto generate only those samples whose density is closeto a target value, and a kind of hard searching problememerges. This aspect is in fact connected with spin-glassand NP-complete problems; we will see in fact that R Γ can be read as a susceptibility of a homogeneous system. Analysis of R Γ . – Given a motif Γ, the density of Γin a graph realization is n Γ = cN (cid:88) i k Γ ( i ) , (5)where k Γ ( i ) counts the number of motifs Γ passing throughthe node i . The coefficient c depends on the definitionof the motif considered and serves to avoid over-countingwhen the motif is symmetric. For example, if the motifΓ is the triangle, we set c = 1 /
3. If the motif is notsymmetric, we can establish to count only those motifsthat pass through a specif node of Γ. For example, ifp-3. Ostilli 1Γ is a triple (two consecutive links), we can set c = 1,but a motif contributes only when the center of the triplecoincides with i . However, since we are interested only inthe relative fluctuations R Γ , we do not need to specify itsince c , as well as any constant, does not play any rolefor R Γ . Let us consider now the numerator of Eq. (4).Note that the hidden variable scheme does not distinguishnodes, therefore we can make use of the fact that nodes areall statistically equivalent. By using this property, fromEq. (5) we get the following susceptibility (cid:10) n Γ (cid:11) − (cid:104) n Γ (cid:105) = c N (cid:110)(cid:10) k Γ ( i ) (cid:11) − (cid:104) k Γ ( i ) (cid:105) (cid:111) + c (cid:110) (cid:104) k Γ ( i ) k Γ ( j ) (cid:105) i (cid:54) = j − (cid:104) k Γ ( i ) (cid:105) (cid:111) , (6)where i and j represent two arbitrary distinct indices. Inthe rhs of Eq. (6) we have a self-term proportional to themotif variance, rescaled by the factor 1 /N , and a mixed-term that accounts for correlations between two motifscentered at two different nodes. Note that, in general,the self-term, despite appears to be order 1 /N , cannot beneglected. In fact, due to exact cancellations in the mixedterm, the mixed- and self-terms give contributions of thesame order of magnitude.For what follows, we find it convenient to introduce an-other symbol for the averages with respect to the PDF ρ ( h ) given by Eq. (2): if f ( · ) is any function of the hiddenvariables h , . . . , h N we define[ f ] = (cid:90) N (cid:89) i =1 dh i ρ ( h i ) f ( · ) . (7)In particular, from Eq. (3) we have (cid:104) a i,j (cid:105) = [ p ( h i , h j )].For γ > k s ∼ N / , therefore Eqs. (1) and (2) imply (cid:104) a i,j (cid:105) = [ p ( h i , h j )] = O (cid:0) N − (cid:1) . (8)Next we analyze Eq. (6) in a few crucial motifs. Link (Γ ). If Γ is the link k Γ ( i ) coincides with the stan-dard definition of degree of the node i . For a given graphrealization, corresponding to a given realization of the h ’s,in terms of adjacency matrix we have k Γ1 ( i ) = (cid:88) l (cid:54) = i a i,l . (9)By using Eq. (3), Eqs. (7)-(9), and the statistical equiva-lence of nodes, we have (cid:104) k Γ1 ( i ) (cid:105) = (cid:88) l (cid:54) = i (cid:104) a i,l (cid:105) = ( N −
1) [ p ( h , h )] . (10)Let us now consider the product k Γ1 ( i ) k Γ1 ( j ). Notice that a i,j = a i,j . We have to distinguish the cases i = j and i (cid:54) = j , see Fig. 3. For i = j we have (cid:10) k Γ1 ( i ) (cid:11) = ( N − N −
2) [ p ( h , h ) p ( h , h )]+( N −
1) [ p ( h , h )] , (11) Fig. 3: Contributions to Eqs. (11) (lower connected motifs)and (12) (upper disconnected motifs and the 3 connected motifslocated in the central part of the figure). Nodes i and j are tobe kept fixed, while the others can vary, provided the topologyis kept fixed. Contributions from disconnected motifs alwayscancel in R Γ . while for i (cid:54) = j we have (cid:104) k Γ1 ( i ) k Γ1 ( j ) (cid:105) = ( N − N −
3) [ p ( h , h )] +3( N −
1) [ p ( h , h ) p ( h , h )] , (12)where the factor 3 comes from the fact that two links em-anating from nodes i and j can share a same node in 3topologically equivalent ways. On plugging Eqs. (10)-(12) into Eq. (4) via Eq. (6) and keeping only terms in N , which cancel exactly, and terms in N , we obtain R Γ1 = 4 N [ p ( h , h ) p ( h , h )][ p ( h , h )] − N + 1 N [ p ( h , h )] . (13)Due to Eq. (8) and its generalizations, for γ >
3, eachterm present in the rhs of Eq. (13) is of order 1 /N ,therefore we have R Γ1 = O (1 /N ) and the network isself-averaging with respect to the link density, while for γ < R Γ1 decays sloweras N − γ . It is interesting however to observe the gen-eral behavior of R Γ1 with respect to γ for finite N . As ageneral rule, [ p ( h , h ) p ( h , h ) . . . p ( h m , h m +1 )] for large γ tends to factorize : [ p ( h , h ) p ( h , h ) . . . p ( h m , h m +1 )] → [ p ( h , h )] [ p ( h , h )] . . . [ p ( h m , h m +1 )]. Therefore, the firsttwo terms in the rhs of Eq. (13) tend to cancel for large γ . However, the last term does not cancel for large γ (thisissue will be discussed elsewhere). Diagrammatic calculus.
From Eq. (13) we see thatthe main term is given by the ratio between a [ · ]-averageof two links sharing a common node, and the square ofthe [ · ]-average of a single link, i.e. , our motif. From thisexample is clear that a correspondence between formu-las and diagrams can be established to avoid unnecessarysimulations and to improve our understanding about themain contributions to R Γ , especially those that can gen-erate non self-averaging. In this sense, we find it con-venient to make use of the compact notation [Γ], whereΓ can be any motif. For example, by referring to Fig.1, we have [Γ ] = [ p ( h , h )], [Γ ] = [ p ( h , h ) p ( h , h )],[Γ ] = [ p ( h , h ) p ( h , h ) p ( h , h )], and so on. The rolep-4 self- vs non-self-averaging phase transition scenario in complex networksplayed by these [ · ]-averages, is similar to the role playedby Green functions in statistical field theory. Moreover,since R Γ is defined in terms of connected correlation func-tions (6), we need to work only with Green functions ofconnected motifs, as the contributions of disconnected mo-tifs always cancel. Next, by using this diagrammatic tool,we evaluate R Γ in the crucial case of k -cliques. We firstanalyze the case k = 3 in detail, and then we look at thegeneral behavior R Γ kc , omitting contributions which arenot essential here. Further details will be given elsewhere. Triangle (Γ ). This is the simplest k -clique. We have R Γ3 = 9 N [Γ × ][Γ ] − N + 14 N [Γ d ][Γ ] + 6 N ] . (14)The factor 9 comes from: 8 ways to build Γ × from themixed term with i (cid:54) = j , and an extra contribution fromthe self-term i = j . Similarly to the last term of Eq. (13),the last two terms of Eq. (14) do not cancel for large γ . k -Clique (Γ kc ). Given Γ kc (( k −
1) is the degree of eachnode), if Γ kc × indicates the motif in which two k -cliquesΓ kc share a common node, we have R Γ kc = b k N [Γ kc × ][Γ kc ] − b k N + O (cid:18) N (cid:19) , (15)where b k is a combinatorial term which depends only on k , and the last term is positive and plays a role similar tothe last two terms of Eq. (14).In Fig. 4 we show the behavior of R Γ1 and R Γ3 vs N for γ = 4 .
2, and show the matching simulations vs di-agrammatic analysis. Notice that the theory (see nextparagraph) predicts R Γ3 → γ >
4, and R Γ3 → ∞ for2 . < γ <
4, however γ = 4 . R Γ3 decays very slowly with N . N -5 -4 -3 -2 -1 R Γ R Γ SimulationsR Γ Diagrammatic AnalysisR Γ SimulationsR Γ Diagrammatic Analysis
Fig. 4: R Γ1 and R Γ3 as functions of N for γ = 4 .
2. Simulations(circles and squares) made with S = 10 samples for N ∈ [10 , · ], S = 25 · for N = 5 · , S = 4 · for N = 10 , and S = 10 for N = 2 · . Diagrammatic analysis(triangles) made by numerical integrations of Eqs. (13) and(14) by using 10 points per integral. Singular terms.
The main message of this Letter is that,if the maximal degree k max of the motif Γ is greater than 1, there are contributions which make R Γ divergent for N → ∞ . Let us analyze the Green function [Γ kc × ] whichappears in Eq. (15). From Eqs. (1)-(2), by enumeratingthe 2 k − kc × with 1 , . . . , k −
1, 2 k − kc × ] = a k − k ( k ) s (cid:90) k − (cid:89) i =1 dh i k − (cid:89) i =1 h k − − γi h k − − γ k − × k − (cid:89) i =1 (cid:18) h k − h i k s (cid:19) − (cid:89) i Hidden variables models provide apowerful tool to investigate complex networks analytically.In this framework we present a first systematic analysisof the fluctuations of the density of motifs. Surprisingly,if the cut-off h max of the model is properly chosen as toreproduce a power-law distribution P ( k ) ∝ k − γ havingsupport k ≤ N , i.e. , h max = O ( N λ ), with λ ≥ 1, aphase transition scenario emerges, with self-averaging re-gions separated from non-self-averaging regions . Underthe simple null-model assumption that the system obeys apower-law, the potential practical consequences of such apicture are dramatic. The existence of non-self-averagingregions implies the instability of the substructures com-posing a network to even small perturbations: patternsand communities structures observed in a given network-realization, will have a totally different configuration inanother network-realization. For the same reason, in sim- We postpone the issue of the minimal value of λ , λ c , above whichsuch a picture still applies. However, for γ > 2, 1 ≥ λ c > / ( γ − p-5. Ostilli 1ulating complex networks, the presence of large fluctua-tions implies a very demanding statistics in order to havea fair evaluation of the averages of the observables of in-terest: given a motif Γ, if we are in a non-self-averagingregion, the number of samples that we need to properlyevaluate the average (cid:104) n Γ (cid:105) will be a growing function of thesystem size N . Hence, since the evaluation of n Γ alreadyrequires N operations per each sample, the total numberof operation to evaluate (cid:104) n Γ (cid:105) scales as N α , with α > γ [27]. In fact,with such small values, even small motifs and commu-nities are guaranteed to be stable as soon as they have k min (cid:38) γ . Whereas, for γ larger, communities for which k min (cid:46) γ (cid:46) k max will be unstable to small perturbations.In practical terms this means that, between two networks,one with say, γ = 2 . 3, and another with say, γ = 2 . 1, thelatter is more stable or, in other words, it has a largerprobability to exist. Moreover, our analysis also showsthat fluctuations of motifs become negligible for γ → γ → vs non-self-averaging phase transition scenario is expected to hold forany hidden-variable model characterized by power laws. Fig. 5: R Γ3 (top with γ (cid:39) . , γ (cid:39) 4) and R Γ4 c (bottom with γ (cid:39) . , γ (cid:39) . 7) as functions of N and γ , from numericalintegrations of Eq. (15) using 10 points per integral. ∗ ∗ ∗ Work supported by DARPA grant No. HR0011-12-1-0012; NSF grants No. CNS-0964236 and CNS-1039646; and by Cisco Systems. We thank D. Krioukov from whomthis research was inspired, and G. Bianconi, Z. Toroczkai,M. Bogu˜n´a, and C. Orsini for useful discussions. REFERENCES[1] B. Bollob´as, Random Graphs , 2nd ed., Cambridge Uni-versity Press (2001).[2] R. Albert, A.L. Barbasi, Rev. Mod. Phys. 47 (2002).[3] S.N. Dorogovtsev, J.F.F. Mendes, Evolution of Networks (University Press: Oxford, 2003).[4] M. E. J. Newman, Networks: An Introduction , (UniversityPress: Oxford, 2010).[5] J. Park, M. E. J. Newman, M. E. J., Phys. Rev. E ,066117 (2004).[6] S. Fortunato, Phys. Rep. , 75-174 (2010).[7] R. Milo, S. Shen-Orr, S. Itzkovitz, N. Kashtan, D.Chklovskii, U. Alon, Science 824 (2002).[8] A. Vazquez, R. Dobrin, D. Sergi, J.-P. Eckmann, Z. N.Oltvai, A.-L. Barab´asi, Proc. Natl Acad. Sci. USA 10 (2004).[10] G. Caldarelli, A. Capocci, P. De Los Rios, M. A. Mu˜noz,Phys. Rev. Lett. , 258702 (2002).[11] M. Bogu˜n´a and R. Pastor-Satorras, Phys. Rev. E ,036112 (2003).[12] J. Park and M. E. J. Newman, Phys. Rev. E , 066146(2004).[13] M. Catanzaro and R. Pastor-Satorras, Eur. Phys. J. B. , 241 (2005).[14] M. Catanzaro, M. Bogu˜n´a and R. Pastor-Satorras, Phys.Rev. E , 027103 (2005).[15] G. Bianconi and M. Marsili, J. Stat. Mech.: Theory Exp.,P06005 (2005).[16] G. Bianconi and M. Marsili, Europhys. Lett. , 740746(2006).[17] G. Bianconi and M. Marsili, Phys. Rev. E , 066127(2006).[18] S. N. Dorogovtsev, J. F. F. Mendes, A. N. SamukhinNucl.Phys. B 396 (2003).[19] G. Bianconi, Eurphys. Lett. E , 28005 (2008).[20] K.-I. Goh, B. Kahng, and D. Kim, Phys. Rev. Lett. ,278701 (2001).[21] S.N. Dorogovtsev, A.V. Goltsev, J.F.F. Mendes, Rev.Mod. Phys. , 1275 (2008).[22] S. N. Dorogovtsev, J. F. F. Mendes, and A. N. Samukhin,Phys. Rev. E ,026120 (2012).[26] C. I. Del Genio, T. Gross, and K. E. Bassler, Phys. Rev.Lett. , 178701 (2011).[27] SIAM Rev., 51(4), 661703 A. Clauset, C. R. Shalizi, andM. E. J. Newman, SIAM Rev. , 661 (2009)., 661 (2009).