High-Dimensional Linear Regression via Implicit Regularization

1 Introduction

• The method uses the Hadamard product over-parameterization β = g∘l and applies gradient descent directly to the quadratic loss rather than an explicit penalty.This contrasts with the prior lasso reformulation, which used a smoothed optimization problem after a change of variables.Source: §1 · p. 2
• With noisy responses, early stopping yields more accurate estimation than explicitly penalized approaches by overcoming their large bias.The rate is optimal when signals are strong, and the estimator also achieves a rate for settings with both strong and weak signals.Source: §1 · p. 2
• For noiseless responses, near-zero initialization and a properly tuned gradient method converge to a least ℓ1-norm solution when the design satisfies the Restricted Isometry Property.The joint objective is highly non-convex and has exponentially many saddle points.Source: §1 · p. 2
• Numerical studies show encouraging improvements over existing explicitly regularization-based approaches.Source: §1 · p. 2
• The paper extends implicit-regularization analysis in first-order methods beyond work focused exclusively on noiseless or perfectly separable data.It presents gradient descent initialized near zero as a way to induce sparsity with provable guarantees.Source: §1 · p. 2
• The proposed estimator suffers less bias than some penalty-induced estimators and has a natural initialization that provably leads to the optimal solution.The paper contrasts this with non-convex-penalty M-estimators that require stringent initialization conditions to avoid poor local minima.Source: §1 · p. 3

2.1 Setup and notations

• The setup treats β* as an unknown s-sparse vector and defines the Hadamard product componentwise as βj = gjlj.It also establishes notation for subvectors, norms, and comparison of sequences used in the analysis.Source: §2.1 · p. 3

2.2 Related literature

• The method optimizes quadratic loss over β = g ˝ l rather than adding an explicit penalty to the lasso objective.The paper describes this Hadamard-product parameterization as shared with Hoff (2017), while distinguishing its unpenalized gradient-descent objective.Source: §2.2 · p. 3, §2.3 · p. 4

2.3 Gradient descent with Hadamard product parametrization

• For the Hadamard-product objective, all local minima are global and all saddle points are strict, while sufficiently small step size yields almost-sure convergence to a global minimum.The stated landscape has no local maxima; the paper also notes exponentially many saddle points can occur.Source: §2.3 · p. 4

2.4 Gradient descent converges to sparse solution

• Under suitable design conditions, gradient descent initialized arbitrarily close to zero converges to a solution of the basis pursuit problem.The discussion motivates this through the minimal ℓ1-norm solution and notes that the formal result uses a Restricted Isometry Property condition.Source: §2.4 · pp. 5–6

2.5 Gradient descent with early stopping

• With noisy responses, the method combines gradient descent with early stopping to prevent overfitting.The paper describes validation-based stopping and reports that the estimation error is significantly lower than lasso’s in its Figure 1 simulation, while the prediction errors are nearly the same.Source: §2.5 · pp. 6, 8, Figure 1 · p. 7

2.6 Early Stopping criterion

• The authors define hold-out, cross-validation, and SURE-based approaches for selecting the gradient-descent stopping time.Cross-validation is the default in the paper’s subsequent work; its simulations report similar performance for SURE-based stopping and hold-out or cross-validation.Source: §2.6 · p. 9

3.1 Assumptions

• The proof relies heavily on restricted isometry, although supplementary numerical results suggest it may not be necessary; the authors conjecture restricted eigenvalue conditions may suffice.The weaker-condition conclusion is a conjecture, not a proved guarantee.Source: §3.1 · p. 10
• The finite-sample guarantee assumes sub-Gaussian noise, an s-sparse bounded signal, a restricted isometry condition, and small initialization and step size.The theorem gives a probability bound and a permitted interval of stopping times under these assumptions.Source: §3.1 · p. 10, §3.2 · p. 11

3.2 Finite sample error bound

• With only strong signals, the estimator attains a parametric rate independent of dimension; with weak signals, its rate depends on the numbers of strong and weak signals and remains nearly minimax-optimal.The weak-signal regime has rate σ^2(s1+s2 log p)/n, as reported in the passages.Source: §3.2 · p. 11

Simulations

• When the weak signals were too weak for any method to select, the proposed method had true-negative error 4 and false detections close to zero.MCP had similar true-negative error but slightly more false positives; lasso and SCAD had higher false-detection risk.Source: §4 · p. 13

5 Discussion

