Extreme risk induced by communities in interdependent networks

Networks in nature not only depend on each other but also have internal community structures, such as infrastructure networks with links within and across geographic regions. The communities play an important role when the networks undergo localized failures in specific regions, for instance when natural disasters or economic sanctions disrupt a local community region and consequently influence the whole system. How a disruption in one community propagates through the entire system is a crucial, but still open, question. Here we find that the community structure embeds extreme risk: weakening the community strength could abruptly drive the system to a precarious state. Examining the business-flight network among cities as a proxy for the world economy, we find this real coupled system evolving towards the extreme vulnerable phase due to ongoing globalization. This shows the community risk indeed exists in real world networks and deserves more attention from the scientific community.The effect of localized disruption or failure in interdependent networks with internal community structure remains an open question. Adopting the generating function approach, the authors are able to uncover rich phase transition behaviours and the associated risks for such system, and by studying real networks under random failures within a community, they find that weakening the community strength could rapidly drive the system to a precarious state.

I nterdependent networks are important since many real work systems have interactions at different levels, and such complexity often leads to new and rich phase transition behaviors not present in single layer networks [1][2][3] . Although much research has focused on the robustness to failure in multilayer interdependent networks 1,[4][5][6][7][8][9][10][11][12] , these studies have assumed that networks are unstructured. This is in sharp contrast to real-world networks that have an internal community structure 6,[13][14][15][16] . For example, transportation networks have more connections within urban regions than between urban regions. A group of countries usually have more economic ties within the group than with countries outside the group. A problem of particular interest is determining how these complex systems with a rich community structure behave under localized disruption. Such community structure has been studied on single layer networks 17 to understand the impact of removing intercommunity links, and rich phase transition has been found when the number of communities changes. Random attack on intercommunity links on interdependent networks have also been found 18,19 to have even more complex transitions and scaling relations. Yet, certain attacks are localized in specific communities, like natural disasters and economic sanctions. The underlying mechanism controlling how such localized disruption in one community of a complex system disseminates throughout the entire system has not been understood, as well as the conditions that the network becomes vulnerable to abrupt collapse. Analytically, it is challenging to study since the network system with both interdependency and communities structure is highly complex.
Here, we present a generalized framework of interdependent multilayer network with community structure 1,6 based on generating functions, and study cascading processes that occur across the entire network initiated by random damage disruption in a single community. Such framework can be generalized to any layers of interdependent networks with arbitrary number of communities. Our analysis reveals that such system has rich phase transition behaviors that are much more complicated than interdependent networks without community structures. In particular, we find both theoretically and empirically that network robustness changes abruptly from safe to vulnerable as the strength of the community changes. The safe region of the system is characterized by the lack of phase transition phenomenon, i.e. the system does not disintegrate even with one whole community fully removed from the system. The vulnerable region is when phase transition is present, i.e. removing enough fraction of one community will disintegrate the whole system. Employing the business-flight network among cities of North America, Asia, and Europe as an example with strong community structure, we confirmed the presence of such risk in this interdependent system. More strikingly, this system is evolving towards the transition point from the safe to the vulnerable phase due to ongoing globalization.

Results
Theoretical analysis. We consider an interdependent multilayer network with a community structure. In a simple network consisting of two network layers A and B, each layer has the same number of communities m in which A i and B i are of size N i , i = 1,2, …...m. Every node in A i has exactly one interdependent node in B i , and vice versa [see Fig. 1a]. We define the generating function of A i as where p iA (k) is the degree distribution of nodes in community i of layer A. We define the community structure using the distribution of both intra-and intercommunity links 14,16 . For a community A i within network layer A and has average degree 〈K iA 〉, we use 〈k ijA 〉 to denote the average number of intercommunity links per node in community A i that connects to community A j . We further let α ij; = 〈k ijA 〉/〈K iA 〉 to be the fraction cross community links in A i . These parameters can also be defined in terms of stub that is a link with one end from a node and the other end not connected yet. Then in the network formation process one make α ii fraction of stubs in community A i to be connected to other stubs in community A i , and connect α ij fraction of stubs in community A i to stubs in community A j . β ii and β ij are defined similarly for community B. The generating functions for communities A i and B i are G iA ðξ 1 ; ξ 2 ; :::; ξ m Þ ¼ P P iA ðk 1 ; k 2 ; :::; k m Þξ k 1 1 :::ξ k m m and G iB ðζ 1 ; ζ 2 ; :::; ζ m Þ ¼ P P iB ðk 1 ; k 2 ; :::; k m Þζ k 1 1 :::ζ k m m , respectively 20 . P iA (k 1 , k 2 , ..., k m ) is the probability of finding a node in community A i with k 1 links connecting to nodes in community 1, k 2 links connecting to nodes in community 2,…, and k m links connecting nodes to community m. P iB (k 1 , k 2 , ..., k m ) is defined analogously for network layer B. These generating functions can be expressed by substituting x in Eq. (1) with (ξ 1 , ξ 2 , ..., ξ m ) (derivation provided in Methods), We next attack community A in one of the network layers (say A i ) by initially removing a fraction of 1−p i nodes in community A i . In interdependent networks, it is usually assumed, based on percolation theory, that nodes become non-functional if they become disconnected to the network giant component 1,5,[21][22][23][24][25] . That is all nodes outside the giant component in network A are removed, and so are their interdependent nodes in network B. This cascading process continues until no more nodes can be removed from the system. We finally obtain the mutually connected giant component of the remaining functional nodes.
As we remove nodes and links during the cascading process, the entire network breaks down into the components connected through inter-and intra-community links. The components in community A i comprises the set of A i nodes belonging to a single percolation component of the entire network. The largest of these components is the giant component of A i , and only nodes of the giant component continue to function. Figure 1b-d shows the various types of giant component in multicommunity networks. Note that some nodes are not linked to nodes in the same community but to nodes in the other community.
At the end of the cascading process, the remaining size of A i and B i is where

