Critical Infrastructure Vulnerability to Spatially ... - Wiley Online Library

0 downloads 0 Views 1MB Size Report
KEY WORDS: Chinese railway system; critical infrastructures; spatially localized failures; vulnerability. 1. ..... can be reached from the ith station within the next.
Risk Analysis

DOI: 10.1111/risa.12708

Critical Infrastructure Vulnerability to Spatially Localized Failures with Applications to Chinese Railway System Min Ouyang,1,2 Hui Tian,1 Zhenghua Wang,3 Liu Hong,1,2,∗ and Zijun Mao4

This article studies a general type of initiating events in critical infrastructures, called spatially localized failures (SLFs), which are defined as the failure of a set of infrastructure components distributed in a spatially localized area due to damage sustained, while other components outside the area do not directly fail. These failures can be regarded as a special type of intentional attack, such as bomb or explosive assault, or a generalized modeling of the impact of localized natural hazards on large-scale systems. This article introduces three SLFs models: node centered SLFs, district-based SLFs, and circle-shaped SLFs, and proposes a SLFs-induced vulnerability analysis method from three aspects: identification of critical locations, comparisons of infrastructure vulnerability to random failures, topologically localized failures and SLFs, and quantification of infrastructure information value. The proposed SLFs-induced vulnerability analysis method is finally applied to the Chinese railway system and can be also easily adapted to analyze other critical infrastructures for valuable protection suggestions. KEY WORDS: Chinese railway system; critical infrastructures; spatially localized failures; vulnerability

1. INTRODUCTION

as the life of its citizens. However, these systems are also subject to many types of hazards, such as unavoidable natural hazards, component aging, sharp demand increase, and climate change. Better protection of these spatially distributed critical infrastructures requires modeling their component fragilities under different types of hazards and then analyzing their system-level vulnerability. Here, the term “vulnerability” has many definitions,(2) without a broadly accepted one. This article interprets that “vulnerability” is associated with a specific initiating event and the vulnerability of a system to a specific initiating event is quantified as its performance drop. Hence, defining the initiating event is the first step for vulnerability analysis. According to the type of the initiating event, infrastructure vulnerability studies in the literature can be grouped into three types. The first type is infrastructure vulnerability analysis under random failures, which is for a series of initiating events, such as equipment failures, vegetation, animals, human errors, and so on, with large variety and uncertainty. These failures can be

Critical infrastructures are defined as systems or assets, whether physical or virtual, so vital to a nation that the incapacity or destruction of such systems and assets would have a debilitating impact on security, national economic security, national public health or safety, or any combination of those matters.(1) These infrastructures usually include electric power, water, telecommunication, transportation systems, and so on, and are playing an important role in supporting the economy of a nation as well 1 School

of Automation, Huazhong University of Science and Technology, Wuhan, China. 2 Key Lab. for Image Processing and Intelligent Control, Huazhong University of Science and Technology, Wuhan, China. 3 Department of Civil and Environmental Engineering, Rice University, TX, USA. 4 College of Public Administration, Huazhong University of Science and Technology, Wuhan, China. ∗ Address correspondence to Liu Hong, School of Automation, Huazhong University of Science and Technology, 1037 Luoyu Road, Wuhan 430074, China; [email protected].

1

C 2017 Society for Risk Analysis 0272-4332/17/0100-0001$22.00/1 

2 modeled by randomly removing a certain fraction or number of infrastructure components;(3–6) or by assigning a failure probability to each component and then comparing this probability to an uniformly distributed random number within [0, 1] to judge the component state (failed or normal);(7–9) or by first selecting the number of failed components according to a given distribution and then randomly removing the selected number of components.(10–12) Note that all the above studies on modeling random failures do not consider component geographical locations. The second type is infrastructure vulnerability analysis under intentional attacks. These attacks usually cause the failures of important components that are determined according to component degree,(13–15) betweenness,(16–18) load levels,(13) or are identified by using optimization techniques to search a set of components whose failures can cause the largest system-level vulnerability.(19–21) In addition to these studies on intentional attacks, which still do not consider component geographical locations, there also exist some other types of intentional attack, such as bomb or explosive assault, which only cause the failures of infrastructure components within a localized area. Berezin et al.(22) and Nicholson et al.(23) modeled the localized area by a circleshaped area with a random attack center in the network region and then analyzed the network vulnerability, but these modeling approaches cannot identify the critical location so that the failures of components inside such a critical location lead to the largest vulnerability. Patterson and Apostolakis(24) modeled the localized areas by dividing the infrastructure map into a generic hexagonal grid with a small radius and each hexagon is a localized area. Johansson and Hassel(25) made a similar analysis by dividing the infrastructure region into a square grid and each square is a localized area. These modeling approaches enable identifying the critical location, but the results depend on how the infrastructure map is partitioned into small localized areas (square, hexagon, etc.). The third type is infrastructure vulnerability analysis under natural hazards, such as earthquakes, hurricanes, floods, and lightening. These hazards can cause infrastructure components distributed within an influence area to fail simultaneously. The impact of these natural hazards on infrastructure components is usually modeled according to fragility curves, which provide the probability of exceeding a certain damage state threshold conditional on a selected hazard intensity measure, such as peak ground

Ouyang et al. acceleration or peak ground velocity or permanent ground deformation for seismic hazards,(26–30) and three-second gust wind speed for hurricane hazards.(31,32) In the literature, natural-hazardsinduced vulnerability studies are usually aimed at relatively small-scale systems, with all their components located in the hazard influence area and subjected to possible damage. Note that the relatively small-scale systems are usually parts of much larger or nationalscale systems, and then the impact of localized natural hazards on large-scale systems can be regarded as localized failures, such as vulnerability analysis of the Chinese railway system under flood hazards, where the flood hazards only occur in some provinces and their induced component failures are local relative to the system scale,(33) and seismic vulnerability analysis of European gas and electricity networks, where seismic hazards only occur in earthquakeprone regions.(34) However, these studies are usually made for specific cases, and a general study of localized failures and their induced vulnerability has been seldom addressed in the field of infrastructure engineering. According to the above literature review, this article will introduce a general type of initiating events in critical infrastructures, called spatially localized failures (SLFs), which are defined as the failure of a set of infrastructure components distributed in a spatially localized area due to damage sustained, while other components outside the area do not directly fail. Based on this definition, SLFs can be regarded as a special type of intentional attack, such as bomb explosion with mass destruction, or a generalized modeling of the impact of localized natural hazards on large-scale systems. Note that SLFs are only initiating events in critical infrastructure systems; how the systems perform under this event, including possible cascading failures, depends on the system operation mechanisms and the spatial correlation or heterogeneity of system parameters, such as supply levels for source nodes and demand levels for load nodes. In addition, different from previous studies modeling the localized areas by regular shapes, such as circle, square, or hexagon, and then identifying the critical components or locations based on the partition of an infrastructure map,(24,25) this article proposes a SLFs-induced vulnerability analysis method that can exactly identify the critical locations without partitioning the infrastructure map. Also, in the literature there exists another type of localized failure, called topologically localized failures (TLFs), which are defined as the failure of a

