Unitary Convolutions for Learning on Graphs and Groups

Authors: bobak t kiani
lukas  fesser
melanie  weber
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 02:35:29

Abstract: Data with geometric structure is ubiquitous in machine learning often arising from fundamental symmetries in a domain, such as permutation-invariance in graphs and translation-invariance in images. Group-convolutional architectures, which encode symmetries as inductive bias, have shown great success in applications, but can suffer from instabilities as their depth increases and often struggle to learn long range dependencies in data. For instance, graph neural networks experience instability due to the convergence of node representations (over-smoothing), which can occur after only a few iterations of message-passing, reducing their effectiveness in downstream tasks. Here, we propose and study unitary group convolutions, which allow for deeper networks that are more stable during training. The main focus of the paper are graph neural networks, where we show that unitary graph convolutions provably avoid over-smoothing. Our experimental results confirm that unitary graph convolutional networks achieve competitive performance on benchmark datasets compared to state-of-the-art graph neural networks. We complement our analysis of the graph domain with the study of general unitary convolutions and analyze their role in enhancing stability in general group convolutional architectures.

Main paper: 
 introduction: 
 In recent years, the design of specialized machine learning architectures for structured data has received a surge of interest. Of particular interest are architectures for data domains with inherent symmetries, such as permutation-invariance in graphs and sets, translation-invariance in images, and other symmetries that arise from fundamental laws of physics in scientific data.Group-convolutional architectures allow for explicitly encoding symmetries as inductive biases, which has led to performance gains in scientific applications [ESM23, WWY20] ; theoretical studies have analyzed the impact of geometric inductive biases on the complexity of the architecture [MMM21, BVB21, KLL + 24]. Graph Neural Networks are among the most popular architectures for graph machine learning and have found impactful applications in a wide range of disciplines, including in chemistry [GRK + 21], drug discovery [ZAL18] , particle physics [SBV20] , and recommender systems [WSZ + 22]. However, despite these successes, several limitations remain. A notable difficulty is the design of stable, deep architectures. Many of the aforementioned applications require accurate learning of long-range dependencies in the data, which necessitates deeper networks. However, it has been widely observed that group-convolutional networks suffer from instabilities as their depths increases. On graph domains, these instabilities have been studied extensively in recent years, notably in the form of over-smoothing effects [RBM23] , which characterizes the fast convergence of the representations of nearby nodes with depth. This effect can often be observed after only a few iterations of message-passing (i.e., small number of layers) and can significantly decrease the utility of the learned representations in downstream tasks. While interventions for mitigating over-smoothing have been proposed, including targeted perturbations of the input graph's connectivity (rewiring) and skip connections, a more principled architectural approach with theoretical guarantees is still lacking. Similar effects, such as exploding or vanishing gradients, have also been studied in more general group-convolutional architectures, specifically in CNNs [TK21, SF21, LJW + 19], for which architectural interventions (e.g., skip connections) have been proposed.In this work, we take a different route. Inspired by a long line of work studying unitary recurrent neural networks [ASB16, HSL16, LCMR19, KBLL22], we propose to replace the standard group convolution operator with a unitary group convolution. By construction, the unitarity ensures that the linear transformations are norm-preserving and invertible, which can significantly enhance the stability of the network and avoid convergence of representations to a fixed point as its depth increases. We introduce two unitary graph convolution operators, which vary in the way the message passing and feature transformation are parameterized. We then generalize this approach to cover more general group-convolutional architectures.Our theoretical analysis of the proposed unitary graph convolutions shows that they enhance stability and prevent over-smoothing effects that decrease the performance of their vanilla counterparts. We further describe how generalized unitary convolutions avoid vanishing and exploding gradients, enhancing the stability of group-convolutional architectures without additional interventions, such as residual connections or batch normalization.We now provide a brief background into work that motivated and inspired this current study, deferring a more complete discussion to App. A. Unitary matrices have a long history of application in neural networks, specifically related to improving stability for deep networks [SMG13, PMB13] enhancing the learning of long-range dependencies in data [BSF94, Hoc91, SMG13] . [ASB16, JSD + 17, HSL16] implemented unitary matrices in recurrent neural networks to address issues with the challenge of vanishing and exploding gradients inherent to learning long sequences of data in RNNs [BSF94, LJH15] . These original algorithms were later improved to be more expressive while still being efficient to implement in practice [HWY18, LCMR19, KBLL22] . For graph convolution, [HSTW20] discuss applications of the exponential map to linear convolutions; here, we use this same exponential map to explicitly apply unitary operators parameterized in the Lie algebra. For image data, various works [SGL18, LHA + 19, TK21, SF21, KBLL22] design and analyze variants of orthogonal or unitary convolution used in CNN layers. These can be viewed as a particular instance of the group convolution we study here over the cyclic group. More recently, proposals for unitary or orthogonal message passing have shown improvements in stability and performance compared to conventional message passing approaches [GZH + 22, AEL + 24, QBY24]. However, in contrast to our work, these methods do not always implement a unitary transformation across the whole input (e.g. only applying it in the feature transformation) and in the case of [QBY24] can be computationally expensive to implement for large graphs (see App. A for more detail). 
background and notation: 
 We denote scalars, vectors, and matrices as 𝑐, w, and M respectively. Given a matrix M , its conjugate transpose and transpose are denoted M † and M ⊤ . Given two matrices, A and B, we denote their tensor product (or Kronecker product) as A⊗B. Given a vector w, its standard Euclidean norm is denoted ∥w∥. For a matrix M , we denote its operator norm as ∥M ∥ and Frobenius norm as ∥M ∥ 𝐹 .Group Theory Basics Symmetries ("invariances") describe transformations, which leave properties of data unchanged ("invariant"), and as such characterize the inherent geometric structure of the data domain. Algebraically, symmetries can be characterized as groups. We say that a group is a matrix Lie group, if it is a differentiable manifold and a subgroup of the set of invertible 𝑛×𝑛 matrices (see App. C). Lie groups are associated with a Lie algebra, a vector space, which is formed by its tangent space at the identity. A comprehensive introduction into Lie groups and Lie algebras can be found in [FH13, Hal15] . Throughout this work we will encounter the 𝑛-dimensional orthogonal 𝑂 (𝑛) and unitary 𝑈 (𝑛) Lie groups, which are defined asEQUATIONLie algebras of 𝑂 (𝑛) and 𝑈 (𝑛) are the set of skew symmetric 𝔬(𝑛) and skew Hermitian 𝔲(𝑛) matrices respectively, i.e.,EQUATIONGiven a matrix M ∈ 𝔬(𝑛) (or 𝔲(𝑛)), the matrix exponential maps the matrix to an element of the Lie group exp(M ) ∈ 𝑂 (𝑛) (or 𝑈 (𝑛)). More details on unitary parametrizations can be found in App. C.2.We denote graphs by G = (𝑉, 𝐸) where 𝑉 and 𝐸 denote the set of nodes and edges respectively. For a graph on 𝑛 nodes, unless otherwise specified, we let 𝑉 = {1, . . . , 𝑛} index the nodes and denote the adjacency matrix by A ∈ R 𝑛×𝑛 and node features x ∈ R 𝑛×𝑑 . We also often use D ∈ R 𝑛×𝑛 to denote the diagonal degree matrix where diagonal entry 𝑖 records the degree of node 𝑖. The normalized adjacency matrix is defined as A = D -1/2 AD -1/2 . Given a node feature matrix X ∈ R 𝑛×𝑑 in where row 𝑖 denotes the 𝑑 in -dimensional feature vector of node 𝑖, graph convolution operators take the general form [KW16, ZK20]EQUATIONwhere W 0 , . . . W 𝑘 ∈ R 𝑑 in ×𝑑 out are trainable parameters. For ease of presentation, we will often omit the adjacency matrix as an explicit input to the operation. Often, only a single "message passing" step is included and the operation takes the simple form 𝑓 conv (X) = AXW . Eq. ( 3) is equivariant under any permutation matrix P 𝜋 ∈ R 𝑛×𝑛 5 since 𝑓 conv (P 𝜋 X; P 𝜋 AP -1 𝜋 ) = P 𝜋 • 𝑓 conv (X; A).(4)We aim to parameterize a particular subset of these operations which preserve unitarity properties.In general, a linear convolution operation conv 𝐺 : C 𝑛×𝑐 → C 𝑛×𝑐 takes the form of a weighted sum over linear transformations that are equivariant to a given group 𝐺. For simplicity, we assume here that we are working with finite groups though this can be generalized to other settings [KT18, CGKW18, CW16] . Given an input X ∈ C 𝑛×𝑐 consisting of 𝑐 channels in a vector space of dimension 𝑛, we study convolutions of the formEQUATIONwhere T 1 , . . . , T 𝑚 ∈ C 𝑛×𝑛 are linear operators equivariant to the group 𝐺 and W 1 , . . . , W 𝑚 ∈ C 𝑐×𝑐 are parameterized weight matrices. The graph setting is recovered by setting T 𝑘 = A 𝑘 . Similarly, for cyclic convolution as in conventional CNNs, one sets T 𝑘 to be the circulant matrices sending basis vector T 𝑘 e 𝑖 = e 𝑖+𝑘 where indexing is taken mod 𝑛.Question: Does the paper discuss the limitations of the work performed by the authors?Answer: [Yes] Justification: We have described limitations of our work at various points. Incorporating unitarity is not a panacea and we mention the challenges in dealing with this (e.g. see discussion section).Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations. 
