Equivariant Machine Learning on Graphs with Nonlinear Spectral Filters

Authors: ya-wei eileen lin
ronen  talmon
ron  levie
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-14 19:36:55

Abstract: Equivariant machine learning is an approach for designing deep learning models that respect the symmetries of the problem, with the aim of reducing model complexity and improving generalization. 
In this paper, we focus on an extension of shift equivariance, which is the basis of convolution networks on images, to general graphs. Unlike images, graphs do not have a natural notion of domain translation. 
Therefore, we consider the graph functional shifts as the symmetry group: the unitary operators that commute with the graph shift operator. 
Notably, such symmetries operate in the signal space rather than directly in the spatial space.
We remark that each linear filter layer of a standard spectral graph neural network (GNN) commutes with graph functional shifts, but the activation function breaks this symmetry. Instead, we propose nonlinear spectral filters (NLSFs) that are fully equivariant to graph functional shifts and show that they have universal approximation properties. 
The proposed NLSFs are based on a new form of spectral domain that is transferable between graphs. 
We demonstrate the superior performance of NLSFs over existing spectral GNNs in node and graph classification benchmarks.

Main paper: 
 introduction: 
 In many fields, such as chemistry [37] , biology [36, 3] , social science [9] , and computer graphics [97] , data can be described by graphs. In recent years, there has been a tremendous interest in the development of machine learning models for graph-structured data [98, 13] . This young field, often called graph machine learning (graph-ML) or graph representation learning, has made significant contributions to the applied sciences, e.g., in protein folding, molecular design, and drug discovery [45, 2, 70, 89] , and has impacted the industry with applications in social media, recommendation systems, traffic prediction, computer graphics, and natural language processing, among others. Geometric Deep Learning (GDL) [13] is a design philosophy for machine learning models where the model is constructed to inherently respect symmetries present in the data, aiming to reduce model complexity and enhance generalization. By incorporating knowledge of these symmetries into the model, it avoids the need to expend parameters and data to learn them. This inherent respect for symmetries is automatically generalized to test data, thereby improving generalization [6, 76, 25] . For instance, convolutional neural networks (CNNs) respect the translation symmetries of 2D images, with weight sharing due to these symmetries contributing significantly to their success [54] . Respecting node re-indexing in a scalable manner revolutionized machine learning on graphs [13, 12] and has placed GNNs as a main general-purpose tool for processing graph-structured data. Moreover, within GNNs, respecting the 3D Euclidean symmetries of the laws of physics (rotations, reflections, and translations) led to state-of-the-art performance in molecule processing [24] .Our Contribution. In this paper, we focus on GDL for graph-ML. We consider extensions of shift symmetries from images to general graphs. Since graphs do not have a natural notion of domain translation, as opposed to images, we propose considering functional translations instead. We model the group of translations on graphs as the group of all unitary operators on signals that commute with the graph shift operator. Such unitary operators are called graph functional shifts. Note that each linear filter layer of a standard spectral GNN commutes with graph functional shifts, but the activation function breaks this symmetry. Instead, we propose non-linear spectral filters (NLSFs) that are fully equivariant to graph functional shifts and have universal approximation properties.In Sec. 3, we introduce our NLSFs based on new notions of analysis and synthesis that map signals between their node-space representations and spectral representations. Our transforms are related to standard graph Fourier and inverse Fourier transforms but differ from them in one important aspect. One key property of our analysis transform is that it is independent of a specific choice of Laplacian eigenvectors. Hence, our spectral representations are transferable between graphs. In comparison, standard graph Fourier transforms are based on an arbitrary choice of eigenvectors, and therefore, the standard frequency domain is not transferable. To achieve transferability, standard graph Fourier methods resort to linear filter operations based on functional calculus. Since we do not have this limitation, we can operate on the frequency coefficients with arbitrary nonlinear functions such as multilayer perceptrons. In Sec. 4 , we present theoretical results of our NLSFs, including the universal approximation and expressivity properties. In Sec. 5, we demonstrate the efficacy of our NLSFs in node and graph classification benchmarks, where our method outperforms existing spectral GNNs. 
background: 
 Notation. For N ∈ N, we denote [N ] = {1, . . . , N }. We denote matrices by boldface uppercase letter B, vectors (assumed to be columns) by lowercase boldface b, and the entries of matrices and vectors are denoted with the same letter in lower case, e.g., B = (b i,j ) i,j∈ [N ] . Let G = ([N ], E, A, X) be an undirected graph with a node set [N ], an edge set E ⊂ [N ] × [N ], an adjacency matrix A ∈ R N ×N representing the edge weights, and a node feature matrix (also called a signal) X ∈ R N ×d containing d-dimensional node attributes. Let D be the diagonal degree matrix of G, where the diagonal element d i,i is the degree of node i. Denote by ∆ any normal graph shift operator (GSO). For example, ∆ could be the combinatorial graph Laplacian or the normalized graph Laplacian given by L = D -A and N = D -1 2 LD -1 2 , respectively. Let ∆ = VΛV ⊤ be the eigendecomposition of ∆, where V is the eigenvector matrix and Λ = diag(λ 1 , . . . , λ N ) is the diagonal matrix with eigenvalues (λ i ) N i=1 ordered by |λ 1 | ≤ . . . ≤ |λ N |. An eigenspace is the span of all eigenvectors corresponding to the same eigenvalue. Let P i = P ∆;i denote the projection upon the i-th eigenspace of ∆ in increasing order of |λ|. We denote the Euclidean norm by ∥X∥ 2 . We define the channel-wise signal norm ∥X∥ sig of a feature matrix X ∈ R N ×d as the vector ∥X∥ sig = (∥X :,j ∥ 2 ) d j=1 ∈ R d , where X :,j is the j-th column on X and 0 ≤ a ≤ 1. We abbreviate multilayer perceptrons by MLP.Spectral GNNs define convolution operators on graphs via the spectral domain. Given a self-adjoint graph shift operator (GSO) ∆, e.g., a graph Laplacian, the Fourier modes of the graph are defined to be the eigenvectors {v i } N i=1 of ∆ and the eigenvalues {λ i } N i=1 are the frequencies. A spectral filter is defined to directly satisfy the "convolution theorem" [11] for graphs. Namely, given a signal X ∈ R N ×d and a function Q : R → R d ′ ×d , the operator Q(∆) : R N ×d → R N ×d ′ defined byEQUATIONis called a filter. Here, d ′ is the number of output channels. Spectral GNNs, e.g., [21, 49, 61, 5] , are graph convolutional networks where convolutions are via Eq. ( 1), with a trainable function Q at each layer, and a nonlinear activation function.Equivariance describes the ability of functions f to respect symmetries. It is expressed as f (H κ x) = H κ f (x), where K ∋ κ → H κ is an action of a symmetry group K on the domain of f . GNNs [83, 14] , including spectral GNNs [22, 49] and subgraph GNNs [31, 4] , are inherently permutation equivariant w.r.t. the ordering of nodes. This means that the network's operations are unaffected by the specific arrangement of nodes, a property stemming from passive symmetries [95] where transformations are applied to both the graph signals and the graph domain. This permutation equivariance is often compared to the translation equivariance in CNNs [56, 57] , which involves active symmetries [44] .The key difference between the two symmetries lies in the domain: while CNNs operate on a fixed domain with signals transforming within it, graphs lack a natural notion of domain translation. To address this, we consider graph functional shifts as the symmetry group, defined by unitary operators that commute with the graph shift operator. This perspective allows for our NLSFs to be interpreted as an extension of active symmetry within the graph context, bridging the gap between the passive and active symmetries inherent to GNNs and CNNs, respectively. 
