Year: 2019

Abstract:
Graph-based semi-supervised learning (SSL) assigns labels to initially unlabelled vertices in a graph. Graph neural networks (GNNs), esp. graph convolutional networks (GCNs), inspired the current-state-of-the art models for graph-based SSL problems. GCNs inherently assume that the labels of interest are numerical or categorical variables. However, in many real-world applications such as coauthorship networks, recommendation networks, etc., vertex labels can be naturally represented by probability distributions or histograms. Moreover, real-world network datasets have complex relationships going beyond pairwise associations. These relationships can be modelled naturally and flexibly by hypergraphs. In this paper, we explore GNNs for graph-based SSL of histograms. Motivated by complex relationships (those going beyond pairwise) in real-world networks, we propose a novel method for directed hypergraphs. Our work builds upon existing works on graph-based SSL of histograms derived from the theory of optimal transportation. A key contribution of this paper is to establish generalisation error bounds for a one-layer GNN within the framework of algorithmic stability. We also demonstrate our proposed methods’ effectiveness through detailed experimentation on real-world data. We have made the code available.

heading: 1 INTRODUCTION
text: In the last decade, deep learning models have been successfully embraced in many different fields and proved to achieve unprecedented performance on a vast range of applications Krizhevsky et al. (2012); Goodfellow et al. (2014); Bahdanau et al. (2015); LeCun et al. (2015). Graph Convolutional Network (GCN) Kipf & Welling (2017) was recently proposed as an adaptation of a particular deep learning model (i.e., convolutional neural networks Lecun et al. (1998)) to enable handling graph-structured data. GCN was shown to be, in particular, effective in semi-supervised learning on attributed graphs. GCNs have inspired the current state-of-the art models for graph-based SSL Wu et al. (2019a); Veličković et al. (2018); Vashishth et al. (2019). GCNs inherently assume that the labels of interest are numerical or categorical variables.
However, in many real-world applications such as co-authorship networks, recommendation networks, etc., vertex labels can be naturally represented by probability distributions or histograms. Moreover, these real-world network datasets have complex relationships going beyond pairwise associations. Such relationships can be modelled naturally and flexibly by hypergraphs. Moreover, hypergraphs can encode additional relationships with directions as illustrated in Figure 1 and these hypergraphs are directed hypergraphs Gallo et al. (1993).
Inspired by a prior work Solomon et al. (2014) that generalised label propagation to graph-based soft SSL setting and motivated by the fact that GNNs have inspired state-of-the-art models for traditional graph-based SSL, we make the following contributions.
• We explore GNNs for soft SSL in which vertex labels are probability distributions. Motivated by real-world applications, we propose DHN (Directed Hypergraph Network), a novel method for directed hypergraphs. DHN can be applied for soft-SSL using existing tools from optimal transportation (Section 3).
• We provide generalisation error bounds for a one-layer GNN within the framework of algorithmic stability. We establish that such models, which use filters with bounded eigenvalues
independent of graph size, can satisfy the strong notion of uniform stability and thus is generalisable. In particular, the algorithmic stability of a one-layer GNN depends on the largest absolute eigenvalue of the graph convolution filter (Section 4).
• We demonstrate DHN’s effectiveness through detailed experimentation on real-world data. In particular, we demonstrate superiority over state-of-the-art hypergraph-based neural networks. We provide new empirical benchmarks for soft-SSL on directed hypergraphs and make the code available to foster reproducible research (Section 5).

heading: 2 RELATED WORK
text: Graph-based deep learning: Geometric deep learning Bronstein et al. (2017) is an umbrella phrase for emerging techniques attempting to generalise (structured) deep neural networks to nonEuclidean domains such as graphs and manifolds. GCN Kipf & Welling (2017) and their various extensions are the current state-of-the art for graph-based SSL Wu et al. (2019a); Veličković et al. (2018); Monti et al. (2018); Vashishth et al. (2019); Ma et al. (2019); Qu et al. (2019) and graphbased unsupervised learning Hamilton et al. (2017a); Veličković et al. (2019) problems. GCNs have also been applied for semi-supervised graph classification Li et al. (2019).
The reader is referred to a comprehensive literature review Bronstein et al. (2017) and extensive surveys Hamilton et al. (2017b); Battaglia et al. (2018); Zhang et al. (2018); Wu et al. (2019b); Sun et al. (2018); Zhou et al. (2018) on this topic. Recently, graph-based deep models (also messagepassing neural networks Gilmer et al. (2017)) have been analysed theoretically Dehmamy et al. (2019); Ying et al. (2019); Knyazev et al. (2019); Chen et al. (2019b); Maron et al. (2019); Morris et al. (2019); Xu et al. (2019b); Chen et al. (2019a); Kawamoto (2018); Xu et al. (2018); Chen et al. (2018). We note a work on stability and generalisation bounds of GCNs Verma & Zhang (2019).
Learning on hypergraphs: Hypergraph is a combinatorial structure consisting of vertices and hyperedges, where each hyperedge is allowed to connect any number of vertices, thus generalizing graphs. This additional flexibility facilitates the capture of higher order interactions among objects; applications have been found in many fields such as computer vision Govindu (2005), network clustering Demir et al. (2008), folksonomies Ghoshal et al. (2009), cellular networks Klamt et al. (2009), and community detection Chien et al. (2018).
The seminal work on hypergraphs Zhou et al. (2007) introduced the popular Agarwal et al. (2006); Feng et al. (2018; 2019) clique expansion of a hypergraph. Hypergraph neural networks (HGNN) Feng et al. (2019) use the clique expansion while HyperGCN Yadati et al. (2019) uses the mediator-
based Laplacian to extend GCNs to hypergraphs. Another line of work uses the mathematically appealing tensor methods Shashua et al. (2006); Bulò & Pelillo (2009); Kolda & Bader (2009) but they are limited to uniform hypergraphs. Recent developments work for arbitrary hypergraphs and fully exploit the hypergraph structure Hein et al. (2013); Zhang et al. (2017); Chan & Liang (2018); Wendler et al. (2019). These developments are motivated from the spectral theory of hypergraphs which is an active area of research Li & Milenkovic (2017); Chien et al. (2019); Chitra & Raphael (2019); Li & Milenkovic (2018b); Li et al. (2018a); Li & Milenkovic (2018a)
Graph-based soft SSL: Researchers have shown that using unlabelled data during training can improve label prediction significantly Chapelle et al. (2010); Zhu et al. (2009); Subramanya & Talukdar (2014); Yang et al. (2016). While most methods assume that labels of interest are numerical or categorical variables, other works “soften” this assumption and handle “soft labels” such as histograms Corduneanu & Jaakkola (2005); Tsuda (2005). One way of propagating histograms is to minimise the Kullback-Leibler (KL) divergence Subramanya & Bilmes (2011). Recent studies have replaced the metric-agnostic KL divergence with metric-aware Wasserstein distance (interactions between histogram bins) for graphs Solomon et al. (2014) and hypergraphs Gao et al. (2019).
Embeddings in Wasserstein space: There exist at least a couple of recent works that embed Gaussian distributions in the Wasserstein space Muzellec & Cuturi (2018); Zhu et al. (2018). Inspired by a recent work Frogner et al. (2019), in this work, we focus on embedding input data as a discrete probability distrirbution on a fixed support set. The Wasserstein distance and its gradient require the solution of a linear program Villani (2008) and are costly to compute Peyré & Cuturi (2019). A popular efficient approximation is the Sinkhorn divergence Cuturi (2013) in which the underlying problem is regularised and is computed efficiently by a fixed-point iteration. Recent works have shown that it is suitable for gradient-based optimisation through automatic differentiation Frogner et al. (2015); Genevay et al. (2018); Frogner et al. (2019). A couple of recent works compute Wasserstein distance between graph pairs Xu et al. (2019a); Vayer et al. (2019).