Critical Infrastructure Vulnerability to Spatially Localized Failures set of the nearest neighboring nodes from a given attack center node. Different from SLFs, the failed nodes in TLFs are the neighboring nodes of an attack center node and could be far from each other, so these failed nodes in TLFs are not necessary located within a spatially localized area that only covers those failed nodes. As studied by Shao et al.,(35) for regular networks (each node has the same number of neighbors) and a peer-to-peer computer network, in the case of only considering system topologies and taking the largest component size as the performance metric, they are less vulnerable to TLFs than random failures. However, these purely topologybased studies usually produce vulnerability results with weak correlations to those results obtained from the case of considering infrastructure flow properties, which is partly because topology-based analysis and flow-based analysis are usually characterized by different performance metrics to express the consequences.(36–38) When the flow properties and the flow-based performance model are considered, whether the above result still holds true and whether the system is further less vulnerable to SLFs than TLFs will be also studied in this article. The rest of the article is organized as follows: Section 2 proposes several SLFs models and a SLFs-induced vulnerability analysis method. Section 3 applies the proposed method to the Chinese railway system. Section 4 discusses the findings and extensions and Section 5 provides conclusions. 2. SPATIALLY LOCALIZED FAILURES AND THEIR INDUCED VULNERABILITY ANALYSIS METHOD To analyze infrastructure vulnerability, there exist many approaches in the literature, such as empirical approaches, agent-based approaches, systemdynamics-based approaches, economic-theory-based approaches, network-based approaches, and others. A detailed review of these approaches is provided in Ref. 39. Modeling SLFs in critical infrastructures requires system topological information, so a network-based approach is used. This article will take the Chinese railway system (CRS) as an example for illustrative purpose. The railway stations are described as nodes and the railway segments between stations are represented by links; thus the physical layout of the CRS is a network. In the network, there are trains departing from and arriving at some stations. The CRS network can be denoted by G = (N, E, F), where N is the set of nodes, E is

3

the set of links, and F represents the train flow information, recording each train ID and its route. Based on this network description, this section will first propose three SLFs models, then briefly introduce a vulnerability assessment model, and finally provide a SLFs-induced vulnerability analysis method. 2.1. Modeling Spatially Localized Failures According to the definition of SLFs, modeling the SLFs needs to first define a spatially localized area and then describe how components within this area fail. Previous studies usually describe a spatially localized area in a two-dimensional plane.(22–25) However, the infrastructure systems are embedded on the earth surface, which is actually the surface of a sphere in a three-dimensional space. Hence, this article describes the localized area in the infrastructure region as a part of the sphere surface in the three-dimensional space. This section mainly introduces three SLFs models, and how to simulate these SLFs and estimate their induced system-level vulnerability will be further introduced in Section 2.3. The first is node-centered spatially localized failures (NCSLFs), where a fraction pnc of nodes (including the attack center node itself) that have the nearest spherical distance from a given attack center node are selected as failed, and the failure area can be regarded as any area encircling those nodes. For the nodes inside the failure area, this article considers a special case that all nodes within the area are totally damaged, which may be unrealistic in some cases but enables at least an upper-bound estimation of the vulnerability. Based on these assumptions, the parameter pnc can be regarded as the failure intensity in NCSLFs. Note that this model is considered mainly for comparisons of SLFs-induced vulnerability to random failures and topologically localized failures, where the latter are also based on a randomly selected attack center node but the failed nodes are the nearest neighboring nodes from the attack center node.(35) The second is district-based spatially localized failures (DBSLFs), where a localized area is determined based on administrative districts, such as country boundaries, province boundaries, or subcompany affiliated areas. For each infrastructure component (node or link) in a district, its failure probability usually depends on the intensity of the hazard parameter used to quantify the magnitude of the hazardous event. This article considers a special case to make DBSLFs comparable with the other two

4 SLFs models: each component has the same failure probability, pdb , that is regarded as the failure intensity in DBSLFs. Note that this model is considered mainly for approximately simulating the impact of localized natural hazards on large-scale systems, like the vulnerability analysis of the Chinese railway system under flood hazards where flood events only occur in some provinces and only railway components in those provinces would fail.(33) The third is circle-shaped spatially localized failures (CSSLFs), where a localized area has a fixed circle shape determined by a center and a radius r. For the components (nodes or links or both) within or passing through the circle area, they are assumed totally damaged. The radius r in this model is regarded as the failure intensity in CSSLFs. This model is considered to simulate the impact of intentional attacks, such as bomb or military attacks with mass destruction. Note that the NCSLFs is not the special case of CSSLFs when the attack centers are only located at the nodes. Different from the NCSLFs, the localized area in CSSLFS is shaped by a circle, but the localized area in NCSLFs is any shape encircling the failed nodes. Also, given a value of failure intensity (pnc for NCSLFs, and r for CSSLFs), the number of failed nodes is fixed for NCSLFs, but variable for CSSLFs, depending the location of the attack center. In addition, the above NCSLFs and CSSLFs models mainly simulate a bomb or explosive assault; in each of these attacks there usually exists an attack center and the attack influence area depends on the direct distance from the attack center. Thus, only the geographical location of each node is used to compute the distances between nodes for defining the localized areas, which is the same with previous modeling approaches on SLFs.(22–25) 2.2. SLFs-Induced Vulnerability Assessment This article interprets that vulnerability is associated with a specific initiating event and the system vulnerability to a specific event is quantified as its performance loss. Given a specific damage scenario generated according to the SLFs models introduced in Section 2.1, how an infrastructure system performs depends on its operation mechanism and the spatial correlation or heterogeneity of some system parameters. For some types of infrastructure systems, such as electric power systems, there also exist cascading failures, which make some initially normal components also subject to damage. In addition, infrastructure systems are usually characterized

