Evolution of Robustness in Growing Random Networks

Networks are widely used to model the interaction between individual dynamic systems. In many instances, the total number of units and interaction coupling are not fixed in time, and instead constantly evolve. In networks, this means that the number of nodes and edges both change over time. Various properties of coupled dynamic systems, such as their robustness against noise, essentially depend on the structure of the interaction network. Therefore, it is of considerable interest to predict how these properties are affected when the network grows as well as their relationship to the growth mechanism. Here, we focus on the time evolution of a network’s Kirchhoff index. We derive closed-form expressions for its variation in various scenarios, including the addition of both edges and nodes. For the latter case, we investigate the evolution where single nodes with one or two edges connecting to existing nodes are added recursively to a network. In both cases, we derive the relations between the properties of the nodes to which the new node connects along with the global evolution of network robustness. In particular, we show how different scalings of the Kirchhoff index can be obtained as a function of the number of nodes. We illustrate and confirm this theory via numerical simulations of randomly growing networks.


INTRODUCTION
Complex networks are broadly used to model interaction within natural and engineered systems [1][2][3].They describe the interaction taking place between individual elements, such as the chemical bonds between atoms that form a molecule, or the communication transmitted between neighboring individuals in flocks of birds or vehicular platoons [4].From their structure important properties of the coupled dynamical systems can be deduced such as their intrinsic natural frequencies or their stability and robustness to external perturbations [5].While in many instances the structure of the coupling network as well as the number of interacting elements composing the system remain constant in time, it is typically not the case in a wide variety of coupled systems such as social networks, vehicular platoon formation, swarming autonomous robots, animal collective behaviors, cells evolution, molecules interacting in chemical reactions [6][7][8][9][10].In all these examples, when an element (commonly represented as a node) or an interaction (represented as an edge) is added to or removed from the system, its overall dynamical properties are modified.In particular, both the steady states and the corresponding transient stability are affected by the evolution of the system.It is therefore an important task to predict how these properties are changing while the network evolves, and be able to anticipate potential development of instabilities.More specifically, if one has to sequentially add agents to a system, how should they interact with the existing units so that stability is preserved, or at least not hindered to importantly.This is the main question that we investigate in this manuscript.Previous works have considered the evolution of some network properties such as the degree distribution in random growing networks with preferential attachment [11,12], or the evolution of the Wienner index in random recursive trees [13].In this manuscript, we investigate the time-evolution of the Kirchhoff index [14][15][16].The latter has proved to be useful in chemistry [14,17,18] and networked dynamical systems [19,20].For coupling networks that are not growing and which are static in time, the robustness of diffusively coupled oscillators have been direclty related to the Kirchhoff index of the effective coupling network [19,21,22].Namely, the larger the Kirchhoff index, the more important are the fluctuations within the dynamical system.Indeed, let us consider a set of N oscillators each with a continuous degree of freedom x i ∈ R , diffusively coupled together and subjected to noise as, where a ij = a ji > 0 are the elements of the adjacency matrix encoding the undirected coupling network, η i 's are white uncorrelated noise inputs, i.e. ⟨η i (t)η j (t ′ )⟩ = η 2 0 δ ij δ(t − t) .Then, the average variance in the long time limit is given by [23], arXiv:2308.05934v3 [nlin.AO] 17 Sep 2023 with Kf 1 the Kirchhoff index of the coupling network (see Sec. for the definition).Similar relations are obtained for deterministic perturbations that have a short correlation time [19].Given this direct connection between a global network index and the fluctuations of the dynamical system supported by the network, it is interesting to investigate how the Kirchhoff index evolve while the network grows.For the evolution of the network, we consider a simple growing algorithm where at each iteration, a single new node is added and connects to existing nodes.We derive analytical expression for the time-evolution of the Kirchhoff.In particular, when connecting the new node to the existing ones, we identify which of their nodal properties influence the scaling of the Kirchhoff index as a function of the total number of nodes.One can then use these properties when adding new nodes to achieve different scalings for the Kirchhoff index and thus, for the fluctuations.The manuscript is organized as follows.In Sec. , we give the definition of the Kirchhoff index and discuss bounds previously derived.In Sec. , we consider growing networks and provide expressions for the time-evolution of the Kirchhoff index when edges and nodes are added.Finally, Sec.gives our conclusions and outlook.

