Cascading Failures in Interdependent Systems: Impact of Degree Variability and Dependence

TL;DR

This paper investigates how degree variability and dependence among nodes in interdependent networks influence the likelihood and extent of cascading failures, revealing that higher variability and positive degree correlations enhance system robustness.

Contribution

It introduces a dependence graph model capturing degree variability and dependence, and analyzes their effects on cascading failure propagation in interdependent systems.

Findings

Higher degree variability reduces cascade size.

Positive degree correlations increase system robustness.

Variability and dependence properties significantly affect failure dynamics.

Abstract

We study cascading failures in a system comprising interdependent networks/systems, in which nodes rely on other nodes both in the same system and in other systems to perform their function. The (inter-)dependence among nodes is modeled using a dependence graph, where the degree vector of a node determines the number of other nodes it can potentially cause to fail in each system through aforementioned dependency. In particular, we examine the impact of the variability and dependence properties of node degrees on the probability of cascading failures. We show that larger variability in node degrees hampers widespread failures in the system, starting with random failures. Similarly, positive correlations in node degrees make it harder to set off an epidemic of failures, thereby rendering the system more robust against random failures.

Tables

Table 1. TABLE I: Degree distributions p 1 ( ℓ ) subscript superscript 𝑝 ℓ 1 p^{(\ell)}_{1} , ℓ = 1 , 2 , 3 ℓ 1 2 3 \ell=1,2,3 .

	$d_1$
$d_2$	0 1 2 3 0 0.40000 0.12000 0.16000 0.01200 1 0.03750 0.01125 0.01500 0.01125 2 0.05000 0.01500 0.02000 0.01500 3 0.01250 0.00375 0.00500 0.00375

Equations

$$q_{i,i}$$

$$h_i(\mathbf{o})$$

$$=$$

$${\bar{d}}_j(\mathbf{o})$$

$$h_{i^{+}}(\mathbf{o})$$

$$=$$

$$\frac{w_{i,i}(d){\wp }_i^{in}(d,1)}{p_{i,i}(d)}\le \frac{w_{i,i}(d+1){\wp }_i^{in}(d+1,1)}{p_{i,i}(d+1)}\text{mbox}foralld\in \mathrm{I}\mathrm{N}.$$

$$\text{mathbbm}M$$

$$f_i(\mathbf{s})=\sum _{\mathbf{o}\in {\text{mathbbm}Z}_{+}^{2N_S}}{h}_i(\mathbf{o})\prod _{j\in {\mathcal{S}}^{+}}{s}_j^{o_j}.$$

$$\mathbf{f}(\mathbf{\mu })$$

$$\text{mathbbm}M$$

$$\xi (\mathbf{x})+\xi (\mathbf{y})\le \xi (\mathbf{x}\vee \mathbf{y})+\xi (\mathbf{x}\wedge \mathbf{y}),$$

$$C_{\mathbf{X}}(F_{X_1}(x_1),\dots ,F_{X_n}(x_n))$$

$$\frac{{\partial }^2}{\partial x_i\partial x_j}{\xi }_{\mathbf{s}}(\mathbf{x})$$

$$f_1(\mathbf{s})$$

$$f_3(\mathbf{s})$$

$${\mu }_1$$

$$\left(\sum _{d_1\in {\text{mathbbm}Z}_{+}}{p}_{1,1}^{(2)}(d_1){\mu }_3^{d_1}\right)\left(\sum _{d_2\in {\text{mathbbm}Z}_{+}}{p}_{1,2}^{(2)}(d_2){\mu }_1^{d_2}\right)$$

$$\left(\sum _{d_1\in {\text{mathbbm}Z}_{+}}{p}_{1,1}^{(2)}(d_1+1){\mu }_3^{d_1}\right)\left(\sum _{d_2\in {\text{mathbbm}Z}_{+}}{p}_{1,2}^{(2)}(d_2){\mu }_1^{d_2}\right)$$

$${\mathbf{f}}^2({\mathbf{\mu }}^1)$$

$$f_i^1({\mathbf{\mu }}^1)$$

$$\frac{{\partial }^2}{\partial x_i\partial x_j}{\psi }_{\mathbf{s}}(\mathbf{x})=\log (s_i^{-1})\log (s_j^{-1}){\psi }_{\mathbf{s}}(\mathbf{x})\ge 0$$

$$f_i^2({\mathbf{\mu }}^1)$$

$$\displaystyle{\mathcal{O}}_{j}:=\big{\{}{\bf o}\in{\mathbbm{Z}}^{2N_{S}}_{+}\ |\ o_{j^{\prime}}=0\mbox{ for all }j^{\prime}\notin{\mathcal{S}}^{+}_{j}\big{\}}.$$

$$\text{mathbbm}E\left[g({\mathbf{O}}_i^{(\ell )})\right]$$

$$\text{mathbbm}E\left[g({\mathbf{O}}_i^{(1)})\right]$$

$$\text{mathbbm}E\left[g({\mathbf{O}}_i^{(\ell )})\right]=\sum _{\mathbf{d}\in {\text{mathbbm}Z}_{+}^{N_S}}{p}_i^{(\ell )}(\mathbf{d})$$

$$\displaystyle\times\Big{[}\underbrace{\sum_{{\bf o}\in{\mathcal{O}}_{i}:\bar{{\bf d}}({\bf o})\leq{\bf d}}g({\bf o})\prod_{j\in{\mathcal{S}}}{{d_{j}}\choose{\bar{d}_{j}({\bf o})}}\left(q^{\ell}_{i,j}\right)^{\bar{d}_{j}({\bf o})}\left(\bar{q}^{\ell}_{i,j}\right)^{d_{j}-\bar{d}_{j}({\bf o})}}_{\text{ \large{$=:g^{*}_{\ell}({\bf d})$}}}\Big{]}.$$

$$\text{mathbbm}E\left[g({\mathbf{O}}_i^{(1)})\right]=\text{mathbbm}E\left[g_1^{\ast }({\mathbf{D}}_i^{(1)})\right]=\sum _{\mathbf{d}\in {\text{mathbbm}Z}_{+}^{N_S}}{p}_i^{(1)}(\mathbf{d})g_1^{\ast }(\mathbf{d})$$

$$(\text{ref}eq:appen2-0)$$

$$g_1^{\ast }(\mathbf{d})=\text{mathbbm}E\left[g(\mathbf{Y}(\mathbf{d},1))\right]\le \text{mathbbm}E\left[g(\mathbf{Y}(\mathbf{d},2))\right]=g_2^{\ast }(\mathbf{d})$$

$$q_{i,i}^1$$

Full text

1 Introduction

Many systems providing critical services to modern societies (e.g., smart grids, manufacturing systems, transportation systems) comprise multiple heterogeneous systems that support each other to enable the functionality of the overall system. In particular, (local) decision makers or subsystems belonging to different constituent or component systems (CSes) rely on each other to perform various functions. For instance, a modern power system not only includes an electrical grid/network, but also depends on an information and communication network (ICN) to monitor the state of the electrical network and to communicate and execute appropriate control actions based on the observed state.

Throughout the paper, we refer to the (local) decision makers or subsystems in various CSes simply as agents. Intricate (inter-)dependence among agents in CSes makes the analysis of these complex systems challenging. Moreover, in some cases, (random or targeted) failures of a small number of agents in one CS have potential to cause unexpected, widespread failures of many agents in multiple CSes.

The 2003 blackout in Italy provides a good example [38]. The onset of the blackout was triggered by an initial failure in the power grid, which caused a disruption to both the power grid and the ICN used for communication between power substations. This secondary failure in ICN further hampered the communication between stations and their regulation, setting off rapid cascading failures throughout a large part of the power grid.

As illustrated by this example, due to increasing reliance of modern societies on such complex systems and interdependence among CSes, there is a growing interest in modeling and understanding the interaction between (agents in) interdependent CSes and the robustness of the overall systems (e.g., [1, 2, 7, 12, 20, 38, 41, 42, 43, 46]). Yet, there is no theory that allows us to answer even a basic question, “Given two different interdependent networks or systems, when can we say that one network or system is more robust than the other?”

The overarching goal of our study, which complements those of existing studies (summarized in Section 2), is to contribute to the emerging theory on complex systems, in particular on the influence of the dependence structure properties between agents on the robustness of the systems with respect to localized, random failures in CSes. Our hope is that the findings will help engineers and researchers identify critical properties of robust systems and incorporate them into design guidelines of complex systems.

To this end, we develop a general model for capturing the propagation of failures from one agent to another both within individual CSes and across multiple CSes. This model is similar to that of [21] and allows us to introduce asymmetric dependence among agents belonging to heterogeneous systems (e.g., electrical network vs. ICN) and to study different ways in which failures can proliferate through diverse CSes.

Some key questions we are interested in are: (a) When is it possible for a localized initial failure in one CS, beginning with one or a small number of randomly chosen agents, to trigger a cascade of failures not only within the CS in which the initial failure originated, but also in other CSes? (b) How does the underlying dependence structure among the agents in various CSes influence the dynamics of failure propagation and the likelihood of such cascading failures? (c) How can we identify susceptible CSes that are more likely to set off widespread failures across many CSes, starting with a few initial failures in the CSes? In this study, we aim to offer partial answers to these important questions.

While we carry out the study in the framework of propagating failures in interdependent systems, we suspect that the basic model and approach as well as some of key findings can be extended to other applications with appropriate changes. These applications include (i) information or rumor propagation or new technology adoption via multiple social networks, (ii) an epidemic of disease across multiples geographic locations (e.g., cities or countries), and (iii) spread of malware in the Internet.

1.1 Summary of main results

We model the (inter-)dependence among agents with the help of what is known in the literature as a degree-based model or Chung-Lu model [10]. A similar model is used in many existing studies, e.g., [6, 7, 22, 44, 46]. In order to capture different manners in which failures can spread both within various CSes and between CSes, we model the dynamics of failure propagation using a multi-type branching process with the assumption that there are no cycles in the set of failed agents in a local neighborhood.

We call other agents which an agent can influence and cause to fail its (dependence) neighbors, and model the number of neighbors that an agent has in each CS, using a random degree vector. The $j$-th element in the degree vector of an agent belonging to the $i$-th CS ($i\neq j$) is the number of agents in the $j$-th CS which are dependent on the CS $i$ agent.

Given fixed distributions of degree vectors for agents in different CSes, the tools from multi-type branching process theory are employed to estimate the probability that a random initial failure in the $i$-th CS will give rise to cascading failures, affecting a large number of agents across the system. Since this probability of suffering an epidemic of failures depends on which CS suffers the initial failure, it also tells us which CSes are more vulnerable to random failures than others.

The primary goal of our study is to investigate how (i) the variability of agents’ degrees (in the aforementioned random degree vectors) and (ii) the dependence structure of the degrees influence the likelihood of a random failure in a CS sparking a chain of failures throughout the system with many agents. To achieve this goal, we adopt well-known (integral) stochastic orders that are partial orders on the set of degree distributions [32, 40]. They allow us to compare degree distributions of interest to us with respect to degree variability and (positive) dependence.

There are many other important properties, such as assortativity [34, 35] and clustering [11, 46] often observed in social networks or engineered systems as well as community structure [29], which influence the dynamics of information or failure propagation. But, as it will be clear, even without these properties, analyzing the role of degree variability and dependence is technically challenging. For this reason, we do not model them here and leave an investigation of their effects for a future study. For instance, clustering is shown to impede global cascades in multiplex networks [46], and thus the probability of cascades we estimate here may serve as an upper bound when there is clustering. However, we suspect that the qualitative findings reported in this study will continue to hold even in the presence of clustering.

