GTST: A Python Package for Graph Two-Sample Testing

Ragnar L. Gudmundarson; Gareth W. Peters

(1) Overview

Introduction

It is increasingly common to obtain sampled data in the form of graph or network realisations. This is either through construction of graph structures as summary statistics for multivariate spatial-temporal data or directly as a graph-network data set from say a social network realisation or the like. In such settings in practice, one often needs to draw inferences about two samples of such graphs or networks. That is, to test whether the collections of graphs/networks in one sample were generated from the same distribution as the collection of graphs/networks in the second sample.

Kernel methods have proven to be useful in pattern recognition tasks such as classification and can be further extended as an inference procedure to two-sample hypothesis testing on structured data. The method embeds the graphs into a reproducing kernel Hilbert space (RKHS) via a feature map which is then extended further to the embedding of a probability distribution. The two-sample null hypothesis is that the generating mechanism behind the two samples is the same and the test statistic, which is called the maximum mean discrepancy (MMD), is the largest distance between the means of the two sample embeddings [].

Graph kernels are already well established and widely used for solving classification tasks on graphs and can further be used to compare samples of graphs and to perform graph screening []. They provide a very flexible way of comparing graphs as they exist for a wide range of different graph structures, for example, weighted, directed, labelled, and attributed graphs. Their performance depends on their expressiveness, that is, their ability to distinguish non-isomorphic graphs. The difficulty of distinguishing two samples of graphs varies strongly based on the type of graphs.

Graph two-sample hypothesis testing is a problem that frequently arises in various disciplines, for example in bio informatics [], community detection [], and risk management []. Graph two-sample hypothesis have mostly been performed by using graph statistics such has the degree centrality and shortest paths. Although these methods can often give good performances they fail to take into account various attributes that are often present in real graphs such as node labels, edge labels, node attributes and edge weights. When the kernel two-sample hypothesis testing was introduced [] a flood gate opened to allow for testing of such attributes and therefore providing a flexible way of performing two-sample hypothesis testing. Luckily, there also exists a vast literature on graph kernels [, ]. Until now, there is no package which allows one to estimate graphs from real valued data matrices and perform hypothesis testing in a flexible manner. The package GTST provides three functionalities

Functions to estimate the two-sample hypothesis test statistic;
Functions to calculate different graph kernels; and
A function to estimate graphs from data matrices.

The package allows other network science researches who may not be programmers to use the MMD testing framework.

There exists a python package called GraKel [] which is dedicated to calculating various graph kernels. The package is very user-friendly so the GTST user can use all graph kernels available in the Grakel package. The GTST package extends the choice of graph kernels available for use in graph two-sample testing, by also developing a collection of kernels not available in GraKel. Such as, the fast random walk kernels which is based on ideas from [] along with an additional fast random walk kernel for edge-labelled graphs; The Wasserstein Weisfeiler-Lehman Graph kernel [] whose original code was adjusted for the package needs. The Deep Graph kernel [], and the graph neural tangent kernel [] whose original code was adjusted for the package needs. The MONK estimator, which is a robust estimator of the MMD, was developed by MONK [] and they do provide the code online and in a packaged environment. However, we have adjusted the code slightly to allow for robust comparing of samples of different sizes. The MMDGraph then estimates the p-value of test by using a bootstrap or a permutation sampling scheme. The package also allows for estimating graphs using sklearn’s graphical lasso []. Additional preprocessing can be done by using the nonparanormal transform []. The best graph is found by using the EBIC criterion []. GTST assumes that the graphs passed are a networkx object []. One can additionally use pre-computed kernels to perform tests.

Implementation and architecture