Definitions
Let us consider a graph (called network in the following) G made of N vertices (called nodes in the following) and M edges.Each edge ϵ (ij) between two nodes i and j has an associated weight a ij > 0 .The network Laplacian matrix is commonly defined as L as L ij = −a ij if i ̸ = j and there exists an edge between nodes i and j , and L ij = 0 otherwise, and L ii = N k=1 a ik for i = 1, ...N .The Kirchhoff index (Kf 1 ) of an undirected network is given by the sum of the effective resistance distances (Ω ij ) between all the nodes, i.e. [14] The resistance distance between node i and j is defined by, where L † denotes the pseudo-inverse of the Laplacian matrix L of the network.Using the eigenvectors u α and eigenvalues λ 1 = 0 < λ 2 < ... < λ N , of L , one conveniently rewrite the resistance distance Ω ij = α>1 (u α,i − u α,j ) 2 /λ α , which yields for the Kirchhoff index [24], Depending on the time-scale of the noise input, the amplitude of the small fluctuations of diffusively coupled oscillators can be expressed in terms of the Kirchhoff index or its generalization which reads [19], For specific network models for which the spectrum is known, the Kirchhoff can be analytically obtained.For example, for the complete, star and cycle network one has Kf 1 ∼ = N, N 2 , N 3 , and Kf 2 ∼ = 1, N 2 , N 5 respectively as the number of nodes N becomes large.These network models will be useful below when we consider the limiting case of random growing networks.In the specific case where the network is a tree, the resistance distance is equal to the shortest path distance in the same network where all the weights on the edges have been replaced by their inverse.In such a situation, the Kirchhoff index is equal to the Wienner index [17], that is defined as the sum of all shortest path distances in the network.In the following we discuss the Kirchhoff index, as the growing networks we consider do not have to be trees.From the resistance distance, one also defines a centrality measure that reads, where C(i) is called the resistance centrality of node i .

Lower bound on Kf1
The Kirchhoff index has been extensively studied and many bounds depending on the number of nodes N , edges n e and other properties have been derived.Relevant for the following, is a lower bound obtained by Zhou and Trinajstić [25] which is: for a connected network with N ≥ 3 , n e edges and a maximum degree ∆ , then From this inequality, one concludes that, as long as n e ∝ N , then Kf 1 /N scales at least as N when the number of nodes becomes large.This is the case in the growing algorithm we investigate below, where a single new node is added at each iteration, such that n e ∝ N .Therefore, the lowest scaling achievable for Kf 1 /N within our growing algorithm is linear in N .