The high-level messages of our analytical findings can be summarized as follows.

F1. Compare two distinct systems with different degree distributions of agents. Then, when there are a small number of random failures in some CS, the system in which agents’ degrees exhibit higher variability is less likely to suffer widespread failures. In particular, suppose that agents have identical or similar average degrees (hence comparable levels of dependence among agents and CSes) in two different systems. In this case, the system with more homogeneous or predictable degrees (thus less variability in agents’ degrees) is more susceptible to extensive failures in the system.

F2. Consider two systems in which agents’ degree distributions have identical marginal distributions. Therefore, loosely speaking, we can say that they display the same level of (inter-)dependence. In this case, the system in which the agents’ degrees are more positively correlated is less likely to experience cascading failures as a result of initially localized, random failures.

The first finding indicates that degree distributions with higher variability, such as power laws, which permit the existence of large degree hubs, are more robust to random failures than more concentrated distributions such as Poisson distributions. Furthermore, it hints that systems in which all agents in each CS have similar degrees are most prone to an outbreak of failures. This observation is consistent with earlier numerical studies (e.g., [8, 9]) that suggest that scale-free networks with power law degree distributions are more resilient to random attacks, but vulnerable to intentional attacks that target high-degree nodes.

The second finding above may be somewhat counter-intuitive at first sight. One might suspect that positive correlations would be helpful to spreading failures because high-degree agents are likely to have even larger aggregate degrees with increasing positive correlations and thus serve as more effective conduits for transmitting failures. However, our finding reveals that stronger positive correlations have similar effects as higher variability in the first finding. We suspect that the reason behind this is that stronger positive correlations increase the variability in the aggregate degree of agents. Consequently, they hinder the proliferation of failures, rendering the system more robust to random failures.

A few words on notation: throughout the paper, we will use boldface letters or symbols to denote (row) vectors or vector functions.111All vectors are assumed to be row vectors. For instance, ${\bf d}$ denotes a vector, and the $j$-th element of ${\bf d}$ is denoted by $d_{j}$. Vectors ${\bf 0}$ and ${\bf 1}$ represent the vectors of zeros and ones, respectively, of appropriate dimensions. The set ${\mathbbm{Z}}_{+}$ (resp. ${\rm I\kern-1.5ptN}$) denotes the set of nonnegative integers $\{0,1,2,\ldots\}$ (resp. positive integers $\{1,2,3,\ldots\}$). Finally, all vector inequalities are assumed componentwise.

2 Related Literature

There is already a large volume of literature on related topics, including cascading failures and robustness of complex systems [1, 4, 5, 43], spread of epidemics and efficient immunization [37, 39, 45], and information or rumor propagation [11, 46]. Given the significant body of studies in related fields, it is not possible to provide a summary of all. For this reason, we limit our discussion to a short list of most pertinent studies in the settings of multiplex or interdependent networks, and do not discuss other relevant studies (e.g., [2, 6, 14, 26, 29, 42]), including many important studies on a single, monolithic network (e.g., [3, 8, 9, 11, 22, 30, 31, 34, 35, 44]), here. We instead refer an interested reader to the references and those therein.

In [7], Buldyrev et al. investigated cascades of failures in two interdependent networks – networks A and B – using numerical studies. In their model, each node in network A (resp. B) depends on a randomly chosen node in network B (resp. A), which is modeled using a directed support link, and the failure of the node on which a node depends causes its own failure. Furthermore, nodes in network A (resp. B) are connected with each other according to a degree distribution $P_{A}$ (resp. $P_{B}$). Initially, a fraction ($1-p$) of network A nodes are removed, which triggers failures of nodes in both networks through connectivity and dependence. They studied the probability that a giant component survives as a function of $p$, and identified a threshold for a first-order phase transition. In addition, their numerical results suggest that broader degree distributions $P_{A}$ and $P_{B}$ make the network more vulnerable to random failures, whereas in a single scale-free network, the opposite has been observed.

The findings in [7] have been extended in a series of follow-up studies. Parshani et al. [36] studied a similar model and demonstrated that, as the coupling between the two networks diminishes, the phase transition changes to a second-order transition (from a first-order transition). In [41], Shao et al. relaxed the assumption that each node is dependent on exactly one node in the other network and modeled the number of support links of nodes using random variables. Huang et al. [18] considered the robustness of the system against targeted attacks by mapping the problem to a previously studied problem with random attacks, and suggested that the presence of high degree nodes in scale-free networks makes it challenging to protect interdependent networks against targeted attacks.

In [46], Zhuang and Ya$\breve{\rm g}$an studied information propagation in a multiplex network with two layers representing an online social network (OSN) and a physical network, both of which have high clustering. Only a subset of vertices in the physical network are assumed to be active in the OSN. Their key findings are: (a) clustering consistently hampers cascades of information to a large number of nodes with respect to both the critical threshold of information epidemics and the mean size of epidemics; and (b) information transmissibility (i.e., average probability of information transmission over a link) has significant impact; when the transmissibility is low, it is easier to trigger a cascade of information propagation with a smaller, densely connected OSN than with a large, loosely connected OSN. However, when the transmissibility is high, the opposite is true.

In another study [17], Hu et al. studied the problem of viral influence spreading, for instance, in adoption of new technologies or scientific ideas. They modeled the spread of adoption using a multiplex network in which there are two different types of links – (i) undirected connectivity links and (ii) directed influence links. Outgoing influence links of a node lead to other nodes whose adoption of a new technology or idea causes the node to adopt it with some fixed probability. Similarly, incoming influence links of a node originate from other nodes that watch the node to see if it adopts the technology or idea first and, if so, follow its trend with fixed probability.

Their key findings include the following: (a) viral cascades are feasible only if there are positive correlations between the connectivity degrees and outgoing influence degrees of nodes. The intuition is that when there are positive correlations, even the adoption of a new technology by some random node would make it easy to influence nodes with high connectivity degrees because they tend to have larger outgoing influence degrees, hence more likely to be influenced by other adopters; and (b) positive correlations between connectivity degrees and incoming influence degrees do not facilitate viral cascades much.

Khamfroush et al. [21] investigated the propagation of phenomena (e.g., failures, infections or rumors) in interdependent networks. By introducing temporal dynamics into the model, they studied how quickly phenomena spread in three different types of networks – scale-free networks, small-world networks, and Erd$\ddot{\rm o}$s-R$\acute{\rm e}$nyi networks. Two of key observations from their simulation studies are: (a) scale-free networks are in general more conducive to spreading phenomena than the other two types; and (b) the choices of initial spreaders can have greater impact than the network type. Based on the latter observation, they proposed a new centrality metric, called path-degree centrality, to better identify more effective initial spreaders.

While many of these studies (e.g., [7, 18, 36, 41]) examine the robustness of multiplex or interdependent networks, their approaches and goals are very different from those of our study. First, most of the aforementioned studies focus on the analysis of the emergence or survival of giant components under the assumption that only the nodes that belong to the giant component can continue to function properly. We, on the other hand, investigate (i) when it is possible to see an epidemic of failures through dependencies among different systems and (ii) how the underlying dependence properties between systems shape the likelihood of such catastrophic failures.

Second, although the existing studies summarized here (and others we are aware of) provide interesting observations and major contributions to the growing understanding of complex systems, to the best of our knowledge, none of these studies aims to present analytical findings that enable us to compare the robustness of different networks on the basis of their dependence properties (which can be partially ordered). In contrast, as summarized in Section 1, our intent is to take another step towards building a comprehensive theory on complex systems which will help us determine when the robustness of such systems improves or deteriorates as a consequence of changes in their dependence properties. Some of our preliminary results are reported in [24].

3 System Model

Let $N_{S}>1$ be the number of CSes in the (global) system we consider. For each $i\in{\mathcal{S}}:=\{1,2,\ldots,N_{S}\}$, let ${\mathcal{A}}_{i}$ be the set of agents in the $i$-th CS. When convenient, we use $a^{i}$ and $a^{i,k}$ to denote a generic agent in CS $i$ and the $k$-th agent in ${\mathcal{A}}_{i}$, respectively.

We model the internal or intra-CS interdependence among agents in the $i$-th CS using a dependence graph ${\mathcal{G}}_{i}=({\mathcal{V}}_{i},{\mathcal{E}}_{i})$: the vertices in ${\mathcal{V}}_{i}$ are the agents in CS $i$. The edges in ${\mathcal{E}}_{i}$ are undirected edges between vertices in ${\mathcal{V}}_{i}$ and indicate mutual dependence relations between the end vertices.222These dependence relations are not necessarily the physical links in a network. For example, in a power system, an overload failure in one part of power grid can cause a failure in another part that is not geographically close or without direct physical connection to the former. An undirected edge $e\in{\mathcal{E}}_{i}$ should be interpreted as a pair of directed edges pointing in the opposite directions. Two agents with an undirected edge between them are said to be (dependence) neighbors.

In addition to the (undirected) edges between agents in the same CS, we model the dependence of an agent in one CS on another agent in a different CS using a directed edge; if there is a directed edge from agent $a^{i}$ to agent $a^{j}$, where $j\neq i$, this means that $a^{j}$ depends on $a^{i}$ and, when $a^{i}$ fails, it could cause $a^{j}$ to crash as well. We do not assume that this dependence is mutual to allow asymmetric dependence among CSes.333Mutual inter-CS dependence can be handled by replacing the directed edges with undirected edges similar to those used to model the intra-CS interdependence, and vice versa, without affecting our key findings. Also, a similar model is employed in [21]. If there is a directed edge $a^{i}\to a^{j}$ ($j\neq i$), we say that $a^{j}$ is a CS $j$ (dependence) neighbor of $a^{i}$ and that $a^{i}$ supports $a^{j}$.

We oftentimes need to distinguish the neighbors in the same CS from those in other CSes. For this reason, we call the neighbors in the same CS (resp. other CSes) internal neighbors (resp. external neighbors). Note that an external neighbor of an agent is another agent in a different CS which it supports. In addition, we refer to the number of internal neighbors of an agent as its internal degree.

3.1 Agent degree distributions

Throughout the paper, we shall use ${\bf D}=(D_{j};\ j\in{\mathcal{S}})$ to denote a $N_{S}$-dimensional (random) vector that describes the number of neighbors that an agent has in each CS $j\in{\mathcal{S}}$. In other words, $D_{j}$, $j\in{\mathcal{S}}$, is the number of CS $j$ neighbors. We call ${\bf D}$ the (dependence) degree vector of the agent. Note that, for a CS $i$ agent, its degrees $D_{j}$, $j\neq i$, represent the number of agents that it supports in other CSes and hence can affect in case of its own failure, but not those that support it. Thus, they denote the outgoing degrees of the agent.

The degrees of a CS $i$ agent are denoted by a random vector ${\bf D}_{i}=(D_{i,j};\ j\in{\mathcal{S}})$ with a distribution $p_{i}$; given ${\bf d}\in{\mathbbm{Z}}_{+}^{N_{S}}$, the probability that a randomly chosen CS $i$ agent has $d_{j}$ neighbors in CS $j$, $j\in{\mathcal{S}}$, is equal to $p_{i}({\bf d})$. We find it convenient to define the marginal distributions $p_{i,j}:{\mathbbm{Z}}_{+}\to[0,1]$, $j\in{\mathcal{S}}$, where p_{i,j}(d_{j})={\mathbbm{P}}\left[\mbox{CS i$agent has$d_{j}$CS$j neighbors}\right].

