Improving Adaptivity via Over-Parameterization in Sequence Models

Authors: yicheng  li
qian  lin
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-11 15:30:52

Abstract: It is well known that eigenfunctions of a kernel play a crucial role in kernel regression.
    Through several examples, we demonstrate that even with the same set of eigenfunctions, the order of these functions significantly impacts regression outcomes.
    Simplifying the model by diagonalizing the kernel, we introduce an over-parameterized gradient descent in the realm of sequence model to capture the effects of various orders of a fixed set of eigen-functions.
    This method is designed to explore the impact of varying eigenfunction orders.
    Our theoretical results show that the over-parameterization gradient flow can adapt to the underlying structure of the signal and significantly outperform the vanilla gradient flow method.
    Moreover, we also demonstrate that deeper over-parameterization can further enhance the generalization capability of the model.
    These results not only provide a new perspective on the benefits of over-parameterization and but also offer insights into the adaptivity and generalization potential of neural networks beyond the kernel regime.

Main paper: 
 introduction: 
 In recent years, the remarkable success of neural networks in a wide array of machine learning applications has spurred a search for theoretical frameworks capable of explaining their efficacy and efficiency. One such framework is the Neural Tangent Kernel (NTK) theory (see, e.g., Jacot et al. [2018] , Allen-Zhu et al. [2019] ), which has emerged as a pivotal tool for understanding the dynamics of neural network training in the infinite-width limit. The NTK theory posits that the training dynamics of wide neural networks can be closely approximated by a kernel gradient descent method with the corresponding NTK, elucidating their convergence behaviors during gradient descent and shedding light on their generalization capabilities. Parallel to this, an extensive literature on kernel regression (see, e.g., Bauer et al. [2007] , Yao et al. [2007] ) has studied its generalization properties, showing its minimax optimality under certain conditions and providing insights into the bias-variance trade-off. Thus, one can almost fully understand the generalization properties of neural networks in the NTK regime by analyzing the kernel regression method.However, the application of NTK theory to analyze neural networks, while invaluable, essentially frames the problem within a traditional statistical method by a fixed kernel. The NTK analysis, by its reliance on the fixed kernel approximation, can not entirely account for the adaptability and flexibility exhibited by neural networks, particularly those of finite width that deviate from the theoretical infinite-width limit [Woodworth et al., 2020] . Moreover, empirical evidence [Wenger et al., 2023, Seleznova and Kutyniok, 2022 ] also suggests that the assumption of a constant kernel during training, a cornerstone of NTK analysis, may not hold in practical scenarios where the network architecture or initialization conditions foster a dynamic evolution of the kernel. Also, Gatmiry et al. [2021] showed the benefits brought by the adaptivity of the kernel on a three-layer neural network. These results underscore the need for a more nuanced understanding of neural network training dynamics, one that considers the intricate interplay between network architecture, initialization, and the optimization process beyond the simplifications of NTK theory.Recently, another branch of research has focused on the over-parameterization nature of neural networks beyond the NTK regime, exploring how over-parameterization can lead to implicit regularization and even improve the generalization. In terms of training dynamics, studies (Hoff [2017] , Gunasekar et al. [2017] , Arora et al. [2019a] , Kolb et al. [2023] , etc.) in this domain have revealed that over-parameterized models, particularly those trained with gradient descent and its variants, exhibit biases towards simpler, more generalizable functions, even in the absence of explicit regularization terms. Moreover, in terms of generalization, recent works [Vaškevičius et al., 2019 , Zhao et al., 2022 , Li et al., 2021a] have shown that in the setting of high dimensional linear regression, over-parameterized models with proper initialization and early stopping can achieve minimax optimal recovery under certain conditions. These results underscore the potential and benefits of over-parameterized models that go beyond the traditional statistical paradigms.In this work, we will incorporate the insights from the kernel regression and the over-parameterization theory to investigate how over-parameterization can improve generalization and also adaptivity under the non-parametric regression framework. As a first step towards this direction, we will focus on the sequence model, which is an approximation of a wide spectrum of non-parametric models including kernel regression. We will show that, by dynamically adapting to the underlying structure of the signal during the training process, over-parameterization method with gradient descent can significantly improve the generalization properties compared with the fixed-eigenvalues method. We believe that our results provide a new perspective on the benefits of over-parameterization and offer insights into the adaptivity and generalization properties of neural networks beyond the NTK regime.Limitations of the (fixed) kernel regression. In this work, we first investigate the limitations of the (fixed) kernel regression method by specific examples, illustrating that the traditional kernel regression method suffers from the misalignment between the kernel and the truth function. We show that even when the eigen-basis of the kernel is fixed, the associated eigenvalues, particularly their alignment with the truth function's coefficients in the eigen-basis, can significantly affect the generalization properties of the method.Advantages of over-parameterized gradient descent. Focusing on the alignment between the kernel's eigenvalues and the truth signal (the truth function's coefficients), we consider the sequence model and introduce an over-parameterization method (8) that can dynamically adjust the eigenvalues during the learning process. We show that with proper early-stopping, the over-parameterization method can achieve nearly the oracle convergence rate regardless of the underlying structure of the signal, significantly outperforming the vanilla fixed-eigenvalues method when the misalignment is severe. In addition, the over-parameterization method is also adaptive by its universal choice of the stopping time, which is independent of the signal's structure.Benefits of deeper parameterization. Moreover, we also consider deeper over-parameterization (14) and explore how depth affects the generalization properties of the over-parameterization method. Our results show that adding depth can further ease the impact of the initial choice of the eigenvalues, thus improving the generalization capability of the model. We also provide numerical experiments to validate our theoretical results in Section C.We denote by ℓ 2 = (a j ) j≥1 | j≥1 a 2 j < ∞ the space of square summable sequences. We write a ≲ b if there exists a constant C > 0 such that a ≤ Cb and a ≍ b if a ≲ b and b ≲ a, where the dependence of the constant C on other parameters is determined by the context. 
limitations of fixed kernel regression: 
 Let us consider the non-parametric regression problem given by y = f * (x) + ε, where ε is the noise with mean zero and variance σ 2 , x ∈ X and X is the input space with µ being a probability measure supported on X . The function f * (x) represents the unknown regression function we aim to learn. Suppose we are given samples {(x i , y i )} n i=1 , drawn i.i.d. from the model. We denote X = (x 1 , . . . , x n ) ⊤ and Y = (y 1 , . . . , y n ) ⊤ .Let k : X × X → R be a continuous positive definite kernel and H k be its associated reproducing kernel Hilbert space (RKHS). The well-known Mercer's decomposition [Steinwart and Scovel, 2012] of the kernel function k gives k(x, y) = ∞ j=1 λ j e j (x)e j (y),where (e j ) j≥1 is an orthonormal basis of L 2 (X , dµ), and (λ j ) j≥1 are the eigenvalues of k in descending order. Moreover, we can introduce the feature map Φ(x) = (λ 1 2 j e j (x)) j≥1 : X → ℓ 2 (as a column vector) such that k(x, x ′ ) = ⟨Φ(x), Φ(x ′ )⟩. With the feature map, a function f ∈ H k can be represented as f (x) = ⟨Φ(x), β⟩ ℓ 2 for some β ∈ ℓ 2 .Defining the empirical loss as L = 1 2n n i=1 (y i -f (x i )) 2 , we can consider an estimator f t = ⟨Φ(x), β t ⟩ ℓ 2 governed by the following gradient flow on the feature spaceEQUATIONThis kernel gradient descent (flow) estimator corresponds to neural networks at infinite width limit by the celebrated neural tangent kernel (NTK) theory [Jacot et al., 2018 , Allen-Zhu et al., 2019 ].An extensive literature [Yao et al., 2007 , Lin et al., 2018 , Li et al., 2024a] has studied the generalization performance of such kernel gradient descent estimator. From the Mercer's decomposition, we can further introduce interpolation spaces for s ≥ 0 asEQUATIONwhich is equipped with the norm ∥f ∥[H k ] s = ∥β∥ ℓ 2 for f = ∞ j=1 β j λ s 2 j e j .Particularly, the interpolation space [H k ] 1 corresponds to the RKHS H k itself. Then, assuming the eigenvalue decay rate λ j ≍ j -γ , the standard results (see, e.g., Yao et al. [2007] , Li et al. [2024a] ) in kernel regression state that the optimal rate of convergence under the source conditionf * ∈ [H k ] s with ∥f * ∥ [H k ] s ≤ 1 is n -sγ sγ+1. However, since the interpolation space [H k ] s is defined via the eigen-decomposition of the kernel, the generalization performance of kernel regression methods is ultimately related to the eigen-decomposition of the kernel and the decomposition of the target function under the basis, so the performance is intrinsically limited by the relation between the target function and the kernel itself. In other words, the choice of the kernel could affect the performance of the method. To demonstrate this quantitatively, let us consider the following examples. Example 2.1 (Eigenfunctions in common order). It is well known that kernels possessing certain symmetries, such as dot-product kernels on the sphere or translation-invariant periodic kernels on the torus, share the same set of eigenfunctions (such as the spherical harmonics or the Fourier basis). If we consider a fixed set of eigenfunctions {e j } j≥1 and a given truth function f * , for two kernels k 1 and k 2 with eigenvalue decay rates λ j,1 ≍ j -γ1 and λ j,2 ≍ j -γ2 respectively, it follows that:f * ∈ [H k1 ] s1 ⇐⇒ f * ∈ [H k2 ] s2 for γ 1 s 1 = γ 2 s 2 .Given that the convergence rate is dependent solely on the product sγ, the convergence rates relative to the two kernels will be identical.Example 2.1 seems to show that when the eigenfunctions are fixed, kernel regression methods yield similar performance across different kernels. However, it's important to note that this similarity is due to both kernels having the same eigenvalue decay order, which aligns with the predetermined order of the basis. In fact, if the eigenvalue decay order of a kernel deviates from that of the true function, even if the eigenfunction basis remain the same, it can lead to significantly different convergence rates. Let us consider the following example to illustrate this point. Example 2.2 (Low-dimensional structure). Consider translation-invariant periodic kernels on the torus T d = [-1, 1) d with the uniform distribution. Then, their eigenfunctions are given by the Fourier basis ϕ m (x) = exp(iπ ⟨m, x⟩), m ∈ Z d . Within this basis, a target function f * (x) can be represented as:f * = m∈Z d f m ϕ m (x).Assuming f * exhibits a low-dimensional structure, specifically f * (x) = g(x 1 , . . . , x d0 ) for some d 0 < d, and considering g belongs to the Sobolev space H t (T d0 ) of order t, the coefficients f m can be shown to simplify to:f m = g m1 , m = (m 1 , 0), m 1 ∈ Z d0 , 0, otherwise.Let us now consider two translation-invariant periodic kernels k 1 and k 2 given in terms of their eigenvalues: k 1 is given by λ m,1 = (1 + ∥m∥ 2 ) -r for some r > d/2, whose RKHS is the fulldimensional Sobolev space H r (T d ); k 2 is given by λ m,2 = (1 + ∥m∥ 2 ) -r for m = (m 1 , 0) and λ m,2 = 0 otherwise. Then, the function f * belongs to both [H k1 ] s and [H k2 ] s for s = t/r. After reordering the eigenvalues in descending order, the decay rates for the two kernels are identified as γ 1 = 2r/d and γ 2 = 2r/d 0 . Thus, the convergence rates with respect to the two kernels are respectively:2t 2t + d and 2t 2t + d 0 .Therefore, we see that when d is significantly larger than d 0 , the convergence rate for the second kernel notably surpasses that of the first.This example illustrates that the eigenvalues can significantly impact the learning rate, even when the eigenfunctions are the same. In the scenario presented, the second kernel benefits from the low-dimensional structure of the target function by focusing only on the relevant dimensions, whereas the first one suffers from the curse of dimensionality since it considers all dimensions. The key point to take away from this example is the alignment between the kernel and the target function. To generalize this example, we can consider the following example where the order of the eigenvalues does not align with the order of the target function's coefficients.Example 2.3 (Misalignment). Let us fix a set of the eigenfunctions (e j ) j≥1 and expand the truth function as f * = j≥1 θ * j e j . Note that by giving (e j ) j≥1 , we already defined an order of the basis in j, but coefficients θ * j of the truth function are not necessarily ordered by j. Suppose that an index sequence ℓ(j) gives the descending order of θ * ℓ(j) . Then we can characterize the misalignment by the difference between ℓ(j) and j. Specifically, we assume thatEQUATIONwhere larger q indicates a more severe misalignment. In terms of eigenvalues, let us consider λ j,1 ≍ j -γ , which is in the order of j, while λ ℓ(j),2 ≍ j -γ , which is in the order of ℓ(j). Then, the convergence rates with the two sequences of coefficients are respectively p p + q and p p + 1 .Therefore, the convergence rates can differ greatly if the misalignment is significant, namely when q is large.From Example 2.2 and Example 2.3, we find that it is beneficial that the eigenvalues of the kernel align with the structure of the target function. However, one can hardly choose the proper kernel a priori, especially when the structure of the target function is unknown, so the fixed kernel regression can be limited by the kernel itself and be unsatisfactory. Motivated by these examples, we would like to explore the idea of an "adaptive kernel approach," where the kernel can be learned from the data.Question: Does the paper discuss the limitations of the work performed by the authors? Answer: [Yes] Justification: We have discussed the limitations in the last section. The assumptions are also explained and justified. Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations. 
adapting the eigenvalues by over-parameterization in the sequence model: 
 Motivated by the examples in the last section, as a first step toward the adaptive kernel approach, we consider adapting the eigenvalues of the kernel with eigenfunctions fixed. To simplify the analysis, we would like to the following sequence model, which captures the essences of many statistical models [Brown et al., 2002 , Johnstone, 2017] .The sequence model Let us consider the sequence model [Johnstone, 2017] z j = θ * j + ξ j , j ≥ 1 (5)where (z j ) j≥1 is the observation, θ * = (θ * j ) j≥1 ∈ ℓ 2 is a sequence of unknown truth parameters and ξ j , j ≥ 1 are (not necessarily independent) ϵ 2 -sub-Gaussian random variables with mean zero and variance at most ϵ 2 . For any estimator θ = ( θj ) j≥1 , the generalization error is measured byR( θ; θ * ) = ∞ j=1 ( θj -θ * j ) 2 .Under the asymptotic framework, we are often interested in the behavior of the generalization error as ϵ → 0. Here, we note that the connection between non-parametric regression and the sequence model yields ϵ 2 ≍ n -1 .To see the connection between the sequence model and the non-parametric regression model, we first write the gradient flow (2) in the RKHS in the matrix form asβt = -∇ β L = - 1 n Φ(X)Φ(X) ⊤ β t + 1 n Φ(X)y,where the feature matrix Φ(X) = (Φ(x 1 ), . . . , Φ(x n )) ∞×n and y = (y 1 , . . . , y n ) ⊤ . Now, since the eigenfunctions (e j ) j≥1 are fixed, intuitively, the gradient flow can be diagonalized in the eigen-basis since 1 n Φ(X)Φ(X) ⊤ ≈ Λ = diag(λ 1 , λ 2 , . . . ) and the noise components are approximately normal with variance σ 2 /n by the central limit theorem. Thus, we reach the sequence model. We refer to Subsection B.1 for a more detailed explanation of the connection between the sequence model and the kernel regression model.Regarding the power series expansion (3) in RKHS, for a sequence (λ j ) j≥1 of descending positive numbers (e.g., λ j = j -γ ), we can consider similarly the parameterization θ j = λ 1 2 j β j , j ≥ 1 in ℓ 2 . Since (λ j ) j≥1 corresponds to the eigenvalues of the kernel in the kernel regression, here we also refer to (λ j ) j≥1 as the "eigenvalues" with a little abuse of terminology.With the component-wise loss function L j (θ j ) = 1 2 (θ j -z j ) 2 , we can apply a gradient descent (gradient flow) with early stopping to derive a component-wise estimator θj . If we directly parameterize θ j = λ 1 2 j β j with only β j trainable, we obtain the vanilla gradient descent method, which is just the diagonalized version of the kernel gradient descent. The estimator is simply given by θj = (1 -e -λj t )z j , where t is the stopping time, and its generalization error is easily computed asEQUATION)We note here that these quantities also correspond to generalization error in the (fixed) kernel regression setting [Li et al., 2024a] . In particular, under the setting of (4) and λ j ≍ j -γ , by choosing t ≍ ϵ -2qγ p+q , we obtain the convergence rate ϵ 2p p+q , which is far from optimal if q is large.By the discussion in the previous section, we find it essential to adjust the eigenvalues beyond the fixed ones (λ j ) j≥1 . Inspired by the over-parameterization nature of neural networks, we can also consider over-parameterization with gradient descent in our sequence model to train the eigenvalues: Replacing λ 1/2 j with trainable parameter a j , let us parameterizeEQUATIONwhere a j aims to learn the eigenvalues and β j aims to learn the signal. We consider the following gradient flow (simultaneously for each component j):EQUATIONHere, (λ j ) j≥1 serves as the initial eigenvalues, while the trainable parameters (a j ) j≥1 are updated to adjust the eigenvalues during the training process.To state our results with the most generality, let us introduce the following quantities on the target parameter sequence θ * :EQUATION)The quantity Φ(δ) measures the number of significant components in the target parameter sequence θ, while Ψ(δ) measures the contribution of the insignificant components, which are commonly considered in the literature on the sequence model [Johnstone, 2017] . For the concrete setting of (4), it is easy to show thatEQUATION)Moreover, we make the following assumption on the span of the significant components. Assumption 1. There exists constants κ ≥ 0 and C sig > 0 such thatEQUATION)Assumption 1 says that the span of the significant components, namely those with θ * j ≥ δ, grows at most polynomially in 1/δ. This assumption is mild and holds for many practical settings, such as cases considered in Example 2.3 (κ = 2q p+1 for the first kernel). In other perspective, it imposes a mild condition on the misalignment between the ordering of the truth signal and the ordering of the eigenvalues, where κ measures the misalignment between the ordering of θ j and the ordering of j itself. Then, the following theorem characterizes the generalization error of the resulting estimator. Theorem 3.1. Consider the sequence model (5) under Assumption 1. Fix λ j ≍ j -γ for some γ > 1 and let θOp be the estimator given by the gradient flow (8) stopped at time t. Then, there exists some constants B 1 , B 2 > 0 such that when B 1 ϵ -1 ≤ t ≤ B 2 ϵ -1 , we haveER( θOp , θ * ) ≲ ϵ 2 Φ(ϵ) + ϵ -1/γ + Ψ (ϵ ln(1/ϵ)) as ϵ → 0.(12)Let us further introduce deeper over-parameterization by adding extra D-layers:θ j = a j b D j β j (13)and consider the gradient flowEQUATIONwhere b 0 is the common initialization of all b j . We remark here if one considers the overparameterization θ j = a j b j,1 • • • b j,D β j with the same initialization b j,k = b 0 , k = 1, . . . , D, then b j,k 's remain to be the same by symmetry, so this is equivalent to our parameterization θ j = a j b D j β j . The following theorem presents an upper bound for the generalization error by this deeper overparameterized gradient flow. Theorem 3.2. Consider the sequence model (5) under Assumption 1. Fix λ j ≍ j -γ for some γ > 1 and let θOp,D be the estimator given by the gradient flow (14) stopped at time t. Then, by choosing b 0 ≍ ϵ 1 D+2 , there exists some constants B 1 , B 2 > 0 such that whenB 1 ϵ -1 ≤ b D 0 t ≤ B 2 ϵ -1 , we have ER( θOp,D , θ * ) ≲ ϵ 2 Φ(ϵ) + ϵ -2 D+2 1 γ + Ψ (ϵ ln(1/ϵ)) as ϵ → 0.(15)Benefits of Over-parameterization Theorem 3.1 and Theorem 3.2 demonstrate the advantage of over-parameterization in the sequence model. Compared with the vanilla fixed-eigenvalues gradient descent method, the over-parameterized gradient descent method can significantly improve the generalization performance by adapting the eigenvalues to the truth signal. For a more concrete example, if we consider the setting of (4), plugging (10) yields the following corollary. Corollary 3.3. Consider the over-parameterized gradient descent in (8) (setting D = 0) or (14). Suppose (4) holds and λ j ≍ j -γ for γ > 1 and γ ≥ p+1 D+2 . Then, by choosing b 0 ≍ ϵEQUATION)In comparison, the vanilla gradient flow method yields the rate ϵ 2p p+q .Ignoring the logarithmic factor, Corollary 3.3 shows that the over-parameterized gradient descent method can achieve a nearly optimal rate ϵ 2p p+1 , while the vanilla gradient descent method only achieves the rate ϵ 2p p+q . The improvement is significant when q is large, which corresponds to the case that the misalignment between the ordering of the truth signal and the ordering of the eigenvalues is severe. Moreover, if we return to the low-dimensional regression function in Example 2.2 with the isotropic kernel k 1 , we can see that while the vanilla gradient descent method suffers from the curse of dimensionality with the rate 2t 2t+d , the over-parameterization leads to the dimension-free rate 2t 2t+d0 . Therefore, the over-parameterization significantly improves the generalization performance.Learning the eigenvalues To further investigate how the eigenvalues are adapted by overparameterized gradient descent, we present the following proposition. Proposition 3.4. Given the same conditions as in Theorem 3.2 or Theorem 3.1 (with D = 0 and b D j = 1 for Theorem 3.1), the term learning the eigenvalues a j (t)b D j (t) is non-decreasing in t for each j. Moreover, letting δ ∈ (0, 1), when ϵ is small enough, the following holds at time t chosen as in Theorem 3.1 or Theorem 3.2:• Signal component: There exist constants C, c > 0 such that for any component satisfying θ * j ≥ Cϵ ln(1/ϵ), it holds with probability at least 1 -δ thatEQUATION• Noise component: There exist constants c, C, C ′ > 0 such that, for any component where θ * j ≤ ϵ and λ j ≤ cϵ 2 D+2 , it holds with probability at least 1 -δ thatEQUATION)From this proposition, we can see that for the signal components, the eigenvalues are learned to be at least a constant times a certain power of the truth signal magnitude. Thus, over-parameterized gradient descent adjusts the eigenvalues to match the truth signal as expected. In the case of noise components, although the eigenvalues are still increasing due to the training process, the eigenvalues do not exceed the initial values by some constant factor, provided that λ j is relatively small. This finding suggests that over-parameterized gradient descent effectively adapts eigenvalues to the truth signal while mitigating overfitting to noise. We remark that when λ j is relatively large, the method still tends to overfit the noise components, contributing an extra ϵ -2 D+2 1 γ term in the generalization error, but this term becomes negligible for large γ. Moreover, we also remark that there is a ln(1/ϵ) gap between the signal and noise components. This is because the signal and the noise can not be distinguished for the components in the middle.Adaptive choice of the stopping time A notable advantage of the over-parameterized gradient descent method is its adaptivity. Consider the scenario described by (4), vanilla gradient descent requires the selection of a stopping time t ≍ ϵ -(2qγ)/(p+q) to achieve the optimal convergence rate. However, this choice of stopping time critically depends on the unknown parameters p and q of the truth parameter, posing a significant challenge in practical applications. In contrast, the over-parameterized gradient descent only need to choose the stopping time as t ≍ ϵ -2D+2 D+2 , which does not rely on the unknown truth parameters, while still achieving the nearly optimal convergence rate. This independence from the truth parameters allows the over-parameterization approach to adaptively accommodate any truth parameter structure by employing a fixed stopping time, without the need for prior knowledge about the truth function's properties.Effect of the depth The results in Theorem 3.2 also show that deeper over-parameterization can further improve the generalization performance. In the two-layer over-parameterization, the extra term ϵ -1/γ in Theorem 3.1 emerges due to the limitation of the adapting large eigenvalues. With the introduction of depth, namely adding extra D layers to the model with proper initialization, this term can be improved to ϵ -2 D+2 1 γ in Theorem 3.2. This improvement suggests that the depth can refine the model's sensitivity to eigenvalue adaptation, enabling a more nuanced adjustment to the underlying signal structure. This finding underscores the importance of model depth in enhancing the learning process, providing also theoretical evidence for the empirical success of deep learning models.Comparison with previous works We will compare our results with the existing literature [Zhao et al., 2022 , Li et al., 2021a , Vaškevičius et al., 2019] on the generalization performance of overparameterized gradient descent in the following aspects:• Problem settings: While the existing literature [Zhao et al., 2022 , Li et al., 2021a , Vaškevičius et al., 2019] investigate the realms of high-dimensional linear regression, focusing on implicit regularization and sparsity, the present study dives into kernel regression and its approximation by Gaussian sequence models, emphasizing the adaptivity of overparameterization to the underlying signal's structure, a leap towards understanding model complexity beyond mere regularization. Moreover, while the literature primarily focuses on the setting of strong signal, weak signal and noise separation, we consider the more general setting of the sequence model with arbitrary signal and noise components. differently for each component. Moreover, we extend the over-parameterization to deeper models by adding extra D layers. Although Vaškevičius et al. [2019] and the subsequent work Li et al. [2021a] also consider the deeper over-parameterization, their over-parameterization is in the form of θ = u ⊙D -v ⊙D with u(0) = v(0) = α1. Unfortunately, though being easy to analysis because of the homogeneous initialization, this setup could not bring insights into the learning of the eigenvalues, which is the key to our results. Furthermore, the analysis for Theorem 3.2 involves the interplay between the differently initialized a j and b j , so our analysis is more involved than the existing works. We also remark that although we only consider the gradient flow in the analysis, the results can be extended to the gradient descent with proper learning rates.• Interpretation of the over-parameterization: The previous works view the overparameterization mainly as a mechanism for implicit regularization, while our work provides a novel perspective that over-parameterization adapts to the structure of the truth signal by learning the eigenvalues. Our theory also aligns with the neural network literature [Yang and Hu, 2022, Ba et al., 2022] , where over-parameterization with gradient descent is known to be beneficial in learning the structure of the target function.• Connection to sparse recovery: Our results can be phrased for the setting of high dimensional regression with sparsity. Taking a sparse signal (θ * j ) j≥1 , e.g., θ * j = 1 for j ∈ S, |S| = s and θ * j = 0 for j / ∈ S, we find that Φ(ϵ) = s and Ψ(ϵ) = 0. Consequently, ignoring the extra error term, the resulting rate obtained by Theorem 3.1 or Theorem 3.2 is Õ(s/n) (ignoring the logarithmic factor). This rate coincides with the minimax rate for sparse recovery in high-dimensional regression.In this subsection, we will provide an outline of the proof of Theorem 3.1 and Theorem 3.2. For the detailed proof, we refer to Section D for the analysis of the gradient flow equation and Section E for the generalization error.The proof of Theorem 3.1 and Theorem 3.2 relies on the analysis of the gradient flow ( 8) and ( 14) for each component j. For notation simplicity, we will suppress the index j in the following discussion. Firstly, the symmetry of the equation allows us to consider only the case z > 0. Then, one can find thatd dt a 2 = d dt β 2 = 1 D d dt b 2 = 2θ(z -θ),so we get the conservation quantities a 2 (t) ≡ β 2 (t) + λ and b 2 (t) ≡ Dβ 2 (t) + b 2 0 . Consequently, for the case D = 0, we can obtain the explicit gradient flow of θ:θ = a 4 0 + 4θ 2 (z -θ), θ(0) = 0.Since a 4 0 + 4θ 2 can be bounded by a multiple of a 2 0 + 2θ, we can consider the other equation θ = (a 2 0 + 2θ)(z -θ), which admits a closed-form solution. For the case D ≥ 1, the equation is more complicated. We will apply a multiple stage analysis concerning both the effect of a(t) and b(t).The generalization error In terms of the generalization error, we first separate the noise case when |ξ j | ≥ θ * j /2 and the signal case when |ξ j | < θ * j /2. For the noise case, we apply the analysis of the equation to show that θ j (t) is bounded roughly by λ j for our choice of t. Moreover, the fact that λ j is summable ensures that error of these noise components does not sum up to infinity. On the other hand, for the signal case, if θ * j ≥ ϵ ln(1/ϵ), we can show that our choice of t allows θ j (t) to exceed z/2 and converge to z close enough, so the error in these components is only caused by the random noise and sum up to ϵ 2 Φ(ϵ). In addition, the remaining signal components contribute to the error term Ψ(ϵ ln(1/ϵ)). Summing up these two terms, we can obtain the desired generalization error bound.Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof? • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