unitary group convolutions: 
 We first describe unitary convolution for data on graphs (equivariant to permutations) which is the main focus of our study and then detail general procedures for performing unitary convolutions equivariant to general finite groups. Implementing these operations often requires special considerations to handle nonlinearities, complex numbers, initialization, etc. which we discuss in App. E.We introduce two variants of unitary graph convolution, which we denote as Separable unitary convolution (short UniConv) and Lie orthogonal/unitary convolution (short Lie UniConv). UniConv is a simple adjustment to standard message passing and treats linear transformations over nodes and features separately. Lie UniConv, in contrast, parameterizes operations in the Lie algebra of the orthogonal/unitary groups. This operation is fully unitary, but does not have the tensor product nature of the separable UniConv.By introducing complex numbers, we can enforce unitarity separately in the message passing and feature transformation. The unitary message passing has a wave-like nature which ensures messages "propagate" through the graph. In contrast, the standard message passing has a unique fixed point corresponding to the all ones vector which inherently causes oversmoothing in the features. Here, 𝑐 is chosen to ensure the operator norm of the matrix I + A is bounded by one.Definition 1 (Separable unitary graph convolution (UniConv)). Given an undirected graph G over 𝑛 nodes with adjacency matrix A ∈ R 𝑛×𝑛 , separable unitary graph convolution (UniConv) 𝑓 Uconv :EQUATION)where U ∈ 𝑈 (𝑑) is a unitary operator and 𝑡 ∈ R controls the magnitude of the convolution.One feature of the complexification of the adjacency matrix is that messages propogate as "waves" as observed for example in Fig. 1 . Since A is a symmetric matrix, exp(𝑖A𝑡) is unitary for all values of 𝑡 ∈ R and corresponds to vanilla message passing up to first order: exp(𝑖A𝑡) ≈ I + 𝑖A𝑡 + 𝑂 (𝑡 2 ). Remark 2. We observe that performance on real-world tasks is usually improved when enforcing unitarity in the node message passing where oversmoothing occurs, but not necessarily when enforcing unitarity in the feature transformation U . Thus, one can choose to leave U in Eq. (6) as a fully parameterized (unconstrained) matrix as we often do in our experiments.More generally, one can parameterize the operation in the Lie algebra by first forming a skew Hermitian convolution operation 𝑔 conv : C 𝑛×𝑑 → C 𝑛×𝑑 and then applying the exponential map. This approach has the benefit that it can be fully implemented using real numbers to obtain an orthogonal operator by enforcing constraints in the real part of the weight matrix W only. Definition 3 (Lie orthogonal/unitary graph convolution (Lie UniConv)). Given an undirected graph G over 𝑛 nodes with adjacency matrix A ∈ R 𝑛×𝑛 , Lie unitary/orthogonal graph convolution (Ortho-Conv) 𝑓 Uconv :EQUATION)The operation of 𝑔 conv can be represented as a vector-matrix operation vec(𝑔 conv (X)) = A ⊗ W ⊤ vec(𝑋) where A ⊗ W ⊤ belongs in the Lie algebra. If W is real-valued, the above returns an orthogonal map since the exponential of a real-valued matrix is real-valued.The exponential map in Definitions 1 and 3 can be performed using accurate approximations with typically constant factor overhead in runtime. We use the simple 𝐾-th order Taylor approximationEQUATION)For experiments, we find that setting 𝐾 = 10 suffices in all cases as the error exponentially decreases with 𝐾. Various other accurate and efficient approximations exist as detailed in App. C.2. We also refer the reader to App. E for other implementation details associated to handling complex numbers, initialization, etc.In the more general setting, we are concerned with parameterizing unitary operations of the formEQUATIONwhere T 1 , . . . , T 𝑚 ∈ C 𝑛×𝑛 are linear operators equivariant to the group 𝐺 and W 1 , . . . , W 𝑚 ∈ C 𝑐×𝑐 are parameterized weight matrices (e.g., set T 𝑘 = A 𝑘-1 to recover graph convolution). One can enforce and parameterize unitary convolutions in Eq. ( 9) in the Lie algebra basis or in the Fourier domain as detailed in App. D.Input: equivariant linear operator L ∈ C 𝑛 → C 𝑛 Input: vector x ∈ C 𝑛 1: L = 1 2 (L -L † ) (skew symmetrize operator) 2: return exp( L)(x) (or approximation thereof)In the Lie algebraic setting (Algorithm 1), one explicitly parameterizes operators in the Lie algebra or orthogonally projects arbitrary linear operators onto this basis by the mapping X ↦ → (X -X † )/2. This parameterization is particularly simple to implement (it is a linear basis) and unitary operators are subsequently implemented by applying the exponential map. This setting covers previous implementations of unitary RNNs and CNNs [LCMR19, SF21] and is detailed in Sec. 3.1 for GNNs. EQUATION)Parameterizing operations on the Lie algebra simply requires that𝑚(𝑔) = 𝑚(𝑔 -1 ) * since R -1 𝑔 = R † 𝑔 .An example implementation of the above for a toy learning task on the dihedral group is in App. F.3. One can also generally implement convolutions in the (block diagonal) Fourier basis of the graph or group (see App. D and Algorithm 3). Here, one employs a Fourier operator which block diagonalizes the input into its irreducible representations or some spectral representation. Fourier representations often have the advantage of being faster to implement due to efficient Fourier transforms. Since we do not use this in our experiments, we defer the details to App. D. 
properties and theoretical guarantees: 
 Unitary operators are now well studied in the context of neural networks and have various properties that are useful in naturally enhancing stability and performance of learning architectures [ASB16, SMG13, LCMR19, JSD + 17, KBLL22]. These properties and their theoretical guarantees are outlined here. We defer all proofs to App. B, many of which follow immediately from the definition of unitarity.Throughout, we will assume that convolution operators act on a vector space V and are built from a basis of linear operators that is equivariant to input and output representation 𝜌(𝑔) of a group 𝐺. We set input/output representations to be equal so that the exponential map of an equivariant operator is itself equivariant.6Typically called cross-correlation in mathematics literature.Fact 1 (Basic properties). Any unitary convolution 𝑓 Uconv : V → V built from Algorithms 1 and 3 meets the basic properties:EQUATIONA simple corollary of the above isometry property leveraged in prior work on unitary CNNs [TK21, SGL18] is that unitary matrices naturally provide robustness guarantees and a provable means to bound the effects of adversarial perturbations (see Corollary 9 in App. B). For graphs, we note that the properties above are generally impossible to obtain with graph convolution operations that perform a single message passing step as we show below. Proposition 4. Let 𝑓 conv : R 𝑛×𝑑 → R 𝑛×𝑑 be a graph convolution layer of the formEQUATION) where W 0 , W 1 ∈ R 𝑑×𝑑 are parameterized matrices. The linear map 𝑓 (•, A) : R 𝑛×𝑑 → R 𝑛×𝑑 is orthogonal for all adjacency matrices A of undirected graphs only if W 1 = 0 and W 0 ∈ 𝑂 (𝑑) is orthogonal. Furthermore, denoting J A ∈ R 𝑛𝑑×𝑛𝑑 as the Jacobian matrix of the map 𝑓 conv (•, A), for any choice of W 0 , W 1 , there always exists a normalized adjacency matrix Â such thatEQUATION)where ∥M ∥ is the operator norm of matrix M .The above shows that one must apply higher powers of A as in the exponential map to achieve a linear operator that is close to orthogonal.Oversmoothing During the training of GNNs we often observe that the features of neighboring nodes become more similar as the depth of the networks (i.e., the number of message-passing iterations) increases. This "oversmoothing" phenomenon has a strong connection to the spectral properties of graphs where convergence of a function on the graph is measured through the Dirichlet form7 or its normalized variant also termed the Rayleigh quotient. Definition 5 (Rayleigh quotient [Chu97] ). Given an undirected graph G = (𝑉, 𝐸) on |𝑉 | = 𝑛 nodes with adjacency matrix A ∈ {0, 1} 𝑛×𝑛 , let D ∈ R 𝑛×𝑛 be a diagonal matrix where the 𝑖-th entry D 𝑖𝑖 = 𝑑 𝑖 and 𝑑 𝑖 is the degree of node 𝑖. Let 𝑓 : 𝑉 → C 𝑑 be a function from nodes to features. Then the Rayleigh quotient 𝑅 G ( 𝑓 ) is equal toEQUATIONwhere A = D -1/2 AD -1/2 is the normalized adjacency matrix and X ∈ C 𝑛×𝑑 is a matrix with the 𝑖-th row set to feature vector 𝑓 (𝑖). We will at times abuse notation and let X be an input to 𝑅 G (X).Given an undirected graph G on 𝑛 nodes with normalized adjacency matrix A = D -1/2 AD -1/2 , we compare the Rayleigh quotient of normalized vanilla and unitary convolution. First, it is straightforward to show that the Rayleigh quotient is invariant to unitary transformations giving a proof that unitary graph convolution avoids oversmoothing. Proposition 6 (Invariance of Rayleigh quotient). Given an undirected graph G on 𝑛 nodes with normalized adjacency matrixA = D -1/2 AD -1/2 , the Rayleigh quotient 𝑅 G (X) = 𝑅 G ( 𝑓 Uconv (X))is invariant under normalized unitary or orthogonal graph convolution (see Definitions 1 and 3).In contrast, oversmoothing commonly occurs with vanilla graph convolution and has been proven to occur in a variety of settings [RBM23, CW20, BDGC + 22, Ker22]. To illustrate this, we exhibit a simple setting below commonly found at initialization where the parameterized matrix is set to be an orthogonal matrix and input features are random. Here, the magnitude of oversmoothing concentrates around its average and grows with the value of Tr( A 3 ) which corresponds to the (weighted) number of triangles in the graph.7This quantity has various equivalent names including the Dirichlet energy, local variance, and Laplacian quadratic form.Proposition 7. Given a simple undirected graph G on 𝑛 nodes with normalized adjacency matrix A = D -1/2 AD -1/2 and node degree bounded by 𝐷, let X ∈ R 𝑛×𝑑 have rows drawn i.i.d. from the uniform distribution on the hypersphere in dimension 𝑑. Let 𝑓 conv (X) = AXW denote convolution with orthogonal feature transformation matrix W ∈ 𝑂 (𝑑). Then, the event below holds with probability 1exp(-Ω( √ 𝑛)):𝑅 G (X) ≥ 1 -𝑂 1 𝑛 1/4and 𝑅 G ( 𝑓 conv (X)) ≤ 1 - 2a ). Node features x 𝑖 ∈ R are a single number set to zero for all but two randomly chosen nodes whose features are set to one. The goal is to predict the distance between these two randomly chosen nodes. For a graph of 𝑛 nodes, conventional message passing architectures require at least 𝑛/2 sequential messages to fully learn this dataset. As shown in Fig. 2b , conventional message passing networks fail to learn this task whereas the unitary convolutional architecture succeeds. We refer the reader to App. F.1 for further details and results for additional architectures.Tr( A 3 ) Tr( A 2 ) + 𝑂 1 𝑛 1/8For complex valued inputs, we additionally require that the functions are holomorphic. For sake of simplicity, we will treat functions 𝑓 : C 𝑛 → C 𝑛 as real-valued functions 𝑓 : R 2𝑛 → R 2𝑛 .9Technically, this holds for almost all x ∈ C 𝑛 due to non-differentiability of the activation at singular points, but this is handled in practice by choosing a gradient in the sub-differential which meets the criteria. Coco and Pascal are node classification tasks.Table 1 shows that the Unitary GCN outperforms standard message passing architectures and is competitive with other state of the art architectures as well, many of which employ global attention mechanisms. This provides evidence that the unitary GCN is nearly as effective at learning long-range signals.For node classification, we consider the Heterophilous Graph Dataset proposed by [PKD + 23]. The dataset contains the heterophilous graphs Roman-empire, Amazon-ratings, Minesweeper, Tolokers, and Questions, which are often considered as a benchmark for evaluating the performance of GNNs on graphs where connected nodes have dissimilar labels and features. Our results in Table 2 show that the unitary GCN outperforms the baseline message passing and graph transformer models. Given the heterophilous nature of this dataset, these findings reinforce the notion that unitary convolution enhances the ability of convolutional networks to capture long-range dependencies. 