Note that the degrees of an agent to different CSes are not assumed to be mutually independent. Put differently, the number of other agents that a CS $i$ agent, say $a^{i}$, supports in different CSes, i.e., ($D_{i,j}$; $j\in{\mathcal{S}}\setminus\{i\}$), could be correlated and depend on its internal degree $D_{i,i}$. This is important because in practice the failure of an important agent in a system may trigger the failure of many other agents across different CSes, suggesting that the degrees of such agents could be correlated. Thus, we wish to study the impact of such degree correlations on system robustness.

Finally, throughout this paper (except for in Section 6), we assume that the internal degree of an agent is at least one with probability one, i.e., $p_{i,i}(0)=0$ for all $i\in{\mathcal{S}}$; otherwise, the agent should not belong to CS $i$.

3.2 Propagation of failures

To study the robustness of a system to failures, we need to model how a failure spreads from one agent to another. Here, we explain the models we employ to approximate the dynamics of failure propagation both within a CS and between agents in different CSes.

P1. Intra-CS failure propagation – We model failure propagation within CS $i$ with the help of a function $\wp^{in}_{i}:{\rm I\kern-1.5ptN}^{2}\to[0,1]$: for fixed $d\in{\rm I\kern-1.5ptN}$ and $1\leq n_{f}\leq d$, $\wp^{in}_{i}(d,n_{f})$ tells us the probability that a CS $i$ agent with an internal degree $d$ will fail when $n_{f}$ internal neighbors collapse.

An example that fits this model is the random threshold model used by Watts in [44], which is also used in other studies (e.g., [3, 6, 21]). In the Watts’ model, every agent $a^{i,k}\in{\mathcal{A}}_{i}$ has some intrinsic value $\xi^{i,k}$. These values $\xi^{i,k}$ of CS $i$ agents are modeled using mutually independent (continuous) random variables (rvs) with some common distribution $F_{i}$. We refer to $\xi^{i,k}$ as their security states.

In his model, a CS $i$ agent, say $a^{i}$, goes down as a consequence of the failures of its internal neighbors when the fraction of its failed internal neighbors exceeds its security state $\xi^{i}$. Therefore, for a given pair $(d,n_{f})$ with $n_{f}\leq d$, $\wp^{in}_{i}(d,n_{f})$ is equal to ${\mathbbm{P}}\left[\xi^{i}<n_{f}/d\right]=F_{i}(n_{f}/d)$.

P2. Inter-CS failure propagation – We model the propagation of a failure from one agent to an external neighbor in a similar manner. Suppose that agent $a^{j}$ is a CS $j$ agent. Denote the number of CS $i^{\prime}$ agents that support $a^{j}$ by $D^{in}_{j,i^{\prime}}$. We call $D^{in}_{j,i^{\prime}}$ the incoming CS $i^{\prime}$ degree of $a^{j}$ (to distinguish it from its outgoing degree to CS $i^{\prime}$).

Although our model can be generalized to allow the failure probability of $a^{j}$ (as a result of failures of CS $i$ agents supporting $a^{j}$) to be a function of all of its incoming degrees ($D^{in}_{j,i^{\prime}}$; $i^{\prime}\in{\mathcal{S}}\setminus\{j\}$) without altering our main results, for the ease of exposition, here we adopt a simpler model in which the failure probability only depends on $D^{in}_{j,i}$: when $a^{j}$ has incoming CS $i$ degree $D^{in}_{j,i}$ and $n^{i}_{f}$ of these supporting CS $i$ agents collapse, the probability that $a^{j}$ also fails as a result is given by some function $\wp^{ex}_{j,i}(D^{in}_{j,i},n^{i}_{f})$. Also, we assume that the external neighbors of a failed agent go down with the prescribed probability independently of each other.

It is clear that this model is general enough to include the one studied in [41], where $a^{j}$ is unaffected by the failure of a supporting CS $i$ agent as long as there is another supporting CS $i$ agent that is still operational.

3.3 Tree-like infection graphs

For our study, we focus on scenarios where collapsed agents do not cause the failures of many neighbors on the average.444When each failed agent triggers many other neighbors to crash as well, cascading failures are likely and should happen often. This may indicate that the system is poorly designed. Instead, we are interested in more realistic scenarios of interest in which cascading failures are possible and do occur, but not too frequently. To make this more precise, we introduce infection graphs: starting with an initial failure of an agent in the system, the infection graph consists of all failed agents and the directed edges used to contribute to the failures of neighbors. In other words, a directed edge from agent $a_{1}$ to agent $a_{2}$ belongs to the infection graph if and only if (iff) agent $a_{2}$ is a neighbor of $a_{1}$ and agent $a_{1}$ failed before agent $a_{2}$ did. We assume that this infection graph can be approximated using a tree-like structure, which we call an infection tree. A similar assumption is introduced in [6, 22, 44, 45].

Although this assumption may not always hold in real systems, it allows us to approximate the dynamics of spreading failures as a multi-type branching process (described in Section 5) and to carry out an analytical study. Moreover, when a graph is sparsely connected in large networks with small average degrees, with high probability there are only few short cycles in the giant component [19]. Thus, it is a reasonable assumption when failed agents do not cause the crash of a large number of neighbors in each CS, which is the scenario of interest to us.

When there is a directed edge from agent $a_{1}$ to agent $a_{2}$ in an infection tree, i.e., the failure of $a_{1}$ causes that of $a_{2}$, we refer to $a_{1}$ (resp. $a_{2}$) as the parent (resp. a child). Also, borrowing from the language of epidemiology, we say that the parent infected the child.

4 Agent Types and Children Distributions

As mentioned in Section 1, we are interested in scenarios where the number of agents in each CS is large. In a system with many agents, the propagation of failures can be approximated with the help of a multi-type branching process under some simplifying assumptions. To be more precise, we shall borrow from the theory of branching processes with finitely many types in order to study the likelihood of a small number of initial failures leading to an epidemic of failures infecting many other agents.

4.1 Agent types

In our model, depending on how an agent is infected, there are two possible types we need to consider for the failed agent. To formalize this, we define the types of agents in varying CSes. Given $N_{S}$ CSes, there are $2N_{S}$ types of interest to us. A CS $i$ agent ($i\in{\mathcal{S}}$) can be either type $i$ or $i+N_{S}$: a type $i$ agent is a CS $i$ agent whose internal neighbors are all functional, i.e., have not crashed. On the other hand, a CS $i$ agent is of type $i+N_{S}$ if it has an internal neighbor that went down. For notational simplicity, we use $i^{+}$ to denote $i+N_{S}$ ($i\in{\mathcal{S}}$) hereafter.

We shall discuss the distribution of the number of children of various types which are produced by a failed agent of type $i\in\{1,2,\ldots,2N_{S}\}=:{\mathcal{S}}^{+}$ shortly. To explain these children (vector) distributions, we first need to describe how we approximate the probability that a neighbor of a failed agent falls victim to infection.

4.2 Infection probability of neighbors

The assumption that the infection graph is tree-like has an important implication: except for the root of the infection graph, each failed agent has exactly one supporting agent whose collapse led to its own failure, namely its parent. This observation helps us simplify the models used to capture the propagation of failures as follows.

Because an agent faces possible infection from at most one failed supporting agent under the assumption, we can approximate the probability that a neighbor of a failed agent, say $a^{i}$, will be infected as explained below.

$\bullet$** Intra-CS infection probability –** Consider an internal neighbor, $\tilde{a}^{i}$, of agent $a^{i}$. Following the explanation in [22, 44], the probability that agent $\tilde{a}^{i}$ has $d$ internal neighbors is proportional to $d\cdot p_{i,i}(d)$. Recall that $\tilde{a}^{i}$ with internal degree $d$ will be infected by the failure of $a^{i}$ with probability $\wp^{in}_{i}(d,1)$. Thus, by conditioning on the internal degree of $\tilde{a}^{i}$, the probability that an internal neighbor of a failed CS $i$ agent will be infected can be approximated using

[TABLE]

where $d_{i,{\rm avg}}:=\sum_{d\in{\rm I\kern-1.5ptN}}d\cdot p_{i,i}(d)$ is the average internal degree of CS $i$ agents, and $w_{i,i}(d):=d\cdot p_{i,i}(d)/d_{i,{\rm avg}}$. In the example of Watts’ model, we have $q_{i,i}=\sum_{d\in{\rm I\kern-1.5ptN}}w_{i,i}(d)\ F_{i}(d^{-1})$.

The average internal degree $d_{i,{\rm avg}}$ in the denominator of (1) serves as a normalizing constant so that ${\bf w}_{i,i}:=(w_{i,i}(d);\ d\in{\rm I\kern-1.5ptN})$ is the internal degree distribution of a randomly picked internal neighbor [22, 44].555This sampling technique is called sampling by random edge selection [27].

$\bullet$** Inter-CS infection probability – Suppose that a CS $j$ agent, $a^{j}$, is an external neighbor of a failed CS $i$ agent $a^{i}$. Recall that the inter-CS dependence is asymmetric. Let $p^{in}_{j,i}$ be the conditional distribution of $D^{in}_{j,i}$ given that $a^{j}$ is a CS $j$ neighbor of a CS $i$ agent. Then, we approximate the probability that $a^{j}$ will be infected by the failure of $a^{i}$ using $q_{i,j}:=\sum_{d\in{\rm I\kern-1.5ptN}}p^{in}_{j,i}(d)\ \wp^{ex}_{j,i}(d,1)$.**

Consider the model studied in [41] where the agent $a^{j}$ survives as long as one other supporting CS $i$ agent is functional. In this case, agent $a^{j}$ will be infected by the failure of $a^{i}$ iff $D^{in}_{j,i}=1$. In other words, $\wp^{ex}_{j,i}(1,1)=1$ and $\wp^{ex}_{j,i}(d,1)=0$ for all $d>1$. Thus, $q_{i,j}=p^{in}_{j,i}(1)$.

The conditional distribution $p^{in}_{j,i}$ is determined by the distribution of incoming degrees $D^{in}_{j,i}$ of CS $j$ agents. In our study, we assume that the incoming degree distributions, hence their conditional distributions and $q_{i,j}$, $j\neq i$, are fixed while we study the influence of the variability and dependence of internal and outgoing degrees of agents.

We believe that many of system parameters used in the model can be estimated in practice, for instance, from historical data, physical laws (e.g., power grid), or simulation studies. These include degree distributions $(p_{i};\ i\in{\mathcal{S}})$ and infection probabilities $(q_{i,j};\ i,j\in{\mathcal{S}})$, and can be used to examine the robustness of the system.

4.3 Distributions of children vectors

Let $a^{i}$ be a randomly picked CS $i$ agent, and assume that it crashes and infects some of its neighbors. Then, the type of its child belongs to $({\mathcal{S}}\cup\{i^{+}\})\setminus\{i\}=:{\mathcal{S}}^{+}_{i}$: if $a^{i}$ triggers the failure of an external neighbor in CS $j$, $j\neq i$, the type of the child is simply $j$. If $a^{i}$ causes an internal neighbor to fail, then the child’s type is $i^{+}$ because agent $a^{i}$ is an infected internal neighbor of the child.

Based on this observation, we can approximate the distribution of the number of children produced by $a^{i}$ by considering its two possible types.