numerical experiments: 
 In this section, we provide some numerical experiments to validate the theoretical results. For more detailed numerical experiments, please refer to Section C. We approximate the gradient flow equation ( 22) and ( 30) by discrete-time gradient descent and truncate the sequence model to the first N terms for some very large N . We consider the settings as in Corollary 3.3 that θ * is given by (4) for some p > 0 and q ≥ 1. We set ϵ 2 = n -1 , where n can be regarded as the sample size, and consider the asymptotic generalization error rates as n grows.We first compare the generalization error rates between vanilla gradient descent and overparameterized gradient descent (OpGD) in Figure 1 on page 10. The results show that the overparameterized gradient descent can achieve the optimal generalization error rate, while the vanilla gradient descent suffers from the misalignment caused by q > 1 and thus has a slower convergence rate. Moreover, with a logarithmic least-squares fitting, we find that the resulting generalization error rates also match the theoretical results in Corollary 3.3 (0.5 for OpGD and 0.33 for vanilla GD).Additionally, we investigate the evolution of the eigenvalue terms a j (t)b D j (t) over time t as discussed in Proposition 3.4. The results are shown in Figure 2 on page 10. We find that the eigenvalue terms can indeed adapt to the underlying structure of the signal: for large signals, the eigenvalue terms approach the signals as the training progresses, while for small signals, the eigenvalue terms do not increase significantly. Moreover, we find that deeper over-parameterization reduces the fluctuations of the eigenvalue terms for the noise components, and thus improves the generalization performance of the model. In summary, the numerical experiments validate our theoretical results and provide insights into the adaptivity and generalization properties of the over-parameterized gradient descent method. Figure 1 : Comparison of the generalization error rates between vanilla gradient descent and overparameterized gradient descent (OpGD). We set p = 1 and q = 2 for the truth parameter θ * , and γ = 1.5 for the left column and γ = 3 for the right column. For each n, we repeat 64 times and plot the mean and the standard deviation.  
