Uncertainty Quantification of Multi-Scale Resilience in Nonlinear Complex Networks using Arbitrary Polynomial Chaos
aa r X i v : . [ n li n . AO ] O c t Uncertainty Quantification of Multi-Scale Resiliencein Nonlinear Complex Networks using ArbitraryPolynomial Chaos
Mengbang Zou, Luca Zanotti Fragonara, Weisi Guo,
Senior Member, IEEE
Abstract —Resilience characterizes a system’s ability to retainits original function when perturbations happen. In the past yearsour attention mainly focused on small-scale resilience, yet ourunderstanding of resilience in large-scale network consideringinteractions between components is limited. Even though, recentresearch in macro and micro resilience pattern has developedanalytical tools to analyze the relationship between topologyand dynamics across network scales. The effect of uncertaintyin a large-scale networked system is not clear, especially whenuncertainties cascade between connected nodes. In order toquantify resilience uncertainty across the network resolutions(macro to micro),an arbitrary polynomial chaos (aPC) expansionmethod is developed in this paper to estimate the resiliencesubject to parameter uncertainties with arbitrary distributions.For the first time and of particular importance, is our abilityto identify the probability of a node in losing its resilience andhow the different model parameters contribute to this risk. Wetest this using a generic networked bi-stable system and this willaid practitioners to both understand macro-scale behaviour andmake micro-scale interventions.
Index Terms —Uncertainty; Resilience; Arbitrary PolynomialChaos Expansion; Dynamic Complex Network
I. I
NTRODUCTION I N a connecting world, local dynamics of components innetworked systems like critical infrastructure (CI) systems,ecosystems, biological systems etc. affect each other throughinteractions among components and all together result ina more sophisticated multi-scale network wide dynamics.Example include a water distribution network couples localpumps and reservoirs to deliver supply via Navier-Stokesdynamics [1], an electric grid that uses power-flow equations,a fully loaded structure that connects beams and joints viathe Ramberg–Osgood equation, a spatially stochastic wirelessnetwork that performs traffic load balancing [2], or a fibreoptic network that connects optic switches via the NonlinearSchrodinger’s dynamic.
A. Network Resilience Modeling
Research on resilience of dynamic networks has beenwidely applied in a wide range of fields from nature to man-made network including blackout in power systems [3] to loss
M. Zou and L. Zanotti Fragonara are with Cranfield University, CranfieldMK43 0AL, U.K.W. Guo is with Cranfield University, Cranfield MK43 0AL, U.K., andalso with the Alan Turing Institute, London, NW1 2DB, U.K. (correspondingauthor e-mail: weisi.guo@cranfield.ac.uk).We acknowledge funding from EPSRC grant CoTRE - Complexity Twinfor Resilient Ecosystems (EP/R041725/1). of biodiversity in ecology [4]. Due to the different researchcontext, the concept of resilience is different in differentfields [5]. Up to now, over 70 definitions of resilience haveappeared in scientific research [6]. In this paper, resilienceis defined as the ability of a system to maintain its originalfunctionality when perturbations happen [7]. This ability is ofgreat importance in reducing risks and mitigating damages [8],[9]. Whilst the past research on resilience gives us insight intohow a few interacting components (small networks) work [9],the loss of resilience in large-scale networked systems (e.g. nodes) is difficult to predict and analyse analytically.These analytical limitations stem from theoretical gaps:most current analytical methods of resilience are appropriatefor analyzing models with a high degree of homogeneity whichenables mean field to be applied [7]. Whilst this has givenus insight into the coupling relationship between topologyand dynamics, it doesn’t enable heterogeneous prediction ofnode level dynamics. Node level is important to make criticalinterventions to specific components whilst preserving ourmulti-scale understanding of general system behaviour. Inorder to precisely identify the node-level resilience function,a sequential heterogeneous mean field estimation approach isproposed recently [10]. B. Uncertainty in Network Resilience
Mathematical models play an important role in representingphysical systems when we analyze the dynamics of networkedsystems. Actually, a lot of uncertainty exists behind the math-ematical models of many practical problems in the real world,which makes these problems more difficult and complex toanalyze. These uncertainty may stem directly from incompleteinformation of the system or measurement noise of initial dataas well as from parameters of models whose values are notknown exactly [11] [12]. In order to know the effect of ar-bitrary parameters uncertainty on the network-level resilience,our previous work introduced a polynomial chaos (PC) method[13] to understand macro-scale network wide resilience lossuncertainty. However, we still do not know the effect ofparameters uncertainty on node-level resilience, especially asparameters uncertainty may cause different effects on nodesin a network. This can paint a different picture to that of theoverall macro-scale network behaviour. That is to say, a macro-scale resilient network may hide non-resilient behaviour at themicro-level, which if not addressed in time can cause long termissues.
1) Uncertainty Modeling Review:
Since the uncertainty iswidespread in practical problems in real world and has aneffect on systems’ performance, how to quantify these uncer-tain factors is the main purpose of research on UncertaintyQuantification (UQ). UQ methods mainly includes: MonteCarlo Methods [14], Perturbation Methods [15], MomentEquation Methods [16], Polynomial approximation methods[17].Due to the high accuracy and computational efficiencycomparing with traditional Uncertainty Quantification (UQ)methods like Monte Carlo method, Polynomial Chaos Ex-pansion (PCE) method has been widely used in dynamicsystems [18]. For example, a Polynomial Chaos Expansionmethod was used to estimate the dynamic response bounds ofnonlinear system with interval uncertainty [19]. A PolynomialChaos Expansion method was applied to analysis the effect ofuncertainty in parameters on the received signal concentrationin molecular signals [20]. The PCE was initially proposed toanalyze stochastic processes based on Hermite polynomials,which were only suitable for random variables (r.v.) followingGaussian distribution [17]. However, uncertainty does notalways obey the Gaussian distribution. Whilst a normal scoretransformation could be used to solve this problem [21], butcan lead to slow convergence [22]. To solve this problem,the generalized polynomial chaos (gPC) has been developed[22] [23]. The gPC extends PCE toward a broader rangeof applications which could be used encompassing the moregeneral Gamma distribution, Beta distribution, and many otherflexible distribution functions. This is further advanced to con-sider stochastic processes represented by r.v. of any probabilitydistributions [24].The methods mentioned above need to know the detailedinformation of the involved probability density functions.While, information about distribution is usually difficult toknow or incomplete in practical applications. In differentmodels or circumstances, statistical information of parametersmay exist in many different forms. They could be discrete,continuous, or discretized continuous, even exist analyticallyin the probability density distribution (PDF) or numericallyas a histogram. There are two main reasons that limit thewidespread use of the above methods. The first reason is thatthere exist strict restrictions in most cases when they are used.The second one is that the information of problems to besolved is not always complete and perfect [25].
2) Arbitrary Polynomial Chaos:C. Novelty and Contribution
While some scientists have begun to research the networkresilience, research in estimating node-level resilience consid-ering uncertainty is rarely. In real world, ignoring uncertaintymay lead to deviation or even error when estimating resilienceof a system as well as a node. Besides, we still don notknow the role that uncertainty plays in the nonlinear dynamicssystem–whether the uncertainty will prevent or acceleratethe loss of resilience of system in macro and micro level.Therefore, under this circumstance, it is nature for us toconsider uncertainty when estimating resilience of each nodein network with nonlinear dynamics. This paper addresses the lack of uncertainty quantification inthe multi-scale resilience of complex networks with nonlineardynamics. The novelty is to enable parameter uncertainty thatfollow arbitrary distributions and estimate the resilience ofthe whole network and each node. To achieve this, a multi-dimensional arbitrary polynomial chaos (aPC) method is de-veloped in this paper to quantify uncertain factors. Besides, wecompare dynamics with certain parameters to dynamics withuncertainty when estimating the micro and macro resilience.By analyzing effects of parameters and network topology withuncertainty on the multi-scale resilience of dynamic network,we will better understand large-scale dynamic networks.II. S
YSTEM S ETUP
A. Node Level Nonlinear Dynamics and Resilience
The resilience behavior of a one-dimensional nonlineardynamic system in ecology [26], engineering [27] etc. couldbe characterized by the equation: ˙ x = f ( β, x ) (1)where f ( β, x ) shows the system’s dynamics, and the parameter β shows environment conditions (show in Figure (1)). One ofthe stable fixed points, x of equation (1) could be found by f ( β, x ) = 0 (2) dfdx (cid:12)(cid:12)(cid:12)(cid:12) x = x < (3)where f is smooth and equation (2) provides the system’ssteady state and equation (3) guarantees its linear stability. Weassume that a stable equilibrium x e > always exists whichis away from the origin. Besides, the bifurcation may occurnear to the origin shown in Figure (2). We define two differentkinds of stable equilibrium: healthy equilibrium and unhealthy equilibrium. The healthy equilibrium is far from the origin andit is a desired state of the system. The unhealthy equilibriumis near to the origin and it is an undesired one. Resilience in this general case is defined by a healthy and an unhealthyequilibrium. Only when does the healthy equilibrium exist inthe system, the system is resilient. While, if in the systemhealthy and unhealthy equilibrium exist at the same time, thesystem will transit from the desired stable fixed point to theundesired one, which indicates the loss of resilience in thesystem.
B. Network Level
Networked system often consists of a large number of com-ponents interacting with each other through the network(showin Figure (1) (b)). ˙ x i = f ( x i , a i ) + n X j =1 a ji g ( x i , x j , b ij ) (4)Node i ’s behavior is characterised by a self-dynamic f ( · ) and a coupling dynamic g ( · ) . A and B both are vectors ofparameters. A = { a , ..., a i } , B = { b , ..., b ij } . So we canrewrite equation (4): ˙ X = F ( X, A , B ) , (5) Fig. 1. It shows the dynamics of one-dimensional system and multi-dimensional system. In a one-dimensional system, resilience is characterised by thefunction x ( β ) . When β > β c , only one stable fixed point (blue) exist in the system. When β < β c , two (or more) stable fixed points will appear, whichmaps to a desirable state (blue) and an undesirable state (red). (b) In a multi-dimensional dynamic system, dynamics of system is captured by the complexweighted network w i instead of β . w i is affected by both environmental conditions and the interaction strength.Fig. 2. Red line describes a system with more than one stable equilibrium(healthy one and unhealthy one both exist). Blue line describes a system withonly one stable equilibrium. where F : R N −→ R N defined by equation (4). Let w i be theweighted in-degree of vertex v i , i.e. w i = N X j =1 a ji , (6) and w av represents the average of all weighted in-degrees. Theweighted out-degree of vertex v j is represented by w out j , i.e. w out j = N X i =1 a ji . (7)Similarly, let d i be the in-degree of v i and d out i be itsout-degree. Generally, the relationship between topology (e.g.properties of M ij ) and resilience of network is still not veryclear. One way to solve this problem is to compress the multi-dimensional dynamics to one-dimensional dynamics and mapthe overall effective dynamics of a networked system to itstopology. Indeed,a common network-level effective dynamicsmay hide different node-level dynamics of different nodes(shown in Figure (3)). In order to understand node-levelresilience and dynamics, a sequential estimation approach isproposed to solve this problem [10]. However, we still do notknow the effect of uncertainty parameters on the resilience atnode-level. III. A PPROACH AND M ETHODOLOGY
To answer this question, an arbitrary polynomial chaosexpansion method is proposed to estimate the resilience atnode-level with uncertainty. We do so by defining arbitraryuncertainty distributions on the network dynamic parameters.
Fig. 3. Similar Network Dynamics Hide Different Node Dynamics. It shows different dynamic response at node level. While the mean dynamic shows thenetwork is resilient, node 1 and node 2 have different dynamics. (a) Node 1 recoveries its resilience eventually, but node 2 loses its resilience. (b)Node 1 onlyhas one healthy equilibrium, but node 2 has a healthy equilibrium and an unhealthy equilibrium.
A. Dynamic network with uncertainty
Uncertainty may exit in self-dynamics, coupling dynamicor topology of networked dynamic system. It is assumedthat parameters in dynamics function could be represented byrandom variables. What’s more, parameters’ value are assumedto be within a range of their real values. Therefore, we get a i = a i (1 + e u i ) , b ij = b ij (1 + e v ij ) , M = M (1 + e r ) ,where u i , v ij , r are r.v. uniform in [ a, b ] and e , e , e areconstants. U = { u , ..., u i } . V = { v , ..., v ij } . The dynamicnetwork’ mathematics model with uncertainty can be writtenas: ˙ x i = f ( x i , a i (1 + e u i ))+ n X j a ji (1 + e r ) g ( x i , x j , b ij (1 + e v ij )) (8) B. Sequential Heterogeneous Mean Field Estimation
A mean field estimation is used to drive sequential substi-tution and estimation of node-level resilience.
Step 0 : Firstly, we could estimate the equilibrium e { } ofthe mean field dynamics of the networked system through ahomogeneous average degree w av = N P Ni =1 P Nj =1 a ij or aweighted average degree [7] w av = < w out w in >< w out > (9) w out = ( w out , w out , ..., w out N ) is a vector of weighted out-degrees, w in = ( w in , w in , ..., w in N ) is a vector of weightedin-degrees, and < w out > = ( N P Ni =1 w out i ) is the averageweighted out-degree. We define := 1 , ... ∈ R N . Themean field approximation of the equilibrium can be calculatedaccording to equation (10). Ξ :=
M ean [ F ( x
1, A, B )] =1 N N X i =1 ( f ( x, a i )) + 1 N N X i =1 w av g ( x, x, b ij ) (10) Ξ( x ) depends on A and B . A and B are both vectors ofr.v., for any x . Therefore, Ξ( x ) is a function depending on the random variable x . And we find x which satisfies the function Ξ( x ) = 0 .For fixed x , f ( x, a i ) is a function depending on iid r.v. a i .We set µ f ( x ) := E [ f ( x, a i )] (11) δ f ( x ) := p Var [ f ( x, a i )] (12)According to Central Limit Theorem (CLT), for big enough n , n P ni =1 f ( x, a i ) can be approximated by a normally dis-tributed random variable with mean µ f ( x ) and standard devi-ation n δ f ( x ) , i.e n n X i =1 f ( x, a i ) ∼ N ( µ f ( x ) , n δ f ( x ) ) (13) g ( x, x, b ij ) is a function depending on random variable x . Weset µ g ( x ) := E [ g ( x, x, b ij )] (14) δ g ( x ) := q Var [ g ( x, x, b ij )] (15)Then we can get n n X i,j =1 M ji g ( x, x, b ij ) ∼ N ( mn µ g ( x ) , mn δ g ( x ) ) (16)From the mentioned above, we know that Ξ( x ) is the sumof 2 normally distributed r.v.. Then we can get Ξ( x ) ∼ N ( µ f ( x ) + mn µ g ( x ) , n δ f ( x ) + mn δ g ( x ) ) (17) Ξ α ( x ) could be realised by drawing ζ α from N (0 , and wecan get Ξ α ( x ) = µ f ( x ) + mn µ g ( x ) + r n δ f ( x ) + mn δ g ( x ) ζ α (18)It is assumed that every realisation of Ξ( x ) has the shapedescribed in Figure (2). The equilibrium e { } could be calcu-lated from equation (18). Since ζ α is a random variable whichis normally distributed, a polynomial chaos expansion (PCE)could be used to approximate the equilibrium e { } . P ( e { } ) . (a) Resilience Bounds and Uncertainty Region with certain parameters(b) Critical Resilience Value at Node Level with certainty parametersFig. 4. Critical Resilience Value Identifies Vulnerable Nodes with Certain Parameters. (a) Resilience Bounds shows the Upper-Bound and Lower-Bound ofequilibrium when links removed. In this figure, it explicitly predicts when the loss of resilience will happen. (b) Critical Resilience shows the relationshipbetween average weight value of network and critical weight value. When w i > w crit , the node is resilient, otherwise it is not. We calculate the smallest positive root ρ { } of Ξ ′ ( x ) and set τ { } = Ξ( ρ { } ) .Since Ξ( x ) is a random variable, ρ { } and τ { } are functionsbased on this random variable. Meanwhile, τ { } is an indicatorfor the saddle-node bifurcation, which means that whether thesystem is resilient could be judged by τ { } . To be specific, fora given realization of ζ α , if τ α > , then only one equilibrium(the healthy one) exists in the system and the system isresilient. Otherwise, if τ α < , three equilibrium exist in thesystem including healthy and unhealthy equilibrium. Then thedynamics is non-resilient. P ( τ > represents the probabilityof the system being resilient. We use PCE to estimate τ ( ζ ) .This PCE is represented by e τ n ( ζ ) and we set pos ( x ) = (cid:26) if x > otherwise (19)Then, the probability of resilience could be calculated √ π + ∞ Z Z −∞ pos ( e τ n ( ζ )) dζ (20) Step 1 : The mean field approximation is used as an initialguess to estimate the dynamics of each node: ˙ x i = f ( x i , a i ) + w i g ( x i , e { } , b ij ) = 0 (21)The solution of equation (21) is a function of w i , i.e. χ { } ( w i ) .Then we have first order approximation e { } i = χ { } ( w i ) .Since parameters a i , b ij does not always belong to commondistribution like Gaussian distribution, Binomial distributionetc. We need to use the arbitrary polynomial chaos (aPC) [25]to approximate e { } i and its distribution. Step 2 : The previous approximation could be used toestimate the effect of the graph on a vertex. The effect ofan in-edge on the dynamics of vertex i is g ( x i , x j ) and theprobability of of a vertex j is on the other side of the in-edge is proportional to its out-degree. So the average effectis P Nj =1 d out j g ( x i , x j , b ij ) / P Nj =1 d out j . To approximate meaneffect of the neighbours, components in g ( · ) are weightedby d out . Therefore, the previous step’s estimation could beused to make further estimation. And we can approximate the (a) The effect of uncertain parameters on network when average weightis 7. U1 −1.00−0.75−0.50−0.250.000.250.500.751.00U2−1.00−0.75−0.50−0.250.000.250.500.751.00e 10.610.811.011.211.411.611.8 equilibrium of a node ith uncertain parameters U ,U (0, 0, 11.13) equilibrium of a node hen uncertain parameters are zero (b) The effect of uncertain parameters on a node’s equilibriumFig. 5. The effect of uncertain parameters on network and node. In system with certain parameters, network and node-level dynamics are certain. Whilein a system with uncertain parameters, network and node-level dynamics are uncertain. In (b) e represents the equilibrium of a node, U , U are uncertainparameters. equilibrium of each node from ˙ x = f ( x i ) + w i P Nj =1 d out j g ( x i , e { } j ) P Nj =1 d out j = 0 (22)Notice that the solution depends on w i . So the second orderapproximation is e { } i = χ { } ( w i ) . Also, e { } and its distri-bution could be approximated by aPC. Step 3 to n : We repeat the above steps with each approxi-mation calculated in the previous step.
C. Arbitrary Polynomial Chaos Expansion1) One-Dimensional aPC: Ξ is defined as a random vari-able with PDF w . Let us consider a stochastic model X = φ (Ξ) . φ is a function that is square integrable on R with aweight function w . For a stochastic analysis of X , the model φ (Ξ) may be expanded as follows: φ (Ξ) = d X i =0 c i P ( i ) (Ξ) (23) c i are expansion coefficients and P ( i ) (Ξ) are the polynomialswhich forms the basis (cid:8) P (0) , P (1) , P (2) , ..., P ( i ) (cid:9) . P ( i ) (Ξ) isorthogonal with respect to w . In aPC, w may have an arbitraryform, which is different general PCE methods. The basis (cid:8) P (0) , P (1) , P (2) , ..., P ( i ) (cid:9) need to be formed specificallyaccording to the statistics information of w .
2) Multi-Dimensional aPC:
In some cases, the number ofinput parameters is more than one, i.e.
Ξ = { Ξ , Ξ , ..., Ξ N } .The model output X could be approximated by a multivariatepolynomial expansion: φ (Ξ , Ξ , ..., Ξ N ) = M X i =1 c i Φ i (Ξ , Ξ , ..., Ξ N ) . (24)In equation (24), the number of input parameters is N andthe number of M is decided by the formula M = ( N + d )! / ( N ! d !) , where d represents the order of expansion. Here, we need to construct the orthogonal polynomial basis Φ i for Ξ , Ξ , ..., Ξ N . Assuming that the input parameters areindependent, the multi-dimensional basis can be constructed asa simple product of the corresponding univariate polynomials Φ i (Ξ , Ξ , ..., Ξ N ) = N Y j =1 P ( α ij ) j (Ξ , Ξ , ..., Ξ N ) , N X j =1 α ij ≤ M, i = 1 , ..., N, (25)In equation (25), α ij is a multivariate indicator with theinformation on how to list all possible products of individualunivariate basis functions. α ij is a multivariate index containinginformation about how to enumerate all possible products ofindividual univariate basis functions. Here, we give an exampleto illustrate it. For example, we assume that d = 2 , N = 2 ,then M = 6 . Therefore, according to equation (25), the basisfunctions are (cid:8) , Ξ , Ξ , Ξ Ξ , Ξ , Ξ (cid:9) .We define the polynomial P ( k ) (Ξ) of degree k in therandom variable Ξ : P ( k ) (Ξ) = k X i =0 P ( k ) i Ξ i , k ∈ [0 , d ] (26)where P ( k ) i are coefficients in P ( k ) (Ξ) .The key of the aPC method is to construct the polynomialsin equation (26) to form an orthonormal basis for arbitrarydistributions which could be discrete, continuous, raw datasets or by their moments. we define the orthonormality forpolynomials P ( k ) and P ( l ) as Z P ( k ) (Ξ) P ( l ) (Ξ) dw (Ξ) = (cid:26) ∀ k = l else (27)Here we assume that the leading coefficients of all polyno-mials: P ( k ) k = 1 ∀ k . The k th raw (crude) moment of therandom variable Ξ is defined as µ k = Z Ξ k dw (Ξ) (28) Fig. 6. Approximate resilience of system by Polynomial Chaos Expansion. We truncate the series to arbitrary order N from 2 to 5. (a) approximate theminimum value of system by PCE. (b) approximate the probability of resilience. It is clear that there is a significant difference in results between N = 2 and N = 3 , , in (a) and (b). The relationship between raw moments of Ξ and theircoefficients can be written in matrix form (the detail processcould be seen in [25]): µ µ ... µ k µ µ ... µ k +1 ... ... ... ... µ k − µ k ... µ k − ... P ( k )0 P ( k )1 ... P ( k ) k − P ( k ) k = ... (29)For multi-dimensional r.v., the polynomial P ( k ) j (Ξ j ) is definedas: P ( k ) j (Ξ j ) = k X i =0 P ( k ) i,j Ξ ij (30)and the unknown polynomial coefficients P ( k ) i,j can be definedfrom the following matrix equation [28]: µ ,j µ ,j ... µ k,j µ ,j µ ,j ... µ k +1 ,j ... ... ... ... µ k − ,j µ k,j ... µ k − ,j ... P ( k )0 ,j P ( k )1 ,j ... P ( k ) k − ,j P ( k ) k,j = ... (31) We now show the results of a real system case study toillustrate how the aPC framework can be used.IV. R ESULTS FOR B I -S TABLE S YSTEMS
Bi-stable dynamical systems are common across social (e.g.population logistic model [29]), ecological (e.g. soil health[30]), climate (e.g. ocean circulation [31]), and human conflictsystems [32]. There exists a stable undesirable state (e.g.population collapse or conflict) and a stable desirable state(e.g. healthy population with collaboration [33]), with anunstable transition brink in between, and this is ideal fordemonstrating the concept of resilience and uncertainty. Net-works that connect such systems represent a wider interactingecosystem and often a mutualistic coupling represents positivereinforcing interactions. Interaction examples include gravity,radiation, or Boltzmann Lotka Volterra (BLV) models [34]frequently use a x i × x j mutualistic attractor component. A. Case Study: Ecological Network
A case of pollinator networks is used to illustrate thedynamics of networked system at micro level and macro level[35]. x i represents the abundance of species i , which is givenby: dx i dt = B i + x i (1 − x i K i )( x i C i −
1) + N X j a ji x i x j D i + E i x i + H j x j (32) B i represents the incoming migration rate of species i fromother ecosystems. The second term on right hand showslogistic growth with carrying capacity K i of the system, andthe Allee effect (low abundance ( x i < C i ) causes negativegrowth) [36]. The third term is a coupling function whichsaturates for large x i or x j ( j ’s positive contribution to x i isbounded).For simplicity, we use homogeneous parameters: B =0 . , C = 1 , K = 5 , D = 5 , E = 0 . , H = 0 . . Besides, itis assumed that some parameters’ value has to be within of its mean. Here, we set C = E [ C ](1+0 . U ) , E = E [ E ](1+0 . U ) , where U , U are random variables uniform in [ − , ( U , U could be r.v. that follow arbitrary distributions. Weonly need to know the statistics information of their rawmoments.). The definition of system resilience in this modelis the ability of the system to recover species abundance fromextinction [10]. To achieve this, the system should be in theregime where only one equilibrium exists. This because if overone equilibrium exist in the system, the system will be trappedin the state with low abundance, which means that the systemcan not recover its species abundance and loses its resilience.In Figure (4), we show what happened when a networkbecomes less connected by removing edges. In this case,parameters are certain and the figure explicitly shows thebounds of equilibrium under different perturbation and theregime where loss of resilience happens. Critical functiondescribes resilience regimes which maps macro (network-level) properties (average weighted degree w av to micro (local-level) properties (critical resilience value w crit )). For each w av ,corresponding w crit could be calculated from equation (33).The critical weight, w crit , is a function of w av since it is afunction of e { } and e { } is a function of w av . In Figure(4) (b), we see the graph of w crit versus w av . Since e { } isdiscontinuous, w crit is also discontinuous. ˙ x i = f ( x i ) + w i g ( x i , e { } ( w av )) (33)In this case, a critical average weight w ∗ is about 7 wherebifurcation will happen. When average weight is greater than7, the system is resilient and almost every node in thissystem is resilient. The critical weight can reveal some basicproperties for the dynamics on a nodal level. For example wesee in Figure (4) (b) that when when w av > w ∗ , w crit is almost0. This implies that if the system on average is in the resilientregion, a vertex will also be in the resilient region even if it isvery weakly connected to the rest of the network. However, inthe case with uncertain parameters, even if the average weightis greater than 7, the system is possibly not resilient. We useaPC to analysis what will happen in the regime where loss ofresilience may exist. In Figure (5), it shows dynamics of thesystem with uncertain parameters when average weight is 7.We use aPC to approximate the minimum value of the functionand whether the system is resilient. B. Analysis on the Effect of Uncertainty
Firstly, we use the method described in
Step 0 to analysethe probability of system to be resilient and approximate theequilibrium. We truncate the series to arbitrary orders N from 2 to 5 (Figure (6)). It is clear that the convergence of thefunction can be improved with the increase of the polynomialorder ( N ). However, with the increasing of the order, muchmore simulation is needed. Therefore, we have to make acompromise between accuracy and computational efficience.In Figure (6), it clearly shows the difference among differentorders especially N = 2 . To calculate the probability ofresilience, a graph of Cumulative Distribution Function (CDF)with different truncation is shown in Figure (6) (b). In Figure(6), when N = 3 , , , the results are almost the same.However, there is an obvious difference for N = 2 .Considering the accuracy and computational efficiency, wechoose N = 3 for the polynomial order. So, we can seethe effect of uncertain parameters on system resilience aswell as node-level resilience. When parameters are certainand average weight is 7, the system is resilient and nodesare resilient. However, when parameters are uncertain in thiscase, the probability of resilience of the system is about 0.798.So according to analysis above, some nodes will also possiblylose resilience.Second, we use the method in Step 1 and
Step 2 to estimatethe equilibrium of each node. The method aPC mentionedabove is used to estimate a node’s equilibrium and we truncatethe series to arbitrary orders N from 2 to 5 shown in Figure(7). In Figure (7) (a), it shows that the node has differentequilibrium when parameters U , U have different values andthe results for N = 2 , , , almost overlap. Meanwhile, it canbe seen that the results of CDF also almost overlap. Therefore,we consider N = 2 as the appropriate choice.In Figure (8), we show the effect of uncertain parameters onthe resilience of whole network and each node. In Figure (8),it is clear that when parameters are certain, network couldmaintain its resilience when average weight is greater than7. However, with the effect of uncertain parameters, networkcould lose resilience even though its average weight is greatthan 7. With the growth of average weight, the network hasmore chance to be resilient. When the average weight isgreater than a critical value, the network is absolute resilient.Similarly, in Figure (8) (d) red part shows that when node’sweight is greater than a critical value under certain averageweight, the node could maintain its resilience. While, withthe effect of uncertainty represented by blue part, a node maylose resilience even though its weight is greater than the criticalvalue in Figure (4) (b). Therefore, the method mentioned abovecould help us understand the effect of uncertainty on network-level and node-level resilience. Also, it help us to predictwhether a node is resilient and the probability of a node tolose resilience. V. C ONCLUSIONS
At present, the research of how to estimate node-levelresilience of dynamic networked system is still limited. Nodelevel is important to make critical interventions to specificcomponents whilst preserving our multi-scale understandingof general system behaviour. In this paper, an arbitrary poly-nomial chaos expansion (aPC) method is used to quantify theuncertainty of arbitrary uncertain distributions. This approach (a) Approximate equilibrium of a node by aPC when we truncate theseries to arbitrary orders N from 2 to 5 (b) The CDF of equilibrium when we truncate the series to arbitraryorders N from 2 to 5Fig. 7. Approximate equilibrium of a node in the networked system by aPC. In (a), the four color surface , blue, green, red, yellow surface, present differenttruncation from 2 to 5. e represents the equilibrium of a node, U , U are uncertain parameters. It can be seen that these surface almost overlap which meansthat their accuracy are similar. In (b), it shows the CFD when we truncate the series from 2 to 5.(a) Probability of resilience when average weight of network isdifferent (b) Critical weight of node with different average weight of network av (c) Difference between certain and uncertain parameters when averageweight of network changes (d) Relationship between critical weight and average weight withcertain and uncertain parametersFig. 8. It shows the effect of uncertainty parameters on resilience of network and each node. (a) (b) show the probability of resilience at network-level andnode-level in system with uncertainty. (c) (d) show the difference between resilience with certain parameters and uncertain parameters at network-level andnode-level (Blue represents system with uncertain parameters and red represents system with certain parameters). can effectively estimate node-level resilience and analyse theeffect of uncertainty on each node. This could help us bettermake a prediction of the probability that a node loses itsresilience and reduce the risk of uncertainty. In the future, wewould like to survey how the community structure of networkaffects network-level and node-level resilience, for example,whether there exists a relationship between modularity ofcommunity in network and resilience.R EFERENCES[1] Z. Wei, A. Pagani, G. Fu, I. Guymer, W. Chen, J. McCann, andW. Guo, “Optimal sampling of water distribution network dynamicsusing graph fourier transform,”
IEEE Transactions on Network Scienceand Engineering , vol. 7, no. 3, pp. 1570–1582, 2020.[2] G. Moutsinas and W. Guo, “Probabilistic stability of traffic load balanc-ing on wireless complex networks,”
IEEE Systems Journal , 2019.[3] R. Arghandeh, A. Von Meier, L. Mehrmanesh, and L. Mili, “On thedefinition of cyber-physical resilience in power systems,”
Renewableand Sustainable Energy Reviews , vol. 58, pp. 1060–1069, 2016.[4] C. N. Kaiser-Bunbury, J. Mougal, A. E. Whittington, T. Valentin,R. Gabriel, J. M. Olesen, and N. Bl¨uthgen, “Ecosystem restorationstrengthens pollination network resilience and function,”
Nature , vol.542, no. 7640, pp. 223–227, 2017.[5] X. Liu, D. Li, M. Ma, B. K. Szymanski, H. E. Stanley, and J. Gao,“Network resilience,” arXiv preprint arXiv:2007.14464 , 2020.[6] L. Fisher, “More than 70 ways to show resilience,”
Nature , vol. 518,no. 7537, pp. 35–35, 2015.[7] J. Gao, B. Barzel, and A.-L. Barab´asi, “Universal resilience patterns incomplex networks,”
Nature , vol. 530, no. 7590, pp. 307–312, 2016.[8] R. Cohen, K. Erez, D. Ben-Avraham, and S. Havlin, “Resilience of theinternet to random breakdowns,”
Physical review letters , vol. 85, no. 21,p. 4626, 2000.[9] R. V. Sole and M. Montoya, “Complexity and fragility in ecologicalnetworks,”
Proceedings of the Royal Society of London. Series B:Biological Sciences , vol. 268, no. 1480, pp. 2039–2045, 2001.[10] G. Moutsinas and W. Guo, “Node-level resilience loss in dynamiccomplex networks,”
Nature Scientific Reports , 2020.[11] R. Pulch, P. Putek, E. J. W. ter Maten, and W. Schoenmaker, “Uncer-tainty quantification: introduction and implementations,” in
Nanoelec-tronic Coupled Problems Solutions . Springer, 2019, pp. 197–221.[12] Z. M. Prince and J. C. Ragusa, “Parametric uncertainty quantificationusing proper generalized decomposition applied to neutron diffusion,”
International Journal for Numerical Methods in Engineering , vol. 119,no. 9, pp. 899–921, 2019.[13] G. Moutsinas, M. Zou, and W. Guo, “Uncertainty of resiliencein complex networks with nonlinear dynamics,” arXiv preprintarXiv:2004.13198 , 2020.[14] G. Fishman,
Monte Carlo: concepts, algorithms, and applications .Springer Science & Business Media, 2013.[15] C. Zhao, S. Xie, X. Chen, M. P. Jensen, and M. Dunn, “Quantifyinguncertainties of cloud microphysical property retrievals with a perturba-tion method,”
Journal of Geophysical Research: Atmospheres , vol. 119,no. 9, pp. 5375–5385, 2014.[16] D. Zhang,
Stochastic methods for flow in porous media: coping withuncertainties . Elsevier, 2001.[17] N. Wiener, “The homogeneous chaos,”
American Journal of Mathemat-ics , vol. 60, no. 4, pp. 897–936, 1938.[18] D. Shen, H. Wu, B. Xia, and D. Gan, “Polynomial chaos expansion forparametric problems in engineering systems: A review,”
IEEE SystemsJournal , 2020.[19] L. Wang, Z. Chen, and G. Yang, “A polynomial chaos expansionapproach for nonlinear dynamic systems with interval uncertainty,”
Nonlinear Dynamics , pp. 1–20, 2020.[20] M. Abbaszadeh, G. Moutsinas, P. J. Thomas, and W. Guo, “Uncertaintyquantification in molecular signals using polynomial chaos expansion,”
IEEE Transactions on Molecular, Biological and Multi-Scale Commu-nications , vol. 4, no. 4, pp. 248–256, 2018.[21] H. Wackernagel,
Multivariate geostatistics: an introduction with appli-cations . Springer Science & Business Media, 2013.[22] D. Xiu and G. E. Karniadakis, “The wiener–askey polynomial chaos forstochastic differential equations,”
SIAM journal on scientific computing ,vol. 24, no. 2, pp. 619–644, 2002. [23] D. Xiu and Karniadakis, “Modeling uncertainty in flow simulations viageneralized polynomial chaos,”
Journal of computational physics , vol.187, no. 1, pp. 137–167, 2003.[24] X. Wan and G. E. Karniadakis, “Multi-element generalized polynomialchaos for arbitrary probability measures,”
SIAM Journal on ScientificComputing , vol. 28, no. 3, pp. 901–928, 2006.[25] S. Oladyshkin and W. Nowak, “Data-driven uncertainty quantificationusing the arbitrary polynomial chaos expansion,”
Reliability Engineering& System Safety , vol. 106, pp. 179–190, 2012.[26] R. M. May, “Thresholds and breakpoints in ecosystems with a multi-plicity of stable states,”
Nature , vol. 269, no. 5628, pp. 471–477, 1977.[27] A. M. Lyapunov, “The general problem of the stability of motion,”
International journal of control , vol. 55, no. 3, pp. 531–534, 1992.[28] S. Oladyshkin, H. Class, R. Helmig, and W. Nowak, “A concept for data-driven uncertainty quantification and its application to carbon dioxidestorage in geological formations,”
Advances in Water Resources , vol. 34,no. 11, pp. 1508–1518, 2011.[29] M. Lundqvist, A. Compte, and A. Lansner, “Bistable, irregular firing andpopulation oscillations in a modular attractor memory network,”
PLoSComput Biol , vol. 6, no. 6, 2010.[30] L. Todman, F. Fraser, R. Corstanje, and et al., “Evidence for functionalstate transitions in intensively-managed soil ecosystems,”
Sci Rep , vol. 8,2018.[31] R. Marsh, A. Yool, T. M. Lenton, M. Y. Gulamali, N. R. Edwards, J. G.Shepherd, M. Krznaric, S. Newhouse, and S. J. Cox, “Bistability of thethermohaline circulation identified through comprehensive 2-parametersweeps of an efficient climate model,”
Climate Dynamics , vol. 23, 2004.[32] G. Aquino, W. Guo, and A. Wilson, “Nonlinear dynamic models of con-flict via multiplexed interaction networks,”
Preprint arXiv 1909.12457 ,2019.[33] J. Ron, I. Pinkoviezky, E. Fonio, O. Feinerman, and N. Gov, “Bi-stabilityin cooperative transport by ants in the presence of obstacles,”
PLOSBiology , 2018.[34] A. Wilson, “Boltzmann, lotka and volterra and spatial structural evo-lution: an integrated methodology for some dynamical systems,”
J. R.Soc. Interface , 2008.[35] J. N. Holland, D. L. DeAngelis, and J. L. Bronstein, “Populationdynamics and mutualism: functional responses of benefits and costs,”
The American Naturalist , vol. 159, no. 3, pp. 231–244, 2002.[36] W. C. Allee, O. Park, A. E. Emerson, T. Park, K. P. Schmidt et al.et al.