Year: 2021

Abstract:
Network-valued data are encountered in a wide range of applications, and pose challenges in learning due to their complex structure and absence of vertex correspondence. Typical examples of such problems include classification or grouping of protein structures and social networks. Various methods, ranging from graph kernels to graph neural networks, have been proposed that achieve some success in graph classification problems. However, most methods have limited theoretical justification, and their applicability beyond classification remains unexplored. In this work, we propose methods for clustering multiple graphs, without vertex correspondence, that are inspired by the recent literature on estimating graphons— symmetric functions corresponding to infinite vertex limit of graphs. We propose a novel graph distance based on sorting-and-smoothing graphon estimators. Using the proposed graph distance, we present two clustering algorithms and show that they achieve state-of-the-art results. We prove the statistical consistency of both algorithms under Lipschitz assumptions on the graph degrees. We further study the applicability of the proposed distance for graph two-sample testing problems.

heading: 1 INTRODUCTION
text: Machine learning on graphs has evolved considerably over the past two decades. The traditional view towards network analysis is limited to modelling interactions among entities of interest, for instance social networks or world wide web, and learning algorithms based on graph theory have been commonly used to solve these problems (Von Luxburg, 2007; Yan et al., 2006). However, recent applications in bioinformatics and other disciplines require a different perspective, where the networks are the quantities of interest. For instance, it is of practical interest to classify protein structures as enzyme or non-enzyme (Dobson & Doig, 2003) or detect topological changes in brain networks caused by Alzheimer’s disease (Stam et al., 2007). In this paper, learning from network-valued data refers to clustering where each network is treated as an entity, as opposed to the traditional network analysis problems that involve a single network of interactions (Newman, 2003).
Machine learning on network-valued data has been an active area of research in recent years, although most works focus on the network classification problem. The generic approach is to convert the network-valued data into a standard representation. Graph neural networks are commonly used for network embedding, that is, finding Euclidean representations of each network that can be further used in standard machine learning models (Narayanan et al., 2017; Xu et al., 2019). In contrast, graph kernels capture similarities between pairs of networks that can be used in kernel based learning algorithms (Shervashidze et al., 2011; Kondor & Pan, 2016; Togninalli et al., 2019). In particular, the graph neural tangent kernel defines a graph kernel that corresponds to infinitely wide graph neural networks, and typically outperforms neural networks in classification tasks (Du et al., 2019). A more classical equivalent for graph kernels is to define metrics that characterise the distances between pairs of graphs (Bunke & Shearer, 1998), but there has been limited research on designing efficient graph distances and developing algorithms for clustering network-valued data.
The motivation for this paper stems from two shortcomings in the literature on network-valued data analysis: first, the efficacy of existing kernels or embeddings have not been studied beyond network classification, and second is the lack of theoretical analysis of these methods, particularly in the small sample setting. Generalisation error bounds for graph kernel based learning exist (Du et al., 2019), but these bounds, based on learning theory, are meaningful only when many networks
are available. However, in many applications, one needs to learn from a small population of large networks and, in such cases, an informative statistical analysis should consider the small sample, large graph regime. To address this issue, we take inspiration from the recent statistics literature on graph two-sample testing—given two (populations of) large graphs, the goal is to decide if they are from same statistical model or not. Although most theoretical studies in graph two-sample testing focus on graph with vertex correspondence (Tang et al., 2017a; Ghoshdastidar & von Luxburg, 2018), some works address the problem of testing graphs on different vertex sets either by defining distances between graphs (Tang et al., 2017b; Agterberg et al., 2020) or by representing networks in terms of pre-specified network statistics (Ghoshdastidar et al., 2017). The use of network statistics for clustering network-valued data is studied in Mukherjee et al. (2017). Another fundamental approach for dealing with graphs of different sizes is graph matching, where the objective is to determine the vertex correspondence. Graph matching is often solved by formulating it as an optimization problem (Zaslavskiy et al., 2008; Guo et al., 2019) or defining graph edit distance between the graphs (Riesen & Bunke, 2009; Gao et al., 2010). Although there is extensive research on graph matching, the efficacy of these methods in clustering network-valued data remains unexplored.
Contribution and organisation. In this work, we follow the approach of defining meaningful graph distances based on statistical models, and use the proposed graph distance in the context of learning from networks without vertex correspondence. In particular, we propose graph distances based on graphons. Graphons are symmetric bivariate functions that represent the limiting structure for a sequence of graphs with increasing number of nodes (Lovász & Szegedy, 2006), but can be also viewed as a nonparametric statistical model for exchangeable random graphs (Diaconis & Janson, 2007; Bickel & Chen, 2009). The latter perspective is useful for the purpose of machine learning since it allows us to view the multiple graphs as random samples drawn from one or more graphon models. This perspective forms the basis of our contributions, which are listed below:
1) In Section 2, we propose a distance between two networks, that do not have vertex correspondence and could have different number of vertices. We view the networks as random samples from (unknown) graphons, and propose a graph distance that estimates theL2-distance between the graphons. The distance is inspired by the sorting-and-smoothing graphon estimator (Chan & Airoldi, 2014).
2) In Section 3, we present two algorithms for clustering network-valued data based on the proposed graph distance: a distance-based spectral clustering algorithm, and a similarity based semi-definite programming (SDP) approach. We derive performance guarantees for both algorithms under the assumption that the networks are sampled from graphons satisfying certain smoothness conditions.
3) We empirically compare the performance of our algorithms with other clustering strategies based on graph kernels, graph matching, network statistics etc. and show that, on both simulated and real data, our graph distance-based spectral clustering algorithm outperforms others while the SDP approach also shows reasonable performance, and they also scale to large networks (Section 3.3).
4) Inspired by the success of the proposed graph distance in clustering, we use the distance for graph two-sample testing. In Section 4, we show that the proposed two-sample test is statistically consistent for large graphs, and also demonstrate the efficacy of the test through numerical simulation.
We provide further discussion in Section 5 and present the proofs of theoretical results in Appendix.