heading: 3 METHOD
text: In this section, we first describe soft SSL on directed hypergraphs and then propose DHN (Directed Hypergraph Network) for the problem.

heading: 3.1 DIRECTED HYPERGRAPH
text: A directed hypergraph Gallo et al. (1993) is an ordered pairH = (V,Ed) where V = {v1, · · · , vn} is a set of n vertices and
Ed = {(t1, h1), · · · , (tm, hm)} ⊆ 2V × 2V
is a set of m directed hyperedges. Each element in Ed is an ordered pair (t, h) where t ⊆ V is the tail and h ⊆ V is the head with t 6= ∅, h 6= ∅. Denote the set of all undirected hyperedges by E.
E = ⋃
(t,h)∈Ed
( t ∪ h ) .
Denote I ∈ {0, 1}|V |×|E| to be the incidence matrix of E i.e. I(v, e) = 1 if v ∈ e and 0 otherwise.

heading: 3.2 SOFT SSL ON DIRECTED HYPERGRAPHS
text: We consider the problem of predicting probability distributions for the vertices inH = (V,Ed) given a typically small subset Vk ⊆ V of vertices with known distributions. In this work, we are concerned with discrete distributions modelled on a metric space i.e. an ordered pair (M,C) in which M is a set and C is the cost function (metric) associated with the set. Furthermore, we assume that we are provided with a feature matrix, XV ∈ Rn×DV , in which each vertex v ∈ V is represented by a DV -dimensional feature vector xv (here n = |V |). We are also provided with a hyperedge feature matrix XE ∈ Rm×DE with xe, e ∈ E as DE-dimensional feature representations (here m = |Ed|).
Our objective is to learn a labelling function Z = h ( H, XV , XE ) that maps each vertex to a probability distribution in the space of discrete probability distributions PF (M) on F atoms (F is number of histogram bins) defined on the metric space (M,C). The cost function C can be represented by a non-negative symmetric matrix of size F × F . Note that each row of Z ∈ [0, 1]n×F maps each vertex v ∈ V to a probability distribution Zv ∈ [0, 1]F .
The function h is going to be trained on a supervised loss, L ,w.r.t to the vertices in Vk so that the trained h can be used to predict distributions of all the vertices in V \ Vk. We now give an example application and then the details of the labelling function h followed by the supervised loss L.
Example application: Predicting topic distributions of authors in co-authorship networks can be posed as a soft SSL problem on directed hypergraphs. V represents the set of authors, E the set of all collaborations (documents), Ed the citation relationships among the documents, F the number of possible research interests of authors (Machine Learning, Theoretical Computer Science, etc.), XV and XE any available features on the authors and documents respectively (e.g. text attributes).

heading: 3.3 DHN (DIRECTED HYPERGRAPH NETWORK)
text: Hypergraphs contain hyperedges in which relationships can go beyond pairwise and hence are challenging to deal with. A flexible way to embed vertices of a hypergraph is to “approximate” the hypergraph by a suitable graph and then apply traditional graph-based methods on the vertices. Two notable candidates of h are Hypergraph neural network (HGNN) Feng et al. (2019) and Hypergraph Convolutional Network (HyperGCN) Yadati et al. (2019). HGNN uses the clique expansion of the hypergraph Zhou et al. (2007) while HyperGCN uses the mediator-based Laplacian Chan & Liang (2018) to approximate the input hypergraph. However, they are restricted to undirected hyperedges and also cannot exploit the hyperedge feature matrix XE .
A key idea of our approach is to treat each hyperedge e ∈ E as a vertex of the graph G = (E,Ed). We then pass G through a graph neural network to obtain HE = fGNN (G, XE) so that the initial features, XE , are refined to HE . We then propose the layer-wise propagation rule of DHN as:
H (t+1) V = σ
([ H
(t) V , I ·H (t) E ·Θ
(t) ]) , t = 0, · · · , τ − 1 (1)
where [·, ·] denotes concatenation, t is the time step, I is the incidence matrix, H(t+1)E = σ1 ( ITH (t) V ) for t = 1, · · · , τ − 1, H(0)E = fGNN (G, XE), σ and σ1 are non-linear activation
functions, and τ is the total number of time steps with H(0)V = XV . Note that the labelling function Z = h(H, XV , XE) = softmax ( HτV ) where softmax is applied row-wise.