conclusion and future work: 
 In this work, we studied the generalization properties of over-parameterized gradient descent in the context of sequence models. We showed that the over-parameterization method can adapt to the underlying structure of the signal and significantly outperform the vanilla fixed-eigenvalues method. These results provide a new perspective on the benefits of over-parameterization and offer insights into the adaptivity and generalization properties of neural networks beyond the kernel regime.However, there are also limitations of this work and many interesting directions for future research. For example, one can directly consider the over-parameterization in the kernel regression by replacing the feature map Φ(x) = (λ 1/2 j e j (x)) j≥1 with the learnable one Φ(x; a) = (a j e j (x)) j≥1 , where a j 's are learnable parameters initialized by a j (0) = λ 1/2 j . However, the analysis would be more challenging since now the components are mutually coupled in the gradient flow dynamics.Perhaps one of the most interesting directions is to study how the over-parameterization method can also learn the eigenfunctions of the kernel during the training process, which leads to the truly "adaptive kernel method". We believe that future studies on this topic will provide a deeper theoretical understanding of the success of neural networks in practice.  
experimental setting/details: 
 Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results? Answer: [Yes] Justification: The related codes are provided in the supplementary material. Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
experiment statistical significance: 
 Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments? Answer: [Yes] Justification: Error bars are reported in Figure 3 on page 21. Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)• The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction). 10. Broader Impacts Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed? Answer: [NA] Justification: It is mainly a theory paper, so there is no societal impact of the work Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.• The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA] Justification: It is mainly a theory paper.Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?Answer: [NA] Justification: The paper does not use existing assets.Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. • If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. 
 
 Appendix: In this section, we will provide additional related works and further discussions.The regression problem with a fixed kernel has been well studied in the literature. It has been shown that with proper regularization, kernel methods can achieve the minimax optimal rates under various conditions [Caponnetto and De Vito, 2007 , Steinwart et al., 2009 , Lin et al., 2018 , Fischer and Steinwart, 2020 , Zhang et al., 2023] . Recently, a sequence of works provided more refined results on the generalization error of kernel methods [Li et al., 2023b , Bordelon et al., 2020 , Cui et al., 2021 , Mallinar et al., 2022 , Li et al., 2023a , 2024a] .The NTK regime of neural networks Over-parameterized neural networks are connected to kernel methods through the neural tangent kernel (NTK) theory proposed by Jacot et al. [2018] , which shows that the dynamics of the neural network at infinite width limited can be approximated by a kernel method with respect to the corresponding NTK. The theory was further developed by many follow-up works [Arora et al., 2019c ,b, Du et al., 2018 , Lee et al., 2019 , Allen-Zhu et al., 2019] . Also, the properties on the corresponding NTK have also been studied [Geifman et al., 2020 , Bietti and Bach, 2020 , Li et al., 2024b] .There has been a surge of interest in understanding the role of over-parameterization in deep learning. One perspective is that over-parameterized models trained by gradient-based methods can expose certain implicit bias towards simple solutions, which include linear models [Hoff, 2017 , Vaškevičius et al., 2019 , Zhao et al., 2022 , Li et al., 2021a] , matrix factorization [Gunasekar et al., 2017 , Arora et al., 2019a , Li et al., 2021b , Razin et al., 2021 , Chou et al., 2023] , linear networks [Yun et al., 2021 , Nacson et al., 2022] and neural networks [Kubo et al., 2019 , Woodworth et al., 2020] . Moreover, variants of gradient descent such as stochastic gradient descent are also shown to have implicit regularization effects [Li et al., 2022 , Vivien et al., 2022] . However, most of these works focus only on the optimization process and the final solution, but the generalization performance is not well understood.Generalization Guarantees for Over-parameterized Models Being the most related to our work, a few works provided generalization guarantees for over-parameterized models, which only include linear models [Zhao et al., 2022 , Li et al., 2021a , Vaškevičius et al., 2019] and single index model [Fan et al., 2021] . In detail, the two parallel works [Zhao et al., 2022 , Vaškevičius et al., 2019] studied the high-dimensional linear regression problem under sparse settings and showed that a two-layer diagonal over-parameterized model with proper initialization and early stopping can achieve minimax optimal recovery under certain conditions. The subsequent work [Li et al., 2021a] obtained similar results for multi-layer diagonal over-parameterized models.The adaptive kernel perspective The idea of an adaptive kernel has appeared in a few recent works in various forms [Chen et al., 2023 , Gatmiry et al., 2021 , LeJeune and Alemohammad, 2023 , Yang and Hu, 2022 , Ba et al., 2022] , which is also known as "feature learning". Notably, Gatmiry et al. [2021] showed the benefits brought by the adaptivity of the kernel on a three-layer neural network, which is similar to our work in the adaptive kernel perspective. However, our work and theirs consider different aspects of the adaptive kernel: while they considered an adaptive kernel space in the form of G ⊙ K ∞ around the NTK space, we consider an eigenvalue-parameterized kernel space with fixed eigen-basis. We believe that these various results, including ours, will contribute to a better understanding of the generalization properties of over-parameterized models as well as neural networks.B.1 The connection between RKHS and the sequence model Diagonalized kernel gradient flow as sequence model Moreover, this connection can also be seen directly from the following. The Mercer's decomposition of the RKHS H associated with the kernel k also provides a series representation of the RKHS:EQUATIONwhere we denote byℓ 2 = (a j ) j≥1 | ∞ j=1 a 2 j < ∞ .Therefore, by introducing the feature mapping Φ : X → ℓ 2 defined byEQUATIONwe establish a one-to-one correspondence between a function f ∈ H in the RKHS and a vector β ∈ ℓ 2 in feature space via f (x) = ⟨Φ(x), β⟩ ℓ 2 for f ∈ H and β ∈ ℓ 2 . Moreover, it is convenient to consider Φ(x) as column vectors and also define the feature matrix Φ(X) = (Φ(x 1 ), . . . , Φ(x n )) for X = (x 1 , . . . , x n ). Then, the gradient flow (2) in the feature space ℓ 2 writesEQUATIONIntuitively, since the j, l-th entry1 n Φ(X)Φ(X) ⊤ j,l = λ 1 2 j λ 1 2 l n n i=1 e j (x i )e l (x i )the law of large numbers implies that 1 n Φ(X)Φ(X) ⊤ converges to the diagonal operator Λ = diag(λ 1 , λ 2 , . . . ); moreover, since j-th entry1 n Φ(X)y j = λ 1 2 j n n i=1 f (x i )e j (x i ) + λ 1 2 j n n i=1e j (x i )ε i , the central limit theorem suggest that it can be approximated by λ 1 2 j z j , where z j = θ j + ξ j and ξ j is a normal random variable with mean zero and variance σ 2 /n. Therefore, with these approximations, the equation can be diagonalized asβj (t) = -λ j β j + λ 1 2 j z j = -λ j (β j (t) -λ -1 2 j z j ). Moreover, we can rewrite the representation f (x) = ⟨Φ(x), β⟩ ℓ 2 into f = ∞ j=1 λ 1 2 j β j e j = ∞ j=1 θ j e j , θ j = λ 1 2 j β j .Then, in terms of θ, we haveβj (t) = λ 1 2 j βj (t) = -λ j (θ j -z j ).This is exactly the vanilla gradient flow in the sequence model in Section 3. Furthermore, we can consider the parameterized feature map Φ a (x) = (a j e j (x)) j≥1 , where a = (a j ) j≥1 . We can consider similar gradient flow in the feature space with both β and a trainable. Then, with similar diagonalizing argument, we can show that the corresponding gradient flow in the sequence model is just the over-parameterized gradient flow in (7). Similar correspondence can be established for the multi-layer models (13).The deduction in Example 2.1 is straightforward. The series expansion (3) can be written as[H k ] s =    f = j≥1 f j e j | j≥1 f 2 j λ -s j < ∞    .We recall that λ j,1 ≍ j -γ1 and λ j,2 ≍ j -γ2 . Let f * = j≥1 f * j e j . Then, for ℓ = 1, 2,f * ∈ [H k ℓ ] s ℓ ⇐⇒ j≥1 (f * j ) 2 λ -s ℓ j,1 < ∞ ⇐⇒ j≥1 (f * j ) 2 j γ ℓ s ℓ < ∞ Consequently, f * ∈ [H k1 ] s1 ⇐⇒ f * ∈ [H k2 ] s2 for γ 1 s 1 = γ 2 s 2 . B.2.2 Example 2.2We justify the claims in Example 2.2. Let us consider the torus T d = [-1, 1) d and the uniform measure µ on T d (namely µ(T d ) = 1). Let us recall that the multidimensional Fourier basis is given by ϕm (x) = exp(iπ ⟨m, x⟩) for m ∈ Z d .The Sobolev space H s (T d ) is defined via the Fourier coefficients asH s (T d ) =    f ∈ L 2 (T d ) | m∈Z d |f m | 2 (1 + ∥m∥ 2 ) s < ∞    ,which is equipped with the inner product (as thus the induced norm) Andreas Christmann, 2008] . Using the boundedness of ϕ m , we have⟨f, g⟩ H s (T d ) = m∈Z d f m g m (1 + ∥m∥ 2 ) s . Now, we briefly show that H s (T d ) is an RKHS when s > d/2. It suffices to show that f → f (x) is bounded for each x ∈ T d [m∈Z d |f m ϕ m | ≤ m∈Z d |f m | ≤   m∈Z d |f m | 2 (1 + ∥m∥ 2 ) s   1 2   m∈Z d (1 + ∥m∥ 2 ) -s   1 2Now,m∈Z d (1 + ∥m∥ 2 ) -s ≲ x∈R d (1 + |x| 2 ) -s dx ≲ ∞ 0 (1 + r 2 ) -s r d-1 dr.Since the last integral is finite when s > d/2, we find that there is a constant C such thatm∈Z d |f m ϕ m | ≤ C∥f ∥ H s (T d ) .Therefore, the series expansion f (x) = m∈Z d f m ϕ m (x) converges absolutely and uniformly, andthus |f (x)| ≤ m∈Z d |f m ϕ m | ≤ C∥f ∥ H s (T d ) , showing that f → f (x) is bounded for each x ∈ T d .Moreover, it is easy to see from the Mercer's decomposition (1) and the power series expansion (3) that the kernel ofH s (T d ) is given by k(x, x ′ ) = m∈Z d (1 + ∥m∥ 2 ) -s ϕ m (x)ϕ m (x ′ ), so its eigenvalues are λ m = (1 + ∥m∥ 2 ) -s .To determine the eigenvalue decay rate of λ m = (1 + ∥m∥ 2 ) -s after reordering them in decreasing order, it suffices to determine the count # {m : λ m > δ}: the eigenvalue decay rate is β if # {m : λ m > δ} ≍ δ -1/β , see, e.g., Proposition A.1 in Li et al. [2024a] . We have# {m : λ m > δ} = # m : (1 + ∥m∥ 2 ) -s > δ ≍ Vol( x ∈ R d : (1 + |x| 2 ) -s > δ ) ≍ Vol( x ∈ R d : |x| < δ -1 2s ) ≍ δ -d 2s .We consider a function f * (x) = g(x 1 , . . . , x d0 ) with low-dimensional structure. Let us denote by x ≤d0 = (x 1 , . . . , x d0 ) and x >d0 = (x d0+1 , . . . , x d ) for simplicity. Then, the Fourier coefficients of f * are given byf m = ⟨f * , ϕ m ⟩ L 2 (T,dµ) = 2 -d T d 0 ×T d-d 0 g(x 1 , . . . , x d0 ) exp(iπ ⟨m ≤d0 , x ≤d0 ⟩) • exp(iπ ⟨m >d0 , x >d0 ⟩)dx ≤d0 dx >d0 = 2 -d0 T d 0 g(x 1 , . . . , x d0 ) exp(iπ ⟨m ≤d0 , x ≤d0 ⟩)dx ≤d0 • 2 -(d-d0) T d-d 0 exp(iπ ⟨m >d0 , x >d0 ⟩)dx >d0 = g m ≤d 0 • 1 {m>d 0 =0} , so we show that f m = g m ≤d 0 , m = (m ≤d0 , 0), m ≤d0 ∈ Z d0 , 0, otherwise.We recall that g ∈ H t (T d0 ), som ≤d 0 ∈Z d 0 g m ≤d 0 2 (1 + ∥m ≤d0 ∥ 2 ) t < ∞.To determine the smoothness of f * on [H k1 ] s and [H k2 ] s , following (3), we computem∈Z d |f m | 2 (1 + ∥m∥ 2 ) r s = m=(m ≤d 0 ,0),m ≤d 0 ∈Z d 0 g m ≤d 0 2 (1 + ∥m ≤d0 ∥ 2 ) r s m ≤d 0 ∈Z d 0 g m ≤d 0 2 (1 + ∥m ≤d0 ∥ 2 ) rs , so f * belongs to [H k1 ] s and [H k2 ] s for s = t/r B.2.3 Example 2.3We recall that f * = j≥1 θ * j e j , θ * l(j) ≍ j -(p+1)/2 and ℓ(j) ≍ j q for p > 0, q ≥ 1, where ℓ(j) gives the descending order of θ * j . To compute the relative smoothness of f * w.r.t. λ j,1 ≍ j -γ , we compute qγs) , so we have s < p/(qγ) (but arbitrarily close) and the corresponding generalization error rate is sγ sγ+1 < p p+q . The generalization error rate w.r.t. λ l(j),2 ≍ j -γ can be computed similarly.j≥1 θ * j 2 λ -s j = j≥1 θ * l(j) 2 λ -s l(j) ≍ j≥1 j -(p+1) (ℓ(j)) γs ≍ j≥1 j -(p+1) j qγs ≍ j≥1 j -1-(p-In this section, we provide numerical experiments to validate the theoretical results. The codes are provided in the supplementary material. We approximate the gradient flow equation ( 22) and ( 30) by discrete-time gradient descent with sufficiently small step size. Moreover, we truncate the sequence model to the first N terms for some very large N . We consider the settings as in Corollary 3.3 that θ * is given by (4) for some p > 0 and q ≥ 1. We set ϵ 2 = n -1 , where n can be regarded as the sample size, and consider the asymptotic performance of the generalization error as n grows. For the stopping time, we choose the oracle one that minimizes the generalization error for each method. We first compare the generalization error rates between vanilla gradient descent and over-parameterized gradient descent (OpGD) in Figure 3 on page 21. The results show that the over-parameterized gradient descent can achieve the optimal generalization error rate, while the vanilla gradient descent suffers from the misalignment caused by q > 1 and thus has a slower convergence rate. Moreover, with a logarithmic least-squares fitting, we find that the resulting generalization error rates are consistent with the theoretical results in Corollary 3.3 (0.5 for OpGD and 0.33 for vanilla GD); the oracle stopping times for the over-parameterized gradient descent also match the theoretical value (0.5). Figure 3 : Comparison of the generalization error rates between vanilla gradient descent and overparameterized gradient descent (OpGD). We set p = 1 and q = 2 for the truth parameter θ * . The left and right columns show respectively the generalization error and the orcale stopping time with respect to n. For the upper row, we set the eigenvalue decay rate γ = 1.5; for the lower row, we set γ = 3. For each n, we repeat 64 times and plot the mean and the standard deviation. Furthermore, we investigate the generalization performance of over-parameterized gradient descent (also with deeper parameterization) for different settings of the truth parameter θ * , the eigenvalue decay rate γ and the depth D. The results are reported in Table C on page 22. We find that the generalization error rates are in general consistent with the theoretical results in Corollary 3.3, while p = 0.6 (r * = 0.37) p = 1 (r * = 0.5) p = 3 (r * = 0.75) γ q = 1 q = 1.5 q = 2 q = 1 q = 1.5 q = 2 q = 1 q = 1.5 q = 2 1. 1 : Convergence rates of the over-parameterized gradient descent (8) under different settings of the truth parameter p, q and the eigenvalue decay rate γ, where r * is the ideal convergence rate. The convergence rate is estimated by the logarithmic least-squares fitting of the generalization error with n ranging from 2000, 2200, . . . , 4000, where the generalization error is the mean of 256 repetitions. p = 0.6 (r * = 0.37) p = 1 (r * = 0.5) p = 3 (r * = 0.75) γ q = 1 q = 1.5 q = 2 q = 1 q = 1.5 q = 2 q = 1 q = 1.5 q = 2 D = 1 there are some fluctuations due to the finite sample size. Comparing the results for γ = 1.1 across depth D = 0, D = 1 and D = 2, we see that the method with D = 0 has the slowest convergence rate, while the method with D = 2 has the fastest convergence rate, justifying that deeper parameterization can improve the generalization performance. In summary, the numerical experiments validate our theoretical results.Additionally, we investigate the evolution of the eigenvalue terms a j (t)b D j (t) over time t as discussed in Proposition 3.4. The results are shown in Figure 4 on page 23. We find that the eigenvalue terms can indeed adapt to the underlying structure of the signal: for large signals, the eigenvalue terms approach the signals as the training progresses, while for small signals, the eigenvalue terms do not increase significantly. Moreover, we find that deeper over-parameterization reduces the fluctuations of the eigenvalue terms for the noise components, and thus improves the generalization performance of the model.We also explore the adaptivity of the over-parameterized gradient descent beyond the sequence model. Let us consider the diagonal adaptive kernel method by parameterizing the feature map with Φ(x; a) = (a j e j (x)) j≥1 introduced in Section 5. We use the setting in Example 2.2 where the eigenfunctions are the trigonometric functions. In particular, we consider the truth function f * (x) = sin(7.5πx 1 ) with x = (x 1 , x 2 ) ∈ R 2 . We present the generalization error curve of a single trial and also the generalization error rates with respect to the sample size n in Figure 5 on page 24. The result also shows the benefit of over-parameterization in adapting to the underlying structure of the signal.In this section, we provide additional experiments to test the eigenvalue misalignment phenomenon on real-world data. Recalling Example 2.2 and Example 2.3, we know that the misalignment happens when the order of the eigenvalues of the kernel mismatches the order of coefficients of the truth function. Therefore, to test the misalignment, we compute the coefficients of the regression function over the eigen-basis of the kernel and examine whether the coefficients decay in the order given by the kernel. For the eigen-basis, we use the multidimensional Fourier basis (the trigonometric functions) considered in Example 2.2, where the order is given by the descending order of λ m = (1 + ∥m∥ 2 ) -r . We consider the two real-world datasets: "California Housing" and "Concrete Compressive Strength". We compute the empirical inner product of the regression function with the Fourier basis functions up to a certain order. Then, we plot the coefficients in the order given by the kernel. The results are shown in Figure 6 on page 25. The figures show that the empirical coefficients exhibit significant spikes. Also, among the coefficients, only very few components have large magnitudes, indicating the sparse structure of the regression function. Together, these results suggest that the eigenvalues of the kernel are misaligned with the truth function in these datasets. Let us consider the gradient flow considered in ( 8), where we remove the subscript j for notational simplicity. Let L = 1 2 (θ -z) 2 and θ = aβ. We are interested in the gradient flow:EQUATION)Symmetry of the solution We can find that the solution of the equation for z < 0 can be obtained by simply a(t), -β(t) for the positive signal case of -z > 0. Therefore, we only need to consider the case of z > 0. In this case, it is obvious that a(t), β(t) and θ(t) are all non-negative and increasing.Gradient flow of θ Now we notice that 1 2EQUATION)Using this conservation quantity, we can prove the following estimations:EQUATION)Moreover, the derivative of θ writesθ = a β + ȧβ = (a 2 + β 2 )(z -θ).Using a 2 + β 2 = (a 2 -β 2 ) 2 + 4a 2 β 2 = √ λ 2 + 4θ 2 , we conclude the follow explicit equation for θ:EQUATION)Then, we have the following approximation of the solution.Lemma D.1. Let us consider the gradient flow (25) andEQUATION)Then we haveEQUATION)Moreover, (26) is solved byEQUATION)Proof. It suffices to consider the case z ≥ 0. It is easy to see from the gradient flow ( 22) that both a(t), β(t) are non-negative. Then, using the elementary inequality1 √ 2 (λ + 2x) ≤ λ 2 + 4x 2 ≤ λ + 2x, we have 1 √ 2 (λ + 2x)(z -x) ≤ λ 2 + 4x 2 (z -x) ≤ (λ + 2x)(z -x).Then, the comparison principal in ordinary differential equation yields ( 27). The verification of ( 28) is straightforward.Now let us consider the deeper parameterization of the formEQUATIONwhere a, b, β are all trainable parameters, and the gradient flowEQUATION)Main idea We provide the main idea of analyzing the gradient flow (30) here. Using the conservation quantities, we can first focus on the equation of β. For the initial stage when t is relatively small, β(t) only grows linearly in t. Next, when β(t) exceeds a certain threshold depending on the initialization (and also the interplay between λ j and b 0 ), β(t) grows exponentially in t, provided that θ(t) ≤ z/2. Thus, we can upper bound the hitting time of θ(t) to z/2. Finally, when θ(t) ≥ z/2, we consider directly the equation of θ and show that θ(t) converges to z exponentially fast.Symmetry of the solution Similar to the two-layer case, we can find that the solution of the equation for z < 0 can be obtained by negating the sign of β(t) for the positive signal case. Therefore, we will focus on the case of z ≥ 0 where a(t), b(t), β(t) and thus θ(t) are all non-negative and non-decreasing.Conservation quantities Similarly, we haveEQUATIONEQUATION)Using these conservation quantities, we can prove the following estimations in terms of β:EQUATION)The evolution of θ It is direct to compute thatEQUATIONAuxiliary notations Let us introduceEQUATIONLemma D.2 (Noise case). For the gradient flow (30), we have the initial estimationEQUATIONwhereEQUATIONEQUATIONwhereEQUATIONProof. It suffices to consider the case z > 0. Recalling (30) and θ ≥ 0, we have β ≤ ab D z.Using the upper bound in (33), when t ≤ T esc , namely whenβ(t) ≤ λ 1 2 and √ Dβ(t) ≤ b 0 , we have β ≤ ( √ 2λ 1 2 )( √ 2b 0 ) D z = 2 D+1 2 λ 1 2 b D 0 z, implying that β(t) ≤ 2 D+1 2 λ 1 2 b D0 zt, for t ≤ T esc . Therefore, we get the lower boundT esc ≥ min(T (1) , T (2) ),where T (1) and T (2) are defined by (37) in the lemma. Combining this again with the upper bound that θ = ab D β ≤ 2 D+1 2 a 0 b D 0 β when t ≤ T esc , we prove (36). For the second part, we consider the case λ1 2 ≤ b 0 / √ D.In this case, we have T (1) ≤ T (2) and thus T (1) ≥ T (1) from the above argument. Now, when t ∈ [T (1) , T (2) ], we turn to the following equation: ) , T (2) ], which yieldsβ ≤ ( √ 2β)( √ 2b 0 ) D z = 2 D+1 2 b D 0 zβ, for t ∈ [T (1β(s + T (1) ) ≤ β(T (1) ) exp 2 D+1 2 b D 0 zs = λ 1 2 exp 2 D+1 2 b D 0 zs , for s ∈ [0, T (2) -T (1) ]. Comparing β(s + T (1) ) with b 0 / √ D gives T (2) -T (1) ≥ 2 D+1 2 b D 0 z -1 ln b 0 √ Dλ 1 2 = T (1) ln b 0 √ Dλ 1 2 , so T (2) ≥ T (1) + T (1) ln b 0 √ Dλ 1 2 ≥ T (1) 1 + ln b 0 √ Dλ 1 2 = T (1,2) .Therefore, the comparison principal yieldsβ(t) ≤ λ 1 2 exp 2 D+1 2 b D 0 |z|(t -T (1) ) + for t ≤ T (1,2) ,where we notice that the bound also holds for t ≤ T (1) ≤ T (1) since at that time ) , T (2) ], while the bound also holds for t ≤ T (1) . Lemma D.3 (Signal case). For the gradient flow (30), we have:β(t) ≤ λ 1 2 . Finally, (38) is obtained by using θ = ab D β ≤ 2 D+1 2 b D 0 β 2 when t ∈ [T (1EQUATIONEQUATIONwhereR = ln (D|z|/2) 1 D+2 b0 , D = 1, 1 D-1 , D > 1.Moreover, we haveEQUATIONProof. It suffices to consider the case z > 0. To provide an upper bound of the signal time T sig , we observe that the lower bound in (33) implies a sufficient condition for θ ≥ z/2 thatEQUATIONWe first consider case that λ1 2 ≤ b 0 / √ D. Let us define T (1) := 2(b D 0 z) -1 and suppose that T sig ≥ 2(b D 0 z) -1, otherwise the statement (40) already holds. Then, we first haveEQUATIONThis implies thatβ(t) ≥ 1 2 λ 1 2 b D 0 zt, for t ≤ T sig ,and thusT (1) ≤ T (1) ≤ T sig .Now, for t ∈ [T (1) , T sig ], we use the alternative bound a ≥ β to obtainβ ≥ 1 2 βb D 0 z, for t ∈ [T (1) , T sig ],giving thatβ(s + T (1) ) ≥ λ 1 2 exp 1 2 b D 0 zs , for s ∈ [0, T sig -T (1) ].Comparing it with (43), we obtainT sig -T (1) ≤ 2(b D 0 z) -1 ln (D -D/2 z/2) 1 D+2 λ 1 2 = T (1) ln (D -D/2 z/2) 1 D+2 λ 1 2, which, together with the upper bound of T (1) , proves (40).The case that λ 2) . Using the estimation (44) again, we getEQUATIONT (2) ≤ T (2) ≤ T sig . Now, when t ∈ [T (2) , T sig ], we use b ≥ √ Dβ to obtain β ≥ 1 2 λ 1 2 D D 2 β D z, for t ∈ [T (2) , T sig ].This implies that for s ≤ T sig -T (2) ,β(s + T (2) ) ≥ b 0 √ D exp 1 2 λ 1 2 D D 2 zs , if D = 1, β(s + T (2) ) ≥ (D -1/2 b 0 ) -(D-1) - D -1 2 λ 1 2 D D 2 zs -1 D-1 , if D > 1.Consequently, comparing it with (43) givesT sig -T (2) ≤ 2 λ 1 2 D D 2 z -1 ln (Dz/2) 1 D+2 b 0 = T (2) ln (Dz/2) 1 D+2 b 0 , if D = 1, T sig -T (2) ≤ 2 D -1 √ Dλ 1 2 b D-1 0 z -1 = 1 D -1 T (2) , if D > 1.Finally, let us consider the convergence stage when t ≥ T sig . Since it is always true that θ ≤ z, using the lower bounds in (32), we havez ≥ θ = ab D β ≥ β • D D 2 β D • β = D D 2 β D+2 , implying that β ≤ D -D 2(D+2) z 1 D+2. Now, plugging this into (34) and noticing that θ ≥ z/2 when t ≥ T sig , we deriveθ = θ 2 (a -2 + Db -2 + β -2 )(z -θ) ≥ θ 2 β -2 (z -θ) ≥ 1 4 z 2 D D D+2 z -2 D+2 (z -θ) = 1 4 D D D+2 z 2D+2 D+2 (z -θ).Therefore, we havez -θ(s + T sig ) ≤ (z -θ(T sig )) exp - 1 4 D D D+2 z 2D+2 D+2 s = 1 2 exp - 1 4 D D D+2 z 2D+2 D+2 s ,Notations For notation simplicity, we will use C, c to denote generic positive constants that may change from line to line.E.1 Proof of Theorem 3.1We recall thatER( θ; θ * ) = E ∞ j=1 ( θj -θ * j ) 2 .Let us define the signal eventEQUATIONThen, on S j we have 1 2 θ * j ≤ |z j | ≤ 3 2 θ * j , while on S ∁ j we have |z j | ≤ 3|ξ j |. Then, we decompose( θj -θ * j ) 2 = ( θj -θ * j ) 2 1 Sj + ( θj -θ * j ) 2 1 S ∁ j .Moreover, when the signal is significant, we use( θj -θ * j ) 2 1 Sj ≤ 2( θj -z j ) 2 1 Sj + 2(z j -θ * j ) 2 1 Sj = 2( θj -z j ) 2 1 Sj + 2ξ 2 j 1 Sj .On the other hand, when the noise is dominating, we apply( θj -θ * j ) 2 1 S ∁ j ≤ 2 θ2 j 1 S ∁ j + 2(θ * j ) 2 1 S ∁ j .Summing over j and taking the expectation, we haveEQUATIONNow, the first term (46), representing the absolute error, is controlled by Proposition E.1 thatE ∞ j=1 ξ 2 j 1 Sj + (θ * j ) 2 1 S ∁ j ≤ 4 ϵ 2 Φ(ϵ) + Ψ(ϵ) .Therefore, we focus on the remaining two terms and obtain the estimations ( 49) and ( 51) in the following.The noise term The term (47) represents the extra error caused by the estimator when the noise dominates. Applying Lemma D.1, we obtainθj = |θ j (t)| ≤ θ(t) = λ j (E j -1) 2|z j | + λ j E j |z j | ≤ 1 2 λ j exp ((6|ξ j | + λ j ) t) on S ∁ j .Let us choose J = min {j ≥ 1 :λ j ≤ ϵ} ≍ ϵ -1/γ . Then, since t ≤ B 2 ϵ -1 , we have (|ξ j | + λ j )t ≤ B 2 (|ξ j |/ϵ + 1) and thus E j≥J θ2 j 1 S ∁ j ≤ C j≥J λ 2 j E exp(C(|ξ j |/ϵ) + C) ≤ C j≥J λ 2 j ≤ CJ -(2γ-1) ≤ Cϵ 2-1/γ ,where we notice that E exp(C(|ξ j |/ϵ) + C) is uniformly bounded by some constant for all j since each |ξ j |/ϵ is 1-sub-Gaussian. On the other hand, using the obvious bound θj ≤ 3|ξ j | on S ∁ j , we obtainE j<J θ2 j 1 S ∁ j ≤ C j<J Eξ 2 j ≤ Cϵ 2 J ≤ Cϵ 2-1/γ .Combining two terms, we conclude thatEQUATION)The signal term The term (47) represents the approximation error when the signal is significant. We apply Lemma D.1 again to deriveθj -z j ≤ θ(t/ √ 2) -z j = 2|z j | + λ j 2|z j | + λ j exp (2|z j | + λ j )t/ √ 2 |z j |.Using the fact that 1 2 θ * j ≤ |z j | ≤ 3 2 θ * j on S j , we derive thatEQUATION)Let us define ν = ϵ ln(1/ϵ) ≥ ϵ. Recalling (9) and using Assumption 1, we have j ≤ max J sig (ϵ) ≤ Cϵ -κ , for j ∈ J sig (ν) ⊆ J sig (ϵ).Now, if we take t ≥ B 1 ϵ -1 for some constant B 1 , since θ * j ≥ ν for j ∈ J sig (ν), we also have1 √ 2 θ * j t ≥ 1 √ 2 tϵ ln(1/ϵ) ≥ cB 1 ln(1/ϵ),and thus when j ∈ J sig (ν),ln λ 2 j exp θ * j t/ √ 2 = 2 ln λ j + 1 √ 2 θ * j t ≥ cB 1 ln(1/ϵ) -C ln j ≥ (cB 1 -C) ln(1/ϵ).Consequently, as long as B 1 is large enough, we haveλ 2 j exp θ * j t/ √ 2 ≥ 1 when j ∈ J sig (ν).Therefore, plugging this into (50), we get( θj -z j ) 2 1 Sj ≤ C ( θ * j + λ j ) 2 θ 2 j λ 2 j exp (2 θ * j + λ j )t/ √ 2 ≤ C exp -( θ * j + λ j )t/( √ 2) ( θ * j + λ j ) 2 θ 2 j ,and hencej∈Jsig(ν) E( θj -z j ) 2 1 Sj ≤ j∈Jsig(ν) exp -( θ * j + λ j )t/( √ 2) ( θ * j + λ j ) 2 θ 2 j ≤ C ∞ j=1 ( θ * j + λ j )t -2 ( θ * j + λ j ) 2 θ 2 j = Ct -2 ∞ j=1 θ 2 j ≤ Cϵ 2 .On the other hand, with the trivial bound θj -z j ≤ |z j |, the remaining terms are bounded byj / ∈Jsig(ν) E( θj -z j ) 2 1 Sj ≤ j / ∈Jsig(ν) (θ * j ) 2 = Ψ(ν).Therefore, we conclude thatEQUATIONE.2 Proof of Theorem 3.2Following the same argument as the proof in the previous section, we introduce the events S j and S ∁ j in ( 45) and decompose the error as in ( 46), ( 47) and ( 48). Then, we will focus on the last two terms and derive the estimations ( 53) and ( 54) in the following. We recall that b 0 ≍ ϵ 1 D+2 and we choose t such that B 1 ϵ -1 ≤ b D 0 t ≤ B 2 ϵ -1 for some constants B 1 , B 2 > 0 that will be determined later. Here, we also note that the in correspondence to the component-wise gradient flow considered in Subsection D.2, the initializations are given by λ = λ j and b 0,j = b 0 in (30). Now, since λ j ≍ j -γ , the indexEQUATION)The noise term To apply the bounds in Lemma D.2, let us denote the eventA j =    ω : 3 • 2 D+1 2 b D 0 |ξ j |t ≤ ln b 0 λ 1 2 j √ D    .Then, since |z j | ≤ 3|ξ j | on S ∁ j , we have t ≤ T(1,2) j on A j , where T(1,2) j is defined via (39). Then, for j > J, applying (38) yieldsθ2 j 1 S ∁ j ∩Aj ≤ 2 D+1 b 2D 0 λ 2 j exp 2 D+5 2 b D 0 |z j |t 1 S ∁ j ∩Aj ≤ 2 D+1 b 2D 0 λ 2 j exp 2 D+5 2 B 2 |z j |/ϵ 1 S ∁ j ∩Aj ≤ Cb 2D 0 λ 2 j exp(C|ξ j |/ϵ)1 S ∁ j ∩Ajwhere we use b D 0 t ≤ B 2 ϵ -1 in the last inequality. Consequently,E j>J θ2 j 1 S ∁ j ∩Aj ≤ Cb 2D 0 j>J λ 2 j E exp(C|ξ j |/ϵ) ≤ Cb 2D 0 j>J λ 2 j ≤ Cb 2D 0 J -(2γ-1) ≤ Cb 2(D+2-1/γ) 0, where in the second inequality we notice that E exp(C|ξ j |/ϵ) is uniformly bounded since each |ξ j |/ϵ is 1-sub-Gaussian.On the other hand, noticing b D 0 t ≤ B 2 ϵ -1 again, we havej ∈ A ∁ j =⇒ Cb D 0 t|ξ j | ≥ ln b 0 λ 1 2 j √ D =⇒ |ξ j |/ϵ ≥ cB -1 2 ln b 0 λ 1 2 j √ D = cB -1 2 ln b 2 0 λ -1 j /D .Hence, using Lemma F.1 with the sub-Gaussian property of ξ j and noticing that b 2 0 λ -1 j /D ≥ 1 when j > J, we obtainE j>J θ2 j 1 S ∁ j ∩A ∁ j ≤ C j>J E ξ 2 j 1 |ξ j |/ϵ ≥ cB -1 2 ln b 2 0 λ -1 j /D ≤ Cϵ 2 j>J exp -c cB -1 2 ln b 2 0 λ -1 j /D 2 ≤ Cϵ 2 j>J exp -c ln b 2 0 j γ 2 ≤ Cϵ 2 ∞ J exp -c ln b 2 0 x γ 2 dx ≤ Cϵ 2 b -2/γ 0 ∞ c exp -c [ln(y γ )] 2 dy, y = b 2/γ 0 x ≤ Cϵ 2 b -2/γ 0 ≤ Cϵ 2 ϵ -2 D+2 1 γ .Finally, using the bound θj ≤ 3|ξ j | again, the remaining terms are bounded byE j≤J θ2 j 1 S ∁ j ≤ ϵ 2 J ≤ Cϵ 2 b -2/γ 0 ≤ Cϵ 2 ϵ -2 D+2 1 γ .In summary, we conclude thatEQUATION)The signal term We will apply the bound in Lemma D.3. Let us denoteJ rec = j : t ≥ 2T sig j ∩ J sig (ϵ).Then, when j ∈ J rec and S j holds, (42) and the fact1 2 θ * j ≤ |z j | ≤ 3 2 θ * j imply |θ j -z j | ≤ 1 2 |z j | exp - 1 4 D D D+2 z 2D+2 D+2 (t -T sig j ) ≤ C θ * j exp -c θ * j 2D+2 D+2 t .Consequently,E j∈Jrec ( θj -z j ) 2 1 Sj ≤ C j∈Jrec (θ * j ) 2 exp -c θ * j 2D+2 D+2 t ≤ C j∈Jrec (θ * j ) 2 θ * j 2D+2 D+2 t -D+2 D+1 = C j∈Jrec t -D+2 D+1 ≤ Ct -D+2 D+1 |J sig (ϵ)| ≤ Cϵ 2 Φ(ϵ),where we use exp(-cx) ≤ Cx -D+2 D+1 in the second inequality.Let us define ν = ϵ ln(1/ϵ) ≥ ϵ. We claim that J ∁ rec ⊆ J sig (ν) ∁ on S j as long as b D 0 t ≥ B 1 ϵ -1 for some large constant B 1 . Then, using the obvious bound θj -z j ≤ |z j | ≤ 3 2 θ * j on S j , we haveEQUATIONTo prove the claim, we show that J sig (ν) ⊆ J rec on S j as long as B 1 is large enough. Recalling (9) and using Assumption 1, for j ∈ J sig (ν) ⊆ J sig (ϵ), we havej ≤ max J sig (ϵ) ≤ Cϵ -κ .Now, we show that t ≥ 2T sig j for j ∈ J sig (ν) on S j for different cases in Lemma D.3.1/2 j ≤ b 0 / √ D, we have (40) and thust ≥ 2T sig j ⇐= b D 0 |z j |t ≥ 1 + ln (D -D/2 z/2) 1 D+2 λ 1/2 j + ⇐= B 1 2 |θ j |ϵ -1 ≥ C ln θ * j λ -1 j + C ⇐= B 1 ϵ ln(1/ϵ)ϵ -1 ≥ Cγ ln j + C ⇐= B 1 ln(1/ϵ) ≥ Cκ ln(1/ϵ) + C. • If λ 1/2 j ≥ b 0 / √ D, (41) gives t ≥ 2T sig j ⇐= √ Dλ 1/2 j b D-1 0 |z j |t ≥ 1 + R + j , where R j = ln (D|zj |/2) 1 D+2 b0 , D = 1, 1 D-1 , D > 1.So for both D = 1 and D > 1, we have similarlyt ≥ 2T sig j ⇐= 1 2 √ Dλ 1/2 j b D-1 0 θ * j t ≥ C ln b -1 0 + C ⇐= 1 2 |θ j |b D 0 t ≥ C ln(1/ϵ) + C ⇐= B 1 ln(1/ϵ) ≥ C ln(1/ϵ) + C.Therefore, for both cases, we have t ≥ 2T sig j as long as B 1 is large enough. This finishes the proof of the claim.The following proposition connect the absolute error term with the ideal risk in Johnstone [2017] . Proposition E.1. For the sequence model ( 5), recalling the signal events (45) and the quantities (9), we haveEQUATION)Proof. It is straightforward to see thatξ 2 j 1 Sj + θ 2 j 1 S ∁ j = ξ 2 j 1 {2|ξj |<|θj |} + θ 2 j 1 {2|ξj |≥|θj |} ≤ min(4ξ 2 j , θ 2 j ), so E ξ 2 j 1 Sj + θ 2 j 1 S ∁ j ≤ E min(4ξ 2 j , θ 2 j ) ≤ 4E min(ξ 2 j , θ 2 j ) ≤ 4 min(ϵ 2 , θ 2 j ).Summing over j yields the inequality. The last equality follows from the definition of Φ(ϵ) and Ψ(ϵ).The fact that a j (t)b D j (t) is non-decreasing follows from the analysis of the gradient flow in Subsection D.1 and Subsection D.2. For δ ∈ (0, 1), let us choose C large enough such that P {|ξ j | ≤ Cϵ} ≥ 1 -δ for any fixed j.The case D = 0 For the signal component where θ * j ≥ 2Cϵ ln(1/ϵ), we have |z j | ≥ 1 2 θ * j with high probability. Then, we follow the analysis of the signal term in Subsection E.1 and obtain that|θ j (t) -z j | 2 ≤ C(θ * j ) 2 t -2 ≤ 1 4 (θ * j ) 2 ,provided that ϵ is small enough. This implies that|θ j (t)| ≥ 1 4 θ * j .Now, the second inequality in (24) implies that |θ j (t)| ≤ a 2 j (t). Consequently, we conclude that a j (t) ≥ 1 2 θ * j 1/2 .For the noise component where θ * j ≤ ϵ, we have |z j | ≤ (C + 1)ϵ with high probability. Moreover, since λ j ≍ j -γ , we have J = min {j ≥ 1 : λ j ≤ ϵ} ≍ ϵ -1/γ . Following the similar analysis of the noise term in Subsection E.1, when j ≥ Cϵ -1/γ for some C > 0, we have|θ j (t)| ≤ λ j exp(C(|z j | + λ j )t) ≤ λ j exp(Cϵt) ≤ Cλ j .Then, the first inequality in (24) gives |θ j (t)| ≥ β 2 j (t), so we have|β j (t)| ≤ Cλ 1/2 j . Finally, a j (t) = λ j + β 2 j (t) ≤ λ j + Cλ j ≤ Cλ 1/2 j .The case D ≥ 1 For the signal component where θ * j ≥ 2Cϵ ln(1/ϵ), we still have we have |z j | ≥ 1 2 θ * j with high probability. Now, from the analysis of the signal term in Subsection E.2, we have t ≥ 2T sig j . Moreover, investigating the proof of Lemma D.3, we see that the analysis in Subsection E.2 actually shows that|β j (t)| ≥ c|z j | 1 D+2 ≥ c θ * j 1 D+2 .Consequently, (32) implies thata j (t)b D j (t) ≥ |β j (t)| D+1 ≥ c θ * j D+1 D+2 .For the noise component where θ * j ≤ ϵ, we also have |z j | ≤ (C + 1)ϵ with high probability. Now, we have|z j |b D 0 t ≤ (C + 1)ϵb D 0 t ≤ C 0 , for some constant C 0 .Following the similar analysis of the noise term in Subsection E.2, we can choosej ≥ Cϵ -2 D+2 1 γ such that 1 + ln b 0 λ 1/2 j √ D ≥ C 0 .Now, this condition guarantees that t ≤ T (1,2) defined in Lemma D.2, so Lemma D.2 gives|β(t)| ≤ λ 1 2 j exp Cb D 0 |z j |(t -T (1) ) + ≤ λ 1 2 j exp Cb D 0 |z j |t ≤ Cλ 1 2 j .Combining it with the upper bound in (32) and noticing t ≤ T (1,2) yielda j (t)b D j (t) ≤ 2 D+1 2 |β j (t)|b D 0 ≤ Cϵ D D+2 λ 1 2 j .Lemma F.1. Suppose X is σ 2 -sub-Gaussian, namely, P {|X| ≥ t} ≤ 2 exp -1 2σ 2 t for t ≥ 0. Then for M ≥ 0, we have the tail boundEQUATION)Proof. Using integration by parts, we haveEX 2 1 {X ≥ M } = 2 ∞ 0 rP {|X| ≥ max(M, r)} dr ≤ 4 ∞ 0 r exp - 1 2σ 2 max(M 2 , r 2 ) dr = M 2 + 2σ 2 exp - M 2 2σ 2 ≤ 4σ 2 exp - 1 4σ 2 M 2 .Lemma F.2. Suppose that (θ j ) j≥1 satisfies θ l(j) ≍ j -(p+1)/2 for some p > 0 and |θ j | = 0 otherwise, where l(j) is a sequence of indices. Defining Φ(δ) and Ψ(δ) as in (9), we haveΦ(δ) ≍ δ -2 p+1 , Ψ(δ) ≍ δ 2p p+1 .Proof. First, from the definition of Φ(δ) and Ψ(δ), we see that they do not depend on ordering of the indices and zero values of θ j . Therefore, we can assume that l(j) = j. Then, assuming thatc 1 j -(p+1)/2 ≤ |θ j | ≤ C 1 j -(p+1)/2 , we have Φ(δ) = |{j : |θ j | ≥ δ}| ≤ j : C 1 j -(p+1)/2 ≥ δ ≤ (δ/C 1 ) -2 p+1 .Moreover,Ψ(δ) = ∞ j=1 |θ j | 2 1 {|θ j | < δ} = j>Φ(δ) |θ j | 2 ≤ C 2 1 j>Φ(δ) j -(p+1) ≤ C 2 1 CΦ(δ) -p ≤ C ′ δ 2p p+1for some constant C ′ > 0. The lower bound of them can be obtained similarly.1. Claims Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope? Answer: [Yes] Justification: The main claims are supported by our main theorems as well as the numerical experiments.Guidelines:• The answer NA means that the abstract and introduction do not include the claims made in the paper. • The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers. • The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings. • It is fine to include aspirational goals as motivation as long as it is clear that these goals are not attained by the paper. 