open access to data and code: 
 Question: Does the paper provide open access to the data and code, with sufficient instructions to faithfully reproduce the main experimental results, as described in supplemental material?Answer Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
discussion: 
 In this paper we introduced unitary graph convolutions for enhancing stability in graph neural networks. We provided theoretical and empirical evidence for the effectiveness of our approach. We further introduce an extension to general groups (generalized unitary convolutions), which can be leveraged in group-convolutional architectures to enhance stability. Limitations Perhaps the biggest challenge in working with unitary convolutions is the overhead associated with maintaining unitarity or orthogonality via approximations of the exponential map or diagonalizations in Fourier or spectral bases. For implementing unitary maps, working with complex numbers also requires different initialization and activation functions. We refer the reader to App. C.2 and E for methods to alleviate these challenges in practice. Separately, there may be target functions or problem instances where unitarity or orthogonality may not be appropriate. For example, one can envision node classification tasks where the target function is neither (approximately) invertible nor isometric. In such instances, non-unitary layers will be required to learn the task. More generally, deciding when to use unitary layers is problem-dependent. In some cases, such as applications with smaller input graphs, simple and more efficient interventions such as adding residual connections or including batch norm will likely suffice for addressing signal propagation problems.In the graph domain, extensions of graph convolution to more advanced methods such as those that better incorporate edge features could widen the range of applications. Exploring hybrid models that combine unitary and non-unitary layers (e.g. global attention mechanisms) could potentially lead to more robust and versatile graph neural networks. Future work can also improve the efficiency of the parameterizations and implementations of the exponential map (see App. C.2). In a similar vein, it is likely that approximately unitary/orthogonal layers suffice in many settings to achieve the performance gains we see in our work. Methods that approximately enforce or regularize layers towards unitarity may be of interest in these instances due to their potential for improved efficiency. In this study, we mainly focused on applications to graph classification and regression tasks; however, the proposed methodology is much more general and could open up a wider range of applications to domains with more general symmetries or different data domains. For example, unitary matrices offer provable guarantees to adversarial attacks (see Corollary 9) and testing this robustness in practice on geometric data has yet to be conducted. AA.1 Further related literature Unitary RNNs Unitary neural networks were initially developed to tackle the challenge of vanishing and exploding gradients encountered in recurrent neural networks (RNNs) processing lengthy sequences of data. Aiming for more efficient information learning compared to established nonunitary architectures like the long-short term memory unit (LSTM) [HS97] , early approaches ensured unitarity through a sequence of parameterized unitary transformations [ASB16, MHRB17, JSD + 17].Other methods like the unitary RNN (uRNN) [WPH + 16], the Cayley parameterization (scoRNN) [HWY18] , and exponential RNN (expRNN) [LCMR19] parameterized the entire unitary space, maintaining unitarity via Cayley transformations or parameterizing the Lie algebra of the unitary group and performing the exponential map. To eliminate the reliance on the matrix inversion or matrix exponential in these prior algorithms which are computationally expensive in higher dimensions, [KBLL22] showed how to fully parameterize unitary operations more efficiently when applying gradient updates in a low rank subspace. We compare here previous GNN architectures which have proposed variants of complex-valued or isometric message passing. To facilitate this comparison, first let us recall some standard notation and implementations. Given an (possibly normalized) adjacency matrix A ∈ R 𝑛×𝑛 acting on 𝑛 nodes and input features X ∈ R 𝑛×𝑑 of dimension 𝑑, a basic linear message passing layer 𝑓 : R 𝑛×𝑑 → R 𝑛×𝑑 ′ takes the formEQUATIONwhere 𝜎 : R → R is a pointwise nonlinearity and W ∈ R 𝑑×𝑑 ′ is a parameterized weight matrix transforming the node features. There are variations of the above which normalize the adjacency matrix or add residual connections, but for sake of comparison, we will study the simplified form above. In contrast, the unitary message passing layer we propose takes the formEQUATIONwhere X and W are now complex-valued and 𝑡 ∈ R is a parameter controlling the strength of message passing.[GZH + 22] propose a version of an orthogonal GNN which enforces either exactly or approximately that the matrix W in Eq. ( 16) is orthogonal, i.e.𝑓 W (X) = 𝜎(AXW ), whereEQUATION)They also discuss regularization of node transformations W to bias towards orthogonality. [AEL + 24] also detail a version of a GNN which enforces orthogonality in the matrix W . We note that our study focuses on orthogonality in the node evolution during message passing. We discuss in the main text how to combine this with orthogonality in the feature transformation W .Our work also has overlap with the recently proposed fractional Graph Laplacian approach of [MPBK24] . There, they analyze a continuous-time message passing scheme called the fractional Schrödinger equation taking the formEQUATIONwhere 𝛼 > 0 is a chosen parameter and L is the normalized graph Laplacian. Unitarity can be strictly enforced in this continuous time format although [MPBK24] do not take this approach.In fact, solving this continuous-time differential equation for some time 𝑡, we obtain vec(X) (𝑡) = exp(𝑖𝑡L 𝛼 ⊗ W ) vec(X) (0) where vec(X) (𝑡) is the vectorized version of X at time 𝑡. Constraining L 𝛼 ⊗ W to be Hermitian obtains a unitary transformation. Beyond enforcing unitarity, our method works over discrete time, and we find it more efficient to parameterize the weight matrix W outside of the exponential map.A recent preprint proposes a unitary message passing algorithm called graph unitary message passing (GUMP) which sends messages through a larger ring graph constructed from the edges of the original graph [QBY24] . Given a graph 𝐺 with nodes and edges 𝑉 and 𝐸 respectively, the steps of this approach loosely are as follows:1. Construct a new directed graph (digraph) 𝐺 ′ with the same nodes and include directed edges (𝑖, 𝑗) and ( 𝑗, 𝑖) for any undirected edge (𝑖, 𝑗) in original graph.2. Convert 𝐺 ′ to its line graph 𝐿(𝐺 ′ ): here, each node is an edge in 𝐺 ′ and two nodes in 𝐿(𝐺 ′ ) share an edge if there is a path through the original edges in 𝐺 ′ (e.g. nodes corresponding to edges (𝑖, 𝑗) and ( 𝑗, 𝑘)).3. Construct new node representations in the line graph where for each edge (𝑖, 𝑗) we append the node representations [x 𝑖 , x 𝑗 ].4. Given the adjacency matrix A 𝐿 for 𝐿 (𝐺 ′ ), find permutation matrices P 1 , P 2 which form a block diagonal matrix P ⊺ 1 A 𝐿 P 2 . 5. Project each block to its closest unitary in Frobenius norm according to the closed form projection operation [Kel75] :EQUATION)In the above, U (𝑛) denotes the set of unitary matrices in C 𝑛×𝑛 .6. Invert the permutations obtaining P 1 A 𝐿 P ⊺ 2 . 7. Perform GNN convolution using the adjacency matrix constructed previously for 𝑘 layers. 8. Attach node representations of 𝐿 (𝐺 ′ ) to corresponding nodes of 𝐺 and append these to the initial representations x to output the final representations.As evident above, the approach here is rather meticulous and many details of the approach are left out. We refer the reader to the preprint [QBY24] for those details. Key to this procedure is that unitarity is roughly enforced by constructing the matrix A 𝐿 which is itself unitary in the vector space of the line graph 𝐿(𝐺 ′ ). Nonetheless, the map from node representations on the original graph 𝐺 to new node representations on 𝐺 will only be approximately unitary as the additional steps here will not guarantee properties of isometry and invertibility in general.For sake of comparison, we should note that the complexity of this approach (in number of nodes Our approach is also noticeably simpler as it only reparameterizes the message passing step to be an exponential map over standard message passing procedures.First, we review the properties shown in Fact 1 that unitary convolution meets the properties of invertibility, isometry, and equivariance. This fact follows virtually immediately from the definition of unitarity and equivariance.Proof of Fact 1. The proofs of isometry and invertibility follow directly from the definition of unitarity. What remains to be shown is the equivariance property. For unitary convolution built using Algorithm 1, equivariance can be checked by noting that the exponential map is a composition of equivariant linear operators. For unitary convolution in the Fourier domain (Algorithm 3), equivariance follows from convolution theorems where linear operators appropriately applied in the block diagonal Fourier basis are equivariant [FH13, KT18] . □As an application, unitary transformations are often used to provide provable robustness guarantees to adversarial perturbations of inputs [TK21, LHA + 19, MMS + 17, SGL18].Corollary 9 (Certified robustness to adversarial perturbations [TK21] ). A simple consequence of the isometry property of the unitary convolution is that ∥ 𝑓 Uconv (x)-𝑓 Uconv (y) ∥ = ∥x-y∥ so the Lipschitz constant of this function is 1. Assume neural network classifier 𝑓 : R 𝑛 → R 𝑚 is composed of unitary transformations and 1-Lipschitz bounded nonlinearities followed by an 𝐿-Lipschitz transformation and for given input x, it has margin M 𝑓 (x) defined asEQUATIONThen, 𝑓 is certifiably robust (i.e. classification is unchanged) to perturbations x + Δ of magnitude∥Δ∥ < M 𝑓 (x) ( √ 2𝐿) -1 .In the main text, we noted in Proposition 4 that the above facts in some way cannot be obtained using the standard graph convolution. We restate this here and prove it. Proposition 4. Let 𝑓 conv : R 𝑛×𝑑 → R 𝑛×𝑑 be a graph convolution layer of the formEQUATION) where W 0 , W 1 ∈ R 𝑑×𝑑 are parameterized matrices. The linear map 𝑓 (•, A) : R 𝑛×𝑑 → R 𝑛×𝑑 is orthogonal for all adjacency matrices A of undirected graphs only if W 1 = 0 and W 0 ∈ 𝑂 (𝑑) is orthogonal. Furthermore, denoting J A ∈ R 𝑛𝑑×𝑛𝑑 as the Jacobian matrix of the map 𝑓 conv (•, A), for any choice of W 0 , W 1 , there always exists a normalized adjacency matrix Â such thatEQUATION)where ∥M ∥ is the operator norm of matrix M .Proof. Note that for (a potentially normalized) adjacency matrix A, the Jacobian J A in the basis of vec(X) is equal toJ A = I ⊗ W ⊤ 0 + A ⊗ W ⊤ 1 .(22) Thus, for the map to be orthogonal, it must hold thatEQUATIONFirst, let A 0 be the empty graph on two nodes so A = 0, thenEQUATIONwhich means W 0 must be orthogonal if J A 0 is orthogonal. A simple instance with the desired structure can be constructed as follows: Consider the graph on two nodes, which is connected with adjacency matrix A 1 . Here, we haveEQUATION) Note that if W 0 is orthogonal, then W 1 = 0 if J A 1 J ⊤ A 1 = I.This proves the first part of the proposition.For the second part, we take traces and note thatEQUATIONIn contrast, Tr(I) = 2𝑑. Therefore, for any choice of the values of ∥W 0 ∥ 2 𝐹 , ∥W 1 ∥ 2 𝐹 , it must hold that eitherTrEQUATIONW.l.o.g. assume the first event holds. Denoting the singular values of a matrix J A 1 J A 1 -I as 𝑠 1 , . . . , 𝑠 2𝑑 , we haveEQUATIONRearranging, we obtain the final result. □One consequence of the above fact is that fully parameterized unitary graph convolution requires higher order powers of A to be implemented. This is essentially one reason why the exponential map (or some approximation thereof) is needed.We now show that the Rayleigh quotient defined in Definition 5 is invariant under unitary transformations. As before, this will follow virtually directly from the unitarity properties. We restate the proposition below. Proposition 6 (Invariance of Rayleigh quotient). Given an undirected graph G on 𝑛 nodes with normalized adjacency matrixA = D -1/2 AD -1/2 , the Rayleigh quotient 𝑅 G (X) = 𝑅 G ( 𝑓 Uconv (X)) is invariantProof. First, consider separable unitary convolution (Definition 1). The Rayleigh quotient takes the formEQUATION)The Frobenius norm is invariant under unitary transformations so the denominator ∥ exp(𝑖 A)XU ∥ 2 𝐹 = ∥X ∥ 2 𝐹 . For the numerator we have by the cyclic property of traceThe matrices exp(-𝑖 A), (I -A), exp(𝑖 A) all share the same eigenbasis so they commute and thusEQUATIONThus, we haveEQUATIONSimilarly for Lie unitary/orthogonal convolution (Definition 3), the isometry property guarantees ∥X ∥ 2 𝐹 = ∥ 𝑓 Uconv (X) ∥ 2 𝐹 . Furthermore, when viewed as a linear map in the basis of vec(X) ∈ C 𝑛𝑑 (i.e. entries of X viewed as a vector), 𝑓 Uconv can be written as a linear map of the form vecEQUATION)where W is the feature transformation matrix in Definition 3. Finally, note that exp(A ⊗ W ⊤ ) commutes with A soEQUATION)Multiplying the above by ∥X ∥ -2 𝐹 recovers 𝑅 G (X). □In contrast, we gave an example (Proposition 7) where the Rayleigh quotient decays with high probability for vanilla convolution, which we restate below.Proposition 7. Given a simple undirected graph G on 𝑛 nodes with normalized adjacency matrix A = D -1/2 AD -1/2 and node degree bounded by 𝐷, let X ∈ R 𝑛×𝑑 have rows drawn i.i.d. from the uniform distribution on the hypersphere in dimension 𝑑. Let 𝑓 conv (X) = AXW denote convolution with orthogonal feature transformation matrix W ∈ 𝑂 (𝑑). Then, the event below holds with probability 1exp(-Ω( √ 𝑛)):𝑅 G (X) ≥ 1 -𝑂 1 𝑛 1/4and 𝑅 G ( 𝑓 conv (X)) ≤ 1 -Tr( A 3 )EQUATION)By symmetry properties of the uniform distribution Unif(𝑆 𝑑-1 ) on the hypersphere,EQUATIONFurthermore, treating 𝑅 G (X) as a function of its node features x 1 , . . . , x 𝑛 where x 𝑖 = X ⊤ 𝑖,: , we have by the Azuma-Hoeffding inequality thatEQUATIONThe above can be shown by noting that the operator norm of I -A is bounded by two and changing the values of any x 𝑖 changes 𝑅 G (X) by at most 4/𝑛. Therefore, setting 𝜖 = 𝜀𝑛 -1/4 , we get that with probability 1exp(-Ω(𝜀 2 √ 𝑛)),EQUATIONFor 𝑅 G ( 𝑓 conv (X)), we haveEQUATION)Again, using the Azuma-Hoeffding argument as before, the numerator above concentrates around its expectation asEQUATIONA similar statement holds for the denominator ∥ AX ∥ 2 𝐹 . Furthermore, we have by symmetry properties of the distribution of X and linearity of expectation:EQUATIONCombining the above facts and applying the union bound, we have that with probability 1exp(Ω(𝜀 2 √ 𝑛)),EQUATIONIn the last equality, we use the fact that A has bounded degree and thus Tr( A 2 ) = Θ(𝑛). □We give a brief overview of matrix Lie groups and Lie algebras here and recommend [Hal15] for detailed exposition. Matrix Lie groups are subsets of invertible matrices that form a differentiable manifold formally defined below. Definition 10 (Matrix Lie groups [Hal15] ). A matrix Lie group is any subgroup of 𝐺 𝐿 (𝑛, C) with the property that any sequence of matrices M 𝑚 ∈ C 𝑛×𝑛 in the subgroup converge to a matrix M that is either an element of the subgroup or not invertible (i.e., not in 𝐺 𝐿(𝑛, C)).The orthogonal and unitary groups whose definitions are copied below meet these criteria.EQUATIONEQUATIONA crucial operation in matrix Lie groups is the exponential map. Given a matrix M : C 𝑛×𝑛 which is an endomorphism End(C 𝑛 ) on the vector space C 𝑛 , the exponential map exp : End(C 𝑛 ) → End(C 𝑛 ) returns another matrix (endomorphism) and is defined asEQUATION)The exponential map, as we will show below, maps from the Lie algebra to the Lie group. It also has an interpretation over Riemannian manifolds which we do not discuss here. We refer the reader to a textbook [Pet06] or prior work in machine learning [LCMR19] for further details of that connection. For compact groups, the exponential map is a smooth map whose image is the connected component to the identity of the Lie group [KJ08, Hal15] . The Lie algebra is the tangent space of a Lie group at the identity element. To see this, note thatEQUATIONandEQUATION)Note that exp(0) = I. The above gives us the Lie algebra to a given group.Definition 11 (Lie algebra [Hal15] ). Given a matrix Lie group 𝐺, the Lie algebra 𝔤 of 𝐺 is the set of matrices X such that 𝑒 𝑡X ∈ 𝐺 for all 𝑡 ∈ R.As an example, consider the unitary group where given a matrix U ∈ 𝑈 (𝑛) and X ∈ 𝔲(𝑛). Here, we haveEQUATIONFirst, let us recall the definition of the exponential map. Given a linear operator L : V → V which is an endomorphism End(V) on a vector space V, the exponential map exp :EQUATIONApplying the exponential map of a linear operator to a given vector in a vector space of dimension 𝑛 to computer precision can scale in practice as 𝑂 (𝑛 3 ) similarly to performing an eigendecomposition. In fact, for matrices M ∈ C 𝑛×𝑛 which have an eigendecomposition M = U DU † (for skew-Hermitian matrices M ∈ 𝔲(𝑛), the spectral theorem guarantees the existence of an eigendecomposition.), exp(M ) = U exp(D)U † where exp(D) can be easily implemented by performing the exponential elementwise on the diagonal entries of D. However, in practice, we can exploit approximations which avoid expensive operations such as eigendecompositions.Taylor approximation Often the simplest and most efficient approximation is a 𝑘-th order truncation to the Taylor series expEQUATIONwhere conv 𝑝 : V → V indicates the composition of the operator conv 𝑝 times. By Taylor's theorem, one can bound the error in this approximation as exp[conv] (X) - 𝑘 ∑︁ 𝑝=0 1 𝑝! conv 𝑝 (X) 2 ≤ 𝑂 ∥conv∥ 𝑘+1 ∥X∥ 2 (𝑘 + 1)! ,where ∥conv∥ denotes the operator norm of conv. Therefore, error decreases exponentially with 𝑘. In practice, we find setting 𝑘 = 12 is sufficiently accurate for the GNNs we train in our experiments. Of course, 𝑘 can be increased in settings with larger graphs or more training time where instabilities are likely to arise. We also by default perform unitary convolution with normalized adjacency matrices D -1/2 AD -1/2 (D is the diagonal degree matrix) so that the operator norm of the linear convolution operation is bounded by one.Padé approximant Padé approximations are optional rational functional approximations of a given function up to a given order. In unitary settings, Padé approximants are often preferred because they return a unitary operator. In contrast, Taylor approximations discussed previously only return a linear operator which can be made arbitrarily close to a unitary operator. Padé approximants of order 𝑘 take the form expEQUATIONwhere 𝑝 𝑘 , 𝑞 𝑘 are degree 𝑘 polynomials. The degree 𝑘 Padé approximant agrees with the Taylor expansion up to order 2𝑘. The first order Padé approximant, also known as the Cayley map, has been used in prior neural network implementations [WPH + 16, HWY18, TK21] and takes the formEQUATION)One practical drawback of Padé approximants are that they typically require implementations of a matrix inverse or approximations thereof which render this more challenging to implement. In our setting, we found that the Taylor approximation sufficed in accuracy and stability so we did not implement this.Finally, we should remark that Padé approximants can be used to pre-compute matrix exponentials to computer precision [Hig09, AMH10, AMH09]. [LCMR19] use these techniques to parameterize unitary layers accurately for RNNs. One could pre-compute matrix exponentials of adjacency matrices in our setting to speed up training, though we have not implemented this in our experiments. We leave this for future work.We briefly mention here for sake of completeness a different approach to approximating the exponential map based on the Baker-Campbell-Hausdorff formula [Hal15] which states that for matrices X, Y ∈ C 𝑛×𝑛 :EQUATION)The Lie-Trotter approximation uses the above fact to implement the matrix exponential of a sum of matrices as a product of matrix exponentials over elements of the sum. This method is often used in high dimensional spaces where the output of the exponential map on each Lie algebra generator is known. The first order expansion takes the form [Tro59, CST + 21]EQUATION)where 𝑘 is a positive integer controlling the level of approximation. The above is accurate to order 𝑂 (1/𝑘), and higher order accurate schemes exist. These methods are commonly used in machine learning, quantum computation, and numerical methods10 [CST + 21, MGL + 23, Llo96, Str68, MQ02]. One advantage of this approach is that the approximation applied to an operator in the Lie algebra always returns an element in the Lie group since the approximation is a product of Lie group elements themselves.In the main text, we described a generalized procedure to implement unitary convolution with parameterized operators in the Lie algebra of a group. We complement that with a mostly informal discussion on how to implement unitary convolution in the Fourier domain. The algorithms are summarized in Algorithm 3 and Algorithm 2.Before detailing the algorithm in Algorithm 3, we provide a brief overview of Fourier-based convolutions over arbitrary groups, which generalizes the classical Fourier transform to that over arbitrary groups. We recommend [S + 77] for a more complete and rigorous background. Throughout this, we will assume that groups are finite, though this can be generalized to other groups, especially those that are compact. Input: equivariant linear operator L ∈ C 𝑛 → C 𝑛 Input: vector x ∈ C 𝑛 1: L = 1 2 (L -L † ) (skew(𝑔) = 𝑄 𝜌(𝑔)𝑄 -1 .The Fourier transform of a function 𝑓 : 𝐺 → C with respect to a set of irreducible representations (irreps) {𝜌 𝑖 } 𝑚 𝑖=1 is given by: fEQUATIONFor abelian groups these irreps are all one dimensional. The set of 𝜌 𝑖 (𝑢) for the cyclic group for example correspond to the entries of the discrete Fourier transform matrix. For non-abelian groups, there always exists an irrep which is at least of dimension 2.In Algorithm 3, we denote the Fourier transform F : C 𝑛 → ⊕ 𝑚 𝑖=1 C 𝑑 𝑖 ×𝑑 𝑖 as a map that takes in the inputs to a function 𝑓 and outputs a direct sum of the Fourier basis of the function f over the irreps. Given two functions 𝑓 , 𝑔 : 𝐺 → C, their convolution outputs another function ( 𝑓 * 𝑔) : 𝐺 → C and is equal toEQUATION)In the main text, we describe these operations as vector operations over a vector space C 𝑛 . Setting 𝑛 = |𝐺 | and taking the so-called regular representation as maps acting on C |𝐺 | can recover the form in the main text. Finally, the Fourier transform of their convolution is the matrix product of their respective Fourier transforms:EQUATIONwhere 𝜌 𝑖 𝑚 𝑖=1 are the irreps of 𝐺. The above assumes that all the functions have input domain over the group. This is not strictly necessary and generalizations exist which map functions on homogeneous spaces to the setting above [KT18] .The general implementation of Fourier convolution is given in Algorithm 3. Here, one employs a Fourier operator which block diagonalizes the input into its irreducible representations or some spectral representation. Then, applying blocks of unitary matrices in this representation and inverting the Fourier transform implements a unitary convolution. The details will depend on the particular form of the Fourier transform and irreducible representations. This method is often preferred when filters are densely supported and efficient implementations of the Fourier transform are obtained. Previous implementations have been designed for CNNs [TK21, KBLL22] . Example 3 (Convolution on regular representation (Fourier basis)). Continuing the previous example, assume the group 𝐺 has unitary irreducible representations 𝜌 1 , . . . , 𝜌 𝑚 (irreps) where 𝜌 𝑖 : 𝐺 → C 𝑑 𝑖 ×𝑑 𝑖 . The group Fourier transform maps input function 𝑥 : 𝐺 → C to its irrep basis as [FH13, KT18] xEQUATIONand group convolution is now block diagonal in this basis(𝑚 ★ 𝑥) (𝜌) = x(𝜌) m(𝜌) † .(60) Implementing unitary convolution requires that m(𝜌) is unitary for all irreps 𝜌.Handling complex numbers and enforcing isometry in neural networks requires changes to some standard practice in training neural networks. We summarize some of the important considerations here and refer the reader to surveys and prior works for further details [BQL21, TK21, LCMR19, KBLL22].Handling different input and output dimensions Unitary and orthogonal transformations are defined on input and output spaces of the same dimension. Maintaining isometry for different input and output dimensions formally requires manipulations of the Stiefel manifold as studied in various prior works [LCMR19, LFT20, NA05, JD15]. When the input dimension is less than that of the output dimension, one simple way to implement semi-unitary or semi-orthogonal convolutions via standard unitary layers is simply to pad the inputs with zeros to match the output dimensionality. We also often will simply use a standard (unconstrained) linear transformation to first embed inputs in the given dimension that is later used for unitary transformations.Nonlinearities For handling complex numbers, one must typically redefine nonlinearities to handle complex inputs. We find that applying standard nonlinearities separately to the real and imaginary parts works well in practice in line with other works [BQL21] . To enforce isometry as well in the nonlinearity, we use the GroupSort : R 2 → R 2 activation [TK21, ALG19], which acts on a pair of numbers as GroupSort(𝑎, 𝑏) = (max(𝑎, 𝑏), min(𝑎, 𝑏)), (61) and is clearly norm-preserving. To apply this nonlinearity to a given layer, we split the channels or feature dimension into two separate parts and apply the nonlinearity across the split.Initialization Given the constraints on unitary matrices and skew-Hermitian matrices, prior work has proposed various forms of initialization that meet the constraints of these matrices. One strategy that has been effective in prior work and proposed in [HSL16, HWY18] is to initialize in 2 × 2 blocks along the diagonal. For skew symmetric matrices, one way of achieving this is to initialize 2 × 2 blocks as 0 𝑠 -𝑠 0 ,where 𝑠 ∼ Unif (-𝜋, 𝜋) for example [HSL16] .Directed graphs Directed graphs present a challenge because their adjacency matrix is not guaranteed to be symmetric. Given an adjacency matrix A ∈ R 𝑛×𝑛 of a directed graph, one simple way to proceed is to split the directed graph into its symmetric and non-symmetric parts A = A sym +A nonsym .Here, we assume that for matrix entry A 𝑖 𝑗 , either A 𝑖 𝑗 = A 𝑗𝑖 or A 𝑖 𝑗 A 𝑗𝑖 = 0 (i.e. one of the transposed entries is zero). Then, we setEQUATIONandEQUATIONFinally, one can then perform graph convolution with the skew-Hermitian matrixEQUATIONWe do not work with directed graphs in our experiments and have not implemented this in our code.Fig. 3 shows additional results for the toy model considered in the main text. As a reminder, this task is to learn the graph distance between pairs of randomly selected nodes in a ring graph of 100 We show here the performance of various models with 5, 10, or 20 layers. The unitary GCN is the only message passing architecture that achieves stable performance with added layers and can learn the task. Apart from message passing architectures, global transformer architectures like GPS can learn the task when given Laplacian positional encoding. The trivial performance corresponding to outputting the average output is shown as a dotted horizontal line.nodes. Message passing architectures need at least 50 sequential messages to fully learn this task.As layers are added to the unitary GCN, the network is able to learn the task better. This is in contrast to other message passing architectures which perform worse with additional layers. Apart from message passing architectures, strong performance is also achieved by transformer architectures like GPS with global attention. The networks we study are the vanilla GCN (Residual GCN includes skip connections) [KW16] , graph attention network [VCC + 17], spectral convolution [ZK20] , transformerbased GPS [RGD + 22], and graph-coupled oscillator network (GraphCON), which is a discretization of a second-order ODE on the graph [RCR + 22]. Hyperparameters and additional network details are reported in App. G.We perform experiments on the ENZYMES, IMDB-BINARY, MUTAG, and PROTEINS tasks from the TU Dataset database [MKB + 20]. As can be seen in table 3 , with the exception of IMDB, a GCN with UniConv layers outperforms all message-passing GNNs tested against on all datasets, by margins of up to 18 percent. We follow the training procedure in [NHN + 23] and report results for GCN and GIN from [NHN + 23]. For GAT and unitary GCN, we tune the dropout and learning rate as detailed in App. G. All results are accumulated over 100 random trials. In addition, in Fig. 4 , we show that the unitary GCN maintains its performance over large network depths in contrast to other message passing layers. The deteriation in performance of conventional message passing layers is a likely a consequence of over-smoothing which we show does not occur in unitary graph convolution [CW20, RBM23] .The dihedral group 𝐷 𝑛 is a group of order 2𝑛 describing the symmetries (rotations and reflections) of a regular polygon. Its elements are generated by a rotation generator 𝑟 and reflection generator 𝑠:EQUATIONIn this task, analaogous to the graph distance task in Sec. 5, the goal is to learn the distance between two group elements 𝑔, 𝑔 ′ in the dihedral group 𝐷 𝑛 . Formally, we aim to learn a target function 𝑓 : R | 𝐷 𝑛 | → Z mapping inputs of dimension 2𝑛 (the vector space of the regular representation) to Table 3: Comparison of Unitary GCN with Lie UniConv layers (Definition 3) with other GNN architectures on the TU datasets. Each complex number is counted as two parameters for our architectures, except for wide Unitary GCN which counts a complex numbers as one parameter so that the width of hidden layers roughly matches that of vanilla GCN.Figure 4 : Test accuracies on Mutag for GCN, GIN, GAT, and a GCN with UniConv layers with increasing number of layers. Except for the unitary network, all other message passing architectures collapse to trivial accuracy levels as the number of layers increases. positive integers. The input space R |𝐷 𝑛 | is indexed by group elements and convolution is performed on the regular representation of the group. Each data pair (x 𝑖 , 𝑦 𝑖 ) 𝑖 is drawn as follows:• Draw two group elements 𝑔, 𝑔 ′ ∼ Unif (𝐷 𝑛 ) from the uniform distribution over the group 𝐷 𝑛 without replacement. • Set x 𝑖 = e 𝑔 + e 𝑔 ′ where e 𝑔 ∈ R | 𝐷 𝑛 | is a unit vector whose entries are all zero except for the entry corresponding to operation 𝑔 set to one. • Set 𝑦 𝑖 = 𝑑 𝐷 𝑛 (𝑔, 𝑔 ′ ) ∈ Z where 𝑑 𝐷 𝑛 (𝑔, 𝑔 ′ ) indicates the number of applications of the generators that need to be applied to go from 𝑔 to 𝑔 ′ , i.e.EQUATIONFig. 5 plots the performance of three different networks in this task. The vanilla network performs standard convolution in each layer where the group convolution is parameterized over elements of the generator. The residual network is the same architecture but includes skip connections after convolutions. The unitary network applies the exponential map to a skew-Hermitian version of the group convolution as outlined earlier. The unitary network is able to learn this task in fewer layers and with added stability. No hyperparameter tuning was performed in these experiments. Adam optimizer parameters were set to the default setting. All networks have a global average pooling after the last convolution layer followed by a 2 layer MLP to output a scalar. The number of channels is set to 32 for each convolution layer.To compare orthogonal (real-valued) and unitary convolutional layers, we repeat our experiments on Peptides-func and Peptides-struct from the main text. We report the results in are not included in this ablation as no edge feature aggregator was included in the architecture. GCN with OrthoConv layers achieves performance levels similar to UniConv without using complex numbers in its weights.For LRGB datasets, we evaluate our GNNs using the GraphGym platform [YYL20] . In Table Our implementation of unitary graph convolution does not take into account edge features. Thus, for datasets with edge features, at times, we included a single initial convolution layer using the GINE [XHLJ18] or Gated GCN [BL17] architecture which aggregates the edge features into the node features. For instances where this is done, we indicate the type of layer as an "edge aggregator" hyperparameter in the tables below.All our experiments were trained on a single GPU (we used either Nvidia Tesla V100 or Nvidia RTX A6000 GPUs number of parameters as the dimension of the Lie Algebra (i.e. the number of parameters to fully parameterize the Lie algebra).Toy model: graph distance All networks consist of the stated number of convolution or transformer layers with feature dimension 128 followed by a global average pooling over nodes. Pooled features are then passed through a single hidden layer MLP with 128 dimension width. For networks with 5, 10 and 20 layers respectively, we set the learning rate to 0.0007, 0.0003 and 0.0001 respectively. Networks are trained with the mean average error (L1) loss. Activation functions for all networks are set to GELU apart from the unitary networks where we choose the norm-preserving activation GroupSort [TK21, ALG19] (see App. E). To ensure unitary layers are also truly normpreserving, we do not include a trainable bias in these layers.Peptides For the Peptides experiments, the training procedure follows that of [TRRG23] . Unless otherwise stated, we use the Adam optimizer with learning rate set to 0.001 and a cosine learning rate scheduler [KB14] . Most hyperparameters were taken from [TRWG21] . For our purposes, we tuned dropout rates in the set {0.1, 0.15, 0.2} and tested networks of layers in the set {6, 8, 10, 12, 14} selecting the hidden width to fit within the parameter budget. For larger levels of dropout, we found that we needed to train the network for more epochs to achieve convergence. Final hyperparameters are listed in Table 6 .COCO-SP and PascalVOC-SP For the COCO-SP and PascalVOC-SP experiments, our training procedure again follows that of [TRRG23] . We use the Adam optimizer with learning rate set to 0.001 and a cosine learning rate scheduler [KB14] , unless otherwise stated. For our purposes, we tuned dropout rates in the set {0.1, 0.15, 0.2} and tested networks of layers in the set {6, 8, 10, 12, 14} selecting the hidden width to fit within the parameter budget. For larger levels of dropout, we found that we needed to train the network for more epochs to achieve convergence. Final hyperparameters are listed in Table 6 .Our experiments are designed as follows: For a given GNN model, we train on a part of the dataset and evaluate performance on a withheld test set using a train/val/test split of 50/25/25 percent. For all models considered, we record the test set accuracy of the settings with the highest validation accuracy. As there is a certain stochasticity involved, especially when training GNNs, we accumulate experimental results across 100 random trials. We report the mean test accuracy, along with the 95% confidence interval. For the TU Datasets, our training procedure and hyperparameter tuning procedure follows that described in [NHN + 23]. Namely, we fix the depth and width of the network as in [NHN + 23] and tune the learning rate and dropout parameters. Final hyperparameters are listed in Table 7 . To ensure fairness and comparability, we conducted the same hyperparameter searches with the baseline models we compare against, but found no improvements over the numbers reported in previous papers. We therefore report their (higher) numbers instead of ours for these models. To output a single number, we use a single convolution layer (SAGE convolution) to map from the higher dimensional space to a single number for each node. [PKD + 23] differs in that they have an MLP in between layers; we opt for a more simple approach. Final hyperparameters are listed in Table 8 . To ensure fairness and comparability, we conducted the same hyperparameter searches with SAGE and GCN, but generally found no significant differences with respect to the numbers reported in previous papers. We therefore report their usually higher numbers instead of ours for these models.We list below the licenses of code and datasets that we use in our experiments.  