heading: 2 GRAPH DISTANCE BASED ON GRAPHONS
text: Clustering or testing of multiple networks requires a notion of distance between the networks. In this section, we present a transformation that converts graphs of different sizes into a fixed size representation, and subsequently, propose a graph distance inspired by the theory of graphons. We first provide some background on graphons and graphon estimation. Graphon has been studied in the literature from two perspectives: as limiting structure for infinite sequence of growing graphs (Lovász & Szegedy, 2006), or as exchangeable random graph model. In this paper, we follow the latter perspective. A random graph is said to be exchangeable if its distribution is invariant under permutation of nodes. Diaconis & Janson (2007) showed that any statistical model that generates exchangeable random graphs can be characterised by graphons, as introduced by Lovász & Szegedy (2006). Formally, a graphon is a symmetric measurable continuous function w : [0, 1]2 → [0, 1] where w(x, y) can be interpreted as the link probability between two nodes of the graph that are assigned values x and y, respectively. This interpretation propounds the following two stage sampling
procedure for graphons. To sample a random graph G with n nodes from a graphon w, in the first stage, one samples n variables U1, . . . , Un uniformly from [0, 1] and constructs a latent mapping between the sampled points and the node labels. In the second stage, edges between any two nodes i, j are randomly added based on the link probability w(Ui, Uj). Mathematically, if we abuse notation to denote the adjacency matrix by G ∈ {0, 1}n×n, we have
U1, . . . , Un iid∼ Uniform[0, 1] and Gij |Ui, Uj ∼ Bernoulli(w(Ui, Uj)) for all i < j.
We consider problems involving multiple networks sampled independently from the same (or different) graphons. We make the following smoothness assumptions on the graphons.
Assumption 1 (Lipschitz continuous) A graphon w is Lipschitz continuous with constant L if |w(u, v)− w(u′, v′)| ≤ L √ (u− u′)2 + (v − v′)2 for every u, v, u′, v′ ∈ [0, 1].
Assumption 2 (Two-sided Lipschitz degree) A graphon w has two-sided Lipschitz degree with constants λ1, λ2 > 0 if its expected degree function g, defined by g(u) = ∫ 1 0 w(u, v)dv, satisfies
λ2|u− u′| ≤ |g(u)− g(u′)| ≤ λ1|u− u′| for every u, u′ ∈ [0, 1].
One of the challenges in graphon estimation is due to the issue of non-identifiability, that is, different graphon functions w can generate the same random graph model. In particular, two graphons w and w′ generate the same random graph model if they are weakly isomorphic—there exist two measure preserving transformations φ, φ′ : [0, 1] → [0, 1] such that w(φ(u), φ(v)) = w′(φ′(u), φ′(v)). Moreover, the converse also holds meaning that such transformations are known to be the only source of non-identifiability (Diaconis & Janson, 2007). This weak isomorphism induces equivalence classes on the space of graphons. Since our goal is only to cluster graphs belonging to random graph models, we simply make the following assumption on our graphons.
Assumption 3 (Equivalence classes) Any reference toK graphons, w1, . . . , wK , assumes that, for every i, j, either wi = wj or wi and wj belong to different equivalence classes. Furthermore, without loss of generality, we assume that every graphon wi is represented such that the corresponding degree function gi is non-decreasing.
Remark on the necessity of Assumptions 1–3. Assumption 1 is standard in graphon estimation literature (Klopp et al., 2017) since it avoids graphons corresponding to inhomogeneous random graph models. It is known that two graphs from widely separated inhomogeneous models (in L2distance) are statistically indistinguishable (Ghoshdastidar et al., 2020), and hence, it is essential to ignore such models to derive meaningful guarantees. Assumption 2 ensures that, under a measurepreserving transformation, the graphon has strictly increasing degree function, which is a canonical representation of an equivalence class of graphons (Bickel & Chen, 2009). Assumption 3 is needed since graphons can only be estimated up to measure-preserving transformation. As noted above, it is inconsequential for all practical purposes but simplifies the theoretical exposition.
Graph transformation. In order to deal with multiple graphs and measure distances among pairs of graphs, we require a transformation that maps all graphs into a common metric space—the space of all n0 × n0 symmetric matrices for some integer n0. While the graphon estimation literature provides several consistent estimators (Klopp et al., 2017; Zhang et al., 2017), only the histogram based sorting-and-smoothing graphon estimator of Chan & Airoldi (2014) can be adapted to meet the above requirement. We use the following graph transformation, inspired by Chan & Airoldi (2014). The adjacency matrix G of size n× n is first reordered based on a non-unique permutation σ, such that the empirical degree based on this permutation is monotonically increasing. The degree sorted adjacency matrix is denoted by Gσ . It is then transformed to a ‘histogram’ A ∈ Rn0×n0 as
Aij = 1
h2 h∑ i1=1 h∑ j1=1 Gσ(i−1)h+i1,(j−1)h+j1 , where h = ⌊ n n0 ⌋ and b·c is the floor function. (1)
Proposed graph distance. Given two directed or undirected graphs G1 and G2 with n1 and n2 nodes, respectively, we apply the transformation (1) to both the graphs with n0 ≤ min{n1, n2}. We propose to use the graph distance
d(G1, G2) = 1
n0 ‖A1 −A2‖F , (2)
where A1 and A2 denote the transformed matrices and ‖ · ‖F denotes the matrix Frobenius norm. Proposition 1 shows that, if G1 and G2 are sampled from two graphons, then the graph distance (2) consistently estimates the L2-distance between the two graphons, which is defined as ‖w1 − w2‖2L2 = ∫ 1 0 ∫ 1 0 (w1(x, y)− w2(x, y))2 dx dy.
Proposition 1 (Graph distance is consistent) Let w1 and w2 satisfy Assumptions 1–3. Let G1 and G2 be random graphs with at least n nodes sampled from the graphons w1 and w2, respectively. If n→∞ and n0 is chosen such that n
2 0 logn n → 0, then with high probability (w.h.p.),∣∣‖w1 − w2 ‖L2 − d(G1, G2)∣∣ = O ( 1n0) . (3)
Proof sketch. We define a novel technique for approximating the graphon. The proof in Appendix A.1 first establishes that the approximation error is bounded using Assumption 1. Consequently, a relation between approximated graphons and transformed graphs is derived using lemmas from Chan & Airoldi (2014). Proposition 1 is subsequently proved using the above two results.
Remark on Proposition 1 for sparse graphs. The defined sampling procedure for graphon generates dense graphs, deviating from the real world sparse networks. To adapt it to sparse graphs, one may modify the sampling procedure to Gij |Ui, Uj ∼ Bernoulli(ρw(Ui, Uj)) where ρ depends only on n (Olhede & Wolfe, 2014). Under this process, the consistency result in Proposition 1 remains unchanged for ρ = Ω( √ log n/n) (proof in Appendix A.2). This bound cannot be improved to the expected real world measure where ρ = Ω(log n/n), because of the degree sorting step in (1). Nevertheless, our analysis allows for relatively sparser graphs with strong consistency result.
Notation. For ease of exposition, Proposition 1 as well as main results are stated asymptotically using the standard O(·) and Ω(·) notations, which subsume absolute and Lipschitz constants. We use “with high probability” (w.h.p.) to state that the probability of an event converges to 1 as n→∞.

heading: 3 GRAPH CLUSTERING
text: We now present the first application of the proposed graph distance (2) in the context of clustering network-valued data. We are particularly interested in the setting where one needs to cluster a small population of large graphs, that is, minimum graph size n grows faster than the sample size m. This scenario is relevant in practice as bioinformatics or neuroscience application often deals with very few graphs (see real datasets in Section 3.3). Theoretically, this perspective complements guarantees for (graph) kernels that are applicable only in supervised setting and large sample regime, m→∞. In contrast, our guarantees are more conclusive for bounded m and large graph size, n→∞.

heading: 3.1 DISTANCE BASED SPECTRAL CLUSTERING (DSC)
text: Given m graphs with adjacency matrices G1, ..., Gm, we propose a distance based clustering algorithm where we apply spectral clustering to an estimated distance matrix. The distance matrix D̂ ∈ Rm×m is computed on all pairs of graphs using the defined estimator function (2), that is D̂ij = d(Gi, Gj). Unlike the standard Laplacian based spectral clustering, which is applicable for adjacency or similarity matrices, we use the method suggested by Mukherjee et al. (2017) that computes the K leading eigenvectors of D̂ (corresponding to the K smallest eigenvalues in magnitude) and applies k-means clustering to the rows of the eigenvector matrix resulting in K number of clusters. We refer to this distance based clustering algorithm as DSC, described in Algorithm 1 of Appendix. To derive the statistical consistency of DSC, we consider the problem of clustering m random graphs of potentially different sizes, each sampled from one ofK graphons. We establish the consistency in Theorem 1 by proving that the number of misclustered graphs→ 0 asymptotically.
Theorem 1 (Consistency of DSC) Consider K graphons satisfying Assumptions 1–3, and m random graphs G1, . . . , Gm, each sampled from one of the K graphons (assume there is at least one graph from each graphon). Define the distance matrix D ∈ Rm×m such that Dij = ‖wi − wj‖L2 where wi and wj are the graphons from which Gi and Gj are generated. Let n be the size of the smallest graph, and γ be the K-th smallest eigenvalue value of D in magnitude. As n→∞, if n0 is chosen such that m
2n20 logn n → 0, then DSC misclusters at most O
( m3
γ2n20
) graphs w.h.p.
Proof sketch. The proof, given in Appendix A.3, uses the Davis-Kahan spectral perturbation theorem to bound the error in terms of ‖D̂ −D‖F , which is further bounded using Proposition 1. While the number of misclustered graphs seem to depend on m3, we note that there is an inverse dependence on γ2 which has dependence onm (see Corollary 1 that illustrates it for a specific case). Moreover, our focus is on the setting where m = O(1) and n, n0 → ∞, in which case, the error asymptotically vanishes. It is natural to wonder whether the dependence on m and n0 is tight in the above bounds. Currently, we do not know the optimal rates, but deriving this would be difficult due to the strong dependency of entries in D̂ and slow rate of convergence of the graph distance in Proposition 1. The presence of γ in the above clustering error bound makes Theorem 1 less interpretable. Hence, we also consider the specific case of K = 2 (two graphons) in the following result, along with the assumption that equal number of graphs are generated from both graphons.
Corollary 1 Let w 6= w′ be two graphons satisfying Assumptions 1–3, and m is a bounded even number. Assume that equal number of graphs are generated from w and w′. For any n0 and large enough constant C such that ‖w − w′‖L2 ≥ C mn0 and m2n20 logn n → 0 as n → ∞, the number of graphs misclustered by Algorithm 1 goes to zero w.h.p.
The corollary implies that given the observed graphs are large enough, and if the choice of n0 is relatively small, n0 √ n/ log n, and the graphons are Ω( 1n0 ) apart in L2-distance, then the clustering is consistent. Intuitively, it can be understood that if we condense large graphs to a small representation (small n0), then the clusters can be identified only if the models are quite dissimilar.