Decision: **Summary:**

This paper investigates how over-parameterization can enhance generalization and adaptivity in non-parametric regression models, particularly focusing on sequence models that approximate kernel regression. The authors introduce an over-parameterized gradient descent method to analyze how varying eigenfunction orders affects performance. Key findings include:

- Limitations of traditional kernel regression methods when eigenvalue alignment is suboptimal
- Over-parameterized gradient flow can adapt to the underlying signal structure and outperform standard methods
- Adding depth to the model further improves generalization by mitigating initialization effects

**Strengths:**

- Novel theoretical analysis providing insights into benefits of over-parameterization
- Clear motivation and presentation of ideas
- Rigorous mathematical framework that could be useful for further research
- Shows limitations of fixed kernel approaches like NTK and provides path forward
- Demonstrates how over-parameterization can lead to adaptive kernel learning

**Weaknesses:**

- Model setting is somewhat simplified/constrained compared to practical neural networks
- Limited experimental validation on real-world datasets and tasks
- Not fully clear how insights translate to more complex architectures used in practice
- Some reviewers felt theoretical claims about optimality of linear decay schedule needed stronger support

**After rebuttal:**

Overall, reviewers found the theoretical contributions and insights valuable, but had some concerns about how directly applicable the results are to practical deep learning. The authors provided reasonable responses explaining the scope of the work as a theoretical foundation to build upon. Most reviewers were satisfied with the responses and some increased their scores, indicating the paper makes a solid contribution to understanding over-parameterization, even if further work is needed to fully bridge to applied settings.