experiment statistical significance: 
 Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?Answer: [Yes] Justification: Error bars included in all experiments and tables. Details are provided in App. G.• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Justification: A single GPU is needed to run an individual experiment; details are provided in App. G.• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines?Justification:Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction). 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?Answer: [NA] Justification: This work presents foundational research on improving learning on graphs and groups. Though this work has impact on the research community in machine learning, we do not see any clear links to questions of privacy, misuse, fairness, or other major societal impacts here.Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA] Justification: Our work is largely foundational and theoretical. All models and datasets are relatively small and scientific in nature.• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort. 
licenses for existing assets: 
 Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?Answer: [Yes] Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained? Answer: [NA] Justification: Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.10In numerical methods, they fall under the umbrella of splitting methods or Strang splitting. 
 
 Appendix:  


Decision: This paper introduces unitary convolutional layers for Graph Neural Networks (GNNs), with the primary aim of mitigating issues like gradient collapse, explosion, and over-smoothing in GNNs. The experimental results demonstrate that these layers can effectively aggregate information across long distances and achieve strong performance on real-world benchmarks.

The paper was well received by the reviewers, who appreciated the formalism and theoretical analysis of the proposed operators. The initial scores generally improved after the authors' rebuttal, which clarified some misunderstandings and provided additional experimental results.