G iB
1 ðζ 1 ; ζ 2 ; :::; ζ m Þ is defined analogously for community B i . For example, for a pair of Erdös-Renyi networks 26 , Eq. (3) reduces to a simple form, the derivation of which is provided in Supplementary Note 1. The internal structure of our network model is a generalization of a model 1 in which networks have no internal community structure.
To demonstrate critical phenomena in our network model, we consider a simple case of two equally sized communities that are symmetrical in each layer of the network 27 , i.e., m = 2. Without loss of generality, we attack community 1 at the initial stage, i.e., p 1 < 1 and p 2 = 1. We also set α 11 = α 22 , Here the set of parameters α 11 and β 11 is sufficient to describe the community structure. A key parameter that quantifies the robustness of the system is the size of the critical value p 1c that describes the threshold below which the entire system disintegrates with no remaining functional giant component. Thus the smaller the p 1c value, the less vulnerable the network, implying that when p 1c = 0 the network is perfectly robust.
Our first main finding is that a stronger community structure does not always increase the robustness of the interdependent networks, a phenomenon significantly different from the one found in single networks. In single networks, the stronger community structure always increases robustness as shown in Fig. 2a (see Supplementary Note 2 for derivation). Differently, Fig. 2b shows the critical point p 1c against the community strength α 11 of network layer A when the community strength β 11 of network layer B is fixed. The average degree is fixed at 〈K A 〉 = 〈K B 〉= 4 in both network layers. When the community structure is strong in layer B-which occurs when β 11 > 0.436-there is a monotonic dependence of p 1c on α 11 , but when the community structure in B is relatively weak (i.e., when β 11 < 0.436) the behavior of p 1c is non-monotonic.
In order to better understand this, we draw the contour of p 1c with respect to the changes of community structures strength in Fig. 2c. Figure 2c shows that when β 11 < 0.436 the contour of p 1c is the bulging of the equipotential line, specified by the constant value of p 1c , cuts the horizontal line characterized by constant β 11 value twice, resulting in non-monotonous changes as the parameter α 11 changes. Besides this, we find that that system robustness falls into one of the three regions (see Fig. 2c): "vulnerable" p 1c >p r 1c À Á , "robust" p 1c <p r 1c À Á , and "safe" (p 1c = 0, which means that even after the removal of all nodes in a community, there still exist a giant component of the system) regions with the increase of community structure strength in both two network layers A and B. Here p r 1c is the critical point of the corresponding interdependent networks without any community structures (α 11 = β 11 = 0.5).
We also note that the contour lines intersect at two interesting symmetrical points: (α 11 , β 11 ) equals to (0.436, 1) and (1, 0.436) (see Fig. 2c and Supplementary Note 3 for deviation), which implies p 1c changes abruptly at these two unusual points. α 11 = 1 or β 11 = 1 means network A or B has two disconnected (localized) communities. This is a good approximation of real-world network segmentations which could be either geographical or political/economic imposed by embargo/sanctions. For instance, the world-wide business-flight network, which we examine below, has β 11 = 0.98. We find that when one network layer has disconnected communities, as the community strength in another network layer weakens, the critical point suddenly jumps from 0 to a finite number, which is our second main finding. This abrupt jump of the critical point value is a first-order phase transition in which a small change in community strength dramatically increases structural risk (see Fig. 2d and Supplementary Note 4 for derivation).
Empirical implications. An example of a multilayered community structure network is the network of global cities. We examine the data from the system of North American, European, and Asian cities-three different communities-in which transportation and business connections among them define two layers of network in the system. As expected, there are more connections among cities located on the same continent than that among cities located on different continents 28 . We collect business and flight data for 145 North American cities, 158 Asian cities, and 334 European cities and for companies across 21 major industrial sectors in 2010 and 2013 (see Supplementary Note 10 for data description) 29,30 . We use the data to construct an interdependent network of business and flight connections among these cities, where cities are nodes and business and flight connections are links between the cities. In the business network layer (Fig. 3a), a connectivity link between two cities is formed when at least 10 pairs of companies have business connections with each other. In the flight network layer (Fig. 3b), a link is formed when there are at least 200,000 passenger trips between the two cities annually. Simulations demonstrate that the results are not sensitive to these two threshold values (see Supplementary Fig. 5).
The business and transportation network layers are interdependent because businessmen must travel to conduct their business, and airport of one city also dependents on the companies of this city. Usually, it is not easy to obtain the interdependency relationships on important infrastructure networks, thus we begin by assuming a single interdependency link 1 between business and transportation in the same city. We then assume the typical failure mechanism for interdependent networks in which the failure of one node leads to the failure of its interdependent node.
In numerical simulations, we first construct three interdependent networks comprising only cities in North America vs Asia, Europe vs Asia, and Europe vs North America. Secondly, altering the business network A but fixing the flight network B (β 11 is close to 1), we tune down community strength α 11 of business network by rewiring some of its links while keeping their node degrees unchanged (see Supplementary Note 5-8 for the technical details). This approach enables us to assess the system robustness, i.e. how p 1c depends on α 11 , and for what range of parameter values does the system reside in the safe phase or the vulnerable phase. As shown in Fig. 3c-e, we see a phenomena similar to Fig. 2d that p 1c changes abruptly from 0 to a finite value, signaling the emergence of the community structure risk in the system. As the trend of globalization continues, we expect the community structure in the business network to become weaker. Indeed this was the case when comparing the value of community strength α 11 value between year 2010 and 2013 in Fig. 3 for each pair of communities North America-Asia, Europe-Asia, and Europe-North America. Our results from Fig. 3c-e show that as globalization drives global economic system to be more integrated, the North American and European economies exhibit larger robustness than the Asian economy. Finally, we construct an interdependent network of business and flight connections comprising all cities in North America, Europe, and Asia and we find that the overall system is still shifting towards the abrupt change through removing nodes within the Europe community (Fig. 3f). Our results and analysis are based on communities divided according to different continents, and use them to serve as proxies to represent the coarse structure of the world economy. Using our mathematical framework, we can estimate the risk towards the unstable region given the trend from 2010 to 2013. In essence, the world economy, according to this mechanism, had been likely heading towards the unstable phase. In Supplementary Note 9 we discuss some additional weaker forms of interdependency between two network layers characterized by probabilistic interdependency links.

