Title: A MECHANISM OF IMPLICIT REGULARIZATION IN DEEP LEARNING

Abstract: Despite a lot of theoretical efforts, very little is known about mechanisms of implicit regularization by which the low complexity contributes to generalization in deep learning. In particular, causality between the generalization performance, implicit regularization and nonlinearity of activation functions is one of the basic mysteries of deep neural networks (DNNs). In this work, we introduce a novel technique for DNNs called "random walk analysis" and reveal a mechanism of the implicit regularization caused by nonlinearity of ReLU activation. Surprisingly, our theoretical results suggest that the learned DNNs interpolate almost linearly between data points, which leads to the low complexity solutions in the over-parameterized regime. As a result, we prove that stochastic gradient descent can learn a class of continuously differentiable functions with generalization bounds of the order of O(n -2 ) (n: the number of samples). Furthermore, our analysis is independent of the kernel methods, including neural tangent kernels.

1. INTRODUCTION:

Deep Neural Networks (DNNs) have demonstrated dominating performance in numerous machine learning tasks, and it shows great generalization performance in the over-parameterized regime. Theoretically, mechanisms of implicit regularization, which is considered as an important factor of such generalization performance in over-parameterized DNNs, still remain unknown (Zhang et al., 2016; Neyshabur, 2017) . In the over-parameterized regime, recent studies report that generalization bounds for DNNs can be obtained by replacing the neural network with its linear approximation model with respect to weight parameters at initialization. Most of these studies rely on the connection between deep learning and neural tangent kernels (NTKs) (Daniely et al., 2016; Jacot et al., 2018; Arora et al., 2019b) , which characterizes the dynamics of network outputs throughout gradient descent training in the infinite width limit. However, a source of implicit regularization in over-parameterized DNNs has not been identified. Recent empirical and theoretical results indicate that generalization performance and implicit regularization of over-parameterized DNNs cannot be captured by NTK analysis (Wei et al., 2018; Allen-Zhu & Li, 2019; Woodworth et al., 2019; Geiger et al., 2019) . Understanding how implicit regularization properly controls the superfluous expressive power of over-parameterized DNNs gives us new insights into the theoretical analysis in deep learning. This leads to the first question: Question 1. What kind of low complexity is caused by implicit regularization in deep learning? Linear networks without activation functions are important subject, and there are a number of theoretical works on the implicit regularization in over-parameterized neural networks mainly focusing on linear models (Ji & Telgarsky, 2018; Gidel et al., 2019; Arora et al., 2019a) . In contrast, whole properties of over-parameterized DNNs that may result from nonlinearity of activation functions cannot be captured by the approximated linear models, specifically, the kernel regression predictor using the NTK. However, some mechanisms of the implicit regularization can depend on nonlinearity. This leads to the next question:

Question 2. Can we identify a mechanism of implicit regularization that depends on nonlinearity of activation functions?

Now for optimization algorithms, the training dynamics of full batch gradient descent (GD) is better understood although GD is too expensive for most applications and one often uses stochastic gra- Figure 1 : Random walks having a step size distributed according to N (0, 2/m) (m = 10 3 ) after α steps. The smaller the number of steps α is, the more straight. For visibility, we transformed the random walk sequence {y i } α i=0 into {z i } α i=0 defined by z i := (y iy 0 ) -(i/α)(y αy 0 ).

dient decent (SGD) instead. In most cases, the authors used GD to derive their results by the NTK analysis. Recent work showed that without any structural assumptions about the data distribution, two-layer or three-layer over-parameterized networks trained by SGD can learn C ∞ -class functions (Allen-Zhu et al., 2018a; Arora et al., 2019c) . However, it is not clear the relation between generalization and implicit regularization for DNNs optimized by SGD. Towards this end, the following question is also unsolved:

Question 3. Is it possible to obtain provable generalization bounds based on the implicit regularization for DNNs optimized by SGD?

To answer these questions, we introduce a novel analysis for DNNs and characterize a mechanism of implicit regularization that caused by nonlinearity of ReLU activation. Our results indicate that the DNNs (trained by SGD) interpolate almost linearly between data points, which leads to the low complexity solutions in the over-parameterized regime. Accordingly, we prove that SGD with random initialization can learn a class of continuously differentiable functions with generalization error bounds of the order of O(n -2 ) (n: the number of samples), which is independent of the NTK analysis.

In order to introduce our analysis for implicit regularization, let us focus on the DNN output on a one-dimensional linear path between training points. We define it as x(s) := (1 -s)x (p) + sx (q) (s ∈ [0, 1]), where x (p) and x (q) are training points. The corresponding function of the DNN with ReLU (ReLU DNN) is continuous piecewise linear on the linear path x(s). In the hidden layer, the corresponding function g l (x) of each unit in the l-th layer is also continuous piecewise linear on the linear path x(s), that is, a composite function (g l • x)(s) is continuous piecewise linear. Since DNN nonlinearity is linked to breakpoints, which is caused by ReLU activation, the set of the breakpoints plays a key role in the network behavior. These kind of breakpoints are also known as kinks or knots (Steinwart, 2019) . Since (g l • x)(s) is continuous piecewise linear, the gradient of (g l • x)(s) has gaps at the breakpoints, which we call gradient gaps. We focus on gradient gaps with respect to the parameter s.1 Our key finding is that the gradient gaps are a constant multiple of independent Gaussian random variables according to N (0, 2/m), where m is the number of units in the hidden layer. Since the gradient of the unit (g l • x)(s) is the sum of the gradient gaps, it is a "Gaussian random walk" 2 (see Step 2 in the proof of Lemma 3). The relation between the step size variance and the number of steps determines the Gaussian random walk behavior. In the corresponding function of the DNN, the step size variance is the variance of weights, and the number of steps is the number of breakpoints. In a regime where ( the number of steps × the step size variance) ≤ O(1), the Gaussian random walk strolls little from the origin with high probability (see Figure 1 ). Our results show that the ReLU DNN is in the same regime and "simplicity" of the gradient depends crucially on the number of breakpoints. Hanin & Rolnick (2019) proves that the average number of breakpoints is linear in the number of hidden units at initialization, and we also give a proof of an upper bound on the number of breakpoints even after training (see §3.1).