The main limitation pointed out by the reviewers concerns the writing and organization of the paper. This is the main reason for the “borderline reject” score from one of the reviewers, who, however, decided to lower his confidence score. The authors have proposed addressing these issues with revisions that should be applied in the camera-ready version. Additionally, the introduction of two unitary convolution operators caused some confusion, and the authors should better motivate their inclusion.

Reviews: review: paper_summary: paper_summary: This paper proposed a new convolution operator UniConv on graphs and general tasks. Both graph convolution and general form convolution form are provided based on the group theory. With the toy experiment, it shows the advantages of the proposed method clearly. Based on both theoretical and experimental analysis, the proposed network is effective in long range graphs and heterophilic graphs.
strength_weakness: strength_weakness: 1. Though the general form of convolution is provided. It is not evaluated in experiment. I will suggest putting these content in the appendix in this case. Otherwise, it is better to provide the corresponding experimental results.

2. The computational cost is higher compared to other methods, but there are no experimental results to show the important details.
questions: questions: See weakness.
limitations: limitations: N/A
paper_summary: paper_summary: This paper introduces convolutional layers for GNNs that are based on group-convolutions. Some of the proposed layers have interesting theoretical properties such as avoiding oversmoothing or vanishing gradients. Experimentally, these layers are able to aggregate information across long distances and achieve strong results on real-world benchmarks.
strength_weakness: strength_weakness: - **(W1 - Clarity and Significance)** The reader is not guided through the paper and it is difficult to understand the flow of ideas. I will give some examples of this:

   - In Section 3.1 (Lie) UniConv layers are introduced which is the main contribution of the pape. However, to me it was not clear that this is the central definition. 

   - In the following Section 3.2, the paper generalizes the setting of unitary convolutions and gives several examples. Yet, to me it is unclear what this section intends to communicate and how it is relevant to the rest of the paper. 

   - Another example of this is the “Vanishing/Exploding gradients”-part of Section 4, here dynamical isometry is defined as a measure as a property to analyze exploding gradients and an example is given that a combination UniConv together with a yet completely undefined activation function (GroupSort) is “perfectly dynamically isometric”. This raises many questions to the reader: Why is an example important enough to warrant the definition of Dynamical isometry? How does this relate to the main contribution of the paper? Are vanishing / exploding gradients a problem with GNNs? What is GroupSort and is it actually used in practice? 