Ouyang et al. by different performance metrics, even for the same type of infrastructure system.(40) This article takes the Chinese railway system as an example for illustrative purposes, and this section will briefly introduce a performance metric and two performance models for railway systems, which have been studied in the authors’ previous work.(20,38) The CRS operates the trains according to a daily timetable; thus the number of trains running in the network is time dependent at different hours of a weekday, but the daily accessibility between stations within a typical weekday, which is quantified as the average fraction of reachable stations from each station by taking one or several trains during the day, keeps constant in each weekday. For illustrative purposes, this article will take this daily accessibility as the performance metric. Also, note that this daily accessibility immediately after the disruptive event is affected by the duration of the event. For simplification, this article assumes that all components’ damage states last longer than 24 hours; thus, the postevent daily accessibility metric and the daily accessibility-based vulnerability metric will be both independent of the event duration. In addition, this model does not consider rescheduling and rerouting in the next day immediately after disruptions. Note that the network congestion may occur after a disruptive event, and modeling network congestion requires more detailed information, such as the event occurrence time and location, each train’s location at the occurrence of the event, and so on. This article simply assumes that all interrupted trains can find some nearby stations or emergency railway segments to stop so that they would not cause congestion and affect other trains. More details of this model are provided in the authors’ previous work.(20,38) Denote n as the number of railway stations, and then in a normal weekday passengers in each station can reach other n – 1 stations. Define nid as the number of railway stations (except the ith station) that can be reached from the ith station within the next 24 hours after an event d, then the daily accessibility in the next day after the event d is defined as follows: 1  nid . n n−1 n

PMd =

(1)

i=1

Note that in a normal day PMn = 1, then the railway vulnerability under the event d is calculated as follows:

Critical Infrastructure Vulnerability to Spatially Localized Failures PMn − PMd 1  nid . = 1 − V = PMn n n−1 n

d

e6

(2)

i=1

Based on these metrics, depending on how much information is used for analysis, two railway performance models are used to compute vulnerability to the event. The first railway performance model is a purely topological model (PTM). An event initially causes the direct failures of some components (stations or railway links or both), which change the railway topology and lead to disconnection of some stations from others. During the next 24 hours immediately after the event, the number of stations nid that are reachable from the ith station is the total number of nodes minus one in the postevent connected cluster that includes the station i, then topology-based vulnerability VTd is computed by Equation (2). The second railway performance model is a train flow model. An event initially causes some railway component failures, then if a train j has its route passing through one of these failed components, the jth train is out of operation. After the event, the new topology together with the remaining normal trains can construct a space-of-changes network,(41) where the nodes still represent the railway stations in the physical network and two nodes are connected by a link if there exists at least one train with its route passing through these two stations. Based on this space-of-changes network, the value nid for the ith station is the number of nodes minus one in the postevent connected cluster that includes the ith station, then the flow-based vulnerability VFd to the event d is computed by Equation (2). 2.3. SLFS-Induced Vulnerability Analysis Based on the SLFs models and their induced vulnerability assessment model, this section will introduce a SLFs-induced vulnerability analysis method from three aspects, including identification of critical locations, comparisons of vulnerability to different failure models, and quantification of infrastructure information value. 2.3.1. Identification of Critical Locations For each SLFs model, this section will introduce the procedure to search the most critical location. For the NCSLFs, given a failure intensity pnc , the

n1

n2 e7 e1 e3 e4 e2 n4 n3 e5

(a)

5 n

e6 n1 e

e1 n3

n2 e7 e3 e4 e5 n4

(b)

n5

e6 n1

n5

n2 e7 e1 e3 e4 e2 n4 n3 e5

(c)

Fig. 1. A simple network and three CSSLFs scenarios. The attack radius is the same in three cases. (a) An attack scenario that causes the failure of node n1. (b) An attack scenario that covers a maximal node set {n1, n2}. (c) An attack scenario that covers a maximal node set and has its circle boundary pass through two nodes n1 and n2.

critical location is determined by the critical attack center node whose failure together with the failures of its nearest nodes make the system the most vulnerable; for the DBSLFs, given a failure probability pdb , the averaged vulnerability under the failures of nodes in each district can be computed over many simulation runs until the 98% confidence interval for the expected vulnerability has a width less than 0.01, and then the critical location is the district that makes the system the most vulnerable. For the CSSLFs, given a circle radius r, the critical location is determined by the critical attack center. Note that the attack center in this model can be in any position of the infrastructure map; thus then the number of candidate attack centers is infinite and how to exactly identify the critical attack center is difficult. Because the number of infrastructure components is limited, as an alternative, the problem of identifying critical location can be converted to first search a set of infrastructure nodes (attack objectives) that can be covered by a circle with radius r and whose failures can make the system the most vulnerable among the failures of all different node sets (each of them can be covered by a circle with radius r). To better introduce the algorithm, this article first takes a simple example, shown in Fig. 1, for demonstration. The network in Fig. 1 has five nodes and seven edges, and the three circle-shaped attack scenarios have the same radius. In Fig. 1(a), the node n1 is covered by a circle with radius r, but this node is not the optimum attack objective because if moving the attack center into some position like in Fig. 1(b), one more node n2 is covered by a circle with radius r. This indicates the attack scenario in Fig. 1(b) can make the system more vulnerable than the scenario in Fig. 1(a), and the attack location in Fig. 1(b) is more critical than that in Fig. 1(a). The

6 next question is whether the attack center in Fig. 1(b) can further move to some position to cover more nodes besides n1 and n2. The answer is negative; thus these two nodes {n1, n2} can be called a maximal node set MNS, which is defined as a set of nodes that can be covered by a circle with radius r and any new node added will make the nodes not able to be covered by a circle with radius r. Given a network position and an attack radius r, if it is possible to identify all its maximal node sets MNSs, then among them, the MNS whose contained node failures make the system the most vulnerable is the optimum attack objective. Hence, the problem of identifying critical location in CSSLFs is further converted into how to identify all MNSs. Note that the nodes in each MNS can be covered by a circle with its center at any position in an area described by a series of inequalities (the spherical distance between the attack center and the location of each node of the attack objectives is less than or equal to r), and then the circle center can move into some position so that it can pass through at least two nodes in the MNS. For example, the attack center in Fig. 1(b) can move to the position shown in Fig. 1(c) that can make the two nodes n1 and n2 still covered and make its boundary pass through the two nodes n1 and n2. Based on this fact, identifying all maximal node sets can be realized according to the following procedures: Step 1: for any two infrastructure nodes, determine the attack center with radius r and its circle boundary passing through these two nodes. The location of an infrastructure node is usually recorded by its latitude and longitude, which actually corresponds to a three-dimensional coordinate system. To make the attack-center-related equations easier to be solved, this article uses a three-dimensional coordinate system for calculation. For any node with longitude g and latitude t, its location can be converted into a point (x, y, z), with x = cos(t) × sin(g), y = cos(t) × cos(g), z = sin(t). For a network with n nodes, there are n(n – 1)/2 pairs of nodes, but searching the attack centers does not need to consider all these node pairs. If the spherical distance between two nodes is larger than 2r, these two nodes cannot be covered by a circle with radius r; thus this node pair is not considered for further analysis as the attack center does not exist. If spherical distance is not larger than 2r, denote the two nodes’ locations by (x1 , y1 , z1 ) and (x2 , y2 , z2 ), respectively, and the earth radius by R ࣈ 6,400 km, then