heading: 3.2 SIMILARITY BASED SEMI-DEFINITE PROGRAMMING (SSDP)
text: We propose another algorithm for clustering m graphs based on similarity between pairs of graphs. The pairwise similarity matrix Ŝ ∈ Rm×m is computed by applying the Gaussian kernel on the distance between the graphs, that is Ŝij = exp ( −d(Gi,Gj)σiσj ) , where σ1, . . . , σn are parameters. For theoretical analysis, we assume σ1 = . . . = σn is fixed, but in experiments, the parameters are chosen adaptively. We use the following semi-definite program (SDP) (Yan et al., 2018; Perrot et al., 2020) to find membership of the observed graphs. Let X ∈ Rm×m be the normalised clustering matrix, that is, Xij = 1/|C| if i and j belong to the same cluster C, and 0 otherwise. Then, the SDP for estimating X is as follows:
max X
trace(ŜX) s.t.X ≥ 0, X 0, X1 = 1, trace(X) = K, (4)
where X ≥ 0, X 0 ensure that X is a non-negative, positive semi-definite matrix, and 1 denotes the vector of all ones. We denote the optimal X from the SDP as X̂ . Once we have X̂ , we apply standard spectral clustering on X̂ to obtain a clustering of the graphs. We refer to this algorithm as SSDP, described in Algorithm 2 of Appendix. We present strong consistency result for SSDP below.
Theorem 2 (Consistency of SSDP) ConsiderK graphons,w1, . . . , wK , satisfying Assumptions 1– 3, andm random graphs, each sampled from one of theK graphons. Let n be the size of the smallest graph. As n→∞, if n0 is chosen such that m
2n20 logn n → 0 and minl 6=l′ ‖wl −wl ′‖L2 = Ω ( m n0 ) , then
the number of graphs misclustered by SSDP is zero w.h.p.
Proof sketch. The proof in Appendix A.4 adapts Perrot et al. (2020, Proposition 1) to the present setting and combines it with Proposition 1 to derive the stated condition for zero error.
Theorem 2 is slightly stronger than Theorem 1, or Corollary 1, since SSDP achieve a zero clustering error for large enough graphs. This theoretical merit of SDP over spectral clustering is known in the statistics literature. Similar to Corollary 1, the choice of n0 is important such that it does not violate the minimum L2-distance condition in the theorem to ensure consistency.
Remark on the knowledge of K. Above discussions assume that the number of clusters K is known, which is not necessarily the case in practice. To tackle this issue, one can estimate K using Elbow method (Thorndike, 1953) or approach from Perrot et al. (2020),and then use it as input in our algorithms, DSC and SSDP. One can modify the SDP (4) and Theorem 2 to the case where K
is adaptively estimated. However, we found the corresponding algorithm, adapted from Perrot et al. (2020), to be empirically unstable in the present context. Hence, the knowledge of K is assumed in the following experiments, which also allows the efficacy of the proposed algorithms and graph distance to be evaluated without the error induced by incorrect estimation of K.

