A Geometric Approach for Real-time Monitoring of Dynamic Large Scale Graphs: AS-level graphs illustrated

by   Loqman Salamatian, et al.
IAE Savoie Mont Blanc

The monitoring of large dynamic networks is a major chal- lenge for a wide range of application. The complexity stems from properties of the underlying graphs, in which slight local changes can lead to sizable variations of global prop- erties, e.g., under certain conditions, a single link cut that may be overlooked during monitoring can result in splitting the graph into two disconnected components. Moreover, it is often difficult to determine whether a change will propagate globally or remain local. Traditional graph theory measure such as the centrality or the assortativity of the graph are not satisfying to characterize global properties of the graph. In this paper, we tackle the problem of real-time monitoring of dynamic large scale graphs by developing a geometric approach that leverages notions of geometric curvature and recent development in graph embeddings using Ollivier-Ricci curvature [47]. We illustrate the use of our method by consid- ering the practical case of monitoring dynamic variations of global Internet using topology changes information provided by combining several BGP feeds. In particular, we use our method to detect major events and changes via the geometry of the embedding of the graph.



There are no comments yet.


page 1

page 2

page 3

page 4


Distributed dynamic modeling and monitoring for large-scale industrial processes under closed-loop control

For large-scale industrial processes under closed-loop control, process ...

Anomaly and Change Detection in Graph Streams through Constant-Curvature Manifold Embeddings

Mapping complex input data into suitable lower dimensional manifolds is ...

Nonparametric Link Prediction in Large Scale Dynamic Networks

We propose a nonparametric approach to link prediction in large-scale dy...

Staged Animation Strategies for Online Dynamic Networks

Dynamic networks – networks that change over time – can be categorized i...

Monitoring dynamic networks: a simulation-based strategy for comparing monitoring methods and a comparative study

Recently there has been a lot of interest in monitoring and identifying ...

Node Alertness-Detecting changes in rapidly evolving graphs

In this article we describe a new approach for detecting changes in rapi...

Growing Graphs with Hyperedge Replacement Graph Grammars

Discovering the underlying structures present in large real world graphs...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1. Introduction

Monitoring large scale dynamic graphs such as the Internet, neuronal connections or social networks is a challenging task, mainly due to their rapidly changing properties. In particular, slight local changes in a graph, such as a link removal, can result in major variations in global properties. For instance, a few disconnected components may yield great changes in the average lengths between vertices. It is often challenging to distinguish changes in a graph that have only a local impact, from those that propagate globally.

The relevance and difficulty of the efficient monitoring of large scale networks is illustrated in the monitoring of global Internet through information provided by the Border Gateway Protocol (BGP). BGP is the standard inter-domain routing protocol used over Internet to primarily exchange reachability information among Border routers in Autonomous Systems (AS). Dynamic variations of the network resulting from local events, e.g., link failures, router reboot, network operator policies changes or malicious events may lead to a continuous stream of prefix updates or removal. The stream that reflects localized changes in the topology of Internet can be observed in the routing tables updates. While the impact of these events remains generally localized on the affected ASes (e.g. an outage in a stub AS might disconnect it from Internet and perturb only direct communications) sometimes the reach is larger and the impact deeper. In August 2017, a BGP leak 111illegitimate advertisement for prefixes generally due to misconfiguations from a single stub AS, (one of Google’s in this case), led to a partial black-out in Japan, resulting into a major and large scope disruption (jap, [n. d.]). These larger scope events might affect multiple ASes, largely beyond the ones directly involved in the source BGP event, in particular with how widespread BGP security risks are (Al-Musawi et al., 2017a).

Several methods for detecting and mitigating the impact of BGP-level network events such as BGP hijacks have been introduced. Proactive mechanisms, like RPKI (Wählisch et al., 2012), rely on cryptography mechanisms to authenticate the source of BGP announcements and the sequence of updates (Butler et al., 2010). Other approaches are reactive and involve first detecting that a suspicious (or abnormal) event has happened and then taking circumvention actions (Kruegel et al., 2003)(Karlin et al., 2008). The detection of BGP related events is often implemented by third-parties that monitor several BGP vantage points and have therefore a broader view of the BGP activity (claffy, 2001)(Orsini et al., 2016)

. These approaches leverage the characteristics of BGP hijacks as seen in practice, and combine the views from several vantage points (up to 100 collection points for some of them) to implement heuristics to detect BGP events that are then made available to AS owners

(dyn, [n. d.]).

While these approaches might be successful identifying the source of the abnormal event, they are unfortunately incapable of tracking the magnitude of an event within the BGP-derived AS graph and of characterizing its importance and global impact. In this paper, our aim is to develop a monitoring system that detects quasi-instantaneously whether important changes have happened in the network and determine how they may have affected the network.

Understanding and characterizing the scope of BGP events is important both for security purposes and for troubleshooting (Li and Brooks, 2011)(Labovitz et al., 1997). We consider two application settings to motivate our method. In the first setting we consider an AS, possibly multi-homed, that would like to monitor large scale changes in the Internet connectivity. In the second setting, we consider a third-party actor, possibly a governmental authority, that is interested in monitoring the global stability of Internet and in this capacity would like to detect major changes on Internet connectivity. We will show later that the approach developed in the paper is relevant to both scenarios.

Overview of the technique: At each snapshot, the network can be defined as a graph with a set of vertices linked through edges. Graph embedding techniques (Cai et al., 2018) can be applied to embed this graph into a space endowed with a metric capturing the interrelation between local and global connectivity. The motivation for embedding graphs in manifolds comes from results in continuous manifold theory, namely the Gauss-Bonnet theorem (Lafontaine, 2015) which relates a local property, the curvature, to global topological properties of the continuous manifold, e.g. number of holes. The Gauss-Bonnet theorem ensures that an increase in curvature in one part of the manifold, has to be compensated by a decrease in another part, unless an important change in the topology has happened, i.e. the manifold acts as a play-dough that when squeezed somewhere will inflate somewhere else.

In this context we use a discretization of the Ricci curvature (Ollivier, 2012). The Ricci curvature has been widely used in continuous manifold theory to obtain volume and propagation inequalities such as the Brunn-Minkowski theorem (Bobkov and Madiman, 2013), and relates the local properties of a manifold with its overall shape. The Ollivier-Ricci curvature defines an extension of the Ricci curvature in discrete settings and is defined through the notions of optimal transport and transport plan. The optimal transport formalizes the problem of finding (if it exists) a function that transforms a distribution of weight over a set of vertices into another desired distribution, while minimizing the product of the amount of transported weight with the distance to which it was transported (Lévy and Schwindt, 2018). We provide further details in Sections 2.2 and 2.3. More in depth descriptions can also be found in (Lafontaine, 2015; Ollivier, 2012).

By projecting every element of a graph with vertices with respect to the Ollivier-Ricci curvature, we obtain a manifold located in where the coordinates of vertex are the Ollivier-Ricci curvatures to all other vertices . A dynamic network can, therefore, be seen as a drifting process where vertex positions are fluctuating over time. This embedding of the graph in is called the Ollivier-Ricci graph Embedding (O.R.E.). In practice, however, embedding the graph in a dimensional space is a daunting task. To reduce the dimensionality, we define landmarks that act as anchors in the space, and we project the graph into a -dimension space. We will present a way of selecting the landmark in a way that the loss of information is minimal.

Assuming we know this embedding, our goal is to detect major changes in the geometry and topology of the BGP graph. When a local change propagates globally within a manifold, it impacts the overall curvature, i.e., it changes a large portion of the shortest paths on the manifold (named geodesics).

Main contributions : The contributions of this paper are as follows:

  1. We propose in section 2.3, a geometric view-point to the analysis of AS-level graphs as measured through multiple BGP collectors. We introduce the background material briefly describing the Riemannian geometry and present metrics applicable to AS-level graphs. We use the Ricci curvature and its extension, the Ollivier-Ricci Curvature to characterize the geometry of discrete graphs, in particular the AS-level graphs. We show that the Ollivier-Ricci (Ollivier, 2009) curvature can be used as a metric to quantify the impact of a local change on the geometry of the network. These changes can be interpreted via the so-called optimal transport, thus giving a detailed picture of how a change will affect the global network traffic. This provides a new way of quantifying the severity of BGP events.

  2. We tackle, in section  3.3, the problem of finding landmarks in highly dynamic settings. We provide and validate a statistical approach to discover optimal landmark points.

  3. We develop, in section 3, a monitoring system that uses the Ollivier-Ricci curvature to evaluate if incoming BGP updates have resulted in geometric changes over the global Internet.

  4. We evaluate, in section  4, this geometry change detection framework over real BGP flows and show that the monitoring system detects all known major BGP-related events of the past years. We also validate that local events can be distinguished from global events. During the evaluation, we compare the Ollivier-Ricci embedding metric with the traditional metrics (the shortest path distance or the spectral distance (Luxburg et al., 2004)). We show that the Ollivier-Ricci distance exhibits better stability in events associated to local changes while showing great sensitivity for events of global changes. We also show that the Ollivier-Ricci metric can uncover relationships between ASes that are invisible to other metrics.