In this work, we prove a priori generalization estimates, which is independent of the posterior data distribution, by analyzing the behavior of DNNs on the linear path x(s). In other words, we show that in Theorem 2, the difference between the network output function and its linear interpolation is evaluated only from the amount of the weight change. We note that in Theorem 2, we dose not use the properties of the trained neural networks. Our key analysis is that the gradient of unit (g l • x)(s) is a "Gaussian random walk" and approximately equal to a straight line on the linear path x(s), which we call random walk analysis. Our technique is based on a priori estimates that variation of the gradients between data points is extremely small and depends essentially on the value of the number of breakpoints times the variance of weights.

Difference from NTK. Interestingly, our analysis idea is different from other previous work based on NTK. Our findings are some novel aspects as follows:

• The NTK is defined using the gradient of the DNN output with respect to weight parameter space. In contrast, the linear approximation (Lemma 3 in this paper) is defined using the gradient of the DNN output with respect to input parameter space. In other words, the variables to be differentiated are different. • The random walk analysis indicates that over-parameterized ReLU DNNs interpolate almost linearly between the data points. For ReLU activation, since the NTK kernel mapping is not Lipschitz but 1/2-Hölder, it is difficult to obtain such a result in the NTK analysis without a tradeoff between smoothness and approximation (Bietti & Mairal, 2019) .

Our Contributions. In this work, we consider an L-layer over-parameterized ReLU neural network with l 2 regression task using SGD from random initialization. We show the constitutive relation between implicit regularization and generalization, which enables us to provide new insights on the role of implicit regularization in deep learning. Our main contributions are as follows:

• Our random walk analysis provides a priori estimates of low complexity in overparameterized deep neural networks, and directly indicates that the unit output between data points is properly controlled by weight initialization and SGD to keep connecting the data points almost straight, which is one of the underlying mechanisms of implicit regularization in the over-parameterized regime. • Our result suggests that implicit regularization is attributed to the nonlinearity of ReLU DNN, which is indicated by the fact that variance of weight and the number of breakpoints determine the Gaussian random walk behavior. • We also prove that in one-dimensional input case, SGD with random initialization can learn C 1 -class functions. Our generalization estimates are based on the implicit regularization.

2. PRELIMINARIES AND NOTATION:

Notation. For n ∈ N, we let [n] = {1, 2, . . . , n}. We use N (µ, σ 2 ) to denote the Gaussian distribution of mean µ and variance σ 2 . We use v 2 to denote the Euclidean norm of a vector v, use v F to denote the Frobenius norm of a vector v. When z is a sub-Gaussian random variable, we let z ψ2 = inf{t > 0 | E[exp(z 2 /t 2 )] ≤ 2} to denote sub-Gaussian norm (Vershynin, 2018) . For a vector v, we denote by [v] i or v i the i-th element of v. For a matrix M , we denote by [M ] i,j or M i,j the entry in the i-th row and j-th column of M , and we denote by [M ] i the i-th row vector of M . We denote by 1 {E} the indicator function for the event E. ReLU activation is given by φ(x) = max{0, x}, and for a vector a ∈ R n , we define φ(a) = φ(a 1 ), φ(a 2 ), . . . , φ(a n ) .

Network structure. In this work, an L-layer fully-connected feed-forward ReLU neural network with m units in each hidden layer is given by f : R d → R c . For input x ∈ R d , unit output g l (x) of each layer l ∈ [L] and network output f (x) are given by the following functions:

EQUATION

where A ∈ R m×d , W l ∈ R m×m (l ∈ L), and B ∈ R c×m are weight matrices. We assume the following Gaussian initialization: A i,j ∼ N (0, 2/m), W l i,j ∼ N (0, 2/m), and B i,j ∼ N (0, 2/c). For input x and weight matrices -→ W := (W 1 , W 2 , . . . , W L ), the network output f (x) is also denoted by f ( -→ W, x). In this work, we only update weights in -→ W and leave A and B at the random initialization. For input x ∈ R d and l ∈ [L] , we denote by G l (x) a diagonal matrix, which represents the activation pattern of the l-th layer, which we call an indicator matrix. More precisely, we define the i-th diagonal element as G l (x) i,i := 1 {g l i (x)≥0} . Since the ReLU activation is positive homogeneous, we obtain the following equality:

g l (x) = W l G l-1 (x)g l-1 (x).