nonlinear spectral graph filters: 
 In this section, we present new concepts of analysis and synthesis under which the spectral domain is transferrable between graphs. Following these concepts, we introduce new GNNs that are equivariant to functional symmetries -symmetries of the Hilbert space of signals rather than symmetries in the domain of definition of the signal [64] .For motivation, we start with the grid graph R with node set [M ] 2 and circular adjacency B, we define the translation operator T m,n by [m, n] as T m,n ∈ R M 2 ×M 2 ; (T m,n x) i,j = x l,k where l = (i-m) mod M and k = (j-n) mod M . Note that any T m,n is a unitary operator that commutes with the grid Laplacian ∆ R , i.e., T m,n ∆ R = ∆ R T m,n , and therefore it belongs to the group of all unitary operators U R that commute with the grid Laplacian ∆ R . In fact, the space of isotropic convolution operators (with 90 o rotation and reflection symmetric filters) can be seen as the space of all normal operators1 that commute with unitary operators from U R [18] . Applying a non-linearity after the convolution retains this equivariance, and hence, we can build multi-layer CNNs that commute with U R . By the universal approximation theorem [19, 34, 58] , this allows us to approximate any continuous function that commutes with U R .Note that such translation equivariance cannot be extended to general graphs. Achieving equivariance to graph functional shifts through linear spectral convolutional layers Q(∆) is straightforward, since these layers commute with the space of all unitary operators U ∆ that commute with ∆. However, introducing non-linearity ρ(Q(∆)X) breaks the symmetry. That is, there exists U ∈ U ∆ such that ρ (Q(∆)UX) ̸ = Uρ (Q(∆)X) , where ρ is any non-linear activation function, e.g., ReLU, Sigmoid, etc. This means that multi-layer spectral GNNs do not commute with U G , and are hence not appropriate as approximators of general continuous functions that commute with U G (see App. A for an example illustrating how non-linear activation functions break the functional symmetry). Instead, we propose in this paper a multi-layer GNN that is fully equivariant to U G , which we show to be universal: it can approximate any continuous graph-signal function (w.r.t. some metric) commuting with U G .We define the symmetry group of graph functional shifts as follows. Definition 1 (Graph Functional Shifts). The space of graph functional shifts is the unitary subgroup U ∆ , where a unitary matrix U is in U ∆ iff it commutes with the GSO ∆, namely, U∆ = ∆U.It is important to note that functional shifts, in general, are not induced from node permutations. Instead, functional shifts are related to the notion of functional maps [73] used in shape correspondence and are general unitary operators that are not permutation matrices in general. The value of the functionally translated signal at a given node can be a mixture of the content of the original signal at many different nodes. For example, the functional shift can be a combination of shifts of different frequencies at different speeds. See App. B for illustrations and examples of functional translations.A fundamental challenge with the symmetry group in Def. 1 is its lack of transferability between different graphs. Hence, we propose to relax this symmetry group. Let g 1 , . . . , g S : R → R be the indicator functions of the intervals {[l s , l s+1 ]} S s=1 , which constitute a partition of the frequency band [l 1 , l S ] ⊂ R. The operators g j (∆), interpreted via functional calculus Eq. ( 1), are projections of the signal space upon band-limited signals. Namely, g j (∆) = i:λi∈[lj ,lj+1] v i v ⊤ i . In our work, we consider filters g j that are supported on the dyadic sub-bands λ N r S-j+1 , λ N r S-j , where 0 < r < 1 is the decay rate. See Fig. 5 in App. F for an illustrated example. Note that for j = 1, the sub-band falls in [0,λ N r S-1 ]. The total band [0, l S ] is [0, λ N ].Definition 2 (Relaxed Functional Shifts). The space of relaxed functional shifts with respect to the filter bank {g j } K j=1 (of indicators) is the unitary subgroup U g ∆ , where a unitary matrix U is in U g ∆ iff it commutes with g j (∆) for all j, namely, Ug j (∆) = g j (∆)U.Similarly, we can relax functional shifts by restricting to the leading eigenspaces. Definition 3 (Leading Functional Shifts). The space of leading functional shifts is the unitary subgroup U J ∆ , where a unitary U is in U J ∆ iff it commutes with the eigenspace projections {P j } J j=1 .We use the terminology of analysis and synthesis, as in signal processing [68] , to describe transformations of signals between their graph and spectral representations. Here, we consider two settings: the eigenspace projections case and the filter bank (of indicators) case. The definition of the frequency domain depends on a given signal X ∈ R N ×d , where its projections to the eigenspaces of ∆ are taken as the Fourier modes. Spectral coefficients are modeled as matrices R that mix the Fourier modes, allowing to synthesize signals of general dimensions.Spectral Index Case. We first define analysis and synthesis using the spectral index up to frequency J. Let P J+1 = I -J j=1 P j be the orthogonal complement to the first J eigenprojections. Let X ∈ R N ×d be the graph signal and R ∈ R (J+1)d×(J+1) d the spectral coefficients to be synthesized, where d represents number of output channels. The analysis and synthesis are defined respectively byA ind (∆, X) = ∥P i X∥ sig J+1 i=1 ∈ R (J+1)d and S ind (R; ∆, X) = P j X ∥P j X∥ a sig + e J+1 j=1 R, (2)where 0 ≤ a ≤ 1 and 0 < e ≪ 1 are parameters that promote stability, and the power ∥P j X∥ a sig as well as the division in Eq. ( 2) are element-wise operations on each entry. Here, P j X/(∥P j X∥ a sig + e) J+1 j=1 ∈ R N ×(J+1)d denotes concatenation. We remark that the orthogonal complement in the (J + 1)-th filter alleviates the loss of information due to projecting to the low-frequency bands, and therefore, the full spectral range of the signal can be captured. This is particularly important for heterophilic graphs, which rely on high-frequency components for accurate label representation. The term index stems from the fact that eigenvalues are treated according to their index when defining the projections P j . Note that the synthesis here differs from classic signal processing as it depends on a given signal on the graph. When treating ∆ and X as fixed, this synthesis operation is denoted by S ind ∆,X (R) := S ind (R; ∆, X). We similarly denote A ind ∆ (X) := A ind (∆, X). Filter Bank Case. Similarly, we define the analysis and synthesis in the filter bank up to band g K as follows. Let g K+1 (∆) = I -K j=1 g j (∆) denote the orthogonal complement to the first K bands. Let X ∈ R N ×d be the graph signal, and let R ∈ R (K+1)d×(K+1) d represent the spectral coefficients to be synthesized, where d refers to the general dimension. The analysis and synthesis in the filter bank case are defined byA val (∆, X) = ∥g i (∆)X∥ sig K+1 i=1 ∈ R (K+1)d and S val (R; ∆, X) = g j (∆)X ∥g j (∆)X∥ a sig + e K+1 j=1(3) respectively, where a, e are as before. Here, g j (∆)X/(∥g j (∆)X∥a sig + e) K+1 j=1 ∈ R N ×(K+1)d .The term value refers to how eigenvalues are used based on their magnitude when defining the projections g j (∆). As before, we denote S val ∆,X and A val ∆ . In App. C.1, we present a special case of diagonal synthesis where d = d. In App. D.3, we show that the diagonal synthesis is stably invertible.We introduce three novel types of non-linear spectral filters (NLSF): Node-level NLSFs, Graph-level NLSFs, and Pooling-NLSFs. Fig. 1 illustrates our NLSFs for equivariant machine learning on graphs.Node-level NLSFs. To be able to transfer NLSFs between different graphs and signals, one key property of NLSF is that they do not depend on the specific basis chosen in each eigenspace. This independence is facilitated by the synthesis process, which relies on the input signal X. Following the spectral index and filter bank cases in Sec. 3.3, we define the Index NLSFs and Value NLSFs byEQUATIONEQUATIONwhere Ψ ind : R (J+1)d → R (J+1)d×(J+1) d and Ψ val : R (K+1)d → R (K+1)d×(K+1) d are called nonlinear frequency responses, and d is the output dimension. To adjust the feature output dimension, we apply an MLP with shared weights to all nodes after the NLSF. In the case when d = d and the filters operator diagonally (i.e., the product and division are element-wise in synthesis), we refer to it as diag-NLSF. See App. C for more details.Graph-level NLSFs. We first introduce the Graph-level NLSFs that are fully spectral, where the NLSFs map a sequence of frequency coefficients to an output vector. Specifically, the Index-based and Value-based Graph-level NLSFs are given byEQUATIONwhereΨ ind : R (J+1)d → R d ′ , Ψ val : R (K+1)d → R d ′, and d ′ is the output dimension.Pooling-NLSFs. We introduce another type of graph-level NLSFs by first representing each graph in a Node-level NLSFs as in Eq. ( 4) and Eq. ( 5). The final graph representation is obtained by applying a nonlinear activation function followed by a readout function to these node-level representations. We consider four commonly used pooling methods, including mean, sum, max, and L p -norm pooling, as the readout function for each graph. We apply an MLP after readout function to obtain a d ′dimensional graph-level representation. We term these graph-level NLSFs as Pooling-NLSFs. To understand which Laplacian and parameterization (index v/s value) of the NLSF are preferable in different settings, we follow the random geometric graph analysis outlined in [74] . Specifically, we consider a setting where random geometric graphs are sampled from a metric-probability space S. In such a case, the graph Laplacian approximates continuous Laplacians on the metric spaces under some conditions. We aim for our NLSFs to produce approximately the same outcome for any two graphs sampled from the same underlying metric space S, ensuring that the NLSF is transferable. In App. D.5, we show that if the nodes of the graph are sampled uniformly from S, then using the graph Laplacian L in Index NLSFs yields a transferable method. Conversely, if the nodes of the graph are sampled non-uniformly, and any two balls of the same radius in S have the same volume, then utilizing the normalized graph Laplacian N in Value NLSFs is a transferable method. Given that graphs may fall between these two boundary cases, we present an architecture that chooses between the Index NLSFs with respect to L and Value NLSFs with respect to N, as illustrated in Fig. 2 . While the above theoretical setting may not be appropriate as a model for every real-life graph dataset, it suggests that index NLSF may be more appropriate with L, value NLSFs with N, and different graphs are more appropriately analyzed by different balances between these two cases.In the Laplacian attention architecture, a soft attention mechanism is employed to dynamically choose between the two parameterizations, given by att (Θ ind (L, X), Θ val (N, X)) = αΘ ind (L, x)∥(1 -α)Θ val (N, X), where 0 ≤ α ≤ 1 is obtained using a softmax function to normalize the scores into attention weights, balancing each NLSFs' contribution. 
theoretical properties of nonlinear spectral filters: 
 We present the desired theoretical properties of our NLSFs at the node-level and graph-level.NLSFs are implemented by computing the eigenvectors of the GSO. Most existing spectral GNNs avoid direct eigendecomposition due to its perceived inefficiency. Instead, they use filters implemented by applying polynomials [49, 22] or rational functions [61, 5] to the GSO in the spatial domain. However, power iteration-based eigendecomposition algorithms, e.g., variants of the Lanczos method, can be highly efficient [80, 55] . For matrices with E non-zero entries, the computational complexity of one iteration for finding J eigenvectors corresponding to the smallest or largest eigenvalues (called leading eigenvectors) is O(JE). In practice, these eigendecomposition algorithms converge quickly due to their super-exponential convergence rate, often requiring only a few iterations, which makes them as efficient as message passing networks of signals with√ J channels.This makes NLSFs applicable to node-level tasks on large sparse graphs, as demonstrated empirically in App. F.5, since they rely solely on the leading eigenvectors. In Sec. 4.4, we show that using the leading eigenvectors can approximate GSOs well in the context of learning on graphs. Note that we can precompute the spectral projections of the signal before training. For node-level tasks, such as semi-supervised node classification, the leading eigenvectors only need to be pre-computed once, with a complexity of O(JE). This During the learning phase, each step of the architecture search and hyperparameter optimization takes O(N Jd) complexity for analysis and synthesis, and O(J2 d 2 )for the MLP in the spectral domain, which is faster than the complexity O(Ed 2 ) of message passing or standard spectral methods if N J < Ed. Empirical studies on runtime analysis are in App. F.For dense matrices, the computational complexity of a full eigendecomposition is O(N b ) per iteration, where N b is the complexity of matrix multiplication. This is practical for graph-level tasks on relatively small and dense graphs, which is typical for many graph classification datasets. In these cases, the eigendecomposition of all graphs in the dataset can be performed as a pre-computation step, significantly reducing the complexity during the learning phase.We demonstrate the node-level equivariance of our NLSFs, ensuring that our method respects the functional shift symmetries. The proof is given in App. D.1.Proposition 1. Index NLSFs in Eq. (4) are equivariant to the graph functional shifts U ∆ , and Value NLSFs in Eq. ( 5) are equivariant to the relaxed graph functional shifts U g ∆ .In this subsection, we discuss the approximation power of NLSFs.Node-Level Universal Approximation. We begin with a setting where a graph is given as a fixed domain, and the data distribution consists of multiple signals defined on this graph. An example of this setup is a spatiotemporal graph [17] , e.g., traffic networks, where a fixed sensor system defines a graph and the different signals represent the sensor readings collected at different times.In App. D.2.1, we prove the following lemma, which shows that linear NLSFs exhaust the space of linear operators that commute with graph functional shifts.Lemma 1. A linear operator R N ×d → R N ×d commutes with U J ∆ (resp. U g ∆ ) iffit is a NLSF based on a linear function Ψ in Eq. (4) (resp. Eq. ( 5)).Lemma 1 shows a close relationship between functions that commute with functional shifts and those defined in the spectral domain. This motivates the following construction of a pseudo-metric on R N ×d . In the case of relaxed functional shifts, we define the standard Euclidean metric dist E in the spectral domain R (K+1)×d . We pull back the Euclidean metric to the spatial domain to define a signal pseudo-metric. Namely, for two signals X and X ′ , their distance is defined by dist ∆ X, X ′ := dist E A(∆, X), A(∆ ′ , X ′ ) . This pseudo metric can be made into a metric by considering each equivalence class of signals with zero distance as a single point in the space. As MLPs Ψ can approximate any continuous function R (K+1)×d → R (K+1)×d (the universal approximation theorem [19, 34, 58] ), node-level NLSFs can approximate any continuous function that maps (equivalence classes of) signals to (equivalence classes of) signals. For details, see App. D.2.2. A similar analysis applies to hard functional shifts.Graph-Level Universal Approximation. The above analysis also motivates the construction of a graph-signal metric for graph-level tasks. For graphs with d-channel signals, we consider again the standard Euclidean metric dist E in the spectral domain R (K+1)×d . We define the distance between any two graphs with GSOs and signals (∆, X) and(∆ ′ , X ′ ) to be dist (∆, X), (∆ ′ , X ′ ) := dist E A(∆, X), A(∆ ′ , X ′) . This definition can be extended into a metric by considering the space G of equivalence classes of graph-signals with distance 0. As before, this distance inherits the universal approximation properties of standard MLPs. Namely, any continuous function G → R d ′ with respect to dist can be approximated by NLSFs based on MLPs. Additional details are in App. D.2.3.Graph-Level Expressivity of Pooling-NLSFs. In App. D.2.4, we show that Pooling-NLSFs are more expressive than graph-level NLSF when any L p norm is used in Eq. (4) and Eq. ( 5) with p ̸ = 2, both in the definition of the NLSF and as the pooling method. Specifically, for every graph-level NLSF, there is a Pooling-NLSF that coincides with it. Additionally, there are graph signals (∆, X) and (∆ ′ , X ′ ) for which a Pooling-NLSF can attain different values, whereas any graph-level NLSF must attain the same value. Hence, Pooling-NLSFs have improved discriminative power compared to graph-level NLSFs. Indeed, as shown in Tab. practice, which can be attributed to their increased expressivity. We refer to App. D.2.5 for additional discussion on graph-level expressivity.Since NLSFs on large graphs are based on the leading eigenvectors of ∆, we justify its low-rank approximation in the following. While approximating matrices with low-rank matrices might lead to a high error in the spectral and Frobenius norms, we show that such an approximation entails a uniformly small error in the cut norm. We define and interpret the cut norm in App. D.4.1, and explain why it is a natural graph similarity measure for graph machine learning.The following theorem is a corollary of the Constructive Weak Regularity Lemma presented in [29] .Its proof is presented in App. D. ∥M -C∥ □ < 3α 2 R J .Note that the bound in Thm. 1 is uniformly small, independently of M and its dimension N . This theorem justifies using the leading eigenvectors when working with the adjacency matrix as the GSO.For a justification when working with other GSOs see App. D.4. 