2. Background

2.1. BGP collection and AS level graphs

Interdomain routing is a collaborative effort among ASes, which interconnects and exchanges routing information using the BGP protocol. Two main types of announcements messages are exchanged between routers running the BGP protocol: updates that contain each an AS path composed of the list of ASes to cross in order to reach a given prefix through the announcing neighboring AS, and withdrawals announcing that a IP address prefix is not anymore reachable through the announcing AS (Caesar and Rexford, 2005).

The decision to announce a path depends on route export policies that derive from contractual relationship between ASes as well as traffic engineering implemented by ASes to manage the traffic exchanges. Generally, both ASes relationships and traffic engineering policies are maintained confidential.

Any BGP router will have to deal with a continuous flow of update or withdrawal messages that reflect every slight change in the network topology. When receiving an update for a prefix, a router check in his routing table if it already knows an AS path to it. If this is the case, it will choose between the known active AS path and the new incoming one (Bonaventure, 2002). This choice can be done using the AS path length, i.e., a shorter AS path length being better, or through more complex mechanism involving local preferences, multi-exit-discriminator path attributes or AS policies. An updated path might be re-announced to neighbors or rather filtered depending on the AS policies. Sometimes AS path are inflated by repeating one of several ASes to make them less likely to be chosen as forwarding path as they become long (Zhang et al., 2004).

Another characteristic of BGP is its lack of authentication for routes announcement. It is therefore possible for an AS able to advertise illegitimately paths for prefixes it does not own. These illegitimate advertisements, called BGP Prefix hijacks, can impact many ASes, or even the whole Internet, affecting service availability and security of communications as packets might follow paths they are not supposed to take, e.g., because surveillance is implemented along this path. BGP hijack might also be caused by router misconfiguration or malicious attacks (Ballani et al., 2007).

When an AS path is announced, with two consecutive AS1 and AS2, we can infer that there is a direct link between AS1 and AS2. Therefore, by looking at the flow of announced paths coming from neighboring ASes one can infer an AS level graph connecting ASes that have appeared at least on one of the announced AS paths (Gao and Rexford, 2001). We leverage this approach to create each minute a snapshot of the AS level graph. For this purpose we combine several BGP feeds coming from several BGP collectors. As announcements are made at prefix level, we can assign to each link in the AS level graph a count attribute that measure the number of IP addresses that can be reached through paths going through the link. We moreover assign to each vertex in the graph some attributes including the number of prefixes issued at the AS represented by the node, and ctime the last time there was an update or an announcement relative to that AS.

Indeed, the AS graphs generated by this approach are known to be incomplete (Chen et al., 2002). In particular, BGP path filtering policies do not expose less preferred paths that would be chosen if the preferred announced path were not available. Moreover, traffic engineering and load balancing can further blur the real AS graph (Knight et al., 2011). This incompleteness adds a dimension to the challenge addressed in this paper. The question we will address in this paper is ”can we detect major changes in Internet structure by just looking at incomplete AS graphs that are gathered through a limited number of collection points”.

2.2. Riemanian Manifolds

Intuitively, a manifold is a surface in a multi-dimensional space that can be locally mapped to an euclidean space, i.e., if we zoom in enough on the surface it resembles an euclidean plane. A Riemannian manifold is one that have attached to it a bilinear distance metric, which then allows for a notion of angles and distances on the manifold. One can look at a manifold from two perspectives: the extrinsic view considers the manifold as an object embedded in a (larger) euclidean space, while the intrinsic view considers properties determined solely by distance within the surface. In other terms, the extrinsic view is the one of an observer that could look at the manifold globally, and the intrinsic view is that of an ant (or an army of ants) moving along the manifold. Relation between these two perspectives has been a major question in geometry, e.g. the ”remarkable theorem” of Gauss, Theorema Egregium, states that the curvature of a surface does not change if one bends the surface without stretching it. Thus the curvature is an intrinsic invariant of a surface (Berger and Gostiaux, 1987). The intrinsic and extrinsic viewpoints are fundamental in our work as well, since we aim at characterizing which local changes will result in major global changes. Through the Theorema Egregium the notion of curvature of a surface becomes a main tool for doing this.

Curvature can be intimately related to the local behavior of geodesics, i.e., shortest paths along the manifold. On a space with (strictly) positive curvature, two distinct geodesics starting at two points close to each other and pointing to the same destination will ultimately converge to the same point (Ollivier, 2012). Inversely with negative curvature, two geodesics will drift apart getting further of each other (see Figure 1).

Figure 1. In flat geometry, parallel lines are equidistant while they tend to get further (resp. closer) in hyperbolic (resp. spherical) space.

Several metrics exist for describing the curvature such as the sectional curvature (Berger and Gostiaux, 1987), Gromov -hyperbolicity (Chen et al., 2012)

, Bakry Émery Tensor

(Bakry and Émery, 1985), Forman-Ricci curvature (Forman, 2003) (Weber et al., 2016) etc. (Berger, 2003). In this paper, we use an extension of the Ricci curvature to discrete objects, and particularly graphs, which we describe in the following section.

2.3. Discrete curvature and the Ollivier-Ricci curvature

Let’s define the discrete extension of Ricci curvature to a graph . This extension uses the concept of transportation cost (Piccoli and Rossi, 2012) over a graph. Let’s assume a graph where a distance metric is defined between each pair of vertices . This distance can be the minimum hop distance, a minimal weight distance or any other arbitrary distance matrix. Consider distributions and over the vertex set , i.e. for and . The optimal transport is defined as follows.

Definition 2.1 (Optimal Transport).

where the minimization is over all transport plans , i.e. such that


The value is referred to as the transportation distance.

Intuitively, the transportation distance evaluates the effort needed to transport a mass distributed following a distribution over the different vertices of the graph to another mass distribution . The optimal transport

always exists and can be computed exactly by solving the linear program (

1), or approximated efficiently (see for example (Cuturi, 2013)).

We can leverage the optimal transport to extend the Ricci curvature to graphs.

Definition 2.2 (Ollivier’s Ricci curvature (Ollivier, 2009)).

The Ollivier Ricci curvature between two vertices and is defined as:


where for is a family of distribution.

The distributions can be general, and we will specify which distribution is suitable for our application later on. It can be shown that the curvature is bounded in general as (Ollivier, 2012). However for distribution such that half of the mass is positioned on the central node and the remaining mass is distributed over the neighbors, one can prove that (Jost and Liu, 2014).

Examples: To illustrate intuitively the Ollivier-Ricci curvature let us look at three extreme cases. For these examples, we assume put weights uniformly among neighbors of , i.e. if is the number of neighbors of vertex , then if is a neighbor of , and 0 otherwise.

  • Let be a clique with vertices. In this case all neighbors of are neighbors of . So the optimal transport plan has not to transport anything from neighbors of to neighbors of and only needs to transport a mass from to . This means that the Ollivier-Ricci curvature equal to . As grows, the curvature tends to its maximal value .

  • Let now be an alignment of nodes, i.e., a line graph, and let and be any two distinct nodes on the graph. It is easy to verify that the optimal transport corresponds to transporting the mass in the left (resp. right) neighbor of to the left (resp. right) neighbor of with a cost of , therefore .

  • Let be two nodes star networks connected by a single link. In this case, the optimal transport plan between the first star center and the second center send all the masses in neighbors of to neighbors of through the link connecting the two stars with a distance cost of 3 and transport the mass in to with a distance cost of 1. This means that the optimal cost will be resulting into an Ollivier-Ricci curvature . It is noteworthy that if the distribution is such that half of the masses are positioned on and respectively with the remaining masses are distributed evenly on their neighbors, the optimal transport cost becomes and the Ollivier-Ricci curvature is

Real graphs like the AS level graphs are not as simple as the above examples. Nonetheless, when the Ollivier-Ricci curvature obtained between two vertices is close to 1 this can be interpreted as the two nodes are in a clique-like structure, i.e., there is a rich network of link connecting the two nodes such that almost all nodes are neighbors of each other. On the other hand, when we have negative curvature close to , all the paths connecting the two vertices are fully relying on a few nodes. This interpretation will be helpful to understand the curvature variation, i.e., increase of curvature means that the network is becoming more connected and decrease of curvature is a sign of the emergence of a bottleneck.