ROBUSTNESS OF GROWING NETWORKS
Networks can grow in two ways: (i) some new nodes are connected to the existing network nodes; (ii) some edges are added within the existing nodes.On one hand, for (i) intuitively, based on the examples given in Sec. and Eq. ( 8), the Kirchhoff index increases at least linearly with N .On the other hand, for (ii) one can show that the Kirchhoff index can only decrease by adding a new edge in the network.Indeed, adding one edge with corresponding weight a kl > 0 between nodes k and l is a rank-1 modification of the Laplacian matrix, i.e.L(t + 1) = L(t) + a kl e kl e ⊤ kl , where [e kl ] i = (δ ik − δ il ) ∈ R Nt , with N t the number of nodes at iteration t .Therefore, if the Kirchhoff index at iteration t is Kf (t) , one can use the Sherman-Morrison-Woodbury formula [26,27] to obtain the Kirchhoff index at step t + 1 , where Ω (2) kl (t) and Ω kl (t) are always positive, Kf 1 can only decrease when an edge is added to the existing network.Below, we discuss how the Kirchhoff index is modified when a single node is added to the existing network together with m new edges.
One new node with a single connection (m = 1) We investigate the evolution of the Kirchhoff index for the growing process depicted on Fig. 1 .When the new node connects to a single existing one, if one starts with a network that is a tree, then the network will remain a tree as it grows.In addition to that, if the selected existing node is uniformly chosen at random, the latter is also called a random recursive tree.In such a situation, one can replace the resistance distance by the geodesic or shortest path distance (and thus the Kirchhoff index by the Wienner index) in the following discussion [13,28].But in general, we do not assume that the starting network is a tree.At iteration t one has the Kirchhoff index Kf If the new node at the iteration t + 1 is connected to node k , one has, where a k(Nt+1) = a new is the weight of the newly added edge between nodes k and N t+1 .In this simple case, one observes that the modification of the Kirchhoff index is given by the sum of the resistance distances from node k to all the other already existing nodes in the network plus N t times the resistance of the newly added edge (see Fig. 1).The less central node k is (in terms of resistance distances with respect to the existing nodes), the more the Kirchhoff index grows.As expected, Kf 1 (t) only increases with the iterations, as no new path within the existing nodes is created.If the node to which the new node connects is uniformly selected at random amongst the existing one at each iteration, then on average the Kirchhoff index will increase as, where N 0 and Kf 1 (0) are respectively, the initial number of nodes and the initial Kirchhoff index, H N = N k=1 k −1 is the N -th harmonic number which can be written as H N = γ + ψ 0 (N + 1) where γ ∼ = 0.577 is the Euler-Mascheroni number and ψ 0 (n) = Γ ′ (n)/Γ(n) is the digamma function.As its integer argument becomes large, the digamma function satifies ψ 0 (n) n→∞ ∝ ln n .Therefore, when the number of iterations becomes large, the last term in Eq. ( 13) will dominate such that, The scaling is confirmed numerically in Fig. 2, where the solid green curves are 20 realizations of the random growing process where one starts from ten connected nodes and then one new node is added at each iteration, that connects uniformly at random to an existing one.One observes that the simulations follow the predicted scaling of Eq. ( 14) given by the dashed black line.Note that it is the same scaling as the Wienner index for random recursive trees [13].This random evolution of the network is bounded by the two limiting cases that we now discuss.Instead of uniformly picking within the existing nodes, one may use some property of the nodes.Let us discuss what happens when one selects the most or least central node in terms of resistance distance.If the least central node k at each iteration t , i.e. with largest Nt l=1 Ω kl (t) , is chosen to be connected to the new node, the network will tend to form a chain.Therefore, assuming that the weights on the edges are of order 1 , when N t becomes large, one has In this case, the Kirchhoff index grows as, which is faster than in the random uniform case Eq. ( 14) .If instead, one selects the most central at each time step, then the network becomes star-like.Indeed, by connecting a new node to the most central existing one, its centrality becomes even more important.This means that all the newly added nodes will connect to the same node.Thus, assuming that the edges weights are of order 1 , one has for large N t , Interestingly, by selecting the most central node, one achieves a scaling for Kf 1 with N t only log N t better than Eq.( 14) where the node is uniformly chosen amongst the existing ones.

FIG. 2. Evolution of the Kirchhoff index divided by the number of nodes
Nt when a new node connects to a single existing one at each iteration.The initial network has 10 nodes and is obtained from a Watts-Strogatz rewiring procedure with nearest neighbors coupling [29].The green curves correspond to 20 realizations where one starts from the initial network, to which nodes are then recursively added, uniformly selecting at random the existing nodes to which they connect.For large Nt , the green curves follow the scaling of Eq. ( 14).The red and blue curves are obtained by selecting at each iteration the least and most central existing node, respectively.When Nt is large, they follow the scalings of Eqs. ( 16), (18) .The dotted, dashed, dotted-dashed black lines give N 2 t , Nt logNt , Nt , respectively.FIG. 3. Evolution of the network from iteration t to t + 1 , where a new node (in red) connecting to two existing ones (in black), k and l , has been added.The label of the new node is Nt+1 = Nt + 1 .In this case, a new path between k and l has been created.