Discussion
In summary, we have studied the robustness of multilayer interdependent networks under community attacks, characteristic of many real system constrained by geographic locations. Through our developed framework, we find that such a system Critical behavior of network(s) with community structure. a Critical transition point p 1c of a single layer random network subjecting to community attack is plotted against its community strength α 11 . As the positive community structure strengthens with increasing α 11 , p 1c decreases monotonically before reaching 0, after which point the two communities both can have a giant component even when the other community is removed from the network. b Critical transition point p 1c of interdependent networks subjecting to community attack is plotted against the community strength α 11 of layer A for different β 11 values of layer B. For large β 11 values (0.5 and 0.8), the behavior is similar to single networks with monotonic decrease of critical threshold p 1c value. Whereas for small β 11 (0 or 0.4), the change of p 1c is non-monotonic and does not reach the safe state of p 1c = 0. The solid red circle represents random failures of an unstructured random network system which is equivalent to α 11 = β 11 = 0.5. c Phase diagram of p 1c with respect to the changes in community structures described by α 11 and β 11 in interdependent networks. The upper right region is the safe region, in which even if we remove a whole community, there still exist a giant component in the whole system. The robust region between boundary p 1c = 0 and p 1c ¼ p r 1c represents a phase that is more robust than unstructured networks without community structures. The vulnerable region below the boundary p 1c ¼ p r 1c represents a phase that is more vulnerable than unstructured networks. In b and c, the two interdependent networks are Erdös-Rényi (ER) networks and both have two communities with equal size and equal average degree 4. d When one layer of the networks is completely localized (β 11 = 1), extreme community structure risk exists when the community strengths α 11 decreases, as p 1c changes discontinuously from 0 to 0.34. In the system, K = 4 for two equally sized communities. The network has a total of N ; = 10,000 nodes, with each community having 5000 nodes. Each data point is an average of 10 simulations and the standard deviation is presented by the error bar ARTICLE COMMUNICATIONS PHYSICS | https://doi.org/10.1038/s42005-019-0144-6 exhibits much more complex behaviors than the ones without clear community structures. Depending on the respective community strengths of the different network layers, the whole system may reside in one of the three different regions of phase diagram, with different phase transition behaviors. In particular, we find that the system embeds extreme risk when one of the layers has strong community structure. In contrary to a single network which is always more robust with enhanced community structure, interdependent networks showing strong community structure sometimes makes the system more vulnerable. In cases of extremely strong community structure with very sparse intercommunity in one layer of network, a small change in the other layer of network's community strength could induce abrupt changes in total systemic resilience if one community is under attack. This new finding adds to the growing knowledge of resilience of interdependent networks with community structure, in particular attacks on specific communities that is representative of realistic events like natural disaster and economic embargo. From the global business-flight networks taken as a proxy for the world economy, we observe that, as globalization weakens the community structure, the entire network is approaching a state in which there is a potential for abrupt disruptions in certain communities. This also leads to the finding that Asia strongly depends on North America and Europe economies but not the other way round. Our result is indicative for a broad range of systems which have community structures defined by geographic location and physical infrastructure. While we mainly focused on the analytical results of two communities on two-layer interdependent networks, our adopted mathematical framework allows for more complex network structures, and their phase transition behaviors could be very interesting for further studies. On the empirical side, our methods can be similarly extended to study other real interdependent networks with community structures, to give better understanding of their dynamical  Table 1 and Supplementary Tables 1 and 2. a Firm and b flight networks between Asia and North America. The red nodes represent cities in Asia, and the blue nodes represent the cities in the North America. The links represent air travel activities in the flight network and business activities in the firm network. We consider the Asian cities as one community and the North America cities as another community in both networks. c Empirical study of the structural risks shows the same abrupt change. The two communities are cities in North America and Asia, and the two layers of networks are business interaction networks and flight networks. The critical transition point p NA c (removing nodes within the North America community) is plotted against community structure strength α 11 of the business layer for the business-flight networks. The original value of α 11 is the empirical value of the system in 2010 and is marked on the horizontal axis (green dot), while the new value of α 11 is measured in 2013. Each data point is an average of 1000 simulations and the shaded region is the bounded by the standard deviation. The result shows that the whole system is shifting towards the abrupt change, where there will be extreme systemic risk. The decrease in α 11 value is possibly due to the effect of globalization, which is weakening the community structure. d, e Similar to c but the results is on the two communities of Europe vs. Asia and North America. From 2010 (green dot) to 2013 (purple dot), the system has shifted from the safe phase to the critical point p E c (removing nodes within the Europe community) close to the vulnerable state. f Similar trend for the systemic risks can be observed in an interdependent network involving all of three communities, and the node removal is done on the Europe community. For all of the simulations in c-f, we only show the attack on one community, as the attack on the other communities does not have the abrupt change on the systemic resilience, i.e. no p c is found.
behaviors for risk assessment, which then leads to the development of efficient risk mitigation strategies in those systems.