3. Ollivier-Ricci Embedding Monitoring System

In this section we will describe how to implement the Ollivier-Ricci Embedding of the AS level graph. Recall that the embedding consists of projecting the AS level graph into a -dimensional Euclidean subspace where each coordinate is the Ricci curvature measured from a node of the graph (an AS) to one of the landmarks. We will define the curvature matrix as the matrix where each line contains the dimensional projection of an AS that have seen a BGP update during the graph collection snapshot.

3.1. Pre-processing steps

The previous section briefly introduced some background of the theoretical tools used by the monitoring system we describe here. We introduced the Ollivier-Ricci curvature as a metric relating local property of a vertex to global properties of the graph.

We assume access to a series of (periodic) snapshots of the AS-level BGP graphs. Each vertex representing an AS in this graph has a set of attributes including a timestamp ctime identifying the last time the AS has seen a prefix update or withdrawal announcement, as well as the number of prefixes it announces. Each link also contains a count attribute measuring the number of prefixes that are announced to be reachable through this link.

(a) Global effect event
(b) Local effect event
(c) Natural drift


Figure 2. Representative sample of the matrix for 3 different types of events. In Fig. (a)a we show the Google leak Incident that impacted Japan the 25th of August 2017. In Fig.(b)b we show a local BGP event, an outage, that impacts only slightly a small part of the matrix (Brazil leak, see table 2). The event pictured in Fig. (c)c corresponds to a period (5 Jan 2018, see in table 2) where there is no either local, outage or hijack, or global event in the network. However, there will still be a natural drift coming from BGP updates and withdrawals.

Choice of mass distribution:

In order to derive the curvature we have to set a source and destination mass distribution for the transport ( and in (1)). In his seminal paper, Ollivier (Ollivier, 2009) used a simplistic distribution that distributes evenly the mass to transport over the neighbors of the vertex. As explained in Section 2.3, we will put a mass in the center node in order to ensure that the Ollivier-Ricci curvature between each two nodes and remains in .

However, this distribution assumes the same importance for all neighbors which might be unrealistic. A good choice for this distribution is to set it according to the amount of traffic flowing from one AS to another. This information is often unavailable for large networks. Alternatively, we need to use a measure that can be readily obtained from BGP feeds.

As we store for each AS to AS link in the network a count attribute that contains the number of IP addresses reachable through this link, we set the source and destination mass distributions following the distribution of count to the neighbors. This distribution, that is recalculated after each update or withdrawal, gives to each neighbor a weight proportional to the number of IP addresses that is accessible through them. Next, we use the Ollivier-Ricci curvature derived from this measure.

Construction of curvature matrix:

In its steady state, the AS graph contains up to 80 K ASes and 200 K links that dynamically change because of permanent updates and withdrawals. It is therefore unfeasible to monitor the variation of curvature between all ASes and landmarks. We therefore resort to use the ctime attribute to detect ASes that have undergone a change and limit the computation to those ASes.

Let’s assume that during the snapshot , ASes have seen a prefix update or withdrawal. We derive for the snapshot a matrix , where is the number of landmarks.

In order to detect changes we generate a similar matrix of curvatures with the same size, but for the last snapshot . We then use as the basic statistics for changes detection. We show three examples of the matrix corresponding to three different BGP events in Figure (f)f.

3.2. Monitoring

A dynamic graph will undergo continuous changes that will alter the graph locally and potentially lead to small modifications at the global level. This results in permanent fluctuations of the curvature measure. Nonetheless, the Ollivier-Ricci variations can become large when an event changes the geometry and impacts the geodesics. Every non-zero value in the , corresponds to either a stretch or a compression of the underlying geometry and indicates a change in the network conditions, possibly induced by an anomaly.

The sign of the spikes characterizes the effect the event had on the connectivity. A positive spike means that the curvature has increased, and that the neighborhoods of the changed AS are now more connected to the rest of the network. Inversely, a negative spike implies that the neighborhoods of the changed AS are more dependent on a reduced number of connections with the landmarks. Generally, a BGP hijacks or leaks generate a positive shift in the Ricci curvature as they involve for the impacted vertex the emergence of new and shorter paths (Sermpezis et al., 2018). Inversely, link (route) failures generally decrease the curvature and the diversity (see for example Fig. (b)b).

Interpretation of Curvature Matrix

As can be seen in Fig. (f)f, the matrix contains vertical and horizontal alignments of decreases or increases in curvature. A column with large variations of curvature means that the geometry toward the landmark attached to this column has experienced a major change. If several landmarks experience important changes, it indicates a global change. Similarly, a horizontal alignment of large variations means that the AS has seen major change toward all landmarks. A global change will induce several horizontal and vertical lines with important changes. We see also in Fig. (c)c a representative case showing the perpetual drift. Generally, the drift will show itself with changes in the curvature that are limited to the impacted ASes (represented by horizontal lines in the heatmap of the curvature matrix). Local events and global events will exhibit vertical lines in the heatmap, meaning that at least one of the landmarks has seen curvature change to numerous ASes. The difference between a local and a global event is in the strength of the change and the number of large curvature changes. While the previous analysis is intuitive, we need to define an automatized detection scheme.

We show in Fig 3 the sum of positive and negative values in the matrix. Following the discussion in section 3.2, we observe a behavior compatible with the Gauss-Bonnet theorem, i.e. zero sum or very small one in period without major BGP events, and a change in period with large event.

Figure 3. Time series of the sum of the positive (resp. negative) elements in along with the global sum.

From the previous discussion one can consider that monitoring the sum of curvature is enough to detect large events. Nonetheless, when the Gauss-Bonnet theorem is valid, the sum of curvature will not change whatever large is the fluctuation. However, the Gauss-Bonnet theorem is not formally valid for Ollivier-Ricci curvature. This means, just looking at the sum of curvature might generate false alarms and misdetections. We have therefore to define another method to translate the matrix into a time-series over which we will detect large events.

Analysis of Curvature MatrixAs said before columns or lines structure in the matrix matrix

are representative of the changes in the topology. Typically, the matrix contains several very small values, relative to small curvature changes, and a limited number of lines or columns with larger curvature changes. In order to analyze the structure of this matrix we will use its singular values

222Singular values are extension of eigenvalues to non-square matrices (sin, 2014), i.e.

, the roots of eigenvalues of the

matrix . The elements on the diagonal of are the column-wise sum of squares of matrix and its trace is the Frobenius norm of , :

The Frobenius norm represent therefore the variance of the values in

, i.e., the strength of the curvature changes.

One can check that if the matrix contains a non-zero column, i.e., its rank is , the matrix will have strictly positive eigenvalues, or equivalently will have non-zero singular values. However, in practice, the curvature matrix is filled with small values, rather than 0, and only some of the lines or columns will have large values. In such case, have some large singular values representing the large values columns, and other eigenvalues close to 0. In such contexts, the stable rank defined as

where is the largest eigenvalue of .

is frequently used as a robust estimator of the rank of a matrix, in particular for generating low rank approximation of noisy matrices

(Cohen et al., 2016). As the rank is directly related to number of large columns in , the stable rank estimates this number, i.e., a stable rank close to 1 means that only a single column of has important changes, while a large stable rank indicates large number of columns, i.e., landmarks, having seen important changes. It is also noteworthy that have the below property :

In other terms, the largest singular value gives the largest norm stretching that the matrix can generate.

The above observations lead into a simple anomaly detection scheme. For each incoming matrix

, we evaluate the normalized Frobenius norm, i.e., , where is the number of ASes having seen changes in the snapshot, and the stable rank of , i.e., . The set defines a time trajectory that will be used to detect large events. We use in place of , as it is bounded, i.e., . A large event will have a large value of and a small , while a local event might have large value of but will have a value of close to 1.

Figure 4. Components of the monitoring system

3.3. Landmarks selection

Our approach relies on constructing the Ollivier-Ricci embedding in . However, monitoring the variations of all distances is unfeasible over a real AS-level graph with over 60k vertices and 3.6 billion distances (cid, [n. d.]). We reduce the dimensionality by choosing a set of landmarks and only monitor the variation of curvature towards them, i.e., we represent the position of vertex

by a vector

where is the curvature from to landmark . The essential property we are seeking to preserve despite the dimensionality reduction, is that the variation of the time series between time and , denoted hereafter by the drift, should be small when there are no major changes in the graph underlying geometry and substantial when such a large change happens. The landmarks should satisfy this property and we will validate this claim in Sec. 4.

Choosing landmarks in a highly dynamic graph is a complex task (Leskovec and Faloutsos, 2006).

We considered multiple landmark selection process: (a) random landmarks selection, (b) highest degree, (c) highest centrality, (d) highest number of triangles, (e) Tier 1 and Tier 2 AS, a random walk starting at (e) a random vertex or (g) at the collector node, and finally (f) a lazy random walk that we introduced and which was developed for this purpose.