Reviews: review: paper_summary: paper_summary: Authors consider the problem of fitting data with given a kernel but with a different estimator than kernel regression, which is based on an overparameterized version of the gradient flow for the squared loss. This estimator is inspired by the behavior of first-order algorithms in the area of deep learning theory under overparameterized models. The interesting regime in their theoretical result is when the order of the coefficients of the target function (that we hope to learn) in the eigenbasis of the kernel does not match the decay in the eigenvalues of the kernel itself; in this regime, while the conventional kernel regression estimator is suboptimal and can have undesirable dimension dependencies, they show that their estimator achieves a better rate, and intuitively is able to exploit the low dimensionality and adapt the eigenvalues of the kernel to the coefficients of the target function, a task which kernel regression with a fixed kernel is unable to do. Their analysis starts with using the Le cam theorem on the equivalence of kernel fitting and sequence models, then instead of running the gradient flow for the L2 loss of the sequence model only on the functions coefficients, they also train on some hyperparameters for the sequence model that resemble the kernel eigenvalues.
strength_weakness: strength_weakness: My major concern is that similar results to this one already exist in the literature and the authors have not compared their results to them and how they differ from the them.In particular, running gradient-based methods on overparameterized nets, such as the quadratic parameterization is known to recover the underlying sparse solutions [3], also [2] section 6.2. In particular, I am curious if your result can imply or be phrased as a sparse recovery type argument for high high dimensional regression (i.e. similar rate to lasso)? And can you quantitatively compare your end results/rates with these papers, independent of the differences between algorithms/initialization? The work of [Zhao et al., 2022] that you have already cited also seems to have such sparse recovery results for high dimensional regression for gradient descent run on overparameterized model. Please also clarify the difference between your setup/rates with them quantitatively.