Discussion
According to Eqs. ( 14), ( 16), (18), the weakest growth in the Kirchhoff index is obtained when the new nodes simply connect to the most central existing one in terms of resistance distance.Using this mechanism to grow a network leads to a very specific structure where a single node is connected to almost all the other ones.While such structure enhance the transient stability of the system by minimizing the growth of the small fluctuations, it also makes the system very vulnerable to any failure of this particular most central node.Indeed, if one removes that node from the system, then most of its components will become disconnected.Given this structural weakness, selecting nodes at random when adding new nodes seems like a more robust option, as the growth of the Kirchhoff index is only log N t worth than selecting the most central one.Moreover, the connection within the network are more uniformly distributed in that situation, reducing the number of disconnected components in case of failures.On the opposite side, if one wants to increase as much as possible the fluctuations in the system, new nodes should be connected to the least central node in terms of resistance distance.
One new node with two connections (m = 2) The case where one node is added together with two edges at each iteration is more complex as an increasing number of loops is introduced into the network.If the new node is connected to the existing nodes k and l , the effective resistance along the new path going from k to l is, The process is depicted on Fig. 3.It is important to remark that, adding node N t+1 is not the same as adding an edge between nodes k and l whose weight would be ω −1 kl .However, replacing the path where the new node is located by an equivalent edge provides a lower bound on the new Kirchhoff index.Indeed, by doing so, one obtains the sum of new resistance distances between the already existing nodes, The variation of the Kirchhoff index is a function of Ω (2) kl (t) and of how central in terms of resistance distance the new node is [see last term in Eq. ( 21)].For the latter term, it is challenging to find a closed form expression.However, one can have an estimate of it based on the resistance distance in the new network.Indeed, once the new node has been added, the resistance distances within the existing nodes at iteration t become where we replaced the new node by an equivalent edge between k and l using the Sherman-Morrison-Woodbury formula similarly as in Eq. ( 9) .Using Eq. ( 22), one can approximate the last term in Eq. ( 21) as the weighted average, We expect this approximation to be valid when the edge weights surrounding the new node, including a k Nt+1 , a l Nt+1 are homogeneous enough, or when a k Nt+1 , a l Nt+1 are much larger than the surrounding edge weights.Indeed, if the a k Nt+1 , a l Nt+1 are weak, the centrality of the new node is expected to be lower than those of k and l .Using this approximation together with Eq. ( 23), yields for Eq. ( 21), where we used the relation between L † ii and the resistance centrality of node i [see Eq. ( 7)].Now this expression gives an approximation of Kf 1 (t + 1) based only on quantities at iteration t .To reduce the increase of Kf 1 , one should therefore find nodes k and l such that Ω (2) kl (t) and Ω kl (t) are large, and which have very different resistance centralities, e.g.k being part of the most central nodes while l belongs to the least central ones.We group together the terms in Eqs.(27) as, Then, one can choose nodes k and l that minimize/maximize the latter quantities.To get more intuition, we investigate numerically Eqs. ( 28) and (29).In particular, we consider the maximization or minimization at each iteration of µ kl (t) , ρ kl (t) , ρ kl (t) − µ kl (t) .This is shown in Fig. 4. We consider edge weights such that a k Nt+1 = a l Nt+1 = 1 , meaning that the last term in µ kl (t) vanishes.As expected, maximizing µ kl (t) + ρ kl (t) (red curve) at each iteration gives the most important increase in Kf 1 (t)/N t which scales as N  [29].The black dashed-dotted and dashed lines give the scalings Nt and N 2 t , respectively.Note that in these simulations we make sure that k ̸ = l but found similar scalings when relaxing this condition.
ρ kl (t) (orange curve).The minimization of µ kl (yellow curve) does not produce such an increase of the Kirchhoff index, which seems to remain linear i.e.Kf 1 (t)/N t ∝ N t as t becomes large.A similar linear scaling is observed for the minimization of ρ kl (t) − µ kl (t) (blue curve), ρ kl(t) (cyan curve) and the maximization of µ kl (t) (green curve).Interestingly, one observes that the maximization of µ kl (t) gives a lower Kirchhoff index than the minimization of ρ kl (t) − µ kl (t) .Therefore, one can tune the increase of the Kirchhoff index by choosing one or another quantity to optimize at each iteration.In Fig. 5, we simulate the case where k and l at uniformly chosen are random at each iteration.One observes that the 20 realizations of the process yield a linear scaling of Kf 1 (t)/N t with N t .