• The theory relies on RIP, and the authors ask whether it remains valid under the weaker restricted eigenvalue condition and whether data-driven complexity measures can select stopping parameters reliably.They identify computationally efficient alternatives to cross-validation for choosing iteration number and step size as an open practical problem.Source: §5 · p. 14

A.1 Dynamical system interpretation for section 2.4

• Under the stated limiting-point condition, the noiseless dynamics yield a condition that, together with Xβ8 = y, coincides with the KKT condition for basis pursuit.Source: §A.1 · p. 15

A.2 Adaptive step size and variable selection

• Coordinate-specific step sizes yield dynamics related to adaptive Lasso; smaller steps act like larger penalties, while larger steps let coordinates move more freely.When p > n, the authors note that preprocessing such as variable screening can be used to construct weights.Source: §A.2 · p. 16
• Hard-thresholding the gradient-descent estimate can remove false detections while retaining strong signals in the weak-signal simulation.The Lasso selected two false components, indexed 6 and 19; hard thresholding the proposed estimate with the same noise-level parameter removed those detections.Source: §A.2 · pp. 16–17

A.4 Extension to Elastic Net

• For the elastic-net reparametrized objective, all local minima are global and there are no local maxima.Source: §A.4 · pp. 17–18

A.3 Empirical Comparisons between different early stopping criteria

• The early-stopping comparison uses the S1–S4 strong-signal simulations, averages estimation errors over 50 trials, and contrasts oracle stopping with five-fold cross-validation.Source: §A.3 · p. 17, Figure 4 · p. 19

A.4 Extension to Elastic Net

• The Hadamard reparametrization extends implicit regularization to an elastic-net objective combining ℓ1 and ℓ2 penalties.The ℓ1 effect is represented in the limiting dynamics, while λ controls the trade-off between the two norms.Source: §A.4 · pp. 17, 19

A.5 Choosing the step size and initial value

• With sufficiently small step size η and initialization α, the optimal estimation rate is attainable with a properly tuned stopping time.The appendix suggests setting α = 1/p and choosing η using a small subset of the training data.Source: §A.5 · p. 19

A.7 Simulations on the effect of error variance on the optimal stopping region

• As σ decreases, the optimal stopping-time region expands, while its starting point remains at the same rate when initialization, condition numbers, and smallest signals are fixed.The simulations fix α and η and vary σ across 0.01, 0.1, 0.5, and 1.Source: §A.7 · p. 20

B.2 Some Useful Results

• Under the (s + 1, δ)-RIP, the Hadamard product (n^-1 X^T X u) ∘ v approximates u ∘ v in infinity norm for s-sparse u, with error at most δ||u||2||v||∞.Source: §B.2 · p. 22
• Under the (2, δ)-RIP, the same infinity-norm approximation holds for arbitrary u and v, with error bounded by δ||u||1||v||∞.Source: §B.2 · p. 22

C.1 Proof of Lemma 2.1

• The reparametrized least-squares objective has no local maxima, and every local minimum is global.Source: §C.1 · pp. 23–24

C.2 Proof of Lemma 2.2

• For twice continuously differentiable functions with only strict saddle points, the cited result states that gradient descent with random initialization and a sufficiently small constant step size almost surely converges to a local minimizer.Source: §C.2 · p. 24

C.3 Gradient descent error dynamics

• During burn-in, strong-signal components grow exponentially until reaching m/2, while weak-signal and error components remain bounded by O(p^-2).Proposition 1 gives this behavior with probability at least 1 − exp(−c1 log p) − exp(−c2 Ms) for t < T1.Source: §C.3 · p. 25
• After entering the good region, the iterates converge linearly toward β* up to statistical error εn and remain within an O(εn) neighborhood through the stated time interval.Source: §C.3 · p. 25
• The two stage propositions jointly identify the theoretical best-solution time interval corresponding to the stabilized region in Figure 1.Source: §C.3 · p. 25

C.4 Proof of Proposition 1

• Under the RIP and theorem assumptions, Proposition 3 controls the error dynamics while maintaining bounds on weak-signal and error components.The induction analysis carries these bounds forward for t < T1.Source: §C.4 · pp. 26, 30
• For strong signals, components move toward the true sign exponentially: positive signals rise and negative signals fall.The proof states this behavior when the iterate is still on the opposite side of half the true signal.Source: §C.4 · p. 28
• After at least Θ(log(m/α^2)/(ηm)) iterations, the proof guarantees each signal component is either within cα^4 of its truth or has magnitude at least half its true magnitude.This is the stated consequence of the signal-dynamics analysis.Source: §C.4 · p. 30