- **(W2 - Clarity)** The experiment section is difficult to read and does not properly explain the experiments or the design decisions. For example, for the graph distance dataset the model is not mentioned in the text. In this case, the model is GCN with the novel Lie UniConv layer (Figure 2). However, it is not explained why the Lie UniConv layers are used and not the UniConv Layers. Similarly, for the experiments on LRGB UniConv layers are used but it is not explained why Lie UniConv layers are not tried.

**To sum up.** To me it seems like this could be a good paper with interesting ideas and solid experiments. However, the writing significantly holds it back and I think more time is required to allow  the authors to rewrite their paper. Furthermore, I think that there is a lot of useful information hidden in the appendix and it might make sense to focus the paper on the main idea (Section 3.1) and move other parts to the appendix to make more space.  As I do believe that this paper requires significant changes to the text I thus vote to reject.
questions: questions: Most of my questions are already part of the weakness section.

- **(Q1)** (Lie) UniConv are implemented as a $k$-th taylor approximation. What is the value for $k$? It seems to me like the ability to model long-range interactions is primarily dependent the value of $k$ and the number of layers $L$, so knowing this value seems important.
limitations: limitations: Yes.
paper_summary: paper_summary: In this paper the authors propose a mathematically consistent treatment of unitary convolutions in the context of neural networks defined on graphs and groups. In the graph neural network context, the main idea is to transform a standard linear convolutional layer $f_\text{conv} = \mathbf{AXW}$ with $\mathbf{A}$ the adjacency matrix, $\mathbf{X}$ the signal input matrix and $\mathbf{W}$ the parameterize into a linear layer $f_\text{Uconv}$ that is unitary (norm preserving), taking into account both the matrices $\mathbf{A}$ and $\mathbf{W}$ (previous literature focus on $\mathbf{W}$ or perform this in a less natural way (with worse computational complexity) [QBY24]. The main idea revolves around the following two steps: (1) project the operation into matrix Lie algebra (e.g., the space of skew-symmetric matrices) (2) apply the exponential map on the projected operation to obtain an unitary operation. In the graph setting the authors do this in two ways: the first is called "Separable unitary graph convolution (UniConv)" and takes the form $f_{\text{Uconv}} (\textbf{X}) = \exp(i\textbf{A}𝑡)\textbf{X}\textbf{U}, \textbf{U} \textbf{U}^{\dagger} = \textbf{I}$. The second is called "Lie orthogonal/unitary graph convolution (Lie UniConv)" and is defined by applying an exponential map on a convolution with skew-symmetric weight matrices in Lie algebra, having the advantage of being parameterizable over the reals only. The authors then discuss how to approximate the exponential map (they resort to simple Taylor truncation).  Subsequently, the authors develop the theory in the case of general groups and show to ways of obtaining unitary convolutions: the first case generalizes the two methods defined for graphs, by first skew-symmetrizing the equivariant map and then applying exp; the second is defined as a (block-diaognal) convolution in the Fourier domain followed by inverse transform. Following the authors show basic properties of the proposed method (invertibility, isometry and equivariance) and prove favorable results which mainly follow by isometry (avoiding oversmoothing, avoiding vanishing/exploding gradients). In the experimental section, the authors perform experiments on toy ring graph (on Lie UniConv) and the standard benchmarks LRGB and Heterophilous Graph Dataset (on UniConv) (also showing group convolution in Appendix on Dihedral group $D_n$).
strength_weakness: strength_weakness: - There aren't many experiments on Lie UniConv. From what I understand, the only experiment with this layer is the graph ring experiment in Figure 2. Also it would have been interesting to see comparison between the two type of convolutions (UniConv vs Lie UniConv), with the more realistic datasets, in order to asses which one performs better in practice.

- The part on unitary maps in Fourier basis feels somehow detached from the general discourse of the paper. I believed the author added it for completeness but (1) the theory to understand it is relatively advanced (had to study externally some representation theory to understand well the definition of Fourier transform on finite groups); a short primer on representation theory in the Appendix, similar to the one for Lie groups / algebras, would have been of great help. (2) There are no experiments in this setting. (3) some parts are a little bit shaky, for example in: ```One can also generally implement convolutions in the (block diagonal) Fourier basis of the graph
or group (Algorithm 2). Here, one employs a Fourier operator which block diagonalizes the input into its irreducible representations or some spectral representation.```, if one considers using the eigendecomposition of some graph Laplacian as a spectral basis, the (graph) Fourier transform is not given by a set of matrices varying over irreps but a signal of the same dimensionality of the input function. I do not think one can use the same approach there.

- Sometimes authors leave facts implicit and the writing could be improved in key places. For example I had difficulty understanding why Lie UniConv returns an unitary map before looking at the tensor product version of $g_{\text{conv}}$ (Eq. (36), Proposition 6). Also the statements in Propositions 4 and 7 are not really clear from a first reading, some more intuition about the statements should be given.
questions: questions: - I don't see how the graph convolutions defined are particular cases of the group convolutions described. Is it the symmetric group S(n) acting on the graph?
- Line 122: Why UniConv has a "tensor product nature"? From what I understand the tensor product is more related to the Lie UniConv because you can write it as $\exp(\mathbf{A} \otimes \mathbf{W}^T)$.
- Line 186: "We set input/output representations to be equal so that the exponential map of an equivariant operator is
itself equivariant": Assuming this does not limit the class of groups or transforms we can use? Not all groups have irreps of same dimensions. Or we say that the transforms are equivariant only in this sense? Also why exponential map is equivariant only if the input output representations are equal?
- Line 212: is (I − A) a type of graph Laplacian?
- Line 272: Why authors don't compare to new architectures as well on this task, like in Table 1? Also I see in Appendix that GPS it is performing pretty well w.r.t. Unitary GCN but it was not put in this graph (maybe to show better improvements?).
- Line 763: I think authors should say that in case of Algo 2 this follows because Fourier transform is unitary. It is not clear if and why this should be the case in finite group setting.
- Line 766: "For unitary convolution in the Fourier domain (Algorithm 2), equivariance follows from convolution theorems where linear operators appropriately applied in the block diagonal Fourier basis are equivariant". This is too succinct and should be more expanded formally.
- Line 821: missing trace on following lines and last line should not be put there ($= R_{\mathcal{G}(\mathbf{X})}$)
limitations: limitations: The authors have discussed limitations of their contribution adeguatly.
paper_summary: paper_summary: The paper proposes to use so-called "unitary group convolutions" in graph neural networks, with main motivation to do steps to overcome gradient collapse or explosion and oversmoothing effects in graph neural networks under equivariance constraints. This fits within a more general line of work in which unitary constraints are applied for regularizing various learning architectures.
strength_weakness: strength_weakness: Edit: after rebuttal, the main issue of originality is 100% solved, and the points in questions 3 and 5 are fully solved. The authors have proposed better experiments. This makes my concern in high part mute, the only small (and out of the scope of this work) point is that it will be great to see more verifications of oversmoothing control of beyond the (already convincing) Rayleigh Quotients preservation, which is proved in the work. As said above, this is for future work, and requires much more space, so the paper is very strong as is.
---------------------------

The originality of this approach is a bit limited, since normalizing message passing is not a very novel idea. 
Doing it with unitary constraints has not been explored before for GNNs. It has been explored in other settings though.

The mathematical setup is sketchy at times, it needs better care. In particular, it is not clear if the proposed convolution is really a unitary operation (see question 3 and the crucial question 5 in "questions" section, especially).

The experiments show that normalizing is a good idea, but they don't fully convince that unitary constraints is the best way to normalize, since comparisons are done mainly against non-normalized GNNs. One has the impression that this is mainly a first exploration of a theme, and it would need more careful benchmarking, since the computational cost for implementing this particular normalization is not always negligible.
questions: questions: What's the relative benefit compared to more naive normalizations? 


Minor correctionsl/doubts:
1) lines 29-31: "it has been widely observed that group-convolutional networks suffer from instabilities as their depth increases" -- if it's widely observed, please indicate a few references where this was observed?