Dataset and loss function. The data are generated from an unknown distribution D over (x, y) ∈ R d × R c , where x is the input data point and y is the label associated with this data point. We assume without loss of generality that for each input x = (x 1 , . . . , x d ), using additional coordinates (x d+1 , x d+2 ), the replacement input x := (x 1 , . . . , x d , x d+1 , x d+2 ) is normalized so that x 2 = 1 and its last coordinate x d+2 = 1/ √ 2.3 we also use x to denote the replacement input x ∈ R d+2 . The training data Z := {(x (1) , y (1) ), . . . , (x (n) , y (n) )} is given as n i.i.d. samples from D. We define the minimum distance of the training data:

δ := min{ x (i) -x (j) 2 : ∀i, j ∈ [n], i = j} > 0.

For the l 2 regression loss ( ŷ, y) := 1 2 ŷy 2 2 and a subset of training data Z (τ ) ⊂ Z, we define our regression objective as follows:

L Z (τ ) ( -→ W) := E (x,y)∼Z (τ ) [ (f ( -→ W, x), y)].

Stochastic gradient descent with Gaussian initialization. We use mini-batch SGD to train the network with a constant learning rate η > 0, a batch size b and iteration number T . Let -→ W (0) := (W 1,(0) , W 2,(0) , . . . , W L,(0) ), A, B be weight matrices generated from the above Gaussian initialization. Suppose we start at

EQUATION

where

Z (t) ⊂ Z is a mini-batch of size b. For input x, each layer l ∈ [L] and each step t ∈ [T ],

we denote by g l,(t) (x) the unit output, f (t) (x) the network output, W l,(t) the weight matrix, and G l,(t) (x) the indicator matrix.

In the above setting, recent paper (Allen-Zhu et al., 2018b) shows that SGD can allow an overparameterized multi-layer network to attain arbitrarily low training error as follows:

Theorem 1 (Convergence of SGD (Allen-Zhu et al., 2018b) ). For any ε ∈

(0, 1], δ ∈ (0, O(1/L)] and b ∈ [n], let m ≥ Ω poly(n,L,δ -1 )d b , η := Θ bδd poly(n,L)m log 2 m , T = Θ poly(n,L) log 2 m bδ 2 log n log m ε ,

and:

-→ W (0) ,A,B are at random initialization. Then, it satisfies with probability at least 1 -e -Ω(log 2 m) over randomness of Z (1) , . . .

Z (T ) : L Z ( -→ W) ≤ ε, W l,(t) -W l,(0) F ≤ O λ log m √ m , [W l,(t) ] i -[W l,(0) ] i 2 ≤ O λ log m m (∀t ∈ [T ]), where λ := n 3.5 √ c δ √ b .

Our results on generalization also crucially depend on this analysis.

3. RANDOM WALK ANALYSIS:

In this work, we consider the behavior of the unit output in each layer on a one-dimensional linear path between two data points, which we denote by x(s) (s ∈ [0, 1]). In this section, for notational simplicity, we drop the superscript with respect to (t), which is the t-th iteration.

Definition 3.1 (One-dimensional Linear Path). For each pair of data points x (p) and x (q) (p = q), We define x(s) := (1 -s) x (p) + s x (q) , s ∈ [0, 1], and denote by v := x (q) -x (p) the direction vector.

Note that the unit output g l i (x(s)) and the network output f i (x(s)) are continuous piecewise linear functions on s ∈ [0, 1], and the network output can be expressed as follows:

f (x(s)) = BG L (x(s))W L G L-1 (x(s))W L-1 • • • G 1 (x(s))W 1 G 0 (x(s))Ax(s).

(3)

Note that g l i (x(s)) and f i (x(s)) have linear approximations at the point s = 0 as follows4 :

EQUATION

We prove that these piecewise linear functions are almost straight, which contributes to the generalization of NN. In other words, the unit output g l i (x(s)) can be approximated by a linear function with high accuracy. Moreover, the small difference between g l i (x(s)) and gl i (s) indicates the low complexity of the unit output. Theorem 2 (A Priori Estimates for Implicit Regularization). Under the same setting as Theorem 1, with probability at least 1 -e -Ω(log 2 m) , for every

EQUATION

)

EQUATION

3.1. INTUITION BEHIND IMPLICIT REGULARIZATION:

To prove Theorem 2, we introduce our key analysis that the gradient of g l,(t) (x(s)) is "Gaussian random walk" and nearly equal to a straight line on the linear path x(s), which we call random walk analysis. For simplicity, we explain the outline of the proof for the initial state of the network (i.e. t = 0). All proofs are given in the Supplementary Material (including 0 ≤ t ≤ T ).

To state our key Lemma, we define the following function: gl (s) := W l φ(g l-1 (s)), which is an analogue of the unit output g l (x(s)) = W l φ(g l-1 (x(s))). In other words, gl (s) is the one in which g l-1 (x(s)) is replaced by gl-1 (s). Next, Lemma 3 shows that each gl i can be well approximated by a linear function gl i (s) with high probability. Lemma 3 (Linear Approximation Analogue). With probability at least 1 -e -Ω(log 2 m) , for every

EQUATION

EQUATION

)

The purpose of this subsection is to give an intuitive explanation of the proof of Lemma 3, which can be divided into two steps. The first step gives estimates of the number of breakpoints. We show that the number of breakpoints of the piecewise linear function gl i (s) is less than or equal to the number of units in the layer (i.e. m).

