Dynamical resilience of networks against targeted attack

Physica A 528 (2019) 121329

Contents lists available at ScienceDirect

Physica A journal homepage: www.elsevier.com/locate/physa

Dynamical resilience of networks against targeted attack ∗

Feifei Xu a , Shubin Si a , Dongli Duan b , , Changchun Lv a , Junlan Xie a a b

School of Mechanical Engineering, Northwestern Polytechnical University, Xi’an 710072, PR China School of Information and Control Engineering, Xi’an University of Architecture and Technology, Xi’an 710311, PR China

highlights • • • •

Quantify the mathematical relationship of network resilience and attack scale f. Get the critical breakdown conditions fc with percolation theory. For BD dynamics, SF and ER networks are both robust with fc ≈ 1 . For R or E dynamics, ER networks are more vulnerable compared to structural robustness.

article

info

Article history: Received 31 October 2018 Received in revised form 3 March 2019 Available online 11 May 2019 Keywords: Resilience Complex networks Dynamical processes Critical thresholds

a b s t r a c t To mitigate the damage caused by external interference, many researchers have studied the resilience changes in network when it is under targeted attack. However, it is still unclear how dynamical processes impact network resilience. Here, with a broad range of steady-state dynamical processes including birth-death processes ( BD), regulatory dynamics (R) and epidemic processes (E ) on Scale-Free (SF ) and Erdős– Rényi (ER) networks, we explore the resilience of complex networks under two attack strategies: from high-degree nodes (HD) and from low-degree nodes (LD). Mapping the multi-dimensional dynamics equation into one-dimensional equation, we quantify the relationship between network resilience and attacked node fraction f and present the critical thresholds fc at which the network loses its resilience. When take dynamical processes into consideration in resilience research, we get some novel conclusions. Compared to structural robustness without dynamical processes, with R or E , ER networks are more vulnerable and the heterogeneity of SF networks has different effects on thresholds fc under HD attack strategy. The theoretical solutions are consistent with the simulation results to some extent, our outcomes are helpful for optimizing networks and enhancing the resilience of networks. © 2019 Elsevier B.V. All rights reserved.

1. Introduction Complex systems and networks have stirred considerable research attention due to their wide applications in sciences and engineering. Network resilience is defined as the ability of network to adjust its activity to retain its basic function when errors, failures or environmental changes occur [1]. The resilience of real-world networks subject to targeted attack is very significant for the study of network safety [2,3], which is a research hot spot in network science nowadays. The resilience can be influenced by the changes of the network structure as well as the dynamics of network components [4]. Liu et al. [5] pointed out that studying the dynamical mechanism of the cascading failure can help us to design more ∗ Corresponding author. E-mail address: [email protected] (D. Duan). https://doi.org/10.1016/j.physa.2019.121329 0378-4371/© 2019 Elsevier B.V. All rights reserved.

2

F. Xu, S. Si, D. Duan et al. / Physica A 528 (2019) 121329