2) line 103: "linear convolution" -- as opposed to what other kind of convolution? isn't convolution always linear?

3) line 144: the following is an important point/issue. 
"When W contains only real-valued entries, the above returns an orthogonal map" -- I don't follow this, why is that? this statement would require a proof, because it far from obvious, and in general it is false. Recall that a product of matrix exponentials is in general not the exponential of the product, so $exp (g_{\mathrm{conv}})(X) \neq e^A X e^{W^*}$ in this case.

4) I don't see what's the importance / relevance of examples 1 and 2 from Section 3. I think these are standard examples, they can be skipped or moved to the appendix. 

5) lines 217-218: this probably the most important theory issue of the paper.
 "the Rayleigh quotient is invariant to unitary transformation giving a proof that unitary graph convolution avoids oversmoothing" 
One can vaguely agree that Rayleigh quotient collapse may in some way be correlated to oversmoothing, but saying that Prop 6 "gives a proof that unitary convolution avoids oversmoothing" is unjustified without further proofs/details. In particular some form of a result saying that "Rayleigh quotients do not collapse implies there is no oversmoothing" is missing. Without this, the main claim of the paper hinges on "vague intuition" about RQ's and the paper is much less valuable at a theory level.

6) Definition 6 (dynamical isometry): what iis the connection of this definition with the definitions from section 3.1? This part on dynamical isometry is at present very weakly (if at all) connected to the core of the paper, and I did not understand how the authors are linking it to the main definitions. See question 3 related to this.

7) The experiments show performance on some benchmarks, but I would have been more convinced if the unitary convolution proposed method would have been compared to other simpler normalization methods which have similar effects on oversmoothing and gradient vanishing/explosion.
limitations: limitations: Main limitations are already included in the "questions" section, esp questions 3,5,7.
paper_summary: paper_summary: The paper introduces two unitary graph convolution operators (UniConv and Lie UniConv), and studies their performance and ability to avoid over-smoothing even in deep Graph Neural Networks. UniConv (short for separable unitary convolution) takes the form f_{Uconv} = exp(iAt)XU (with A the adjacency matrix, i the imaginary unity, and U a unitary operator), while Lie UniConv (short for Lie Unitary Convolution) is defined as f_{Uconv} = exp(g_conv)(X) (with g_conv(X)=AXW, and W a skew-symmetric matrix). In both cases, the exponential map is approximated with a truncated Taylor series to contain complexity and avoid any form of eigen-decompositon. In experimental evaluation, GNNs implemented with UniConv achieved good performance on a variety of real-world datasets, while Lie UniConv was the only model to succeed in a toy example exhibiting long-range dependencies when compared to classic GNNs (such as GCN, GAT and Residual GCN).
strength_weakness: strength_weakness: While I generally appreciated the paper (hence my slightly positive score), I’m not fully convinced by the Lie UniConv approach. The authors propose the two methods I highlighted above as a solution to avoid overfitting. As overfitting is caused by a contraction of all frequencies associated to all but the lowest frequency of the Laplace operator, converting the adjacency matrix of the provided graph into a unitary operator (as it is done in UniConv) makes sense to avoid the issue (as all eigenvalues of the obtained diffusion operator now have norm equal to 1). However, the Lie UniConv approach doesn’t go in that direction, but rather applies an exponential map on operator g_conv to obtain an orthogonal map. In this case, why the exponential of such operator is an orthogonal map is not specified in the paper, and I struggle to understand the rationale behind this approach. In this direction, in the experimental evaluation of the paper, I didn’t see any result using Lie Uniconv besides the experiment of the synthetic long rage benchmark I highlighted above. On top to this, on such benchmark, UniConv is not applied in the paper. If this was achieving comparable, or even superior, results to Lie Uniconv, then there would not be any result in the paper showing the benefit of using the Lie Uniconv method.

I would kindly ask the authors if they could clarify my doubts above, and present additional results with Lie Uniconv in the rebuttal (unless I missed them somehow in the paper).

In addition to this, I would like to point out that the results presented in Table 2 are incomplete, and while the proposed Uniconv outperforms the reported methods, it is possible to see better performing solutions on Roman E in “A critical look at the evaluation of gnns under heterophily: Are we really making progress?”. I would thus ask the authors to report all the relevant results from the paper, to ensure a fair comparison.
questions: questions: See Weaknesses
limitations: limitations: Yes

Comments: title: All doubts clarified, thanks
comment: Dear authors, I am sorry for the delay, I submitted a reply with the wrong reader-group. You should now be able to see my reply above.
I thank the authors for their thorough rebuttal. The additional experiments do strengthen the paper and including the results for both architectures does in my opinion simplify the experiment section.

> Can the reviewer expand on what they mean by central definition? From our understanding, the definitions of UniConv operations are central to the paper. 