C1. Type $i$ agent $a^{i}$ ($i\in{\mathcal{S}}$) – Let $D_{i,i}$ be its internal degree. Some of these internal neighbors, however, may not be affected by the failure of agent $a^{i}$ and remain uninfected. Similarly, an external neighbor of agent $a^{i}$ in CS $j$ will go down (as a consequence of $a^{i}$’s failure) with probability $q_{i,j}$. As a result, the actual number of CS $j$ neighbors infected by $a^{i}$, which we denote by $O_{i,j}$, can differ from $D_{i,j}$, $j\in{\mathcal{S}}\setminus\{i\}$.

Summarizing this argument, we approximate the probability that $a^{i}$ will produce children given by a children vector ${\bf o}=(o_{j};\ j\in{\mathcal{S}}^{+})\in{\mathbbm{Z}}_{+}^{2N_{S}}$, where $o_{j}$ is the number of type $j$ children, using the following children distribution:

[TABLE]

where $q_{i,j}$ are as defined earlier, $\bar{q}_{i,j}:=1-q_{i,j}$ for all $j\in{\mathcal{S}}$, and $\bar{{\bf d}}:{\mathbbm{Z}}_{+}^{2N_{S}}\to{\mathbbm{Z}}_{+}^{N_{S}}$ with

[TABLE]

C2. Type $i^{+}$ agent $a^{i}$ – For type $i^{+}$ agents ($i\in{\mathcal{S}}$), the children distribution is closely related to that of type $i$ agents with a minor difference: for ${\bf o}\in{\mathbbm{Z}}_{+}^{2N_{S}}$,

[TABLE]

where ${\bf e}_{i}$ is a zero-one vector whose only nonzero entry is the $i$-th entry.

The only difference between (2) and (10) is that, for type $i^{+}$ agents, we first remove one of internal neighbors before counting the neighbors that can be infected by the agents. The reason for this is that a type $i^{+}$ agent has a parent in CS $i$ and the number of remaining internal neighbors that it can potentially infect is its internal degree minus one. We denote the marginal distribution of the number of type $j$ children of a type $i$ agent by $h_{i,j}$, $i,j\in{\mathcal{S}}^{+}$.

Before we proceed, let us comment on a simplifying assumption we implicitly introduced in Section 4.2 and our approximations in (2) and (10). Suppose that $a^{i}$ is a CS $i$ agent that is infected by an internal neighbor. Given that $a^{i}$’s failure is caused by another internal neighbor, its conditional degree distribution will likely be different from $p_{i}$ we assumed in (2) and (10).

For instance, if $d\cdot\wp^{in}_{i}(d,1)$ is increasing (resp. decreasing), its internal degree $D_{i,i}$ conditional on the event that it is infected by an internal neighbor, is larger (resp. smaller) than the unconditional internal degree with respect to the usual stochastic order [40]. The reason for this is that the internal degree distribution of an internal neighbor of $a^{i}$ is given by ${\bf w}_{i,i}$, where $w_{i,i}(d)\propto d\cdot p_{i,i}(d)$. Therefore, if $d\cdot\wp^{in}_{i}(d,1)$ is increasing in $d$, we have

[TABLE]

Moreover, this conditional degree distribution could differ from the conditional degree distribution we would see provided that $a^{i}$ is infected by a parent in a different CS, which would also depend on the CS to which the parent belongs. Therefore, it is clear that the conditional degree distribution of an infected agent will likely differ from $p_{i}$ and depend on how it was infected.

Unfortunately, computing and adopting accurate conditional degree distributions for the analysis is quite challenging for several reasons. For example, in order to compute the necessary conditional probabilities, for each failed agent $a^{i}$, we need to know exactly how it was infected. More precisely, we must take into account the history or the sequence of agents that collapsed and led to the infection of $a^{i}$, as well as the joint distributions of agents’ incoming and outgoing degrees. The reason for this is that the conditional degree distribution of the parent of $a^{i}$ in turn depends on its parent and so on.

Iteratively computing the conditional degree distributions of all infected agents while accounting for the history and relevant joint degree distributions quickly becomes intractable. For this reason, in order to maintain mathematical tractability of the model, we make a simplifying assumption that the (conditional) degree distribution of infected CS $i$ agents can be approximated by $p_{i}$ for all $i\in{\mathcal{S}}$. However, we believe that this is a reasonable assumption, especially when we compare two systems with similar degree distributions for local comparison (as the assumption or its failure would affect both of them alike).

5 Multi-type Branching Process for Modeling

the Spread of Failures

We approximate the propagation of failures, using a multi-type branching process: when a type $i$ agent ($i\in{\mathcal{S}}^{+}$) fails, it produces children of various types in accordance with the children distribution described in the previous section, independently of other infected agents in the system.

5.1 Infection tree

Suppose that a CS $i$ agent ($i\in{\mathcal{S}}$), say $a^{i}$, is the first agent to experience a random failure. As mentioned in Section 1, we are interested in determining: (i) if it is possible for this single failure to lead to widespread infection of many other agents through dependence among agents in multiple CSes, and (ii) if so, how likely it is for the system to suffer such a cascade of failures.

To answer these questions, we consider a (directed) infection tree that captures the spread of failures, which is rooted at agent $a^{i}$. We denote the tree by ${\mathcal{T}}:=({\mathcal{V}}_{{\mathcal{T}}},{\mathcal{E}}_{{\mathcal{T}}})$, where ${\mathcal{V}}_{\mathcal{T}}$ is the set of failed agents and ${\mathcal{E}}_{\mathcal{T}}$ is the set of directed edges via which failures transmitted. For each $k\in{\mathbbm{Z}}_{+}$, let ${\mathcal{N}}(k)=({\mathcal{N}}_{j}(k);\ j\in{\mathcal{S}}^{+})$ denote the set of $k$-hop neighbors of $a^{i}$ in ${\mathcal{T}}$, where ${\mathcal{N}}_{j}(k)$ is the set of type $j$ $k$-hop neighbors.

Note that ${\mathcal{N}}(k)$, $k\in{\rm I\kern-1.5ptN}$, are random sets, and we are mostly interested in the cardinalities $N_{j}(k)=|{\mathcal{N}}_{j}(k)|$. In the example with three CSes ($N_{S}=3$) shown in Fig. 1, the initial failure occurs in CS 2 (Root shown as a filled red circle). The dotted arrows indicate how the failures transmitted between agents. Here, $N_{j}(1)=1$ for $j\in\{1,5\}$ and $N_{j}(1)=0$ for $j\in\{2,3,4,6\}$. Similarly, $N_{4}(2)=2$ and $N_{j}(2)=0$ for $j\in\{1,2,3,5,6\}$. The tree ${\cal T}$ consists of the root and infected agents (filled orange circles) along with the dotted arrows.

Regrettably, computing the exact distribution of the total number of failed agents (i.e., $1+\sum_{k\in{\rm I\kern-1.5ptN}}\left(\sum_{j\in{\mathcal{S}}^{+}}N_{j}(k)\right)$) is challenging, if possible at all, for large systems. For this reason, we follow a similar approach employed in [6, 22, 44, 45] and, rather than analyzing a finite system, consider an infinite system in which the degree vector of each CS $i$ agent is given by a random vector with a common distribution $p_{i}$, independently of each other. In other words, the degree vectors of CS $i$ agents are given by independent and identically distributed (i.i.d.) random vectors with the distribution $p_{i}$. Moreover, the degree vectors of agents in different CSes are assumed mutually independent. This degree-based model is also known as the Chung-Lu model in the literature [10].

By the strong law of large numbers, the fraction of CS $i$ agents with degree vector ${\bf d}$ converges to $p_{i}({\bf d})$ almost surely for all ${\bf d}\in{\mathbbm{Z}}_{+}^{N_{S}}$. Using this model, we will first look for a condition under which ${\mathbbm{P}}\left[\limsup_{k\to\infty}N_{j}(k)=0\mbox{ for all }j\in{\mathcal{S}}^{+}\right]<1$. Put differently, there is positive probability that the failures will continue to propagate forever in an infinite system. We shall use this probability of cascading failures (PoCF) in an infinite system to approximate the probability that a large system would experience an epidemic of failures.

The answer to this question can be obtained by studying a multi-type branching process with $2N_{S}$ types. Let ${\bf N}=\{{\bf N}(k),\ k\in{\mathbbm{Z}}_{+}\}$, where ${\bf N}(k)=(N_{i}(k);\ i\in{\mathcal{S}}^{+})$ and $N_{i}(k)$, $i\in{\mathcal{S}}^{+}$, is the number of type $i$ agents in the $k$-th generation. Recall that, for $i\in{\mathcal{S}}$, a type $i$ agent (resp. $i^{+}$ agent) is a CS $i$ agent with no infected internal neighbor (resp. with a failed internal neighbor).

5.2 Probability of extinction

The probability ${\mathbbm{P}}\left[\limsup_{k\to\infty}N_{j}(k)=0\mbox{ for all }j\in{\mathcal{S}}^{+}\right]$ is called the probability of extinction (PoE) [15]. Obviously, the PoCF is equal to one minus the PoE. Since the initial failure can originate in any of $N_{S}$ CSes, we denote the PoE starting with a random failure in CS $i$ by $\mu_{i}$.

For each $i\in{\mathcal{S}}^{+}$, let ${\mathbbm{E}}\left[h_{i}\right]=({\mathbbm{E}}\left[h_{i,j}\right];\ j\in{\mathcal{S}}^{+})$ be a $1\times 2N_{S}$ row vector, whose $j$-th element is the expected number of type $j$ children from a failed type $i$ agent. Define ${\mathbbm{M}}=[M_{i,j}]$ to be a $(2N_{S})\times(2N_{S})$ matrix, whose $i$-th row is ${\mathbbm{E}}\left[h_{i}\right]$, i.e., $M_{i,j}={\mathbbm{E}}\left[h_{i,j}\right]$ for all $i,j\in{\mathcal{S}}^{+}$.

An example of ${\mathbbm{M}}$ with $N_{S}=2$ is shown below.

[TABLE]

Note that $M_{i,i}=M_{i^{+},i}=0$ for all $i\in{\mathcal{S}}$ because a CS $i$ agent infected by another CS $i$ agent will be type $i^{+}$ ($i\in{\mathcal{S}}$). Similarly, $M_{i,j^{+}}=0$ for all $i,j\in{\mathcal{S}}$ and $j\neq i$ under the assumption of no cycles in the infection tree.

Definition 1

A square matrix ${\bf A}$ is said to be (positively) regular if there exists $k\in{\rm I\kern-1.5ptN}$ such that ${\bf A}^{k}$ is positive, i.e., all entries are positive.

Assumption 1

*We assume that ${\mathbbm{M}}$ is (positively) regular.

One can show that a sufficient condition for the (positive) regularity of ${\mathbbm{M}}$ is that (i) $p_{i,i}(0)+p_{i,i}(1)<1$ for all $i\in{\mathcal{S}}$ and (ii) ${\mathbbm{M}}$ is irreducible.666Let ${\bf A}=[A_{i,j}]$ be a $2N_{S}\times 2N_{S}$ matrix, where $A_{i,j}={\bf{1}}\left\{M_{i,j}>0\right\}$ for all $i,j\in{\mathcal{S}}^{+}$. This matrix ${\bf A}$ can be viewed as an adjacency matrix for a directed graph ${\cal G}$ with $2N_{S}$ vertices that represent the $2N_{S}$ types. There is a directed edge from vertex $i$ to vertex $j$ if a type $i$ agent produces a type $j$ agent with positive probability. The matrix ${\mathbbm{M}}$ is irreducible if and only if the directed graph ${\cal G}$ is strongly connected [13]. The first condition simply means that when a CS agent $i$ fails, there is positive probability that it will infect another internal neighbor. The second condition ensures that a random initial failure in any CS can eventually cause some agents in every other CS $j$, $j\neq i$, to go down with positive probability, following a sequence of infections.