resilient networks. Yuan et al. [6] studied the dynamical tolerance of oscillator networks against targeted attacks, which provides a reference for studying the dynamical resilience of other types of dynamical systems. The complex networks with dynamical behavior behave universal resilience patterns which allow us to predict the system’s resilience [1] and tipping points [7], thus guiding us for reasonably deploying networks and enhancing security of the networks. In this paper, we focus on the network resilience against targeted attack with dynamical processes. To capture the resilience in networks, one needs a mathematical scaffold which able to seize two basic aspects: the structural robustness of the networks and the failure spreading behaviors within the networks. So far, the first issue was dealt with by examining the changes of the structure connectivity caused by node removal [8], the critical properties for the disintegration of networks were unveiled with percolation theory [9]. In particular, Crucitti et al. [10] proved that SF networks were robust to random removal and vulnerable to intentional attack. Wu et al. [11] proposed the natural connectivity to accurately quantify the structural robustness of networks. Increasing evidence showed that many actual networks interact with each other, triggering a research hotspot for the structural robustness of interdependent networks [12,13]. As for the second issue, many attempts were made to model the failure spreading process. Motter et al. [14] introduced the overload mechanism of cascading failure to study the cascades triggered by random breakdown and intentional attacks. Wang et al. [15] found the abnormal fluctuation of cascading dynamics. Subsequently, many theoretical models were raised to describe the emergence of the cascading failure, such as the threshold model [16–18], the epidemic model [19–21], the cascade model in interdependent networks [22–24] and so on. The researches above mainly focused on the changes of structural connectivity after initial failure and attack in the networks. The spreading mechanisms of dynamical behavior in real-world complex systems often can be represented by dynamical processes, including birth-death processes (BD), regulatory dynamics (R), epidemic processes (E ) and so on. In the dynamical equations which describe the dynamical processes, the state of a node is affected by the state of nodes connected to it, the network topology, and some inherent parameters, dynamically changing with time. Barzel et al. [25] developed a mathematical framework that uncovered the universality of the interaction between the topology and the dynamics of complex systems. Recently, the resilience research of networks began to be combined with the dynamical processes. Gao et al. [1] derived effective one-dimensional dynamics which can be used to predict the system’s resilience. Jiang Junjie et al. [7] proposed a mechanism to predict the critical points of the ecological dynamical systems. Jiang Jian et al. [26] analyzed an epidemic dynamical process on a two-layer network to explore how time delay impacts epidemic spreading on multilayer networks under limited resources. We can see that it remains unexplored how dynamical processes impact network resilience after a fraction f of nodes attacked. In this paper, the main contribution is to quantify the relationship between network resilience and attacked node fraction f with BD, R and E dynamics. We also obtain the critical thresholds fc at which the network loses its resilience, which is helpful to understand the dynamical mechanism under critical conditions in real-world networks. We find that network resilience with dynamical processes has different properties compared to structural robustness without dynamics. Our researches provide a new framework for the analysis of network resilience with dynamical processes. This paper is arranged as follows: with average activity ⟨x⟩ characterizing network resilience, Section 2 maps the multidimensional dynamics equations into one-dimensional rate equations, then the theoretical solutions of ⟨x⟩ for BD, R and E dynamics are obtained. Applying these three dynamical models on SF and ER networks under two attack strategies: from high-degree nodes (HD) and from low-degree nodes (LD), Section 3 presents the simulation results of average activity ⟨x⟩, and shows the comparison of simulation results and theoretical solutions. In Section 4, the critical thresholds fc of attacked node fraction are obtained; we get several novel conclusions through threshold analysis. Section 5 concludes and discusses this paper. 2. Dynamical models and theory solutions of network resilience 2.1. Dynamical models Barzel et al. proposed a general description of systems governed by pairwise interactions [25] dxi dt

= W (xi ) +

N ∑

Aij Q xi , xj ,

(

)

(1)

j=1

where N is the number of nodes in the system, xi represents the activity of node i. The first term on the right of the equation describes the self-dynamics of each (node,)while the second term shows the interactions between node i and its neighbors. Aij is the adjacency matrix and Q (xi , xj ) represents the dynamical laws governing the pairwise interactions. With the appropriate choice of W (xi ) and Q xi , xj , Eq. (1) can be used to model several kinds of different dynamical systems, such as cellular, ecological and social systems. Eq. (1) is a multi-dimensional model over the complex parameter space which cannot be treated analytically. If we can reduce the multi-dimensional equation to a one-dimensional equation, then we can analyze it analytically [1]. First, we need to condense the multi-dimensional variables in the equation into one-dimensional variables. Considering that ∑N the network average activity can characterize network resilience, we average the activity of all the nodes as ⟨x⟩ = 1/N i=1 xi , so, the N-dimensional variables xi are condensed into a one-dimensional variable ⟨x⟩. Then for the second term on the

F. Xu, S. Si, D. Duan et al. / Physica A 528 (2019) 121329

3

Table 1 Dynamical models and resilience breakdown conditions. Dynamics

Rate equation

BD

dxi dt

R

dxi dt

E

dxi dt

∑ = −Bxi b + Nj=1 Aij Rxaj ∑ xh = −Bxi + Nj=1 Aij R 1+jxh j ∑ = −Bxi + Nj=1 Aij R(1 − xi )xj

One-dimensional equation

⟨x⟩c

d⟨x⟩ dt

= −B⟨x⟩b + ⟨s⟩ R⟨x⟩a

0

d⟨x⟩ dt

⟨x⟩ = −B⟨x⟩ + ⟨s⟩R 1+⟨ x⟩h

(h − 1) h

d⟨x⟩ dt

= −B⟨x⟩ + ⟨s⟩R(1 − ⟨x⟩)⟨x⟩

0

h

1