Let G(V, S) denote a graph with vertex set V and edge set S. In the two-sample testing of graph-valued data, we assume we are given two sets of samples/observations that comprise collections of graph-valued data {G₁, …, G_n} and ${G ′ 1, …, G ′ n ′}$ M1 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \{{G^{\prime}_1},\cdots,{G^{\prime}_{n^{\prime}}}\} \] \end{document} where $G i, G ′ j ∈ Ω, ∀ i, j$ M2 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {G_i},{G^{\prime}_j} \in \Omega,\quad \forall i,j \] \end{document} . The graphs in the two samples are all generated independently from two probability spaces $(Ω, ℱ, ℙ)$ M3 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ (\Omega,{\mathcal F},\mathbb{P}) \] \end{document} and $(Ω, ℱ, ℚ)$ M4 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ (\Omega,{\mathcal F},\mathbb{Q}) \] \end{document} , and the goal is to infer whether ℙ = ℚ. We note that in this general testing framework the vertex sets and edge-sets do not have to be equal, but they can be common if it is desirable for the application. This is therefore a very general testing framework. In the simplest case, where we assume both sets of sample graphs come from a common set of vertices, then the sample space Ω contains all possible edges that can occur in a graph G, that is $Ω = {(v 1, v 2), …, (v 1, v | V |), (v 2, v 1), …, (v | V |, v | V | − 1)}$ M5 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \Omega = \{ ({v_1},{v_2}), \ldots,({v_1},{v_{|V|}}),({v_2},{v_1}), \ldots,({v_{|V|}},{v_{|V| - 1}})\} \] \end{document} . As the sample space is discrete we can define the σ-algebra as the power set of Ω, namely, $ℱ = (Ω)$ M6 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\mathcal F} = {\mathcal P}(\Omega) \] \end{document} . The probability function $ℙ : ℱ ↦ [0, 1]$ M7 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \mathbb{P}:{\mathcal F} \mapsto [0,1] \] \end{document} then defines the probability of obtaining a certain graph in the sample set of graph-valued data. As an example we can define for instance a population distribution to be uniform

ℙ (G (V, S)) = 1 (M | S |)

M8 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \mathbb{P}(G(V,S)) = \frac{1}{{\left(\begin{array}{c} M\\ |S| \end{array} \right)}} \] \end{document}

where $M = (2 | V |)$ M9 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ M = {{\left(\begin{array}{c} |V|\\ 2 \end{array}\right)}} \] \end{document} is the total number of possible edges and G(V, S) is a graph with |V| vertices and |S| number of edges. This setting is illustrated in Figure 1.

Figure 1

The graph two sample testing scenario. Here we have observed 4 graphs from ℙ and 3 graphs from ℚ. The sample space of ℙ and ℚ is the same ( $(25)$ M33 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \left({_2^5} \right) \] \end{document} possible edges).

Now, returning to the concept of two-sample testing for graph-valued data. The goal is to infer whether the two samples of graphs are generated according to the same distribution. This involves developing a statistical test $T ({G} i = 1 n, {G ′} i = 1 n ′)$ M10 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ T(\{ G\} _{i = 1}^n,\{ G^{\prime}\} _{i = 1}^{n^{\prime}}) \] \end{document} to determine from the population samples whether there is sufficient evidence to reject the null hypothesis that both population distributions generating the two samples of graphs are equivalent, where $T ({G} i = 1 n, {G ′} i = 1 n ′) : {G} i = 1 n × {G ′} i = 1 n ′ ↦ {0, 1}$ M11 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ T(\{ G\} _{i = 1}^n,\{ G^{\prime}\} _{i = 1}^{n^{\prime}}):\{ G\} _{i = 1}^n \times \{ G^{\prime}\} _{i = 1}^{n^{\prime}} \mapsto \{ 0,1\} \] \end{document} is a function that distinguishes between the null hypothesis and the alternative hypothesis:

H 0 : ℙ = ℚ H 1 : ℙ ≠ ℚ .

M12 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \begin{array}{*{20}{c}} {{H_0}:}&{\mathbb{P} = \mathbb{Q}}\\ {{H_1}:}&{\mathbb{P} \ne \mathbb{Q}.} \end{array} \] \end{document}