heading: 3.4 THE SUPERVISED LOSS L
text: A crucial observation here is that because of the softmax layer, the output of h is (already) inherently a probability distribution. For each vertex v ∈ Vk, the predicted distribution Zv and the (known) true distribution Yv must be “close” to each other. A natural way to compare probability distributions is to use the KL-divergence between Yv and Zv . However, KL-divergence cannot exploit the metric space (M,C) and suffers from stability issues Chen et al. (2016). In this work, we use the more stable Wasserstein distance to exploit the metric space Gao et al. (2019).
L = ∑ v∈Vk Wp ( Zv, Yv ) where Wp(µ, ν) = ( inf π∈Π(µ,ν) ∫ M×M C(x1, x2) pdπ(x1, x2) ) 1 p .
(2)
For discrete distributions,Wp is the solution of a linear program. For practical purposes, we compute the regularised distance using the Sinkhorn algorithm. Please see Appendix A.1 for more details.
Optimisation: We call DHN optimised with the Wasserstein loss as Soft-DHN. All parameters are learned using stochastic gradient descent (SGD). Please see Appendix for time complexity.

heading: 4 THEORETICAL ANALYSIS: GENERALISATION ERROR BOUND
text: In this section, we establish generalisation error bounds for a one-layer GNN by extending the results of a traditional GCN Verma & Zhang (2019) to the soft SSL setting with Wasserstein loss. The main novelty is to generalise the error bounds to the learning problem “valued in the Wasserstein space”. The main challenge is that the Wasserstein space is an abstract metric space without linear structure.
The section is organised as follows. We first introduce all the notations needed (ego-graph view, semi-supervised learning setting, etc.). We then give single layer and SGD bounds using the nota-
tions. We finally give the main result (Theorem 1) which states that a GNN trained with Wasserstein loss has the same generalisation error bound as the traditional GCN (trained with cross entropy).
Let G = (V,E) be a connected graph with |V | = n vertices. We consider GNN of a single layer f (X,Θ) = σ (KXΘ) (3)
where X ∈ Rn×d is the feature matrix (n is the number of vertices in a graph, d is the dimension of the feature vectors), K = g (LG) is a graph filter (typically symmetrically normalised adjacency with self loops, and LG ∈ Rn×n is the graph Laplacian), and Θ ∈ Rd×F is the set of parameters. We note that our proposed DHN falls under this formulation in special circumstances. Specifically if the non-linearity σ1 in Equation 1 is removed we get the kernel K = IIT (also known as the clique expansion of the hypergraph Zhou et al. (2007).) The non-linearity σ in Equation 3 is the softmax function acting on each row of the product g (LG)XΘ ∈ Rn×F ; the output is of dimension n× F , where each output row is a discrete probability distribution, i.e.,
f (X,Θ) ≥ 0 and f (X,Θ)1F = 1n where 1F = (1, . . . , 1)
F ∈ RF , and similarly for 1n. Without loss of generality, we assume d = 1. Note that in order for the output to be nontrivial probability distributions, we must assume F > 1.
We adopt an ego-graph view (Verma & Zhang, 2019, formula (2)) to simplify our discussion for local behavior of the soft GCN at a particular vertex. Whenever no confusion arises, we identify a vertices x and χ in the graphG with their respectiveD-dimensional feature vectors. Thus the output of f at x ∈ V is
f (x,Θ) = σ  ∑ χ∈N (x) KxχχΘ  = σ  ∑ χ∈N (x) Kxχχ  ·Θ 
where N (x) denotes for the one-hop neighborhood of x with respect to the adjacency relation defined by matrix K, and Kxχ ∈ R stands for the entry in K ∈ Rn×n that describes the adjacency relation between vertices x and χ. Let Ex := ∑ χ∈N (x)Kxχχ ∈ R so that f (x,Θ) = σ (Ex ·Θ) .
We consider the supervised learning setting, and learn GNN from the training set {zi = (xi, yi) , i = 1, . . . ,m} sampled i.i.d. from the product space V × PF with respect to probability distribution D on this product space, where PF is the space of discrete probability distributions on F atoms. The output of softmax lies in PF , which is a convex cone. For any new data z = (x, y) ∼ D, we evaluate the performance of GNN f using a Wasserstein cost
` (f (·,Θ) , z) = ` (f (·,Θ) , (x, y)) = W (f (x,Θ) , y) . Here the Wasserstein cost is defined with respect to a cost function penalizing moving masses across bins. Since we are working only with histograms in GNN, we shall use a cost function C ∈ RF×F that is defined for pairs of histogram bins. The transport problem is a linear program with z = f(x,Θ) as in Equation 10. Please see Appendix Sections A.3.1 and A.3.3 for the assumptions made and the notations used for Algorithmic stability. We now derive single layer and SGD bounds.
4.0.1 SINGLE LAYER BOUND
By the triangle inequality,
|W1 (f (x,ΘS) , y)−W1 (f (x,ΘS′) , y)| ≤W1 (f (x,ΘS) , f (x,ΘS′)) . By (Villani, 2008, Theorem 6.13), recall that the diameter of the support is D, we have
W1 (f (x,ΘS) , f (x,ΘS′)) ≤ D ‖f (x,ΘS)− f (x,ΘS′)‖TV where ‖·‖TV is the total variation distance, which by definition is
‖f (x,ΘS)− f (x,ΘS′)‖TV = 1
2 F∑ i=1 |[f (x,ΘS)]i − [f (x,ΘS′)]i| = 1 2 ‖σ (Ex ·ΘS)− σ (Ex ·ΘS′)‖1
where ‖·‖1 is the L1-distance on RF . Since the softmax function is Lipschitz continuous, we have
‖f (x,ΘS)− f (x,ΘS′)‖TV ≤ Lσ 2 |Ex| · ‖ΘS −ΘS′‖1
and thus |EA[W1(f(x,ΘS),y)−W1(f(x,ΘS′ ),y)]|≤ LσD 2 supx∈V |Ex|·EA‖ΘS−ΘS′‖1= LσD 2 gλEA‖ΘS−ΘS′‖1 (4) where we used notation gλ := supx∈V |Ex| as defined in Verma & Zhang (2019), which is known to be upper bounded by λmaxG , the spectrum of the graph Laplacian LG with largest absolute value. 4.0.2 SGD BOUND
It now remains to bound EA ‖ΘS −ΘS′‖1 resulting from the SGD iterations. The main technical challenge, as noted before, is to generalise the results to the Wasserstein space which is an abstract metric space without linear structure. Specifically, we have to modify the “gradient” in the Wasserstein space as the straightforward version Verma & Zhang (2019) does not satisfy the Lipschitz condition required in the algorithmic stability framework. To the best of our knowledge this modification is not seen in existing literature and can be thought of as a generalisation of the “gradient clipping” operation Hardt et al. (2016). The entire proof is in the Appendix A.3.4. We state the main result here.
EA [ ‖ΘS,T −ΘS′,T ‖1 ] ≤ 2ηFgλD
m
T∑ t=1 ( 1 + 3 2 ηDLσg 2 λ )t−1 . (5)
Combining equation 4 and equation 23 gives us
|EA [W1 (f (x,ΘS) , y)−W1 (f (x,ΘS′) , y)]| ≤ ηFLσg
2 λD 2
m
T∑ t=1 ( 1 + 3 2 ηDLσg 2 λ )t−1 . (6)
Therefore, we actually have the uniform algorithmic stability equation 14 holds with
βm = ηFLσg
2 λD 2
2m
T∑ t=1 ( 1 + 3 2 ηDLσg 2 λ )t−1 . (7)
Note that here βm = O( 1m ) (needed to obtain a tight generalisation bound).