The second step gives estimates of gradient gaps of gl i (s). We show that the gradient gaps are independent Gaussian random variables and the gradient of gl (s) is the sum of the gradient gaps, which indicates that the gradient is a "Gaussian random walk". Note that in this setting, the number of breakpoints is equal to the total number of steps of the random walk.

Step 1. By the definition of gl (s), we may write

EQUATION

Note that the input to the above φ, that is g l-1 j (s), is a linear function on s ∈ [0, 1]. For each j ∈ [m], the linear equation g l-1 j (s) = 0 has at most one solution (s = s * ). If the solution s = s * satisfies 0 < s * < 1, then g l i (s) has a breakpoint at s = s * . This shows that the number of breakpoints of g l i (s) is bounded by the number of the linear equations {g l-1

j (s) = 0} j∈[m] , which is equal to m. We denote by 0 < s 1 < s 2 < • • • < s α < 1 breakpoints of g l i (s)

, where α is the number of breakpoints. For β ∈ [α], we set an open interval of breakpoints I β := (s β , s β+1 ).

Step 2. To give an intuition that the gradient of gl i (s) -g l i (s) is "Gaussian random walk", we fix some notation. For β ∈ [α], we define the gradient ∇ β :

EQUATION

)

Note that g l i (s) is a linear function, and there is no breakpoints for gl i (s) on I β . Note also that since gl i (s) -g l i (s) is linear on I β , the gradient is constant on I β . For β ∈ [α], we define a gradient gap:

x β := ∇ β -∇ β-1

, and we have the following estimate with probability at least 1 -e -Ω(log 2 m) :

EQUATION

Note that ω β is an element of the weight matrix W l . Thus, {x β } are independent Gaussian random variables, and ∇ β is the sum of x β , that is

EQUATION

)

This shows that the gradient ∇ β is a "Gaussian random walk" and depends essentially on the number of breakpoints α and the variance of weights 2/m. Using the randomness of W l and general Hoeffding's inequality (Vershynin, 2018) , we have

EQUATION

)

According to Step 1, the number of breakpoints α is less than or equal to m. Thus, ∇ β is bounded by log 2 m √ m v 2 with probability at least 1 -e -Ω(log 2 m) . Therefore, with the same probability, we

have d ds gl i (s) -gl i (s) ≤ log m √ m v 2 , (∀s ∈ [0, 1]

). Thus, integrating this inequality from s = 0 to s = 1 shows the estimates of (7) in Lemma 3. Note that the proof of (8) in Lemma 3 is idential to the above proof, except the fact that each entry of B follows from N (0, 2/c) insted of N (0, 2/m).

Note that as illustrated in Figure 1 , for sufficiently large m, the gradient ∇ β can be made arbitrarily small, which means that the networks interpolate almost linearly between the data points.

3.2. PROOF THEOREM 2:

The purpose of this subsection is to give the proof of Theorem 2. We proceed by induction on the layer l. Note that since g 0 (x(s)) = Ax(s) is linear and g0 (s) = g 0 (x(s)), g1 (s) is identically equal to g 1 (x(s)) = W 1 φ(Ax(s)). Thus, the case l = 1 is true. Assume the theorem holds for any layer l = k, and let us prove it for l = k + 1. Using the triangle inequality, we have

EQUATION

By Proposition 11.3 of (Allen-Zhu et al., 2018b) , for any y, ỹ ∈ R m , there exists a diagonal matrix D ∈ R m×m such that |D j,j | ≤ 1 and φ(y) -φ( ỹ) = D(yỹ). Thus, for the first term on the RHS, we may write

EQUATION

In the last inequality, we use general Hoeffding's inequality (Vershynin, 2018) , and this inequality holds with probability at least 1 -e -Ω(log 2 m) . We can now apply our induction hypothesis to the above estimate. Thus, applying (7) in Lemma 3 to the second term on the RHS of ( 14) shows the estimates of (5) in Theorem 2. Next, Using the triangle inequality for f (t) j , we have

EQUATION

In a similar way, using the randomness of B (recall each entry of B follows from N (0, 2/c)), with probability at least 1 -e -Ω(log 2 m) , we have

EQUATION

In the last inequality, we use (7) in Lemma 3. Thus, applying (8) in Lemma 3 to the second term on the RHS of ( 18) shows the estimates of ( 6) in Theorem 2. We complete the proof.

4. GENERALIZATION:

In this section, considering the mechanism of implicit regularization revealed by random walk analysis, we provide a priori estimates for the generalization performance of over-parameterized deep neural networks, in an l 2 regression task on a one-dimensional input-space. Since the network is over-parameterized, the expressive power of the network is rich enough to considerably overfit the data. Nevertheless, it is known empirically that properly initialized over-parameterized deep neural networks can achieve the good generalization performance while fitting all training data. We estimate the low complexity of over-parameterized deep networks and show that the error between a trained over-parameterized neural network and the target function can be uniformly bounded by an arbitrarily small positive number.

Setting. We propose a new type of generalization bounds (Theorem 4), and prove this theorem in models with one dimensional input: d = 1. In the following, we restrict ourselves to a onedimensional regression task on an interval

[0, ν] ⊂ R. The training dataset { x (i) , y (i) } i∈[n]