The test statistic used in this case is the largest distance between the expectation of some function w.r.t. to the two probability distributions. Let $ℋ B$ M13 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {{\mathcal H}_B} \] \end{document} be a class of functions $f : Ω → R$ M14 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ f:\Omega \to R \] \end{document} . The maximum mean discrepancy (MMD) is defined as:

MMD [ℋ B, ℙ, ℚ] : = sup f ∈ ℋ B (E G ~ ℙ [f (G)] − E G ′ ~ ℚ [f (G ′)]) .

M15 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\rm{MMD}}[{{\mathcal H}_B},\mathbb{P},\mathbb{Q}]: = \mathop {\sup }\limits_{f \in {{\mathcal H}_B}} ({E_{G\sim\mathbb{P}}}[f(G)] - {E_{G^{\prime}\sim\mathbb{Q}}}[f(G^{\prime})]). \] \end{document}

Where the class of function, $ℋ B$ M16 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {{\mathcal H}_B} \] \end{document} , is the unit ball in a RKHS, then the squared population MMD becomes kernalized:

MMD 2 [ℋ 1, ℙ, ℚ] = E ℙ, ℙ [k (G, G)] − 2 E ℙ, ℚ [k (G, G ′)] + E ℚ, ℚ [k (G ′, G ′)],

M17 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\rm{MM}}{{\rm{D}}^2}[{{\mathcal H}_1},\mathbb{P},\mathbb{Q}] = {E_{\mathbb{P},\mathbb{P}}}[k(G,G)] - 2{E_{\mathbb{P},\mathbb{Q}}}[k(G,G^{\prime})] + {E_{\mathbb{Q},\mathbb{Q}}}[k(G^{\prime},G^{\prime})], \] \end{document}

where k is some graph kernel. The kernel matrix plays a central role and is defined as:

K i j = k (G i, G j),

M18 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {{\bf K}_{ij}} = k({G_i},{G_j}), \] \end{document}

where ${G i} i = 1 2 n$ M19 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \{{G_i}\}_{i = 1}^{2n} \] \end{document} is the data. Note, both sampels are included in the data. In the context of two sample graph testing, we have two data sources, sample 1 of graphs assumed drawn from P and sample 2 of graphs assumed drawn from ℚ. It can be good to order the kernel matrix such that K has a block structured as follows:

(1)

K = [K ℙ ℙ K ℙ ℚ K ℚ ℙ K ℚ ℚ],

M20 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\bf K} = \left[{\begin{array}{*{20}{c}} {{{\bf K}_{\mathbb{PP}}}}&{{{\bf K}_{\mathbb{PQ}}}}\\ {{{\bf K}_{\mathbb{QP}}}}&{{{\bf K}_{\mathbb{QQ}}}} \end{array}} \right], \] \end{document}

where K_ℙℚ is the Kernel function evaluated at data points within the sample coming from the unknown distributions ℙ and ℚ. K_ℙℙ, K_ℚℙ and K_ℚℚ are defined analogously. Note $K ℙ ℚ = K ℚ ℙ T$ M21 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {{\bf K}_{\mathbb{PQ}}} = {\bf K}_{\mathbb{QP}}^T \] \end{document} .

Different estimates of MMD2 are:

Unbiased estimate

An unbiased estimate of MMD² for n and n’ samples of graphs is given by:

(2)

MMD^u 2 [ℋ B, ℙ, ℚ] = 1 n (n − 1) ∑ i = 1 n ∑ j ≠ i n k (G i, G j) + 1 n ′ (n ′ − 1) ∑ i = 1 n ′ ∑ j ≠ i n ′ k (G ′ i, G ′ j) − 2 n n ′ ∑ i = 1 n ∑ j = 1 n ′ k (G i, G ′ j) .