If ${\mathbbm{M}}$ is regular, there exists $k^{\star}\in{\rm I\kern-1.5ptN}$ such that, starting with any random failure, regardless of which CS experiences the initial failure, there is strictly positive probability that ${\bf N}(k^{\star})\geq{\bf 1}$, i.e., there is a failure in every CS.

Let $\bm{\mu}=(\mu_{i};\ i\in{\mathcal{S}}^{+})$. Although the PoEs of interest to us are ($\mu_{i}$; $i\in{\mathcal{S}}$), we will compute ($\mu_{i^{+}}$; $i\in{\mathcal{S}}$) as well. Under the (positive) regularity assumption, $\bm{\mu}={\bf 1}$ if (i) $\rho({\mathbbm{M}})<1$ or (ii) $\rho({\mathbbm{M}})=1$ and there is at least one type for which the probability that it produces exactly one child is not equal to one, where $\rho({\mathbbm{M}})$ is the spectral radius of ${\mathbbm{M}}$ [16]. Similarly, if $\rho({\mathbbm{M}})>1$, then $\bm{\mu}<{\bf 1}$ and there is strictly positive probability that the cascading failures continue forever in an infinite system, suggesting that there could be an epidemic of failures in a large system.

It is noteworthy that whether or not there could be cascading failures in the infinite system depends only on the mean number of children of varying types that each type produces. However, the exact PoEs $\bm{\mu}$ vary from one set of children distributions $\{h_{i};i\in{\mathcal{S}}^{+}\}$ to another set $\{\tilde{h}_{i};i\in{\mathcal{S}}^{+}\}$ even when the matrix ${\mathbbm{M}}$ remains the same.

For each $i\in{\mathcal{S}}^{+}$, define a generating function $f_{i}:[0,1]^{2N_{S}}\to[0,1]$, where

[TABLE]

Then, the PoE vector $\bm{\mu}$ is given as a fixed point that satisfies

[TABLE]

where ${\bf f}({\bf s})=(f_{1}({\bf s}),\ldots,f_{2N_{S}}({\bf s}))$. When $\rho({\mathbbm{M}})>1$, there exists a unique $\bm{\mu}$ that satisfies (16) with strict inequality, i.e., $\bm{\mu}<{\bf 1}$ [15].

A key question of interest to us is how the degree distributions ($p_{i}$; $i\in{\mathcal{S}}$) affect the PoE vector $\bm{\mu}$, especially with fixed ${\mathbbm{M}}$. To be more precise, we will investigate how the variability and dependence structure of agents’ degree vectors shape the PoEs. To this end, we introduce several stochastic and dependence orders that we employ to compare the degree distributions. Using these orders, we first examine a simple scenario consisting of two symmetric interdependent CSes in the subsequent section, followed by more general settings in Section 7.

6 Two symmetric interdependent systems

The goal of this section is, by studying simpler scenarios first, to highlight some insights on how (i) the variability of degrees of agents (i.e., the number of neighbors in two different CSes) and (ii) the dependence of the two degrees influence the PoEs, even when the mean degrees remain fixed.

Consider a system comprising two interdependent CSes ($N_{S}=2$), and suppose that the degree distributions $p_{1}$ and $p_{2}$ are symmetric, i.e., $p_{1}(d_{1},d_{2})=p_{2}(d_{2},d_{1})$ for all $(d_{1},d_{2})\in{\mathbbm{Z}}_{+}^{2}$. Moreover, in order to simplify the analysis and shed some light on our main findings in general settings to follow, we set $q_{i,j}=1$ for all $i,j\in{\mathcal{S}}$. This assumption will be relaxed in the subsequent section. To correctly interpret this assumption and the findings in this section, a reader should view the symmetric degree distributions $p_{i}$, $i=1,2,$ as children distributions $h_{i}$, $i=1,2,$ instead; otherwise, if all internal dependence graphs ${\cal G}_{i}$, $i\in{\mathcal{S}}$, are connected and $q_{i,j}=1$ for all $i,j\in{\mathcal{S}}$, under positive regularity assumption of ${\mathbbm{M}}$, every agent will eventually be infected, starting with any failure. For this reason, we remove the assumption that the internal degree is at least one in this section to allow for the possibility that some agents do not produce any children. This assumption will be reintroduced in the following section.

Throughout this and following sections, we assume that $\rho({\mathbbm{M}})>1$ under all considered degree distributions so that it is possible for a random failure to trigger cascading failures. In the case of two interdependent CSes, the assumed symmetry of the degree distributions and the uniqueness of the fixed point $\bm{\mu}$ satisfying ${\bf f}(\bm{\mu})=\bm{\mu}<{\bf 1}$ tell us $\mu_{1}=\mu_{2}$ and $\mu_{3}=\mu_{4}$.

6.1 The effects of degree variability

We first study the variability of degrees. One common way to compare the variability of two rvs is the second-order stochastic dominance (SSD) [40]. Loosely speaking, if rv $Y$ dominates rv $X$ with respect to SSD ($X\leq_{SSD}Y$) and ${\mathbbm{E}}\left[X\right]={\mathbbm{E}}\left[Y\right]$, it means that $Y$ is more predictable than $X$. It turns out $X\leq_{SSD}Y$ is equivalent to $X$ being smaller than $Y$ with respect to increasing concave (ICV) order ($X\leq_{icv}Y$); for all increasing, concave functions $\xi:{\rm I\kern-2.5ptR}\to{\rm I\kern-2.5ptR}$, ${\mathbbm{E}}\left[\xi(X)\right]\leq{\mathbbm{E}}\left[\xi(Y)\right]$.777These inequalities in the definitions of various stochastic and dependence orders are required to hold only for the functions for which the expectations are well defined. Since $\xi(x)=x$ is concave and increasing, $X\leq_{SSD}Y$ implies ${\mathbbm{E}}\left[X\right]\leq{\mathbbm{E}}\left[Y\right]$.

In order to eliminate the effects of the correlations between two degrees and focus on the role of their variability on PoEs, we assume that the two degrees of an agent are independent in this subsection.

Assumption 2

$p_{1}(d_{1},d_{2})=p_{1,1}(d_{1})\cdot p_{1,2}(d_{2})$* for all $d_{1},d_{2}\in{\mathbbm{Z}}_{+}$.

The following lemma illustrates how the variability in degrees affects the PoEs when the degrees of an agent are independent.

Lemma 1

Consider two degree distributions $p^{(\ell)}_{1}$, $\ell=1,2$. Let ${\bf D}^{(\ell)}_{1}$, $\ell=1,2,$ be a random vector with distribution $p^{(\ell)}_{1}$. Suppose that Assumption 2 holds for $p_{1}^{(\ell)}$, $\ell=1,2$, and $D^{(1)}_{1,j}\leq_{SSD}D^{(2)}_{1,j}$ for $j=1,2$. Then, $\bm{\mu}^{2}\leq\bm{\mu}^{1}$, where $\bm{\mu}^{\ell}<{\bf 1}$, $\ell=1,2$, is the PoE vector under degree distribution $p^{(\ell)}_{1}$.

Proof:

A proof is provided in Section 8.1. ∎

The lemma tells us that $\mu^{2}_{1}=\mu^{2}_{2}\leq\mu^{1}_{1}=\mu^{1}_{2}$. Thus, an implication of Lemma 1 is that even when the mean degrees of agents are fixed, the PoEs tend to increase as the degrees of agents become more spread out, i.e., have greater variability, suggesting that widespread failures would be less likely as the degrees of agents vary more widely.

We say that $X$ is smaller than $Y$ with respect to first-order stochastic dominance (FSD) or usual stochastic order if, for all increasing functions $\xi:{\rm I\kern-2.5ptR}\to{\rm I\kern-2.5ptR}$, ${\mathbbm{E}}\left[\xi(X)\right]\leq{\mathbbm{E}}\left[\xi(Y)\right]$ [32, 40]. This is equivalent to $F_{X}(t)\geq F_{Y}(t)$ for all $t\in{\rm I\kern-2.5ptR}$. Clearly, by definition, FSD implies SSD. Hence, $D^{(1)}_{1,j}\leq_{FSD}D^{(2)}_{1,j}$, $j=1,2$, is a sufficient condition for Lemma 1 to hold and, as one would expect, when agents’ degrees become larger, an outbreak of failures is more likely.

Example 1: Consider two degree distributions $p^{(1)}_{1}$ and $p^{(2)}_{1}$ shown in Table I. Even though these distributions may not be realistic or representative, we use them to illustrate our findings with numerical examples. One can easily verify that (i) $D^{(\ell)}_{1,1}$ and $D^{(\ell)}_{1,2}$, $\ell=1,2$, are independent and (ii) $D^{(1)}_{1,j}\leq_{SSD}D^{(2)}_{1,j}$, $j=1,2$. In addition, both distributions yield

[TABLE]

with $\rho({\mathbbm{M}})=1.021>1$. The entropies of $D^{(1)}_{1}$ and $D^{(1)}_{2}$ (resp. $D^{(2)}_{1}$ and $D^{(2)}_{2}$) are 1.786 and 1.003 (resp. 1.461 and 0.884), respectively, suggesting that $D^{(1)}_{j}$, $j=1,2$, are more unpredictable than $D^{(2)}_{j}$, $j=1,2$.

The PoE vector satisfying (16) for $p^{(1)}_{1}$ (resp. $p^{(2)}_{1}$) is $\bm{\mu}^{1}=(0.9646\ 0.9646\ 0.9761\ 0.9761)$ (resp. $\bm{\mu}^{2}=(0.9586\ 0.9586\ 0.9720\ 0.9720)$). Thus, although the two degree distributions yield the same matrix ${\mathbbm{M}}$, the PoE is larger under distribution $p^{(1)}_{1}$. Equivalently, the PoCF, beginning with a random failure in either CS, is 0.0354 (resp. 0.0414) under $p^{(1)}_{1}$ (resp. $p^{(2)}_{1}$), which represents roughly a 17 percent difference in PoCF.

6.2 The effects of degree dependence

We now turn our attention to the role of dependence between the two degrees of an agent. To this end, we adopt a well-known dependence order, called concordance order (CO) [32]: suppose that ${\bf X}=(X_{1},X_{2})$ and ${\bf Y}=(Y_{1},Y_{2})$ are two bivariate rvs with identical marginal distributions. This means that the variability of each rv remains fixed. Then, ${\bf X}$ is smaller than ${\bf Y}$ in CO $({\bf X}\leq_{C}{\bf Y})$ if Cov($\xi_{1}(X_{1}),\xi_{2}(X_{2})$) $\leq$ Cov($\xi_{1}(Y_{1}),\xi_{2}(Y_{2})$) for all increasing functions $\xi_{j}$, $j=1,2$. Note that this implies Cov($X_{1}$, $X_{2}$) $\leq$ Cov($Y_{1},Y_{2}$).