Methods
Multivariate generating function. In this subsection, we harness multivariate generating function to formalize the cascading process of interdependent networks with multiple communities 20,30,31 . We first construct the multivariate generating function of the interdependent networks with communities. For conciseness, we illustrate the generating function associated with the bi-community structure (a network with two communities). The result can be easily generalized to account for multicommunity systems. For bi-community system, the probability normalization requires that α 12 = 1 − α 11 and α 21 = 1 − α 22 and the inter-link identity requires that N 1 ⋅〈K 1A 〉⋅α 12 = N 2 ⋅〈K 2A 〉⋅α 21 . Therefore one can describe the structure of network A by one parameter α 11 , because the other three parameters can be calculated from the above three relations. By the definition of α 11 , the probability of finding a node in community 1 with intra-degree k 1 and inter-degree k 2 should be k 1 þ k 2 k 1 α k 1 11 ð1 À α 11 Þ k 2 . Such node with degree k 1 + k 2 has a prior probability p iA (k 1 + k 2 ). Thus one can write the generating function of community 1 in terms of ξ 1 and ξ 2 : The prefactor of ξ k 1 1 ξ k 2 2 represents the probability of finding a node in community 1 of network A with intra-degree k 1 and inter-degree k 2 . By a similar derivation, one can write out the generating function for community i in multicommunity network A and B as G iA ðξ 1 ; ξ 2 ; ::; Percolation equation. Here we show the percolation equation in the form of generating functions 32 . Initially a fraction (1 − p i ) of nodes from community i are randomly removed. The initial removal is followed by a cascading process of percolation failures of the rest nodes. Recall a single complex network with generating function G(x) and branching process G 1 (x) = G′(x)/G′(1) (ref. 33 ), it is known that after a random removal of 1 − p nodes, the generating function of the degree distribution of the remaining nodes can be written by G remain (x) = G(1 − p (1 − x)). At step n, say the size of the giant component of community i in layer A and B due to percolation failures are g n i and h n i . Analogously, we introduce u n i and v n i as the probability of finding a non-functional node by following a link of a randomly chosen node in community i of network A and B at step n after the initial removing of 1−p i fraction of nodes from each community, respectively. We define f n i by The physical meaning of f n i is the probability of finding a node along a randomly selected link in community i of the original network to be non-functional at step n. So the fraction of nodes that are still functional at step n + 1 is :::; f n m Þ; ð8Þ or equivalently À v n 2 Þ; :::; 1 À p m ð1 À u n m Þð1 À v n m ÞÞ: ð9Þ in which G iA 1 ðu n 1 ; u n 2 ; :::; u n m Þ is the generating function associated with the branching process: Similarly, the fraction of functional nodes h nþ1 i in community i in network B at step n ; + 1 is in which v n i is the probability of finding functional node by following a link of a randomly chosen node in community i of network B. It satisfies the following transcendental equation: where G iB 1 ðv n 1 ; v n 2 ; :::; v n m Þ is the generating function related to the branching process: At equilibrium when no more nodes can be removed due to percolation and interdependency failures, the above important parameters u 1 i and v 1 i should obey the following self-consistent equations: and the remaining size of i'th community μ 1 i becomes i ¼ 1; 2; :::; m

Data Availability
The data that support the findings of this study are available from the corresponding author upon reasonable request.
Received: 13 September 2018 Accepted: 11 March 2019 Community size, average degree, and community structure index are indicated α r ii ¼
Competing interests: The authors declare no competing interests.
Reprints and permission information is available online at http://npg.nature.com/ reprintsandpermissions/ Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/ licenses/by/4.0/.

Supplementary Figures
Universal boundary between Safe region and Robust region in terms of (hK 1 i↵ 11 , hK 2 i 11 ); Dashed are their asymptotic lines; Under this boundary, when we remove a full community, the other community will be fragmented. Above the boundary, when we remove a full community, the other community can survival with a mutual giant component at the equilibrium state. (b) Optimum ↵ opt 11 vs 11 .The optimized ↵ opt 11 means that it can maximize the percolation robustness (minimizing p 1c ). When 11 < 0.436, ↵ opt 11 increases from 0.5 to 1 as 11 increases. When 11 0.436, ↵ opt 11 remains 1. To increase ↵ 11 , each time one intra link from each community is randomly selected (denoted by the dashed lines). Then they are replaced by two inter-links denoted by the thick (red online) lines. The process repeats until ↵ 11 reaches the target value. To decrease ↵ 11 . Each time two inter links between these two communities are selected (denoted by the dashed lines). They are then replaced by two intra links (denoted by thick lines (red online)) connecting nodes of the same community. The above process repeats until ↵ 11 reaches expected value. This way can keep the original network structure as much as possible. Note that, if the two communities are of different sizes or average degrees, the maximum ↵ 11 can be 1, but the minimum ↵ 11 may not reach 0.  Fig.3 in the main text. We only present the attack on the first community which is North America of (a), Europe for (b) and (c). For the attack on the second community in respective plots, the system does not disintegrate with when the whole of second community is removed, leading to p 2c = 0 regardless of ↵ 11 values, and it is not plotted here.

Supplementary Note 1 Interdependent Erdős-Rényi networks
Here we show the percolation theory for Erdős-Rényi (ER) networks. For simplicity we assume N 1 = N 2 . The generating functions for ER networks with two equal communities can be analytically expressed by the average 4 degrees: By substituting the above equations into Eq. (14) in the main text, the self-consistent equations reduce to and the final size of the network becomes Because each layer of networks have two equal communities, thus we can use the fact that ↵ 11 = ↵ 22 and 11 = For unbiased attack, p 1 = p 2 = p, the above expressions for µ 1 1 and µ 1 2 can be written as the following form: One can see this equation is independent of the community structure, which implies that the percolation failures of interdependent ER networks with two equal communities is independent of the community structure.

