The Implicit Bias of Minima Stability in Multivariate Shallow ReLU Networks
Authors: Mor Shpigel Nacson
Rotem Mulayoff
Greg Ongie
Tomer Michaeli
Daniel Soudry
Venue: ICLR 2023 poster (2023)
Keywords: Implicit bias
implicit regularization
stability
Hessian
dynamical systems
depth separation
approximation

Abstract: We study the type of solutions to which stochastic gradient descent converges when used to train a single hidden-layer multivariate ReLU network with the quadratic loss. Our results are based on a dynamical stability analysis. In the univariate case, it was shown that linearly stable minima correspond to network functions (predictors), whose second derivative has a bounded weighted $L^1$ norm. Notably, the bound gets smaller as the step size increases, implying that training with a large step size leads to `smoother' predictors. Here we generalize this result to the multivariate case, showing that a similar result applies to the Laplacian of the predictor. We demonstrate the tightness of our bound on the MNIST dataset, and show that it accurately captures the behavior of the solutions as a function of the step size. Additionally, we prove a depth separation result on the approximation power of ReLU networks corresponding to stable minima of the loss. Specifically, although shallow ReLU networks are universal approximators, we prove that stable shallow networks are not. Namely, there is a function that cannot be well-approximated by stable single hidden-layer ReLU networks trained with a non-vanishing step size. This is while the same function can be realized as a stable two hidden-layer ReLU network. Finally, we prove that if a function is sufficiently smooth (in a Sobolev sense) then it can be approximated arbitrarily well using shallow ReLU networks that correspond to stable solutions of gradient descent.

Published as a conference paper at ICLR 2023
THE IMPLICIT BIAS OF MINIMA STABILITY IN MULTIVARIATE SHALLOW RELU NETWORKS
Mor Shpigel Nacson' & Rotem Mulayoff Electrical & Computer Engineering, Technion
Greg Ongie Mathematical and Statistical Sciences Marquette University
Tomer Michaeli & Daniel Soudry Electrical & Computer Engineering, Technion
ABSTRACT
We study the type of solutions to which stochastic gradient descent converges when used to train single hidden-laye multivariate ReLU network with the quadratic loss. Our results are based on a dynamical stability analysis. In the univariate case, was shown that linearly stable minima correspond to network functions (predictors), whose second derivative has a bounded weighted 7,7 norm. Notably, the bound gets smaller as the step size increases, implying that training with a large step size leads to 'smoother' predictors. Here we generalize this result the multivariate case, showing that a similar result applies to the Laplacian of the predictor. We demonstrate the tightness of our bound the MNIST dataset, and show that it accurately captures the behavior of the solutions as a function of the step size. Additionally, we prove depth separation result the approximation power of ReLU networks corresponding to stable minima of the loss. Specifically, although shallow ReLU networks are universal approximators, we prove that stable shallow networks are not. Namely, there is a function that cannot be well- approximated by stable single hidden-layer ReL.U networks trained with non- vanishing step size. This is while the same function can be realized as stable two hidden-layer ReLUU network. Finally, we prove that a function is sufficiently smooth (in a Sobolev sense) then it can be approximated arbitrarily well using shallow ReLU networks that correspond to stable solutions of gradient descent.
INTRODUCTION
Neural networks (NNs) have been demonstrating phenomenal performance in a wide array of fields, from computer vision and speech processing medical sciences. Modern networks are typically taken to be highly overparameterized. In such setting, the training loss usually has multiple global minima, which correspond models that perfectly fit the training data. Some of those models clearly sub-optimal in terms of generalization. Yet, the training process seems to consistently avoid those bad global minima, and somehow steer the model towards global minima that generalize well A long line of works attributed this behavior to "implicit biases" of the training algorithms, e.g., (Zhang al., 2017; Gunasekar et al., 2017; Soudry et al., 2018: Arora al., 2019).
Recently, it has been recognized that dominant factor affecting the implicit bias of gradient descent (GD) and stochastic gradient descent (SGD), is associated with dynamical stability. Roughly speaking, the dynamical stability of minimum point refers to the ability of the optimizer to stably converge to that point. Particular research efforts have been devoted to understanding linear stability, namely the dynamical stability of the optimizer's 's linearized dynamics around the minimum (Wu et al., 2018; Nar & Sastry, 2018; Mulayoff et al., 2021; Ma & Ying, 2021). For GD and SGD, it well known that a minimum is linearly stable if the loss terrain is sufficiently flat w.r.t. the step size 1)-
Concretely, a necessary condition for a minimum to be linearly stable for GD and SGD is that the top eigenvalue of the Hessian at that minimum point be smaller than 2/n (see Sec. 2). Although this
Indicates equal contribution. Correspondence: (mor. shor.shpigo shppielll..... gel, rotem. mu ayo @ggm i com.

---

Published as a conference paper at ICLR 2023
20
0-
-20
-20
-2
-2
272
(a) Small step size (1) = 0.001) (b) Medium step size (1) = 0.01) (c) Large step size (1) = 0.1)
Figure I: Larger step size leads to smoother prediction function. We train a single hidden-layer ReLU network on a regression task with two-dimensional data, depicted by red points. The different panels show the predictor function f obtained when training with different step sizes.
condition only characterizes the linearized dynamics, it has been empirically shown to hold in real- world neurll training (Cohen al.. 2020; Gilmer et al., 2022). The linear stability condition turns out to have a strong effect on the nature of the network that is obtained upon convergence, both in terms of the end-to-end predictor function (Mulayoff al., 2021), ), and in terms of the way this function is implemented by the network (Mulayoff & Michaeli, 2020).
Mulayoff' et al. (2021) studied how linear stability affects single hidden-layer univariate ReLU network, when trained with the quadratic loss. They showed that in this setting, stable solutions of SGD with step size ny correspond to functions f satisfying
f (1)
1 1
7 2
where denotes the network input-output function, and g is weight function that depends only on the training data. This result implies that for univariate shallow ReLU networks, SGD is biased towards 'smooth' solutions'. Moreover, the larger the step size 17, the smoother the solution becomes.
In this paper, we study the stable solutions of single hidden-layer ReLU networks with multidimen- sional inputs, trained using SGD and the quadratic loss. Particularly, in Sec. we generalize the result of Mulayoff et al. (2021) to the multivariate setting. As it turns out, the natural extension of (1) involves the Radon transform of the Laplacian of the predictor function, Af (see Thm. 1). However. we show this result can also be interpreted in primal space as
(A(()))((2)dd <
(2)
JRd
2
where p is some weighting function. Thus, stable solutions of SGD in the multivariate case also correspond to smooth predictors (i.e., functions whose Laplacian has small weighted L' norm). The larger the step size, the smoother the function becomes. Figure illustrates this phenomenon.
Additionally, study the approximation power of single hidden-layer ReLU networks corresponding to stable minima. It well known that shallow ReI.U networks can approximate any continuous function over a compact set (Pinkus, 1999). However this does not imply that SGD can stably converge to such approximations. If there exist functions whose approximations are all unstable, then this property may of limited practical interest. In Sec. we prove that every convergent sequence of stable networks has a limit function that also satisfies the stability condition (Thm. 1). Building on this, we prove a depth separation result. Specifically, we show that there exists function that does not satisfy the stability condition for any positive step size. Namely, it cannot be stably approximated by a single hidden-layer ReLU network trained with non-vanishing step size. Yet, the same function can be realized as two hidden- layer ReLU network corresponding to stable minimum. Moreover, in Sec. 5 we show that a function sufficiently smooth (Sobolev) then it can be approximated arbitrarily well using single hidden-layer ReL.U networks that correspond to stable solutions of GD.
Finally, in Sec. 3.3 and 6 we demonstrate our results. Particularly, we illustrate how our stable minima characterization (Thm. 1) can be used to predict certain properties of the solution. For example, for
In slight abuse of terms, in this paper we say a function is 'smooth' if some weighted T,1 norm of its second derivative is bounded.
2

---

Published as a conference paper at ICLR 2023
certain isotropic data (e.g., Gaussian), show that a large step size tends to increase the biases of all neurons. We also demonstrate on the MNIST dataset the tightness of our stability bound, and that predicts well the dependence of the stability and generalization performance on the step size.
2 BACKGROUND: MINIMA STABILITY OF SGD
In this section we give brief survey on minima stability. Let us consider the problem of minimizing empirical loss using SGD. We are interested in the typical regime of overparameterized models. In this setting, there exist multiple global minimizers of the loss. Yet, SGD cannot stably converge to any minimum. The stability of minimum is associated with the dynamics of SGD in its vicinity. Specifically, a minimum said to be stable if once SGD arrives near the minimum, it stays in its vicinity. If SGD repels from the minimum, then we say that it is unstable.
Formally, let l R4 - I be differentiable almost everywhere for all j E [n]. Here we consider a loss function L and its stochastic analogue,
1 EE (0)
-(0)= (0) and (.( - -
JEB,
where B, is batch of size B sampled iteration 1. We assume that the batches [B.} are drawn uniformly from the dataset, independently across iterations. SGD's update rule is given by
0 = 0, - v()((,,))
where 1j is the step size. Analyzing the full dynamics of this system is intractable most cases. Therefore, several works studied the behavior of this system near minima using linearized dynamics (Wu et al., 2018; Ma & Ying, 2021; Nar & Sastry, 2018; Mulayoff et al., 2021), which is common practice for characterizing the stability of nonlinear systems.
Definition (Linear stability). Let 0 be a twice differentiable minimum of Consider the linearized stochastic dynamical system
O2+1 - 0. On-n -n V2 ("",(0")) 00)).
(5)
Then 0* is e linearly stable for any 0y in the s-ball B, (f"), we have lim B[ - 0- <
t+o00
Namely, minimum is E linearly stable if once θ, enters an E e-ball around the minimum, it ends up at distance no greater than E from in expectation. Under mild conditions, any stable minimum of the nonlinear system is also linearly stable (Vidyasagar, 2002, p. 268). We have the following condition. Lemma 1 (Necessary condition for linear stability (Mulayoff et al., 2021, Lemma 1)), Consider SGD with step size g,, where batches are drawn uniformly from the training set, independently across iterations. If 0" is an linearly stable minimum of L, then
Amax (V2 17
This condition states that stable minima of SGD are flat w.r.f. the step size. Although this result was proved for the linearized dynamics, was observed to hold also in practice, where the full nonlinear dynamics apply. Particularly, much empirical evidence real-world neural-network training (Cohen al., 2020; Gilmer et al., 2022) points out that GD and SGD converge only to linearly stable minima, i.e., minima satisfying (6). More dynamical stability and its interaction with common practices (e.g., learning rate decay, absence of & and B in Lemma result, etc.) in App. A.
LARGE STEP SIZE BIASES TO SMOOTH FUNCTIONS
Consider the set of multivariate functions over 1d that can be implemented by single hidden-layer ReLU network with & neurons,
Fr S: J f(x)
(2) wi
x NTTw!!
(1) -1(2)) + ((2)
i=1
where (-) denotes the ReLU activation function. Each f E F& is piecewise linear function with at most & knots?, Given some training set {æ y3 -=> are interested in functions that globally
A 'knot' a boundary between two pieces (i.e., intersection between hyperplanes) See Fig. for illustration.
3

---

Published as a conference paper at ICLR 2023
minimize the quadratic loss"
1 2n j=1
C()
Σ
(f(x)) - -1j)
(8)
Definition (Solution). A function f E F is 'solution' iflution'' if(()) = i.e,, [(x,) - M/ Vj E [m].
We focus on the overparameterized regime (kd > n) in which there exist multiple solutions. We want to study the properties of solutions which correspond to stable minima of SGD. However, key challenge is that any solution Fs typically has infinitely many different parameterizations, In other words, there are various parameter vectors
T ]"
0 a W (I)T
(1)T W
b (U
u (2)
by
k
that can implement the same function f. Different parameterizations correspond to different minima, which may have different Hessian eigenvalues. Therefore, for a given step size 7), parameteriza- tions of f may be stable while others may not. Thus, to determine whether SGD can stably converge solution f. we need to check whether there exists some stable minimum 0, which corresponds to a parametrization of f. We therefore use the following definition.
Definition (Stable solution). A solution f E Fk said to be stable for step size n if there exists minimum 0 of the loss that corresponds to f. where 0" is linearly stable for SGD with step size r,.
The next theorem characterizes stable solutions using the Radon transform R (see App. C) and the Laplace operator A. Particularly, we use the inverse of the dual Radon transform, (R*)-1 and interpret A in the weak sense, i.e., sum of weighted Dirac delta functions (see App. D).
Theorem 1 (Properties of stable solutions). Let f be a linearly stable solution for SGD with step size 17. Assume that the knots of f do not coincide with any training point. Then
1 1_1 1
11 S
(10)
17 2
where I|-1RR.g is the stability norm, defined as
| - [ - [[(R))+AA)] (0,b)l((e b)ds())d,,
(11)
and g(v, b) e min (ğ(v, b), 0( -b) is non-negative weighting function, with ğ given by
g(v,b) b x t.s b) x [
+1.
(12)
Here X is random vector drawn from the dataset's 's distribution (i.e., sampled uniformly from (œ, )
This theorem, whose proof is provided in App. E, shows that the step size constrains the stability norm of the solution. Notably, the constraint becomes stricter the step size increases. Before interpreting this result, let us note that although it depends only on the step size, other hyper-parament parameerss (e.g., batch size, initialization) may potentially improve the bound. Yet, as we discuss App. A, the effect of other understood in primal space and in Radon space. In the following, we discuss both interpretations and give examples.
seems secondary in practical settings. The implications of Thm. can be
3.1 PRIMAL SPACE INTERPRETATION
Theorem stated in Radon space, which may be difficult interpret. However, in some cases can also be interpreted in primal space, by deriving an alternative form for the stability norm || R Specifically, in App. G we show that if is piecewise continuous and L initeggra,lle" then for all Fk and p = R R--g we have
1 = ( |AS()))o())d))
(13)
'We focus MSE loss for simplicity, but the results be extended to other loss functions, see App. J.
"Which is true, for example, when the training set is finite.
"This integral should he interpreted in the distrihutional sense (see App. for details).
4

---

Published as a conference paper at ICLR 2023
0.5
x10-1
2
34
0.4
2-
0.3
0
1-
0.
21
-2
0.1
-20
2
-2
2
-2
-2
22
211
(a) Surface of f
(b) Laplacian magnitude [Af]
(c) Weight function
Figure 2: Illustration of the stability norm. Panel (a) depicts an interpolating function f. Panel (b) displays the absolute value of the Laplacian of f. i.e., |a). Here the color codes the amplitude of the delta functions. Panel (c) presents the weight function p. The stability norm is the weighted sum of line integrals of p, according to |A
In this presentation of the stability norm, P is not necessarily non-negative. Nevertheless, all its hyper-plane integrals are non- negative,e,, since Rp = g 2 0. Thus, the stability norm can interpreted as a non-negative linear combination of hyper-plane integrals of p along the knots of f. This is visualized in Fig. 2. Hence, Thm. combined with (13) implies that the larger the step size v,, smoother the solution becomes.
3.2 RADON SPACE INTERPRETATION
Another interesting interpretation of Thm. 1 can be derived in Radon space. First, let us examine how the weight function g(v. behaves as function of b. For every fixed v. the function g(v,)) has finite support, [min; [ v), max; [2:] v)]. Moreover, g(v, -) typically has most of its mass concentrated around the center of the distribution of the projected data points {ar] v). and it decays towards the endpoints (see e.g., Fig. 3).
Next, let us interpret how the term (R*)-1A) Af behaves. For single hidden-layer ReLU network, (R is sum of Dirac deltas. Specifically, shown in (Ongie et al., 2020), if f is a function of the form f(x) Ein a..(U, x b.) + with ||e:|22 1 for all i [k]. then (see App. F.3)
R (
a (v.))))))
(14)
where Af is the (distributional) Laplacian of f, and ((m,b) denotes Dirac delta centered at (v, gd-1 IR. We can thus define parameter space representation for the stability norm as (see App. F.3)
Sg
aa y (v., bz).
(15)
Generally, this parametric representation of the stability norm satisfies || S Sg. where equality happens whenever the ReLU knots of the representation do not coincide (i.e., there is one Dirac function for each ReLU unit). Yet, this parametric view of the stability norm also obeys (see App. F.3)
So
(16)
7 2
Hence, larger step sizes 17 push Sg be smaller, and from (15) we see that Ja,| will tend be small. Also, since gvv.-) typically decays towards the boundary of its support, this pushes the neurons' biases, b,, away from the center of the distribution. The resulting effect that the predictor function becomes flatter, especially near the center of the distribution. This is illustrated in Figs. and 4(b).
3.3 EXAMPLES
Earlier we introduced two interpretations for the stability norm I1-1RR.9" one in primal space, which uses the weight function p. and one in Radon space, which uses the weight function g. In this section, we compute and p for two toy examples, for which (13) holds.
5

---

Published as a conference paper at ICLR 2023
x10
x10-1
0.3
Offset b
0.2
101
Offset
5.
6.1
0.1
-2
Angle (deg.)
Angle (dog.)
(a) Two data points: g
(b) Two data points:
(c) Gaussian data: g
(d) Gaussian data: P
Figure 3: Visualization of of Thm. and p of (13) for two toy examples. (a), (b) Two data points
21 (1,0) and 22 (-1, 0). (c)(d) Two dimensional Gaussian data, i.e. X ~ N (0, f).
Example 1: Two data points in R. Assume the dataset contains two points: æ1 (1, 0) and 222 = (-1,0). In this case we can analytically calculate g and (see App. H.1). Figure depicts these functions. Here, has singularities 1 22 and at the origin. Yet, despite these singularities, all line integrals of are finite and thus the expression fgu (x) p() p(z)dx well-defined. Moreover, while p takes negative values, all its line integrals take positive values.
Example Isotropic distribution. Suppose the data isotropically distributed, i.e., P(X 6) does not depend on the direction of v. Thus g is independent of v, which implies that R a radial function. In App. H.2 we give an analytic expression for g for any isotropic distribution. In the special case of X ~ (0, I). ggv) -) decays monotonically with b, and thus, discussed in Sec. 3.2, large step sizes will tend to increase the biases all neurons. Additionally, we show in App. H.2 that for 2D data, p is positive and strietly decreasing in ||=]|] and it satisfies the asymptoties
p(x) = ([ogg(|==))) as |x|| - O, and p(x) 0[ [|| as ||>0|| - co. Figure 3 visualizes and p for two dimensional Gaussian data.
4 STABILITY LEADS TO DEPTH SEPARATION
Single hidden- -ayee neural networks are universal approximators. i.e., they can approximate arbitrarily well any continuous function over compact sets (Pinkus, 1999). However, some of these approxi- mations may correspond to unstable minima that are virtually unreachable by training via SGD. To understand what is the effective approximation power of neural networks, we need to identify the class of functions that have stable approximations. We have the following (see proof in App. I.I).
Proposition 1. Let X be the interior of the convex hull of the training points, and f X - R. be any function. Suppose there exists sequence of single hidden-layer ReLU networks {fk] with bounded
stability norm that converges to f in L.' over X. Then [|] . is finite, and lim [ R.g = Iff) R.
k=+00
Let } be convergent sequence of stable solutions with a growing number of knots, i.e., Vfx E FA |lfellR., S 1/m 1/2 (see Thm. 1). Then, by the proposition above we have that the limit function f also satisfies this inequality. Therefore, the effective class of functions that can be approximated arbitrarily well by single hidden-layer ReLU networks includes only continuous functions / that satisfy the stability condition |I/l| S 1/1 - 1/2. As the step size decreases, more functions satisfy this condition, suggesting that more functions can be stably approximated by single hidden-layer ReLU networks. Surprisingly. there exists at least one continuous function that has ||p|| R. oc and therefore does not satisfy the stability condition for any positive step size (see proof App. I.2). Therefore from Prop. and Thm. this function cannot be approximated arbitrarily well by single layde ReLU networks trained with non-vanishing step size.
Proposition 2. Assume the input dimension d 2 2, and p(æ) ([1-||x0] Suppose the support of is contained in the interior of the convex hull of the training points. Then || ||
Intriguingly, this function does have an implementation as a finite-width two hidden-layer network, p(x) = =o(1- (o(;) +
(-)x)))) which is stable solution for a fixed step size. Indeed, in App. 1.3 we demonstrate that for an appropriate choice of 7,, GD is able to converge to this implementation. Thus, have depth separation result: the function p cannot be approximated by stable minima of one hidden-layer networks trained with a non-vanishing step size, yet with two hidden-layers, GD can converge to this function with fixed step size.

---

Published as a conference paper at ICLR 2023
5 SHALLOW NETWORK APPROXIMATIONS OF SMOOTH FUNCTIONS
In Sec. 4 we showed that stable single-hidden layer ReLU networks are not universal approximators. In this section we give an approximation guarantee under smoothness assumptions. That is, we show that a function is sufficiently smooth, then it can be approximated arbitrarily well using single hidden-layer networks that correspond to stable solutions of GD.
Let w W denote the weighted Sobolev space of all functions whose weak partial derivatives up to order d+1 1 are bounded in weighted L' nnorm ||- I, with weight function w(x) R² [1+1b][(x).
Let (l denote the corresponding Sobolev norm
Wlllww (ad) 11 Σ 1a1=
||00 I
(17)
where 9 is only. Our results use the "R-norm" || | introduced by Ongie et al. (2020) (see Sec. 7 for details), and the stability norm I| ||R. with a different weight function 9 defined below (see proof in App. L.).
For technical convenience, we restrict ourselves odd input dimensions d
Proposition 3. Assume d is odd and let f Wd+1, (R.). Then, there exists sequence of single hidden-layer ReLU network functions {fc} such that ft E F& converges to f in L' over any compact subset K R.a i.e. lim JK (fx((=) f(x) dd O, and satisfies the bounds k4000
[| + 1.ffllRR S ca.illf 11 for all k, where
b(v,b) - (x P(0,6)=P(x'>>b) /. (xTo-))) (x'0--))'''>0)11+ +
[|X]] x x x'
(18)
and ca.g is a constant depending on d and but independent of f. Here X is drawn uniformly at random from the dataset.
This proposition shows that for any f E Wd+1.1 (R2) there exists a sequence of single hidden-layer ReLU network approximations { for which [[] and {Il (l.fll/...)))) are bounded. To prove that these functions can have stable parameterizations for GD, we need show that if both the stability norm and R-norm are bounded (a function space property). then there exists a corresponding minimum with bounded sharpness in parameter space. To this end, we derive an upper bound on the minimal sharpness of a solution f among its different parameterizations in terms of the stability norm and the R-norm (see proof in App. K).
Lemma 2. Let f E Fx be a solution for which the knots do not coincide with any training point. such that
Then there exists implementation 0" corresponding
Amaxx ((2)())) s 1+ 211flle + (11f111e + inf [Vf(()I)) V Vrmaax (xx) 1 EE
X
wERd
(19)
Here X is drawn uniformly at random from the dataset, and Σx is the covariance matrix of X. Combining Prop. 3 and Lemma 2 get that any ʃ Wdd111 can be approximated arbitrarily well by sequence of stable solutions for GD with a fixed step size
Theorem 2. Suppose the input dimension d is odd, and let f Wd+1.1 (RR.). Then, there exist rj > 0 and a sequence of single hidden-layer ReLU network functions { such that fs. E F& converges to I in 1 over any compact subset K C IRd and every Sv. is stable for GD with step size 1).
fetworrfunctions
This theorem state that any sufficiently smooth function can be stably approximated in the limit of infinitely many neurons. We can also use Lemma to guarantee the stability of solutions in the finite case. Since Xmax < 2/1j sufficient condition for stability in GD, we have the following.
Theorem Let f E Fs be solution for which the knots do not coincide with any training point. If
Il + [|5
inf IVf(=)I)) VAmAn (ixxx<Ex)) (Ex) (Ex)
E
X 2
<
1
(20)
xERA
2
then f is a stable solution for GD with step size 17.
Theorem complements Thm. 1, as gives a sufficient condition for stability in function space.
Note that for GD, 1) 2/Amax is necessary and sufficient condition for linear stability.
7

---

Published as a conference paper at ICLR 2023
10
2/1 Aese 21lfilss
6-
Sharpness
[[b]]
10a
2-
109
0
10
10
10
10
10
10
10-2
10
7
7
(a) Sharpness versus step size
(b) Bias versus step size
Figure 4: Validating the bounds on synthetic data. We trained a two-layer ReLU network on regression task with synthetic data using GD (see Sec. 6). Panel (a) depicts the sharpness of the minima to which GD converged, as a function of the step size As 17 increases, the minima get flatter in parameter space (yellow curve), which translates to smoother predictors function space (purple curve). Panel (b) shows the norm of the bias vector b as function of the step size. Here we see that the bias vector grows with the step size, as the predictor function gets smoother.
EXPERIMENTS
We now demonstrate our theoretical results. We start with regression task synthetic data. Here, we drew = 100 pairs (zz, Vy) in R x R from the standard normal distribution to serve as our training set. We fit a single hidden-layer ReLU network with & 40 neurons to the data using GD with various step sizes (runs were stopped when the loss dropped bellow 10 10-3). ). For each minimum 0" which GD converged, we computed the loss' sharpness, Amax ( (V²)) and our lower bound on the sharpness, + 2 || (Lemma 3 in the appendix). Additionally, we numerically determined the sharpness of the flattest implementation for every solution. Figure 4(a) depicts the results for this experiment. The red line marks the border of the stable region, which 2/1. Namely, (S)GD cannot stably converge to minimum whose sharpness is above this line. The dashed yellow line shows the sharpness of the minima to which GD converged practice. As can seen, here GD converged at the edge of stability (the two lines coincide), phenomenon discussed in (Cohen et al., 2020). The blue curve the sharpness the flattest implementation each solution (see App. A), while the the purple curve is lower bound. We see that our bound is quite tight (blue purple). Furthermore, as the step size increases the minima get flatter in parameter space (yellow curve), which translates to smoother predictors in function space (purple curve). Additionally, see from Fig. 4(b) that the norm of the bias vector b increases with the step size, as our theory predicts (Sec. 3.2).
Next, we present experiment with binary classification on MNIST (LeCun, 1998) using SGD. In this experiment we used n = 512 samples from two MNIST classes, '0' and '1', The classes were labeled as y and y -1, respectively. For the validation set we used 4000 images from the remaining samples in each class. We trained single hidden-layer ReLU network with k = 200 neurons using SGD with batch size B 16, and the quadratic loss. To perform classification at inference time, we thresholded the net's output at O. We ran SGD until the loss dropped below 10-8 for 2000 consecutive epochs. Figure 5(a) shows the same indices as in the previous experiment. Here we see again that as the step size increases, the minima get flatter in parameter space (yellow curve), which translates to smoother predictors in function space (purple curve). Figure 5(b) shows the classification accuracy the validation set, where we see that the network generalizes better as the step size increases, as past work showed, egg., (Keskar et al., 2017). More experiments in App. N.
7 RELATED WORK
Dynamical stability analysis was applied to neural network training in several works. particular, Nar & Sastry (2018) analyzed Lyapunov stability of one hidden-layer ReLU networks without bias. They proved a bound on the network's output which depends on the step size and training data and implies that the network's 's output should be smaller for training samples with larger magnitude. Mulayoff & Michaeli (2020) characterized the flattest minima for linear nets and showed that these minima have unique properties. Yet, in their setting all minima implement the same onput-output
8

---

Published as a conference paper at ICLR 2023
100
10
Accuracy [%6]
90
Sharpness
80
2/11
10a
70
Aux
min. Xmna
60
10
50
10
10-3
10
10
10
10
10-2
10-1
7
7
(a) Sharpness versus step size
(b) Val. accuracy vs. step size
Figure 5: Validating the bounds on MNIST. We trained a single hidden-layer ReLU network for binary classification on two classes from MNIST using SGD (see Sec. 6). Panel (a) depicts the sharpness versus the step size 1j. Here as increases, the minima get flatter in parameter space (yellow curve), which translates to smoother predictors in function space (purple curve). Panel (b) shows the performance on the validation set. Here the trained model generalizes better as the step size increases.
function. Thus, their results only show that SGD is biased toward certain implementations of the same function, whereas our result shows that SGD is biased toward certain functions.
Wu et al. (2018) proved a sufficient condition for dynamical stability of SGD in expectation using second moment. Ma & Ying (2021) extended their result by showing a necessary and sufficient condition for dynamical stability in higher moments. In addition, the authors combined this condition with the multiplicative structure of neural nets to prove an upper bound on the Sobolev seminorm of the model's input-output function at stable interpolating solutions. Their upper bound extends to deep nets, yet it depends on the norm of the first layer of the network, which in general can be large.
Mulayoff et al. (2021) characterized the stable solutions of SGD for univariate single hidden-layer ReLU networks with the square loss. Our Thm. is the natural extension of (Mulayoff et al., 2021, Thm. 1) to the multivariate case. To prove it, we combine the proof technique of lower bounding the top eigenvalue of the Hessian, used by Mulayoff et al. (2021), with the Radon transforms analysis used by Ongie et al. (2020). Combining these techniques is not a priori trivial, since Radon transform was not used before for Hessian analysis. Also, it required several subtle steps that are not encountered in the univariate setting nor in (Ongie et al., 2020) (e.g., working with the inverse of the dual Radon transform to obtain the primal space representation).
Ongie et al. (2020) studied the space of functions realizable as infinite-width single hidden-layer ReLU nets with bounded weights norm. Their settings assumes explicit regularization, i.e., min-norm solution, whereas here we derived our results for SGD withoul regularization, via implicit bias. On the technical level, they introduced the "R-norm" || that is closely related to the stability norm. Particularly, ||- || = || |R.. for 1 the constant 1 function. They proved similar results of depth separation and approximation guarantees, shown here Sees. 4-5. More related work in App. B.
8 CONCLUSION
Large step sizes often used to improve generalization (Li al., 2019). This work suggests an explanation to this practice. Specifically, we showed that large step sizes lead to smaller stability norm and thus can bias towards smooth predictors in shallow multivariate ReLU networks. We find the smoothness measure depends on the data via specific functions or p, and exemplify their properties. Moreover, we studied the approximation power of ReLU networks that correspond to stable solutions. Although shallow networks are universal approximators, we proved that stable solutions of these networks are not. Namely, there is a function that cannot be well-approximated by stable single hidden-layer ReLU networks trained with non-vanishing step size. Yet we showed that the same function can be realized as a stable two hidden-layer network, leading to a depth separation result. This result explain the success of deep models over shallow ones. Finally, we gave approximation guarantees for stable shallow ReLU networks. In particular, we proved that any Sobolev function can be approximated arbitrarily well using GD with single hidden-layer ReL.U networks.

---

Published as a conference paper at ICLR 2023
ACKNOWLEDGMENTS
The research of RM was supported by the Planning and Budgeting Committee of the Isracli Council for Higher Education, and by the Andrew and Erna Finci Viterbi Graduate Fellowship. GO was supported by NSF CRII award CCF-2153371. The research of DS was funded by the European Union (ERC, A-B-C-Deep. 101039436). Views and opinions expressed are however those of the author only and not necessarily reflect those of the European Union or the European Research Council Executive Agency (ERCEA). Neither the European Union nor the granting authority can be held responsible for them. DS also acknowledges the support of Schmidt Career Advancement Chair in AI. TM was supported by grant 2318/22 from the Israel Science Foundation and by the Ollendorff' Center of the Viterbi Faculty Electrical and Computer Engineering at the Technion.
10

---

Published as a conference paper at ICLR 2023
REFERENCES
Sanjeev Arora, Nadav Cohen, and Elad Hazan. On the optimization of deep networks: Implicit accelcration by overparameterization. In International Conference on Machine Learning, pp. 244-253. PMLR, 2018.
Sanjeev Arora, Nadav Cohen, Wei Hu, and Yuping Luo. Implicit regularization in deep matrix factorization, Advances in Neural Information Processing Systems, 32:7413-7424, 2019.
Shahar Azulay, Edward Moroshko, Mor Shpigel Nacson, Blake Woodworth, Nathan Srebro, Amir Globerson, and Daniel Soudry. On the implicit bias of iitiiaiiztiinn shape: Beyond infinitesimal mirror descent. International Conference on Machine Learning. 2021.
David Barrett and Benoit Dherin. Implicit gradient regularization. In International Conference Learning Representations, 2021.
Mohamed Ali Belabbas. On implicit regularization: Morse functions and applications to matrix factorization. arXiv preprint arXiv:2001.04264, 2020.
Emmanuel Candès. Harmonic analysis of neural networks. Applied and Computational Harmonic Analysis, 6(2):197-218, 1999.
Emmanuel J Candès and David L Donoho. Ridgelets: A key higher-dimensional intermittency? Philosophical Transactions of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences, 357(1760):2495 24002 1999.
Sean M. Carroll and Bradley W. Dickinson. Construction of neural nets using the Radon transform.
In International Joint Conference on Neural Networks, volume 1, pp. 607-611, 1989.
Lenaic Chizat, Edouard Oyallon, and Francis Bach. On lazy training in differentiable programming. Advances in Neural Information Processing Systems, 32, 2019.
Jeremy Cohen, Simran Kaur, Yuanzhi Li, J Zico Kolter, and Ameet Talwalkar. Gradient descent on neural networks typically occurs the edge of stability. International Conference on Learning Representations, 2020.
Laurent Dinh, Razvan Pascanu, Samy Bengio, and Yoshua Bengio. Sharp minima can generalize for deep nets. In International Conference on Machine Learning, pp. 1019-1028. PMLR, 2017.
Armin Eftekhari and Konstantinos Zygalakis. Implicit regularization in matrix sensing: Initialization rank matters. arXiv preprint urXiv:2008.12000 2021.
Charles L. Epstein. Introduction to the mathematics of medical imaging. Society for Industrial and Applied Mathematics (SIAM), 2007.
Gauthier Gidel, Francis Bach, and Simon Lacoste-Julies Implicit regularization of discrete gradient dynamics linear neural networks. Advances in Neural Information Processing Systems, 32, 2019.
Justin Gilmer, Behrooz Ghorbani. Ankush Garg, Sneha Kudugunta, Behnam Neyshabur, David Cardoze, George Edward Dahl, Zachary Nado, and Orhan Firat. A loss curvature perspective on training instabilities of deep learning models. In International Conference on Learning Represen- tations, 2022.
Suriya Gunasekar, Blake E Woodworth, Srinadh Bhojanapalli, Behnam Neyshabur, and Nati Srebro. Implicit regularization in matrix factorization. Advances in Neural Information Processing Systems, 30, 2017.
Suriya Gunasekar, Jason Lee, Daniel Soudry, and Nathan Srebro. Characterizing implicit bias in terms of optimization geometry. In International Conference on Machine Learning, pp. 1832-1841. PMLR, 2018a.
Suriya Gunasekar, Jason Lee, Daniel Soudry, and Nati Srebro. Implicit bias of gradient descent on linear convolutional networks. Advances in neural information processing systems, 31, 2018b.
11

---

Published as a conference paper at ICLR 2023
Sigurdur Helgason. The Radon transform. In The Radon Transform. Springer, 1999.
Yoshifusa Ito. Representation of functions by superpositions of step or sigmoid function and their applications to neural network theory. Neural Networks, ((3):3853394, 1991.
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.
Ziwei Ji and Matus Telgarsky. Gradient descent aligns the layers of deep linear networks. In International Conference on Learning Representations, 2019a.
Ziwei Ji and Matus Telgarsky. The implicit bias of gradient descent nonseparable data. In Conference on Learning Theory, pp. 1772-1798. PMLR, 2019b.
Ziwei Ji and Matus Telgarsky. Characterizing the implicit bias via primal-dual analysis. In Algorithmic Learning l'heory, pp. 772-804. PMLR, 2021.
Ziwei Ji, Miroslav Dudík, Robert E Schapire, and Matus Telgarsky. Gradient descent follows the regularization path for general losses. In Conference on Learning Theory, pp. 2109-2136. PMLR, 2020.
Hui Jin and Guido Montúfar. Implicit bias of gradient descent for mean squared error regression with wide neural networks. arXiv preprint arXiv:2006.0733: 2020.
Nitish Shirish Keskar, Dheevatsa Mudigere, Jorge Nocedal, Mikhail Smelyanskiy, and Ping Tak Peter Tang. On large-batch training for deep learning: Generalization gap and sharp minima. In International Conference on Learning Representations, 2017.
Yann LeCun. The MNIST database of handwritten digits. http://aann.eec..coom//xdddddeee
1998.
Yuanzhi Li, Tengyu Ma, and Hongyang Zhang. Algorithmic regularization in parammetrateeeed matrix sensing and neural networks with quadratic activations. In Conference On Learning Theory, pp. 2-47. PMLR, 2018.
Yuanzhi Li, Colin Wei, and Tengyu Ma. Towards explaining the regularization effect of initial large learning rate in training neural networks. Advances in Neural Information Processing Systems, 32, 2019.
Donald Ludwig. The Radon transform Euclidean space. Communications on pure and applied mathematics, 19(1):44 (111))44-8... 1966.
Kaifeng Lyu and Jian Li. Gradient descent maximizes the margin of homogeneous neural networks. International Conference on Learning Representations, 2020.
Chao Ma and I Lexing Ying. On linear stability of SGD and input-smoothness of neural networks. Advances in Neural Information Processing Systems, 34, 2021.
Cong Ma, Kaizheng Wang, Yuejie Chi, and Yuxin Chen. Implicit regularization in nonconvex statistical estimation: Gradient descent converges linearly for phase retrieval and matrix completion. In International Conference on Machine Learning, pp. 3345-3354. PMLR, 2018a.
Siyuan Ma, Raef Bassily, and Mikhail Belkin. effectiveness of SGD in modern over-parametrized learning. In International Conference on Machine Learning, pp. 3325-3334. PMLR, 2018b.
The power of interpolation: Understanding the
Paul Malliavin, Hélène Airault, Leslie Kay, and Gérard Letac Integration and probability, volume 157. Springer Science & Business Media, 1995.
Rotem Mulayoff and Tomer Michaeli. Unique properties of flat minima in deep networks. In International Conference on Machine Learning, pp. 7108-7118, PMLR, 2020.
Rotem Mulayoff, Tomer Michaeli, and Daniel Soudry. The implicit bias of minima stability: A view from function space. Advances in Neural Information Processing Systems, 34:17749-17761, 2021.
12

---

Published as a conference paper at ICLR 2023
Mor Shpigel Nacson, Suriya Gunasekar, Jason Lee, Nathan Srebro, and Daniel Soudry. Lexicographic and depth-sensitive margins in homogeneous and non-homogeneous deep models. In International Conference Machine Learning, pp.
PMLR, 2019a.
Mor Shpigel Nacson, Jason Lee, Suriya Gunasekar, Pedro Henrique Pamplona Savarese, Nathan Srebro, and Daniel Soudry. Convergence of gradient descent on separable data. In The 22nd International Conference on Artificial Intelligence and Statistics, pp. 3420-3428. PMLR, 2019b.
Mor Shpigel Nacson, Nathan Srebro, and Daniel Soudry. Stochastic gradient descent on separable data: Exact convergence with a fixed learning rate. In The 22nd International Conference Artificial Intelligence and Statistics, pp. 3051-3059. PMLR, 2019c.
Kamil Nar and Shankar Sastry. Step size matters in deep learning. Advances in Neural Information Processing Systems, 31, 2018.
Greg Ongie, Rebecca Willett, Daniel Soudry, and Nathan Srebro. A function space view of bounded norm infinite width relu nets: The multivariate case. In International Conference on Learning Representations, 2020.
Rahul Parhi and Robert Nowak. Banach space representer theorems for neural networks and ridge splines. Journal of Machine Learning Research, 22(43):1-40, 2021.
Allan Pinkus. Approximation theory of the MLP model in neural networks. Acta nunerica, 8: 143-195, 1999.
Noam Razin and Nadav Cohen. Implicit regularization in deep learning may not be explainable by norms, Advances in neural information processing systems, 33:21174-21187 2020.
Michael Shub. Global stability of dynamical systems. Springer Science Business Media, 2013.
Samuel L Smith, Picter-Jan Kindermans, Chris Ying, and Quoe Le. Don't decay the learning rate, increase the batch size. In International Conference on Learning Representations, 2018.
Samuel L Smith, Benoit Dherin, David Barrett, and Soham De. On the origin of implicit regularization in stochastic gradient descent. In International Conference on Learning Representations, 2021.
Donald C Solmon. Asymptotic formulas for the dual radon transform and applications. Mathematische Zeitschrift, 195(3):321-343, 1987.
Sho Sonoda and Noboru Murata. Neural network with unbounded activation functions is universal approximator. Applied and Computational Harmonic Analysis, 43(2):233-22) 268, 2017.
Daniel Soudry, Elad Hoffer, Mor Shpigel Nacson, Suriya Gunasekar, and Nathan Srebro. The implicit bias of gradient descent on separable data. The Journal of Machine Learming Research, 19(1): 2822-2878, 2018
Gal Vardi and Ohad Shamir. Implicit regularization in ReLU networks with the square loss. In Conference on Learning Theory, pp. 3224-4255 PMLR, 2021.
Mathukumalli Vidyasagar. Nonlinear systems analysis. Society for Industrial and Applied Mathemat- ics (SIAM), 2002.
Blake Woodworth, Suriya Gunasekar, Jason D Lee, Edward Moroshko, Pedro Savarese, Itay Golan, Daniel Soudry, and Nathan Srebro. Kernel and rich regimes in overparametrized models. In Conference on Learning Theory, pp. 3635-3673. PMLR, 2020.
Lei Wu, Chao Ma, et al. How SGD selects the global minima in over parameterize learning: A dynamical stability perspective. Advances in Neural Information Processing Systems, 31: 8279-8288, 2018.
Chiyuan Zhang, Samy Bengio, Moritz Hardt, Benjamin Recht, and Oriol Vinyals. Understand- ing deep learning requires rethinking generalization. In International Conference on Learning Representations, 2017.
13

---

Published as a conference paper at ICLR 2023
A ADDITIONAL DISCUSSION
The independence of Lemma 1 and Theorem 1 on the batch size B. Theorem 1 relies on Lemma (see App. E), which was proved in (Mulayoff et al., 2021). This lemma states that if a minimum 0 is linearly stable for SGD with batch-size B, then the Hessian H of the loss at O" must satisfy Amax ( H)) < 2/1,, where 1 is the step-size. Importantly, this necessary condition holds true for any batch size B. and thus Theorem 1 ís independent of B. We note, however, that the precise stability threshold of SGD might depend on B. yet the important points to notice are: (1) Here all we need is necessary condition, and Lemma provides a simple bound that holds for any B. (2) Empirical evidence shows that there is not much room for improvement upon this batch- balch-size-indepeneen bound in real-world settings. Specifically, for practical batch sizes, the gap between 2/1 and the stability threshold of SGD is often very small (see Fig. 5 in our paper and Figures 2-3 in (Gilmer et al., 2022)).
The proof of Lemma is actually quite short and easy to follow (see (Mulayoff et al., 2021, App. II).
The idea is that if 0 is an linearly stable minimum, then by definition we have
lim sup |[ 0. 0* s,
(21)
(-
where {0,) are governed by the linearized stochastic dynamics given in Eq. (5). Using Jensen's 's inequality, for all & we get [[E[O.)] 0"| s [E[[0 0" [I] Thus,
sup
E [. 0"l S lim sumEs 0. - Se. Se.
(22)
Note that under the linearized dynamics, {E 0 o are precisely GD steps. Therefore, we have that if 0" is linearly stable for SGD, then it must be linearly stable also for GD. Now, well-known fact is that 0" is linearly stable for GD if and only if Amax (H) < 2/1. This is how we get the necessary condition in Lemma 1, which does not depend on the batch size.
The independence of Lemma and Theorem 1 on e. Lemma states that a necessary condition for a twice differentiable minimum to be e linearly stable is >max (V²L(0") < 2/1). That is, the condition does not depend on which might seem not intuitive. However, the reason Lemma I does not depend on E is because refers to linear stability, opposed to non-linear dynamical stability. In linear stability for twice-diffi diicediifferentablle minima, all we care about is the second-order Taylor approximation of the loss at the minimum. In see previous paragraph we explained that Lemma gives a necessary condition through reduction to GD. Now, when applying GD on a quadratic loss (with a PSD matrix), for any E only one of two things can happen:
1. Either 30, E (0")) lim sup 100, - 0*|| = +00 (unstable for any & > 0),
t-+00
2. or VO, E (0"): lim sup |00 0'| < 0o 0"|l < (stable for any >> 0).
(+220
In any outcome, the result does not depend on E, and therefore does not appear in the result of Lemma 1. Note that for con-differentiable minima, does affect linear stability in GD, however Lemma 1 only refers to twice-differentiable minima.
Theorem I is based on Lemma and therefore does not depend on e. Yet, beyond this technical reasoning, it is important to note that here we consider interpolating solutions. For those solutions, the global minimum of the loss is also a global minimum w.r.t. each data sample (æx My;) separately. Therefore, despite the stochasticity of SGD, every step points towards a global minimum. This implies that if the stability criterion is satisfied, then SGD converges to the minimum (lim sup E[] 0, - 00
t-
0) and if it is not satisfied, then SGD repels from the minimum (lim sup E[ 0" oo for the (-
linearized dynamics). This is also seen in simulations where models are overfit to training data using SGD, egg., (Ma et al., 2018b). Particularly, in our simulations the loss always converged to 0 when it converged (we arbitrarily decided to stop each run when the loss dropped below 10-*). In the general case of non-interpolating solutions, the expected final distance the minimum in mini-batch SGD (lim sup [[]] 0* II]) can be strictly positive finite number, and therefore in those cases e does to
play role.
14

---

Published as a conference paper at ICLR 2023
Large step size training and warmup. While Theorem 1 applies to any positive step size, it is most interesting when considering large step sizes. High learning rates are standard practice, as they are associated with good generalization (Li et al., 2019). However, there cases, e.g., large initialization, in which high learning rate might cause the training to diverge. In these cases, a learning rate warmup is applied, enabling training with large step sizes.
Learning rate decay. Practitioners often work with learning rate schedule, which typically reduces the step size toward the end of training. this scenario, 2/n can be quite high at the end, making Theorem 1 loose. Here, empirical evidence shows that when reducing the step size at a late stage, the sharpness of the obtained minimum is often still controlled by the initial step size, e.g. Figs. 1 and 3 in (Gilmer et al., 2022). Moreover, although learning rate decay is popular technique, there are other popular training schemes in which the learning rate is not reduced, e.g., (Smith et al., 2018). Lastly, as for the depth separation results (Sec. 4) and the approximation results (Sec. 5), they apply for any fixed positive step size. In other words, the learning rate decay does not affect these results.
Initialization independent results. Our results are independent of the initialization. In the past was shown that initialization can have large effect which minimum GD converges under certain conditions. However, this do not contradict our results, as explained below.
For very small step sizes, the GD trajectory follows that of gradient flow (GF). Under certain conditions, e.g., infinite width or vanishing initialization, it was shown that the network does not change much the evolution trajectory of GF. In this case, the initialization dominates the properties of the obtained solution. This is known as kernel regime or Neural Tangent Kernel (NTK) regime (Jacot et al., 2018; Chizat al., 2019). However, for practical step sizes and standard initialization, recent work (Cohen et al.. 2020) showed that GD typically deviates from the GF trajectory, while entering the Edge of Stability regime. This occurs when the stability threshold achieved during training, i.e., Amax ( (V)) > 2/n for some & > O. In this case, GD converges to a different minimum than GF, i.e., GD escapes the NTK regime. Similar behavior was shown also for SGD (Gilmer et al., 2022).
Theorem 1 Proof idea. Theorem 1 is a result of two properties of twiicediffereeniablle minima. First, we show in Lemma 3 in the appendix that these minima satisfy Amaxx (222) 2 1+2 || R Second, we know from Lemma 1 that stable minima satisfy Amax (22C) I 2/n. Together, these properties imply that 2/17 2 1++2 2 ||fll from which we get the result of the theorem. |l S 1/n 1/2. Note that twice-differentiable minima correspond to functions whose knots do not coincide with any training point. Although we prove our result only for such functions, we observed that in practice the condition Amax (V2 I 2/1 is always met around minima to which SGD converges (see Sec. 6 and the next paragraph). For full derivation of the theorem see App. E.
Theorem assumption and minima that are not tiiee-differenttable. Theorem assumes that the knots of f do not coincide with any training point. This assumption is done for technical simplicity, that is to ensure the minimum is twice-differentiable. In practice, we observed that usually only small fraction of the knots coincide with training points. For example, in our MNIST experiment we had only 30 training points coinciding with knots of ʃ (out of n 512 training points and k = 200 neurons). Note that the differriffeettiabblitt assumption in Theorem 1 is required for Lemma to hold as Theorem invokes Lemma 1. However, we observed that in practical settings, Lemma I still applies in settings where the assumption is violated. This can be appreciated by the red curve upper bounding the orange curve in Figs. 4(a) and 5(a). Namely, despite the fact that the minima not always wwice-differentiabll in our experiments, the stability criterion for SGD is still observed to be upper bounded by 2/n.
It is possible to extend the analysis to minima that are not twice-differentiable by using the same method as in (Mulayoff et al., 2021). However, this makes the analysis much more complex.
The meaning of min Xmax in figures 4(a) and 5(a). The curve min Amax shows the sharpness of the flattest implementation of each solution f. In more detail, as discussed earlier, Theorem is a result of two properties of twice-differentiabll minima. On the one hand, we know from Lemma 1 that stable minima satisfy Amaxe (V2C) C) < 2/n. On the other hand, Lemma 3 in App. E asserts that (V3
these minima satisfy Xmax (V) 2 1 + 2 Note that each function ʃ E Ft has multiple
15

---

Published as a conference paper at ICLR 2023
implementations, i.e., different minima in parameter space which all correspond to f. These minima can have different sharpness. Here, we can look at the best implementation, i.e. a solution to min Amax, where the minimum is taken over all loss' minima {0} that implement f. Overall, given minimum A with a corresponding function f. we have
1 + 2 I min Amax V2C) S Xmax V22)= S 2/n.
(23)
These inequalities give us the result of Theorem 1, 1 + Rl S 1/17 1/2. To understand the tightness of each part of our analysis, we added to the plots Amax and min Amax- In both figures, Amaxx V C) equals just below 2/n, phenomenon known edge of stability (Cohen et al., 2020). Additionally. min Amax (V)C) is close to 1 + 2|lfRR.. in these experiments. Yet, Figure shows that Amax C) can be quite larger than min A (V)C), meaning that there exists a far flatter minimum that implements the same function. This fact was used by Dinh et al. (2017) to show that sharp minima can generalize.
B ADDITIONAL RELATED WORK
Implicit bias. A long line of works studied the implicit bias of the training procedure in attempt to better understand generalization in overparameterized models. For the classification setting, in the case of linear prediction function, linearly separable data, and exponentially tailed loss functions (e.g. logistic and exponential), Soudry et al. (2018) 8) showed that GD converges in the direction of the SVM solution. This result was later extended to linear fully connected and convolutional neural networks (Gunasekar et al., 201 8b; Ji & Telgarsky, 2019a), more loss functions (Nacson et al., 2019b; Ji & Telgarsky, 2021), SGD optimization algorithm (Nacson et al., 2019c), other generic optimization methods (Gunasekar et al., 201 8a), sen-repe data (Ji & Telgarsky, 2019b), and homogeneous prediction functions (Nacson et al., 2019a; Lyu & Li, 2020; Ji et al., 2020). However, all those results do not depend on the step size, except for the requirement that it be sufficiently small.
Another line of works studied the implicit bias in the context of linear models with quadratic loss such matrix factorization (Gunasekar et al., 2017: Li et al., 2018; Arora et al., 2018; 18; 2019; Belabbas, 2020; Eftekhari & Zygalakis, 2021; Gidel et al., 2019; Ma et al., 201 2018a; Woodworth et al., 2020; Azulay et al., 2021). However, all of these works relied on either small infinitesimal step size (i.e.. gradient flow). Thus, they do not capture how the step size affects the implicit bias. Moreover, they assumed a manifold property (Azulay et al., 2021). As pointed out by Razin & Cohen (2020) and Vardi & Shamir (2021) 1). these ussumptions do not always apply. In contrast, our result is based on a stability condition of SGD, which depends on the step size and does not require the manifold assumption.
How the step size affects the implicit bias. To investigate the implicit bias of the step size, Barrett & Dherin (2021) and Smith et al. (2021) suggested using a modified loss. Under this modified loss, gradient flow approximates the trajectory of (S)GD on the original loss. However, the step size should be sufficiently small for the approximation to hold true. Moreover, the induced regularization term this method yields is expressed in terms of the model's parameters. Additionally, this term increases linearly with the step size and vanishes at any stationary point.
Radon transform analysis of shallow networks. Radon transform analysis has previously been used in studies of the approximation capabilities of single hidden-layer neural networks with bounded activation functions (Carroll & Dickinson, 1989; Ito, 1991), and more general activation functions in the ridgelet framework (Candès & Donoho, 1999; Candès, 1999). More recently, Sonoda & Murata (2017) used ridgelet transform analysis to study the approximation properties of two-layer neural networks with unbounded activation functions, including the ReI.U.
Parts of this work extend results by Ongie et al. (2020), which defined a similar Radon-domain seminorm (the "R-norm") to determine the space of functions realizable as an infinite-width single hidden-layer ReLU networks with square- squaaee weights Parhi & Nowak (2021) proved representer theorem for single hidden-layer ReLU networks using the R-norm. Finally, an L version of the R-norm is used by Jin & Montúfar (2020) describe the function space implicit bias of training a single hidden- layer ReLU network using gradient descent in the neural tangent kernel regime.
16

---

Published as a conference paper at ICLR 2023
C THE RADON TRANSFORM
For a function f: TR + R. the d-dimensional Radon transform Rf is the collection of all integrals of f over (d - (d= diddinnniin affine hyperplanes in Rd. Every hyperplane pair (v, b) E gd-1 x I. where v is unit normal to the hyperplane and b E R is its distance from the origin. Therefore, the Radon transform R. is the function over (v, b) E Sd I given by
be parametrized
R. (v, b) A [ TT==
f(x))d(()))
(24)
where ds(z) represents integration with respect to the (d 1) dimeensional surface measure on the hyperplane v x = b. Note that the Radon transform is an even function, i.e, R.f R,,(- b), since (v, b) and (-v, bb) describe the same hyperplane.
The dual Radon transform RD maps functions defined on S IR to functions on R2 by
{ 1p (u, v ds(v)
R" (²)((x)
(25)
ga-1
for all æ R4, where ds(v) represents integration with respect the (d d-ddienssion surface measure on the unit sphere Sd-1
The Radon transform and its dual invertible over spaces of smooth functions via the inversion formulas:
R 7a(-A) = R
(26)
(R") yaR( ( A) ,
(27)
where the fractional Laplacian operator ( (-A) + is defined by application of ramp function in Fourier domain (i.e., multiplication by || in Fourier domain), and
ad 2(2m)d-1
(28)
is a dimension dependent constant.
These transforms may be extended to spaces of distributions (e.g., Dirac deltas) in standard way (Ludwig, 1966; Helgason, 1999), which we summarize in App. D. Important for this work is the distributional dual inverse Radon transform (R*)-1 which maps a distribution defined over Euclidean space Rd to distribution in Radon domain gd-1 x R.
D DISTRIBUTIONAL FRAMEWORK
f IR + IR be any locally integrable function. Then its Laplacian Af can be interpreted tempered distribution, meaning that A defined via the duality pairing
(AS. v) A (.. (p) = [(a [(x) Axp(x) da.
(29)
where is any Schwartz test function on Rd i.e., smooth function such that the function and its partial derivatives of all orders have sufficiently fast decay at infinity; denote this space of functions by S(R4). For example, if f consists of single ReLU unit, i.e., f(x) a(v b) such that
|| = 1, then is easy to show Aʃ 5(v æ b), meaning (Af.4) f( == s=x)) da = R). In other words, Af is the distribution given by evaluation of the Radon transform of test function at the point (v, 6) gd-1 1 x
Next, we describe how to understand the operator (R*)-1 in distributional sense. Let Sy (Sd-1 R) denote the image of Schwartz functions S(Rd) under the classical Radon transform R. The space SH (S x RR is characterized in (Ludwig, 1966; Helgason, 1999); it is the space of all even Schwartz functions defined on gd-1 x & that additionally satisfy some moment conditions It is also shown
'Specifically, for all positive integers k; the function v -+ Sica 6(v, b)b* db needs to be homogeneous polynomial of degree k.
17

---

Published as a conference paper at ICLR 2023
by Ludwig (1966) that the classical inverse Radon transform R-1 is linear homeomorphism of Su (S R) onto S(R"). Therefore, we may define its distributional transpose (R-1)* = -R*)
applied to any tempered distribution h by
((R")- o) (h, (h,))
(30)
for all @ E SH(S" x R). Note that (R*)~1h h is distribution belonging to S/ (Sa-1 x R), the topoloical dual SH(S x R).
Returning to the example where f is single ReLU unit, i.e., f(x) a x b) with ||l| = 1,
then for any test function E Sy (Sd-1
x R we have
) Φ) = (AS, R o) = (RR | = 06( b).
(31)
This shows (R b). i.e., Dirac delta centered at b). If f(x)) - (0 x
b,) + is any single hidden-layer ReLU network such that ||l/|| = 1 for all 1, ,,,. k, then by linearity we have
RR") (R*)
Σ a,o (,o(v..)))
(32)
Finally, we may define the total variation
||
||TV for any distribution E (S4-1 xR by
||a|lr
sup
Ka, d)
(33)
S (5- xx3]
If [|y]|TT finite, then a a distribution of order-O. this case, since Sy (Sd-1 the space of even and continuous functions on gd-1 x R that vanish at infinity, can be extended uniquely to an even signed measure on 3d-1 x R, and ||a||ry is equal to the total variation norm of a.
R) dense in
E PROOF OF THEOREM 1
In the proof of the theorem we use the following lemma (for the proof of this lemma Appendix F). Lemma (Top eigenvalue lower bound). Let f E Fx be twice-differentiiabl minimizer of the loss function, then
Amnax (V)C) 2 2
(34)
where I|-IIR.. denotes the stability norm.
Let f F& be a stable solution of the loss function. Then, according to Definition 3, there exists a linearly stable minimum point 0 E (dd+))k++ such that the network at this minimum implements f. Due to the fact that the knots of f do not contain any training point, we have that 0 is a twice- differentiable minimum. From Lemma 1, since 0 is twice differentiable stable minimum then
Amax (V2C) < -
(35)
On the other hand, from Lemma 3 we have that
Amax (VQC) 2 112 1f
(36)
Using (35) and (36) we get
1 |l
-- ___ 1
(37)
1) 2
F PROOF OF LEMMA 3
The proof of the Lemma consists of the following steps:
1. Calculating V&CC. and showing that at global minimum takes the form V3C (Appendix F.1).
-
18

---

Published as a conference paper at ICLR 2023
2. Lower bounding Amax \ by using
2
mmax (V²c)=
max #ES(++2)
(ViC)
e=
max DEB(++2) 72
+
max u
[
72
(38)
and lower bounding the right hand side (Appendix F.2).
3. Simplifying the lower bound to obtain more interpretable version which does not depend on the specific implementation of f (Appendix F.3).
F.1
HESSIAN COMPUTATION
Recall that
1 " (())
Σ ==1
(s (x,) - uz)²
(39)
2n
where
((x) Σ
(2)
x w
(1)
+ +1(1)) (22)
(40)
We denote
(1))
(1)
E BRdk
6(1) = [0 11 (
E ERE
w ;
w
w
(2)
w
w
RR 6(2) E R
(41)
and
(w()) vec
0
(1)
(d(42)4+1
(42)
(2(
b(2)
Using these notations, assuming that 0* is twice differentiable global minimum of L, we have that the gradient is
n
VoL
(f(x)) 9)) Vof (x,)
(43)
The Hessian is given by
T + 1 n
Vof (x,) Vol (x₂)
Σ j=1
(f (2,) 8) (z,)
n j=1
Vol (x,) Vel (æ,)
(44)
where in the last transition we used Vj E |n] f(x,) = y (see Def. 2). From direct calculation we obtain 1=1
(w (2) o I (x: 0)) & æ
0f owar) f V mtan f Bea f
vec
w 0I
Vof [x)
((MI)) x ++6
T
(45)
(x
where C denotes the Hadamard product, & represents the Kronecker product and Rd
R(d+2)++1 {0, 1 is the activation pattern of all neurons for input æ. namely [[(;; 0)]; if
x T + 6 > 0 and [I[[;; 0)}, = 0 otherwise. Let us denote the tangent features matrix by
(1)
+- oo (x1) Vof (xz)
Vof (æ,)] E +2++
(46)
Then the Hessian can be expressed as V2L = ΦΦ T /n, and its maximal eigenvalue can be written as
Amax (2
max #ES(d+2)
V2L
max +C
T T
max u
(47)
19

---

Published as a conference paper at ICLR 2023
F.2
LOWER BOUNDING THE TOP EIGENVALUE
Continuing from the previous section's 's calculation, if we take u - J/2 1, we obtain
1
[
MAx MER
$1 ||
2
2
14
+
R
w
(2)
Σ
n
mi
3
+
12
x
w
(1)
+ (00))
j=1
j=1
2
2
11
(wf())
Σ
Fzgy 5ز
Σ
j
Σ
23 w
T (1) b5"))
2
2
1
2
k
(2)
Σ
d
Σ
+
Σ I
Σ
2
T
w
+
(48)
1=1
where in (*) we used a² + B² 2 2 |o:91 Let C. S {23} } be the set of training points for which the ith neuron is active, and denote nu = [c:], that is
n2
j.
(49)
j=l
Then,
Amax
(V2)) 21
2
w
(2)
Σ
+ 2
Σ
2 w
+ ((2)
t=1
Ec
(2)
n2
Σ x =EC;
1
Σ
(1) + b
= 11 + 2
W
+1
T4i
i=1
k
1 + 2
w
(2)
(P E 2 V | E CIIl +
IE
w + bg x E C.
(50)
where X is random sample from the dataset under uniform distribution. Next, we define
(1
w
w
(51)
w
w
Using these notations we obtain
Amax (V²L) 2 +1++
w
w
(1)
(P
E c.))² VIEEXIX E C,ll +1
x E X ( x EC
1+2 2
w
1=1
(x4)
1+2
g
(52)
20

---

Published as a conference paper at ICLR 2023
where in (+) and (**) we defined, respectively,
9 w. b) IP x > 8)
x
b
x
w
[
1,
9 (w b) min (ā io, b) (-w,
(53)
From derivation far, obtain that
Amax (V)2) 2 1++2
w
(2)
w
(1)
w
J.
(54)
(2) We denote a, = |10 (1) and the representation dependent stability norm
So A
la:
i 1 .b)
(55)
F.3 IMPLEMENTATION FREE LOWER BOUND
In this section, our goal is to give a simpler lower bound of the multivariate stability norm Sg. that does not depend on the specific representation of f. Let be the signed measure over ga-1 x E given by
Ei=1 a,d (w." (w.").))) ) and whose total variation measure |a| given by |a]
(w." ((0) Recall that Ets l0:10 (001
& Σ
f(x)
w
(1) (2)
α w0 (1) -1, (10)) + ((2)
aga x w
-) + t600) (56)
1=1
1=1
and thus
82 (x))
Af(x) Σ
1=1
k
Σ
1 T (1)
((T)
w
i=1
/
a (w, 8 (w To b) ds d
88-113
(w, w x) ds(w)
33-1
RRaa.
(57)
Namely, Af f is weighted sum of Diracs supported on hyperplanes. From the last equation obtain (R* (R*)-IAA. Combining these results we get that
So
/
gdla| ([al.9) >
[(R)) (v, (v. ds(v)db = ||
R.
33d 1xR
J3d 1xR
(58)
where the inequality in the third step is due to g being non-negative and the scenarios in which multiple deltas become active at the same location, namely 3i 전 j w. = w(" and 8(1) = 8
Note that if the deltas do not align, then Sg = ||
RR." Overall we have that
Xanax ViC) 2 1 2Sg 2 1 2
(59)
G DERIVATION OF THE STABILITY NORM IN PRIMAL SPACE
We defined the multivariate stability norm as || = (((R))-AAf A.f..gse 1x2, where (5)) xx denotes the integral inner-product Sd-1 R. Supposing that the inverse Radon transform of
21

---

Published as a conference paper at ICLR 2023
exists, then by purely formal reasoning we ought to have
l = ([(R)) 9)ga-1x8
= (1())) Aʃ| R.R. 99)8--1
(R' (R)) ,R (2-99) pd
(60)
Therefore, making the (formal) definitions | R (R((R))-1 and = R- 9. we may also interpret the stability norm [|/l| R. as the quantity
(
(61)
In the event that f and g are smooth, and is in the range of the classical Radon transform, then the above expression equal [..l However, this is not generally the case in our setting, and below we show how give more precise interpretation of this integral formula using distributional theory. In particular, we show that when f is a finite-width ReLU network | is equal to the total variation measure of Af (i.e., the measure- theooreet analog of the absolute value of a function). Additionally, in the event that g does not have classically defined Radon inverse, we show how the integral in (61) can be interpreted using a smoothing approach.
Let f a finite width single hidden-layer ReLU network, i.e., f E FR for some finite k. Recall that (R*)-IAf is finite weighted sum of Diracs in gd-1 RR Let ((R*)-1AA) be the associated total variation measure, and define A = R (-' Afl. where R is the distributional dual Radon transform. Here | is a tempered distribution given by a (positive) weighted sum of Diracs supported on hyperplanes. For example, if f a single ReLU unit of the form f(x) = a (a: b) with || la 1, then | fl is the distribution |a|8(x v - b), that is, for any test function é we have
[ARR D)R= lal
c( = aa [(
(62)
wubb
On the other hand, treating Aʃ measure defined over a compact subset of R2 its total variation measure |All is also equal to |a (a(6(x)) To - b). The following result shows that, more generally, when f is any finite width ReLU net then AflR and |Afl are equal measures.
Proposition 4. Let| J E Fs, then |s| R. = [ as measures defined over any compact set of TR..
Proof. Since f E FR. there exists a representation of as f(x)= = (,o(v)) x by) + x q + c where ||w;|| = 1 for all i, # 0 for all i, and (v,, b.) * =(V,, by) for all * (i.e., the knots of all ReLU units are distinct"), and where K' < k. Therefore, each ReLU unit in this representation of f maps a distinct Dirac ((.),b) in Radon space after applying the operator (R*)-IA and
=
(R))~~ = Σ a ((r.))))
(63)
Let X be any compact subset of TRdd and let ø be any continuous test function defined over X. Then,
(LASR Φ) (((, Σ |a (Rp(,, bz).
(64)
Now, we show the same equality holds with | in place of | Let 1+ denote the set of indices i such that a > 0 and 7- denote the set of indices such that o, < O. Define the measures 44 = =ierr á, -b,) 6,) and μ = Σier- a; ó(v. -- b²). Observe that μ4 and μ are both positive measures whose supports only possibly intersect set of measure zero, and A p4 -L-.
This implics the total variation measure of A is given by A = u+ + a_ = E 1 |a (a,] b,). Hence,
(Af), Φ)
Ro)
(65)
as claimed.
'If (v,, b:) = (wy, by) then one of the neurons is redundant. If (v., b.) = -w,, by) then these units can be combined into an affine function, which we "absorb" into the term x q
22

---

Published as a conference paper at ICLR 2023
Note, however, that when f corresponds to finite-widh ReLU network, | does not have finite total variation considered as a measure defined over all R Due to this technicality, when p is not compactly supported, we need understand the integral in (61) as being with respect to the distribution f(x in place of the measure |
Now show that when p R does exist in classical sense, the integral in (61) can still be interpreted using a smoothing approach.
Proposition 5. Let J E FR and suppose g is an even, piecewise continuous T,2 function on 2d-1 x IR and let p = R~ does not intersect the set of points where g is discontinuous. Then [lfll is finite and
distibutional Radon inverse. Further, assume the support of (R*)
ll -fe- (Af(=))o(())dd,
(66)
JR"
where the integral above is understood the finite limit
lim c-0
(|As (R, Pc)
(67)
where Pe is smooth approximation of p defined independently of f whose classical Radon transform
Ic = R.p exists for all e > O, and Ge - 9 uniformly as e -y on any closed subset of TRd over which g is continuous.
Proof. For any e >> 0 let o. S(gd-1 x R) be a compactly supported even function acting smooth upproximation of the identity, i.e., for any continuous, even function vanishing at infinity we have o, * f + h uniformly as e + 0. where * denotes convolution of functions on S x R. Define Ge = ( g) Xc. where Xe (11 b) is smooth cutoff' function that is equal to one if [b] S 1/c and rapidly decays to zero for bb| 2 1/c. Observe that J. is an even Schwartz function by construction. Furthermore, since g is piecewise continuous and L (and in particular, vanishes at infinity), this implies 9e g uniformly over any closed set that does not intersect the set of points where is discontinuous. Therefore, if we let U C Sá-1 x R be any closed set containing the support of I(R))- A that does not intersect the set of points where g is discontinuous (which is guaranteed to exist since the support of (RR*)=1 A is finite set), then we have Je + g uniformly over U. Therefore
|| =
9/ = lim c>
-
(68)
where the limit guaranteed to exist since the finite measure ((R')-'Af Afl is a continuous linear functional over CUUU), the space of continuous functions over U vanishing infinity. Finally, since ge is Schwartz, (Solmon, 1987, Thm. 7.7) guarantees P. = R RR exists as C -moothh function on Rd that is also integrable along hyperplanes and for which the classical Radon inversion formula holds: Rp. = J.. Therefore, we have
|l = lim (I(R*)-1 Af1, R.fl,RReess Rec)sa--xE
-+0
lim (R" (RR) Pe)z2 c-0
lim
(69)
c++0
as claimed.
The assumption made above that the support of (R*)-1Af does not intersect the set of points where g is not overly restrictive. For example, this assumption holds when f corresponds to differenniable minimizer of the squared loss defined in terms of a finite set of training points and is the data dependent weighting function defined in (12). In this case, the discontinuity set of g b) corresponds to the set of hyperplanes {æ Rd x v = b} that intersect one or more the training points. And f corresponds a differentiable minimizer if and only if the hyperplanes defined by the knots of the ReLU units making up f (i.e., the support of (R*)-11 do not intersect any training points.
H EXAMPLES OF g AND P
H.1 TWO DATAPOINTS
For this example, it is easy to calculate that g(v,b) b)
(lua) 1b1),
(70)
23

---

Published as a conference paper at ICLR 2023
where a is positive constant. We compute R RR'9 by first determining its Laplacian, Ap, then inverting and Laplacian to recover p. First, the intertwining property of the Laplacian and the Radon transform gives
0 RAp= = ōb
8
Rp
Ob29.
(71)
Therefore, by the Fourier slice theorem (Helgason, 1999), for all (0, sd-1 x IR we have
a ag (1,
((A))(sn) = 86
(72)
where 3{-} is the 2-D Fourier transform, and Bo is the Fourier transform in the b variable. For fixed v, the function b -> g(v,,)) is continuous and piecewise linear with knots 0 and +/11], and it easy to see that
0 )e b) a(8(b- (0 + 8(6 loi1) 20(b))
(73)
which implies
{ { . (v, =a J2 + denlvsl -
(s) 30
(74)
we restrict the unit-norm vector v (Ui, 02) to such that en 2 0 and define & = so then we see that a & = s|vg and a. & -8/11/ Therefore, we have
((AP)(E) a (e-)2maTe - -
(75)
and inverting the Fourier transform gives
Ap(x) (5(x x1) + 8(x x2)) 22(x)))
(76)
Finally, since (xx) =
log(1|=]]) is the fundamental solution of Poisson's equation in 2D (i.e.,
AØ 6, where ó a Dirac centered origin), we have
p(x)
[log[[ x1 Il) log( x2 log([[=))))
(77)
2T
While each term in the sum above not absolutely integrable along lines, their sum is absolutely integrable along lines. This is because the function + log(It)) is absolutely integrable over any
neighborhood of the origin, and by multipole expansion we may show that p(xr) = O [||ar] 122 as + oc, which is also absolutely integrable along lines.
H.2 ISOTROPIC DATA DISTRIBUTION
For an isotropic data distribution, i.e., P (xx (x = =M (b) for any that satisfies ||a|| = 1, we will have that
2 @ (v, -M(b) ( M(s)dey/ (o+ M(b) / M(s)dz)²
+1
(78)
Then, from symmetry and assuming that ğ is decreasing we obtain
g(e, b) ((66 / M((s)d (11+ M( / M())dz + 1.
2
1
(79)
Note that only depends on |b| and is decreasing in (b). Considering the parameter space representa- tion for the stability norm
b So la g (v, bz)
(80)
we can see that solutions with larger [b,], i.e., solutions which are more flat function space, will have smaller stability norm.
24

---

Published as a conference paper at ICLR 2023
We now characterize p R. g. For simplicity, we focus the two-dimensional setting (d = 2). Since g b) does not depend on v. we drop this dependence and simply write g(b). Note that this implies is radial function. Let p(r) denote the radial profile of p. i.e,, p(x) =( We additionally make the following assumptions; g(b) is twice continuously differentiable away from the origin, and both g and its weak derivative g' are bounded and absolutely integrable. In this case, P has the integral formula"
g'(b) p(r) - - VB2 - p2 db.
(81)
The assumptions on above are sufficient to show the integrand in (81) is absolutely integrable over [r. 00) for > 0. Since ((b)) is assumed to be decreasing in (b), we have (g(b) 2 for all >
which shows p(r) 2 0 for all > O. However, if is not smooth at the origin, then 9 g'(b) = (1)
as b + 0+, and elementary analysis shows p(r) = O(log(r) as r - 0+ and p(r) (11r)) as r +
Finally, if we additionally assume g' (bb is non-increasing for > O, then p(r) is strictly decreasing for 0. To see this, fix any r' > т, and define 8 = r' -r. Using the change of variables b + b 8, we may show
g((b+ db.
p(r))
(82)
(b+ p2
Since g' is assumed to be non-increasing, we have g'bb+ 5) S (b) and it is elementary to show ((b+ + 8)2 - - < (b2 - -)p- -1/2 for all > r, which shows the integrand in (82) is pointwise strictly bounded above by the integrand in (81) for all b > r. hence p(r') < p(r).
In the case of 2D Gaussian distributed data X ~ NO, I). then M()) is the complementary of the CDF of normal random variable: (b) v2x fh e e- db. It is easy to verify that the resulting g function satisfies the above assumptions (g(b) decreasing, twice continuously differentiable away from the origin, both and its weak derivative g' are bounded and absolutely integrable, and g' is non-increasing). Therefore, the resulting p has the all the properties outlined above.
I DEPTH SEPARATION PROOFS
Before giving the proofs in this section we introduce some additional notation. Let X denote the closed convex hull of the training points and X its open interior. Additionally, let Y = {(v, b) E gd-1 x R : o² x: > b for some x E X), and let Y denote its closure. Note that for any smooth function á with support contained in X, the Radon transform Re has support contained in Y. Finally, for any distribution lo and open set U. we let hlp denote its restriction U.
1.1 PROOF OF PROPOSITION
First, we show that the convergence of a sequence of functions fo to f in L' -norm over X implies that the sequence of distributions Afslx converges weakly to the distribution A For all test functions @ E S(X) we have
KAA o) = Iffe f. (0) s 1l (x))) r=
(83)
where we used Holder's inequality to achieve the final bound. Therefore, we have limp-soc (Ase, Φ) + (A,, ()), which proves the weak convergence.
Next, we show (R*) A fy converges weakly to (R*)-1Afly For all test functions & E (Y) we have
(84) Since & vanishes outside Y. by the support theorem (Helgason, 1999, Corollary 2.8) we ensured that R R has support contained in X, hence R R11 E S(X). The desired result now follows immediately by weak convergence of Afxlx Afxx.
((R*)-'Afi - (R')-'Af, = (Afr - Af.R-(p)
Since y has support contained in Y, this further implies that the distribution g- (R*)-1AA converges weakly to the distribution g (R*)= ff. Finally, since [lfklR = 1lg (R*) -1 ffx Il is bounded
See Proposition 3.5.1 in (Epstein, 2007).
25

---

Published as a conference paper at ICLR 2023
by assumption, this implies each g- (("))"AAf f is measure having finite total variation, hence their weak limit g -R**=1 Af is also a measure with finite total variation (ie., the weak limit of order-0 distributions that are bounded in TV-norm is also order-O distribution). Therefore, [|l = lig (R*)-IAflIrv is finite, as claimed.
1.2 PROOF OF PROPOSITION
To show the pyramid function p has infinite stability norm, we prove that g-(Rt)-1AAp distribution of order > 0, which implies
|P ||I|R
sup #ESH (S4-1 x3)
(g (R*)
Ap. Ap, () = +00.
(85)
First, observe that the Laplacian Ap is an order-O distribution whose support is contained in the unit ball... This implies the distribution (R*)-1As) Ap, is supported on a compact set K Radon domain. By our assumption Since g piecewise continuous and K compact, this implies there exists constants C1 02 > such that S g( b) S 02 for all (v, b) K. Therefore, we see that (R*) Ap is order-0 distribution if and only (R*) -1 Ap is an order-O distribution. However, by result in (Ongie et al., 2020), we know (R-)- Ap has order >> 0 (i.e., in the terminology of (Ongie et al., 2020), p has infinite R-norm). Additionally, we show this by direct calculation in App. M in the case of input dimension d = 2. Therefore, 9 (R*)-1 App must be distribution of order > 0 and hence IpllR. +oo as claimed.
the convex hull of the training points, g(v, b) > 0 for all (v. b) E K.
1.3 STABILITY OF THE TWO HIDDEN-LAYER IMPLEMENTATION OF p(æ)
Let us focus on the under-parametern predictor p(æ) that globally minimizes the loss. In this case, the set of all global minima corresponds to different implementations of p(x). Under this setting, we will prove that there exists a set of nonzero Lebesgue measure such that for any initialization inside this set, GD necessarily converges to p(æ).
setting, in which there exists a single optimal onput-output
To do so, we will first prove that for any minimum point 0- E I corresponding to implementation of p(æ), there exists a nonzero step size v) with which 0* is linearly stable. Furthermore, we will show that there exists a set The (0") embedded in subspace of dimension m MNull in which any initialization converges to 0" Here m is the number of parameters in our two hidden-layer network, and Mulll the number of zero eigenvalues of V C at 0*, Next, we will show that there connected set of minima O" around 0", such that the union UoER- The (0) has a nonzero Lebesgue measure.
Let us start with some minimum point 0", which corresponds to an implementation of p(x). GD's update rule is
i 4+1 = ())
(86)
Define the mapping
T(0) = 0 - (())
(87)
Then (86) can be written
0,11 T(O)).
(88)
This equation describes the full dynamics of GD using the nonlinear mapping T. Note that this representation, 0* is equilibrium point of T, i.e.. T(00) 0" We would like to show that it is possible to converge to 0" Assume there is a finite number of training samples n, and none of them coincide with the knots of p. Then T is differentiable in small neighborhood of 0" The Jacobian matrix of 7 is
i ---I-nv²'C0
(89)
00
and its eigenvalues are
A
T = x; -I-
=1 -1
(90)
00
26

---

Published as a conference paper at ICLR 2023
In this setting, the loss' Hessian at 0" has non-negative and bounded eigenvalues. Particularly, there exists" 10 a sufficiently small step size 1 that satisfies
0 < - 7,,, (V ((0))) 1,
(91)
for all Using the Center and Stable Manifold Theorem (Shub, 2013, Th. III.7), there exists a bounded set Too (0") such that
The (0") (0=) and
44 (0'): [T 0"1
S a 0 0*
(92)
for some 0 1. Here Tug. (0") is tangent the hyperplane that contains 0 and is spanned by the nonzero eigenvectors of V =(0"). Thus, assume that the initial point 0g E T then
(T(O,) 1''1 a' (TT(0)) 0* 3 0,
(93)
which shows that with any initialization in (0"), GD's 's iterations converge to e"
Next, note that there is a neighborhood of 0 within which the set of global minima of the loss form a smooth WNaat-dmmnnsion
manifold Let us denote this set of global minima around 0" by 0*, and 17 < 2/ maxøge- {Amax ( [)} and limit the set O* such that Vi 7 * 1/X, Then, for each minimum 0 E O according to the first part of the proof, there exists (m - mNul) -NNn) ddimensiooa The (0). Now, since each The (0) is contained in hyperplane that is orthogonal
the tangent of @" at 0, and the dimension of the tangent of @* at 0' is MNull we have that the dimension of the union of these sets, Us (0), is m. Thus, the set UsEe. The (0) is of nonzero Lebesgue measure within Th and for any initialization in Uoce- The GD converges p(ar).
J GENERAL LOSS FUNCTIONS
In this section, we discuss how our results can be extended to general loss function with a unique finite root. Assume some general loss function
1 e (f(x;),,;)
c(f)
(94)
where É(a, b) is twice differentiable w.r.t. a and is minimized when a = b, i.e.,
① E(a, b) = va= =
E' (a, b) a
(95)
da
Then, we can calculate the loss' gradient
VoL = - Σ
e' v)
Vof(x,)
(96)
j=l
and Hessian matrix
1 VOL
e" (f((z)))) Vol (z,) Vol (x;) + 1 n
e (((z)) (x,)
T1
j=1
Vol 8d (f(m))))) Vol (2,) Vol
(97)
n 1=1
Bo2 l(a, b). and in the last transition used f(z,) yj and ('(a, a) 0 for all E R. If (z))).. y;) = C > for all training points, then we can generalize our results by simply multiplying the RHS of (47) by C. If not, the analysis can still be used but need to add weigtning term to the stability norm which depends on the value of L" (f(æ;), Vj) for each data point.
where (a) A 8²
10 "Specifically, any 7 such that 17 S 2/ Amas: ( (VL(())) and : 17 # 1/A,, satisfies this condition. This set corresponds to multiplying the weights of corresponding neurons within different layers by positive factors whose product is 1.
27

---

Published as a conference paper at ICLR 2023
K PROOF OF LEMMA 2
In this section, our goal is to upper bound the top eigenvalue the flattest implementation of predictor function f. Here right singular vector of ②, then
notation and some derivations form App. F.1. Let q denote the top
mmax
1
|
71
2
2
1
Σ
n
4
Σ
(2) a
+
Σ
99 w
(2)
1i.t
71
1=1
27
Σ
[1) w! b(4))
9jo
21 +
Σ
2 3 11
+ Σ
2 j.c
i=1
1=1
Σ
w
(1)
σ
2 w
(1)
(())))
(98)
where in the inequality we used
EE=
Uz)
S n for all u E S and substituted
w
(1)
w
(99)
w (1)
w
Let O(f) be the set of all implementations corresponding to f. Since substituting
(1) + -
b
(2) - c.w,
(100)
does not affect the network's 's functionality f, we have
min Xmax (V)C) (E
< min (EE(()
1
n
+
Σ
k
(2(2))
(2))
d
Σ
+
22 Σ
1=1
[=]
==
27
Σ
w
(1)
[
w!!
B
2
d Σ
"
"
min n
w1
(2)
Σ
+
Σ
(ce())
1=1
j=1
27
+0p
wi
(1) 2
Σ
1io
x
T
w
(1)
((1)
(101)
28

---

Published as a conference paper at ICLR 2023
A necessary condition for optimality is that the derivative of the objective with respect to e; is equal to zero:
12))2
2
20
wi
3
+
Σ
/=1
2
2c; 3
200
Σ 9io [0) 20 !1 8,01
(1)
Σラ Qio (2. (1! - b(:l))
w;
(102)
(2)
(Σzz
م 2 + (Σn g;L1js)
w
V
It is easy to verify that these solutions for } are indeed global minima. Plugging this in, get
min (Ee(S)
Amax (Voc)
min (Ee(S)
1 n
n +
2
Σ
w
(1)
w
(2)
Σ
n
40
T
w
(1)
)l4)
j=1
2
2
Σ
Σ
(103)
Now, by Sauchy-Sehwarz inequality three times we get
Σ
9
2 T (1)
w
-
$ |191,
Σ
x
w
11)
i=1
j=1
n
Σ
||q11 Σ
x 2
/=1
Σ j=1
< ||ll2
Σ
(104)
Since || = 1. the right hand sides of these inequalities are independent of Thus, using these inequalities further upper bound the top eigenvalue we have
min Amax (V)L) eee())
min (Ee(f)
1
+
2
W
W
(2)
7
2
2
w
(1)
-
(0))
Σ
Σ x 19l
+
Σ Lii
π
1=1
1=1
++ min n oee(f)
w
(1)
(2)
2
x
W (1)
Σ
||ر5مll + 1)
Is..
(105)
To continue upper bounding the sharpness (Amax (V))) of the flattest implementation, we can consider some implementation of Specifically, since ʃ E Fs, it can be represented as
f(x) = Σ ai (v. x b.) + 3x h + c,
(106)
i=1
29

---

Published as a conference paper at ICLR 2023
where ||w/|| = for all i. ||h|| 1, a, + 0 for all i, and (v,, b.) + = by) for all + j (i.e., the knots of all ReLU units are distinct ), and where k' < k. Thus, we use the following implementation of f
m (1)
vi 0,00 bis
mi (2)
ais b(2== c BT,
(107)
for E [k]], where
n h Fjj
(108)
Additionally, if needed, we add two ReLU neurons to implement the linear component.
WK' 11+1 = h, 6/41 = -,, 4'4 1 = 3, w (1)
(2)
4'42 = --hh -h, b142 = T, w 2 = - -3.
(2)
(109)
Thus,
1+
min nt (E
2
w
(1)
(2)
Σ
1
02
1 w
((.))
Σ ||e,| 112 +1) Uus
j=1
$ 11
+
2
Σ
02
x 0i bi)
Σ j=1
a
(1 ||=,||| - 1) U3,
+ 72
181
Σ
02 [2. h - -)
Σ
(102211 + +1) T
/=1
1=1
+181
Σ
2 [-xe T h +T)
Σ
||0r|1 + 1) 11:
/=1
j=1
2 K'
S 1++ +
las
(x) 0i bi)
Σ
||æ, - =1) م4.
77
j=1
5=1
4|3]
n
+
Σ [a] h - 77²²
Σ j=1
||,, 11 +1).
(110)
n
where in the last inequality we used
max {o² (x.]h-1) ,), 0² (-x,n++ < (z] h=)²;
(111)
and ||= 2| +1 > O and thus removing the indicator term only increases the RHS of (110). Recall that
we denoted C. S { be the set of training points for which the ith neuron is active, and n, = [C,]. Then.
k''
lasl,
Σ 22 [x. U;
- bi)
Σ
||225 112 +
1) T
/=1
k'
|0:
n4
(x. vs - -b))
1
Σ
(112,)11 + +1)
n
n; jECt
ni
i=1
jEC:
K'
la IP
(x T v > b.)
v E (0 v, b,)²
x To, > b.] v E 1 + |
?
x T uu
> by]
(112)
i=1
If (v., b:) = (v,, by) then one of the neurons is redundant. If (w., (vv.b:) b.) = - (u,, by) then these units can be combined into an affine function, which we "absorb" into the term oz: h +
30

---

Published as a conference paper at ICLR 2023
Additionally,
4|3
[ h - 7)²
Σ
by / + 1) = 4131 / Var (x h)) 11+ E
x |?
(113)
j=1
j=1
Define
g(v = (x (xJ > b) B 1 - 6)2120 To b] 1 + ||
(114)
Then, we have
W
min Amax (V2c) & +1+2 (E3(1)
Σ
|alā(0,, 6,) + 4113 V Var (x: h 1 + ||==|
i=1
1 +1+2
[[(R)) (z. 0(v. ()ds
ga-1xR
+ 4 181 Var Var (z)Th))/
/ 11 E ||0||
Var Var (zo 1 2 [[ + 4131 /
1)) (1) E ||0||
(115)
Note that,
Var (a: hxTh)= S V Amax Ez)
(116)
where En is the covariance matrix of æ. Additionally,
[Vf(z)]
Σ a 0,00 3h
2 |0 |/|
a; 1 u. x 6=
1=1
|8|
Σ
lail ||
|8|
lail
1=1
|8| [
(R+) 11 Af ds(v)db
Jg
|8
||
(117)
Therefore, for any x E Rd
|8| < [f(x)] + +1fllm
(118)
Taking the tightest bound obtain
| < lffl + inf [[V](=)] #eRd
(119)
Overall, combining (115), (116), and (119) we obtain
min Xmaa (V2C) (Ee(()
<1+2 +4 (15fl1 + inf mERe
[1T)(())) V Ammox(Ia (Ex) v
1 -
(120)
31

---

Published as a conference paper at ICLR 2023
L PROOF OF PROPOSITION 3
Let f E W.t... (R4). First, we show that this implies both l and [I are finite. Since we assume d is odd (-A) (- (d is an integral power of the negative Laplacian applied to f. hence can be expanded as a linear combination of d d order partial derivatives, and so
I[-A) (d+ (d+1)/2 fll.. S add Σ 19l=d+1
|10 fllı, I aalll Callffllwsttttee))
(121)
where 0a is constant depending on d but independent of ʃ. Therefore, |[(-A)(+1)/2 f'l.... is finite. In particular, this shows ( (-A)(d+1)/² E L (BR), and so ((-A))(d4)/)) exists in classical sense. This implies have the formulas I|
fllı and llflRR yallå R( (d+11//2 fll1 where Va 2(2))=r (see (Ongie al., 2020, Prop. 1)).
Recall that v((x) = [1 + [a)) = cs + Cal|=r|| with ca = Sor 1 du and Ca = saa 1o ldv. Therefore we have that
(I ( () f11
(- f(x) w(x:)da
R (- //} (1, 6) (1+ (b)ds(v)db
1 XR
R {
(-
(0. (1 +1+ [bb))ds(b))db
34-1 XR
1 [8]
1 1x-
R.
(-))
(1 + g(v. ())ds
1 g(v. b)
2 C9
R {-
(1 - g(v, ())dd
3d--1xR
% Callfflle + (1.l).))
(122)
1+ where Cg = inf (v. I+)(w,b) is finite and non-zero because for all E gd-1 we have g(u,b) b) O(b)) as |6| + where the implied constant is independent of v. Therefore, we have shown |f and [ are finite as claimed.
Let a = (R*) AA = R function, which can be identified with finite signed measure. Since | bl9 aili is also finite, we see that a(v, b) == (1 ÿ(v ()))) b) is also an L function which can be identified with finite signed measure. By (Malliavin et al., 1995, Thm. 6.9), this implies there exists sequence of finite atomic measures {āx ), such that each â& consists of a sum of at most k Diracs, converging narrowly 13 to à with | S ||al|]. Define (1 b) = (v b) (1 + ((v.))), which is also an atomic measure. Then it is easy show ag a narrowly, well. By Lemma 5 of (Ongie et al., 2020), this implies there exists a sequence of single hidden-layer ReL.U networks fx E Fk converging f pointwise. Therefore, for all we have
Then | = ||a]] is finite, and so a is an L
[lfull + fa a = 1lfellR,144 = || || $ ||ā]]| = Ilfll = llfllR + llflRR
(123)
Combining this inequality with the bound on |l + |R. given above, we see that
[ + [[fa|R S || + [| S ca.o fllll
(124)
where Cd.g C is constant defined independently of f.
Finally, if K is any compact subset, the pointwise convergence of fú to on K can be upgraded to L -eooverrennce using Lebesgue's dominated convergence theorem: by the bounds on the Lipschitz constant of function given in terms of the R-norm in Proposition 8 of (Ongie et al., 2020), have ( (x)] < |[=]|[C + ll < B|| for some constants C, B 2 O, and since x + B B|[=|| is L'-integrable over any compact subset, the hypotheses of Lesbesgue's dominated convergence theorem hold.
Namelyy (ax.sp) (a, sp) for all continuous and bounded functions y gd-1 R R.
32

---

Published as a conference paper at ICLR 2023
M STABILITY NORM OF "PYRAMID" FUNCTION
Here, to provide a better understanding the depth separation result in Proposition 2, we show by direct calculation that the "pyramid" function in d 2 dimensions, given by
p(æ) p(21, x2) [1 |i [
(125)
fails have finite stability norm. particular, we explicitly compute (R*)-1Ap as tempered distribution and show it is not a finite measure (i.e., must be distribution of order > 0), which implies cannot have finite stability norm under the assumptions Proposition 2.
First, observe that Ap is linear combination of Diracs supported on finite line segments 4z defining the "edges" the pyramid:
Ap(x) = Ecudew
(126)
This means that if Φ is any Schwartz class test function, then
(Ap(a)) = Σ ck / (x))ds(x))
(127)
Note that Ap(x) is finite measure (i.e., distribution of order zero), since
(A)()))))) S Σ [callal] l1
(128)
where [2 is the length of the line segment l.. See Fig. 6 below for illustration.
07
45
45
04
05
-15
-05
05
1.5
(a) Pyramid function p
(b) Approximation of Ap
Figure Visualizations the pyramid function p and Laplacian - Ap.
Now we compute RAp, the Radon transform of Ap. which exists as tempered distribution. First, we show how to compute Ro; where & denotes general line segment. Let to is any Schwartz class test function defined in Radon domain, then
(Roe, :p) = ( ()z))s(())-
/./ / p(0, v xx)ddds())= p(0,00 wdds(
(129)
where the exchange of integrals is justified by Fubini's theorem since d, is finite measure.
33

---

Published as a conference paper at ICLR 2023
Suppose l is vertical line segment l ((0,t) :t E [c, d. Assuming is such that 02 # 0, then the inner integral above with place of Eg above simplifies as
J. p(o v z)ds(x))
((v,, u2t)dt
lould
so b)db
D2
Jezle
[ /2/ (12)) (b)p(1,))
(130)
02
In the event that 02 = 0 we have
J (v, ((),))))
/- o(v, ld- clo
-ld-cl0.(0)) = lddddd -d (n))
(131)
Therefore, we have shown
102 (12alc//21d d clo(b)
(b)
if 02 + 0,
Rōg(U, b)
(132)
if 0/2 = 0.
Now, consider one of line segments E coinciding with the edges of the pyramid. This can be parameterized as & {bb + top : E [ck, dx]] where us E St is unit vector, by E R is constant, and vł is orthogonal to us. Let 02 is the angle such that ((cos())),si((00)) then lx a rotation of the vertical line segment through the angle Os, and translation by bg Therefore, by properties of Radon transforms,
[cos(Ox),, sin(O))],
R5, (v(0),b) - (5 0x), by cos - 62)),
(133)
where we set v00) = [cos(0) (), sin( for all 0 E T). More concretely, we can express every slice RRe (v. as either a weighted indicator function when v * t which is non-zero when b such 6} intersects the line segment lg, or weighted Dirac when
that the line Lu,b {a 10k, i.e.,
[sin(0 0k)l if (o),, O & is singleton, |exl0(b bx) if parallel vk: else.
R5 (v(0),b)
(134)
For fixed v(f) E gd-1 the set of b E R for which L.()) O to 0 is always closed interval [ore, Bn]. Therefore, for 0 + Og can write Ro, (0,-) sin 0% 0x)1-11a oa (0)) for some ag (0) and 3x (0) that vary continuously with A.
Finally, by linearity, we obtain RAp- = Σ. c. Rótz See Figure for an approximate plot RAp.
Now we compute (R')-1Ap. Recall that (R*)-1 KR where K. = Ha, filtering step with H being the Hilbert transform applied separably in the b-variable (Helgason, 1999). For smooth function g.
1 Hg(bb) p.v.
g(b') db'
(135)
π
where p.v. indicates principle value integral. Therefore, for any + Or, we have
1
(((),) = sin(( 0+)l (Hoar -70 (H)).(0))
1
* sin (0 0x)l p.v. b (0
P.V. b (2)
(136)
and for 0 = 0g we have
14x p.v. (b by
1
(RRr, (v(02),6) = 14x1746 -x)
(137)
34

---

Published as a conference paper at ICLR 2023
1000
-15
1.5
200
0.0
40.5
200
offetb
05
0.5
460
-00
1.5
2
0.5
angle (n radians)
2,5
0.5
angle (n racians)
1.5
2.5
(a) Approximation of RApp
(b) Approximation of (R*) 1110)_1 Ap
Figure 7: Visualizations of RAp and (R*)-1Ap.
Finally, by linearity of the operator KRR, we have
(R*)~'Ap=KRRA
Σ CuRRRE
(138)
Thus.
[(R')~App(<(0))bb
Σ ck sin(0 -an(0-0441-- 0k)1-1 P.v. b-v.b- (x) -..V b-x(0)
-1
Σ cullxl0(0 02) p.v. (b bz.2
(139)
See Figure for an approximate plot of has singularities along 1-D manifold in Radon domain. This set corresponds to all lines in the primal domain passing through the corners of the pyramid.
As evidenced by the plot, this density
Finally, show that (R*)-IAp not finite measure (i.e., it is not order zero distribution). Intuitively, this is because the "density" a(v, 6) is not absolutely integrable, since every 1-D angular slice has singularities like 1/16]. Below we prove this more formally.
To prove a cannot be an order zero distribution, we construct a family of uniformly bounded test functions (pc)eso such that (a, sp,) 2 p(e)llpeeloo. where p(e) function such that p(e) - +o0 as + +
be a small fixed constant less than one. For every e < 7, consider the "rainbow-shaped subset of Radon domain Ω2 defined by the inequalities - < 0 < 7/2 and cos(f) - < < cos ( where. In primal domain, the set corresponds a collection lines that nearly intersect the corner point (),
Let y >>
Only three terms in the sum making up (R*)-LAA in (139) are dominant the region Ω, corre- sponding to the three line segments in the support of Ap that arise from the right-most corner of the pyramid. Elementary calculations show these three terms are specified by the parameters:
G = 2, 0 = =/22, 3 (0) = cos( ())
C2 V2, 0g r/4, 82( cos (
c3 /2, 03 //4, (0 = cos(())
(140)
Therefore, a well-approximated on S2 by
v2
2
2
1
ā
P.v.
(141)
sin T/4)|
si ग /
sin(0 π. 2/2)
cos(0) bb'
35

---

Published as a conference paper at ICLR 2023
where omit the terms p.v. 6. and p.v. P.V.(6-b,)= [2. since points in Ω2 far from their singularity set. In particular, can show terms are locally smooth and bounded). Let g(0) be the function of 0 in front of the principle value independent of
an order zero distribution when restricted to Q2 (i.e. all other
in ã. Note that g > B > 0 for all 0 E f2 where the constant
Let (0, b) be smooth function supported in S2, such that 0 < (0, b) I 1 and pe (8. b) on 7/2 and cos c = S cos e² Then for
the region defined by the inequalities -7/2 (9.(0.)) db is bounded below by I job log(s any fixed 0 E (-7/2, the integral p.v. cos(O)- Therefore, we have
((a, (())).) 2 1(ā, v.)1 - -lāā- aa,
7/2
2
(0o, (00, )):00 - C1lspello
7/2
2 (y ( lllleloo..
(142)
Since loo = 1 and B log(e-1) C + +00 as - 0+, this shows cannot be distribution of order zero, i.e., it cannot be identified with a finite measure.
N ADDITIONAL EXPERIMENTS
The experiments in Sec. 6 are designed to demonstrate Theorem a diverse range of step sizes (1O-+,0.1] Since flat minima of the loss landscape are concentrated near the origin in parameter space, and training with small step size near flat minima is inefficient, we used large initialization (about 10 times larger than standard methods). Here repeat the MNIST experiment using various initialization scales on a higher range of step sizes (110 3 0.2).. Figure 8 presents the sharpness curves for the different scales. For large initialization, 10 and x15, get the same behavior as depicted in Sec. 6. For small initialization, x and x5, the sharpness of the obtained solutions is fixed for small learning rates up to a critical step size 17" At this threshold, the sharpness equals 2/1²". and any increment in the step size makes the minimum unstable. This pushes SGD to flatter minima for larger step sizes, ones that satisfy the stability criterion. Here is important note that for standard initialization, shown in Fig. 8(a), the threshold is well before the standard step size of 1 = 0.1. Namely, this phenomenon happens using standard initialization and standard learning rate.
36

---

Published as a conference paper at ICLR 2023
Sharpnoes
Sharpnees
10a
10
10
10
(a) Initialization factor (standard)
(b) Initialization factor 5
R/1 Aass
R/4
Ave
Auea Alll
Sharpoess
Sharpness
10
10º
10
10
10
10-1
10-2
10-1
(c) Initialization factor
(d) Initialization factor
Figure 8: Sharpness vs. step size for different initialization scales. We trained a single hidden-layer ReLU network for binary classification on two classes from MNIST using SGD (see Sec. for details). Specifically, we initialized the network using different scales, and for each scale trained the network using multiple step sizes. We see that as r) increases, the minima get flatter in parameter space (yellow curve), which translates smoother predictors in function space (purple curve).
37


Reviews:
[{"review_id": "wtYDZUrepOm", "reviewer": "Reviewer_7NLR", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper studies the implicit bias of minima stability for two-layer neural network with multiple inputs. Based on previous works on the same topic for two-layer neural networks with single input, the authors extended the analysis to include multi-dimensional input, and showed that stability controls a weighted norm of the Laplacian of the function represented by the neural network. Discussions are made on this regularization effect in both the primal and Radon spaces. Then, the approximation capability of neural networks under stability is studied. For two-layer neural networks, it is shown that functions in a weighted Sobolev space and be approximated by networks with increasing width at stable minima. A depth separation result is proven showing that three-layer neural networks (two hidden-layer NNs) can represent a function that cannot be approximated by two-layer neural networks under any stability requirement. ", "summary_of_the_review": "This paper studies the implicit bias of minima stability and the approximation of neural networks at stable minima. The theoretical results are limited (as detailed in the \"Strength and Weaknesses\" part), but still interesting because (1) it accurately describes an implicit regularization effect of GD, and (2) it discusses the approximation of neural network under a more realistic setting. ", "strength_and_weaknesses": "Strength: \n\n1. The paper extends previous results on the implicit bias of minima stability for single-input two-layer neural networks to multi-input networks, making the results slightly more general.\n2. The approximation capability of neural networks is studied under stability requirement. This is a more realistic setting to study the approximation properties of neural networks. \n\nWeaknesses:\n\n1. In the analysis of stability, only the learning rate is considered, while the batch size is ignored. Hence the stability criterion is tight only for GD. For SGD, it is just a necessary condition. The difference between GD and SGD cannot be seen in the analysis.\n2. The approximation results in Section 4 and 5 are training data dependent, while in practice the training data are sampled randomly. Is it possible to develop similar results that hold with high probability over the sampling of training data? ", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "Not applicable"}, "extra_text": {}}, {"review_id": "gYxC9W9hIL", "reviewer": "Reviewer_iAv3", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "The paper gives a necessary condition on a function $f$ for the existence of a local minimum (of a shallow ReLU network which represents $f$) that is stable with a fixed learning rate $\\eta$. The stability condition on $f$ requires a stability norm of $f$ to be bounded. This norm can be interpreted as a form of weighted L1 norm of the Laplacian of $f$.\n\nThe authors then give a depth separation result: an example function $f$ which has infinite stability norm (and hance cannot be approximated with a fixed learning rate in shallow ReLU networks), but such that there is a deeper network representing $f$ in stable manner. \n\nFinally the authors show (under some conditions) that any function $f$ in a weighted Sobolev space can be approximated by networks of increasing width $k$ in a stable manner. ", "summary_of_the_review": "The paper gives a rather complete answer to the question of minima stability in shallow ReLU networks, using a combination of techniques from two previous papers.", "strength_and_weaknesses": "I think that the question of minima stability is very interesting and crucial to understand what type of functions are learned by DNNs. The analysis proposed here seems quite complete for shallow ReLU networks, though it requires a number of technical assumptions and it gives only a necessary condition for stability. \n\nThe depth separability example function is both surprisingly simple and interesting. Tt suggests that stability conditions in deeper networks might be very different from those of shallow networks.", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "TBxPoTdm5", "reviewer": "Reviewer_EjYW", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper studies the properties of solutions that are stable for SGD, i.e., the ones that satisfy $\\lambda_1(H)\\leq 2$/step size. Specifically, the authors define a stability norm and prove that the stability norm of a stable solution must be no greater than $1/\\eta-1/2$. Then, the authors study the property of stability norm in both the primal space and dual space. In primal space, the stability norm can be explained as a certain weighted Laplacian norm, suggesting that training with a larger step size tends to find smoother models. The explanation in dual space is also discussed but no clear implication is made. Lastly, it is shown that this stability class can approximate functions with certain high-order smoothness. \n\nThe authors also come up with a target function, showing that this function cannot be approximated by stable two-layer neural networks but can be exactly represented by stable three-layer neural networks. \n", "summary_of_the_review": "This paper provides a new analysis of how the dynamical stability of SGD can implicitly regularize the \"complexity\" of the learned model. The result is well supported both theoretically and empirically. I generally feel this is a good work. ", "strength_and_weaknesses": "### Strong point\n\nThe paper is overall well-written and I can follow most results smoothly. \nThe theoretical result suggests that the stability of SGD implicitly regularizes the smoothness of learned model whereby SGD with large learning rate  generalizes well. Preivous similar works all are limited to linear networks and univariate two-layer neural networks. This work considers the standard two-layer neural networks, which is much more interesting. Given recent progresses in understanding the interplay between dynamical stability and implicit bias of SGD, I think this work is of great importance. \n\n### Weak point\n\n1. (Ma et al., 2021) also established the connection between stability and the input smoothness of implemented neural networks, and their analysis is also applicable to deep nets. I think the authors should provide a careful discussion of the similarity and differences. \n\n2. In Figure 4, what is the meaning of $\\min. \\lambda_{\\max}$?\n\n3. I find the paragraph above Lemma 2 (the one used to motivate Lemma 2) quite confusing and the writing there can be improved. \n\n4. In Figure 5, it is shown that the accuracy of MNIST has only around 55% for small LR but nearly 100% for large LR. This dramatical difference between small LR and large LR is very unexpected since this task is extremely simple, for which a simple logistic regression might have near 90% accuracy.  I do not understand why neural networks have only 55% accuracy even if it is trained with small LR. Is it because a large initialization is used? ", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "Not applicable"}, "extra_text": {}}, {"review_id": "OxqiNL0aOd", "reviewer": "Reviewer_KQ4p", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper generalizes the results of Mulayoff et al. 21 on the implicit bias of minima stability to functions that have multivariate output dimension. The paper considers training ReLU 2-layer networks with squared loss. The main contributions are:\n\n* [minima stability implies smoothness] Show that linearly-stable minimizers of the loss correspond to smooth functions (where a certain weighted integral of the laplacian is not too large)\n\n* [2-layer vs. 3-layer separation] Show that there are functions which 2-layer networks cannot represent at stable minima; but on the other hand, 3-layer networks can.\n\n* [universal representation of smooth functions] Show that sufficiently smooth functions (in a Sobolev sense) can be approximated by stable minima.\n", "summary_of_the_review": "Understanding implicit bias of neural networks is an important topic, and this paper presents significant original results. I recommend acceptance.", "strength_and_weaknesses": "### Strengths\n* The proofs are clearly written and easy to follow. (I appreciated the Section C providing background on the Radon transform.)\n* The analysis is novel, and significantly extends the previously known case of functions f : R \\to R to functions f : R^d \\to R\n* The results are of interest to the neural networks community, since they give new guidance on how the step size affects the implicit bias of training.\n\n### Weaknesses\n* In Proposition 2 statement, a distribution on x should be specified?\n\nTypos: \"globaly\", \"can be quiet high\", \"this method yield\", \"in the space even\"", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works.", "empirical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works."}, "extra_text": {}}]

Author rebuttal:

Decision: Accept: poster