Roughly speaking, ${\bf X}\leq_{C}{\bf Y}$ means that $Y_{i}$, $i=1,2$, are more positively correlated than $X_{i}$, $i=1,2$. In addition, as explained in [32, p. 109], CO is the only integral stochastic order that satisfies natural properties that one would expect of a stochastic order for comparing dependence.

The second lemma examines how the (positive) dependence of degrees influences the PoEs. Its proof is omitted here due to a space constraint and can be found in [25].

Lemma 2

Consider two degree distributions $p^{(\ell)}_{1}$, $\ell=1,2$, with identical marginal distributions. Let ${\bf D}^{(\ell)}_{1}$, $\ell=1,2,$ be a random vector with distribution $p^{(\ell)}_{1}$. Suppose ${\bf D}^{(1)}_{1}\leq_{C}{\bf D}^{(2)}_{1}$. Then, $\bm{\mu}^{1}\leq\bm{\mu}^{2}$, where $\bm{\mu}^{\ell}$, $\ell=1,2$, is the PoE vector under the degree distribution $p^{(\ell)}_{1}$.

A key finding of Lemma 2 is that as the two degrees of agents become more positively correlated, it becomes more difficult to set off cascading failures. One possible way to interpret this finding is that as the degrees become more positively correlated, the variability in the total degree of an agent, i.e., the sum of two degrees, also grows. Hence, Lemma 1 suggests that the PoEs should increase.

**Example 2: ** For the second example, consider degree distributions $p^{(2)}_{1}$ and $p^{(3)}_{1}$ given in Table I. To obtain $p^{(3)}_{1}$, we modified $p^{(2)}_{1}$ in order to introduce weak positive correlations between the two degrees by (i) adding 0.007 to $p^{(2)}_{1}(1,0)$ and $p^{(2)}_{1}(3,2)$ and (ii) subtracting 0.007 from $p^{(2)}_{1}(3,0)$ and $p^{(2)}_{1}(1,2)$. The correlation coefficient of $D^{(3)}_{1,1}$ and $D^{(3)}_{1,2}$ is 0.0368, indicating weak positive correlations.

One can show that (i) $p^{(2)}_{1}$ and $p^{(3)}_{1}$ have the same marginal distributions and (ii) ${\bf D}^{(2)}_{1}\leq_{C}{\bf D}^{(3)}_{1}$. Therefore, Lemma 2 tells us that $\bm{\mu}^{2}\leq\bm{\mu}^{3}$. Indeed, $\bm{\mu}^{3}=(0.9604\ 0.9604\ 0.9732\ 0.9732)$, which is larger than $\bm{\mu}^{2}$ from the previous example. Accordingly, the PoCF decreases from 0.0414 to 0.0396, which represents approximately a 4.5 percent reduction in PoCF, despite very weak correlations in $p^{(3)}_{1}$; although $p^{(2)}_{1}$ and $p^{(3)}_{1}$ are close (with Kullback-Leibler divergence $D_{{\rm KL}}(p^{(2)}_{1}||p^{(3)}_{1})=$ 0.0094), the likelihood of experiencing widespread failures changes somewhat noticeably. This points to possible sensitivity of PoEs to the degree distributions, including their dependence structure, in some cases.

7 General Settings

In Section 6, we considered scenarios with two CSes and deterministic transmission of infections, and studied how the variability of degrees and dependence between the two degrees of agents alter the PoEs. In this section, we return to the general settings described in Section 3 and examine how the degree distributions $(p_{i};i\in{\mathcal{S}})$ shape the PoEs.

General settings pose additional challenges that we did not have to cope with in the simpler two-CS scenarios. First, unlike in univariate or bivariate cases, choosing a suitable stochastic order for comparing degree distributions becomes more problematic. The reason for this is that there are several different stochastic orders one can consider, which can be viewed as extensions of a single stochastic order for univariate rvs to random vectors. Second, perhaps more importantly, if the infection probabilities $q_{i,j}$, $i,j\in{\mathcal{S}}$, are not equal to one, even when two different sets of degree distributions $(p^{(\ell)}_{i};\ i\in{\mathcal{S}})$, $\ell=1,2$, can be ordered using some stochastic order, the associated children distributions $(h^{(\ell)}_{i};\ i\in{\mathcal{S}}^{+})$ are in general not guaranteed to preserve the ordering with respect to the same stochastic order.

Consider two sets of degree distributions ${\bf p}^{(\ell)}=(p^{(\ell)}_{i};\ i\in{\mathcal{S}})$, $\ell=1,2$. Let ${\bf D}^{(\ell)}_{i}=(D^{(\ell)}_{i,j};\ j\in{\mathcal{S}})$, $\ell=1,2$, and $i\in{\mathcal{S}}$, be a random vector with distribution $p^{(\ell)}_{i}$. In order to make progress, we introduce the following assumption on internal failure probability functions $\wp^{in}_{i}$, $i\in{\mathcal{S}}$.

Assumption 3

Assume that $\phi_{i}^{\star}:{\rm I\kern-1.5ptN}\to[0,\infty)$, where $\phi_{i}^{\star}(d):=d\cdot\wp^{in}_{i}(d,1)$, is non-decreasing and satisfies $\alpha\cdot\phi^{\star}_{i}(d_{1})+(1-\alpha)\phi^{\star}_{i}(d_{3})\leq\phi^{\star}_{i}(d_{2})$ for all $d_{1}\leq d_{2}\leq d_{3}$, where $\alpha=(d_{3}-d_{2})/(d_{3}-d_{1})$.

Roughly speaking, Assumption 3 states that an agent has a higher total risk of experiencing a failure (due to the failure of one of the internal neighbors) with an increasing internal degree. An example of $\wp^{in}_{i}$ that satisfies Assumption 3 is $\wp^{in}_{i}(d,1)=d^{-\beta}$ with $\beta\in[0,1]$. Obviously, $\wp^{in}_{i}(d,1)=\sum_{k}c_{k}\cdot d^{-\beta_{k}}$, where $c_{k}>0$ and $\beta_{k}\in[0,1]$ for all $k$, also satisfies the assumption. We point out that, as illustrated by the example, this assumption captures the heightened aggregate risk of failure seen by higher degree agents, even though a single neighbor might pose less risk (i.e., smaller $\wp^{in}_{i}(d,1)$).

7.1 The effects of dependence

The dependence order we adopted in the previous section, namely CO, can be generalized to random vectors consisting of more than two rvs: suppose that ${\bf X}$ and ${\bf Y}$ are $n$-dimensional random vectors with $n>2$. Then, ${\bf X}$ is said to be smaller than ${\bf Y}$ in CO if $F_{{\bf X}}({\bf t})\leq F_{{\bf Y}}({\bf t})$ and $\bar{F}_{{\bf X}}({\bf t})\leq\bar{F}_{{\bf Y}}({\bf t})$ for all ${\bf t}\in{\rm I\kern-2.5ptR}^{n}$, where $F_{{\bf X}}={\mathbbm{P}}\left[{\bf X}\leq{\bf t}\right]$ and $\bar{F}_{{\bf X}}={\mathbbm{P}}\left[{\bf X}>{\bf t}\right]$. One can verify that these conditions imply that the marginal distributions are identical.

7.1.1 Supermodular order

In our study, we instead consider a dependence order that is somewhat stronger than CO. This is called supermodular order (SMO) [32]: a function $\xi:{\rm I\kern-2.5ptR}^{n}\to{\rm I\kern-2.5ptR}$ is called supermodular if, for all ${\bf x},{\bf y}\in{\rm I\kern-2.5ptR}^{n}$,

[TABLE]

where ${\bf x}\wedge{\bf y}$ and ${\bf x}\vee{\bf y}$ are $(\min(x_{i},y_{i});\ i=1,2,\ldots,n)$ and $(\max(x_{i},y_{i});\ i=1,2,\ldots,n)$, respectively. If the inequality in (18) goes the other way, the function $\xi$ is called submodular.

A random vector ${\bf X}$ is smaller than a random vector ${\bf Y}$ in SMO (${\bf X}\leq_{sm}{\bf Y}$) if ${\mathbbm{E}}\left[\xi({\bf X})\right]\leq{\mathbbm{E}}\left[\xi({\bf Y})\right]$ for all supermodular functions $\xi$ [32]. SMO is a multivariate positive dependence order that satisfies the nine natural properties discussed in [32, pp. 110-111]. Furthermore, for bivariate cases, CO and SMO are equivalent. For $n>2$, however, SMO implies CO, but they are no longer equivalent. Finally, if ${\bf X}\leq_{sm}{\bf Y}$, we have Cov($X_{i},X_{j}$) $\leq$ Cov($Y_{i},Y_{j}$) for all $i\neq j$.

The following theorem generalizes Lemma 2. We will defer the proof of the theorem till after Theorem 2 in the subsequent subsection.

Theorem 1

*Suppose that Assumption 3 holds and ${\bf D}^{(1)}_{i}\leq_{sm}{\bf D}^{(2)}_{i}$ for all $i\in{\mathcal{S}}$. Let $\bm{\mu}^{\ell}$, $\ell=1,2$, be the PoE vector under the degree distributions $(p^{(\ell)}_{i};\ i\in{\mathcal{S}})$. Then, we have $\bm{\mu}^{1}\leq\bm{\mu}^{2}$.

7.2 The effects of variability

In Section 6 with bivariate degrees, we assumed independence of two degrees and studied the influence of variability of each degree on the PoEs with the help of SSD (or ICV order). In this section, we adopt a stochastic order that allows us to examine the impact of variability with a common dependence structure captured by what is called copula [33], without having to assume independence.

7.2.1 Copulas

Suppose that ${\bf X}$ is an $n$-dimensional random vector with a joint distribution function $F_{{\bf X}}$. A copula of ${\bf X}$ (or associated with $F_{{\bf X}}$) is a function $C_{{\bf X}}:[0,1]^{n}\to[0,1]$ satisfying

[TABLE]

where $F_{X_{i}}$ is the marginal distribution of $X_{i}$, $i=1,2\ldots,n$. For instance, a copula of mutually independent $X_{i}$, $i=1,2\ldots,n$, is a product function, i.e., for all ${\bf u}\in[0,1]^{n}$, we have $C_{{\bf X}}({\bf u})=\prod_{i=1}^{n}u_{i}$. Also, if $F_{X_{i}}$, $i=1,2,\ldots,n$, are continuous, there is a unique copula that satisfies (19).

It is clear from (19) that a copula of a random vector captures most of the dependence structure properties that do not depend on the marginal distributions. In this sense, two random vectors with a common copula have similar dependence structure among the comprising rvs. For a more detailed discussion of copulas, we refer an interested reader to a manuscript by Nelson [33].

7.2.2 Increasing directionally concave order

A function $\xi:{\rm I\kern-2.5ptR}^{n}\to{\rm I\kern-2.5ptR}$ is said to directionally concave (DCV) if, for all ${\bf x}_{i}\in{\rm I\kern-2.5ptR}^{n}$, $i=1,2,3,4$, with ${\bf x}_{1}\leq{\bf x}_{2},{\bf x}_{3}\leq{\bf x}_{4}$ and ${\bf x}_{1}+{\bf x}_{4}={\bf x}_{2}+{\bf x}_{3}$, we have $\xi({\bf x}_{1})+\xi({\bf x}_{4})\leq\xi({\bf x}_{2})+\xi({\bf x}_{3})$. If the inequality goes the other way, the function is called directionally convex (DCX). Clearly, $\xi$ is DCV iff $-\xi$ is DCX. It turns out that a function $\xi$ is DCV (resp. DCX) iff it is submodular and componentwise concave (resp. supermodular and componentwise convex) [40, p. 335].

