On The Expressive Power of GNN Derivatives
Authors: Yam Eitan
Moshe Eliasof
Yoav Gelberg
Fabrizio Frasca
Guy Bar-Shalom
Haggai Maron
Venue: ICLR 2026 Poster (2026)
Keywords: Graph Neural Networks
GNNs
Expressivity
Message Passing
Geometric deep learning
Differential geometry
Symmetry
Primary area: learning on graphs and other geometries & topologies
TL;DR: We present High Order Derivative-GNN (HOD-GNN), a method that leverages the input derivatives of a base MPNN to enhance the expressive power of a downstream GNN by injecting  derivative-informed node features.

Abstract: Despite significant advances in Graph Neural Networks (GNNs), their limited expressivity remains a fundamental challenge. Research on GNN expressivity has produced many expressive architectures, leading to architecture hierarchies with models of increasing expressive power. Separately, derivatives of GNNs with respect to node features have been widely studied in the context of the oversquashing and over-smoothing phenomena, GNN explainability, and more. To date, these derivatives remain unexplored as a means to enhance GNN expressivity. In this paper, we show that these derivatives provide a natural way to enhance the expressivity of GNNs. We introduce High-Order Derivative GNN (HOD-GNN), a novel method that enhances the expressivity of Message Passing Neural Networks (MPNNs) by leveraging high-order node derivatives of the base model. These derivatives generate expressive structure-aware node embeddings processed by a second GNN in an end-to-end trainable architecture. Theoretically, we show that the resulting architecture family's expressive power aligns with the WL hierarchy. We also draw deep connections between HOD-GNN, Subgraph GNNs, and popular structural encoding schemes. For computational efficiency, we develop a message-passing algorithm for computing high-order derivatives of MPNNs that exploits graph sparsity and parallelism. Evaluations on popular graph learning benchmarks demonstrate HOD-GNN’s strong performance on popular graph learning tasks.

Published as a conference paper at ICLR 2026
ON THE EXPRESSIVE POWER OF GNN DERIVATIVES
Yam Eitan
Moshe Eliasof
Technion - Israel Institute of Technology
Ben-Gurion University of the Negev
University of Cambridge
Yoav Gelberg
Fabrizio Frasca
University of Oxford
Technion Israel Institute of Technology
Guy Bar-Shalom Technion Israel Institute of Technology
Haggai Maron
Technion Israel Institute of Tochnology NVIDIA Research
ABSTRACT
Despite significant advances in Graph Neural Networks (GNNs), their limited expressivity remains a fundamental challenge. Research on GNN expressivity has produced many expressive archilectures, leading to archi- tecture hierarchies with models of increasing expressive power. Separately, derivatives of GNNs with respect to node features have been widely studied in the context of the oversquashing and smothhng phenomena, GNN explainability, and more. To date, these derivatives remain unexplored as a means to enhance GNN expressivity. In this paper, we show that these derivatives provide a natural way to enhance the expressivity of GNNs. We introduce High-Order Derivative GNN (HÖD-GNN), a novel method that enhances the expressivity of Message Passing Neural Networks (MPNNs) by leveraging high-order node derivatives of the base model. These derivatives generate expressive structure-arss aware node embeddings processed by a second GNN in an end-to-end trainable architecture. Theoretically, we show that the resulting architecture family's expressive power aligns with the WL. hierarchy. We also draw deep connections between HOD-GNN, Subgraph GNNs, and popular structural encoding schemes. For computational efficiency, we develop a messago- pagggggganinggggggo algorithm for computing high-order derivatives of MPNNs that exploits graph sparsity and parallelism. Eval- uations on multiple graph learning benchmarks demonstrate HOD-GNN 's excellent performance on popular graph learning tasks.
INTRODUCTION
Graph Neural Networks (GNNs) have become foundational tools in geometric deep learning, with widespread applications in domains such life sciences (Wong et al., 2024h, social sciences (Monti et al. 2019). optimization (Cappart et al. 2023). and more. Despite their empirical success, many GNNs suffer from fundamental limitation: their expressive power is inherently bounded. In particular, the widely used family of Message Passing Neural Net- works (MPNNs) at most as expressive as the Weisfeiler Lehman (I-WL) graph isomorphism test Morris et al. (2019); Xu et al. (2018), limiting their ability to distinguish between even simple non-isomorphic graphs and capture intricate structural patterns Chen et al. (2020) To address this shorteoming, growing body of work has proposed more expressive GNN architectures, typically organized into expressivity hierarchies that balance computational
cost. with representational power Maron et al. (209hh Morris al. 2019 2021

---

Published as a conference paper at ICLR 2026
Concurrently to advances in GNN expressivity, the derivatives of the final node represen-
tations h. and the graph-level output hou with respect to the initial features X, have playod a key role in several research directions. For over-squashing analysis (Di Giovanni et al.)
On!") 2023a b), both first-order derivatives 8x, and mixed partial derivatives 8x.0x. quantify
inter-node inffuence and communication capacity. In wver-smoothigg studies like Arroyo et al.
2hnm 2025), derivatives a.X." are used to analyze vanishing gradients, connecting over-smoothing
to diminished gradient flow. GNN gradient-based explainability methods (Baldassarre & Azizpour 2019 Pope el al.; 2019) also use these derivatives to identify inlluential nodes and leatures. The prevalence of these derivatives across diverse contexts in GNN research suggests they encode valuable information that may be informative for graph learning tasks.
Our approach. In this work, we reveal a surprising connection between these two lines of research. We show that incorporating derivatives of a base MPNN with respect to initial node features as additional inputs to downstream MPNN enhances the expressivity of the base components. One intuitive way to understand this connection is through the mechanism by which GNNs with marking Papp & Wattenhofer 2022 Pellizzoni et al.) (0244) improve expressivity: they choose node from the input graph and add to it. a unique identilier before processing it through an MPNN. While these identifiers are often implemented through an explicit, often discrete, perturbations to the node features, our approach instead computes derivatives of the MPNN output, capturing the effect of infinitesimal perturbations. Thus, giving the model access to derivative information both leveragos quantities which frequently arisc in thcoretical analyses and may thus encode valuable structural information, and yiclds expressivity gains equivalent to GNNs with marking. See Section 3 for details.
We introduce High-Order Derivative GNN (HOD-GNN), a novel expressive GNN family that leverages the derivatives of base MPNN to improve its expressive power. We first introduce 1-HOD-GNN, which consists of three components: a base MPNN, derivative encoder network, and a downstream GNN. I-HOD-GNN computes high-order derivatives
("A(7) and of the base MPNN with respect to the features of a single node at time, i.e., 8X g
8x: These derivatives are then encoded into new derivativo-aware node features via the encoder werw which are then passed to the downstream GNN. Theoretically, We show that 1-HOD-GNN models are more expressive than standard GNNs, can compute popular structural encodings, and are tightly related to Subgraph GNNs Cotta et al. (2021); ;; Bevilacqua al. (2021). Empirically, we demonstrate several desirable properties of our model: it achieves strong performance across a range of standard graph benchmarks, scales to larger graphs that remain out of reach for other expressive GNNs, and can accurately count graph substructuros, providing direct empirical evidence of its expressive power.
We then extend 1-HOD-GNN to H-HOD-GNN, which supports mixed derivatives with respect to l distinct node features (i.e. g011-- 1- ((T) A hout, Like 1-HOD-GNN,
"x.......""" the k-HOD-GNN forward pass begins by computing higher-order mixed derivatives of base MPNN, which form a l--nnexed derivative tensor. This tensor is then used to construct new node features using higher-order encoder network (as in Maron et al. (2018) ; Morris let al. (2019)), ), which are subsequentlyy passed to a downstream GNN for final prediction. We theoretically analyze k-HOD-GNN. showing that it can distinguish between graphs that are indistinguishable to the k-WL test" resulting in model that is more expressive than any of its individual components alone. Furthermore, we leverage results from Zhang et al. 2024b) to analyze K-HOD-GNN's ability to compute homomorphism counts, demonstrating its capacity to capture fine-grainod structure.
Efficiently computing high-order node derivatives is a core component of HOD-GNN. To this end, we develop a novel algorithm for computing these derivatives via an analytic, mossage-passing-like procedure. This approach yields two key benefits. First, being fully analytic, it enables differentiation through the derivative computation itself, allowing
1 (t) is the representation of the node v after the t-th GNN layer. We refer here to the folklore WT. test rather than the oblivious variant; see Morris et al. 2023 for detailed discussion of the differences.
2

---

Published as a conference paper at ICLR 2026
4
oout
Base MPNN
UiD A
1
4
(
h
(
Downetream GNN
h fin
D
Dout
X2
X3
dkj?
db)
PLE
U upi 8
Lppout Uppr%, 4
ax;
6X;
ABPLE
6X
6x,
68.
ox;
Efficient slg
261"
06
6x
AX,
Node-wise Deriv
Out. Deriv
Figure 1: The HOD-GNN pipeline. Given an input graph, we compute the outputs and derivatives of base MPNN. The derivatives are processed by two encoders (denoted U) to produce features that are concatenated with the base MPNN outputs and passed to downstream GNN for final prediction.
HOD-GNN to be trained end-to-end. Second, the message- message-passing- passing-lii like structure exploits the sparsity of graph data, improving scalability (see Section for detailed complexity analysis). Combined with the empirical observation that HOD-GNN remains effective even when using base MPNNs with small hidden dimensions, our method scales to benchmarks containing larger graphs that are often out of reach for other expressive GNN architectures.
Our contributions. (1) We introduce k-HOD-GNN, novel expressive GNN family that integrates derivative-baser embeddings; (2) We provide theoretical analysis of its expressivity and computational properties; (3) We propose an algorithm for efficient derivative computation graphs; (4) We demonstrate consistently high empirical performance across seven standard graph classification and regression benchmarks. Additionally, we show that HOD-GNN scales to benchmarks containing larger graphs that are typically out of reach for many expressive architectures on standard hardware.
2 PRELIMINARIES AND PREVIOUS WORK
Notation. The size of a set is denoted by |S]. @ denotes concatenation. We denote graphs by G = (A, X), where A E gax= is the adjacency matrix and X Raxd is the node feature matrix, with n nodes and d-dimensional features per node. The node set of graph is denoted by V(g).
MPNNs and GNN expressivity, MPNNs(GGimms et al. 2017) are a widcly used class of GNNs that update node representations through iterative aggregation of local neighborhood information. At each layer t, the representation !0 of node v is updated via:
h MLP(() ((((-)) (AG(( [nlt"))uNN(o))]
(1)
where (v) denotes the neighbors of node in the graph and AGG are permutation- preserving aggregation function. After T cassage-passingg layers, graph-level representation typically obtained by applying a global pooling operation over all node embeddings:
h out AGG ( uE v(
|
MPNNs have inherent expressivity limitations Morris et al. 2019 Xu et al. 2018 Weisfeiler Leman 1968). as they cannot distinguish graphs that are indistinguishabbl by the WL test. To address this, wide range of more expressive GNN architectures have been
proposed (Morris et al.] 2021
Maron et al.] 2018 Puny et al. 2023 Cotta et al. 2021 2023) (Sec Appendix A for details or
Rieck et al. 2019 [Sato et al.
2021 Dwivedí al.
Sato 2020 Morris et al. 2021
Jogelkn 2022 Li &
Leskovec 2022:
Zhang ct al. 2024a) for
comprehensive surveys.)
Subgraph GNNs. Subgraph GNNs Zhang Li 2021
Cotta et al.
2021 Bevilacqua et al.
2021 Frasca et. al. 2022 Zhang et al. 2023b
Bar-S
Bar-Shalon
al.
2024b)
expressive
3

---

Published as a conference paper at ICLR 2026
GNNs that operate set of subgraphs Bg = {S, v E (G)), where (g) denotes the set of all k-tuples of nodes in the imput graph G, and each subgraph S, corresponds to one such tuple. In this work, we focus on the widely adopted node-marking DS-GNN (See e.g.
Cotta et al. (2021); Bevilacqua et al. [0021 Papp & Wattenhofer (2022) and their higher- order gencralization, k-OSAN (Qian al. 2022), though note that many other variants
of Subgraph GNNs cxist. For precise defnniion of DS-GNN and K-OSAN, see Appendix E.1
Derivatives of MPNNs.
Derivatives frequently appear in the analysis of GNNs. In the study of oversquashing the failure of information to propagale through graph struc-
tures AAlon & Yahav 2020 Topping al.; play a key role (For a comprehensive overview, see kknnha (2023)) Node derivatives are
2021
Di Giovanni et. al. 2023a b) derivatives
also used in GNN explainability (Ying al.
2019 Luo et al. 2020 Baldassarre & AAiiipour
2019 Pope ct al. Backpropagation
2019). Gradient-based approaches such as Sensitivity Analysis, Guided
Baldassarre & Azizpour 2019), and Grad-CAM (Pope et al., 2019) rely on
derivative magnitudes Finally, several standalone works make use of mode-baser derivatives.
E.g., Arroyo et al. (2025) use node derivatives to draw a connection between vanishing gra-
dients, and over-smoothing, and Keren Taraday et al. (2024) propose aggregation functions, designed to induce non-zero mixed node derivatives. See Appendix A for further discussion.
3 METHOD
We begin this section with a discussion motivating the use of MPNN derivatives and their contribution to improving expressivily. We then introduce k-HOD-GNN, an expressive GNN architecture that enhances representational power by leveraging derivatives of base MPNN. We first present the full details of the I-HOD-GNN model, followed by an overview of its higher-order generalization. A comprehensive treatment of the general k-order case is provided in Appendix C.1. We emphasize that in k-HOD-GNN, the parameter k refers to the number of distinct nodes with respect to which derivatives are taken, not the total hh
but not derivative order. For instance, 1-HOD-GNN uses derivatives of the form 8=X.'
g01+ which involve multiple nodes. mixed derivatives such as Jot Xu, 0 X
3.1 MOTIVATION
Beyond being widely used and informative quantity in GNN analysis, MPNN derivatives can enhance expressivity. To build intuition to why that is, we begin with a simple example, showing that first-order derivatives allow us to count triangles, a task that standard MPNNs cannot perform. Consider the model M(A. = A which can be implemented by three-layer GCN with identity activation. For any node v, the derivative of its final feature vector (T) with respect to its own input feature vector X, is exactly A B.e. aggregating these derivatives, we can compute E. A3 A which is exactly the number of triangles in the graph.
To illustrate how higher-order derivatives further enhance expressivity, we recall that GNNs with marking (Papp & Wattenhofer 2022) improve expressive power over standard MPNNs by selecting a node v in the input graph G = (A, X) and attaching a unique identifier to it, yielding the modified input X tee, for some E R. The oulput. is then h" = _MA, X + ee,.). If M employs an analytic activation function a (See definition E.1 in the appendix), then M itself is analytic. Consequently, its output can be approximated by the Taylor expansion:
m 2' MPAA X + Ie,)
m out
M(A, X + tee))
(3)
Oxt
==0
deg
the output of M without marking. This shows that by leveraging the higher-order derivatives of an MPNN, one can approximate the output of a GNN with marking to arbitrary precision. As a result, derivatives strictly extend the expressive power of MPNNs. An expanded intuitive discussion of these expressivity gains, along with the natural connection between HOD-GNN and subgraph GNNs, is provided in Appendix B
where hout - = X)
While GNNs with marking can select multiple nodes, for clarity we focus the single-node
case. The general case is discussed in Appendix B
4

---

Published as a conference paper at ICLR 2026
3.2 THE 1-HOD-GNN ARCHITECTURE
Overview. A I-HOD-GNN model, denoted v consists of two GNNs: a base MPNN M and a downstream GNN T as well as two derivative encoder networks Unode and U Given an input graph 9, the computation of (()) proceeds in four steps (see Figure D: (1) Compute the final node representations and output of M; (2) Compute derivative tensors of the output with respect to the input node features (defined below); (3) Use Unode and U to extract new derivative informed node features from the derivative tensors. (4) Apply T to the input graph enriched with derivalive-inlormeed features. Importantly, we develop an efficient algorithm for step (2) that. enables backpropagation through the derivative computation itself, making all four above steps differentiable (see Section 3.2.1 Consequently, the entire HOD-GNN model can be trained end-to-end, a strategy we adopt in all our experiments.
Steps 1 & 2. In the first two stages, we compute the final node representations tT)) and the output vector hout using the base MPNN M, along with their corresponding derivative tensors, defined below:
Definition 3.1. Given a graph g g= (A, X) with n nodes, an MPNN M and an intermediate node feature representation matrix h E Raxd", the derivative tensor of D(h) E R R2xaxd xdxm is defined by:
8th ax%
(h))vv,,dd,
where v, u E Vg) are nodes, e [d]].j [d] spocify the feature dimensions of the node feature vectors h. X, respectively, and a E [m] where m E N is a hyperparameter specifying the maximum order of derivatives to be considered. Similarly, given a graph-level prediction vector hout E Rdd the derivative tensor D(hout) E Traxd xdxm is defined by:
(u, = axa a.j
In 1-HOD-GNN, we compute the output derivative tensor D = D(hout), which captures how the output of the base MPNN M responds to perturbations in the input node features. In parallel, we compute the node-wise derivative tensor D(T), where ((t) = D(hl)) for I=,,..., T. These tensors characterize how each node's representation at layer & changes in response to variations in the input fcatures. Derivative tensors are computed using Algorithm [T] described in Section 3.2.1 and claborated on in Appendix [D] The algorithm leverages the sparsity of the input graph to enable efficient computation of high-order derivatives. Crucially, Algorithm [ is fully differentiable with respect to the weights of M, enabling end-to-end training of 4.
Step 3. In the third stage of our method, we extract new node features from the derivative tensors Duut and (() using the encoder networks Unode and Uoue_ First, as Dout e Raxdd xdxm is a tensor indexed by a single node, it can be directly interpreted as a node feature matrix by flattening the remaining dimensions. We thus define the encoder network Tout to be a DeepSets (Zaheer et al. 2017) update
yout (D out ). MLP out (D.
(6)
Secondly, since (T) is indexed by pairs of nodes in G, it shares the structure of the adjacency matrix A, which is also pairwise-indexed. We can thus define the encoder network Unode Rn= xd' xdxm - Rnd to be any GNN architecture which maps adjacency matrices with edge features to node feature matrices.
To enhance sensilivity to global interactions, we select Uuode to be a 2-Invariant Graph Network (IGN)Marrn Maron et al.) 2018 optionally using searsity-preservin simplified variant
formally defined in Appendix C.1
Empirically we observe that setting Uout 0 produces similar results; we nevertheless keep the module in our formulation for completeness.
5

---

Published as a conference paper at ICLR 2026
We construct the derivative-infommee node features h by combining information from the base MPNN M, the pooled intermediate derivatives, and the output derivatives:
h der
hT) @ U ( @ Unode (node(D(T)) (D()).
(7)
Step 4. In the final stage, replace the original node features of with the derivative- informed features dder and apply a downstream GNN T to produce a graph-level prediction. For the remainder of this work, assume T is an MPNN, though our approach is compatible with any GNN architecture.
3.2.1 EFFICIENT DERIVATIVE TENSOR COMPUTATION
We now describe efficient algorithm for computing the derivative lensors in 1-HOD- GNN model 4 with base MPNN M. For clarity, we focus on the case where M is a GIN Xu et al. 2018). in which case the message- saseaage-pssingg and readout functions are given by:
h = MLP ((1 + e)h +
Σ 1 !))
h MLP
Σ
(8)
uEN(()
NEV(S)
An extension of this algorithm to general MPNNs and higher-order mixed derivatives is provided in Appendix [D] For convenience, we decompose the node update in Equation into two parts: an aggregation step and DeepSets-based update, given respectively by:
agg. update
DeepSets update
hg-1) -1 eh
Σ wEN(()
h (t-1))
h MLP(h(-"))
(9)
The algorithm is based on the following two observations: First since h (t-1) is linear combination of h((-1) and its neighboring node features (h (h(-1) [Nv)), the derivatives of (t-1) are likewise linear combinations of the derivatives of h (t-1) and (h! (t-1) MuuA ((o)).
More explicitly:
Σ D D(--1) [
(10)
mEN(e)
This computation mirrors the GIN aggregation update in Equation leveraging the sparsity of the graph. Second, since the DeepSets update applies an MLP independently to each node feature h (t-1),
apply Faá di Bruno's formula (see e.g. (Hardy 2006)) to compute the derivatives of cach n based on the derivatives of ((11)) This results in a "DeepSets-like" derivative update, allowing us to compute (( directly from ((-((-))) Iteratively applying these two steps yields the final node-wise derivative tensor ((T) through a differentiable, pressagepasssig--lee procedure. A similar approach allows efficient computation of Dout from (T), See Appendix D for full details of the algorithm.
Computational Complexity, An important property of the above algorithm is that, for sparse graphs or relatively shallow base MPNNs, it is computationally efficient. To see this, first notice that since (0) X, the lensor D (0) is extremely sparse, satislying:
(0) [0, u, v.j.j [1 if i i= a a=
(11)
otherwise.
Thus, it can be stored efficiently using sparse matrices. At each layer t, the derivative aggregation step (Equation 10) increases the number of non-zero entries only in proportion to the number of node pairs that exchange messages for the first time. Thus, for small
values of I or for sparse graphs G, the derivative tensor (t) remains sparse. Moreover, the algorithm's gressagepaassing-aee structure ensures runtime efficiency as well. For full complexity analysis of our algorithm, see Section 4
3.3 K-HOD-GNN VIA MIXED DERIVATIVES
We generalize 1-HOD-GNN, which operates on single-node derivatives (i.e., derivatives (( Aou,
of the form 8xa
0x to k-HOD-GNN, which extracts information from mixed

---

Published as a conference paper at ICLR 2026
partial derivatives across k nodes (i.e.; 8X " T or "x "x k-HOD-GNN offers increased expressive power at the cost of greater computational complexity. We begin by formally defining the k-indexed derivative tensors. For simplicity, we assume the input node features are I-dimensional, handling the more general case in Appendix C.I]
2 2
Definition 3.2. Given a graph g (A, X) with nodes, and an MPNN M and intermediate node feature malrix h E Baxd", the k-indexed derivative tensor of (h) Raxn" xd xm is defined by:
go1t
Du(h)(u, (, a, 0x 0x X
(12)
where E V(g), (u un) E (*) , i E [d"] and aa),..,, E [m] DR is
defined similarly for graph-level prediction vectors.
The k-indexed derivative tensors capture how the output and node representations of M change under joint perturbations to the fcatures nodes, thereby encoding rich higher-order structural interactions within the graph. We compute these lensors using an extension of the derivative computation process described in Section 3.2. (See Appendix D for more details).
As in the 1-HOD-GNN case, k-HOD-GNN model t consists of base MPNN M and a downstream network T. Given an input graph G, the computation of ((G) proceeds in the same four stages established earlier. First, we compute the output and final node representations of M along with the k-indexed derivative tensors D and D We then use h D( and D to extract new derivative informed node features. These node features are computed through:
hdder MLP (N(T) ( Unnde (D @ U (D*tt)
(13)
where Uuode and Uout are learned (k: + 1)-IGN and k-IGN encoders respectively. Finally, we substitute the original node features of graph & with deer, and input the resulting graph into T to generate the final graph-level prediction. For more details, see Appendix C.1
4 THEORETICAL ANALYSIS
In this section, we analyze the expressive power and computational complexity of k-HOD-GNN. Formal statements and complete proofs of all results in this section are provided in Appendix E.2
Expressive power. We begin by formally relating k-HOD-GNN to both K-OSAN subgraph GNNs as well as (k + 2)-IGNs, revealing new insights into HOD-GNNs' expressive power, and their position in the WT, hierarchy.
Theorem 4.1 (informal). Any OSON model can be approzimated by a k- HOD-GNE model using an analytic activation function, to any precision. Additionally, any k- HOD-GNN model can be approzimated by a (k 2)-TGN model.
Corollary 4.2. There exist non-isomo graphs that are indistinguishable by the folklore k- WL (k-FWL) test but are distinguishable by k-HOD--NN GNN.. Additionally, any pair of graphs that is indistinguishable by the (k + F)L test is also indistinguishable by k k-HOD-GNNA
The proof of Theorem 4.1 relics on the analyticity of the activation functions used by our base MPNN and the use of higher-order derivatives. However, in what follows, we show that even when restriccted to first-order derivatives and using the commonly employed ReL.U activation, HOD-GNN remains strictly more expressive than a widely used technique for enhancing GNN expressivity: incorporating Random Walk Structural Encodings (EWSEs)(Dwiver 2021 into base MPNN.
et al.
Theorem 4.3 (Informal). Even when lémiled to firsl-order derivalives and ReLU aclivalions,
1 I-HOD-GNN is strictly more expressive than MPNNs enhanced with random walk structural encodings.
The first part of Theorem 4.3 is constructive: it shows that a simple initialization of the base MPNN's weights yields derivatives equal to RWSEs. In our experiments, we use slightly mod-
7

---

Published as a conference paper at ICLR 2026
ified version of this initialization (see Appendix G allowing HOD-GNN to serve learnable extension of RWSE. Further analysis the expressive power of HOD-GNN when using edge- feature derivatives, or when only using output-level derivatives presented in Appendix F
Space and time complexity. To conclude this section, we analyze the computational complexity of k-HOD-GNN and compare it to other expressive architectures, namely, (k+1)- IGN and k-OSAN. We show that k-HOD-GNN achieves better complexity when using relatively shallow base MPNNs, while maintaining comparable complexity with deeper ones. The primary source of computational overhead in k-HOD-GNN lies in the derivative tensor computation, and the encoder network forward pass. We now analyze each of these components.
First, while k-HOD- GNN computes derivative tensors D (4) with O( potential entries, these tensors are sparse for relatively shallow base MPNNs. Moreover, each D can be efficiently computed from (("")) This is formalized in the following proposition:
Proposition 4.4. In a k-HOODGNT model applied to a graph with nodes and marimum degree d The number of non-2ero entries in (()) at most (n minfn d² Additionally, each D (t) can be computed from (t-1) in time O(d n min (t=11-))).
Focusing next on the encoder networks, we show that they can be designed to exploit derivative sparsity for improved efficiency, while retaining the full expressivity of the k-OSAN architecture:
Proposition 4.5. In a k-HOD-GNN model, the encoder functions Umode and yous can be chosen such that the model retains the expressive power of k-OSAN, while the computation of Umode (D (T) and Dou (Dott) has both time and space complerity O (n .min(n², db-T))-
Propositions 1.4 and 4.5 suggest that K-HOD-GNN model with base MPNN of depth T has space complexity O O(n min(n", dT) and time complexity O (d n min(n*, d=-(T-1)j)) In comparison, K-OSAN has space complexity O(n and time complexity O(d **),, while (K+1)-IGN incurs both time and space complexity of O(n++1). Assuming the input graph is sparse (i.e., d & n), k-HOD-GNN is more efficient then k-OSAN and (k + 1)-IGN when the base MPNN is shallow (d" n), while all three models have comparable complexity when the base MPNN is deep (dd > >n).
5 EXPERIMENTS
Our experimental study is designed validate the thcorctical arguments from the previous section and to address the following guiding questions: (Q1) How does HOD-GNN perform on real-world datasets when compared against strong GNN baselines? (Q2) Can HOD-GNN scale to larger graphs that are beyond the reach of Subgraph GNNs, and how does it perform in this regime? (Q3) How does the expressive power of HOD-GNN compare with natural and widely used GNN baselines? We evaluate HOD-GNN across eight benchmarks, with additional experimental details provided in Appendix
Basclines. We compare HOD-GNN against strong representatives from three natural families of GNNs. First, motivated by its connection to positionaal/struutural encodings (PSEs, Section 4). we consider encoding-augmentted MPNNs, including Laplacian PEs (Dwivedi el al. 2023), RWSEs (Dwivedi el. al.) (021), SiguNet. (Lim et al. 2022). random node identifiers (Abboud et al. 2020 Sato et al. 2021) as well as recent methods such as GPSE (Cantürk et al. 2024] and MOSE (Bao et al. 2024h. Second, since HOD-GNN theoretically related to Subgraph GNNs, we compare with representative models like 2022), and Subgraphormer (Bar-Shalom GNN-AK (Zhao et al.] 2022), SUN (Frasca et al.
et al. 2023). Because such models often struggle to scale, we also include sampling- based 2024). HyMN (Southern et al. 2025 and Subgraphormer with random sampling. Finally, we benchmark against widely used and mod- general-purpose GNNs, including GIN (Xu et al. 2018). GCN (Kipf & Welling 2016). GatedGCN [Bresson & Laurent 2017). GPS (Rampášek et al. 2022), and GraphViT (He
variants such as Policy-Learn (Bevilacqua et al.
et. al. 2023 Across experiments we include representatives from each family, while the
8

---

Published as a conference paper at ICLR 2026
Table 1: Performance on OGB and ZINC datasets (4 sceds). First and sccond best scores are highlighted. Scores sharing a color are not statistically distinguishable based on Welch's I-test with a relaxed Uhreshold of p < 0.2. a denoles results not previously reported, and "x" indicates that digits beyond this point were not provided.
Mechod + Dataset
ZINC-12K (MAE J)
MOLIIV (ROC-AUC t) (ROC-AUC f) (ROC-AUC 1)
MOLTOX21
MOLBACE
Common Baselines
GCN
Kipf& & Welling 2016
0.32110.000
75.2910.69
79.151144 72.97+4.00
76.06:00.97
GIN
Xu al 2018
0. 0.163+00001
74.9110.51
75.58+1.40
PNA
Corso ot al. 2020
0.761=0.002
13.30=
19.05+11.22
GPS Rampaack
0.070=0.004
05.70=00.40
18.80+1001
GraphviT Hi et
0.085+0.005
78.51+0.77
17.92+1.49
Subgraph GNNs
Reconstr. GNN (Cotta et al. 2021
75.15=1.40
76.32+11.0
GNN-AK+ (7hao et al. 2022
0.091+0.011
79.61+11.11
SUN (EGO+F Frasca et al 2022
0.084=0.002
80.03
Full [Bevilacqua el al. 2024
0.087=0.003
76.25#1.12
78.41+1.94
76.54+1.37
OSAN [Qian al. 2122
0.177+0.016
6..30+6.60
Random BBeylarqua
2024
0.102:00.03
76.62-10.63
78.14:12.36
77.30-12.66
Policy-Learn Beviacquu ot a. 2024
0.097-10005
77.36:00.60
78.39112.28
78.49:11.01
Subgraphormer Bar-Shalom of al
2024a
0.063-10.001
84.35:10.65
09.58:0.335
HyMN (outhern et al. 2026
0.080=00003
77.82=0.59
81.16+1.21
81.01+1.17
PSEs
GIN Laplscian PE Dwivedi et al.
2023
0.162=0.014
76 0.
80.40+1.5x
75.60=1.1x
GIN RWSE (Dwivedr et al. 2021
0.128=0.005
76.30=0.55x
29.60=2.8x
18.
SignNet [Lim el 2022
0.102=00002
RNT Abboud dt al
2020
0.136:+0.0070
61. 61.944-2.51
77.74+00.98
GSN
Bouriisas cb al 2022
0.101100010
76.08.10.79
77.40
80.39.10.00
ENGNN (Wang Zhang
7025
0.11410.005
78.51:10.86
GPSE
Canturk et 2024
0.065=0.003
77.40=0.8x
80.80+3.1x
1..13
MOSE
Bao
0.062
Ours HOD-GNN
0.0066:0000035
77.9910.71
82.10+1.45
0..8
specific choice of baselines in each task reflects relevance and standard practice in prior work. Throughout the paper, we report results directly from prior work and include any relevant baseline, even if values for some of the benchmarks were not reported. This allows for a broad and fair comparison rather than excluding useful basclincs. Missing entrics are markod by "L",
OGB and ZINC To evaluate HOD-GNN's real-world performance (Q1), we benchmark it on standard graph property prediction datasets: ZINC (Irwin et al.) 2012) for regression, and three molecular classification tasks from the OGB suite (Hu et al.j 20206] molhiv, mo bbace, and mo1tox21. These benchmarks provide standardized splits and are the de facto choice for assessing GNN performance. As shown in Table [ HOD-GNN delivers excellent results across all tasks, standing out as the only model that consistently ranks within the top two tiers.
Peptides. Section 4 established that HOD-GNN has improved computational complexity compared to Subgraph GNNs. To demonstrate its scalability in practice (Q2) and to further assess its performance on real-world data (Q1), we evaluate HOD-GNN on the Peptides datasets from the LRGB benchmark (Dwivedi et al.) 2022), where the goal is to predict global structural and functional properties of peptides represented as graphs. As stated in prior work (Southern et al.) 2025; Bar-Shalom et al. 2023). full-bag Subgraph GNNs cannot process these graphs directly using standard hardware, requiring the use of subsampling strategies that can reduce expressivity and introduce optimization challenges due to randomness. In contrast, HOD-GNN handles these graphs directly without subsampling. As shown in Table 2 HOD-GNN surpasses all sampling- based Subgraph GNNs and is the only model that consistently ranks within the top two tiers, underscoring both its scalability and elfectiveness on challenging real-world molecular tasks.
Synthetic experiments. To evaluate the realized expressiveness of HOD-GNN (Q3) and empirically support Theorems 4.1 and we conduct two synthetic studies. First, following the protocol of Huang et al. (2022). we assess the ability of I-HOD-GNN to learn to count small substructures, a standard proxy for practical GNN expressivity (Bouritsas ct al. 2022 Arvind et al.) 2020). Thcorems and 4.3 predict that: (i) with analytic
activations, 1-HOD-GNN matches the power of certain Subgraph GNNs, and (ii) with ReLU, it is strictly stronger than MPNNs with RWSEs Table (Appendix confirms both predictions. Additional details are provided in Appendix
9

---

Published as a conference paper at ICLR 2026
We additionally test 1- and 2-HOD-GNN on the regular graph pairs in the BREC dataset (Wang & Zhang 2024], which inchude 50 pairs separable by 3-WL but not 2-WL and 90 pairs Indistinguisaable even by 3-WL, Table 7 shows that HODD-NN separates 34/90 of the -WW--indstinguuihhbbll pairs, placing it among the strongest models and empirically validating the theoretical advantage predicted by Theorem 4.1. Additionally, 1-HOD-GNN performs similarly to DS-GNN, consistent with Theorem 4.1. Additional details are provided in Appendix H.2
Table 2: Performance on PEPTIDES (4 seeds). First and second best scores are highlighted. Same color scores are not statistically distinguishable based on Welch's t-test with a relaxed threshold of p < 0.2.
Additional ablations and empirical insights. Appendix H includes further ablations and anal- ysis. We first evaluate how the hyperparameter m from Definition 3.1 which sets the maximum derivative order used in HOD-GNN, affects ex- pressive power. Results show that increasing m consistently strengthens expressivity, until the per- formance gains eventually plateau. This suggests that small values such effective in practice. We also analyze the stability of HOD-GNN, showing that for derivative orders m = =I,.... 1,, the training curves remain stable and the norms of the derivative tensors stay well- behaved relative to the final node-feature norms produced by the base MPNN. Finally, Inspired by recent work linking overly expressive GNN expres- sivity to poor generalization (Franks et al. 2024 Maskey et al. 2025 Carrasco et al. |0225) we evalu- ate HOD-GNN 's train Lest gaps on OGB datasels. HOD-GNN consistently exhibits smaller gaps than both less expressive (e.g., GIN, GCN) and more expressive (e.g., DSS-GNN (FD)) baselines, indi- cating strong generalization without overfitting. See Appendix H.5 for details.
Poptides-tume fune Peptides-struct AP
Model
MAE
Common Baoulimus
E {2, (243,4) are already
GCN
50.30=0.22 0.23
0.3456=00000
GINE:
56.98-0.77 70
0.3547::00000 0045
GONI
55.43=0.77 78
0.34771=0.00111
GatedGIC
58.64_00.7
0.342001:0 0613
DIGL MPNN LapPE
68.30 26
0.26116. 0018
MisHop- GCN. LapPE
68.40=0.00
0.2614:50.00 0023
Dite..G LapPE
71.50:10.00 44
0.2536
LaN-LapPE
63.84=1.21
0,2683=00.00 00443
GoubGGS LapPE
65.36=0.4
0.2500=0.0000
Exphormer GraphNIT
03.17:0000
0.24812::
09.19:20.85
0.247420.0016
Subgraph GNNS
Policy- Leatn
04.50=0018
0.2476=20.001
Salagrahhommer 50%
64.16m0.55 52
0.2494:20.0020
HyMN
68.57=0.56
0.2464 0015
PSE.
GCN- Laplacian PE
62.18200.0 50
0.249222: 0019
GCN + RWSR
60.87:0.00 0..9
0.257440..00 0.0000
SignNet
63.14:0.30 0.50
GPSE GCN
63.16:0..8 0.85
0.24877:00.00
GPSE GPS
66.88=1.51
0.24644:0.00m
MOSE
63 63.5act:l11
0.318=00010
Ours
HOD-GNN
60.68100.50 55
0.24802000011
Summary. Across all experiments, we find consistent evidence supporting the guiding questions outlined above. (A1) Across the ZINC, OGB, and Peptides datasets, HOD-GNN is the only architecture that consistently ranks within the top two liers. (A2), On the challenging Peptides datasets, HOD-GNN scales to larger graphs that full-bag Subgraph GNNs cannot handle, while maintaining strong predictive performance. (A.) On several synthetic experiments, 1-HOD-GNN matches the expressivity of Subgraph GNNs while surpassing encoding-augmenttd MPNNs, showing that HOD-GNN is as expressive yet more scalable than subgraph GNNs. Moreover, 2-HOD-GNN exhibits even stronger expressive power, demonstrating the benefits of increasing derivative order. Together, these findings establish HOD-GNN as a scalable, expressive, and broadly effective GNN framework across synthetic, molecular, and large-scale real-world tasks
6 CONCLUSION
We introduce HOD-GNN, a GNN that enhances the expressivity of base MPNN by leveraging its high-order derivatives. We provide a theoretical analysis of HOD-GNN's expressive power, connecting it to the k-OSAN framework and RWSE encodings, and show that. it can offer better scalability than comparable expressive models. Empirically, HOD-GNN has strong performance across a range of benchmarks, matching or surpassing encoding based methods, Subgraph GNNs, and other common baselines.
Limitations and Future Work. Although HOD-GNN has favorable theoretical complexity (Section g its practical efficiency depends on sparse matrix operations, and more optimized implementations could significantly improve scalability. Additionally, the connection between MPNN derivatives and oversquashing or oversmoothing, alongside HOD-GNN's strong perfor- mance with deep base MPNNs and small hidden dimensions, suggests a deeper link not fully addressed in this work. Exploring this connection is promising direction for future research.
10

---

Published as a conference paper at ICLR 2026
ACKNOWLEDGEMENTS
HM is supported by the Israel Science Foundation through a personal grant (ISF 264 264/23) and an equipment grant (ISF 532/23), and by the Career Advancement Chairs in Artificial Intelligence Schmidt Futures. ME acknowledges support from the Isracli Ministry of Innovation, Science and Technology. YG is supported by the UKRI Engineering and Physical Sciences Research Council (EPSRC) CDT in Autonomous and Intelligent Machines and Systems (grant reference EP (P/SS22050/1) F.F. conducted this work supported by Aly Kaufman Post-Doctoral Fellowship. G.B. is supported by the Jacobs Qualcomm PhD Fellowship.
REFERENCES
Ralph Abboud, Ismail Ilkan Ceylan, Martin Grohe, and Thomas Lukasiewicz. The sur- prising power of graph neural networks with random node initialization. orXiv preprint 011::200000119,, 2020.
Singh Akansha. Over-squashing in graph neural networks: A comprehensive survey, arXiv preprint arXiv:2308. 15t.:888881156,, 2023.
Marjan Albooych, Daniele Bertolini, and Siamak Ravanbakhsh. Incidence networks for geometric deep learning. arXie preprint ar arxiv:1905. 11460, 2019.
Uri Alon and Eran Yahav. On the bottleneck of graph neural networks and its practical implications. arXiv preprint ArXiv:2006. 05205, 2020.
Álvaro Arroyo, Alessio Gravina, Benjamin Gutteridge, Federico Barbero, Claudio Gallicchio, Xiaowen Dong, Michael Bronstein, and Pierre Vandergheynst. On vanishing gradients, over-smoothing, and over-squahhing in gnns: Bridging recurrent and graph learning. arXiv preprint orxiv:2002. 10818, 2025.
Vikraman Arvind, Frank Fuhlbrück, Johannes Köbler, and Oleg Verbitsky. On weisfeiler- leman invariance: Subgraph counts and related graph properties. Journal of Computer and System Sciences, 113:42 59, 2020.
Waiss Azizian and Marc Lelarge. Exprcssive power of invariant and equivariant graph neural networks. arXie preprint arXiv:2006.15646, 2020.
Federico Baldassarre and Hossein Azizpour. Explainability techniques for graph convolutional networks. arXi preprint arXiv:1905.13666 2019.
Limus Bao, Emily Jin, Michael Bronstein, ismail ilkan Ceylan, and Matthias Lanzinger. Homomorphism counts structural encodings for graph learning. arXiv preprint arXiv:2410. 186766 2024.
Guy Bar-Shalom, Beatrice Bevilacqua, and Haggai Maron. Subgraphormer: Subgraph GNNs moet graph transformers. In NewrIPS 2023 Workshop: New Frontiers in Groph Learning, 2023. URL https:// wwt.://wwwtt.oeeeeoeomee fd
Guy Bar-Shalom, Beatrice Bevilacqua, and Haggai Maron. Subgraphormer: Unifying subgraph GNNs and graph transformers via graph products. In Forty-first International net/forumiia
Conference on Machine Learning, 2024a. URL
(VVJDW
Guy ChalSma Yam Eitan, Fabrizio Frasca, and Haggai Maron. A flexible, equivariant framework for subgraph gnns via graph products and graph coarsening, arXiv preprint or Xix::22406 09291, 2024b.
Beatrice Bevilacqua, Fabrizio Frasca, Derek Lim, Balasubramaniam Srinivasan, Chen Cai, Gopinath Balamurugan, Michael M Bronstein, and Haggai Maron. Equivariant subgraph aggregation networks. arXiv preprind ar 2xxm:2110.02111 2021.
11

---

Published as a conference paper at ICLR 2026
Beatrice Bevilacqua, Moshe Eliasof, Eli Mcirom, Bruno Ribeiro, and Haggai Maron. Efficient subgraph gnns by learning effective selection policies. International Conference on Learning Representations, 2024.
Cristian Bodnar, Fabrizio Frasca, Nina Otter, Yuguang Wang, Pictro Lio, Guido F Montufar, and Michael Bronstein. Weisfeiler and lehman go cellular: Cw networks. Advances in Newral Informalion Processing Systems, 34:2625 2-22222-640, 2021.
Giorgos Bouritsas, Fabrizio Frasca, Stefanos Zafeiriou, and Michael M Bronstein. Improving graph neural network expressivity via subgraph isomorphism counting. IEEE Transactions on Pallern Analysis and Machine Intelligence, 45(())6557-6) -66657666666666 2022.
Xavier Bresson and Thomas Laurent. Residual gated graph convnets. arXiv preprint ar Xix: 7111.11111111111 1.005755, 2017.
Semih Cantürk, Renming Liu, Olivier Lapointe-Gagné, Vincent Létourncau, Guy Wolf, Dominique Beaini, and Ladislav Rampášck. Graph positional and structural encoder. In Forty-first International Conference on Machine Learning. 2024. URL https:// openreview neti/eewweeeeeeeeeeeeeeeeeeeeeee tttts/forummiddddue d
Quentin Cappart, Didier Chétclat, Elias Khalil, Andrea Lodi, Christopher Morris, and Petar Velićković. Combinatorial optimization and reasoning with graph neural networks. Journal of Machine Learning Research, 24(130):1 61, 2023.
Martin Carrasco, Caio Netto, Vahan Martirosyan, Aneeqa Mchrab, Ehimare Okoyomon, and Caterina Graziani. Rademacher meets colors: More expressivity, but at what cost? In New Perspectives in Graph Machine Learning, 2025.
Zhengdao Chen, Lei Chen, Soledad Villar, and Joan Bruna. Can graph neural networks count substructures? Advances in neural information processing systems, 33:10383 10395, 2020.
Gabriele Corso, Luca Cavalleri, Dominique Beaini, Pietro Liò, and Petar Veličković. Principal neighbourhood aggregation for graph nets. Advances in Neural Information Processing Systems, 33:13260 13271, 2020.
Leonardo Cotta, Christopher Morris, and Bruno Ribeiro. Reconstruction for powerful graph representations, In Advances in Neural Information Processing Systems, volume 34, 2021.
Francesco Di Giovanni, Lorenzo Giusti, Federico Barbero, Giulia Luise, Pietro Lio, and Michael M Bronstein. On over-squashing in message passing neural networks: The impact of width, depth, and topology. In International Conference on Machine Learning, pp. 7865 7885. PMLR, 2023a.
Francesco Di Giovanni, T Konstantin Rusch, Michael Bronstein, Andreea Deac, Mare Lackenby, Siddhartha Mishra, and Petar Veličković. How does over-squashing affect the power of gnns? Transactions Machine Learning Reseorch, 2023b.
Vijay Prakash Dwivedi, Anh Tuan Luu, Thomas Laurent, Yoshua Bengio, and Xavier Bresson. Graph neural networks with learnable structural and positional representations. arXiv 10. 07875, 2021.
preprint arXiv:1 10.
Vijay Prakash Dwivedi, Ladislav Rampášek, Michael Galkin, Ali Parviz, Guy Wolf, Anh Tuan Luu, and Dominique Beaini. Long range graph benchmark. Advonces in Neural Information Processing Systems, 35:22326 15:22326-22340, 22340, 2022.
Vijay Prakash Dwivedi, Chaitanya K Joshi, Anh Tuan Luu, Thomas Laurent, Yoshua Bengio, and Xavier Bresson. Benchmarking graph neural networks. Journal of Machine Learning Research, 24(43):1-48) 2023.
Yam Eitan, Yoav Gelberg, Guy Bar-Shalom, Fabrizio Frasca, Michael Bronstein, and Haggai Maron. Topological blind spots: Understanding and extending topological deep learning through the lens of expressivity. arXiv preprinl xxXiv:2408.054 05486, 2024.
12

---

Published as a conference paper at ICLR 2026
Moshe Eliasof, Fabrizio Frasca, Beatrice Bevilacqua, Eran Treister, Gal Chechik, and Haggai Maron. Graph positional encoding via random feature propagation. In International Conference on Machine Learning, pp. 9202 9223. PMLR, 2023.
Matthias Fey and Jan Eric Lenssen. Fast graph representation learning with pytorch geometric. arXiv preprint arx: 1et.19990..228,, 2019.
Billy J Franks, Christopher Morris, Ameya Velingker, and Floris Geerts, Weisfeiler-leman at the margin: When more expressivity matters. arXwv preprint urXiv:2422.0055 2024.
Fabrizio Frasca, Beatrice Bevilacqua, Michael Bronstein, and Haggai Maron. Understand- ing and extending subgraph gnns by rethinking their symmetries. Advances in Neural Information Processing Systems, 01313766-1390 2022.
Floris Geerts. The expressive power of kth-order invariant graph networks. arXiv preprint arXtv:2007. 12035, 2020.
Yoav Gelberg, Yam Eitan, Aviv Navon, Aviy Shamsian, Theo Putterman, and Haggai Maron. Gradmetanet: An cquivariant architecture for lcarning gradients. In Workshop Newral Network Weights as a New Data Modality, 2025.
Justin Gilmer, Samuel S Schoenholz, Patrick Riley, Oriol Vinyals, and George E Dahl. Neural message passing for quantum chemistry. In International conference on machine learning, pp. 1263-1272. PMLR, 2017.
Michael Hardy. Combinatorics of partial derivatives. arXiv preprint moth /0t//060114, 2006.
Xiaoxin He, Bryan Hooi, Thomas Laurent, Adamn Perold, Yann LeCun, and Xavier Bresson. A generalization of vit /mlp-mixerr to graphs. In International conference machine learning, pp. 12724-12745. PMLR, 2023.
Nicholas J Higham. Accuracy and stability of numerical algorithms. SIAM, 2002.
Shengding Hu, Zheng Xiong, Meng Qu, Xingdi Yuan, Marc-Alexandre Côté, Zhiyuan Liu, and Jian Tang. Graph policy network for transferable active learning graphs. Advances in Neural Information Processing Systems, 33:10174-10151 2020a.
Weihua Hu, Matthias Fey, Marinka Zitnik, Yuxiao Dong, Hongyu Ren, Bowen Liu, Michele Catasta, and Jure Leskovec. Open graph benchmark: Datasets for machine learning graphs. Advances in neural information processing systems, 33:22118 122113, 2020b.
Yinan Huang, Xingang Peng, Jianzhu Ma, and Muhan Zhang, Boosting the cycle counting power of graph neural networks with i--gnns In The Eleventh International Conference on Learning Representations, 2022.
John J Irwin, Teague Sterling, Michael M Mysinger, Eric S Bolstad, and Ryan Coleman. Zine-a free database of commercially available compounds for virtual screening. Journal of Chemical Information and Modeling, 52(7):1157-1768; 2012.
Stefanie Jegelka. Theory of graph neural networks: Representation and learning. In The International Congress of Mathemattiiians pp. 1 1-23, 2022.
Mitchell Keren Taraday, Almog David, and Chaim Baskin. Sequential signal mixing ag- gregation for message passing graph neural networks. Advances in Neural Information Processing Systems, 37:93985-990 94021, 2024.
Thomas N Kipf and Max Welling, Semi-supervised classification with graph convolutional networks. International Conference on Learning Representations, 2016.
Pan Li and Jure Leskovec. The expressive power graph neural networks. Graph Neural Networks: Foundations, Frontiers, and Applications, pp. 63-98, 2022.
Derek Lim, Joshua Robinson, Lingxiao Zhao, Tess Smidt, Suvrit Sra, Haggai Maron, and Stefanie Jegelka. Sign and basis invariant networks for spectral graph representation learning. Xiv preprint. arXio:2202.1 19013, 2022.
13

---

Published as a conference paper at ICLR 2026
Ilya Loshchilov and Frank Hutter. Decoupled weight decay regularization. arXiv preprint arXiv: 1ttt:1711.005 .51001, 2017.
Dongsheng Luo, Wei Cheng, Dongkuan Xu, Wenchao Yu, Bo Zong, Haileng Chen, and Xiang Zhang. Parameterized explainer for graph neural network. Advances in neural information processing systems, 93:19620-1990 2020.
Haggai Maron, Heli Ben-Hamu, Nadav Shamir, and Yaron Lipman. Invariant and equivariant graph networks. arXiv preprint arXi:: 09902, 2018.
Haggai Maron, Heli Ben-Hamu, Hadar Serviansky, and Yaron Lipman. Provably powerful graph networks. Advances in neural information processing systems, 32, 2019.
Sohir Maskey, Raffaele Paolino, Fabian Jogl, Gitta Kutyniok, and Johannes F Lutzeyer. Graph representational learning: When does more expressivity hurt generalization? arXiv preprint arx:v:2505. 11055555.119999 11298, 2025.
Eric Mitchell, Charles Lin, Antoine Bossclut, Chelsea Finn, and Christopher D Manning. Fast model editing at sealc. arXiv preprint arXiv:2110.11 110009 2021.
Fedorico Monti, Fabrizio Frasca, Davide Eynard, Damon Mannion, and Michael M Bronstein. Fake news detection on social media using geometric deep learning. arXiv preprint ar xrv: 06073, 2019.
Christopher Morris, Martin Ritzert, Matthias Fey, William L Hamilton, Jan Eric Lenssen, Gaurav Rattan, and Martin Grohe Weisfeiler and leman go neural: Higher-order graph neural networks. In Proceedings of the A AAAT conference on artificial intelligence, volume 33, pp. 4602-4604 2019.
Christopher
Morris, Yaron Lipman, Haggai Maron, Bastian Rieck, Nils M Kriege, Martin Grohe, Matthias Fey, and Karsten Borgwardt. Weisfeiler and leman go machine learning: The story so far. arXi preprint arXiu:21 12.09992, 2021.
Christopher Morris, Yaron Lipman, Haggai Maron, Bastian Rieck, Nils M Kriege, Martin Grohe, Matthias Fey, and Karsten Borgwardt. Weisfeiler and leman go machine learning: The story so far. The Journal of Machine Learning Research, 24(1):1580 15923, 2023.
Pál András Papp and Roger Wattenhofer. A theoretical comparison of graph neural network extensions. In International Conference on Machine Learning, pp. 17823-17345. PMLR, 2022.
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al. Pytorch: An imperative style, high-performance deep learning library. Advances in neural information processing systems, 32, 2019.
Paolo Pellizzoni, Till Hendrik Schulz, Dexiong Chen, and Karsten Borgwardt. On the expressivity and sample complexity of node-individualized graph neural networks. Advances in Neural Information Processing Systems, 17:120221-120255 2024.
Phillip E Pope, Soheil Kolouri, Mohammad Rostami, Charles Martin, and Heiko Hoffmann. Explainability methods for graph convolutional neural networks. In Proccedings of the IEEE/ CVF conference on computer vision and paltern recognition, pp. 10772 10781, 2019.
Omri Puny, Derek Lim, Bobak Kiani, Haggai Maron, and Yaron Lipman. Equivariant polynomials for graph neural networks. In International Conference on Machine Learning, pp. 28191 28222. PMLR, 2023.
Chendi Qian, Gaurav Rattan, Floris Geerts, Mathias Niepert, and Christopher Morris. Ordered subgraph aggregation networks. Advances in Neural Information Processing Systems, 35:21030-2100 2022.
14

---

Published as a conference paper at ICLR 2026
Ladislav Rampášck, Michael Galkin, Vijay Prakash Dwivedi, Anh Tuan Luu, Guy Wolf, and Dominique Beaini. Recipe for a general, powerful, scalable graph transformer. Advances in Neural Information Processing Systems, 35:14501-14515, 2022.
Ronald C Read and Robin J Wilson. An allas of gruphs. Oxford University Press, 1998.
Bastian Rieck, Christian Bock, and Karsten Borgwardt, A persistent weisfeiler-lehman procedure for graph classification. In International Conference on Machine Learning, pp. 5448-5458. PMLR, 2019.
Ryoma Sato. A survey on the expressive power of graph neural networks. arXiv preprint arXiv:2003.0400 2020.
Ryoma Sato, Makoto Yamada, and Hisashi Kashima. Random features strengthen graph neural networks. In Proceedings of the 2021 SIAM international conference on data mining (SDM), pp. 333 341. SIAM, 2021.
Nimrod Segol and Yaron Lipman. On universal equivariant set networks. ar preprint arXi::1910.022. 02421, 2019.
Joshua Southern, Yam Fitan, Guy Bar-Shalom, Michael M. Bronstein, Haggai Maron, and Fabrizio Frasca. Balancing efficiency and expressiveness: Subgraph GNNs with walk--aased centrality, 2025. URL attpss//opeerevate
Jan Tönshoff, Martin Ritzert, Eran Rosenbluth, and Martin Grohe, Where did the gap go? reassessing the long-rang graph benchmark. arXiv preprint Xri::22309. 00367, 2023.
Jake Topping, Francesco Di Giovanni, Benjamin Paul Chamberlain, Xiaowen Dong, and Michael M. Bronstein. Understanding over-squashing and bottlenecks on graphs via curvature. arXiv preprint or arXw:2111.14422 2021.
Xiyuan Wang and Muhan Zhang. Using random noise equivariantly to boost graph neural networks universally. arXie preprint arXiv:2502. 0244,, 2025.
Yanbo Wang and Muhan Zhang. An empirical study of realized gnn expressiveness. In Forty-first International Conference on Machine Learning, 2024.
Boris Weisfeiler and Andrei Leman. The reduction of graph canonical form and the algebra which appears therein. nli, Series, 2(9):12-16, 1968.
Felix Wong, Erica Zheng, Jacqueline A Valeri, Nina M Donghia, Melis N Anahtar, Satotaka Omori, Alicia Li, Andres Cubillos-Ruiz, Aarti Krishnan, Wengong Jin, et al. Discovery of structural class of antibiotics with explainable deep learning. Nature, 12677997) 2024.
Dejia Xu, Peihao Wang, Yifan Jiang, Zhiwen Fan, and Zhangyang Wang. Signal processing for implicit neural representations. Advances in Neural Information Processing Systems, 35:13404-13411 2022.
Keyulu Xu, Weihua Hu, Jure Leskovec, and Stefanie Jegelka. How powerful are graph neural networks? International Conference Learning Representations, 2018.
Zuoyu Yan, Junru Zhou, Liangeai Gao, Zhi Tang, and Muhan Zhang. An efficient subgraph gnn with provable substructure counting power. In Proceedings of the 30th ACM SIGKDD Conference on Knowledgee Discovery and Data Mining, pp. 3702-3713, 2024.
Zhitao Ying, Dylan Bourgeois, Jiaxuan You, Marinka Zitnik, and Jure Leskovec. Gu- nexplainer: Generating explanations for graph neural networks. Advances in neural information processing systems, 32, 2019.
Jiaxuan You, Jonathan M Gomes-Selman, Rex Ying, and Jure Leskovec. Identity-aware graph neural networks. In Proceedings of the AAAI conference on artificial intelligence, volume 35, pp. 10737-10745 2021.
15

---

Published as a conference paper at ICLR 2026
Manzil Zaheer, Satwik Kottur, Siamak Ravanbakhsh, Barnabas Poczos, Russ R Salakhutdinov,, and Alexander Smola. Deep sets. Advances in neural information processing systems, 30, 2017.
Bingxu Zhang, Changiun Fan, Shixuan Liu, Kuihua Huang, Xiang Zhao, Jincai Huang, and Zhong Liu. The expressive power of graph neural networks: A survey. IEEE Transactions on Knowledge and Data Engineering, 2024a.
Bohang Zhang, Guhao Feng, Yiheng Du, Di Hc, and Liwei Wang. A complete expressivo- ness hierarchy for subgraph gnns via subgraph weisfeiler-lehman tests. In International Conference on Machine Learning, pp. 41019 41077. PMLR, 2023a.
Bohang Zhang, Shengjic Luo, Liwei Wang, and Di Hc. Rethinking the expressive power of gnns via graph biconnectivity. arXiv preprint or Xiw:2301.09500, 2023b.
Bohang Zhang, Jingchu Gai, Yiheng Du, Qiwei Ye, Di He, and Liwei Wang. Beyond weisfeiler- lehman: A quantitative framework for gum expressiveness. arXiv preprint 240100001,
2024b.
Muhan Zhang and Pan Li. Nested graph neural networks. In Advances Neural Information Processing Systems, volume 34, 2021.
Lingxiao Zhao, Wei Jin, Leman Akoglu, and Neil Shah. From stars to subgraphs: Uplifting any GNN with local structure awareness. In International Conference on Learning Represenlations, 2022.
16

---

Published as a conference paper at ICLR 2026
A PREVIOUS WORK
Expressive Power and Hierarchies in GNNs. The expressive power of GNNs is oflen
measured by their ability to distinguish non-isomorphic graphs. Foundational results (Morris et al. 2019 Xu et al. 2018) show that standard MPNNs are bounded by the 1- I-Weisteiler- Lehman [L) test [Weisfeiler & Leman 1968), motivating the development of more expressive architectures, see (Sato) 2020 Morris et al., 2021 Jegelka 2022 Li & Leskovec 2022 Zhang et al. 2024a) for comprchensive surveys. For instance, Morris et al. (2019) and
Maron al. (2018) introduced GNN hicrarchics matching the expressivity of the k-WL Lest at computational cost of Onn*) in both time and memory. Other approaches include
equivariant polynomial-bused models (Maron el. al. 2019 Puny et al. 2023 Subgraph
GNNs (Zhang & Li 2021: Cotta et al.
2021 Bevilacqua et al. 2021 Frasca et al. 2022
Zhang et al. 2023b a Sar-Shalom et al.
2024b topologically enhanced GNNs (Rieck et al.
2019 Bodnar et al., 2021 Fitan et al.
2024h and more. A complementary line of work improves expressivity by cnriching node features with informative structural descriptors,
such substructure and homomorphism counts (Bouritsas et al.
2022
Bao et al. 2024
random node features (Abboud et al. 2020 Sato et al. 202 D Caaso et al. 2023 or spectral
methods (Dwivedi al. 2023 Lim al. 2022
Derivatives of GNNs. Derivatives frequently appear in the analysis of GNNs. A prominent
example is the study of oversquashing-the failure of information propagate through graph structures (Alon & Yahav 2020 Topping et al. 2021 Di Giovanni et al. 2023a bj. A central tool in analyzing oversquashing the use of derivatives: specifically, the gradients of final node representations with respect to initial features (sce Di Giovanni ct al. 2023a and
mixed output derivatives with respect to pairs of input nodes (see Di Giovanni et al. (2023b). For a comprehensive overview of oversquashing, play a key role in GNN explainability (Ying et al. 2019 Luo el al. 2020 Baldassarre & Azizpour 2019 Pope et al. (011)). Gradient based approaches such as Sensitivity Analysis,
Akansha (2023). Node derivatives also
Guided Backpropagation [Baldassarre & Azizpour 2019), and Grad-CAM (Pope et al. 2019)
rely derivative magnitudes to compute importance scores. Finally, several standalone works make use of node-based derivatives. For instance, AArooo et al. (2025) use node derivatives to draw a connection between vanishing gradients, and over-smoothing, In a different direction, Keren Taraday et al. (2024) propose new aggregation functions for MPNNs, designed specifically induce non-zero mixed node derivatives.
Learning over derivative input. Beyond GNNs, several recent works have explored learning directly from derivative- based inputs. Xu
(2022) propose framework that processes spatial derivatives of implicit neural representations (INRs) to modify them without explicit decoding. [Mitchell al. (2021) introduce a learned method for fact editing in LLMs using their gradients. Gelberg et al. (2025) present general architecture for learning over sets of gradients, with applications in meta-optimization, domain adaptation, and curvature estination.
B MOTIVATION
Beyond being widely used and informative quantity GNN analysis, MPNN derivatives enhance expressivity. We now provide intuition for why this occurs by drawing connection to Subgraph GNNs.
Consider a DS-GNN model V composed of base MPNN M followed by a downstream MPNN T. We assume for simplicity that M outputs graph-level scalars, and that the activation function used in M is analytic with infinite convergence radiu That is, for every x E IR, a satisfies:
(a)
(
=
(14)
a!
Many commonly used functions, including sin, cos, and exp, are analytic with infinite convergence radius. See Appendix E.2 for discussion of the case where is not analytic.
17

---

Published as a conference paper at ICLR 2026
For each node of given imput graph 9, we define function fo :RRRA R by:
fu(z) M(A, @ a e,),
(15)
which corresponds to the ouptput of M obtained by scaling the node marking feature at node v by z, Observe that fo is analyti since M is composed of analytic functions, and that by definition, fo(1) - h i.c., the representation of the graph augmented with mark for the node v. Also observe that by expanding fa around x = O, we can approximate h to any desired precision using a finite number of derivatives, up to order m.
since
m
(()) (0)
f0(1)
~
(16)
a!
a=0
Morcover, each derivative f (0) corresponds to a partial derivative of the output of M with respect to the v-th coordinate of the augmented input:
(0) MM
hhh
(17)
00
3ne
This suggests that by using an encoder network to extract node features from the first m derivatives of M with respect each node feature, we can effectively reconstruct hiub, By passing these derivative-baset features to the downstream GNN 'T, we can approximate the behavior of V, and therefore be at least as expressive.
C HOD-GNN VARIANTS
C.1 K-HOD-GNN
In this section, elaborate on the k-HOD-GNN architecture, discussed in Section 3.3 Similar to 1-HOD-GNN, k-HOD-GNN model, denoted &, consists a base MPNN M, two derivative encoders Unode Uout and a downstream network T. Given an input graph g= (A, X); the computation of (G) proceeds in four stages: (1) compute the output and final node representation the input graph using M; (2) compute the k-indexed derivative Lensors D (T) D (3) extract new derivative informed node features form the derivative tensor; (4) Use these new node features for downstream processing.
Steps 1 & 2. In the first two stages, we compute the final node representations h²)) and the output vector hou using the base MPNN M, along with their corresponding derivative lensors defined below (the k-indexed derivative tensor was defined in Section 3.3 for he case where the input features are 1-dimensional).
Definition C.1. Given a graph = (A, Xo) with n nodes, Xo E Raxd an MPNN M and some intermediate node feature representation matrix h E Raxd' the k-indexed derivative tensor of D (h)) E IRaxn"xe xd' dx is defined by:
glal
D (h)(,,, u, h (Xo)
(18)
0x
9in X=X,
uidi
where E V(g), (ui,..., uk) E V*(g). i [d'), a (ag f0....,m. m 1 1]axk and Ja|- = E
Similarly, given a geaph prediction vector hout E TRdd the derivative tensor D (hout) ga xd' x xm is defined by:
D2 (a(hluu, af 8thh h;
(19)
The derivative tensors D % passing- like procedure detailed in Appendix |
((h(T)) ) and D D(h are computed using a message- which enables both efficient computation and allows us to backpropagate through it, supporting end-to-end training of k-HOD-GNN.
For definition of multii-dimensionl analytic functions Appendix F.1
18

---

Published as a conference paper at ICLR 2026
Step 3. In the third stage of our method, extract new node features informed by the derivative tensors D ( D& using the encoder networks Unode R.ax"" xd' xmxxxx dx + R
and Uout 1= xd' xm ggad"" As ( (T) is a tensor indexed by (k 1) nodes and u tensor indexed by k nodes, natural choices for Uuode and C are (k + 1)-IGN and k-IGN (Maron et al. 2018) designed specifically to process such data.
Furthermore, Proposition 1.4 shows that for sparse graphs or relatively shallow base MPNNs, the derivative tensor D (?) itself becomes sparse, with space complexity O (n min d), where d is the maximum degree of the input graph. This motivates the choice of node encoder Uhnde that preserves this sparsity structure and exploits it for improved computational
efficiency.
To enable this, we can define the encoder Uuode as a subclass of the general (k+1)-IGN architecture, implemented as a DeepSel (Zaheer et al. 2017) operating independently on the derivative entries associated with each node feature h. Specifically, we define:
Uds-node D(T)) = DeepSet((( [v, u. t(u,
u(u, u, a)) u (9),iE E [m]] (m]
(20)
Here we assume nodes are given in index form, that is U E (I,...,n) u E fl,..., n)t, and t(c, a) is function that encodes the derivative pattern associated with the index tuple (v, u. a); that is t(01, U1 01) f(u2, U2 (22)
a1 Q2 and Eo Sn such that e ((2)), in o (uz))
(21)
where Sn denotes the symmetric group on n elements. In other words, maps each derivative index tuple to a comonical identifier that is invariant under permutations of the node indices but sensitive to the derivative multi-index.
This design improves both space and time complexity, as the sets over which the DecpSet operates are typically sparse. The proof of Theorem 4.1 provided in appendix E.2 shows that this encoder architecture is enough to be as expressive as K-OSAN.
Finally, we proceed the same way 1-HOD-GNN constructing the derivative-informed node features deer by combining information from the base MPNN, the pooled intermediate derivatives, and the output derivatives:
h der h
+ Uout (Datt) @ Uuode (T)
(22)
Step 4. This step is identical to that of I-HOD-GNN described in Section 3.2
C.2 EDGE-HOD-GNN
In this section, we extend the 1-HOD-GNN formulation to incorporate edge-feature derivatives rather than node-feature derivatives. We refer to this variant as edge-HOD-GGNI The construction of this architecture closely parallels that of Section 3.2 and is detailed below.
Steps 1 & 2. Similarly to steps 1 & 2 in Section 3.2 in the first two stages, compute the final node representations hT) and the output vector nut using the base MPNN M. The key difference is that now compute the edge derivative tensors, defined below:
Definition C.2. Given graph G E) with n nodes, I edges, an edge feature matrix E, an MPNN M and an intermediate node feature representation matrix E Ruxdd the edgo-derivative tensor D(h) E gaxixd xdxm is defined by:
hh Dedge (h)le, e,,,,, OE,,"
(23)
where E V(G), E(G), e [d'),j [d] specify the feature dimensions of the feature vectors h, E. respectively, and a E [m] where m E N is a hyperparameter specifying the maximum order of derivatives to be considered. Similarly, given graph-level prediction
19

---

Published as a conference paper at ICLR 2026
vector nou E Rd the derivative tensor (ho l xdxm is defined by:
Dedge (hou) j, h
(24)
d
Step 3. Similarly to step 3 in Section 3.2 in the third stage of our method, we extract new edge features from the edge derivative tensors Dddge and ( odge using encoder networks edge and Uadge First, Dedge gxxd' xdxm is tensor indexed by single edge, it can be directly interpreted as an edge feature matrix by flattening the remaining dimensions. We thus define the encoder network U to be a DeepSets (Zaheer et al. 2017 update:
Uedge (Deuge edgeee MLP Deut le.ele.....)).
(25)
Secondly, since ( edge is indexed by node-edge pairs, it shares the structure of the incidence matrix of 9 denoted by B. We can thus define the encoder network Usage Baxiix xdxm + Rixd to be any GNN architect which maps incidence matrices
to edge feature matrices. An example of such architecture is Albooyeh et al. 2019))
Step In the final stage, we replace the original edge features of with the derivative- informed features Eder, and apply a downstream GNN T to produce graph-level prediction.
D DERIVATIVE COMPUTATION
We now extend the derivative computation algorithm presenled in Section 3.2. to account for k-mixed derivatives as well as more general MPNN architectures. Similarly Section 3.2.1 we first split the node update procedure of an MPNN into two parts: an aggregation step:
( = (((-1) e AGG
(
(26)
and DeepSet step:
00 MLP((
(27)
Our algorithm begins by computing the initial derivative tensor (0) (we abuse notation
and ommit the subscript k:), and then recursively constructs cach D (t-1) (t-h(t))) from
(t-1) and then (t)) from D (t-1) Finally, the output derivative tensor Dout is obtained from (T), See Algorithm for the full procedure.
To analyze the time and memory complexity of each step, we define the derivative sparsity of a node of an input graph G at layer t of the base MPNN, as the number of non-zero derivatives corresponding lo h That is
Se, {(u, i, ax)
D [ u, aJ) # +07|-
(28)
min Se.t vEV(g)
(29)
When is small, the tensor (t) can be stored efficiently in memory using sparse representa- tions, requiring only O(n st) space. The quantities 8,,, and s are leveraged in Appendix E.2
to derive concrete asymptotic bounds for the complexity of Algorithm I
Computing (0), Since (0) X, the derivatives straightforward to compute:
1 if 3s s.t. v u. 0g,, = 1, Σ. Ms.jE(d] avs = 0
(0) b. u, a
(30)
0 otherwise.
Although D (0) is high dimensiaan,,nnl,, it is extremely sparse with so = O(1).
Computing D (t-1) from ((-1), Abusing notation, and assuming that for every node
in the input graph, both h (t-1) and h (t-l).aak = (GG(t) (Ih!-1 :u ((E)()))) lie in Rdd,
20

---

Published as a conference paper at ICLR 2026
we observe that since h 0,,0''-11 h (t-1) the derivalives of the first. d' coordinates of h are precisely D (t-1), In contrast, the derivatives with respect to the last d' coordinates depend heavily on the choice of aggregation function AGG.
However, when the aggregation function is linear, i.e., of the form AGG ((h!- ((())))
Σ v'E((")
b (t-1)
(t-1)
(31)
for some coefficients bet that depend only the adjacency matrix A of the imput graph (this is the case in most widely used MPNNs), the computation of derivatives simplifies. Since differentiation commutes with linear operations, the computation in Equation 31 carries over, yielding
D(h (t 1),85 |u,...]==
Σ vEE((1)
be (,
(32)
Ag Aggregating the derivatives of neighboring nodes mirrors the structure of message passing, which endows it with several beneficial properties. First, using Equation 32 the Lensor
D ((-1) [v....] can be computed in time O(d ((-) where d is the maximal degree of the input graph[] for each node, we aggregate neighbor derivative vectors, each containing at most 8t-1 non-zero entries. This means that. the total computation time of this step is Od n 8t-))
The argument above also implics that the total number of non-zero elements in (f- 1) bounded by (min(d- 8t-1, n²)). Thus, the above update leverages the sparsity of the graph to achieve elficiency in both space and time complexity.
[v....]
Computing D (t) from D (t-1)
We begin by assuming that the MLP in Equation 27 has depth 1, i.e.,
MLP (x) o (w b())).
(33)
In the general case where the MLP has depth 1, the procedure described below is applied recursively I times.
We define the intermediate linear activation as
h (t-1),Lin ( h000 + b())
(34)
and describe how compute ((h from (t-1) followed by the computation of (t) from (1=(tt)).
First, since the update in Equation 34 is affine, we can drop the bias term and get:
D D(h (t-1)in (o, u. i, = w(t) D (tt"1)o, [.. u, a,
(35)
Second, notice that
h = o ( We can use the Faá di Bruno's formula (see e.g. (Hardy 2006) which generalizes the chain rule higer-order derivatives, for our last step. Faà di Bruno's 's formula states that for a pair of functions g RR + R, f R + R, y ( 2n,) the following holds, regardless of whether the variables Iı, In are all distinct, identical, or grouped into distinguishable categories of indistinguishable variables:
(t- 1).Lin
(36)
0x ... Ox f(v)
8
Σ (((n)()) zEnl
II BEz
B
y Oz,
(37)
Jee
where:
We slightly abuse notation by using d to denote both the input feature dimension and other quantities; the meaning should be clear from context.
21

---

Published as a conference paper at ICLR 2026
n)) The notation B E indicates that B is one of the subsets (or "blocks") in the partition
II denotes the collection of all partitions of the index set fl....,
For any set A, the notation |A| represents its cardinality. In particular, | is the number of blocks in the partition, and [B] the number of elements in the block B.
Equations 36 and 37 Imply that are able to compute ( u....] based only ((h [1, u, u.... and the derivatives of at the point he ( Combining this
update with Equation 35 results way to compute ( from D ((--1)
Importantly, the update above computes each entry of [ u....] using only the cor- responding entries of D (t,
[1, u.... i.e. those associated with the same node tuple (v, u). Moreover, from Equations 35 and 37 it follows that if " [v, u....] = 0, then ( [e, u u..] = as well. This implies that like D (t-1) the number non-zero entries in each D ( [u.... is also bounded by O [min(d- st-1, n²)). Consequently, the update performed only over the non-zero entries of D (t-1) has runtime complexity of O (n. min(d 81-1, n²)).
Computing Dout from D (T))
Recall that
out AGG (( DE v
(38)
where AGGaa denotes the final aggregation over node embeddings. This operation can be treated analogously to the node update step: For most common MPNNs, it decomposcs into a linear aggregation followed by an MLP. Consequently, the derivative tensor D out can be computed from ( using the same primitives described above.
Algorithm 1 Efficiont Computation of Derivative Tensors
Require: Graph G = (A, X), base GIN M with T layers
1: (0) < X
D node feature init. deriv. init through Eq 30
2: D (0) D(X)
3: for t = 1 to T do
h (t-1)
+- h (t-1)
+ EEVEN(O) bor boi
4:
b lincar node agg.
5:
D (t-1)
lo... .....)(((t-))) (v, 1 1 D(..---
deriv. agg.
6:
h
MLP (t- (h (1-1))
D DeepSet update deriv. DeepSet update. extract output deriv.
( 1= gegggee- gel- (c-1) [v, [v..... h
J.....,h!!! (
7:
8: end for
Dout
get-out-der (2))
9: return out
E EXTENDED THEORETICAL ANALYSIS
E.1 DEFINITIONS
Before delving into the proofs, we begin by formally defining several key concepts used throughout the analysis:
Definition E.1 (Analytic Function). A function f R + is said to be analytic at 2c E R if for some R R holds that for all [X| R:
f(x) = Σ
ga f(xo)
(39)
a!!
|a|=0
22

---

Published as a conference paper at ICLR 2026
where E No and we use the following notation:
a 01 + 02 + an:
a! an! ! or!. an!!
(a) = (x1))1 )ar (xz)
(..))
... (In)
((a) ōx, 1 Dz2 ...
pal
x=a
The largest such R. is called the radius of convergence. A function f: R² + R is analytic all functions fi...., fm are analytic.
Definition E.2 (k-OSAN). k-OSAN model v consists of a base MPNN M that produces updated node features (as opposed to directly outputting a graph-level prediction). followed by a downstream MPNN T that aggregates these features to produce a final graph-level output. Given an input graph G = (A, the output (G) is computed in four stages:
Step 1: Construct bag of subgraphs Bg = (S.. u E (²() each of the form S, = (A, 3 e"). Here e" guxk is "node marking" 5 feature matrix assigning unique identifier to each node u1,... ug. That is:
[ 0 else.
Uj
(40)
Step 2:
Compute the (k: + 1)-node indexed tensor: Hlu, u] M(Su)
(41)
Step 3:
set aggregation function to produce new node features:
h = AGG((H(u, ul u E (*))).
(42)
Step 4: Compute the final output through:
(g) -((A,h²ub)
(43)
For k = 1, A-OSAN are also reffered to as DS-GNNs.
Definition E.3 (Random Walk Structural Encoding). For a graphh = (A, the Random Walk Structural Encoding (RWSE) with number of steps defined as
h rwse A diag(AA))
(44)
where = A Diag(deg(ur)-1 deg(u, (1) is the row-normalized adjacency matrix, diag(-) denotes the vector of diagonal entries of matrix, and Diag(v) denotes the diagonal matrix with vector v its diagonal.
E.2 PROOFS
Theorem 4.1 We begin by formally stating and proving Theorem 4.1 splitting it into two Thoorey pairs.
Theorem E.4 (N-HOD-GNN is as expressive k-OSAN). Let (G' (A', i [/} be fonite set of graphs, and Y be a k-OSAN model. for any > O, there erists k-HOD-GNN model & such that for each i E [!]
(VG))) ((G)) <E.
(45)
Alhhugh alternative methods for initializing node features in subgraphs have been proposed, they offer the same expressive power. We therefore focus on the simple approach used here.
23

---

Published as a conference paper at ICLR 2026
Corollary E.5. There erist non-isomorphicc graphs that are indistinguishable by the folklore k- WL test but are distinguishable by k-HOD-GNN, Additionally, k-HOD-GNN able to compute the homomorphism count of k-apez forest graphs.
Proof. We begin the proof by making a few simplifying assumptions on V, which we can do without loss of generality. We begin by assuming that all activation functions used in M, the base MPNN of V, are analytic with infinite radius of convergence (e.g., exp(x), sin(r)). This assumption can be made without loss of generality: Since MLPss with non- pon-polynomial analytic activations are universal approximators, each MLP in M can be replaced with one using analytic activation function that approximates the original to arbitrary precision. Furthermore, since the composition of analytic functions with infinite convergence radius remains analytic with infinite convergence radius, it follows that M-as composition of affine transformations and activation functions is itself analytic.
Secondly, we assume that the final node representations produced by M are one-dimensional. This assumption can be made without loss of generality: we can append a final MLP to M that compresses each node's feature vector downstream network T that reconstructs the original feature dimension. This effectively amounts to inserting an autocncoder, where the encoder is a final pointwise update in M and the decoder is an initial pointwise update in T. Since this architecture can approximate the original architecture to arbitrary precision, we may assume without loss of generality that M produces 1-dimensiioal node embeddings.
scalar, and prepend an MLP to the
Additionally, we can assume the input node feature matrices X [1] are also all dimensional. This follows from the same argument as above.
Finally, consider k-HOD-GNN models that only use the node derivative tensor ( which we use to extract node features through an IGN encoder U, and disregard Dout
We prove the theorem in three steps
1. We show that an intermediate representation of (T) can encode the tensor Ha,. that is produced at stage (2) of the forward pass of v (see Definition E.2
2. we show that U(D(T)) can approximate the node feature matrix huu produced al stage (3) of the forward pass of V.
3. We show o can approximate T to finite precision.
Step 1 For each k-tuple of nodes u E V*()") and matrix Y E R**k, we define the node feature matrix broad, (Y) E Bax (here we abuse notation and not include the graph index i) by
{o U5.r v ==
broadu (Y)
(46)
else.
Additionally, for each node v E G', we define fú IR xxk : R
(Y) M(A, X broad,, (Y))-
(47)
Finally, define Yu R by
Y
e "aij'
(48)
where " is the node marking node feature matrix introduced in Definition E.2 (That is,
broad, (Y e e")
First, it is easy to that
(Y") = M(S)). H' u).
(49)
where St. and H'(v, u/ are introduced Definition E.2
9In cases where u has repeated entries uj Miz only consider values of Y for which
V. Ya.t
24

---

Published as a conference paper at ICLR 2026
Second, M is analytic with infinite convergence radius, the functions fū are all analytic with infinite convergence radi, and so
fe (Y")=
Σ 0a (() |a|=0
(50)
Here, a E {O...., 1] and (Y")a= = II 22 Since we are concerned with a finite number of graphs, for any <> 0 we can choose an integer 7 such that for all graphs g and all e E V(g'), u E (g') it holds that:
fi.j2
1 (Y")
Σ
80
(0))())
<6.
(51)
a
Defining
H[
Σ |a|-0
a!!
0a
Se (O)
(52)
we get that
H'
(53)
Additionally, from the definition fi it holds that the derivatives of fi at zero correspond to the entries of the k-indexed derivative tensor of g denoted by D (7) that is
0 fe. (0) D (T) a.,
(54)
where we only take derivatives with respect the last feature dimensions, which correspond to the "marking vectors"
Moreover, let & be function such that l(v, a:) encodes the derivative pattern associated with the index tuple (v, u, a)). That is, I satisfies:
((11, u1, 0n) (U2,, U22 02)
ai 02; and ; 3o Sn such that iı ( 4 (us))
(55)
where Sn is the symmetric group on n elements. The value of ai (Y")a is determined entirely by MV, u, ax), and thus can be recovered from it.
This implics that the DeepSet encoder Uds-node defined in Equation 20 in Appendix C.1 have an intermediate layer L such that
can
L(D A'.
(56)
Here L simply multiples each entry D (T) 0 u, ai) by 1 (Y")a, followed by summing over the ar indices. Thus, by Equations 53 and 56 choosing M as the base MPNN of k-HOD-GNN model allows us to approximate each H to arbitrary precision using its derivatives. Note that since the DeepSet encoder is restricted instance of a (k + 1)-IGN encoder, it can achieve the same effect.
Step Equation 42 shows that the node features hanh, are constructed by applying set-wise aggregation function over the set (H* [u. u] u E V*(g))) Any continuous set-wise aggregation function can be approximated to arbitrary precision by DeepSet architecture (Zaheer et al., 2017) (see Segol Lipman (2019) for proof). Moreover, any DeepSet model applied parallel over the first index of k-indexed node tensor can be exactly implemented by a k-IGN, since each layer in such a model consists a linear equivariant transformation followed by a pointwise nonlinearity. Thus, we construct our encoder
U to first approximate the mapping D - A' as an intermediate representation, and then
here we abuse notation and omit the subscrit in D
25

---

Published as a conference paper at ICLR 2026
approximate the subsequent mapping H - hauba, Since both approximations can be made to arbitrary precision, this completes the proof of the claim in Step 2.
Step The final forward step applying a downstream GNN to the updated node features to produce graph-lovel representation- identical in both OSAN and HOD-GNN. Therefore, by choosing the downstream GNN the k-HOD-GNN model to match that of V the proof is complete.
Using Theorem we now prove corollary
Proof of corollary E.5 It was shown by Qian et al. 2022 that k-OSAN models can dis- tinguish between graphs that are indistinguishable by the k- WL test. As Theorem E.4 establishes that method can approximate any k-OSAN model to arbitrary precision, it follows that k-HOD-GNN can do the same. Similarly, Zhang et al. (2024b) showed that k-OSAN models compute homomorphism counts of k-apex lorests graphs in which the removal of al most & nodes yields a forest. Therefore, by Theorem E.4 Corollary E.5 follows.
Theorem E.6 ((k + 2)-IGNs are expressive as k-HOD--NNN) Let d be H-HHODGNN model and let G,g' g' a pair of graphs such that
4 # (G').
(57)
There exists a (k (k+22 2)) model such that: V()) (G')).
(58)
Corollary E.7. k-HOD-GNN is unable to distinguish any pair of (k + 1)-FWL indistin- guishable graphs.
Proof. First, recall that the k-HOD GNN 1OD-GNA is composed of a base MPNN M, downstream MPNN T, (k + I)- encoder Unode and a k-IGN encoder Uout all of which are less expressive than the (k + 2)-IGN architecture. Thus, is enough to show that (k + 2)-IGN is able to simulate the ellicient derivative algorithm presented in Appendix D to compute the k-order derivative tensors (T) and D (out) We now show how (k: 2)-IGN is able to simulate each step of this algorithm.
Before we begin, recall that (k + 2)-IGN operates on tensors T E R2=+22d where d is called "the feature dimension". As the derivative tenors are of the form D I xd xm we slightly change the notation for the tensors T which (k + 2)-IGN operates on to T e R=* 1x allowing multiple feature dimensions. We stress this is only a notation convenience, as we can transform T back to single feature dimension simply by "flattening" the different feature dimensions. Thus, for the rest of the proof, similarly to definition we assume T' is indexed by T(u, i, a) where u = (U1...., Ux+2) VV)) (arj)sadd ( 1) dx and lo| Σ aa
[d]
Throughout the proof, tensor representations corresponding to G, G' will be denoted by T, T respectively.
Step 1: Computing (()) For an input graph (A, Recall that:
Since h (0) X, we have:
[I if 3s s.t. == 0s, 1,,:=1,E. sa ass =
D (0) (0, =
(59)
otherwise.
The initial tensor (()) used by (k: + 2I- is such that Tlu,... T 1 if and only the map i, + u, is a graph isomorphism on the subgraphs of g induced by o and u respectively.
26

---

Published as a conference paper at ICLR 2026
In addition, one of the core operations of (A + 2)-GNNs allows it to apply "pointwise" linear
layer followed by an activation on any entry .... of the tensor simultancously, That is, using (k + 2)-IGN layers Maron et al. (019); Frasca et al. (2022) for more details).
update the tensor T via MLP(( + T(u.... (See
Define
S1 = (T [u. lu.... ... u2) U (T [u....]Iu =22
(60)
52 (T (u. (l,...|/u E (a( U (T (u,...]luu .uu1 No (u2))
(61)
S3 (T( [u, .u. & (Na(u2) U (u2))) U (T Here (-) is the neighborhood of a node.
(0) [u, [u,...]lu & (Wer (uz) U {u2))}
(62)
From the definition of the initial tensor ( above, S1, S3 all have pairwise emply intersections. Thus, we can define an MLP such that Vn E S1
MLP MLP(=))i,a)= 1 if 3s s.t. a.i = 1, 1, Zajeed asd = 0
(63)
otherwise.
Var S2:
MLP(=)[i,a) =1
(64)
and Var Sz:
ML,(x)[t,al
(65)
by updating Tlu.... (LP(T [u, T f'u,.... MLP(T [u....]) we now get:
T(u .... Uk+1, Mk+l.
( U2, U4+2....
(66)
and for 11 # u2
T(u,, . 2..2,...) 12.....,H42....1 Aur
(67)
Thus T now stores both the information of the tensor D (0) and the adjacency A.
Step 2: Computing (t-1) from (t-1).
For simplicity assume the base MPNN used in k-HOD-GNN is GIN architecture (the proof can be casily generalized to the general case)
Recall from Appendix D that
[u....] D (...
Σ D (u',.) e'EN(a)
(68)
This amounts to constructing a "flattened" " node feature vector D E x defined by
Danr [M] = (t=") J.l..ttn))
(69)
and performing standard message passing on it.
We can follow similar path, defining "flattened" matrix T E 1 defined by Taat [U1, Ug] = T(u1, u2. J.fatten)))
(70)
and then use (A: + 2I-IGN layer to perform 2-IGN updates on THat- From equation the tensor retains the adjacency information of the input graph. Following arguments
67
27

---

Published as a conference paper at ICLR 2026
presented in Maron et al. 2019 2-IGN layer simulate message passing. Thus,
able compute
completing the step.
Step 3: Computing D from ((-11) As shown in Appendix D D (c) can be computed
from D (t-1) by a point-wise update. That is, there exists continuous function f which depends the choice of the activation and weights of the base MPNN M such that
D Juiia f(D (1 [u, al).
(71)
We can thus choose an MLP which approximates f to finite precision and update our tensors T, T' according to T[u,...] MLP T' [..] MLP((T(u,...)) This finishes
the current step.
Iteratively updating the tensors and T' according to steps 1-3, we reach final tensors such that
T [22, U2 Uk42 = D (2)) [u,,..., 2....
(72)
[u2, U2 u2.....
uk+2 .... D D(22) [u2. ......... 1k..22.....
(73)
(Ic + 2)-IGN can then apply projection P. 1=" *+2 xxe k xxd,
defined by
((T) U2, MK++2----. - u2, ...., uk+2. ...2....1
.I= = (T [12,-.. ... uk42..... (74)
Thus (k: + 2)-IGN can recover (T))
Step 4: Computing Dout from (T)) For most common MPNNs, this step decomposes
into a linear aggregation (which can be handled exactly like step 2) followed by an MLP (which is equivalent to the update of step 3). Consequently, the derivative tensor Dout can be computed from D (T) using the same primitives described above.
Proof of corollary
E.7 Recall that the (k + 1)-FWL has equivalent expressive power the
(k + 2b-oblivious WL test (see Morris et al. (2023) ). In addition, as shown in Geerts (2020):
Azizian & Lelarge (2020); Maron et al. 2019) (T 2I-IGNs have the same expressive power
as the [K + 2) oblivious WL test. This together with Theorem E.6 completes the proof.
Theorem 4.3 We now formally state and prove Theorem
4.3
Theorem E.8 [HOD-GNN is strictly more expressive than RWSE HOWSE+MPPN MPNN). For any MPNN
T augmented with random walk structural encodings (see Definition E.3). there exists 1-HOD- GNN model & that uses ReLU activations and only first-order derivatives such that, for every graph G, holds that
T(g) = (g)).
(75)
Moreover, there crist pair of graphs G and g2 such that for every (WWSE-augmentee MPNN
T(g = T(g?),
(76)
yet there exists a 1-HOD-GNN model 1, using ReLU activations and only first-order deriva- lives, such. that.
((G') + ((g2).
(77)
Proof. To prove the first part the theorem for RSWE with L mumber of steps, we begin with simple preprocessing step. For each input graph G = (A, X) with node feature matrix
28

---

Published as a conference paper at ICLR 2026
QU15
QU19
Figure 2: Two quartic vertex-transitive graphs that cannot be distinguished by MPNNs augmented with RWSE, but can be distinguished using HOD-GNN.
X E Raxd we define an extended feature matrix X + 1, E (xx(d+L) by padding X with constant vecotr of lenght L and value 1. We then pass X to the HOD-GNN model instead of the original
It is now sulficient lo construct base MPNN M with T layers such that:
1. The first d coordinates of each final node embedding satisfy (T) e,0:d-1 Xus, i.e, the original features are preserved.
2. The first-order derivatives of the remaining L coordinates (i.e., for indices I = d, d.....d d satisfy:
0h (I) v
Á
(78)
ox,,1
where À is the row-normalized adjacency matrix (see Definition E.3
If the above conditions hold, choose MLP such that
h MLP (T) U (D()) e (Dou)) X, $ h rwse
(79)
and choose our downstream model to be exactly T thus satisfying Equation 76
We now construct MPNN M that satisfies both of the conditions above.
For each layer i t= 0,...., . T 1, we define the update rule of M to act separately on two parts of the node feature vector: the first d + coordinates and the remaining T coordinates. Specifically:
h (I+1)
(1) (()
(80)
"v,0d4t-1 1000::+t-11
h (H+1)
sdd+t--11 ReLU
Σ
(t) a
(81)
deg(v) aEA((o)
29

---

Published as a conference paper at ICLR 2026
We set the number of layers in M to be T L, where L is the number of random walk steps used in the original RWSE encoding.
First, the proposed update rule is straightforward to implement within the MPNN framework
defined by Equation !
as follows a standard message- passiing structure.
Second, Equation 80 guarantees that the first d coordinates of cach node's feature vector
remain unchanged throughout the layers, thereby satisfying condition (1).
Finally, Equations
%0 and
81 together imply that for all t O....,T T - 1, it holds that
(3) = A' 1r.
(82)
This implies that
0h 0h u,f A An
(T)
(83)
ax 7,7
and so condition (2) holds, completing the first part of the proof.
For the second part of the proof, following the nolation of (Read & Wilson 1998), let 91 and G2 denote the quartic vertex-transitive graphs Q115 and QUT9, respectively (see Figure ( with all imitial node features equal to 12 R2, Here, trantix-tanssitvee" means that for any pair of nodes, there exists a graph automorphism mapping one to the other, and "quartic" indicates that the graphs are 4-regular. This pair of graphs was shown in (Southern et al. 2025) to be indistinguishable by MPNNs augmented with RWSE. Therefore, to conclude the proof, it suffices to construct an HOD-GNN model that can distinguish between them.
We define the base MPNN M to consist of two layers, specified as follows. The first layer performs standard neighbor aggregation:
h ReLU Σ
(0)
(84)
(EN(())
In the second layer, the first coordinate each node feature is preserved, while the second coordinate is updated via another aggregation step. That is
h
(2) v.0
= -
(85)
ReLU
Σ
(86)
r.1
MEN(
from the same argument the first. part of the proof, gel that for M, it holds that
0h A u,u: 0x n.
(87)
Note that for pair of nodes u, v, we have Au, = 1 and A 0 if and only if E (v) and, for all u' e WVv), it holds that u & (N(')) that is, there is no path of length exactly 2 from u to v. A straightforward computation shows that no such node pairs exist in g!. whereas several such pairs appear in g². Thus, the set of off-diagonal derivative vectors {D (2).2 [u, v 1,:] u A (())) contains values that do not appear in the corresponding set {D (2)1 [ 1.:] u u V)
30

---

Published as a conference paper at ICLR 2026
We define the node derivative encoder Unode 2-IGN model that operates in two steps. First, it constructs filtered tensor D by zeroing out all entries of D that are not off- diagonal or are farther than c = 4 from the vector (1,0). Then, performs row-wise summation over D to produce node features:
Unode (D (2)
D u].
(88)
(EV(
It follows that Unode (D (2),1 while Unode (D (D (2) + O. Therefore, when these vectors are passed to the downstream MPNN, it will be able to distinguish between the two graphs, completing the proof.
Computational complexity.
We now prove Propositions 4.4 and 4.5 which analyze the time and space complexity of HOD-GNN.
Proof- Recall that in Appendix D we have defined the derivative sparsity of a node o of an input graph G at layer of the base MPNN, as
Se, { itu.i,a) D D u. i, ai + +0}l
(89)
and the maximal derivative sparsity layer is
min 80,4-
(90)
EVV(G)
the tensor (() has memory complexity of O(n 8t) as it can be stored using sparse matrix representations. Algorithm for elficient derivative computation iteratively computes ( using (t-1) and finally computes D Dut from (T)) We now prove by induction that for an input graph G with maximal degree maximal number of derivatives, the input node fcature dimension, and the dimension of ((), all of which are hyperparamcters, such that
there exists a constant which depends on the
< C min (n*)
(91)
We additionally show that D be constructed from (t-1) in time complexity of
Odd n min minfdl (t-1), nn))
To begin, as saw in Appendix D
if Es s.t. v= v= Jos 2s, 0g = 1, I, ag =0
D [v. u. i,
(92)
otherwise,
and so Se,0 < C = 1k+1 where m is the maximal derivative degree, and | is the dimension of the imput node features. Now assuming Equation 91 holds for 1, in Appendix D we saw that (t)) can be computed from (t-1) in time complexity of O(d $1-1) = (dd n (t-(dd)((t1), ,n*)), secondly, we saw that 8; = O (min(d st-1, n²]) O (min n²)) completing the induction. As the memory complexity of storing (()) is O(n st) O (n. m n**) this proves Proposition 4.4 Recall now that have seen in the proof of thcorem E.4 that k-HOD-GNN able to achieve the same expressivity as k-OSAN by disregarding
out and choosing the encoder Umode to be
uds-node (D (D DeeeSett(DD (D u. al u E (9), i [d'], jd] a E [m]*)).
(93)
This encoder can leverage the sparsity of ( having runtime complexity of O(n ST) =
O(n min(n db-T)), completing the proof of proposition 4.5
We slightly abuse notation by using d to denote the maximal degree as it usually denotes the input feature dimension, we denote this quantity I for this discussion.
31

---

Published as a conference paper at ICLR 2026
F ADDITIONAL PROPOSITIONS
In this section, we present three additional propositions which further explore the exact. expressivity gains obtained from using derivative signals.
The first proposition shows that even when using only output-lovel derivatives, 1-HOD-GNN is still fairly expressive.
Proposition F.1. Even when selling Umode 0, 1-HOD-GNN still as expressive as the DS- WL for node based policies Lest. defined in Bevilacqua el al. 2021).
Proof. We begin by recalling the definition of the DS-WL test. Consider two graphs 91 = X1) and G2 = (A2, X2). For each graph, the DS-WL test constructs a bag of node-marked subgraphs, where each subgraph is obtained by appending unique mark to the feature of exactly one node. Formally, for i E (1,2},
B (
V(G.))
(94)
where each marked subgraph is defined as
S (A,, X, e,),
(95)
with e, denoting a one-hot indicator applied to node v.
The DS-WL test then applies the WT. refinement procedure independently to each marked subgraph S until convergence, and assigns it a final color
cs WL(S).
(96)
This produces, for each graph, a multiset of resulting WL colors:
C ccs | B
i (1,2).
(97)
Finally, DS-WT declares 91 and G2 to be non-isomorphic if and only if C1 & Cz. Tf the two multisets coincide, the test does not distinguish the graphs.
As shown in Morris et al. (2019), for any finite collection of graphs there exists an MPNN M whose output equivalent the WL coloring. In particular,
M(S) = M(S') O cs',
(98)
for any marked subgraphs S and S'.
Consequently, Gi and G2 are separated by the DS-WL test if and only if there exists an MPNN such that
(M(S) S E B:} # (M(S) S B2).
(99)
As previously discussed, we may assume without loss generality that the activation function used by M is analytic. Indeed, if this is not the case, the MLP components of M can be approximated to arbitrary precision by MLPs employing an analytic activation function.
For every node V(G),, define the scalar function
fi(x) M(A,, X, Oxe,)
(100)
where x E R.
We note the following properties:
1. f((1)= = M(S,,)) by definition of the marked subgraph S.,
2. For every j E N,
= a M(
Ozj f((()
(101)
ax!
==
since the perturbation only modifies the feat ure al node v.
3. Each f is analytic, as M is a composition of linear transformations and analytic activation functions.
32

---

Published as a conference paper at ICLR 2026
Since analytic functions uniquely determined by their derivatives at point, obtain that for any pair of nodes u E V(1) and u E V(G2),
0 0x
0 M(S,,u for all E EEN) 0x.
o(
(102)
Therefore, if M(S # M(S2,u) analyticity implies that there exists some integer such that
nT +
8m
M(S
M(Sz,u),
(103)
0
0xm
and consequently
Dout (S).)) + Dout
(104)
Hence,
(M(S) | S E B.} + {M(() S B²)
(105)
can hold only if
(S) S E Bi} # { (S) S E B2).
(106)
Finally, by choosing the encoder - to be a fully expressive multiset function (e.g., DeepSets Zaheer et al. (2017)), 1-HOD-GNN is able to distinguish 91 from G2 using only the derivative tensor Doutt This completes the proof.
The following proposition analyses the expressive power of edge-HOD-GNN.
Proposition F.2. edge-HOO Odge-HOO-GNA GNN is eapressive as the DS- WL for edge based policies lest defined in Bevilacqua et al. (022)) and can thus seporate pairs of 2-FWL indistinguishable graphs.
Proof. The proof of this proposition follows the same argument as the previous one; we include it here for completeness. For brevity, for the rest of this prove we refer to the odge-marking WL test simply as the WL test, and the edgo-derivative tensor simply as the derivative tensor, denoted by Dout We begin by recalling the definition of the DS- WL test with edge marking. Consider two graphs Gi = (Ai, X1, Ei) and G2 (A2, X2, E2) where E,, Eg represent edge feature matrices. For each graph, the DS-WL test constructs a bag of edge-marked subgraphs, where each subgraph is obtained by appending unique mark to the feature of exactly one edge. Formally, for i {1,2),
B, (S,,e le E(G.))
(107)
where each marked subgraph is defined as
Sec (A,, X,, E, R @e.),
(108)
with e, denoting a one-hot indicator applied to edge e.
The DS test then applies the WL refinement procedure independently to each marked subgraph S.,, until convergence, and assigns it a final color
cs = WL(S).
(109)
This produces, for each graph, a multi-sct of resulting WL colors:
C {cs E ()) i (1,2).
(110)
Finally, DS-WI declares 91 and G2 to be non-isomorphic if and only C₁ * Cz. If the two multi-sets coincide, the test does not distinguish the graphs.
As shown in Morris et al. (2019), for any finite collection of graphs there exists an MPNN M whose output equivalent to the WL coloring. In particular,
M(S) MM(S') + cs CS
(111)
for any marked subgraphs S and S'.
33

---

Published as a conference paper at ICLR 2026
Consequently, 91 and G2 are separated by the DS-WL test if and only if there exists MPNN such that
(M(S) S E & (M(S) B2).
(112)
As previously discussed, we may assume without loss of generality that the activation function used by M analytic. Indeed, if this not the case, the MLP components of M can be approximated to arbitrary precision by MLPs employing an analytic activation function.
For every edge e E(G,), define the scalar function
fx( M(A, x, E, 0x.
(113)
where E R.
We note the following properties:
1. fd(1) M(S..e) by definition of the marked subgraph St.e-
2. For every j E N,
dx fo (x) = O M(S..c)
(114)
dE
since the perturbation only modifies the feature at edge e.
3. Each f is analytic, as M is a composition of linear transformations and analytic activation functions.
Since analytic functions are uniquely determined by their derivatives at point, we obtain that for any pair edges e, E E(g)) i [2]. [2].
Or M(S2,c2) for all j jEN] = M(S M(S2.r OE
M(S cl
(115)
OE
Therefore, if M(SL,,)) + M(S,,e.), analyticity implics that there exists some integer m such that.
am M(S1.e) * 3m
dmm M(S,.e)).
(116)
dmm
and consequently
Dou # (S2.ez)))
(117)
Hence,
( (M(S) S SEB. B.} + (M(S) 8 E B2)
(118)
can hold only if
(S) S S Bi} (S) S B2).
(119)
Finally, by choosing the encoder U to be fully expressive multiset function (e.g., DeepSets Zaheer el al. (2017)), edge-HOD-GNN is able to distinguish G1 from G2 using only the edge derivative tensor Dout the edge marking WI. test was shown in Bevilacqua et al. 2021 to be able to separate 2-FWL indistinguishable graphs. This completes the proof.
Finally, the next proposition shows that for a fixed base MPNN M, increasing the hyperpa- rameter m (Definition 3.1), which specifics the highest derivative order used by HOD-GNN, strictly improves expressivity.
Proposition F.3. For any m E N, there erist choices of the base MPNN M such that a H.HD-GGN using M with saper-parameter m is unable count the number of triangles of given input graph, while a &-HOD-GNN using M with m + can.
Proof. Let E N. We define the base messago-pasee passing network Mm as follows.
Consider an input graph g = XA,X), where X E R is a node- feature matrix with scalar features.
34

---

Published as a conference paper at ICLR 2026
Initialization. We define the first hidden representation by node-wise MLP update (which ignores all neighbors):
X m++
(1)
(120)
m 1)!
This update depends only on X, and is therefore realizable by a local node transformation.
Message passing. For E (2 4}, we define
10 = Ah (-11)
i.e.,
h (4)
Σ
(1-1)
(121)
(EN(1)
Each update is realizable by standard MPNN aggregation step without applying a post- aggregation MLP.
Readout. The global output Mm is defined by summing the final node representations:
h out
Σ
(122)
DEV
Thus, Mm is a message-passing neural network consisting of one node-wise MLP layer followed by three pure aggregation layers, and final sum readout.
As a consequence the chain rule, all partial derivatives both (4 and h of order at most vanish. Therefore, 1-HOD-GNN model restricted to the hyperparameter has expressive power equivalent to a standard MPNN, and in particular cannot count. triangles.
Morcover, have
m (4) h
(A)),
(123)
0x, m+ 0.0
Hence, a 1-HOD-GNN model with hyperparameter 1 can exploit this derivative to compute
(A3) N.v
Σ
(124)
which equals the number of triangles in the input graph G. This completes the proof.
G EXPERIMENTAL DETAILS
In this soction, we provide details on the experimental validation doscribed and discussed in Section
G.1 ARCHITECTURE
In all our experiments, we have used a 1-HOD-GNN architecture, We now describe its components in detail.
Base MPNN. all experiments, we use a GIN Xu et al. 2018) architecture as the base OPPN-described in Equation g (Section 3 due its simplicity and maximal expressivity. All MLPs used in the node updates are two layers deep and employ ReLU activations. We extract first-order derivatives from the base MPNN, which by Theorem 4.3 yields greater expressivity than MPNNs augmented with RWSE. After computing the final node representations h (1), we apply a residual connection by aggregating all intermediate layers:
T
h0
(T)
(125)
t
1=1
where normalizing in !! helped stabalize training.
Derivative Encoding.
35

---

Published as a conference paper at ICLR 2026
For simplicity, we ignore the output derivatives of the base MPNN and use only the final node-wise derivative tensor ( The encoder Uhode is implemented lightweight, efficient module that applies pointwise MLP to the diagonal entries of D(T;
noddo (())) = MLP (u.
(126)
Downstrcam GNN.
For simplicity and to isolate the effects of HOD-GNN, we restrict our experiments to MPNNs as downstream architectures We use GIN for all experiments but peptieds fune in which we used a GCN (Kipf & Welling, 2016 All MLPs used in the node updates are two layers deep with ReLU activations. The final MT.P head also uses ReL. activations and consists of 1 to 3 layers, following the default settings from [Southenn et al. (2025), without tuning. In experiments with parameter budgets, we adjust only the hidden dimension to within the limit, selecting the largest value that satisfies the constraint.
Initialization. We initialize the base MPNN such that all calegorical feature embeddings are set to the constant vector 1, all MLP weights are initialized to the identity, and the parameters in the GIN update (Equation 9 are set to -1. This initialization is motivated by the proof of Theorem 4.3 following a similar line of reasoning. it implies that the diagonal of
the derivative tensor
[e, U...] corresponds exactly to the centrality encoding proposed
in Southern et al. (0225
Optimization. We use separate learning rates for the base and downstream MPNNs. For the downstream MPNN, we adopt the learning rates used in SSuthern et al. (2025), while for the base MPNN, we fine-tune by selecting either the same rate or one-tenth of it.
G.2 EXPERIMENTS
We provide below the details of the datasets and hyperparameter configurations used in our experiments. Our method is implemented using PyTorch (Paszke et al. 2019) and Py Torch Geometric (Fey & Lenssen 2019), and is based on code provided in Southern et al. (2025) and Rampášek al. (2022). :. Test performance is evaluated at the epoch achieving the best validation score and is averaged over four runs with different random seeds. We optimize all models using AdamW (Loshchilov & Hutter 2017), with a linear learning rate warm-up followed by cosine decay. We track experiments and perform hyperparameter optimization using the Weights and Biases platform. All experiments were conducted on a single NVIDIA A100-SXM4-40GB GPU.
OGB datascts. We evaluate on three molccular property prediction benchmarks from the OGB suite (Hu et al. 2020b): MOLHIV, MOLBACE, and MOLTOX21. These datasets share a standardized node and edge featurization capturing chemophysical properties. We adopt the challenging scaffold split proposed in (Hu et al. 2020a). To prevent memory issues, we use a batch size of 128 for MOLHIV and 32 for the remaining datasets. All rewm models use a hidden dimension of 300, consistent, with prior work (Hu et al. 2020a Bevilacqua et al. 2024). We sweep over several architectural choices, including the hidden dimension of the Base MPNN k 8, 10...., 30, 32, the initial learining rate of the ..... Base MPNN k = 0.001, 0.0001, and the dropout rate k = 0.1, 0.2, 0.3, 0.4, 0.5 The number of layers in the base MPNN was selected to match the positional encoding step used in the corresponding experiment from SSuthern et al. (025)). Hyperparameter tuning was performed on the validation set using four random seeds. Results are reported at the test epoch corresponding to the best. validation performance. All models were trained for 100 epochs. The final parameters used for each experiment are reported table
Zine The ZINC dataset (Dwivedi et al. 2023) includes 12k molccular graphs of commercially available chemical compounds, with the task of predicting molecular solubility. We follow the predefined dataset splits and report the Mean Absolute Error (MAE) both the loss and evaluation metric. Our downstream MPNN for this task includes 6 message-passing layers and 3 readout layers, with hidden size of 120 and no dropout. We use batch size of 32 and train for 2000 epochs. We performed a small sweep over the depth of the base MNPNN k 10, 12, 14, 16, 18, 20 and the hidden dimension k = 30, 35. ... 80. The hidden dimension
36

---

Published as a conference paper at ICLR 2026
of the downstream MPNN was chose as 120 to meet the 500k parameter constraint. The final hyperparameters are listed in Table
Peptides Peptides-fune and Peptides-struce introduced by Dwivedi et al 2022 consist of graphs representing atomic peptides. Peptides-fune is a multi-label classification benchmark with 10 nonexclusive peptide function labels, while Peptides-struct is regression task involving 11 different structural attributes derived from 3D conformations.
For both datasets, adopt the hyperparameter setup proposed by Tönshoff et al. (2023 for the rewnstream GNN, which has a parameter budget under 500k and where they use 250 epochs. We set the number of message-passing layers in our base MPNN with the positional encoding steps to be 20, aligned with the number of steps used for the random maak walk structural encoding. The only tuned component is the learning configurations are summarized in Table 5
of the base MPNN. The final
Key empirical findings. Across all benchmarks, HOD-GNN is highly competitive and the only architecture that consistently ranks within the top two model tiers. Additionally, the strong performance of HOD-GNN on the large-scale Peptides datasets where full-bag Subgraph GNNs are generally unable to run shows ils ability to scale effectively. Notably, the base MPNNs in HOD-GNN are often significantly deeper and narrower than those in typical GNNs. For instance, the OGB datasets, the base MPNNs use 17-20 layers with hidden dimensions as low as 16 32. Despite their compact size, these base MPNNs yield notable performance gains standard GINE, suggesting potential robustness to oversquashing. Morcover, the fact that they are significantly deeper than typical MPNNs further suggests that HOD-GNN may help mitigale oversmoothing, hypothesis leave for fulure work.
Table 3: Best-performing hyperparameters for each OGB dataset.
Hyperparameter MOLHIV MOLBACE MOLTOX21
Downstream model
# Layers
2
8
10
# #Readout Layers
Hidden Dimension
300
300
300
Dropout
0.0
0.5
0.3
Learning Rate
0.0001
0.0001
0.001
Base MPNN
# #Layers
16
20
20
Hidden Dimension
16
16
Dropout
0.2
0.5
0.2
Learning Rate
0.0001
0.0001
0.0001
Parameters
450,848
1,723,448
2,165,782
Table 4: Best-performing hyperparameters for the ZINC dataset.
Hyperparameter
ZINC
Downstream model
4 Layers
6
4 Reedoutt Layers
3
Hidden Dimension
120
Dropout
0.0
Learning Rate
0.001
Base MPNN
#Layers
12
Hidden Dimension
75
Dropout
0.0
Learning Rate
0.0001
Parametees
498,144
37

---

Published as a conference paper at ICLR 2026
Table 5: Best-performing hyperparameters for the Peptides-fune and Peptides-struct datasets.
Hyperparamoter
Peptides-func
Poptides-struct
Downstream model
#Layers
10
#Readout Layers
3
Hidden Dimension
234
143
Dropout
0.1
0.2
Learning Rate
0.001
0,001
Base MPNN
#Layers
20
20
Hidden Dimension
8
Dropout
0.1
0.2
Learning Rate
0.0001
0.001
# Parameters
498,806
493,849
Table Normalized MAE results on the counting subgraphs dataset. Cells below 0.01 are highlighted in yellow.
Method
3-Cycle 4-Cycle
5-Cycle 6-Cycle
Tailed Tri.
Chordal Cycle
4-Clique
4-Path
MPNN
0.3515
0.2742
0.2088
0,1555
0.3631
0.3114
0,1645
0.1592
MPNN + RWSE
0.0645
0.0264
0.0746
0.0578
0.0505
0.1008
0.0905
0.0217
GPS RWSE
0.0185
0.0433
0.0172
0.0551
0.0446
0.0974
0.0836
0.0284
HyMN
0.0384
0.0933
0.1350
0.0936
0.0084
0.0746
0.0680
0.0120
GNN-AK+
0.0004
0.0040
0.0133
0.0238
0.0043
0,0112
0.0049
0.0075
HOD-GNN ReLU
0.0012
0.0046
0.0210
0.0380
0.0083
0.0510
0.0293
0,0081
HOD-GNN SILU
0.0008
0.0042
0.0068
0.0222
0.0066
0.0195
0.0055
0.0069
H ADDITIONAL EXPERIMENTS
H.1 SUBSTRUCTURE COUNTING
We adopt the synthetic node-level subgraph counting experiment used in Huang et al. 2002 Yan al. (2024). The dataset consists of 5,000 graphs generated from a mixture of distributions (see Zhao et al. (2022) for more details), with a lainn/vaidaatin, ietiatiioo//ttttt test split of 0.3 0.3/0.2/0. The task is node-level regression: predicting the number of substructures such as 3-cycles, 4-cycles, 5-cycles, 6-cycles, tailed triangles, chordal cycles, 4-cliques and 4-paths, where continuous outputs approximate discrete counts. We report the normalized MAE for each bascline, and highlight cases where the error falls below 0.001, since in these instances rounding the predictions yields exact counts.
For training, use the AdamW optimizer with an initial learning rate of 0.001, a cosine scheduler with warmup, and train for 5,00 epochs. The batch size set to 128. We compare HOD-GNY GNN against both positional/ positionnal/structuaal structural encoding methods (MPNN RWSE (D wiivedi
et al. 2021,, , RWSS (Rampášek et al. 2022), HyMN (Southern et al. 2025 and
node- based subgraph GNNs (GNN-AK (Zhao et al. 2022
Nested GNN (Zhang & Li
2021). ID-GNN (You et al. 2021 HyMN (Southern et al. 2025 ).
To examine the role of activation functions, evaluate HOD-GNN with the non-analytic
ReLU and the analytic SiLU. The results, summarized in Table B are consistent with
theory. In line with Theorem 4.1 HOD-GNN achieves performance comparable to surpassing other subgraph GNNs. Furthermore, consistent with Theorem 4.3 we observe that analytic activations enhance expressivity, while even with con-analytic ReLU, HOD-GNN still outperforms encoding-beas based methods.
H.2
GRAPH SEPARATION ABILITY OF - k-HOD-GNN
To further evaluate the expressive power of &-HOD-GNN and to empirically validate Theo- rem 4.1 experimented with both 1-HOD-GNN and 2-HD-GNT
the family of regular graph pairs from the BREC benchmark (Wang & Zhang 2024h. This dataset contains 140 pairs of regular graphs: 50 pairs that are distinguishable by 3-W WT. but not by 2-WT, and
38

---

Published as a conference paper at ICLR 2026
Table 7: Results on separation of pairs of regular graphs from the BREC dataset.
Model
Regular Graphs Number
Regular Graphs Accuracy
3-WL
50
35.7%
NGNN
48
34.3%
DE+NGNN
50
35.7%
DS-GNN
34.3%
DSS-GNN
48
34.3%
SUN
50
35.7%
SSWL P
50
35.7%
GNN-AK
50
35.7%
KP-GNN
106
75.7%
P-GNN
100
71.4%
OSAN
5.7%
I-HOD GNN
47
33.5%
2-HOD GNN
60.0%
90 pairs that remain indistinguishable even under 3-WL. We follow the exact training and evaluation procedures proposed in Wang & Zhang (2024).
Table 7 reports the separation performance our models alongside several subgraph-based GNNs. As shown, 2-HOD-GNN successfully distinguishes 34/90 of the 3- B-WL-indistingushhbb pairs, placing it among the top-performing models and providing strong empirical support for the theoretical advantage predicted in Theorem
Furthermore, 1-HOD-GNN obtains accuracy comparable to DS-GNN on the same benchmark, aligning with Theorem 4.1 for the k case and reinforcing that, without accessing higher- order derivatives, the model behaves similarly to standard MPNNs.
H.3 ABLATION ON DERIVATIVE ORDER
To assess the impact of higher-order derivatives on HOD-GNN's expressive power, we conduct an ablation study on the most challenging substructure counting task: 6-cycle prediction, which consistently yields the highest MAE. Specifically, we evaluate models restricted to derivatives of order at most k. for h = 0,....44 The results, presented in Table 8 demonstrate that increasing the maximal derivative order consistently improves MAE, with performance saturating at k 4. All experiments use the analytic SiLU activation function. We hypothesize that the observed saturation arises because higher-order derivatives of SiLU rapidly diminish toward zero, limiting the additional expressive gain.
Table 8: Effect of maximal derivative order on 6-Cycle MAE.
Maximal derivative order
6-Cycle MAE
0
0.1555
1
0,0275
2
0.0223
0.0221
0.0231
H.4 STABILITY ANALYSIS
We evaluate the stability of HOD-GNN both terms of training dynamics and the behavior of the derivative tensor norm. Across all tasks, we consistently observe smooth and stable loss curves, even when incorporating higher-order derivatives. To further assess stability, we conduct two dedicated experiments.
Training Loss Dynamics. We examine convergence behavior on the 6-cycle subgraph counting task by varying the maximum derivative order E 1,2, 3, 4. Figure 3 shows the training loss curves, showing consistent and reliable convergence across all derivative orders.
39

---

Published as a conference paper at ICLR 2026
MAE over Epochs for HOD-GNN with Different
HOD-GNN d=1
0.8
HOD GNN d-2
HOD-GNN d=3
0.7
HOD-GNN d-4
0.6
0.5
MAE
0.3
0.2
0.0
100
300
400
500
Epoch
Figure 3: The training loss of HOD-GNN with maximal derivative d {1,2,3, 4] on the 6-cycle counting task.
Relative Norm over Epochs
Relative (d=1) Relative Norm (d=2) 2)
0.09
0.08
Relative Norm
0.07
0.06
0.05
0.04
0.03
20
40
60
80
100
Epoch
Figure 4: The relative norm of the derivative tensor with maximal derivative d E (1,2) with respect to the norm of the final node feature matrix of the base MPNN on the MOLBACE dataset.
Norm of the Derivative Tensor. We additionally assess the stability of the derivative tensors used by HOD-GNN. We compute the norm of the derivative tensor during training on the MOLBACE dataset. Following standard practice (Higham 2002), we report the relative norm: the ratio between the derivative tensor norm and the final node feature norm of the base MPNN. Figure shows that the relative norm consistently remains significantly lower than that of the node features throughout training, confirming that HOD-GNN operates in well-conditioned regime.
H.5 ANALYSIS GENERALIZATION BEHAVIOR HOD-GNN
Recent studies Franks et al. 2024 Maskey et al. 2025 Carrasco et al. 2025 have raised concerns about the generalization ability of highly expressive GNNs, showing that performance
40

---

Published as a conference paper at ICLR 2026
Table 9: Comparison of model train-test performance gap on MOLHIV and MOLTOX21.
Metlod
MOLRIV Train MOLHIV Teol
Gup-MOLLHY
MOLTOXI
MOLTOXI
CapMOLTOOXI
GNN-GNN (BD)
01.21
76.43
2.12
15.28
92:38 1.57
15.12
17:26
DSS-GNN (ND)
89.70 5.20
2.15
75.34
16
HOD ONN
88.15 1 0.51
80.36
0.52
1.29
90.73 1: 1.00
77:50 1 0.71
16.74
can degrade when model becomes "too expressive." " These findings suggest that GNNs often gencralize best when they have just the right capacity, a balance influenced both by architectural expressivity and by paramcter scale (c.g., overly wide hidden dimensions can also harm generalization).
Motivated by these observations, we examine train-test, gaps for HOD-GNN the OGB datasets MOLHIV and MOLTOX21 (Table g). HOD-GNN exhibits moderate gaps that are smaller than those of (i) less expressive architectures such as GCN and GIN, (ii) models of comparable expressivity such as DSS-GNN (ND), and (iii) more expressive variants such as DSS-GNN (ED). These findings indicate that HOD-GNN displays strong generalization. Bascline values are taken from Bevilacqua et al. b022
We propose two hypotheses as to why HOD-GNN generalizes well, informed by recent literature:
Compact parameterization enabled by higher-order derivatives. Across our experiments, the base MPNN within HOD-GNN consistently performs best with small hidden dimensions (8-32), resulting in models with relatively few parametors, thus aligning well with the "just the right capacity" approach mentioned above. We hypothesize that higher-order derivatives enrich these compact represenlations, effectively increasing expressive power without increasing model size, thereby supporting beller generalization.
Principled initialization that recovers RWSE. The constructive proof following Theorem 4.2 gives an explicit initialization for HOD-GNN's derivative features that exactly matches Random Walk Structural Encodings (RWSE) a widely used structural prior known to improve expressivity without overfitting. This suggests an alternative view of HOD-GNN:
At the beginming of training, the derivative-informee features remain close to RWSE, maintaining well known baseline behavior.
When the task benefits from additional capacity, the model can naturally move beyond RWSE and leverage higher-order derivative information.
This yields a natural "just expressive enough" bias: HOD-GNN begins from a well- understood RWSE bascline and increases expressivity only when the data requires it.
We view further investigation of HOD-GNN's generalization properties as a promising direction for future work.
H.6 COMPARISON OF RUNTIME AND MEMORY USAGE AGAINST SUBGRAPH GNNS
In section show that HOD-GNN is able scale to the Peptides datasets, which are unreachable for full-bag subgraph GNNs on standard hardware, as stated in Southern et al. (2025); (Bar-Shllmm et al. (2023). To further demonstrate the scalability of HOD-GNN, we benchmark its runtime and memory usage against subgraph-basee GNNs on the MOLHIV datasct. Since the choice of subgraph selection policy is the primary factor determining the asymptotic complexity of subgraph GNNs (Bevilacqua et al. 2021), we evaluate range of policies using the default parameters from the original work. For subgraph GNNss we adopt the hyperparameters reported in [evilacqua et al. 2021 for MOLHIV. For HOD-GNN, we use the same hyperparameters as in Section 5 detailed in Appendix
reports GPU memory usage and per-epoch training and test runtimes. HOD-GNN achieves improvements in memory efficiency, requiring less than half the GPU memory compared to subgraph GNNs, In terms of runtime, HOD-GNN is faster than edge- deletion and node-deletion policies, while achieving comparable training time to ego-nets
Results. Table 10
41

---

Published as a conference paper at ICLR 2026
Table 10: Runtime and memory comparison on the MOLHIV dataset. HOD-GNN demon-
strates both improved memory efficiency and competitive runtime.
GNN
GPU Memory (MiB) Training Time / Epoch (s)
Test Time / Epoch (s)
odge-dclction
32,944
62.13
5.70
node-deletion
29,826
58.80
4.01
ego-nets
25,104
53.16
3.07
ego-nets-t
25,211
54.19
3.20
HOD-GNN
12,964
53.34
2.28
Table 11: MAE comparison of different GNN architectures on ZINC-12K.
GNN
MAE
GCN HOD-GNN + GCN
0.321 + 0.009
0.080 + 0.006
GIN
0.163 + 0.004 0.066 0.003
HOD-GNN GIN
GPS
0.070 + 0.004
HOD-GNN GPS
0.064 + 0.002
and ego-nets+- Notably, ego-net. policies are known to be less expressive (Bevilacqua et al. 2021h, highlighting that HOD-GNN achieves both efficiency and expressivily.
Discussion. The improvements primarily stem from the analytic computation of higher- order derivatives in HOD-GNN, which avoids the costly enumeration of subgraphs. While these advantages already translate to lower memory usage and faster runtimes in practice, emphasize that HOD-GNN's scalability potential is not yet fully realized. Tn particular, it relies on efficient sparse matrix multiplications, which are currently suboptimally implemented in popular GNN libraries such Py'oorch Geometric. We therefore anticipate that further optimization of sparse kernels would amplify the scalability benefits of HOD-GNN.
H.7 ABLATION ON BACKBONE GNN
In most experiments, employ GIN [u et al. 2018 as the backbone of HOD-GNN. A key advantage of our approach, however, is its compatibility with any message-passingg backbone. To assess the effect of backbone choice, evaluate HOD-GNN on the ZINC dataset using GCN, GIN, and GPS as base architectures. The results, summarized in Table demonstrate that HOD-GNN consistently improves upon its backbone across all settings. Among the tested architectures, GPS achieves the strongest performance, followed by GIN, with GCN ranking third reflecting its lower expressivity relative the other backbones.
ETHICS STATEMENT
This paper advances the theoretical and empirical study of graph neural networks through derivativo-basod architectures. Our experiments are conducted exclusively on widely used, publicly available benchmark datasets (ZINC, OGB, and Peptides molecular tasks), which contain no personally identifiable or sensilive information. We believe our work does not raise immediate ethical concerns. Nevertheless, representation learning may be applied to sensitive domains (e.g., biological or social network data). We encourage responsible use of our methods in accordance with the ICLR Code of Ethics.
acknowledge that improvements in graph
REPRODUCIBILITY STATEMENT
We have taken several steps to ensure reproducibility. All datasets used are standard bench- marks with clearly defined splits. A detailed description of architectures, hyperparameters, and training settings as well as full proofs of the theoretical results appear in the Appendix. Additionally, anonymized source code is provided supplementary material and will be
42

---

Published as a conference paper at ICLR 2026
released publicly upon publication. Together, these resources should allow rescarchers to fully reproduce our experimental and theoretical results.
43


Reviews:
[{"review_id": "ClSj03WfrX", "reviewer": "Reviewer_WWrV", "summary": "This paper introduces new class of Graph Neural Networks that aims to enhance expressivity by leveraging derivatives of MPNNs with respect to node features. The authors argue that derivatives naturally encode structural information and can thus extend the expressive power of base GNNs, drawing parallels with Weisfeiler–Lehman hierarchies, subgraph GNNs, and positional encodings. The paper includes theoretical expressivity proofs, an efficient message-passing algorithm for derivative computation, and experiments.", "questions": "See weaknesses.", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1) The paper is well written and easy to follow. The idea is clear, simple but not trivial, and the conceptual move of translating the use of derivatives, which are widely explored in other contexts, into the domain of expressivity is interesting. \n\n2) It is particularly interesting that k-HOD-GNN can approximate any k-OSAN model to arbitrary precision while remaining computationally more efficient. This means HOD-GNN advances the field in a desirable direction, i.e., toward architectures that are highly expressive yet computationally tractable.\n\n3) Generally, the theoretical exposition is coherent and convincing, and it is always clear what is being done and why. Overall, the contribution feels clean: the narrative is consistent, and the motivation for using derivatives is well justified.", "weaknesses": "**Minor issues:**\n\n1) There are a few small notational inconsistencies:\n\n-Line 150: $\\mathbf{X}$ should be bold, and h should carry a subscript (v for node-level, out for global as the authors define), since it’s unclear otherwise which representation is being referred to.\n\n-Line 158: when mentioning the final node representation $h_v$, it would be more consistent to add the superscript ($T$) as done earlier.\n\n-Line 125: you use $V^k$, which is understandable but never formally defined. Similarly, what exactly a “model” $\\mathcal{M}$ is could be defined more formally, though the meaning is still clear.\n\n2) The overview paragraph (lines 191–199) could be made clearer. In particular, point (4) says that $\\mathcal{T}$ is applied to “derivative-informed features”, but since $\\mathcal{T}$ is itself a GNN, the phrasing might mislead the reader to think that it acts only on those features, rather than on the full graph. Moreover, earlier in the paragraph $\\mathcal{T}$ is only called “downstream network,” which increases potential confusion. For such a crucial overview section, absolute clarity would help.\n\n**Conceptual concerns:**\n\nI think the paper could engage more critically with the current stage of the field. The literature on GNN expressivity went through a strong wave of works proposing models with expressive power going beyond 1-WL. However, in recent times the community has started to question whether such additional expressivity is truly useful in practice [1]—especially on feature-rich datasets, where node attributes already break most of the symmetries that 1-WL cannot distinguish. One could argue that even if 1-WL failures are rare in real-world datasets, developing models that go beyond the WL limit is still valuable, as it guarantees the ability to distinguish as many input as possible, regardless of whether such cases are frequent or not. However, several recent works have shown that higher expressivity often comes at the cost of worse generalization [2,3,4]. In this sense, proposing a model that only aims to surpass the WL boundary, without discussing why this added expressivity is beneficial or how it impacts generalization, feels a bit outdated relative to how the field has evolved. It would strengthen the paper to address questions like: 1) Do you have any hints about the generalization behavior of your model? 2) Why might HOD-GNN preserve generalization compared to other expressive models? 3) Can you provide an intuition for how derivative-based features may induce better inductive bias? A discussion about this would make the contribution feel more in tune with the current research conversation.\n\n\n**Experimental completeness:**\n\nIf results for certain baselines were not reported in prior papers, they were not rerun here. This leads to some missing entries, particularly for two datasets. It’s not necessarily a major issue, but it would be useful to explain the reason for the omission explicitly. For instance, were those models prohibitively expensive to run?\n\n\n[1] Jogl, Fabian, Pascal Welke, and Thomas Gärtner. \"Is Expressivity Essential for the Predictive Performance of Graph Neural Networks?.\" NeurIPS 2024 Workshop on Scientific Methods for Understanding Deep Learning. 2024.\n\n[2] Franks, Billy Joe, et al. \"Weisfeiler-Leman at the margin: When more expressivity matters.\" Forty-first International Conference on Machine Learning.\n\n[3] Maskey, Sohir, et al. \"Graph Representational Learning: When Does More Expressivity Hurt Generalization?.\" arXiv preprint arXiv:2505.11298 (2025).\n\n[4] Carrasco, Martin, et al. \"Rademacher Meets Colors: More Expressivity, but at What Cost?.\" arXiv preprint arXiv:2510.10101 (2025).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "RaOsio2LaF", "reviewer": "Reviewer_S9Dt", "summary": "The authors consider the expressive power of GNN from a derivative perspective. Combining high-order derivatives, the authors enhance the expressivity of GNN.", "questions": "1.I am interested in a more detailed characterization of the relationship between the k-WL test and GNN derivatives. For example, on a high level, can the k-WL test be equivalently expressed using derivatives of order less than k? Such a connection could directly elucidate the source of enhanced expressive power.\n\n2.I would like to clarify the impact of different choices for $\\alpha$ on expression ability in k-HOD-GNN? Theorem 4.1 only states existence, but it remains unclear: whether any possible choice of $\\alpha$ in k-HOD-GNN has stronger expression power than the k-WL test, when k and m are fixed. Then, it would be valuable to clarify the conditions on $\\alpha$ that guarantee stronger expression power. In addition, a systematic method for selecting the optimal $\\alpha$ is also an important theoretical justification.", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "In this paper, the authors consider the expressive power for the derivatives of GNN, an innovative and interesting perspective to analyze expressive ability. By leveraging to high-order node derivatives, the authors enhance expressive power of GNN. The analysis seems reasonable and the conclusion seems reliable.", "weaknesses": "While this is an interesting work, this finding is not surprising. Normally, the node feature vectors can be viewed as the zero-order information, i.e., the function value $h_v$. High-order derivatives can be seen as the high-order information for the function value. Combining much more high-order information can distinguish many graphs doesn't seem surprising to me. Instead, I am more interested in uncovering deeper insights from the derivative perspective, as such insights could provide a more fundamental understanding of the expressive power of GNNs.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "2zOBGHQUdW", "reviewer": "Reviewer_986R", "summary": "This paper introduces HOD-GNN, an architecture that enhances GNN expressivity by computing derivatives of a base MPNN with respect to input node features and using these derivatives as additional structural information to input to the second GNN.  The paper proves that k-HOD-GNN can approximate k-OSAN subgraph GNNs, which is strictly more expressive than MPNNs with random walk structural encodings, and achieves favorable computational complexity for sparse graphs with shallow base MPNNs by exploiting an efficient message-passing-like derivative computation algorithm. Empirically, HOD-GNN consistently ranks in the top two tiers across seven benchmarks.", "questions": "See weaknesses.", "limitations": "", "rating": 6, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "1. The paper presents a fresh perspective on enhancing GNN expressivity through derivatives. The intuition is clearly explained through the connection to marking-based GNNs and Taylor expansion. I think the presented approach is conceptually elegant. \n\n2. The theoretical analysis provides a complete expressivity characterization by proving k-HOD-GNN approximates k-OSAN subgraph GNNs while offering constructive separation examples that precisely identify when the method is efficient.\n\n3. The three-component pipeline from base MPNN, intermediate derivative computation, and the downstream GNN is clear and easy to follow.", "weaknesses": "1. The claimed \"inductive bias toward derivative-aware representations\" is mentioned but never explained.\n\n2. The triangle counting example in the Motivation section is thin as a driver of the whole method. This motivating toy showing first-order derivatives recover triangle counts uses a very special linear stack, i.e., identity activations, and does not show whether derivatives remain informative for typical non-linear, normalized, or regularized MPNN on real tasks. Thus, it’s illustrative but not compelling evidence that derivatives are broadly the right signal. Also, in the theoretical analysis, the separation beyond k-WL depends on analytic activations that  iss fine for softplus/tanh, but not ReLU.\n\n3. Flattening $D^{out}$ could be large in practice. \n\n4. For $U^{node}$, a 2-IGN is powerful but can be heavy. The paper claims scalability via sparsity, but constants can bite (derivative channels $\\times$ feature dims $\\times$ orders). \n\n5. I didn’t see normalization for $D^{out}$ or $D^{(T)}$. Some stabilization may be necessary because  derivatives can vary wildly with feature scaling.\n\n6. I’m confused about the rationale for computing two derivative tensors that target different functions. Specifically, $D^{(T)}$ computes derivatives of the final node embeddings, while $D^{out}$ computes that based on the graph-level output. Why is it necessary to encode both, rather than using only node-level derivatives with an invariant head, or using only graph-level derivatives if the target is graph classification?", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "DAKi1Qhfsf", "reviewer": "Reviewer_nnEj", "summary": "This paper proposes HOD-GNN, a framework that enhances GNN expressiveness by computing high-order derivatives of a base MPNN with respect to node features and using them as structural encodings. Experiments on graph classification and regression benchmarks demonstrate competitive performance, with particular success on larger graphs where full-bag subgraph GNNs fail.", "questions": "1. Can you provide experiments validating k-HOD-GNN for $k \\ge 2$? Even a simple synthetic task distinguishing graph pairs requiring 2-node marking would substantially strengthen the main contribution. Without this, the paper's scope should be revised to focus on 1-HOD-GNN.", "limitations": "", "rating": 4, "confidence": 3, "soundness": 2, "presentation": 3, "contribution": 3, "strengths": "1. The paper introduces a novel approach to GNN expressivity. Unlike prior work that uses derivatives for analysis, this paper uses derivatives as input features to enhance expressiveness. The connection between Taylor expansion of marking and derivatives provides clear intuition.\n\n2. The efficient derivative computation algorithm is a notable technical contribution. Exploiting graph sparsity through message-passing-like updates makes the approach practical, as demonstrated by favorable memory usage and successful scaling to Peptides datasets where subgraph GNNs fail.\n\n3. The experimental evaluation is thorough for k=1 case, covering diverse benchmarks, multiple baseline families.", "weaknesses": "1. Despite the paper's title and emphasis on \"k-HOD-GNN\", all experiments use k=1. The core theoretical claim that `k-HOD-GNN matches k-OSAN expressivity for general k` is  unvalidated. No experiments test k=2 or higher, no comparison with k-OSAN for k>1, and no demonstration of distinguishing graph pairs that require k>1. This makes the paper's main contribution more claimed than demonstrated.\n\n2. Theorem 4.2 only proves superiority over RWSE+MPNN for ReLU, not k-WL equivalence. The paper fails to adequately address this gap between theoretical claims and practical implementation. The k-WL distinguishability claim in the introduction is therefore misleading for the models actually tested.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Poster)