Title: CONSTANT CURVATURE GRAPH CONVOLUTIONAL NETWORKS

Abstract: Interest has been rising lately towards methods representing data in non-Euclidean spaces, e.g. hyperbolic or spherical, that provide specific inductive biases useful for certain real-world data properties, e.g. scale-free, hierarchical or cyclical. However, the popular graph neural networks are currently limited in modeling data only via Euclidean geometry and associated vector space operations. Here, we bridge this gap by proposing mathematically grounded generalizations of graph convolutional networks (GCN) to (products of) constant curvature spaces. We do this by i) introducing a unified formalism that can interpolate smoothly between all geometries of constant curvature, ii) leveraging gyro-barycentric coordinates that generalize the classic Euclidean concept of the center of mass. Our class of models smoothly recover their Euclidean counterparts when the curvature goes to zero from either side. Empirically, we outperform Euclidean GCNs in the tasks of node classification and distortion minimization for symbolic data exhibiting non-Euclidean behavior, according to their discrete curvature.

1. INTRODUCTION:

Graph Convolutional Networks. The success of convolutional networks and deep learning for image data has inspired generalizations for graphs for which sharing parameters is consistent with the graph geometry. Bruna et al. (2014) ; Henaff et al. (2015) are the pioneers of spectral graph convolutional neural networks in the graph Fourier space using localized spectral filters on graphs. However, in order to reduce the graph-dependency on the Laplacian eigenmodes, Defferrard et al. (2016) approximate the convolutional filters using Chebyshev polynomials leveraging a result of Hammond et al. (2011) . The resulting method (discussed in appendix A) is computationally efficient and superior in terms of accuracy and complexity. Further, Kipf & Welling (2017) simplify this approach by considering first-order polynomials approximations obtaining high scalability. The proposed graph convolutional networks (GCN) is interpolating node embeddings via a symmetrically normalized adjacency matrix, while this weight sharing can be understood as an efficient diffusion-like regularizer. Recent works extend GCNs to achieve state of the art results for link prediction (Zhang & Chen, 2018) , graph classification (Hamilton et al., 2017; Xu et al., 2018) and node classification (Klicpera et al., 2019; Veličković et al., 2018) .

Euclidean geometry in ML. In machine learning (ML), data is most often represented in a Euclidean space for various reasons. First, some data is intrinsically Euclidean, such as positions in 3D space in classical mechanics. Second, intuition is easier in such spaces, as they possess an appealing vectorial structure allowing basic arithmetic and a rich theory of linear algebra. Finally, a lot of quantities of interest such as distances and inner-products are known in closed-form formulae and can be computed very efficiently on the existing hardware. These operations are the basic building blocks for most of today's popular machine learning models. Thus, the powerful simplicity and efficiency of Euclidean geometry has led to numerous methods achieving state-of-the-art on tasks as diverse as machine translation (Bahdanau et al., 2014; Vaswani et al., 2017) , speech recognition (Graves et al., 2013) , image classification (He et al., 2016) or recommender systems (He et al., 2017) .

Riemannian ML. In spite of this success, certain types of data (e.g. hierarchical, scale-free or spherical data) have been shown to be better represented by non-Euclidean geometries (Defferrard et al., 2019; Bronstein et al., 2017; Nickel & Kiela, 2017; Gu et al., 2019) , leading in particular to the rich theories of manifold learning (Roweis & Saul, 2000; Tenenbaum et al., 2000) and information geometry (Amari & Nagaoka, 2007) . The mathematical framework in vigor to manipulate non-Euclidean geometries is known as Riemannian geometry (Spivak, 1979) . Although its theory leads to many strong and elegant results, some of its basic quantities such as the distance function d(•, •) are in general not available in closed-form, which can be prohibitive to many computational methods.

Representational Advantages of Geometries of Constant Curvature. An interesting trade-off between general Riemannian manifolds and the Euclidean space is given by manifolds of constant sectional curvature. They define together what are called hyperbolic (negative curvature), elliptic (positive curvature) and Euclidean (zero curvature) geometries. As discussed below and in appendix B, Euclidean spaces have limitations and suffer from large distortion when embedding certain types of data such as trees, e.g. fig. 1 . In these cases, the hyperbolic and spherical spaces have representational advantages providing a better inductive bias for the respective data.

Figure 1 : Euclidean space quickly "runs out of space" when fitting exponentially volume growing data such as trees. The embedding distance between the red and the green node keeps decreasing as the number of tree nodes increases, meaning that the graph distance (shortest path length) is no longer accurately represented. Details in appendix B.

The hyperbolic space can be intuitively understood as a continuous tree: the volume of a ball grows exponentially with its radius, similarly as how the number of nodes in a binary tree grows exponentially with its depth. Its tree-likeness properties have long been studied mathematically (Gromov, 1987; Hamann, 2017; Ungar, 2008) and it was proven to better embed complex networks (Krioukov et al., 2010) , scale-free graphs and hierarchical data compared to the Euclidean geometry (Cho et al., 2019; Sala et al., 2018; Ganea et al., 2018b; Gu et al., 2019; Nickel & Kiela, 2018; 2017; Tifrea et al., 2019) . Several important tools or methods found their hyperbolic counterparts, such as variational autoencoders (Mathieu et al., 2019; Ovinnikov, 2019) , attention mechanisms (Gulcehre et al., 2018) , matrix multiplications, recurrent units and multinomial logistic regression (Ganea et al., 2018a) .

Similarly, spherical geometry provides benefits for modeling spherical or cyclical data (Defferrard et al., 2019; Matousek, 2013; Davidson et al., 2018; Xu & Durrett, 2018; Gu et al., 2019; Grattarola et al., 2018; Wilson et al., 2014) .

Computational Efficiency of Constant Curvature Spaces (CCS). CCS are some of the few Riemannian manifolds to possess closed-form formulae for geometric quantities of interest in computational methods, i.e. distance, geodesics, exponential map, parallel transport and their gradients. We also leverage here the closed expressions for weighted centroids.

"Linear Algebra" of CCS: Gyrovector Spaces. In order to study the geometry of constant negative curvature in analogy with the Euclidean geometry, Ungar (1999; 2005; 2008; 2016) proposed the elegant non-associative algebraic formalism of gyrovector spaces. Recently, Ganea et al. (2018a) have linked this framework to the Riemannian geometry of the space, also generalizing the building blocks for non-Euclidean deep learning models operating with hyperbolic data representations.

However, it remains unclear how to extend in a principled manner the connection between Riemannian geometry and gyrovector space operations for spaces of constant positive curvature (spherical). By leveraging Euler's formula and complex analysis, we present to our knowledge the first unified gyro framework that smoothly interpolates between geometries of constant curvatures irrespective of their signs. This is possible when working with the Poincaré ball and stereographic spherical projection models of respectively hyperbolic and spherical spaces.

How should one adapt graph neural networks to non-flat geometries of constant curvature?