Let ${\bf X}$ and ${\bf Y}$ be two $n$-dimensional random vectors. Random vector ${\bf X}$ precedes ${\bf Y}$ in DCV order (${\bf X}\leq_{dcv}{\bf Y}$) if ${\mathbbm{E}}\left[\xi({\bf X})\right]\leq{\mathbbm{E}}\left[\xi({\bf Y})\right]$ for all DCV functions $\xi$. Note that ${\bf X}\leq_{dcv}{\bf Y}$ iff ${\bf Y}\leq_{dcx}{\bf X}$. If the inequality is required only for increasing DCV functions, we say that ${\bf X}$ precedes ${\bf Y}$ in increasing DCV (IDCV) order. As expected, ${\bf X}\leq_{idcv}{\bf Y}$ iff ${\bf Y}$ precedes ${\bf X}$ in decreasing DCX (DDCX) order (${\bf Y}\leq_{ddcx}{\bf X}$).

If ${\bf X}\leq_{idcv}{\bf Y}$, then $X_{i}\leq_{icv}Y_{i}$ for all $i=1,2,\ldots,n$. Therefore, $Y_{i}$, $i=1,2,\ldots,n$, are in a way more predictable than $X_{i}$, $i=1,2,\ldots,n$. As pointed out in [32, p. 135], the DCV (or DCX) order goes one step further and allows us to compare random vectors with a common copula, but with different variability in the marginals. Moreover, an example is provided to illustrate that convex order is not suitable for this purpose.

Utilizing the IDCV order (or, equivalently, DDCX order), the following theorem sheds some light on how the variability of agents’ degrees influences the PoEs, even when the mean degrees stay fixed.

Theorem 2

Suppose that Assumption 3 holds and ${\bf D}^{(1)}_{i}\leq_{idcv}{\bf D}^{(2)}_{i}$ with ${\mathbbm{E}}\left[{\bf D}^{(1)}_{i}\right]={\mathbbm{E}}\left[{\bf D}^{(2)}_{i}\right]$ for all $i\in{\mathcal{S}}$. Let $\bm{\mu}^{\ell}$, $\ell=1,2$, be the PoE vector under the distributions $(p^{(\ell)}_{i};\ i\in{\mathcal{S}})$. Then, we have $\bm{\mu}^{2}\leq\bm{\mu}^{1}$.

Proof:

A proof of Theorem 2 is given in Section 8.2. ∎

As explained before, when ${\bf D}^{(1)}_{i}\leq_{idcv}{\bf D}^{(2)}_{i}$ for some $i\in{\mathcal{S}}$, the degrees of CS $i$ agents with the distribution $p^{(1)}_{i}$ have greater variability than with the distribution $p^{(2)}_{i}$. As a result, Theorem 2 can be viewed as a generalization of Lemma 1.

We are now ready to provide the proof of Theorem 1.

Proof of Theorem 1: Recall that a function is DCX iff it is both supermodular and componentwise convex. Thus, it is obvious that SMO implies DCX order, hence DDCX order. For this reason, if ${\bf D}^{(2)}_{i}\leq_{sm}{\bf D}^{(1)}_{i}$, then ${\bf D}^{(2)}_{i}\leq_{ddcx}{\bf D}^{(1)}_{i}$ or, equivalently, ${\bf D}^{(1)}_{i}\leq_{idcv}{\bf D}^{(2)}_{i}$. Now Theorem 1 follows from Theorem 2.

7.3 *Comparison on the basis of

children distributions*

In Theorems 1 and 2, the inequalities in stochastic orders are imposed on the degree distributions $(p^{(\ell)}_{i};\ i\in{\mathcal{S}})$. In some cases, however, it may be possible to estimate the children distribution. If we could directly compare the children distributions (as we implicitly did in Section 6 for two symmetric CSes), we can prove a stronger result than Theorem 2.

Suppose that ${\bf X}$ and ${\bf Y}$ are two $n$-dimensional random vectors. We say that ${\bf X}$ is smaller than ${\bf Y}$ in Laplace transform (LT) order (${\bf X}\leq_{LT}{\bf Y}$) if ${\mathbbm{E}}\left[\exp\left(-{\bf s}{\bf X}^{T}\right)\right]\geq{\mathbbm{E}}\left[\exp\left(-{\bf s}{\bf Y}^{T}\right)\right]\ \mbox{ for all }{\bf s}>{\bf 0}$.

One can easily show that, for every ${\bf s}>{\bf 0}$, the function $\xi_{{\bf s}}:{\rm I\kern-2.5ptR}^{n}\to{\rm I\kern-2.5ptR}$, where $\xi_{{\bf s}}({\bf x})=\exp(-{\bf s}{\bf x}^{T})$ for ${\bf x}\in{\rm I\kern-2.5ptR}^{n}$, is DDCX: since $\xi_{{\bf s}}$ is twice differentiable, it is DCX iff $\frac{\partial^{2}}{\partial x_{i}\partial x_{j}}\xi_{{\bf s}}({\bf x})\geq 0$ for all ${\bf x}\in{\rm I\kern-2.5ptR}^{n}$ [32, Theorem 3.12.2, p. 132]. For all $i,j\in\{1,2,\ldots,n\}$, we have

[TABLE]

Clearly, $\xi_{{\bf s}}({\bf x})$ is decreasing in ${\bf x}$ and, hence, is DDCX. This tells us that if ${\bf X}\leq_{ddcx}{\bf Y}$, then ${\bf X}\leq_{LT}{\bf Y}$.

Suppose that $(h^{(\ell)}_{i};\ i\in{\mathcal{S}}^{+})$, $\ell=1,2$, are the children distributions under degree distributions $(p^{(\ell)}_{i};\ i\in{\mathcal{S}})$. Let ${\bf O}^{(\ell)}_{i}=(O^{(\ell)}_{i,j};\ j\in{\mathcal{S}}^{+})$ be the random children vector with distribution $h^{(\ell)}_{i}$, $\ell=1,2$, and $i\in{\mathcal{S}}^{+}$. The following theorem holds under a much weaker condition (namely, LT order) than IDCV order without Assumption 3.

Theorem 3

*Suppose that ${\bf O}^{(1)}_{i}\leq_{LT}{\bf O}^{(2)}_{i}$ for all $i\in{\mathcal{S}}^{+}$. Let $\bm{\mu}^{\ell}$, $\ell=1,2$, be the PoE vector under the children distributions $(h^{(\ell)}_{i};\ i\in{\mathcal{S}}^{+})$. Then, $\bm{\mu}^{2}\leq\bm{\mu}^{1}$.

Proof:

Please see Section 8.4 for a proof. ∎

7.4 Discussion

In this subsection, we briefly discuss some of modeling assumptions and their roles.

$\bullet$ Asymmetric inter-CS dependence and symmetric intra-CS dependence – In our model outlined in Section 3, we purposely assumed somewhat different intra-CS dependence and inter-CS dependence among agents. More specifically, we assumed that intra-CS dependence is symmetric, while inter-CS dependence is asymmetric.

Suppose that the inter-CS dependence is symmetric instead. Then, the inter-CS infection probabilities $q_{i,j}$, $i\neq j$, can be computed in a manner similar to that of $q_{i,i}$, $i\in{\mathcal{S}}$, as outlined in Section 4.2. Since Theorems 1 - 3 hold with symmetric internal dependence, it is not surprising that the same results hold with symmetric inter-CS dependence if $\wp^{ex}_{j,i}(d,1)$, $j\neq i$, satisfy Assumption 3. Similarly, if the intra-CS dependence is asymmetric, under the same assumption that the conditional degree distributions of failed agents are similar to the prior degree distributions $p_{i}$, $i\in{\mathcal{S}}$, the results still hold. These observations suggest that our findings are true under more general settings and are not sensitive to specific modeling assumptions introduced in the study.

$\bullet$ Effects of failure probability functions – Recall that the functions $\wp^{in}_{i}$ and $\wp^{ex}_{j,i}$ are used to model the vulnerability of agents to the failures of those that support them. It is clear from Section 4.2 that, except for Assumption 3 on $\wp^{in}_{i}$, our main findings in Theorems 1 - 3 do not impose any other conditions on $\wp^{in}_{i}$ and $\wp^{ex}_{j,i}$. In particular, although it is reasonable to expect agents with larger incoming degrees to be less susceptible to the failure of a single supporting agent, the functions $\wp^{ex}_{j,i}$ need not be monotonic. This suggests that our reported findings are insensitive to the exact choices of these functions as long as Assumption 3 is met by $\wp^{in}_{i}$, $i\in{\mathcal{S}}$.

But, we also note that this may no longer be true if the dependence graph exhibits assortativity, i.e., the degrees of neighbors are correlated. We do not study the impact of assortativity here, and refer an interested reader to [23] for a study of its influence in a single homogeneous network with strategic agents.

8 Proofs of Main results

This section provides the proofs of our main results. A reader who is not interested in the proofs can safely skip the section.

8.1 Proof of Lemma 1

Let ${\bf f}^{\ell}$, $\ell=1,2$, be the generating functions corresponding to the degree distributions $p^{(\ell)}_{1}$, $\ell=1,2$. Using the definition of the generating function, we obtain

[TABLE]

The expressions for $f_{2}({\bf s})$ and $f_{4}({\bf s})$ can be easily obtained using the assumed symmetry of degree distributions $p_{1}$ and $p_{2}$. Recall that, from the assumed symmetry of the degree distributions, we know $\mu_{1}=\mu_{2}$ and $\mu_{3}=\mu_{4}$.

Substituting $\mu_{1}$ for $\mu_{2}$ in (20) for the PoEs and using the assumed independence of the two degrees of an agent,

[TABLE]

Clearly, for any fixed $s\in(0,1)$, $s^{x}$ is a decreasing, convex function of $x$ or, equivalently, $-s^{x}$ is an increasing, concave function of $x$. Because $D^{(1)}_{1,j}\leq_{SSD}D^{(2)}_{1,j}$, $j=1,2,$ we have, for all $\bm{\mu}<{\bf 1}$,

[TABLE]

Repeating the same steps, starting with (21), yields

[TABLE]

Since these two inequalities hold for all $\bm{\mu}<{\bf 1}$, we get

[TABLE]

Corollary 2 [15, p. 42] tells us that the inequality ${\bf f}^{2}(\bm{\mu}^{1})\leq\bm{\mu}^{1}$ in (22) means $\bm{\mu}^{2}\leq\bm{\mu}^{1}$. This completes the proof.

8.2 Proof of Theorem 2

In order to prove the theorem, we will make use of the following lemma. Its proof is provided in Section 8.3.

Lemma 3

*Under Assumption 3, we have ${\bf O}^{(1)}_{i}\leq_{idcv}{\bf O}^{(2)}_{i}$ for all $i\in{\mathcal{S}}^{+}$, where ${\bf O}^{(\ell)}_{i}$, $\ell=1,2$, is a random vector with the children distribution $h^{(\ell)}_{i}$.

We know from Corollary 1 [15, p. 42] that the only solutions of (16) in the unit cube are $\bm{\mu}$ and ${\bf 1}$. However, when $\rho({\mathbbm{M}})>1$, we have $\bm{\mu}<{\bf 1}$.