experimental results: 
 We evaluate the NLSFs on node and graph classification tasks. Additional implementation details are in App. E, and additional experiments, including runtime analysis and ablation studies, are in App. F.We first demonstrate the main advantage of the proposed Node-level NLSFs over existing GNNs with convolution design on semi-supervised node classification tasks. We test three citation networks [84, 100] : Cora, Citeseer, and Pubmed. In addition, we explore three heterophilic graphs: Chameleon, Squirrel, and Actor [79, 90] . For more comprehensive descriptions of these datasets, see App. E.We compare the Node-level NLSFs using Laplacian attention with existing spectral GNNs for nodelevel predictions, including GCN [49] , ChebNet [22] , ChebNetII [41] , CayleyNet [61] , APPNP [50] , GPRGNN [16] , ARMA [5] , JacobiConv [96] , BernNet [42] , Specformer [7] , and OptBasisGNN [38] .We also consider GAT [93] and SAGE [39] . For datasets splitting on citation graphs (Cora, Citeseer, and Pubmed), we apply the standard splits following [100] , using 20 nodes per class for training, 500 nodes for validation, and 1000 nodes for testing. For heterophilic graphs (Chameleon, Squirrel, and Actor), we use the sparse splitting as in [16, 41] , allocating 2.5% of samples for training, 2.5% for validation, and 95% for testing. We measure the classification quality by computing the average classification accuracy with a 95% confidence interval over 10 random splits. We utilize the source code released by the authors for the baseline algorithms and optimize their hyperparameters using Optuna [1] . Each model's hyperparameters are fine-tuned to achieve the highest possible accuracy. Detailed hyperparameter settings are provided in App. E.Tab. 1 presents the node classification accuracy of our NLSFs using Laplacian attention and the various competing baselines. We see that att-Node-level NLSFs outperform the competing models on the Cora, Citeseer, and Chameleon datasets. Notably, it shows remarkable performance on the densely connected Squirrel graph, outperforming the baselines by a large margin. This can be explained by the sparse version in Eq. ( 21) of Thm. 1, which shows that the denser the graphs is, the better its rank-J approximation. For the Pubmed and Actor datasets, att-Node-level NLSFs yield the second-best results, which are comparable to the best results obtained by APPNP. Recently, the work in [77] identified the presence of many duplicate nodes across the train, validation, and test splits in the dense split setting of the Chameleon and Squirrel [75] . This results in train-test data leakage, causing GNNs to inadvertently fit the test splits during training, thereby making performance results on Chameleon and Squirrel less reliable. To further validate the performance of our Node-level NLSFs on these datasets in the dense split setting, we use both the original and filtered versions of Chameleon and Squirrel, which do not contain duplicate nodes, as suggested in [77] . We use the same random splits as in [77] , dividing the datasets into 48% for training, 32% for validation, and 20% for testing.Tab. 2 presents the classification performance comparison between the original and filtered Chameleon and Squirrel. The baseline results are taken from [77] , and we include the following competitive models: ResNet+SGC [77] , ResNet+adj [77] , GCN [49] , GPRGNN [16] , FSGNN [69] , GloGNN [62] , and FAGCN [8] . The detailed comparisons are in App. F. We see in the table that the att-Node-level NLSFs consistently outperform the competing baselines on both the original and filtered datasets. att-Node-level NLSFs demonstrate less sensitivity to node duplicates and exhibit stronger generalization ability, further validating the reliability of the Chameleon and Squirrel datasets in the dense split setting. We note that compared to the dense split setting, the sparse split setting in Tab. 1 is more challenging, resulting in lower classification performance. A similar trend of significant performance difference between the two settings of Chameleon and Squirrel is also observed in [41] .We further illustrate the power of NLSFs on eight graph classification benchmarks [47] . Specifically, we consider five bioinformatics datasets [10, 53, 86] : MUTAG, PTC, NCI1, ENZYMES, and PROTEINS, where MUTAG, PTC, and NCI1 are characterized by discrete node features, while ENZYMES and PROTEINS have continuous node features. Additionally, we examine three social network datasets: IMDB-B, IMDB-M, and COLLAB. The unattributed graphs are augmented by adding node degree features following [26] . For more details of these datasets, see App. E.In graph classification tasks, a readout function is used to summarize node representations for each graph. The final graph-level representation is obtained by aggregating these node-level summaries and is then fed into an MLP with a (log)softmax layer to perform the graph classification task. We compare our NLSFs with two kernel-based approaches: GK [87] and WL [86] , as well as nine GNNs: GCN [49] , GAT [93] , SAGE [39] , ChebNet [22] , ChebNetII [41] , CayleyNet [61] , APPNP [50] , GPRGNN [16] , and ARMA [5] . Additionally, we consider the hierarchical graph pooling model DiffPool [101] . For dataset splitting, we apply the random split following [94, 101, 67, 104] , using 80% for training, 10% for validation, and 10% for testing. This random splitting process is repeated 10 times, and we report the average performance along with the standard deviation. For the baseline algorithms, we use the source code released by the authors and optimize their hyperparameters using Optuna [1] . We fine-tune the hyperparameters of each model to achieve the highest possible accuracy. Detailed hyperparameter settings for both the baselines and our method are provided in App. E.Tab. 3 presents the graph classification accuracy. Notably, the att-Graph-level NLSFs (i.e., NLSFs without the synthesis and pooling process) perform the second best on the ENZYMES and PROTEINS. Additionally, att-Pooling-NLSFs consistently outperform att-Graph-level NLSFs, indicating that the node features learned in our Node-level NLSFs representation are more expressive, corroborating our theoretical findings in Sec. 4.3. Furthermore, our att-Pooling-NLSFs consistently outperform all baselines on the MUTAG, PTC, ENZYMES, PROTEINS, IMDB-M, and COLLAB datasets. For NCI1 and IMDB-B, att-Pooling-NLSFs rank second and are comparable to ChebNetII. 