is given as n i.i.d. samples from some unknown distribution D. We assume that the corresponding target function for the regression task f

* : [0, ν] → R c is a C 1 -class function. This implies that y (i) = f * (x (i)

). Without loss of generality, we may assume that the input data {x (i) } i∈[n] follow a uniform distribution on [0, ν], and after relabeling, we may assume that the data points {x (1) , . . . , x (n) } are ordered by index:

0 < x (1) < • • • < x (n) < ν.

We define (x (0) , y (0) ) := (0, f * (x (0) )), (x (n+1) , y (n+1) ) := (ν, f * (ν)) and T := sup k∈[n+1] x (k) -x (k-1) . We denote by f : [0, ν] → R the linear interpolation of the data points

x (i) , f * (x (i) ) i∈[n] .

A priori generalization bounds. Now, we introduce a novel approach for a priori generalization bounds, which is based on random walk analysis in the over-parameterized regime. This can be interpreted as the significant expressive power of over-parameterized neural networks is controlled by implicit regularization.

Theorem 4. Suppose f * (x) is a C 1 -class function on [0, ν]. Under the same setting as Theorem 1, for δ ∈ (0, 1/2], then with probability at least 1 -(δ + e -Ω(log 2 m) ), we have

EQUATION

)

Proof sketch of generalization. The purpose of this paragraph is to give an intuitive explanation of the proof of Theorem 4 for one-dimensional output c = 1. All proofs are given in the Supplementary Material.

To estimate E (x,y)∼D f (T ) (x), y , we evaluate f (T ) (x) -f * (x) 2 , which we may write as

EQUATION

The first term of the RHS of ( 21) represents an error between the trained network output f (T ) (x) and the piecewise linear function f (x). We can use the results of random walk analysis to evaluate this critical term as follows. Theorem 2 provides that the linear approximation error of the network output, which we denote by ε f (T ) , is small in each interval [x (k-1) , x (k) ], and hence in each interval [x (k-1) , x (k) ], the difference between f (T ) (x) and f (x) falls within the error ε f (T ) . This suggests that the network output between adjacent training points is properly controlled by weight initialization and SGD to keep connecting the points almost straight, which results in low complexity of over-parameterized neural networks. This statement can be extended to the interval [0, ν]. Fix ε > 0, and suppose m ≥ O T 6 /ε 3 . Then with high probability at least 1 -e -Ω(log 2 m) , we have

sup x∈[0,ν] f (T ) (x) -f (x) ≤ ε 2 (22)

The second term of RHS of ( 21) is the error of the piecewise linear approximation of f * (x) by f (x) on the interval [0, ν]. Note that the error can be reduced by increasing n. Hence, for any fixed ε > 0, there exists δ > 0 such that if n ≥ O ν √ εδ , with probability at least 1 -δ, we have

sup x∈[0,ν] | f (x) -f * (x)| ≤ ε

2 . Thus, from ( 21), for m and n sufficiently large, the error between the trained network output and the target function is uniformly bounded by ε on the interval [0, ν] with probability at least 1 -e -Ω(log 2 m) + δ . This gives a priori estimates for the generalization performance of over-parameterized neural networks.

5. RELATED WORK:

Implicit regularization in neural networks has recently become an active area of research in machine learning. A number of works have focused on the behavior of gradient descent on overparameterized neural networks (Neyshabur et al., 2014; Lin et al., 2016; Zhang et al., 2016; Soudry et al., 2018; Rahaman et al., 2018) . In order to get a handle on implicit regularization in deep neural networks, the majority of theoretical attention has been devoted to linear neural networks (Ji & Telgarsky, 2018; Gidel et al., 2019; Arora et al., 2019a) .

Many works try to explain generalization of over-parameterized neural networks. Recent works have shown that on sufficiently over-parameterized neural networks, the learning dynamics of gradient descent are governed by the NTK (Daniely et al., 2016; Jacot et al., 2018; Du et al., 2018; Allen-Zhu et al., 2018b; Lee et al., 2019; Arora et al., 2019b) . In these settings, the implicit regularization and the generalization error of the resulting network can be analyzed via NTK and the reproducing kernel Hilbert space (RKHS) (Bietti & Mairal, 2019; Nakkiran et al., 2019) . To extend it to SGD, Hayou et al. (2019) introduce a stochastic differential equation dependent on the NTK. Closely related work Cao & Gu (2019) showed that the expected 0-1 loss of a wide enough ReLU network trained with SGD and random initialization can be bounded by the training loss of a random feature model induced by the network gradient at initialization. This contrasts with other recent results that show a provable separation between the generalization error obtained by neural networks and kernel methods (Wei et al., 2018; Allen-Zhu et al., 2018a; Allen-Zhu & Li, 2019) . Several papers suggest that training deep models with gradient descent can behave differently from kernel methods, and have much richer implicit regularization (Chizat et al., 2019; Woodworth et al., 2019; Yehudai & Shamir, 2019; Geiger et al., 2019) .

6. CONCLUSION:

In this work, probability estimates for the network output behavior (i.e. random walk analysis) provide a priori generalization estimates for l 2 regression problems. We prove that even after training, network gradients between the data points are approximately Gaussian random walks, and the variation of the gradients between the data points is extremely small and depends essentially on the number of breakpoints and the variance of weights. To the best of our knowledge, this paper is the first to show a mechanism of implicit regularization and to prove the generalization bounds by using the implicit regularization for deep (three or more layer) neural networks with ReLU activation. As a result, we also show that over-parameterized deep neural networks can learn C 1 -class functions. Importantly, our analysis is independent of the kernel generalization analysis, and the generalization bounds are different from the NTK inductive bias of the RKHS norm.

A PROOFS:

A.1 PROOF OF LEMMA 3 Lemma 3 (reshown, see §3.1). With probability at least 1 -e -Ω(log 2 m) , for every x (p) , x (q) (p, q ∈ [n], p = q), t ∈ [T ] and l ∈ [L], we have

EQUATION

EQUATION

Proof. We prove Lemma 3 for a fixed l ∈ [L] and t ∈ [T ] every pair of data points x (p) , x (q) (p, q ∈ [n], = q), because we can apply union bound at the end.

Recall the definition g l,(t) (x) (the unit output), G l,(t) (x) (the indicator matrix function) and x(s) = (1 -s)x (p) + sx (q) , s ∈ [0, 1] (a linearly interpolating path between data points x (p) , x (q) ), we have

EQUATION

Using the definition of g l,(t) (s) and g l,(t) (s), we also have

EQUATION

EQUATION

where G l-1,(t) (s) is an indicator diagonal matrix for s ∈ [0, 1] as follows:

EQUATION

We begin by proving that the number of breakpoints of the function

EQUATION

Note that the input to the above φ, that is g

EQUATION

(s) = 0 has at most one solution s = s * .5 If the solution s = s * satisfies inequality 0 < s * < 1, then g l,(t) i (s) has a breakpoint at s = s * . This shows that the number of breakpoints of g l,(t) i (s) is bounded by the number of the linear equations

EQUATION

which is clearly equal to m.

We denote by 0

EQUATION

, where α is the number of breakpoints. For notational simplicity, we set s 0 = 0 and s α+1 = 1. Note that, with probability 1, breakpoints are all distinct from each other. Each breakpoint s β (β ∈ [α]) corresponds to a linear equation g l-1,(t) j (s) = 0 for some j ∈ [m]. In other words, for each breakpoint s β (β ∈ [α]), there exists a unique element j ∈ [m] such that g l-1,(t) j (s) = 0. Therefore, for s = s β , we set j = j β

EQUATION

EQUATION

This implies,

EQUATION

)

Note that from the definition of the indicator matrix,

[G l-1,(t) (x(0))] i,i = 1 {[g l-1,(t) (x(0))]i≥0} ,

and from the definition of the linear approximation, g l-1,(t) (0) = g l-1,(t) (x(0)), which says that G l-1,(t) (0) equals G l-1,(t) (x(0)).

For each i ∈ [m], we may write

EQUATION

where M (t) is a matrix:

M (t) := W l-1,(t) G l-2,(t) (x(0))W l-2,(t) • • • W 1,(t) G 0,(t) (x(0))A (40)

and λ j is the difference of indicator functions:

λ j := 1 g l-1,(t) j (s)≥0 -1 g l-1,(t) j (0)≥0 . (41) For β ∈ [α],

we set an open interval of breakpoints I β := (s β , s β+1 ), and recall the notation of the gradient ∇ β

EQUATION

and the gradient gap x β := ∇ β -∇ β-1 . Note that since there is not a breakpoint in I β , the gradient

∇ β is constant in I β .

Hence, for β ∈ [α] and s ∈ I β , we may write

EQUATION

)

Note that the gradient gap x γ (γ ∈ [α]) can be rewritten as

EQUATION

where ζ γ is a difference of indicator functions:

EQUATION

)

least 1 -e -Ω(log 2 m) , we have

EQUATION

)

Therefore, we have the estimate of x (0)

γ .

Next, we will estimate the 2nd term x

(1) γ . By definition, we may write

EQUATION

)

Note that we use Theorem 1 (the convergence theorem) to write

EQUATION

Thus, applying the estimate of h γ to eq. ( 58), we have

EQUATION

EQUATION

)

The last inequality uses λ ≤ o( √ m), which is indicated by the result of Theorem 1, therefore we have the estimate of x

(1) γ .

Next, we will estimate the 3rd term x

(2) γ . By definition, we may write

EQUATION

)

We set u

(2) γ := ζ γ dh γ and u (2) := (u

EQUATION

β ). Using the randomness of W l and general Hoeffding's inequality (Vershynin, 2018) , for each fixed vector v, we have

EQUATION

and choose τ := log m √ m u (2) 2 , with probability 1 -e -Ω(log 2 m) , we have

EQUATION

)

The last inequality uses λ ≤ o( √ m) and the result of Theorem 1:

EQUATION

Therefore, we have the estimate of x

(2)

γ .

Next, we will estimate the last term x

(3) γ . Note that since ζ γ = ±1, we have | ζ γ | = 1, and conditioning on this event (56) happens, we may write

EQUATION

)

Putting this all together, with probability 1 -e -Ω(log 2 m) , we have

EQUATION

Therefore, we have

EQUATION

This finishes the proof of (7) in Lemma 3. Finally, the proof of (8) in Lemma 3 is idential to the above proof, except the fact that each entry of B follows from N (0, 2/c) insted of N (0, 2/m). We complete the proof.

