We thank the reviewer for the positive feedback and also happy to hear that you enjoyed reading the article. We address the specific questions in detail below:

• The method is just regression with a specific GP, Theorem 4.1 seems (at least intuitively) to be clear. While there is nothing wrong with the paper, I did miss some deeper insights.

We thank the reviewer for this comment. To emphasize the strengths of our approach and contrast it with traditional GP regression on graphs, we briefly review prior graph‑GP methods. Most of these methods place a full prior

$f \sim \mathcal{GP}\bigl(0,K\bigr),\quad \text{where} \quad K\in\mathbb{R}^{N\times N},$

where $K$ may be a graph aware Matérn kernel [1] or a wavelet‐inspired kernel [2] to name a few. Inference is done by obtaining the posterior distribution and it involves inverting $K + I$ matrix, incurring $\mathcal{O}(N^3)$ time and $\mathcal{O}(N^2)$ memory. While they can quantify uncertainty—unlike GCN or GAT, which yield only point estimates. GP based models still exhibit some key limitations: 1. Lack of feature‑wise interpretability, since every prediction mixes all features and nodes through the dense inverse; 2. Scalability due to matrix inversion operation; 3. Cannot be linearized, because the most of these graph aware kernels do not satisfy the conditions stated by Bochners theorem. Importantly, these black box approaches to be interpretable/explainable requires a post-hoc explainers that adds complexity and lacks transparency. In contrast, our current framework is self interpretable by design.

Our framework is inspired by additive models and achieves feature‑level interpretability. We treat each feature independently (also check below for discussion on the correlated features), place a GP prior on its graph‑filtered values, and choose a kernel that satisfies Bochner’s conditions for shift invariance. More specifically, our kernel is RBF that operates on the filtered feature vector $\tilde x_i$, it remains shift‑invariant —unlike conventional graph kernels. This design admits a RFF linearization, yielding a lightweight, linear‑in‑parameters model with $\mathcal{O}(NMD)$ forward‑pass complexity instead of costly matrix inversion.

Importantly this is the first work that introduces learnable FIR filters into the GP framework where each coefficient $\alpha_h$ directly quantifies the contribution of the $h$-hop neighborhood, providing structural interpretability alongside feature‑wise insights. To the best of our knowledge, this is the first work to provide a theoretical justification for embeddings that are not only permutation equivariant but also stable under graph perturbations, owing to the use of finite impulse response (FIR) filters.

Finally, Theorem 4.1 states permutation equivariance property that captures the graph symmetry. We emphasize that this property is non-trivial due to non-linearity involved (contrary to any GP based regression algorithms). The novelty again lies in integration of FIR filters and RFF based approximations, which induces the structure dependencies while maintaining the permutation equivariance.

In summary, G‑NAMRFF transforms a black‑box, cubic‑cost GP regression into a glass‑box, linear‑cost, interpretable model, going well beyond standard GP approaches.

[1] Borovitskiy et al., “Matérn Gaussian processes on graphs,” AISTATS ’21.

[2] Fang et al., “Wavelet kernels for graph Gaussian processes,” ICML ’21.

• Comment on expressivity

We thank the reviewer for the insightful comment. As pointed out correctly, G-NAMRFF models the output as a sum of univariate functions, each acting on a single feature independently. This design avoids feature mixing basically for two reasons: (1) To reduce the parameter complexity and design a very light-weight model (2) To position the model in line with recent interpretable models such as NAMs and GNAN. We acknowledge that this imposes a limitation on expressivity, as it does not explicitly capture inter-feature correlations. Despite this constraint, G-NAMRFF being light-weight matches or outperforms several black-box GNNs across benchmarks, demonstrating that the core design is effective for many practical tasks.

To further enhance expressivity while preserving interpretability, we highlight that G-NAMRFF can be readily extended using ideas from NODE-GAM [1] and GAMI-Net [2], which incorporate structured pairwise interactions. Specifically, G-NAMRFF can model both independent and pairwise feature effects as:

