Implicit Graph Neural Networks (IGNNs) define node representations as fixed points of graph neural operators, enabling effectively infinite-depth propagation, iteration-independent parameterization, and flexible test-time computation. Yet these benefits depend on the equilibrium being unique and reached by fixed-point iteration. Existing constructions often impose constraints on recurrent updates to obtain these guarantees, limiting the transformations available at equilibrium. This raises a central question: can IGNNs gain expressiveness through richer, edge-dependent transformations while retaining the inherent strengths of their equilibrium formulation? We introduce SheafDEQ, a subhomogeneous deep-equilibrium architecture with adaptive neural-sheaf propagation. Its learned, matrix-valued sheaf restriction maps can align, mix, or reverse neighbouring representations. Under mild regularity conditions, we prove that SheafDEQ admits a unique equilibrium reached globally by fixed-point iteration from any positive initial state. Contractivity further guarantees convergence under bounded communication staleness. We evaluate SheafDEQ on distributed-inference tasks requiring repeated nonlocal aggregation and on community detection whose rewiring increasingly favours cross-community interactions. SheafDEQ improves over implicit baselines on Sums, MNIST Terrain, and Coordinates, and on community detection as connectivity becomes increasingly heterophilic. SheafDEQ remains competitive with recurrent models without equilibrium guarantees. Continued-iteration diagnostics show decreasing residuals and retained test accuracy through 100 iterations, unlike the finite-horizon control, while delayed-update experiments show low sensitivity to bounded communication staleness.
Figures & tables
Architecture
Propagation primitive
Well-posedness mechanism
Gradient computation
Inference procedure
Finite-horizon recurrent models
APPNP ( Gasteiger et al., 2019 )
graph convolution
not required
unrolled differentiation
K -step recurrence
Loop Transformer (control)
graph attention
not required
unrolled differentiation
K -step recurrence
GraphAttnDEQ (ablation)
graph attention
not provided
unrolled differentiation
K -step recurrence
Energy-defined equilibrium model
EnergyGNN ( Solodova et al., 2025 )
edge energy
strongly convex energy
implicit differentiation
optimization
Table 1: Comparison of the evaluated graph architectures.
Model
Lhetero
1
3
5
7
9
10
Non-adaptive models
MLP
33.6 [33.3,34.8]
33.3 [33.3,37.4]
33.3 [32.8,33.4]
33.3 [31.5,33.5]
33.3 [33.3,34.0]
33.3 [33.3,36.3]
APPNP
40.0 [39.4,40.5]
40.8 [39.4,40.8]
40.7 [40.6,41.0]
40.3 [40.0,40.8]
40.4 [40.3,40.7]
40.0 [39.5,40.2]
IGNN
42.2 [40.9,42.8]
43.6 [43.0,44.6]
44.6 [43.9,45.0]
45.0 [44.8,45.4]
45.8 [45.5,46.7]
52.7 [51.9,54.2]
Adaptive finite-depth/horizon models
Table 2: Community detection. Test accuracy (%) across rewiring levels Lhetero∈{1,3,5,7,9,10} , reported as median [min,max] over five parameter-initialization seeds. For each model family, we select the configuration with the highest mean validation accuracy at each level; the two SheafDEQ normalization variants are selected separately. Pre-norm =✓ denotes the proposed model and pre-norm =× its ablation. 1st , 2nd , and 3rd mark the highest test medians. Full results appear in Appendix Table 11 .
Model
Counting
Sums
MNIST
Coordinates
−MSE↑
−MSE↑
Acc. ↑
−Ldist↑
Non-adaptive and energy-based models
MLP
-0.367 [-0.574,-0.347]
-0.154 [-0.999,-0.145]
53.3 [53.3,53.3]
-0.105 [-0.106,-0.100]
APPNP
-0.335 [-0.348,-0.309]
-0.136 [-0.140,-0.134]
74.8 [74.3,75.3]
-0.030 [-0.032,-0.029]
IGNN
-0.336 [-0.681,-0.307]
-0.130 [-0.998,-0.128]
76.1 [75.6,79.9]
-0.028 [-0.028,-0.027]
EnergyGNN
-0.705 [-0.895,-0.305]
-0.184 [-0.589,-0.140]
60.4 [53.3,68.2]
— †
Table 3: Synthetic distributed tasks. Test performance reported as median [min,max] over five parameter-initialization seeds; higher is better. For each model family, we select the configuration with the highest mean validation performance for each task; the two SheafDEQ normalization variants are selected separately. The † denotes that the EnergyGNN optimization solver did not converge for any seed; no predictive score is reported. 1st , 2nd , and 3rd mark the highest medians; ties at the displayed precision share a rank. Full results and baseline-coverage details appear in Appendix Table 12 and Appendix E.9 .
Figure 1: Finite-iteration recurrence diagnostics. Panels show (a) the median relative update residual and (b) test accuracy under continued recurrence. For the Loop Transformer, panel (a) uses one-step logit updates. For SheafDEQ, solid curves and shading in panel (a) show medians and interquartile ranges across six initial-state conditions. The Loop Transformer does not accept an externally specified recurrent initial state, so its curve is unshaded; the dashed line in panel (b) marks its training horizon Ktrain=20 .
Model
dout
SheafDEQ, pre-norm =✓
7.08×10−4[7.08,7.09]×10−4
SheafDEQ, pre-norm =×
4.12×10−2[1.57,5.32]×10−2
Loop Transformer (20 loops)
7.09×10−2[6.32,7.93]×10−2
Table 4: Bounded-staleness sensitivity at (supd,τmax)=(5,2) . Median [min,max] relative logit discrepancy dout across 20 schedules; lower is better. Complete results appear in Appendix E.11 .
Appendix figures & tables12 assets
Supplementary material from the paper’s appendix.
Appendix
Symbol
Shape/domain
Meaning
G=(V,E)
∣V∣=n
Undirected graph with node set V and edge set E .
A
Rn×n
Raw adjacency matrix, without added self-loops.
DG
Rn×n
Scalar node-degree matrix DG=diag(A1) .
A+,D+
Rn×n
Self-loop quantities A+=A+I and D+=diag(A+1) .
A
Rn×n
Symmetric normalized adjacency DG−1/2ADG−1/2 .
A+
Rn×n
Symmetric normalized adjacency with self-loops (D+)−1/2A+(D+)−1/2 .
Appendix
Table 5: Core graph, sheaf, and equilibrium notation.
Symbol
Domain
Meaning
k
N
Iteration index.
K
N
Finite iteration or layer budget.
Ktrain
N
Iteration horizon used during training.
Kref
N
Iteration used to define the numerical reference state and prediction.
kNN
N
Number of nearest neighbours used to construct a nearest-neighbour graph.
Lhetero
{0,…,10}
Community-detection perturbation level.
Appendix
Table 6: Iteration, benchmark, and diagnostic notation.
Model
Forward solve
Backward gradient
MLP
feedforward (no iteration)
autograd
APPNP family
dataset-dependent K Picard steps
unrolled autograd
IGNN
Picard to tolerance
implicit (adjoint Picard)
EIGNN
closed-form (spectral)
autograd through eigh(M(F))
EnergyGNN (edge)
L-BFGS to tolerance
implicit (CG + Hessian–vector products)
GraphAttnDEQ
dataset-dependent K weight-tied steps
unrolled autograd
Appendix
Table 7: Forward/backward mechanism per comparison model. IGNN, EIGNN, and EnergyGNN use implicit or closed-form solutions specifically to avoid backpropagating through an unrolled iterate trajectory; the APPNP family, GraphAttnDEQ, and Loop Transformer are finite-horizon, weight-tied models trained by ordinary unrolled backpropagation.
Dataset
Total nodes
Train
Validation
Test
Counting
1,275
765
255
255
Sums
100,000
60,000
20,000
20,000
MNIST Terrain
147,800
88,680
29,560
29,560
Coordinates
30,000
18,000
6,000
6,000
Appendix
Table 8: Train, validation, and test split sizes for the distributed benchmarks. The first three tasks use node-level masks; Coordinates uses a graph-level split.
Dataset
c
K
d
q
Metric
Counting
6
50
3
2
MSE
Sums
6
50
3
2
MSE
MNIST Terrain
16
10
2
8
Acc.
Coordinates
16
10
2
8
MSE
Community (all Lhetero )
48
20
3
16
Acc.
Appendix
Table 9: Per-dataset configuration used by the recurrent comparison models; DNSD uses separate 12-layer configurations for the community benchmark (Appendix D.6 ). c : total node representation dimension. K : iteration/layer budget. d : stalk dimension for SheafDEQ; d also sets the GraphAttnDEQ and Loop Transformer head counts. q=c/d : number of sheaf feature channels.
Model
Community
MLP (depth parity)
42,675
APPNP
2,691
IGNN
2,643
EIGNN
2,643
EnergyGNN
340
GraphAttnDEQ
9,555
Appendix
Table 10: Audited trainable-parameter counts from the current implementation for the community configuration. Experimental checkpoints may differ by at most one c -dimensional bias vector (at most 48 parameters) because of a minor code revision; this bookkeeping difference does not affect any reported result. APPNP shift variants, DNSD activation variants of a fixed map type, and the two SheafDEQ normalization variants have identical counts.
Model
L1
L3
L5
L7
L9
L10
MLP
33.6 [33.3,34.8]
33.3 [33.3,37.4]
33.3 [32.8,33.4]
33.3 [31.5,33.5]
33.3 [33.3,34.0]
33.3 [33.3,36.3]
APPNP
39.4 [39.3,39.8]
38.6 [38.5,39.1]
37.8 [36.3,38.9]
37.9 [37.7,37.9]
38.4 [37.1,38.7]
38.0 [37.8,38.4]
APPNP- βshift=1.25
40.0 [39.4,40.5]
40.8 [39.4,40.8]
40.7 [40.6,41.0]
40.1 [39.8,40.7]
40.2 [39.8,40.7]
40.0 [39.5,40.2]
APPNP- βshift=1.65
40.6 [40.5,40.9]
40.8 [40.6,40.9]
40.3 [39.0,40.8]
40.3 [40.0,40.8]
40.4 [40.3,40.7]
40.1 [39.3,40.8]
IGNN
42.2 [40.9,42.8]
43.6 [43.0,44.6]
44.6 [43.9,45.0]
45.0 [44.8,45.4]
45.8 [45.5,46.7]
52.7 [51.9,54.2]
Loop Transformer
40.4 [40.0,43.4]
41.0 [38.7,42.9]
41.8 [40.7,45.5]
98.8 [97.0,99.2]
90.5 [86.2,92.6]
99.1 [98.8,99.1]
Appendix
Table 11: Community detection: test accuracy (%) across rewiring levels L1–L10. Each entry reports the median [min,max] over five parameter-initialization seeds. GraphAttnDEQ rows show all βshift configurations. SheafDEQ rows show the same configurations ( βshift∈{1.2,1.3,1.4,1.6} ) for both the pre-restriction-normalized variant (covered by the conditional contraction guarantee) and its ablation. Wide ranges at intermediate levels indicate substantial variation across parameter-initialization seeds, but the aggregate statistics do not distinguish optimization variability from differences in the finite-iteration states reached at inference. 1st , 2nd , 3rd best per column.
Figure 2: Seed-dependent training trajectories. Validation and test accuracy for pre-restriction-normalized SheafDEQ at Lhetero=3 with βshift=1.2 . Colours denote parameter-initialization seeds; dashed lines mark the best validation epoch for each seed.
Model
Counting
Sums
MNIST
Coordinates
−MSE↑
−MSE↑
Acc. ↑
−Ldist↑
MLP
-0.367 [-0.574,-0.347]
-0.154 [-0.999,-0.145]
53.3 [53.3,53.3]
-0.105 [-0.106,-0.100]
APPNP
-0.335 [-0.348,-0.309]
-0.136 [-0.140,-0.134]
74.8 [74.3,75.3]
-0.030 [-0.032,-0.029]
APPNP- βshift=1.25
-0.342 [-0.360,-0.318]
-0.135 [-0.142,-0.134]
53.3 [53.3,53.3]
-0.044 [-0.044,-0.043]
APPNP- βshift=1.65
-0.337 [-0.358,-0.320]
-0.139 [-0.142,-0.138]
61.3 [60.8,62.2]
-0.051 [-0.064,-0.041]
IGNN
-0.336 [-0.681,-0.307]
-0.130 [-0.998,-0.128]
76.1 [75.6,79.9]
-0.028 [-0.028,-0.027]
Appendix
Table 12: Complete distributed benchmark results. Each entry reports the median [min,max] over five independent parameter-initialization seeds. MNIST reports accuracy (%); Counting and Sums report −MSE ; Coordinates reports −Ldist (all higher is better). GraphAttnDEQ rows show all βshift∈{1.2,1.3,1.4,1.6} configurations. SheafDEQ rows show the same configurations for both normalization variants. The † denotes that the EnergyGNN optimization solver did not converge for any seed; no predictive score is reported. Other “—” entries denote unavailable evaluations. 1st , 2nd , 3rd best per column; ties at the displayed precision share a rank.
Figure 3: Complete finite-iteration recurrence diagnostics. Panels show (a) the median relative update residual, (b) the median state discrepancy from the numerical reference state at Kref=100 , (c) the median prediction discrepancy, and (d) test accuracy under continued recurrence. For the Loop Transformer, panel (a) uses one-step logit updates and panel (b) is unavailable. Solid curves and shading show medians and interquartile ranges across six initial-state conditions; the dashed line marks its training horizon Ktrain=20 . The references are defined at the finite iteration Kref , and the default prediction discrepancy is zero at Kref by construction.
Without pre-restriction normalization K=500 numerical reference; residual 1.04×10−2
No normalization
(2,1)
4.00×10−2;[7.02×10−3,4.76×10−2]
Appendix
Table 13: Complete delayed-schedule results. Each entry is the median [min,max] relative Frobenius logit discrepancy across 20 schedules. Both SheafDEQ models use u=500 mean local updates per node. Loop Transformer is evaluated for its native 20 loops and is reported only as a finite-horizon control.
Figure 4: Post-slot asynchronous inference trajectories. Relative Frobenius logit discrepancy from the synchronous reference versus mean local updates per node. Thin curves denote the 20 sampled schedules; thick curves show their median. All points are recorded after complete logical-slot commits, at their observed work coordinates, without interpolation. For pre-restriction-normalized SheafDEQ, the delayed conditions (2,1) , (5,2) , and (8,4) reach nearly identical discrepancies at the evaluation budget. The ablation without pre-restriction normalization displays substantially larger discrepancies and schedule variation. The (1,0) synchronous control is shown only through iteration 420: later discrepancies reflect comparison with a finite numerical reference, rather than delayed communication.
Deep Graph Neural Networks (GNNs) are essential for capturing complex dependencies in graph-structured data. However, scaling GNNs to depth remains challenging, as stacking layers leads to representation collapse and diminishing sensitivity due to repeated aggregation. While Neural Sheaf Diffusion (NSD) provides strong theoretical guarantees against such collapse, these guarantees do not translate to practice: as depth increases, the disagreement signal of the sheaf Laplacian vanishes, limiting the contribution of deeper layers. We identify mechanisms that hinder NSD effectiveness at depth and propose Deep Neural Sheaf Diffusion (DNSD), which replaces the sheaf Laplacian with a sheaf adjacency operator to maintain informative signals across layers. This is complemented by normalization, odd nonlinearities, and gating. To provide a principled explanation of the expected performance improvement, we contrast sheaf diffusion to graph attention mechanisms, highlighting that DNSD replaces scalar attention scores with matrix-valued edge functions and normalizes node representations rather than attention scores. We demonstrate empirically that DNSD effectively utilizes deep aggregation in graph tasks, outperforming GNN and NSD baselines with up to 30pp accuracy on synthetic long-range datasets, and consistently outperforming them on real-world benchmarks. These results position sheaf-based architectures as a promising building block for graph foundation models by supporting effective deep architectures.
Sheaf Neural Networks (SNNs) were introduced as an extension of Graph Convolutional Networks to address oversmoothing on heterophilous graphs by attaching a sheaf to the input graph and replacing the adjacency-based operator with a sheaf Laplacian defined by (learnable) restriction maps. Prior work motivates this design through theoretical properties of sheaf diffusion and the kernel of the sheaf Laplacian, suggesting that suitable non-identity restriction maps can avoid representations converging to constants across connected components. Since oversmoothing can also be mitigated through residual connections and normalization, we revisit a trivial sheaf construction to ask whether the additional complexity of learning restriction maps is necessary. We introduce an Identity Sheaf Network baseline, where all restriction maps are fixed to the identity, and use it to ablate the empirical improvements reported by sheaf-learning architectures. Across five popular heterophilic benchmarks, the identity baseline achieves comparable performance to a range of SNN variants. Finally, we introduce the Rayleigh quotient as a normalized measure for comparing oversmoothing across models and show that, in trained networks, the behavior predicted by the diffusion-based analysis of SNNs is not reflected empirically. In particular, Identity Sheaf Networks do not appear to suffer more significant oversmoothing than their SNN counterparts.
Ferran Hernandez Caralt, Mar Gonzàlez i Català, Adrián Bazaga +1
Department of Computer Science and Technology, University of Cambridge, Cambridge, United Kingdom · Microsoft
Sheaf Neural Networks (SNNs) generalize message passing by replacing scalar edge weights of standard Graph Neural Networks (GNNs) with learnable, edge-dependent restriction maps between node stalks. Despite their strong theoretical foundations and promising transductive results, SNNs have been evaluated almost exclusively on transductive node classification, leaving their behaviour under inductive protocols unknown. We address this gap through the first systematic benchmark of the sheaf design space, evaluating three diffusion mechanisms (neural sheaf diffusion, sheaf attention, and sheaf attention with Graph Attention Network v2), three restriction-map parameterizations, three stalk dimensions, and six modern GNN architectural components, within a message-passing reformulation that never assembles the heavy sheaf Laplacian, making the full design space trainable under cross-graph batching. Across 1,890 controlled experiments on 14 inductive datasets, multiple insights emerge: restriction maps are the dominant design choice and general maps are preferable, larger stalks add capacity but not long-range reach, architectural components explain more performance variation than the entire sheaf-specific design space itself. Under a matched protocol, SNNs transfer to inductive settings but do not reach the strongest baselines, with gaps being dataset-dependent. Practically, a single sheaf configuration can generalize across datasets, so effort is better spent tuning the surrounding architectural recipe than the sheaf operator itself.
Stefano Fiorini, Edoardo Coppola, Pietro Liò
1Independent Researcher · Department of Computer Science and Technology, University of Cambridge