Supplementary Note 2 Single networks with community structure
Our percolation model can be applied to single networks with community structure by removing terms that are associated with network B. The self consistency equations are simplified to and the remaining size of i'th community µ 1 i becomes with i = 1, 2, ..., m. If we consider ER networks, the above self consistency equations can be written as and the remaining giant component sizes reduce to with i = 1, 2, ..., m, corresponding to the m communities in the network A. For community attack one has p 1 = p and p 2 = 1. For community attack, it can be shown that the transition point p 1c satisfies We can see that, if p 1c = 0, it implies that hk 22 i  1. This allow us got a 22  1 K2A .

Supplementary Note 3 Safe and robust boundaries on ER network systems
For community attack, say we only attack community 1, i.e. p 1 = p and p 2 = 1. Now we consider the safe-unsafe transition. At the transition point, the complete failure of community 1 leads to the failure of community 2. Within the safe region, the removal of community 1 does not lead to the failure of community 2. In other words, whether or not community 2 can self survive determines the safe-unsafe transition boundary. Thus we set µ 1 1 = 0, the second self-consistent equation becomes We have hk 22A i = ↵ 11 hK 2A i and hk 22B i = 11 hK 2B i, and let x = µ 1 2 , thus x = 0 is the solution to the critical point. If the system is at transition from safe (which has nonzero solution of x) to unsafe region (only zero solution of x exists), one expects the derivatives of the two sides of the above equation should equal: i .
By solving the above equations, one can obtain the universal boundary between safe and unsafe regions. The physical meaning of the above equations are as follow: If community 2 can survive from the full failures of community 1, the intra-links in community 2 can prevent the system from spontaneous collapse. The symmetric solution of the above equations is hk 22A i = hk 22B i = 2.4554, which is the critical average degree (two layers have equal average degree) for an unstructured interdependent network to be stable [1].
To obtain robust-vulnerable boundary, we let p 1 = p r 1c and p 2 = 1 in Eq. 2. p r 1c is the critical transition point for a unstructured network under random attack. For conciseness, let x = µ 1 1 and y = µ 1 2 . To find solution to Eq. 2 graphically, one can plot the two equations in terms of two curves x(y) and y(x), the cross points of which are the solutions. At the critical point p 1 = p r 1c and p 2 = 1, we have nontrivial (nonzero) solutions to x and y.
When p 1 < p r 1c , the system's giant component size drops to 0, or equivalently we have only trivial solutions x = 0 and y = 0. Therefore we should expect that at the critical point the two curves are tangential to each other, i.e. @x @y @y @x = 1 [2]: @x @y Therefore we have 3 equations with 4 variables x, y, ↵ 11 and 11 . By eliminating x and y from the above 3 equations, we obtain an equation that describes the robust-vulnerable boundary curve in terms of ↵ 11 and 11 . The Consider community attack (to community 1) and network B being localized, one has p 1 < 1, p 2 = 1 and 11 = 1.
We look for the transition point p 1c when the whole system fails. Eqs. 2 reduce to in which we use x and y to denote µ 1 1 and µ 1 2 respectively, for conciseness' sake. When community 1 completely fails, one has x = 0 and the second equation above further reduces to critical p 1c requires that the derivatives of the two sides in the above equation should equal, therefore one has The above p 1c is finite. If hK 1A i = hK 1B i = 4 and ↵ 11 = ↵ 11c , p 1c = 0.34. Eq. 6 is similar to Eq.

Supplementary Note 5 Simulations
In this subsection we describe simulation details. To generate networks with community structure, we first generate a Poisson degree distribution for each network. Next we assign each node with a probability of being connected by a randomly chosen link. The probability is proportional to the node's degree. Last we loop over all links and assign each of them with two ends according to nodes' probabilities. For higher probability nodes, they are more likely to be connected by the links. The actual number of links they are connected are proportional to their preassigned degrees. Hence the final network recovers the degree distribution. In order to obtain community structure, we use the same protocol for both intra and inter links. To control the community structure, the number of intra and inter links are evaluated in advance. We generate random networks with size N 1 = N 2 = 10, 000. To calculate giant component size after random or community attacks, we perform at least 10 realizations for each parameter set and take the average.
To find the giant components in our interdependent system, we first look for the largest percolation cluster in the networks without considering the community structure. We later categorize the nodes in the cluster by their associated communities. The giant component in each community is therefore obtained. Specifically, if the communities are disconnected (localized), we look for the largest cluster in each disconnected communities.
Therefore we allow multiple components in our final giant component. As in the real system, localized community can be self supportive without relying on external resources.

Supplementary Note 7 How to change community strength ↵ 11
Supplementary Figure 3 shows a typical case of how we vary community structure to interdependent networks.
After the reconnecting the links, the nodes' degrees preserve and we also make the reconnected network as close as possible to the original one.

Supplementary Note 8 Empirical studies
To manually change community structure of an empirical network, we use the same protocol described in Supplementary Note 7 to regenerate the links in the network. To control the community structure, we evaluate in advance the number of inter-links and intra-links so that the average degree and community index would be obeyed in the following process. Lastly, we assign two nodes to each link to form the network. We do so according to the link property. If the link is an inter-link, we assign nodes of different communities to it; while if the link is an intra-link of community 1, we restrict the selection of nodes to community 1 during the process. To minimize statistical error during the manipulation of empirical data, we perform 100 independent simulations by randomly carrying out the above process, and take the average of any measurement we need.

Supplementary Note 9 Partial dependency case
Here we have also empirically studied the partial interdependent case [2]. For the North America-Asia Firm-Flight interdependent networks, we setup the interdependent ratio between the nodes in two layer to be Q < 1, which 10 is more realistic. It means that a failed node in one layer leads to the failure of the same node in other layer with probability Q. We also detect the abrupt change of p 1c when Q = 0.9 shown in Supplementary Figure 4.    Retail sale in non-specialised stores Retail sale in non-specialised stores with food, beverages or tobacco predominating Other retail sale in non-specialised stores Retail sale of food, beverages and tobacco in specialised stores Retail sale of fruit and vegetables in specialised stores Retail sale of meat and meat products in specialised stores Retail sale of fish, crustaceans and molluscs in specialised stores Retail sale of bread, cakes, flour confectionery and sugar confectionery in specialised stores