A.2 PROOF OF GENERALIZATION THEOREM:

Theorem 2 (reshown, see §3). Under the same setting as Theorem 1, with probability at least 1 -e -Ω(log 2 m) , for every

EQUATION

EQUATION

In order to give provable guarantees for the generalization performance of over-parameterized deep networks in the l 2 regression task on the one-dimensional input-space, we apply the above Theorem 2 to bound the error between the trained network output f (T ) (x) and the piecewise linear function f (x). Note that if we apply Theorem 1 to this case (one-dimensional models), simple preprocessing is needed: It just consider a mapping [0, ν] → S 1 + (a half circle).

S 1 + = {x = (x 1 , x 2 ) ∈ R 2 : x 2 = 1, x 2 ≥ 1} It is apparent that there exists a C ∞ diffeomorphism from [0, ν] to S 1

+ . Under this setting, the input from S 1 + is the two-dimension. More precisely, the input needs another axis of coordinates to express the biases term. This procedure is written in the subsection of Dataset and loss function in the paper p.4. As a result, we deal with not the case of d = 1 but the case of d = 3 in Theorem 1. Lemma 5. Under the same setting as Theorem 1 and 2, for δ ∈ (0, 1/2], then with probability at least 1 -(δ + e -Ω(log 2 m) ), we have

EQUATION

Proof. Without loss of generality, we may assume that the entry of the weight matrix B i,j ∼ N (0, 2/c), where c = Θ m 1/3 .6 

For k ∈ [n + 1], we denote by |I k | the length of the interval I k = x (k-1) , x (k) . We set c = c = m 1/3 in Theorem 2, with probability at least 1 -e -Ω(log 2 m) , we have

sup x∈I k f (T ) j (x) - f (T ) j (x) ≤ C j log L+1 m m 1/6 sup k∈[n+1] x (k) -x (k-1) = C j log L+1 m m 1/6 T (85)

where C j is a constant, and T := sup k∈[n+1] x (k) -x (k-1) . Here we set ε : k) ), by eq. ( 85) we have

EQUATION

EQUATION

Note that the piecewise linear function fj linearly connects x (k-1) , fj (x (k-1) ) and k) ) . Thus, the following inequality holds.

EQUATION

EQUATION

We use eq. ( 85) and eq. ( 87) to see that

EQUATION

Note that ε is not dependent on k. Thus, we have

EQUATION

Next, we evaluate T = sup k∈[n+1] x (k) -x (k-1) . Recall that the training dataset {x (k) } k∈[n] is generated from the uniform distribution in the interval [0, ν]. We denote by µ := n/ν the density of training samples. In this setting, it is well known that for every r > 0 (r < ν), the number of training samples in the interval [0, r] follows the Poisson distribution with mean µr, then the sequence of inter-sample distances is independent and identically distributed exponential random variables having mean 1/µ (= ν n ) (Daley & Vere-Jones, 2007) . Define ∆ to be the above exponential distribution having mean ν/n, which is the inter-sample distance. Then, using Markov's inequality, for τ > 0, we have

EQUATION

Thus, we can choose τ := ν nτ , which implies that with probability at least 1 -δ, we have T ≤ ν nδ .

(91)

Combining this with eq. ( 89), with probability at least 1 -(δ + e -Ω(log 2 m) ), we have

EQUATION

)

We complete the proof of the Lemma 5.

Lemma 6. Suppose f * j (x) is a C 1 -class function on [0, ν] and fj (x) is the linear interpolation of the data points x (i) , f * j (x (i) ) i∈ [n] . For δ ∈ (0, 1/2], with probability at least 1 -δ, we have

EQUATION

Proof. We use f * j (x) to denote the first derivative of f * j (x). Since f * j is a C 1 -class function on the compact set [0, ν], f * j is a bounded function. Thus, we can define S (k) and S as follows: S (k) := sup S (k) .

(94)

Since f * j (x) and fj (x) are continuous functions, there exists h (k) (0 ≤ h (k) ≤ x (k) -x (k-1) ) such that sup x∈I k f * j (x) -fj (x) = f * j (x (k-1) + h k ) -fj (x (k-1) + h k ) .

(95)

By the mean value theorem, there exists a real number c (k) ∈ I k such that, fj (x (k-1) + h (k) ) = f * j (c (k) )h (k) + f * j (x (k-1) ).

(96)

Note that by the definition of fj , we have f * j (x (k-1) ) = fj (x (k-1) ), therefore, we also have that

EQUATION

≤

x (k-1) +h (k)

x (k-1) Now we apply the estimate of T from eq. ( 91). We may therefore rewrite eq. ( 100) as

EQUATION

EQUATION

We complete the proof of Lemma 6.

Proof of Theorem 4 (Generalization):

Proof. From eq. ( 89) and eq. ( 93), with probability at least 1 -(δ + e -Ω(log 2 m) ), we have for every j = 1, 2, • • • , c, f (T ) j (x) -f * j (x) ≤ f Using the triangle inequality, this gives

f (T ) j (x), y = 1 2 f (T ) (x) -f * (x) 2 2 = 1 2 c j=1 f (T ) j (x) -f * j (x) 2 ≤ 1 2 C 2 c ν 2 n 2 δ 2 1 + log L+1 m m 1/6 2 (∀x ∈ [0, ν]).

