Neuron Pursuit for Deep Neural Networks

• Neuron Pursuit is a greedy, stagewise algorithm for training deep homogeneous neural networks by iteratively adding neurons based on dominant KKT directions.
• It utilizes constrained maximization of a Neural Correlation Function to select neurons that improve performance and promote balanced sparsity in the network.
• The algorithm alternates neuron addition with gradient descent optimization, demonstrating descent properties and connections to methods like OMP and Frank–Wolfe.

Neuron Pursuit (NP) is a greedy, stagewise algorithm for training deep homogeneous neural networks. It alternates between expanding the model by adding neuron(s) whose incoming and outgoing weights are selected through constrained maximization of a Neural Correlation Function (NCF), and minimizing the training loss over the augmented network. The method is motivated by a gradient-flow analysis for homogeneous networks in the small-initialization regime, especially the emergence of sparsity at the origin, the structure of the first saddle point encountered after escape from the origin, and empirical saddle-to-saddle dynamics (Kumar et al., 15 Sep 2025).

1. Homogeneous-network setting

The framework is formulated for neural networks that are $L$-positively homogeneous in their parameter vector $\mathbf{w}\in\mathbb{R}^k$, meaning

$$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$

Examples explicitly identified in the formulation include feed-forward networks with activations of the form $\sigma(x)=\max(x,\alpha x)^p$, including ReLU ($p=1$), leaky ReLU ($p=1$), squared ReLU ($p=2$), and linear networks (Kumar et al., 15 Sep 2025).

For feed-forward $L$-layer networks without explicit biases,

$$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$

with $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$, $\mathbf{w}\in\mathbb{R}^k$0, and $\mathbf{w}\in\mathbb{R}^k$1. Biasing can be implemented by appending a constant input. Positive homogeneity yields a scaling law for gradients, $\mathbf{w}\in\mathbb{R}^k$2, which underlies the time-rescaling identities used throughout the analysis.

Training is posed with squared loss on data $\mathbf{w}\in\mathbb{R}^k$3. Writing $\mathbf{w}\in\mathbb{R}^k$4, $\mathbf{w}\in\mathbb{R}^k$5, and $\mathbf{w}\in\mathbb{R}^k$6, the objective is

$\mathbf{w}\in\mathbb{R}^k$7

and the continuous-time gradient-flow dynamics are

$\mathbf{w}\in\mathbb{R}^k$8

where $\mathbf{w}\in\mathbb{R}^k$9 is the Jacobian of $$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$0 with respect to $$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$1.

A central analytic object is the NCF,

$$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$2

together with the constrained NCF maximization

$$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$3

In NP, the residual vector and layer-specific homogeneous responses of currently inactive neurons determine which neuron should be added next.

2. Gradient flow from the origin to sparse saddles

The theoretical motivation begins with the behavior of gradient flow under small initialization. For $$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$4, gradient flow on the NCF either converges to the origin or diverges in norm while converging in direction to a non-negative KKT point of the constrained NCF problem. If initialization is $$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$5 with $$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$6 small, gradient flow remains near the origin for time $$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$7, and the weights align toward KKT directions when the initial direction lies in their stable set (Kumar et al., 15 Sep 2025).

At positive KKT points of the constrained NCF for feed-forward homogeneous networks, the layers are balanced. For each hidden neuron $$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$8 at layer $$\mathcal{H}(x;c\mathbf{w})=c^L\,\mathcal{H}(x;\mathbf{w}),\qquad c\ge 0.$$9,

$\sigma(x)=\max(x,\alpha x)^p$0

This proportionality implies that neurons with approximately zero incoming weights also have approximately zero outgoing weights. The consequence is emergent sparsity and the presence of inactive neurons early in training.

The analysis then extends beyond the origin. After the initial phase, gradient flow escapes the origin along a path determined by the dominant KKT direction and reaches a saddle point $\sigma(x)=\max(x,\alpha x)^p$1 of the training loss. The first saddle encountered has a specific sparsity pattern: the incoming and outgoing weights of a subset of hidden neurons are exactly zero. Along the limiting gradient-flow path beyond the origin, the corresponding zero rows and columns remain zero, so sparsity persists through the escape-to-saddle phase (Kumar et al., 15 Sep 2025).