M22 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \begin{array}{l} \widehat {{\rm{MMD}}}_u^2[{{\mathcal H}_B},\mathbb{P},\mathbb{Q}] = \frac{1}{{n(n - 1)}}\sum\limits_{i = 1}^n {\sum\limits_{j \ne i}^n k } ({G_i},{G_j}) + \frac{1}{{n^{\prime}(n^{\prime} - 1)}}\sum\limits_{i = 1}^{n^{\prime}} {\sum\limits_{j \ne i}^{n^{\prime}} k } ({{G^{\prime}_i}},{{G^{\prime}_j}})\\ \quad \quad \quad \quad \quad \quad \quad \;\,{\kern 1pt} - \frac{2}{{nn^{\prime}}}\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^{n^{\prime}} k } ({G_i},{{G^{\prime}_j}}). \end{array} \] \end{document}

Biased estimate

A biased estimate of MMD² for n and n’ samples of graphs is given by:

(3)

MMD^b 2 [ℋ B, ℙ, ℚ] = 1 n 2 ∑ i = 1 n ∑ j = 1 n k (G i, G j) + 1 n ′ 2 ∑ i = 1 n ′ ∑ j = 1 n ′ k (G ′ i, G ′ j) − 2 n n ′ ∑ i = 1 n ∑ j = 1 n ′ k (G i, G ′ j) .

M23 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \begin{array}{l} \widehat {{\rm{MMD}}}_b^2[{{\mathcal H}_B},\mathbb{P},\mathbb{Q}] = \frac{1}{{{n^2}}}\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n k } ({G_i},{G_j}) + \frac{1}{{{{n^{\prime}}^2}}}\sum\limits_{i = 1}^{n^{\prime}} {\sum\limits_{j = 1}^{n^{\prime}} k } ({{G^{\prime}_i}},{{G^{\prime}_j}})\\ \quad \quad \quad \quad \quad \quad \quad \;\,{\kern 1pt} - \frac{2}{{nn^{\prime}}}\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^{n^{\prime}} k } ({G_i},{{G^{\prime}_j}}). \end{array} \] \end{document}

Unbiased linear time estimate

An unbiased can be computed in O(n) time. Assume that n = n’ and define n₂ = ⌊n/2⌋, where ⌊∙⌋ is the floor function, then the linear estimate is computed as:

MMD^l 2 [ℋ B, ℙ, ℚ] = 1 n 2 ∑ i = 1 n 2 (k (G 2 i − 1, G 2 i) + k (G ′ 2 i − 1, G ′ 2 i) − k (G 2 i − 1, G ′ 2 i) − k (G 2 i, G ′ 2 i − 1)),

M24 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \begin{array}{l} \widehat {{\rm{MMD}}}_l^2[{{\mathcal H}_B},\mathbb{P},\mathbb{Q}] = \frac{1}{{{n_2}}}\sum\limits_{i = 1}^{{n_2}} (k({G_{2i - 1}},{G_{2i}}) + k({{G^{\prime}_{2i - 1}}},{{G^{\prime}_{2i}}}) - k({G_{2i - 1}},{{G^{\prime}_{2i}}})\\ \quad \quad \quad \quad \quad \quad \quad \;\; - k({G_{2i}},{{G^{\prime}_{2i - 1}}})), \end{array} \] \end{document}

while $MMD^l 2$ M25 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \widehat {{\rm{MMD}}}_l^2 \] \end{document} has higher variance than $MMD^u 2$ M26 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \widehat {{\rm{MMD}}}_u^2 \] \end{document} , it is computationally much more appealing.

Robust Estimate

Consider the partition ${S q} q ∈ 1 : Q$ M27 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\{{S_q}\}_{q \in 1:Q}} \] \end{document} of $| S q | = n / Q$ M28 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ |{S_q}| = n/Q \] \end{document} partitions. By the representer theorem [] we can express a function f within the RKHS space as:

f = ∑ i a i k (⋅, x i) + ∑ i b i k (⋅, y i) .