In order to evaluate the performance of different landmark selection methods, we define two metrics that are calculated over a landmark set . A first metric, , evaluates the diversity arising from , by assessing the proportion of distinct neighbors of vertices in among all possible neighbors, and the average distance in hops between landmarks, i.e.:


High diversity is a necessary condition for a landmark to have high curvature. So we use the metric as a low complexity proxy to detecting potentially high curvatures. A good set of landmark should therefore at the same time achieve a good diversity, large value of , and be far apart, large value of . We consider 30 of our collected AS graphs and run over them the different landmark selection algorithm and see the outcomes. For random walk based approaches (lazy random walk, fully random and random walk), we do 10 runs each with 35k random steps (10k for 10 landmarks) and for deterministic approach (triangle, degree, tier 1& 2), we execute the code over the 30 graphs. We compute for each case the score and of resulting landmarks obtained by different methods. We report the average of in 1 and look at the distribution of the measure for 20 landmarks in 5.

Figure 5. Distribution of for 20 landmarks for different landmarks selection approaches
of landmarks 10 (10K) 20 (35K) 30 (35k)
Lazy Rand Walk 5.447742 7.07428 5.2616129
Fully random 3.670968 4.118387 3.841935
Highest degree 1.123871 1.265484 1.106452
Highest centrality 1.16788 1.278206 1.125779
Highest triangle 1.144516 1.209032 0.984516
Mix Tier 1 and Tier 2 1.565161 1.548065 1.620645
Rand. Walk 1.859355 2.0933387 1.830323
Rand. Walk Ripe 1.741935 1.720806 1.823971
Table 1. Average of for different landmarks selection approaches

We can see that the lazy random walk achieves a better performance both in term of diversity and of average distance, we will therefore use the 20 landmarks coming from this method in the forthcoming.

4. Evaluation and validation

The previous sections presented several components of a large-scale graph monitoring system and its application to the particular case of AS level BGP graph. We present in Fig. 4 the components of the monitoring system. In this section, we will assert the performance of our approach and compare it with other approaches.

4.1. Datasets

The first stage of the monitoring system pipeline is the BGP feeds gathering.

For the results shown in this paper we used the python version of libBGPstream, pybgpstream (lib, [n. d.]) and gathered BGP feeds coming from all 13 active RIPE Routing Information Service open feeds (RIP, [n. d.]). The BGP gathering tools use the approach described in section 2.1 to update an AS level graph that is saved each minute into graph snapshot in a GML formatted file. Each trace is initiated with an empty graph, meaning that there is a transient phase during which the AS level graph size increases according to arriving BGP updates. Generally after up to one hour, the graph reaches a steady state with around 68 k ASes. It is noteworthy that even in steady state the graph is permanently updated with new paths that appear and one can never assume that the AS level graph is becoming static. Because of the permanent dynamic of the AS level graph, we are making the analysis over all snapshots regardless of being in transient phase or steady state. As we will see this can result in observing relatively large geometry changes that are arising from gaining visibility for a large part of the topology.

As historical data are available going back for example to 1999 for the rrc00, we have reconstructed the BGP feeds as they were announced during 8 different major BGP even of the past years (see table 2 for a list of monitored events). In addition to these events, we also run the BGP gathering components over long period without known major events.

Date Name Trace duration #ASes #links Extend of the event
23 Jan 2008 Middle East cable cut 7h00 K Medium
25 Aug 2017 Google leakage 3h30 K K Large
21 Oct 2017 Brazil leakage 3h30 K K Large
22 Apr 2016 Innofield incident 1h10 K K Large
12 Dec 2017 Russia hijack 2h00 K K Large
21 Oct 2017 Before Brazil leakage 2h00 K K No event
5 Jan 2018 NA 20h00 K K No event
23 Apr 2016 NA 1h45 K K Medium
22 Apr 2016 After Innofield incident 1h45 K K No event
21 Oct 2017 Before Brazil Leakage 8h30 K K No event
Table 2. List of major BGP events studied in this paper

4.2. Complexity

The different steps of the monitoring system have different complexities. We are currently executing the BGP gathering step on a single computer using a multithreaded (4 threads) python code. We are observing over different trace an average of 1200 BGP announcements per second, that can peak to up to 10 000 during period of high intensity. The current non-optimized python code is running in almost real-time, i.e., one second in the BGP feed is processed in less than one sec, beyond when BGP updates peak times when the current system lags. The memory usage of the tool is relatively large, in particular the routes database might contain up to 12 millions BGP route occupies around 6 GBytes of memory.

The most costly step of the monitoring system is the derivation of the matrix. Let’s remind that this matrix contains curvature elements. As explained before, we have used throughout the paper landmarks. The number of ASes, , that have seen an update or withdrawal during an AS graph snapshot interval varies accordingly to the dynamic of the network from 20 to more than 10000. This means that we need to solve up to 200000 linear programming optimization in the form given in Eq. 1. The main issue is that each one of these optimizations needs a distance matrix between all neighbors of source and destination vertices. Even if the calculating the shortest distance have a linear complexity calculating this matrix 200000 times over a graph with up to 60K nodes is challenging. Currently, using a single processor and with a python code without major optimization, we are able to derive the matrix and the event detection step, in average, in 3 minutes for each snapshot, i.e., we have a 3 minutes of detection delay. We are confident that this time can be reduced with a better implementation in a more suitable language and over a multiprocessor machine.

4.3. Large-scale event detection validation

In this section we will see how the proposed monitoring system behaves in real operational settings. However, before dealing with our proposed system we have to compare it with alternative monitoring systems. Most of BGP monitoring systems target detection of localized event like outage and more importantly BGP hijacks. While these issues are indeed important, they generally have not a global impact. In (Wang et al., 2002) and (Chen et al., 2015) the BGP update volumes were used to detect and characterize large-scale BGP events. Both these papers showed that although the update volume is correlated to the scale of BGP events, this was not enough to construct a reliable detector, because of it is hard to predict and detect BGP sessions resets that can also result in large volume of updates.

Other alternative embedding metrics that could be used in order to monitor changes toward landmarks are the shortest path distances and the spectral distances. The spectral distances are obtained through a spectral projection of the AS level graph (Luxburg et al., 2004; Yehuda et al., 2017). We compare in Fig. 6 the distribution of the variations of these embedding metrics toward some landmarks and compare it with the Ollivier-Ricci embedding we have developed in this paper. We would like embedding metrics to exhibit small variations when strictly local events happen and large variations for large one. This means that the most of the time a suitable embedding metric should be around 0, with large value being less frequent. However, both shortest path distance and spectral distances suffer from high sensitivity to local changes. Any local change on the shortest distance might result at least to +1 or -1 on the embedding metric. Spectral distance is also know to be very sensitive to local neighborhood change, i.e., node with same neighborhood are projected to the same spectral coordinates, and slight change in neighborhood can result in large variation of the spectral distance. The sensitivity of the two alternative metrics can be seen in Fig. 6. The shortest distance metric show large discrete variation, and the spectral distance exhibit also frequently large variation, while the Ollivier-RRicci embedding show a higher concentration of values around 0 drift with less frequent large values. This indicates that the two alternative embedding metrics are likely to not have good performances for graph monitoring applications. We will show in the forthcoming that despite the two others Ollivier-Ricci embedding attain very good performance for monitoring large-scale graphs.

(a) Shortest path embedding
(b) Spectral embedding
(c) Ricci embedding
Figure 6. Distribution of the variations of 3 possible embedding metrics

We illustrate the method for large-xfscale even detection by showing in Fig. (a)a the time evolution of obtained over AS level snapshots gathered around Aug. 25, 2017. As can be seen from this diagram at 03:23 GMT the normalized energy exhibits a large peak. At 03:22 GMT, it has been reported that in Chicago, Illinois, Google generated a massive BGP leaks, that had major consequences for the Internet in Japan, and therefore the peak observed in (a)a is related to this leak. The peak is well separated from other points of the phase diagram showing that Google BGP leak changed severely overall Internet geometry. A more in-depth analysis of the matrix shows a large increase of curvature, i.e., the network becomes more similar to a clique and the overall shortest distances decrease. This is characteristic of a BGP leak that provide a lot of new paths to the network. We show in Fig. (b)b, the phase diagram, a 2 dimensional diagram with in the -axis the normalized energy, , and on the y-axis the inverse of the stable rank, relative to the previous time plot.As can be seen from phase diagram there are several points that all happen after 3:23 GMT with large . All these points have overall decreasing curvature, i.e., the new shortest paths leaked by the Google leakage are withdrawn gradually. However, for all these point , meaning that only a small number of landmarks are affected by the change, making the change local. In other term the global event at 03:22 GMT is resolved by a sequence of more limited local events.