This first saddle is the prototype for NP. The paper studies saddle points of the form $\sigma(x)=\max(x,\alpha x)^p$2, where $\sigma(x)=\max(x,\alpha x)^p$3 denotes active weights and $\sigma(x)=\max(x,\alpha x)^p$4 denotes incoming and outgoing weights of currently inactive neurons. These saddles exhibit a sparse block structure: active blocks carry the current representation, while inactive blocks define a low-norm subspace in which future neuron activation will occur.

3. Directional convergence near saddles

Near a saddle of the form $\sigma(x)=\max(x,\alpha x)^p$5, the analysis assumes a decomposition of the network output into the contribution of the active subnetwork and a leading homogeneous response of the inactive block. For fixed $\sigma(x)=\max(x,\alpha x)^p$6 and small $\sigma(x)=\max(x,\alpha x)^p$7,

$\sigma(x)=\max(x,\alpha x)^p$8

together with corresponding gradient expansions in $\sigma(x)=\max(x,\alpha x)^p$9 and $p=1$0, where $p=1$1 is $p=1$2-homogeneous in $p=1$3, $p=1$4, and has locally Lipschitz gradients. For feed-forward networks with $p=1$5, this decomposition is proved explicitly; in the leaky ReLU case with $p=1$6, smoothness for the Taylor expansions requires $p=1$7, whereas for $p=1$8 the decomposition is polynomial (Kumar et al., 15 Sep 2025).

A second ingredient is a Lojasiewicz-type condition for the restricted loss

$p=1$9

namely

$p=1$0

for $p=1$1 within a neighborhood of $p=1$2. This condition is used to keep the active block near the saddle for long times and to bound the path length.

The time-rescaling identity for $p=1$3-homogeneous gradient flows is

$p=1$4

Thus, small initialization slows the dynamics by a factor $p=1$5 without changing the limiting direction. Inactive weights therefore evolve slowly enough for their directional behavior to be characterized while the active block remains close to the saddle.

Let $p=1$6 be the residual at the saddle, and define the NCF associated with the leading inactive-neuron response by

$p=1$7

Under the decomposition and Lojasiewicz assumptions, gradient flow initialized at $p=1$8 satisfies, for sufficiently small $p=1$9,

$p=2$0

for all $p=2$1. If the initial inactive direction lies in the stable set of a KKT point of the constrained NCF, then $p=2$2 remains small in norm but converges in direction to that KKT point. If not, its norm remains only $p=2$3. This is the direct theoretical rationale for adding neurons aligned with dominant KKT directions rather than relying on arbitrary random growth (Kumar et al., 15 Sep 2025).

For feed-forward networks, the leading homogeneous response of the inactive block is explicitly decomposed into neuron-wise terms. When $p=2$4, each term is $p=2$5-homogeneous in the inactive weights, and the NCF splits into a sum of neuron-wise homogeneous components that can be optimized in parallel with independent initializations. Only components associated with dominant KKT points survive in the joint limiting direction. A plausible implication is that overparameterization is advantageous because multiple candidate inactive-neuron directions can be explored simultaneously.

The same proportionality law reappears at KKT points near a saddle: $p=2$6 This balancedness motivates the balanced small initialization used when NP adds a new neuron.

4. Neuron Pursuit algorithm

NP operationalizes the saddle analysis as an explicit training procedure. The objective remains

$p=2$7

but instead of fixing the architecture width in advance, NP grows the network iteratively (Kumar et al., 15 Sep 2025).

The initialization phase begins with one neuron per layer. It solves

$p=2$8

by projected gradient ascent with multiple random initializations. The initial network is set to $p=2$9, and gradient descent is run on $L$0 until a stationary point $L$1 is reached.

At iteration $L$2, the residual is

$L$3

For each hidden layer $L$4, NP solves the layer-specific constrained NCF

$L$5

where

$L$6

and

$L$7

Here $L$8, $L$9, and $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$0 is the downstream subnetwork determined by the current weights.