Minor comments
This work [1] also shows that training overparameterized three-layer neural nets with a specific architecture leads to the algorithm adaptively picking the kernel (instead of NTK where the kernel is fixed), which seems similar in flavor to your result, so it is good mention/compare.
questions: questions: - Can you compare your rate with the results in sparse recovery, if you use lasso instead of your estimator? specifically do you think the low dimensional example that you mention can be observed as some kind of sparsity that can be exploited by lasso?
limitations: limitations: The literature review/comparison with current result still needs work, upon addressing this I am happy to increase my score.
paper_summary: paper_summary: This paper investigates how over-parameterization can enhance generalization and adaptivity within the non-parametric regression framework. Drawing on insights from kernel regression and over-parameterization theory, the authors focus on the sequence model, which can approximate a wide range of non-parametric models, including kernel regression.

The authors highlight three key findings through mathematically stylized examples. First, they identify limitations of traditional kernel regression, including neural tangent kernel (NTK) theory, by demonstrating that the alignment (ordering) of eigenvalues can significantly impact generalization properties, even when eigenfunctions are correctly specified (Section 2). Second, by analyzing the gradient flow for a two-layer Hadamard neural network corresponding to the sequence model, they claim that this method can dynamically adjust eigenvalues to adapt to the underlying structure of the signal, achieving nearly the oracle rate with suitable regularization (Section 3.1). Finally, they show that adding depth can mitigate the impact of eigenvalue initialization, thereby improving generalization capability (Section 3.2).
strength_weakness: strength_weakness: While this work makes significant theoretical contributions and provides a fresh perspective on the benefits of over-parameterization, there are areas where it could be improved:

**1. Narrow/Restrictive Settings:** The sequence model and the parameterization in Eq. (8) and Eq. (14), viz., the Hadamard neural network, might be too simplistic. I am curious if the insights in this paper can be extended to more realistic parameterizations, such as fully connected networks.

**2. Dynamically Evolving Kernels:** This paper addresses issues with misaligned eigenvalues and the benefits of over-parameterized gradient flow to mitigate them. However, during training, the kernel itself (and thus its eigenfunctions) may evolve, as empirically evidenced by [Seleznova and Kutyniok 2022] and [Wenger et al., 2023]. While the authors mention in Section 4 that this is left for future work, it would be helpful to understand if the eigenvalue alignment discussed in this paper can shed light on over-parameterization in practical settings beyond the stylized examples provided.
questions: questions: 1. The meaning of the sentence in lines 172-175 is unclear. Could the authors elaborate on this? Additionally, it would be helpful if the sequence model were described earlier to enhance its comprehension. Perhaps the authors could include a background section or rearrange the description of the sequence model in lines 183-186 to be presented upfront, with appropriate references.

2. Regarding the adaptive choice of the stopping time discussed in lines 279-288, I can roughly follow the authors' point. However, it appears that there is an unspecified constant, which is not known a priori, needed to choose the stopping time $t \asymp \epsilon^{\frac{2D+2}{D+2}}$, making it not practically feasible either.  Thus, I am wondering if the advantage is more aligned with universality, rather than practicality.  Could the authors clarify this?

3. Including the numerical experiments from Appendix A into the main text would be beneficial, if possible. Additionally, could the authors design a simple experiment to support the practical relevance of eigenvalue misalignment and the discussions in this paper? For instance, would it be possible to demonstrate that it is common in real-world datasets that kernels are well-identified but eigenvalues are easily misaligned?

