Topological Vulnerability Evaluation Model Based on Fractal Dimension of Complex Networks

With an increasing emphasis on network security, much more attentions have been attracted to the vulnerability of complex networks. In this paper, the fractal dimension, which can reflect space-filling capacity of networks, is redefined as the origin moment of the edge betweenness to obtain a more reasonable evaluation of vulnerability. The proposed model combining multiple evaluation indexes not only overcomes the shortage of average edge betweenness’s failing to evaluate vulnerability of some special networks, but also characterizes the topological structure and highlights the space-filling capacity of networks. The applications to six US airline networks illustrate the practicality and effectiveness of our proposed method, and the comparisons with three other commonly used methods further validate the superiority of our proposed method.


Introduction
Complex network is widely used to model the structure of many complex systems in nature and society [1,2]. Some network models are used to solve natural problems, such as climate issues [3]. While some are used to analyze social and practical problems, such as supply chain management [4], transportation networks [5], power grid networks [6][7][8][9], water distribution networks [10], network optimization [11][12][13], as well as game theory in operation research [14,15] and etc. An open issue is how to assess the vulnerability of complex networks [16][17][18], whose main objective is to understand, predict, and even control the behavior of a networked system under vicious attacks or any types of dysfunctions [19].
Different approaches to characterize network vulnerability and robustness have recently been proposed, which can be grouped into two types broadly [9,20,21]. The first type of the approaches is structural robustness-topological properties of networks [22,23]. Such as, the average edge betweenness [19], the network connectivity level, the size of largest component [16] and the average geodesic length [6] and etc., are directly used to define network vulnerability. The second one concerns dynamical robustness [24][25][26]. The dysfunctions of a node or link will cause the redistribution on the the load of other nodes or links, with the risk that some other nodes or links may be overloaded, which will further cause a sequence of failures and even threaten the global stability [6,20,27]. Such behavior is called cascading failures [28][29][30]. In addition, the evaluation of network vulnerability is of uncertainty, and many methods have promising aspect to address such problem. For example, fuzzy set theory is efficient to model linguistic information [31] wihle Dempster-Shafer evidence theory can combine different information in an efficient manner [32][33][34]. The vulnerability analysis of complex systems, such as physic protect systems, is modelled by these mathematical tools [35].
One of the mostly used methods is proposed by Boccaletti et.al [19]. They construct a multi-scale evaluation model of vulnerability, which makes use of the average edge betweenness and introduces a key coefficient p. Due to its effectiveness, this method has been heavily studied [20]. One limitation of the original average edge betweenness model is that it cannot differentiate two different networks in some situations. To solve this problem, a key coefficient p is introduced to improve the original model in their work. However, a straight problem is how to determine the coefficient p. The way used in Boccaletti et.al's work is lack of physical significance.
The main motivation of our work is that we believe this coefficient p should be determined by the network itself but not just geometrically. To address this issue, we take the fractal dimension of complex network into consideration.
Dimension is one of the fundamental properties of complex networks characterizing not only its topological properties but also dynamic characteristics [36][37][38], which are two key aspects to determine network vulnerability. Analyses of a variety of real complex networks show that self-similar characteristic exists on all length scales, that is, many complex networks exhibit fractal properties [23,39,40]. Then researchers found that this self-similar characteristic can used to characterize dimension, i.e., fractal dimension, which gives a good definition to dimension of many geometric images and the network with fractal properties, such as Koch curve, Sierpinski triangle, the Coast of Britain, social networks, power grid networks, airline networks and so on [41][42][43]. In addition, fractal dimension has also been confirmed to be able to measure the space-filling capacity of a pattern reflecting the complexity of network directly [44], which has relationships with network vulnerability. Based on this idea, we propose that fractal dimension is a promising alternative to be the key coefficient in determining network vulnerability in this paper. This paper is organized as follows. Section 2 introduces the preliminaries. Section 3 presents details of the proposed method and the steps about its application. Section 4 compares the proposed method with the existing methods listed in other papers by calculating the vulnerability of them. Finally, we summarize our results in Section 5.

Preliminaries
In this section, we introduce Boccaletti et.al's vulnerability model [19] and three other commonly used methods of vulnerability evaluation [6,16,20]. In general, the complex networks can be represented by an undirected and unweighted graph G = (V, E), where V is the set of nodes and E is the set of edges. Each edge connects exactly one pair of nodes, and a vertex-pair can be connected by maximally one edge, i.e. loop is not allowed.