Birth-death processes BD [27,28]: Birth-death processes provide a framework for modeling large variety of biological processes, including genome evolution, population dynamics, and speciation. This equation describes the population density at site i, emerging in queuing theory, population dynamics and biology. Here we set b = 2 to represent pairwise depletion, and a = 1 to describe linear flow. Regulatory dynamics R [29,30]: this dynamics of gene regulation is described by the Michaelis–Menten (MM) equation. Here, the first term on the right side of dynamical equation describes degradation, where xi is the expression level of a gene; the second term captures the genetic activation, where the parameter h is the Hill coefficient which describes the level of cooperativity in gene regulation, we set h = 1. Epidemic dynamics E [31,32]: in this equation, the spread of diseases is described by the susceptibleinfected-susceptible (SIS) model. In this model, each of the individuals in the network is assumed to occupy one of two states: susceptible (S) or infected (I). Each node can be seen as an individual and each link is a connection along which the diseases can spread to other individuals. For all the three dynamics, B and R are rate constants, which we set to B = R = 1.

right of Eq. (1), the adjacency matrix Aij includes N 2 parameters whose changes capture different kinds of perturbations, such as node removal, link removal, weight reduction or any combination. If not the type of ∑ perturbations, ∑Nconsidering ∑N N we can condense the N 2 parameter of Aij into a single parameter as ⟨s⟩ = 1/N i=1 j=1 Aij , where si = j=1 Aij is the degree of node i, so ⟨s⟩ is the average degree of all the nodes. Therefore, Eq. (1) can be mapped to a one-dimensional equation [1] d ⟨x ⟩

= W (⟨x⟩) + ⟨s⟩ Q (⟨x⟩ , ⟨x⟩). dt The regulatory dynamics (R) captured by the Michaelis–Menten (MM) equation can be described as dxi dt

= −Bxi +

N ∑

Aij R

j=1

xhj 1 + xhj

,

(2)

(3)

where xi is the expression level of a gene, B and R are rate constants, h is the Hill coefficient which describes the level of cooperativity in gene regulation. Mapping Eq. (3) to the form of Eq. (2), we arrive at the one-dimensional equation d ⟨x ⟩ dt

= −B ⟨x⟩ + ⟨s⟩ R

⟨x ⟩h . 1 + ⟨x ⟩h

In order to provide the system’s steady state, we set f (⟨x⟩ , ⟨s⟩) =

) B 1 + ⟨x⟩h

(4) d⟨x⟩ dt

= 0, then we can get ⟨s⟩ in function of ⟨x⟩ as