Typos/minor suggestions:
- Line 89: add transpose to $X$
- LIne 149: defines -> defined
- Line 205: Inspired -> Inspired by
- Line 235: shows the generalization error -> presents an upper bound for the generalization error
- Line 260: have -> present?
limitations: limitations: This is primarily a theoretical work, and the authors discussed the potential limitations of the work and potential future research directions.
paper_summary: paper_summary: This paper proposed an overparameterized gradient descent method. Its benefits on generalization has been verified both theoretically and empirically.
strength_weakness: strength_weakness: 1.The insightful explanations for why the proposed method can improve generalization are lacked.

2.While the experiments in appendix have verified the effectiveness of overparameterized gradient descent, it is recommended to compare the performance of this new method with other algorithms, such as SGD or SAM.
questions: questions: See weakness.
limitations: limitations: No negative social impact.
paper_summary: paper_summary: The authors analyze the generalization error for kernel regression sequence models in the overparameterized regime.
They rigorously show that due the overparameterization, the learned parameters are better adapted to the underlying structure.
They also provide some numerical experiments corroborating their findings (in the appendix).
strength_weakness: strength_weakness: The main concern here, is that the setting is too constrained, in that the kernel map is assumed to be fixed.
questions: questions: - The NTK is mentioned multiple times in the introduction, could you elaborate a bit more on the connection of your work to NTK?
- How would you go about approximating the (infinite sequence ) $\lambda_i$ in practice? How does this affect performance and guarantees?
limitations: limitations: Limitations have been addressed.
paper_summary: paper_summary: This paper investigates the influence of over-parameterization on the adaptivity and generalization of sequence models. The work highlights the significance of eigenfunctions in kernel regression and introduces an over-parameterized gradient descent method to analyze the effects of varying eigenfunction orders. Theoretical results demonstrate that over-parameterized models can adapt to the underlying signal structure and outperform traditional methods. The research shows that deeper over-parameterization further enhances model generalization. This approach provides insights into improving the flexibility and performance of neural networks, especially in practical scenarios where network architecture and initialization conditions dynamically evolve.
strength_weakness: strength_weakness: One significant weakness of the paper is the limited experimental validation. While the theoretical results are compelling, the empirical evidence provided is not sufficient to fully support the claims made about the superiority of the over-parameterized models in practical scenarios. The paper would benefit from a more extensive set of experiments across various datasets and model architectures to demonstrate the robustness and generalizability of the proposed methods. Additionally, there is a need for more practical examples that illustrate how the approach can be applied in real-world settings, especially given the dynamic nature of network architecture and initialization conditions mentioned. Including these aspects would enhance the paper's contributions and practical relevance.
questions: questions: 1. Can you provide more detailed theoretical justification or insights into why deeper layers contribute to better performance? Specifically, how does the introduction of additional layers influence the adaptation of eigenvalues, and are there any diminishing returns or optimal depth considerations that should be taken into account?

2. The theoretical results are promising, but the empirical validation is limited. Have you considered testing the proposed over-parameterized gradient descent method on a wider range of datasets and model architectures? How does the method perform in different domains, such as natural language processing or computer vision, compared to traditional methods and other over-parameterization techniques?
limitations: limitations: None

Comments: title: Response
Re:
Experiment details
Response to the Authors' Rebuttal
Start discussions with authors, please
Dear Reviewers
Remaining response for Question 3
References
comment: Thank you for your insightful comments.

**Regarding the first point**, we apologize for any miscommunication. 
Our intention was to emphasize that while most research on over-parameterized models focuses on structural
results, there has been less exploration into the generalization benefits of over-parameterization. We appreciate your
reference to [1], which is indeed a significant work on the generalization capabilities of over-parameterized neural
networks that we had missed.

Our work shares some conceptual similarities with [1], particularly in viewing over-parameterization through an adaptive
kernel perspective. However, there are still differences in our approach: while [1] explores an adaptive kernel space
in the form of $G\odot K^\infty$ around the NTK space, our study examines an eigenvalue-parameterized kernel space with
a fixed eigen-basis. We believe that the diversity of approaches, including those in [1], will contribute to a deeper
understanding of the generalization properties in over-parameterized models. We will ensure to include [1] in our
revised introduction and add a detailed comparison within the paper.

**Regarding the second point**, there is indeed a fundamental difference between their approach and our approach in the weighted
infinite case. Specifically:

1. We have to track the differently weighted initialization across components;
2. We provide new insights on the benefits of over-parameterization by learning the eigenvalues, see Proposition 3.4;
3. We extend our study to deeper layers, demonstrating additional benefits, whereas Zhao et al., 2022 is limited to the
   two-layer setting.

We appreciate your thorough review and the opportunity to improve our work.
Thank you for your response.

regarding your claims in the new response that previous work in NN theory only show structural result about convergence of gradient based methods and does not consider generalization, that is not true I think, e.g. reference [1] that I sent above also consider the generalization capability of the final NN weights (which achieves similar rates to what you get in classical sparse recovery). So I think I don't understand that part of your claim.

Regarding your comparison with Zhao et al., 2022, you mentioned they consider the l-2 case without weights so it cannot be used for the infinite weighted case, do you mean there is a fundamental barrier in generalizing their approach to cover the weighted infinite case, or is it that they just did not investigate this case in that paper?
I thank the authors for responding to my comments. I'd like to increase my score based on the authors' response.
I thank the authors for answering my questions. I keep my score.
We would like to thank you for raising the score and providing further valuable comments on our paper. 
In the new revision, we will follow your advice to improve the clarity.

For the experiments on the real-world dataset, it seems that we missed the opportunity to post a global rebuttal where we can show figures and tables. 
We apologize for this. We will include the results in the new revision of the paper.
Nevertheless, we would like to describe in detail the experiments and results here.

First, as discussed in Example 2.2 (Low-dimensional example) or Example 2.3, the misalignment happens when the order of the eigenvalues of the kernel mismatches the order of coefficients of the truth function. 
For the experiment setup, we consider the multidimensional Fourier basis (the trigonometric functions) as the eigen-basis of the kernel,
which can be given by $e_{\mathbf{m}}(x) = \exp(2\pi i \langle{\mathbf{m},x}\rangle)$, where $\mathbf{m} \in \mathbb{Z}^d$.
For the commonly considered kernel that correspond to Sobolev space (see also Example 2.2), the kernel and also its eigenvalues are isotropic in the sense that
 $\lambda_{\mathbf{m}} \asymp (1+||{\mathbf{m}}||_2^2)^{-r}$ for some $r > 0$.
Therefore, it gives an order of the eigen-basis. 

For the real-world dataset, we compute the coefficients of the regression function over the eigen-basis of the kernel by the empirical inner product.
For "California Housing" and "Concrete Compressive Strength", the dimension of $x$ is 8, so we choose $\mathbf{m}$ up to $|m_i| \leq 2$, resuling in $5^8 = 390625$ coefficients.
Then, we plot the magnitudes of the coefficients in the order given by the kernel.
If the kernel is well aligned with the truth function, the coefficients should exhibit a smooth decay in this order.
However, the resulting plot show **multiple significant spikes**.
Also, among the coefficients, only very few components have large magnitudes, indicating the **sparse structure** of the regression function.
Together, these results suggest that the eigenvalues of the kernel are misaligned with the truth function in these datasets.

Please feel free to ask if you have any further questions or need more details on the experiments. Thank you.
I thank the authors for addressing my questions and concerns. Since most of my inquiries were for clarification, I trust that the authors will incorporate the additional explanations and reorganization into their revision. However, I was unable to locate the real-world dataset experiments mentioned in response to my Question 3;  if the experiments have already been conducted, could the authors attach them to the global rebuttal? It would be also beneficial to include these (along with the authors' response to my Q3) in the revision.

Overall, I believe this work makes a solid theoretical contribution, though I suggest some revisions to further improve clarity regarding the scope of contributions and limitations, as well as their relation to existing work (e.g., following the authors' discussion with Reviewer EWHq and Reviewer kYgJ) and future directions. Specifically, I recommend the following among others, some of which the authors have already committed to:

  (i) Properly introduce the sequence model, including a concise yet self-contained exposition of the Le Cam equivalence, potentially in the Appendix if necessary.

  (ii) Further clarify the dynamic adaptivity of the kernel. While I agree with the authors that a dynamic update of eigenvalues can lead to adaptive kernel choices, even with a fixed eigenbasis, this offers limited adaptivity compared to a more general adaptive kernel/feature model, as the authors also noted. It would be helpful if the authors could elaborate on the extent to which the current approach is effective and whether they believe an extension to a more general kernel model is needed.

With that understanding, I am inclined to raise my rating to a 6.
Dear Reviewers,

Thank you for your valuable contributions to the review process. As we enter the discussion phase (August 7-13), I kindly request your active participation in addressing the authors' rebuttals and engaging in constructive dialogue.

Please:

- Carefully read the authors' global rebuttal and individual responses to each review.

- Respond to specific questions or points raised by the authors, especially those requiring further clarification from you.

- Engage in open discussions about the paper's strengths, weaknesses, and potential improvements.

- Be prompt in your responses to facilitate a meaningful exchange within the given timeframe.

- Maintain objectivity and professionalism in all communications.

If you have any concerns or need guidance during this process, please don't hesitate to reach out to me.

Your continued engagement is crucial for ensuring a fair and thorough evaluation process. 
Thank you for your dedication to maintaining the high standards of NeurIPS.

Best regards,

Area Chair
We would like to thank the reviewers for their comments. We’ve taken the time to carefully go through each review. If anything is still unclear, please just let us know and we’ll get back to you quickly. Thanks again for your support!
3. > Including the numerical experiments from Appendix A into the main text would be beneficial, if possible. ...

   Thank you for your advice. With the additional content page if our paper is accepted, we will be able to include the
   numerical experiments from Appendix A into the main text to enhance the readability of the paper.

   For real-world datasets, we believe that the eigenvalue misalignment is a common phenomenon in practice.

    * From the theory side, as discussed in Example 2.2 (Low-dimensional example), if the regression function of the
      data exhibits an unknown low dimensional structure on the covariates, then using a fixed kernel that takes all
      covariates into account would lead to misalignment of the eigenvalues.
      Since the low-dimensional structure is very common in real-world datasets, the eigenvalues would be misaligned in
      practice.

    * To conduct experiment on real-world datasets,
      we can estimate the regression function's coefficients over the eigen-basis of the kernel and check whether they
      exhibit a sparse structure.
      We have conducted new experiments on datasets like "California Housing" and "Concrete Compressive Strength", where
      we
      choose the eigen-basis to be the multi-dimensional Fourier basis.
      The results show that the coefficients over the eigen-basis concentrate over a few components, indicating the
      sparse structure of the regression function,
      so misalignment of the eigenvalues would occur.
      Therefore, the eigenvalue misalignment might be very common in real-world datasets.
[Hoff 2017] Peter D. Hoff. Lasso, fractional norm and structured sparse estimation using a Hadamard product
parametrization. Computational Statistics & Data Analysis, 115:186–198, November 2017. ISSN 0167-9473. doi:
10.1016/j.csda.2017.06.007

[Gunasekar et al. 2017] S. Gunasekar, B. Woodworth, S. Bhojanapalli, B. Neyshabur, and N. Srebro. Implicit
regularization in matrix factorization. In Advances in Neural Information Processing Systems, volume 2017-December,
pages 6152–6160, 2017.

[Vaškevicius et al., 2019] Tomas Vaškeviˇcius, Varun Kanade, and Patrick Rebeschini. Implicit regularization for optimal
sparse recovery, September 2019. URL http://arxiv.org/abs/1909.05122.

[Zhao et al., 2022] Peng Zhao, Yun Yang, and Qiao-Chu He. High-dimensional linear regression via implicit
regularization. Biometrika, 109(4):1033–1046, November 2022. ISSN 0006-3444, 1464-3510. doi: 10.1093/biomet/asac010.

[Johnstone 2017] Iain M. Johnstone. Gaussian estimation: Sequence and wavelet models. 2017.