Graph Comparison via the Non-backtracking Spectrum

TL;DR

This paper introduces a new graph comparison metric based on the non-backtracking spectrum, enabling reliable comparison of graphs of different sizes and types, with applications in classification tasks.

Contribution

It proposes a novel spectral metric using the non-backtracking operator, addressing limitations of existing methods in graph comparison and scalability.

Findings

The metric effectively compares graphs of different sizes.

Watts-Strogatz graphs form a manifold in the spectral embedding.

The metric improves graph classification performance.

Abstract

The comparison of graphs is a vitally important, yet difficult task which arises across a number of diverse research areas including biological and social networks. There have been a number of approaches to define graph distance however often these are not metrics (rendering standard data-mining techniques infeasible), or are computationally infeasible for large graphs. In this work we define a new metric based on the spectrum of the non-backtracking graph operator and show that it can not only be used to compare graphs generated through different mechanisms, but can reliably compare graphs of varying size. We observe that the family of Watts-Strogatz graphs lie on a manifold in the non-backtracking spectral embedding and show how this metric can be used in a standard classification problem of empirical graphs.

Tables

Table 1. Table 1: Empirical datasets consider in this section, in addition to three synthetic graph models for comparison. Each dataset contains a collection of graphs of varying size, however within each collection all graphs are assumed to be generated from the same mechanism. For synthetic graphs we consider graphs with mean size 200 200 200 and variance 40 40 40 . ER graphs have rewiring parameter p = 0.1 𝑝 0.1 p=0.1 ., the WS graphs have k = 20 𝑘 20 k=20 and p = 0.1 𝑝 0.1 p=0.1 , and the BA graphs have incoming vertices connect to m = 2 𝑚 2 m=2 vertices with the system initialised with a single vertex. All graphs are taken to be simple and undirected.

Dataset		Count	$\mathrm{𝐦𝐚𝐱}\left(𝒏\right)$	$\mathrm{𝐦𝐢𝐧}\left(𝒏\right)$
Facebook100 (FB)	[41]	26	4563	769
World Subways (SW)	[36]	15	433	82
Autonomous Systems (AS)	[27]	30	3144	2948
Metabolic Networks (MN)	[22]	10	1593	505
Erdős-Rényi (ER)	[32]	30	244	151
Watts-Strogatz (WS)	[32]	30	244	151
Barabási-Albert (BA)	[32]	30	244	151

Table 2. Table 2: The clustering performance of the distributional non-backtracking spectral distance (d-NBD), the non-backtracking spectral distance (NBD), and the Laplacian (LAP) distance. Results are averaged over a 10-fold cross-validation. For the graph collections considered, the d-NBD performs without fault achieving perfect accuracy, precision, and recall. Both the NBD and LAP clusterings also show relatively high accuracy but a number of graphs are consistently labelled incorrectly.

	Train			Test
	Rec.	Prec.	Acc.	Rec.	Prec.	Acc.
d-NBD	1.00	1.00	1.00	1.00	1.00	1.00
NBD	0.96	0.92	0.94	0.97	0.95	0.95
LAP	0.92	0.80	0.84	0.92	0.79	0.84

Equations