(

⟨s ⟩ =

R⟨x⟩h−1

In addition, we set

−B + ⟨s⟩ R (

.

(5)

∂ f (⟨x⟩,⟨s⟩) ∂⟨x⟩

h⟨x⟩h−1

1 + ⟨x ⟩h

= 0 to obtain the critical threshold, we get

)2 = 0.

(6)

By solving the simultaneous equations of Eqs. (5) and (6), we get the threshold 1

⟨x⟩c = (h − 1) h .

(7)

The threshold ⟨x⟩c represents the critical point at which the networks transform from active state to inactive state. With the method described above, we calculate the ⟨x⟩c of three dynamical models BD, R and E , which are shown in Table 1. The critical thresholds ⟨x⟩c are only related to different dynamics, and independent of the network topology. Characterizing the network resilience by average activity ⟨x⟩, we can predict the resilience of complex networks effectively and accurately by exploring the changes of ⟨x⟩. 2.2. Theory solutions of network resilience We extend the study of network resilience to targeted attack on a fraction f of nodes in the networks with above dynamical processes. To explore the analytic relationship between average activity ⟨x⟩ and attacked node fraction f , we

4

F. Xu, S. Si, D. Duan et al. / Physica A 528 (2019) 121329 Table 2 Theory solutions of average activity ⟨x⟩. Dynamics

Simplified one-dimensional equation after being attacked

BD R E

d⟨x⟩ dt d⟨x⟩ dt

= −⟨x⟩2 + ⟨s⟩ ⟨x⟩ (1 − f )

d⟨x⟩ dt

= − ⟨x⟩ + ⟨s⟩ (1 − ⟨x⟩) ⟨x⟩ (1 − f )

= − ⟨x⟩ + ⟨s⟩

⟨x⟩ 1+⟨x⟩

Theory solution

⟨x⟩ = ⟨s⟩ − f ⟨s⟩

(1 − f )

⟨x⟩ = ⟨s⟩ − f ⟨s⟩ − 1 ⟨x⟩ =

⟨s⟩−f ⟨s⟩−1 ⟨s⟩−f ⟨s⟩

divide Q (⟨x⟩ , ⟨x⟩) in Eq. (2) into two parts: attacked nodes and intact nodes. We set attack intensity q = 1 for targeted attack, then the activity of each attacked node changes to zero. Setting specific values for the coefficients and exponents of the three one-dimensional equations in Table 1, we obtain simplified one-dimensional equations after a fraction f of nodes attacked. Then the theory solutions as the expression of average activity ⟨x⟩ are shown in Table 2. ⟨s⟩ is the average degree of the network after being attacked. If we attack a fraction f of nodes randomly, the new average degree will become ⟨s⟩ = ⟨k⟩ · (1 − f ) [8], where ⟨k⟩ is the average degree before attack. For SF network whose degree distribution follows P (k) ∼ k−γ , we can get [33]

∫ Kmax γ −1 ⟨km ⟩ = (γ − 1) Kmin km−γ dk Kmin ( ) γ −1 m−γ +1 m−γ +1 γ −1 = m−γ K K − K , max min +1 min

(8)

here, γ is the power exponent of degree distribution, Kmax and Kmin are the maximum and minimum node degree 1

respectively, Kmax = Kmin N γ −1 [33], N is the number of nodes in the network. If we attack hubs of fraction ∫ K fd in the network ′ ′ with HD, the maximum degree Kmax will change into Kmax , and Kmax ≤ Kmax , it can be written as K ′max P (k) dk = fd [8]. γ −1

We can get fd = Kmin

( fd =

(

)γ −1

Kmin

1−γ

1−γ

max

)

K ′ max − Kmax , then substituting Kmax , we get 1

−

K ′ max

N

.

(9)

Intuitively, if we attack the network randomly, the fraction f of nodes attacked can be written as

∫ Kmax f =

K ′ max

kP (k) dk

⟨k⟩

=

1 γ −1

⟨k⟩ 2 − γ

γ −1

Kmin

(

2−γ

2−γ Kmax − K ′ max

)

.

(10)

Substituting Kmax and solving the simultaneous equations of Eqs. (9) and (10), we can get f =

1 γ −1

⟨k⟩ 2 − γ

⎡ Kmin ⎣N

2−γ γ −1

( − fd +

1 N

−γ ) 12−γ

⎤ ⎦.

Setting m = 1 in Eq. (8) and substituting Kmax into Eq. (8), we get ⟨k⟩ =

(11)

γ −1 K 2−γ min

( N

2−γ γ −1

) − 1 . Finally, the fraction f

of nodes attacked randomly, equivalent to the condition of attacking fraction fd of hubs with HD, can be expressed as 2−γ

f =

(

N γ −1 − f d +

1 N

) 12−γ −γ .

2−γ

(12)

N γ −1 − 1 Similar with the above way, for attacking fraction fa of nodes with LD on SF networks, we get 2−γ

f =

(1 − fa ) 1−γ − 1 2−γ

.

(13)

N γ −1 − 1 Substituting Eqs. (12) and (13) into ⟨s⟩ = ⟨k⟩·(1 − f ), then substituting them into the theory solutions in Table 2 again, we can get the analytic relationship between average activity ⟨x⟩ and attacked node fraction fd and fa for SF networks. Then for ER networks, the fraction f ′ of attacking nodes intentionally( can be) considered approximately equal to the fraction f of attacking nodes randomly. Hence, by substituting ⟨s⟩ = ⟨k⟩ · 1 − f ′ into the theory solutions in Table 2, we can get the analytic relationship between average activity ⟨x⟩ and attacked node fraction f ′ for ER networks. The theory solutions of network average activity for BD, R and E are following

⎧ 2 ⎨⟨x⟩BD = ⟨k⟩ · α , 2 ⟨x⟩R = ⟨k⟩ · α − 1, 2 −1 ⎩ ⟨x⟩E = ⟨k⟨⟩·α . k⟩·α 2

(14)

F. Xu, S. Si, D. Duan et al. / Physica A 528 (2019) 121329

here, α =

2−γ (fd +1/N ) 1−γ 2−γ N γ −1 −1

−1

for SF networks with HD, α =

2−γ 2−γ N γ −1 −(1−fa ) 1−γ 2−γ N γ −1 − 1

5

for SF networks with LD, α = 1 − f ′ for ER networks.

We can see that α is only related to the attacked node fraction for ER networks, however, α is also related to network size and power exponent for SF networks. The average activity ⟨x⟩ of R is 1 less than that of BD, and the average activity ⟨x⟩ of E is the quotient of R and BD. Eq. (14) shows that average activity ⟨x⟩ is related to average degree ⟨k⟩ for both SF and ER networks, we will explore the relationship between ⟨x⟩ and ⟨k⟩ in Fig. 3. The theory solutions in Eq. (14) provide insight into the resilience of many complex systems and illustrate the important impact of dynamical processes on network resilience. Our outcomes are helpful for optimizing networks and enhancing the resilience of networks. 3. Simulation for dynamical resilience To illustrate the resilience change trend for dynamical processes, we simulate the behaviors of BD, R and E on SF and ER networks under two attack strategies (HD and LD) with four different attack intensity. SF networks are heterogeneous networks with power-law degree distribution [34], ER networks are homogeneous networks with Poisson degree distribution [35], SF and ER networks are selected as simulation objects to explore the impact of network structure on resilience. The nodes in network are attacked in strict order of their degree, targeted attacks on the networks generally start from nodes with high degree, however, nodes with low degree can also play an important role in some cases [36,37]. So we explore the impact of attack strategy on resilience with two attack strategies HD and LD. In most cases, the nodes in the network do not completely fail after being attacked, but lose a portion of their function which means they are attacked by certain intensity. Therefore, we set different attack intensity to explore their impact on resilience. The process is following: initially, the network is at a steady state with an activity vector x˙ = {x1 , x2 , . . . , xN }, then a fraction f of nodes is attacked with intensity q in the order from high to low node degree or in reverse order, the activity of attacked node changes to be xi (t) = (1 − q)xi permanently. Subsequently, the network will come to a new steady state with a new activity vector x˙ (t) = {x1 (t), x2 (t), . . . , xN (t)} following the rate equations in Table 1. We examine the average activity of all the nodes in the steady state, and use it to characterize the network resilience. The differential equation group with dynamical behaviors is simulated with Runge–Kutta method [38]. The simulation results for BD, R and E dynamics as f increases are shown in Fig. 1. We can see from Fig. 1 that the larger the attack intensity q, the greater the change of ⟨x⟩, and when q = 1, ⟨x⟩ will reduce to 0. In Fig. 1(a), when nodes are attacked with HD, ⟨x⟩ decreases as f increases, and it decreases dramatically at the beginning. In Fig. 1(b), when nodes are attacked with LD, ⟨x⟩ tends to decrease linearly as f increases for the three dynamical models. By comparing the attack effects of the two attack strategies HD and LD, we can get three influencing factors for the selection of attack strategy: (i) network structure: for SF networks, HD is better than LD, but for ER networks, the difference between the attack effects of HD and LD is small. (ii) attack intensity: when the attack intensity is small, the attack effects of HD and LD are not much different. The greater the attack intensity, the greater the difference between the attack effects of HD and LD. (iii) attacked node fraction: when the attacked node fraction f is small, HD is better than LD. As f increases, the difference between the attack effects of HD and LD becomes smaller. Simulation experiments are conducted to explore the changes of the average activity during the attack process. This can help us to understand how the resilience changes when the network is attacked, and can be used to predict the changes in network resilience for dynamical processes. We can get four factors that can affect the resilience of complex networks from Fig. 1: (i) attacked node fraction: ⟨x⟩ decreases as f increases, the thresholds of ⟨x⟩ are zero for three dynamical models BD, R and E (set h = 1 for R). If the attack intensity is large enough, ⟨x⟩ will reach its threshold when f increases to a certain value, which means the network loses its resilience. (ii) attack intensity: the larger q is, the faster ⟨x⟩ changes with f . When q is small, the network will not lose its resilience completely. (iii) network structure: when the nodes are attacked with HD, with large attack intensity, SF networks lose resilience at f = 0.5 approximately, however, this proportion for ER networks is almost 0.9. (iv) attack strategy: when the nodes are attacked with LD, the results are different from HD. The trends of ⟨x⟩ with f are approximately linear, and for both SF networks and ER networks, ⟨x⟩ reaches its threshold only when both q and f are large enough. According to the analytic relationship between average activity ⟨x⟩ and attacked node fraction fd , fa and f ′ in Eq. (14), we make the comparison of theory solutions and simulation results, shown in Fig. 2. We simulate 30 networks and calculate the mean of ⟨x⟩ to ensure the universality of simulation results. We set the negative values in the theory solutions as zero. For BD and R, the theory solutions agree with the simulation results, but for E , they are not matched exactly. The reason may be that we did some approximations in the analysis. In the process of dimension reduction, we assume the dynamical equations which describe the dynamical processes are linear and the adjacency matrix has little degree correlations, indeed, these dynamical equations are nonlinear and the networks have a rather high level of degree correlations. In addition, in the process of theoretical analysis, we assume that SF networks still obey the original power law distribution after being attacked. In fact, after some nodes in the network are removed, the network topology may change and the degree distribution of the remaining nodes no longer obeys the original distribution law. These reasons lead to the deviation between theory solutions and simulation results, especially for E .

6

F. Xu, S. Si, D. Duan et al. / Physica A 528 (2019) 121329

Fig. 1. Simulation results of average activity as f increases for BD, R and E on SF and ER networks with four attack intensity and two attack strategies (N = 1000, ⟨k⟩ = 4 for SF networks and N = 1000, ⟨k⟩ = 10 for ER networks). (a) HD -attack; (b) LD-attack.

4. Critical breakdown conditions of dynamical resilience In this section we mainly focus on the critical conditions when networks lose resilience. Average activity ⟨x⟩ as average degree ⟨k⟩ increases for different dynamics on ER networks is shown in Fig. 3 (It can be seen from Eq. (14) that when α is fixed, the relationship between ⟨x⟩ and ⟨k⟩ for SF networks is the same as ER networks). We can see that when f is fixed, network average activity ⟨x⟩ shows a growth trend as average degree ⟨k⟩ increases. ⟨x⟩ increases continually for BD and R, while ⟨x⟩ increases firstly and then tends to 1 for E . We can get a conclusion that the resilience for BD and R can be enhanced by increasing average degree ⟨k⟩, but not feasible for E whose average activity is always less than 1. With ⟨k⟩ decreasing, which is equivalent to edge attacks, the average activity ⟨x⟩ becomes smaller and it reaches the threshold at the same time for R and E , earlier than BD. By setting ⟨x⟩ = 0 in Eq. (14), we can get the critical thresholds fc of attacked node fraction at which the network loses its resilience. We can find from Eq. (14) that the thresholds for R and E are equal. BD • For ER = 1. (ii) the threshold for R and E is only related to average degree ⟨k⟩: fc R = fc E = √ networks, (i) fc 1 − ⟨1k⟩ .

• For SF networks, (i) fc BD = 1 − N1 for both HD and LD. The size of the actual network is generally large so we can assume N1 ≈ 0, then fc BD ≈ 1. (ii) for R and E , the thresholds fc are also related to degree distribution power ( ) 1−γ √ exponent γ and network size N. (iii) assuming N is infinite, we get fc R = fc E = 1 − 1/⟨k⟩ 2−γ for HD and (√ ) 1−γ fc R = fc E = 1 − 1/⟨k⟩ 2−γ for LD. √ We can see that the thresholds fc BD for both SF and ER networks are 1. The threshold fc R = fc E = 1 − ⟨1k⟩ for ER networks is smaller than structural robustness threshold 1 − ⟨1k⟩ . In resilience research, when take the dynamical processes into consideration, we can get some novel conclusions. (1) Both SF and ER networks with BD are very robust to targeted attack. (2) Compared to structural robustness without dynamical processes, ER networks with R or E are more vulnerable. For SF networks, the thresholds fc R and fc E under two attack strategies HD and LD are shown in Fig. 4.

F. Xu, S. Si, D. Duan et al. / Physica A 528 (2019) 121329

7

Fig. 2. Comparison of theory solutions and simulation results for BD, R and E on SF and ER networks under two attack strategies (N = 1000,

⟨k⟩ = 4 for SF networks with γ = 3 and N = 1000, ⟨k⟩ = 10 for ER networks). (a) HD-attack; (b) LD-attack.

Fig. 3. Analytic average activity ⟨x⟩ versus average degree ⟨k⟩ at fixed attacked node fraction f for three dynamical processes on ER networks. (a) Attacked node fraction f = 0.2; (b) Attacked node fraction f = 0.5; (c) Attacked node fraction f = 0.8.

Fig. 4. Thresholds plot of attacked node fraction in the parameter (γ , ⟨k⟩)–plane for R and E on SF networks under two attack strategies. (a) HD-attack; (b) LD-attack.

8

F. Xu, S. Si, D. Duan et al. / Physica A 528 (2019) 121329

Fig. 4 shows the pseudo color map of thresholds fc changing with power exponent γ and average degree ⟨k⟩ for two attack strategies, we can see demonstratively how thresholds change as γ and ⟨k⟩ change. In Fig. 4(a), when network is attacked with HD, the thresholds fc show growth trends as power exponent γ and average degree ⟨k⟩ increase. For structural robustness without dynamics, SF networks are vulnerable to targeted attack, the thresholds fc are small, and they increase first and then decrease to 0 as the increase of γ , see Ref. [9] for detailed explanation. Power exponent γ represents the heterogeneity of networks, so we can get a conclusion that the heterogeneity of SF networks with R or E has different effects on thresholds fc compared to structural robustness without dynamics. In Fig. 4(b), when network is attacked with LD, the thresholds fc still show a growth trend as average degree ⟨k⟩ increases and eventually fc approach 1, while show a downward trend as power exponent γ increases, which is opposite to HD attack strategy. The larger γ is, the lower the heterogeneity of the network, we can get a conclusion that the heterogeneity of SF networks has different effects on thresholds fc under different attack strategies. 5. Summary and outlook We have proposed a framework for exploring the resilience of complex networks when a fraction f of nodes is attacked, and shown that dynamical processes can affect network resilience. We solve the network resilience both analytically and experimentally. In our resilience research, we take dynamical processes into consideration and focus on the difference between network resilience and structural robustness. We obtain the critical thresholds fc of attacked node fraction at which the network loses its resilience, getting some novel conclusions: (a) The threshold for BD is always 1, both SF and ER networks with BD are very robust to targeted attack. (b) The thresholds for R and E are equal, with these two dynamics, compared to structural robustness, (i) ER networks are more vulnerable; (ii) the heterogeneity of SF networks has different effects on thresholds fc under HD attack strategy. In this paper, we give the analytic relationship between average activity ⟨x⟩ and attacked node fraction f . In real networks, the situation is not always that the entire node is removed; instead, the node is attacked by certain intensity. The future work is to extend our approach to the analytic relationship between average activity and attacked node fraction f as well as attack intensity q. In Fig. 2, the theory solutions are not matched exactly with the simulation results, especially for E . We will consider the change of network degree distribution law in the theoretical analysis process in future work. Our methods are applicable not only to the three dynamical models studied in this paper, but also to many other dynamical models. Acknowledgments Authors acknowledge support from the China Postdoctoral Science Foundation (project n. 2017M613336), of the Natural Science Foundation of Shaanxi Province (project n. 2017JM7011), and of the National Natural Science Foundation of China (projects n. 71401178, 71771186, 71631001). References [1] J. Gao, B. Barzel, A.-L. Barabási, Universal resilience patterns in complex networks, Nature 530 (7590) (2016) 307–312. [2] S. Auer, K. Kleis, P. Schultz, J. Kurths, F. Hellmann, The impact of model detail on power grid resilience measures, Eur. Phys. J. Spec. Top. 225 (3) (2016) 609–625. [3] J. Wang, L. Rong, Cascade-based attack vulnerability on the US power grid, Saf. Sci. 47 (2009) 1332–1336. [4] J. Gao, X. Liu, D. Li, S. Havlin, Recent progress on the resilience of complex networks, Energies 8 (10) (2015) 12187–12210. [5] X. Liu, H. Peng, J. Gao, Vulnerability and controllability of networks of networks, Chaos Solitons Fractals 80 (2015) 125–138. [6] T. Yuan, K. Aihara, G. Tanaka, Robustness and fragility in coupled oscillator networks under targeted attacks, Phys. Rev. E 95 (1) (2017) 012315. [7] J. Jiang, Z.G. Huang, T.P. Seager, W. Lin, C. Grebogi, A. Hastings, Y.C. Lai, Predicting tipping points in mutualistic networks through dimension reduction, Proc. Natl. Acad. Sci. USA 115 (4) (2018) 639–647. [8] R. Cohen, K. Erez, D. ben Avraham, S. Havlin, Resilience of the internet to random breakdowns, Phys. Rev. Lett. 85 (21) (2000) 4626–4628. [9] R. Cohen, K. Erez, D. ben Avraham, S. Havlin, Breakdown of the internet under intentional attack, Phys. Rev. Lett. 86 (16) (2001) 3682–3685. [10] P. Crucitti, V. Latora, M. Marchiori, A. Rapisarda, Efficiency of scale-free networks: error and attack tolerance, Physica A 320 (2003) 622–642. [11] J. Wu, M. Barahona, Y. Tan, H. Deng, Spectral measure of robustness in complex networks, IEEE Trans. Syst. Man Cybern. A 41 (6) (2011) 1244–1252. [12] G.J. Baxter, S.N. Dorogovtsev, A.V. Goltsev, J.F.F. Mendes, Avalanche collapse of interdependent networks, Phys. Rev. Lett. 109 (24) (2012) 248701. [13] D. Cellai, E. López, J. Zhou, J.P. Gleeson, G. Bianconi, Percolation in multiplex networks with overlap, Phys. Rev. E 88 (5) (2013) 052811. [14] A.E. Motter, Y.C. Lai, Cascade-based attacks on complex networks, Phys. Rev. E 66 (2) (2002) 065102. [15] J. Wang, L. Cai, B. Xu, E. Sun, Z. Zhiguo, Out of control: Fluctuation of cascading dynamics in networks, Physica A 462 (2016) 1231–1243. [16] D.E. Whitney, Dynamic theory of cascades on finite clustered random networks with a threshold rule, Phys. Rev. E 82 (2010) 066110. [17] J.P. Gleeson, D.J. Cahalane, Seed size strongly affects cascades on random networks, Phys. Rev. E 75 (2) (2007) 056103. [18] P. Singh, S. Sreenivasan, B.K. Szymanski, G. Korniss, Threshold-limited spreading in social networks with multiple initiators, Sci. Rep. 3 (7) (2013) 2330. [19] W. Wang, M. Tang, H.F. Zhang, Y.C. Lai, Dynamics of social contagions with memory of nonredundant information, Phys. Rev. E 92 (1) (2015) 012820. [20] I. Pozzana, K. Sun, N. Perra, Epidemic spreading on activity-driven networks with attractiveness, Phys. Rev. E 96 (4–1) (2017) 042310. [21] W. Cota, A.S. Mata, S.C. Ferreira, Robustness and fragility of the susceptible-infected-susceptible epidemic models on complex networks, Phys. Rev. E 98 (2018) 012310. [22] J. Gao, D. Li, S. Havlin, From a single network to a network of networks, Nat. Sci. Rev. 1 (3) (2014) 346–356. [23] J. Wang, Y. Li, Q. Zheng, Cascading load model in interdependent networks with coupled strength, Physica A 430 (2015) 242–253.

F. Xu, S. Si, D. Duan et al. / Physica A 528 (2019) 121329 [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38]

9

Y. Zhang, N. Yang, Vulnerability analysis of interdependent R & D networks under risk cascading propagation, Physica A 505 (2018) 1056–1068. B. Barzel, A.-L. Barabási, Universality in network dynamics, Nat. Phys. 9 (10) (2013) 673–681. J. Jiang, T. Zhou, The influence of time delay on epidemic spreading under limited resources, Physica A 508 (2018) 414–423. A.S. Novozhilov, G.P. Karev, E.V. Koonin, Biological applications of the theory of birth-and-death processes, Brief. Bioinform. 7 (1) (2006) 70–85. F.W. Crawford, M.A. Suchard, Transition probabilities for general birth–death processes with applications in ecology, genetics, and evolution, J. Math. Biol. 65 (3) (2012) 553–580. G. Karlebach, R. Shamir, Modelling and analysis of gene regulatory networks, Nat. Rev. 9 (10) (2008) 770–780. M. Heckera, S. Lambecka, S. Toepferb, E. Somerenc, R. Guthkea, Gene regulatory network inference: Data integration in dynamic models - a review, BioSystems 96 (2009) 86–103. P.S. Dodds, D.J. Watts, A generalized model of social and biological contagion, J. Theoret. Biol. 232 (4) (2005) 587. X. Meng, S. Zhao, T. Feng, T. Zhang, Dynamics of a novel nonlinear stochastic SIS epidemic model with double epidemic hypothesis, J. Math. Anal. Appl. 433 (2016) 227–242. A.-L. Barabási, R. Albert, Emergence of scaling in random networks, Science 286 (15) (1999) 509. A.-L. Barabási, Scale-free networks: A decade and beyond, Science 325 (412) (2009) 412–413. P.L. Krapivsky, S. Redner, Organization of growing random networks, Phys. Rev. E 63 (2) (2001) 066123. G. Tanaka, K. Morino, A. Kazuyuki, Dynamical robustness in complex networks: the crucial role of low-degree nodes, Sci. Rep. 2 (1) (2012) 232. Z. He, S. Liu, M. Zhan, Dynamical robustness analysis of weighted complex networks, Physica A 392 (18) (2013) 4181–4191. M. Riza, H. Aktöre, The Runge–Kutta method in geometric multiplicative calculus, LMS J. Comput. Math. 18 (1) (2015) 539–554.