Multi-scale Evaluation of Vulnerability
In Boccaletti et.al's work [19], the original method to evaluate the vulnerability is represented by the average edge betweenness, which is defined as: Where |E| is the number of the edges, and b l is the edge betweenness of the edge l, define as: Where n jk is the total number of geodesics (shortest path) from node j to k, and n jk (l) is the number of geodesics from j to k that contain the link l. However, this evaluation method of b 1 (G) gives no relevant new information about the vulnerability of some special networks. For example, two networks referred in [19] shown in Fig 1 can't be distinguished using this method. By evaluating the vulnerability according to Eq 1, one gets b 1 (G) = b 1 (G 0 ) = 43/13. It's absolute that the "bat" graph G is more vulnerable than the "umbrella" graph G 0 , but Eq 1 gives the same evaluation results about them.
In order to overcome the original method's limitation of failing to distinguish some networks, a key coefficient p was introduced by Boccaletti et.al in the improved model, which is called multi-scale evaluation of vulnerability [19] and shown as below: If we want to compare two networks G and G 0 , first computes b 1 .
To get the coefficient p, Boccaletti et.al define a relative function of p like: The coefficient p is obtained when the function has a maximal value. For more detailed information to determine the coefficient p, refer [19]. It's clear that, the definition of coefficient p is based on geometrical definition and lack of physical significance. In our opinion, the coefficient p should be defined by the complex network itself, we take the fractal dimension as consideration.

Fractal Dimension
One of the most typical ways to calculate the fractal dimension is Box covering algorithm, a power-law relation between the number of boxes needed to cover the network and the size of the box [45][46][47]. For a given network G and box size l B , a box is a set of nodes where all distances l ij between any two nodes i and j in the box are smaller than l B . The minimum number of boxes required to cover the entire network is denoted by N B . The detailed illustration of the calculation of the fractal dimension referred in [46] is given in Fig 2. The fractal dimension or box dimension d B calculated with the box covering algorithm is given as follows [40,46]:

Comparison Preparation
For the sake of comparison, three other commonly used methods to calculate vulnerability are described as follows. The first method is the average inverse geodesic length l −1 [6]: Where d(v, w) is the length of the geodesic between v and w (v, w 2 V), and the factor N(N − 1) is the number of pairs of nodes. The larger l −1 is, the better the network functions. The second method is the size of largest component LCS (0 < LCS < 1) [16], which quantifies the number of nodes in the largest connected subgraph and is defined as follows: Where N s is the size of the largest connected subgraph. The larger LCS is, the more robust the network is.
And the third method is the normalized average edge betweenness b nor (G) [20], which is on the base of the Eq 3 while p = 1 and is defined as: Where G complete is a complete graph and G path is a path graph. The larger LCS is, the more vulnerable the network is.

Methods
In this section, the proposed method is detailed. As mentioned in introduction section, average edge betweenness can't highlight difference between some special networks, which is caused by the average process. Because for a network whose geodesic distribution concentrates strongly around a single link or node, the potential risk of the failure in this critical link/node will be hidden by the average process with the rest of minor links. For example, edge betweenness of some edges in "bat" graph are much higher than "umbrella" graph, which represents the "bat" graph is more vulnerable than "umbrella" graph, but the average process hides their differences and gives the same vulnerability of them. The improved multi-scale evaluation model is equivalent to the p origin moment of edge betweenness. Since p origin moment can highlight the difference between data to eliminate the effect of average process. But the problem is how to determine the key coefficient p to make effective and reasonable evaluation, which is the motivation of our work. In our opinion, the fractal dimension of the complex network is a promising alternative to redefine the coefficient p. It is well-known that the fractal dimension can characterize the network structure and basic physical properties including space-filling capacity. For a given network with fixed number of nodes, the higher the space-filling capacity is, the more edges connecting the nodes in this network, and then the smaller of the diameter of the network, thus the less boxes will be required to cover the whole network, i.e., the smaller the fractal dimension is. So the fractal dimension decreases with the space-filling capacity of network. We also know that given fixed number of nodes of a network, the more edges, the more connected and robust of this network. That is, the fractal dimension has direct relationships with network vulnerability. So using the fractal dimension to redefine p can not only highlight difference between networks properly, but also ensure practical significance of p, which is more reasonable than geometric coefficient p in multi-scale model. Therefore, we use the fractal dimension to redefine p.
Definition. Given a set of edge betweenness B = {b l 1 , b l 2 , . . ., b l n } and fractal dimension d B of network G, the proposed network vulnerability model is defined as: The larger the V d B , the more vulnerable the network. Any network that exhibits fractal properties can be applied to this proposed method to evaluate its vulnerability.
To demonstrate the usage of the proposed method, we apply it to two synthetic networks: Erdős-Rényi(ER) random networks [48] and Barabási-Albert(BA) model of scale-free networks [49]. The vulnerability of a network should be calculated as follow steps: Step 1 calculate the fractal dimension d B of the network using box-covering algorithm [40,46], i.e. Eq 5. The results of the two networks are illustrated in Fig 3. Step 2 Calculate the average edge betweenness according to Eq 2, and normalized by N ðN À1Þ 2 .
Step 3 Calculate the vulnerability V d B in accordance with Eq 9. Table 1 shows the results: V d B (BA)>V d B (ER), which mean that the ER network is more robust than the BA network. One point should be noted is that, the two networks are synthesized randomly, the results will vary with different network structure.