$$B_{(u\to v),(x\to y)}=\{\begin{array}{ll}1& \text{mbox}ifx=v\text{mbox}andy\ne u\\ 0& \text{mbox}otherwise.\end{array}$$

$$\mathrm{tr}{A}^k=\sum _{i=1}^n{\lambda }_i^k$$

$$\mathrm{tr}{B}^k=\sum _{i=1}^{2m}{\lambda }_i^k,$$

$$\det (B-\lambda I_n)=({\lambda }^2-1{)}^{|E|-|V|}\det (Q_{\lambda })$$

$$\displaystyle B^{\prime}=\left(\begin{array}[]{cc}A&I-D\\ I&0\end{array}\right)$$

$$\det (Q_{\lambda })=\det (B^{\prime }-\lambda I_n).$$

$$\vec{v}=[{\alpha }_1,\dots ,{\alpha }_k,{\beta }_1,\dots ,{\beta }_k],$$

$$\displaystyle B^{\prime}=\left(\begin{array}[]{cc}A&I-D\\ I&0\end{array}\right)$$

$$|{\lambda }_i|\to {\log }_{\rho (B)}|{\lambda }_i|.$$

$$F(r,\theta )=\frac{1}{2n}\sum _{i=1}^{2n}\text{mathds}{1}_{\{|{\hat{\lambda }}_i|\le r\}}\text{mathds}{1}_{\{0\le \arg ({\hat{\lambda }}_i)\le \theta \}}$$

$$d(G_1,G_2)=\frac{1}{\pi }{\left({\int }_0^{\pi }{\int }_0^1|F_1(r,\theta )-F_2(r,\theta ){|}^{2}drd\theta \right)}^{\frac{1}{2}},$$

$$\vec{v}=[$$

$$F(r_1,{\theta }_2),\dots ,F(r_k,{\theta }_2),$$

$$\dots$$

$$F(r_k,{\theta }_1)\dots ,F(r_k,{\theta }_k)]$$

Full text

1 Introduction

Comparing graph structures is a fundamental task in graph theory. In particular the need to identify similar graph structure arises in order to determine common function, behaviour, or generative process. This need is universal across disciplines and applications range from image processing [11], chemistry [1, 26], and social network analysis [24, 29].

For small graphs, identifying two identical graph structures (up to an isomorphism, or relabelling of vertices) is a trivial and computationally feasible task due to the limited number of possible configurations. For graphs with four vertices there are six possible graph configurations. However, for graphs with ten vertices there are $11,716,571$ possible configurations. Although there is not a closed formula for the number of graph configurations with $n$ vertices, it is clear that it soon becomes intractable to enumerate or compare large graphs through brute-force.

The problem of determining whether two graphs are isomorphic is equivalent to finding a permutation matrix $P$ such that $AP=PB$ where $A$ and $B$ are the adjacency matrices of the two graphs. Despite the difficulty in the task111 It is currently unknown as to whether this problem lies in the class of P or NP-hard problems [4]

there have been a number of attempts to define graph distance by minimising $||AP-PB||$ for a suitably defined matrix norm $||\cdot||$. The matrices $A$ and $B$ can also take different forms to shift from a local perspective (adjacency matrix) to consider global properties, where $A$ and $B$ capture the number of paths between vertices. Examples of such distances include the chemical distance [26], CKS distance [8], and edit distance [18, 39]. Recent work has focused on reducing the space over which to find a permutation matrix [6] which makes the problem tractable, although still computationally restrictive for large graphs.

Beyond scalability issues, one of the drawbacks of these methods is that they require graphs to be of the same order (i.e. the same number of vertices). For many applications this may be desirable, especially where the addition of a single vertex or edge can have drastic consequences on graph function (in graphs of chemical compounds for example). However for other applications we may be interested in the structure at a coarser level. For example, a ring graph of $100$ vertices is structurally similar to a ring graph of $200$ vertices, however the two cannot be compared using traditional isomorphism arguments.

Contrasting approaches use a number of graph properties, or features, to characterise graph structure [43, 17, 23, 35, 33]. The most promising of these make use of the graph spectrum, either of the adjacency matrix directly, or of a version of the graph Laplacian [13, 44, 46, 16, 28]. The spectra of graphs or graph operators are useful since the eigenvalues and eigenvectors characterise the topological structure of a graph in a way which can be interpreted physically. Eigenvalues of the graph Laplacian describe how a quantity (information, heat, people, etc.) localised at a vertex can spread across the graph. These eigenvalues also dictate the stability of dynamics acting across the graph [32]. Beyond the Laplacian, recent work has investigated the spectral properties of higher-order operators such as the non-backtracking graph operator (described in the next section) [40]. For a more detailed (but not complete) examination of graph distances we refer the reader to a recent survey [15].

In this article we present a new method to compare graph structure using the distribution of eigenvalues of the non-backtracking operator in the complex plane. In Section 2 we describe the non-backtracking operator and investigate some of its spectral properties before introducing the distributional non-backtracking spectral distance (d-NBD) in Section 3. We show empirical results for both synthetic and real graphs in Section 4 before discussing the significance of the results and future research in Section 5.

2 The Non-backtracking Operator and Spectrum

Consider a simple undirected binary graph $G=(V,E)$ where $V\subset\mathbb{N}$ is a set of vertices and $E\subset V^{2}$ is a set of edges. Here we use the word simple to mean that the graph contains no self-loops or multiple edges between vertices. Let $n=|V|$ be the number of vertices and $m=|E|$ be the number of edges.

Spectral analysis of graphs is typically conducted using the Laplacian matrix $L=A-D$ where $A$ is the adjacency matrix which encodes $G$, and $D$ is a diagonal matrix of vertex degrees with $D_{ii}=\sum_{j}A_{ij}$ and zeros elsewhere. In this work we chose instead to use a different linear graph operator, namely the non-backtracking operator. The non-backtracking (or Hashimoto) matrix $B$ is a $2m\times 2m$ matrix defined on the set of directed edges of $G$. Here each undirected edge $u\leftrightarrow v$ is replaced by two directed edges $u\to v$ and $v\to u$. The non-backtracking matrix is then given by

[TABLE]

The non-backtracking matrix is asymmetric and the entries of $B^{k}$ describe the number of non-backtracking walks of length $k$ across $G$. Of special note is $(B^{k})_{ii}$ which counts the non-backtracking cycles of length $k$ which start and end at edge $i$222Non-backtracking does not mean that a walk only visits distinct vertices, only that it does not immediately return to a vertex it has just arrived from.. The total number of non-backtracking cycles of length $k$ is therefore captured in the trace of $B^{k}$.

2.1 Spectral Properties

The spectral properties of the non-backtracking matrix are well understood, especially for regular graphs [31, 3, 38]. As $B$ is asymmetric the eigenvalues can also take complex values. In Figure 1 we show the spectrum in the complex plane of three real-world graphs and three synthetic graphs of varying size.

The first observation is that the bulk of the eigenvalues $\lambda_{k}$ lie within the disk $|\lambda_{k}|<\sqrt{\rho(B)}$ where $\sqrt{\rho(B)}=\max_{k}|\lambda_{k}|$ is the spectral radius of $B$. This is commonly used as a heuristic for the bulk of the spectrum, however in particular cases more is known. For random regular graphs of degree $d$ there are $n$ pairs of conjugate eigenvalues which lie exactly on a circle of radius $\sqrt{d-1}$, and for a stochastic block model the expected magnitude of the eigenvalues can be shown to be less than or equal to $\sqrt{\rho(B)}$ in the limit $n\to\infty$ [25].

Another property of the spectrum of $B$ is that any tree structure (whether connected to the graph, or disconnected) contributes zero eigenvalues to the spectrum. This is a consequence of any non-backtracking walk becoming ‘stuck’ at the leaves of the tree. Furthermore any unicyclic component gives rise to eigenvalues in $\{0\}\cup\{\lambda:|\lambda|=1\}$. These two properties are easily explained as we make the connection between the number of non-backtracking cycles in the graph and the eigenvalues of $B$.

For any matrix $A$ the trace of $A^{k}$ can be calculated as

[TABLE]

where $(\lambda_{i})_{i=1}^{n}$ are the eigenvalues of $A$. This result is trivial for symmetric matrices which can be diagonalised, however the result also holds for asymmetric matrices. We can therefore write

[TABLE]

noting that $\operatorname{tr}{B^{k}}$ captures a count of all non-backtracking cycles of length $k$. This means that the eigenvalues of $B$ encode all information regarding the number of non-backtracking cycles in $A$. Using this relation it is then clear that tree structures contribute zero eigenvalues since they cannot be part of a non-backtracking cycle.

The last property of the spectrum we detail is derived from the Ihara determinant formula [5, 20, 21]. For any finite and undirected graph

[TABLE]

where $Q_{\lambda}=(D-I_{n})-\lambda A+\lambda^{2}I_{n}$ and $D$ is the matrix with vertex degrees on the diagonal with zeros elsewhere. Using a linearisation of the quadratic polynomial $Q_{\lambda}$333This is in the form of a quadratic eigenvalue problem (QEP). one can show that the $2n\times 2n$ matrix

[TABLE]

shares the same eigenvalues as $Q_{\lambda}$, i.e.,

[TABLE]

Therefore all eigenvalues of $B$ save for $|E|-|V|$ eigenvalues of $\pm 1$ are eigenvalues of the smaller matrix $B^{\prime}$. This result has very practical implications for calculating the spectrum of $B$ given that $n$ can often be substantially smaller than $m$. The matrix $B^{\prime}$ provides a means to calculate a significant proportion of the spectrum of $B$ in a fashion that scales linearly with the total number of vertices.

3 Distributional Non-backtracking Spectral Distance

Much like the spectrum of the graph Laplacian has found use in community detection [14], dynamics on graphs [32] and graph clustering [34] our aim is to exploit the non-backtracking spectrum to differentiate between graph structures. Recent work [10, 40] has shown that the two-cores extracted from two graphs are isomorphic when the set of all non-backtracking cycles in each graph are equal (this is referred to as the length spectrum). This means that should we be able to enumerate all possible non-backtracking cycles we can effectively compare two graph structures. In practice enumerating all possible cycles of all possible lengths is infeasible, and even so it remains unclear how two such sets should be compared.

Torres et al. [40] address this issue by considering a ‘relaxed’ length spectrum, calculating the number of cycles of length $k$ as opposed to the set of cycles themselves. As shown in (1) the number of these cycles is captured completely by the spectrum of $B$. Based on this argument the distance between the number of cycles (and therefore eigenvalues) should provide a reasonable approximation to the graph isomorphism.

In [34] the number of non-backtracking cycles of length $k$, $c_{k}$ are used directly as an embedding of the graph. Each graph is embedded as a vector $\vec{v}_{1}=[c_{3},c_{4},c_{5},c_{6},c_{7},\ln(c_{2m})]$ with the graph distance subsequently defined as the Euclidean distance between vectors. These distances correlate well with the graph edit distance and perform marginally better than a truncated Laplacian spectrum embedding when applied to a computer vision task.

In contrast, the approach in [40] uses the largest $k$ eigenvalues of $B$ (in order of magnitude) to construct a feature vector. In particular for an ordered sequence of eigenvalues $(\lambda_{i})_{i=1}^{2m}$ the graph is embedded as

[TABLE]

where $\lambda_{s}=\alpha_{s}+\beta_{s}i$ and again defining graph distance using the Euclidean distance between embeddings. Taking the largest eigenvalues is an intuitive step as these contribute the most to the cycle counts. There are however many instances where the magnitude of the eigenvalues are approximately equal (seen in Figure 1) and so the method is heavily sensitive to the choice of $k$ and the method of eigenvalue calculation. This is especially problematic for a $d$-regular random graph where all non-real eigenvalues of $B$ are conjugate pairs on the circle of radius $\sqrt{d-1}$.

The two previous approaches to non-backtracking distance are limited in their ability compare graphs of different size. Consider two graphs of sizes $n_{1}$ and $n_{2}$. Naturally as we increase the number of vertices and edges the number of possible cycles will increase and therefore any cycle counting distance will diverge as ${|n_{1}-n_{2}|\to\infty}$. Furthermore, in order to compare the eigenvalues of the non-backtracking matrices directly we are restricted to choosing the largest $\min\{2n_{1},2n_{2}\}$ eigenvalues. We are therefore neglecting a potentially significant part of the spectrum for the larger graph, i.e. we are using a truncated spectrum. The issues surrounding the truncated spectrum become clear for $d$-regular random graphs. Since the bulk of all eigenvalues like on the circle $|\lambda|=\sqrt{d-1}$ we cannot order them by magnitude. We therefore need a secondary ordering, say in descending order of $\Re(\lambda)$. The full spectrum will cover all eigenvalues on the circle, however the truncated spectrum will lie on an arc of $|\lambda|=\sqrt{d-1}$. An effective comparison is therefore not possible.

This motivates us to consider a property of the spectrum that can be compared across graphs of different size, namely the spectral distribution. In this sense we hope to cluster graphs which share the same properties of generating mechanism, regardless of size. For computational efficiency we consider the eigenvalues of

[TABLE]

where now $A$ is the two-core of the original graph (vertices of degree one are iteratively removed), and $D$ is the corresponding degree matrix. This means that we are operating on a reduced spectrum with no zero eigenvalues and $|E|-|V|$ eigenvalues of $1$ or $-1$ omitted. The spectral radius of $B^{\prime}$ scales with the largest vertex degree, however the corresponding eigenvalues are not localised around high degree vertices to the same extent as they are with the graph Laplacian. This is due to walks on $B^{\prime}$ not being permitted to return to high degree vertices immediately after leaving them. This suggests a rescaling of the eigenvalues

[TABLE]

The rescaled eigenvalues $(\hat{\lambda})_{i=1}^{2n}$ then lie exclusively in the disk $|\hat{\lambda}|\leq 1$ with the bulk of the eigenvalues distributed in $|\hat{\lambda}|\leq\frac{1}{2}$.

Figure 2 shows the rescaled spectrum for Erdős-Rényi graphs of varying size with edge connection probabilities $p=0.2$ (top) and $p=0.8$ (bottom).

Here we can clearly distinguish between the two parameter regimes and we see that these distributions are consistent across graphs irrespective of size.

We capture the distribution of the eigenvalues through the empirical cumulative spectral density given by

[TABLE]

This is defined on $r\in[0,1]$, $\theta\in[0,\pi]$ and for $n>0$. A special case to mention is when the two-core of a graph is empty. In this case we set $F(r,\theta)=0$ $\forall r,\theta$ by convention444 Alternatively we could set $F(r,\theta)=\infty$ so that all distances to the empty graph are infinite. . Since the complex eigenvalues of $B$ (and therefore $B^{\prime}$) are conjugate pairs we gain no further information by considering the lower half of the complex plane and is hence omitted. This naturally leads us to define the distributional non-backtracking spectral distance (d-NBD) as the distance between empirical spectral densities. For two graphs $G_{1}$ and $G_{2}$ with rescaled spectral densities $F_{1}$ and $F_{2}$ respectively, the d-NBD is given by

[TABLE]

i.e. $d(G_{1},G_{2})$ is proportional to the $p$-norm $||F_{1}-F_{2}||_{p}$ with $p=2$.

The d-NBD satisfies many desirable properties. It is non-negative, symmetric, and satisfies the triangle inequality (via the Minkowski inequality). While it is evident that $d(G,G)=0$ it is not the case that $d(G_{1},G_{2})=0$ implies that $G_{1}=G_{2}$, only that they share the same two-core. The d-NBD is therefore a pseudometric.

Finally we define a $k^{2}$-dimensional embedding of a graph through a discretization of the empirical spectral density. Here a graph is represented by

[TABLE]

where $(r_{1},\theta_{1})=(0,0)$ and $(r_{k},\theta_{k})=(1,\pi)$. This embedding is useful in the next section to visualise the d-NBD by plotting the graph embeddings in low-dimensional space.

4 Results

In this section we give a number of results from both synthetic and real-world graph examples.

4.1 Synthetic Graphs

In order to investigate the properties of the d-NBD we consider three random graph models; the Erdős-Rényi graph (ER), the Watts-Strogatz graph (WS), and the $d$-regular random graph (RR). One of the main advantages of the d-NBD over other graph distances is that ability to compare graphs of varying size. Figure 3 illustrates the stability of the d-NBD for graphs of size $50$ to $1000$.

Considering a benchmark graph with $100$ vertices we subsequently calculate the distances to an ensemble of graphs all generated from the same model. That is, we generate a reference graph $G^{*}_{100}~{}\rm{Model}$ and for each $n\in[50,100,\dots,1000]$ we sample $200$ graphs $(G_{n}^{(i)})_{i=1}^{200}$ from that model and of that size before computing $d(G^{*}_{1}00,G_{n}^{(i)})$ for each graph. Due to randomness in each model the distances are small but non-zero (Figure 3(a)) however this distance remains constant as $n$ grows. The distances do however become noisier for small $n$ as the graph properties for smaller graphs are more susceptible to stochastic fluctuations.

In contrast, taking the largest $k$ eigenvalues of $B^{\prime}$ does not provide a stable distance as $n$ varies (Figure 3(b)). For both ER and WS graphs the distance from the reference graph is monotonically increasing with $n$. The distance for the regular graph is however stable. This can be attributed to the fact that the largest degree, which correlates with the largest eigenvalue, is fixed. The d-NBD therefore is able to compare graphs generated from the same model, even if the degree distribution changes.

A further test we can consider on synthetic graphs is whether the d-NBD has a continuous dependence on the model parameters. For this purpose we consider the Watt-Strogatz model since we can consider dependency is both the nearest neighbour connection parameter $k$ and the edge rewiring parameter $p$. We can interpolate between regularised ring structure ($p=0,k\geq 2$) to completely randomised connections. Similarly we can interpolate between sparsely connected graphs ($k=2$) to the complete graph ($k=n-1$).

By considering a $10,000$-dimensional embedding of these graphs, via (8), we can see that the space of Watts-Strogatz networks lies on a manifold (Figure 4.)

The embedding has been reduced to three dimensions using principle component analysis and each point is taken as the average position of $100$ graphs with $200$ vertices. This observation gives further credit to the distance as changes in the underlying graph model give predictable changes in the underlying embedding. It must be stressed that the changes in the embedding are non-trivial functions of the parameter changes. This is clear by the uneven spacing of points in Figure 4 which were generated by an even spacing of parameter space.

These results on synthetic graphs have shown that the d-NBD captures, but does not necessarily inform about, the generating process for these graphs and is able to detect subtle changes in the underlying model parameters. By their very nature, these synthetic graphs are generated from simple mechanisms. Is the d-NBD is able to capture more complex generative mechanisms?

4.2 Empirical Graphs

Graph distances are useful in practice to be able to determine either the function of a graph or the (closely related) mechanisms for graph formation. To illustrate how the d-NBD can distinguish between these mechanisms we consider a number of graph collections (detailed in Table 1) drawn from different fields.

For example, we consider a sample of $26$ graphs from the Facebook100 dataset which describes the social connections within universities in the US. These networks are human-made (although technologically assisted) and capture a number of human behaviours such as assortativity and clustering [41]. In contrast the metabolic networks dataset captures the reactions between metabolites within a number of different organisms. These graphs are generated based on chemical and biological interaction, although arguably there has been human interference in the design and curation of such graphs. A useful distance measure should therefore be able to distinguish between these different graph types by having small intra-type distances and relatively large inter-type distances.

Figure 5 shows a two-dimensional projection of the spectral embedding (8) for each graph, each of which is coloured by the graph collection it came from.

Visually it is evident that the embedding performs well as each data collection forms its own isolated cluster. This is even the case for network collections, such as the Facebook100 graphs where the largest graph is approximately six times as large as the smallest graph in the collection.

To assess the equality of the clustering more formally we conduct a $k$-nearest neighbour classification [2] with the number of nearest neighbours $k=10$. We compare the clustering of the d-NBD against two benchmarks. First, the non-backtracking spectral distance (NBD) [40] which is simply the Euclidean distance between the largest eigenvalues. Since the smallest graph we consider has $82$ vertices we can only consider up to $164$ eigenvalues. We also compare against the spectrum of the graph Laplacian, again taking up to $82$ of the largest eigenvalues. Similarly to the non-backtracking matrix, the graph Laplacian (in particular the spectrum) has been previously used to characterise a graph [34].

We assess the quality of the clustering using both the training and testing accuracy, precision, and recall. The results are presented in Table 2.

To prevent overfitting we average the results over a $10$-fold cross validation. The d-NBD achieves perfect performance for both the training and test data. Arguably this may not be a difficult classification task and so these results should be be interpreted as a flawless measure. However even in this ‘simple’ task both the NBD and Laplacian methods struggle to classify a number of graphs correctly. This is likely due to the methods not being able to utilise the full spectrum when comparing graphs of dissimilar size.

We can alternatively use more sophisticated clustering algorithms such as a random forest classifier. These algorithms provide a slight improvement in accuracy for both NBD and the Laplacian method however neither achieve the accuracy of the d-NBD clustering.

5 Conclusions

In this article we have presented a method of graph comparison which can be used to compare graphs of different sizes. Through application to both synthetic and real graphs we have shown that the d-NBD can successfully identify the differences between graphs generated by different mechanisms. This distance performs better than distances defined on the truncated spectra of both the Laplacian and non-backtracking matrix, however this comes at a cost of having to compute the entire graph spectrum. This is not an issue for graphs of moderate size, but can become computationally taxing for graphs formed of millions of vertices.

One possible remedy to this may be to consider the pseudospectra of such matrices. The pseudospectra gives some information on the location of eigenvalues and can be calculated to a prescribed precision within the complex plane [42]. The pseudospectrum also has the added property that it gives a measure of the stability of the eigenvalues to small perturbations in the adjacency matrix. It should therefore give an indication to the effect of adding or removing an edge from the graph.

Another unexplored avenue is the use of graph distances to infer the graph generating mechanism rather than using the distance to distinguish between potentially different mechanisms. Approaches in this area are predominantly statistical, such as the work of Chen et al. [9] who use approximate Bayesian computation to select a likely graph generation mechanism. These approaches are heavily computational especially if the number of model choices is large. The d-NBD could be used as a precursor to such analysis by filtering possible models based on the graph distance to representative samples from these models.

The d-NBD has yet to be tested on a larger set of data which includes a wider range of graph categories (social, biological, etc.) but also with a diversity of graph collections within each category (Twitter, Facebook, and Instagram social graphs for instance). We anticipate that the d-NBD will continue to outperform other truncated spectral methods on the basis that more information is captured in the complete spectrum.

There also lies the more fundamental question of why spectral methods work in clustering graphs (both for the Laplacian and non-backtracking operators). The non-backtracking operator is closely related to the deformed Laplacian [19], also know as the Bethe Hessian [37]. This suggests there may be a physical interpretation to elucidate the underpinnings of such operators. What does a family of Laplacian-like operators tell us about the graph structure, and how many such operators are sufficient to characterise a graph? Connecting these two areas could lead to possible generalisations which could provide significant improvements in graph comparison.

Acknowledgments

The authors would like to thank Renaud Lambiotte for useful discussion and both Renaud and Ebrahim Patel for a careful reading of early versions of this article.