Organizations: LIX, Ecole Polytechnique, IP Paris · Université Paris Saclay, Université Paris Cité, ENS Paris Saclay, CNRS, SSA, INSERM, Centre Borelli, F-91190, Gif-sur-Yvette, France · Centre d’Analyse et de Mathématique Sociales, EHESS, CNRS, 75006 Paris, France..
Graph Neural Networks (GNNs) learn node representations through iterative network-based message-passing. While powerful, deep GNNs suffer from oversmoothing, where node features converge to a homogeneous, non-informative state. We re-frame this problem of representational collapse from a bifurcation theory perspective, characterizing oversmoothing as convergence to a stable ``homogeneous fixed point.'' Our central contribution is the theoretical discovery that this undesired stability can be broken by replacing standard monotone activations (e.g., ReLU) with a class of functions. Using Lyapunov-Schmidt reduction, we analytically prove that this substitution induces a bifurcation that destabilizes the homogeneous state and creates a new pair of stable, non-homogeneous patterns that provably resist oversmoothing. Our theory predicts a precise, nontrivial scaling law for the amplitude of these emergent patterns, which we quantitatively validate in experiments. Finally, we demonstrate the practical utility of our theory by deriving a closed-form, bifurcation-aware initialization and showing its utility in real benchmark experiments.
Figures & tables
Figure 1: Effective potential landscape across bifurcation. Left: Supercritical Pitchfork bifurcation. Right: Transcritical bifurcation. For both bifurcation types, below the critical coupling w<1 , the system is at the subcritical regime – the homogeneous state a=0 is the unique stable minimum, corresponding to oversmoothing. Above w>1 , the origin becomes unstable (hollow circle) and stable minima emerge (filled circles), representing non-homogeneous pattern formation that resists oversmoothing.
Figure 2: Phase diagram validating Theorem 3.2 . Left : Dirichlet energy ED as a function of activation slope α and coupling w (dashed red) precisely separates the oversmoothing regime ED≈0 for the pattern forming regime ED>0 . Right : Node representations at marked points. In the subcritical regime all representations collapse to zero (homogeneous fixed point). In the supercritical regime, a structured pattern emerges.
Figure 3: Empirical validation of theoretical predictions. Dirichlet energy, amplitude and scaling laws versus normalized coupling wkw on Erdos-Renyi, Barabasi-Albert, Watts-Strogatz and Random-Regular graph topologies
Figure 4: Spectral filtering controls mode selection. Top : Bandpass filters P(λ) centered at low, mid and high frequencies. The dashed line is the critical threshold αw1 , only modes that exceed this threshold can bifurcate. Bottom : Stable patterns on a 10×10 grid graph obtained by iterating x(ℓ+1)=ϕ(wP(A)x(ℓ)) to convergence from random initializations. Each filter selects distinct eigenmode: smooth (low frequency), stripped (mid-frequency), or checkerboard (high-frequency), which illustrates topological prior selection without collapse.
Figure 5: Depth Robustness and Phase Transitions on CORA. (a, b) ED0 and accuracy vs. depth. Activations with stabilizing cubic terms (sin, tanh) resist oversmoothing up to 64 layers, whereas ReLU collapses. (c, d) Phase transition at depth 64. Varying the bifurcation parameter δ reveals a transition near δ=0 ; only supercritical initialization ( δ>0 ) enables pattern formation and successful learning.
Method
Texas
Wisconsin
Squirrel
Chameleon
Cornell
Computer
CiteSeer
PubMed
Cora
CoauthorCS
Photo
ogbn-ArXiv
Homophily
0.11
0.21
0.22
0.23
0.30
0.78
0.74
0.80
0.81
0.81
0.83
0.83
GCN
55.14 ± 5.16
51.76 ± 3.06
53.43 ± 2.01
64.82 ± 2.24
60.54 ± 5.30
82.6 ± 2.4
76.50 ± 1.36
88.42 ± 0.50
86.98 ± 1.27
91.1 ± 0.5
91.2 ± 1.2
72.17 ± 0.33
GCN + PairNorm
60.27 ± 4.34
48.43 ± 6.14
50.44 ± 2.04
62.74 ± 2.82
58.92 ± 3.15
73.59 ± 1.47
–
87.53 ± 0.44
85.79 ± 1.01
–
–
–
GraphSAGE
82.43 ± 6.14
81.18 ± 5.56
41.61 ± 0.74
58.73 ± 1.68
75.95 ± 5.01
82.4 ± 1.8
76.04 ± 1.30
88.45 ± 0.50
86.90 ± 1.04
91.3 ± 2.8
91.4 ± 1.3
–
GAT
52.16 ± 6.63
49.41 ± 4.09
40.72 ± 1.55
60.26 ± 2.50
61.89 ± 5.05
78.0 ± 19.0
76.55 ± 1.23
87.30 ± 1.10
86.33 ± 0.48
90.5 ± 0.6
85.7 ± 20.3
73.65 ± 0.11
CGNN
71.35 ± 4.05
74.31 ± 7.26
29.24 ± 1.09
46.89 ± 1.66
66.22 ± 7.69
80.29 ± 2.0
76.91 ± 1.81
87.70 ± 0.49
87.10 ± 1.35
92.3 ± 0.2
91.39 ± 1.5
57.80 ± 2.5
Table 1: Node classification accuracy (%) on benchmark datasets. Results are reported as mean ± standard deviation across multiple runs. Datasets are ordered by increasing homophily levels, where higher values indicate more homophilic graphs. Gold , silver , and bronze indicate 1st, 2nd, and 3rd best results. Results for baseline methods on ogbn-ArXiv are taken from the corresponding papers.
Method
Mean Rank
Sine-Poly(A)-sh
3.00 ± 2.13
Sine-Poly(A)-pl
3.18 ± 3.31
GELU-Poly(A)-pl
4.82 ± 1.78
GELU-Poly(A)-sh
5.32 ± 2.56
ReLU-Poly(A)-pl
6.64 ± 2.41
ReLU-Poly(A)-sh
6.91 ± 3.56
Table 2: Ablation by average rank ( ± std) across benchmark datasets. Methods are ranked by accuracy on each dataset (1 = best), and the average rank is computed across all datasets except large-scale OGBN-arXiv. Lower is better.
Appendix figures & tables4 assets
Supplementary material from the paper’s appendix.
Appendix
Dataset
#Nodes
#Edges
#Features
#Classes
h%
Avg. N.D.
Diameter
Class Type
Cora
2,708
5,278
1,433
7
80.4
3.9
∼ 19
Citation networks
Citeseer
3,327
4,676
3,703
6
73.5
2.8
∼ 16
Citation networks
Pubmed
19,717
44,327
500
3
80.2
4.5
∼ 21
Citation networks
Computers
13,752
245,861
767
10
77.7
35.8
∼ 14
Co-purchase networks
Photo
7,650
119,081
745
8
82.7
31.1
∼ 12
Co-purchase networks
CoauthorCS
18,333
81,894
6,805
15
80.0
8.9
∼ 15
Co-authorship networks
Appendix
Table 3: Statistics of the benchmark datasets.
Method
Texas
Wisconsin
Squirrel
Chameleon
Cornell
Computer
CiteSeer
PubMed
Cora
CoauthorCS
Photo
Sine-Poly(A)-sh
(3,2)
(3,2)
(4,1)
(4,1)
(3,2)
(4,2)
(2,2)
(4,2)
(3,1)
(2,2)
(3,2)
Sine-Poly(A)-pl
(4,2)
(3,2)
(4,1)
(4,2)
(4,2)
(4,1)
(2,2)
(3,2)
(3,2)
(3,2)
(3,2)
ReLU-Poly(A)-sh
(4,2)
(2,1)
(2,1)
(2,1)
(2,2)
(3,1)
(2,2)
(3,2)
(3,2)
(3,1)
(3,2)
ReLU-Poly(A)-pl
(4,1)
(3,2)
(3,1)
(3,1)
(2,2)
(2,2)
(3,2)
(3,2)
(4,2)
(4,2)
(3,1)
Appendix
Table 4: Best hyperparameters (depth,K) for polynomial- A models.
Method
Mean Rank
Sine-Poly(A)-sh
3.00 ± 2.13
Sine-Poly(A)-pl
3.18 ± 3.31
GELU-Poly(A)-pl
4.82 ± 1.78
GELU-Poly(A)-sh
5.32 ± 2.56
ReLU-Poly(A)-pl
6.64 ± 2.41
ReLU-Poly(A)-sh
6.91 ± 3.56
Appendix
Table 5: Average rank ( ± std) across all benchmark datasets. Methods are ranked by accuracy on each dataset (1 = best), and the average rank is computed across all datasets except large-scale OGBN-arXiv. Lower is better.
Method
Texas
Wisconsin
Squirrel
Chameleon
Cornell
Citeseer
Computers
Pubmed
Cora
CoauthorCS
Photo
Hom. level
0.11
0.21
0.22
0.23
0.30
0.74
0.78
0.80
0.81
0.81
0.83
Sine-Poly(A)-sh
90.09 ± 0.1
85.62 ± 0.0
60.81 ± 0.56
70.81 ± 0.24
76.58 ± 0.1
78.18 ± 1.01
91.56 ± 0.32
89.71 ± 0.09
89.78 ± 0.20
95.84 ± 0.10
95.69 ± 0.17
Sine-Poly(A)-pl
93.69 ± 0.2
80.39 ± 0.0
61.79 ± 0.40
70.85 ± 0.92
76.87 ± 0.1
78.15 ± 0.75
91.94 ± 0.07
88.67 ± 0.27
89.89 ± 0.50
95.42 ± 0.07
95.55 ± 0.21
Sine GCN
52.68 ± 1.19
50.57 ± 2.77
29.25 ± 0.55
44.27 ± 1.47
40.50 ± 2.92
78.12 ± 0.38
91.25 ± 0.10
88.38 ± 0.18
89.85 ± 0.19
94.15 ± 0.12
94.50 ± 0.35
ReLU-Poly(A)-sh
92.79 ± 0.1
86.27 ± 0.1
53.57 ± 0.82
69.63 ± 0.35
67.57 ± 0.2
77.56 ± 0.55
91.80 ± 0.43
88.36 ± 0.14
89.36 ± 0.18
93.42 ± 0.14
94.37 ± 0.25
ReLU-Poly(A)-pl
90.09 ± 0.1
84.31 ± 0.0
53.67 ± 0.46
69.89 ± 0.49
76.58 ± 0.1
77.53 ± 0.55
91.85 ± 0.26
88.46 ± 0.19
89.47 ± 0.21
93.33 ± 0.19
94.09 ± 0.27
Appendix
Table 6: Activation ablation. Node classification accuracy (%) on benchmark datasets without polynomial filtering and for different base activation functions. GELU and Swish alternatives are added. Results are reported as mean ± standard deviation across multiple runs. Datasets are ordered by increasing homophily levels (Hom.), where higher values indicate more homophilic graphs.
In deep graph neural networks, increasing depth enlarges the receptive field but often leads to over-smoothing, where node representations tend to align. We develop a unified, mode-wise stability framework for deep GNN propagation that provides a principled characterization of over-smoothing. By interpreting layer depth as time and layer updates as graph-coupled dynamics, over-smoothing can be understood as an undesirable dynamical synchronization of features, for which the master stability curve provides a theoretical tool to assess the stability of synchrony. Guided by this theory, we further propose Stability-Shaped Deep Graph Learning (SDGL) to mitigate over-smoothing in deep GNNs. SDGL has two complementary instantiations: one induces controlled Turing instability to replace synchronization with spatial pattern formation, and the other maintains stable near-critical propagation. Experiments on diverse node- and graph-level benchmarks demonstrate the improved depth scaling and consistent accuracy gains over strong baselines, including graphs exhibiting long-range dependencies.
Junyou Zhu, Langzhou He, Fenying Cai +7
Department of Complexity Science, Potsdam Institute for Climate Impact Research, 14473 Potsdam, Germany · Machine Learning Group, Technical University of Berlin, 10587 Berlin, Germany · Department of Computer Science, University of Illinois at Chicago, Chicago, IL 60607, USA +6
Oversmoothing is a fundamental limitation of deep graph neural networks (GNNs), where repeated message passing causes node representations to become increasingly similar, eventually collapsing toward a low-dimensional subspace. This phenomenon limits the effective depth of message-passing architectures and motivates the search for mechanisms that preserve representation diversity. In this paper, we study a recurrent graph neural network in which independent Gaussian noise is injected after every propagation step and analyze the resulting architecture as a stochastic dynamical system. Under a standard global contraction assumption on the deterministic update, we prove that the hidden representations form a geometrically ergodic Markov chain admitting a unique invariant probability measure. Our main theoretical result establishes an explicit positive lower bound on the expected stationary Dirichlet energy, proportional to both the noise variance and the spectral gap of the underlying graph. Consequently, the stationary representations cannot collapse onto the constant manifold, providing a rigorous guarantee that asymptotic oversmoothing is prevented in the sense of non-vanishing Dirichlet energy. Our analysis reveals persistent stochastic perturbations as a fundamentally different mechanism for combating oversmoothing, complementing existing deterministic approaches based on residual connections, normalization, and graph rewiring. Finally, numerical experiments on both linear and nonlinear recurrent graph neural networks closely match the theoretical predictions, illustrating the emergence of a stationary distribution and the predicted dependence of the limiting Dirichlet energy on the noise intensity.
Mostafa Haghir Chehreghani
Department of Computer Engineering Amirkabir University of Technology (Tehran Polytechnic) Tehran, Iran
Graph Neural Networks (GNNs) face two fundamental challenges when scaled to deep architectures: oversmoothing, where node representations converge to indistinguishable vectors, and oversquashing, where information from distant nodes fails to propagate through bottlenecks. Both phenomena are intimately tied to the underlying graph structure, raising a natural question: can we optimize the graph topology to mitigate these issues? This paper provides a theoretical investigation of the computational complexity of such graph structure optimization. We formulate oversmoothing and oversquashing mitigation as graph optimization problems based on spectral gap and conductance, respectively. We prove that exact optimization for either problem is NP-hard through reductions from Minimum Bisection, establishing NP-completeness of the decision versions. Our results provide theoretical foundations for understanding the fundamental limits of graph rewiring for GNN optimization and justify the use of approximation algorithms and heuristic methods in practice.
Mostafa Haghir Chehreghani
Department of Computer Engineering, Amirkabir University of Technology (Tehran Polytechnic), Tehran, Iran