Ouyang et al. the attack center (x0 , y0 , z0 ) satisfies the following equations:

⎧ x02 + y02 + z02 = R2 ⎪ ⎪ ⎪ ⎪   r 2 ⎨ (x0 −x1 )2 +(y0 − y1 )2 +(z0 −z1 )2 = 2R sin .(3) 2R ⎪ ⎪   r 2 ⎪ ⎪ ⎩(x0 −x2 )2 +(y0 − y2 )2 +(z0 −z2 )2 = 2R sin 2R

The solution (x0 , y0 , z0 ) to Equations (3) can be obtained by the formulas in the Appendix. Step 2: for each of the attack centers computed in Step 1, compute the spherical distance between the attack center and all infrastructure nodes to identify the nodes that can be covered by a circle with radius r, and then these nodes are a candidate maximal node set. For a node (xi , yi , zi ), whether it can be covered by a circle can be judged by the following inequality: (xi − x0 )2 + (yi − y0 )2 + (zi − z0 )2   r 2 ≤ 2R sin . 2R

(4)

Step 3: refine all candidate maximal node sets so that each MNS is not the proper subset of any other maximal node set, then the refined candidate maximal node sets are just all maximal node sets. For each MNS, the system vulnerability under its contained nodes’ failures can be computed. The optimum MNS that makes the system the most vulnerable can be identified, and the critical location is an area that can cover all nodes in the optimum MNS. Note that the above algorithm only considers node failures in a circle area to identify critical locations. This assumption does not affect the results in most scenarios compared with the case that considers both node and link failures in the attack area because the failure of a node will also cause all its attached links’ failures. But for some special link whose endpoints are outside the attack area and a part of which is inside the area, this link should be also regarded as failed and needs to be additionally addressed, which will be further discussed in the Section 4. 2.3.2. Comparisons of Vulnerability to Different Failure Models This section introduces the vulnerability comparison method under different failure models, including random failures (RFs) where a certain fraction prf of randomly selected nodes are failed,(3,5) topologically localized failures (TLFs) where a certain fraction ptl

Critical Infrastructure Vulnerability to Spatially Localized Failures of nodes that are the nearest neighboring nodes of a given node are failed,(35) and SLFs. This comparison addresses two questions. One is whether the system is less vulnerable to SLFs than TLFs or RFs. As studied by Shao et al.,(35) for regular networks (each node has the same number of neighbors) and a peer-to-peer computer network, in the case of only considering system topologies and taking the largest component size as the performance metric, the networks are less vulnerable to TLFs than RFs. When the flow properties and the flow-based performance model are considered, whether the above result still holds true and whether the infrastructure system is less vulnerable to SLFs than both TLFs and RFs will be studied. The other question is whether the system is less vulnerable to all SLFs scenarios than RFs with the same failure intensity. To answer these questions and make different failure models comparable, only node failures are considered. For the RFs, the failure intensity prf is set to range from 0 to 1 with a step of 1/n (n is the total number of nodes), and for each failure intensity, a certain number of failure scenarios are generated through randomly selecting prf n nodes from the network (similar to the n-k scenarios in the literature, where k is the number of failed components) to compute system-level vulnerability until the 98% confidence interval for the expected vulnerability has a width less than 0.01. The averaged vulnerability is computed and denoted by E[VT ,RFs (prf )] for topology-based vulnerability and E[VF ,RFs (prf )] for flow-based vulnerability. For the TLFs, each node together with its ptl n nearest neighboring nodes are selected and removed from the network to compute the system-level vulnerability. To make SLFs comparable with the above two failure models, only the NCSLFs model is considered to address the first question. The NCSLFs are simulated by removing each of n nodes and its pnc n – 1 nearest nodes in terms of the spherical distance, then the average vulnerability among all n cases can be computed and denoted by E[VT,NCSLFs (pnc )] for topologybased vulnerability and by E[VF,NCSLFs (pnc )] for flow-based vulnerability. According to the above procedures, the RFs-, TLFs-, and NCSLFs-induced vulnerability can be compared for different failure intensities. For the second question, the comparison analysis will be made for all three SLFs models. For the NCSLFs, comparing the vulnerability to each node i centered failures with E[VRFs (prf )] can check

7

whether the system is more vulnerable to this node centered failure than to random failures. For the DBSLFs, denote the number of all nodes in the ith district by di , and comparing the vulnerability to the failures of all these nodes with E[VRFs (di /n)] can check whether the system is more vulnerable to this district failure than to random failures. For the CSSLFs, given a circle radius r, for each maximal node set i with a size si , comparing the vulnerability to the failures of all nodes in this MNS with E[VRFs (si /n)] can check whether the system is more vulnerable to this MNS failure than to random failures. 2.3.3. Quantification of Infrastructure Information Value This article defines a piece of infrastructure information as the differences of attack consequences when it is exposed to the attackers or not. Based on this definition, the attack consequences are measured from the infrastructure company’s point of view and then a flow-based vulnerability metric mainly concerned by the company is used, such as the flow-based accessibility vulnerability metric for the Chinese railway system. This article assumes that attackers will use the bomb assault, which is modeled by CSSLFs. In addition, the attackers’ strategies are modeled according to various rules, depending on different levels of infrastructure information they have (1) When the attackers have no information, the attack centers could be located anywhere on the infrastructure map. (2) When the attackers only have infrastructure nodes’ geographical locations, the attack center is assumed to cover as many nodes as possible and then the maximal node set with the largest size is the attack objective. (3) When the attackers further have the infrastructure topological information, the attack center is assumed to cause the largest topology-based vulnerability and then the maximal node set whose failures can result in the largest topology-based vulnerability will be the attack objective. (4) When the attackers further have the infrastructure flow information, the attack center is assumed to cause the largest flow-based vulnerability and then the maximal node set whose failures can result in the largest flow-based vulnerability is the attack objective. Based on the above assumptions, this article considers three types of infrastructure information: node geographical information (NGI), system topological information in the case of knowing NGI (STI | NGI), and system flow information in the case of knowing

8

Ouyang et al.

Fig. 2. A geographical representation of the Chinese railway system. The railway stations’ located provinces and affiliated railway bureaus (corporations) are also shown on the map. The number behind a province name or a bureau index letter is the province or the bureau ID.

NGI and STI (SFI | NGI, STI). Their information values INGI , ISTI | NGI , ISFI | NGI,STI can be respectively quantified as follows: RAC LMNS INGI (r ) = VF,CSSLFs (r ) − E VF,CSSLFs (r ) , (5) LT BV LMNS IST I|NGI (r ) = VF,CSSLFs (r ) − VF,CSSLFs (r ),

