Initialization Schemes for Kolmogorov–Arnold Networks: An Empirical Study
Authors: Spyros Rigas
Dhruv Verma
Georgios Alexandridis
Yixuan Wang
Venue: ICLR 2026 Poster (2026)
Keywords: Kolmogorov-Arnold networks
weight initialization
NTK
Function Fitting
PDEs
Primary area: other topics in machine learning (i.e., none of the above)
TL;DR: We propose new initialization schemes for KANs, showing that our Glorot-inspired method outperforms the default on large architectures, while our empirical power-law scheme consistently achieves the best results.

Abstract: Kolmogorov–Arnold Networks (KANs) are a recently introduced neural architecture that replace fixed nonlinearities with trainable activation functions, offering enhanced flexibility and interpretability. While KANs have been applied successfully across scientific and machine learning tasks, their initialization strategies remain largely unexplored. In this work, we study initialization schemes for spline-based KANs, proposing two theory-driven approaches inspired by LeCun and Glorot, as well as an empirical power-law family with tunable exponents. Our evaluation combines large-scale grid searches on function fitting and forward PDE benchmarks, an analysis of training dynamics through the lens of the Neural Tangent Kernel, and evaluations on a subset of the Feynman dataset. Our findings indicate that the Glorot-inspired initialization significantly outperforms the baseline in parameter-rich models, while power-law initialization achieves the strongest performance overall, both across tasks and for architectures of varying size. This work underscores initialization as a key factor in KAN performance and introduces practical strategies to improve it.

Published as a conference paper at ICLR 2026
INITIALIZATION SCHEMES FOR KOLMOGOROV- ARNOLD NETWORKS: AN EMPIRICAL STUDY
Spyros Rigas"
Dhruv Verma
Department of Digital Industry Technologies National and Kapodistrian University of Athens 344 00 Psachna, Greece spyrigas@uoa.
Applied and Computational Mathematics California Institute of Technology
Pasadena, CA 91125, United States of America dverma@ www.ech..eccccccoccccoooooeceeee edu
Georgios Alexandridis
Yixuan Wang
Department of Digital Industry Technologies National and Kapodistrian University of Athens 344 00 Psachna, Greece
Applied and Computational Mathematics California Institute Technology
Pasadena, CA 91125, United States of America
gealexandri@uoa.cc gr
zoywang@caltech. edu
ABSTRACT
Kolmogorov-Arold Networks (KANs) are recently introduced neural architec- ture that replace fixed nonlinearities with trainable activation functions, offering enhanced flexibility and interpretability. While KANs have been applied suc- cessfully across scientific and machine learning tasks, their initialization strate- gies remain largely unexplored. In this work, study initialization schemes for spline-based KANs, proposing two theory-driven approaches inspired by LeCun and Glorot, well
empirical power-law family with tunable exponents. Our evaluation combines large-scale grid searches on function fitting and forward PDE benchmarks, an analysis of training dynamics through the lens of the Neu- ral Tangent Kernel, and evaluations on a subset of the Feynman dataset. Our findings indicate that the Glopot-inspired initialization significantly outperforms the baseline in parameter-rich models, while power-law initialization achieves the strongest performance overall, both across tasks and for architectures of varying size. This work underscores initialization as key factor in KAN performance and introduces practical strategies to improve it.
INTRODUCTION
Kolmogorov-Ammol Networks (KANs) (Liu et al. 2025) have recently emerged as an alter- native backbone architecture to Multilayer Perceptrons (MLPs), drawing inspiration from the Kolmogorov-Amot representation theorem (Kolmogorov 1957) in manner analogous to how the learning of MLPs relies on universal approximation theorems. Unlike MLPs, which use fixed nonlinear activation functions and trainable synaptic weights, KANs comprise grid-dependent train- able activation functions. This provides them with flexibility in modeling complex nonlinear re- lationships, while requiring fewer and smaller layers. Since their introduction, KANs have found numerous applications, often surpassing the performance of their MLP-based counterparts Yu ct al.
2024 Poeta et al.) 2024). There have been many notable results in scientific problem-solving do-
mains, including function fitting and symbolic regression (Liu et al.) 2024 Shukla et al. 2024), partial differential equations (PDEs) (Shukla al.
2024 Rigas al. 2024 Wang et al. 2025b)
and operator learning (Abueidda et al., 2025 Shukka et al. 2024; Lee et al. 2025). among other
applications (Howard et al. 2024 Kundu al. 2024 Kashen 2025).
Beyond these benchmarks, there has also been significant progress in the theoretical understanding of KANs (Zhang & Zhou 2025: Alter et al. 2025 Wang et al.] (202aa). However, one important
theoretical and practical aspect that remains understudied pertains to their initialization strategies.
"Corresponding author.

---

Published as a conference paper at ICLR 2026
Current literature mainly relies the standard initialization method proposed in the introductory KAN paper (Liu et al. 2025). or explores alternative KAN variants such as Chebyshev-bard for- mulations (Rigas al.
2026) This highlights clear gap and motivates investigation into more
effective initialization approaches for the standard spline-based architecture. Effective initialization is erucial, as good "initial guess" for the network weights can significantly accelerate training (Mishkin & Matas 2016; Skorski et al. (222)) and prevent early saturation of hidden layers (Glorot & Bengio 2010). However, despite extensive research into initialization methods for MLP-based ar- chitectures, these results cannot be directly applied to KANs. Furthermore, even within MLP-baser architectures, initialization methods often require separate consideration depending on the specific architecture design (Huang et al. 2020), activation function (He et al. 2015). or even on a complete case-by-case basis (Skorski et al. 2021).
In response to this research gap, this work explores initialization strategies for the standard, spline- based KAN architecture. Drawing parallels with MLPs, we propose serrance-preserving schemes inspired by LeCun (LeCun et al. 1998) and Glorot (Glorot & Bengio 2010) initializations, includ- ing a variant that employs normalize zed spline basis functions. addition, recognizing that theoretical frameworks may not always align with empirical performance (Mishkin & Matas 2016), we further propose an empirical family of power-law initializations parameterized by two exponents. We per- form grid searches to identify suitable exponent choices for the power-law method on function fitting and forward PDE benchmarks and then evaluate all initialization schemes on these tasks. We sub- sequently fix the exponents to values that lie within the identified well-performing range and select two representative architectures (in terms of parameter count) to analyze training dynamics through the evolution of training loss curves and the Neural Tangent Kernel (NTK) spectrum (Jacot al. 2018 Wang et al., (2022). Finally, we evaluate said architectures on subset of the Feynman dataset (Udrescu Tegmark| [2020), which, although widely used for symbolic regression, is formulated here as a function fitting benchmark as in Liu et al. (2025).
BACKGROUND
2.1 KOLMOGOROV-ARNOLD NETWORKS
Within the standard formalism, the output, y E R"",,, of KAN KANlayer is related to its input, x e R",
via:
Flm Σ
G+
Σ bjim Bn (2:)) (ct)) = I,..., Douts
Tzi Rx,) + Cje
m=1
where Tjs. Cjs and fjim are the layer's trainable parameters, R(x) corresponds to residual function, typically chosen as the SiLU, i.c., R(x) = x (1+e and Bm (x) denotes a univariate spline basis function of order k defined a grid with G intervals. For each of the layer's trainable parameters, the original KAN formulation initializes the scaling weights as Cji = 1, the residual weights ry using Glorot initialization (Glorot & Bengio, 2010), and the basis weights bjim from normal distribution with zero mean and small standard deviation, typically set to = 0.1. We will henceforth refer to this configuration as the "baseline initialization".
2.2 RELATED WORK
In the existing KAN literature, initialization strategies have only been explored in certain KAN vari- ants (see, e.g., Guilhoto & Perdikaris (2025), while the standard spline-basee architecture has not yet received dedicated attention this regard. A natural starting point for studying initialization is to follow the historical developments in MLP-based architectures, beginning with variance- preservin schemes such as those proposed by LeCun (LeCun et al., 1998) and Glorot (Glorot & Bengio 2010), which stabilize activation variance across layers and mitigate progressive vanishing or explosion. Within the KAN family, Glorot-inspired initialization has been applied successfully to Chebyshev- based variants (Rigas al. 2026). though this setting differs substantially from the spline-based case studied here, since it removes the residual term of Eq. (T) and employs completely different basis function. Consequently, it remains unclear whether such strategies directly transfer the standard
2

---

Published as a conference paper at ICLR 2026
KAN formulation, motivating the investigation presented in this work. To the best of our knowl- edge, the present work provides the first systematic study of initialization strategies for spline-based KANs.
3 METHODOLOGY
3.1 PROPOSED INITIALIZATIONS
Since the three weight types in KAN layer independent, may initialize the scaling weights Cjz to I and focus exclusively on the initialization of the residual weights ry and basis weights bjon- We assume that these weights are drawn from zero-mean distributions with standard deviations or and ob, respectively. To determine suitable values for a, and Øs we follow the principle of variance preservation proposed by LeCun (LeCun et al.) 1998), which stipulates that the variance of the outputs of each layer should match that of its inputs, thereby avoiding amplification or attemuation of the signal across layers. Assuming statistical independence among terms and an equal contribution to the variance from each of the (G+A+t terms in the summand of Eq. (I), derive the following expressions for the standard deviations
Var(æ,)
Var(x;)
or
(in(G k 1) (0)
Tin G k MB
where
PR = E (x,))
(0) #B = Bm (x:)²]
(3)
with Pg denoting the expectation over both the input distribution and all spline basis indices, and PR denoting the expectation over the input distribution alone.
(0)
If we further assume that each component of is drawn from a given distribution (e.g., the uniform distribution 4 (-1, 1), as often the case tasks like function fitting or PDE solving), then all statistical quantities in Eq. ( can be evaluated directly, except for pn
(0) Due to the dependence of
the baline-basis functions on the underlying grid, no general analytic expression exists for age This
leads to two practical alternatives: one may either estimate 11g numerically by sampling large number of input points from the assumed distribution initialization, or set the expectation value to unity by modifying the architecture of the KAN layer to use normalized spline basis functions, defined as
(0)
Bm (xx)
E [m (x:))
Bm (x:)
(4)
PB
E2 Bm (x+)
where the expectation values are computed over the layer inputs during the forward pass. We will refer to the former alternative as "LoCun-numerical" initialization, while the latter is referred as norrmaleorma initialization.
While these LeCun-inspired schemes focus on preserving the variance of forward activations, they do not explicitly account for the propagation of gradients. To address this, we also consider Glorot- inspired initialization, which aims to balance forward- and backward- pasds variance by maintaining stable variance for both activations and gradients across layers. Under the same assumptions as before, we derive the following expressions for the standard deviations?
2
o, V G+k+1
(0) Thn PR
11)
9b V G+++1 Thm HB nou HB
(0)
())
(5)
+ Pl M
See Appendix for detailed derivations.
"ee Appendix
B for detailed derivations.

---

Published as a conference paper at ICLR 2026
where
(1) PR = R' (
(1) = E B', (x;)²
(6)
with the expectations defined analogously to PR and #B computed using automatic differentiation of the spline basis functions, together with the numerical sampling strategy discussed for the
(I) (0) in Eq. 3 In practice, PB can be
case, while ( can be evaluated analytically for standard choices of R(x) such the SiLU.
In addition to these theory-driven initialization strategies, we also investigate an empirical approach based on power-law scaling of the KAN layer's architectural parameters. Specifically, initialize the weights such that their standard deviations follow the form
Thin (G
k)
1iin (G l:
where and 8 are tunable exponents selected from the set {0.0, 0.25. ..., 1.75, 2.0}. The motivation behind this empirical scheme is to perform a grid search over (a. (3) configurations in order to identify trends specific exponent pairs that consistently improve training speed and convergence. Such searches can be carried out on per-domain basis (e.g., function fitting, forward PDEs), after which the resulting well-performing ranges may serve reusable heuristics for future problems of the same type.
3.2 EXPERIMENTAL SETUP
We evaluate initialization strategies on two benchmark families: function fitting tasks and forward PDE problems. For function fitting, iterations (epochs), while for PDEs we consider the Allen-Cahn equation, Burgers' equation, and the two-dimensional Helmholtz equation, using KANs trained for 5,000 iterations within the Physics- Informed Machine Learning (PIML) framework (Raissi et al. 2019). Across both benchmarks, performance is measured using the final training loss and the relative L error with respect to the reference solution. For the purposes of the initial grid search, we test architectures with 1-4 hidden layers, widths equal to 2 for i ...... 6, and grid sizes G (5, 10, 20, 40) for function fitting, while for PDEs we restrict to G E (5, 10, 20). All experiments presented herein are repeated with five random seeds, except in the power-law grid search where we use three seeds reduce computational cost, and we report the median outcome across runs. Further implementation details, including the explicit formulas of the target functions and the PDE setups, are provided in Appendix C All experiments are implemented in JAX (Bradbury et al.) 2018), with KANs trained using the JaxKAN framework (Rigas & Papachristou 2025). Training is performed on a single NVIDIA GeForce RTX 4090 GPU.
five wwo-dimensional target functions and train for 2,000
EXPERIMENTS & DISCUSSION
4.1 GRID-SEARCH RESULTS
The grid search over (a, (3) configurations and the architectural variations described in Section 3.2 resulted 126,240 trained KAN model instances for function fitting. After aggregating the repeated runs by their median outcome, this number reduces to 40,800 representative results. From these, we retain only the best-performing (a,13) configuration per setting, yielding 2,400 final entries. Table 1 reports, for each target function and initialization scheme, the percentage of runs that outperform the bascline initialization described Section 2.1 Results are compared with respect to final training loss and relative I error, and additionally report the percentage of runs where both metrics improve simultaneously.
The LeCun-insping schemes rarely outperform the baseline the smaller architectures, but their effectiveness increases with depth, width and grid size. In some of the larger settings, the nor- malized variant in particular achieves improvements exceeding two orders of magnitude relative
4

---

Published as a conference paper at ICLR 2026
the baseline. However, in terms of absolute frequency, Table [ clearly shows that the baseline still outperforms both LeCun-basee variants: for more than 70% of configurations, the resulting rela- tive L error under LeCun initialization is higher than under the baseline scheme. Between the two LeCun variants, the normalized version consistently performs better than the numerical one, which is consistent with its design, variance preservation is enforced by construction. On the other hand, the Glorot-inspired initialization performs more robustly. Across all five functions, it yields success rates of approximately 60-75% for both loss and relative L2 error, indicating that simultaneously balancing forward- and backward- -aaaaddpasssssa variances is considerably more effective than preservation alone. The few cases where the baseline performs better occur pre- dominantly for the smaller architectures.
As far as the power-law initialization is concerned, it exhibits the strongest overall performance. Table [ shows that, for nearly every architecture and target function, there exists at least one (a, 8) pair that outperforms the baseline, often by substantial margin. The most favorable region concentrated around small residual exponents, i.e., E {0.25, 0.5}. combined with larger basis- function exponents, namely 3 2 1.0. Full grid-search results illustrating these trends are provided in Appendix D and the supplementary material. Notably, even when fixing single configuration within this region, the method remains highly robust. For example, with (a, B) = (0.25, 1.0), the initialization yields simultaneous improvements in both loss and 1,2 error over the baseline in 96.88% of runs for Si (z. y). 95.83% for f2 [2, y). 97.92% for fa y), 87.50% for sa y), and 89.58% for fz(2.)) y). This indicates that once a suitable exponent range is identified for given problem type, a fixed choice within that region can systematically outperform both the baseline and the Glorot initialization.
Following the same procedure for the PDE benchmarks, we trained 56,882 models, number which reduced to 18,360 representative results after aggregation and 1, 1,080 final entries after selecting the best (a, 3) per setting. Table 2 summarizes the outcomes in terms of final training loss, relative L error, and their joint improvement over the baseline initialization.
When comparing the two LeCun- schemes to each other, the observed behavior is essentially the inverse of what we observed in function fitting. The normalized variant, which performed better than the numerical one for function fitting, fails almost entirely in PDE problems: in nearly all configurations, it produces no improvement over the baseline, with success rates effectively equal to zero. This discrepancy can be attributed to the fact that PDE losses involve not only the network output but also its derivatives and nonlinear combinations thereof. Normalizing the spline basis propagates all multiplicative constants (the standard deviation in the denominator of Eq. 44) in this case) into all derivatives, altering the stiffness of the residuals. While the numerical variant avoids this issue, the results are still poor: although it occasionally outperforms the baseline on larger architectures, it remains generally ineffective. The Glorot-inspired initialization again shows more consistent improvements. As in function fitting, it performs significantly better than the baseline on parameter-rich architectures and the cases where it underperforms correspond to smaller models.
Table 1: Percentage of runs that outperform the baseline initialization on function fitting bench- marks. Columns correspond to target functions, while rows correspond to initialization schemes and evaluation metries. Best results per function are shown in bold.
Initialization
Metric
f (a, M)
fa (x. v)
fa (x., y)
fx; y)
fo (x., M)
Loss
18.75%
14.58%
12.50%
25.00%
26.04%
LeCun-numerical
L²
6.25%
4.17%
5.21%
14.58%
2.08%
Both
1.04%
0.00%
0.00%
8.33%
0.00%
Loss
19.79%
28.13%
19.79%
41.67%
31.25%
LeCun nommalized
L?
11.46%
9.38%
11.46%
26.04%
6.25%
Both
2.08%
5.21%
5.21%
16.67%
1.04%
Loss
78.13%
76.04%
78.13%
63.54%
72.92%
Glorot
r,²
78.13%
75.00%
78.13%
64.58%
72.92%
Both
78.13%
75.00%
78.13%
60.41%
64.59%
Loss
100.00%
100.00%
100.00%
100.00%
98.96
Power Law
L2
100.00
100.00% %
100.00%
94.79%
96.88
Both
100.00 %
100.00% %
100.00%
94.79%
95.83
5

---

Published as a conference paper at ICLR 2026
Table 2: Percentage of runs that outperform the baseline initialization forward PDE benchmarks. Columns correspond to the three PDEs considered, while rows correspond to initialization schemes and evaluation metries. Best results per PDE are shown in bold.
Infitalization
Metric
Allen-Cahn
Burgers
Hclmholtz
Loss
11.11%
11.11%
8.13%
LeCun-mumerical
1,2
16.67%
22.22%
15.28%
Both
8.33%
6.94%
2.78%
Loss
2.78%
0.00%
0.00%
LeCun-normalized
L²
0.00%
0.00%
0.00%
Both
0.00%
0.00%
0.00%
Loss
55.56%
50.00%
76.39%
Glorot
L
51.39%
54.17%
72.22%
Both
36.11%
Loss
98.61 %
100.00 %
Power-Law
L
94.44%
73.61
87.50%
Both
94.44%
73.61% %
87.50%
The power-law initialization remains the strongest of all approaches, though the advantage is less pronounced in the PDE case. For more than 70% of configurations (up to 94.44% in the Allen- Cahn equation) there exists at least one (ar,3) pair that outperforms the baseline simultaneously in both loss and relative L error. The preferred exponent region differs slightly from the function fitting case: while small residual exponents remain strongly favored, in this case lower values of 3 (typically 0.75 < B < 1.25) yield the best results. The complete grid-search results supporting these observations are provided in Appendix DameedixDDandd and the supplementary material. For the previously discussed configuration (a,, () = (0.25, 1.0), the power-law initialization outperforms the baseline on both metrics in 83.33% of runs for Alhn, 54.93% for Burgers' and 59. .2% for Helmholtz, with even higher success rates when each metric is considered individually.
4.2 TRAINING DYNAMICS ANALYSIS
The previous experiments established that the Glorot- and power-law-based schemes provide strong alternatives the baseline initialization, with the latter yielding the most consistent improvements overall. In contrast, the AcCun-based variants exhibited substantially weaker performance. To better understand the mechanisms behind these trends, we next examine training dynamics greater detail for each initialization scheme, excluding only the LeCun-normalized variant due to its complete breakdown in PDE problems for the aforementioned reasons. We begin with the function fitting benchmarks: Figure [ shows the evolution of the training loss for two representative settings, "small" architecture (G = 5, two hidden layers with 8 neurons each) and a "large" architecture (G' = 20, three hidden layers with 32 neurons each). For consistency across experiments, we fix the power-law parameters to a = 0.25 and 8 = 1.75, which lie within the range identified as favorable in the grid search for function fitting. The curves averaged over five seeds, with shaded regions indicating the standard error.
Across all settings, the loss curves in Figure 1 reinforce the conclusions drawn from the previous experiments: the power-law initialization consistently outperforms all other schemes, converging both faster and to lower final losses. For the small architecture, the baseline and Glorot initial- izations behave similarly, whereas for the larger architecture Glorot has a clear advantage over the baseline. The LeCun initialization, while stable and not prone divergence, persistently underper- forms the remaining schemes in both architectures. The oscillatory behavior observed for Glorot and pawer-law in the farge-architecture setting is a consequence of using a fixed learning rate throughout training. We intentionally avoided learning-rate scheduling because initialization and learning-rate adaptability are known to interact in subtle ways (e.g.. YYang et al.) (02T)) and our goal was to isolate the effect of initialization alone. In Appendix E we provide results for the same experiments using a learning-rue rate scheduler; the resulting curves are significantly smoother, however the relative performance ranking of the initialization schemes remains unchanged.
The corresponding analysis for the PDE benchmarks is shown in Figure 2 where again we use the power-law law parameters region for 3, in order to demonstrate that the method remains robust even when not tuned specifically for PDEs. For the small architecture, all initialization schemes ultimately perform comparably, al-
0.25 and 1.75, despite the grid search identifying different optimal

---

Published as a conference paper at ICLR 2026
fiX. y)
104
5ix.
10y
f(X.y)
10a
((0.y)
UOPM 2 1pdpp S-D
Training Loss
10-1
10a
10a4
10
10-1
10"4
10:
10-4
10"2
10
10a4
10-2
10
10-3
10.
UPM E uphpp 122 -D
Training Loss
10-3
10-
10"
10-2
10-4
10-
10a
10-
10-1
2000
Training Iteration
Training Iteration
Training Reration
Training Iteration
Training Iteration
Beseline
LeOun
Glorot
Power-Lav
Figure Training loss curves for function fitting benchmarks under baseline, LeCun-numerical, Glorot and power-law (a 0.25, = 1.75) initializations. Results averaged over five seeds, with shaded regions indicating the standard error. Top row: "small" architecture (G = 5, two hidden layers with 8 neurons each). Bottom row: "large" architecture (G 20, three hidden layers with 32 neurons each).
beit different rates. Notably, the power-law initialization exhibits markedly lower variance across seeds and, in the Allen-Cahn case in particular, reaches a minimum significantly faster than the alter- natives. In the larger architecture, however, the trends closely mirror those observed in the function fitting experiments: the LeCun initialization lags behind by several orders of magnitude, the base- line also underperforms, and both Glorot and power-law provide substantial improvements. These two consistently achieve the lowest losses, with the power-law scheme displaying advantage in Allen-Cale and Burgers equations. As with the function fitting results, the oscillatory behavior observed for Glorot and power-law arises from the fixed learning rate used during training. The corresponding smoothed curves obtained with learning-rate scheduling are provided in Appendix
Allen: Cahn
Burgers
Helmholtz
10t
2000 -2 W00
Training
1.7
0
10:1
30a
1000
2000
3000
4000
5000
1000
2000
3000
4000
5000
1000
2000
3000
4000
5000
300
ZE GRDW E udop
3.0
10
10a
Training
107
107
10a
101
1000
2000
4300
5000
1000
2000
3000
4000
5000
1000
2000
3000
4000
5003
Trainin Iteration
Training Beration
Training Iteraiio
Baseline
LeCun
Glorot
PowerLaw
Figure Training loss curves for forward PDE benchmarks under baseline, LeCun-numericee Glorot, and power-law (a 0.25, 3 1.75) initializations. Results are averaged over five seeds, with shaded regions indicating the standard emoor. Top row: "small" architecture (G = 5, two hidden layers with 8 neurons each). Bottom row: "large" architecture (G 20, three hidden layers with 32 neurons each).
To gain further insight into why the baseline and LeCun initializations underperform while Glorot and plw consistently succeed, we complement the loss-curve analysis with study NTK dynamics. Since the discrepancies between initialization strategies most pronounced in larger

---

Published as a conference paper at ICLR 2026
models, focus here the "large" architecture. For the power-law method keep 0.25 and B = 1.75, however we also provide the spectra for other configurations in Appendix F.2 Fig- ure shows the NTK eigenvalue spectra initialization, at intermediate training iterations, and at convergence, for the task of fitting (a y) (results for other targets are provided in Appendix F.3
Baseline
LeCun
Gloret
Power-Law
104
164
Eigmmaliss
10-1
104
155
10-
30
10a
30
Indices
Indices
Indices
Intiaaization
Intermediate Iterations
Final Iteration
Figure 3: Eigenvalue spectra of the NTK matrix at initialization (solid blue), intermediate iterations (dashed teal), and final iteration (dashed green) for function fitting benchmark fa y) under dif- ferent initialization strategies. Results correspond to the "large" architecture (G = 20, three hidden layers with 32 neurons each). The power-law initialization uses a = 0.25, 3 1.75.
The spectra reveal several notable differences. The baseline initialization exhibits a spectrum dom- inated by very large leading eigenvalues and steep decay, which collapses further during training, indicating poor conditioning and a low effective rank. LeCun behaves similarly but with even more extreme magnitudes. In contrast, the Glorot initialization produces well-spread spectrum with stable leading and trailing cigenvalues throughout optimization. The power-law scheme yields an even better spectrum, closely following a power-law decay initialization and remaining perfectly stable during training, suggesting balanced sensitivity across modes. These observations align with the performance trends reported so far: Glorot and power-law initializations induce stable, well- conditioned NTK spectra and therefore correspond to faster optimization and lower final error, while the baseline and LeCun produce highly skewed or collapsing spectra and thus consistently under- perform during training.
We also extend the NTK analysis PDE benchmarks, focusing on the Allen-Cahn equation as a representative case (results for Burgers and Helmholtz are provided in Appendix F3 To this end, we adopt the NTK formalism developed for PIML (Wang et al.) 2022) and adapt it to account for Residual-Basse Attention (RBA) weights (Anagnostopoulos et al. 22024) which are applied in the loss functions studied herein (see Appendix C22 for details). The resulting kernel is identical to the standard PINN NTK, except that incorporates the corresponding RBA weights (see Appendix [.]] for the full derivation). Figure Isows the NTK eigenvalue spectra separately for the PDE residual term (top row) and the boundary/initial conditions (bottom row).
The PDE residual spectra largely mirror the function fitting case: baseline and LeCun initializations yield poorly conditioned cigenvalues that collapse over training, while both Glorot and power-law maintain stability. The key difference lies in the boundary/initial condition terms, where Glorot shows some irregularities, though far less severe than the baseline and LeCun. The plwer-law scheme stands out as the most consistent, providing stable and well-structured spectra across both PDE residual- and boundaryfinitial-eem
4.3 FEYNMAN DATASET BENCHMARKS
As final benehmark, we turn to a subset of the Feynman dataset, restricted to dimensionless equa- tions. The explicit formulas of the target functions and implementation details are provided in Ap- pendix C.3 Having already established via empirical evidence and NTK analysis that the LeCun- based schemes are consistently underperforming, we benchmark here only the two competitive ini- tializations (Glorot and power-law) alongside the baseline. We again fix the power-law exponents to (a, 3) (0.25, 1.75) and evaluate all methods using the same "small" and "large" architectures defined in Section 4.2 Tables and report the results in terms of final training loss and relative L error with respect to the reference solutions, for the small and large settings, respectively.
8

---

Published as a conference paper at ICLR 2026
Eigenvaluss (MDE)
1000
12a4
10a
10
1244
10t
1311
10
10
160
207
122
201
13a
307
1.2
20'
107
337
1.0
101
Elganwatus
40
10-
Indices
Indices
inttalizatinn Intermediste iterations
Final Iteration
Figure 4: NTK eigenvalue spectra for the Allen-Cahn PDE benchmark under bascline, LeCun numerical, Glorot, and pawer-law (a 0.25, 3 1.75) initializations. Top row: spectra corre- sponding to the PDE residual term. Bottom row: spectra for the boundary/initial condition terms. Solid blue lines show the initialization, dashed teal lines show intermediate iterations, and dashed green lines show the final iteration. Results correspond to the "large" architecture (G 20, three hidden layers with 32 neurons each).
Table 3: Results on the Feynman benchmark for the "small" architecture (G = 5, two hidden layers with neurons each). Reported values correspond to median final training loss and relative L error with respect to the reference solution. Best results per equation are shown in bold. The power-law initialization uses a 0.25, 8 1.75.
Baseline
Glorat
Power- Law
Function
Loss
1.2
Loss
1.2
Loss
L
1.62
5.17 -3
9.86 10
4.22 10
1.88
4.14 10
1.6.2
3.58 10
4,28
8.69 10
5.01
10
1.97
4.37 10
1.12111
1.40
3.67
1.30
3.75
1.12
1.07
113.12
2.24 10a
1,86 10
3,51 10a
7.65 10
1.75 10
2.36 10
L16.6
2.62 10-4
3.55 10-2
2.94 10
3.63
10
1.19 10-4
2.92 10
118.4
1.39 10a
1.00 10
2.31 10a
1.00 10
1.04 10a
1.00 10)
1.2.2
5.00 10
7.21 10
1.40 10
1.19
10
9.99 10
3.13 10
127.6
1,87 10
1.00 10
1.24 10
1.00 10
1.77
1.00 10
129.16
1,05 10
1.14 10
1.22 10
1.24 10
3.14 10
6.83 10-3
L30.3
4,00 10-6
4.62 10
9,00 10
6.84 10
4.88 10
1.73 10-3
L40.1
10
4,76 10
3.90 10
10
1.74 10
1.81
L50.26
1.40 10-6
10
10
3.47
1.17 10
1.20
11...42
10
10
2.50 10
7.74 10
1.44 10
16
-0
10-3
1.16 10
4,97
10
1..1
10
-0
10-5
10
10
1..1
6.20
1.21
11.35.18
1.10 10
1.38
1.84
1.48
1136.38
3.50 10
1.17 10
6.50
1.57 10
2.71 10
3.43 10
11.10.19
1.40
10
3.14 10
1.50
2.90 10
8.26 10
2.24 10
11.13.37
2.60
1.05
5.30
1.1
4.35
4.37
The results confirm the same overall trends observed in the earlier benchmarks. both settings, power-law initialization achieves the best performance on the majority of equations, often by large margins in terms both final training loss and relative L error. Glorot initialization also provides substantial improvements over the baseline, particularly in the large architecture, where it consis- tently narrows the gap to pawer-law. A comparison between Tables and 4 further highlights the role of initialization: with Glorot and power-law. the richer architecture is able to drive the loss down by several orders magnitude and simultaneously reduce the L error, whereas under the baseline initialization, performance often degrades when moving from the small to the large setting.

---

Published as a conference paper at ICLR 2026
Table 4: Results the Feynman benchmark for the "large" architecture (G 20, three hidden layers with 32 neurons each). Reported values correspond to median final training loss and relative L error with respect the reference solution. Best results per equation are shown in bold. The power-law initialization uses 0.25, 3 1.75.
Function
Loss
L
Loss
L
Loss
L
16.2
1,09
1.51 10
4,80
4.19
5.20 10
3,85 10
16.2b
1.36 10-a
1.64 10
7.60 10-3
5.80 10
2.18 10-6
4.59
10
112.11
1.64 10-4
3.77 10
3,00 10-6
1.77 10
2.16 10
1.66 10
L13.12
2.70- 10
3.08 10
10
10
2.53
5,49-10"
L16.6
1.63 10
10
10-6
L488 10
L18.4
10a
1.00-10" 10
10
10"
10 10-2
100
10
10-5
10-3
10-7
10-3
127.6
10
1074
10"
10
129.16
2.01
2.06
130.3
1,18 10
7.72
10
2.92 10
2.17
4.17
140.1
2.26 10
6,70
5,00 10
3.39 10
1.41
6.17 10
1.50.26
2.03
4,38
2.00
1.50
3.70
2.25
11.2.42
1.52
6,86
4,00
2.62
8.54 10
4,98
4
11.. 15a
6.60 10
7.60 10
2,00 10
5.47 10
8,13 10
4,40 10
11117
1.75 10
9.78 10
1.10 10-3
1.01
10
1,80 10
3.00 10-3
1L11227
8.80 10
1.04 10
1,00 10
3.76 10
1.54 10 10-7
1.95 10
11335.1
7.40 10
1.19 10
6,00 10
1,18 10
2.95 10
7.77 10
11.6638
1,93
9.48 10
8.00
1.11 10
3 10
4,92 10
11.10.19
1.81 10
2.74 10
1.00 10-6
9.89 10-4
9.87 10-2
8.70 10
11.77337
1.45 10
9.10 10
4,90 10-5
1.31 10
6.45 10-6
5.14 10
CONCLUSION
In this work, we proposed and systematically evaluated new initialization strategies for spline-based KANs. Specifically, we introduced two theory-driven schemes inspired by LeCun and Glorot, in- cluding a variant with normalized basis functions, as well as an empirical family of power-law initializations. Through large-scale grid searches, identified favorable exponent ranges for the power-law method for both function fitting and forward PDE problems. Across all evaluations, in- cluding loss curve analysis and NTK dynamics, we showed that initialization plays critical role in KAN performance. In particular, our results demonstrate that while LeCun-inspire schemes of- fer limited benefits, Glorot-inspired initialization emerges as a strong candidate for parameter- rich architectures and the empirical power-law family provides the most robust improvements overall, achieving faster convergence and lower errors across benchmarks. These findings highlight ini- tialization key component of training KANs and identify effective practical strategies for the process.
5.1 LIMITATIONS AND FUTURE WORK
While study establishes the importance of initialization in spline-based KANs, it also comes with limitations. Our power-law scheme, although empirically effective, currently lacks a rigor- ous theoretical foundation, and understanding why specific exponent ranges perform well remains an open question. Moreover, although we considered both supervised function fitting and physics- informed PDE benchmarks, further exploration in other domains such as reinforcement learning generative modeling could provide additional insights. Addressing these limitations offers natural directions for future work, including deriving principled theory for power-law initialization, inves- tigating transfcrability across KAN variants (c.g.. Chebyshev- based or residual-free forms), and exploring initialization strategies in conjunction with adaptive optimization techniques.
REPRODUCIBILITY STATEMENT
The full code (including selected seeds for each experiment) and the processed data from the grid-search experiments are publicly available at the following GitHub repository: https: / / thub. sri gotasb/g////rr////aa/ /KAN Initialization Schemes.
10

---

Published as a conference paper at ICLR 2026
ACKNOWLEDGMENTS
S. R. and G. A. supported by the Innovative Health Initiative Joint Undertaking (IHI JU) under grant agreement No 101253520. The JU receives support from the European Union's Horizon Eu- rope research and innovation programme and COCIR, EFPIA, Europa Bío, MedTech Europe, and Vaceines Europe.
LLM USAGE
Large Language Models (LLMs) were used during peer review for grammar and syntax refinement only; all ideas, technical content, analyses and conclusions remain the authors' work.
REFERENCES
Diab W. Abueidda, Panos Pantidis, and Mostafa E. Mobasher. DeepOKAN: Deep operator net-
work based on Kolmogorov Arnold networks for mechanics problems. Comput. Methods Appl. Mech. Eng., 436: 17699, 2025. doi: hanp::/oo..ogggoo.........m..mm.
URL https:
/ 7/www. www.sciieeess. rect. wtt.///ccccccce.c http://ccc..cccccooccccoooooa/aa
Tal Alter, Raz Lapid, and Moshe Sipper. On the robustness of Kolmogorov-Amold networks:
An adversarial perspective. Transactions on Machine Learning Research, 2025. URL https:
//openzeview
Sokratis J. Anagnostopoulos, Juan Diego Toscano, Nikolaos Stergiopulos, and George Em Karni- adakis. Residual- based attention in physics- informeed neural networks. Comput. Methods Appl. Mech. Eng., 421: 116805 2024. doi: 111016/.ccm..2000 16805.
James Bradbury, Roy Frostig, Peter Hawkins, Matthew James Johnson, Chris Leary, Dougal Maclaurin, George Necula, Adam Paszke, Jake VanderPlas, Skye Wanderman- and Qiao Zhang. JAX: composable transformations of Python+NumPy programs, 2018. URL http:
//github. com jaapmm/jaammmmaaaaa
Xavier Glorot and Yoshua Bengio. Understanding the difficulty of training deep feedforward neural networks. In Proceedings of the Thirteenth International Conference Artificial Intelligence and Statistics, volume 9, pp. 249 256, 2010. URL httttts::/ppcceeednittt mlr. press/ v9 09/glorot: (a. htm
Leonardo Ferreira Guilhoto and Paris Perdikaris. Deep learning alternatives of the Kolmogorov superposition theorem. The Thirteenth International Conference on Learning Representations, 2025. URL https:://peeneeeee...ttt wtwpeerwwwwwwweeeeeeeeeeeeeeee http://eet.tt.ottt.t..m/mm/m
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Delving deep into rectifiers: Surpassing human-level performance on ImageNet classification. In 2015 IEEE International Conference on Computer Vision (ICCV), pp. 1026-1034, 2015. doi: 10. 00.1109//CCV.2015.123
Amanda A. Howard, Bruno Jacob, Sarah H. Murphy, Alexander Heinlein, and Panos Stinis. Fi- nite basis Kolmogorov-Amoo networks: domain decomposition for data-driven and physics- informed problems, 2024. htt.htttt////xlle.e
abprg//6bb/224.... 9662
Xiao Shi Huang, Felipe Perez, Jimmy Ba, and Maksims Volkovs. Improving transformer optimiza- tion through better initialization. In Proceedings of the 37th International Conference on Machine Learning, volume 119, pp. 4483-4483, 2020. URL httth:://troceeeeeeet htt.pg/........../...
.ml.
v 2119/huang200 h
Arthur Jacot, Franck Gabriel, and Clement Hongler. Neural Tangent Kernel: Convergence and generalization in neural networks. In Advances in Neural Information Processing Systems, vol- ume 31, 2018. URL ott.r//ooooooooor///oooooooooooo cappae fawerrffeeeeepppppee
neur

---

Published as a conference paper at ICLR 2026
Ali Kashefi. Kolmogorov-Amold PointNet: Deep learning for prediction of fluid fields on irregular geometries. Comput. Methods Appl. Mech. Eng., 439:117888, 2025. doi: https//ddi..org/11.0016
j.cma.2025.1 117888. URL https: wwtp///www.... www.sacciecedireec.oo www./aciens///tttii.....
p p11/80045782225001 6044
Andrey Kolmogorov. On the representation of continuous functions of several variables by superpo- sition of continuous functions one variable and addition. Doklady Akademii Nauk SSSR, 114:
369-373, 1957.
Akash Kundu, Aritra Sarkar, and Abhishek Sadhu. KANQAS: Kolmogorov-Aroll network for quantum architecture search. EPJ Quantum Technol., 1 11:76, 2024. doi: https://doi.org/10.. 1140/ 140/ c2ig//s4005070024-00022----
Yann LeCun, Leon Bottou, Genevieve B. OIT, and Klaus-Robert Müller. Efficient backprop. In Genevieve B. OIT and Klaus-Robert Müller (eds.), Neural Networks: Tricks of the Trade, pp. 9-50. Springer, 1998. ISBN 978-3-540-494300 doi: 10. 10.100773-540-444300 2.
Jin Lee, Ziming Liu, Xinling Yu, Yixuan Wang, Haewon Jeong, Murphy Yuezhen Niu, and Zheng Zhang. KANO: Kolmogorov-Amool neural operator, 2025. URL. https:://axiis org/
abs 2os/2509. 16825
Ziming Liu, Pingchuan Ma, Yixuan Wang, Wojeiech Matusik, and Max Tegmark. Kan 2.0:
Kolmogorov-Amold networks meet science, 2024. URL htttps://rxxiv.orge
2408. 10205
Ziming Liu, Yixuan Wang, Sachin Vaidya, Fabian Ruchle, James Halverson, Marin Soljacic, Thomas Hou, and Max Tegmark. KAN: Kolmogorov-Arooe networks. The Thirteenth International Conference Learning Representations, 2025. URL/ttts:: URLLatps::/pppneeeeee /oopnnevvew. net ftte/forum:id=0z00qp000
Dmytro Mishkin and Jiri Matas. All you need is a good init. In 4th International Conference on Learning Representations, ICLR 2016, 2016. URL
06422
Eleonora Poeta, Flavio Giobergia, Eliana Pastor, Tania Cerquitelli, and Elena Baralis. bench- marking study of Kolmogorov-Amole networks on tabular data. In 2024 IEEE 18th Interna- tional Conference on Application of Information and Communication Technologies (AICT), pp- 1000///ACCTT1888222.00477744
1-6, 2024. doi:
Maziar Raissi, Paris Perdikaris, and George Em Karniadakis. Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. Comput. Phys., 378:686-7001 2019. doi: https///bo.orr//1.0000/////
2018. 10.045. URL https:/ hwtp///www... www.ssscecnceiiirce...t.... con////ceenne///riiill//..... 500219991 8307125
Spyros Rigas and Michalis Papachristou. jaxKAN: A unified JAX framework for Kolmogorov- Arnold networks. Journal of Open Source Sofwware 10(108):7830, 2025. doi: 10.21 105/joss. 07830. URL
211.5 joss. 07830
Spyros Rigas, Michalis Papachristou, Theofilos Papadopoulos, Fotios Anagnostopoulos, and Geor- gios Alexandridis. Adaptive training of grid-dependent physics-informed Kolmogorov-Arool 200.CEEESSS.0223349922222
networks. IEEE Access, 12: 129172922-776998 2024. doi: 10.
Spyros Rigas, Fotios Anagnostopoulos, Michalis Papachristou, and Georgios Alexandridis. Training deep physics- innfommedd Kolmogoroy-Amold networks. Comput. Methods Appl. Mech. Eng., 452: 118761, 2026. doi: 10. 2.016//.cma.2266... 18866.
Khemraj Shukla, Juan Diego Toscano, Zhicheng Wang, Zongren Zou, and George Em Karniadakis. A comprehensive and FAIR comparison between MLP and KAN representations for differential equations and operator networks. Comput. Methods Appl. Mech. Eng., 431:117290,, 2024. doi: http:://bo.oay///l.lo.oo...ammap URL https: www. sciencedirect. com/
12

---

Published as a conference paper at ICLR 2026
Maciej Skorski, Alessandro Temperoni, and Martin Theobald. Revisiting weight initialization of deep neural networks. In Proceedings of The 13th Asian Conference on Machine Learning, volume 157, pp. 1192-1207, 2021. URL hktps://proceedlinn /
skorski21a ht.m
Silviu-Marian Udrescu and Max Tegmark. AI Feynman: A physics-inspired method for symbolic regression. Science Advances, 6116))aayy6631 2020. doi: 10.1 (0.1116//siadd..aay26333
Sifan Wang, Xinling Yu, and Paris Perdikaris. When and why pinns fail train: A neural tangent doi: 10.
kernel perspective. Journal of Computational Physics, 449:110768,, 2022. 2021.1 10768.
Sifan Wang, Bowen Li, Yuhan Chen, and Paris Perdikaris. Piratenets: Physics-informed deep learn- ing with residual adaptive networks. Journal of Machine Learning Research, 25 2024.
Yixuan Wang, Jonathan W. Siegel, Ziming Liu, and Thomas Y. Hou. On the expressiveness and spectral bias of KANs. In The Thirteenth International Conference on Learning Representations, 2025a. URL https: Dtttps://operrevieee wttp://wwwwweeeeeooeeeeeeeeeeee "oorummiddyyllDUuuuGm
Yizheng Wang, Jia Sun, Jinshuai Bai, Cosmin Anitescu, Mohammad Sadegh Eshaghi, Xiaoying Zhuang. Timon Rabczuk, and Yinghua Liu. Battogooro-Amooi-iinfoooe
neural network: A physics-informed deep learning framework for solving forward and inverse problems based on Kolmogorov-AAmold networks. Comput. Methods Appl. Mech. Eng., 433:1 2025b. doi: https://bbboog//0...0////..mmmmp http://w///aattea/aa///aaaaaaaaa Scieene//rrttceee/1//05557330000
URL https: /www wwtp://wwienndddeee..... com/
Ge Yang, Edward Hu, Igor Babuschkin, Szymon Sidor, Xiaodong Liu, David Farhi, Nick Ry- der, Jakub Pachocki, Weizhu Chen, and Jianfeng Gao. Tuning large neural networks via zero- shot hyperparameter transfer. M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and J. Wortman Vaughan (eds.), Advances in Neural Information Processing Systems, volume 34, pp. 17084-17097, 2021. URL htttps::/proceeddi newrseisssennnnenrsss......... wttp:/rrppp.ccpppprrreeeele paper/2021/11 e / httta///ccccccaaaaaa..ee....a./.
Rumpeng Yu, Weihao Yu, and Xinchao Wang. KAN MLP: A fairer comparison, 2024. URL
Xianyang Zhang and Huijuan Zhou. Generalization bounds and model complexity for Kolmogorov- Amool networks. In The Thirteenth International Conference on Learning Representations, 2025. URL https::/ net /eteeffoo forum? 1d=g52MyAUUGG
13

---

Published as a conference paper at ICLR 2026
A DERIVATION OF LECUN-INSPIRED INITIALIZATION SCHEME
In this appendix, we provide a derivation of Eqs. from the main text. Assuming statistical independence between each term in the outer sum Eq. and requiring the output variance to match the input variance, one finds
G+
Var (x:) 7inn Var Túi R (2;) + cjs Σ
bjmm Bm (aii)
where the right-hand side contains the variance of a sum of G +++ terms: one residual term and G++ k spline basis terms. We adopt simplifying assumption that the total variance is approximately equipartitioned across all components SaaennenttSalwwwwe Flloowng us bypass pairwise covariance terms. This leads to the following expressions for the residual and spline basis terms, respectively;
Var (2:4
Var (2:)
Tin Var (rye R (2+)),
G+k+1 ninn Var [byem Bm (a4))
c+k+1
Since the trainable weights rys are independent of the residual function R (zt). the variance of their product becomes
Var [r, R (xx)) = E² (rys) Var [R (x;)] E² [ (2x)) Var (1ys) + Var (ryx) Var [R (x;))
Var (rzs) Var (xx)) E² [R()) of E [1e (a:)]
(10)
and, in completely analogous manner, we find
Var (ijm Bm (2x4) of [B (x:)]
(11)
Substitution of the expressions of Eqs. 10
into Eqs.
yields Eq.
from Section 3.1
'This assumption does not necessarily hold in general. For example, one could consider a 50%-509, split between the residual and hasis function terms. We experimented with this alternative and found that it yielded poorer results compared the variance partitioning that leads Eqs. 2
14

---

Published as a conference paper at ICLR 2026
B DERIVATION OF GLOROT-INSPIRED INITIALIZATION SCHEME
In this appendix, we derive Eqs. 5 from the main text. Unlike the LeCun-inspired scheme, which focuses solely on variance preservation in the forward pass, the Glorot principle (Glorot & Bengio, 2010) requires that the variance of both activations and backpropagated gradients remain constant across layers. For analytical tractability. and following the standard assumption this setting, further approximate these constant values by unity.
Var(æ;;) = Var(y:) = 1, Var(dz;; = Var(óy;) = 1,
(12)
where dz, = aC/Oz; and ōyj OC/Oyj, with L denoting the loss function. This approximation consistent with the common assumption of i.idd inputs with zero mean and unit variance. In practice, when this assumption does not hold, an additional gain factor can be introduced to rescale the initialization,
standard in frameworks such as PyTorch.
Using the result from Appendix A together with the first condition of Eq. the constraints for
variance preservation in the forward pass can be written as
= (G + k + 1) Thnn of PR
(0)
I= Rin of PB
(0) (13)
where PR { and PB (0)
= B (x:)²² as defined Eq.
For the backward pass, differentiating Eq. with respect to z; gives
θy: = rs R' (x;) + Cjt
G+
Σ
bjam Byn (x,);
(14)
dr,
m=1
Setting Czi = the chain rule yields
Oyi 8y3
G4+bk Beat
dx;
Ta Rt (x;)
(aa))og dys
jjim By (2x gys
(15)
Do'
residual contribation
m-th spline contribution
and applying the second condition of Eq.
gives
G+k
= Nout Var Ty R'(;) + Σ byrm By, (xx)
(16)
m=1
where we have adopted the standard Glorot assumptions: the dy, are zero-mean, mutually inde- pendent, and independent of weights and inputs. At this point may again equipartition the total variance across the (G+++1) components (one residual term and G+ spline terms), exactly mir- roring the forward-pass treatment. This leads to
1 = (G + + 1) Nout Var [rys R' (x:))s
1= (G + k 1) Nout Var (Dyim B'r (xx)),
(17)
and, following the same arguments in Appendix A
find
(G l: 1) ouu o P
(1)
1 1 (G Nout nour of PB
(1)
(18)
(1) (xx,) and f (1) where P
E[B", (z:)²1 defined in Eq. 6
Equations 13 and (77) are the forward- and backward-pass constraints, respectively. Balancing
them in the Glorot manner (i.e., by harmonic averaging) yields the standard deviations in Eq. ( of
15

---

Published as a conference paper at ICLR 2026
the main text. As sanity check, consider an MLP: the residual term is absent, and the linear layer followed by a nonlincarity can be viewed as a single basis function. For the common hyperbolic langent activation, (0) (1)
PB PB R 1 (Gloro & Bengio, 2010), so our scheme reduces to
2
2
(19)
G
(0) B
+
Nout
PB
1iin Flout
which recovers the classical Glorot initialization.
16

---

Published as a conference paper at ICLR 2026
C IMPLEMENTATION DETAILS
This appendix provides the full specifications of the benchmarks used in our experiments, including the functional forms of the target problems, training setups, and data generation procedures. We separate the discussion into three parts: function fitting, forward PDE problems and the Feynman dataset.
C.I FUNCTION FITTING
For the function fitting experiments of Section 4.1 and Section 4.2 we study five wwo-dimensional functions ranging from simple expressions to more complex, nonlinear, or piecewise-de delined forms.
Specifically, we consider the following functions the |-1, x [-1] 1] domain:
f (x. = xy
f2 (z.v) = exp (sin(Tra:) + y²)
fa (x..) = I1 (x) + exp exp bl) I1 (3) + sin (zy)
fa v) -s s(fs (z, erf (v) c fs (x, v) erf (v)
Js [s v) = s0. - x) + erf(s) min ryy. my
where T1 (x)) is the modified Bessel function of first order, sgn (ar) is the sign function, erf (x) is the emor function and S (a:), C (x) are the Fresnel integral functions defined as
/ sin
= fo
(x) J
dl,
C (xx)
cos
dt.
(20)
The reference surfaces for these functions are shown in Figure
f1(x, v)
fx(x
1.0
1.0
1.0
1.0
2.5
2.0
1.5
0.0
0.0
0.0
0.0
1.0
2
0.5
0.0
-1.
-1.0
-1.0
-1.0
-1.0
-0.5
0.0
1.0
-1.0
-0.5
0.0
0.5
1.0
-1.0
-0.5
0.0
1.0
falx. y)
fg
1.0
1.0
0.4
0.5
0.5
0.3
0.0
02
0.0
0.1
-1.0
0.0
-1.0
Figure 5: Reference surfaces for the five two-dimensional target functions fi through I's used in the function fitting experiments.
The KAN models used to fit these functions utilize spline basis functions of order & 3. defined over an augmented, uniform grid within the [-1, 1] domain (Liu et al. 2025). Training is performed using the Adam optimizer with fixed learning rate of 10 33 with the objective of minimizing the mean squared error between the predicted and reference function values. For each target function fa (xr, y), with i 1,...; generate 4,000 random input samples uniformly distributed over the
17

---

Published as a conference paper at ICLR 2026
60x
Burgers
Uix
Helmholtz
utx, V)
1.0
1.0
1.0
0.5
0.5
0.5
0.5
0.5
0.5
0.0
0.0
0.0
0.0
0.0
-0.5
00.5
-1.8.00
0.25
0.50
0.75
1.00
-..00
0.25
0.50
0.75
1.00
1.0
-0.5
0.5
Figure 6: Reference solutions for the three PDE problems considered.
domain [-1, x [],,]], and calculate the corresponding outputs to serve as ground truth during training. compute the relative L error between the model predictions and reference solutions, we evaluate all trained models on a uniform 200 200 grid covering the same domain.
C.2 FORWARD PDE PROBLEMS
In addition to function fitting, we consider three representative forward PDEs commonly used as PIML benchmarks. For each case, we specify the governing equation, domain and boundary/initial conditions.
Allen-Cahi equation.
We solve the Allen-Cahn equation on (f,x) E [0, x [-1,1]:
(1, - D (1, - c(u(l,,)) u((,xx))) 0,
(21)
with diffusion coefficient = 10 and reaction strength c = The initial and boundary condi- tions
u(0,x)) x² cos((Tx),
x rE(-1,1]
(22)
u(t, 1) cut, 1) 11,
te [0, 1].
(23)
Since the Allen Cahn equation has no analytic closed-form solution, we use the reference solution used in Wang et al. (2024), which is depicted in the left plot of Fig.
Burgers' equation.
We solve the viscous Burgers' equation on (t, E [0, [-1,1]:
(t, + u(t, x) ue (t, - Unx (t, = 0,
(24)
for v/ = 0.01/T, with initial and boundary conditions
x(0, -
(-1,111),
(25)
u(t,-1 1(t,
t [0, 1].
(26)
Similar to the Alhen-Cahn equation, Burger's equation has analytic closed-fome solution, there- fore we use the reference solution used in Wang et al. (2024), which depicted the middle plot of Fig. 6
18

---

Published as a conference paper at ICLR 2026
Helmholtz equation. We solve two-dimensional Helmholtz problem on (z. [-11, 1]2 with unit wavenumber and a separable sinusoidal source:
(2r y) + (2, y) + u(x,y) = f(x.8). y).
(27)
where
fz, v)
xa +2)) sin(Taa sin(Ta
(28)
and 01 = and 02 4. We consider homogeneous Dirichlet boundary conditions:
u(x, v)
for (x, E ((-1,113).
(29)
The analytic solution to this PDE problem is
Ure (I,, v) = sin (πx) sin (4my),
(30)
and is depicted in the right plot of Fig. ooor z, sampled on a uniform 512x512 grid.
The PDE problems are solved using the Residual-Bass -asse Attention (RBA) weighting scheme (Anag- nostopoulos et al. 2024) within the PIML framework (Raissi al., 2019), where the training objec-
tive is defined as sum of weighted residuals associated with the PDE differential operator and the boundary/initial condition operators. Specifically, we minimize
No.
(0)
(pde) (pdo) a
y,
(0)
(bc)
(bc) r
(0)
(31)
N
pdc
Noc
(pde) where |1-1/2 denotes the 1,2 norm. Here, r represents the residual of the governing PDE evaluated (bc) denotes the residual of the boundary or initial condition (both
at the i-th collocation point, while :
are included in the second summation). The weights a (E) (E E (pde, bc}) are initialized to 1 and updated after each training iteration according to
(((()
(E), (new) a,
(E), (old)
-
(32)
Ng max, (1r 1
with hyperparameters = 0.999 and = 0.01. This formulation ensures that collocation points with larger relative residuals assigned greater importance during optimization
We minimize the loss function Eq. (31) using the Adam optimizer with a fixed learning rate of 10-3, operating in full-batch mode. For each PDE, we sample Npdc 2 collocation points uniformly from 26 28 grid, while for boundary and initial conditions points per condition, sampled uniformly along the corresponding axis. The spline basis functions defined as in the function fitting case (see Appendix C.1
2 collocation
C.3 FEYNMAN DATASET
As a third benchmark, we consider the subset of the Feynman dataset used in Section 4.3 The
implementation details are identical those of the function fitting benchmarks in Appendix C.I with the exception of sampling. In this case, we generate 4,000 random input samples uniformly distributed over the domain (-1,0) (0, 1), explicitly excluding the points 11, 0, and 1 to avoid singularities in certain formulas.
To compute the relative L error between model predictions and reference solutions, we evaluate all trained models uniform 200 x 200 grid for two-dimensional functions and a uniform 30 30 x
Withouut RBA, the models trained to solve the Allen-Cahn equation would yield highly inaccurate solu- tions, preventing meaningful comparison of initialization schemes.
19

---

Published as a conference paper at ICLR 2026
30 grid for three-dimensional functions. Table 5 lists the indices of the functions included in this benchmark, together with their explicit expressions for reference.
Table 5: Dimensionless formulas from the Feynman dataset used the benchmark. Each entry shows the dataset index and the corresponding explicit expression.
Index
Formula
I6.2
fi (x1, = exp
(-2-7 (2T=²) 1/22
225
L6.2b
22, exp
(111-22)
(2m233 1/22
203
L12.11
[2( = 1 (12()
L13.12
f4(21. = (1/2² 1)
L16.6
f6(21,02) = (I1 + /22/ + x1
L18.4
x2) = (1 + x .2 / (1 + x1)
126.2
fr(21,02) = aarssin و sin(
127.6
f6(21, 22) = 1/(1 + 21
129.16
02, x2) vI+ + 2x cos x3)
L.003
f10(21,22) sin sin (x122/ sin (x2/2)
L40.1
(11( mn exp(-12)
L50.26
(2((11)2)) cos(æ1) cos? (I1)
11.2.42
f1s(11,22) (æn 1)12
I1.6.55
f1l21, ma) = = /I² + 23
I1.11.7
f1s(11, 22, x3) = x1 (1 + x2 cos(ca))
I1.11.27
f16(51, = (2)02) (1 -
11.35.18
fir(11, =1/(exp(02) + exp( -2))
II.36.38
(2, =1 1233
10.10..1
f19(11,22) /1 . + x
III.1..37
22,1 43) = x2 (1 æ (o3)
20

---

Published as a conference paper at ICLR 2026
D INDICATIVE RESULTS FOR POWER-LAW GRID-SEARCH
To illustrate the performance landscape of the power-law initialization, we present heatmaps over (a, (3) configurations for representative cases. Specifically, Figures 7 show results for the func- tion fa y) across the four grid sizes, while Figures I IT-L3 show the corresponding results for the Allen- Cahn PDE. In each heatmap. the horizontal axis corresponds to and the vertical axis S, with rows and columns indicating different network widths and depths, respectively. These visu- alizations highlight the regions where power-law initialization provides the greatest improvements, and help motivate the choice of (a, 3) = (0.25, 1.75) used for the architectures studied in Sections [42 and 4.33 Complete heatmaps for all benchmarks included the supplementary material (see Reproducibility Statement).
21

---

Published as a conference paper at ICLR 2026
Depth
Depth
Depth
Depth
2.0
10-1
10-4
1.5
10
100
10
10
10
0.5
10-3
10
0.0
2.0
10-4
1.5
10-2
10
10-
= 1PPM
1.0
10
10-3
0.5
10
10-
10
10
0.0
10-5
10-4
10
2.0
10-4
10
1.5
10-
10-1
HPM
1.0
10
10-3
20-1
10
10
0.5
0.0
10-4
10
10-5
2.0
10-1
1.5
10-1
10
1.8
10
20
0.5
10-
10-1
10
10-3
10-4
10-4
0.0
2.0
10-2
101
10-2
10-1
1.5
10-3
10-
10-1
28 MER W
1.0
10-
10
10
0.5
10-1
10-5
10
00
10-4
2.0
10
10-1
10
10-1
15
10-2
10-3
10-
10-4
10
10
10
69
10-3
10-9
10
10-5
10-1
1.0
15
0.0
1.5
0.0
0.5
1.0
2.0
1.5
Figure 7: Grid search for the power-law initialization applied fit function fs (z, for G 5. Each heatmap corresponds to an architecture, with the horizontal and vertical axis representing and 8, respectively, and color denoting final training loss.
22

---

Published as a conference paper at ICLR 2026
Depth
Depth
Depth
Depth
2.0
10
10
10-1
1041
1.5
10
104
1.
10
10
0.5
10-3
10
104
0.
10-4
2.0
10-1
10
10-1
1041
1.5
10a4
1.
10a4
10
10-3
10
10
0.5
10
10-5
10-3
3044
0.0
2.0
10-1
10-1
1.5
10-2
10-2
10-
1.0
10-
10
10
1044
10
10-
0.5
10-5
10-5
104
00
10
10
10-4
2.8
10-1
10-2
10-2
10
1.5
10-2
10
10-3
10-
3.0
10-
10-4
10-
10-
104
10-
0.5
10
10-5
10
10
0.0
10-5
10
10
2.0
10-1
10-2
10-2
10-
101
1.0
1.5
10-4
10-4
10
10-2
10
10
0.5
10-
10
10
10
10-4
0.0
10-7
2.0
10
10-1
1.5
10-2
10-3
10"
10-4
10
10
1.
10-4
10-4
0.5
10
10
10-4
0.0
10
10-4
10-6
10-7
10-7
155
1.5
2.0
0.0
0.5
1.0
2.0
1.5
2.0
Figure 8: Grid search for the power-law initialization applied to fit function fs (z. y) for G 10. Each heatmap corresponds to an architecture, with the horizontal and vertical axis representing a and 8, respectively, and color denoting final training loss.
23

---

Published as a conference paper at ICLR 2026
Depth
Depth
Depth
Depth
2.0
10-3
10
10
1.5
10
10-2
104
10
10
10
6.5
10
10
0.
10-4
10-3
10
2.0
10-
10-1
1.5
10a2
10
10
10-3
1.0
10
10-4
10-4
10
0.5
10-4
10-5
10-5
10-1
0.0
10-4
10-
10-4
2.0
10
10-1
10-1
1.5
10-2
1044
10
10-1
10
1.0
10
10
10-4
05
10-5
10
0.0
10
10-
104
2.0
10-1
10a1
10"
10-2
10-
20
1.5
10
10-3
10-1
1.0
10-5
10-4
0.5
10-4
10-1
0.0
20-7
2.0
10-4
10
102
10-1
100
1.5
10-3
10'
10-2
26 19PM
1.0
10
10
10-1
0.5
10-5
0.0
10-1
10-7
10-4
10-4
2.0
10-1
10-2
10a
10-2
15
10-1
10-3
10-2
10'
69 MPP
10
10-4
10
10-5
10-t
10-4
10-1
10
1.0
155
1.5
2.0
0.0
0.5
1.0
2.0
1.5
Figure 9: Grid search for the power-law initialization applied to fit function fs (z. y) for G 20. Each heatmap corresponds to an architecture, with the horizontal and vertical axis representing o and 8, respectively, and color denoting final training loss.
24

---

Published as a conference paper at ICLR 2026
Depth
Depth
Depth
2.0
1044
10:
10-1
1.5
104
10
10
10
6.5
10
10
10
0.0
10
10-3
10
2.0
10-1
10-1
1.5
10
10
10
1PPM
10
10-3
20
1.0
10
0.5
10
10-5
10-3
1044
0.0
10
1044
2.0
10-1
10
10
1.5
20
10
10
10
1.0
10-4
10-
10
0.5
10-4
10
6a
10-t
10-4
2.0
10-1
10
100
10
1.5
10-1
10
10
20 M9PM
1.
10
10-
10-2
10
0.5
10
10-
0.0
10
10-1
2.0
10
10-2
10-
10a
1.5
10
10-
ZE VAaA
10-4
1.0
10
10-1
10
0.5
10-
0.8
10-7
10-7
10-4
10-5
2.0
100
10-t
10*
1.5
10-5
10
10
100
10
10-4
10-4
10
10-7
10-9
10-4
10-4
1.0
1.5
0.0
1.5
2.0
0.0
0.5
1.0
2.0
1.5
Figure 10: Grid search for the pawer-law initialization applied fit function fs (æc, y) for G 40. Each heatmap corresponds to an architecture, with the horizontal and vertical axis representing and 8, respectively, and color denoting final training loss.
25

---

Published as a conference paper at ICLR 2026
Depth
Depth
Depth
Depth
2.0
10
10
1.5
20c
10a
10
10
10
0.5
10-1
20-1
0.0
2.0
101
1.5
100
102
1.0
201
0.5
10-4
10-1
0.0
2.0
101
10
1.5
20'
8 Midii
1.8
100
20*
10
0.5
10-3
10
10-2
0.0
2.8
10*
10a
10
1.5
10r
97 19PM
1.0
100
10:
10
20'
10P
10
0.5
10-1
10
10-2
10-2
0.8
2.0
104
2011
10a
1.5
104
309
10a
10
25 LADWA
10F
1.
100
10
10a
0.5
10
0.0
10-2
10-1
2.0
10*
10a
2012
102
1.5
10a
101
10
1.0
109
10'
10a
10*
p9
05
10a
10:4
20a
1.0
15
1.5
2.0
0.0
0.5
1.0
2.0
1.5
Figure 11: Grid search for the power-law initialization applied for the solution of the Allen-Cahn equation for G = 5. Each heatmap corresponds to an architecture, with the horizontal and vertical axis representing and 3, respectively, and color denoting final training loss.
26

---

Published as a conference paper at ICLR 2026
Depth
Depth
Depth
Depth
2.0
10
1.5
10
10c
UPPM
10
1.0
10t
0.5
20
10
20-1
0.0
2.0
101
10
10*
1.5
10
100
10a
1.0
10
0.5
10-2
10-1
0.0
2.0
10a
10*
104
1.5
10:
102
101
= MAPN
1.0
10c
10
10
10a
20
0.5
10-2
10-2
10
6a
2.0
20
1.5
10a
101
102
10*
1.
10
102
0.5
10-1
10a
0.0
10-2
2.0
104
10a
1.5
104
10*
1.
10
16F
10a
6.5
10-1
10-1
0.0
10-1
2.0
10c
10a
1.5
10a
10
10a
2012
10a
20a
200
1.0
1o4
10*
05
100
10
10*
p9
10-3
10
10-1
20a
1.0
15
1.5
2.0
0.0
0.5
1.0
2.0
1.5
Figure 12: Grid search for the power-law initialization applied for the solution of the Allen-Cahn equation for G = 10. Each heatmap corresponds to an architecture, with the horizontal and vertical axis representing and 3, respectively, and color denoting final training loss.
27

---

Published as a conference paper at ICLR 2026
Depth
Depth
Depth
Depth
2.0
101
20le
1.5
10
10
104
104
1.0
10*
0.5
10-4
10-1
10
20<
109
0.0
10
2.0
101
10
10:
1.5
10
10*
20*
10
19PM
1.0
10a
10*
10-1
0.5
10"
10*
10-4
20-1
0.0
10
2.0
10a
204
10a
1.5
16-
1.0
10c
10a
10*
10-1
10-1
10
0.5
10
00
10-3
10-2
10-1
2.0
10a
10*
10a
10a
1017
1.5
10
10a
30*
9t 100
1.0
10
102
104
10-1
0.5
20*
10-2
10-1
10P
00
10-3
2.0
10ah
10
2014
1.5
101
10*
10'
104
1PM
1.0
10a
20
10*
10-1
0.5
102
10-2
10
10-2
0.0
10-4
2.0
10c
1032
2017
10a
15
10'
10a
100
104
107
10F
109
10
10º
10-1
10a
10-2
10
10-1
1.0 155
2.0
1.5
2.0
0.0
0.5
1.0
2.0
1.0
1.5
Figure 13: Grid search for the power-law initialization applied for the solution of the Allen-Cahn equation for G = 20. Each heatmap corresponds to an architecture, with the horizontal and vertical axis representing and 3, respectively, and color denoting final training loss.
28

---

Published as a conference paper at ICLR 2026
E TRAINING CURVES WITH LEARNING RATE SCHEDULING
In Section 4.2 of the main text, the training curves shown in Figures 1 and were obtained using a fixed learning rate in order isolate the effect of initialization, as initialization and learning- rate adaptability are known to interact (c.g.. (Yang al., 2021)). However, the fixed learning rate induces oscillations the loss curves, particularly for the Glorot and power-law schemes in the larger architectures. To verify that these oscillations purely an artifact of the constant learning rate, we repeat the same training experiments using a learning-rate scheduler.
Figure 14 shows the results for the function fitting benchmarks, where all settings identical to those in the main text except for the use a learning-rate scheduler: training begins with a learning rate of 10 followed by exponential decay with decay factor 0.9 every 50 iterations. Similarly, Figure 15 shows the results for the PDE benchmarks, again using all the same hyperparameters as
fd0.y)
Mxx.y)
4220 2 gpdp S
10
Training Loss
30"
10-2
10-
10-1
10-
10-3
10a4
10-4
10-4
10-4
10-4
10-3
10
10
1000
2000
1000
2000
1000
100
101
10
1.0
UODAM E uphpp 12
10-
10-2
Training Loss
10--
10-1
10-4
10
10-1
10-1
10"2
10-4
10:4
10-4
10-5
1000
2000
1000
2000
1000
Tranning Iteration
Training Iteration
Training Iteration
Training Iteration
Training Iteration
Baseline
LoCun
Glorst
Power-Law
Figure 14: Training loss curves for function fitting benchmarks under baseline, JeCunnnnumeiian Glorot and power-law (a 0.25, 8 1.75) initializations when using learning-rate scheduler. Results are averaged over five seeds, with shaded regions indicating the standard error. Top row: "small" architecture (G 5, hidden layers with 8 neurons each). Bottom "large" archi- tecture (G = 20, three hidden layers with 32 neurons each).
Allen: Cahn
Burgers
Helmholtz
8 1900 2. 12509
Training Loss
1.0
10-4
10-1
1.022
1000
2000
3000
4000
5000
1300
2000
3000
4000
5000
1300
2000
3000
4000
3000
1.00
33 upm. E- uppp
10a
0z 0
Training
16a
10a
10-1
10-3
10a4
1000
4300
5000
1000
2000
4000
1000
2000
3000
4000
5000
Taining Iteration
Training Reration
Training Beration
Baseline
Powerl aw
Training loss curves for forward PDE benchmarks under baseline, LeCun-numerical Glorot and power-law (a 0.25, 8 1.75) initializations when using learning-rai -aae scheduler. Results are averaged over five seeds, with shaded regions indicating the standard error. Top row: "small" architecture (G' 5, two hidden layers with 8 neurons each). Bottom row: "large" archi- tecture (G = 20, three hidden layers with 32 neurons each).
Figure 15;
29

---

Published as a conference paper at ICLR 2026
in the main text except for the scheduler: training begins with learning rate of 10 3, followed by exponential decay with decay factor 0.85 every 100 iterations.
10-3,
Across both sets of benchmarks, the learning-rate schedulers eliminate the oscillatory behavior ob- served under fixed learning rate, yielding smoother training curves. Quantitatively, the final losses slightly higher than those reported in the main text, due the learning rate decaying even in gions where larger step size would allow for further progress. Nonetheless, this influences only the numerical values: the qualitative picture remains unchanged, and the relative performance ordering of the initialization schemes is consistent with the fixed learning rate setting.
30

---

Published as a conference paper at ICLR 2026
F NEURAL TANGENT KERNEL ANALYSIS
In this work, we use NTK analysis (Jacot et al. 2018) to better understand the effect of initialization schemes on function fitting and PDE benchmarks, both in terms of stability and conditioning.
F.1 NTK FOR PIML WITH RBA WEIGHTS
In this subsection, we derive the NTK formalism used in our PDE experiments. Specifically, we extend the standard NTK framework for PIML (Wang et al. 2022) to cover the RBA-weighted loss function of Eq. 31).
(pde)
We denote the PDE and boundary/initiall condition residuals at the i-th collocation point by r and r (be) respectively, as in Appendix C.2
We may re-weight the loss function of Eq. (3) to follow Wang et al. (2022) and subsequently write it vector form as
(()) = 2 || (0) +
1 _ =|
(0)1
F(() A (5)
(33)
are the RBA weights and (E) (E) Throughout single gradient step we treat (5) constants, as they are updated only between steps by Eq. (2). outside of the gradient descent scheme.
where p(() stacks the residuals of type & E (pde, bc), (E) (a
(E) No
Let (0) E Nex be the Jacobian of the residuals with respect to the parameters, i.e., its i-th
row is ((()(() (0) (r(()((0)))) (0)/00T For parameter update AO -1 V (0), first-order expansion around 0 yields.
A (0) = A (r) (0) ~ (55)
(34)
Using Eq.
33) and the chain rule, the full-batch gradient
N pde
Noc
VoC(O)
(pd))
(0) Ve P.
(pde) (0) +
(hc)
0 Vo Vor!
((bc) (0).
(35)
Since (()) a, r, and (5) is held fixed within the step.
Vo ((() (0) = a ( (0) = a (I)
(36)
Substituting A0 = -n ((( into Eq.
34 and grouping terms gives the linear dynamics
(F())(0) -1 (A (15) (0)) (A (pde) ( (ppe) (((dee
K (S.pde)
A (() (()) (A (bc) ((bc)
(37)
(E.bc)
As mentioned the main text, Eq. 37 shows that the weighted residual vectors ((() evolve under a weighted NTK with blocks
(E,0) (A(
&. e {pde, bc).
(38)
31

---

Published as a conference paper at ICLR 2026
F.2 NTK SPECTRA FOR VARYING POWER-LAN LAW EXPONENTS
To complement the heatmaps of Appendix[Daand to further illustrate the robustness of the power-law initialization, examine in this Appendix how the NTK spectrum varies across different (a, .8) configurations. The goal of this analysis is twofold. First, it provides an NTK-based view of the "good regions" identified in the grid search, showing how favorable exponent choices correspond to well-conditioned and stable kernels. Second, it demonstrates that the power-law scheme is not sensitive to a single finely tuned pair of exponents: more than one (ax, 5) combinations within the identified range yield well-behaved spectra. This supports the idea that one may tune the exponents
once per problem domain and thereafter select any configuration from the favorable region.
Figure 16 displays the NTK eigenvalue spectra for all exponent pairs considered in the grid search, using the "large" architecture (G 20, three hidden layers with 32 neurons each) and the function fitting target fI 3). Figure 17 shows the corresponding results for the PDE residual term of Burgers' equation. both cases, well well-conditioned spectra concentrate in the same regions suggested by the grid-search results.
0.0
025
0.5
- 0.75
# 1.4
9 1.25
1
1.75
$ 2,0
S E
154
Siperoall (log MAMM
Index
Index
Index
Index
Index
Index
Index
Index
Index
Ligenvalue index flog scale) Intermediate Iterations
Inttillization
Iteration
Figure 16: NTK eigenvalue spectra for the large architecture -G 20, three hidden layers with 32 neurons each) on function fitting target fa (x. y). shown for all (a, (8) configurations considered in the grid search. Each panel corresponds to one exponent pair and displays spectra at initialization, mid-training, and convergence.
32

---

Published as a conference paper at ICLR 2026
0.25
0m7.75
1m1..55
1=15
1.1775
2,5
SIRE RERT 5
ETTEE
SERR R
Eigenvaluns (log scale)
BEER
5 1 E
Index
Index
Index
Index
Index
Index
Index
Index
Eigenvalue index flog scale)
Inttialization
Figure 17: NTK eigenvalue spectra for the large architecture (G = 20, three hidden layers with 32 neurons each) on the PDE residual term of Burgers' equation, shown for all (a, 3) configurations considered in the grid search. Fach panel corresponds to exponent pair and displays spectra initialization, mid-training, and convergence.
33

---

Published as a conference paper at ICLR 2026
F.3 ADDITIONAL NTK SPECTRA
For completeness, we report additional NTK spectra not included in the main text. Figures 18 show the results for the remaining function fitting benchmarks (fi, ʃ2, f4. fs), while Figures correspond to the Burgers' and Helmholtz PDEs. All results are obtained using the "large"
21
22 23
architecture (G 20, three hidden layers with 32 neurons each) and values a 0.25, B 1.75 for the power-law initialization, consistent with the setting analyzed in Section
Baseline
LeCun
Gloret
Power-Law
10'
10a
104
Eigenvaluss
155
101
10a
10a
Indices
Indices
Indices
Indices
Intiaaization
Intermediate Itarations
Iteration
Figure 18: Eigenvalue spectra of the NTK matrix initialization (solid blue), intermediate iterations (dashed teal), and final iteration (dashed green) for function fitting benchmark (a, y) under dif- ferent initialization strategies. Results correspond the "large" architecture (G 20, three hidden layers with 32 neurons each). The power-law initialization uses a = 0.25, B 1.75.
Baseline
LeCun
Glorot
Power-Law
Eigenvaouss
164
12-
105
10-
104
10-7
25
Indioes
Indices
Intialization
Intermediate Iterations
Final Iteration
Figure 19: Eigenvalue spectra of the NTK matrix at initialization (solid blue), intermediate iterations
(dashed teal), and final iteration (dashed green) for function fitting benchmark f2(x, y) under dif- ferent initialization strategies. Results correspond the "large" architecture (G 20, three hidden layers with 32 neurons each). The power-law initialization uses a = 0.25, B 1.75.
Baseline
LeCun
Gloret
Power-Law
10'
12P
12a
10.
Figenvaluns
155
20a
101
10-
10a
Indices
Indices
Indices
Inzialization
Intermediate Iterations
Final Iteration
Figure 20: Eigenvalue spectra of the NTK matrix at initialization (solid blue), intermediate iterations
(dashed teal), and final iteration (dashed green) for function fitting benchmark fa y) under dif- ferent initialization strategies. Results correspond to the "large" architecture (G = 20, three hidden layers with 32 neurons each). The plwer-law initialization uses a = 0.25, B 1.75.
34

---

Published as a conference paper at ICLR 2026
Baseline
Gloret
10a
161
Eigenvalioss
307
10
10a
Indices
Indices
Indices
Indices
Figure 21: Eigenvalue spectra of the NTK matrix initialization (solid blue), intermediate iterations
(dashed teal), and final iteration (dashed green) for function fitting benchmark fs y) under dif- ferent initialization strategies. Results correspond the "large" architecture (G 20, three hidden layers with 32 neurons each). The power-law initialization uses = 0.25, B 1.75.
Glarot
1832
127
Eiganvakues (PPEE)
122
104
4200
162
1234
10
10-0
30
107
101
1.0
16a
101
3.7
12a
107
1.2
100
10a
107
Eignnaauuss 190
29a
132
134
20
10th
10
157
10a
122
1.7
10'
Indices
Indices
Indices
Indices
Initia liatton
Intermediate iterations
Final Iteration
Figure 22: NTK eigenvalue spectra for the Burgers' PDE benchmark under baseline, LeCun numerical, Glorot, and power-law (a 0.25, 3 1.75) initializations. Top row: spectra corre- sponding to the PDE residual term. Bottom row: spectra for the boundary/initial condition terms. Solid blue lines show the initialization, dashed teal lines show intermediate iterations, and dashed green lines show the final iteration.
35

---

Published as a conference paper at ICLR 2026
Baseline
LeCun
Giorot
Power-Law
Eigenvalues (POI)
1017
1022
19a
10
2014
100
199
2013
100
2014
134
157
2013
2014
10a
Eigievauuss (BC)
12a
16a
10'
1.05
10.4
10
10-4
10
100
10mh
10
Indices
Indices
Indices
Indices
Initaaliatton
Intermediate Iterations
Final Iteration
Figure 23: NTK eigenvalue spectra for the Helmholtz PDE benchmark under bascline, LeCun- numerical, Glorot, and plwer-law (a = 0.25, 3 1.75) initializations. Top row: spectra corre- sponding the PDE residual term. Bottom spectra for the boundary/initial condition terms. Solid blue lines show the initialization, dashed teal lines show intermediate iterations, and dashed green lines show the final iteration.
36


Reviews:
[{"review_id": "J7ZdPU86eI", "reviewer": "Reviewer_bFkh", "summary": "The paper proposes three families of initialization for spline‑based KANs: (i) LeCun‑inspired schemes (one ``numerical,'' one using batch‑normalized spline bases), (ii) a Glorot‑inspired scheme that balances forward/backward variances, and (iii) an empirical power‑law that scales residual/basis weights via exponents $\\alpha, \\beta$. Experiments include 2‑D function fitting, three forward PDEs (Allen–Cahn, Burgers, 2‑D Helmholtz), and a subset of the Feynman dataset; training‑dynamics are visualized with both loss curve analysis and NTK spectra. The numerics show that both Glorot and tuned power‑law tends to beat the original KAN baseline, and tuned power‑law is the strongest overall.", "questions": "The author's best region favors small $\\alpha$ and larger $\\beta$. Could the authors provide a heuristic/NTK‑based rationale and show how eigenvalue conditioning varies across $(\\alpha,\\beta)$?", "limitations": "", "rating": 6, "confidence": 3, "soundness": 2, "presentation": 3, "contribution": 3, "strengths": "1. Overall the paper is clearly written with clear problem statement and clean formulation for the KAN layer. The LeCun/Glorot initializations are derived explicitly.\n2. The paper include broad and thorough experimentation. The NTK spectra are a nice lens to visualize conditioning differences for different initialization.\n3. Even though the paper is empirical, better initialization recipes in different scenarios for KANs are of practical relevance.", "weaknesses": "1. Although the authors claim that the power-law initialization has the best overall performance, its performance is based on a grid search over $\\alpha,\\beta$ and only the best‑performing configuration is retained. In practice, such tuning can be costly, and the paper offers no low‑budget recipe to make the scheme usable out of the box.\n2. The LeCun/Glorot initialization are derived by applying the standard variance‑preserving idea to the KAN layer and rely on a strong equipartition assumption across the \n$G+k+1$ terms. The paper itself notes this simplification and that an alternative split performed worse. Meanwhile, the best results come from the empirical power‑law, for which the paper offers no mechanistic explanation beyond trends in heatmaps. The authors acknowledge this limitation, but as it stands, the conceptual contribution feels incremental and the empirical win might due to a tuned heuristic.", "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": "UkrpKdWNis", "reviewer": "Reviewer_tHHq", "summary": "The paper explores three new initialization schemes for KANs, identifying a current gap that KAN implementations rely on the standard initialization presented in the original KAN paper. The first two initialization methods explored LeCun and Glorot, are based in theory, the first driven by variance preservation, in which the variance of the inputs should match the variance of the outputs (LeCun) and the other driven by stabilizing variance for activations and gradients across layers. Additionally, the paper explores a third, empirical initialization based on the power-law scaling of KAN structure, that includes tunable exponents on the KAN's architecture (termed \"power-law\" in the paper). The paper goes on to evaluate the accuracy of these initialization methods on two main tasks, function fitting on five 2D target functions, and forward PDE problems. On both of these tasks, power-law initialization significantly outperforms other methods, with Glorot initialization as a distant second method. Additionally, the authors explore training stability, comparing loss curves and eigenvalue spectra of the NTK matrix.", "questions": "I would love to understand better what is driving the dynamics of the LeCun and Glorot initialization methods. I would also like to better understand the decision to add tunable parameters to an initialization method, whether there is a precedent here, and why we believe that would not result in overfitting.", "limitations": "", "rating": 4, "confidence": 3, "soundness": 2, "presentation": 3, "contribution": 2, "strengths": "The paper addresses an open question in KAN training and successfully demonstrates different initialization variants and their relative benefits. The paper motivates well the use of the two theoretically-driven initialization methods, and selects a good set of tasks and demonstrations to illustrate the relative benefits of each method. This is very fundamental work that is important to enabling the community to practically use KANs.", "weaknesses": "I appreciated the theoretical motivation of LeCun and Glorot initializations but was disappointed that there wasn't much discussions on why these methods behaved in the way they did. Table 1 for example shows that only about 15% of runs with LeCun initialization outperform the baseline, but there isn't an exploration into why this might be the case, or an exploration of the dynamics of this method of initialization. The same holds true for the results in Table 2. \n\nI'm not quite sure I understand what is being demonstrated in Figure 1 for the large architecture, as training instability for the power law initialization dominates the dynamics. The authors say \"The more pronounced oscillations observed in its curves stem from the fact that the model already reaches very low loss values under a fixed learning rate, leaving no smoothing effect from adaptive adjustments.\" but the learning rate in these plots is fixed to 10e-3 so I'm not sure what is meant by \"no smoothing effect from adaptive adjustments\" unless the learning rate is in fact changing. The same instability concern holds true for the PDE case.\n\nOverall, I have concerns with the power law initialization, as I believe this method is prone to overfitting, and that results of this method are artificially inflated. I'm not quite sure that comparing the performance of this method with the theoretically motivated initializations.", "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": "3wszlqpBPL", "reviewer": "Reviewer_c2rF", "summary": "This paper investigates initialization strategies for spline-based Kolmogorov-Arnold Networks (KANs), a relatively new neural network architecture. The authors argue that proper initialization is crucial for KAN performance but remains under-explored, and methods from traditional MLPs are not directly applicable. The paper proposes and evaluates several initialization schemes: two theoretically-motivated approaches inspired by LeCun and Glorot/Xavier variance preservation principles (including a variant using batch-normalized spline bases), and an empirical power-law family with tunable exponents. These schemes are evaluated through large-scale grid searches on function fitting and forward PDE approximation tasks, analysis of training dynamics (loss curves, Neural Tangent Kernel spectrum), and testing on a subset of the Feynman dataset. Key findings suggest that the Glorot-inspired scheme performs well, especially for larger models, while the empirical power-law initialization consistently achieves the best performance across tasks and architectures.", "questions": "1. Could you comment on how the performance of KANs with the proposed initializations compares to standard, well-initialized MLPs on some of the benchmark tasks?\n\n2.Could the sensitivity to the variance equipartition assumption in the LeCun/Glorot derivations be investigated? What happens if different partitions are assumed?\n\n3.Can the authors provide more insight into why the LeCun-normalized scheme failed to perform well? Is there an interaction between batch normalization and the spline basis learning?", "limitations": "", "rating": 6, "confidence": 2, "soundness": 2, "presentation": 3, "contribution": 3, "strengths": "1: Addresses a critical and timely aspect of training KANs – initialization,  and provides the first systematic comparison of different initialization strategies specifically tailored for spline-based KANs. \n\n2: The proposed power-law initialization demonstrably leads to superior performance across a wide range of settings in the experiments, offering a practical and effective strategy. The finding that small $\\alpha$ and large $\\beta$ are generally preferred provides a useful heuristic.", "weaknesses": "1. The paper is purely empirical, lacking a theoretical justification for the effectiveness of methods. For example, why do small $\\alpha$ and large $\\beta$ work well for power-law initialization? If relying solely on empirical search, there may be a suspicion of parametric-innovation.\n\n2. The data in experimental is relatively simple, considering only function fitting and PDE, which limits the generalizability of the author's conclusions.\n\n3. Lack of in-depth analysis on the reasons for poor initialization performance, such as for LeCun-normalized.", "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)