In this work, we propose constant curvature GCNs to model non-Euclidean data. Node embeddings lie in spaces of constant curvature or product of those instead of a Euclidean space, thus leveraging both the representational power of these geometries and the effectiveness of GCNs.

Concurrent to our work, Chami et al. (2019) ; Liu et al. (2019) propose hyperbolic graph neural networks using tangent space aggregation. 

2. THE GEOMETRY OF CONSTANT CURVATURE SPACES:

Riemannian Geometry. A manifold M of dimension d is a generalization to higher dimensions of the notion of surface, and is a space that locally looks like R d . At each point x ∈ M, M can be associated a tangent space T x M, which is a vector space of dimension d that can be understood as a first order approximation of M around x. A riemannian metric g is given by an inner-product g x (•, •) at each tangent space T x M, g x varying smoothly with x. A given g defines the geometry of M, because it can be used to define the distance between x and y as the infimum of the lengths of smooth paths γ : [0, 1] → M from x to y, where the length is defined as (γ) :=

1 0 g γ(t) ( γ(t), γ(t))dt .

Under certain assumptions, a given g also defines a curvature at each point. Unifying all curvatures κ. There exist several models of respectively constant positive and negative curvatures. For positive curvature, we choose the stereographic projection of the sphere, while for negative curvature we choose the Poincaré model which is the stereographic projection of the Lorentz model. As explained below, this choice allows us to generalize the gyrovector space framework and unify spaces of both positive and negative curvature κ into a single model which we call the κ-stereographic model. The κ-stereographic model. For a curvature κ ∈ R and a dimension d ≥ 2, it is defined as

st d κ = {x ∈ R d | -κ x 2 2 < 1} equipped with its Riemannian metric g κ x =

4 (1+κ||x|| 2 ) 2 I =: (λ κ x ) 2 I. Note in particular that when κ ≥ 0, st d κ is R d , while when κ < 0 it is the open ball of radius 1/ √ -κ. Gyrovector spaces & Riemannian geometry. As discussed in section 1, the gyrovector space formalism is used to generalize vector spaces to the Poincaré model of hyperbolic geometry (Ungar, 2005; 2008) . In addition, important quantities from Riemannian geometry can be rewritten in terms of the Möbius vector addition and scalar-vector multiplication (Ganea et al., 2018a) . We here extend gyrovector spaces to the κ-stereographic model, i.e. allowing positive curvature.

For κ > 0 and any point x ∈ st d κ , we will denote by x the unique point of the sphere of radius κ -1 2 in R d+1 whose stereographic projection is x. As detailed in appendix C.2.2, it is given by

EQUATION

For x, y ∈ st d κ , we define the κ-addition, in the κ-stereographic model by:

x ⊕ κ y = (1 -2κx T y -κ||y|| 2 )x + (1 + κ||x|| 2 )y 1 -2κx T y + κ 2 ||x|| 2 ||y|| 2 ∈ st d κ .

(2)

The κ-addition is defined in all the cases except for spherical geometry and x = y/(κ y 2 ) as stated by the following theorem proved in Appendix C.2.1.

Theorem 1 (Definiteness of κ-addition). We have 1 -2κx T y + κ 2 ||x|| 2 ||y|| 2 = 0 if and only if κ > 0 and x = y/(κ y 2 ).

For s ∈ R and x ∈ st d κ (and |s tan -1

κ x | < κ 1 2 π/2 if κ > 0)

, the κ-scaling in the κ-stereographic model is given by:

EQUATION

where tan κ equals κ -1/2 tan if κ > 0 and (-κ) -1/2 tanh if κ < 0. This formalism yields simple closed-forms for various quantities including the distance function inherited from the Riemannian manifold (st d κ , g κ ), the exp and log maps, and geodesics, as shown by the following theorem.

Theorem 2 (Extending gyrovector spaces to positive curvature). For x, y ∈ st d κ , x = y, v = 0, (and x = -y/(κ y 2 ) if κ > 0), the distance function is given by a :

EQUATION

the unit-speed geodesic from x to y is unique and given by

EQUATION

and finally the exponential and logarithmic maps are described as:

EQUATION

a We write -x ⊕ y for (-x) ⊕ y and not -(x ⊕ y).

Proof sketch:

The case κ ≤ 0 was already taken care of by (Ganea et al., 2018a) . For κ > 0, we provide a detailed proof in Appendix C.2.2. The exponential map and unit-speed geodesics are obtained using the Egregium theorem and the known formulas in the standard spherical model. The distance then follows from the formula d κ (x, y) = log κ x (y) x which holds in any Riemannian manifold.

Around κ = 0. One notably observes that choosing κ = 0 yields all corresponding Euclidean quantities, which guarantees a continuous interpolation between κ-stereographic models of different curvatures, via Euler's formula tan(x) = -i tanh(ix) where i := √ -1. But is this interpolation differentiable with respect to κ? It is as shown by the following theorem, proved in Appendix C.2.3.

Theorem 3 (Smoothness of st d κ w.r.t. κ around 0). Let v = 0 and x, y ∈ R d , such that x = y (and x = -y/(κ y 2 ) if κ > 0). Quantities in Eqs. (4,5,6) are well-defined for |κ| < 1/ min( x 2 , y 2 ), i.e. for κ small enough. Their first order derivatives at 0 -and 0 + exist and are equal. Moreover, for the distance we have:

EQUATION

Note that for x T y ≥ 0, this tells us that an infinitesimal change of curvature from zero to small negative, i.e. towards 0 -, while keeping x, y fixed, has the effect of increasing their distance.

As a consequence, we have a unified formalism that interpolates smoothly between all three geometries of constant curvature.

3. κ-GCNS:

We start by introducing the methods upon which we build. We present our models for spaces of constant sectional curvature, in the κ-stereographic model. However, the generalization to cartesian products of such spaces (Gu et al., 2019) follows naturally from these tools.

3.1. GRAPH CONVOLUTIONAL NETWORKS:

The problem of node classification on a graph has long been tackled with explicit regularization using the graph Laplacian (Weston et al., 2012) . Namely, for a directed graph with adjacency matrix A, by adding the following term to the loss:

i,j A ij f (x i ) -f (x j ) 2 = f (X) T Lf (X)

, where L = D -A is the (unnormalized) graph Laplacian, D ii := k A ik defines the (diagonal) degree matrix, f contains the trainable parameters of the model and X = (x j i ) ij the node features of the model. Such a regularization is expected to improve generalization if connected nodes in the graph tend to share labels; node i with feature vector x i is represented as f (x i ) in a Euclidean space.

With the aim to obtain more scalable models, Defferrard et al. (2016) ; Kipf & Welling (2017) propose to make this regularization implicit by incorporating it into what they call graph convolutional networks (GCN), which they motivate as a first order approximation of spectral graph convolutions, yielding the following scalable layer architecture (detailed in appendix A):

EQUATION

where Ã = A + I has added self-connections, Dii = k Ãik defines its diagonal degree matrix, σ is a non-linearity such as sigmoid, tanh or ReLU = max(0, •), and W (t) and H (t) are the parameter and activation matrices of layer t respectively, with H (0) = X the input feature matrix.