Comparison and Discussion
In this section, to testify the correctness and effectiveness of the proposed method, three commonly used methods presented in section 2, that is, the average inverse geodesic length l −1 , the largest component size LCS and the normalized average edge betweenness b nor (G), are applied to some real situations to compare with the proposed method. In order to get vulnerability embodying dynamic characteristics of networks to make a more reasonable comparison, we apply the RB attack strategy [6] to networks when using the three methods reflecting just static  topological properties. RB attack strategy means that one should remove the node with highest betweenness value and recalculate the betweenness of the network. Similarly, removing other nodes until all required number of nodes have been removed. Finally, we can get the vulnerability of the network based on the remaining network. In this paper, l −1 , LCS and b nor (G) are computed on attack strategy of removing 1% nodes from the original networks. In this paper, six unweighed US airline networks, categorized by year, are used to do analysis. That is, US airline network in 2005, 2007, 2010, 2013. These data are downloaded from the Bureau of Transportation Statistics (BTS) Transtats site with the following filters [50]: Geography = all; Year = 2011; Months = all; and columns: Passengers, Origin, Dest. Based on this table, the airport codes are converted into id numbers; if there are flights between two airports, an edge between two nodes will be connected correspondingly. Also ties with a weight of 0 are removed (only cargo), self-loops and small subgraph are removed, constructing connected and unweighed networks.
Firstly, l −1 , LCS and b nor (G) are applied to these six airline network to calculate vulnerability, and the results are illustrated in Table 2. Then the fractal dimension of these networks are calculated and illustrated in Fig 4. It's absolute that all these networks exhibit fractal properties, which reaches the prerequisite of the proposed method to use fractal dimension to characterize network vulnerability. Table 2 and Table 3 show the evaluation results of all these method, which are concluded as follows: 1. The average edge betweenness of these networks are very small, while the evaluation results of our method is larger than b 1 at magnitude level, which indicates that the fractal dimension is appropriate to be the p origin moment of edge betweenness to highlight the difference between networks effectively.
2. The vulnerability order given by other three commonly used methods are unreasonable to some extent. l −1 and LCS show that UAN2013 is the most vulnerable network and b nor (G) illustrates a high robustness of UAN2003 and UAN2005.
3. Our method give the most reasonable results among these methods, and the vulnerability order obtained by it is: UAN2003 > UAN2005 > UAN2007 > UAN2009 > UAN2011 > UAN2013. The robustness of these airline networks increases over years, which is consistent with the real situations in the country.
The obvious advantage of the proposed method is due to its multiple evaluation index, that is, edge betweenness, p origin moment and fractal dimension, which are all crucial to the vulnerability evaluation as described in section 3. While the evaluation index of other three methods is obvious inferior to ours, in detail, f l À1 just utilizes the geodesic length of network, g LCS Table 2. General characteristics of these networks and results illustrations. The proposed vulnerability evaluation method V d B is calculated based on fractal dimension d B . The normalized average inverse geodesic length f l À1 , normalized largest component size g LCS and the normalized average edge betweenness b nor (G) are computed after 1% of vertices are removed and are normalized by the corresponding values of the initial networks.  just considers the simple structure properties-connectivity, and b nor (G) only involves the average betweeness.
In conclusion, p origin moment is an effective approach to improve the average process of edge betweenness and the application of fractal dimension can highlight the difference between networks efficiently. So the proposed method has its geometrical and physical significance, and the reasonable application results also illustrate its effectiveness and practicality.

Conclusions
In this paper, the fractal dimension is redefined as p origin moment of edge betweenness to improve the vulnerability evaluation of multi-scale model. The experiments indicate that the fractal dimension indeed has relationships with network vulnerability, because a higher spacefilling capacity means that less boxes will be required to cover the whole network, thus with smaller fractal dimension; while a network with higher space-filling capacity will contribute to more connected and robust structure. The applications to several real airline networks and comparison with other commonly used methods also illustrate the effectiveness and practicability of using the fractal dimension to evaluate network vulnerability. One of our ongoing works is to explore the specific relationship between vulnerability and topological properties, such as connectivity, cluster degrees and so on. doi:10.1371/journal.pone.0146896.t003