Let ${\bf f}^{\ell}$ be the generating function under the children distributions ${\bf h}^{(\ell)}=(h^{(\ell)}_{i};\ i\in{\mathcal{S}}^{+})$, $\ell=1,2$. Then, ${\bf f}^{1}(\bm{\mu}^{1})=\bm{\mu}^{1}<{\bf 1}$. From the definition of the generating function,

[TABLE]

Since ${\bf O}^{(1)}_{i}\leq_{idcv}{\bf O}^{(2)}_{i}$ or, equivalently, ${\bf O}^{(2)}_{i}\leq_{ddcx}{\bf O}^{(1)}_{i}$ for all $i\in{\mathcal{S}}^{+}$, we can show ${\bf f}^{2}({\bf s})\leq{\bf f}^{1}({\bf s})$ for all ${\bf s}<{\bf 1}$ as follows: for any ${\bf 0}<{\bf s}<{\bf 1}$, define $\psi_{{\bf s}}:{\rm I\kern-2.5ptR}^{2N_{S}}\to{\rm I\kern-2.5ptR}$, where $\psi_{{\bf s}}({\bf x})=\exp\left(-\sum_{i\in{\mathcal{S}}^{+}}\log(s_{i}^{-1})\ x_{i}\right).$ According to Theorem 3.12.2 [32, p. 132], $\psi_{{\bf s}}$ is DCX iff $\partial^{2}\psi_{{\bf s}}({\bf x})/\partial x_{i}\partial x_{j}\geq 0$ for all ${\bf x}\in{\rm I\kern-2.5ptR}^{2N_{S}}$ and $i,j\in{\mathcal{S}}^{+}$. Clearly,

[TABLE]

for all ${\bf x}\in{\rm I\kern-2.5ptR}^{2N_{S}}$ and $i,j\in{\mathcal{S}}^{+}$, and $\psi_{{\bf s}}$ is DCX.

Together with ${\bf O}^{(2)}_{i}\leq_{ddcx}{\bf O}^{(1)}_{i}$ and $\bm{\mu}^{1}<{\bf 1}$, this tells us that, for all $i\in{\mathcal{S}}^{+}$,

[TABLE]

and ${\bf f}^{2}(\bm{\mu}^{1})\leq\bm{\mu}^{1}$. Corollary 2 [15, p. 42] states that, if ${\bf 0}\leq{\bf f}^{2}(\bm{\mu}^{1})\leq\bm{\mu}^{1}\leq{\bf 1}$, then $\bm{\mu}^{2}\leq\bm{\mu}^{1}$. This completes the proof of the theorem.

8.3 Proof of Lemma 3

We will first prove the lemma for $i\in{\mathcal{S}}$ and then for $i\in{\mathcal{S}}^{+}\setminus{\mathcal{S}}=\{N_{S}+1,\ldots,2N_{S}\}$.

$\bullet$ ${\bf O}^{(1)}_{i}\leq_{idcv}{\bf O}^{(2)}_{i}$, $i\in{\mathcal{S}}$: We prove the claim using the definition of IDCV order. Suppose that $g:{\rm I\kern-2.5ptR}^{2N_{S}}\to{\rm I\kern-2.5ptR}$ is IDCV. Recall that this is equivalent to saying that $-g$ is DDCX. For notational convenience, for each $j\in{\mathcal{S}}$, define

[TABLE]

For $\ell=1,2$,

[TABLE]

where the mapping $\bar{{\bf d}}$ was defined in (9), $q^{\ell}_{i,j}=q_{i,j}$ for $j\neq i$, and $q^{\ell}_{i,i}$ is the probability that an internal neighbor of a CS $i$ agent will fall victim to infection, under the internal degree distribution $p^{(\ell)}_{i,i}$, which is given in (1).

We shall prove $\mathbbm{E}\left[g({\bf O}^{(1)}_{i})\right]\leq\mathbbm{E}\left[g({\bf O}^{(2)}_{i})\right]$ in two steps. First, we will show

[TABLE]

Then, we will demonstrate that $(\ref{eq:appen2-0})\leq\mathbbm{E}\left[g({\bf O}^{(2)}_{i})\right]$. If this is true for all IDCV functions $g$, it implies ${\bf O}^{(1)}_{i}\leq_{idcv}{\bf O}^{(2)}_{i}$ by definition.

After interchanging the order of the two summations in (8.3), we get

[TABLE]

Note that $g^{*}_{\ell}({\bf d})$ is $\mathbbm{E}\left[g({\bf Y}({\bf d},\ell))\right]$, where ${\bf Y}({\bf d},\ell)$ is a vector consisting of $N_{S}$ mutually independent binomial rvs. In other words, with a little abuse of notation, the $j$-th element $Y_{j}(d_{j},\ell)$ is a Binomial($d_{j},q^{\ell}_{i,j}$) rv and, hence, can be viewed as a sum of $d_{j}$ i.i.d. Bernoulli($q^{\ell}_{i,j}$) rvs.

In order to finish the proof of the first step, we make use of the following lemma.

Lemma 2.17 [28]: Let $\left(S^{(i)}_{j},j=1,2,...\right)$ be independent sequences of i.i.d. nonnegative rvs. Suppose $\theta:{\rm I\kern-2.5ptR}^{m}\to{\rm I\kern-2.5ptR}$ is increasing DCX (resp. IDCV). Then, $\zeta:{\rm I\kern-1.5ptN}^{m}\to{\rm I\kern-2.5ptR}$ defined by $\zeta({\bf n})={\mathbbm{E}}\left[\theta\left(\sum_{j=1}^{n_{1}}S^{(1)}_{j},\ldots,\sum_{j=1}^{n_{m}}S^{(m)}_{j}\right)\right]$ is increasing DCX (resp. IDCV).

Since the function $g$ was assumed IDCV and ${\bf Y}({\bf d},1)$ consists of mutually independent binomial rvs (each of which is a sum of i.i.d. Bernoulli rvs), the above lemma tells us that $g^{*}_{1}$ is IDCV. Since $g^{*}_{1}$ is IDCV and, with a little abuse of notation, $p^{(1)}_{i}\leq_{idcv}p^{(2)}_{i}$, we have

[TABLE]

This proves the first step.

In order to prove the second step

[TABLE]

it suffices to show that $q^{1}_{i,j}\leq q^{2}_{i,j}$ for all $j\in{\mathcal{S}}$; if this is true, for all ${\bf d}\in{\mathbbm{Z}}_{+}^{N_{S}}$, $Y_{j}(d_{j},1)\leq_{st}Y_{j}(d_{j},2)$ for all $j\in{\mathcal{S}}$, where $\leq_{st}$ denotes the inequality with respect to usual stochastic order. Since $Y_{j}(d_{j},\ell)$, $j\in{\mathcal{S}}$, are mutually independent, Theorem 3.3.8 [32, p. 93] tells us ${\bf Y}({\bf d},1)\leq_{st}{\bf Y}({\bf d},2)$. Because the function $g$ is assumed increasing, this implies

[TABLE]

from the definition of usual stochastic order [32].

First, recall that $q^{1}_{i,j}=q^{2}_{i,j}=q_{i,j}$ for all $j\neq i$ because they are assumed fixed. Thus, we only need to show $q^{1}_{i,i}\leq q^{2}_{i,i}$. Using the definition in (1),

[TABLE]

Because $\phi_{i}^{*}$ satisfies Assumption 3 and $p^{(1)}_{i}\leq_{idcv}p^{(2)}_{i}$, which implies $p^{(1)}_{i,i}\leq_{icv}p^{(2)}_{i,i}$, we obtain

[TABLE]

Recall that $d^{(1)}_{i,{\rm avg}}=d^{(2)}_{i,{\rm avg}}$ from the assumption in the theorem. This proves $\mathbbm{E}\left[g({\bf O}^{(1)}_{i})\right]\leq\mathbbm{E}\left[g({\bf O}^{(2)}_{i})\right]$. Since this inequality holds for every IDCV function $g$ (with well defined expectations), we have ${\bf O}^{(1)}_{i}\leq_{idcv}{\bf O}^{(2)}_{i}$.

$\bullet$ ${\bf O}^{(1)}_{i}\leq_{idcv}{\bf O}^{(2)}_{i}$, $i\in{\mathcal{S}}^{+}\setminus{\mathcal{S}}$: First, note that if $\varphi:{\rm I\kern-2.5ptR}^{2N_{S}}\to{\rm I\kern-2.5ptR}$ is IDCV, then so is $\varphi^{*}_{{\bf e}}:{\rm I\kern-2.5ptR}^{2N_{S}}\to{\rm I\kern-2.5ptR}$ with $\varphi^{*}_{{\bf e}}({\bf x})=\varphi({\bf x}+{\bf e})$ for all ${\bf e}\in{\mathbbm{Z}}^{2N_{S}}$. This follows directly from the observation that $\varphi^{*}_{{\bf e}}$ satisfies the characterization (ii) and (iii) of DCV functions in Theorem 3.12.2 [32, p.132].

This tells us that $\tilde{p}^{(\ell)}_{i}$, $\ell=1,2$, with $\tilde{p}^{(\ell)}_{i}({\bf d}):=p^{(\ell)}_{i}({\bf d}+{\bf e}_{i})$ satisfy $\tilde{p}^{(1)}_{i}\leq_{idcv}\tilde{p}^{(2)}_{i}$. The claim that ${\bf O}^{(1)}_{i}\leq_{idcv}{\bf O}^{(2)}_{i}$ now follows from the proof of the previous case.

8.4 Proof of Theorem 3

From the definition of $\bm{\mu}^{2}$, we have ${\bf f}^{2}(\bm{\mu}^{2})=\bm{\mu}^{2}<{\bf 1}$. Using the given children distributions, for all $i\in{\mathcal{S}}$,

[TABLE]

where $\pi^{2}_{j}=\log(1/\mu^{2}_{j})>0$, $j\in{\mathcal{S}}^{+}$. Similarly, for all $i^{+}\in{\mathcal{S}}^{+}\setminus{\mathcal{S}}$,

[TABLE]

Recall that ${\bf O}^{(1)}_{i}\leq_{LT}{\bf O}^{(2)}_{i}$ for all $i\in{\mathcal{S}}^{+}$. Thus, because $\pi^{2}_{j}$, $j\in{\mathcal{S}}^{+}_{i}$, are all positive, we get

[TABLE]

which yields $\mu^{2}_{i}\leq f^{1}_{i}(\bm{\mu}^{2})$ for all $i\in{\mathcal{S}}$. Following the same argument starting with (26) gives $\mu^{2}_{i^{+}}\leq f^{1}_{i^{+}}(\bm{\mu}^{2})$ for all $i^{+}\in{\mathcal{S}}^{+}\setminus{\mathcal{S}}$. Together, we obtain $\bm{\mu}^{2}\leq{\bf f}^{1}(\bm{\mu}^{2})$. Corollary 2 [15, p. 42] now tells us $\bm{\mu}^{2}\leq\bm{\mu}^{1}$, completing the proof of the theorem.

9 Conclusion

We investigated the impact of variability and correlations in degrees of agents on the robustness of interdependent systems. Our findings suggest that they both can have significant influence on the likelihoods of having catastrophic failures in complex systems comprising multiple heterogeneous systems via dependency among the agents. In particular, our results revealed that both increasing variability and positive dependence render the system more robust against random failures.

We are currently working to incorporate other graph properties displayed by both natural and engineered systems, such as assortativity and clustering, and to understand their role in the robustness of interdependent systems. Our goal is to identify a suitable way of imposing a partial order on the underlying dependence graphs and compare their resilience against both random and targeted attacks.