C.5 Proof of Proposition 2

• The Proposition 2 proof expands the squared estimation-error update and decomposes the remaining bounds around errors induced by strong signals.It uses the update formula and orthogonality between subvectors to organize the error terms.Source: §C.5 · pp. 31–32

• In the burn-in stage, the smallest strong-signal component grows exponentially until it reaches m/2.The proof places this stage within t ≤ Θ(log(√m/α)/(ηm)).Source: §C.6 · p. 33
• After βt enters ΘLR, convergence is geometric up to a high-order error term.Proposition 1 places entry into ΘLR within Θ(log(1/α)/(ηm)) iterations.Source: §C.6 · pp. 33–34

• The first experiment varies initialization across independent and correlated Gaussian designs, with n = 200, s = 4, and p ∈ {500, 2000, 10000}.The independent design satisfies RIP with high probability in the stated setting; the correlated design does not.Source: D Simulations Beyond RIP Conditions · pp. 34–35
• When RIP fails in the correlated design, the algorithm still recovers β˚ as α → 0, but more slowly.The authors attribute the slower convergence to a smaller allowable step size and a larger condition number of X, and leave the conjecture without RIP proof open.Source: D Simulations Beyond RIP Conditions · p. 35
• In the null-space-property-violating example, the algorithm converges to the least ℓ1-norm solution β: rather than the sparsest solution β˚.The estimation error decreases as the initialization level α decreases.Source: D Simulations Beyond RIP Conditions · p. 35

• Algorithm 2 achieves comparable median standardized RMSE to Algorithm 1 across 50 repeated simulations.Algorithm 2 uses MCP as the pilot estimator and gradient descent for post-debiasing.Source: E Simulation results for Algorithm 2 · p. 35, Table 4 · p. 36

F Point-to point Comparisons with Vaskevicius et al. (2019)

• The theoretical contributions include β = a2 − b2, separating strong and weak signals by recovery iterations, and adaptively improving the final error rate through early stopping.The paper uses β = g ∘ l as its formulation and describes β = a2 − b2 as a proof technique.Source: F Point-to point Comparisons with… · p. 36
• The target rate is s2 log p/n, reaching the parametric root-n rate when s2 = 0; the compared rate has slightly worse dependence on s1.The paper says its stronger conditions support this target, while Vaskevicius et al. (2019) developed a less stringent analysis.Source: F Point-to point Comparisons with… · p. 36
• The paper’s optimal stopping region has different rates for a and b; a smaller noise level expands the region through a larger b.The authors say this matches their supplementary simulations; the compared theorem gives a and b the same rate.Source: F Point-to point Comparisons with… · pp. 36–37
• The vanilla algorithm’s computational aspects are not elaborated, while the compared adaptive algorithm has O(pn) convergence iterations; preprocessing yields a revised algorithm with the same efficiency.The authors say the comparator’s time-varying component-specific step sizes alter the limiting dynamics, leaving convergence to the minimal ℓ1 solution unclear.Source: F Point-to point Comparisons with… · p. 37
• The methodology extends to coordinate-specific step sizes and elastic-net penalties, compares stopping criteria, and proposes two-stage debiasing to refine explicit estimators.The refinement is reported to boost convergence rate when there is no weak signal.Source: F Point-to point Comparisons with… · p. 37
• Simulations consider a relaxed restricted eigenvalue condition; without the null space property, the algorithm cannot reach the sparsest solution even in the noiseless case.The authors connect this observation to sparsity arising from the ℓ1 norm rather than an ℓ0-type method.Source: F Point-to point Comparisons with… · pp. 37–38
• In the noiseless regime, the analyzed non-convex objective has no local maxima and all local minima are global; the paper also studies variable selection, unlike the compared work.The paper reports simulations and real-data analysis of variable selection and empirical rules using cross-validated lasso λ.Source: F Point-to point Comparisons with… · p. 38

G Discussion

• The discussion identifies two open questions: whether the theory holds without RIP and whether data-driven complexity measures can select stopping and step-size parameters reliably.The paper characterizes RIP as stronger than the restricted eigenvalue condition it calls the minimal possible assumption in the literature.Source: G Discussion · p. 38