Probabilistic measures of edge criticality in graphs: a study in water distribution networks

The issue of vulnerability and robustness in networks have been addressed by several methods. The goal is to identify which are the critical components (i.e., nodes/edges) whose failure impairs the functioning of the network and how much this impacts the ensuing increase in vulnerability. In this paper we consider the drop in the network robustness as measured by the increase in vulnerability of the perturbed network and compare it with the original one. Traditional robustness metrics are based on centrality measures, the loss of efficiency and spectral analysis. The approach proposed in this paper sees the graph as a set of probability distributions and computes, specifically the probability distribution of its node to node distances and computes an index of vulnerability through the distance between the node-to-node distributions associated to original network and the one obtained by the removal of nodes and edges. Two such distances are proposed for this analysis: Jensen–Shannon and Wasserstein, based respectively on information theory and optimal transport theory, which are shown to offer a different characterization of vulnerability. Extensive computational results, including two real-world water distribution networks, are reported comparing the new approach to the traditional metrics. This modelling and algorithmic framework can also support the analysis of other networked infrastructures among which power grids, gas distribution and transit networks.

the assessment of how much this impacts the ensuing increase in vulnerability. The resilient operation of complex networks depends on their structural connectivity i.e., the existence of redundant paths between pair of nodes.
The most used analysis of vulnerability is carried out by studying different measures of the connectivity of the graph as expressed by centrality indices. According to a widely used metric (Latora and Marchiori 2007) an increase in vulnerability is the loss of efficiency as a consequence of the failure of a set of nodes/edges and their removal from the network. The structure and functions of the network strongly rely on the existence of paths between pair of nodes: when nodes and/or links are removed the length of such paths will increase and eventually some couples of nodes will become disconnected. A closely related analysis can be carried out using spectral graph theory: algebraic connectivity, given by the smallest positive eigenvalue of the Laplacian matrix of the network, was introduced in Fiedler (1973). The larger is the algebraic connectivity and the more difficult is to cut the network into disconnected components.
The work presented in this paper focuses on analyzing the criticality of links towards the overall vulnerability of the network by combining notions from graph theory with probabilistic analysis of dissimilarity between networks.
While previously quoted approaches are based on average values of shortest paths the approach proposed in this paper is based on the node-to-node distance probability distributions. These node related distributions are then aggregated into a network wise distribution. At this point the similarity between graphs can be associated to a distance between distributions. The 2 measures considered are the Jensen-Shannon (JS) divergence, based on Kullback-Leibler divergence, and the Wasserstein-1 (WST-1, aka the Earth-Mover) distance based on optimal transport theory. Extensive computational results, including two real-world water distribution networks, are reported comparing the distance-based approach with the traditional robustness metrics based on centrality measures, the loss of efficiency and spectral analysis.