3.2. TOOLS FOR A κ-GCN:

Learning a parametrized function f θ that respects hyperbolic geometry has been studied in (Ganea et al., 2018a) : neural layers and hyperbolic softmax. We generalize their definitions into the κ-stereographic model, unifying operations in positive and negative curvature. We explain how curvature introduces a fundamental difference between left and right matrix multiplications, depicting the Möbius matrix multiplication of (Ganea et al., 2018a) as a right multiplication, independent for each embedding. We then introduce a left multiplication by extension of gyromidpoints which ties the embeddings, which is essential for graph neural networks.

3.3. κ-RIGHT-MATRIX-MULTIPLICATION:

Let X ∈ R n×d denote a matrix whose n rows are d-dimensional embeddings in st d κ , and let W ∈ R d×e denote a weight matrix. Let us first understand what a right matrix multiplication is in Euclidean space: the Euclidean right multiplication can be written row-wise as (XW) i• = X i• W. Hence each d-dimensional Euclidean embedding is modified independently by a right matrix multiplication. A natural adaptation of this operation to the κ-stereographic model yields the following definition.

Definition 1. Given a matrix X ∈ R n×d holding κ-stereographic embeddings in its rows and weights W ∈ R d×e , the κ-right-matrix-multiplication is defined row-wise as

(X ⊗ κ W) i• = exp κ 0 ((log κ 0 (X)W) i• ) = tan κ ||(XW) i• || ||X i• || tan -1 κ (||X • i||) (XW) i• ||(XW) i• ||

where exp κ 0 and log κ 0 denote the exponential and logarithmic map in the κ-stereographic model.

This definition is in perfect agreement with the hyperbolic scalar multiplication for κ < 0, which can also be written as r ⊗ κ x = exp κ 0 (r log κ 0 (x)). This operation is known to have desirable properties such as associativity (Ganea et al., 2018a) .

3.4. κ-LEFT-MATRIX-MULTIPLICATION AS A MIDPOINT EXTENSION:

For graph neural networks we also need the notion of message passing among neighboring nodes, i.e. an operation that combines / aggregates the respective embeddings together. In Euclidean space such an operation is given by the left multiplication of the embeddings matrix with the (preprocessed) adjacency Â: l) . Let us consider this left multiplication. For A ∈ R n×n , the matrix product is given row-wise by:

H (l+1) = σ( ÂZ (l) ) where Z (l) = H (l) W (

(AX) i• = A i1 X 1• + • • • + A in X n•

This means that the new representation of node i is obtained by calculating the linear combination of all the other node embeddings, weighted by the i-th row of A. An adaptation to the κ-stereographic model hence requires a notion of weighted linear combination. We propose such an operation in st d κ by performing a κ-scaling of a gyromidpointwhose definition is reminded below. Indeed, in Euclidean space, the weighted linear combination αx+βy can be re-written as (α+β)m E (x, y; α, β) with Euclidean midpoint m E (x, y; α, β) := α α+β x + β α+β y. This motivates generalizing the above operation to st d κ as follows. (Ungar, 2008) and alternatively, here in eq. ( 10).

Definition 2. Given a matrix X ∈ R n×d holding κ-stereographic embeddings in its rows and weights A ∈ R n×n , the κ-left-matrix-multiplication is defined row-wise as

EQUATION

The κ-scaling is motivated by the fact that

d κ (0, r ⊗ κ x) = |r|d κ (0, x) for all r ∈ R, x ∈ st d κ .

We remind that the gyromidpoint is defined when κ ≤ 0 in the κ-stereographic model as (Ungar, 2010) :

EQUATION

with λ κ x = 2/(1 + κ x 2 ). Whenever κ > 0, we have to further require the following condition:

EQUATION

)

For two points, one can calculate that (λ κ x -1) + (λ κ y -1) = 0 is equivalent to κ x y = 1, which holds in particular whenever x = -y/(κ y 2 ). See fig. 3 for illustrations of gyromidpoints.

Our operation κ satisfies interesting properties, proved in Appendix C.2.4: Theorem 4 (Neuter element & κ-scalar-associativity). We have I n κ X = X, and for r ∈ R,

r ⊗ κ (A κ X) = (rA) κ X.

The matrix A. In most graph neural networks, the matrix A is intented to be a preprocessed adjacency matrix, i.e. renormalized by the diagonal degreee matrix D ii = k A ik . This normalization is often taken either (i) to the left:

D -1 A, (ii) symmetric: D -1 2 AD -1 2 or (iii)

to the right: AD -1 . Note that the latter case makes the matrix right-stochastic 1 , which is a property that is preserved by matrix product and exponentiation. For this case, we prove the following result in Appendix C.2.5:

Theorem 5 (κ-left-multiplication by right-stochastic matrices is intrinsic). If A, B are rightstochastic, φ is a isometry of st d κ and X,Y are two matrices holding κ-stereographic embeddings:

EQUATION

The above result means that A can easily be preprocessed as to make its κ-left-multiplication intrinsic to the metric space (st d κ , d κ ). At this point, one could wonder: does there exist other ways to take weighted centroids on a Riemannian manifold? We comment on two plausible alternatives.

Fréchet/Karcher means. They are obtained as arg min x i α i d κ (x, x i )2 ; note that although they are also intrinsic, they usually require solving an optimization problem which can be prohibitively expensive, especially when one requires gradients to flow through the solutionmoreover, for the space st d κ , it is known that the minimizer is unique if and only if κ ≥ 0. Tangential aggregations. They are defined by lifting the points in a chosen tangent space via the logarithmic map, performing a linear combination and then projecting back via the exponential map, and were in particular used in the recent works of Chami et al. (2019) and Liu et al. (2019) . The below theorem describes that for the κ-stereographic model, this operation is also intrinsic, i.e. commutes with isometries. We prove it in Appendix C.2.6.

Theorem 6 (Tangential aggregation is intrinsic). Define the tangential aggregation of x 1 , . . . ,

EQUATION

)

For any isometry φ of st d κ , we have

EQUATION

3.5 LOGITS Finally, we need the logit and softmax layer, a neccessity for any classification task. We here use the model of (Ganea et al., 2018a) , which was obtained in a principled manner for the case of negative curvature. We leave for future work the adaptation of their analysis to positive curvature and use in our experiments the straightforwardly adaptated formula to positive curvature, which we detail appendix D.

3.6. κ-GCN:

We are now ready to introduce our κ-stereographic GCN (Kipf & Welling, 2017), denoted by κ-GCN 2 . Assume we are given a graph with node level features G = (V, A, X) where X ∈ R n×d with each row X i• ∈ st d κ and adjacency A ∈ R n×n . We first perform a preprocessing step by mapping the Euclidean features to st d κ via the projection X → X/(2 |κ|||X|| max ), where ||X|| max denotes the maximal Euclidean norm among all stereographic embeddings in X. For l ∈ {0, . . . , L -2}, the (l + 1)-th layer of κ-GCN is given by:

EQUATION

where

H (0) = X, σ ⊗κ (x) := exp κ 0 (σ(log κ 0 (x)))

is the Möbius version (Ganea et al., 2018a ) of a pointwise non-linearity σ and Â = D-1 2 Ã D-1 2 . The final layer is a κ-logit layer (appendix D):

EQUATION

where W (L-1) contains the parameters a k and p k of the κ-logits layer. A very important property of κ-GCN is that its architecture recovers the Euclidean GCN when we let curvature go to zero:

κ-GCN κ-→0 ---→ GCN.

4. EXPERIMENTS:

We evaluate the architectures introduced in the previous sections on the tasks of node classification and minimizing embedding distortion for several synthetic as well as real datasets. We detail the training setup and model architecture choices to appendix E. Minimizing Distortion Our first goal is to evaluate the graph embeddings learned by our GCN models on the representation task of fitting the graph metric in the embedding space. We desire to minimize the average distortion, i.e. defined similarly as in (Gu et al., 2019) :

1 n 2 i,j d(xi,xj ) d G (i,j) 2 -1 2 ,

where d(x i , x j ) is the distance between the embeddings of nodes i and j, while d G (i, j) is their graph distance (shortest path length).

We create three synthetic datasets that best reflect the different geometries of interest: i) "Tree'": a balanced tree of depth 5 and branching factor 4 consisting of 1365 nodes and 1364 edges. ii) "Torus":

We sample points (nodes) from the (planar) torus, i.e. from the unit connected square; two nodes are connected by an edge iff their toroidal distance (the warped distance) is smaller than a fixed R = 0.01; this gives 1000 nodes and 30626 edges. iii) "Spherical Graph": we sample points (nodes) from S 2 , connecting nodes iff their distance is smaller than 0.2, leading to 1000 nodes and 17640 edges.

For the GCN models, we use 1-hot initial node features. We use two GCN layers with dimensions 16 and 10. The non-Euclidean models do not use additional non-linearities between layers. All euclidean parameters are updated using the ADAM optimizer with learning rate 0.01. Curvatures are learned using (stochastic) gradient descent and learning rate of 0.0001. All models are trained for 10000 epochs and we report the minimal achieved distortion. The results shown in table 1 reveal the benefit of our models. One can notice that estimated curvatures correspond to our geometric knowledge about these specific datasets.

4.1. NODE CLASSIFICATION:

We consider the popular node classification datasets Citeseer (Sen et al., 2008) , Cora-ML (McCallum et al., 2000) and Pubmed (Namata et al., 2012) . Node labels correspond to the particular subfield the published document is associated with. Dataset statistics and splitting details are deferred to the appendix E due to the lack of space.

Curvature Estimations of Datasets To understand how far are the real graphs of the above datasets from the Euclidean geometry, we first estimate the graph curvature of the four studied datasets using the deviation from the Parallelogram Law (Gu et al., 2019) Training Details We trained the Euclidean models with the hyperparameters chosen as reported in (Klicpera et al., 2019) . Namely, for GCN we use one hidden layer of size 64, dropout on the embeddings and the adjacency of rate 0.5 as well as L 2 -regularization for the weights of the first layer with λ = 0.02. Only for Cora-ML we had to adjust the regularization factor λ to 0.002 to ensure similar scores as achieved in (Klicpera et al., 2019) .

All Non-Euclidean models use biased-L2 regularization with α = 10 and λ = 2e -2. Euclidean models used L2 regularization with the same parameter λ. We used a combination of dropout and dropconnect for the non-Euclidean models. All models have the same number of parameters. We use 2 GCN layers, hidden dimension 64. Product models split hidden dimension into [32, 32] and also input features equally. Non-Euclidean models do not use additional non-linearities. Euclidean parameters use a learning rate of 0.01 for all models using ADAM. The curvatures are learned using gradient descent with a learning rate of 0.01. We show the values of the learned curvatures in appendix E. We use early stopping: we first train for a maximum of 2000 epochs, then we check every 200 epochs for improvement in the validation cross entropy loss; if that is not observed, we stop.

Node classification results. These are shown in table 2. It can be seen that our models are competitive with the two Euclidean GCN considered (with or without non-linearities), showcasing the benefit of our proposed architecture.

A GCN -A BRIEF SURVEY:

A.1 CONVOLUTIONAL NEURAL NETWORKS ON GRAPHS One of the pioneering works on neural networks in non-Euclidean domains was done by (Defferrard et al., 2016) . Their idea was to extend convolutional neural networks for graphs using tools from graph signal processing.

Given a graph G = (V, A), where A is the adjacency matrix and V is a set of nodes, we define a signal on the nodes of a graph to be a vector x ∈ R n where x i is the value of the signal at node i. Consider the diagonalization of the symmetrized graph Laplacian L = UΛU T , where Λ = diag(λ 1 , . . . , λ n ).

The eigenbasis U allows to define the graph Fourier transform x = U T x ∈ R n .

In order to define a convolution for graphs, we shift from the vertex domain to the Fourier domain:

x G y = U U T x U T y

Note that x = U T x and ŷ = U T y are the graph Fourier representations and we use the element-wise product since convolutions become products in the Fourier domain. The left multiplication with U maps the Fourier representation back to a vertex representation. As a consequence, a signal x filtered by g θ becomes y = Ug θ (Λ)U T x where g θ = diag(θ) with θ ∈ R n constitutes a filter with all parameters free to vary. In order to avoid the resulting complexity O(n), (Defferrard et al., 2016) replace the non-parametric filter by a polynomial filter:

g θ (Λ) = K-1 k=0 θ k Λ k

where θ ∈ R K resulting in a complexity O(K). Filtering a signal is unfortunately still expensive since y = Ug θ (Λ)U T x requires the multiplication with the Fourier basis U, thus resulting in complexity O(n 2 ). As a consequence, (Defferrard et al., 2016) circumvent this problem by choosing the Chebyshev polynomials T k as a polynomial basis, g θ (Λ) = K k=0 θ k T k ( Λ) where Λ = 2Λ λmax -I. As a consequence, the filter operation becomes y = K k=0 θ k T k ( L)x where L = 2L λmax -I. This led to a K-localized filter since it depended on the K-th power of the Laplacian. The recursive nature of these polynomials allows for an efficient filtering of complexity O(K|E|), thus leading to an computationally appealing definition of convolution for graphs. The model can also be built in an analogous way to CNNs, by stacking multiple convolutional layers, each layer followed by a non-linearity.

A.2 GRAPH CONVOLUTIONAL NETWORKS (Kipf & Welling, 2017) extended the work of (Defferrard et al., 2016) and inspired many follow-up architectures (Chen et al., 2018; Hamilton et al., 2017; Abu-El-Haija et al., 2018; Wu et al., 2019) . The core idea of (Kipf & Welling, 2017) is to limit each filter to 1-hop neighbours by setting K = 1, leading to a convolution that is linear in the Laplacian L:

g θ x = θ 0 x + θ 1 Lx

They further assume λ max ≈ 2, resulting in the expression

g θ x = θ 0 x -θ 1 D -1 2 AD -1 2 x

To additionally alleviate overfitting, (Kipf & Welling, 2017) constrain the parameters as θ 0 = -θ 1 = θ, leading to the convolution formula

g θ x = θ(I + D -1 2 AD -1 2 )x

Since I + D -1 2 AD -1 2 has its eigenvalues in the range [0, 2], they further employ a reparametrization trick to stop their model from suffering from numerical instabilities:

g θ x = θ D-1 2 Ã D-1 2 x

where Ã = A + I and Dii = n j=1 Ãij .

Rewriting the architecture for multiple features X ∈ R n×d1 and parameters

Θ ∈ R d1×d2 instead of x ∈ R n and θ ∈ R, gives Z = D-1 2 Ã D-1 2 XΘ ∈ R n×d2

The final model consists of multiple stacks of convolutions, interleaved by a non-linearity σ:

H (k+1) = σ D-1 2 Ã D-1 2 H (k) Θ (k)

where H (0) = X and Θ ∈ R n×d k .

The final output H (K) ∈ R n×d K represents the embedding of each node i as h i = H i• ∈ R d K and can be used to perform node classification:

Ŷ = softmax D-1 2 Ã D-1 2 H (K) W ∈ R n×L

where W ∈ R d K ×L , with L denoting the number of classes.

In order to illustrate how embeddings of neighbouring nodes interact, it is easier to view the architecture on the node level. Denote by N (i) the neighbours of node i. One can write the embedding of node i at layer k + 1 as follows:

h (k+1) i = σ   Θ (l) j∈Ni∪{i} h (k) j |N (j)||N (i)|  

Notice that there is no dependence of the weight matrices Θ (l) on the node i, in fact the same parameters are shared across all nodes.

In order to obtain the new embedding h (k+1) i of node i, we average over all embeddings of the neighbouring nodes. This Message Passing mechanism gives rise to a very broad class of graph neural networks (Kipf & Welling, 2017; Veličković et al., 2018; Hamilton et al., 2017; Gilmer et al., 2017; Chen et al., 2018; Klicpera et al., 2019; Abu-El-Haija et al., 2018) .

To be more precise, GCN falls into the more general category of models of the form

z (k+1) i = AGGREGATE (k) ({h (k) j : j ∈ N (i)}; W (k) ) h (k+1) i = COMBINE (k) (h (k) i , z (k+1) i ; V (k) )

Models of the above form are deemed Message Passing Graph Neural Networks and many choices for AGGREGATE and COMBINE have been suggested in the literature (Kipf & Welling, 2017; Hamilton et al., 2017; Chen et al., 2018) .

B GRAPH EMBEDDINGS IN NON-EUCLIDEAN GEOMETRIES:

In this section we will motivate non-Euclidean embeddings of graphs and show why the underlying geometry of the embedding space can be very beneficial for its representation. We first introduce a measure of how well a graph is represented by some embedding f :

V - → X , i → f (i): Definition 3. Given an embedding f : V - → X , i → f (i) of a graph G = (V, A) in some metric space X , we call f a D-embedding for D ≥ 1 if there exists r > 0 such that r • d G (i, j) ≤ d X (f (i), f (j)) ≤ D • r • d G (i, j)

The infimum over all such D is called the distortion of f . The r in the definition of distortion allows for scaling of all distances. Note further that a perfect embedding is achieved when D = 1.

B.1 TREES AND HYPERBOLIC SPACE:

Trees are graphs that do not allow for a cycle, in other words there is no node i ∈ V for which there exists a path starting from i and returning back to i without passing through any node twice. The Figure 5 : Euclidean embeddings of trees of different depths. All the four most inner circles are identical. Ideal node embeddings should match in distance the graph metric, e.g. the distance between the pink and green nodes should be the same as their shortest path length. Notice how we quickly run out of space, e.g. the pink and green nodes get closer as opposed to farther. This issue is resolved when embedding trees in hyperbolic spaces. number of nodes increases exponentially with the depth of the tree. This is a property that prohibits Euclidean space from representing a tree accurately. What intuitively happens is that "we run out of space". Consider the trees depicted in fig. 5 . Here the yellow nodes represent the roots of each tree. Notice how rapidly we struggle to find appropriate places for nodes in the embedding space because their number increases just too fast.

Moreover, graph distances get extremely distorted towards the leaves of the tree. Take for instance the green and the pink node. In graph distance they are very far apart as one has to travel up all the way to the root node and back to the border. In Euclidean space however, they are very closely embedded in a L 2 -sense, hence introducing a big error in the embedding. This problem can be very nicely illustrated by the following theorem: Theorem 7. Consider the tree K 1,3 (also called 3-star) consisting of a root node with three children. Then every embedding {x 1 , . . . , x 4 } with x i ∈ R k achieves at least distortion 2 √ 3 for any k ∈ N.

Proof. We will prove this statement by using a special case of the so called Poincaré-type inequalities (Deza & Laurent, 1996) :

For any b 1 , . . . , b k ∈ R with k i=1 b i = 0 and points x 1 , . . . , x k ∈ R n it holds that k i,j=1 b i b j ||x i -x j || 2 ≤ 0

Consider now an embedding of the tree x 1 , . . . , x 4 where x 1 represents the root node. Choosing b 1 = -3 and b i = 1 for i = 1 leads to the inequality

||x 2 -x 3 || 2 + ||x 2 -x 4 || 2 + ||x 3 -x 4 || 2 ≤ 3||x 1 -x 2 || 2 + 3||x 1 -x 3 || 2 + 3||x 1 -x 4 || 2

The left-hand side of this inequality in terms of the graph distance is

d G (2, 3) 2 + d G (2, 4) 2 + d G (3, 4) 2 = 2 2 + 2 2 + 2 2 = 12

and the right-hand side is

3 • d G (1, 2) 2 + 3 • d G (1, 3) 2 + 3 • d G (1, 4) 2 = 3 + 3 + 3 = 9

As a result, we always have that the distortion is lower-bounded by 12 9 = 2 √ 3

Euclidean space thus already fails to capture the geometric structure of a very simple tree. This problem can be remedied by replacing the underlying Euclidean space by hyperbolic space.

Consider again the distance function in the Poincaré model, for simplicity with c = 1:

d P (x, y) = cosh -1 1 + 2 ||x -y|| 2 (1 -||x|| 2 )(1 -||y|| 2 )

Assume that the tree is embedded in the same way as in fig. 5 , just restricted to lie in the disk of radius 1 √ c = 1. Notice that as soon as points move closer to the boundary (||x|| -→ 1), the fraction explodes and the resulting distance goes to infinity. As a result, the further you move points to the border, the more their distance increases, exactly as nodes on different branches are more distant to each other the further down they are in the tree. We can express this advantage in geometry in terms of distortion: Theorem 8. There exists an embedding x 1 , . . . , x 4 ∈ P 2 for K 1,3 achieving distortion 1 + for > 0 arbitrary small. Proof. Since the Poincaré distance is invariant under Möbius translations we can again assume that x 1 = 0. Let us place the other nodes on a circle of radius r. Their distance to the root is now given as