The dominant KKT pair is then selected across layers: $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$1 A single neuron is added in layer $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$2 with small KKT-aligned incoming and outgoing weights,

$$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$3

while the other layers are left unchanged. Gradient descent is then run again from the augmented network to obtain the next stationary point. The procedure repeats for $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$4 iterations or until a target level of training loss, training error, or test error is reached.

The algorithm is therefore stagewise in a precise sense: each stage identifies the currently most correlated homogeneous response of an inactive-neuron template to the residual, inserts the corresponding neuron, and reoptimizes the enlarged network. This emulates the empirically observed saddle-to-saddle training dynamics rather than attempting to learn all neurons simultaneously.

5. Descent properties, related methods, and empirical behavior

The analysis provides a per-addition descent statement. If the dominant KKT score is positive, then for sufficiently small $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$5, adding the selected neuron strictly decreases the loss before any subsequent gradient-descent refinement: $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$6 For the first stage this compares $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$7 with the loss after adding the initial small network; at later stages it compares the loss at the current stationary point with the loss after appending the KKT-aligned neuron. If gradient descent thereafter decreases the loss, then NP has a monotonic descent property (Kumar et al., 15 Sep 2025).

The method is explicitly related to several classical greedy or stagewise procedures. In diagonal two-layer linear networks, NP coincides with Orthogonal Matching Pursuit (OMP) under mild assumptions, because the constrained NCF reduces to selecting the column of $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$8 most correlated with the residual, followed by least-squares on the selected support. NP also has an analogy with functional gradient boosting, since it adds at each step the basis function most aligned with the residual, and with Frank–Wolfe, since neuron addition is obtained by constrained maximization over a unit-norm set of candidate atoms. Relative to growth and pruning methods, NP constructs sparsity from scratch through staged neuron addition guided by directional-convergence theory.

Implementation uses full-batch gradient descent in the minimization phase and projected gradient ascent for the constrained NCF problems with $$\mathcal{H}(x;W_1,\dots,W_L)=W_L\,\sigma\big(W_{L-1}\sigma(\cdots \sigma(W_1x)\cdots)\big),$$9 random restarts. The scalar $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$0 controls the initial norm of added weights and is reduced multiplicatively if the loss does not decrease; the example given is a factor $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$1 for a few trials. For ReLU pointer value retrieval (PVR), weights are rebalanced after each iteration so that incoming and outgoing norms per neuron are equal, described as Saul’s balancing for ReLU. The per-iteration cost comprises $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$2 restarts of projected gradient ascent for each hidden layer, with cost similar to a few forward and backward passes, plus full-batch gradient descent to convergence for the augmented network.

The experimental study covers synthetic sparse functions and algorithmic tasks. On the hypersphere, with $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$3 uniform on $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$4 and $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$5, the targets $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$6 and $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$7 are 3-layer ReLU networks with $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$8 first-layer neurons and $W_l\in\mathbb{R}^{k_l\times k_{l-1}}$9 second-layer neurons, respectively. On the hypercube, with $\mathbf{w}\in\mathbb{R}^k$00 uniform on $\mathbf{w}\in\mathbb{R}^k$01, the targets are $\mathbf{w}\in\mathbb{R}^k$02 and $\mathbf{w}\in\mathbb{R}^k$03, learned using 3-layer leaky-ReLU networks with $\mathbf{w}\in\mathbb{R}^k$04. On Gaussian inputs $\mathbf{w}\in\mathbb{R}^k$05, the targets are $\mathbf{w}\in\mathbb{R}^k$06 and

$\mathbf{w}\in\mathbb{R}^k$07

learned with 4-layer leaky-ReLU networks. The metrics are relative training and test errors $\mathbf{w}\in\mathbb{R}^k$08. The baseline is gradient descent with small initialization, $\mathbf{w}\in\mathbb{R}^k$09 neurons per layer, and up to $\mathbf{w}\in\mathbb{R}^k$10 iterations; NP is run up to $\mathbf{w}\in\mathbb{R}^k$11 iterations.