summary: 
 We presented an approach for defining non-linear filters in the spectral domain of graphs in a transferable and equivariant way. Transferability between different graphs is achieved by using the input signal as a basis for the synthesis operator, making the NLSF depend only on the eigenspaces of the GSO and not on an arbitrary choice of the eigenvectors. While different graph-signals may be of different sizes, the spectral domain is a fixed Euclidean space independent of the size and topology of the graph. Hence, our spectral approach represents graphs as vectors. We note that standard spectral methods do not have this property since the coefficients of the signal in the frequency domain depend on an arbitrary choice of eigenvectors, while our representation depends only on the eigenspaces. We analyzed the universal approximation and expressivity power of NLSFs through metrics that are pulled back from this Euclidean vector space. From a geometric point of view, our NLSFs are motivated by respecting graph functional shift symmetries, making them related to Euclidean CNNs.Limitation and Future Work. One limitation of NLSFs is that, when deployed on large graphs, they only depend on the leading eigenvectors of the GSO and their orthogonal complement. However, important information can also lie within any other band. In future work, we plan to explore NLSFs that are sensitive to eigenvalues that lie within a set of bands of interest, which can be adaptive to the graph.We offer a construction for complex-valued signals and first order derivative as the graph Laplacian, which is easy to extend to the real case and second order derivative. Consider a circular 1D grid with 100 nodes and the first order central difference as the Laplacian. Here, the Laplacian eigenvectors are the standard Fourier modes. In this case, one can see that a graph functional shift is any operator that is diagonal in the frequency domain and multiplies each frequency by any complex number with the unit norm. Consider the nonlinearity that takes the real part of the signal and then applies ReLU, which we denote in short by ReLU. Consider a graph signal x = (e iπ10n/100 + e iπ20n/100 ) 99 n=0 . We consider a graph functional shift S that shifts frequencies 10 and 20 at a distance of 5, and every other frequency is not shifted. Namely, frequency 10 is multiplied by e -iπ50/100 = -i, frequency 20 by e -iπ=-iπ100/100 = -1, and every other frequency is multiplied by 1. Consider also the classical shift D that translates the whole signal by 5 uniformly. Since x consists only of the frequencies 10 and 20, it is easy to see that Sx = Dx. Hence, ReLU(Sx) = ReLU(Dx) = D(ReLU(x)). Conversely, if we apply ReLU(x) and only then shift, note that ReLU(x) consists of many frequencies in addition to 10 and 20. For example, by nonnegativity of ReLU, ReLU(x) has a nonzero DC component (zeroth frequency). Now, S(ReLU(x)) only translates the 10 and 20 frequencies, so we have S(ReLU(x)) ̸ = D(ReLU(x)) = ReLU(S(x)).Here we present an additional discussion and illustrations on the new notions of symmetry.To illustrate the concepts of relaxed symmetry and functional translation, we present a toy example involving the classical translation and a functional translation of a Gaussian signal on a 2D grid with a standard deviation of 1.The classical translation involves moving the entire signal uniformly across the grid. This uniform movement can be represented as I ′ (x, y) = I(x -t x , y -t y ), where t x and t y are the translation amounts in the x and y directions, respectively, and x -t x and y -t y are circular translations, namely, subtractions modulo the size of the circular grid. For simplicity, we consider t x = t y = t. For instance, if t = 5, the Gaussian is shifted by 5 units in both the x and y directions. Fig. 3 top-row shows the classical translation for t = 0, t = 5, t = 10, and t = 15. We see that every part of the signal shifts at the same rate and direction, preserving the overall shape of the Gaussian signal. Note that classical translation is equivalent to modulation of the frequency domain, i.e. utx+vty) , where I(u, v) is the Fourier transform of I(x, y). Therefore, we can view the classical translation as a specific type of functional translation.I ′ (u, v) = I(u, v)e -i2π(Next, we illustrate a functional translation that is based on a frequency-dependent movement. Specifically, the Gaussian signal is decomposed into low and high-frequency components based on two indicator band-pass filters. In our example, the translation parameter t is different for the low and high-frequency components: low frequencies are shifted by t low while high frequencies are shifted by t high . This functional translation is defined via modulations in the frequency domain given by I ′ low (u, v) = I low (u, v)e -i2π(ut x,low +vt y,low ) for low-frequency components and I ′ high (u, v) = I high (u, v)e -i2π(ut x,high +vt y,high ) for high-frequency components. The combined functionally translated signal in the frequency domain is then I ′ = ( I ′ low , I ′ high ). Fig. 3 bottom-row demonstrates the functional translation for t = 0, t = 5, t = 10, and t = 15. We observe that the low-frequency components (smooth parts) of the signal move at one speed, while high-frequency components move at another. This demonstrates that functional symmetries are typically more rich than domain symmetries. We present another toy example to illustrate that functional translations are more stable than hard symmetries of the graph, namely, graph automorphisms (isomorphisms from the graph to itself). Automorphisms are another analog to translations on general graphs, which competes with our notion of functional shifts. Consider as an example the standard 2D circular grid (discretizing the torus). The automorphisms of the grid include all translations. However, this notion of hard symmetry is very sensitive to graph perturbations. If one adds or deletes even one edge in the graph, the automorphism group becomes much smaller and does not contain the translations anymore.In contrast, we claim that functional shifts are not so sensitive to graph perturbations. To demonstrate this empirically, we conducted the following toy experiment. We add a Gaussian noise to the edge weights of the 2D grid to create a perturbed graph. Given a standard domain shift, we optimize the coefficients of a functional shift so it is as close as possible to the classical domain shift of the clean grid in Frobenius error. Fig. 4 presents an example of classically and functionally shifted image of the digit 5. The original digit (left), a classical domain translation of the clean grid (middle), and a functional translation constructed from the perturbed graph (right) to match the classical translation. We see that standard translations can be approximated by a graph functional shift on a perturbed grid.We further demonstrate the robustness to graph perturbations of NLSFs trained on MNIST classification. The experimental setup follows previous work in [22, 71, 61] . We compare the NLSF on the clean grid (i.e., without Gaussian noise) to the NLSF on the perturbed grid. Tab. 4 presents the classification accuracy of the NLSF. We see that the classification performance is almost unaffected by the perturbation, indicating that NLSFs on the perturbed grid can roughly express CNN-like operations (translation equivariant functions).Original Translated Func. translated We present two special cases of NLSFs as follows.In Sec. 3, we introduced the NLSFs with the output dimension d, which is a tunable hyperparameter.Here, we present a special case when d = d such that the multiplication and division in synthesis are operated diagonally. Specifically, the diagonal analysis and synthesis in the spectral index case and in the filter bank case are defined respectively byA ind, diag (∆, X) = ∥P j X∥ sig J+1 j=1∈ R (J+1)×d and ( 7)EQUATIONandEQUATIONEQUATION)where the product and division are element-wise along the channel direction. That is,r j ⊙ Pj X ∥Pj X∥ a sig +e = r j 1 Pj X:,1 ∥Pj X:,1∥ a 2 +e , . . . , r j d Pj X :,d ∥Pj X :,d ∥ a 2 +ein the spectral index case (resp.r j ⊙ gj (∆)X ∥gj (∆)X∥ a sig +e = r j 1 gj (∆)X:,1 ∥gj (∆)X:,1∥ a 2 +e , . . . , r j d gj (∆)X :,d ∥gj (∆)X :,d ∥ a 2 +ein the filter bank case). Here, R = (r j ) J+1 j=1 ∈ R (J+1)×d in the spectral index case (resp. R = (r j ) K+1 j=1 ∈ R (K+1)×d in the filter bank case) are the spectral coefficients to be synthesized and a, e are as before.For Node-level diag-NLSFs, we define the Index diag-NLSFs and Value diag-NLSFs as followsEQUATIONEQUATIONwhere Ψ : R (J+1)d → R (J+1)d in Θ ind, diag and Ψ : R (K+1)d → R (K+1)d in Θ val, diag . To adjust the feature output dimension at each node, we apply an MLP with shared weights to all nodes after the NLSF. We present the empirical study of diag-NLSFs in App. F.7.In Sec. 3.3, we introduce the NLSFs with the orthogonal complement. Specifically, the (J + 1)-th filter in the Index NLSFs is given by P J+1 = I -J j=1 P j , and the (K + 1)-th filter in Vale NLSFs is defined as g K+1 (∆) = I -K j=1 g j (∆). To explore the effects of the orthogonal complement, we focus on the leading NLSFs in both the Index NLSFs and Vale NLSFs, considering only the first J and K filters without including their orthogonal complements. In this case, the analysis and synthesis in the spectral index case and in the filter bank case are respectively given byEQUATIONandS ind, lead (R; ∆, X) = P j X ∥P j X∥ a sig + e J j=1R andEQUATIONwhere P j and g j (∆) are defined as in Sec. 3.2. The lead-NLSFs are then defined byEQUATIONEQUATIONwhere Ψ ind : R Jd → R Jd×J d , Ψ val : R Kd → R Kd×K d , and d is the output dimension. To adjust the feature output dimension, we apply an MLP with shared weights to all nodes after the NLSF. Similar to App. C.1, when considering d = d such that the multiplication and division in synthesis are operated diagonally, we have the lead-diag-NLSFs defined byEQUATIONEQUATION)where Ψ : R Jd → R Jd in Θ ind, lead, diag and Ψ : R Kd → R Kd in Θ val, lead, diag . We present the empirical study of lead-NLSFs and lead-diag-NLSFs in App. F.7.We note that the numbering of the statements corresponds to the numbering used in the paper, and we have also included several separately numbered propositions and lemmas that are used in supporting the proofs presented. We restate the claim of each statement for convenience.We start with two simple lemmas that characterize the graph functional shifts. The lemmas directly follow the fact that two normal operators commute iff the projections upon their eigenspaces commute.Projection to Eigenspaces Case. Note thatR N =   J+1 j=1 P j R N   ,where denotes direct products of linear spaces and the (J +1)-th filter is the orthogonal complement, given by P J+1 = I -J j=1 P j . Denote by ⊕ the direct sum of operators.Lemma D.1. U is in U J ∆ iff it has the form U = J+1 j=1 U j, where U j is a unitary operator in P j R N .Proof. Let U ∈ U J ∆ . Since U commute with P j for j = 1, . . . , J, and with I, it also commutes with P J+1 = I -J j=1 P j . Therefore, we can writeU = J+1 j=1 P j U = J+1 j=1 P 2 j U = J+1 j=1 P j UP j .Now, when restricting P j UP j to an operator in P j R N , it is unitary. Indeed, for every v = P j v ∈ P j R N and u = P j u ∈ P j R N , since U is unitary and P j is self-adjoint and satisfies P 2 j = P j , we have ⟨P j UP j v, u⟩ = v, P j U -1 P j u , and P j U -1 P j is the inverse of P j UP j in P j R N , since for every v ∈ P j R NP j U -1 P j P j UP j v = P j U -1 UP j v = P j v = v,and similarly P j UP j P j U -1 P j v = v. Here, an invertible normal operator commutes with an orthogonal projection if and only if its inverse does.The other direction is trivial.Projection to Bands Case. Note thatR N =   K+1 j=1 g j (∆)R N   ,where the (K + 1)-th filter is the orthogonal complement, given by g K+1 (∆) = I -K j=1 g j (∆). Lemma D.2. U is in U g ∆ iff it has the form U = K+1 j=1 U j, where U j is a unitary operator in g j (∆)R N .The proof is analogous to the proof of Lemma D.1.Proposition 1. Index NLSFs in Eq. (4) are equivariant to the graph functional shifts U ∆ , and Value NLSFs in Eq. (5) are equivariant to the relaxed graph functional shifts U g ∆ .Proof. We start with Index-NLSFs. Consider the Index NLSF defined as in Eq. ( 4)Θ ind (∆, X) = P j X ∥P j X∥ a sig + e J+1 j=1 Ψ ind ∥P i X∥ sig J+1 i=1 .We need to show that for any unitary operator U ∈ U ∆ , Θ ind (∆, UX) = UΘ ind (∆, X). Consider U ∈ U ∆ and apply it to the graph signal X. The Index NLSF with the transformed input is given byΘ ind (∆, UX) = P j UX ∥P j UX∥ a sig + e J+1 j=1 Ψ ind ∥P i UX∥ sig J+1 i=1 .Since U ∈ U ∆ , it commutes with P j for all j, i.e., UP j = P j U. Using this commutation property, we can rewrite the norm and the projections as ∥P j UX∥ sig = ∥UP j X∥ sig .In addition, since U is a unitary matrix, it preserves the norm. Therefore, we have ∥UP j X∥ sig = ∥P j X∥ sig .Substituting these expressions back into the definition of the Index NLSFs givesΘ ind (∆, UX) = UP j X ∥P j X∥ a sig + e J+1 j=1 Ψ ind ∥P i X∥ sig J+1 i=1 .Notice that U appears linearly in the numerator of the fraction. Hence, we can factor it out byΘ ind (∆, UX) = U P j X ∥P j X∥ a sig + e J+1 j=1 Ψ ind ∥P i X∥ sig J+1 i=1 .The expression inside the summation is exactly the original Index NLSFs applied to X, so Θ ind (∆, UX) = UΘ ind (∆, X). Therefore, we have shown that applying the unitary operator U to the graph signal X results in the Index NLSFs being transformed by the same unitary operator U, proving the equivariance property.The proof for value-NLSF follows the same steps.In this section, we focus on value parameterized NLSFs. The analysis for index-NLSF is equivalent.Commute with Functional ShiftsWe next show that node-level linear NLSFs exhaust the space of linear operators that commute with graph functional shifts (on the fixed graph).Lemma 1. A linear operator R N ×d → R N ×d commute with U J ∆ (resp. U g ∆ ) iffit is a NLSF based on a linear function Ψ in Eq. (4) (resp. Eq. ( 5)).Proof. For simplicity, we restrict the analysis to the case of 1D signals (d = 1). The extension to a general dimension d is natural. By Lemma D.2, the unitary operators U in U g ∆ are exhausted by the operators of the formU =   K+1 j=1 U j   ,where U j is any unitary operator in g j (∆)R N . Hence, since the unitary representation T → T of the group of unitary operators in R m (for any m ∈ N) is irreducible, by Schur's lemma [85] any linear operator B that commutes with all operators of U g ∆ must be a scalar times the identity when projected to g j (∆)R N . Namely, B has the formB =   K+1 j=1 (b j I j )   ,where I j is the identity operator in g j (∆)R N . This means that linear node-level NLSFs exhaust the space of linear operators that commute with U g ∆ . The case of hard graph functional shifts is treated similarly.The above analysis motivates the following construction of a metric on R N . Given the filter bank {g j } K+1 j=1 , on graph with d-dimensional signals, we define the standard Euclidean metric dist E in the spectral domain R (K+1)×d . We pull back the Euclidean metric to the spatial domain to define a graph-signal metric. Namely, for two signals X and X ′ , their distance is defined to bedist ∆ X, X ′ := dist E A(∆, X), A(∆ ′ , X ′ ) .This defines a pseudometric on the space of graph-signals. To obtain a metric, we consider each equivalence class of signals with zero distance as a single point. Namely, we define the space N ∆ := R N /dist ∆ to be the space of signals with 1D features modulo dist. In N ∆ , the pseudometric dist becomes a metric, and A ∆ /dist ∆ an isometry of metric spaces2 . Now, since MLPs Ψ can approximate any continuous function R (K+1)×d → R (K+1)×d ′ (by the universal approximation theorem), and by the fact that A ∆ /dist ∆ is an isometry, node-level NLSFs based on MLPs in the spectral domain can approximate any continuous function from signals with d features to signal with d ′ features (N ∆ ) d → (N ∆ ) d ′ .The above analysis also motivates the construction of a graph-signal metric for graph-level tasks. For graphs with d-channel signals, we define the standard Euclidean metric dist E in the spectral domain R (K+1)×d . We pull back the Euclidean metric to the spatial domain to define a graph-signal metric.Namely, for two graphs with Laplacians and signals (∆, X) and (∆ ′ , X ′ ), their distance is defined to bedist (∆, X), (∆ ′ , X ′ ) := dist E A(∆, X), A(∆ ′ , X ′ )This defines a pseudometric on the space of graph-signals. To obtain a metric, we consider equivalence classes of graph-signals with zero distance as a single point. Namely, we define the space G to be the space of graph-signal modulo dist. In G, the function dist becomes a metric, and A an isometry.By the isometry property, G inherits any approximation property from R (K+1)×d . For example, since MLPs can approximate any continuous function R (K+1)×d → R d ′ , the space of NLSFs based on MLPs Ψ has a universality property: any continuous function G → R d ′ with respect to dist can be approximated by a NLSF based on an MLP.We now show that pooling-NLSF are more expressive than graph-level NLSF if the norm in Eq. ( 4) and Eq. ( 5) is L p with p ̸ = 2.First, we show that there are graph-signals that graph-level-NLSFs cannot separate and Pooling-NLSFs can. Consider an index NLSF with norm L 1 normalize by 1/N . Namely, for a ∈ R N ,∥a∥ 1 = 1 N N n=1 |a n | .The general case is similar.For the two graphs, take the graph LaplacianEQUATIONwith eigenvalues 0, 1, and corresponding eigenvectors (1, 1) and (1, -1). Take the graph LaplacianEQUATIONwith eigenvalues 0, 1, 3. The first two eigenvectors are (1, 1, 1) and (1, 0, -1) respectively.Consider an Index-NLSF based on two eigenprojections P 1 , P 2 . As the signal of the first graph take (1, 1) + (1, -1), and for the second graph take (1, 1, 1) + 3 2 (1, 0, -1). Both graph-signals have the same spectral coefficients (1, 1), so graph-level NLSF cannot separate them. Suppose that the NLSF isΘ val (∆, X) = (∥P 1 X∥ a sig + e) P 1 X ∥P 1 X∥ a sig + e + (∥P 2 X∥ a sig + e) P 2 X ∥P 2 X∥ a sig + e .The outcome of the corresponding Pooling NLSF on the two graphs is1 = ∥(1 + 1, 1 -1)∥ 1 ̸ = ∥(1 + 3/2, 1, 1 -3/2)∥ 1 = 4 3 .Hence, Pooling-NLSFs separate these inputs, while graph-level-NLSFs do not.Next, we show that Pooling-NLSFs are at least as expressive as graph-level NLSFs. In particular, any graph-level NLSF can be expressed as a pooling NLSF.Let Ψ be a graph-level NLSF. Define the node-level NLSF that chooses one spectral index j with a nonzero value (e.g., the band with the largest coefficient) and projects upon it the value Ψ(∆, X)(∥P j X∥ a sig + e)/ ∥P j X∥ sig . Hence, before pooling, the NLSF givesΘ(∆, X) = Ψ(∆, X)(∥P j X∥ a sig + e) P j X ∥P j X∥ a sig + e / ∥P j X∥ sig ,where j depends on the spectral coefficients (e.g., it is the index of the largest spectral coefficient). Hence, after pooling, the Pooling NLSF returns Ψ(∆, X), which coincides with the output of the graph-level NLSF. Now, let us show that for p = 2, they have the same expressivity. For any NLSF,Θ(∆, X) = J j=1 Ψ ∥P i X∥ sig J i=1 j P j X ∥P j X∥ a sig + e 2 sig = J j=1 Ψ ∥P i X∥ sig J i=1 2 j ∥P j X∥ 2 sig ∥P j X∥ a sig + e .This is a generic fully spectral NLSF.We offer additional discussion on the expressivity of graph-level NLSFs. Note that graph-level NLSFs are bounded by the expressive power of MPNNs with random positional encodings. For example, in [81] , the authors showed that random features improve GNN expressivity, distinguishing certain structures that deterministic GNNs cannot. The Lanczos algorithm for computing the leading J eigenvectors can be viewed as a message-passing algorithm with a randomly initialized J-channel signal.The example in App. D.2.4 can be written as a standard linear graph spectral filter, followed by the pooling (readout) function. This shows that graph-level NLSFs (without synthesis and pooling) are not more expressive than standard spectral GNNs with pooling. In the other direction, there are no graph-signals that can be separated by a graph-level NLSF but not by a spectral GNN. Given two graph-signals that an NLSF can separate, we can build a specific standard spectral GNN that gives the same output as the NLSF specifically for the two given graph-signals. However, several considerations should be noted:1. the feature dimension of this GNN would is as large as the combined number of eigenvalues of the two graphs times the number of features, 2. this GNN is designed to fit these two specific graphs, and will not work for other graphs, and, 3. if implemented via polynomial filters, the polynomial order would have to be the combined number of eigenvalues of the two graphs, which makes it very unstable.Let us explain the architecture next. Given the two graphs (G 1 , f 1 ) with N 1 vertices and (G 2 , f 2 ) with N 2 vertices, with node features of dimension D, define band-pass filters that separate the combined set of eigenvalues of the two given graphs. Concatenate these filters to build a single multi-channel filter. Namely, this filter maps the signal to the sequence of band-pass filtered signal -each band at a different channel, for a total of D(N 1 + N 2 ) channels. If the band-pass filters are based on polynomials, each polynomial would be chosen to be zero on all eigenvalues except for one. Apply L 2 -norm pooling on each channel to obtain a single feature of dimension D(N 1 + N 2 ). This gives the sequence of norms of the signal at the bands, where, for the two given graphs, the two pooled signals of the two graphs are supported on disjoint sets of indices. Hence, the linear spectral filter contains the frequency information of the NLSF. Then, one can apply the MLP of the NLSF to the channels corresponding to the first graph, and to the channels corresponding to the second graph. This means that the linear spectral GNN gives the exact same outputs on (G 1 , f 1 ) and (G 2 , f 2 ) as the NLSF.Regarding pooling-NLSF, we do not offer an answer to the question whether standard spectral GNNs are more powerful (assuming an unlimited budget) than pooling-NLSFs, or vice-versa. More practically meaningful questions would be:1. Compare the expressivity of standard spectral GNNs to NLSFs for a given budget of parameters.2. How many graph-signal can a single spectral GNN separate, vs a single NLSF, of the same budget.We leave these questions as open problems for future research. We note that, in practice, NLSFs with the same budget as standard spectral GNNs perform better.We present the analysis for diagonal synthesis defined from App. C.1. In the fixed graph setting, given any 1D signal X such that g j (∆)X ̸ = 0 for every j ∈ [K + 1], we show that the synthesis operator S val, diag ∆,X is stably invertible.Since we consider 1D signal X, the diagonal synthesis in filter bank case in Eq. ( 10) can be written asS val, diag ∆,X = Hr ′ ,whereEQUATION∆)X ∥gj (∆)X∥ a 2 +e ∈ R N , and r ′ = [r 1 , r 2 , . . . , r K+1 ] ∈ R K+1 . Since (g j1 (∆)X) ⊤ (g j2 (∆)X) = 0 for all j 1 , j 2 ∈ [K + 1], the matrix H ⊤ H is a diagonal matrix with entriesg1(∆)X ∥g1(∆)X∥ a 2 +e 2 2. Therefore, the singular values of H areσ j = ∥g 1 (∆)X∥ 2 ∥g 1 (∆)X∥ a2 + e the right singular vectors are the standard basis elements e j in R K+1 , and the left singular vectors are g j (∆)X ∥g j (∆)X∥ 2 .for j ∈ [K + 1]. Hence, we haveS val, diag ∆,X 2 = max j σ j and S val, diag ∆,X -1 2 = 1 min j σ j .Suppose that a = 1 and e = 0. In this case, S val, diag ∆,X is an isometry from the spectral domain to a subspace of the signal space R N . Analysis is the adjoint of synthesis. This analysis can be extended to the index parametrization case for diagonal synthesis, and to higher dimensional signals.In this section, we develop the setting under which the low-rank approximation of GSOs with their leading eigenvectors can be interpreted as a uniform approximation (Sec. 4.4).The cut norm of M ∈ R N ×N is defined to beEQUATIONThe distance between two matrices in cut norm is defined to be ∥M -M ′ ∥ □ .The cut norm has the following interpretation, which has precise formulation in terms of the weak regularity lemma [32, 65] . Any pair of (deterministic) graphs are close to each other in cut norm if and only if they can be described as pseudo-random graphs sampled from the same stochastic block model. Hence, the cut norm is a meaningful notion of graph similarity for practical graph machine learning, where graphs are noisy and can represent the same underlying phenomenon even if they have different sizes and topologies. In addition, the distance between non-isomorphic graphons is always positive in cut norm [66] . In this context, the work in [59] showed that GNNs with normalized sum aggregation cannot separate graphs that have zero distance in the cut norm. This means that the cut norm is sufficiently discriminative for practical machine learning on graphs.The following lemma, called the constructive weak regularity lemma in Hilbert spaces, was proven in [29] . It is an extension of the classical respective result from [65] .We define the Frobenius norm normalized by 1/N as∥M∥ F = 1 N N i,j=1 |m i,j | 2 .Lemma D.3. [ [29] ] Let {K j } j∈N be a sequence of nonempty subsets of a real Hilbert space H and let δ ≥ 0. Let J > 0, let R ≥ 1 such that J/R ∈ N, and let g ∈ H. Let m be randomly uniformly sampled from [J]. Then, in probability 1 -1 R (with respect to the choice of m), any vector of the formg * = m j=1 γ j f j such that γ = (γ j ) m j=1 ∈ R m and f = (f j ) m j=1 ∈ K 1 × . . . × K mthat gives a close-to-best Hilbert space approximation of g in the sense thatEQUATION)where the infimum is over γ ∈ R m and f ∈ K 1 × . . . × K m , also satisfies ∀w ∈ K m+1 , |⟨w, g -g * ⟩| ≤ ∥w∥ ∥g∥ R J + δ.∥M -C∥ □ < 3α 2 R J .Proof. Let us use Lemma D.3, with H = R N ×N , and K j = K the set of symmetric rank one matrices of the form vv ⊤ where v ∈ R N is a column vector. Denote by Y m the space of linear combinations of m elements of K, which is the space of symmetric matrices of rank bounded by m.For the Hilbert space norm, we take the Frobenius norm. In the setting of the lemma, we take g = M, and g * ∈ Y m , and δ = 0. By the lemma, with probability 1 -1/R, any Frobenius minimizer C, namely, that satisfies ∥M -C∥ F = min C ′ ∈Ym ∥M -C ′ ∥ F , also satisfies ⟨Y, M -C⟩ ≤ ∥Y∥ F ∥M∥ F R J ≤ α ∥Y∥ F R J for every Y ∈ Y m . Hence, for every choice of subset S, T ⊂ [N ], we have i∈S j∈T (m i,j -c i,j ) = 1 2 i∈S j∈T (m i,j -c i,j ) + i∈T j∈S (m i,j -c i,j ) = 1 2 i∈S∪T j∈S∪T (m i,j -c i,j ) - i∈S j∈S (m i,j -c i,j ) - i∈T j∈T (m i,j -c i,j ) ≤ 1 2 1 S∪T 1 ⊤ S∪T , M -C F + 1 2 1 S 1 ⊤ S , M -C F + 1 2 1 T 1 ⊤ T , M -C F ≤ 3α 2 R J ,where for a set S ⊂ [N ], denote by 1 S ∈ R N the column vector with 1 at coordinates in S and 0 otherwise.Hence, we also have∥M -C∥ □ ≤ 3α 2 R J .Lastly, note that by the best rank-m approximation theorem (Eckart-Young-Mirsky Theorem [92, Thm. 5.9]), any Frobenius minimizer C is the projection upon the m leading eigenvectors of M (or some choice of these eigenvectors in case of multiplicity higher than 1). Now, one can use the adjacency matrix A as M in Thm. 1. When working with sparse matrices of E ≪ N 2 edges, to achieve a meaningful scaling of cut distance, we re-normalize the cut norm and define∥M∥ (E) □ = N 2 E ∥M∥ □ .With this norm, Thm. 1 givesEQUATION)While this bound is not uniformly small in N , it is still independent of M. In contrast, the error bounds for spectral and Frobenius norms do depend on the specific properties of M.Now, if we want to apply Thm. 1 to other GSOs ∆, we need to make some assumptions. Note that when the GSO is a Laplacian, we take as the leading eigenvectors the ones correspoding to the smallest eigenvalues, not the largest ones. To make the theorem applicable, we need to reorder the eigenvectors of ∆. This can be achieved by applying a decreasing function h to ∆, such as h(∆) = ∆ -1 . The role of h is to amplify the eigenspaces of ∆ in which most of the energy of signals interest (the ones that often appear in the dataset) is concentrated. Under the assumption that h(∆) has entries bounded by some not-too-high α > 0, one can now justify approximating GSOs by low-rank approximations based on the smallest eigenvectors.projection upon the Paley-Wiener space of L. R N : L 2 (S) → L 2 [N ] be the sampling operator, defined by (R N f ) j = f (x j ).Let f : S → R be a bounded signal over the metric space and suppose that f ∈ PW(J). Denote f = (f (x n )) N n=1 = R N f . By Monte Carlo theory, we have(Lf ) j ≈ w(x j )N Br(xj ) (f (x) -f (y))dy = w(x j )N Lf (x j ),So, in case the sampling is uniform. i.e., w(x) = 1, we have L ≈ N L, pointwise. To make this precise, let us recall Hoeffding's Inequality [43] . Theorem D.1 (Hoeffding's Inequality [43] ). Let Y 1 , . . . , Y N be independent random variables such that a ≤ Y n ≤ b almost surely. Then, for every k > 0,P 1 N N n=1 (Y n -E[Y n ]) ≥ k ≤ 2 exp - 2k 2 N (b -a) 2 .Using Hoeffding's Inequality, one can now show that there is an event of probability more than 1 -p in which for every j ∈ [N ], the error between ( 1 N LR N f ) j and Lf satisfiesR N Lf - 1 N LR N f = O J log(1/p) + log(N ) + log(2) N .The fact that different graphs of different sizes N approximate L with different scaling factors means that L is not value transferable. Let us show that L is index transferable.We now evoke the transferability theorem -a slight reformulation of Section 3 of Theorem 4 in [60] . Theorem D.2. Let L be a compact operator on L 2 (S) and L an operator on C N . Let J ∈ N and let R N be the sampling operator defined above. Let PW(J) be the space spanned by the J leading eigenvectors of L, and let f ∈ PW(J). Let g : R → R be Lipschitz with Lipschitz constant D. Then,R N g(L)f - 1 N g(L)R N f ≤ D √ J R N Lf - 1 N LR N f .For every leading eigenvalue λ j of L let γ ij be the leading eigenvalue of L closest to λ, where, if there are repeated eigenvalues, i j is chosen arbitrarily from the eigenvalues of L that best approximate λ i . Let δ = min{α, β}, where α = minj∈[J] min min ij ̸ =i∈[N ] |λ j -γ i | , min j̸ =m∈[J] |λ j -λ m | , min j̸ =m∈[J] γ ij -λ m and β = min i∈[J],n∈[N ] |γ i -γ n |and suppose that δ > 0. For each j ∈ [N ] there exists a Lipschitz continuous function g j with Lipschitz constant δ -1 such that g j (λ j ) = g j (γ ij ) = 1 and γ is zero on all other leading eigenvalues of L and all other eigenvalues of L. Hence,g j (L)f = p L j f, g j (L)R N f = p L ij f, where p L j is the projection upon the space spanned by the j-th eigenvector of L, and p L ij is the projection upon the space spanned by the i j -th eigenvector of L . Now, by the transferability theorem,R N p L j f -p L ij f ≤ D √ J R N Lf - 1 N Lf = O( J/N ).By induction over j, with base j = 1, we must have i j = j for every j ∈ [J]. Now, note that by standard Monte Carlo theory (evoking Hoeffding's inequality again and intersecting events), we have in high probability. Hence, by the fact thatR N p L j f 2 2 -p L j f 2 2 = O(J -1/2 )R N p L j f 2 2 -p L j f 2 2 = R N p L j f 2 -p L j f 2 R N p L j f 2 + p L j f 2 and R N p L j f 2 + p L j f 2 is bounded from below by the constant p L j f 2 , we have R N p L j f 2 -p L j f 2 = O(J -1/2 ).This shows thatp L j f 2 ≈ p L j f 2 ,which shows index transferability. Namely, for two graphs of N and N ′ nodes sampled from S, with corresponding Laplacians L and L ′ , by the triangle inequality,p L j f 2 ≈ p L ′ j f 2 .Next, we show value transferability for N. Here,Nf (x j ) ≈ 1 V (x j ) Lf (x j ).Hence, if V (x) = V is constant, we have N ≈ L up to a constant that does not depend on N . In this case, by a similar analysis to the above, N is value transferable. We note that N is also index transferable, but value transferability is guaranteed in a more general case, where we need not assume a separable spectrum.We describe the experimental setups and additional details of our experiments in Sec. 5 . The experiments are performed on NVIDIA DGX A100.We provide a detailed overview of the experimental settings for semi-supervised node classification tasks, along with the validated hyperparameters used in our benchmarks.Datasets. We consider six datasets in node classification tasks, including Cora, Citeseer, Pubmed, Chameleon, Squirrel, and Actor. The detailed statistics of the node classification benchmarks are summarized in Tab. 5. For datasets splitting on citation graphs (Cora, Citeseer, and Pubmed), we follow the standard splits from [100] , using 20 nodes per class for training, 500 nodes for validation, and 1000 nodes for testing. For heterophilic graphs (Chameleon, Squirrel, and Actor), we adopt the sparse splitting method from [16, 41] , allocating 2.5% of samples for training, 2.5% for validation, and 95% for testing. The classification quality is assessed by computing the average classification accuracy across 10 random splits, along with a 95% confidence interval.Baselines. For GCN, GAT, SAGE, ChebNet, ARMA, and APPNP, we use the implementations from the PyTorch Geometric library [28] . For other baselines, we use the implementations released by the respective authors. ], and the number of the bands in Value NLSFs within K ≤ S. The graph-level representation is then fed into an MLP with a (log)softmax classifier, using a cross-entropy loss function for predictions over the labels. Specifically, the MLP consists of three fully connected layers with 256, 128, and 64 neurons, respectively, followed by a (log)softmax classifier. We conduct a grid search on the hyperparameters for each dataset using Optuna [1] . An early stopping criterion is employed during training, stopping the process if the validation loss does not decrease for 100 consecutive epochs.Here, we present additional experiments on node and graph classification benchmarks, ablation studies, runtime analysis, and uniform sub-bands.We present additional experimental results for semi-supervised node classification using random splits, adhering to the protocol established by [41] . The results, summarized in Tab. 7, demonstrate the classification accuracy across six benchmark datasets: Cora, Citeseer, Pubmed, Chameleon, Squirrel, and Actor. Our att-Node-level NLSFs achieve the highest accuracy on four out of the six datasets, outperforming other models significantly. On Pubmed, it records a close second-highest accuracy, slightly behind APPNP. Our method achieves a competitive second place in the Actor dataset. att-Node-level NLSFs demonstrate substantial improvements, particularly in challenging datasets like Chameleon and Squirrel. The comparison models, including GCN, GAT, SAGE, ChebNet, ChebNetII, CayleyNet, APPNP, GPRGNN, and ARMA, varied in their effectiveness, with ChebNetII and APPNP also showing strong results on several datasets. These findings highlight the efficacy of the att-Node-level NLSFs in semi-supervised node classification tasks.Following [77] , we conduct additional experiments on the original and filtered Chameleon and Squirrel datasets in the dense split setting. We use the same random splits as in [77] , dividing the datasets into 48% for training, 32% for validation, and 20% for testing. We compare the Nodelevel NLSFs using Laplacian attention with GCN [49] , SAGE [39] , GAT [93] , GT [88] , H 2 GCN [103] , CPGNN [102] , GPRGNN [16] , FSGNN [69] , GloGNN [62] , FAGCN [8] , GBKGNN [23] , JacobiConv [96] , and ResNet [40] with GNN models [77] . The study by [103] demonstrates the benefits of separating ego-and neighbor-embeddings in the GNN aggregation step when dealing with heterophily. Therefore, [77] also adopts this approach for the GNN aggregation step in GAT and GT models, denoted as "sep." The baseline results used for comparison are taken from [77] . Tab. 8 presents the full performance comparison on the original and filtered datasets. Note that Tab. 8 is the same as Tab. 2 but with more baseline methods. att-Node-level NLSFs achieve the highest accuracy on both the filtered Chameleon and Squirrel datasets. Additionally, att-Node-level NLSFs demonstrate strong performance on the original Chameleon dataset, achieving the highest In Sec. 3.5, we note that using the graph Laplacian L in Index NLSFs and the normalized graph Laplacian N in Value NLSFs is transferable. Since real-world graphs often fall between these two boundary cases, we present the Laplacian attention NLSFs that operate between them at both the node-level and graph-level. Indeed, as demonstrated in Sec. 5, the proposed att-Node-level NLSFs, att-Graph-level NLSFs, and att-Pooling-NLSFs outperform existing spectral GNNs.We conduct an ablation study to evaluate the contribution and effectiveness of different components within the att-Node-level NLSFs, att-Graph-level-NLSFs, and att-Pooling-NLSFs on node and graph classification tasks. Specifically, we compare the Index NLSFs and Value NLSFs using both the graph Laplacian L and the normalized graph Laplacian N to understand their individual and Laplacian attention impact on these tasks.The ablation study of att-Node-level NLSPs for node classification is summarized in Tab. 9. We investigate the performance on six node classification benchmarks as in Sec. 5, including Cora, Citeseer, Pubmed, Chameleon, Squirrel, and Actor. Tab. 9 shows that using the graph Laplacian L in Index NLSFs (denoted as Θ ind (L, •)) and the normalized graph Laplacian N in Value NLSFs (denoted as Θ val (N, •) ) has superior node classification accuracy compared to using the normalized graph Laplacian N in Index NLSFs (denoted as Θ ind (N, •) ) and the graph Laplacian L in Value NLSFs (denoted as Θ val (L, •)). This is in line with our theoretical findings in App. D.5. Moreover, the att-Node-level NLSFs using the Laplacian attention att(Θ ind (L, •), Θ val (N, •)) yield the highest accuracies across all datasets, corroborating the findings in Sec. 3.5. We also note that without Laplacian attention, the Node-level NLSFs (Θ ind (L, •) and Θ val (N, •)) alone still achieve more effective classification performance compared to existing baselines, as shown in Tab. 1.Tab. 10 demonstrates the ablation study of Graph-level NLSFs for graph classification tasks. We examine the eight graph datasets as in Sec. 5, including MUTAG, PTC, ENZYMES, PROTEINS, NCI1, IMDB-B, IMDB-M, and COLLAB. Similar to the above, we investigate the Index and Value settings using graph Laplacian L and normalized graph Laplacian N, including Φ ind (L, •), Φ ind (N, •), Φ val (L, •), and Φ val (N, •), along with their variants using Laplacian attention. Here, we see that att-Graph-level NLSFs do not show significant improvement over the standard Graph-level NLSFs. Notably, Φ val (N, •) outperforms other models in social network datasets (IMDB-B, IMDB-M, and COLLAB), where the node features are augmented by the node degree. We plan to investigate the limited graph attribution in future work.The ablation study of att-Pooling-NLSFs for graph classification is reported in Tab. 11. Similar to Tab. 10, we consider eight graph classification benchmarks: MUTAG, PTC, ENZYMES, PROTEINS, NCI1, IMDB-B, IMDB-M, and COLLAB. Tabl. 11 demonstrates that using the graph Laplacian L in Index Pooling-NLSFs (denoted as Θ P ind (L, •)) and the normalized graph Laplacian N in Value Pooling-NLSFs (denoted as Θ P val (N, •)) achieves superior graph classification accuracy compared to using the normalized graph Laplacian N in Index Pooling-NLSFs (denoted as Θ P ind (N, •)) and the graph Laplacian L in Value Pooling-NLSFs (denoted as Θ P val (L, •)). This finding aligns with our theoretical results in App. D.5. Moreover, the att-Pooling-NLSFs using the Laplacian attention att P (Θ ind (L, •), Θ val (N, •)) yield the highest accuracies across all datasets, corroborating the findings in Sec. 3.5. In addition, att-Pooling-NLSFs consistently outperform att-Graph-level NLSFs as shown in Tab. 10, indicating that the node features learned in our Node-level NLSFs representation are more expressive. This supports our theoretical findings in Sec. 4.3.The ablation study demonstrates that the Laplacian attention between the Index and Value NLSFs significantly enhances classification accuracy for both node and graph classification tasks across various datasets, outperforming existing baselines.Table 13 : Experimental results on large heterophilic graphs. The results for BernNet, ChebNet, ChebNetII, and GPRGNN are taken from [41] , while the results for OptBasisGNN are taken from [38] . All other competing results are taken from [63] . In our NLSFs, the eigendecomposition can be calculated once for each graph and reused in the training process. This step is essential as the cost of the forward pass during model training often surpasses the initial eigendecomposition preprocessing cost. Note that the computation time for eigendecomposition is considerably less than the time needed for model training. For medium and large graphs, the computation time is further reduced when partial eigendecomposition is utilized, making it more efficient than the training times of competing baselines. To evaluate the computational complexity of our NLSFs compared to baseline models, we report the empirical training time in Tab. 12. Our Node-level NLSFs showcase competitive running times, with moderate values per epoch and total running times that are comparable to the most efficient models like GCN, GPRGNN, and APPNP. Notably, Node-level NLSFs are particularly efficient for the Cora, Citeseer, and Pubmed datasets. For the dense Squirrel graph, our Node-level NLSFs exhibit efficient performance with a moderate running time, outperforming several models that struggle with significantly higher times.To demonstrate the scalability of our method, we conduct additional tests on five large heterophilic graphs: Penn94, Pokec, Genius, Twitch-Gamers, and Wiki datasets from [63] . The experimental setup is in line with previous work by [16, 63, 41] . We use the same hyperparameters for our NLSFs as reported in App. E. Tab. 13 presents the classification accuracy. We see that NLSFs outperform the competing methods on the Penn94, Pokec, Genius, and Wiki datasets. For the Twitch-Gamers dataset, NLSFs yield the second-best results. Our additional experiments show that our method could indeed scale to handle large-scale graphs effectively.F.6 Index-by-Index Index-NLSFs and Band-by-Band Value-NLSFs Our primary objective in this work is to introduce new GNNs that are equivariant to functional symmetries, based on a novel spectral domain that is transferable between graphs. We emphasize the unique aspects of our method rather than providing an exhaustive list of operators in this group, which, while important, is a key direction for future research.Here, we present a new type of NLSFs: index-by-index Index-NLSFs and band-by-band Value-NLSFs. We denote them as follows:Γ ind (∆, X) = J j=1 γ j ∥P j X∥ sig P j X ∥P j X∥ a sig + e and Γ val (∆, X) = K j=1 γ j ∥g j (∆)X∥ sig g j (L)X ∥g p (∆)X∥ a sig + e ,where a, e are as before in Sec. 3.3. Note that these operators are also equivariant to our group actions.We investigate index-by-index Index-NLSFs and band-by-band Value-NLSFs in graph classification tasks as described in Sec. 5, including MUTAG, PTC, ENZYMES, PROTEINS, NCI1, IMDB-B, IMDB-M, and COLLAB datasets. Unlike Graph-NLSFs, which are fully spectral and map a sequence of frequency coefficients to an output vector, index-by-index Index-NLSFs and band-by-band Value-NLSFs do not possess such a sequential spectral form. In index-by-index Index-NLSFs and bandby-band Value-NLSFs, the index-by-index (or band-by-band) frequency response is projected back to the graph domain. Consequently, for graph classification tasks, we apply the readout function as defined for Pooling-NLSFs in Sec. In App. C, we presented the diag-NLSFs, considering d = d, and lead-NLSFs for leading filters that do not include orthogonal complements. We explore the leading filters, orthogonal complement, and diagonal operation in our NLSFs. The results of these investigations are summarized in Tab. 15 , which shows the node classification accuracy achieved using various configurations, including NLSFs, diag-NLSFs, lead-NLSFs, and lead-diag-NLSFs. We see that incorporating both the orthogonal complement and the multi-channel approach yields the highest classification accuracy.Our primary objective in this work is to introduce new GNNs that are equivariant to functional symmetries, based on a novel spectral domain transferable between graphs using analysis and synthesis. Our NLSFs in Sec. 3.4 consider filters g j that are supported on the dyadic sub-bands λ N r S-j+1 , λ N r S-j . The closer to the low-frequency range, the denser the sub-bands become. We illustrate an example of filters g j supported on dyadic sub-bands with r = 1 2 and S = 4 in Fig. 5 , showing that sub-bands become denser as they approach the low-frequency range. Our primary goal is to highlight the unique aspects of our method rather than sub-band separation, which, while crucial, is an important consideration across spectral GNNs. Therefore, we also present uniform sub-bands, where filters g j are supported on the uniform sub-bands (j -1) λ N S , j λ N S . Note that the modifications required for our NLSFs are minimal, and most steps can be seamlessly applied with filters supported on the uniform sub-bands. Equivariant GNNs are designed to handle graph data with symmetries. The output of an equivariant GNN respects the same transformation applied to the input. Depending on the application, the transformation could involve translations, reflections, or permutations, to name but a few [82, 46, 78] . Due to respecting the symmetries, equivariant GNNs can reduce model complexity and improve generalization [35, 27] , which can be applied to test data and produce more interpretable representations [99, 51] . When symmetry plays a critical role, such as in physical simulations, molecular modeling, and protein folding, equivariant GNNs have been demonstrated to be effective [20] . For example, [91, 33] employ spherical harmonics to handle 3D molecular structures, and [82] simplify computations for explicit rotation and translation matrices by focusing on relative distances between nodes. In addition, [30] embeds symmetry transformations by parameterizing convolution operations over Lie groups, such as rotations, translations, and scalings, through Lie algebra.• 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 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. 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 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: The paper poses no such risks. Guidelines:• 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. 12. 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] Justification: We provide credit to the original papers of all external data, code, and models used in this work, ensuring that their contributions are acknowledged and their licensing terms are followed. Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL.• The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. • assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. • If this information is not available online, the authors are encouraged to reach out to the asset's creators. 
new assets: 
 Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets? Answer: [Yes] Justification: We provide detailed descriptions of our method's implementation within the paper and report the hyper-parameters in App. E. Guidelines:• The answer NA means that the paper does not release new assets.• Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc. • The paper should discuss whether and how consent was obtained from people whose asset is used. • 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 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.The operator B is normal iff B * B = BB * . Equivalently, iff B has an orthogonal eigendecomposition.A ∆ /dist ∆ operates on an equivalence class [X] of signals by applying A ∆ on an arbitrary element Y of [X]. The output of A ∆ on Y does not depend on the specific representative Y, but only on [X]. 
 
 Appendix:  


Decision: The paper introduces Nonlinear Spectrum Filters (NLSFs), a novel filter operation for Graph Neural Networks (GNNs) that is equivariant to functional symmetries. The authors propose node-level, graph-level, and pooling non-linear spectral filters, demonstrating that these filters outperform standard convolutional GNNs in both semi-supervised node classification and graph classification tasks. The paper is well-supported by rigorous theoretical analysis and extensive empirical validation across various benchmarks.

While there were minor concerns about the mathematical density of the paper, which may limit accessibility, the reviewers acknowledged that the authors effectively addressed the most relevant concerns. There was limited discussion among reviewers, with the primary concern from the "borderline accept" reviewer focusing on the scalability of the method. The authors addressed this by providing results on larger datasets.

Given the overall positive feedback from the reviewers, particularly regarding the paper’s novelty, theoretical rigor, and comprehensive experimental validation, I recommend accepting the paper.

Reviews: review: paper_summary: paper_summary: The paper tackles the task of Network design for graph neural networks. The suggested approach is based on spectral properties of graphs. So far in the literature spectral methods were limited in assuming that the graph domain is fixed. To address this, a relaxed version of symmetry is proposed based on band-limited projections. In addition, a nonlinear spectral filter design is suggested, suggesting node-level, graph-level, and pooling operations. The method is evaluated on several graph learning tasks, demonstrating improvement in generalization over existing spectral methods.
strength_weakness: strength_weakness: Writing Quality: Some sections of the manuscript could benefit from revision. For instance, reordering the paragraphs in the introduction could improve readability. Specifically, mentioning what was missing from previous works earlier rather than at the end would help.

More examples can be found in the method section: i) The discussion on problem the activation functions is missing some details, e.g., what rho is exactly? ii) The paper states that "It is important to note that functional shifts, in general, are not induced from node permutations. In stead, functional shifts are related to the notion of functional maps...". This sentence is too vague. Consider adding more details to make it clearer.