Related works
The analysis of WDNs using measures and analytical tools from complex network theory has been, at the authors' knowledge, first suggested in the seminal paper (Yazdani and Jeffrey 2011) and further developed in Shuang et al. (2014), Archetti et al. (2015), Soldi et al. (2015). A closely related approach, spectral analysis, has been also used in Candelieri et al. (2017), Diao et al. (2016), Herrera et al. (2016), Maiolo et al. (2018) for assessing the resilience in WDNs.
Di Nardo et al. (2018) focuses on spectral techniques and shows how spectral metrics and algorithms support critical tasks of WDN management by just using topological and geometric information. More recently Ulusoy et al. (2018) proposed a hydraulically informed measure of criticality called water flow edge betweenness centrality (WFEBC), built upon on shortest-path and random walks betweenness measures (Herrera et al. 2016(Herrera et al. , 2015. Shuang et al. (2019) is a wide survey of quantitative resilience methods of WDN including network base approaches. Diao (2020) extends this analysis to multiscale resilience in water distribution and drainage systems.
In more general terms the very issue of differential robustness hinges on a notion of distance or similarity between networks or graphs. The distributional approach to network dissimilarity has been suggested in Schieber et al. (2017).
Entropy based methods like Kullback-Leibler, and Jensen-Shannon are most widely used but can run into problems when the distributions have different supports. This limitation might be relevant when a network has been derived deleting nodes/edges. For this reason, we have introduced in this paper the use of the Wasserstein distance. Wasserstein distance in particular has become important in image analysis and text analysis (Bonneel et al. 2016). More references will be given in "Conclusions and perspectives" section. Here we only quote Deza and Deza (2009) for a general survey of distances, Cover and Thomas (2006) as a general introduction to information theoretic arguments and Villani (2008) for a mathematical framework of optimal transport problems which the Wasserstein distance belongs to.

The contributions of this paper
The key element in this paper is that nodes and graphs are represented as probability distributions: this allows a richer representation of the structural properties of the network than that enabled by the average values of shortest paths and other connectivity metrics. The probability distribution associated to each node, of the node to node distance, is a signal of the relative importance of that node in the graph. The node related distributions can then be aggregated into a distributional representation of the whole network.
The main contribution of this paper is the introduction of a new set of vulnerability metrics given by the distance between the probability distributions of node-node distances between the original network and that resulting from the removal of nodes/ edges. These metrics allow the formulation of measures of criticality of single edges in the network.
Two distances have been analysed: Jensen-Shannon divergence based on information theoretic arguments and Wasserstein (WST) distance (aka the Earth-Mover) distance, based on optimal transport theory. The computational results confirm that the values of the distances JS and WST are strongly related to the criticality of the removed nodes/ edges and can capture numerically the impact of the deletion of a specific node/edge. A key advantage of the Wasserstein distance is that, under quite general conditions specified in Peyré and Cuturi (2019), it is a differentiable function of the parameters of the distributions which makes possible its use to assess the sensitivity of the network robustness to distributional perturbations.

Organization of the paper
The structure of the paper is as follows: "Entropy based measures of distance between networks" section describes entropy-based measures of distance, in particular JS, between networks. "Wasserstein distance" section describes distances, in particular WST, based on transport theory. "Experimental setting and computational results" section describes the modeling and algorithmic structure of the analysis framework proposed and the computational results on the networks used in this study. Finally, in "Conclusions and perspectives" section some conclusions and perspectives are provided.
To get a better focus on the main arguments, the background notions are contained in the appendices (A. graph basic definitions, network measures and spectral analysis, B. efficiency and vulnerability measures, and C. for the data and software resources).

Entropy based measures of distance between networks
Given a graph G(V , E) we associate to each node i = 1, . . . , n a discrete probability distribution as the fraction of nodes which are connected to i at a distance k. The distribution is displayed in Fig. 1. Figure 1a displays Anytown, used in the literature as a benchmark (Farmani et al. 2005): the associated graph G consists of 25 nodes and 44 edges.
The support of this distribution is 1, . . . , D(G) where D(G) is the diameter of G. The distance distribution over the whole network is represented by the global histogram in Fig. 2. (1) P i (k) =  The node-node distance distributions at network level of G and G ′ Let G′ be the graph without the red edge and p and p′ the distributions P G (k) and P G ′ (k) . For Anytown D(G) = 8 and the two histograms are displayed in Fig. 2.
The most widely used distance measure is the Kullback-Leibler (KL) divergence which has the drawback of being asymmetric and possibly infinite when there are points such that p = 0 and p′ > 0.
The Jensen-Shannon divergence is built on KL and is symmetric and always definite.
The use of Jensen-Shannon divergence in computing the dissimilarity between networks has been considered in (Schieber et al. 2017) along with the distance D JS :

Wasserstein distance
Wasserstein distance is a measure of the distance between two probability distributions. It is also called Earth Mover's (EM) distance from its informal interpretation as the minimum cost of moving and transforming a pile of sand in the shape of one probability distribution to the shape of the other distribution. The cost is quantified by the amount of sand moved times the moving distance. If the distribution domain is continuous the formula for the Earth-Mover (EM) distance is: where �(p, p′) denotes the set of all joint distributions γ x, y whose marginals are respectively p and p′ . One joint distribution γ x, y ∈ �(p, p′) describes one transport plan: intuitively γ x, y indicates how much mass must be transported from x to y in order to transform the distribution p into the distribution p′ . Therefore, the marginal distribution of over x adds up to p x γ x, y = p ′ y and analogously y γ x, y = p(x) . If x is the starting point and y the destination the total amount of sand moved is γ x, y and the traveling distance is x − y and thus the cost is γ x, y x − y .
The expected cost averaged over all the x, y pairs can be computed as: Finally, we take the minimum among the costs of all sand moving solutions as the EM distance: the EM distance is the cost of the optimal transport plan.
This paper is concerned with the Wasserstein metric for discrete distributionsspecifically one-dimensional histograms-in the Euclidean space. We focus on concepts related to the practical calculation of this metric. The Wasserstein metric is motivated by the classical optimal transportation problem first proposed in Monge (1781).
This problem received its modern linear programming formulation by Kantorovich (1942). There are two sets of points x i with i = 1, . . . , m and y j with j = 1, . . . , n . The cost of transporting one unit from x i to y j is given by c ij . The transport plan can be expressed in the form of an m × n matrix where the element γ ij represent the amount (of sand) transported from x i to y j . The central problem is to find the transport plan that minimizes the cost of total transportation cost which is the sum of the cost on the available roots m i=1 n j=1 c ij γ ij . Recently the Wasserstein distance has become a key tool in image processing (Bonneel et al. 2016) and machine learning (Frogner et al. 2015) and it has been used also for the generation of adversarial networks (Arjovsky et al. 2017;Weng 2019).
The formulation, of the WST distance for general probability measures requires sophisticated mathematical models (Villani 2008). Its computation raises challenging mathematical and computational problems (Peyré and Cuturi 2019).
In particular when the cost matrix is 0 on the diagonal and 1 elsewhere, Earth Mover distance is given by the l 1 norm (Manhattan distance) of the difference between histograms.
In the case of water distribution networks, the distributions of node-node distances are discrete and 1-dimensional, as shown in Fig. 2, and the computation of WST distance reduces to the comparison of two 1-D histograms: the corresponding Wasserstein distance can be computed very efficiently by using a simple sorting.
There are two key advantages of WST over JS: also in the cases when the distributions are supported in different spaces, even without overlaps, WST can still provide a meaningful representation of the distance between distributions. Furthermore, a key advantage of WST is its differentiability. Both points are illustrated in the following example (Fig. 3).
Let's consider Z = U (0, 1) the uniform distribution on the unit interval. Let P be the distribution of (0, Z) (0 on the x-axis and the random varable Z on the y axis and P θ = (θ , Z).
Therefore, Wasserstein provides a smooth measure which is useful for any optimization and learning process using gradient descent (Arjovsky et al. 2017).
When G′ is derived from G removing some edges, we have D(G′) ≥ D(G). Since we are in the particular case where the distributions are represented by histograms one can extend to G the support of G′ setting µ G (k) = 0 for k = D(G) + 1, . . . , D(G′).
The informal interpretation of WST, from which its name Earth Mover is derived, is captured by the formula.

Experimental setting and computational results
In our problem the network space is given by the basic network G and the subgraphs obtained by the removal of one or more edges. The elements of this space are represented as probability distributions of node-2-node distances (Eq. 1) and aggregated into a distributional representations of the basic network and of the subgraphs (Eq. 2). In this space it takes place the computation of the Wasserstein distance. The result of this computation is mapped back into the network space as measures of network dissimilarity and labels of criticality of individual components. In 4.1 the experimental setting is described.

Experimental setting
The networks considered are: • Neptun is the WDN of the Romanian city of Timisoara, with an associated graph of 333 nodes and 339 edges. • Abbiategrasso refers to a pressure management zone in Milan (namely, Abbiategrasso) with an associated graph consisting of 1213 nodes and 1391 edges. The description of these networks is also given in "Appendix C". As far as efficiency and algebraic connectivity measures are concerned, we have used the standard measures reported in the "Appendix B".
The first step in the analysis is clustering in order to identify the specific edges whose removal induces a disconnection of the network. The number of clusters K is set according to context information about the districtualization adopted by the water utility. Failures affecting also only one pipe may imply a reduction in the efficiency of the network and an increase in vulnerability. The software used is described in "Appendix C.2".
The two real-world WDNs analyzed are very sparse (with density q lower or equal to 0.006). The two real-world WDNs (given in "Appendix C.1") are effectively planar and "almost" regular. This fact can be due to the fact that their structure is strongly constrained by spatial characteristics making a classification based on nodal degree distribution less meaningful.
In particular their degree distribution does not follow a power law and measures, given in Table 1, for Neptun and Abbiategrasso, really set them apart from other kinds of networks like transportation, communications and social. The usefulness of measures based on centrality or spectral measures in assessing the increase of vulnerability or the loss of robustness is quite limited, according to our results, in WDN of realistic size.

Computational results
The computational results, which are quite unique in the literature given the size of the networks analysed, demonstrate that probabilistic distance measures show better capacity to discriminate between different networks not only globally but also edge-wise. The results are organized around figures and tables. Figures 4 and 5 display (respectively for Neptun and Abbiategrasso) the output of clustering, that is critical edges whose removal generates a disconnection; these edges are highlighted in red. Table 1 displays centrality measures. Table 2 efficiency, vulnerability, and algebraic connectivity. Table 3 the values of distributional distances compared with the loss of efficiency. The most relevant result is that using discrete probability distributions to represent the graph associated to the WSN, as well as all its elements, provides a novel and more effective analytical framework to assess similarity between two networks (Table 1), as well as the same network before and after some modifications, such as expansion, disruptions, etc. Especially in the case of disruptive events, the novel WST based analytical framework enable the definition and the computation of new and more discriminant criticality indices (comparison between Tables 2 and 3), ideally applicable to any critical networked infrastructure (transportation networks, energy grids, communications networks, etc.). Another relevant benefit provided by the novel WST based analysis is that it can deal with distributions with different supports as well as sparsity, and it is not affected by different binning schemes. This is a relevant advance with respect to previous approaches based on probability distributions, for instance those based on the Kullback-Leibler or the Jenses-Shannon divergences, which require to adopt some drawbacks for dealing with the two mentioned issues. Even more important, the proposed WST based analysis allows for an intuitive visualization of the criticality of each edge, matching the typical requirement for a "priority ranking" of critical elements, driving the preparation of a budget-constrained rehabilitation plan. For the WST distance it has been computed also the heat map of edge wise criticality for Neptun and Abbiategrasso. In Fig. 6, the colour associated to each edge i, j is an index of criticality given by the Wasserstein distance between G and G\ i, j .

Conclusions and perspectives
The main result of this paper is that probabilistic measures based on the probability distribution of node-node distances, yield a distance, between the original network and that resulting from the removal of some nodes/edges, which can provide a set of indicators of the increase in vulnerability. The computational results confirm that the value of the distances Jensen-Shannon and Wasserstein confirm the results of clustering.
Albeit the measures of vulnerability analyzed in this paper are based on topological arguments they are anyway important to narrow the set of critical components before moving to the hydraulic verification.  Significantly the differentiability of WST allows to include in the differential robustness analysis edges weighted by parameters derived from the flow dynamics in the network correlating connectivity analysis to hydraulic performance indicators.
This analysis framework supports decision making at design stage, to simulate alternative network layouts of different robustness, and also at operational stage where the decision to be taken can be, which nodes/edges are to temporarily be removed for maintenance and rehabilitation. Indeed, critical tasks of WDN management can be supported by just using topological and geometric information. The analysis framework also helps for the efficient and automatic definition of district metered areas and to facilitate the localization of water losses through the definition of an optimal network partitioning.
The modelling and algorithmic framework platform developed can be straightforwardly translated to many networked infrastructures among which power grids, transit networks but also global supply chains network whose vulnerability has been exposed in the recent COVID crisis.

Appendix A
Graph theory is the mathematical basis to provide a unifying language for the study of networks: with this in mind it is useful to give some basic definitions which will be used in the sequel. For a wide-ranging analysis of the role of graph theory in the analysis of networks the reader is advised to look at (Newman 2010).

A.1 Basic definitions
Let denote a graph with G = (V , E) , where V is the set of nodes and E is the set of edges. Each edge of G is represented by a pair of nodes i, j with i = j , and i, j ∈ V and with n = |V | and m = |E| . If i, j ∈ E, i and j are called adjacent nodes. A graph G is undirected if i, j and j, i represent the same edge. A graph G is simple if no self-loops are admitted (edges starting from a node and ending on the same node) and only one edge can exist between each pair of nodes i, j , with i = j . The adjacency relationship between the nodes of G can be represented through a non-negative n × n matrix A (i.e., the adjacency matrix of G ). The entry a i,j of the adjacency matrix A is 1 if i and j are adjacent Fig. 6 Heatmap of edge wise criticality for Neptun (left) and Abbiategrasso (right). The colour associated to each edge (i, j) is an index of criticality given by the WST distance between G and G\(i, j) nodes (i.e., i, j ∈ E ), and 0 otherwise. Furthermore, a ij = a ji if G is undirected and a ii (entries on the diagonal) are 0 if G is simple.
• The degree of the node i , k i is the number of edges having i as one of the two nodes on the edge: k i = n j=1 a i,j . Anyone of the edges having i as one of its nodes is called incident on i.
• When G is directed, meaning that the order of the two nodes of an edge is relevant for its definition, the k i can be split into out-degree (number of edges having i as first node) and in-degree (number of edges having i as second node). • A path in a graph is a sequence of nodes connected by edges the length of the path is the number of edges. A connected component is a maximal subgraph when all nodes can be reached from every other. • The shortest path between i and j is the path with the smallest length. This length is called the distance between i and j d i,j . The largest distance among each possible pair of nodes in G is named diameter D(G). • The characteristic path length is the average distance for every possible pair of nodes i, j .
A useful representation is to arrange the distances in the distance matrix D = d i,j i, j = 1, . . . , n . The maximum entry of row i max j=1,...,n d i,j is also known as the eccentricity of node i . The maximum eccentricity among the nodes is equal to D(G).
A subgraph G′ = (V ′, E′) of G is a graph such that V ′ ⊆ V and E′ ⊆ E ; a connected component of G is maximal if is the largest possible subgraph for which you could not find another node in the graph that could be added to the graph with all the nodes be still connected.
The core concept is centrality which addresses the question "which are the most important nodes in a network?". There are many centrality measures from the simplest like node degree, which can anyway be illuminating, to eigenvector-based measures like Page Rank.

A.2 Basic measures
The density of the network is the fraction of edges which are present in the network: The number of edges m = 1 2 n i=1 k i .
If c is the mean node degree, c = 1 n n i=1 k i we get c = 2m n and q = c n−1 .
The density is in the range (0, 1).
A cut-set, specifically a node cut-set, is a set of nodes whose removal disconnects i and j . A minimum cut-set is the smallest cut-set. Analogously for edge cut-set.
The centrality measures address the issue of the relative importance of nodes/edges. The most widely used measures are: • Closeness centrality C = n d ij is based on the mean distance from i to j averaged on all nodes. • Betweenness centrality: let be η i st = 1 if node i lies on the shortest path from s to t and 0 otherwise. Then, betweenness centrality is given by b i = 1 n 2 n s,t=1 η i st . It measures the extent to which a node lies on the paths between other nodes. We can similarly define an edge betweenness that counts the number of shorter paths that run along the edge.
• Link-per-node ratio (e), as the number of edges of a graph with respect to the number of its nodes. e = m n . • Central point dominance c ′ b , based on betweenness centrality is a measure for characterizing the organization of a network according to its path-related connectivity; where b i is the betweenness centrality of the node i and b max is the maximum value of betweenness centrality over all the n nodes of the network. • The clustering coefficient (CC) is the number of triangles with respect to the overall number of possible connected triples, where a triple consists of three nodes connected at least by two edges while a triangle consists of three nodes connected exactly by three edges: There are other definitions of CC for which the reader is addressed to Newman (2010).
To compute the centrality indices in this paper, the open-source software Cytoscape (http:// www. cytos cape. org/) has been adopted (Morris et al. 2011). This point will be further discussed in "Wasserstein distance" section.

A.3 Spectral measures
Spectral graph theory studies the eigenvalues of matrices that embody the graph structure. One of the main objectives in spectral graph theory is to deduce structural characteristics of a graph from such eigenvalue spectra.
In case of undirected graphs, the adjacency matrix A(G) is symmetric and all its eigenvalues are real. The eigenvalues µ 1 (G) ≤ µ 2 (G) ≤…µ n (G) of A(G) are called the spectrum of G. The largest eigenvalue of the adjacency matrix µ n (G) is called spectral radius of G and is denoted by ρ(G) . An important property is given by the following inequality.
where �(G) = max k i : i = 1, . . . , n that relates the spectral radius with �(G), the maximum degree of the nodes.
The spectrum of A(G) allows to define the Eigenvector centrality of the node i , is x i = ρ(G) −1 n j=1,j� =i a ij x j . Katz centrality and Page Rank algorithm are just parametrized version of eigenvector centrality (Newman 2010).
The difference s(G) = ρ(G) − µ n−1 (G) between the spectral radius of G and the second largest eigenvalue of the adjacency matrix A(G) is called the spectral gap of G. A small value of s(G) is usually observed through low connectivity, and the presence of bottlenecks and bridges whose removal cuts the graph into disconnected parts.
The Laplacian matrix of G is an n × n matrix L(G) = D(G) − A(G) , where D(G) = diag(k i ) . The matrix L(G) is positive semi-definite in case of simple graph. The eigenvalues of L(G) are called the Laplacian eigenvalues of G. The Laplacian eigenvalues 1 (G) = 0 ≤ 2 (G), ... ≤ n (G) are all real and nonnegative. The smallest eigenvalue is always equal to 0 with multiplicity equals to the number of connected components of G . The second smallest eigenvalue is called the algebraic connectivity of G which is one of the most widely used measures of connectivity. Larger values 2 (G ) represent higher robustness against efforts to disconnect the graph, so the larger it is, the more difficult it is to cut a graph into independent components. An important inequality for the algebraic connectivity is given by that relates it with the minimum degree of the nodes δ(G) = min i=1,...,n k i . In case of a connected graph, also the following inequality can be proved that relates the algebraic connectivity with the diameter of the graph. Another spectral distance is based on the analysis of the eigenvectors of the Laplacian.
where G\{v} denotes the network G without the node i . The loss of efficiency of G is defined as The larger the value V E (G) , the larger is the vulnerability. Analogous formulas can be written removing the edges.

B.2 Spectral vulnerability measures
There is no specific formula, contrary to those reported in the previous subsections, linking spectral analysis to a measure of vulnerability related to the removal of a node. However, both algebraic connectivity 2 , the second smallest eigenvalue of the Laplacian L(G) and spectral gap s(G) , the difference between the spectral radius of G and the second eigenvalue of the adjacency matrix A(G) , are indicators of difficulty to split the graph. The larger the algebraic connectivity, the more difficult it is to disconnect the graph. It is also related to the min-cut problem in spectral clustering. A large value of the spectral gap, together with a uniform degree distribution, results in higher structural sturdiness and robustness against node and link failures. On the contrary, low values of spectral gap indicate a lack of good expansion properties usually represented by bridges, and network bottlenecks. The larger the spectral gap the more robust is the network.

C.1 Data resources
• Neptun is the WDN of the Romanian city of Timisoara, with an associated graph of 333 nodes and 339 edges. This network has been a pilot in the European project Icewater. • Abbiategrasso refers to a pressure management zone in Milan (namely, Abbiategrasso) with an associated graph consisting of 1213 nodes and 1391 edges. This network has also been a pilot in the European project Icewater.
In analyzing WDNs one must consider that most of the end-users are supplied by single connections. A preliminary preprocessing has been performed in order to identify the simple core of the network by removing all nodes with degrees smaller than two (also called leaves of the network) along with the parent edge connecting them to their associated root-node.)

C.2 Software resources
The tools for the analysis of connectivity and efficiency have been coded in Pythton.
For the network visualization it has been used the Python package WINTR, based on EPANET 2 in WDN's. which enables network editing and hydraulic simulation.
The software used for the computation of JS is the function JensenShannonDivergence from the Java library Mallet. The software used for the computation of the Wasserstein distance is the function EarthMoversDistance from the Java library Apache Commons Math. For the clustering it has been used an effective computational scheme which uses a data representation in the lower dimensional space spanned by the most relevant eigenvectors of the normalized Laplacian matrix. The basic steps are: