A simple mean field model of feature learning

## Abstract

Abstract Feature learning (FL), where neural networks adapt their internal representations during training, remains poorly understood. Using methods from statistical physics, we derive a tractable, self-consistent mean-field (MF) theory for the Bayesian posterior of two-layer non-linear networks trained with stochastic gradient Langevin dynamics (SGLD). At infinite width, this theory reduces to kernel ridge regression, but at finite width it predicts a symmetry breaking phase transition where networks abruptly align with target functions. While the basic MF theory provides theoretical insight into the emergence of FL in the finite-width regime, semi-quantitatively predicting the onset of FL with noise or sample size, it substantially underestimates the improvements in generalisation after the transition. We trace this discrepancy to a key mechanism absent from the plain MF description: self-reinforcing input feature selection . Incorporating this mechanism into the MF theory allows us to quantitatively match the learning curves of SGLD-trained networks and provides mechanistic insight into FL.

## 1 Introduction

The ability of deep neural networks to automatically learn useful features during training is widely regarded as a central factor in their success (Bengio et al., 2014; LeCun et al., 2015). The best understood theories of generalisation work in the infinite-width limit (Jacot et al., 2018; Lee et al., 2018). These yield valuable insights but cannot account for key finite-width phenomena. In particular, finite networks exhibit feature learning (FL) effects such as improved sample complexity (e.g., sparse parity (Damian et al., 2022; Daniely & Malach, 2020) and (multi-)index functions (Bietti et al., 2022)) and task-aligned changes in hidden representations (Papyan et al., 2020; Chizat et al., 2019; Corti et al., 2025; Nam et al., 2024). This gap highlights the need for theoretical frameworks that capture the mechanisms underlying FL.

### 1.1 Related work

Attempts to bridge this theoretical gap fall into two main categories. (1) Dynamical theories describe the evolution of network properties during training via integro-differential equations (Bordelon & Pehlevan, 2022; 2023; Mei et al., 2018; Montanari & Urbani, 2025; Celentano et al., 2025; Lauditi et al., 2025; Shi et al., 2022). These approaches, often centered on the transition from lazy to rich regimes, can yield accurate empirical predictions but are mathematically complex. (2) Static theories analyze ensembles of neural networks trained with stochastic gradient Langevin dynamics (SGLD), a process equivalent to sampling from the Bayesian posterior that simplifies the analysis of the final state after training (Welling & Teh, 2011; Teh et al., 2015). By focusing on the stationary distribution of the posterior, these theories capture feature learning (FL) through data-dependent kernel adaptations, either in double asymptotic limits of data and width (Pacelli et al., 2023; Baglioni et al., 2024; Fischer et al., 2024; Rubin et al., 2025) or via kernel renormalization (Howard et al., 2025; Aiudi et al., 2025). While both categories offer valuable insights, their reliance on complex theoretical machinery often obscures the core mechanisms driving FL.