No qualitative results are provided. Is it possible to visualize learned features as in the illustration in figure 2? Is it possible to design a toy experiment showcasing the suggested notion of relaxed symmetry, for which the suggested network design generalizes adequately?
questions: questions: No question other than the weakeness stated above.
limitations: limitations: Yes
paper_summary: paper_summary: The authors propose nonlinear spectral filters (NLSFs) that achieve full equivariance to graph functional shifts, demonstrating that these filters have universal approximation properties. These NLSFs are designed based on transferable spectral domain, potentially improving GNN performance in node and graph classification tasks across diverse graph structures.
strength_weakness: strength_weakness: Despite these merits, I have the following concerns about the paper.

1- While the paper presents a compelling method with potential applications in graph analysis, one significant limitation is its scalability, particularly concerning large-scale graphs. The reliance on specific spectral properties, such as the leading eigenvectors, may not only limit the method's capacity to capture diverse graph dynamics but also result in computational inefficiencies when applied to extensive graph datasets.

2- The datasets used in the paper predominantly consist of mid-sized, homophilic graphs, which may not fully represent the diverse range of real-world applications, particularly in contexts involving heterophilic graphs.

3- The efficiency of the proposed models in terms of computation and resource utilization is not adequately discussed.
questions: questions: (i) Does your theory adapt differently when applied to heterophilous graphs compared to homophilous graphs, and if so, how are these differences addressed within your methodology?