d P (x i , 0) = cosh -1 1 + 2 ||x i || 2 1 -||x i || 2 = cosh -1 1 + 2 r 2 1 -r 2

By invariance of the distance under centered rotations we can assume w.l.o.g. x 2 = (r, 0). We further embed

• x 3 = r cos( 2 3 π), r sin( 2 3 π) = -r 2 , √ 3 2 r • x 4 = r cos( 4 3 π), r sin( 4 3 π) = -r 2 , √ 3 2 r .

This procedure gives:

d P (x 2 , x 3 ) = cosh -1   1 + 2 || 3r 2 , - √ 3 2 r || 2 (1 -r 2 ) 2   = cosh -1 1 + 2 3r 2 (1 -r 2 ) 2

If we let the points now move to the border of the disk we observe that

cosh -1 1 + 2 3r 2 (1-r 2 ) 2 cosh -1 1 + 2 r 2 1-r 2 r-→1 ---→ 2

But this means in turn that we can achieve distortion 1 + for > 0 arbitrary small. QED.

The tree-likeliness of hyperbolic space has been investigated on a deeper mathematical level. (Sarkar, 2011) show that a similar statement as in theorem 8 holds for all weighted or unweighted trees. The interested reader is referred to (Hamann, 2017; Sarkar, 2011) for a more in-depth treatment of the subject.

Cycles are the subclasses of graphs that are not allowed in a tree. They consist of one path that reconnects the first and the last node: (v 1 , . . . , v n , v 1 ). Again there is a very simple example of a cycle, hinting at the limits Euclidean space incurs when trying to preserve the geometry of these objects (Matousek, 2013) . Theorem 9. Consider the cycle G = (V, E) of length four. Then any embedding (x 1 , . . . , x 4 ) where

x i ∈ R k achieves at least distortion √ 2.

Proof. Denote by x 1 , x 2 , x 3 , x 4 the embeddings in Euclidean space where x 1 , x 3 and x 2 , x 4 are the pairs without an edge. Again using the Poincaré-type inequality with b 1 = b 3 = 1 and b 2 = b 4 = -1 leads to the short diagonal theorem (Matousek, 2013) :

||x 1 -x 3 || 2 + ||x 2 -x 4 || 2 ≤ ||x 1 -x 2 || 2 + ||x 2 -x 3 || 2 + ||x 3 -x 4 || 2 + ||x 4 -x 1 || 2

The left hand side of this inequality in terms of the graph distance is d G (1, 3) 2 + d G (2, 4) 2 = 2 2 + 2 2 = 8 and the right hand side is 1 2 + 1 2 + 1 2 + 1 2 = 4. Therefore any embedding has to shorten one diagonal by at least a factor √ 2.

It turns out that in spherical space, this problem can be solved perfectly in one dimension for any cycle.

Theorem 10. Given a cycle G = (V, E) of length n, there exists an embedding {x 1 , . . . , x n } achieving distortion 1.

Proof. We model the one dimension spherical space as the circle S 1 . Placing the points at angles 2πi n and using the arclength on the circle as the distance measure leads to an embedding of distortion 1 as all pairwise distances are perfectly preserved.

Notice that we could also use the exact same embedding in the two dimensional stereographic projection model with c = 1 and we would also obtain distortion 1. The difference to the Poincaré disk is that spherical space is finite and the border does not correspond to infinitely distant points. We therefore have no since we do not have to pass to a limit.

C SPHERICAL SPACE AND ITS GYROSTRUCTURE:

Contrarily to hyperbolic geometry, spherical geometry is not only in violation with the fifth postulate of Euclid but also with the first. Notice that, shortest paths are not unique as for antipodal (oppositely situated) points, we have infinitely many geodesics connecting the two. Hence the first axiom does not hold. Notice that the third postulate holds as we stated it but it is sometimes also phrased as: "A circle of any center and radius can be constructed". Due to the finiteness of space we cannot have arbitrary large circles and hence phrased that way, the third postulate would not hold. Finally, we replace the fifth postulate by:

• Given any straight line l and a point p not on l, there exists no shortest line g passing through p but never intersecting l.

The standard model of spherical geometry suffers from the fact that its underlying space depends directly on the curvature κ through a hard constraint -κ x, x = 1 (similarly to the Lorentz model of hyperbolic geometry). Indeed, when κ → 0, the domain diverges to a sphere of infinite radius which is not well defined.

For hyperbolic geometry, we could circumvent the problem by moving to the Poincaré model, which is the stereographic projection of the Lorentz model, relaxing the hard constraint to an inequality. A similar solution is also possible for the spherical model.

C.1 STEREOGRAPHIC PROJECTION MODEL OF THE SPHERE:

In the following we construct a model in perfect duality to the construction of the Poincaré model. Fix the south pole z = (0, -1 √ κ ) of the sphere of curvature κ > 0, i.e. of radius R := κ -1 2 . The stereographic projection is the map:

Φ : S n R - → R n , x → x = 1 1 + √ κx n+1

x 1:n with the inverse given by

Φ -1 : R n - → S n R , x → x = λ κ x x, 1 √ κ (λ κ x -1)

where we define λ κ x = 2 1+κ||x|| 2 . Again we take the image of the sphere S n R under the extended projection Φ((0, . . . , 0, -1 κ )) = 0, leading to the stereographic model of the sphere. The metric tensor transforms as:

g κ ij = (λ κ x ) 2 δ ij C.2 GYROVECTOR SPACE IN THE STEREOGRAPHIC MODEL C.2.1 PROOF OF THEOREM 1

Using Cauchy-Schwarz's inequality, we have A :

= 1 -2κx T y + κ 2 ||x|| 2 ||y|| 2 ≥ 1 -2|κ| x y + κ 2 ||x|| 2 ||y|| 2 = (1 -|κ| x y ) 2 ≥ 0.

Since equality in the Cauchy-Schwarz inequality is only reached for colinear vectors, we have that A = 0 is equivalent to κ > 0 and x = y/(κ y 2 ).

C.2.2 PROOF OF THEOREM 2:

Let us start by proving that for x ∈ R n and v ∈ T x R n the exponential map is given by

EQUATION

Indeed, take a unit speed geodesic γ x,v (t) starting from x with direction v. Notice that the unit speed geodesic on the sphere starting from x ∈ S n-1 is given by Γ

x ,v (t) = x cos κ (t) + 1 √ κ sin κ (t)v

. By the Egregium theorem, we know that Φ(γ x,v (t)) is again a unit speed geodesic in the sphere where Φ

-1 : x → x = λ κ x x, 1 √ κ (λ κ x -1) . Hence Φ(γ x,v (t)

) is of the form of Γ for some x and v . We can determine those by