TRANSPORTATION AND STORAGE
Land transport and transport via pipelines Passenger rail transport, interurban Freight rail transport Other passenger land transport Urban and suburban passenger land transport Taxi operation Other passenger land transport n.e.c. Freight transport by road and removal services Freight transport by road Removal services Transport via pipeline Water transport Sea and coastal freight water transport 21 Inland freight water transport Air transport Passenger air transport Freight air transport and space transport Freight air transport Space transport Warehousing and support activities for transportation Warehousing and storage Support activities for transportation Service activities incidental to land transportation Service activities incidental to water transportation Service activities incidental to air transportation Cargo handling Other transportation support activities Postal and courier activities Postal activities under universal service obligation Other postal and courier activities 9. ACCOMMODATION AND FOOD SERVICE ACTIVITIES Accommodation Hotels and similar accommodation Holiday and other short-stay accommodation Camping grounds, recreational vehicle parks and trailer parks Other accommodation Food and beverage service activities Restaurants and mobile food service activities Event catering and other food service activities Event catering activities Other food service activities Beverage serving activities

INFORMATION AND COMMUNICATION
Publishing activities Publishing of books, periodicals and other publishing activities Book publishing Publishing of directories and mailing lists Publishing of newspapers Publishing of journals and periodicals Other publishing activities Software publishing Publishing of computer games Other software publishing Motion picture, video and television programme production and music publishing activities Motion picture, video and television programme activities Motion picture, video and television programme production activities Motion picture, video and television programme post-production activities Motion picture, video and television programme distribution activities Motion picture projection activities Sound recording and music publishing activities Programming and broadcasting activities Radio broadcasting Television programming and broadcasting activities 22 Telecommunications Wired telecommunications activities Wireless telecommunications activities Satellite telecommunications activities Other telecommunications activities Computer programming, consultancy and related activities Computer programming activities Computer consultancy activities Computer facilities management activities Other information technology and computer service activities Information service activities Data processing, hosting and related activities; web portals Data processing, hosting and related activities Web portals Other information service activities News agency activities Other information service activities n.e.c.

FINANCIAL AND INSURANCE ACTIVITIES Financial service activities, except insurance and pension funding Monetary intermediation Central banking
Other monetary intermediation Activities of holding companies Trusts, funds and similar financial entities Other financial service activities, except insurance and pension funding Financial leasing Other credit granting Other financial service activities, except insurance and pension funding n.e.c. Insurance, reinsurance and pension funding, except compulsory social security Insurance Life insurance Non-life insurance Reinsurance Pension funding Activities auxiliary to financial services and insurance activities Activities auxiliary to financial services, except insurance and pension funding Administration of financial markets Security and commodity contracts brokerage Other activities auxiliary to financial services, except insurance and pension funding Activities auxiliary to insurance and pension funding Risk and damage evaluation Activities of insurance agents and brokers Other activities auxiliary to insurance and pension funding Fund management activities

REAL ESTATE ACTIVITIES
Real estate activities Buying and selling of own real estate Renting and operating of own or leased real estate Real estate activities on a fee or contract basis Real estate agencies Temporary employment agency activities Other human resources provision Travel agency, tour operator and other reservation service and related activities Travel agency and tour operator activities Travel agency activities Tour operator activities Other reservation service and related activities Security and investigation activities Private security activities Security systems service activities Investigation activities Services to buildings and landscape activities Combined facilities support activities Cleaning activities General cleaning of buildings Other building and industrial cleaning activities Other cleaning activities Landscape service activities Office administrative, office support and other business support activities Office administrative and support activities Combined office administrative service activities Photocopying, document preparation and other specialised office support activities Activities of call centres Organisation of conventions and trade shows Business support service activities n.e.c. Activities of collection agencies and credit bureaus Packaging activities Other business support service activities n.e.c.

PUBLIC ADMINISTRATION AND DEFENCE; COMPULSORY SOCIAL SECURITY
Public administration and defence; compulsory social security Administration of the State and the economic and social policy of the community General public administration activities Regulation of the activities of providing health care and other social services Regulation of and contribution to more efficient operation of businesses Provision of services to the community as a whole Foreign affairs Defence activities Justice and judicial activities Public order and safety activities Fire service activities Compulsory social security activities

EDUCATION Education
Pre-primary education Pre-primary education Primary education Primary education Secondary education General secondary education Technical and vocational secondary education Higher education Post-secondary non-tertiary education Tertiary education Other education Sports and recreation education Cultural education Driving school activities Other education n.e.c. Educational support activities

HUMAN HEALTH AND SOCIAL WORK ACTIVITIES
Human health activities Hospital activities Medical and dental practice activities General medical practice activities Specialist medical practice activities Dental practice activities Other human health activities Residential care activities Residential nursing care activities Residential care activities for mental retardation, mental health and substance abuse Residential care activities for the elderly and disabled Other residential care activities Social work activities without accommodation Social work activities without accommodation for the elderly and disabled Other social work activities without accommodation Child day-care activities Other social work activities without accommodation n.e.c.