(6)

LF BV LT BV ISF I|NGI,ST I (r ) = VF,CSSLFs (r ) − VF,CSSLFs (r ), (7) RAC LMCS LT BV where E[VF,CSSLFs (r )], VF,CSSLFs (r ), VF,CSSLFs (r ), LF BV and VF,CSSLFs (r ) are, respectively, the average flowbased vulnerability under random attack centers (RAC), the flow-based vulnerability under the failures of all nodes in the largest maximal node set (LMNS), the flow-based vulnerability under the failures of all nodes in the maximal node set

that can result in the largest topology-based vulnerability (LTBV), and the flow-based vulnerability under the failures of all nodes in the maximal node set that can result in the largest flowbased vulnerability (LFBV). In addition, a certain number of random attack centers in the infrastructure map are generated by the ARCGIS software to make the 98% confidence interval RAC (r )] have a width less than 0.01 for difof E[VF,CSSLFs ferent r. 3. APPLICATIONS AND RESULTS The proposed approach is illustrated with a Chinese railway system example in this section. This system plays a crucial role in the economy of China and supports the well-being of its citizens. In 2012,

Critical Infrastructure Vulnerability to Spatially Localized Failures it transported around 1.89 billion passengers and approximately 3.89 billion metric tons of cargo. This article mainly focuses on passenger transport for vulnerability analysis. According to the recent timetable “Chinese Railway Passenger Train Timetable” published in 2010,(42) there are more than 2,000 passenger stations in total, among them only stations in the cities at prefecture level and above in China as well as some other stations that are railway hubs or have trains departing from or arriving at are selected, while multiple stations in a city are combined into one station, and finally the coarse-grained railway system has n = 399 stations connected together by e = 500 railway links. These stations and railway links are separately managed by 18 railway bureaus (or corporations) and are distributed in 32 provinces (autonomous regions and municipalities), as shown in Fig. 2. Based on this coarse-grained railway data, the Shenyang railway bureau is responsible for the largest number of stations (57 stations), while the Heilongjiang province has the largest number of stations (33 stations). In the rail layout, on a typical weekday there are 4,196 trains transporting passengers between cities. The initial departure and final arrival stations as well as the detailed routes of each train are obtained from the timetable. For this system, the SLFs-induced vulnerability analysis is performed as follows. First, the critical locations are analyzed for different SLFs models when the flow-based performance model is used. For the NCSLFs, when the failure intensity pnc is varied from 1/399 to 1.0 with a step of 1/399, Fig. 3(a) shows the corresponding ID(s) of the critical center node(s). When pnc ࣙ 92.7%, there exists a sharp increase of the number of critical center nodes, which all can cause the largest vulnerability, but when pnc < 92.7%, there are only one or two critical center nodes, and it can be found that some critical nodes repeat at different pnc . Fig. 3(b) uses the star-shaped nodes to emphasize the nodes that are identified as critical for at least one value of pnc < 92.7% and uses the star size to display how frequently these critical nodes repeat in different pnc . Note that the tree-like structure on the left of Fig. 3(b) makes the network vulnerable, which is validated by the proposed method as some root nodes in this tree-like structure are emphasized as critical. But these root nodes (with small star sizes) are only critical for one or two values of pnc , while the critical center nodes with large star size (repeated as critical a lot) are usually centralized in the north of

9

China. This is partly due to the railway topology and if the attack center nodes are those stations, both the upper- left part and the upper-right part of the system will be disconnected from the remaining network, making the system more vulnerable than other node centered failures. The train flow is another reason for the locations of critical center nodes; the Linyi Station in Shandong province, as shown in the Fig. 3(b), repeats as critical the most (18 times), this station is not special from the topological perspective, but when checking the train flow, there are a large number of trains passing it. For the DBSLFs, this article divides the districts according to the provinces or the railway bureaus. When the failure probability pdb in a district ranges from 0.05 to 1.0 with a step of 0.05, Fig. 4(a) shows the critical province(s) and bureau ID (s). When the districts are divided according to the railway bureaus, the critical bureau is always the Shenyang railway bureau, with its location shown in Fig. 4(b). Note that node failures in this bureau can disconnect many railway stations in the upper-right-most side of the system and then result in relatively larger vulnerability. When the districts are divided according to the provinces, the critical province is Gansu province if pdb ࣘ 0.25, and Heilongjiang province if 0.25 < pdb ࣘ 0.425, and Jilin province when pdb > 0.425, with their locations also shown in Fig. 4(b). Similar to the results in NCSLFs, these districts are critical partly due to their topological positions and partly due to the train flow passing through them. For the CSSLFs, Fig. 5(a) shows the largest and average vulnerability as a function of circle radius r ranging from 10 to 500 km with a step of 10 km when the nodes in each of the maximal node sets are subject to attack. The inset figure shows the ratio of the largest and average vulnerability to r to simply reflect the attack efficiency. From the figure, with the increase of r, the average vulnerability increases smoothly, but the largest vulnerability has a sharp increase when r is between 100 km (0.143) and 140 km (0.409), and it shows a peak of the attack efficiency at r = 140 km. Fig. 5(b) shows the optimum attack centers under different r, and the optimum attack areas for r = 90, 170, and 500 km are also shown on the map. When r ࣘ 90 km, the optimum attack centers are almost all located around the center of the critical circle for 90 km; when 90 km < r ࣘ 170 km, the circle centers are around the center of the critical circle for the 170 km; when r keeps increasing, the circle centers are still in the north part of the system. This is because the failures of stations in these critical

10

Ouyang et al.

Fig. 3. (a) Critical center node ID(s) under different failure intensities in the NCSLFs. (b) Display of critical center nodes in the railway map when the failure intensity is less than 0.927. The star-shaped nodes are critical nodes, and the star size is proportional to the repeat frequency of the critical nodes under different failure intensities.

13 12

Critical District ID

11 10 9 8 7 6 5 4

Railway Bureau Province Province 0

0.2 0.4 0.6 0.8 Failure Probability of Nodes in Each District

1

(a)

(b)

Fig. 4. (a) Critical district (railway bureau or province) IDs as a function of the failure probability of nodes in each district; the districts are divided according to railway bureaus and provinces. (b) Display of critical districts and their related nodes in the railway map.

areas are able to disconnect the railway topology into separated parts, which indicates the critical locations are largely affected by infrastructure topology. Second, the flow-based vulnerability results under random failures (RFs), topologically localized failures (TLFs), and NCSLFs are shown in Fig. 6(a). The topology-based vulnerability curves show similar trends and are not shown in this article. From the figure, when the flow properties are considered, similar to the results found by Shao et al.,(35) the system is less vulnerable to TLFs than RFs, but the system is more vulnerable to TLFs than SLFs. The second question in Section 2.3.2 is whether the system is less

vulnerable to all SLFs scenarios than random failures with the same failure intensity. The results are shown in Figs. 6(b)–6(d) for the three SLFs models. Fig. 6(b) shows the fraction of NCSLFs scenarios more harmful (making the system more vulnerable) than random failures under different fractions pnc of failed nodes when two railway models are considered, respectively. For both models, it can be found that when pnc is small or close to 1.0, there exist many NCSLFs scenarios more harmful than random failures. When 0.2 < pnc ࣘ 0.9, none of the NCSLFs scenarios make the system more vulnerable than random failures.

Critical Infrastructure Vulnerability to Spatially Localized Failures

11

Fig. 5. (a) The largest and average vulnerability as a function of the circle radius in CSSLFs. The inset figure shows the ratio of the largest and average vulnerability to the attack radius to reflect the attack efficiency. (b) Display of critical attack centers and several critical attack areas.

Fig. 6(c) shows the topology-based and flowbased vulnerability under DBSLFs when the districts are divided according to railway bureaus and provinces with all nodes failed in each district (pdb = 1). Similar results are also found for other values of pdb . It can be found that when the topologybased performance model is used, there exists one railway bureau (Lanzhou railway bureau, 18 railway stations, VT ,DBSLFs = 0.194 and E[VT,RFs (18/399)] = 0.157) and two provinces (Gansu province, 13 railway stations, VT,DBSLFs = 0.171 and E[VT,RFs (13/399)] = 0.116; Jilin province, 25 railway stations, VT, DBSLFs = 0.298 and E[VT,RFs (25/399)] = 0.218) that the failures of all stations inside each of these districts make the system more vulnerable than random failures with the same failure intensity (the number of randomly selected nodes to fail equals to the number of nodes in the corresponding district); when the flowbased performance model is used, there is only one province (Gansu province, 13 railway stations, VF, DBSLFs = 0.180 and E[VF,RFs (13/399)] = 0.175) that the failures of all stations inside this province are more harmful than random failures with the same failure intensity. Also, it can be found that for the DBSLFs, the failures of all stations in a district with a larger number of nodes do not necessarily make the system more vulnerable, specially for the Jilin province; it has 25 railway stations, less than the number of stations in Heilongjiang province (33 stations), but the failures of all stations in Jilin province result in both the largest topology-based and the

flow-based vulnerability (VT ,DBSLFs = VF ,DBSLFs = 0.298), which are much larger than those for the Heilongjiang province (VT ,DBSLFs = 0.191, VF ,DBSLFs = 0.195). Fig. 6(d) shows the number of maximal node sets in CSSLFs as a function of r and also displays the fraction of MNSs failures more harmful than random failures with the same intensity. Similar trends are found for two performance models. From the figures, with the increase of the attack radius, the number of MNSs increases almost linearly, but the fraction of MNSs failures more harmful than random failures fluctuates; particularly when r is around 170 km, there exists an increase for both curves. This is because at r around 170 km, there exist many MNSs located near the critical circle for 90 and 170 km in Fig. 5(b) and the failures of stations in these MNSs can easily disconnect the CRS topology and then make the system more vulnerable than random failures with the same intensity. When r becomes larger and larger, the MNSs located in the north part of the CRS can still disconnect the system topology, but the number of nodes inside some of these MNSs also increases and the random failures with a corresponding higher intensity make the CRS vulnerability increase faster than those MNSs failures, then less and less MNSs failures are more harmful than random failures. When r is 500 km, the fraction of MNSs failures more harmful than random failures is close to zero. The above results indicate that the CRS is less vulnerable to most of SLFs scenarios than random

12

Ouyang et al.

Fig. 6. (a) Flow-based accessibility vulnerability of the Chinese railway system under random failures (RFs), topologically localized failures (TLFs), node-centered spatially localized failures (NCSLFs). (b) Fraction of NCSLFs scenarios more harmful than random failures when two performance models are respectively considered. (c) Topology-based and flow-based accessibility vulnerability in DBSLFs when the districts are divided according to railway bureaus and provinces. (d) Number of maximal node sets in CSSLFs and the fraction of MNSs failures more harmful than random failures as a function of circle radius.

failure with the same intensity, but there still exist a few SLFs scenarios that should be paid attention to. Third, the information values for node geographical information (NGI), system topological information in the case of knowing NGI (STI | NGI), and system flow (train flow) information in the case of knowing NGI and STI (SFI | NGI,STI) are all quantified, with results shown in Fig. 7. From the figure, when an attacker only has the node geographical information NGI and the maximal node set with the largest size is the attack objective, LMCS (r ) the resulting flow-based vulnerability VF,CSSLFs is always strictly larger than the flow-based average RAC (r )] under random attack vulnerability E[VF,CSSLFs centers, and then the NGI information value esti-

mated according to Equation (5) is always larger than zero; when an attacker has the system topological information STI plus the NGI information and the maximal node set that can cause the largest topology -based vulnerability is the attack objective, the LT BV (r ) is resulting flow-based vulnerability VF,CSSLFs LMCS larger than VF,CSSLFs (r ) at most values of r, but at LT BV LMCS (r ) is less than VF,CSSLFs (r ), r = 160 km, VF,CSSLFs LT BV LMCS (r ), while at r = 170 km, VF,CSSLFs (r ) = VF,CSSLFs which makes the STI | NGI information have negative or zero value. This indicates that exposing more system information to attackers does not necessarily make the system more vulnerable; when an attacker has the train flow information plus the NGI and STI

Critical Infrastructure Vulnerability to Spatially Localized Failures

Fig. 7. Information value as a function of circle radius in the CSSLFs for three types of railway information: node geographical information (NGI), system topological information in the case of knowing NGI (STI | NGI), and railway system flow information in the case of knowing NGI and STI (SFI | NGI,STI).

information and the maximal node set whose failures can cause the largest flow-based vulnerability is the attack objective, the resulting flow-based vulneraLF BV LT BV (r ) is never less than VF,CSSLFs (r ), bility VF,CSSLFs LF BV LT BV (r ) but at many values of r, VF,CSSLFs (r ) = VF,CSSLFs and the SFI | NGI information have zero value. When comparing the values of different types of information, it can be found that if r ࣙ 190 km, the STI | NGI information is the most valuable, while the SFI |NGI, STI information is the least valuable with its value less than 0.04 and equal to zero at many r; if 150 km ࣘ r < 190 km, the NGI information is the most valuable. Also, the STI | NGI information is always more valuable than the SFI | NGI,STI information. These results indicate that to protect the CRS for circle-shaped spatially localized attacks, the railway topological information and station geographical information is much more valuable than the train flow information. 4. DISCUSSION The SLFs-induced vulnerability analysis on the Chinese railway system in Section 3 shows that (1) for different SLFs models, the critical locations are concentrated in the north part of the CRS map, which is partly due to the railway topology and partly due to the train flow; (2) the SLFs make the CRS less vulnerable than both random failures and topologically localized failures in most cases, but there still exist a few SLFs scenarios that are more harm-

13

ful than random failures, and these scenarios should be given more attention; (3) information about the railway topology and station geography is much more valuable than information about the train flow, with the railway topology being most valuable for the majority of attack radius. These results provide useful suggestions for Chinese railway protection. However, the above results are based on several assumptions, and this section will further discuss these assumptions and consider possible extensions. (1) For the flow-based railway performance model, the daily accessibility metric only quantifies the possibility of travel among cities, without considering the travel feasibility due to flow-related issues, such as overload-induced unavailability to change the travel vehicles or unavailability to get travel tickets. When more flow-based data are available, by assigning all network links with capacity data, the railway performance model can be adapted to make a more practical flow-based vulnerability analysis. In addition, the spatially localized attacks might cause large-scale panic issues, which would affect the need to travel, and then some flow-based performance metric, such as the number of transported customers, will be further affected. In this case, additional psychological analysis should be considered and combined into the framework to provide a more comprehensive vulnerability analysis. (2) This article only considers node failures in a spatially localized area. Despite that node failures can bring all their attached links’ failures and then almost all links passing through the failure area will also fail, there also exist some cases where two endpoints of a link are outside the failure area, while a part of the link is inside the area; in this case, this link failure should be also considered as failed. To see how the vulnerability analysis results will be affected in this case, this article further studies the railway vulnerability when all nodes inside a province and all links passing through the province are failed. The results show that among 32 provinces, when using the topology-based performance model, two provinces have increased topology-based vulnerability; when using the flow-based performance model, six provinces have increased flow-based vulnerability. Most of provinces have the vulnerability increased less than 30%, but there exists one province, Hebei province, whose both topology-based and flow-based vulnerability are four to five times the vulnerability under only node failures. This province then becomes the most critical province, and the most critical Jilin province under only node failures becomes the

14 second critical under both node and link failures. These results indicate that for spatially localized attacks that consider both node and link failures, some links with endpoints outside the attack area and link segments inside the area should be paid much more attention for critical location identification. However, in this case, as an infrastructure link is usually not straight line and cannot be expressed by a simple mathematical expression, checking whether some part of a link is covered by a failure area is difficult by purely mathematical method. As an alternative, it can be realized through dividing a link into a lot of short segments and labeling each segment by a point and then checking whether there exists one point covered by the failure area, or by manually checking in the infrastructure map based on the results under only node failures with the help of an ARCGISbased decision support system. Hence, the proposed approach for only node failures provides an initial step to identify critical locations in the case of considering both node and link failures, and enables at least a lower-bound estimation of the vulnerability. 5. CONCLUSIONS This article proposes several models of SLFs, and introduces a SLFs-induced vulnerability analysis method from three aspects: identification of critical locations, comparisons of vulnerability to random failures, topologically localized failures and SLFs, and quantification of infrastructure information value. Taking the Chinese railway system as an example for application study, some results valuable for its protection are found and discussed. This article only takes the Chinese railway system as an example for illustrative purposes, but the proposed SLFs models and their induced vulnerability analysis method can be easily adapted to analyze other critical infrastructures. Note that this article does not consider a fully multiobjective approach for determining a Pareto front that maximizes vulnerability and minimizes the number of nodes failures, which should be related to protection analysis against SLFs and is a direction for future research. Also, this article only considers a single system for analysis. When considering multiple systems with interdependencies, whether some different results can be found will be a direction for future research. Moreover, the information value quantification is based on a strong assumption that the attackers will use the strategies introduced in Section 2.3.3; relaxing those assumptions and modeling the interactions between

Ouyang et al. attacker and defender(21,43) will be more practical for the value quantification of different infrastructure information. Finally, this work is based on an assumption that the failure of a set of components will cause the vulnerability not smaller than the failure of any subset of these components, which holds true for the purely topological model and the train flow model studied in this article. But this assumption might not hold true for some other flow models, such as direct current power flow model. When a power system operates close to its capacity limit, there may exist Braess’s paradox so that more component failures can cause less vulnerability. In this case, the proposed algorithm will be invalid, and then using the concept of maximal component set to design a new algorithm is an interesting direction for future research. ACKNOWLEDGMENTS This material is based upon work supported in part by the National Science Foundation of China under Grants 71671074, 51208223, 61572212, 71303085, and the Fundamental Research Funds for the Central Universities under Grant 2014QN166. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the sponsors. APPENDIX The solution (x0 , y0 , z0 ) to Equations (3) is determined by the following formulas: x0 =

−(c1 d1 + e1 f1 ) ±

(c1 d1 + e1 f1 )2 − (1 + c12 + e12 )(d12 + f12 − R2 ) 1 + c12 + e12

,

y0 = c1 x0 + d1 , z0 = e1 x0 + f1 , c1 = d1 =

x1 z2 − x2 z1 , y2 z1 − y1 z2

(z1 − z2 ) R2 − 2R2 sin2

 r 2R ,

y2 z1 − y1 z2 −x1 − y1 c1 , z1  r −y1 d1 + R2 − 2R2 sin2 2R . f1 = z1

e1 =

REFERENCES 1. Title 42 US Code, Sec 5195c et seq. 2006 Supp IV, 2011. Critical Infrastructures Protection. Available at: http://www.gpo.gov/, Accessed December 24, 2015. 2. Aven T. On some recent definitions and analysis frameworks for risk, vulnerability and resilience. Risk Analysis, 2011; 31:515–522.

Critical Infrastructure Vulnerability to Spatially Localized Failures 3. Motter AE, Lai YC. Cascade-based attacks on complex networks. Physical Review E, 2002; 66:065102-–065102-4. 4. Crucitti P, Latora V, Marchiori M. Model for cascading failures in complex networks. Physical Review E, 2004; 69:0451041–045104-4. 5. Albert R, Albert I, Nakarado GL. Structural vulnerability of the North American power grid. Physical Review E, 2004; 69:025103-1–025103-4. 6. Yazdani A, Jeffrey P. Complex network analysis of water distribution systems. Chaos, 2011; 21:016111-1–016111-10. 7. Dobson I, Carreras BA, Lynch VE, Newman DE. Complex systems analysis of series of blackouts: Cascading failure, critical points and self-organization. Chaos, 2007; 17:026103-1– 026103-13. ˜ 8. Ouyang M, Duenas-Osorio L. Time-dependent resilience assessment and improvement of urban infrastructure systems. Chaos, 2012; 22:033122-1–033122-11. 9. Carreraa BA, Lynch VE, Dobson I, Newman DE. Complex dynamics of blackouts in power transmission systems. Chaos, 2004; 14:643–652. 10. Dobson I, Carreras BA, Newman DE. A branching process approximation to cascading load-dependent system failure. Thirty-Seventh Hawaii International Conference on System Sciences, Big Island, Hawaii, 2004. 11. Ren H, Dobson I. Using transmission line outage data to estimate cascading failure propagation in an electric power system. IEEE Transactions on Circuits and Systems-II: Express Brief, 2008; 55:927–931. ˜ 12. Ouyang M, Duenas-Osorio L, Min X. A three-stage resilience analysis framework for urban infrastructure systems. Structural Safety, 2012; 36–37:23–31. 13. Hines P, Cotilla-Sanchez E, Blumsack S. Do topological models provide good information about electricity infrastructure vulnerability? Chaos, 2010; 20:033122-1–033122-5. 14. Holmgren AJ. Using graph models to analyze the vulnerability of electric power networks. Risk Analysis, 2006; 26:955– 969. 15. Rosas-Casals M, Valverde S, Sole RV. Topological vulnerability of the European power grid under errors and attacks. International Journal of Bifurcation and Chaos, 2007; 17:2465– 2475. 16. Holme P, Kim BJ, Chang NY, Seung KH. Attack vulnerability of complex networks. Physical Review E, 2005; 65:056109. 17. Bompard E, Wu D, Xue F. The concept of betweenness in the analysis of power grid vulnerability. Complexity in Engineering, 2010; 5:52–54. 18. Mishkovski I, Biey M, Kocarev L. Vulnerability of complex networks. Communications in Nonlinear Science and Numerical Simulation, 2011; 16:341–349. 19. Zio E, Golea LR, Rocco CM. Identifying groups of critical edges in a realistic electrical network by multi-objective genetic algorithms. Reliability Engineering and System Safety, 2012; 99:172–177. 20. Ouyang M, Pan ZZ, Hong L, He Y. Vulnerability analysis of complementary transportation systems with applications to railway and airline systems in China. Reliability Engineering and System Safety, 2015; 142:248–257. 21. Alderson DL, Brown GG, Carlyle WM. Operational models of infrastructure resilience. Risk Analysis, 2015; 35:562–586. 22. Berezin Y, Bashan A, Danziger MM, Li D, Havlin S. Localized attacks on spatially embedded networks with dependencies. Scientific Reports, 2015; 5:8934-1–8934-5. 23. Nicholson C, Barker K, Ramirez-Marquez J. Flow-based vulnerability measures for network component importance: Experimentation with preparedness planning. Reliability Engineering and System Safety, 2016; 145:62–73. 24. Patterson SA, Apostolakis GE.Identification of critical locations across multiple infrastructures for terrorist actions. Reliability Engineering and System Safety, 2007; 92:1183–1203.

15

25. Johansson J, Hassel H. An approach for modeling interdependent infrastructures in the context of vulnerability analysis. Reliability Engineering and System Safety, 2010; 95(12):1335– 1344. 26. Duenas-Osorio L, Craig JI, Goodno BJ. Seismic response of critical interdependent networks. Earthquake Engineering and Structural Dynamics, 2007; 36:285–306. 27. Adachi T, Ellingwood BR. Serviceability of earthquakedamaged water systems: Effects of electrical power availability and power backup systems on system vulnerability. Reliability Engineering and System Safety, 2008; 93:78–88. 28. Esposito S, Iervolino L, d’Onofrio A, Santo A, Cavalieri F, Franchin P. Simulation-based seismic risk assessment of gas distribution networks. Computer-Aided Civil and Infrastructure Engineering, 2015; 30:508–523. 29. Franchin P, Cavalieri F. Probabilistic assessment of civil infrastructure resilience to earthquakes. Computer-Aided Civil and Infrastructure Engineering, 2015; 30:583–600. 30. Cavallaro M, Asprone D, Latora V, Manfredi G, Nicosia V. Assessment of urban ecosystem resilience through hybrid social–physical complex networks. Computer-Aided Civil and Infrastructure Engineering, 2014; 29:608–625. 31. Federal Emergency Management Agency. Hazards U.S. Multi-Hazard (HAZUS-MH) Assessment Tool v1.4. Available at: www.fema.gov/plan/prevent/hazus/index.shtm, 2015. ˜ 32. Ouyang M, Duenas-Osorio L. Multi-dimensional hurricane resilience assessment of electric power systems. Structural Safety, 2014; 48:15–24. 33. Hong L, Ouyang M, Peeta S, He XZ, Yan YZ. Vulnerability assessment and mitigation for the Chinese railway system under floods. Reliability Engineering and System Safety, 2015; 137:58–68. 34. Poljansek K, Bono F, Gutierrez E. Seismic risk assessment of interdependent critical infrastructure systems: The case of European gas and electricity networks. Earthquake Engineering and Structural Dynamics, 2012; 41:61–79. 35. Shao S, Huang X, Stanley HE, Havlin S. Percolation of localized attack on complex networks. New Journal of Physics, 2015; 17:023049-1–023049-5. 36. Cavalieri F, Franchin P. Models for seismic vulnerability analysis of power networks: Comparative assessment. ComputerAided Civil and Infrastructure Engineering, 2014; 29:590–607. 37. Ouyang M. Comparisons of purely topological model, betweenness based model and direct current power flow model to analyze power grid vulnerability. Chaos, 2013; 23: 023114-1–023114-11. 38. Ouyang M, Zhao LJ, Hong L, Pan ZZ. Comparisons of complex network based models and real train flow model to analyze Chinese railway vulnerability. Reliability Engineering and System Safety, 2014; 123:38–46. 39. Ouyang M. Review on modeling and simulation of interdependent critical infrastructure systems. Reliability Engineering and System Safety, 2014; 121:43–46. 40. Pitilakis K, Franchin P, Khazai B, Wenzel H. SYNER-G: Systemic Seismic Vulnerability and Risk Assessment of Complex Urban Utility, Lifeline Systems and Critical Facilities: Methodology and Applications, Vol. 31. Geotechnical, Geological and Earthquake Engineering, Springer Netherlands, 2014. 41. Kurant M, Thiran P. Extraction and analysis of traffic and topologies of transportation networks. Physical Review E, 2006; 74:036114-1–036114-10. 42. Shipping Office of Ministry of Railways. Chinese Railway Passenger Train Timetable. Beijing: China Railway Publishing House, 2010. 43. Murray-Tuite PM. A methodology for assessing transportation network terrorism risk with attacker and defender interactions. Computer-Aided Civil and Infrastructure Engineering, 2010; 25:396–410.