x = Φ -1 (γ(0)) = Φ -1 (x) = λ κ x x, 1 √ κ (λ κ x -1) v = Γ(0) = ∂Φ -1 (y) ∂y γ(0) γ(0) Notice that ∇ x λ κ x = -κ(λ κ x ) 2 x and we thus get v = -2κ(λ κ x ) 2 x T vx + λ κ x v - √ κ(λ κ x ) 2 x T v

We can obtain γ x,v again by inverting back by calculating γ x,v (t) = Φ(Γ x ,v (t)), resulting in

γ x,v (t) = (λ κ x cos κ (t) - √ κ(λ κ x ) 2 x T v sin κ (t))x + 1 √ κ λ κ x sin κ (t)v 1 + (λ κ x -1) cos κ (t) - √ κ(λ κ x ) 2 x T v sin κ (t) Denoting g κ x (v, v) = ||v|| 2 λ κ x we have that exp κ x (v) = γ x, 1 √ g κ x (v,v) v g κ x (v, v)

which concludes the proof of the above formula of the exponential map. One then notices that it can be re-written in terms of the κ-addition. The formula for the logarithmic map is easily checked by verifying that it is indeed the inverse of the exponential map. Finally, the distance formula is obtained via the well-known identity d κ (x, y) = log κ x (y) x where v x = g κ x (v, v). Note that as expected, exp κ

x (v) → κ→0 x + v, converging to the Euclidean exponential map.

C.2.3 PROOF OF THEOREM 3:

We first compute a Taylor development of the κ-addition w.r.t κ around zero:

EQUATION

We then notice that using the Taylor of • 2 , given by x

EQUATION