ARTS, ENTERTAINMENT AND RECREATION
Creative, arts and entertainment activities Performing arts Support activities to performing arts Artistic creation Operation of arts facilities Libraries, archives, museums and other cultural activities Library and archives activities Museums activities Operation of historical sites and buildings and similar visitor attractions Botanical and zoological gardens and nature reserves activities Gambling and betting activities Sports activities and amusement and recreation activities Sports activities Operation of sports facilities Activities of sports clubs Fitness facilities Other sports activities Amusement and recreation activities Activities of amusement parks and theme parks Other amusement and recreation activities To increase ↵ 11 , each time one intra link from each community is randomly selected (denoted by the dashed lines). Then they are replaced by two inter-links denoted by the thick (red online) lines. The process repeats until ↵ 11 reaches the target value. To decrease ↵ 11 . Each time two inter links between these two communities are selected (denoted by the dashed lines). They are then replaced by two intra links (denoted by thick lines (red online)) connecting nodes of the same community. The above process repeats until ↵ 11 reaches expected value. This way can keep the original network structure as much as possible. Note that, if the two communities are of different sizes or average degrees, the maximum ↵ 11 can be 1, but the minimum ↵ 11 may not reach 0.  Fig.3 in the main text. We only present the attack on the first community which is North America of (a), Europe for (b) and (c). For the attack on the second community in respective plots, the system does not disintegrate with when the whole of second community is removed, leading to p 2c = 0 regardless of ↵ 11 values, and it is not plotted here.

Supplementary Note 1 Interdependent Erdős-Rényi networks
Here we show the percolation theory for Erdős-Rényi (ER) networks. For simplicity we assume N 1 = N 2 . The generating functions for ER networks with two equal communities can be analytically expressed by the average 4 degrees: By substituting the above equations into Eq. (14) in the main text, the self-consistent equations reduce to and the final size of the network becomes i , therefore substituting Eq. 2 into the above equation yields Because each layer of networks have two equal communities, thus we can use the fact that ↵ 11 = ↵ 22 and 11 = For unbiased attack, p 1 = p 2 = p, the above expressions for µ 1 1 and µ 1 2 can be written as the following form: One can see this equation is independent of the community structure, which implies that the percolation failures of interdependent ER networks with two equal communities is independent of the community structure.

Supplementary Note 2 Single networks with community structure
Our percolation model can be applied to single networks with community structure by removing terms that are associated with network B. The self consistency equations are simplified to , ..., 1 p m (1 u 1 m )) 5 and the remaining size of i'th community µ 1 i becomes with i = 1, 2, ..., m. If we consider ER networks, the above self consistency equations can be written as and the remaining giant component sizes reduce to We can see that, if p 1c = 0, it implies that hk 22 i  1. This allow us got a 22  1 K2A .

Supplementary Note 3 Safe and robust boundaries on ER network systems
For community attack, say we only attack community 1, i.e. p 1 = p and p 2 = 1. Now we consider the safe-unsafe transition. At the transition point, the complete failure of community 1 leads to the failure of community 2. Within the safe region, the removal of community 1 does not lead to the failure of community 2. In other words, whether or not community 2 can self survive determines the safe-unsafe transition boundary. Thus we set µ 1 1 = 0, the second self-consistent equation becomes We have hk 22A i = ↵ 11 hK 2A i and hk 22B i = 11 hK 2B i, and let x = µ 1 2 , thus x = (1 e hk22Aix )(1 e hk22B ix ).
x = 0 is the solution to the critical point. If the system is at transition from safe (which has nonzero solution of x) to unsafe region (only zero solution of x exists), one expects the derivatives of the two sides of the above equation should equal: i .
By solving the above equations, one can obtain the universal boundary between safe and unsafe regions. The physical meaning of the above equations are as follow: If community 2 can survive from the full failures of community 1, the intra-links in community 2 can prevent the system from spontaneous collapse. The symmetric solution of the above equations is hk 22A i = hk 22B i = 2.4554, which is the critical average degree (two layers have equal average degree) for an unstructured interdependent network to be stable [1].
To obtain robust-vulnerable boundary, we let p 1 = p r 1c and p 2 = 1 in Eq. 2. p r 1c is the critical transition point for a unstructured network under random attack. For conciseness, let x = µ 1 1 and y = µ 1 2 . To find solution to Eq. 2 graphically, one can plot the two equations in terms of two curves x(y) and y(x), the cross points of which are the solutions. At the critical point p 1 = p r 1c and p 2 = 1, we have nontrivial (nonzero) solutions to x and y.
When p 1 < p r 1c , the system's giant component size drops to 0, or equivalently we have only trivial solutions x = 0 and y = 0. Therefore we should expect that at the critical point the two curves are tangential to each other, i.e. @x @y @y @x = 1 [2]: @x @y Therefore we have 3 equations with 4 variables x, y, ↵ 11 and 11 . By eliminating x and y from the above 3 equations, we obtain an equation that describes the robust-vulnerable boundary curve in terms of ↵ 11 and 11 . The boundary curve is numerically solved and shown in Supplementary Figure 1.

7
Supplementary Note 4 Discontinuous transition from safe to unsafe region Consider community attack (to community 1) and network B being localized, one has p 1 < 1, p 2 = 1 and 11 = 1.
We look for the transition point p 1c when the whole system fails. Eqs. 2 reduce to in which we use x and y to denote µ 1 1 and µ 1 2 respectively, for conciseness' sake. When community 1 completely fails, one has x = 0 and the second equation above further reduces to  Figure   2, one can visualize how the nontrivial solution emerges with ↵ 11 increases. At safe-unsafe boundary, the failures of the two communities occur simultaneously. So we assume y = 0, Eq. 4 becomes critical p 1c requires that the derivatives of the two sides in the above equation should equal, therefore one has The above p 1c is finite. If hK 1A i = hK 1B i = 4 and ↵ 11 = ↵ 11c , p 1c = 0.34. Eq. 6 is similar to Eq. 5 except there is a prefactor p 1 . From Supplementary Figure 2 one can tell that with sufficiently small p 1 , one can always eliminate nontrivial solution. So if ↵ 11 < ↵ 11c , Eq. 5 has no nontrivial solutions, while with sufficiently small but 8 finite p 1 one can make Eq. 6's solution trivial, thus p 1c is finite; If ↵ 11 > ↵ 11c , the second equation of Eq. 4 with x = 0 always has nontrivial solutions, therefore p 1c becomes indefinite; If ↵ 11 = ↵ 11c , one finds p 1c = 0.34. The change of ↵ 11 leads to an abrupt change of the transition point p 1c from finiteness to indefiniteness. However if 11 6 = 1, we no longer have Eq. 5. The equations for y also involves x, ↵ 11 can continuously drops to 0 before nontrivial solution appears.