M29 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ f = \sum\limits_i {{a_i}} k(\cdot,{{\rm{x}}_i}) + \sum\limits_i {{b_i}} k(\cdot,{{\rm{y}}_i}). \] \end{document}

A robust MMD estimator is found by solving:

max c : c T K c ≤ 1 m e d q ∈ 1 : Q {1 | S q | [1 q, − 1 q] K c},

M30 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ \mathop {\max}\limits_{{\boldsymbol c}:{{\boldsymbol c}^T}{\bf K}{\boldsymbol c} \le 1} me{d_{q \in 1:Q}}\left\{{\frac{1}{{|{S_q}|}}[{{\bf 1}_q}, - {{\bf 1}_q}]{\bf K}{\boldsymbol c}} \right\}, \] \end{document}

where $c = [a, b] ∈ ℝ 2 n, K = [K ℙ ℙ, K ℙ ℚ; K ℚ ℙ, K ℚ ℚ] ∈ ℝ 2 n × 2 n$ M31 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {\boldsymbol c} = [{\boldsymbol a},{\boldsymbol b}] \in {\mathbb{R}^{2n}}, {\bf K} = [{{\bf K}_{\mathbb{PP}}},{{\bf K}_{\mathbb{PQ}}};{{\bf K}_{\mathbb{QP}}},{{\bf K}_{\mathbb{PP}}}] \in {\mathbb{R}^{2n \times 2n}} \] \end{document} , where K_ℙℚ is the Kernel function evaluated at data points within the sample coming from the unknown distributions ℙ and ℚ, and $1 q ∈ ℝ n$ M32 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[ {1_q} \in {\mathbb{R}^n} \] \end{document} is an indicator vector of block q. For more details see [, ].

The package includes graph kernels such as:

Random walk kernels which can be used on weighted, directed, undirected, and bipartite graphs with node labels, node attributes, edge labels, and edge attributes.
Deep graph kernels can be used on node labelled binary graphs.
Graph neural tangent kernels can be used on binary graphs with node labels/attributes.
Wasserstein Weisfeiler-Lehman graph kernels can be used on undirected graphs with node labels and can be used on node-attributed graphs with the right embedding scheme as well.

Other graph kernels can be utilized via the graph kernel library GraKel [] which is dedicated to calculating various graph kernels. The package is very user-friendly so the GTST user can use all graph kernels available in the Grakel package. For more details on various graph kernels see [, , ].

Example

The workflow is as follows: 1) Use two data arrays to estimate two sequences/samples of graphs using the graphical lasso []. This step can be skipped if the practitioner already has the samples of graphs in a networkx format []. 2) Select a graph kernel. 3) Select an estimator of the MMD or try multiple estimators to obtain a p-value. A quick example for the random walk kernel is:

import numpy as np
import networkx as nx
import GTST
 
# Generate random samples
n1 = n2 = 50 # sample sizes
# generate two samples
 
g1 = [nx.fast_gnp_random_graph (30,0.3) for _ in range (n1)]
g2 = [nx.fast_gnp_random_graph (30,0.4) for _ in range (n2)]
 
# Random Walk, r is number of eigen-pairs, c is the discount constant
MMD_out = GTST.MMD()
MMD_out.fit(G1 = g1,
        G2 = g2,
        kernel = 'RW_ARKU_PLUS',
        mmd_estimators = 'MMD_u',
        r = 6,
        c = 0.001)
print (f" RW_ARKU_plus {MMD_out.p_values}")

For more exampels please see https://github.com/ragnarlevi/GTST.

Quality control

The requirements for running test for GTST is listed in the requirements_dev.txt in the Github page. The tests involve testing if kernels matrices are positive semi definite and whether the kernels, test statistic, and permutation method for the p-value are able to reject the null when the null is “extremely” false using known random data sets. Once the required testing packages have been installed, the tests can be performed by running the command:

pytest

in the root folder, which will take around 15 minutes. This will generate a coverage report which can be found in the htmlcov directory. To view it run:

cd htmlcov python - m http.server

and open the localhost link (something like http://localhost:8000/) in a browser.

To test GTST in a clean environment for all python versions from 3.7–3.10, we use Tox. This can be achieved by running

tox

in the root directory. Note that this takes significantly longer to run, so is best performed as a final check. Whenever the code is pushed to the remote repository, the Tox test suite is automatically run using GitHub actions. To investigate this process, consult the file found at .github/workflows/tests.yml on the github page. The output of the report can be found in the github actions tab.

A coverage report was made by using coveralls locally, but can be found by clicking the coverage shield on the github page.

(2) Availability

Operating system

This package can be run on any operating system where python can be run (GNU/Linux, Mac OSX, Windows).

Programming language

python 3.7+

Additional system requirements

None

Dependencies

Packages that are required are scipy>=1.6.1, scikit-learn>=1.0.2, tqdm>=4.59.0, numpy>=1.20.1, networkx>=2.5, POT>= 0.8.2.

The Grakel package is optional but is recommended to be installed so an extra suite of graph kernels can be used.

Software location

Archive

Name: Zenodo

Persistent identifier: https://doi.org/10.5281/zenodo.8037197

Licence: MIT

Publisher: Ragnar Leví Guðmundarson

Version published: 0.0.6

Date published: 14/06/23

Code repository

Name: Github

Persistent identifier: https://github.com/ragnarlevi/GTST

Licence: MIT

Date published: 11/11/22

Language

English

(3) Reuse potential

The code was originally used in the paper [] to compare pairwise asset return process relationships to study and understand risk and return in portfolio management practice. This allows one to statistically test for significance of any detected differences in portfolio diversification between any portfolio investment strategy when applying differing investment screening criteria or optimal investment strategies. We further remark that this package can further be used in other fields other than portfolio comparison. For example, the package allows, in a straight forwards manner, to compare different community structures, to detect changes in communities, to detect change-point events, to test for differences in traffic networks, and to compare ego-networks of entities.

Every function has docstrings to ensure clarity. Examples for all graph kernels, some Grakel kernels and estimators can be found on the Github page, along with an example of the graph estimation functionality. GTST is released under the MIT license and welcomes any contributions. We encourage users to submit feedback using GitHub issue tracker, or by emailing Ragnar.

List of Contributors

RLG wrote the software package, developed the methodology, wrote the paper, and tested the software. GWP guided the problem formulations and solutions conceptually.

[B1] Bassett, DS, Bullmore E, Verchinski BA, Mattay VS, Weinberger DR, Meyer-Lindenberg A. Hierarchical organization of human cortical networks in health and schizophrenia. The Journal of Neuroscience. 2008; 28(37): 9239–9248. DOI: https://doi.org/10.1523/JNEUROSCI.1929-08.2008

[B2] Buitinck L, Louppe G, Blondel M, Pedregosa F, Mueller A, Grisel O, Niculae V, Prettenhofer P, Gramfort A, Grobler J, Layton R, VanderPlas J, Joly A, Holt B, Varoquaux G. API design for machine learning software: experiences from the scikit-learn project. In ECML PKDD Workshop: Languages for Data Mining and Machine Learning. 2013; 108–122. URL: https://scikit-learn.org/.

[B3] Carreno JG, Cifuentes R. Identifying complex core-periphery structures in the interbank market. Journal Of Network Theory In Finance; 2017. DOI: https://doi.org/10.21314/JNTF.2017.035

[B4] Du SS, Hou K, Salakhutdinov RR, Poczos B, Wang R, Xu K. Graph neural tangent kernel: Fusing graph neural networks with graph kernels. In: Wallach H, Larochelle H, Beygelzimer A, dAlché-Buc F, Fox E, Garnett R, (eds.), Advances in Neural Information Processing Systems. 2019; 32: Curran Associates, Inc. URL: https://proceedings.neurips.cc/paper_files/paper/2019/file/663fd3c5144fd10bd5ca6611a9a5b92d-Paper.pdf.

[B5] Friedman J, Hastie T, Tibshirani R. Sparse inverse covariance estimation with the graphical lasso. Biostatistics. 2007; 9(3): 432–441. DOI: https://doi.org/10.1093/biostatistics/kxm045

[B6] Girvan M, Newman MEJ. Community structure in social and biological networks. Proceedings of the National Academy of Sciences. 2002; 99(12): 7821–7826. DOI: https://doi.org/10.1073/pnas.122653799

[B7] Gretton A, Borgwardt KM, Rasch MJ, Schölkopf B, Smola A. A kernel two-sample test. Journal of Machine Learning Research. 2012; 13(25): 723–773. URL: http://jmlr.org/papers/v13/gretton12a.html.

[B8] Gudmundarson RL, Peters G. Assessing portfolio diversification via two-sample graph kernel inference. a case study on the inuence of ESG screening. SSRN Electronic Journal; 2023. DOI: https://doi.org/10.2139/ssrn.4348306

[B9] Hagberg A, Conway D. Networkx: Network analysis with python; 2020. URL: https://networkx.github.io.

[B10] Kang U, Tong H, Sun J. Fast random walk graph kernel. In Proceedings of the 2012 SIAM International Conference on Data Mining. Society for Industrial and Applied Mathematics; 2012. DOI: https://doi.org/10.1137/1.9781611972825.71

[B11] Kriege NM, Johansson FD, Morris C. A survey on graph kernels. 2020; 5. DOI: https://doi.org/10.1007/s41109-019-0195-3

[B12] Lerasle M, Szabo Z, Mathieu T, Lecue G. {MONK} outlier-robust mean embedding estimation by median-of-means. PMLR. 2019; 97: 3782–37. URL: http://proceedings.mlr.press/v97/lerasle19a.html.

[B13] Liu H, Lafferty J, Wasserman L. The nonparanormal: Semiparametric estimation of high dimensional undirected graphs. Journal of Machine Learning Research. 2009; 10(80): 2295–2328. URL: http://jmlr.org/papers/v10/liu09a.html.

[B14] Nikolentzos G, Siglidis G, Vazirgiannis M. Graph kernels: A survey. J. Artif. Int. Res. 2022; 72: 943–1027. DOI: https://doi.org/10.1613/jair.1.13225

[B15] Orzechowski P, Moore JH. Ebic: A scalable biclustering method for large scale data analysis. GECCO 2019 Companion – Proceedings of the 2019 Genetic and Evolutionary Computation Conference Companion. 2019; 31–32. DOI: https://doi.org/10.1145/3319619.3326762

[B16] Schölkopf B, Herbrich R, Smola AJ. A generalized representer theorem. In: Helmbold D, Williamson B (eds.), Computational Learning Theory. Springer Berlin Heidelberg, Berlin, Heidelberg. 2001; 416–426. DOI: https://doi.org/10.1007/3-540-44581-1_27

[B17] Siglidis G, Nikolentzos G, Limnios S, Giatsidis C, Skianis K, Vazirgiannis M. Grakel: A graph kernel library in python. Journal of Machine Learning Research. 2020; 21(54): 1–5. URL: https://github.com/ysig/GraKeL.

[B18] Togninalli M, Ghisu E, Llinares-López F, Rieck B, Borgwardt K. Wasserstein Weisfeiler-Lehman Graph Kernels. NIPS ‘19, Curran Associates Inc., Red Hook, NY, USA; 2019. URL: https://proceedings.neurips.cc/paper/2019/file/73fed7fd472e502d8908794430511f4d-Paper.pdf.

[B19] Yanardag P, Vishwanathan SVN. Deep graph kernels. Association for Computing Machinery. 2015; 1365–1374. DOI: https://doi.org/10.1145/2783258.278341.

Journal of Open Research Software

Software Metapapers