Discussion
In Sec. , we saw that selecting uniformly at random the existing node to which the new one connects, is log N t worth for the growth of the Kirchhoff index compared to selecting the most central one.Quite interestingly here, when the new nodes connect to two existing ones, selecting them uniformly at random produces that same scaling as the one obtained by minimizing the relevant quantity ρ kl (t) − µ kl (t) .Therefore, when growing a network by connecting the new node to two existing ones, one should only make sure that the nodes are selected uniformly at random to achieve the best scaling.Of course, the latter is true as long as the approxiamtion in Eq. ( 23) holds.If instead one wants to disrupt the system, a scaling of the Kirchhoff index N t times worth compared to the previous situation is obtained by maximizing either ρ kl (t) − µ kl (t) or ρ kl (t) .The latter can be achieved by choosing nodes that are close in terms of Ω kl (t) and Ω (2) kl (t) and rather peripheral in the network i.e. small C(k, t) , C(l, t) .

Remark
One has to be careful when interpreting Eqs. ( 21), (24), and notice that the Kirchhoff index increases at least linearly with N t on average.This can be seen from Eq. ( 8).More intuitively, in the case m = N t , which means that the number of edges added at each iteration grows with the system size, assuming an initial all-to-all network one has, which monotonically increases.In this situation, one adds as many new paths between the existing nodes as possible when including a single new node.Therefore, for any m < N t , the Kirchhoff index must increase.In Eq. ( 21), one might however reduce the amplitude of the increase or sometimes even decrease Kf 1 by carefully selecting k and l .But the latter can only be true for a few iterations.

CONCLUSION
We considered the evolution of random networks where at each iteration, a new node is added and connected to one or two existing ones.When the new nodes only connects to one existing node, the scaling of Kf 1 (t)/N t with the number of nodes is between N t and N 2 t as the number of iterations becomes large.When the existing node is randomly uniformly chosen, the scaling is only logarithmically worse than the lower bound, i.e.Kf 1 (t)/N t t→∞ ∝ N t log N t .In the more complex situation where the new nodes connects to two existing ones, we derived a recursive expression for the evolution of Kf 1 (t) .The latter is essentially given by ρ kl (t) − µ kl (t) which can be expressed in terms of the resistance distances and centralities i.e.Ω kl , Ω kl , C −1 (k, t) , C −1 (l, t) , [see Eqs.( 28), (29)].We showed that, by introducing a bias in the selection of k and l toward the minimum/maximum of these quantities, one can tune the increase of Kf 1 (t)/N t from linear to quadratic in N t .For m > 2 , it is much more challenging to obtain analytical expression for the evolution of Kf 1 .The same applies to the case where m is a function of the number of nodes.Using the lower bound in Eq. ( 8) yet allows one to obtain the minimal scaling of the Kirchhoff index by correctly choosing n e (N ) .
The scenario we consider here applies to evolving systems where a single new node is added at each iteration and connects to existing nodes.This would represent the case where a new molecule forms bonds with another group of molecules, or an autonomous vehicle that joins a platoon by interacting with one or many of its members.With the results presented here, one can anticipate the scaling of the Kirchhoff index based on how the new units connect to the existing ones.Thus, they also give insights on the evolution of fluctuations in networked systems such as consensus dynamics and synchronized systems that are diffusively coupled.

FIG. 1 .
FIG.1.Evolution of the network from iteration t to t + 1 , where a new node (in red) connecting to a single existing one (in black), k has been added.The label of the new node is Nt+1 = Nt + 1 .No new path is create within the existing nodes.

FIG. 5 .
FIG.5.Evolution of the Kirchhoff index divided by the number of nodes Nt when a new node connects to two existing ones at each iteration.The two nodes are selected uniformly at random amongst the existing ones at each iteration.Each grey line (20 in total) is one realization of the process.The initial network has 10 nodes and is obtained from a Watts-Strogatz rewiring procedure with nearest neighbors coupling[29].The black dashed line gives the linear scaling with Nt .
2 t .The same scaling is obtained if one maximizes only FIG. 4. Evolution of the Kirchhoff index divided by the number of nodes Nt when a new node connects to two existing ones (k and l) at each iteration.Two nodes are selected by minimizing/maximizing µ kl (t) , ρ kl (t) , ρ kl (t) − µ kl (t) at each iteration.The meaning of each curve is given in the legend in insets.The initial network has 10 nodes and is obtained from a Watts-Strogatz rewiring procedure with nearest neighbors coupling