Supplementary Note 5 Simulations
In this subsection we describe simulation details. To generate networks with community structure, we first generate a Poisson degree distribution for each network. Next we assign each node with a probability of being connected by a randomly chosen link. The probability is proportional to the node's degree. Last we loop over all links and assign each of them with two ends according to nodes' probabilities. For higher probability nodes, they are more likely to be connected by the links. The actual number of links they are connected are proportional to their preassigned degrees. Hence the final network recovers the degree distribution. In order to obtain community structure, we use the same protocol for both intra and inter links. To control the community structure, the number of intra and inter links are evaluated in advance. We generate random networks with size N 1 = N 2 = 10, 000. To calculate giant component size after random or community attacks, we perform at least 10 realizations for each parameter set and take the average.
To find the giant components in our interdependent system, we first look for the largest percolation cluster in the networks without considering the community structure. We later categorize the nodes in the cluster by their associated communities. The giant component in each community is therefore obtained. Specifically, if the communities are disconnected (localized), we look for the largest cluster in each disconnected communities.
Therefore we allow multiple components in our final giant component. As in the real system, localized community can be self supportive without relying on external resources.

Supplementary Note 7 How to change community strength ↵ 11
Supplementary Figure 3 shows a typical case of how we vary community structure to interdependent networks.
After the reconnecting the links, the nodes' degrees preserve and we also make the reconnected network as close as possible to the original one.

Supplementary Note 8 Empirical studies
To manually change community structure of an empirical network, we use the same protocol described in Supplementary Note 7 to regenerate the links in the network. To control the community structure, we evaluate in advance the number of inter-links and intra-links so that the average degree and community index would be obeyed in the following process. Lastly, we assign two nodes to each link to form the network. We do so according to the link property. If the link is an inter-link, we assign nodes of different communities to it; while if the link is an intra-link of community 1, we restrict the selection of nodes to community 1 during the process. To minimize statistical error during the manipulation of empirical data, we perform 100 independent simulations by randomly carrying out the above process, and take the average of any measurement we need.

Supplementary Note 9 Partial dependency case
Here we have also empirically studied the partial interdependent case [2]. For the North America-Asia Firm-Flight interdependent networks, we setup the interdependent ratio between the nodes in two layer to be Q < 1, which Construction of roads and railways Construction of roads and motorways Construction of railways and underground railways Construction of bridges and tunnels Construction of utility projects Construction of utility projects for fluids Construction of utility projects for electricity and telecommunications Construction of other civil engineering projects Construction of water projects Construction of other civil engineering projects n.e.c. Specialised construction activities Demolition and site preparation Demolition Site preparation Test drilling and boring Electrical, plumbing and other construction installation activities Electrical installation Plumbing, heat and air-conditioning installation Other construction installation Building completion and finishing Plastering Joinery installation Floor and wall covering Painting and glazing Other building completion and finishing Other specialised construction activities Roofing activities Other specialised construction activities n.e.c.

WHOLESALE AND RETAIL TRADE; REPAIR OF MOTOR VEHICLES
Wholesale and retail trade and repair of motor vehicles and motorcycles Sale of motor vehicles Sale of cars and light motor vehicles Sale of other motor vehicles Maintenance and repair of motor vehicles Sale of motor vehicle parts and accessories Wholesale trade of motor vehicle parts and accessories Retail trade of motor vehicle parts and accessories Sale, maintenance and repair of motorcycles and related parts and accessories Wholesale trade, except of motor vehicles and motorcycles Wholesale on a fee or contract basis Agents involved in the sale of agricultural raw materials, live animals, textile raw materials Agents involved in the sale of fuels, ores, metals and industrial chemicals Agents involved in the sale of timber and building materials Agents involved in the sale of machinery, industrial equipment, ships and aircraft Agents involved in the sale of furniture, household goods, hardware and ironmongery Agents involved in the sale of textiles, clothing, fur, footwear and leather goods Agents involved in the sale of food, beverages and tobacco Agents specialised in the sale of other particular products Agents involved in the sale of a variety of goods Wholesale of agricultural raw materials and live animals Wholesale of grain, unmanufactured tobacco, seeds and animal feeds Temporary employment agency activities Other human resources provision Travel agency, tour operator and other reservation service and related activities Travel agency and tour operator activities Travel agency activities Tour operator activities Other reservation service and related activities Security and investigation activities Private security activities Security systems service activities Investigation activities Services to buildings and landscape activities Combined facilities support activities Cleaning activities General cleaning of buildings Other building and industrial cleaning activities Other cleaning activities Landscape service activities Office administrative, office support and other business support activities Office administrative and support activities Combined office administrative service activities Photocopying, document preparation and other specialised office support activities Activities of call centres Organisation of conventions and trade shows Business support service activities n.e.c. Activities of collection agencies and credit bureaus Packaging activities Other business support service activities n.e.c.

PUBLIC ADMINISTRATION AND DEFENCE; COMPULSORY SOCIAL SECURITY
Public administration and defence; compulsory social security Administration of the State and the economic and social policy of the community General public administration activities Regulation of the activities of providing health care and other social services Regulation of and contribution to more efficient operation of businesses Provision of services to the community as a whole Foreign affairs Defence activities Justice and judicial activities Public order and safety activities Fire service activities Compulsory social security activities

EDUCATION Education
Pre-primary education Pre-primary education Primary education Primary education Secondary education General secondary education Technical and vocational secondary education