Similar behavior is observed for other large-scale event. We show in Fig. (a)a another event that happened in Oct. 21, 2017, when a major BGP leakage happened in Brazil at 10:09 GMT. Here also we see a large peak of with .

(a) Google Leakage
(b) Phase diagram of the Google leakage
(c) Brazil Leakage
Figure 7. Validation of the large-scale events detection

One typical question that arise is related to the dynamic building of the AS level graph, i.e. new ASes and link are appearing progressively in the graph. One can wonder if the addition of new nodes is not perturbing the detection process, resulting into spurious events that are caused by the collection process rather than real events over Internet. We show in Fig. 8 the time evolution of the number of AS nodes and links in the AS graph along with the evolution of the energy in matrix around a detected large-scale event. As can be seen from the picture at time of large increase in the number of nodes and links there is a slight increase in the , that is normal and expected, as overall by adding new nodes and links we are changing the structure of the network. However, there is no co-occurrence of large changes in with change in the number of nodes or links, meaning that large-scale events are not impacted by the dynamic process of building the AS graph on the fly. This is a remarkable property that show the robustness of the large-scale event detector.

Figure 8. Time evolution the AS level graph AS nodes and links number

4.4. Analysis in the wild

(a) Phase diagram of all datasets
(b) Zoom of the all datasets phase diagram
Figure 9. Phase diagram obtained with all datasets described in Table 2

In this section we will describe outcomes of applying the AS graph monitoring approach over traces containing well-documented large-scale events. We have selected a corpus of diverse incidents whose influences were variable (from a local event as a local prefix hijacking to global events such as an important black-out). We pick 8 events summarized in table 2. There have been 4 major BGP events in 2016 and 2017 that get reported in major blogs monitoring BGP events like the BGPMON blog (Yan et al., 2009). All these events are reported in table 2 as ”large” events. We add to the table some datasets representing periods without any major reported events.

We show in Fig. 9 the phase diagram resulting from combining all datasets in Table 2. We show the full phase diagram in Fig. (a)a and a zoom on the values with normalized energy larger than 1 in Fig. (b)b. As can be seen all major events that have been reported in the major BGP blogs appear in the zoom region. Meaning that all such events generate large spike in their normalized energy. This observation confirms that the large-scale event detector spots correctly all large-scale events. A more detailed analysis show that all events are detected within 2 one minute AS level graph snapshot, except from the ”2008 cable cut” dataset which does not appear in the large-scale event detection region. Nonetheless, we saw for this dataset, several points with large energy (around 0.85) but high (around 0.8) value of that were compatible with large local changes in the underlying geometry. A deeper analysis of the matrix at that times shows a decrease in curvature along with a balancing curvature increase relative to only one or two landmarks. These observations are consistent with reports (see (cab, [n. d.])) that stated that the cable cut was not reported until several day later as its effect was not really felt globally.

4.5. Calibrating the large event detector

In Fig. 9 one can see a clear threshold for detecting large-scale event at . However, setting this threshold needs a little more motivation. The

matrix can be considered as a random matrix with values distributed between 1 and -1. In

(Rudelson and Vershynin, 2009), Rudelson et al. derived, through an estimate of the smallest and largest singular value, an asymptotic bound for the normalized Frobenius norm of a random sub-gaussian matrix. Through this estimate for large enough , an matrix with independent and identically distributed sub-gaussian entries with mean 0 and variance 1, the normalized Frobenius norm will be distributed following a

distribution with 1 degree of freedom. In our case, the curvature values in

are bounded, and therefore there are sub-gaussian and the theoretical conditions are valid.

Now if the entries of the matrix have a mean 0 and a variance equal to , one can expect that will be distributed following a distribution with 1 degree of freedom. We tested if the empirical distribution of is following a distribution with 1 degree of freedom. It is a well known property of distribution with

degree of freedom that it is equivalent with a Gamma distribution with parameters

. We fitted a Gamma distribution to the all values gathered over all datasets and we estimated the parameter to be that is fully compatible with a distribution with 1 degree of freedom showing that the theoretical asymptotic distribution seems to hold in practice.

This asymptotic distribution gives us a way to derive the decision threshold for a large event. By setting a decision threshold at , i.e., if

a large event is detected, the probability of false alarm will be

, where

is the complementary cumulative distribution function of the

distribution with 1 degree of freedom.

Empirically over all matrices we gathered, the elements in the matrix have mean and variance , meaning that with a threshold the probability of false alarm, i.e. wrongly deciding that an observed spike is a large-scale event while it was not, is less than which gives in average a false alarm each 13,5 years (assuming one snapshot by minute). Interestingly even with such a high threshold, we were able to detect all large-scale events reported during the past 3 years without any false alarm. As we achieved a detection rate of 1 without false alarms, we did not show the ROC curve (Fawcett, 2006) as the result would be trivial.

4.6. Interpretation of large curvature changes

In previous sections we described how to monitor and detect large change in the geometry of AS level graphs. When such an event happens, it is desirable to have some ideas about the source of the problem and to implement a root cause analysis. Interestingly the optimal transport plan that is derived during the calculation of Ollivier-Ricci curvature (see Eq. 1), gives a lot of insights about the causes of the curvature changes. The optimal transport plan can be represented as a matrix where each row is relative to a neighbor of the source and each column is relative to a neighbor of the destination . One can therefore compare for 2 given nodes and this matrix before and after an large-scale event. We show in Fig. 10 such a comparison using the difference between two optimal transport plans for the node that have seen the largest curvature change during a large-scale event. The Figure shows clearly that the mass at ASN 23106 that used to be transferred mainly through ASNs 4637 and 8928 (the dark blue elements in Fig. 10) are transferred after the large-scale event by ASNs 6461 and 40620 (the dark red elements in Fig. 10)

Figure 10. Difference between the transport plan

5. Related Works

Anomaly detection of BGP feeds Many attempts have been made to detect anomalous Internet events through dissecting BGP updates and tables. A complete and recent survey of BGP Anomaly detection techniques is provided in (Al-Musawi et al., 2017b). Among all these work (Labovitz et al., 1997) used Fourier analysis tools to analyze routing update rates and to relate them to BGP anomalies. Mei et al. developed BAlet, a wavelet based approach that cluster fast changing BGP feeds to identify possible source of large-scale anomalies (Mai et al., 2008). Wavelet Transforms were also used in the BGP-lens (Prakash et al., 2009) to analyze the rate of updates of per prefix. Other statistical approaches like PCA have also been used (Huang et al., 2007) over BGP update rate.

Several approaches used machine learning techniques to detect link failure or other direct anomalies like BGP Hijacks or indirect one caused by routing overload generated by work propagation for example

(Li and Brooks, 2011; Lutu et al., 2014). However, almost all of them leverage the rate of BGP updates (see for more details (Al-Musawi et al., 2017b)).

Gao et al. (Gao and Rexford, 2001) infers the non-disclosed policies between ASes by investigating their preferred routes. An ”anomaly” is detected when the chosen path does not correspond to the expected policy. This approach is very localized and cannot be extended to the whole graph as it needs one of the monitored ASes to be affected by the ”anomaly” for the detection to happen. Moreover, inferring non-disclosed policies might be hard specially with malicious actors.

While the above approaches target detecting BGP anomalies, they are not dealing with the topic of this paper that is graph monitoring and detecting large-scale event. Moreover, all the above approaches are suffering from the self-similar and heavy tailed structure of update rates caused by BGP session resets that are known to result into large BGP exchanges and inflated update rates. They also do not provide a way to detect further the location and to interpret the anomalies, i.e., to infer the cause relative to AS graph change to the observed surge in update rates.

Large-scale events detection: The previous described approaches were not targeting specifically large-scale events. Comarela et al. (Comarela and Crovella, 2014) similarly to us targeted large-scale. They aggregate next-hop changes into a metric tensor. By looking at the concentration of elements within the tensor, they proposed a large-scale events detection system. However, they use a daily granularity, which is not satisfactory for small time scale detection. The analysis also heavily relies on hop-count changes that is not a stable enough metric to understand the real evolution of the network as we discussed above (see section 4.4). In comparison, we provide an almost real-time event detection that is detecting a larger set of events that the one that hop-count could detect. Chen et al. (Chen et al., 2015), extended this paper and developed the LBE, a metric that aims into locating large-scale events. However, LBE depends mainly on the number of BGP updates happening between two detection time and not on the nature of the update. In (Xu et al., 2005), PCA is used to separate updates triggered by distinct underlying events. This separation simplifies the root cause analysis but remains highly sensitive to false alarms.