Yes, they are. However, for some reason this was not immediately clear for me when reading this section. Maybe a bit more guiding of the reader is required.


> If we understand the reviewer correctly, they have requested that we focus more on section 3.1, defer parts of section 3.2 to appendix, and make more space for clarifying some points such as when to use Lie UniConv vs. UniConv. We will make these changes in the revised manuscript. Are there other parts of the paper that the reviewer feels should be moved from the appendix to the main text or vice versa?

You understand me correctly, I think that moving 3.2 to the appendix is a good idea. This  should allow you to focus more on your essential results (key definitions + theorems/propositions) and intuitively explain what they mean and why they matter. 

> We are hesitant to make significant changes to the writing and format given this concern was not shared by all the reviewers, but are of course open to more feedback.

I agree with this sentiment. While I do think that my issues with the presentation / clarity persist (and this is somewhat mirrored by feJe), I cannot ignore that this does not seem to be an issue for the other reviewers. Furthermore, upon a further look I have to concede that this paper does contain many non-trivial theoretical insights. I will slightly increase my score but adjust my confidence downward.
I thank the authors for addressing different points such as additional experiments on Lie UniConv and improving the clarity of the paper where I have pointed out different misunderstandings. The papers is a solid contribution and should definitely be accepted.
Thank you once again for your valuable feedback and the constructive discussion. We appreciate your decision to raise your score to an 8.

The score change does not appear on our side. Would you mind verifying that the updated score has been recorded in the system?
Thank you again for your response. My comment was mostly dictated from the fact that I find somehow weird to have two methods proposed in a paper, when one produces consistently better performance than the other, and we have no empirical evidence of any advantage proposed by the method performing the worst. That is why I was asking for running times comparing the two methods you proposed, I was looking for proof of a good reason why one would consider implementing Lie UniConv over Uniconv. 

This said, you are right pointing out that in any case the paper proposes some methodology that outperforms prior art. I'll raise my score to 6 in view of that.
I understand now. I was interpreting $g_{conv}^k(X)$ as $(g_{conv}(X))^k$ and now that I noticed the difference everything is clear.

My main concerns have been all lifted this is the case for concern for novelty (which is now "reverted": I believe the paper has strong novelty), and the topics of questions 3 and 5 (technical points). 

The concern that Rayleigh Quotients may be only one of several possible measures for oversmoothing is still valid for me, but I think it can/should be considered in future work, and due to understanding the difference between normalizing 2-norm and imposing unitary transformations, I am confident in the result being "fortifiable". What I mean is mainly to prove nondegeneracy according to more of the used metrics for oversmoothing. 

Given the above, I think it's fair to raise my score to 8.
No problem; we are happy to clarify. Let us answer your questions in order.


1. Formally, unitary transformations must be invertible so the input and output dimension must be equal. Thus the row and column dimensions of $W \in \mathbb{C}^{b\times b}$ must be equal and similarly for $A \in \mathbb{C}^{a \times a}$. However, this does not mean that $X$ has equal row and column dimension. The map $f_{Uconv}$ maps matrices of size $a \times b$ to matrices also of size $a \times b$ so dimensionality is preserved in this sense.
As noted in the draft though, this does not limit implementation: for $W$ with arbitrary row and column dimension, one can also work within the Stiefel manifold or even more simply pad inputs and outputs so that dimensionalities are equal to actually implement things in practice (see appendix D). 
2. It does not need to hold that $AXW =- (AXW)^\dagger$.  In fact, this will not hold true in general since the row and column dimensions of $X$ are not equal so this equation formally may not even make sense. What does hold true however is that the linear operator $g_{conv}$ is in the Lie algebra. As a reminder, exponential map is equal to
$$\exp(g_{conv})(X) = \sum_{k=0}^\infty \frac{g_{conv}^{(k)}(X)}{k!} = \sum_{k=0}^\infty \frac{A^kXW^k}{k!}. $$
Viewed in the vectorized form from our previous comment this is equivalent to:
$$vec[\exp(g_{conv})(X)] = \sum_{k=0}^\infty \frac{(A \otimes W^T)^k}{k!} vec[X],$$
so as stated earlier, one can check that $A \otimes W^T$ is in the Lie algebra to ensure unitarity.


Please let us know if we can clarify anything further.
Ok, thanks for the detailed reply. I hope it won't bother you if I continue a bit with the clarifications on the same topic (sorry if it does):

1) (minor point) So if $A=A^T$ and $W=-W^\dagger$ then these are square matrices, and then $X$ must be a square matrix, right? Or how can I compute $\exp(AXW)$ if $AXW$ is not a square matrix? (I'm confused by the notation $X\in \mathbb C^{a\times b}$, is it that necessarily $a=b$?) 

2) Assume now that $A=Id$, and $X,W$ are real-valued square matrices of equal dimension $n$, of which $W=-W^T$, to simplify things. Say $A=Id$, and $W=\left[\begin{array}{cc}0&-1\\\\ 1&0\end{array}\right]$ and $X=\left[\begin{array}{cc}0&1\\\\ 0&0\end{array}\right]$. Then $AXW=\left[\begin{array}{cc}1&0\\\\ 0&0\end{array}\right]\neq -(AXW)^T$ right? Then $\exp(AWX)$ is not unitary or orthogonal. I am confused, so maybe clarifying what I misunderstood with this example can help.
Thank you for your quick response and especially for recognizing the novelty of our work. We would be happy to hear about your thoughts on how we can prove stronger results on oversmoothing if you are willing to share it. 

___
Regarding your question, yes we are referring to $f_{Uconv}$. Let us go through this in more detail. 

As a reminder, we have that $f_{Uconv}(X)=\exp(g_{conv})(X)$ and $g_{conv}(X) = AXW$ and for an undirected graph $A = A^\dagger$ (for directed graphs, see Appendix D for proper handling of this). Furthermore, we constrain $W$ so that $W + W^\dagger=0$. Treating inputs and outputs as vectors, we can equivalently view $g_{conv}$ as acting on a vector rather than a matrix:

$$ g_{conv}(X) = AXW \leftrightarrow vec( g_{conv}(X) ) = (A \otimes W^T) vec(X) , $$
where we vectorize $X \in \mathbb{C}^{a \times b}$ with the command $vec(X) \in \mathbb{C}^{ab}$ that stacks columns on top of each other. Now, note that

$$ (A \otimes W^T)  + (A \otimes W^T) ^\dagger = A \otimes (W + W^\dagger)^T = 0, $$

and thus $(A \otimes W^T)$ is in the Lie Algebra of the Unitary/Orthogonal group. Finally, we now have
$$vec(\exp(g_{conv}(X))) = \exp(A \otimes W^T) vec(X),$$ 
and thus $\exp(A \otimes W^T)$ is in the Unitary/Orthogonal Lie group. 

Hopefully, viewing it this way in the "vector sense" clarifies things. Let us know if additional details are requested. We will clarify this in the updated manuscript.
I'll slowly go through the points of your rebuttal, so this may not be my only comment (sorry).

First of all, thanks for pointing out the difference to unitary layers, I was wrong on the novelty part, I apologize and I am grateful that you replied in detail. Reflecting upon this error of mine and your clarification, I now think that actually you could prove even stronger results about oversmoothing prevention by your approach, than the one on Rayleigh Quotients, but due to heavy review overload in 1 week time, I don't have a proposal right now, and that's anyway material for future work. 

About my correction/doubt number 3: when $W$ is an antisymmetric real matrix, $\exp(W)$ is orthogonal, we agree this far.

But in line 144 I'm concerned that you say "the above outputs an orthogonal map". In that comment, you refer to $exp(W)$ or to $f_{Uconv}(X)$? Given that $f_{Uconv}$ has a "U" in it, and given this paper's title, I assume it's $f_{Uconv}$. 

How does "exponential of real-valued antisymmetric implies orthogonal" apply to $f_{Uconv}(X)=\exp(g_{conv}(X))=\exp(AXW)$ ? are you using the fact that $AXW$ is antisymmetric real-valued too if $W$ is? I don't follow.. I'd like to see a proof of the map being orthogonal.

Thanks for the detailed rebuttal so far.
Thank you for your comment and quick feedback. 

Can we ask for a point of clarification? When you state `there is still no evidence on the benefit of using Lie Uniconv in the paper`, are you asking for evidence that this would outperform UniConv? As we state in the draft and in our additional experiments, we do not expect these to perform very differently. These are simply two different means of parameterizing unitary maps. They both nonetheless outperform other non-unitary message passing methods in our experiments.

&nbsp;

___
We are happy to provide empirical runtimes validating that unitary GCN is simply a constant factor overhead over vanilla GCN. See below.

**Time in seconds for a single forward pass on a single batch of data:**
```markdown
| Number of Nodes | Vanilla GCN | Spectral GCN | Gated GCN |   GPS  | Unitary GCN |
|-----------------|:-----------:|:------------:|:---------:|:------:|:-----------:|
| 200             |    0.0024   |    0.0054    |   0.0028  | 0.0047 |    0.0061   |
| 1000            |    0.0036   |    0.0106    |   0.0064  | 0.0124 |    0.0122   |
| 5000            |    0.0084   |    0.0383    |   0.0331  | 0.3009 |    0.0479   |
| 25000           |    0.0733   |    0.2640    |   0.1541  | 7.1890 |    0.1894   |
| 100000          |    0.2383   |    1.1645    |   1.1916  |   oom  |    1.3550   |
```

The table above shows runtimes for a forward pass of a model with two convolution layers of the specified form of width 128. As expected, unitary GCN has constant factor overhead over vanilla GCN scaling similarly to other enhanced message passing schemes like spectral GCN and Gated GCN. GPS scales poorly to large graphs and runs out of memory (oom) for the largest graph of 100,000 nodes. Unitary GCN uses UniConv layers though Lie UniConv performed virtually the same in terms of runtime. All experiments are on synthetic data using the layers as implemented in pytorch geometric.

Finally, we should note that we did not implement some of the efficiency improving techniques that one can implement in this line of work. For example, there exist approximations other than the Taylor approximation that may perform faster, implementations in the spectral basis can speed up unitary convolution in some instances, and there are other tricks as well (e.g. see projUNN for very fast unitary matrix updates using low rank updates). These are described in detail in the draft in appendix A and C.
I thank the reviewers for their response (which clarifies some of my doubts) and for the additional results listed in their attached PDF. Unfortunately, there is still no evidence on the benefit of using Lie Uniconv in the paper, which holds me back from increasing my score. Do the authors have any result comparing the speed of the two approaches perhaps?