where C := max j=1,••• ,c C j .

Thus, we have the following estimate: With probability at least 1 -(δ + e -Ω(log 2 m) ), we have

EQUATION

)

We complete the proof.

:

It is worth pointing out that in this paper, each gradient is a derivative not with respect to weight parameters but with respect to input parameters.

If {Xi} α i=1 are independent Gaussian random variables, then {Yi} α i=1 , where Yi := X1 + X2 + • • • + Xi, is called a Gaussian random walk.

Without loss of generality, one can rescale and assume x 2 ≤ 1/ √ 2 for every input x. Again, without loss of generality, one can pad each x by an additional coordinate x d+1 to ensure x 2 = 1/ √ 2. Finally, without loss of generality, one can pad each x by an additional coordinate x d+2 to ensure x 2 = 1. This last coordinate x d+2 = 1/ √ 2 is equivalent to introducing a (random) bias term and this procedure is described inAllen-Zhu et al. (2018a).

Although the components of g l (s), f (s) are non-differentiable at breakpoints, with probability 1, the start point s = 0 is not a breakpoint.

Otherwise, the linear function is identically zero (i.e. g l-1,(t) j (s) ≡ 0), which does not affect the number of breakpoints.

For each output y = (y1, y2, • • • , yc), an additional coordinates (yc+1, yc+2, . . . , yc) = (0, 0, . . . , 0) can always be padded to the output. We did not try hard to improve the exponent 1/3.

Figures and Tables:

Peer review (iclr2020):
This paper suggests a new technique to analyze the implicit regularization caused by ReLU activations. They bound the generalization error by two terms: 1) one term that represents the distance between the trained network output and a piecewise linear function built based on the set of training points and 2) another term that represents the distance between the piecewise linear approximation and the desired target. The first term is bounded using a random walk type of analysis, which to the best of my knowledge is novel. I find this technique rather interesting and technically sound, although I do have a number of concerns and I*m at the moment more on the reject side, although I will re-consider my score if the authors can provide satisfactory answers. Generalization to more complex activation functions If I understand correctly, the interpolation technique between two points only works for ReLU functions. If one were to try to generalize the analysis to more complex non-linear functions by using a more complex interpolation schemes, wouldn’t you then have a random walk in high-dimensions? If so, wouldn’t that be a problem given the different properties of Brownian motion in high-dimensions? Generalization to smooth activation functions Another question related to the previous one is whether one could hope to generalize the analysis to smooth activation functions. I believe this is also a drawback of combinatorial techniques such as Hanin and Rolnick which have to rely on the discrete nature of the breakpoints. Generalization bound is only derived for 1-d functions Theorem 2 is derived for each dimension independently while the generalization results in Theorem 4 are for 1-dimensional inputs. Where is the difficulty in generalizing these results to higher dimensions? Prior work on generalization of SGD I was really expecting a discussion about how the generalization bound derived in this paper compares to prior work, e.g. Hardt, Moritz, Benjamin Recht, and Yoram Singer. *Train faster, generalize better: Stability of stochastic gradient descent.* arXiv preprint arXiv:1509.01240 (2015). Kuzborskij, Ilja, and Christoph H. Lampert. *Data-dependent stability of stochastic gradient descent.* arXiv preprint arXiv:1703.01678 (2017). Brutzkus, Alon, et al. *Sgd learns over-parameterized networks that provably generalize on linearly separable data.* arXiv preprint arXiv:1710.10174 (2017). And many others… For instance the bound derived in Hardt et al. is also of the order O(n^-2). The bounds in Kuzborskij are also data-dependent and so are yours since your generalization bound depends on the density of the training points. Can you comment on this? What specific insights do we gain your analysis? Noise SGD My understanding is that the authors assume that the noise of SGD is Gaussian. Although this is commonly used when analyzing SGD, there is evidence that the noise is actually not Gaussian, see e.g. Daneshmand, Hadi, et al. *Escaping saddles with stochastic gradients.* arXiv preprint arXiv:1803.05999 (2018). Simsekli, Umut, Levent Sagun, and Mert Gurbuzbalaban. *A tail-index analysis of stochastic gradient noise in deep neural networks.* arXiv preprint arXiv:1901.06053 (2019). I feel this is worth pointing out and one could perhaps also extend this analysis to Heavy-tail noise. I would expect that the results would still hold in expectation but perhaps with a slightly worse probability. Influence step size SGD Using larger step sizes in the SGD updates increase the influence of the noise. I was expecting this to somehow be captured in your analysis but I fail to see where it appears. Can you comment on this? Proof Lemma 3 The derivation of Eq. 10 does not seem completely justified in the proof in the appendix. The authors essentially prove that the length of the gradient gap is bounded by |S| but why is the coefficient omega distributed according to a normal distribution. It seems to me that you need the noise of SGD to be Gaussian for such statement to hold. Can you confirm this? If so, I think this needs to be clearly stated as an assumption since -- as pointed out above -- this is not necessarily true in practice. Minor: proof Theorem 2 It seems rather trivial but for completeness, you should write the proof of Eq. (6) in Theorem 2. “A priori estimates” This is a terminology that is often used in the paper but never defined. What do you mean by “a priori” in this context? Minor comment I would move footnote 3 directly in the main step. I think it is important to point out that the steps of the random walk correspond to the breakpoints.