$y_{i} = \sum_{k=1}^{D} f_{k}(x_{i,k}) + \sum_{k=1}^{D} \sum_{k'>k}^{D} f_{k,k'}(x_{i,k},x_{i,k'})$
where $f_{k,k'}$ is modeled as GP is modeled as a Gaussian Process with a 2D input. These functions can be efficiently approximated using 2D random Fourier features (RFFs) as

$f_{k,k'}(x_{i,k},x_{i,k'}) = \mathbf{\Phi}_ {a}^{T}([\tilde{x}_{i,k},\tilde{x} _{i,k'}])\mathbf{w} _{k,k'}.$

The overall prediction of node $i$ then becomes

$y_{i} = \sum_{k=1}^{D} \mathbf{\Phi}_ {a}^{T}(\tilde{x}_{i,k}) \mathbf{w} _{k} + \sum _{k=1}^{D} \sum _{k'>k}^{D} \mathbf{\Phi} _{a}^{T}([\tilde{x} _{i,k},\tilde{x} _{i,k'}])\mathbf{w} _{k,k'}$

As discussed earlier, this further increases the complexity due to increase in the number of learnable parameters. We again emphasize that the main motivation involved in design of G-NAMRFF is the low-complexity without compromising the performance compared to the black-box models.

In the final draft, we will add a remark stating how G-NAMRFF can be extended to handle correlated features.

[1] Chang, Chun-Hao, Rich Caruana, and Anna Goldenberg. "NODE-GAM: Neural Generalized Additive Model for Interpretable Deep Learning." International Conference on Learning Representations 2022.

[2] Yang, Zebin, Aijun Zhang, and Agus Sudjianto. "GAMI-Net: An explainable neural network based on generalized additive models with structured interactions." Pattern Recognition 120 (2021): 108192.

• Could you make this deep, by staking G-NAMRFF layers on top of each other?

Yes — in principle, G-NAMRFF layers can be stacked to increase model capacity. However, doing so introduces significant interpretability overhead. Each additional layer brings its own set of FIR filter coefficients $\{\alpha_h^{(\ell)}\}$ and RFF weight vectors $\{w_k^{(\ell)}\}$, making it harder to attribute predictions to specific input features or graph structure components.

Unlike message-passing GNNs such as GCNs, which expand their receptive field across layers via

$X^{(\ell+1)} = \sigma\bigl(AX^{(\ell)}W^{(\ell)}\bigr),$

our model uses explicit FIR graph filtering, which already aggregates $R$-hop neighborhood information in a single layer. Therefore, deeper stacking is not necessary purely for increasing receptive field — a single G-NAMRFF layer with properly tuned filter order $R$ and RFF dimension $M$ is often sufficient.

Empirically, we evaluate stacking multiple G-NAMRFF layers on node classification task on Cora dataset and observe only marginal gains or even performance degradation. We attribute this to two competing effects: (1) oversmoothing, where repeated graph filtering causes node features to become indistinguishable, and (2) overfitting from stacking multiple layers of high-dimensional RFF projections.

In summary, while stacking G-NAMRFF is technically feasible and may increase function class capacity, it compromises interpretability and offers limited performance benefits. Our results suggest that a carefully designed single-layer G-NAMRFF already achieves competitive performance with far greater transparency.

• I did not understand Theorem 2.1 (for instance what do you mean with approximated?) and the statement could be wrongly formulated. Could you provide a reference or a proof?

We thank the reviewer for the comment. The statement of Theorem 2.1 is a standard result: By Bochner’s theorem, a shift‑invariant kernel can be expressed as $K(x,x') = \mathbb{E}_{a,b}(\mathbf{\Phi} _{a}^{T} (x) \mathbf{\Phi} _{a}(x))$. The finite dimensional inner product is therefore an unbiased Monte carlo estimate of true kernel approximating it to the sampling error. Complete proof along with the convergence and the concentration bounds can be found in Theorem 1 in [1].

We will add this reference in the revised manuscript.

[1] Rahimi, Ali, and Benjamin Recht. "Random features for large-scale kernel machines." Advances in neural information processing systems 20 (2007).

• Comment on limitations

We thank the reviewer for the suggestion. Following your suggestion, we will include a section in the revised draft mentioning limitations and future directions as follows:

\section{Limitations and Future Directions}

G-NAMRFF, as discussed, is primarily developed under the assumption of a static graph structure and static node features. While this is consistent with most existing interpretable graph learning frameworks, it limits applicability to dynamic real-world data from social or communication networks, where both the graph topology and node attributes evolve over time.

Extending G-NAMRFF to dynamic graphs is a promising direction for future work. In particular, both continuous-time dynamic graphs (CTDGs) and discrete-time dynamic graphs (DTDGs) are widely used to model temporal graph data. Developing such lightweight, interpretable models for dynamic graphs remains an exciting and important research challenge.

We thank the reviewer for their time and for the positive feedback on the paper.

We thank the reviewer for the positive feedback on the model strengths and the presentation of the paper. We address the specific questions in detail below.

• Comment on typos

We thank the reviewer for the thorough reading of the manuscript and pointing out the typos. We will include all the suggestions in the revised manuscript. More specifically, we will add the details such as label for X-axis, references for the Theorem 2 and carefully check for the other typos as well.

• Comment on time complexity

The complexity of G-NAMRFF is dominated by two key steps as outlined in Algorithm 1 in supplementary material, namely (1) Filtering (2) RFF approximation. For a graph with $N$ nodes and $D$ dimensional features, we have Laplacian as $\mathbf{L} \in \mathbb{R}^{N \times N}$, graph data as $\mathbf{X} \in \mathbb{R}^{N \times D}$. Then the graph filtering operation costs $\mathcal{O}(N^{2}D)$. Whereas the RFF approximation costs $\mathcal{O}(NMD)$. Hence the total time complexity for dense graph i.e., worst case complexity of the proposed algorithm is $\mathcal{O}(N^{2}D + NMD)$. Exploiting the sparsity in the Laplacian, let say if $\mathbf{L}$ has atmost $d$ (aka degree) non-zero entries the graph filtering operation costs $\mathcal{O}(dND)$. Total time complexity is $\mathcal{O}(dND+NMD)$.

For dense graphs, the time complexity of GCN with $K$ layers is $\mathcal{O}(KN^2F + K N F^2)$. Exploiting the sparsity the time complexity for GCN is given by $\mathcal{O}(KEF+ KNF^2)$, where $E$ is the total number of edges in the graph. So the time complexity is comparable to the black-box model.

• The authors mention that GNAN requires O(N+E) distance computations per node. Why is this considered 'quadratic scaling'?

As explained in the introduction, GNAN requires computing the distance of a node $i$ to all $N$ nodes in the graph, leading to a per-node complexity of $\mathcal{O}(N + E)$. Aggregating over all $N$ nodes results in an overall complexity of $\mathcal{O}(N(N + E))$. In the case of a complete graph, where $E = \mathcal{O}(N^2)$, the total complexity becomes $\mathcal{O}(N^3)$, making GNAN computationally expensive. This quadratic (or cubic) scaling limits its scalability even if sparsity can be exploited in more realistic sparse graphs.

We will make the code public including all the suggestions in the final draft.

We thank the reviewer for their time and the very positive feedback on the paper.

We thank the reviewer for the very positive feedback on the strengths of the work and the presentation of the paper. We address the specific questions in detail below.

• Comment on expressivity

We thank the reviewer for the insightful comment. As pointed out correctly, G-NAMRFF models the output as a sum of univariate functions, each acting on a single feature independently. This design avoids feature mixing basically for two reasons: (1) To reduce the parameter complexity and design a very light-weight model (2) To position the model in line with recent interpretable models such as NAMs and GNAN. We acknowledge that this imposes a limitation on expressivity, as it does not explicitly capture inter-feature correlations. Despite this constraint, G-NAMRFF being light-weight matches or outperforms several black-box GNNs across benchmarks, demonstrating that the core design is effective for many practical tasks.

To further enhance expressivity while preserving interpretability, we highlight that G-NAMRFF can be readily extended using ideas from NODE-GAM [1] and GAMI-Net [2], which incorporate structured pairwise interactions. Specifically, G-NAMRFF can model both independent and pairwise feature effects as:

$y_{i} = \sum_{k=1}^{D} f_{k}(x_{i,k}) + \sum_{k=1}^{D} \sum_{k'>k}^{D} f_{k,k'}(x_{i,k},x_{i,k'})$, where $f_{k,k'}$ is modeled as a Gaussian Process with a 2D input. These functions can be efficiently approximated using 2D random Fourier features (RFFs) as

$f_{k,k'}(x_{i,k},x_{i,k'}) = \mathbf{\Phi}_ {a}^{T}([\tilde{x}_{i,k},\tilde{x} _{i,k'}])\mathbf{w} _{k,k'}.$

The overall prediction of node $i$ then becomes

$y_{i} = \sum_{k=1}^{D} \mathbf{\Phi}_ {a}^{T}(\tilde{x}_{i,k}) \mathbf{w} _{k} + \sum _{k=1}^{D} \sum _{k'>k}^{D} \mathbf{\Phi} _{a}^{T}([\tilde{x} _{i,k},\tilde{x} _{i,k'}])\mathbf{w} _{k,k'}$

As discussed earlier, this further increases the complexity due to increase in the number of learnable parameters .

In the final draft, we plan to add this as a remark.

[1] Chang, Chun-Hao, Rich Caruana, and Anna Goldenberg. "NODE-GAM: Neural Generalized Additive Model for Interpretable Deep Learning." International Conference on Learning Representations 2022.

[2] Yang, Zebin, Aijun Zhang, and Agus Sudjianto. "GAMI-Net: An explainable neural network based on generalized additive models with structured interactions." Pattern Recognition 120 (2021): 108192.

• Comparison of number of parameters with other baselines

We thank the reviewer for the comment. Although in the main paper we compared the parameter efficiency against GNAN, here we present the comparison against the other black box and glass box models. In particular, we present the ratio of number of parameters of black/glass box models to G-NAMRFF. Considering $H_u$ as size of the hidden layer, $D$ as input feature dimension, $M$ as number of RFFs, $R$ as filter order, $L$ as number of GCN layers or MLP layers and $K$ as number of heads in GAT, we present the parameter efficiency in the below table where it is clear that proposed model is parameter efficient.

• Comment on equation 1

We thank the reviewer for pointing out the typo. We will correct this while submitting the final draft.

• Comment on limitations

We thank the reviewer for the suggestion. Following your suggestion, we will include a section in the revised draft mentioning limitations and future directions as follows:

\section{Limitations and Future Directions} G-NAMRFF, as discussed, is primarily developed under the assumption of a static graph structure and static node features. While this is consistent with most existing interpretable graph learning frameworks, it limits applicability to dynamic real-world data from social or communication networks, where both the graph topology and node attributes evolve over time.

Extending G-NAMRFF to dynamic graphs is a promising direction for future work. In particular, both continuous-time dynamic graphs (CTDGs) and discrete-time dynamic graphs (DTDGs) offer formal models for temporal graph data. The modular nature of G-NAMRFF—particularly its use of FIR graph filters and additive modeling—suggests that it can be adapted to these settings by incorporating time-varying graph filters and temporally evolving RFF features. Developing such lightweight, interpretable models for dynamic graphs remains an exciting and important research challenge.

We thank the reviewer for the detailed review and also for the positive feedback on the strengths of the proposed model. We address specific questions in detail below:

• The effectiveness of the model seems to depend on careful fine-tuning of a large set of hyperparameters.

The proposed G-NAMRFF model includes a limited and interpretable set of hyperparameters: the filter order ($R$), number of random Fourier features ($M$) and kernel bandwidth ($\theta$). Among these $R$ controls the receptive field depth and is typically constrained to small values (e.g., 1–5), as higher values lead to oversmoothing, a well-established limitation in graph representation learning, $M$ and $\theta$ are kernel approximation parameters which are common in any Random Fourier Feature (RFF) or Gaussian Process-based methods and are not specific to our model.

We further highlight that compared to GNAN, our model has significantly less hyperparameters. GNAN requires separate deep neural networks for each feature and for neighborhood weighting, each with tunable number of layers, hidden units, and additional regularization mechanisms (e.g., dropout), increasing tuning complexity.

We also emphasize that the performance of G-NAMRFF is relatively robust across different settings of hyperparameters as demonstrated in the supplementary material. This clarifies that G-NAMRFF is both parameter-efficient and more manageable in terms of tuning compared to prior interpretable baselines.

• The interpretability aspect is limited to static data

G-NAMRFF, as correctly pointed out, is a self-interpretable and lightweight model developed primarily for static graphs. Our study explicitly focuses on static graph settings, as detailed in the paper, given the early stage of research in interpretable graph models—where, to the best of our knowledge, only one prior self-interpretable method (GNAN) has been proposed. The goal of this work is to develop a first-of-its-kind interpretable and parameter-efficient architecture that performs on par with black-box GNNs, while offering real-time interpretability.

That said, the modular design of G-NAMRFF readily allows for extensions to streaming or temporal graph data. In particular, the model’s core strength lies in its use of graph-aware FIR filters, which are functions of the graph Laplacian. In the dynamic setting—commonly modeled via Discrete-Time Dynamic Graphs (DTDGs)—both the graph Laplacian $\mathbf{L}$ and node features evolve with time. Existing works [1, 2] in graph signal processing (GSP) has extensively studied time-varying graph filters which can be directly used to generate time-varying filtered features. These filtered signals can then be passed through the same RFF-based univariate modeling framework, preserving interpretability while capturing temporal dynamics. While our current implementation is for static graphs, the architectural components are well-positioned for principled extension to dynamic graphs which we consider an exciting direction for future work. We will add this in a future direction section while submitting the final draft.

[1] Isufi, Elvin, et al. "Filtering random graph processes over random time-varying graphs." IEEE Transactions on Signal Processing 65.16 (2017): 4406-4421.

[2] Bohannon, Addison W., Brian M. Sadler, and Radu V. Balan. "A filtering framework for time-varying graph signals." Vertex-Frequency Analysis of Graph Signals. Cham: Springer International Publishing, 2018. 341-376.

• The effectiveness of G-NAMRFF beyond graph or node classification tasks is not clear.

We thank the reviewer for the comment. Following prior works in interpretable and explainable graph learning we focused on node and graph classification tasks to enable a fair and direct comparison. However, the core strength of G-NAMRFF lies in its design: graph signal filtering followed by additive modeling via RFFs. This structure is task-agnostic and readily extends to other downstream tasks, such as link prediction and graph regression. To validate this, we conducted additional experiments with downstream task as link prediction. Given a node pair (say $u,v$), we compute their embeddings using G-NAMRFF and predict the likelihood of an edge via a dot product decoder followed by a sigmoid activation. Importantly, since node embeddings are formed through additive combinations of filtered, RFF-transformed features, the model retains feature-level interpretability in the link prediction setting. Specifically, we can decompose the link score into contributions from individual features and hop depths for each node, allowing us to quantify which features most strongly influence the prediction. We quantitatively present the results on benchmark datasets with AUC as a metric. From the results we emphasize that the G-NAMRFF is task agnostic and offers interpretability without compromising the accuracy.

• Comment on limitations

We thank the reviewer for the suggestion. Following your suggestion, we will include a section in the revised draft mentioning limitations and future directions as follows:

\section{Limitations and Future Directions}

G-NAMRFF, as discussed, is primarily developed under the assumption of a static graph structure and static node features. While this is consistent with most existing interpretable graph learning frameworks, it limits applicability to dynamic real-world data from social or communication networks, where both the graph topology and node attributes evolve over time.

Extending G-NAMRFF to dynamic graphs is a promising direction for future work. In particular, both continuous-time dynamic graphs (CTDGs) and discrete-time dynamic graphs (DTDGs) are widely used to model temporal graph data. The modular nature of G-NAMRFF—particularly its use of FIR graph filters and additive modeling—suggests that it can be adapted to these settings by incorporating time-varying graph filters and temporally evolving RFF features. Developing such lightweight, interpretable models for dynamic graphs remains an exciting and important research challenge.

Dear Area Chair,

I have considered all aspects of the authors' rebuttal (which is in line with the essence of the reviewing process) and chosen to retain my original rating. The concerns raised by me in the review are corroborated by authors' response, which limits my recommendation for this work.

Thank you.