Traffic volumes anomaly detection There exists a large literature on anomaly detection in traffic volumes, e.g., (Soule et al., 2005) (Lakhina, 2007). Time series analysis approaches are the most widely spread analytic tool used in the context of anomaly detection. Those approaches rely on the definition of a Linear Time independent (LTI) state space synthetic stochastic model capturing observed temporal and spatial correlations between observed traffic volumes during ”normal” periods. The tracking of this normal model can be done with a Kalmann Filter (Soule et al., 2005) or without it (Lakhina, 2007).

Our approach in this paper differs mainly from the above method by the fact that we are analyzing the graph, a discrete object, rather than a quantitative metric like traffic volume. One can indeed extend method applicable to volume anomaly detection to the normalized energy time series we have developed in this paper. But the difference between large-scale event and local drift is so large in our observations that a simple thresholding was shown to be sufficient.

Ricci curvature on graphs Recent works in Riemannian Geometry from Gromov (Gromov et al., 1999) and Perelman (Perelman, 2002), were extended in the past decade to discrete objects and graphs. This extension use intensively the optimal transport approach developed by Lott and Viliani (Lott and Villani, 2004). Ollivier developed in(Ollivier, 2009) a tractable formula to compute the Ricci curvature over graphs. This attracted a lot of interest in the research community in machine learning and bioinformatics (Whidden and Matsen, 2017)(Pouryahya et al., 2018)(G. Ache and Warren, 2014). In (Weber et al., 2016), Weber et al. presented an approach to obtain a geometric description of highly dynamic networks. They use a different approach of discretization of the Ricci curvature elaborated by Forman (Forman, 2003). The relationship between the two approaches have been explained in (Samal et al., 2017)

. Forman’s curvature depends on a topological description of the Ricci and is defined with respect to polygonal meshes while Ollivier’s curvature is centered around the propagation of Markov Chains and optimal transports characterization of the space. Our assumption is that the AS graph evolves following a latent stochastic process and therefore it is natural for us to favor Ollivier’s description.

6. Conclusion

We have presented a new graph embedding technique for monitoring large scale networks. Our method relies on the curvature, a geometric property linking local scale to global properties. The Ricci curvature that has been initially defined in continuous manifold is extended to graph through the Ollivier-Ricci curvature. We use the Ollivier-Ricci curvature to embed the AS level graph into a -dimensional space. In the embedded space we evaluate the variations of the projection of nodes that have seen a BGP update. This enables us to define a new monitoring system for tracking and detecting large scale events that are changing the overall topology of Internet. We validate the monitoring system over BGP feeds gathered at 20 collectors and we were able to detect in less than 2 minutes all major BGP events that have made the headlines in BGP monitoring blogs.

The described methodology can be easily extended to any type of complex networks and is attractive for purposes beyond the one mentioned in this paper, such as social network analysis, intra and inter protein networks, etc.

One of the main challenge remains to improve the time required for the computation of the optimal transport. This is a very active topic in the machine learning community, see (Seguy et al., 2018) (Solomon et al., 2015). This would allow the methodology to be deployed live and will provide a complementary monitoring tool for Network Operating Centers.

7. appendix

7.1. Landmark Selection algorithms

One obvious way of selecting landmarks is to use the BGP collecting points. But these trivial landmarks will be biased. Nevertheless, collection points have a major advantage as they are accessing an actual BGP feed while for other choice of landmark the AS level topology should be inferred from any other point.

In order to avoid the bias of collection points and to choose suitable landmarks, several approaches are possible. The random choice is appealing because of its simplicity, with the added benefits of security as an attacker may not know which parts of the network he should avoid perturbing to remain undetected. However, random choice in graphs is known to be biased and can result into an over-representation of certain parts of the network.

Another approach leverages the vertex position in the graph. Vertices with high centrality are potentially good landmarks as a large number of shortest paths go through them. Similarly, vertices with the highest degrees are also good candidates for being landmarks.

Finally, in the context of AS level graphs, tier-1 and tier-2 ASes are another set of candidates. Indeed these vertices are also closely related to the high centrality and high degree vertices.

All these above approaches besides the random approach have a concentration bias, i.e., only central part of the graph where most paths converge is accounted (Lee et al., 2006). Moreover, picking highly connected vertices implies more expensive computations for calculating the optimal transport distance. Curvature to these central points are also less sensitive to the effects of events.

In order to deal with these issues, we develop a sampling method to pick good landmarks. To understand intuitively properties of a good set of landmarks, we should go back to the definition of curvature. Let’s first look at the sphere that has a positive curvature. By choosing only two landmarks at opposites poles of the sphere and by monitoring the geodesics from these two poles one can detect any changes happening on the sphere surface. While for negative and flat curvature shapes we need an infinity of landmarks to cover all the surface. In general, the geometry resulting from a real graph cannot be simply interpreted as a sphere. One can assume the graph and its underlying geometry is divided into flat valleys or districts separated by high mountains. Good landmarks should still be positioned in such a way that they have large positive curvature in between them in order to have diverging geodesics and not be positioned in the same valleys or district. We leverage this intuition and inspired by the connection made in (Ollivier, 2009) between Markov Chain, and curvature, we use a random walk approach to find a satisfying set of landmarks, positioned in with high positive curvature place, far apart from each other to ensure they are not in the same region of the graph and randomly disseminated through the graph to avoid a collusion with an attacker .

However, deriving curvature to all other vertices of the network is too complex to be implemented in practice. A necessary condition for a vertex in the graph to have high curvature is to have a relatively large number of neighbors. One can therefore use the landmark set diversity as defined in Eq. 4 as a simple metric to assess possibility of high curvature.

2:Initial Sample
3:Quality Score
4:number of iterations
6:Result Sample
8:procedure Sampling
10:     for  := 1 to  do
11:         Select vertex randomly
13:         Select vertex randomly
15:         Select randomly from [0,1]
22:         if  then
24:              if  then
26:              end if
27:         end if
28:     end for
29:end procedure
Algorithm 1 MCMC landmark selection algorithm

Based on the above considerations, we have developed a Metropolis-Hasting landmarks selection mechanism, we name lazy random walk, that will select landmarks that validate the above constraints. The Metropolis-Hasting algorithm consists of randomly perturbing a sample set according to a proposal distribution . The distribution gives the probability of changing a sample S to a sample S’. We accepts or rejects a new sample set based on a quality score and a variable . The sampling process consists of first choosing randomly among candidate landmarks, one vertex as a candidate for being exchanged with another vertex , i.e., . The vertex is chosen through a random walk among the neighbor of that are not in the set and itself. The sampling probabilities are and , where is the set of neighbors of that are not in plus the vertex and is the set of neighbors of beside that are not in . This sampling mechanism ensures that landmark candidate have not too similar neighborhoods and therefore are in different districts of the graph.

Now let’s define the score function . As explained before, the landmarks should be chosen in such a way that they have large positive curvature in between each other. It is however infeasible to derive for each step of the Hasting-Metropolis the curvatures values between landmarks. We therefore resort to a simpler heuristic. We define as being the set of non-intersecting shortest path from candidate landmarks in to the AS identified as Tier 1 and Tier 2. By non-intersecting, we mean that there exist no points in common between these twi shortest path. The Tier-1 and Tier-2 vertices are chosen because we already know that these vertices will seat in the center of the network, i.e., around the equator for a spherical geometry. Now, by comparing the size of the two sets and , one can assess which of the landmark sets or have the largest positive curvature to the center of the network. This means that we can choose the scoring function as :