heading: 3.3 EXPERIMENTAL ANALYSIS
text: In this section, we evaluate the performance of our algorithms DSC and SSDP, both in terms of accuracy and computational efficacy. We measure the performance of the algorithms in terms of error rate, that is, the fraction of misclustered graphs by using the source of the graphs as labels. Since clustering provides labels up to permutation, we use the Hungarian method (Kuhn, 1955) to match the labels. The performance can also be measured in terms of Adjusted Rand Index (results in Appendix C.3). We use both simulated and real datasets for evaluation and obtain all our experimental results using Tesla K80 GPU instance with 12GB memory from Google Colab.
Simulated data. We generate graphs of varied sizes from four graphons, W1(u, v) = uv, W2(u, v) = exp { −max(u, v)0.75 } , W3(u, v) = exp { −0.5 ∗ (min(u, v) + u0.5 + v0.5) } and W4(u, v) = |u − v|. The simulated graphs are dense and the graph sizes are controlled to study how algorithms scale. Their corresponding L2 distances between pairs of graphons is shown later in Figure 3 and the heatmap of the graphons are visualised in Figure 4 in Appendix.
Real data. We analyse the performance of algorithms using datasets from two contrasting domains: molecule datasets from Bioinformatics and network datasets from Social Networks. The Bioinformatics networks are smaller whereas the latter has relatively larger graphs. We use Proteins (Borgwardt et al., 2005), KKI (Pan et al., 2016), OHSU (Pan et al., 2016) and Peking 1 (Pan et al., 2016) datasets from Bioinformatics, and Facebook Ct1 (Oettershagen et al., 2020), Github Stargazers (Rozemberczki et al., 2020), Deezer Ego Nets (Rozemberczki et al., 2020) and Reddit Binary (Yanardag & Vishwanathan, 2015) datasets from Social Networks. We sub-sample a few graphs from each dataset by setting a minimum number of nodes to validate the case of clustering small number of large graphs (small m, large n). The number and size of the graphs sampled from each dataset are listed in ‘#graphs’ and ‘#nodes’ columns of tables in Figure 1. We evaluate the clustering performance on all combinations of the datasets for three and four clusters in both the domains separately.
Choice of n0 and σi. As noted in our algorithms DSC and SSDP, n0 is an input parameter. Theorems 1 and 2 show that the performance of both DSC and SSDP depend on the choice of n0 = O (√ n/log n ) . In the experiments, we set n0 = √ n/ log n where n is the minimum number of nodes. In Appendix C.2, we use simulated data to show that the above choice of n0 is reasonable (if not the best) for both DSC and SSDP. Furthermore, the similarity matrix Ŝ in SSDP is computed using parameters σ1, . . . , σn. In the experiments, we set σi = d(Gi, G5nn) where G5nn is the fifth nearest neighbour of Gi. Hence, apart from knowledge of K, our algorithms are parameter-free.
Performance comparison with existing methods. We compare our algorithms with a range of approaches for measuring similarity or distance among multiple networks. Most methods discussed below provide a kernel or distance matrix to which we apply spectral clustering to obtain the clusters:
1) Network Clustering based on Log-Moments (NCLM) is an embedding based clustering strategy for different sized graphs (Mukherjee et al., 2017) that is based on network statistics called log moments. Log moments for a graph with adjacency matrix A and number of nodes n is obtained by (log (m1(A)) , log (m2(A)) , . . . , log (mJ(A))) wheremi(A) = trace(A/n)i and J is a parameter.
2) Wasserstein Weisfeiler-Lehman Graph Kernels (WWLGK) is a recent graph kernel that is based on the Wasserstein distance between the node feature vector distributions of two graphs proposed by Togninalli et al. (2019).
3) Graph Neural Tangent Kernel (GNTK) is another graph kernel that describes infinitely wide graph neural networks derived by Du et al. (2019). Both WWLGK and GNTK provide state-of-the-art performance in graph classification with GNTK outperforming most graph neural networks.
4) Network Clustering algorithm based on Maximum Mean Discrepancy (NCMMD) considers a graph metric (MMD) to cluster the graphs. MMD distance between random graphs is proposed as an efficient test statistic for random dot product graphs (Agterberg et al., 2020). We compute
MMD between the graphs that are represented by latent finite dimensional embedding called spectral adjacency embedding with the dimension r as a parameter.
5) In Network Clustering algorithm based on Graph Matching Metric (NCGMM), we match two graphs of different sizes by appending null nodes to the small graph as described in Guo et al. (2019) and compute Frobenius norm between the matched graphs as their distance. Although both the considered graph metrics (MMD and graph matching) are for different purposes, we evaluate their efficacy in the context of clustering.
The different parameters to tune in the algorithms include n0 and σi in our algorithms DSC and SSDP, J in NCLM, number of iterations (#itr) to perform in WWLGK, number of layers (#layer) in graph neural networks for GNTK, r in NCMMD and none in NCGMM. We fix n0 in our algorithms using the theoretical bound and σi is set adaptively as discussed, whereas we tune the parameters for other algorithms by grid search over a set of values.
Evaluation on simulated data. We sample 10 graphs of varied sizes between 50 and 100 nodes from each of the four graphons in Figure 4, and perform the experiments by considering all combinations of three and four clusters of the chosen graphons. Based on the theoretical bound, n0 is fixed to 5 as n = 50. We report the performance for J = 8, r = 3, #itr = 1 and #layer = 2 as these produce the best results. The first row of Figure 1 shows the average performance of the algorithms computed over 5 independent runs. We observe that our algorithm DSC outperforms all the other algorithms, achieving nearly zero error in all cases, and SSDP also performs competitively by standing second or third best. The graph kernels, WWLGK and GNTK, and the graph metric based method NCGMM typically do not perform well. NCMMD either performs very well or quite poorly. We sample small graphs since otherwise GNTK cannot run due to memory requirement for dense large graphs and NCGMM has high computation time. Appendix C.4 includes evaluation of the algorithms except GNTK and NCGMM on larger graphs, where we observe similar behaviour.
Evaluation on real data. We consider all combinations of three and four clusters of both Bioinformatics and Social Networks separately. The second and third rows of Figure 1 show the performance with n0 = 30, J = 8, r = 3, #itr = 1 and #layer = 2, and the upper limit of 7200 seconds (2 hours) as running time of algorithms. We observe DSC outperforms other algorithms by a large margin in majority of the combinations, while in the other combinations like {Proteins,KKI,Peking 1}, DSC performs well with a very small margin to the best performing one. Although NCLM and GNTK compare favorably in Social Networks datasets, they typically have high error rate in Bioinformatics or simulated datasets, suggesting that they could be well suited for large networks, whereas DSC is more versatile and suitable for all networks. SSDP performs moderately on real data, but it achieves the smallest error in some cases, implying that it is suited for certain types of networks.
Computation time comparison. Figure 2 shows the time (measured in seconds) taken by each algorithm for four clusters case, plotted in log scale. Similar behavior is observed in three clusters case
also and the result can be found in Appendix C.5. Our algorithms, DSC and SSDP, perform competitively with respect to time as well. In addition, it scales effectively for large graphs unlike other
algorithms. It is worth noting that although NCLM takes lesser time than DSC and SSDP for small graphs, it takes longer for large social networks datasets, thus favoring our methods in terms of both accuracy and scalability. Graph matching based algorithm, NCGMM, has severe scalability issue demonstrating the inapplicability of such methods to learning problems. We also evaluate the scalability of all the algorithms by measuring the time taken for clustering different sets of varied sized graphs from graphons W1,W2,W3 and W4 and detailed discussion is provided in Appendix C.6. The experimental results also illustrate the high scalability of DSC and SSDP compared to the other algorithms.

heading: 4 GRAPH TWO-SAMPLE TESTING
text: Inspired by the remarkable performance of the proposed graph distance (2) in clustering, we analyse the applicability of the distance for graph two-sample testing. Two-sample testing is usually studied in the large sample case m → ∞, and several nonparametric tests are known that could also be applied to graphs. However, in the context of graphs, it is relevant to study the small sample setting, particularly m = 2, that is, the problem of deciding if two large graphs are statistically identical or not (Ghoshdastidar et al., 2020; Agterberg et al., 2020).
We consider the following formulation of the graph two-sample problem, stated under the assumption that the graphs are sampled from graphons. Given two random graphs, G1 sampled from some model (here, graphon w1), and G2 from another model w2, the goal is to determine which of the following hypothesis is true: H0 : {w1 = w2} or Ha : { w1 6= w2 : ‖w1 − w2‖L2 ≥ φ } for some φ > 0. Existing works consider alternative random graph model, such as inhomogeneous ErdősRényi models or random dot product graph models, which are more restrictive. The condition φ > 0 is necessary if one only has access to finitely many independent samples (Ghoshdastidar et al., 2020). A two-sample test T is a binary function of the given samples such that T = 1 denotes that the test rejects the null hypothesis H0 and T = 0 implies that the test rejects the alternate hypothesis Ha. The goodness of a two-sample test is measured in terms of the Type-I and Type-II errors, which denote the probabilities of incorrectly rejecting the null and alternate hypotheses, respectively.
The goal of this section is to show that one can construct a test T that has arbitrarily small Type-I and Type-II errors. For this purpose, we consider the test
T : I {d(G1, G2) ≥ ξ} (5) for some ξ > 0, where I{·} is the indicator function and d(G1, G2) is the proposed graph distance for some choice of integer n0. We state the following theoretical guarantee for the test T .
Theorem 3 Assume that the graphons w1, w2 satisfy Assumptions 1–3, and let the graphsG1 ∼ w1 and G2 ∼ w2 have at least n nodes. As n → ∞, there is a choice of ξ such that the Type-I and Type-II errors of the test T in (5) go to 0 if n
2 0 logn n → 0 and φ ≥ C n0
, where the constant C depends only on the Lipschitz constants.
Theorem 3 shows that the test T in (5) can distinguish between any pair of graphons that have separation ‖w1 − w2‖L2 = Ω(1/n0) with arbitrarily small error, if the graphs are large enough. Empirical analysis. We empirically validate the consistency result in Theorem 3 by computing power of the proposed two-sample test T , which measures the probability of rejecting the null hypothesis H0. Intuitively, power of the test for graphs sampled from same graphons should be small (close to a pre-specified significance level) since H0 must not be rejected, whereas, it should be close to 1 for graphs sampled from different graphons. As known in the testing literature, theoretical threshold, ξ in (5), is typically conservative in practice and the rejection/acceptance is decided based on p-values, computed using bootstrap samples. To this end, we follow the bootstrapping strategy in Ghoshdastidar & von Luxburg (2018, Boot-ASE algorithm). In addition, we compare the proposed test T by replacing d(G1, G2) in 5 with two other statistics, log moments from Mukherjee
et al. (2017) and MMD, an efficient test statistics for random dot product graphs (Agterberg et al., 2020). We perform the experiment by sampling two graphsG1 ∼ w1 andG2 ∼ w2 of size n and 2n, respectively, where w1 and w2 are chosen from the graphons W1,W2,W3,W4. Suppose n = 100 and thus n0 = 10 from the theoretical bound for evaluating T . The power of T is computed for the significance level 0.05, averaged over 500 trials of bootstrapping 100 samples generated from all pairs of graphons. The plots in Figure 3 show the average power of T with our proposed distance, log moments and MMD as d(G1, G2), respectively. From the results, it is clear that T using our proposed distance can distinguish between graphons that are quite close too (smallest L2 distance for W2 and W3), whereas, other test statistics are weak as log moments statistic accepts H0 even when it is wrong (see W1 and W4) and MMD based test rejects it strongly almost always (see diagonal). Similar result for small and large n, and evaluation on real datasets are presented in Appendix C.

heading: 5 CONCLUSION
text: There has been significant progress in learning on complex data, including network-valued data. However, much of the theoretical and algorithmic development have been in large sample problems, where one has access to m → ∞ independent samples. Practical applications of network-valued data analysis often leads to small sample, large graph problems—a setting where the machine learning literature is quite limited. Inspired by graph limits and high-dimensional statistics, this paper proposes a simple graph distance (2) based on non-parametric graph models (graphons).
Sections 3–4 demonstrate that the proposed graph distance leads to provable and practically effective algorithms for clustering (DSC and SSDP) as well as two-sample testing (5). Extensive empirical studies on simulated and real data show that the clustering based on the graph distance (2) outperforms methods based on more complex graph similarities or metrics, both in terms of accuracy and scalability. Figures 1–2 show that DSC achieves best performance for both small dense graphs (simulated graphons) as well as large sparse graphs (social networks). On the other hand, popular machine learning approaches—graph kernels or graph matching—can be computationally expensive in large graphs and their performance may not improve as n→∞, see WWLGK in Figure 7. Statistical approaches, such as the proposed clustering algorithms and two-sample test, show better performance on large graphs (Figures 1, 3 and Appendix C.4). Theorems 1–3 theoretically support this observation by showing consistency of the clustering and testing methods in the limit of n→∞. The theoretical results, however, hinge on Assumptions 1–3. We remark that such smoothness and equivalence assumptions could be necessary for meaningful non-parametric approaches, which is also supported by the graph testing and graphon estimation literature. Further insights about the necessity of smoothness assumptions would aid in theoretical and algorithmic development.
The poor performance of graph kernels and graph matching in clustering and small sample problems calls for further studies on these methods, which have shown success in network classification. Fundamental research, combining graphon based approaches and kernels, could lead to improved techniques. Instead of sorting-and-smoothing graphon estimator, other histogram based methods that reduce graphons to a fixed size by identifying block structures and averaging (Olhede & Wolfe, 2014) can be explored. However, the distance defined on such reduced graphons require alignment of the blocks between the graphons. While such methods do not require Assumption 2 and allow sparse graphs since there is no sorting, the computationally inefficient alignment step pose a practical challenge. Algorithmic modifications, such as estimation of K, would be also useful in practice.

heading: 6 ETHICS STATEMENT
text: The proposed clustering algorithms for network-valued data can naturally be evaluated on realworld data, but no social interpretation about the results will be drawn. Moreover, development of fair clustering algorithms is not the focus of this project.

heading: 7 REPRODUCIBILITY STATEMENT
text: The assumptions for the theory are stated clearly in Assumptions 1–3 and all the theoretical results, Proposition 1, Theorems 1-3, Corollary 1, are proved in detail in Appendix A. The implementation of the considered algorithms are provided in graph clustering.zip with the filename as the corresponding algorithm. Datasets used in the experiments are public and the links are provided in download datasets function of utils.py. The experimental results can be reproduced by following graph clustering.ipynb.

heading: A PROOFS OF THEORETICAL RESULTS
text: We discuss the proofs of Proposition 1 and Theorems 1–3 with supporting lemmas in this section.

heading: A.1 PROPOSITION 1
text: The distance function defined in (2) estimates the L2-distance between graphons that are continuous. To prove this, we introduce a method to discretize the continuous graphons in the following so that it is comparable with the transformed graphs described in the graph distance estimator (2).
Graphon discretization. We discretize the graphon by applying piece-wise constant function approximation that is inspired from Chan & Airoldi (2014), similar to the graph transformation. More precisely, any continuous graphon w is discretized to a matrix W of size n0 × n0 with
Wij = 1
1/n20
∫ 1 n0
0
∫ 1 n0
0
w(x+ i
n0 , y +
j
n0 ) dxdy (6)
We make Assumptions 1–3 to derive the concentration bound stated in Proposition 1. The proof structure is as follows:
1. We bound the point-wise deviation of graphon with its discretized correspondence using Lipschitzness assumption (Lemma 1).
2. We derive the error bound between the L2-distance and Frobenius norm of the discretized graphons using Lemma 1 (Lemma 2).
3. We establish a relation between Frobenius norm of the histograms of graphons and graphs (Lemma 3).
4. Finally, we prove Proposition 1 by combining Lemmas 2 and 3.
Lemma 1 (Lipschitz condition) For any graphon w and corresponding discretization W , define a piecewise constant function w(x, y) = Wij where x ∈ [ in0 , i n0 + 1n0 ] and y ∈ [ j n0 , jn0 + 1 n0
]. Using the Lipschitz continuous assumption, we have∣∣∣w(x, y)− w(x, y)∣∣∣ ≤ 2√2L
n0 , (7)
where L is the Lipschitz constant in Assumption 1.
Proof of Lemma 1 (Lipschitz condition). We use Assumption 1 on Lipschitzness to prove this lemma. The following holds for a graphon w with Lipschitz constant L,∣∣∣∣w(x+ in0 , y + jn0 )− w(x, y) ∣∣∣∣ ≤ L √ i2 n20 + j2 n20 ≤ L √ 2 n0 (8)
We prove Lemma 1 using (8) and the definition of w(x, y) = Wij where x ∈ [ in0 , i n0 + 1n0 ] and y ∈ [ jn0 , j n0
+ 1n0 ],∣∣∣w(x, y)− w(x, y)∣∣∣ = ∣∣∣∣w(x, y)− w(x, y)± w(x+ in0 , y + jn0 ) ∣∣∣∣
≤ ∣∣∣∣w(x, y)− w(x+ in0 , y + jn0 ) ∣∣∣∣+ ∣∣∣∣w(x+ in0 , y + jn0 )− w(x, y) ∣∣∣∣ = 2√2Ln0
Lemma 2 (Error bound of discretization) For two graphons w1 and w2, the error bound between the L2-distance and the Frobenius norm of the corresponding discretized graphons W1 and W2 satisfies ∣∣∣∣‖w1 − w2 ‖L2 − 1n0 ‖W1 −W2‖F ∣∣∣∣ ≤ 4√2Ln0 . (9)
Proof of Lemma 2 (Error bound of the approximation). Lemma 1 is used to prove this lemma. Let L1 and L2 be the Lipschitz constants of w1 and w2.
‖w1 − w2 ‖2L2 = ∫ 1
0 ∫ 1 0 (w1(x, y)− w2(x, y))2 dx dy
±w1(x, y) and ±w2(x, y) within the square, expand and apply Lipschitzness condition from Lemma 1
≤ 8L 2 1
n20 + 8L22 n20 + 16L1L2 n20 + n0∑ k=0 n0∑ l=0 1 n20 ( (W1)kl − (W2)kl )2 + 4 √ 2(L1 + L2)
n0
n0∑ k=0 n0∑ l=0 1 n20 ∣∣∣(W1)kl − (W2)kl∣∣∣ ≤ (2√2(L1 + L2)
n0
)2 + 1
n20 ‖W1 −W2‖2F +
4 √
2(L1 + L2)
n30 ‖W1 −W2‖1
(a) ≤ (2√2(L1 + L2)
n0 +
1
n0 ‖W1 −W2‖F
)2 (10)
(a) : ‖x‖1 ≤ √ n ‖x‖F
Similarly, by ±w1(x, y) and ±w2(x, y) to 1
n20 ‖W1 −W2‖2F and applying Lipschitzness condition
from Lemma 1 we get,
1
n20 ‖W1 −W2‖2F ≤ (2√2(L1 + L2) n0 + ‖w1 − w2‖L2 )2
(11)
Combining (10) and (11), we prove∣∣∣∣‖w1 − w2 ‖L2 − 1n0 ‖W1 −W2‖F ∣∣∣∣ ≤ 2√2(L1 + L2)n0 (b)≤ 4 √ 2L n0 , (b): L = max{L1, L2}
We derive the following relation between histogram of graphs and graphons by adapting lemmas from Chan & Airoldi (2014) for our problem and the error bound of discretization (9).
Lemma 3 Let G1 ∼ w1 and G2 ∼ w2 have respective graph transformations A1 and A2. Let W1 and W2 be the corresponding discretized graphons of w1 and w2, respectively. As n → ∞ and n20 log n
n → 0, then for any > 0,
|‖A1 −A2 ‖F − ‖W1 −W2 ‖F | ≤ 4 (12) with probability converging to 1.
Proof of Lemma 3. The proof of this lemma is inspired from Chan & Airoldi (2014). For i = {1, 2}, let matrix Wi of size n0 × n0 be the discretized graphon wi and let matrix Âi of size n0 × n0 be another transformation of graphGi based on the true permutation σ̂i, that is, σ̂i denotes the ordering of the graph Gi based on the corresponding graphon wi. In other words, the discretized graphon Wi is the expectation of Âi. The reordered graph is denoted by Gσ̂i and then Âi is obtained by,
(Âi)kl = 1
h2 h∑ k1=0 h∑ l1=0 (Gσ̂i)kh+k1,lh+l1 where h = b n n0 c and b·c is the floor function.
We bound ‖A1 −A2‖F using Â1, Â2, W1 and W2 as,
‖A1 −A2‖F ≤ ‖A1 − Â1‖F + ‖A2 − Â2‖F + ‖Â1 − Â2‖F ≤ 2 max
i={1,2} ‖Ai − Âi‖F + ‖Â1 − Â2‖F
≤ 2 max i={1,2} ‖Ai − Âi‖F + 2 max i={1,2} ‖Âi −Wi‖F + ‖W1 −W2‖F (13)
We have the following for all i using Assumption 2 on Lipschitzness of the expected degree function g(u)i of graphon wi and Lemma 3 of Chan & Airoldi (2014),
E[‖Ai − Âi‖2F ] ≤ n40 n2i
( 2 + 4C2i L 2 i
log ni ni
) + n20 ( 4C2i L 2 i
log ni ni
) ; Ci depends on Lipschitz constants of g(u)i
≤ 2n 4 0 n2 + n40 n2 4C2i L 2 i log n n + n204C 2 i L 2 i log n n = O
( n20 log n
n
) (14)
Applying Markov’s inequality to bound the probability P (
max i={1,2}
‖Ai − Âi‖2F ≥ 21 ) using (14),
P (
max i={1,2}
‖Ai − Âi‖2F ≥ 21 ) (a) ≤ ∑
i={1,2}
E [ ‖Ai − Âi‖2F ] 21 (14) = O ( n20 log n 21n ) ; (a) : Union bound
(15)
Thus asymptotically, as n→∞ and n 2 0 log n
n → 0, for all i and any 1 > 0, ‖Ai − Âi‖F < 1 with
probability converging to 1.
From Lemma 4 of Chan & Airoldi (2014), we have the following for all i,
E[‖Âi −Wi‖2F ] ≤ n40 n2i ≤ n 4 0 n2 (16)
Applying Markov’s inequality to bound the probability P (
max i={1,2}
‖Âi −Wi‖2F ≥ 22 ) using (16),
P (
max i={1,2}
‖Âi −Wi‖2F ≥ 22 ) (a) ≤ ∑
i={1,2}
E [ ‖Âi −Wi‖2F ] 22 (16) ≤ 2n 4 0 22n 2 ; (a) : Union bound
(17)
Again asymptotically as n → ∞, for all i and any 2 > 0, ‖Âi − Wi‖F < 2 with probability converging to 1. Lets assume 2 = 1 = . Substituting (15) and (17) in (13),
‖A1 −A2‖F ≤ 4 + ‖W1 −W2‖F (18)
with probability converging to 1 as n→∞ and n 2 0 log n
n → 0.
The lower bound can similarly be obtained,
‖A1 −A2‖F ≥ ‖W1 −W2‖F − 2 max i={1,2} ‖Âi −Wi‖F − 2 max i={1,2} ‖Âi −Ai‖F
≥ ‖W1 −W2‖F − 4 (19)
with probability converging to 1 as n→∞ and n 2 0 log n
n → 0.
Thus, |‖A1 −A2‖F − ‖W1 −W2‖F | ≤ 4 satisfy for any > 0 with probability converging to 1 as n→∞ and n 2 0 log n
n → 0 from equations (18) and (19).
Proof of Proposition 1 (Graph distance is consistent). Proposition 1 immediately follows from Lemmas 2 and 3 after a simple decomposition step as shown below.
|‖w1 − w2 ‖L2 − d(G1, G2)| = ∣∣∣∣‖w1 − w2 ‖L2 − d(G1, G2)± 1n0 ‖W1 −W2‖F ∣∣∣∣ ≤ ∣∣∣∣‖w1 − w2 ‖L2 − 1n0 ‖W1 −W2‖F ∣∣∣∣+ ∣∣∣∣ 1n0 ‖W1 −W2‖F − d(G1, G2) ∣∣∣∣
2,3 ≤ 4 √ 2L
n0 +
4 n0 = O
( 1
n0
) holds for any > 0 as n→∞ and n 2 0 log n
n → 0 .

heading: A.2 PROPOSITION 1 FOR SPARSE GRAPHS
text: To adapt the dense graph model, graphon, to sparse graphs, the sampling procedure is modified to Gij |Ui, Uj ∼ Bernoulli (ρw(Ui, Uj)) where ρ depends only on n. With this modification, equation 14 will be changed to
E[‖Ai − Âi‖2F ] ≤ 2 n40 n2 + n40 n2 4C2i L 2 i ρ 2 log n n + n204C 2 i L 2 i ρ log n n = O
( n20ρ log n
n
) (20)
and the proof involves modification to Lemma 2 of Chan & Airoldi (2014). We state the modified lemma and briefly sketch the proof below.
Modified Lemma 2 of Chan & Airoldi (2014). Let σ(i) be the oracle permutation such that Uσ(1) < Uσ(2) < . . . < Uσ(n). Then, if ∣∣∣∣σ(i)n − σ(j)n ∣∣∣∣ < 16L1ρ √ log n n then
∣∣dσ(i) − dσ(j)∣∣ <√ log n n
(21)
with probability at least 1 − 8 exp ( − 1
18L21
log n
ρ2
) . Conversely, if (21) holds with probability at
least 1− 8 exp ( − 1
18L21
log n
ρ2
) then
∣∣∣∣σ(i)n − σ(j)n ∣∣∣∣ < 1ρ √ log n n ( 1 3L1 + 1 3L1L2 + 1 L2 ) (22)
with probability at least 1− 40 exp ( − 1
18L21
log n
ρ2
) .
Changes to Proof in Chan & Airoldi (2014). Suppose ∣∣∣σ(i)n − σ(j)n ∣∣∣ < δρ for δ > 0 and δ < ρ.
Then, P (∣∣Uσ(i) − Uσ(j)∣∣ > 3 δρ) ≤ 4 exp(−2n δ2ρ2). Consequently,
P (∣∣g(Uσ(i))− g(Uσ(j))∣∣ > 3L1δ) ≤ P(∣∣Uσ(i) − Uσ(j)∣∣ > 3 δ
ρ
) ≤ 4 exp ( −2nδ 2
ρ2
) .
Following their proof, we have P (∣∣dσ(i) − dσ(j)∣∣ > 6L1δ ∣∣ Uσ(i), Uσ(j)) ≤ 4 exp ( −9
2 nL21δ 2
) + 4 exp ( −2nδ 2
ρ2
) ≤ 8 exp ( −2nδ 2
ρ2
) , (23)
when ρ > 23L1 . Note that there is a small mistake in their proof when Hoeffding’s inequality is
applied, where a factor of n2 is written in place of n. Putting δ = 16L1 √ logn n and considering
ρ = Ω (√ logn n ) since ρ depends on n and δ < ρ, we get
P (∣∣dσ(i) − dσ(j)∣∣ >√ log n n ∣∣ Uσ(i), Uσ(j) ) ≤ 8 exp ( − 1
18L21
log n
ρ2 ) Converse can similarly be proved.
(20) can be derived by substituting the above changed lemma in their proof. Consequently, Proposition 1 still holds under this formulation, with slight change to the condition as n → ∞ and n20ρ log n
n → 0 for ρ = Ω
(√ log n
n
) . We will not consider sparsity ρ for analysing consistency of
our algorithms in the next sections.

heading: A.3 DISTANCE BASED SPECTRAL CLUSTERING (DSC)
text: We make Assumptions 1–3 on the K graphons to analyse the Algorithm 1. We establish the consistency of this algorithm by deriving the number of misclustered graphs |M| through the following steps.
1. We establish deviation bound between the estimated distance matrix D̂ and the ideal distance matrix D (Lemma 4).
2. We formulate Davis-Kahan theorem in terms of the deviation bound using the result from Mukherjee et al. (2017) (Lemma 5).
3. We derive the number of misclustered graphs from Lemma 5.
As stated previously, D̂ in Algorithm 1 is an estimate of D ∈ Rm×m, where we define Dij = ‖wi − wj‖L2 . Note that D is a block matrix with rank K, since Dij = 0 for all Gi, Gj generated from same graphon, and equals the distance between the graphons i and j otherwise.
We derive the deviation bound for the distance matrix using Lemma 3 and the result is as follows.
Lemma 4 (Distance deviation bound) As n→∞ and n 2 0 log n
n → 0, we establish∥∥∥D̂ −D∥∥∥
F = O
( m
n0
) (24)
with probability converging to 1.
Proof of Lemma 4 (Distance deviation bound). From Proposition 1 and the definitions of D̂ij and Dij , it is easy to see that ∣∣∣D̂ij −Dij∣∣∣ = O( 1
n0
) with probability converging to 1 as n → ∞ and
n20 log n
n → 0.
Using Lemmas 2 and 3, and the definitions of D̂ij and Dij , we have∣∣∣D̂ij −Dij∣∣∣ = ∣∣∣∣ 1n0 ‖Ai −Aj‖F − ‖wi − wj‖L2 ∣∣∣∣
= ∣∣∣∣ 1n0 ‖Ai −Aj‖F ± 1n0 ‖Wi −Wj‖F − ‖wi − wj‖L2 ∣∣∣∣
≤ 1 n0 |‖Ai −Aj‖F − ‖Wi −Wj‖F |+ ∣∣∣∣‖wi − wj‖L2 − 1n0 ‖Wi −Wj‖F ∣∣∣∣
(12)
≤ 4 n0
+ 4 √ 2L
n0 with prob. → 1 asymptotically Thus asymptotically, ∣∣∣D̂ij −Dij∣∣∣ = O( 1
n0
) with probability converging to 1.
Hence, ‖D̂ −D‖F = O ( m
n0
) with probability converging to 1 as n→∞ and m 2n20 log n
n → 0.
A variant of Davis-Kahan theorem (Mukherjee et al., 2017) and the derived deviation bound (24) are used to prove the following lemma.
Lemma 5 (Davis-Kahan theorem) Let V and V̂ be them×K matrices whose columns correspond to the leadingK eigenvectors ofD and D̂, respectively. Let γ be theK-th smallest eigenvalue value of D in magnitude. As n→∞ and n 2 0 log n
n → 0, there exists an orthogonal matrix Ô such that,∥∥∥V̂ Ô − V ∥∥∥ F = O ( m
γn0
) (25)
with probability converging to 1.
Proof of Lemma 5 (Davis-Kahan theorem). A variant of Davis Kahan theorem from Proposition A.2 of Mukherjee et al. (2017) states the following for matrix D of rank K. Let V̂ and V be m ∗K matrices whose columns correspond to the leading K eigenvectors of D̂ and D, respectively, and γ be the K-th smallest eignenvalue of D in magnitude, then there exists an orthogonal matrix Ô of size K ∗K such that,∥∥∥V̂ Ô − V ∥∥∥
F ≤
4 ∥∥∥D̂ −D∥∥∥
F
γ
4 = O
( m
γn0
) as n→∞ and m 2n20 log n
n → 0.
The number of misclustered graphs is |M| ≤ 8mT ‖V̂ Ô−V ‖2F where mT is the maximum number of graphs generated from a single graphon (Mukherjee et al., 2017). Since mT = O(m), |M| = O ( m3
γ2n20
) by substituting (25) in |M|. Hence proving Theorem 1.
Proof of Theorem 1. The number of misclustered graphs |M| ≤ 8mT ∥∥∥V̂ Ô − V ∥∥∥2
F from
Mukherjee et al. (2017). Thus, we prove the theorem using Lemma 5. That is, as n → ∞ and m2n20 log n
n → 0,
|M| ≤ 8mT ∥∥∥V̂ Ô − V ∥∥∥2
F
5 = O
( m3
γ2n20
)
Proof of Corollary 1. This corollary deals with a special case where K = 2 and equal number of graphs are generated from the two graphons w and w′. Therefore,mT in the number of misclustered graphs |M| is m/2. The ideal distance matrix D will be of size m×m with 0 and ‖w − w′‖L2 as entries depending on whether the samples are generated from the same graphon or not. For such a block matrix D, the two non zero eigenvalues are ±m
2 ‖w−w′‖L2 . Therefore, γ is
m 2 ‖w−w′‖L2 .
Corollary 1 can be derived by substituting the derived γ in the number of misclustered graphs |M| in Theorem 1 as shown below.
|M| = O (
m
‖w − w′‖2L2n 2 0 ) Let us assume ‖w −w′‖L2 ≥ C m
n0 where C is a large constant, then as n→∞, m
2n20 log n
n → 0,
|M| → 0.

heading: A.4 SIMILARITY BASED SEMI-DEFINITE PROGRAMMING (SSDP)
text: We make Assumptions 1–3 on the K graphons to study the recovery of clusters from Algorithm 2. The proof structure for cluster recovery stated in Theorem 2 is as follows:
1. We establish deviation bound between the estimated similarity matrix Ŝ and the ideal similarity matrix S (Lemma 6).
2. We derive the recoverability condition by adapting Proposition 1 of Perrot et al. (2020) and the obtained deviation bound (Lemma 7).
The ideal similarity matrix S ∈ Rm×m is symmetric with K ×K block structure, and S = ZΣZT where Z ∈ {0, 1}m×K be the clustering membership matrix and Σ ∈ RK×K such that Σll′ represents ideal pairwise similarity between graphs from clusters Cl and Cl′ . From the definition of Sij ,
Σll′ = exp
( −‖wl − wl‖L2
σlσl′
) where wl and wl′ are graphons corresponding to clusters Cl and Cl′ ,
respectively. Ŝ is the estimated similarity matrix of S as mentioned earlier. Since X̂ ∈ RK×K is the normalised clustering matrix, X̂ = ZN−1ZT where N is a diagonal matrix with 1
|C1| , . . . ,
1
|CK | .
We derive the deviation bound for the similarity matrices using Lemma 3 and the result is as follows.
Lemma 6 (Similarity deviation bound) As n→∞, n 2 0 log n
n → 0, we establish
|Ŝij − Sij | = O ( 1
n0
) (26)
with probability converging to 1. Hence, from the result ‖Ŝ − S‖F = O ( m
n0
) with probability
converging to 1.
Proof of Lemma 6 (Similarity deviation bound). We derive the bound using Lemmas 2 and 3, and the definitions of Ŝij and Sij .
Ŝij = exp
( −‖Ai −Aj‖F
n0σiσj
) Consider σi = σj = σ
3 ≥ exp ( −‖Wi −Wj‖F + 4
n0σ2
) with probability→ 1 asymptotically
2 ≥ exp ( −‖wi − wj‖L2
σ2
) exp ( −4 + 4 √ 2L
n0σ2
) exp(−x) ≥ 1− 2x for x > 0
≥ Sij ( 1− 8( + √ 2L)
n0σ2
)
≥ Sij − Sij 8( +
√ 2L)
n0σ2 Sij ∈ [0, 1]
≥ Sij − 8( +
√ 2L)
n0σ2 (27)
Ŝij = exp
( −‖Ai −Aj‖F
n0σ2 ) 3 ≤ exp ( −‖Wi −Wj‖F − 4
n0σ2
) with probability→ 1 asymptotically
2 ≤ exp ( −‖wi − wj‖L2
σ2
) exp ( 4 + 4 √ 2L
n0σ2
) exp(x) ≤ 1 + 2x for x > 0
≤ Sij ( 1 + 8( + √ 2L)
n0σ2
)
≤ Sij + Sij 8( +
√ 2L)
n0σ2 Sij ∈ [0, 1]
≤ Sij + 8( +
√ 2L)
n0σ2 (28)
Thus, from (27) and (28), we get |Ŝij −Sij | ≤ 8( +
√ 2L)
n0σ2 = O
( 1
n0
) for any , with probability
converging to 1 as n → ∞ and n 2 0 log n
n → 0. Hence, ‖Ŝ − S‖F ≤
8m( + √ 2L)
n0σ2 = O
( m
n0
) ,
with probability converging to 1 as n→∞ and m 2n20 log n
n → 0.
The condition for exact recovery of clusters is derived by adapting Proposition 1 of Perrot et al. (2020). The proposition states the recoverability condition for such an SDP defined in (4) in terms of the similarity deviation bound. Thus, we use the derived bound in Lemma 6 and establish condition on the L2-distance to satisfy the proposition from Perrot et al. (2020). First, we state the adapted proposition.
We define ∆1 and ∆2 as,
∆1 = min l 6=l′ (1− Σll′) and ∆2 = max ij |Ŝij − Sij |.
Then, the following should be satisfied for X̂ to be the unique optimal solution of the SDP in (4):
‖Ŝ − S‖F ≤ min l |Cl|min { ∆1 2 ,∆1 − 6∆2 } .
The minimum cluster size min l |Cl| in our case is 1. Consequently, the recoverability condition is derived and is as follows.
Lemma 7 (Recoverability of clusters) As n→∞, m 2n20 log n
n → 0, the min l 6=l′ ‖wl−w′l‖L2 should be Ω ( m
n0
) so that X̂ is the unique optimal solution of the SDP (4).
Proof of Lemma 7 (Recoverability of clusters). We derive the condition to satisfy the stated proposition.
∆1 = 1−maxl 6=l′ Σll′ and ∆2 = min ij |Ŝij − Sij |. The minimum cluster size in our case can be 1.
The analyses of the two cases of the Proposition is as follows. Case 1. Let us assume ∆2 ≤ ∆1 12 , then min { ∆1 2 ,∆1 − 6∆2 } will be ∆12 . Therefore,
8m( + √ 2L)
n0σ2 ≤ 1 2 − 1 2 max l 6=l′ exp
( −‖wl − wl
′‖L2 σ2 ) exp ( −min
l 6=l′ ‖wl − wl′‖L2 σ2
) ≤ 1− 16m( + √ 2L)
n0σ2
min l 6=l′ ‖wl − wl′‖L2 σ2 ≥ − log
( 1− 16m( + √ 2L)
n0σ2
)
min l 6=l′ ‖wl − wl′‖L2 ≥ σ2 ∞∑ k=1
( 16m( + √ 2L)
n0σ2
)k 1
k
min l 6=l′ ‖wl − wl′‖L2 = Ω
( m
n0
) (29)
Case 2. Let us assume ∆2 > ∆1 12
, then min {
∆1 2 ,∆1 − 6∆2 } will be ∆1 − 6∆2. Therefore,
8m( + √ 2L)
n0σ2 ≤ 1−max l 6=l′ exp
( −‖wl − wl
′‖L2 σ2
) − 6 ∗ 8( + √ 2L)
n0σ2 with probability→ 1 aymptotically
min l 6=l′ ‖wl − wl′‖L2 σ2 ≥ − log
( 1− 8(m+ 6)( + √ 2L)
n0σ2
)
min l 6=l′ ‖wl − wl′‖L2 ≥ σ2 ∞∑ k=1
( 8(m+ 6)( + √ 2L)
n0σ2
)k 1
k
min l 6=l′ ‖wl − wl′‖L2 = Ω
( m
n0
) (30)
Thus, from (29) and (30), we must satisfy minl 6=l′ ‖wl − wl′‖L2 = Ω ( m
n0
) for the Proposition to
hold. Consequently, Theorem 2 is the direct reflection of this lemma.

heading: A.5 GRAPH TWO-SAMPLE TESTING
text: Theorem 3 of two-sample testing is proved by deriving the probability of Type-1 and Type-2 errors. We make Assumptions 1–3 for this case. Let W1 and W2 be the n0×n0 discretized graphons of w1 and w2, respectively, obtained using (6). Then, the alternate hypothesis Ha can be rewritten using Lemma 2 in the following way,
1
n0 ‖W1 −W2‖F +
4 √ 2L
n0
(9) ≥ φ
‖W1 −W2‖F ≥ n0φ− 4 √ 2L = ρ (31)
We derive the probability of the errors using Lemma 3 and is stated in the following lemmas.
Lemma 8 (Probability of Type-1 error) The probability of Type-1 error, i.e. rejecting the null hypothesis when it is actually true, is
P(T = 1|H0 : True) ≤ C ξ2 log n n (32)
where C depends only on the Lipschitz constants.
Proof of Lemma 8 (Probability of Type-1 error). The Type-1 error is rejecting H0 when it is true. Therefore, in this scenario, ‖W1 −W2‖F = 0. Thus, from Lemma 3, (15) and (17), we have ‖A1 −A2‖F ≤ 4 with 1− C 2 n20 log n
n probability. Therefore, the probability of Type-1 error is,
P(T = 1|H0 : True) = P(d(G1, G2) ≥ ξ) = P(‖A1 −A2‖F ≥ n0ξ) err only when n0ξ ≤ 4 (15) ≤ C
ξ2n20
n20 log n
n
Lemma 9 (Probability of Type-2 error) The probability of Type-2 error, i.e. accepting the null hypothesis when the alternate hypothesis is actually true, is
P(T = 0|Ha : True) ≤ C1(
φ− 4 √
2C2 n0 )2 log nn (33) where C1 and C2 depend only on the Lipschitz constants.
Proof of Lemma 9 (Probability of Type-2 error). The Type-2 error is evaluating to null hypothesis when the alternate hypothesis is true. Therefore, from (31) ‖W1 −W2‖F ≥ ρ. From Lemma 3, (15) and (17),
‖A1 −A2‖F ≥ ‖W1 −W2‖F − 4 with probability 1− C1 2 n20 log n
n ≥ ρ− 4
The probability of Type-2 error is,
P(T = 0|Ha : True) = P(d(G1, G2) < ξ) = P(‖A1 −A2‖F < n0ξ) err only when n0ξ ≤ ρ− 4 ; let n0ξ = 4
≤ P(‖A1 −A2‖F < ρ
2 )
≤ C1 ρ2
n20 log n
n
We get the probability by substituting ρ n0 = φ− 4
√ 2L
n0 from (31) in the above equation. Theorem
3 can be proved by asymptotic analysis of Lemmas 8 and 9.

heading: B DSC AND SSDP ALGORITHMS
text: The proposed algorithms DSC and SSDP are described as follows:
Algorithm 1: Distance based Spectral Clustering (DSC) input : Adjacency matrices G1, ..., Gm, histogram size n0 output: K clusters C1, ..., CK Construct distance matrix Compute D̂ ∈ Rm×m, where D̂ij = d(Gi, Gj)
Clustering Apply spectral clustering to D̂ with K number of clusters resulting in C1, ..., CK
Algorithm 2: Similarity based Spectral Clustering (SSDP) input : Adjacency matrices G1, ..., Gm, histogram size n0 output: K clusters C1, ..., CK Construct similarity matrix Compute Ŝ ∈ Rm×m, where Ŝij = exp ( −d(Gi,Gj)σiσj
) with σ1 = . . . = σn
Clustering Find X̂ using (4) and apply standard spectral clustering to X̂ resulting in C1, ..., CK .

heading: C EXPERIMENTAL DETAILS
text: In this section, we present experimental details and additional experiments.

heading: C.1 SIMULATED DATA - HEATMAP OF GRAPHONS
text: Figure 4 shows the heatmap of the considered four graphons W1,W2,W3 and W4. We sample graphs from these graphons for the experiments.
C.2 CHOICE OF n0 We validate the theoretically deduced bound for n0 = O( √ n/ log n) by sampling 5 graphs with a fixed number of nodes n from each of the four graphons, in total 20 graphs, and measuring the performance of DSC and SSDP for different n0 = {5, 10, 15, 20, 25, 30}. We perform three simulations with n = {50, 100, 500} and fix neighbourhood of one in SSDP. Figure 5 shows the average error of both the algorithms over 5 independent trials. Based on the theoretical considerations for n0 ( √ n/ log n), we evaluate n0 = {5, 7, 15} for n = {50, 100, 500}, respectively. The experimental results show that the derived bound for n0 serves as a reasonable choice (if not the best) for both DSC and SSDP irrespective of n. Hence, the choice of n0 can be deterministic and adaptive with respect to n, thus making our algorithms parameter-free.

heading: C.3 EXPERIMENTAL RESULTS USING ADJUSTED RAND INDEX (ARI)
text: In this section, we provide the results for evaluation of algorithms on simulated and real data under the same setting as described in Section 3.3. Figure 6 shows the evaluation of all the discussed algorithms using ARI, where the observations made from error rate hold.

heading: C.4 EVALUATION ON LARGE SIMULATED DATA
text: As mentioned in Section 3.3, we evaluate algorithms except GNTK and NCGMM on large graphs sampled from the four graphons W1,W2,W3 and W4 with nodes between 100 and 1000. Figure 7 shows the results measured using average error rate and average ARI, respectively. The proposed algorithm DSC outperforms the others in all the case and SSDP stands second or third best, as observed in simulated data with small graphs.

heading: C.5 COMPUTATION TIME OF ALGORITHMS
text: Table 1 shows the time (measured in seconds) taken for the algorithms on the considered dataset combinations. Three clusters also show similar behavior to four clusters.

heading: C.6 SCALABILITY EXPERIMENT
text: We evaluate the scalability of the considered algorithms using simulated data by measuring the time taken for clustering 40 random graphs, 10 sampled from each of the graphons W1,W2,W3 and W4. We did 7 experiments in which the size of the sampled graphs are varied as [50,max size] where max size = {100, 200, 300, 400, 500, 600, 700}. Figure 8 shows the experimental results which illustrates high scalability of DSC, SSDP and NCLM over other algorithms. Note that the experiment shows NCLM as scalable as DSC and SSDP since the sampled graphs are small.

heading: C.7 TWO-SAMPLE TESTING
text: In this section, we evaluate the efficacy of the proposed test T with different d(G1, G2) by varying the graph sizes n. We consider n = {50, 100, 150} and fix n0 = 10 from the theoretical bound for evaluating the test T . The power is computed using the test T for the significance level 0.05, and the plots in Figure 9 show the average power computed over 500 trials of bootstrapping 100 samples generated from all pairs of graphons for d(G1, G2) as our proposed distance, log moments
and MMD, respectively. From the result for graph sizes (50, 100), we observe that the graphon pair (W2,W3) is not easily distinguishable (low H0 rejection probability), which can be explained by their respective L2 distance that is shown in the left plot of Figure 3. This issue does not arise in testing larger graphs as the result shows for graph sizes (100, 200) and (150, 300). Therefore, test T with the proposed distance can distinguish between pairs of graphons that are quite close provided that the observed graphs are sufficiently large, thus proving to be consistent. On the other hand, log moments and MMD based tests show weakness in distinguishing the graphons, where log moments based test T accepts the null hypothesis in most cases even when the graphons are different for all graph sizes. For instance, the result for graphon pairW1 andW4 is indistinguishable using log moments statistic for any graph size. On the contrary, MMD based test T rejects the null hypothesis almost always for larger graphs (diagonal values in all graph cases). Thus, we conclude that the proposed test T in 5 is consistent and this experiment illustrates the efficiency of the test T compared to other plausible test statistics.
Subsequently, we evaluate the efficacy of the above tests on the discussed real datasets – Bioinformatics and Social Networks. We consider graphs from a dataset to belong to a population and hence the objective of the test statistic is to distinguish graphs from different populations, that is, graphs from two different datasets. Since the populations are not known and the real graphs are treated as representatives of the population, we compute p-value of the test instead of power to measure the efficacy. The p-value is the evidence for rejecting the null hypothesis which implies that the smaller the p-value, the stronger the evidence that the null hypothesis should be rejected. Therefore, the pvalue should be high (greater than the significance level) for graphs from same population and low (' 0) for graphs from different populations. Figure 10 shows the result for both the dataset cases and different tests. From the results, it is clear that the log moments and MMD based tests are poor
and inefficient on real datasets as log moments based test has high acceptance of null hypothesis for almost all the pair of graphs from any population and MMD based test rejects the null hypothesis always except when the graphs are the same. Whereas, the test using our proposed distance perform well on large graphs from Social Networks datasets, for instance, S4 with other datasets and within itself. This test performs well even for small graphs when compared to the other two tests.

Review summaries:
Motivated by the sort-and-smooth graphon estimator of Chan and Airoldi, 2014, this paper proposes two new clustering algorithms (graph distance based spectral clustering and similarity based semidefinite programming) for multiple graphs observed without vertex correspondence. The idea is to use the graphon approximation method to obtain a histogram estimator with the same number of bins for each graph and subsequently apply an existing spectral clustering approach on the graphon based distance matrix. Theoretical properties of their approach are studied (under reasonable smoothness assumptions on the generating process i.e. graphon) and consistency is established. The graph distance metric is applied to test for similarity of two collections of graphs.  

The paper proposes a graph distance based on graphons that can operate on the small sample size regime; this is possible by exploiting a large number of nodes. The proposed paper is meant to address two shortcomings that the authors have found in the literature:  - a lack of study of current graph processing methods beyond graph-level classification, and  - a lack of theoretically sound methods.    The paper also shows theoretical guarantees associated with two existing clustering methods and a two-sample statistical tests when operating with the proposed distance.