The reported outcomes are mixed but structured. On hypersphere $\mathbf{w}\in\mathbb{R}^k$12, NP achieves low test error with fewer samples than gradient descent; on $\mathbf{w}\in\mathbb{R}^k$13, both methods require similar sample sizes. On hypercube data, NP matches gradient descent on $\mathbf{w}\in\mathbb{R}^k$14 and outperforms it on $\mathbf{w}\in\mathbb{R}^k$15, described as consistent with “leap complexity” phenomena. On Gaussian $\mathbf{w}\in\mathbb{R}^k$16, NP slightly outperforms gradient descent; on $\mathbf{w}\in\mathbb{R}^k$17, NP requires more samples than gradient descent. The number of NP iterations often peaks at intermediate sample sizes, where the test error transitions from high to low, and then stabilizes.

Two additional tasks highlight algorithmic structure. For modular addition, a two-layer network with squared activation $\mathbf{w}\in\mathbb{R}^k$18 is trained on $\mathbf{w}\in\mathbb{R}^k$19 classes using one-hot inputs with an extra constant $\mathbf{w}\in\mathbb{R}^k$20 for a bias in the first layer. NP is trained until zero training error, yields small test classification error, and learns weights with sinusoidal structure in the DFT, described as consistent with mechanistic interpretability results on grokking modular arithmetic. For PVR, a three-layer ReLU network takes inputs containing a pointer encoded in symmetric binary form and a bias coordinate, and predicts $\mathbf{w}\in\mathbb{R}^k$21. NP converges in approximately $\mathbf{w}\in\mathbb{R}^k$22 iterations with low test error; the learned first-layer weights tied to data coordinates are sparse and show localized dominant entries in each row (Kumar et al., 15 Sep 2025).

Training curves exhibit plateaus followed by sharp drops, matching the saddle-to-saddle picture. Small-weight blocks converge in direction while remaining small in norm, incoming and outgoing norms remain proportional, and sparsity patterns persist across successive saddles. For some ReLU settings, the leading constrained NCF can be zero, including the examples of deep linear networks or specific data configurations; in such cases, higher-order homogeneous terms govern directional convergence, and NP can still progress by maximizing the next nonzero homogeneous component.

6. Limitations, sensitivity, and acronym ambiguity

Several limitations are explicit. The Lojasiewicz inequality near $\mathbf{w}\in\mathbb{R}^k$23 with $\mathbf{w}\in\mathbb{R}^k$24 is used to keep the active block near the saddle for times $\mathbf{w}\in\mathbb{R}^k$25, but broader validity is left open. The dominant KKT point depends on the initial direction of the small inactive weights near a saddle; NP approximates this dependence using random restarts, so discrepancies from true gradient-flow behavior can arise when the relevant KKT direction is not found. Layer rescaling can alter selection decisions, which is why balancing is emphasized, especially for ReLU. NP also introduces overhead per iteration because each stage requires both constrained NCF optimization and gradient descent to stationarity, although the number of added neurons remains small and the model stays relatively sparse (Kumar et al., 15 Sep 2025).

A common source of confusion is the acronym “NP.” In the literature it also denotes node perturbation, a forward-gradient-style learning rule that estimates credit assignment by injecting perturbations into unit activities, performing two forward passes, and correlating the induced loss change with the perturbations (Dalm et al., 2023). That usage is unrelated to Neuron Pursuit. Neuron Pursuit is a greedy network-growth method for deep homogeneous networks, driven by constrained NCF maximization and saddle-aligned neuron addition; node perturbation is a perturbation-based alternative to backpropagation that uses no backpropagated derivatives and relies on a global scalar performance signal. The shared acronym therefore names two technically distinct research programs.

Taken together, NP presents a theoretically guided view of training in homogeneous networks: early gradient-flow dynamics generate sparse and balanced structures, successive saddles identify which inactive-neuron directions are most correlated with the residual, and explicit addition of neurons along dominant KKT directions converts that analysis into a stagewise optimization procedure. This suggests a training paradigm in which architecture growth, implicit sparsity, and optimization dynamics are treated as a single coupled process rather than as separate design choices.