With the above-defined function, we can use the classical Hasting-Metropolis rejection scheme described in algorithm 1 and implement the sampling.


  • (1)
  • cab ([n. d.]) [n. d.]. 2008 submarine cable disruption. ([n. d.]). https://en.wikipedia.org/wiki/2008_submarine_cable_disruption
  • jap ([n. d.]) [n. d.]. BGP leak causing Internet outages in Japan and beyond. https://bgpmon.net/bgp-leak-causing-internet-outages-in-japan-and-beyond/. ([n. d.]).
  • cid ([n. d.]) [n. d.]. CIDR. https://www.cidr-report.org/as2.0/. ([n. d.]).
  • dyn ([n. d.]) [n. d.]. Dyn. ([n. d.]). https://dyn.com/dns/
  • lib ([n. d.]) [n. d.]. Libbgpstream. ([n. d.]). https://bgpstream.caida.org/docs/api/libbgpstream
  • RIP ([n. d.]) [n. d.]. Ripe. ([n. d.]). https://www.ripe.net/manage-ips-and-asns/db
  • sin (2014) 2014.

    The Singular Value Decomposition

    Wiley-Blackwell, Chapter 14, 229–241. https://doi.org/10.1002/9781118835647.ch14 arXiv:https://onlinelibrary.wiley.com/doi/pdf/10.1002/9781118835647.ch14
  • Al-Musawi et al. (2017a) B. Al-Musawi, P. Branch, and G. Armitage. 2017a. BGP Anomaly Detection Techniques: A Survey. IEEE Communications Surveys Tutorials 19, 1 (Firstquarter 2017), 377–396. https://doi.org/10.1109/COMST.2016.2622240
  • Al-Musawi et al. (2017b) B. Al-Musawi, P. Branch, and G. Armitage. 2017b. BGP Anomaly Detection Techniques: A Survey. IEEE Communications Surveys Tutorials 19, 1 (Firstquarter 2017), 377–396. https://doi.org/10.1109/COMST.2016.2622240
  • Bakry and Émery (1985) Dominique Bakry and M Émery. 1985. Diffusions hypercontractives. 84 (01 1985).
  • Ballani et al. (2007) Hitesh Ballani, Paul Francis, and Xinyang Zhang. 2007. A Study of Prefix Hijacking and Interception in the Internet. SIGCOMM Comput. Commun. Rev. 37, 4 (Aug. 2007), 265–276. https://doi.org/10.1145/1282427.1282411
  • Berger (2003) Marcel Berger. 2003. From Curvature to Topology. Springer Berlin Heidelberg, Berlin, Heidelberg, 543–635. https://doi.org/10.1007/978-3-642-18245-7_12
  • Berger and Gostiaux (1987) M. Berger and B. Gostiaux. 1987. Geometrie differentielle: varietes, courbes et surfaces. Presses universitaires de France. https://books.google.com.au/books?id=bkPvAAAAMAAJ
  • Bobkov and Madiman (2013) Sergey G Bobkov and Mokshay M Madiman. 2013. On the problem of reversibility of the entropy power inequality. In Limit theorems in probability, statistics and number theory. Springer, 61–74.
  • Bonaventure (2002) Olivier Bonaventure. 2002. Interdomain routing with BGP: Issues and challenges. (01 2002).
  • Butler et al. (2010) K. Butler, T. R. Farley, P. McDaniel, and J. Rexford. 2010. A Survey of BGP Security Issues and Solutions. Proc. IEEE 98, 1 (Jan 2010), 100–122. https://doi.org/10.1109/JPROC.2009.2034031
  • Caesar and Rexford (2005) M. Caesar and J. Rexford. 2005. BGP routing policies in ISP networks. IEEE Network 19, 6 (Nov 2005), 5–11. https://doi.org/10.1109/MNET.2005.1541715
  • Cai et al. (2018) H. Cai, V. W. Zheng, and K. Chang. 2018. A Comprehensive Survey of Graph Embedding: Problems, Techniques and Applications. IEEE Transactions on Knowledge and Data Engineering (2018), 1–1. https://doi.org/10.1109/TKDE.2018.2807452
  • Chen et al. (2015) Meng Chen, Mingwei Xu, Qing Li, Xirui Song, and Yuan Yang. 2015. Detect and analyze Large-scale BGP events by bi-clustering Update Visibility Matrix. In 2015 IEEE 34th International Performance Computing and Communications Conference (IPCCC). 1–8. https://doi.org/10.1109/PCCC.2015.7410284
  • Chen et al. (2002) Qian Chen, Hyunseok Chang, R. Govindan, and S. Jamin. 2002. The origin of power laws in Internet topologies revisited. In Proceedings.Twenty-First Annual Joint Conference of the IEEE Computer and Communications Societies, Vol. 2. 608–617 vol.2. https://doi.org/10.1109/INFCOM.2002.1019306
  • Chen et al. (2012) Wei Chen, Wenjie Fang, Guangda Hu, and Michael W. Mahoney. 2012. On the Hyperbolicity of Small-World and Tree-Like Random Graphs. (2012). arXiv:arXiv:1201.1717
  • claffy (2001) k. claffy. 2001. CAIDA: Visualizing the Internet. Internet Computing Online (Jan 2001), 88.
  • Cohen et al. (2016) Michael B. Cohen, Jelani Nelson, and David P. Woodruff. 2016. Optimal Approximate Matrix Product in Terms of Stable Rank. In 43rd International Colloquium on Automata, Languages, and Programming (ICALP 2016) (Leibniz International Proceedings in Informatics (LIPIcs)), Ioannis Chatzigiannakis, Michael Mitzenmacher, Yuval Rabani, and Davide Sangiorgi (Eds.), Vol. 55. Schloss Dagstuhl–Leibniz-Zentrum fuer Informatik, Dagstuhl, Germany, 11:1–11:14. https://doi.org/10.4230/LIPIcs.ICALP.2016.11
  • Comarela and Crovella (2014) Giovanni Comarela and Mark Crovella. 2014. Identifying and Analyzing High Impact Routing Events with PathMiner. In Proceedings of the 2014 Conference on Internet Measurement Conference (IMC ’14). ACM, New York, NY, USA, 421–434. https://doi.org/10.1145/2663716.2663754
  • Cuturi (2013) Marco Cuturi. 2013. Sinkhorn Distances: Lightspeed Computation of Optimal Transport. In Advances in Neural Information Processing Systems 26, C. J. C. Burges, L. Bottou, M. Welling, Z. Ghahramani, and K. Q. Weinberger (Eds.). Curran Associates, Inc., 2292–2300. http://papers.nips.cc/paper/4927-sinkhorn-distances-lightspeed-computation-of-optimal-transport.pdf
  • Fawcett (2006) Tom Fawcett. 2006. An Introduction to ROC Analysis. Pattern Recogn. Lett. 27, 8 (June 2006), 861–874. https://doi.org/10.1016/j.patrec.2005.10.010
  • Forman (2003) Forman. 2003. Bochner’s Method for Cell Complexes and Combinatorial Ricci Curvature. Discrete & Computational Geometry 29, 3 (01 Feb 2003), 323–374. https://doi.org/10.1007/s00454-002-0743-x
  • G. Ache and Warren (2014) Antonio G. Ache and Micah Warren. 2014. Coarse Ricci curvature and the manifold learning problem. (10 2014).
  • Gao and Rexford (2001) Lixin Gao and Jennifer Rexford. 2001. Stable Internet Routing Without Global Coordination. IEEE/ACM Trans. Netw. 9, 6 (Dec. 2001), 681–692. https://doi.org/10.1109/90.974523
  • Gromov et al. (1999) M. Gromov, J. LaFontaine, and P. Pansu. 1999. Metric Structures for Riemannian and Non-Riemannian Spaces. Birkhäuser. https://books.google.com.au/books?id=XUfvAAAAMAAJ
  • Huang et al. (2007) Yiyi Huang, Nick Feamster, Anukool Lakhina, and Jim (Jun) Xu. 2007. Diagnosing Network Disruptions with Network-wide Analysis. SIGMETRICS Perform. Eval. Rev. 35, 1 (June 2007), 61–72. https://doi.org/10.1145/1269899.1254890
  • Jost and Liu (2014) Jürgen Jost and Shiping Liu. 2014. Ollivier’s Ricci Curvature, Local Clustering and Curvature-Dimension Inequalities on Graphs. Discrete Comput. Geom. 51, 2 (March 2014), 300–322. https://doi.org/10.1007/s00454-013-9558-1
  • Karlin et al. (2008) Josh Karlin, Stephanie Forrest, and Jennifer Rexford. 2008. Autonomous security for autonomous systems. Computer Networks 52, 15 (2008), 2908 – 2923. https://doi.org/10.1016/j.comnet.2008.06.012 Complex Computer and Communication Networks.
  • Knight et al. (2011) S. Knight, H. X. Nguyen, N. Falkner, R. Bowden, and M. Roughan. 2011. The Internet Topology Zoo. IEEE Journal on Selected Areas in Communications 29, 9 (October 2011), 1765–1775. https://doi.org/10.1109/JSAC.2011.111002
  • Kruegel et al. (2003) Christopher Kruegel, Darren Mutz, William Robertson, and Fredrik Valeur. 2003. Topology-Based Detection of Anomalous BGP Messages. In Recent Advances in Intrusion Detection, Giovanni Vigna, Christopher Kruegel, and Erland Jonsson (Eds.). Springer Berlin Heidelberg, Berlin, Heidelberg, 17–35.
  • Labovitz et al. (1997) Craig Labovitz, G. Robert Malan, and Farnam Jahanian. 1997. Internet Routing Instability. SIGCOMM Comput. Commun. Rev. 27, 4 (Oct. 1997), 115–126. https://doi.org/10.1145/263109.263151
  • Lafontaine (2015) Jacques Lafontaine. 2015. The Euler-Poincaré Characteristic and the Gauss-Bonnet Theorem. Springer International Publishing, Cham, 323–348. https://doi.org/10.1007/978-3-319-20735-3_8
  • Lakhina (2007) Anukool Lakhina. 2007. Network-wide Traffic Analysis: Methods and Applications. Ph.D. Dissertation. Boston, MA, USA. Advisor(s) Crovella, Mark E. AAI3232904.
  • Lee et al. (2006) Sang Hoon Lee, Pan-Jun Kim, and Hawoong Jeong. 2006. Statistical properties of sampled networks. Phys. Rev. E 73 (Jan 2006), 016102. Issue 1. https://doi.org/10.1103/PhysRevE.73.016102
  • Leskovec and Faloutsos (2006) Jure Leskovec and Christos Faloutsos. 2006. Sampling from Large Graphs. In Proceedings of the 12th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (KDD ’06). ACM, New York, NY, USA, 631–636. https://doi.org/10.1145/1150402.1150479
  • Lévy and Schwindt (2018) Bruno Lévy and Erica L. Schwindt. 2018. Notions of optimal transport theory and how to implement them on a computer. Computers & Graphics 72 (2018), 135 – 148. https://doi.org/10.1016/j.cag.2018.01.009
  • Li and Brooks (2011) J. Li and S. Brooks. 2011. I-seismograph: Observing and measuring Internet earthquakes. In 2011 Proceedings IEEE INFOCOM. 2624–2632. https://doi.org/10.1109/INFCOM.2011.5935090
  • Lott and Villani (2004) John Lott and Cedric Villani. 2004. Ricci curvature for metric-measure spaces via optimal transport. (2004). arXiv:arXiv:math/0412127
  • Lutu et al. (2014) A. Lutu, M. Bagnulo, J. Cid-Sueiro, and O. Maennel. 2014. Separating wheat from chaff: Winnowing unintended prefixes using machine learning. In IEEE INFOCOM 2014 - IEEE Conference on Computer Communications. 943–951. https://doi.org/10.1109/INFOCOM.2014.6848023
  • Luxburg et al. (2004) Ulrike von Luxburg, Olivier Bousquet, and Mikhail Belkin. 2004.

    Limits of Spectral Clustering. In

    Proceedings of the 17th International Conference on Neural Information Processing Systems (NIPS’04). MIT Press, Cambridge, MA, USA, 857–864.
  • Mai et al. (2008) J. Mai, L. Yuan, and C. N. Chuah. 2008. Detecting BGP anomalies with wavelet. In NOMS 2008 - 2008 IEEE Network Operations and Management Symposium. 465–472. https://doi.org/10.1109/NOMS.2008.4575169
  • Ollivier (2009) Yann Ollivier. 2009. Ricci curvature of Markov chains on metric spaces. Journal of Functional Analysis 256, 3 (2009), 810–864.
  • Ollivier (2012) Yann Ollivier. 2012. A visual introduction to Riemannian curvatures and some discrete generalizations.
  • Orsini et al. (2016) Chiara Orsini, Alistair King, Danilo Giordano, Vasileios Giotsas, and Alberto Dainotti. 2016. BGPStream: A Software Framework for Live and Historical BGP Data Analysis. In Proceedings of the 2016 Internet Measurement Conference (IMC ’16). ACM, New York, NY, USA, 429–444. https://doi.org/10.1145/2987443.2987482
  • Perelman (2002) Grisha Perelman. 2002. The entropy formula for the Ricci flow and its geometric applications. (2002). arXiv:arXiv:math/0211159
  • Piccoli and Rossi (2012) Benedetto Piccoli and Francesco Rossi. 2012. Generalized Wasserstein distance and its application to transport equations with source. (2012). https://doi.org/10.1007/s00205-013-0669-x arXiv:arXiv:1206.3219
  • Pouryahya et al. (2018) Maryam Pouryahya, Jung Hun Oh, James C. Mathews, Joseph O. Deasy, and Allen R. Tannenbaum. 2018. Characterizing Cancer Drug Response and Biological Correlates: A Geometric Network Approach. In Scientific Reports.
  • Prakash et al. (2009) B. Aditya Prakash, Nicholas Valler, David Andersen, Michalis Faloutsos, and Christos Faloutsos. 2009. BGP-lens: Patterns and Anomalies in Internet Routing Updates. In Proceedings of the 15th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (KDD ’09). ACM, New York, NY, USA, 1315–1324. https://doi.org/10.1145/1557019.1557160
  • Rudelson and Vershynin (2009) Mark Rudelson and Roman Vershynin. 2009. Smallest singular value of a random rectangular matrix. Comm. Pure Appl. Math (2009), 1707–1739.
  • Samal et al. (2017) Areejit Samal, R. P. Sreejith, Jiao Gu, Shiping Liu, Emil Saucan, and Jürgen Jost. 2017. Comparative analysis of two discretizations of Ricci curvature for complex networks. (2017). arXiv:arXiv:1712.07600
  • Seguy et al. (2018) Vivien Seguy, Bharath Bhushan Damodaran, Remi Flamary, Nicolas Courty, Antoine Rolet, and Mathieu Blondel. 2018. Large Scale Optimal Transport and Mapping Estimation. In International Conference on Learning Representations. https://openreview.net/forum?id=B1zlp1bRW
  • Sermpezis et al. (2018) Pavlos Sermpezis, Vasileios Kotronis, Alberto Dainotti, and Xenofontas Dimitropoulos. 2018. A Survey among Network Operators on BGP Prefix Hijacking. (2018). arXiv:arXiv:18
  • Solomon et al. (2015) Justin Solomon, Fernando de Goes, Gabriel Peyré, Marco Cuturi, Adrian Butscher, Andy Nguyen, Tao Du, and Leonidas Guibas. 2015. Convolutional Wasserstein Distances: Efficient Optimal Transportation on Geometric Domains. ACM Trans. Graph. 34, 4, Article 66 (July 2015), 11 pages. https://doi.org/10.1145/2766963
  • Soule et al. (2005) Augustin Soule, Kavé Salamatian, and Nina Taft. 2005. Combining Filtering and Statistical Methods for Anomaly Detection. In Proceedings of the 5th ACM SIGCOMM Conference on Internet Measurement (IMC ’05). USENIX Association, Berkeley, CA, USA, 31–31. http://dl.acm.org/citation.cfm?id=1251086.1251117
  • Wählisch et al. (2012) Matthias Wählisch, Olaf Maennel, and Thomas C. Schmidt. 2012. Towards detecting BGP route hijacking using the RPKI. In SIGCOMM.
  • Wang et al. (2002) Lan Wang, Xiaoliang Zhao, Dan Pei, Randy Bush, Daniel Massey, Allison Mankin, S. Felix Wu, and Lixia Zhang. 2002. Observation and Analysis of BGP Behavior Under Stress. In Proceedings of the 2Nd ACM SIGCOMM Workshop on Internet Measurment (IMW ’02). ACM, New York, NY, USA, 183–195. https://doi.org/10.1145/637201.637231
  • Weber et al. (2016) Melanie Weber, Emil Saucan, and Jürgen Jost. 2016. Can one see the shape of a network? CoRR abs/1608.07838 (2016).
  • Whidden and Matsen (2017) Chris Whidden and Frederick A. Matsen. 2017. Ricci–Ollivier curvature of the rooted phylogenetic subtree–prune–regraft graph. Theoretical Computer Science 699 (2017), 1 – 20. https://doi.org/10.1016/j.tcs.2017.02.006 Special issue on Analytic Algorithmics and Combinatorics.
  • Xu et al. (2005) Kuai Xu, Jaideep Chandrashekar, and Zhi Li Zhang. 2005. A first step toward understanding inter-domain routing dynamics. 207–212.
  • Yan et al. (2009) H. Yan, R. Oliveira, K. Burnett, D. Matthews, L. Zhang, and D. Massey. 2009. BGPmon: A Real-Time, Scalable, Extensible Monitoring System. In 2009 Cybersecurity Applications Technology Conference for Homeland Security. 212–223. https://doi.org/10.1109/CATCH.2009.28
  • Yehuda et al. (2017) G. Yehuda, D. Keren, and I. Akaria. 2017. Monitoring Properties of Large, Distributed, Dynamic Graphs. In 2017 IEEE International Parallel and Distributed Processing Symposium (IPDPS). 2–11. https://doi.org/10.1109/IPDPS.2017.123
  • Zhang et al. (2004) Ke Zhang, Xiaoliang Zhao, and S. F. Wu. 2004. An analysis on selective dropping attack in BGP. In IEEE International Conference on Performance, Computing, and Communications, 2004. 593–599. https://doi.org/10.1109/PCCC.2004.1395106