Finally Taylor developments of tan κ (|κ|

EQUATION

For κ → 0 -, tan κ (|κ|

EQUATION

The left and right derivatives match, hence even though κ → |κ|

1 2 is not differentiable at κ = 0, the function κ → tan κ (|κ| 1 2 u) is. A similar analysis yields the same conclusion for κ → |κ| -1 2 tan -1 κ (u) yielding For κ → 0, |κ| -1 2 tan -1 κ (u) = u -κu 3 /3 + O(κ 2 ).

(31) Since a composition of differentiable functions is differentiable, we consequently obtain that ⊗ κ , exp κ , log κ and d κ are differentiable functions of κ, under the assumptions on x, y, v stated in Theorem 3. Finally, the Taylor development of d κ follows by composition of Taylor developments:

d κ (x, y) = 2 κ -1 2 tan -1 κ ( (-x) ⊕ κ y ) = 2( x -y + κ((-x) T y) x -y 2 )(1 -(κ/3)( x -y + O(κ)) 2 ) + O(κ 2 )

= 2( xy + κ((-x) T y) xy 2 )(1 -(κ/3) xy 2 ) + O(κ 2 ) = 2 xy -2κ (x T y) xy 2 + xy 3 /3 + O(κ 2 ).

C.2.4 PROOF OF THEOREM 4:

If A = I n then for all i we have j A ij = 1, hence

(I n X) i• = 1 2 ⊗ κ   j δ ij λ κ xj k δ ik (λ κ x k -1) x j   (32) = 1 2 ⊗ κ λ κ xi (λ κ xi -1) x i (33) = 1 2 ⊗ κ (2 ⊗ κ x i ) (34) = x i (35) = (X) i• .

(36) For associativity, we first note that the gyromidpoint is unchanged by a scalar rescaling of A. The property then follows by scalar associativity of the κ-scaling.

C.2.5 PROOF OF THEOREM 5:

It is proved in (Ungar, 2005) that the gyromidpoint commutes with isometries. The exact same proof holds for positive curvature, with the same algebraic manipulations. Moreover, when the matrix A is right-stochastic, for each row, the sum over columns gives 1, hence our operation κ reduces to a gyromidpoint. As a consequence, our κ commutes with isometries in this case. Since isometries preserve distance, we have proved the theorem. and, thus, using the formula of the exp map from theorem 2 we obtain: tg φ(x) ({φ(x i )}; {α i }) = φ(x) ⊕ κ gyr[z, Rx]R(-x ⊕ κ exp κ

x (v)) (47) Using eq. ( 39), we get that gyr[z, Rx]Rw = -φ(x) ⊕ κ (z ⊕ κ (Rx ⊕ κ Rw)), ∀w ∈ st d κ (48) giving the desired tg φ(x) ({φ(x i )}; {α i }) = z ⊕ κ R exp κ

x (v) = φ (tg x ({x i }; {α i }))

D LOGITS:

The final element missing in the κ-GCN is the logit layer, a necessity for any classification task. We here use the formulation of (Ganea et al., 2018a) . Denote by {1, . . . , K} the possible labels and let

a k ∈ R d , b k ∈ R and x ∈ R d .

The output of a feed forward neural network for classification tasks is usually of the form p(y = k|x) = softmax( a k , x -b k )

In order to generalize this expression to hyperbolic space, the authors of (Ganea et al., 2018a) realized that the term in the softmax can be rewritten as Using this equation, they were able to obtain the following expression for the logit layer:

EQUATION

where a k ∈ T 0 st d κ ∼ = R d , x ∈ st d κ and p k ∈ st d κ . Combining all these operations leads to the definition of a hyperbolic feed forward neural network. Notice that the weight matrices W and the normal vectors a k live in Euclidean space and hence can be optimized by standard methods such as ADAM (Kingma & Ba, 2015) .

For positive curvature κ > 0 we use in our experiments the following formula for the softmax layer:

EQUATION

which is inspired from the formula i sin(x) = sinh(ix) where i := √ -1. However, we leave for future work the rigorous proof that the distance to geodesic hyperplanes in the positive curvature setting is given by this formula.

E MORE EXPERIMENTAL DETAILS:

We here present training details for the node classification experiments.

We closely follow the training and evaluation scheme from previous work, e.g. (Klicpera et al., 2019) . We split the data into training, early stopping, validation and test set. Namely we first split the dataset into a known subset of size n known and an unknown subset consisting of the rest of the nodes. For

:

M is right-stochastic if for all i, j Mij = 1.

To be pronounced "kappa" GCN; the greek letter κ being commonly used to denote sectional curvature

CONCLUSIONIn this paper, we introduced a natural extension of graph convolutional networks to the stereographic models of both positive and negative curvatures in a unified manner. We show how this choice of models permits to smoothly interpolate between positive and negative curvature, allowing the curvature of the model to be trained independent of an initial sign choice. We hope that our models will open new exciting directions into non-Euclidean graph neural networks.

Figures and Tables:

Figure 2: Geodesics in the Poincaré disk (left) and the stereographic projection of the sphere (right).

Figure 3: Left: Euclidean Linear combination αx + βy. Middle: Poincaré gyromidpoints (red dots) of two points for different weights on the left and Poincaré gyromidpoints in a hyperbolic triangle on the right with two equal and one free weight. Right: Möbius gyromidpoint in the Poincaré model defined by(Ungar, 2008) and alternatively, here in eq. (10).

Figure 4: Histogram of Curvatures from "Deviation of Parallogram Law"

a k , x -b k = sign( a k , x -b k )||a k ||d(x, H a k ,b k ) where H a,b = {x ∈ R d : x, a -b = 0} = {x ∈ R d : -p + x, a = 0} = Ha,p with p ∈ R d .As a first step, they define the hyperbolic hyperplane asHκ a,p = {x ∈ st d κ -p ⊕ κ x, a = 0}where now a ∈ T p st d κ and p ∈ st d κ . They then proceed proving the following formula: p ⊕ κ x, a | (1 + κ|| -p ⊕ κ x|| 2 )||a|| (50)

Minimum achieved average distortion of the different models. H and S denote hyperbolic and spherical models respectively.

Table content: <table><tr><td>Model</td><td>Tree</td><td colspan="2">Toroidal Graph Spherical Graph</td></tr><tr><td>GCN (Linear)</td><td>0.045</td><td>0.0607</td><td>0.0415</td></tr><tr><td>GCN (ReLU)</td><td>0.0502</td><td>0.0603</td><td>0.0409</td></tr><tr><td>H 10 -GCN</td><td>0.0029</td><td>0.272</td><td>0.267</td></tr><tr><td>S 10 -GCN</td><td>0.473</td><td>0.0485</td><td>0.0337</td></tr><tr><td colspan="2">H 5 × H 5 -GCN 0.0048</td><td>0.112</td><td>0.152</td></tr><tr><td>S 5 × S 5 -GCN</td><td>0.51</td><td>0.0464</td><td>0.0359</td></tr><tr><td>H 2 4 -GCN</td><td>0.025</td><td>0.084</td><td>0.062</td></tr><tr><td>S 2 4 -GCN</td><td>0.312</td><td>0.0481</td><td>0.0378</td></tr></table>

Node classification: Average accuracy across 10 splits with estimated uncertainties at 95 percent confidence level via bootstrapping on our datasplits. H and S denote hyperbolic (Poincaré ball model) and spherical (stereographic projection) models respectively.

Table content: <table><tr><td>Model</td><td>Citeseer</td><td>Cora-ML</td><td>Pubmed</td><td>MS Academic</td></tr><tr><td>GCN (ReLU)</td><td>75.7 ± 0.36</td><td colspan="3">83.31 ± 0.36 79.05 ± 0.52 92.14 ± 0.25</td></tr><tr><td colspan="4">GCN (Linear) 76.28 ± 0.30 83.81 ± 0.35 78.94 ± 0.50</td><td>92.3 ± 0.21</td></tr><tr><td>H 64 1 -GCN S 64 1 -GCN</td><td colspan="3">76.29 ± 0.3 76.18 ± 0.37 83.97 ± 0.31 79.04 ± 0.5 83.6 ± 0.34 79.01 ± 0.58</td><td>92.06±0.22 92.1±0.31</td></tr><tr><td>Prod-GCN</td><td>75.91 ± 0.34</td><td>82.9 ± 0.6</td><td>78.7 ± 0.53</td><td>91.9 ± 0.40</td></tr></table>

Peer review (iclr2020):
Summary: The authors propose using non-Euclidean spaces for GCNs. This is inspired by the recent work into non-Euclidean, and especially hyperbolic, embeddings. A few papers have recently tried to go past embeddings into building non-Euclidean models, requiring the lifting of standard operations in Euclidean space to non-Euclidean settings. This has been done in particular in hyperbolic space, but some datasets benefit from more complex spaces. The authors combine the mixed-curvature product formalism that uses products of Euclidean, hyperbolic, and spherical spaces for embeddings, but use these for GCN operations. Doing this requires, in particular, developing a reasonable way to perform these operations in spherical space (since Euclidean is trivial and hyperbolic has been recently worked on). The authors do a nice lifting via complex operations, and both the hyperbolic and spherical spaces can devolve into the flat Euclidean space when their curvature goes to 0. The authors implement these GCNs, train the curvatures, and demonstrate performance improvements over Euclidean only versions on node classification on benchmark datasets. They also give a fairly nice introduction to all of these ideas in an extended appendix. Strengths, Weaknesses, Recommendation: This paper is reasonably interesting---it joins an effort to produce non-Euclidean models in a tractable way, which is fairly challenging, but could have a good impact. On the plus side, it*s great that the authors added the nice development for the spherical operations, since that will come in handy for many models. The experiments are also good. On the downside, everything here is an extension of existing work, and the body of the paper is hard to read (though this may be inevitable, there*s a lot of background to go over here). Overall I recommend accepting it; I think it*s a solid contribution. Comments: - I don*t understand why the authors say that their space *interpolates smoothly* just because the limit in the curvature is the same from the left and right side. For example, the absolute value function has the same limit from the left and the right at 0, but it*s not differentiable there. Is it actually true that if we take the derivatives of the piecewise hyperbolic/spherical distance function that it*s differentiable at c=0? - There are a couple of recent papers that also consider hyperbolic GCNs, and in fact use similar ideas for the aggregation and update steps (i.e., same lift to hyperbolic space). However, these were recently NeurIPS papers, and the text is not yet out, so I don*t think this should affect the authors* independent work (and also the product part is new). I do recommend that the authors compare against those results in a future update of this work. The papers are *Hyperbolic Graph Convolutional Neural Networks* by Chami et al and Hyperbolic Graph Neural Networks by Liu et al. - One thing that I didn*t see discussed by the authors is that there are subtle difference between hyperbolic and spherical spaces. For example, the weighted midpoint of Def. 3.2 doesn*t immediately extend to spherical space (or at least won*t be unique). As an example, consider S^2 and the mean of two antipodal points on it---there*s many choices for the midpoint. You probably have to limit the operation to a half-sphere (there*s some ideas for this in Gu et al). - For the synthetic tree, why is the number of edges 2(|V|-1) rather than |V|-1? - Are the curvatures the same for each layer for the GCNs? This is an interesting point to discuss (some of the NeurIPS papers I mentioned train the curvature for each layer). Also, how do you select the number of factors of each type? - Minor, but some of these citations can be updated. The *De Sa* et al 2018 arxiv citation is really Sala et al and is an ICML *18 paper. Similarly, Gulcehre et al is a 2019 ICLR paper, and so on. It*s always good to get these right. - Is there any actual empirical importance from recovering the Euclidean case exactly for 0 curvature? The reason I ask is that my experience is that the hyperboloid is typically easier to work with. - One useful thing to point out in B.3.3 is that in general, it need not be a diffeomorphism for all of M for any manifold, which leads to non-uniqueness. In differential geometry, the *cut locus* is the region beyond which there is this non-uniqueness. - In the appendix, the statement *Sarkar (2011) show that a similar statement as in Theorem 2 holds for a very general class of trees* is confusing to me. The *general class*, as far as I know, is actually *all* trees, weighted or unweighted.