(ii) Given that your Nonlinear Spectral Filters (NLSFs) are motivated by respecting graph functional shift symmetries, similar to Euclidean CNNs, do you have any claims or observations regarding how NLSFs fit within or potentially extend the Weisfeiler-Lehman hierarchy of expressivity? Additionally, could you elaborate on how the expressivity of NLSFs, as informed by metrics from the Euclidean vector space, compares to traditional graph neural network models?
limitations: limitations: Also, the efficiency of the proposed models in terms of computation and resource utilization is not adequately addressed. For practical applications, especially in resource-constrained environments, understanding the computational overhead is essential.
paper_summary: paper_summary: The authors introduce spectral GNNs that are equivariant to functional symmetries. Specifically, they introduce node-level, graph-level and pooling non-linear spectral filters and show that these are able to outperform standard convolutional GNNs on (semi-supervised) node classification and graph classification tasks.
strength_weakness: strength_weakness: - While the authors did a good job in trying to introduce all the relevant concepts, the paper is quite dense with mathematical details and notions that will likely be unfamiliar to many GNN researchers and may therefore hinder the accessibility of the manuscript.
questions: questions: - While it's very intuitive to understand what is meant by "shift" in the context of images and CNNs, this doesn't come across very clear in the paper in the context of graphs: what is the rationale behind the decision to "model the group of translations on graphs as the group of all unitary operators on signals that commute with the graph shift operator"? If more space in the paper can be used to make the underlying concepts more accessible (perhaps moving some of the material on the theoretical properties of the NLSFs to the appendix) I think the paper would greatly gain in accessibility, potentially increasing its impact beyond the graph signal processing community.
limitations: limitations: The limitations are clearly discussed.
paper_summary: paper_summary: This paper proposes a spectral GNN called non-linear spectral filters (NLSF), which aims to enhance GNNs with nonlinear functions. Since general GNNs with nonlinear functions do not commute with unitary operators, this paper defines Graph Functional Shifts, which is a set of unitary matrices commuting with a normal graph shift operator (GSO). It then formulates two functions for spectral index and filter bank, respectively, and concatenates these two functions as graph attention. In the experiment section, NLSF is compared with GAT, SAGE, and other spectral-like GNNs. In the node-classification task, att-Node-level NLSF shows outstanding performance among these models. In the graph classification task, att-Graph-level NLSF achieves comparable results with these models, and att-Pooling-NLSF performs better than other models in the graph classification task.
strength_weakness: strength_weakness: 1. The method proposed in this paper requires the use of eigenvalues, hence it necessitates eigen decomposition of the GSO. The time complexity of eigendecomposition is relatively high, especially for very large graphs.
questions: questions: None
limitations: limitations: N/A
paper_summary: paper_summary: This paper addresses the issue of inadequate modeling of graph equivariance in existing spectral GNNs due to nonlinear operations. The authors investigate the concept of domain translation in graph space as functional translations, drawing from the convolutional operations defined on images. Based on a series of in-depth analyses, they propose nonlinear spectral filters (NLSFs) that are fully equivariant to graph functional shifts and demonstrate universal approximation capabilities.
strength_weakness: strength_weakness: 1. The paper is missing important spectral GNN models such as [1,2,3,4].
2. There is a lack of discussion on equivariant GNNs, such as [5,6,7,8,9,10,11,12]. While the focus is on spectral GNNs, it is also essential to discuss works in the spatial GNN domain, especially given that the paper’s title starts with "Equivariant Machine Learning on Graphs." Positioning your work within the broader GNN field can further elucidate the significance and unique advantages of your contributions.
3. The experiments are somewhat lacking in comprehensiveness. While you have demonstrated the effectiveness of your proposed model in common node and graph classification tasks, your core contribution is enhancing equivariant representation learning in spectral GNNs. I suggest adding experiments that specifically show how the improved performance is due to enhanced equivariant learning, such as graph isomorphism tests.