heading: 4.1 PUTTING EVERYTHING TOGETHER
text: Lemma 1 (Verma & Zhang, 2019, Lemma 4): Let λmaxG be the maximum absolute eigenvalue of LG. Let Gx be the ego-graph of a vertex x ∈ V with corresponding maximum absolute eigenvalue λmaxGx . Then the following eigenvalue (singular value) bound holds ∀x ∈ V ,
λmaxGx ≤ λ max G (8)
Lemma 2 (Verma & Zhang, 2019, Theorem 2) A uniform stable randomised algorithm (AS , βm) with a bounded loss function 0 ≤ `(AS ,y) ≤ B, satisfies the following generalisation bound with probability at least 1− δ, over the random draw of S,y with δ ∈ (0, 1),
EA [ R(AS)−Remp(AS) ] ≤ 2βm + (4mβm +B) √ log 1δ 2m
(9)
where R(AS) is the generalisation error/risk and Remp(AS) is the empirical error. Please see Appendix A.3.2 for definitions. Finally, by Equation 8, and Equation 9, our result i.e. Equation 7 immediately gives the following theorem:
Theorem 1 Let AS be a one-layer GNN algorithm (of 3) equipped with the graph convolutional filter g(LG) and trained on a dataset S for T iterations. Let the loss and activation functions be Lipschitz-continuous and smooth. Then the following expected generalisation gap holds with probability at least 1− δ, δ ∈ {0, 1}:
ESGD [ R(AS)−Remp(AS) ] ≤ 1 m O ( (λmaxG ) 2T ) + ( O ( (λmaxG ) 2T ) +B )√ log 1δ 2m
where the expectation ESGD is taken over the randomness inherent in SGD, m is the no. training samples, and B is a constant which depends on the loss function. Our theorem states that GNN trained withe Wasserstein loss enjoys the same generalisation error bound of the traditional GCN (trained with cross entropy). We now discuss experiments.

heading: 5 EXPERIMENTS
text: To demonstrate the effectiveness of our proposed DHN, we conducted experiments on 5 real-world directed hypergraphs. Four of them are Co-authorship datasets and one of them is a recommendation dataset. Table 1 shows the statistics of the datasets. For more details on the construction of the datasets, please see Appendix A.6.

heading: 5.1 EXPERIMENTAL SETUP
text: Inspired by the experimental setups of prior related works Kipf & Welling (2017); Liao et al. (2019), we tune hyperparameters using the Cora citation network dataset alone and use the optimal hyperparameters for all the other datasets. We hyperparameterise the cost matrix (base metric of the Wasserstein distance) as follows:
C =  1 η η . . . η η η 1 η . . . η η ... ... ... . . . ...
... η η η . . . η 1  The cost matrix C is an F ×F matrix (F is the number of histogram bins) with ones on the diagonal and a hyperparameter η elsewhere. We could have used a matrix of all ηs. But it is no different from a matrix of all ones from the optimisation perspective and so we used the above more general matrix. Details of hyperparameter tuning and optimal hyperparameters are in Appendix A.5.

heading: 5.2 BASELINES
text: We used both Wasserstein distance and KL divergence to train different models. As already noted, we used the Sinkhorn algorithm to compute the (regularised) Wasserstein distance. Please see Appendix A.1 for more details. We compared DHN with the following baselines:
• KL-MLP: We used a simple multi-layer perceptron (MLP) on the features of the vertices and trained it using KL-divergence
• OT-MLP: We trained another MLP with the Wasserstein distance as the loss function. Note that this baseline and the previous baseline do not use the structure (graph / hypergraph)
• KLR-MLP: We regularised an MLP with explicit KL-divergence-based regularisation that uses the structure (graph / hypergraph) Subramanya & Bilmes (2011).
• OTR-MLP: We regularised an MLP with explicit Wasserestein-distance-based regularisation that uses the structure (graph / hypergraph) Solomon et al. (2014). For hypergraphs we used the clique expansion of the hypergraph Gao et al. (2019).
• KL-HGNN / KL-HyperGCN: We trained the different GCN-based methods on hypergraphs with KL divergence loss function on the labelled vertices.
• Soft-HGNN / Soft-HyperGCN: We trained the different GCN-based methods on hypergraphs with the Wasserstein distance as the loss function.
Metric for comparison: We use the mean squared error (MSE) between true and predicted distributions on the test set of vertices. Table 2 shows MSEs on the test split for all the three datasets.

heading: 5.3 DISCUSSION
text: We used a simple one-layer architecture for our proposed DHN and a 2-hop simplified GCN Wu et al. (2019a) as the GNN model on the graph G = (E,Ed) i.e.
Z = softmax(I ·HE ·Θ1), HE = A2XEΘ2
whereA is the symmetically normalised adjacency (with self loops) of the graph G. We demonstrate that this simple model is effective enough through an ablation study in Appendix Table 4. Our results demonstrate strong performances across all the datasets esp. on the co-authorships networks.
Specifically, we observe that Soft models (that use the Wasserstein loss) are almost always superior to their counterparts that use the KL divergence as the loss function. This is because the Soft models can exploit the distance matrix C while KL-divergence does not.
Moreover, our proposed DHN outperforms two strong hypergraph baselines viz. HGNN Feng et al. (2019) and HyperGCN Yadati et al. (2019). We believe this is because they do not exploit the rich structural information in the directed hyperedges (connections among hyperedges) while our proposed DHN does exploit them.
We also experimented on standard graph benchmark node-classification datasets such as Cora, Citeseer, and Pubmed by treating the class label as one-hot probability distribution. We used the Soft variants of GCN Kipf & Welling (2017), Simple GCN Wu et al. (2019a), and GAT Veličković et al. (2018). We achieved competitive results as shown in Appendix Table 3

heading: 6 CONCLUSION
text: We have proposed DHN, a novel method for soft SSL on directed hypergraphs. DHN can effectively propagate histograms to unknown vertices by integrating vertex features, directed hyperedges and undirected hypergraph structure. As a key contribution, we have established generalisation bounds for DHN within the framework of algorithmic stability. We have also demonstrated DHN’s effectiveness through detailed experimentation on real-world hypergraph datasets.

heading: A APPENDIX
text: The appendix is organised as follows:
1. Regularised Wasserstein and the Sinkhorn algorithm 2. Time complexity of the proposed DHN 3. Proofs and more notations used for theoretical analysis 4. Additional experiments (Graph-based soft SSL and Ablation Study) 5. Details of hyperparameters 6. Sources of real-world datasets

heading: A.1 REGULARISED WASSERSTEIN AND THE SINKHORN ALGORITHM
text: For discrete distributions, Wp of Equation 2 is the solution of a linear program:
Wp (z, y) p = min F∑ i=1 F∑ j=1 Cpijπij
s.t. πij ≥ 0, F∑ j=1 πij = [z]i , F∑ i=1 πij = [y]j ∀1 ≤ i, j ≤ F
(10)
where for an arbitrary element q ∈ PF , [q]i stands for the probability mass in the ith bin.

heading: A.1.1 SINKHORN DIVERGENCE
text: The expensive linear program 10 can be efficiently solved by entropic regularisation Cuturi (2013):
Wλp (z, y) p = min F∑ i=1 F∑ j=1 Cpijπij + λ · Tr ( π ( log π − 11T )T)
(11)
s.t. πij ≥ 0, F∑ j=1 πij = [z]i , F∑ i=1 πij = [y]j ∀1 ≤ i, j ≤ F
where log(·) is applied element-wise and λ ≥ 0 is a hyperparameter. The optimal solution π∗ for λ > 0 takes the following form.
π∗ = diag(r) · exp ( − C p
λ
) · diag(c)
where diag(z) is a diagonal matrix with the components of z in the diagonal places.
Sinkhorn algorithm Cuturi (2013): We optimise Equation 11 for r and c via matrix balancing, i.e., start with an initialK := exp(−C p
λ ) and alternately ensure the marginal constraints are satisfied until convergence:
r ← z./(Kc) c← y./(KT r)
where ./ is element-wise division. We use the above efficient algorithm in our experiments.

heading: A.2 TIME COMPLEXITY OF DHN
text: We consider the problem of predicting probability distributions for the vertices in H = (V,Ed) given a typically small subset Vk ⊆ V of vertices with known distributions. In this work, we are concerned with discrete distributions modelled on a complete separable metric space (M,C). Furthermore, we assume that we are provided with a feature matrix, XV ∈ Rn×DV , in which each vertex v ∈ V is represented by a Dv-dimensional feature vector xv (here n = |V |). We are also provided with a hyperedge feature matrixXE ∈ Rm×DE with xe, e ∈ E asDE-dimensional feature representations. Let T be the time taken for the Sinkhorn algorithm on all the vertices with known distributions. Note that T can be approximated in near-linear time Altschuler et al. (2017). Further, let T be the total number of epochs of training. Define
N := ∑ e∈E |e|, Nc := ∑ e∈E |e|C2
The time complexity oa one-layer DHN isO ( |Ed| ·DE ·h1+ |E| ·DV ·h2 ) where h1 is the number
of hidden units of the GNN layer and h2 is the number of output channels.

heading: A.3 THEORETICAL ANALYSIS: PROOFS
text: This section is organised as follows:
1. Assumptions/notations
2. Definitions of generalisation and empirical errors
3. Framework of Algorithmic Stability
4. Proof of SGD Bound

heading: A.3.1 ASSUMPTIONS/NOTATIONS
text: To avoid unnecessary technical complications, assume the histogram admits a geometric realisation over the one-dimensional Euclidean space, such that the ith bin is placed at location bi ∈ R, and set
Cij := |bi − bj | , ∀1 ≤ i, j ≤ F.
Without loss of generality, we assume b1 ≤ b2 ≤ · · · ≤ bF , and write hi := bi+1 − bi ≥ 0 for all i = 1, . . . , F −1. Denote the diameter of the support by D := max1≤i,j≤F |bi − bj | = bF − b1.We take the Wasserstein cost as the Wasserstein-1 distance: W (µ, ν) := W1 (µ, ν) . In this particular one-dimensional setting, we have a particularly simple form for the cost function:
W1 (µ, ν) = ∫ 1 0 ∣∣F−1µ (s)− F−1ν (s)∣∣ds = ∫ ∞ −∞ |Fµ (t)− Fν (t)| dt (12)
where Fµ : R → [0, 1], Fν : R → [0, 1] are the cumulative distribution functions of µ, ν, respectively; F−1µ , F −1 ν are the generalized inverses of Fµ and Fν , respectively, defined as (similar for F−1ν ) F−1µ (t) := inf {b ∈ R : Fµ (b) > t} , ∀t ∈ [0, 1] . (13)
This characterisation is seen in any standard literature on optimal transport, e.g., (Villani, 2003, Theorem 2.18).

heading: A.3.2 DEFINITIONS: GENERALISATION AND EMPIRICAL ERRROS
text: Let the learning algorithm AS on a dataset S be a function from ζm to (Y)X . where X is the input Hilber space, Y is the output Hilbert space, ζ = X × Y . The training set of datapoints, labels is S = {z1 = (x1, y1), · · · , zk = (xk, yk))}. Let the loss function be ` : ζm × ζ → R. Then the generalisation error or risk R(AS) is deined as
R(AS) := E [ `(AS , z) ] = ∫ `(AS , z)p(z)dz
where p(z) is the probability of seeing the sample z ∈ S. The empirical error, on the other hand, is defined as
Remp(AS) := 1
k k∑ j=1 `(AS , zj)

heading: A.3.3 ALGORITHMIC STABILITY
text: Denote S ⊂ V × PF for an arbitrary training data set sampled with respect to distribution D. A learning algorithm for GNN, denoted as A, maps a training set S to a trained GNN f (·,ΘS) : x 7→ f (x,ΘS). Let S′ be another training data set that differs from S by exactly one data. Our goal in this section is to establish a uniform stability
sup S⊂V×PF
(x,y)∈V×PF
|EA [` (f (·,ΘS) , (x, y))]− EA [` (f (·,ΘS′) , (x, y))]|
= sup S⊂V×PF
(x,y)∈V×PF
|EA [W1 (f (x,ΘS) , y)]− EA [W1 (f (x,ΘS′) , y)]|
= sup S⊂V×PF
(x,y)∈V×PF
|EA [W1 (f (x,ΘS) , y)−W1 (f (x,ΘS′) , y)]| ≤ 2βm (14)
with βm = O (1/m), which in turn can be used to establish the generalisation error bound of the form (Verma & Zhang, 2019, Theorem 1), following the framework Hardt et al. (2016) for SGDbased learning algorithms. Following the arguments in Verma & Zhang (2019), this boils down to checking a few Lipschitz properties for the cost function.

heading: A.3.4 SGD BOUND: PROOF
text: Note that
W1 (µ, ν) = ∫ 1 0 ∣∣F−1µ (s)− F−1ν (s)∣∣ds = ∫ ∞ −∞ |Fµ (t)− Fν (t)| dt (15)
and from the single-layer bound
|EA[W1(f(x,ΘS),y)−W1(f(x,ΘS′ ),y)]|≤ LσD 2 supx∈V |Ex|·EA‖ΘS−ΘS′‖1= LσD 2 gλEA‖ΘS−ΘS′‖1 (16)
Uniform stability:
sup S⊂V×PF
(x,y)∈V×PF
|EA [` (f (·,ΘS) , (x, y))]− EA [` (f (·,ΘS′) , (x, y))]|
= sup S⊂V×PF
(x,y)∈V×PF
|EA [W1 (f (x,ΘS) , y)]− EA [W1 (f (x,ΘS′) , y)]|
= sup S⊂V×PF
(x,y)∈V×PF
|EA [W1 (f (x,ΘS) , y)−W1 (f (x,ΘS′) , y)]| ≤ 2βm (17)
It now remains to bound EA ‖ΘS −ΘS′‖1 resulting from the SGD iterations. Given training set S, applying SGD to GNN amounts to performing the updates
ΘS,t+1 = ΘS,t − η∇Θ` (f (·,Θ) , (xit , yit)) = ΘS,t − η∇ΘW1 (f (xit ,ΘS,t) , yit)
where η > 0 is the learning rate and zit = (xit , yit) are random data i.i.d. uniformly sampled from the training set. By the simple formulae equation 15 for one-dimensional optimal transport, we can explicitly write out for any parameter set Θ and data z = (x, y)
W1 (f (x,Θ) , y) = ∫ ∞ −∞ ∣∣Ff(x,Θ) (t)− Fy (t)∣∣ dt = F−1∑ i=1 (xi+1 − xi) ∣∣∣∣∣∣ i∑ j=1 ([f (x,Θ)]i − [y]i) ∣∣∣∣∣∣ =
F−1∑ i=1 hi ∣∣∣∣∣∣ i∑ j=1 ([σ (Ex ·Θ)]i − [y]i) ∣∣∣∣∣∣ where again we used notation [y]i to denote the probability mass of y ∈ PF in the ith bin, for all i = 1, . . . , F . Thus
∂
∂Θk W1 (f (x,Θ) , y) = F−1∑ i=1 hi · sgn  i∑ j=1 ( [σ (Ex ·Θ)]j − [y]j ) · i∑ j=1 ∂ ∂Θk [σ (Ex ·Θ)]j
= F−1∑ i=1 hi · sgn  i∑ j=1 ( [σ (Ex ·Θ)]j − [y]j ) · Ex i∑ j=1 [σ (Ex ·Θ)]j (δjk − [σ (Ex ·Θ)]k) = Ex · [σ (Ex ·Θ)]k F−1∑ i=1 hi · sgn  i∑ j=1 ( [σ (Ex ·Θ)]j − [y]j ) 1− j∑ j=i [σ (Ex ·Θ)]j
 =Ex·[σ(Ex·Θ)]k ∑F−1 i=1 hi·sgn{ ∑i j=1([σ(Ex·Θ)]j−[y]j)}( ∑j j=i([y]j−[σ(Ex·Θ)]j)+ ∑F j=i+1[y]j)
=−Ex·[σ(Ex·Θ)]k·W1(σ(Ex·Θ),y)+Ex·[σ(Ex·Θ)]k ∑F−1 i=1 hi·sgn{ ∑i j=1([σ(Ex·Θ)]j−[y]j)} ∑F j=i+1[y]j
where sgn {·} is the sign function, and δik is the Kronecker delta notation. The second equality used the specific form of the derivative of the softmax function. Unfortunately, this gradient is not Lipschitz continuous due to the sign function in the second term. Nevertheless, if we use a “modified gradient” that drops the second term, i.e., choose to update the parameter Θk in the direction
∂̃
∂Θk W1 (f (x,Θ) , y) := −Ex · [σ (Ex ·Θ)]k ·W1 (σ (Ex ·Θ) , y) (18)
then obviously this new choice of “descent” direction is certainly Lipschitz continuous, as∣∣∣∣∣ ∂̃∂ΘkW1 (f (x,ΘS) , y)− ∂̃∂ΘkW1 (f (x,ΘS′) , y) ∣∣∣∣∣
= |Ex| · |[σ (Ex ·ΘS)]k ·W1 (σ (Ex ·ΘS) , y)− [σ (Ex ·ΘS′)]k ·W1 (σ (Ex ·ΘS′) , y)|
≤ |Ex| · [ |[σ (Ex ·ΘS)]k − [σ (Ex ·ΘS′)]k| ·W1 (σ (Ex ·ΘS) , y)
+ [σ (Ex ·ΘS′)]k · |W1 (σ (Ex ·ΘS) , y)−W1 (σ (Ex ·ΘS′) , y)| ] ≤ |Ex| · (LσD |Ex| · ‖ΘS −ΘS′‖1 +W1 (σ (Ex ·ΘS) , σ (Ex ·ΘS′)))
≤ |Ex| · 3
2 LσD · |Ex| · ‖ΘS −ΘS′‖1 ≤
3 2 LσD · g2λ ‖ΘS −ΘS′‖1 .
Therefore, if we define
∇̃ΘW1 (f (x,Θ) , y) :=
( ∂̃
∂Θ1 W1 (f (x,Θ) , y) , · · · ,
∂̃
∂ΘF W1 (f (x,Θ) , y)
)> (19)
then the stochastic update algorithm1
ΘS,t+1 = ΘS,t − η∇̃ΘW1 (f (xit ,ΘS,t) , yit) (20) 1Note that this is not even a stochastic gradient descent algorithm! the “gradient” involved is a “fake”
gradient — this is the counter-intuitive part.
will still satisfy the generalization bound, due to the Lipschitz continuity∥∥∥∇̃ΘW1 (f (x,ΘS) , y)− ∇̃ΘW1 (f (x,ΘS′) , y)∥∥∥ 1 ≤ 3D 2 Lσg 2 λ ‖ΘS −ΘS′‖1 . (21)
In fact, equation 21 is exactly the GNN analogy of (Verma & Zhang, 2019, Lemma 1) which establishes the stability for the “same sample loss” case. The “different sample loss” analogy, or the lemma in (Verma & Zhang, 2019, Lemma 2), can be trivially obtained by the definition equation 18. In fact, noting that∣∣∣∣∣ ∂̃∂ΘkW1 (f (x,Θ) , y)
∣∣∣∣∣ = |Ex · [σ (Ex ·Θ)]k ·W1 (σ (Ex ·Θ) , y)| ≤ gλD, we easily obtain ∣∣∣∣∣ ∂̃∂ΘkW1 (f (x,ΘS) , y)− ∂̃∂ΘkW1 (f (x′,ΘS′) , y′)
∣∣∣∣∣ ≤ 2gλD and it follows that∥∥∥∇̃ΘW1 (f (x,ΘS) , y)− ∇̃ΘW1 (f (x,ΘS′) , y)∥∥∥
1 ≤ 2FgλD. (22)
Putting together equation 21 and equation 22, we obtain the following analogy of (Verma & Zhang, 2019, Lemma 3): Starting with two training data sets S, S′ that differs by exactly one sample, after each iteration t we have
EA [ ‖ΘS,t+1 −ΘS′,t+1‖1 ] = EA [∥∥∥ΘS,t − η∇̃ΘW1 (f (xt,ΘS,t) , yt)−ΘS′,t + η∇̃ΘW1 (f (x′t,ΘS′,t) , y′t)∥∥∥] ≤EA[‖ΘS,t+1−ΘS′,t+1‖1]+(1− 1m )·η· 3D2 Lσg2λEA[‖ΘS,t−ΘS′,t‖1]+ 1m ·η·2FgλD
≤ ( 1 + 3
2 ηDLσg
2 λ ) EA [ ‖ΘS,t −ΘS′,t‖1 ] + 2ηFgλD
m .
Solving this first-order recursion gives the stability after T random update steps:
EA [ ‖ΘS,T −ΘS′,T ‖1 ] ≤ 2ηFgλD
m
T∑ t=1 ( 1 + 3 2 ηDLσg 2 λ )t−1 . (23)
Combining equation 16 and equation 23 gives us
|EA [W1 (f (x,ΘS) , y)−W1 (f (x,ΘS′) , y)]| ≤ ηFLσg
2 λD 2
m
T∑ t=1 ( 1 + 3 2 ηDLσg 2 λ )t−1 . (24)
Therefore, we actually have the uniform algorithmic stability equation 17 holds with
βm = ηFLσg
2 λD 2
2m
T∑ t=1 ( 1 + 3 2 ηDLσg 2 λ )t−1 . (25)

heading: A.4 ADDITIONAL EXPERIMENTS

heading: A.4.1 ABLATION STUDY
text: In this section, we compared our proposed one-layer DHN against deeper layers and without exploiting directed hyperedges.

heading: A.5 DETAILS OF HYPERPARAMETERS
text: Inspired by the experimental setups of prior related works Kipf & Welling (2017); Liao et al. (2019), we tune hyperparameters using the Cora co-authorship network dataset alone. The optimal hyperparameters are fixed and then used for all the other datasets. Table 5 shows the list of hyperparameters used in the datasets. Prior works Kipf & Welling (2017); Liao et al. (2019) have extensively performed tuning of hyperparameters such as hidden size, learning rate, etc and we fixed their reported optimal hyperparameters. It should be noted that self training and co-training methods Li et al. (2018b) can also be used in case of absence of validation data. We hyperparameterise the cost matrix (base metric of the Wasserstein distance) as follows:
C =  1 η η . . . η η η 1 η . . . η η ... ... ... . . . ...
... η η η . . . η 1  Table 6 shows the best results on the validation split of Cora (with optimal hyperaparameters). The training set had 140 vertices, the validation set 1000 vertices and the rest of the vertices were used to test the models. The results reported are after 200 epochs of training with a seed value of 598.

heading: A.6 SOURCES OF THE REAL-WORLD DATASETS
text: Co-authorship data: All authors co-authoring a paper are in one hyperedge. We used the author data2to get the co-authorship hypergraph for cora. We manually constructed the DBLP dataset from Arnetminer3.
2https://people.cs.umass.edu/ mccallum/data.html 3https://aminer.org/lab-datasets/citation/DBLP-citation-Jan8.tar.bz

heading: A.6.1 CONSTRUCTION OF THE DBLP DATASET
text: We downloaded the entire dblp data from https://aminer.org/lab-datasets/ citation/DBLP-citation-Jan8.tar.bz Tang et al. (2008). The steps for constructing the dblp dataset used in the paper are as follows:
• We defined a set of 5 conference categories (histograms for the SSL task) as “algorithms”, “database”, “datamining”, “intelligence”, and “vision” • For a total of 4304 venues in the entire dblp dataset we took papers from only a sub-
set of venues from https://en.wikipedia.org/wiki/List_of_computer_ science_conferences corresponding to the above 5 conferences • From the venues of the above 5 conference categories, we got 22535 authors publishing at
least two documents for a total of 43413 • We took the abstracts of all these 43413 documents, constructed a dictionary of the most
frequent words (words with frequency more than 100) and this gave us a dictionary size of 1425
• We then extracted the 117215 citation links among these documents

heading: A.6.2 CONSTRUCTION OF THE AMAZON OFFICE PRODUCT DATASET
text: We downloaded the entire Amazon data He & McAuley (2016); McAuley et al. (2015). The steps for constructing the dataset used in the paper are as follows:
• We downloaded the office product ratings subset from the entire dataset • We constructed a hypergraph of items with each hyperedge representing a user connecting
all the items that they bought • We removed hyperedges of size 1 • We connected a pair of hyperedges (bi-directional) if they had more than 1 item in common

heading: A.6.3 CONSTRUCTION OF THE ACM DATASET
text: We downloaded the entire ACM data from https://lfs.aminer.org/lab-datasets/ citation/acm.v9.zip Tang et al. (2008). The steps for curating the dataset in the paper are as follows:
• Based on the number of papers published, we identified the six most popular venues: “Journal of Computational Physics”, “IEEE Transactions on Pattern Analysis and Machine Intelligence”, “Automatica (Journal of IFAC)”, “IEEE Transactions on Information Theory”, “Expert Systems with Apllications: An International Journal”, and “IEEE Transactions on Computers” • We then listed the set of all authors published in these venues (we got a total of 67057
authors). • We finally obtained the citation relationships of all the documents co-authored by these
authors (total number of documents is 25511 and total number of citations is 59884)

heading: A.6.4 DETAILS OF THE ARXIV DATASET
text: We downloaded the entire arXiv dataset Clement et al. (2019) from https://github.com/ mattbierbaum/arxiv-public-datasets/releases/tag/v0.2.0. The steps for curating the dataset in the paper are as follows:
• We removed papers without any authors and got a total of 13, 54, 752 edges • We extracted 67, 28, 683 citation edges among these papers • The total number of authors in these papers is 7, 90, 790

Review summaries:
In this paper, the authors propose a soft semi-supervised learning approach on a hypergraph. On the one hand, the vertex labels should be not only numerical or categorical variables but also probability distributions. On the other hand, hypergraphs provide a much more flexible mean to encode the real-world complicated relationship compared to the essential pairwise association. Specifically, the authors carefully obtain the generalization error bounds for a one-layer graph neural network. The experiments support the theory and provide the empirical verification of the proposed method. The appendix offers plenty of details that are helpful for the reader to understand both the theoretical and practical perspectives of this paper. Also, the code looks good to me. The structure of the provided codebase is clear and well documented.