Figure: Figure 1: a) Generalisation error vs. training set size $P$ for 2-layer $N=512$ width RELU networks trained with SGLD on $k$-sparse parity compared to the predicted error of the 3 theories presented in the paper (training details in [Appendix˜F](https://arxiv.org/html/2510.15174v1#A6)). b) Hierarchy of theories presented in this paper.
Refer to caption: x1.png

Here, we offer a complementary perspective, prioritizing simplicity and mechanistic interpretability. We begin by employing standard methods from statistical physics to derive a self-consistent mean-field (MF) theory for the posterior of a two-layer non-linear network trained with SGLD (building on Rubin et al. (2024)). This MF model can be viewed as a minimal extension beyond the fixed-kernel NNGP limit (see [Figure˜1](https://arxiv.org/html/2510.15174v1#S1.F1) b)). It predicts the onset of FL with increasing dataset size $P$ or decreasing noise strength $\kappa$ as a symmetry breaking phase transition where the initial isotropy of the weights is broken towards task-relevant directions given by the data.
While this MF theory predicts a mechanism for the onset of FL, it substantially underestimates the generalization gains of SGLD-trained networks after FL kicks in (see [Figure˜1](https://arxiv.org/html/2510.15174v1#S1.F1) a). We trace this discrepancy to a missing mechanism eliminated by our simple MF approximation: input feature selection (IFS), where networks dynamically amplify weights for relevant input dimensions, specializing subsets of neurons while leaving others inactive (see [Figure˜3](https://arxiv.org/html/2510.15174v1#S3.F3) b)). The homogeneity assumption in the basic MF model cannot capture these effects.

Our key contribution is to show that incorporating this IFS mechanism requires only a minimal, principled modification to the MF theory: endowing the weight prior with a learnable, coordinate-wise variance. The resulting model, which we call MF-ARD (Automatic Relevance Determination), preserves the simplicity and tractability of the MF framework while capturing the essence of FL.

Our contributions are:

- •
Interpretable MF theory: We derive a simple, self-consistent MF theory for the posterior of a non-linear two-layer NN. We show that this model captures the onset of FL as a phase transition but fails to predict the full learning curve of SGLD-trained networks.
- •
Identification of a core FL mechanism: We identify IFS as the key missing mechanism in standard MF theory and introduce MF-ARD, which captures this mechanism while preserving the simplicity of the MF form. We prove ([Theorem˜4.1](https://arxiv.org/html/2510.15174v1#S4.Thmtheorem1)) that this ARD extension eliminates the $\mathcal{O}(d)$ penalty in input dimension $d$ inherent in standard MF theory, providing a mechanistic understanding of how FL can overcome the curse of dimensionality.
- •
Quantitative prediction of generalisation error: We demonstrate that MF-ARD quantitatively predicts the generalisation error of SGLD-trained networks across varying dataset sizes and noise levels (see [Figure˜1](https://arxiv.org/html/2510.15174v1#S1.F1) a), [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5)).

## 2 Theory: A hierarchy of models (SGLD → \rightarrow MF → \rightarrow NNGP)

SGLD is the limit of full-batch GD with weight decay and injected Gaussian noise. It links optimisation and Bayesian inference by viewing parameter trajectories as samples from a posterior distribution (Welling & Teh, 2011; Teh et al., 2015).
Given a neural network described by $f_{{\bm{\theta}}}({\bm{x}}_{\mu})$, a dataset $\mathcal{D}=\{({\bm{x}}_{\mu},{\bm{y}}_{\mu})\}_{\mu=1}^{P}$, drawn from an input distribution $q({\bm{x}},{\bm{y}})$, the SGLD update equation is:

$$ $\displaystyle\Delta\theta_{i,t}:=-\,\eta\!\left[\gamma_{i}\,\theta_{i,t}+\nabla_{\theta_{i}}\!\left(\frac{1}{P}\sum_{\mu=1}^{P}\bigl(f_{{\bm{\theta}}}(\mathbf{x}_{\mu})-{\bm{y}}_{\mu}\bigr)^{2}\right)\right]+\sqrt{2T\,\eta}\;\xi_{i,t},\quad\xi_{i,t}\sim\mathcal{N}(0,1),$ (1) $$

with full-batch gradients on parameters $\theta_{i}$, batch noise replaced by Gaussian noise $\xi$, learning rate $\eta$, weight decay $\gamma_{i}$, and noise strength set by $T$. The stationary posterior $p_{\textrm{GD}}$ of [Equation˜1](https://arxiv.org/html/2510.15174v1#S2.E1) as $\eta\to 0$ for an $L$-layer network (widths $N_{l}$) is

$$ $\displaystyle-\ln p_{\mathrm{GD}}({\bm{\theta}}|\mathcal{D})\propto\underbrace{\sum_{l=1}^{L}\frac{1}{2\sigma_{i}^{2}}\sum_{i=1}^{N_{l}}\|{\bm{\theta}}_{i}\|^{2}}_{-\ln p_{\text{prior}}}+\underbrace{\frac{1}{2\kappa^{2}P}\sum_{\mu=1}^{P}(f_{{\bm{\theta}}}({\bm{x}}_{\mu})-{\bm{y}}_{\mu})^{2}}_{-\ln p_{\text{L}}},$ (2) $$

where ${\bm{\theta}}_{i}$ are weight matrices, and $\sigma_{i}^{2}=T/\gamma_{i}$ is related to noise $\kappa^{2}=T/2$. Training with SGLD can be viewed as sampling from the posterior in [Equation˜2](https://arxiv.org/html/2510.15174v1#S2.E2) (see [Section˜A.1](https://arxiv.org/html/2510.15174v1#A1.SS1) for more details).

### 2.1 SGLD posterior: Fully interacting theory

In this paper, we focus on two-layer fully-connected networks with input dimension $d$, hidden width $N$, input weights ${\bm{w}}_{i}\in\mathbb{R}^{d}$, output weights $a_{i}$, output dimension 1, and nonlinearity $\phi$:

$$ $\displaystyle f_{{\bm{\theta}}}({\bm{x}})=\frac{1}{N^{\gamma}}\sum_{i=1}^{N}a_{i}\,\phi({\bm{w}}_{i}^{\!\top}{\bm{x}}),\qquad w_{ij}\sim\mathcal{N}\!\left(0,\frac{\sigma_{w}^{2}}{d}\right),\quad a_{i}\sim\mathcal{N}(0,\sigma_{a}^{2}).$ (3) $$

The initial parameters are set with Gaussians.
The scale factor $N^{-\gamma}$ on the output layer enables interpolation between $\gamma=1/2$ (NTK scaling) and $\gamma=1$ (mean field scaling, not to be confused with MF theory, see e.g. Mei et al. (2019) for details).
Using the explicit form of $f_{{\bm{\theta}}}$ in [Equation˜3](https://arxiv.org/html/2510.15174v1#S2.E3) to rewrite [Equation˜2](https://arxiv.org/html/2510.15174v1#S2.E2) by expanding the square inside $-\ln p_{\mathrm{L}}$ yields (see [Section˜A.2](https://arxiv.org/html/2510.15174v1#A1.SS2) for the algebra):

#### FL mechanism as complex neuron-neuron interaction

The following competing forces shape the posterior: The self-energy $\Sigma$ acts as a regularizer, keeping neurons from becoming too large. The data coupling term $J_{\mathcal{Y}}$ is the learning signal: it rewards neurons when they align with the target function. Finally, the interaction kernel $G$ captures how neurons influence each other. As $P$ increases (or $\kappa$ decreases), the likelihood tilts the posterior away from the isotropic Gaussian prior, creating anisotropy along task-relevant directions and cooperative alignment across neurons via $G$. This can yield a transition to a non-Gaussian posterior concentrated on low-rank, task-aligned structures.

This posterior is not analytically tractable due to the neuron-neuron coupling via $G$: every neuron’s optimal weight is a function of every other neuron (see [Figure˜2](https://arxiv.org/html/2510.15174v1#S2.F2) for an illustration). From a Bayesian perspective, this means we cannot write the posterior distribution as a simple product of independent terms; it is instead a complex, high-dimensional distribution.

Figure: Figure 2: SGLD Posterior: The model is a fully interacting system where every neuron (row ${\bm{w}}_{i}$) is coupled to every other neuron (${\bm{w}}_{i^{\prime}}$) through the interaction kernel $G$. This results in a complex, high-dimensional posterior that is computationally intractable. MF theory: The pairwise interactions are replaced by an effective field that represents the average influence of all other neurons. Each neuron is now treated as an independent sample from a single, shared, data-dependent distribution, $p_{\textrm{MF}}({\bm{w}})$. NNGP Limit: All interactions are removed. Each neuron is an independent and identically distributed sample drawn from the fixed prior distribution.
Refer to caption: x2.png

### 2.2 MF Theory: Removing off-diagonal couplings

One of the simplest ways to make the posterior in [Equation˜4](https://arxiv.org/html/2510.15174v1#S2.E4) tractable is inspired by statistical physics and called self-consistent mean-field (MF) theory. The main idea is to replace the highly-correlated posterior over the entire weight matrix, $p(\bm{W}|\mathcal{D})$, with a fully factorized approximation where each of the $N$ neurons is independent: $p(\bm{W}|\mathcal{D})\approx\prod_{i=1}^{N}p_{\text{MF}}({\bm{w}}_{i})$ (see [Figure˜2](https://arxiv.org/html/2510.15174v1#S2.F2) for an illustration). Instead of interacting with all other neurons, each interacts with the average behavior of all other neurons. Take any neuron ${\bm{w}}_{i}$, replace the interaction term $\sum_{i^{\prime}\neq i}a_{j}G({\bm{w}}_{i},{\bm{w}}_{i}^{\prime})$, which couples neuron ${\bm{w}}_{i}$ to other neurons ${\bm{w}}_{i^{\prime}}$, with their collective average effect, the ‘mean field’ $\langle f\rangle$:

$$ $\displaystyle\sum_{i^{\prime}\neq i}a_{i^{\prime}}\,G({\bm{w}}_{i},{\bm{w}}_{i^{\prime}})$ $\displaystyle=\frac{1}{P}\sum_{\mu}\phi({\bm{w}}_{i}^{\top}{\bm{x}}_{\mu})\,\overbracket{\sum_{{i^{\prime}}\neq i}a_{i^{\prime}}\,\phi({\bm{w}}_{i^{\prime}}^{\top}{\bm{x}}_{\mu})\;\xrightarrow{\ \mathrm{MF}\ }\ N^{\gamma}\langle f({\bm{x}})\rangle}^{\text{MF closure}}.$ (5) $$

For the approximation to be valid, it must be self-consistent: the average behavior of a neuron drawn from our approximate posterior $p_{\text{MF}}({\bm{w}})$ must exactly reproduce the mean field that we assumed in the first place, giving rise to the following self-consistency equation: $\langle f({\bm{x}})\rangle=N^{1-\gamma}\cdot\mathbb{E}_{({\bm{w}},a)\sim p({\bm{w}},a)}\bigl[a\,\phi\bigl({\bm{w}}^{\top}{\bm{x}}\bigr)\bigr]$.

The posterior is best understood in terms of an orthonormal basis $\{\chi_{A}\}$ w.r.t. the input distribution $q({\bm{x}},{\bm{y}})$, where $A$ indexes the basis. Expanding the MF $\langle f(\mathbf{x})\rangle=\sum_{A}m_{A}\chi_{A}(\mathbf{x})$ in said basis, where the feature coefficients are $m_{A}=\mathbb{E}_{{\bm{x}}}\bigl[\langle f({\bm{x}})\rangle\,\chi_{A}({\bm{x}})\bigr]$, gives the following theory:

The core simplification of the MF theory is that each neuron is treated as an independent sample drawn from the same distribution, which is a function of the prior and the data. This strongly simplifies the problem. We replace the intractable problem over the entire interacting $N\times d$ weight matrix with a much simpler, self-consistent problem for a single $d$-dimensional weight vector.

These fixed-point equations can easily be solved by iterating to self-consistency, enabling the computation of all statistics $\Sigma({\bm{w}}),J_{A}({\bm{w}}),J_{\mathcal{Y}}({\bm{w}})$ directly from the finite training dataset of size $P$. In this way, we naturally capture the effects of training with limited samples (see [Appendix˜E](https://arxiv.org/html/2510.15174v1#A5) for details).

### 2.3 NNGP limit: No interactions

The MF model can be simplified even further in the infinite-width limit with $\gamma=1/2$ (NTK-scaling). For this scaling the limit is well-behaved, the self-energy term in [Equation˜6](https://arxiv.org/html/2510.15174v1#S2.E6) vanishes as well as any $m_{A}$-dependent tilt, resulting in $p({\bm{w}})$ collapsing to its prior (see [Section˜D.1](https://arxiv.org/html/2510.15174v1#A4.SS1) for a proof).

This NNGP limit represents the most restrictive limit in our approximation hierarchy. While MF theory eliminates inter-neuron interactions ($w_{ij}\leftrightarrow w_{i^{\prime}j}$), the NNGP limit additionally removes intra-neuron coordinate coupling ($w_{ij}\leftrightarrow w_{ij^{\prime}}$), see [Figure˜2](https://arxiv.org/html/2510.15174v1#S2.F2) for an illustration. Every weight is sampled from the same distribution (the prior).

#### No FL mechanism

In the infinite-width limit, there is *no* FL in the sense above. With no data-dependent term ($J_{\mathcal{Y}}$) to break symmetry, all neurons remain frozen at their prior. The feature coefficients $m_{A}^{\infty}$ are now determined purely by the fixed kernel $K_{AB}=\mathbb{E}_{{\bm{w}}\sim p_{\infty}}[J_{A}({\bm{w}})J_{B}({\bm{w}})]$ (see [Section˜A.5](https://arxiv.org/html/2510.15174v1#A1.SS5) how this connects to the usual kernel formulation of the NNGP limit).

## 3 Problem: Simple MF does not capture a central FL mechanism

In the following section, we will argue that FL in 2-layer finite-width networks is a two-stage process:

- (i)
Onset (phase transition): The signal from the data is strong enough to overcome the isotropic prior/self-energy, so the feature coefficients $m_{A}=\mathbb{E}_{\bm{x}}[\langle f({\bm{x}})\rangle_{p}\chi_{A}({\bm{x}})]$ turn nonzero.
- (ii)
Specialization through self-reinforcing input feature selection (IFS): After the symmetry breaking, the neurons aligning with the target receive disproportionately larger updates, producing heavy-tailed weight marginals on task-relevant coordinates and neuron-wise sparsification.

Plain MF captures *(i)* but largely misses *(ii)*, which is why for a target $y({\bm{x}})=\chi_{S}({\bm{x}})$, it underestimates the growth of $m_{S}$ and the post-transition generalisation (e.g., [Figure˜1](https://arxiv.org/html/2510.15174v1#S1.F1) a), [Figure˜3](https://arxiv.org/html/2510.15174v1#S3.F3) a)).

### 3.1 Stage 1: The Onset of Feature Learning as a Phase Transition

The simplicity of the MF theory allows us to interpret FL as a *self-consistent symmetry breaking* in the posterior $p_{\textrm{MF}}$, where learning emerges from a competition between two opposing forces. Below a critical $(P_{c},\kappa_{c})$ the only stable fixed point solution of [Equations˜6](https://arxiv.org/html/2510.15174v1#S2.E6) and [7](https://arxiv.org/html/2510.15174v1#S2.E7) is $m_{A}=0\;\forall A$ and $p_{\textrm{MF}}$ corresponds to the Gaussian prior. The Gaussian weight prior and the self-energy term $\Sigma({\bm{w}})$ act as regularizing forces, favouring an isotropic, high-entropy state that penalizes large weights and encourages neurons to remain small, unaligned and near their initialization. Upon crossing the critical values $(P_{c},\kappa_{c})$ when the signal from the data (controlled by dataset size $P$ and noise $\kappa$) becomes strong enough, the $m_{A}=0$ fixed point becomes unstable (see [Theorem˜D.7](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem7) for a proof of the phase transition in $\kappa$ and an explicit formula for $\kappa_{c}$). $p_{\textrm{MF}}$ becomes non-Gaussian, as the neuron-data coupling term ($J_{\mathcal{Y}}$), which acts like an external field, rewards neurons whose activations $\phi({\bm{w}}^{\top}{\bm{x}})$ correlate with the target function $y({\bm{x}})$, breaking the prior’s symmetry and pulling the weights toward a task-aligned configuration. For the case where the target function is a single basis function $y({\bm{x}})=\chi_{S}({\bm{x}})$, the relevant order parameter is the feature coefficient $m_{S}$ of the target function. In this case, equation [7](https://arxiv.org/html/2510.15174v1#S2.E7) becomes (see [Section˜D.5](https://arxiv.org/html/2510.15174v1#A4.SS5) for the derivation):

$$ $\displaystyle m_{S}\;=\;\frac{N^{1-2\gamma}}{\kappa^{2}}\,(1-m_{S})\;\underbrace{{\Biggl\langle\frac{J_{S}({\bm{w}})^{2}}{\displaystyle\sigma_{a}^{-2}+\frac{\Sigma({\bm{w}})}{\kappa^{2}N^{2\gamma}}}\Biggr\rangle_{{\bm{w}}\sim p({\bm{w}}\mid m_{S})}}}_{\omega_{0}}.$ (11) $$

The phase transition occurs when $\omega_{0}>1$, as shown in Figure [3](https://arxiv.org/html/2510.15174v1#S3.F3)a): for the SGLD-trained network and MF theory, the order parameter $m_{S}=0$ until the critical dataset size $P_{c}$, when it increases abruptly. The NNGP model, lacking this non-linear feedback, shows only a smooth, continuous rise. In [Section˜A.4](https://arxiv.org/html/2510.15174v1#A1.SS4) we link this FL phase transition to the emergence of an outlier eigenvalue in the Gram kernel $G({\bm{w}}_{i},{\bm{w}}_{i’})$.

However, Figure [1](https://arxiv.org/html/2510.15174v1#S1.F1) a) reveals the limitation: after the transition, the SGLD-trained network’s generalisation error decreases rapidly, while the MF model’s error decreases much more slowly. This discrepancy reveals a second powerful mechanism at play that MF theory misses.

### 3.2 Stage 2: Specialization via IFS (what simple MF misses)

What enables the SGLD-trained network to improve so rapidly post-transition? The answer lies in a self-reinforcing dynamic that is lost in the MF approximation.

A consequence of IFS is sparsification: only a small number of neurons strongly align with the target while others remain near initialization. We quantify this through neuron specialization $\tilde{m}_{S}=\langle aJ_{S}({\bm{w}})\rangle_{p_{\mathrm{GD}}}$, measuring how strongly individual neurons couple to the target function (see [Figure˜4](https://arxiv.org/html/2510.15174v1#S3.F4) b)). Before the phase transition $P=100<P_{c}$, $p(\tilde{m})$ is sharply peaked at zero. After the transition $P=1000>P_{c}$, $p(\tilde{m})$ develops pronounced tails, most neurons remain unspecialized while a few become strongly task-aligned (Also, see the plot of $\bm{W}\bm{W}^{\top}$ in [Figure˜6](https://arxiv.org/html/2510.15174v1#A2.F6) where task-relevant coordinates ($j=0,1,2,3$) have much higher norm.).

Figure: Figure 3: $k$-sparse parity target $y({\bm{x}})=\chi_{S}({\bm{x}})=\Pi_{j\in S}x_{j}$ with index $S=\{0,1,2,3\}$ in $d=35$ with a ReLU network ($N=512$, see [Section˜A.3](https://arxiv.org/html/2510.15174v1#A1.SS3) and [Appendix˜F](https://arxiv.org/html/2510.15174v1#A6) for more details.). a) $m_{S}$ vs. $P$: SGLD, MF and MF-ARD exhibit a phase transition at $P_{c}$, NNGP grows smoothly (no FL). b) Weight matrix plotted for the SGLD-trained NN: After the phase transition, only a very small set of neurons (9 out of 512) pick up a high correlation with the target ($m_{S}$ is color coded) and only for these neurons there is a strong symmetry breaking where the norm of the first $j=0,1,2,3$ coordinates is much larger than for the rest $d-k$ coordinates (coded as the height in the plot).
Refer to caption: x3.png

### 3.3 Why plain MF misses specialization

MF replaces the neuron–neuron interaction by a single effective field. As a result, *all* neurons are sampled from the *same* single-neuron distribution $p({\bm{w}},a)$ and are pushed to align in the *same average* way. This homogeneity suppresses the greedy dynamics that drive IFS and sparsification: The marginals $p({\bm{w}}_{j})$ for MF theory in [Figure˜4](https://arxiv.org/html/2510.15174v1#S3.F4) b) show very weak tails compared to SGLD ([Figure˜4](https://arxiv.org/html/2510.15174v1#S3.F4) a)) $\Rightarrow$ weak IFS. Similarly, the lack of heavy-tailed specialization $\Rightarrow$ weak sparsification.
Hence, the simplification of removing neuron-neuron interaction weakens IFS. We hypothesize that this explains the persistent generalization error gap we observed relative to SGLD in many settings.

Figure: Figure 4: Setup as in [Figure˜3](https://arxiv.org/html/2510.15174v1#S3.F3). a) Coordinate-wise weight marginals $p(w_{j})$: For SGLD, $p(w_{j})$ has high variance and strong non-Gaussianity on input relevant coordinates $j=0,1,2,3$ (symmetry breaking with IFS). Plain MF also shows symmetry breaking for coordinates $j=0,1,2,3$, but much weaker, while MF-ARD restores the strong anisotropy. b) Distribution of specialization $\tilde{m}{=}aJ_{S}({\bm{w}})$ before vs. after phase transition: SGLD develops heavy tails after the transition (sparsification), MF remains narrow, MF-ARD produces the heavy-tailed specialization. Together, MF-ARD qualitatively recovers the two mechanisms, IFS and sparsification, that MF misses.
Refer to caption: x4.png

## 4 A minimal MF model of feature learning

How can we capture IFS within the tractable MF framework given that it emerges from neuron-neuron interactions that MF explicitly removes? Fully restoring these interactions would make the theory intractable again. The key insight is that IFS’s essential signature is the heavy-tailed coordinate distributions $p(w_{j})$, not the specific interactions that create them. The isotropic prior $p({\bm{w}})=\prod_{j}\mathcal{N}(0,\sigma_{w}^{2}/d)$ prevents this by imposing a uniform norm penalty across all coordinates. Coordinate-dependent variances $p({\bm{w}})=\prod_{j}\mathcal{N}(0,\rho_{j}^{2})$ can selectively reduce penalties on task-relevant coordinates, enabling the heavy-tailed distributions that characterize IFS, without requiring explicit neuron-neuron coupling.

### 4.1 Model - ARD as the smallest change that enables IFS

Having introduced coordinate-dependent variances, we must specify an appropriate prior $p(\rho_{j})$. This approach is well-established in Bayesian neural networks under the name: Automatic Relevance Determination (ARD) (MacKay, 1996; Tipping, 2001). We follow the standard practice and use a gamma distribution prior, which ensures positivity of the variances and provides conjugacy with the Gaussian likelihood which enables computational tractability:

$$ $\displaystyle p({\bm{w}}\mid\bm{\rho})$ $\displaystyle=\prod_{j=1}^{d}\mathcal{N}\!\big(w_{j}\mid 0,\ \rho_{j}^{-1}\big),\qquad p(\bm{\rho})=\prod_{j=1}^{d}\Gamma\!\big(\rho_{j}\mid\alpha_{0},\beta_{0}\big).$ (13) $$

We set $\beta_{0}$ so that $\mathbb{E}[\rho_{j}]=\alpha_{0}/\beta_{0}=d/\sigma_{w}^{2}$, i.e. the ARD prior matches the isotropic initialization in expectation and the model reduces to plain MF at the beginning of the fixed point iteration,
leaving only $\alpha_{0}$ as a free hyperparameter.
(^1^11For widths used here, $\alpha_{0}\ll N$, so results are insensitive to the exact value of $\alpha_{0}$.)
Integrating out $\rho_{j}$ yields heavy–tailed marginals
$p(w_{j})\propto(\beta_{0}+\tfrac{1}{2}w_{j}^{2})^{-(\alpha_{0}+1/2)}$:
strong shrinkage near zero (irrelevant coordinates) and weak shrinkage in the tails (relevant coordinates).
The ARD fixed point is given by stationarity of the negative-log evidence w.r.t. $\rho_{j}$ (see [Section˜A.6](https://arxiv.org/html/2510.15174v1#A1.SS6) for the explicit derivation):

The MF-ARD theory is still a MF theory, as neurons remain independent and interact only through the self-consistent averages $\{m_{A}\}$ and $\{\rho_{i}\}$. The crucial difference is the introduction of the new order parameters $\{\rho_{j}\}$, which allow the model to learn the relevance of each input feature.

#### How ARD induces IFS

We assume a target $y({\bm{x}})=\chi_{S}({\bm{x}})$. The map $\langle w_{j}^{2}\rangle_{p_{\mathrm{ARD}}}\mapsto\rho_{j}\mapsto p_{\mathrm{ARD}}({\bm{w}})$ realises IFS: if coordinate $j$ aligns with the target ($j\in S$), the term $J_{\mathcal{Y}}({\bm{w}})-\sum_{A}m_{A}J_{A}({\bm{w}})$ increases, making $\langle w_{j\in S}^{2}\rangle_{p_{\mathrm{ARD}}}\uparrow$ larger, which *decreases* $\rho_{j\in S}\downarrow$ via [Equation˜15](https://arxiv.org/html/2510.15174v1#S4.E15), thereby *reducing* shrinkage on $w_{j\in S}$ and further increasing $\langle w_{j\in S}^{2}\rangle_{p_{\mathrm{ARD}}}\uparrow$.
Conversely, for non-aligned coordinates $j\notin S,\langle w_{j}\rangle^{2}$ remains small, which increases $\rho_{j}$ and strengthens shrinkage. This positive feedback realises IFS and leads to sparsification of the posterior (see [Figure˜4](https://arxiv.org/html/2510.15174v1#S3.F4) a) and b)). Accordingly, in ${\bm{w}}=(w_{1},...,w_{k},w_{k+1},...,w_{d})$ the first $k$ weights are much larger than the remaining $d-k$. This increases $J_{S}({\bm{w}})$ while decreasing $J_{A}({\bm{w}})\ \forall A\neq S$, driving $m_{S}$ up strongly after the phase transition while decreasing $m_{A}\ \forall A\neq S$. This explains the strong generalisation error drop of MF-ARD after the phase transition (unlike MF).
Note that MF-ARD is still a *static* theory. The mechanism above reflects how the fixed point map contracts to its solution (via [Equations˜14](https://arxiv.org/html/2510.15174v1#S4.E14) and [15](https://arxiv.org/html/2510.15174v1#S4.E15)), not a model of the actual training dynamics. See [Section˜A.7](https://arxiv.org/html/2510.15174v1#A1.SS7) for a discussion of the infinite width limit of the MF-ARD model.

### 4.2 How MF-ARD beats the curse of dimensionality

Let $y({\bm{x}})=\chi_{S}({\bm{x}})$ with $|S|=k$ in $d$ dimensions and $\phi=\mathrm{ReLU}$.
Write $\kappa_{c}^{2}$ for the critical noise at onset of the phase transition. We assume the infinite data limit and focus on the phase transition in the critical noise as it is analytically more tractable than critical data size.

###### Theorem 4.1 .

There exists an outer iteration $t_{0}=\mathcal{O}(1)$ of the fixed point algorithm and a constant $\varepsilon_{0}>0$ (independent of $d$ such that we get $\varepsilon_{0}$-level symmetry breaking towards $S$: $\min_{j\in S}\langle w_{j}^{2}\rangle_{p_{\textrm{ARD}}}^{t_{0}}\;-\;\max_{j\notin S}\langle w_{j}^{2}\rangle_{p_{\textrm{ARD}}}^{t_{0}}\;\geq\;c\cdot\varepsilon_{0}$ (see [Section˜D.6](https://arxiv.org/html/2510.15174v1#A4.SS6) for details). Then, the critical noise scales as:

$$ $\kappa_{c}^{2}\;\asymp\;\begin{cases}\Theta\!\big(\sqrt{1/(dk)}\big)&\text{(plain MF)},\[2.0pt] \Theta\!\big(\sqrt{1/k}\big)&\text{(MF--ARD)}.\end{cases}$ (16) $$

###### Proof.

See [Section˜D.6](https://arxiv.org/html/2510.15174v1#A4.SS6) for a proof and exact constants. Here we summarize the scaling.
∎

$\varepsilon$-symmetry breaking assumes that there is already a very small difference in the variance for weights that are relevant vs. irrelevant for learning the target. Given this asymmetry, MF-ARD can amplify it far more effectively than plain MF. For plain MF, the noise threshold required to transition from high to low generalisation error (equivalently $m_{S}=0\rightarrow m_{S}\neq 0$) scales with the ambient dimension $d$, requiring noise levels $1/\sqrt{d}$ times smaller than MF-ARD. This dimensional dependence is eliminated in MF-ARD. Consequently, the phase boundary scales with the intrinsic problem dimension $k$ rather than the ambient dimension $d$.

###### Theorem 4.1 .

###### Proof.

### 4.3 Results- numerical prediction of learning curves

We evaluate on $k$-sparse parity, a setting where fixed kernels are known to require super-polynomial samples in $d,k$, while finite-width networks can have lower sample complexity due to FL (Daniely & Malach, 2020; Damian et al., 2022; Barak et al., 2022). [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5) shows test-MSE heatmaps over dataset size $P$ and noise $\kappa$ for
(a) SGLD-trained networks, (b) plain MF, and (c) MF–ARD. All three display a transition from high
to low error as $P$ increases or $\kappa$ decreases. Plain MF detects *when* learning
starts but yields a diffuse, shifted boundary and overestimates post-transition error. MF–ARD, with a
single additional set of order parameters $\{\rho_{j}\}$, closely tracks both the *location* and the
*sharpness* of the SGLD phase boundary and reproduces the “helpful noise” regime in which
moderate $\kappa$ lowers the critical sample size $P_{c}$ (the kink around $\kappa=0.05$).

In [Figure˜7](https://arxiv.org/html/2510.15174v1#A2.F7) in the appendix we present the same phase diagram analysis for a single index-model task (Bietti et al., 2022), where the target function takes the form $y({\bm{x}})=g(\mathbf{v}^{T}{\bm{x}})$. Here, $\bm{v},{\bm{x}}\in\mathbb{R}^{d}$, the input ${\bm{x}}$ follows a standard multivariate normal distribution $\mathcal{N}(0,\mathbbm{1})$, $v_{j\in S}=1/\sqrt{k},\ v_{j\notin S}=0$ for an index set with $|S|=k$ and $g$ is a non-linear function (Hermitian polynomial in our case).

Figure: Figure 5: Generalisation phase diagram (Test MSE) vs. dataset size $P$ and noise $\kappa$ for $k$-sparse parity target $y({\bm{x}})=\chi_{S}({\bm{x}})$ with $S=\{0,1,2,3\}$ in $d=35$. a) SGLD, b) MF, c) MF-ARD: All panels show a transition from high to low test error as $P$ increases or $\kappa$ decreases. However, the plain MF theory strongly underestimates the sharpness of the transition while the MF-ARD theory largely reproduces the shape/location of this boundary. Training details are in [Appendix˜F](https://arxiv.org/html/2510.15174v1#A6).
Refer to caption: x5.png

## 5 Conclusion

We present a MF theory for the posterior of SGLD-trained two-layer networks that is simple to interpret, with a clear infinite-width NNGP limit. Extending it via Automatic Relevance Determination (MF-ARD) to include coordinate-dependent weight variances preserves this simplicity while quantitatively matching network performance as a function of data set size and noise on tasks such as $k$-sparse parity and index models, where feature learning is essential.

Our framework characterises feature learning as a two-stage process: a data-driven phase transition mechanism that initiates feature learning, followed by a self-reinforing input feature selection mechanism that leads to sparsification and drives improved generalisation. This explains why finite-width networks outperform kernel methods: they reshape their feature space in a self-reinforcing way that amplifies task-relevant signals and dramatically improves sample complexity. Standard MF theories capture only the transition, but MF-ARD captures both mechanisms.

Limitations of our approach include the focus on two-layer networks and tasks with sparse structures for which the ARD mechanism may be best suited. Extending to deeper architectures (e.g., convolutional or attention layers), and settings requiring distributed or smooth representations remains an open challenge (Petrini et al., 2022), as does connecting the static posterior view to training dynamics.

## References

- Aiudi et al. (2025)
R Aiudi, R Pacelli, P Baglioni, A Vezzani, R Burioni, and P Rotondo.
Local kernel renormalization as a mechanism for feature learning in overparametrized convolutional neural networks.
*Nature Communications*, 16(1):568, 2025.
- Baglioni et al. (2024)
P Baglioni, R Pacelli, R Aiudi, F Di Renzo, A Vezzani, R Burioni, and P Rotondo.
Predictive power of a bayesian effective action for fully connected one hidden layer neural networks in the proportional limit.
*Physical Review Letters*, 133(2):027301, 2024.
- Barak et al. (2022)
Boaz Barak, Benjamin Edelman, Surbhi Goel, Sham Kakade, Eran Malach, and Cyril Zhang.
Hidden progress in deep learning: Sgd learns parities near the computational limit.
*Advances in Neural Information Processing Systems*, 35:21750–21764, 2022.
- Bengio et al. (2014)
Yoshua Bengio, Aaron Courville, and Pascal Vincent.
Representation learning: A review and new perspectives, 2014.
URL [https://arxiv.org/abs/1206.5538](https://arxiv.org/abs/1206.5538).
- Bietti et al. (2022)
Alberto Bietti, Joan Bruna, Clayton Sanford, and Min Jae Song.
Learning single-index models with shallow neural networks, 2022.
URL [http://arxiv.org/abs/2210.15651](http://arxiv.org/abs/2210.15651).
- Bordelon & Pehlevan (2022)
Blake Bordelon and Cengiz Pehlevan.
Self-consistent dynamical field theory of kernel evolution in wide neural networks.
*Advances in Neural Information Processing Systems*, 35:32240–32256, 2022.
- Bordelon & Pehlevan (2023)
Blake Bordelon and Cengiz Pehlevan.
Dynamics of finite width kernel and prediction fluctuations in mean field neural networks.
*Advances in Neural Information Processing Systems*, 36:9707–9750, 2023.
- Celentano et al. (2025)
Michael Celentano, Chen Cheng, and Andrea Montanari.
The high-dimensional asymptotics of first order methods with random data.
*arXiv preprint arXiv:2112.07572*, June 2025.
doi: 10.48550/arXiv.2112.07572.
URL [https://arxiv.org/abs/2112.07572](https://arxiv.org/abs/2112.07572).
- Chizat et al. (2019)
Lenaic Chizat, Edouard Oyallon, and Francis Bach.
On lazy training in differentiable programming.
*Advances in neural information processing systems*, 32, 2019.
- Corti et al. (2025)
Andrea Corti, Rosalba Pacelli, Pietro Rotondo, and Marco Gherardi.
Microscopic and collective signatures of feature learning in neural networks.
*arXiv preprint arXiv:2508.20989*, August 2025.
URL [https://arxiv.org/abs/2508.20989](https://arxiv.org/abs/2508.20989).
- Damian et al. (2022)
Alexandru Damian, Jason Lee, and Mahdi Soltanolkotabi.
Neural networks can learn representations with gradient descent.
In *Proceedings of Thirty Fifth Conference on Learning Theory*, pp. 5413–5452. PMLR, 2022.
URL [https://proceedings.mlr.press/v178/damian22a.html](https://proceedings.mlr.press/v178/damian22a.html).
ISSN: 2640-3498.
- Daniely & Malach (2020)
Amit Daniely and Eran Malach.
Learning parities with neural networks.
*Advances in Neural Information Processing Systems*, 33:20356–20365, 2020.
- Fischer et al. (2024)
Kirsten Fischer, Javed Lindner, David Dahmen, Zohar Ringel, Michael Krämer, and Moritz Helias.
Critical feature learning in deep neural networks.
*arXiv preprint arXiv:2405.10761*, 2024.
- Howard et al. (2025)
Jessica N Howard, Ro Jefferson, Anindita Maiti, and Zohar Ringel.
Wilsonian renormalization of neural network gaussian processes.
*Machine Learning: Science and Technology*, 6(2):025038, 2025.
- Jacot et al. (2018)
Arthur Jacot, Franck Gabriel, and Clément Hongler.
Neural tangent kernel: Convergence and generalization in neural networks.
*Advances in neural information processing systems*, 31, 2018.
- Lauditi et al. (2025)
Clarissa Lauditi, Blake Bordelon, and Cengiz Pehlevan.
Adaptive kernel predictors from feature-learning infinite limits of neural networks.
*arXiv preprint arXiv:2502.07998*, 2025.
- LeCun et al. (2015)
Yann LeCun, Yoshua Bengio, and Geoffrey Hinton.
Deep learning.
*Nature*, 521(7553):436–444, May 2015.
doi: 10.1038/nature14539.
URL [https://www.nature.com/articles/nature14539](https://www.nature.com/articles/nature14539).
- Lee et al. (2018)
Jaehoon Lee, Yasaman Bahri, Roman Novak, Samuel S. Schoenholz, Jeffrey Pennington, and Jascha Sohl-Dickstein.
Deep neural networks as gaussian processes, 2018.
URL [https://arxiv.org/abs/1711.00165](https://arxiv.org/abs/1711.00165).
- MacKay (1996)
David JC MacKay.
Bayesian non-linear modeling for the prediction competition.
In *Maximum Entropy and Bayesian Methods: Santa Barbara, California, USA, 1993*, pp. 221–234. Springer, 1996.
- Mei et al. (2018)
Song Mei, Andrea Montanari, and Phan-Minh Nguyen.
A mean field view of the landscape of two-layer neural networks.
*Proceedings of the National Academy of Sciences*, 115(33):E7665–E7671, 2018.
- Mei et al. (2019)
Song Mei, Theodor Misiakiewicz, and Andrea Montanari.
Mean-field theory of two-layers neural networks: dimension-free bounds and kernel limit.
In *Conference on learning theory*, pp. 2388–2464. PMLR, 2019.
- Montanari & Urbani (2025)
Andrea Montanari and Pierfrancesco Urbani.
Dynamical decoupling of generalization and overfitting in large two-layer networks.
*arXiv preprint arXiv:2502.21269*, 2025.
- Nam et al. (2024)
Yoonsoo Nam, Nayara Fonseca, Seok Hyeong Lee, Chris Mingard, Niclas Goring, Ouns El Harzli, Abdurrahman Hadi Erturk, Soufiane Hayou, and Ard A Louis.
Disentangling rich dynamics from feature learning: A framework for independent measurements.
*arXiv preprint arXiv:2410.04264*, 2024.
- Pacelli et al. (2023)
Rosalba Pacelli, Sebastiano Ariosto, Mauro Pastore, Francesco Ginelli, Marco Gherardi, and Pietro Rotondo.
A statistical mechanics framework for bayesian deep neural networks beyond the infinite-width limit.
*Nature Machine Intelligence*, 5(12):1497–1507, 2023.
- Papyan et al. (2020)
Vardan Papyan, XY Han, and David L Donoho.
Prevalence of neural collapse during the terminal phase of deep learning training.
*Proceedings of the National Academy of Sciences*, 117(40):24652–24663, 2020.
- Petrini et al. (2022)
Leonardo Petrini, Francesco Cagnetta, Eric Vanden-Eijnden, and Matthieu Wyart.
Learning sparse features can lead to overfitting in neural networks.
*Advances in Neural Information Processing Systems*, 35:9403–9416, 2022.
- Rubin et al. (2024)
Noa Rubin, Inbar Seroussi, and Zohar Ringel.
Grokking as a first order phase transition in two layer networks.
*arXiv preprint arXiv:2310.03789*, May 2024.
URL [https://arxiv.org/abs/2310.03789](https://arxiv.org/abs/2310.03789).
Preprint.
- Rubin et al. (2025)
Noa Rubin, Kirsten Fischer, Javed Lindner, David Dahmen, Inbar Seroussi, Zohar Ringel, Michael Krämer, and Moritz Helias.
From kernels to features: A multi-scale adaptive theory of feature learning.
*arXiv preprint arXiv:2502.03210*, February 2025.
URL [https://arxiv.org/abs/2502.03210](https://arxiv.org/abs/2502.03210).
- Shi et al. (2022)
Zhenmei Shi, Junyi Wei, and Yingyu Liang.
A theoretical analysis on feature learning in neural networks: Emergence from inputs and advantage over fixed features.
*arXiv preprint arXiv:2206.01717*, 2022.
- Teh et al. (2015)
Yee Whye Teh, Alexandre Thiéry, and Sebastian Vollmer.
Consistency and fluctuations for stochastic gradient langevin dynamics, 2015.
URL [https://arxiv.org/abs/1409.0578](https://arxiv.org/abs/1409.0578).
- Tipping (2001)
Michael E Tipping.
Sparse bayesian learning and the relevance vector machine.
*Journal of machine learning research*, 1(Jun):211–244, 2001.
- Welling & Teh (2011)
Max Welling and Yee Whye Teh.
Bayesian learning via stochastic gradient langevin dynamics.
In *Proceedings of the 28th International Conference on Machine Learning (ICML 2011)*, ICML’11, pp. 681–688, Bellevue, Washington, USA, June 2011. Omnipress.
ISBN 978-1-4503-0619-5.
doi: 10.5555/3104482.3104568.
URL [https://dl.acm.org/doi/10.5555/3104482.3104568](https://dl.acm.org/doi/10.5555/3104482.3104568).

## Appendix A Additional Background

### A.1 Equivalence of Langevin Dynamics and Bayesian Inference

In this section, we provide a detailed derivation for the equivalence between the stationary distribution of a network trained with Stochastic Gradient Langevin Dynamics (SGLD) and the posterior distribution of a corresponding Bayesian model.

#### Network and Priors

We analyze a two-layer neural network with the functional form:

$$ $f({\bm{x}})=\frac{1}{N^{\gamma}}\sum_{i=1}^{N}a_{i}\phi({\bm{w}}_{i}^{\top}{\bm{x}}).$ (17) $$

The parameters are drawn from independent Gaussian priors, which corresponds to choosing a specific prior in a Bayesian model. The initializations are given by:

- •
Weights: $w_{ij}\sim\mathcal{N}(0,g_{w}^{2})$ with $g_{w}^{2}=\frac{\sigma_{w}}{d}$. This implies a prior probability $p({\bm{w}}_{i})\propto\exp(-\frac{1}{2g_{w}^{2}}\|{\bm{w}}_{i}\|^{2})$.
- •
Amplitudes: $a_{i}\sim\mathcal{N}(0,g_{a}^{2})$ with $g_{a}^{2}=\sigma_{a}^{2}$. This implies a prior probability $p(a_{i})\propto\exp(-\frac{1}{2g_{a}^{2}}\|a_{i}\|^{2})$.

Let ${\bm{\theta}}_{i}=({\bm{w}}_{i},a_{i})$ denote the parameters for the $i$-th neuron.

#### Stochastic Gradient Langevin Dynamics

We consider a full-batch Gradient Descent (GD) update for a parameter set $\theta_{i}$, which includes a weight decay term with coefficient $\gamma$ and injected isotropic Gaussian noise $\xi_{t}\sim\mathcal{N}(0,I)$. This algorithm is known as Stochastic Gradient Langevin Dynamics (SGLD):

$$ $\displaystyle\Delta\theta_{i,t}$ $\displaystyle:=\theta_{i,t+1}-\theta_{i,t}$ (18) $\displaystyle=-\eta\bigl(\gamma\,\theta_{i,t}+\nabla_{\theta_{i}}\ell(f_{{\bm{\theta}}})\bigr)\,+\sqrt{2T\eta}\ \,\xi_{t}\,.$ (19) $$

Here, $\eta$ is the learning rate, $T$ is a scalar temperature that controls the noise magnitude, and $\ell(f_{{\bm{\theta}}})=\frac{1}{P}\sum_{\mu=1}^{P}(y_{\mu}-f({\bm{x}}_{\mu}))^{2}$ is the mean squared error loss over the dataset of size $P$.

In the continuous-time limit ($\eta\rightarrow 0$), this discrete update equation corresponds to a Langevin stochastic differential equation. The stationary distribution of this process, reached as $t\rightarrow\infty$, is given by the Gibbs-Boltzmann distribution:

$$ $\displaystyle p({\bm{\theta}}_{i})\propto\exp\!\Bigl(-\frac{1}{T}\Bigl(\frac{\gamma}{2}\|\theta_{i}\|^{2}+\ell(f_{{\bm{\theta}}})\Bigr)\Bigr)$ (20) $$

We can separate the weight decay terms $\gamma$ into two separate parameters for weights and amplitudes, $\gamma_{w}$ and $\gamma_{a}$, respectively. The stationary distribution for the parameters of a single neuron $({\bm{w}}_{i},a_{i})$ is then:

$$ $\displaystyle p({\bm{w}}_{i},a_{i})\propto\exp\left(-\frac{1}{T}\left(\frac{\gamma_{w}}{2}\|{\bm{w}}_{i}\|^{2}+\frac{\gamma_{a}}{2}\|a_{i}\|^{2}+\frac{1}{P}\sum_{\mu=1}^{P}(y_{\mu}-f({\bm{x}}_{\mu}))^{2}\right)\right).$ (21) $$

#### Bayesian Posterior Distribution

From a Bayesian perspective, we aim to find the posterior distribution of the parameters given the data $\mathcal{D}=\{({\bm{x}}_{\mu},y_{\mu})\}_{\mu=1}^{P}$. The posterior is given by Bayes’ theorem: $p({\bm{\theta}}|\mathcal{D})\propto p(\mathcal{D}|{\bm{\theta}})p({\bm{\theta}})$.

- •
The prior $p({\bm{\theta}})$ is defined by our choice of initialization: $p({\bm{\theta}})=\prod_{i}p({\bm{w}}_{i})p(a_{i})\propto\exp\left(-\sum_{i}\left(\frac{1}{2g_{w}^{2}}\|{\bm{w}}_{i}\|^{2}+\frac{1}{2g_{a}^{2}}\|a_{i}\|^{2}\right)\right)$.
- •
The likelihood $p(\mathcal{D}|{\bm{\theta}})$ is chosen to be a Gaussian distribution with variance $\kappa^{2}$, corresponding to the mean squared error loss: $p(\mathcal{D}|{\bm{\theta}})\propto\exp\left(-\frac{1}{2\kappa^{2}P}\sum_{\mu=1}^{P}(y_{\mu}-f({\bm{x}}_{\mu}))^{2}\right)$.

Combining these, the log-posterior is proportional to the negative of an energy function $\mathcal{E}({\bm{\theta}},\mathcal{D})$. The posterior distribution for a single neuron’s parameters is:

$$ $\displaystyle p({\bm{w}}_{i},a_{i}|\mathcal{D})\propto\exp\left(-\left(\frac{1}{2g_{w}^{2}}\|{\bm{w}}_{i}\|^{2}+\frac{1}{2g_{a}^{2}}\|a_{i}\|^{2}+\frac{1}{2\kappa^{2}P}\sum_{\mu=1}^{P}(y_{\mu}-f({\bm{x}}_{\mu}))^{2}\right)\right).$ (22) $$

#### Identifying the Distributions

To establish the equivalence, we equate the functional forms of the SGLD stationary distribution and the Bayesian posterior. By comparing the exponents term-by-term, we find the following correspondences:

- •
Weights: $\frac{\gamma_{w}}{2T}\|{\bm{w}}_{i}\|^{2}=\frac{1}{2g_{w}^{2}}\|{\bm{w}}_{i}\|^{2}\implies\frac{\gamma_{w}}{T}=\frac{1}{g_{w}^{2}}$
- •
Amplitudes: $\frac{\gamma_{a}}{2T}\|a_{i}\|^{2}=\frac{1}{2g_{a}^{2}}\|a_{i}\|^{2}\implies\frac{\gamma_{a}}{T}=\frac{1}{g_{a}^{2}}$
- •
Loss Term: $\frac{1}{TP}\sum_{\mu}\ell_{\mu}=\frac{1}{2\kappa^{2}P}\sum_{\mu}\ell_{\mu}\implies T=2\kappa^{2}$

This implies that the SGLD algorithm with temperature $T$ effectively samples from a Bayesian posterior with data noise variance $\kappa^{2}=T/2$, and with prior variances $g_{w}^{2}=T/\gamma_{w}$ and $g_{a}^{2}=T/\gamma_{a}$.

#### Final Update Equations

By substituting these relations back into the SGLD update rules, we obtain the dynamics for sampling from the desired posterior:

$$ $\displaystyle\Delta a_{t,i}$ $\displaystyle=-\eta\left(\frac{T}{g_{a}^{2}}a_{t,i}+\nabla_{a_{i}}\left(\frac{1}{P}\sum_{\mu}\ell_{\mu}\right)\right)+\sqrt{2T\eta}\ \,\xi_{t,a}$ (23) $\displaystyle\Delta{\bm{w}}_{t,i}$ $\displaystyle=-\eta\left(\frac{T}{g_{w}^{2}}{\bm{w}}_{t,i}+\nabla_{{\bm{w}}_{i}}\left(\frac{1}{P}\sum_{\mu}\ell_{\mu}\right)\right)+\sqrt{2T\eta}\ \,\xi_{t,w}.$ (24) $$

For training, this means we must set the temperature to $T=2\kappa^{2}$.

### A.2 Derivation of equation 4

We start from the negative log posterior in Eq. equation [2](https://arxiv.org/html/2510.15174v1#S2.E2) and the two–layer model in Eq. equation [3](https://arxiv.org/html/2510.15174v1#S2.E3):

$$ $-\ln p_{\mathrm{GD}}({\bm{\theta}}\mid\mathcal{D})=\underbrace{\sum_{l}\frac{1}{2\sigma_{l}^{2}}\sum_{j}\|{\bm{\theta}}_{l,j}\|^{2}}_{-\ln p_{\text{prior}}}\;+\;\underbrace{\frac{1}{2\kappa^{2}P}\sum_{\mu=1}^{P}\bigl(f_{\bm{\theta}}({\bm{x}}_{\mu})-y_{\mu}\bigr)^{2}}_{-\ln p_{\text{L}}}.$ $$

For our two–layer network,
$f_{{\bm{\theta}}}({\bm{x}})=\frac{1}{N^{\gamma}}\sum_{i=1}^{N}a_{i}\,\phi({\bm{w}}_{i}^{\top}{\bm{x}})$.
Define the shorthand $\phi_{i\mu}:=\phi({\bm{w}}_{i}^{\top}{\bm{x}}_{\mu})$.
Then the data term expands as

$$ $\displaystyle-\ln p_{\text{L}}$ $\displaystyle=\frac{1}{2\kappa^{2}P}\sum_{\mu=1}^{P}\!\left(\frac{1}{N^{\gamma}}\sum_{i=1}^{N}a_{i}\phi_{i\mu}-y_{\mu}\right)^{2}$ $\displaystyle=\frac{1}{2\kappa^{2}P}\sum_{\mu=1}^{P}\!\left[\frac{1}{N^{2\gamma}}\sum_{i=1}^{N}\sum_{j=1}^{N}a_{i}a_{j}\,\phi_{i\mu}\phi_{j\mu}-\frac{2}{N^{\gamma}}\sum_{i=1}^{N}a_{i}\,\phi_{i\mu}y_{\mu}+y_{\mu}^{2}\right].$ (25) $$

Introduce the dataset–averaged quantities

$$ $\displaystyle\Sigma({\bm{w}})$ $\displaystyle:=\frac{1}{P}\sum_{\mu=1}^{P}\phi({\bm{w}}^{\top}{\bm{x}}_{\mu})^{2},$ $\displaystyle G({\bm{w}},{\bm{w}}^{\prime})$ $\displaystyle:=\frac{1}{P}\sum_{\mu=1}^{P}\phi({\bm{w}}^{\top}{\bm{x}}_{\mu})\,\phi({\bm{w}}^{\prime\top}{\bm{x}}_{\mu}),$ $\displaystyle J_{\mathcal{Y}}({\bm{w}})$ $\displaystyle:=\frac{1}{P}\sum_{\mu=1}^{P}\phi({\bm{w}}^{\top}{\bm{x}}_{\mu})\,y_{\mu}.$ (26) $$

Using $\sum_{i,j}=\sum_{i=j}+\sum_{i\neq j}$ in equation [25](https://arxiv.org/html/2510.15174v1#A1.E25), we obtain

$$ $\displaystyle-\ln p_{\text{L}}$ $\displaystyle=\frac{1}{2\kappa^{2}N^{2\gamma}}\left[\sum_{i=1}^{N}a_{i}^{2}\,\underbrace{\frac{1}{P}\sum_{\mu}\phi_{i\mu}^{2}}_{\Sigma({\bm{w}}_{i})}+\sum_{i\neq j}a_{i}a_{j}\,\underbrace{\frac{1}{P}\sum_{\mu}\phi_{i\mu}\phi_{j\mu}}_{G({\bm{w}}_{i},{\bm{w}}_{j})}\right]$ $\displaystyle\qquad-\frac{1}{\kappa^{2}N^{\gamma}}\sum_{i=1}^{N}a_{i}\,\underbrace{\frac{1}{P}\sum_{\mu}\phi_{i\mu}y_{\mu}}_{J_{\mathcal{Y}}({\bm{w}}_{i})}\;+\;\underbrace{\frac{1}{2\kappa^{2}P}\sum_{\mu=1}^{P}y_{\mu}^{2}}_{\text{const.}}.$ (27) $$

The prior part for our parameterization $w_{ij}\sim\mathcal{N}(0,\sigma_{w}^{2}/d)$ and $a_{i}\sim\mathcal{N}(0,\sigma_{a}^{2})$ is

$$ $-\ln p_{\text{prior}}=\frac{1}{2\sigma_{a}^{2}}\sum_{i=1}^{N}a_{i}^{2}+\frac{d}{2\sigma_{w}^{2}}\sum_{i=1}^{N}\|{\bm{w}}_{i}\|^{2}.$ (28) $$

Combining equation [27](https://arxiv.org/html/2510.15174v1#A1.E27) and equation [28](https://arxiv.org/html/2510.15174v1#A1.E28), and discarding the ${\bm{\theta}}$-independent constant, yields Eq. equation [4](https://arxiv.org/html/2510.15174v1#S2.E4):

$$ $\displaystyle-\ln\ p_{\mathrm{GD}}(\bm{W},{\bm{a}}\mid\mathcal{D})$ $\displaystyle=\frac{1}{2\sigma_{a}^{2}}\sum_{i=1}^{N}a_{i}^{2}+\frac{d}{2\sigma_{w}^{2}}\sum_{i=1}^{N}\|{\bm{w}}_{i}\|^{2}+\frac{1}{2\kappa^{2}N^{2\gamma}}\sum_{i=1}^{N}a_{i}^{2}\,\Sigma({\bm{w}}_{i})$ $\displaystyle\quad+\frac{1}{2\kappa^{2}N^{2\gamma}}\sum_{i\neq j}a_{i}a_{j}\,G({\bm{w}}_{i},{\bm{w}}_{j})-\frac{1}{\kappa^{2}N^{\gamma}}\sum_{i=1}^{N}a_{i}\,J_{\mathcal{Y}}({\bm{w}}_{i})+\text{const}.$ $$

This derivation holds for any nonlinearity $\phi$.

### A.3 k k -sparse parity target function

A $k$-sparse parity target function teacher is a single Walsh basis function on the Boolean hypercube.
Let $S\subseteq[d]$ with $|S|=k$. The Walsh function indexed by $S$ is

$$ $\chi_{S}({\bm{x}})\;=\;\prod_{j\in S}x_{j},\qquad{\bm{x}}\in\{\pm 1\}^{d},$ $$

(and, if ${\bm{x}}\in[-1,1]^{d}$, one may use $\chi_{S}({\bm{x}})=\prod_{j\in S}\mathrm{sign}(x_{j})$).
The family $\{\chi_{S}\}_{S\subseteq[d]}$ forms an orthonormal basis under the uniform product measure, i.e.
$\mathbb{E}[\chi_{S}({\bm{x}})\chi_{T}({\bm{x}})]=\mathbbm{1}\{S=T\}$.
In our experiments the teacher is $y({\bm{x}})=\chi_{S}({\bm{x}})$, when $|S|=k$ this is exactly the $k$-parity problem.

### A.4 Relation to kernel eigenvalue outliers

The transition from $m_{S}=0$ to $m_{S}>0$ manifests as an outlier eigenvalue in the learned kernel. Recall that $m_{S}=N^{1-\gamma}\langle aJ_{S}({\bm{w}})\rangle$ measures the mean alignment between neurons and the target mode, where $J_{S}({\bm{w}})=\mathbb{E}[\phi({\bm{w}}^{\top}{\bm{x}})\chi_{S}({\bm{x}})]$ quantifies how well a single neuron with weights ${\bm{w}}$ correlates with $\chi_{S}$. When $m_{S}$ becomes non-zero, it indicates that neurons have collectively aligned their weights to capture the target structure. Their $J_{S}({\bm{w}}_{i})$ values have grown large. This alignment directly impacts the empirical kernel through

$$ $\widehat{R}_{S}^{(a)}=\frac{\chi_{S}^{\top}K\chi_{S}}{\operatorname{tr}(K)}=\frac{\sum_{i=1}^{N}a_{i}^{2}\,J_{S}({\bm{w}}_{i})^{2}}{\sum_{i=1}^{N}a_{i}^{2}\,\Sigma({\bm{w}}_{i})},$ (29) $$

where $K_{\mu\nu}=\frac{1}{N}\sum_{i=1}^{N}\phi({\bm{w}}_{i}^{\top}{\bm{x}}_{\mu})\,\phi({\bm{w}}_{i}^{\top}{\bm{x}}_{\nu})$ is the learned kernel. The numerator $\sum_{i}a_{i}^{2}J_{S}({\bm{w}}_{i})^{2}$ grows quadratically with the alignment strengths $J_{S}({\bm{w}}_{i})$, while the denominator (trace) remains roughly constant. In the kernel regime ($m_{S}=0$), neurons remain randomly oriented so $J_{S}({\bm{w}}_{i})\sim O(N^{-1/2})$ and the ratio stays $O(1/P)$. However, when FL occurs ($m_{S}>0$), the enhanced $J_{S}({\bm{w}}_{i})$ values cause $\chi_{S}^{\top}K\chi_{S}$ to grow substantially, making $\chi_{S}$ an outlying eigendirection of $K$. This anisotropic deformation, from the isotropic NNGP to a low-rank, plus isotropic structure, provides a direct spectral signature of the FL phase transition.

### A.5 Connection to the standard NNGP formulation

The kernel representation in the function basis $\{\chi_{A}\}$ can be transformed to recover the standard NNGP formulation in input space. We start with the function-basis kernel

$$ $K_{AB}=\mathbb{E}_{{\bm{w}}\sim p_{\infty}}[J_{A}({\bm{w}})J_{B}({\bm{w}})],$ (30) $$

where $J_{A}({\bm{w}})=\mathbb{E}_{{\bm{x}}}[\phi({\bm{w}}^{\top}{\bm{x}})\chi_{A}({\bm{x}})]$ projects the neuron’s output onto basis function $\chi_{A}$. Expanding this definition:

$$ $\displaystyle K_{AB}$ $\displaystyle=\mathbb{E}_{{\bm{w}}}\left[\mathbb{E}_{{\bm{x}}}[\phi({\bm{w}}^{\top}{\bm{x}})\chi_{A}({\bm{x}})]\cdot\mathbb{E}_{{\bm{x}}^{\prime}}[\phi({\bm{w}}^{\top}{\bm{x}}^{\prime})\chi_{B}({\bm{x}}^{\prime})]\right]$ (31) $\displaystyle=\mathbb{E}_{{\bm{w}}}\mathbb{E}_{{\bm{x}},{\bm{x}}^{\prime}}\left[\phi({\bm{w}}^{\top}{\bm{x}})\phi({\bm{w}}^{\top}{\bm{x}}^{\prime})\chi_{A}({\bm{x}})\chi_{B}({\bm{x}}^{\prime})\right]$ (32) $\displaystyle=\mathbb{E}_{{\bm{x}},{\bm{x}}^{\prime}}\left[\chi_{A}({\bm{x}})\chi_{B}({\bm{x}}^{\prime})\cdot\underbrace{\sigma_{a}^{2}\mathbb{E}_{{\bm{w}}}[\phi({\bm{w}}^{\top}{\bm{x}})\phi({\bm{w}}^{\top}{\bm{x}}^{\prime})]}_{=:K({\bm{x}},{\bm{x}}^{\prime})}\right],$ (33) $$

where $K({\bm{x}},{\bm{x}}^{\prime})$ is the standard NNGP kernel in input space. The fixed-point equation $m_{A}=\frac{\sigma_{a}^{2}}{\kappa^{2}}(\Xi_{A}-\sum_{B}K_{AB}m_{B})$ with $\Xi_{A}=\sum_{S}y_{S}K_{AS}$ becomes

$$ $m=K(K+\tau I)^{-1}y,\qquad\tau=\kappa^{2}/\sigma_{a}^{2}.$ (34) $$

To see this explicitly, consider data points $\{{\bm{x}}_{\mu}\}_{\mu=1}^{P}$ with labels $y_{\mu}$. The predictor in the function basis is $f({\bm{x}})=\sum_{A}m_{A}\chi_{A}({\bm{x}})$, while in the input basis it becomes $f({\bm{x}}_{\mu})=\sum_{\nu=1}^{P}\alpha_{\nu}K({\bm{x}}_{\mu},{\bm{x}}_{\nu})$ where $\alpha=(K+\tau I)^{-1}y$. The equivalence follows from the change of basis: if $\chi_{A}({\bm{x}})=\delta_{{\bm{x}},{\bm{x}}_{A}}$ (point evaluation basis), then $K_{AB}=K({\bm{x}}_{A},{\bm{x}}_{B})$ directly recovers the Gram matrix. For general orthogonal bases, the kernel ridge regression solution remains invariant under this transformation.

### A.6 Derivation of the ARD precision fixed point

Recall the ARD prior on per-coordinate precisions $\bm{\rho}=(\rho_{1},\ldots,\rho_{d})$ and weights:

$$ $\displaystyle p({\bm{w}}\mid\bm{\rho})=\prod_{j=1}^{d}\mathcal{N}\!\big(w_{j}\mid 0,\rho_{j}^{-1}\big),\qquad p(\bm{\rho})=\prod_{j=1}^{d}\Gamma(\rho_{j}\mid\alpha_{0},\beta_{0}),$ (35) $$

and the single-neuron MF–ARD action ([Equation˜14](https://arxiv.org/html/2510.15174v1#S4.E14) in the main text)

$$ $\displaystyle-\ln p_{\mathrm{ARD}}({\bm{w}},a\mid\{m_{A}\},\bm{\rho},\mathcal{D})=\frac{a^{2}}{2\sigma_{a}^{2}}+\frac{1}{2}\sum_{j=1}^{d}\rho_{j}\,w_{j}^{2}+\frac{a^{2}}{2\kappa^{2}N^{2\gamma}}\Sigma({\bm{w}})-\frac{a}{\kappa^{2}N^{\gamma}}\!\left(J_{\mathcal{Y}}({\bm{w}})-\sum_{A}m_{A}J_{A}({\bm{w}})\right).$ (36) $$

For a width-$N$ two-layer network, the joint action is the sum over neurons $i=1,\dots,N$ of [Equation˜36](https://arxiv.org/html/2510.15174v1#A1.E36), and the (negative) log-evidence (free energy) for fixed $\{m_{A}\}$ is

$$ $\displaystyle\mathcal{F}(\bm{\rho};\{m_{A}\})=-\ln Z(\bm{\rho};\{m_{A}\})-\ln p(\bm{\rho})\quad\text{with}\quad Z(\bm{\rho};\{m_{A}\})=\!\!\int\!\!\prod_{i=1}^{N}\!d{\bm{w}}_{i}\,da_{i}\;e^{-\sum_{i=1}^{N}\!S_{\mathrm{ARD}}({\bm{w}}_{i},a_{i})}.$ (37) $$

Stationarity of the negative log-evidence w.r.t. $\rho_{j}$ gives the ARD FP:

$$ $\displaystyle 0=\frac{\partial\mathcal{F}}{\partial\rho_{j}}=\underbrace{\Big\langle\frac{\partial}{\partial\rho_{j}}\sum_{i=1}^{N}\!S_{\mathrm{ARD}}({\bm{w}}_{i},a_{i})\Big\rangle_{p_{\mathrm{ARD}}}}_{\text{energy term}}\;\underbrace{-\;\frac{N}{2}\,\frac{1}{\rho_{j}}}_{\text{Gaussian normalizer}}\;-\;\underbrace{\frac{\partial\ln p(\bm{\rho})}{\partial\rho_{j}}}_{\text{prior term}}=\frac{1}{2}\sum_{i=1}^{N}\big\langle w_{ij}^{2}\big\rangle_{p_{\mathrm{ARD}}}-\frac{N}{2}\frac{1}{\rho_{j}}-\frac{\alpha_{0}-1}{\rho_{j}}+\beta_{0},$ (38) $$

where we used $\partial S_{\mathrm{ARD}}/\partial\rho_{j}=\tfrac{1}{2}\sum_{i}w_{ij}^{2}$ and $\partial\ln\Gamma(\rho_{j}\mid\alpha_{0},\beta_{0})/\partial\rho_{j}=(\alpha_{0}-1)/\rho_{j}-\beta_{0}$. By MF symmetry, all neurons are i.i.d. under $p_{\mathrm{ARD}}$, so
$\sum_{i=1}^{N}\langle w_{ij}^{2}\rangle_{p_{\mathrm{ARD}}}=N\,\langle w_{j}^{2}\rangle_{p_{\mathrm{ARD}}}$, and [Equation˜38](https://arxiv.org/html/2510.15174v1#A1.E38) yields the closed-form FP update

$$ $\displaystyle\rho_{j}^{\star}\;=\;\frac{\alpha_{0}-1+\tfrac{N}{2}}{\beta_{0}+\tfrac{N}{2}\,\langle w_{j}^{2}\rangle_{p_{\mathrm{ARD}}}}\;\;\stackrel{{\scriptstyle\text{MAP corr.}}}{{\approx}}\;\;\frac{\alpha_{0}+\tfrac{N}{2}}{\beta_{0}+\tfrac{N}{2}\,\langle w_{j}^{2}\rangle_{p_{\mathrm{ARD}}}},$ (39) $$

where the MAP conjugacy correction (absorbing the $-1$ in $\alpha_{0}$) gives the form used in the main text ([Equation˜15](https://arxiv.org/html/2510.15174v1#S4.E15)). With the scale-matching choice $\beta_{0}=\alpha_{0}/d$, [Equation˜39](https://arxiv.org/html/2510.15174v1#A1.E39) becomes

$$ $\displaystyle\rho_{j}^{\star}\;=\;\frac{\alpha_{0}+\tfrac{N}{2}}{\tfrac{\alpha_{0}}{d}+\tfrac{N}{2}\,\langle w_{j}^{2}\rangle_{p_{\mathrm{ARD}}}}.$ (40) $$

### A.7 NNGP limit of the MF-ARD model

#### Usual (kernel) form

Let $C_{\rho}=\operatorname{diag}(\rho)^{-1}$ and $K_{\rho}(\mathbf{x},\mathbf{x}^{\prime})=\sigma_{a}^{2}\,\mathbb{E}_{\mathbf{w}\sim\mathcal{N}(0,C_{\rho})}[\phi(\mathbf{w}^{\top}\mathbf{x})\phi(\mathbf{w}^{\top}\mathbf{x}^{\prime})]$.
Then

$$ $\displaystyle m^{\infty}(\rho)=K_{\rho}(K_{\rho}+\tau I)^{-1}y,\quad K_{\rho}\chi_{A}=\lambda_{A}(\rho)\chi_{A}\ \Rightarrow\ m_{A}^{\infty}(\rho)=\frac{\lambda_{A}(\rho)}{\lambda_{A}(\rho)+\tau}\,y_{A}.$ (45) $$

This shows that $\rho$ rescales coordinates before the nonlinearity, so $K_{\rho}=\mathbb{E}_{\mathbf{z}}[J_{A}(C_{\rho}^{1/2}\mathbf{z})\,J_{B}(C_{\rho}^{1/2}\mathbf{z})]$ is a nonlinear deformation of the isotropic kernel. Any nontrivial change in $\rho$ therefore produces an $O(1)$ change in the Gaussian measure over $\mathbf{w}$, hence an $O(1)$ change in $K$ (it does not vanish with width). ARD improves spectral alignment by increasing $\lambda_{A}(\rho)$ along task-aligned directions. When $\lambda_{A}(\rho)$ crosses the scale $\tau$, $m_{A}^{\infty}$ exhibits a jump, yielding the leading-order FL effect at infinite width.

## Appendix B Additional figures

Figure: Figure 6: Plot of $\bm{W}\bm{W}^{\top}$ for a ReLU $N=512$ network, trained with SGLD on the same setting as [Figure˜3](https://arxiv.org/html/2510.15174v1#S3.F3).
Refer to caption: WWT.png

Figure: Figure 7: Same plot as [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5) a) SGLD, b) MF-ARD, but for a single-index model with Gaussian inputs ${\bm{x}}!\sim!\mathcal{N}(0,\mathbbm{1})$ and teacher $y=\mathrm{He}_{p}({\bm{w}}^{\top}{\bm{x}})$, where ${\bm{w}}_{j}=1/\sqrt{k}$ on a $k$-sized support and 0 otherwise, here $d=18$, $k=2$, $p=4$.
Refer to caption: x6.png

## Appendix C Related work: further details

Motivated by the observation that existing theories of deep learning are often neither simple nor mutually consistent, recent work seeks a unified account that (i) retains analytical simplicity, (ii) captures the essential structure of learning curves , and (iii) explains when and how finite networks escape the curse of dimensionality. The objective across these efforts is to develop theories that consistently cover both “lazy” kernel regimes and genuine FL. A recurring result is that a small number of control parameters, most notably an *effective learning rate*, organize the transition between regimes while preserving MF tractability. The emerging conclusion is that appropriately reduced MF or Bayesian descriptions can remain simple yet accurately capture kernel evolution, alignment to target structure, and the finite-width mechanisms by which deep networks surpass their infinite-width kernel limits.

As mentioned in the main text, attempts to bridge the theoretical gap between infinite-width kernel limits and finite-width FL largely fall into two categories. (1) Dynamical theories derive integro–differential or MF equations that track how representations, kernels, and losses evolve *during training* (Bordelon & Pehlevan, 2022; 2023; Mei et al., 2018; Montanari & Urbani, 2025; Celentano et al., 2025; Lauditi et al., 2025). (2) Static theories analyze the *learned posterior* after training, connecting finite-width learning to data-dependent kernel adaptation (Pacelli et al., 2023; Baglioni et al., 2024; Fischer et al., 2024; Rubin et al., 2025), or to explicit renormalization of the kernel (Howard et al., 2025; Aiudi et al., 2025).

#### Dynamical (training-time) theories

A central line of work develops dynamical MF theory (DMFT) for deep networks, showing that a single control parameter, the effective learning rate, governs whether training remains in the lazy NNGP/NTK regime or enters a FL regime with evolving kernels. Lauditi et al. (2025) formalize this for deep architectures, proving exact analytic solutions in the deep *linear* case and introducing numerical solvers for *nonlinear* activations. Comparisons on synthetic Gaussian data and real networks (including CIFAR tasks) demonstrate close agreement between DMFT predictions and observed training, confirming that the MF reduction captures genuine FL even beyond linear models.

Earlier, Bordelon & Pehlevan (2022) derived DMFT saddle-point equations intended for broad classes of infinite-width networks. These equations are solvable in closed form for deep linear networks, while nonlinear cases require numerical approximation; experiments on deep linear models and wide CNNs trained on two-class CIFAR tasks showed that both loss dynamics and feature-kernel evolution match theory across widths and scaling laws. The same line of work identifies the effective learning rate as the parameter driving faster training, larger kernel movement, and stronger alignment to target functions, with Bordelon & Pehlevan (2023) extending these predictions to *finite-width* networks. Complementing this, Fischer et al. (2024) provide a DMFT description for two-layer nonlinear networks trained by SGD, yielding closed equations for kernels and fields (exact for linear networks; Monte-Carlo DMFT for nonlinear activations). Their experiments on Gaussian-mixture single-index tasks and binary CIFAR subsets again validate the MF predictions and the regime split controlled by effective learning rate. Rubin et al. (2025) generalize the framework to *arbitrary depth*, retaining exact solutions for deep linear networks and approximate numerical solutions for nonlinear activations; experiments on Gaussian-mixture and single-index tasks and CIFAR subsets show that DMFT tracks real training behavior across regimes.

Orthogonal to DMFT but still dynamical, Mei et al. (2018) show that, in the infinite-width MF limit, SGD on two-layer networks follows a distributional-dynamics PDE in Wasserstein space. With noise or regularization this becomes a free-energy flow that guarantees global convergence. Unlike later DMFT work, they do not analyze NTK limits or a kernel–feature phase transition, instead they establish propagation-of-chaos, convergence of noisy SGD, and stability of fixed points in the noiseless case. Experiments on Gaussian classification and ReLU models show the PDE tracks SGD closely and diagnose when poor activations cause failure. In closely related simplified settings, Montanari & Urbani (2025) analyze *single-index* models with targets such as $h(z)=0.9z+z^{3}/6$. Here neurons become exchangeable samples from a common evolving distribution that aligns with the latent direction, the DMFT equations admit exact analytic solutions (as in linear networks), and test error converges to the Bayes limit, with experiments confirming the theory’s accuracy. Additional high-dimensional analyses further develop these dynamics (Celentano et al., 2025).

#### Static (posterior-time) theories

A complementary direction studies the learned posterior and reframes FL as *kernel adaptation*. Pacelli et al. (2023) reduce FL to *layerwise rescaling* of the kernel. Because the Bayesian predictor involves a matrix inverse, these rescalings act as mode-specific changes to the effective regularization of every kernel eigenmode, without rotating or performing relative scaling of feature directions. They derive analytic saddle-point equations for a pair of scalar descriptors per layer and, by solving them, obtain test-error predictions. On MNIST and CIFAR-10, where fully connected networks are not believed to do in strong representation learning, this coarse description already matches the observed bias–variance trade-off and the generalization error of SGLD-trained networks. Pushing this simplification further, Baglioni et al. (2024) collapse all finite-width effects into a *single global rescaling* of the infinite-width kernel, the resulting Bayesian effective action closely matches numerical experiments across MNIST and CIFAR-10.

Static renormalization viewpoints make the same idea explicit. Howard et al. (2025) interpret learning as a Wilsonian renormalization of the kernel, while Aiudi et al. (2025) show that architectural inductive bias matters: fully connected networks can be captured by a scalar kernel rescaling, but convolutional networks induce a *spatially indexed local kernel*. Training then learns a matrix $Q^{*}_{ij}$ that reweights patch–patch correlations, enabling selective emphasis or suppression of local interactions, this richer, weight-sharing–enabled mechanism provides a genuine route to FL that can outperform infinite-width kernel counterparts.

#### Summary

Across both dynamical and static perspectives, a coherent picture emerges. DMFT and related dynamical theories identify an effective learning-rate control parameter that cleanly separates lazy NTK behavior from FL with evolving kernels, with exact solutions in linear (and certain single-index) settings and accurate numerical solvers for nonlinear networks (Bordelon & Pehlevan, 2022; 2023; Lauditi et al., 2025; Fischer et al., 2024; Rubin et al., 2025; Mei et al., 2018; Montanari & Urbani, 2025; Celentano et al., 2025). Static Bayesian and renormalization accounts show that much of finite-width generalization can be captured by low-dimensional kernel *rescalings*, global, layerwise, or spatially local, without tracking full feature rotations (Pacelli et al., 2023; Baglioni et al., 2024; Howard et al., 2025; Aiudi et al., 2025). Together, these results support a simple, consistent theory that reproduces the characteristic shape and offsets of learning curves and clarifies how finite networks can break the apparent curse of dimensionality via structured kernel adaptation and alignment.

## Appendix D Proofs

### D.1 Kernel limit

###### Theorem D.1 .

With $\gamma=1/2$ the infinite width limit has the following FP equations:

$$ $\displaystyle m_{A}^{\infty}=\frac{\sigma_{a}^{2}}{\kappa^{2}}\left(\mathbb{E}_{{\bm{w}}\sim p_{\infty}}[J_{A}({\bm{w}})J_{\mathcal{Y}}({\bm{w}})]-\sum_{B}\mathbb{E}_{{\bm{w}}\sim p_{\infty}}[J_{A}({\bm{w}})J_{B}({\bm{w}})]m_{B}^{\infty}\right).$ (46) $$

###### Proof.

The proof proceeds by taking the $N\to\infty$ limit of the self-consistency equations derived from the cavity method’s free energy.

###### Theorem D.1 .

###### Proof.

#### 1. Self-Consistency from Free Energy.

We begin with the MF free energy functional derived from the cavity method:

$$ $\mathcal{F}(\{m_{A}\})=\frac{1}{2\kappa^{2}}\sum_{A}(y_{A}-m_{A})^{2}-N\ln\mathcal{Z}_{1}(\{m_{A}\})$ (47) $$

where $\mathcal{Z}_{1}(\{m_{A}\})=\int d{\bm{w}}\,e^{-\bar{S}_{\text{eff}}^{(\infty)}({\bm{w}};\{m_{A}\})}$
is the single-neuron partition function. The effective action for a single weight vector ${\bm{w}}$, after integrating out the amplitude $a$, is:

$$ $\bar{S}_{\text{eff}}^{(\infty)}({\bm{w}};\{m_{A}\})=\frac{d}{2\sigma_{w}^{2}}\|{\bm{w}}\|^{2}+\frac{1}{2}\ln\left(\frac{1}{\sigma_{a}^{2}}+\frac{\Sigma({\bm{w}})}{N^{2\gamma}\kappa^{2}}\right)-\frac{\left(J_{\mathcal{Y}}({\bm{w}})-\sum_{A}m_{A}J_{A}({\bm{w}})\right)^{2}}{2\kappa^{4}N^{2\gamma}\left(\frac{1}{\sigma_{a}^{2}}+\frac{\Sigma({\bm{w}})}{N^{2\gamma}\kappa^{2}}\right)}$ (48) $$

The equilibrium state is found by the stationarity condition $\partial\mathcal{F}/\partial m_{A}=0$, which yields the self-consistency equation:

$$ $m_{A}=N^{1-\gamma}\left\langle\mu({\bm{w}})J_{A}({\bm{w}})\right\rangle_{p({\bm{w}}|\{m_{B}\})}$ (49) $$

where $p({\bm{w}}|\{m_{B}\})=\frac{1}{\mathcal{Z}_{1}}e^{-\bar{S}_{\text{eff}}^{(\infty)}}$ is the posterior distribution on a single neuron’s weights, and $\mu({\bm{w}})$ is the posterior mean of its amplitude $a$:

$$ $\mu({\bm{w}})=\frac{\frac{1}{\kappa^{2}N^{\gamma}}\left(J_{\mathcal{Y}}({\bm{w}})-\sum_{B}m_{B}J_{B}({\bm{w}})\right)}{\frac{1}{\sigma_{a}^{2}}+\frac{\Sigma({\bm{w}})}{N^{2\gamma}\kappa^{2}}}$ (50) $$

#### 2. The Infinite-Width Limit with Critical Scaling.

To obtain a non-trivial limit as $N\to\infty$, we must set the scaling to the critical value $\gamma=1/2$. For $\gamma>1/2$, the prefactor $N^{1-\gamma}\to 0$, decoupling the neurons and leading to $m_{A}\to 0$. For $\gamma<1/2$, the interaction term diverges. With $\gamma=1/2$, the terms of order $1/N$ in the effective action $\bar{S}_{\text{eff}}^{(\infty)}$ vanish. Specifically:

$$ $\frac{\Sigma({\bm{w}})}{N^{2\gamma}\kappa^{2}}=\frac{\Sigma({\bm{w}})}{N\kappa^{2}}\xrightarrow{N\to\infty}0$ (51) $$

As a result, the data-dependent and $m_{A}$-dependent terms in the exponent of $p({\bm{w}}|\{m_{A}\})$ vanish. The distribution over weights collapses to its prior:

$$ $p({\bm{w}}|\{m_{A}\})\xrightarrow{N\to\infty}p_{\infty}({\bm{w}})\propto\exp\left(-\frac{d}{2\sigma_{w}^{2}}\|{\bm{w}}\|^{2}\right)$ (52) $$

Simultaneously, the denominator in $\mu({\bm{w}})$ simplifies to $1/\sigma_{a}^{2}$. The expression for the posterior mean amplitude becomes:

$$ $\mu({\bm{w}})\xrightarrow{N\to\infty}\frac{\sigma_{a}^{2}}{\kappa^{2}N^{1/2}}\left(J_{\mathcal{Y}}({\bm{w}})-\sum_{B}m_{B}J_{B}({\bm{w}})\right)$ (53) $$

#### 3. Deriving the Linear System.

We now substitute these limiting forms back into the self-consistency equation, using $m_{A}^{\infty}$ to denote the solution in this limit:

$$ $\displaystyle m_{A}^{\infty}$ $\displaystyle=N^{1-1/2}\left\langle\left[\frac{\sigma_{a}^{2}}{\kappa^{2}N^{1/2}}\left(J_{\mathcal{Y}}({\bm{w}})-\sum_{B}m_{B}^{\infty}J_{B}({\bm{w}})\right)\right]J_{A}({\bm{w}})\right\rangle_{p_{\infty}({\bm{w}})}$ (54) $\displaystyle m_{A}^{\infty}$ $\displaystyle=\frac{\sigma_{a}^{2}}{\kappa^{2}}\left\langle\left(J_{\mathcal{Y}}({\bm{w}})-\sum_{B}m_{B}^{\infty}J_{B}({\bm{w}})\right)J_{A}({\bm{w}})\right\rangle_{p_{\infty}({\bm{w}})}$ (55) $\displaystyle m_{A}^{\infty}$ $\displaystyle=\frac{\sigma_{a}^{2}}{\kappa^{2}}\left(\mathbb{E}_{{\bm{w}}\sim p_{\infty}}[J_{A}({\bm{w}})J_{\mathcal{Y}}({\bm{w}})]-\sum_{B}\mathbb{E}_{{\bm{w}}\sim p_{\infty}}[J_{A}({\bm{w}})J_{B}({\bm{w}})]m_{B}^{\infty}\right)$ (56) $$

Let us define the NNGP kernel matrix $K$ with elements $K_{AB}:=\mathbb{E}_{{\bm{w}}\sim p_{\infty}}[J_{A}({\bm{w}})J_{B}({\bm{w}})]$ and the data-kernel coupling vector $\Xi$ with elements $\Xi_{A}:=\mathbb{E}_{{\bm{w}}\sim p_{\infty}}[J_{A}({\bm{w}})J_{\mathcal{Y}}({\bm{w}})]=\sum_{S}y_{S}K_{AS}$. The equation becomes a linear system for the vector $m^{\infty}$:

$$ $m^{\infty}=\frac{\sigma_{a}^{2}}{\kappa^{2}}(\Xi-Km^{\infty})=\frac{\sigma_{a}^{2}}{\kappa^{2}}(Ky-Km^{\infty})$ (57) $$

#### 4. Solution as Kernel Ridge Regression.

Rearranging the linear system, we get:

$$ $\displaystyle\left(I+\frac{\sigma_{a}^{2}}{\kappa^{2}}K\right)m^{\infty}$ $\displaystyle=\frac{\sigma_{a}^{2}}{\kappa^{2}}Ky$ (58) $\displaystyle\left(\frac{\kappa^{2}}{2\sigma_{a}^{2}}I+\frac{1}{2}K\right)m^{\infty}$ $\displaystyle=\frac{1}{2}Ky$ (59) $$

Let the ridge be $\tau=\kappa^{2}/(2\sigma_{a}^{2})$. The equation is $(\tau I+K)m^{\infty}=Ky$. Solving for $m^{\infty}$ gives the final KRR solution:

$$ $m^{\infty}=K(K+\tau I)^{-1}y$ (60) $$

This completes the proof. The solution depends only on the NNGP kernel $K$, which is fixed by the network architecture and priors, not the training data labels. This demonstrates the absence of FL in this specific infinite-width limit. ∎

### D.2 Integrating out a a

As the action is quadratic in $a$ we can integrate it out.

###### Theorem D.2 .

The distribution is given by

$$ $\displaystyle S_{\textnormal{MF}}({\bm{w}},\{m_{A}\})$ $\displaystyle=const.+\frac{1}{2}\ln\left(\frac{1}{\sigma_{a}^{2}}+\frac{\Sigma({\bm{w}})}{N^{2\gamma}\kappa^{2}}\right)-\frac{\left(J_{\mathcal{Y}}(\mathbf{w})-\sum_{A}m_{A}J_{A}(\mathbf{w})\right)^{2}}{2N^{2\gamma}\kappa^{4}\left(\frac{1}{\sigma_{a}^{2}}+\frac{\Sigma({\bm{w}})}{N^{2\gamma}\kappa^{2}}\right)}$ (61) $\displaystyle+\frac{d}{2\sigma_{w}^{2}}\|{\bm{w}}\|^{2}$ (62) $$

$$ $\displaystyle p(a|{\bm{w}})$ $\displaystyle=\mathcal{N}(\mu({\bm{w}}),\sigma^{2}({\bm{w}}))$ (63) $\displaystyle\sigma({\bm{w}})^{2}$ $\displaystyle=\frac{1}{2\alpha}=\left(\frac{1}{\sigma_{a}^{2}}+\frac{\Sigma({\bm{w}})}{N^{2\gamma}\kappa^{2}}\right)^{-1}$ (64) $\displaystyle\mu({\bm{w}})$ $\displaystyle=\sigma^{2}\beta=\frac{\beta}{2\alpha}=\frac{\frac{1}{N^{\gamma}\kappa^{2}}\left(J_{\mathcal{Y}}(\mathbf{w})-\sum_{A}m_{A}J_{A}(\mathbf{w})\right)}{\left(\frac{1}{\sigma_{a}^{2}}+\frac{\Sigma({\bm{w}})}{N^{2\gamma}\kappa^{2}}\right)}$ (65) $$

###### Proof.

We start with the standard Gaussian integral

$$ $\displaystyle\int dad{\bm{w}}e^{-[\alpha\cdot a^{2}-\beta\cdot a+c]}=\int d{\bm{w}}\sqrt{\frac{\pi}{\alpha}}e^{\frac{\beta^{2}}{4\alpha}-c}$ (66) $$

We have

$$ $\displaystyle\mathcal{S}_{\text{MF}}({\bm{w}},a,\{m_{A}\})$ $\displaystyle=\frac{1}{2\sigma_{a}^{2}}\sum_{i=1}^{N}a_{i}^{2}+\frac{d}{2\sigma_{w}^{2}}\sum_{i=1}^{d}\|{\bm{w}}_{i}\|^{2}+\frac{a^{2}}{2\kappa^{2}N^{2\gamma}}\Sigma({\bm{w}})$ (67) $\displaystyle-\frac{a}{\kappa^{2}N^{\gamma}}\left(J_{\mathcal{Y}}(\mathbf{w})-\sum_{A}m_{A}J_{A}(\mathbf{w})\right)$ (68) $$

For 1 neuron we get

$$ $\displaystyle\mathcal{S}_{\text{MF}}({\bm{w}},a,\{m_{A}\})=a^{2}(\frac{1}{2\sigma_{a}^{2}}+\frac{\Sigma({\bm{w}})}{2\kappa^{2}N^{2\gamma}})-\frac{a}{\kappa^{2}N^{\gamma}}\left(J_{\mathcal{Y}}(\mathbf{w})-\sum_{A}m_{A}J_{A}(\mathbf{w})\right)+\frac{d}{2\sigma_{w}^{2}}\|{\bm{w}}\|^{2}$ (69) $$

with

$$ $\displaystyle\alpha$ $\displaystyle:=\frac{1}{2\sigma_{a}^{2}}+\frac{\Sigma(w)}{2N^{2\gamma}\kappa^{2}},$ (70) $\displaystyle\beta$ $\displaystyle:=\frac{J_{\mathcal{Y}}(w)-\sum_{A}m_{A}J_{A}(w)}{N^{\gamma}\kappa^{2}}.$ (71) $\displaystyle c$ $\displaystyle=\frac{d}{2\sigma_{w}^{2}}\|{\bm{w}}\|^{2}$ (72) $$

and $\alpha(a-\frac{\beta}{2\alpha})^{2}-\frac{\beta^{2}}{4\alpha}+c$ comparing to $\frac{-(a-\mu)^{2}}{2\sigma^{2}}$ we get $\alpha=\frac{1}{2\sigma^{2}}$ and $\mu=\frac{\beta}{2\alpha}$ for the Gaussian identification.
This gives

$$ $\displaystyle\int d{\bm{w}}\sqrt{\frac{\pi}{\alpha}}e^{\frac{\beta^{2}}{4\alpha}-c}=\int d{\bm{w}}e^{-[\frac{1}{2}\ln(\pi)-\frac{1}{2}\ln(\alpha)+\frac{\beta^{2}}{4\alpha}-c]}$ (73) $$

where the exponent is identified as the effective action.
∎

###### Theorem D.2 .

###### Proof.

### D.3 MF- FP equation

###### Lemma D.3 .

Consider a single target mode $y({\bm{x}})=\chi_{S}({\bm{x}})$ and assume other overlaps vanish. Using the self–consistency $m_{A}=N^{1-\gamma}\,\langle a\,J_{A}({\bm{w}})\rangle$ and the conditional mean $\mu({\bm{w}})$ above, we obtain

$$ $\displaystyle m_{S}\;=\;N^{1-\gamma}\,\Big\langle\mu({\bm{w}})\,J_{S}({\bm{w}})\Big\rangle\;=\;\frac{N^{1-2\gamma}}{\kappa^{2}}\,(1-m_{S})\;\Biggl\langle\frac{J_{S}({\bm{w}})^{2}}{\displaystyle\sigma_{a}^{-2}+\frac{\Sigma({\bm{w}})}{\kappa^{2}N^{2\gamma}}}\Biggr\rangle_{{\bm{w}}\sim p({\bm{w}}\mid m_{S})}.$ (74) $$

and in the $\kappa\rightarrow 0$ limit it is

$$ $\displaystyle m_{S}\;\approx\;(1-m_{S})\,N\,\Big\langle\frac{J_{S}({\bm{w}})^{2}}{\Sigma({\bm{w}})}\Big\rangle_{{\bm{w}}\sim p({\bm{w}}\mid m_{S})}$ (75) $$

###### Proof.

Consider a single target mode $y({\bm{x}})=\chi_{S}({\bm{x}})$ and assume other overlaps vanish. Using the self–consistency $m_{A}=N^{1-\gamma}\,\langle a\,J_{A}({\bm{w}})\rangle$ and the conditional mean $\mu({\bm{w}})$ above, we obtain

$$ $m_{S}\;=\;N^{1-\gamma}\,\Big\langle\mu({\bm{w}})\,J_{S}({\bm{w}})\Big\rangle\;=\;\frac{N^{1-2\gamma}}{\kappa^{2}}\,(1-m_{S})\;\Biggl\langle\frac{J_{S}({\bm{w}})^{2}}{\displaystyle\sigma_{a}^{-2}+\frac{\Sigma({\bm{w}})}{\kappa^{2}N^{2\gamma}}}\Biggr\rangle_{{\bm{w}}\sim p({\bm{w}}\mid m_{S})}.$ $$

In the small–noise regime $\kappa\to 0$ (with $\Sigma({\bm{w}})>0$), the denominator simplifies as
$\sigma_{a}^{-2}+\frac{\Sigma({\bm{w}})}{\kappa^{2}N^{2\gamma}}\approx\frac{\Sigma({\bm{w}})}{\kappa^{2}N^{2\gamma}}$,
so that

$$ $m_{S}\;\approx\;(1-m_{S})\,N\,\Big\langle\frac{J_{S}({\bm{w}})^{2}}{\Sigma({\bm{w}})}\Big\rangle_{{\bm{w}}\sim p({\bm{w}}\mid m_{S})}\quad\Longrightarrow\quad m_{S}\;\approx\;\frac{N\Big\langle J_{S}({\bm{w}})^{2}/\Sigma({\bm{w}})\Big\rangle}{1+N\Big\langle J_{S}({\bm{w}})^{2}/\Sigma({\bm{w}})\Big\rangle}\;$ $$

By Cauchy–Schwarz, $J_{S}({\bm{w}})^{2}\leq\Sigma({\bm{w}})$, so $0\leq m_{S}<1$ unless all mass concentrates on perfectly aligned ${\bm{w}}$.
∎

###### Lemma D.3 .

###### Proof.

### D.4 Solution to the fixed point equation

###### Lemma D.4 .

Consider $\phi=\textnormal{ReLU}$. Let the inputs $x_{j}\in\{\pm 1\}$ be i.i.d. and $y({\bm{x}})=\chi_{S}(\mathbf{x})=\prod_{j\in S}x_{j}$ with $S=\{0,1,...,k-1\}$.
We define $R_{k}({\bm{w}})=\frac{J_{S}({\bm{w}})^{2}}{\Sigma({\bm{w}})}$. It holds that

$$ $\displaystyle{\bm{w}}^{*}=\max_{{\bm{w}}}R_{k}({\bm{w}})$ (76) $$

is given by ${\bm{w}}^{*}=(\overbrace{\alpha,\ldots,\alpha}^{k},0,...,0)$.

###### Proof.

Decompose $w_{S}=\alpha\frac{\mathbbm{1}_{S}}{\sqrt{k}}+u,\quad u\perp\mathbbm{1}_{S}$ and keep arbitrary ${\bm{w}}^{C}$. Because the data distribution is invariant to permutations of the $k$ coordinates in $S$ the only $S$-dependent statistic that survives in $\chi_{S}$-weighted expectations is the sum $s({\bm{x}})=\sum_{j\in S}x_{j}$. Any component orthogonal to $\mathbbm{1}_{S}$ averages out in$J_{S}$ but increases $\Sigma({\bm{w}})$ (by convexity of $x\mapsto\phi(x)^{2}$ ). Likewise, weights in $S^{C}$ contribute variance to $\Sigma$ but contribute nothing to $J_{S}$
(they are independent of $\chi_{S}$.
and average to zero under the sign symmetry). Thus the maximizer of $R_{k}({\bm{w}})$ lives in the span of $\mathbbm{1}_{S}$ and ${\bm{w}}^{C}$.
∎

###### Lemma D.5 .

Consider $\phi=\textnormal{ReLU}$. Let the inputs $x_{j}\in\{\pm 1\}$ be i.i.d. and $y({\bm{x}})=\chi_{S}(\mathbf{x})=\prod_{j\in S}x_{j}$ with $|S|=k$. We assume ${\bm{w}}^{*}$ from [Lemma˜D.4](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem4). Then, the ratio of the squared neuron-target coupling $J_{S}(\mathbf{w})^{2}$ to the neuron’s self-energy $\Sigma(\mathbf{w})$ is a constant $R_{k}$ independent of the scale $\alpha$.

###### Proof.

Let $r\sim\text{Binomial}(k,1/2)$ be the number of components $x_{j}=-1$ for $j\in S$. The inner product is $s(\mathbf{x})=\mathbf{w}^{\top}\mathbf{x}=\alpha\sum_{j\in S}x_{j}=\alpha(k-2r)$, and the target function is $\chi_{S}(\mathbf{x})=(-1)^{r}$.

The neuron-target coupling $J_{S}(\mathbf{w})$ is the expectation $\mathbb{E}[\phi(\mathbf{w}^{\top}\mathbf{x})\chi_{S}(\mathbf{x})]$.

$$ $\displaystyle J_{S}(\mathbf{w})$ $\displaystyle=\mathbb{E}\left[\alpha\cdot[k-2r]_{+}\cdot(-1)^{r}\right]$ (77) $\displaystyle=\alpha\sum_{r=0}^{k}P(r)\cdot[k-2r]_{+}\cdot(-1)^{r}$ (78) $\displaystyle=\alpha\cdot 2^{-k}\sum_{r=0}^{\lfloor(k-1)/2\rfloor}\binom{k}{r}(k-2r)(-1)^{r}$ (79) $$

where $[z]_{+}=\max(0,z)$, and the sum is restricted to the terms where $k-2r>0$.

The neuron’s self-energy $\Sigma(\mathbf{w})$ is the expectation $\mathbb{E}[\phi(\mathbf{w}^{\top}\mathbf{x})^{2}]$.

$$ $\displaystyle\Sigma(\mathbf{w})$ $\displaystyle=\mathbb{E}\left[(\alpha\cdot[k-2r]_{+})^{2}\right]$ (80) $\displaystyle=\alpha^{2}\sum_{r=0}^{k}P(r)\cdot[k-2r]_{+}^{2}$ (81) $\displaystyle=\alpha^{2}\cdot 2^{-k}\sum_{r=0}^{\lfloor(k-1)/2\rfloor}\binom{k}{r}(k-2r)^{2}$ (82) $$

We define the scale-independent constants $D_{k}$ and $C_{k}$:

$$ $\displaystyle D_{k}$ $\displaystyle:=2^{-k}\sum_{r=0}^{\lfloor(k-1)/2\rfloor}\binom{k}{r}(k-2r)(-1)^{r}$ (83) $\displaystyle C_{k}$ $\displaystyle:=2^{-k}\sum_{r=0}^{\lfloor(k-1)/2\rfloor}\binom{k}{r}(k-2r)^{2}$ (84) $$

such that $J_{S}(\mathbf{w})=\alpha D_{k}$ and $\Sigma(\mathbf{w})=\alpha^{2}C_{k}$. The ratio is then

$$ $R_{k}:=\frac{J_{S}(\mathbf{w})^{2}}{\Sigma(\mathbf{w})}=\frac{(\alpha D_{k})^{2}}{\alpha^{2}C_{k}}=\frac{D_{k}^{2}}{C_{k}}$ $$

which is independent of $\alpha$, thus proving the proposition.

∎

###### Lemma D.4 .

###### Proof.

###### Lemma D.5 .

###### Proof.

### D.5 Exact solution of the FP equation

###### Theorem D.6 .

Using the setup from [Lemma˜D.5](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem5), espeically the proposed ${\bm{w}}^{*}$, the FP is given by

$$ $\displaystyle m_{S}=(1-m_{S})NR_{k}\left(1-\frac{\sigma_{a}^{-2}}{A^{\star}(m_{S})}\right),\quad A^{\star}(m_{S})=\frac{-C_{k}+\sqrt{C_{k}^{2}+4\left(\frac{dk}{\sigma_{w}^{2}}\right)N^{2\gamma}\,(1-m_{S})^{2}D_{k}^{2}\,\sigma_{a}^{-2}}}{2\left(\frac{dk}{\sigma_{w}^{2}}\right)\kappa^{2}N^{2\gamma}}$ (85) $$

###### Proof.

From [Lemma˜D.3](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem3), we know:

$$ $m_{S}=\frac{N^{1-2\gamma}}{\kappa^{2}}\,(1-m_{S})\;\left\langle\frac{J_{S}(\mathbf{w})^{2}}{\sigma_{a}^{-2}+\Sigma(\mathbf{w})/(\kappa^{2}N^{2\gamma})}\right\rangle_{\mathbf{w}|m_{S}}$ $$

Assuming the posterior distribution of weights is sharply peaked around the optimal direction given by the ansatz in [Lemma˜D.5](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem5), the expectation $\langle\cdot\rangle_{\mathbf{w}|m_{S}}$ reduces to an evaluation at the optimal weight scale $\alpha^{\star}$. The value of $\alpha^{\star}$ is determined by minimizing the single-neuron effective action $\mathcal{S}_{\text{eff}}(\alpha)$:

$$ $\mathcal{S}_{\text{eff}}(\alpha)=\underbrace{\frac{dk}{2\sigma_{w}^{2}}\,\alpha^{2}}_{\text{Weight Prior}}+\underbrace{\frac{1}{2}\ln A(\alpha)}_{\text{Normalizer}}-\underbrace{\frac{1}{2}\frac{\left[(1-m_{S})J_{S}(\alpha)\right]^{2}/(\kappa^{2}N^{2\gamma})}{A(\alpha)}}_{\text{Data Gain}}$ $$

where $A(\alpha)=\sigma_{a}^{-2}+\Sigma(\alpha)/(\kappa^{2}N^{2\gamma})=\sigma_{a}^{-2}+\alpha^{2}C_{k}/(\kappa^{2}N^{2\gamma})$.

Setting $\frac{\partial\mathcal{S}_{\text{eff}}}{\partial\alpha}=0$ yields a quadratic equation for the optimal value $A^{\star}=A(\alpha^{\star})$:

$$ $\left(\frac{dk}{\sigma_{w}^{2}}\right)\kappa^{2}N^{2\gamma}\,(A^{\star})^{2}+C_{k}\,A^{\star}-\frac{(1-m_{S})^{2}D_{k}^{2}\,\sigma_{a}^{-2}}{\kappa^{2}}=0$ $$

The positive root of this equation gives the solution for $A^{\star}(m_{S})$:

$$ $A^{\star}(m_{S})=\frac{-C_{k}+\sqrt{C_{k}^{2}+4\left(\frac{dk}{\sigma_{w}^{2}}\right)N^{2\gamma}\,(1-m_{S})^{2}D_{k}^{2}\,\sigma_{a}^{-2}}}{2\left(\frac{dk}{\sigma_{w}^{2}}\right)\kappa^{2}N^{2\gamma}}$ $$

We now substitute this back into the self-consistency equation. Using the definitions of $J_{S}$, $\Sigma$, and $A^{\star}$, we have $(\alpha^{\star})^{2}=\frac{\kappa^{2}N^{2\gamma}}{C_{k}}(A^{\star}-\sigma_{a}^{-2})$.

$$ $\displaystyle m_{S}$ $\displaystyle=(1-m_{S})\frac{N^{1-2\gamma}}{\kappa^{2}}\frac{J_{S}(\alpha^{\star})^{2}}{A^{\star}}$ (86) $\displaystyle=(1-m_{S})\frac{N^{1-2\gamma}}{\kappa^{2}}\frac{(\alpha^{\star})^{2}D_{k}^{2}}{A^{\star}}$ (87) $\displaystyle=(1-m_{S})\frac{N^{1-2\gamma}}{\kappa^{2}}\frac{1}{A^{\star}}\left[\frac{\kappa^{2}N^{2\gamma}}{C_{k}}(A^{\star}-\sigma_{a}^{-2})\right]D_{k}^{2}$ (88) $\displaystyle=(1-m_{S})N\frac{D_{k}^{2}}{C_{k}}\frac{A^{\star}-\sigma_{a}^{-2}}{A^{\star}}$ (89) $\displaystyle=(1-m_{S})NR_{k}\left(1-\frac{\sigma_{a}^{-2}}{A^{\star}(m_{S})}\right)$ (90) $$

This gives the final fixed-point equation for the order parameter $m_{S}$.
∎

###### Theorem D.7 .

Consider the fixed–point equation in the infinite $P$-limit

$$ $\displaystyle m_{S}\;=\;(1-m_{S})\,NR_{k}\!\left(1-\frac{\sigma_{a}^{-2}}{A^{\star}(m_{S})}\right),$ (91) $$

with $A^{\star}(m_{S})$ given implicitly as the positive root of

$$ $\Bigl(\tfrac{dk}{\sigma_{w}^{2}}\Bigr)\,\kappa^{2}N^{2\gamma}\,\big(A^{\star}\big)^{2}\;+\;C_{k}\,A^{\star}\;-\;(1-m_{S})^{2}D_{k}^{2}\,\sigma_{a}^{-2}\;=\;0,$ (92) $$

and define $C:=\tfrac{dk}{\sigma_{w}^{2}}$. Then the *critical noise* level

$$ $\displaystyle\quad\kappa_{c}^{2}\;=\;\frac{\sqrt{\,C_{k}^{2}+4\,C\,N^{2\gamma}D_{k}^{2}\sigma_{a}^{-2}\,}\,-\,C_{k}}{2\,C\,\sigma_{a}^{-2}\,N^{2\gamma}}$ (93) $$

marks a phase transition:

- (i)
If $\kappa^{2}\geq\kappa_{c}^{2}$, then $A^{\star}(0)\leq\sigma_{a}^{-2}$ and $m_{S}=0$ is a fixed point; in particular for $\kappa^{2}=\kappa_{c}^{2}$ we have $A^{\star}(0)=\sigma_{a}^{-2}$ and the only solution is $m_{S}=0$.
- (ii)
If $\kappa^{2}<\kappa_{c}^{2}$, then $A^{\star}(0)>\sigma_{a}^{-2}$ and there exists a unique nontrivial solution $m_{S}\in(0,1)$, thus the system exhibits symmetry breaking $m_{S}:0\to m_{S}>0$ as $\kappa$ crosses $\kappa_{c}$ from above.

###### Proof.

Set $m_{S}=0$ in equation [92](https://arxiv.org/html/2510.15174v1#A4.E92) to get

$$ $\displaystyle C\,\kappa^{2}N^{2\gamma}\,(A^{\star})^{2}+C_{k}A^{\star}-D_{k}^{2}\sigma_{a}^{-2}=0.$ (94) $$

At the onset of FL the trivial fixed point $m_{S}=0$ changes stability precisely when the right–hand side of the FP map ceases to vanish at $m_{S}=0$, i.e. when

$$ $\displaystyle NR_{k}\!\left(1-\frac{\sigma_{a}^{-2}}{A^{\star}(0)}\right)=0\quad\Longleftrightarrow\quad A^{\star}(0)=\sigma_{a}^{-2}.$ (95) $$

Plugging $A^{\star}=\sigma_{a}^{-2}$ into the quadratic and solving for $\kappa^{2}$ gives

$$ $\displaystyle C\,\kappa_{c}^{2}N^{2\gamma}\sigma_{a}^{-4}+C_{k}\sigma_{a}^{-2}-D_{k}^{2}\sigma_{a}^{-2}=0,$ (96) $$

which, after rearrangement, yields

$$ $\displaystyle\kappa_{c}^{2}=\frac{\sqrt{\,C_{k}^{2}+4\,C\,N^{2\gamma}D_{k}^{2}\sigma_{a}^{-2}\,}\,-\,C_{k}}{2\,C\,\sigma_{a}^{-2}\,N^{2\gamma}},$ (97) $$

the stated expression.

For $\kappa^{2}>\kappa_{c}^{2}$ the quadratic gives $A^{\star}(0)\leq\sigma_{a}^{-2}$, hence
$1-\sigma_{a}^{-2}/A^{\star}(0)\leq 0$ and the FP map evaluates to 0 at $m_{S}=0$, so $m_{S}=0$ is a (and in fact the only) solution. For $\kappa^{2}<\kappa_{c}^{2}$ we have $A^{\star}(0)>\sigma_{a}^{-2}$ so the FP map at $m_{S}=0$ is strictly positive, continuity and the fact that the map is strictly decreasing in $m_{S}$ (because $(1-m_{S})$ and $A^{\star}(m_{S})$ both decrease with $m_{S}$) imply a unique intersection with the diagonal in $(0,1)$, i.e. a unique $m_{S}\in(0,1)$ solves the FP equation.
∎

###### Theorem D.6 .

###### Proof.

###### Theorem D.7 .

###### Proof.

### D.6 How MF-ARD beats the curse of dimensionality

We use the following notation. Let $S\subset[d]$ with $|S|=k$ denote the (unknown) support. For a weight coordinate $w_{j}$ at a given outer iterate, write

$$ $\displaystyle v_{j}:=\langle w_{j}^{2}\rangle_{p},\qquad v^{\prime}_{j}(\rho):=\mathbb{E}_{p_{\rho}}[w_{j}^{2}]$ (98) $$

for the current and next second moments, respectively. For a vector of ARD precisions $\rho\in\mathbb{R}_{+}^{d}$ define the explicit ARD map

$$ $\rho_{j}(v_{j})\;=\;\frac{\alpha_{0}+\frac{N}{2}}{\frac{\alpha_{0}}{d}+\frac{N}{2}\,v_{j}}\;=:\;\frac{A}{B(v_{j})},\qquad A:=\alpha_{0}+\frac{N}{2},\ \ B(v):=\frac{\alpha_{0}}{d}+\frac{N}{2}v.$ (99) $$

We write $w_{-j}:=(w_{1},\dots,w_{j-1},w_{j+1},\dots,w_{d})$ for all coordinates except $j$ and write $p_{\rho}$ for the posterior $p_{\textnormal{ARD}}$ to make the $\rho$ dependence explicit.

#### Assumption

Here, we will state the assumption that is needed to prove the theorem: $\varepsilon$ symmetry breaking towards $S$. There exists an outer iterate $t_{0}=\mathcal{O}(1)$ and a constant $\varepsilon_{0}>0$ (independent of $d$) such that

$$ $\displaystyle\min_{j\in S}v_{j}^{\,t_{0}}\;-\;\max_{j\notin S}v_{j}^{\,t_{0}}\ \geq\ c\,\varepsilon_{0},\qquad v_{j}^{\,t}:=\langle w_{j}^{2}\rangle_{p}\text{ at outer time }t.$ (100) $$

We need to establish the following global bound.

###### Lemma D.8 .

Fix $j\in[d]$ and $w_{-j}$. Assume $\|x\|_{\infty}\leq 1$ (e.g. $x\in\{\pm 1\}^{d}$) and $\phi(z)=\max(0,z)$, as well as $m_{S}\leq 1$ Let

$$ $\displaystyle g(\bm{w})=\frac{1}{2}\ln\!\Big(\sigma_{a}^{-2}+\frac{\Sigma(\bm{w})}{\kappa^{2}N^{2\gamma}}\Big)-\frac{\big(J_{\mathcal{Y}}(\bm{w})-m_{S}J_{S}(\bm{w})\big)^{2}}{2\kappa^{4}N^{2\gamma}\big(\sigma_{a}^{-2}+\frac{\Sigma(\bm{w})}{\kappa^{2}N^{2\gamma}}\big)}.$ (101) $$

Then there exists $L_{\star}>0$, independent of $d$ and $w_{-j}$, such that the map $t\mapsto g(w_{-j},t)$ satisfies

$$ $\displaystyle\big|\partial_{t}^{2}g(w_{-j},t)\big|\ \leq\ L_{\star}\qquad\text{for a.e.\ }t\in\mathbb{R}.$ (102) $$

Consequently, for all $t\in\mathbb{R}$,

$$ $\displaystyle g(w_{-j},t)\ \geq\ g(w_{-j},0)-\frac{L_{\star}}{2}\,t^{2}.$ (103) $$

###### Proof.

Fix $j\in[d]$ and $w_{-j}$. Write $t:=w_{j}$ and, for each input $x$, set
$z(t,x):=\bm{w}^{\top}x=t\,x_{j}+c_{x}$ with $c_{x}:=w_{-j}^{\top}x_{-j}$.
Throughout, $\phi(z)=\max(0,z)$ and $\|x\|_{\infty}\leq 1$.

1. Derivative of $\Sigma$ and $J_{A}$

For $\Sigma(\bm{w})=\mathbb{E}_{x}[\phi(z)^{2}]=\mathbb{E}_{x}[z(t,x)^{2}\,\mathbf{1}\{z(t,x)>0\}]$ we have, for a.e. $t$,

$$ $\displaystyle\partial_{t}\Sigma\;=\;2\,\mathbb{E}_{x}\big[z\mathbf{1}_{\{z>0\}}\,x_{j}\big],\qquad\partial_{t}^{2}\Sigma\;=\;2\,\mathbb{E}_{x}\big[x_{j}^{2}\,\mathbf{1}_{\{z>0\}}\big]\;\leq\;2,$ (104) $$

because $|x_{j}|\leq 1$. By Cauchy–Schwarz,

$$ $|\partial_{t}\Sigma|\;\leq\;2\big(\mathbb{E}_{x}[z^{2}\mathbf{1}_{\{z>0\}}]\big)^{1/2}\big(\mathbb{E}_{x}[x_{j}^{2}\mathbf{1}_{\{z>0\}}]\big)^{1/2}\;\leq\;2\,\sqrt{\Sigma}.$ (105) $$

For $J_{A}(\bm{w})=\mathbb{E}_{x}[\phi(z)\,\chi_{A}(x)]$ (with $|\chi_{A}|\leq 1$), we get for a.e. $t$:

$$ $\displaystyle\partial_{t}J_{A}\;=\;\mathbb{E}_{x}[\mathbf{1}_{\{z>0\}}\,x_{j}\,\chi_{A}(x)]\,,\quad\partial_{t}^{2}J_{A}\;=\;0\ \text{ (a.e.)},\qquad|\partial_{t}J_{A}|\;\leq\;\mathbb{E}|x_{j}|\;\leq\;1,$ (106) $$

and by Cauchy–Schwarz again

$$ $|J_{A}|\;\leq\;\big(\mathbb{E}_{x}[\phi(z)^{2}]\big)^{1/2}\big(\mathbb{E}_{x}[\chi_{A}(x)^{2}]\big)^{1/2}\;=\;\sqrt{\Sigma}.$ (107) $$

2. Decompose $g$ and bound each second derivative

Let

$$ $\displaystyle a:=\sigma_{a}^{-2},\qquad c:=\frac{1}{\kappa^{2}N^{2\gamma}},\qquad J_{\Delta}:=J_{\mathcal{Y}}-m_{S}J_{S},\qquad q:=J_{\Delta}^{2},$ (108) $$

so

$$ $\displaystyle g(\bm{w})\;=\;\underbrace{\tfrac{1}{2}\log\big(a+c\Sigma\big)}_{=:g_{1}(\Sigma)}\;-\;\underbrace{\frac{q}{2\kappa^{4}N^{2\gamma}\,(a+c\Sigma)}}_{=:g_{2}(\Sigma,J_{\Delta})}.$ (109) $$

*Term 1:*
Set $h_{1}(s):=\tfrac{1}{2}\log(a+cs)$, so $h_{1}^{\prime}(s)=\frac{c}{2(a+cs)}$ and
$h_{1}^{\prime\prime}(s)=-\frac{c^{2}}{2(a+cs)^{2}}$. By the chain rule,

$$ $\displaystyle\partial_{t}^{2}g_{1}\;=\;h_{1}^{\prime\prime}(\Sigma)\,(\partial_{t}\Sigma)^{2}\;+\;h_{1}^{\prime}(\Sigma)\,\partial_{t}^{2}\Sigma.$ (110) $$

Using equation [105](https://arxiv.org/html/2510.15174v1#A4.E105) and $\partial_{t}^{2}\Sigma\leq 2$,

$$ $\displaystyle\Big|h_{1}^{\prime\prime}(\Sigma)\,(\partial_{t}\Sigma)^{2}\Big|\;\leq\;\frac{c^{2}}{2(a+c\Sigma)^{2}}\cdot 4\Sigma\;=\;2c^{2}\,\frac{\Sigma}{(a+c\Sigma)^{2}}\;\leq\;\frac{c\,\sigma_{a}^{2}}{2},$ (111) $$

where the last inequality follows by maximizing $u\mapsto\frac{u}{(a+cu)^{2}}$ at $u=a/c$.
Also $\big|h_{1}^{\prime}(\Sigma)\,\partial_{t}^{2}\Sigma\big|\leq\frac{c}{2(a+c\Sigma)}\cdot 2\leq c\,\sigma_{a}^{2}$.
Hence

$$ $\big|\partial_{t}^{2}g_{1}\big|\ \leq\ \frac{3}{2}\,c\,\sigma_{a}^{2}\,.$ (112) $$

*Term 2:*
Write $g_{2}=-\alpha\,\frac{q}{a+c\Sigma}$ with $\alpha:=\frac{1}{2\kappa^{4}N^{2\gamma}}$.
Differentiating twice and grouping terms gives (for a.e. $t$):

$$ $\displaystyle\partial_{t}^{2}g_{2}$ $\displaystyle=-\alpha\Bigg[\frac{q^{\prime\prime}}{a+c\Sigma}-\frac{2\,q^{\prime}\,c\,\Sigma^{\prime}}{(a+c\Sigma)^{2}}-\frac{q\,c\,\Sigma^{\prime\prime}}{(a+c\Sigma)^{2}}+\frac{2\,q\,(c\Sigma^{\prime})^{2}}{(a+c\Sigma)^{3}}\Bigg],$ $$

where primes denote $\partial_{t}$. We now bound $q,q^{\prime},q^{\prime\prime}$ with step 1 and using
equation [107](https://arxiv.org/html/2510.15174v1#A4.E107) and $|\partial_{t}J_{A}|\leq 1$.

- •
$q$: Using $|J_{\Delta}|\ \leq\ |J_{\mathcal{Y}}|+|m_{S}|\,|J_{S}|\ \leq\ (1+|m_{S}|)\sqrt{\Sigma}=:C_{\Delta}\sqrt{\Sigma}$ we get $q=J_{\Delta}^{2}\leq C_{\Delta}^{2}\Sigma$.
- •
$|q^{\prime}|$: We get
$\displaystyle|q^{\prime}|$
$\displaystyle=2|J_{\Delta}||\partial_{t}J_{\Delta}|\leq 2C_{\Delta}\sqrt{\Sigma}\cdot\big(|\partial_{t}J_{\mathcal{Y}}|+|m_{S}||\partial_{t}J_{S}|\big)\leq 2C_{\Delta}\sqrt{\Sigma}(1+|m_{S}|)$
$\displaystyle\leq 2C_{\Delta}^{2}\sqrt{\Sigma}$
- •
$|q^{\prime\prime}|$: We get $q^{\prime\prime}=2(\partial_{t}J_{\Delta})^{2}\leq 2C_{\Delta}^{2}$

Together with $\Sigma^{\prime}\leq 2\sqrt{\Sigma}$ and $\Sigma^{\prime\prime}\leq 2$, we get

$$ $\displaystyle\Big|\tfrac{q^{\prime\prime}}{a+c\Sigma}\Big|$ $\displaystyle\leq\frac{2C_{\Delta}^{2}}{a},$ (113) $\displaystyle\Big|\tfrac{2\,q^{\prime}\,c\,\Sigma^{\prime}}{(a+c\Sigma)^{2}}\Big|$ $\displaystyle\leq\frac{2\cdot 2C_{\Delta}^{2}\sqrt{\Sigma}\cdot c\cdot 2\sqrt{\Sigma}}{(a+c\Sigma)^{2}}=8C_{\Delta}^{2}c\,\frac{\Sigma}{(a+c\Sigma)^{2}}\ \leq\ \frac{2C_{\Delta}^{2}}{a},$ (114) $\displaystyle\Big|\tfrac{q\,c\,\Sigma^{\prime\prime}}{(a+c\Sigma)^{2}}\Big|$ $\displaystyle\leq\frac{C_{\Delta}^{2}\Sigma\cdot c\cdot 2}{(a+c\Sigma)^{2}}\ \leq\ \frac{C_{\Delta}^{2}}{a},$ (115) $\displaystyle\Big|\tfrac{2\,q\,(c\Sigma^{\prime})^{2}}{(a+c\Sigma)^{3}}\Big|$ $\displaystyle\leq\frac{2\,C_{\Delta}^{2}\Sigma\cdot c^{2}\cdot 4\Sigma}{(a+c\Sigma)^{3}}=8C_{\Delta}^{2}c^{2}\,\frac{\Sigma^{2}}{(a+c\Sigma)^{3}}\ \leq\ \frac{32C_{\Delta}^{2}}{27\,a},$ (116) $$

where in the last two lines we used that
$u\mapsto\frac{u}{(a+cu)^{2}}$ and $u\mapsto\frac{u^{2}}{(a+cu)^{3}}$ are maximized at
$u=\frac{a}{c}$ and $u=\frac{2a}{c}$, respectively, with finite maxima depending only on $a,c$.
Therefore $\big|\partial_{t}^{2}g_{2}\big|\leq\alpha\,K_{2}$ for a constant
$K_{2}$ depending only on $(a,c,m_{S})$, and independent of $d$, $w_{-j}$ and $t$.

3: Uniform bound on $\partial_{t}^{2}g$

Combining equation [112](https://arxiv.org/html/2510.15174v1#A4.E112) and the bound for $\partial_{t}^{2}g_{2}$, there exists

$$ $\displaystyle L_{\star}\;:=\;\frac{3}{2}\,c\,\sigma_{a}^{2}\;+\;\alpha\,K_{2}$ (117) $$

such that $|\partial_{t}^{2}g(w_{-j},t)|\leq L_{\star}$ for a.e. $t\in\mathbb{R}$. This proves the first claim.

4: Standard taylor expansion

For any twice–differentiable $h$ with $|h^{\prime\prime}|\leq L_{\star}$ a.e., the 1D Taylor inequality gives

$$ $h(t)\ \geq\ h(0)\ +\ h^{\prime}(0)\,t\ -\ \tfrac{L_{\star}}{2}\,t^{2}\qquad(\forall t\in\mathbb{R}).$ $$

Applying this to $h(t)=g(w_{-j},t)$ yields

$$ $g(w_{-j},t)\ \geq\ g(w_{-j},0)\ +\ \partial_{t}g(w_{-j},0)\,t\ -\ \frac{L_{\star}}{2}\,t^{2}.$ $$

In the parity setting and for $j\notin S$ (off–support), symmetry implies
$\partial_{t}g(w_{-j},0)=0$, giving the stated global quadratic lower bound
$g(w_{-j},t)\geq g(w_{-j},0)-\frac{L_{\star}}{2}t^{2}$.
∎

This lemma provides a precise bound on the conditional second moment of off-support $w_{j}$. The key insight is that when the precision parameter $\rho_{j}$ is large enough to dominate the coupling term, the second moment behaves essentially like that of a Gaussian with precision $\rho_{j}$, up to small corrections.

###### Lemma D.9 .

Let $j\notin S$ be an off-support coordinate and condition on $w_{-j}$. Given Lemma [D.8](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem8), there exists $L_{\star}>0$ (independent of $d$ and $w_{-j}$) such that the coupling function $t\mapsto g(w_{-j},t)$ satisfies:

$$ $|\partial_{t}^{2}g(w_{-j},t)|\leq L_{\star}\quad\text{for a.e. }t\in\mathbb{R}.$ $$

Assume furthermore the off-support symmetry condition:

$$ $\partial_{j}g(w_{-j},0)=0.$ $$

Consider the conditional distribution for $w_{j}$:

$$ $p(w_{j}\mid w_{-j},\rho)\propto\exp\big(-U_{j}(w_{j})\big),$ $$

where the potential is given by:

$$ $U_{j}(w_{j})=\frac{1}{2}\rho_{j}w_{j}^{2}+g(w_{j},w_{-j}).$ $$

If $\rho_{j}>L_{\star}$, then for every fixed $\theta\in(0,1]$, there exist constants $C,c>0$ (depending only on $L_{\star}$ and $\theta$, but independent of $d$ and $w_{-j}$) such that:

$$ $\mathbb{E}\big[w_{j}^{2}\mid w_{-j}\big]\leq\frac{e^{L_{\star}\theta^{2}}}{\rho_{j}-L_{\star}}+C\,e^{-c\,(\rho_{j}-L_{\star})\,\theta^{2}}.$ $$

In particular, for any $\varepsilon>0$, one can choose $\theta\in(0,1]$ small enough so that $e^{L_{\star}\theta^{2}}\leq 1+\varepsilon$, yielding:

$$ $\mathbb{E}[w_{j}^{2}\mid w_{-j}]\leq\frac{1+\varepsilon}{\rho_{j}-L_{\star}}+C\,e^{-c\,(\rho_{j}-L_{\star})}.$ $$

Hence, if $\rho_{j}=\Omega(d)$, then:

$$ $\mathbb{E}[w_{j}^{2}\mid w_{-j}]\leq\frac{1+o_{d}(1)}{\rho_{j}}$ $$

with the $o_{d}(1)$ term uniform in $w_{-j}$.

###### Proof.

The proof uses a careful splitting argument to handle the near-Gaussian behavior in the core region while controlling exponential tails. The key is to leverage the smoothness bound to create uniform upper and lower envelopes for the potential function.

Fix $j\notin S$ and condition on $w_{-j}$. For notational simplicity, write $w:=w_{j}$, $\rho:=\rho_{j}$, and define:

$$ $U(w)=\frac{1}{2}\rho w^{2}+g(w),$ $$

where $g(w):=g(w,w_{-j})$ is the coupling term.

Step 1: Establishing potential envelopes

By Lemma [D.8](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem8), we have $|g^{\prime\prime}(t)|\leq L_{\star}$ for almost every $t$, and the off-support symmetry gives $g^{\prime}(0)=0$.

Using Taylor expansion with integral remainder, for all $w\in\mathbb{R}$:

$$ $-\frac{L_{\star}}{2}w^{2}\leq g(w)-g(0)\leq\frac{L_{\star}}{2}w^{2}.$ (118) $$

Define the effective precision $m:=\rho-L_{\star}>0$. From equation [118](https://arxiv.org/html/2510.15174v1#A4.E118), we obtain two crucial envelopes:

Global lower envelope:

$$ $U(w)\geq U(0)+\frac{m}{2}w^{2}\quad\text{for all }w\in\mathbb{R}.$ (119) $$

Local upper envelope: For any $\theta>0$ and $|w|\leq\theta$:

$$ $U(w)\leq U(0)+\frac{\rho}{2}w^{2}+\frac{L_{\star}}{2}w^{2}\leq U(0)+\frac{m}{2}w^{2}+L_{\star}\theta^{2}.$ (120) $$

The global lower envelope ensures integrability, while the local upper envelope allows us to approximate the distribution by a Gaussian in the core region.

Step 2: Setting up the moment computation

Let $\mu$ be the conditional measure with density proportional to $e^{-U(w)}$. Define the normalization and second moment integrals:

$$ $Z:=\int_{\mathbb{R}}e^{-U(w)}\,dw,\qquad N:=\int_{\mathbb{R}}w^{2}e^{-U(w)}\,dw,$ $$

so that $\mathbb{E}_{\mu}[w^{2}]=N/Z$.

We partition $\mathbb{R}$ into the inside region $I:=\{|w|\leq\theta\}$ and the outside region $O:=\{|w|>\theta\}$, writing $Z=Z_{I}+Z_{O}$ and $N=N_{I}+N_{O}$.

Step 3: Bounding the inside region contribution

Using the local upper envelope equation [120](https://arxiv.org/html/2510.15174v1#A4.E120):

$$ $Z_{I}=\int_{|w|\leq\theta}e^{-U(w)}\,dw\geq e^{-U(0)}e^{-L_{\star}\theta^{2}}\int_{|w|\leq\theta}e^{-\frac{m}{2}w^{2}}\,dw.$ $$

Using the global lower envelope equation [119](https://arxiv.org/html/2510.15174v1#A4.E119):

$$ $N_{I}=\int_{|w|\leq\theta}w^{2}e^{-U(w)}\,dw\leq e^{-U(0)}\int_{|w|\leq\theta}w^{2}e^{-\frac{m}{2}w^{2}}\,dw.$ $$

Combining these bounds:

$$ $\frac{N_{I}}{Z}\leq\frac{N_{I}}{Z_{I}}\leq e^{L_{\star}\theta^{2}}\cdot\frac{\int_{|w|\leq\theta}w^{2}e^{-\frac{m}{2}w^{2}}\,dw}{\int_{|w|\leq\theta}e^{-\frac{m}{2}w^{2}}\,dw}\leq\frac{e^{L_{\star}\theta^{2}}}{m}.$ (121) $$

The final inequality follows because the ratio represents the truncated second moment of a zero-mean Gaussian with precision $m$, which is bounded by the untruncated value $1/m$.

Step 4: Controlling the outside region contribution

For the tail contribution, we use the global lower envelope equation [119](https://arxiv.org/html/2510.15174v1#A4.E119):

$$ $N_{O}=\int_{|w|>\theta}w^{2}e^{-U(w)}\,dw\leq e^{-U(0)}\int_{|w|>\theta}w^{2}e^{-\frac{m}{2}w^{2}}\,dw.$ $$

A standard Gaussian tail bound (obtained by integration by parts) gives, for $a>0$ and $\theta>0$:

$$ $\int_{\theta}^{\infty}t^{2}e^{-at^{2}}\,dt\leq\left(\frac{\theta}{2a}+\frac{1}{4a^{2}\theta}\right)e^{-a\theta^{2}}.$ $$

Applying this with $a=m/2$ and doubling for both tails:

$$ $\int_{|w|>\theta}w^{2}e^{-\frac{m}{2}w^{2}}\,dw\leq\frac{2}{m}\left(\theta+\frac{1}{m\theta}\right)e^{-\frac{m}{2}\theta^{2}}.$ (122) $$

Using the lower bound on $Z_{I}$ from Step 3 and equation [122](https://arxiv.org/html/2510.15174v1#A4.E122):

$$ $\frac{N_{O}}{Z}\leq\frac{N_{O}}{Z_{I}}\leq\frac{e^{L_{\star}\theta^{2}}}{c_{0}}\cdot\frac{2}{m}\left(\theta+\frac{1}{m\theta}\right)e^{-\frac{m}{2}\theta^{2}}\sqrt{m}\leq C_{1}e^{-c_{1}m\theta^{2}},$ $$

where $c_{0},C_{1},c_{1}>0$ are absolute constants depending only on $L_{\star}$ and $\theta$. The polynomial factors in $m$ are absorbed into the constant since the exponential term dominates for $m>0$.

Step 5: Combining the contributions

From equation [121](https://arxiv.org/html/2510.15174v1#A4.E121) and the outside region bound:

$$ $\mathbb{E}_{\mu}[w^{2}]=\frac{N}{Z}=\frac{N_{I}}{Z}+\frac{N_{O}}{Z}\leq\frac{e^{L_{\star}\theta^{2}}}{m}+Ce^{-cm\theta^{2}},$ $$

which establishes the main inequality with $m=\rho_{j}-L_{\star}$.

For the refined bound, choosing $\theta>0$ small enough so that $e^{L_{\star}\theta^{2}}\leq 1+\varepsilon$ gives the result.

Finally, if $\rho_{j}=\Omega(d)$, then $m=\rho_{j}-L_{\star}=\Omega(d)$, so the exponential tail term is $e^{-\Omega(d)}$ uniformly in $w_{-j}$, and:

$$ $\frac{1+\varepsilon}{\rho_{j}-L_{\star}}=\frac{1+\varepsilon}{\rho_{j}\left(1-\frac{L_{\star}}{\rho_{j}}\right)}=\frac{1+o_{d}(1)}{\rho_{j}}.$ $$

∎

The next lemma establishes that off-support $w_{j}$ contract towards equilibrium values of order $O(1/d)$ when initialized in a suitable bootstrap region. The key insight is that by controlling the initialization within $O(1/d)$, we can ensure the dynamics remain stable and converge.

###### Lemma D.10 .

Let $r>0$ be a radius parameter and let $L_{\star}$ be the global smoothness constant from Lemma [D.8](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem8). We define the following key quantities:

$$ $\displaystyle A$ $\displaystyle:=\alpha_{0}+\frac{N}{2},\qquad$ (total prior mass) (123) $\displaystyle B(v)$ $\displaystyle:=\frac{\alpha_{0}}{d}+\frac{N}{2}\cdot v,\qquad$ (effective local prior) (124) $$

For any threshold $K_{\star}>0$ and dimension $d\in\mathbb{N}$, we introduce the bootstrap parameters:

$$ $\displaystyle U_{d}(K_{\star})$ $\displaystyle:=\frac{\alpha_{0}+\frac{N}{2}K_{\star}}{d},\qquad$ (maximum local prior) (125) $\displaystyle a_{d}(r,K_{\star})$ $\displaystyle:=\frac{e^{L_{\star}r^{2}}\cdot\alpha_{0}}{A-L_{\star}\cdot U_{d}(K_{\star})},\qquad$ (affine drift term) (126) $\displaystyle b_{d}(r,K_{\star})$ $\displaystyle:=\frac{e^{L_{\star}r^{2}}\cdot(N/2)}{A-L_{\star}\cdot U_{d}(K_{\star})}.\qquad$ (contraction coefficient) (127) $$

Note that as $d\to\infty$, we have the asymptotic behavior:

$$ $b_{d}(r,K_{\star})\to b_{\infty}(r):=\frac{e^{L_{\star}r^{2}}\cdot(N/2)}{A}.$ $$

Choose $r>0$ sufficiently small so that $b_{\infty}(r)<1$. This ensures asymptotic contraction. Then there exist $d_{0}\in\mathbb{N}$, a finite threshold $K_{\star}>0$ (both independent of $d$), and positive constants $C,c,c^{\prime}>0$ such that for all $d\geq d_{0}$:

- 1.
One-step contraction bound: If $v_{\mathrm{off}}^{t}\leq K_{\star}/d$, then with $\rho_{\mathrm{off}}^{t}=A/B(v_{\mathrm{off}}^{t})$,
$\displaystyle v_{\mathrm{off}}^{t+1}$
$\displaystyle\leq\frac{e^{L_{\star}r^{2}}\cdot B(v_{\mathrm{off}}^{t})}{A-L_{\star}B(v_{\mathrm{off}}^{t})}+C\,e^{-c\,(\rho_{\mathrm{off}}^{t}-L_{\star})\,r^{2}}$
(128)
$\displaystyle\leq\frac{a_{d}(r,K_{\star})}{d}+b_{d}(r,K_{\star})\cdot v_{\mathrm{off}}^{t}+\eta_{d},$
(129)
where the exponential tail satisfies $\eta_{d}\leq Ce^{-c^{\prime}d}$.
- 2.
Bootstrap invariance: There exists a finite $K_{\star}$ such that
$v_{\mathrm{off}}^{t}\leq\frac{K_{\star}}{d}\quad\Longrightarrow\quad v_{\mathrm{off}}^{t+1}\leq\frac{K_{\star}}{d}\quad\text{for all }t\geq 0.$
In particular, any initialization with $v_{\mathrm{off}}^{0}\leq K_{\star}/d$ remains within the bootstrap region for all time.
- 3.
Equilibrium bound: Along any trajectory that remains in the bootstrap region,
$v_{\mathrm{off}}^{\star}\leq\frac{a_{d}(r,K_{\star})}{(1-b_{d}(r,K_{\star}))\cdot d}+O(e^{-c^{\prime}d})=\Theta\left(\frac{1}{d}\right).$

###### Proof.

The proof proceeds in four main steps: establishing a one-step bound, deriving an affine upper bound, proving bootstrap invariance, and analyzing convergence.

Fix $r>0$ and let $L_{\star}$ be the smoothness constant from Lemma [D.8](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem8). Consider an off-support coordinate at outer iteration $t$. We write $v:=v_{\mathrm{off}}^{t}$, $B(v)=\frac{\alpha_{0}}{d}+\frac{N}{2}v$, and $\rho:=\rho_{\mathrm{off}}^{t}=\frac{A}{B(v)}$ with $A=\alpha_{0}+\frac{N}{2}$ (see equation [99](https://arxiv.org/html/2510.15174v1#A4.E99)).

Throughout, we assume $v\leq K_{\star}/d$ for some threshold $K_{\star}>0$ to be determined.

Step 1: Establishing the one-step bound

We begin by applying the second moment bound from Lemma [D.9](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem9). For any $r>0$ and whenever $\rho>L_{\star}$, this yields:

$$ $v_{\mathrm{off}}^{t+1}=\mathbb{E}[w_{j}^{2}]\leq\frac{e^{L_{\star}r^{2}}}{\rho-L_{\star}}+C\,e^{-c\,(\rho-L_{\star})\,r^{2}}.$ (130) $$

Substituting $\rho=\frac{A}{B(v)}$, the main term becomes:

$$ $\frac{e^{L_{\star}r^{2}}}{\rho-L_{\star}}=e^{L_{\star}r^{2}}\cdot\frac{B(v)}{A-L_{\star}B(v)}.$ $$

Since we assume $v\leq K_{\star}/d$, we have the upper bound:

$$ $B(v)\leq U_{d}(K_{\star}):=\frac{\alpha_{0}+\frac{N}{2}K_{\star}}{d}=:\frac{C_{B}}{d},$ $$

which implies $\rho=\frac{A}{B(v)}\geq\frac{A\cdot d}{C_{B}}$. This lower bound on $\rho$ allows us to control the exponential tail:

$$ $C\,e^{-c\,(\rho-L_{\star})\,r^{2}}\leq C\,e^{-c^{\prime}d}=:\eta_{d}$ (131) $$

for some $c^{\prime}=c^{\prime}(A,C_{B},L_{\star},r)>0$ that depends only on the fixed parameters.

Step 2: Deriving the affine upper bound

To obtain a tractable recursion, we upper bound the main term using convexity. Define the function:

$$ $f(u):=\frac{u}{A-L_{\star}u},$ $$

which is convex on the interval $[0,A/L_{\star})$.

For all $u\in[0,U_{d}(K_{\star})]$, monotonicity of the denominator gives:

$$ $f(u)=\frac{u}{A-L_{\star}u}\leq\frac{u}{A-L_{\star}U_{d}(K_{\star})}.$ $$

Applying this with $u=B(v)$, we obtain the key affine upper bound:

$$ $e^{L_{\star}r^{2}}\cdot\frac{B(v)}{A-L_{\star}B(v)}\leq\frac{e^{L_{\star}r^{2}}}{A-L_{\star}U_{d}(K_{\star})}\left(\frac{\alpha_{0}}{d}+\frac{N}{2}v\right)=:\frac{a_{d}(r,K_{\star})}{d}+b_{d}(r,K_{\star})\cdot v,$ (132) $$

where we have defined:

$$ $\displaystyle a_{d}(r,K_{\star})$ $\displaystyle:=\frac{e^{L_{\star}r^{2}}\cdot\alpha_{0}}{A-\frac{L_{\star}C_{B}}{d}},$ (133) $\displaystyle b_{d}(r,K_{\star})$ $\displaystyle:=\frac{e^{L_{\star}r^{2}}\cdot(N/2)}{A-\frac{L_{\star}C_{B}}{d}}.$ (134) $$

Note that $a_{d}(r,K_{\star})=\Theta(1)$ and $b_{d}(r,K_{\star})\to b_{\infty}(r)=\frac{e^{L_{\star}r^{2}}\cdot(N/2)}{A}$ as $d\to\infty$.

Step 3: Establishing bootstrap invariance

The key is to choose parameters ensuring contraction. First, pick $r>0$ small enough so that:

$$ $b_{\infty}(r)=\frac{e^{L_{\star}r^{2}}\cdot(N/2)}{A}<1.$ $$

This is possible because $\frac{N/2}{A}<1$ (since $\alpha_{0}>0$) and $e^{L_{\star}r^{2}}\to 1$ as $r\to 0$.

With $r$ fixed, there exists $d_{0}$ such that for all $d\geq d_{0}$:

$$ $b_{d}(r,K_{\star})\leq\frac{b_{\infty}(r)+1}{2}<1.$ $$

The bound is uniform in $K_{\star}$ because $U_{d}(K_{\star})=O(1/d)$.

Combining equations equation [130](https://arxiv.org/html/2510.15174v1#A4.E130), equation [132](https://arxiv.org/html/2510.15174v1#A4.E132), and equation [131](https://arxiv.org/html/2510.15174v1#A4.E131), for $d\geq d_{0}$ and $v\leq K_{\star}/d$:

$$ $v_{\mathrm{off}}^{t+1}\leq\frac{a_{d}(r,K_{\star})}{d}+b_{d}(r,K_{\star})\cdot v+\eta_{d}.$ $$

Now choose $K_{\star}$ large enough so that:

$$ $K_{\star}\geq\sup_{d\geq d_{0}}\frac{a_{d}(r,K_{\star})}{1-b_{d}(r,K_{\star})}<\infty.$ $$

The supremum is finite because $a_{d}(r,K_{\star})=\Theta(1)$ and $1-b_{d}(r,K_{\star})$ is bounded away from zero for $d\geq d_{0}$.

Finally, enlarge $d_{0}$ if necessary so that $\eta_{d}\leq\frac{1}{2}(1-b_{d}(r,K_{\star}))\frac{K_{\star}}{d}$ for all $d\geq d_{0}$. This ensures bootstrap invariance:

$$ $v_{\mathrm{off}}^{t}\leq\frac{K_{\star}}{d}\quad\Longrightarrow\quad v_{\mathrm{off}}^{t+1}\leq\frac{K_{\star}}{d}.$ $$

Step 4: Convergence analysis

Within the invariant bootstrap region, the affine recursion:

$$ $v_{\mathrm{off}}^{t+1}\leq\frac{a_{d}(r,K_{\star})}{d}+b_{d}(r,K_{\star})\cdot v_{\mathrm{off}}^{t}+\eta_{d}$ $$

contracts towards its fixed point. Since $b_{d}(r,K_{\star})<1$, the equilibrium value satisfies:

$$ $v_{\mathrm{off}}^{\star}\leq\frac{a_{d}(r,K_{\star})}{(1-b_{d}(r,K_{\star}))\cdot d}+O(e^{-c^{\prime}d})=\Theta\left(\frac{1}{d}\right).$ $$

∎

###### Lemma D.11 .

Assume [Lemma˜D.10](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem10) 3. For all sufficiently large $d$,

$$ $\min_{j\in S}v_{j}^{\,t_{0}}\ \geq\ c\,\varepsilon_{0}+\max_{j\notin S}v_{j}^{\,t_{0}}\ \geq\ c\,\varepsilon_{0}-O(d^{-1})\ \geq\ \tfrac{c}{2}\varepsilon_{0}\ =\ \Theta(1).$ $$

Consequently, at (and after) time $t_{0}$, the ARD precisions satisfy $\rho_{\textrm{on}}=\Theta(1)$ while $\rho_{\textrm{off}}=\Theta(d)$.

###### Proof.

The first display follows directly from the assumption [D.6](https://arxiv.org/html/2510.15174v1#A4.SS6.SSS0.Px1) and $v_{\textrm{off}}^{\,t_{0}}=O(d^{-1})$ (bootstrap entry).
For $j\in S$, $B(v_{j})=\frac{\alpha_{0}}{d}+\frac{N}{2}v_{j}\geq\frac{N}{2}\cdot\tfrac{c}{2}\varepsilon_{0}$ for large $d$,
so $\rho_{j}=A/B(v_{j})=\Theta(1)$ uniformly in $d$. For $j\notin S$, by [Lemma˜D.10](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem10) 3), $v_{j}=\Theta(d^{-1})$,
hence $B(v_{j})=\Theta(d^{-1})$ and $\rho_{j}=\Theta(d)$.
∎

###### Lemma D.12 .

For a neuron with equal weights on $S$:

$$ $C_{\textnormal{MF}}=\frac{dk}{\sigma_{w}^{2}}\qquad\text{and}\qquad C_{\textnormal{ARD}}(t)=\sum_{j\in S}\rho_{j}^{\,t}=\Theta(k)\ \text{ for all }t\geq t_{0}\ \text{under \lx@cref{creftypecap~refnum}{lem:von-O1-weak}.}$ $$

###### Proof.

*MF:* With prior penalty $\frac{d}{2\sigma_{w}^{2}}\|\bm{w}\|^{2}$, along $\bm{w}=\alpha\,\mathbf{1}_{S}$ we get
$\frac{dk}{2\sigma_{w}^{2}}\alpha^{2}$, hence $C_{\textnormal{MF}}=\frac{dk}{\sigma_{w}^{2}}$. *ARD:* By [Lemma˜D.11](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem11),
$\rho_{j}^{\,t}=\Theta(1)$ for $j\in S$ and $t\geq t_{0}$, so $C_{\textnormal{ARD}}(t)=\Theta(k)$.
∎

###### Lemma D.13 .

At the symmetric initialization point $m_{S}=0$, the FL threshold satisfies:

$$ $\kappa_{c}^{2}=\frac{\sqrt{C_{k}^{2}+4CN^{2\gamma}D_{k}^{2}\sigma_{a}^{-2}}-C_{k}}{2C\sigma_{a}^{-2}N^{2\gamma}},$ (135) $$

where:

- •
$C$ is the quadratic curvature along the equal-weights $S$-direction,
- •
$C_{k},D_{k}$ are geometric constants depending only on the support size $k$,
- •
$N$ is the number of training samples, $\gamma$ is a scaling exponent,
- •
$\sigma_{a}$ controls the output noise level.

For large curvature $C$ (corresponding to strong regularization), the threshold simplifies to:

$$ $\kappa_{c}^{2}\sim\frac{|D_{k}|\sigma_{a}}{N^{\gamma}\sqrt{C}}=:\frac{\Lambda_{k}}{\sqrt{C}},$ (136) $$

where $\Lambda_{k}:=|D_{k}|\sigma_{a}/N^{\gamma}$.

###### Proof.

By Lemma [D.4](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem4), the optimal weight configuration is $\mathbf{w}^{*}=\alpha\mathbf{1}_{S}$. By Lemma [D.5](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem5), this gives scale-independent constants $C_{k}=\mathbb{E}[Z_{+}^{2}]$ and $D_{k}=\mathbb{E}[Z_{+}\chi_{S}(x)]$ such that:

$$ $\Sigma(\alpha)=C_{k}\alpha^{2},\qquad J_{S}(\alpha)=D_{k}\alpha,\qquad R_{k}=\frac{D_{k}^{2}}{C_{k}}.$ $$

From Theorem [D.7](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem7), the critical threshold occurs when $A^{*}(0)=\sigma_{a}^{-2}$, where $A^{*}$ solves the quadratic equation at $m_{S}=0$.

Step 2: Taylor expansion analysis

We analyze the one-dimensional effective potential along $\mathbf{w}=\alpha\mathbf{1}_{S}$:

$$ $U(\alpha)=\frac{C}{2}\alpha^{2}+g(\alpha),$ $$

where $C$ is the prior curvature (Lemma [D.12](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem12)) and $g(\alpha)$ is the coupling term from Lemma [D.8](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem8).

At onset ($m_{S}=0$), we have $J_{\Delta}=J_{S}=D_{k}\alpha$. The coupling term becomes:

$$ $g(\alpha)=\frac{1}{2}\log\big(a+c\Sigma(\alpha)\big)-\frac{J_{S}(\alpha)^{2}}{2\kappa^{4}N^{2\gamma}(a+c\Sigma(\alpha))},$ $$

with $a=\sigma_{a}^{-2}$, $c=\frac{1}{\kappa^{2}N^{2\gamma}}$, $\Sigma(\alpha)=C_{k}\alpha^{2}$, and $J_{S}(\alpha)=D_{k}\alpha$.

Step 3: Computing the Taylor expansion

Expanding around $\alpha=0$:

For the logarithmic term:

$$ $\frac{1}{2}\log\big(a+cC_{k}\alpha^{2}\big)=\frac{1}{2}\log a+\frac{cC_{k}}{2a}\alpha^{2}+O(\alpha^{4}).$ $$

For the rational term:

$$ $\frac{D_{k}^{2}\alpha^{2}}{2\kappa^{4}N^{2\gamma}(a+cC_{k}\alpha^{2})}=\frac{D_{k}^{2}\alpha^{2}}{2\kappa^{4}N^{2\gamma}a}+O(\alpha^{4})=\frac{\sigma_{a}^{2}D_{k}^{2}}{2\kappa^{4}N^{2\gamma}}\alpha^{2}+O(\alpha^{4}).$ $$

Therefore:

$$ $U(\alpha)=\text{const}+\frac{1}{2}\left[C+\frac{\sigma_{a}^{2}C_{k}}{\kappa^{2}N^{2\gamma}}-\frac{\sigma_{a}^{2}D_{k}^{2}}{\kappa^{4}N^{2\gamma}}\right]\alpha^{2}+O(\alpha^{4}).$ $$

Step 4: Onset condition and threshold

FL onset occurs when the quadratic coefficient vanishes:

$$ $C+\frac{\sigma_{a}^{2}C_{k}}{\kappa^{2}N^{2\gamma}}-\frac{\sigma_{a}^{2}D_{k}^{2}}{\kappa^{4}N^{2\gamma}}=0.$ $$

Multiplying by $\kappa^{4}N^{2\gamma}/\sigma_{a}^{2}$ and rearranging:

$$ $CN^{2\gamma}\sigma_{a}^{-2}\kappa^{4}+C_{k}\kappa^{2}-D_{k}^{2}=0.$ $$

Solving this quadratic in $\kappa^{2}$ gives:

$$ $\kappa_{c}^{2}=\frac{\sqrt{C_{k}^{2}+4CN^{2\gamma}D_{k}^{2}\sigma_{a}^{-2}}-C_{k}}{2C\sigma_{a}^{-2}N^{2\gamma}}.$ $$

Step 5: Asymptotic approximation

For large $C$, the square root is dominated by the term $4CN^{2\gamma}D_{k}^{2}\sigma_{a}^{-2}$, so:

$$ $\sqrt{C_{k}^{2}+4CN^{2\gamma}D_{k}^{2}\sigma_{a}^{-2}}\approx 2|D_{k}|\sigma_{a}^{-1}\sqrt{CN^{2\gamma}}.$ $$

This yields:

$$ $\displaystyle\kappa_{c}^{2}\sim\frac{2|D_{k}|\sigma_{a}^{-1}\sqrt{CN^{2\gamma}}}{2C\sigma_{a}^{-2}N^{2\gamma}}=\frac{|D_{k}|\sigma_{a}}{N^{\gamma}\sqrt{C}}=\frac{\Lambda_{k}}{\sqrt{C}}.$ (137) $$

∎

Now the proof of [Theorem˜4.1](https://arxiv.org/html/2510.15174v1#S4.Thmtheorem1) follows by the lemmata above.

###### Proof.

*MF:* Combine $C_{\textnormal{MF}}=\frac{dk}{\sigma_{w}^{2}}$ ([Lemma˜D.12](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem12)) with
equation [137](https://arxiv.org/html/2510.15174v1#A4.E137) to get $\kappa_{c,\textnormal{MF}}^{2}\asymp\Lambda_{k}\sqrt{\sigma_{w}^{2}/(dk)}=\Theta(\sqrt{1/(dk)})$.
*MF–ARD:* By [Lemma˜D.11](https://arxiv.org/html/2510.15174v1#A4.Thmtheorem11), for $t\geq t_{0}$ we have $C_{\textnormal{ARD}}(t)=\Theta(k)$, hence
$\kappa_{c,\textnormal{ARD}}^{2}\asymp\Lambda_{k}/\sqrt{k}=\Theta(\sqrt{1/k})$ from equation [137](https://arxiv.org/html/2510.15174v1#A4.E137).
∎

###### Lemma D.8 .

###### Proof.

###### Lemma D.9 .

###### Proof.

###### Lemma D.10 .

###### Proof.

###### Lemma D.11 .

###### Proof.

###### Lemma D.12 .

###### Proof.

###### Lemma D.13 .

###### Proof.

###### Proof.

## Appendix E Algorithms

### E.1 FP algorithm

#### Model.

Given $(X,y)$ with $X\in\{\pm 1\}^{P\times d}$ and $y\in\mathbb{R}^{P\times 1}$, we approximate the predictor by a particle ensemble,

$$ $\displaystyle f({\bm{x}})=s_{f}\sum_{b=1}^{B}a_{b}\,\phi({\bm{w}}_{b}^{\top}{\bm{x}}),\qquad s_{f}=\frac{N^{1-\gamma}}{B}.$ (138) $$

We *draw a single dataset of size $P$ once at initialization and keep it fixed for the entire run*. On the training set let $f_{p}=f({\bm{x}}_{p})$ and $r_{p}=y_{p}-f_{p}$. The Langevin temperature is fixed by the likelihood noise as $T\;=\;2\,\kappa^{2}$ ([Section˜A.1](https://arxiv.org/html/2510.15174v1#A1.SS1)). This choice makes SGLD asymptotically sample from the Bayesian posterior. Because all objectives and gradients are computed as empirical averages over this fixed sample (via $1/P$ factors), the dynamics naturally exhibit finite-$P$ fluctuations.

#### Sufficient statistics

We define the following low-rank statistics per particle $b$, with $z_{pb}={\bm{x}}_{p}^{\top}{\bm{w}}_{b}$:

$$ $\displaystyle C_{1,b}$ $\displaystyle=\sum_{p=1}^{P}\phi(z_{pb})\,r_{p},\qquad C_{2,b}=\sum_{p=1}^{P}\phi(z_{pb})^{2},$ (139) $\displaystyle G^{\text{data}}_{b}$ $\displaystyle=\;-\,\frac{2}{P}\sum_{p=1}^{P}\Big(r_{p}-s_{f}a_{b}\,\phi(z_{pb})\Big)\,\phi^{\prime}(z_{pb})\,\big(s_{f}a_{b}\big)\,x_{p}\ \in\mathbb{R}^{d}.$ (140) $$

These quantities are the only per-pass summaries we need to form gradients for $a_{b}$ and $w_{b}$.

#### Prior and SGLD potential

We impose a diagonal (ARD) Gaussian prior on the weights and an i.i.d. Gaussian prior on amplitudes:

$$ $\displaystyle E_{\text{prior}}=\tfrac{1}{2}\sum_{b=1}^{B}\rho^{\top}({\bm{w}}_{b}\odot{\bm{w}}_{b})+\tfrac{1}{2\sigma_{a}^{2}}\sum_{b=1}^{B}a_{b}^{2},\quad\mathcal{L}_{\text{data}}=\frac{1}{P}\sum_{p=1}^{P}(y_{p}-f_{p})^{2},$ (141) $$

and update parameters by Langevin dynamics on the potential

$$ $\displaystyle U(W,a)\;=\;T\,E_{\text{prior}}+\mathcal{L}_{\text{data}}.$ (142) $$

#### Gradients used by SGLD

From the streamed statistics we get closed-form gradients:

$$ $\displaystyle\nabla_{a_{b}}U=\frac{T}{\sigma_{a}^{2}}\,a_{b}\;-\;\frac{2\,s_{f}}{P}\,C_{1,b}\;+\;\frac{2\,s_{f}^{2}}{P}\,C_{2,b}\,a_{b},\qquad\nabla_{w_{b}}U=G^{\text{data}}_{b}\;+\;T\,(\rho\odot{\bm{w}}_{b}).$ (143) $$

These are the only quantities used inside the inner SGLD loop.

#### SGLD updates

We apply Euler–Maruyama steps with isotropic Gaussian noise:

$$ $\displaystyle w_{b}\leftarrow w_{b}-\eta\,\nabla_{w_{b}}U+\sqrt{2T\eta}\,\xi_{w_{b}},\qquad a_{b}\leftarrow a_{b}-\eta\,\nabla_{a_{b}}U+\sqrt{2T\eta}\,\xi_{a_{b}},$ (144) $$

where $\xi_{w_{b}}\!\sim\!\mathcal{N}(0,I_{d})$ and $\xi_{a_{b}}\!\sim\!\mathcal{N}(0,1)$. We use polynomial decay on $\eta$ always matching the SGLD-trained full NNs.

#### ARD update

The ARD update is:

$$ $\displaystyle\alpha_{\text{post}}=\alpha_{0}+\tfrac{B}{2},\qquad\beta_{\text{post}}=\beta_{0}+\tfrac{1}{2}\sum_{b=1}^{B}\|w_{b}\|_{2}^{2},\qquad\rho\leftarrow(1-\lambda)\,\rho+\lambda\,\frac{\alpha_{\text{post}}}{\beta_{\text{post}}}.$ (145) $$

Figure: Algorithm 1 Simple streaming SGLD for $(a,w)$ with optional ARD

#### Fixed-point view and the K K inner steps

The outer loop implements a fixed-point iteration on the network $f$. Writing $r=y-f$, define the map $\mathcal{G}_{K}$ as: (i) run $K$ inner SGLD steps on $({\bm{w}}_{b},a)$ using the current residual $r$; (ii) recompute $f^{\text{new}}(x_{p})=s_{f}\sum_{b}a_{b}\,\phi(w_{b}^{\top}x_{p})$. As $K\!\to\!\infty$ and $\eta\!\to\!0$ the inner Markov chain approaches its stationary law, and the iteration solves the MF FP equations described in the theory sections.

#### Empirical calculation of m S m_{S} and generalization error

Let the teacher be a single Walsh mode $\chi_{S}$, so $y({\bm{x}})=\chi_{S}({\bm{x}})$. On a held-out set $\{x_{\mu}\}_{\mu=1}^{P_{\mathrm{eval}}}$ define the vector $c\in\mathbb{R}^{P_{\mathrm{eval}}}$ by $c_{\mu}=\chi_{S}({\bm{x}}_{\mu})$, the empirical Gram scalar

$$ $g\;=\;\frac{1}{P_{\mathrm{eval}}}\,c^{\top}c,$ (146) $$

and the (scalar) empirical overlap

$$ $m_{S}\;=\;\frac{1}{P_{\mathrm{eval}}}\sum_{\mu=1}^{P_{\mathrm{eval}}}\chi_{S}({\bm{x}}_{\mu})\,f(x_{\mu})\;=\;\frac{1}{P_{\mathrm{eval}}}\,c^{\top}f.$ (147) $$

Let $\bar{f}^{2}=\tfrac{1}{P_{\mathrm{eval}}}\sum_{\mu}f({\bm{x}}_{\mu})^{2}$. Then the empirical test MSE decomposes as

$$ $\displaystyle\widehat{\mathcal{E}}_{\mathrm{test}}$ $\displaystyle=\frac{1}{2P_{\mathrm{eval}}}\sum_{\mu=1}^{P_{\mathrm{eval}}}\big(f({\bm{x}}_{\mu})-\chi_{S}({\bm{x}}_{\mu})\big)^{2}=\underbrace{\tfrac{1}{2}(1-m_{S})^{2}}_{\text{mode term}}\;+\;\underbrace{\tfrac{1}{2}\big(\bar{f}^{2}-2\,m_{S}^{2}+g\,m_{S}^{2}\big)}_{\text{noise / orthogonal term}}$ (148) $\displaystyle=\frac{1}{2}\Big(\bar{f}^{2}-2\,m_{S}+g\Big),$ (149) $$

where the second line is the direct empirical expression. When the Walsh basis is orthonormal on the evaluation set ($g=1$), this reduces to $\widehat{\mathcal{E}}_{\mathrm{test}}=\tfrac{1}{2}(\bar{f}^{2}-2m_{S}+1)$.

## Appendix F Training details (hyperparameters)

Here, we present the hyperparameters for SGLD and MF-ARD for the different figures. The hyperparameters for [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5) are specified below.

- •
Data in [Figure˜1](https://arxiv.org/html/2510.15174v1#S1.F1) is a slice of [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5) for $\kappa=5\cdot 10^{-3}$.
- •
Data in [Figure˜3](https://arxiv.org/html/2510.15174v1#S3.F3) is a slice of [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5) for $\kappa=5\cdot 10^{-3}$.
- •
[Figure˜4](https://arxiv.org/html/2510.15174v1#S3.F4) are averages over 3 trained models from the data of [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5) for $\kappa=5\cdot 10^{-3}$.

**Table 1: Algorithm outlined in [Section˜A.1](https://arxiv.org/html/2510.15174v1#A1.SS1). Hyperparamters for [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5) a).**
| Hyperparameter | SGLD |
| --- | --- |
| $d$ (input dimension) | 35 |
| $P$ (train set sizes) | $\{10,\,100,\,500,\,750,\,1000,\,2133,\,3666,\,5000,\,7500,\,10000\}$ |
| $\kappa$ values | $\{5{\times}10^{-4},\,10^{-3},\,5{\times}10^{-3},\,7.5{\times}10^{-3},\,10^{-2},\,5{\times}10^{-2},\,10^{-1}\}$ |
| $E$ (experiments / config) | 3 |
| teacher set S | $\{0,1,2,3\}$ |
| data distribution | $X\!\in\!\{-1,+1\}^{d}$, $y=\prod_{j\in S}X_{j}$ |
| activation | ReLU |
| $N$ (hidden units) | 512 |
| $\gamma$ (scaling exponent) | $0.5$ |
| $g_{w},\ g_{a}$ (prior variances) | $1.0,\ 1.0$ |
| initialization | $w\!\sim\!\mathcal{N}(0,g_{w}/d),\ \ a\!\sim\!\mathcal{N}(0,g_{a})$ |
| temperature $T$ | $2\,\kappa^{2}$ |
| epochs (max) | $7{,}500{,}000$ |
| batch size | full-batch |
| loss | mean MSE |
| learning rate $\eta$ (final) | $5{\times}10^{-4}$ |
| start LR $\eta_{\text{start}}$ | $5{\times}10^{-3}$ |
| LR scheduler | polynomial (power 2): $\eta_{\text{start}}\!\to\!\eta$ over $2{\times}10^{6}$ steps |

**Table 2: Algorithm outlined in [Section˜E.1](https://arxiv.org/html/2510.15174v1#A5.SS1). Hyperparamters for [Figure˜5](https://arxiv.org/html/2510.15174v1#S4.F5) b) with ARD disabled and c) with ARD enabled.**
| Hyperparameter | MF-ARD |
| --- | --- |
| $d$ (input dimension) | 35 |
| $P$ (train set sizes) | $\{10,100,500,750,1000,2133,3666,5000,7500,10000\}$ |
| $\kappa$ values | $\{5{\times}10^{-4},\,10^{-3},\,5{\times}10^{-3},\,7.5{\times}10^{-3},\,10^{-2},\,5{\times}10^{-2},\,10^{-1}\}$ |
| $E$ (experiments / config) | 3 |
| teacher set S | $\{0,1,2,3\}$ |
| data distribution | $X\!\in\!\{-1,+1\}^{d}$, $y=\prod_{j\in S}X_{j}$ |
| activation | ReLU |
| $B$ (particles) | 512 |
| $N$ | 512 |
| $\gamma$ | 0.5 |
| $\sigma_{a},\ \sigma_{w}$ | $1.0,\ 1.0$ |
| initialization | $w\!\sim\!\mathcal{N}(0,\sigma_{w}^{2}/d),\ a\!\sim\!\mathcal{N}(0,\sigma_{a}^{2})$ |
| outer steps | $7{,}500{,}000$ |
| learning rate scheduler | poly-2 decay: $5\times 10^{-3}\!\rightarrow\!5{\times}10^{-4}$ over $2{\times}10^{6}$ steps |
| SGLD inner steps $K$ | $K_{0}{=}12\ \to\ K_{\min}{=}2$ (decay over $6{\times}10^{5}$ steps) |
| temperature $T$ | $2\kappa^{2}$ |
| ARD | on: $\alpha_{0}{=}4.0$, $\beta_{0}{=}4/35$, EMA $0.5$, $\rho\!\in\![0,10^{18}]$ |

For the single index model we introduced a bias vector in the architecture.

**Table 3: Algorithm outlined in [Section˜A.1](https://arxiv.org/html/2510.15174v1#A1.SS1). Hyperparamters for [Figure˜7](https://arxiv.org/html/2510.15174v1#A2.F7) a).**
| Hyperparameter | SGLD (Hermite single-index) |
| --- | --- |
| $d$ (input dimension) | 18 |
| $P$ (train set sizes) | $\{75{,}000,50{,}000,25{,}000,10{,}000,5{,}000,1{,}000,100,50\}$ |
| $\kappa$ values | $\{10^{-5},10^{-4},10^{-3}10^{-2},10^{-1}\}$ |
| $E$ (experiments / config) | 4 |
| teacher type | single-index Hermite |
| Hermite degree $p$ | $4$ |
| support size $k$ | $2$ (random per experiment) |
| data distribution | $X\sim\mathcal{N}(0,I_{d})$ |
| teacher $w$ | $w_{i}=\frac{1}{\sqrt{k}}$ on support, else 0 |
| labels | $y=\mathrm{He}_{p}(Xw)$ |
| activation | ReLU |
| $N$ | 1024 |
| $\gamma$ | 0.5 |
| $\sigma_{a},\ \sigma_{w},\ \sigma_{b}$ | $1.0,\ 0.5,\ 1.0$ |
| initialization | $w\!\sim\!\mathcal{N}(0,\sigma_{w}^{2}/d),\ a\!\sim\!\mathcal{N}(0,\sigma_{a}^{2}),\ b\!\sim\!\mathcal{N}(0,\sigma_{b}^{2})$ |
| temperature $T$ | $2\,\kappa^{2}$ |
| epochs (max) | $4{,}000{,}000$ |
| batch size | full-batch |
| loss | mean MSE |
| learning rate $\eta$ (final) | $5{\times}10^{-4}$ |
| start LR $\eta_{\text{start}}$ | $2{\times}10^{-3}$ |
| LR scheduler | polynomial (power 2): $\eta_{\text{start}}\!\to\!\eta$ over $2{\times}10^{6}$ steps |

**Table 4: Algorithm outlined in [Section˜E.1](https://arxiv.org/html/2510.15174v1#A5.SS1). Hyperparamters for [Figure˜7](https://arxiv.org/html/2510.15174v1#A2.F7) b) with ARD enabled.**
| Hyperparameter | MF-ARD (Hermite single-index) |
| --- | --- |
| $d$ (input dimension) | 18 |
| $P$ (train set sizes) | $\{75{,}000,50{,}000,25{,}000,10{,}000,5{,}000,1{,}000,100,50\}$ |
| $\kappa$ values | $\{10^{-5},10^{-4},10^{-3}10^{-2},10^{-1}\}$ |
| $E$ (experiments / config) | 4 |
| teacher type | single-index Hermite |
| Hermite degree $p$ | $4$ |
| support size $k$ | $2$ (random per experiment) |
| data distribution | $X\sim\mathcal{N}(0,I_{d})$ |
| teacher $w$ | $w_{i}=\frac{1}{\sqrt{k}}$ on support, else 0 |
| labels | $y=\mathrm{He}_{p}(Xw)$ |
| activation | ReLU |
| $B$ (particles) | 1024 |
| $N$ | 1024 |
| $\gamma$ | 0.5 |
| $\sigma_{a},\ \sigma_{w},\ \sigma_{b}$ | $1.0,\ 0.5,\ 1.0$ |
| initialization | $w\!\sim\!\mathcal{N}(0,\sigma_{w}^{2}/d),\ a\!\sim\!\mathcal{N}(0,\sigma_{a}^{2}),\ b\!\sim\!\mathcal{N}(0,\sigma_{b}^{2})$ |
| outer steps | $4{,}000{,}000$ |
| learning rate scheduler | poly-2 decay: $2\times 10^{-3}\!\rightarrow\!5{\times}10^{-4}$ over $2.5{\times}10^{6}$ steps |
| SGLD inner steps $K$ | $K_{0}{=}12\ \to\ K_{\min}{=}2$ (decay over $6{\times}10^{5}$ steps) |
| temperature $T$ | $2\kappa^{2}$ |
| ARD | on: $\alpha_{0}{=}0.1$, $\beta_{0}{=}\alpha_{0}/d{=}0.1/18$, EMA $0.5$, $\rho\!\in\![0,10^{18}]$ |

## Peer Reviews

confidence: 4
contribution: 2
initial_rating: 4
presentation: 4
questions: Given my concern on the novelty (see "Contribution"), my first concrete
question to the authors is to check whether the previous theory by van
Meegen & Sompolinsky is able to explain their numerical findings or
whether there is still a gap between feature learning theory and
numerics that their theory is able to close.

My second question is more conceptual. I understand the heuristic motivation
of the approach. What I fail to see is why the proposed log likelihood (Eq. 14)
is the log of the stationary distribution of the Langevin training in (Eq. 1).
To me this seems like the attempt to capture heuristically the observations
from numerics rather than deriving this result from first principles.
soundness: 2
strengths: I very much like the form of presentation chosen by the authors, summarizing
main equations in boxes and providing explanations and intuition for the appearing order parameters. Also the authors throughout provide careful comparisons between theory and numerics.
summary: Summary

The manuscripts proposes a novel form of a heuristically
defined mean-field theory for feature learning in single
hidden layer networks in a setting of Bayesian inference.
The aim is to describe feature learning beyond previous works
of kernel adaptation and kernel renormalization.

The main novelty is the proposal of a heterogeneous prior distribution
on the first layer's weights. Deriving mean-field equations for a set of
order parameters, including the variances of the weights' priors, allows
the interpetation of the latter as a measure of relevance of input
features and describes settings where individual neurons specialize
on distinct input features.
Numerically trained networks with Langevin gradient descent with
weight decay show good agreement with the theory.

Soundness
Most of the manuscript is sound, in particular the first part of the
manuscript, until including section 3, which is equivalent to earlier
work on feature learning in a setting of Bayesian inference. The authors
nicely explain the mean-field approach also in intuitive terms, yet provide
all detail in the appendix.

The novel part starting in section 4. Here the main puzzling point to me
is that the heterogeneous variances and the hyper-prior for the variances
of individual neurons' weights are introduced ad-hoc. This is in contrast to
the first part, where the startionary distribution (Eq. (2)) is a direct
consequence of the applied Langevin training dynamics.

Presentation
The presentation of the work is very good and provides excellent explanations
that can also be followed intuitively, even for readers who do not follow
all mathematical details.

Contribution
The main contribution the work is to provide a simple ad-hoc (not strictly
derived from first principles) mean-field theory that explains specialization
of individual neurons to explain feature learning.
A point that needs to be clarified in the rebuttal is the novelty.

In particular in the light of previous work by van Meegen & Sompolinsky 2024
"Coding schemes in neural networks learning classification tasks" Nat Comm.
The latter work considers the same problem as the authors and derives a mean-field theory to explain the specialization of neurons on particular features.
The authors should, for example, look at Figure 8 of said paper. The theory
by van Meegen and Sompolinsky performs a mean-field approach where the readout
weights (variable "a" in the current manuscript) concentrate.
weaknesses: The missing mention of van Meegen & Sompolinksy 2024 as a tightly related
work casts some doubts about the novelty and contribution. In case this
point can be clarified in the rebuttal I am happy to raise my score.
(see section contribution in summary box).

---

confidence: 4
contribution: 1
initial_rating: 2
presentation: 2
questions: 1.  Why did the authors choose to introduce anisotropy explicitly through the ARD prior rather than studying finite-width ($1/N$ or $P/N$) corrections to the isotropic mean-field theory that might be still tractable?
These corrections are known to produce neuron correlations and spontaneous anisotropy without changing the generative model.

2. Can the authors clarify how they are taking the asymptotic limit (of $P$ data and $N$ width)? 
Of course these limits do not commute in practice, and to get a mean-field theory one usually takes $N\to \infty$ at fixed $P$, but Eqs. (6) and (7) have a population average. 

3. Fig.4 (b) shows that the theoretical model at small $P$ over-estimates the anisotropy of hidden neurons. Could the authors comment on that?

4. In Fig. 1a, the test MSE stays flat even with increased samples $P$ before dropping sharply. Why does learning only start at such large $P$? Could the authors comment on what would happen at larger width?  

Minor annotations:
- Lines 44-45: kernel renormalization theories are in the asymptotic limits of data and width so the sentence and citation is imprecise.
- Eq. (4) is missing a $+$ sign (typo).
- Line 160: fully factorized posterior for mean field theory is not an approximation, it’s exact in the limit $N\to \infty$.
- Line 249: $m_s$ mentioned before it has been introduced
soundness: 2
strengths: The paper - like many recent works at the intersection of statistical physics and machine learning - addresses a relevant theoretical question: how mean-field theory can be extended beyond the infinite-width, small-data limit to capture internal feature specialization and neuron-neuron correlations in finite-width networks. The proposed MF-ARD formulation is conceptually clear and analytically tractable. The manuscript is clearly written, and the authors provide numerical comparisons between their theoretical predictions and SGLD simulations, showing qualitative agreement.
summary: The paper proposes a *Mean-Field Automatic Relevance Determination* (MF-ARD) framework for two-layer Bayesian neural networks trained with SGLD. The goal is to model the “specialization transition” in fully connected networks at the interpolation regime, where just a finite subset of neurons align with the target function. Since standard mean-field theory (exact at infinite width $N\to \infty$ when sample $P=\mathcal{O}(1)$) yields a fully factorized posterior, the authors introduce a coordinate-dependent Gaussian prior with Gamma-distributed precisions (inverse variances controlling per-weight decay strength), enabling anisotropic shrinkage across weight coordinates. They derive fixed-point equations for these precisions and show, through comparisons with SGLD simulations, that the resulting model resemble the qualitative emergence of neuron specialization and sparse, task-aligned representations observed in finite-width networks.
weaknesses: **Theoretical inconsistency between MF-ARD and SGLD dynamics.**

This paper is meant to describe the stationary distribution of SGLD dynamics for two-layer neural networks through a mean-field approximation and **uniform** weight decay, but the proposed MF-ARD formulation fundamentally changes the underlying Bayesian model. By introducing coordinate-dependent Gaussian priors with Gamma-distributed precisions $\rho_j$ the authors define a *different* posterior from the isotropic one actually sampled by SGLD. As a consequence, I don’t see how their “MF-ARD theory’’ can be interpreted as a faithful large $P$ limit of the true SGLD dynamics: they are effectively comparing two different models.

**Imposed rather than spontaneous specialization transition in the theory.**

Since the MF-ARD model introduces anisotropy explicitly through the prior rather than deriving it as a spontaneous effect of finite-width fluctuations (i.e., $\mathcal{O}(P/N)$ corrections to mean-field theory), the “specialization transition’’ is imposed by construction and has nothing to do with the original posterior of Eq.(4) they intend to study. 

**Unrealistic neuron-wise sparsification.**

The transition described by the theory corresponds to a *hard sparsification*, where only a finite number of neurons acquire non-zero overlap $m_s>0$ with the teacher while the rest remain near initialization. This is fundamentally different to what happens in finite width neural networks during feature learning. In realistic dynamics, all neurons move substantially, but their updates become strongly correlated, so that the learned representation effectively spans a low-dimensional subspace.

**Restrictive assumptions on data and teacher.**

The theoretical analysis relies on a sparse teacher model, where the number of nonzero coefficients $m_A$ directly determines the number of “relevant” features or basis functions contributing to the target. In practice, this means the teacher vector is assumed to have support only on a small subset of input dimensions ($k$-sparse vector). While this setup is analytically convenient, it is not representative of realistic data or high-dimensional learning tasks. In real datasets, the relevant subspace is *not known a priori* and typically spans a continuum of directions with a decaying or heavy-tailed spectrum of features.  Given that the theoretical setup strongly relies on having $k$ finite, I don't see how it can be extended to practical scenarios. 

**Theory describes a numerically unstable training setup.**

The proposed theory corresponds to an impractical actual training of a neural network with SGLD. Indeed, running SGLD or SGD with neuron-dependent (or coordinate-dependent) weight decay would require maintaining and updating a distinct regularization coefficient for each neuron or parameter, effectively introducing $\mathcal{O}(N) additional hyperparameters. This would break permutation symmetry across neurons, complicates optimization dynamics, and would make training highly unstable. 

In summary, the paper starts from a standard SGLD framework but then shifts to analyzing a modified Bayesian model (MF-ARD) whose dynamics and posterior differ from those of the system it aims to describe. The resulting formulation is not a faithful large $P$ limit of a mean-field theory, nor does it capture the genuine feature-learning phenomenology observed in finite-width networks. For these reasons, I believe the paper requires major conceptual revision before it could be considered for publication. Its current scope and assumptions are limited to shallow networks and restricted to sparse targets, and the theoretical insights it provides are not clearly connected to modern learning dynamics in realistic settings.

---

confidence: 4
contribution: 2
initial_rating: 2
presentation: 3
questions: 1. As far as I can tell, the MF theory provided here is identical to the one appearing in Seroussi et al. (2023), Rubin et al. (2025), and Ringel et al. (2025). However, the authors claim that one of the main contributions of this work is "Interpretable MF theory", which seems unjustified. In this context, I have the following questions/suggestions, and I would be happy to raise my score if they were to be answered as well:
(a) Could the authors explicitly relate their derived mean-field equations to the MF formulations in other recent static theories (e.g., Seroussi et al., 2023; Fischer et al., 2024; Rubin et al., 2024)? Or at least some of them?
(b) If, as I suspect, these result in the same distribution (possibly up to negligible factors of $1/N^2$), then positioning "interpretable FL theory" as a primary contribution is an overstatement and should be clarified in the text. 
(c) Otherwise, the explicit difference between existing approaches should be highlighted.

2. The MF-ARD posterior for a single neuron (Eq. 14) has a very similar form to the plain MF posterior (Eq. 6), with the main difference being the introduction of $d$ new order parameters, the precisions $\rho_j$. Could the authors expand on in what sense does this modification constitute a more interpretable theory of feature learning as posited in the introduction?

3. Could the authors clarify the scaling of the noise parameter $\kappa$ with the network width $N$? In some mean-field treatments, $\kappa$ also acquires a width-dependent scaling, such as $\kappa \sim \mathcal{O}(N^{\gamma-1})$. Is such a scaling assumed here, and if not, how does its absence affect the results?

4. What are the limitations\assumptions on the MF approximation? This should be detailed explicitly in the main text.

5. The paper claims that the "homogeneity assumption" of the MF model prevents it from capturing the effect of specializing neurons. However, MF theory is seemingly capable of capturing distributions where only a small number of neurons specialize through multimodal Gaussian distributions. Essentially, taking a distribution where there is a $1/N$ probability of being in a certain nonzero mode is equivalent to saying that one neuron is drawn from a distribution that is some nonzero Gaussian, and the remaining neurons are drawn from the other modes of the distribution. How does this picture differ from the picture provided here? Contrasting this specialization with multi-modal distributions in the paper itself would be very helpful. 

6. The paper's “related works” section implies that many static MF theories are restricted to the proportional limit where $P \propto N$, where $P,N$ are sample size and network width respectively. This is not a general requirement, as shown in other works (e.g., Ringel et al. (2025) - not cited in this paper). There, other scalings of $P$ were considered, specifically taking $N\propto d$ (where $d$ is the input dimension), but $P$ scales as $d^{3/4}$. Is taking $N\propto d\mapsto\infty$ what the authors were refering to in the statement: “either in double asymptotic limits of data and width”? If so, this should be clarified in the main text.
soundness: 3
strengths: **Clarity**: The paper is exceptionally well-written and structured. The theoretical difficulty hierarchy of SGLD $\rightarrow$ MF $\rightarrow$ NNGP provides a clear and effective pedagogical framework for situating the paper's contributions and understanding the different levels of theoretical approximation. The core argument is easy to follow, and the visualizations are very helpful.

**Quality**: The empirical validation for the specific setting under consideration is of high quality. The quantitative agreement between the MF-ARD theory and SGLD experiments, as shown in Figures 1 and 5, is impressive and demonstrates that the proposed model accurately captures the learning dynamics in this regime. 

**Originality & Significance**: The use of the IFS mechanism as a simple extension to the common MF approaches discussed in the paper is a sound and valuable scientific contribution. The authors have successfully predicted an emergent mechanism of FL, which has a simple intuitive explanation.
summary: This paper presents a theoretical analysis of FL in two-layer non-linear networks trained with SGLD. The authors argue that standard MF approaches to this kind of problem, while correctly predicting the onset of FL as a phase transition, fail to capture the significant generalization improvements that follow. They attribute this failure to a FL mechanism they term Input Feature Selection (IFS). Relying on simple intuitive reasoning, the authors argue that SGLD provides a self-reinforcing dynamic where the network amplifies weights corresponding to relevant input dimensions, leading to neuronal sparsification. The authors show that such a process cannot be described by standard MF approach. To remedy this, they propose a minimal extension, MF-ARD (Automatic Relevance Determination), which incorporates a learnable, coordinate-wise variance into the weight prior to the MF derivation. This model is shown to quantitatively reproduce the learning curves of SGLD-trained ReLU networks on sparse-parity tasks.
weaknesses: Despite its strengths, the paper suffers from several significant weaknesses related to its framing, the generality of its claims, and the novelty of its contribution. My primary concerns are detailed below. Other than these critical points, I think that this is a good paper, and I would be happy to raise my rating if the concerns here are appropriately addressed in the revised text. 

**Critical reference omission**: The paper's analysis is critically hampered by its failure to engage with the concurrent work of Van Meegen & Sompolinsky (2024) (VM&S), which provides a general framework for understanding learned representations in the feature-learning regime. Crucially, VM&S. demonstrate that the weight distribution is determined by the neuronal nonlinearity. Specifically, they show that ReLU networks produce sparse distributions as in this work, while other nonlinearities like sigmoidal, produce qualitatively different distributions. 

**Overstated novelties**: 
(a) It seems that the omitted reference may have led the authors to misinterpret their central finding. The paper identifies sparsification (driven by IFS) as the key mechanism that standard MF theory fails to capture, whereas following VM&S, this seems to be only an outcome of the ReLU activation function used in this work. Considering the previous point, the scope of this paper may not be quite as wide as implied by the authors. 
(b) A central claim of the authors is that they developed an interpretable MF theory, where other works involve “... complex theoretical machinery [that] often obscures the core mechanisms driving FL”.  This is not an accurate representation of existing works, where the “machinery” is solving an equation in one (Li & Sompolinsky, (2021); Pacelli et al. (2023; Ringel et al. (2025)) or two (Rubin et al. (2024)) variables. Moreover, existing frameworks produce simple, intuitive, and interpretable descriptions of FL, and make accurate predictions in various regimes. For instance, the works of Li & Sompolinsky (2021) and Pacelli et al. (2023) offer a simple rescaling explanation, Fischer et al. (2024) find interpretable structures emerging in kernels, and Rubin et al. (2024) show that weights learn to prefer relevant directions. The authors even claim that "Our key contribution is to show that incorporating this IFS mechanism requires only a minimal,
principled modification to the MF theory". By this reasoning, the interpretability of their framework cannot be a significant improvement over MF approaches.

---

confidence: 4
contribution: 2
initial_rating: 4
presentation: 2
questions: 1.	Contextualization and Novelty: please properly cite provide work and provide detail comparison highlighting the novelly of the current work, and how your results go beyond the established phase diagrams and insights provided in the paper by Meegen-Sompolinsky.

2.	Introduction and MF Theory Definition: The Introduction needs complete rewriting to clearly state the problem and the specific theoretical gap. 

3. Line 107: Please remove the statement in Line 107 regarding the comparison between SGLD noise and SGD "batch noise," as this is misleading and conceptually inaccurate.

4.	Theorem Clarity (Theorem 4.1, Eq. (16)):
        (a)  Equation (16) uses both asymptotic notation ($\asymp$) and asymptotic bounds ($\Theta$). Please use only one type of notation for consistency.

        (b)  Please explicitly state in the main text whether the theorem result is specifically derived for the parity model or holds more generally.

       (c)  Please explicitly specify all assumptions required for the theorem to hold, including the assumptions on $\chi_S$ and any extra assumptions currently hidden in the Appendix.

5.	Generalization of $\phi$: The generalization results focus heavily on the overlap $\phi$ with the teacher function. All results both theoretical and numerical are stated for RELU. Can the authors provide more explicit results or discussion on other activation functions, and how do the phase transition and generalization error change? Is the same prior valid for all activation functions? How does it compare with existing MF theories. 

6.	Can you provide examples where the sparse priors would be inappropriate?
7.    Extension to depth: Can the ARD-MF idea apply layer-wise in deeper networks? Do you have any preliminary evidence?

Van Meegen, A., & Sompolinsky, H. (2025). Coding schemes in neural networks learning classification tasks. Nature Communications, 16(1), 3354.
soundness: 2
strengths: 1.	Analytical Tractability and Completeness: The paper provides a self-consistent Mean-Field framework for the Bayesian posterior of a two-layer network. This level of analytical derivation, together with the appropriate approximation, leads to closed-form solutions, which is a valuable result.
2.	Explicit Verification: The theoretical predictions are explicitly validated against SGLD simulations on the teacher-student setup, demonstrating a strong match between the Mean-Field approximation and the numerical results, particularly regarding the phase transition phenomenon.
summary: This paper develops a tractable mean-field (MF) theory to understand feature learning in two-layer neural networks trained with SGLD. The authors identify a two-stage feature learning process: (1) a phase transition where networks align with target functions, and (2) self-reinforcing input feature selection (IFS). This provides an extension over the recently developed MF theory. Experiments on k-sparse parity and single-index models demonstrate quantitative agreement of the theory with SGLD-trained networks.
weaknesses: 1.	Limited Novelty: The main results, including the phase transition and the sparsification of neurons, are not genuinely novel. Furthermore, the phenomena of feature learning and the splitting of network behavior to different activation function were thoroughly explored by Van-Meegen-Sompolinsky Nat. Comm 2025 in earlier, comprehensive studies, which the current work fails to cite and compare there results to. The current paper does not provide key new insights beyond these existing frameworks, especially regarding the challenging question of extending these ideas to deeper neural networks.

2.	Narrow Scope and Lack of Generalization: The analysis is restricted to a single, synthetic teacher-student task (Hermite polynomials) and demonstrates results limited to this toy setting. The paper fails to provide any insights or validation on real-world data, which limits the scope and relevance of the conclusions compared to related works.

3.	The Introduction is insufficient, failing to clearly articulate the existing literature and the specific gap the paper aims to fill, especially given the extensive prior work on mean-field feature learning.


Decision: Reject