[1] How powerful are spectral graph neural networks

[2] Bernnet: Learning arbitrary graph spectral filters via bernstein approximation

[3] Specformer: Spectral graph neural networks meet transformers

[4] Graph neural networks with learnable and optimal polynomial bases

[5] Universal Invariant and Equivariant Graph Neural Networks

[6] E(n) Equivariant Graph Neural Networks

[7] On the Generalization of Equivariant Graph Neural Networks

[8] Expressive Power of Invariant and Equivariant Graph Neural Networks

[9] Approximately Equivariant Graph Networks

[10] Equivariant Polynomials for Graph Neural Networks

[11] Graph Neural Networks for Learning Equivariant Representations of Neural Networks

[12] Sign and Basis Invariant Networks for Spectral Graph Representation Learning
questions: questions: see weaknesses
limitations: limitations: yes, the authors discuss the limitations in Section 6

Comments: title: reply to authors
comment: We warmly thank the reviewer for their comments and for acknowledging our rebuttal.
I appreciate the clarifications and support the proposed updates in the camera ready version of the paper.
We warmly thank the reviewer for their comments and for acknowledging our rebuttal.
We warmly thank the reviewer for their comments and for acknowledging our rebuttal.
We thank the reviewer for their comments and for acknowledging our rebuttal. We would be happy to address any remaining concerns if the reviewer has any.
Thank the authors for their extensive rebuttal which has answered some of my concerns. After reading the rebuttal and other reviewers' comments, I decided to keep my current rating.
Thanks for the authors' detailed rebuttal. As all of my concerns have been addressed, I have raised my score.
I appreciate the author’s rebuttal and believe that their suggestions would greatly benefit the paper. I have no further concerns.