Sparse Polyak: an adaptive step size rule for high-dimensional M-estimation
Authors: Tianqi Qiao
Marie Maros
Venue: NeurIPS 2025 poster ()
Keywords: high-dimensional statistics
non-convex optimization
adaptive step size
Primary area: optimization

Abstract: We propose and study Sparse Polyak, a variant of Polyak’s adaptive step size, designed to solve high-dimensional statistical estimation problems where the problem dimension is allowed to grow much faster than the sample size. In such settings, the standard Polyak step size performs poorly, requiring an increasing number of iterations to achieve optimal statistical precision-even when, the problem remains well conditioned and/or the achievable precision itself does not degrade with problem size. We trace this limitation to a mismatch in how smoothness is measured: in high dimensions, it is no longer effective to estimate the Lipschitz smoothness constant. Instead, it is more appropriate to estimate the smoothness restricted to specific directions relevant to the problem (restricted Lipschitz smoothnes constant). Sparse Polyak overcomes this issue by modifying the step size to estimate the restricted Lipschitz smoothness constant. We support our approach with both theoretical analysis and numerical experiments, demonstrating its improved performance.

Sparse Polyak: an adaptive step size rule for high-dimensional M-estimation
Tianqi Qiao Texas A&M University College Station, TX, USA tianqi qiaoStamu.edu
Marie Maros Texas A&M University College Station, TX, USA
mmaros@tamu.ed edu
Abstract
We propose and study Sparse Polyak, variant of Polyak's adaptive step size, de- signed to solve high-dimensional statistical estimation problems where the problem dimension is allowed grow much faster than the sample size. In such settings, the standard Polyak step size performs poorly, requiring an increasing number itera- tions achieve optimal statistical precision-even when, the problem remains well conditioned and/or the achievable precision itself does not degrade with problem mismatch in how smoothness measured: in high
size. We trace this limitation dimensions, it is no longer effective to estimate the Lipschitz smoothness constant. Instead, it more appropriate estimate the smoothness restricted to specific di- rections relevant to the problem (restricted Lipschitz smoothness constant). Sparse Polyak overcomes this issue by modifying the step size to estimate the restricted Lipschitz smoothness constant. We support our approach with both theoretical analysis and numerical experiments, demonstrating improved performance.
Introduction
Consider the high-dimensional statistical estimation problem
min f(0) = R0301000
((z;, 0), 1=1
(1)
with data points 2, E Rd, - .... We focus on the regime in which the dimensionality grows much faster than the sample size, i.e. & co. To obtain consistent estimates in this setting, is necessary to assume that the true solution exhibits some low-dimensional structure-such as sparsity, In (I) sparsity is enforced through the to constraint, which guarantees that will have at most non-zero elements. This constraint renders the problem in (1) non-convex and, in general, NP-hard, regardless of the objective function f Nataarajan (995).. Nevertheless, under certain assumptions on the data, various algorithms have been developed to efficiently find approximate solutions to (I). Notably, under suitable assumptions that hold for a variety of statistical models, the Iterative Hard Thresholding (IHT) algorithm has been shown efficiently find sufficiently accurate solutions (1) The IHT algorithm results from applying projected gradient descent to T and reads
θ = HT, (0: - ((.))
where HT, denotes the hard thresholding operator. HT, retains the s largest-magnituce components of its input and sets the remaining to zero. Here, O denotes the step-size which ought to be chosen as 2/ 2/(3L) Jain et al. (2014), where L denotes the restricted Lipschitz smoothness (RSS) constant, and be interpreted the Lipschitz smoothness constant of f when restricted sparse directions.
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

For variety of statistical models, such as Generalized Linear Models (GLMs). L remains constant long as #log(d) remains constant, even if oc. This insight underpins the rate invariance of IHT: under suitable conditions, the number of iterations required to achieve (near) optimal statistical precision remains constant
both d and grow.
Analogous to the Lipschitz smoothness constant, the RSS constant must be estimated practice. Thus, the natural question in this context, and the starting point to the work in this paper is: (i) Do already existing approaches to adaptively tune Y via the estimation of the Lipschitz smoothness constant work in the high-dimennional setting? Our criteria to determine whether step-size rule works the high dimension, additional to convergence will be determined by the answers to the following questions: (ii) Can they achieve the same better guarantees than by choosing the optimal fixed step-size? (iii) Can they guarantee rate invariance 4 oc?
1.1
Related works
Over the past three decades, numer- ous methods have been proposed to solve the problem in (T) including
Logistic Regression
d 5000
Matching Pursuit Mallat and Zhang
(1993), Orthogonal Matching Pur- suit Pati et al. (1993), and CoSaMP Needell and Tropp (2009). Iterative Hard Thresholding (IHT) was first in- troduced in Blumensath and Davies (2009) with many variants proposed since. While initial convergence guar- antees for IHT required the Restricted Isometry Property (RIP)-a condition often too stringent in practice-Jar al. (2014) extended IHT's conver- gence guarantees to problems fulfill- ing the restricted strong convexity (RSC) and RSS; establishing linear convergence to near optimal statisti- cal precision. Observe that in the M- estimation context, convergence ar- bitrary precision is unnecessary and convergence to the best achievable sta- tistical precision is preferred. How- ever, for the results ja et al. (2014)
0.6
10000
20000
f(O)) f(0)
0.4
0.2
0
10
20
30
40
50
Iteration t
Figure Performance of Polyak's step size (dashed) and Sparse Polyak (solid) on logistic regression problems with increasing d and n. The quantities 8, s*, and log(d) remain constant. With Polyak's s step-size the performance degrades as increases whereas Sparse Polyak exhibits rate invari- ance, i.e. the number of iterations to achieve (near) optimal statistical precision does not change.
to hold, if the optimal parameter to recover is s sparse, O (R²s*) is required, where R is the restricted condition number. Khanna and Kyrillidis (2018) proposed an accelerated version of IHT which was extended to the stochastic setting in Zhou al. (2018) Zhou et al. (2018) establishes that with a sufficiently large mini- size, acceleration can be achieved and the faster rate requires only O(Rs") Axiotis and Sviridenko (2022) proposed variant of IHT with adaptively chosen weighted Ez penalty for which only s 2 O (Rs" is required. They further establish that for IHT > (Rss" is in fact necessary condition achieve near optimal statistical precision. Li et al. (2016) and Shen and Li (2018) introduced variants of IHTT that incorporate variance reduction. Shen and Li (201 7b) and Yuan al. (2018) propose Partial Hard Thresholding and Gradient Hard Thresholding pursuit respectively, with a focus on support recovery under high SNR assumptions.
Yuan and Li (202T) establishes generalization bounds for solutions found via IHT. Further, Zhang et al. (2025) establishes that IHT, for a range of K, finds solutions that can be shown to achieve the oracle estimation rate under a high SNR condition.
With no exception, the discussed works establishing convergence under the RSS require knowledge of the RSS constant L.. Knowledge of L. is crucial for the requirement > O (R²s*) being sufficient for convergence to optimal statistical precision. More generally, convergence can be established
with > S Consequently, overestimating L. forces smaller step size 7, which in turn slows down convergence and leads to denser solutions.
whenever 2 0 (20g2)

---

Several recent works have proposed adaptive step-size schemes based on estimating the local Lipschitz constant Malitsky and Mishchenko (2000 2024); Latafat et al. (2024). With the goal to jointly exploit the function's s local regularity and the algorithm's s trajectory, Mishkin et al. (2024) derive new convergence results for gradient descent as function of the objective function 's local Lipschitz smoothness and strong convexity parameters. They additionally establish that Polyak's step-size obtains fast and path-dependent rates. However, their results do not naturally extend to constrained problems, particularly non-convex ones like
Polyak's 's original step-size rule, first proposed in Polyak (1969), J, has gathered renewed attention in the machine learning community Ren et al. (2022) Hazan and Kakade (2019); Loizou et al. (202T); Wang et al. (2023); Zamani and Glineur (2024), but remains ill-suited for non- non-convex constrained settings such as ① unless additional assumptions are imposed. Some efforts have been made to adapt Polyak's step size to constrained problems. For instance, Cheng and (2012) addresses box constraints, while Devanathan and Boyd (2024) consider convex constraints.
To this day, except for Li al. (2024), and works that employ inexact line-search strategies Xiao and Zhang (2013); Wang et al. (2014). there has been no study of adaptive step-size schemes in the high-dimensional context. The work in Li al. 2024 handles (T) in the stochastic setting and proposes the use of Polyak's step-size with no modification. The results in Li et al. (2024) are limited even with no stochasticity, as they imply bounds on the restricted condition number, i.e. R RS n-1
where 17 (v5 1)2/4. Further, it is worth mentioning that the notion of RSS assumed in Li et al. (204)) is more restrictive than that assumed in the present paper and Dain et al. (2014). Further, as shown in Fig. [ Polyak's step-size with no modification presents performance that degrades as the size of the problem increases even if R and the optimal statistical precision (of the order of O(* log d) for this particular statistical model) remain constant. This effect is highly undesirable in the nigh-dimensional setting and as shown in the present work, can be avoided with a suitable modification of Polyak's 's step-size.
1.2
Major contributions
Our main contribution is the first adaptive step-size rule that performs well in high-dimensions and preserves the rate invariance property. To develop this scheme, we address question (i) posed in the introduction. We observe that adaptive step-size rules that estimate the Lipschitz smoothness constant do not necessarily work well in high-dimensions. This is because, in many common statistical models, the Lipschitz smoothness constant scales as (d) with high probability. Consequently, we answer questions (i) through (iii) in the negative, by demonstrating empirically (c.f. Fig. I that estimating the Lipschitz smoothness constant via Polyak's 's step-size (dashed line) does not yield rate invariance as d grows even if log(d) is held constant.
To overcome this limitation, we design an adaptive step-size rule that estimates the restricted Lipschitz smoothness constant instead. With this modification, we answer questions (ii) and (iii) in the affirmative. This is captured in Theorem [ and its Corollaries, which particularize the results of Theorem to relevant statistical models. In Sections 3 6 and in Appendix D we provide theoretical and empinical evidence, both on synthetic and real data, that our proposed method outperforms Polyak's 's step-size for high- ditennsinaii M-estimation tasks and achieves rate invariance. Moreover, we thcoretically and empirically show that Sparse Polyak converges to optimal statistical precision at least as fast IHT with the optimal fixed step-size 7 Ol) We also establish sufficient conditions under which we can guarantee support recovery, and particularize our results for specific statistical learning models in Section 4
Our guarantees are derived under standard assumptions, identical to those in Jain et al. (2014); [hhen and Li (b017b); Yuan et al. (2018). To guarantee convergence to optimal statistical precision we require the knowledge of f(0") where denotes the true parameter. While f(0") is known to be of the order O logd) its exact value is typically unknown. Consequently, we provide in Appendix B
double loop method that estimates surrogate to f(0") and allows convergence to optimal statistical precision as long as lower bound to f((") is known. Observe that in our context f f(0)) 2 0 making valid lower bound.
In addition, we prove linear convergence for statistical models that do not satisfy the regularity conditions in Jain et al. (2014); Axiotis and Sviridenko (2022), Yuan and Li (2021) but instead fulfill
3

---

the weaker condition in Loh and Wainwright (2015). This extends the results in Jain et al. (2014); Yuan and Li (202T) cover additional GLMs with an adaptive step-size rule. We provide these additional results in Appendix C.2
Our proof technique may be of independent interest it provides clear pathway establishing convergence of the IHT algorithm. We believe this is key in extending theoretical guarantees to adaptive step-size schemes solve ( We provide sketch of the proof of our main results in Section and the formal proof in Appendix
analysis applies only to Polyak's step-size, conjecture that other adaptive step-size rules such as Barzilai and Borwein (1988), /Zhou et al. (2025), Malitsky and Mishchenko (2020) may also suffer from similar performance degradation in high dimensions, and would therefore require analogous re-engineering to be effective.
Finally, although
Notation. Throughout this paper, we adopt the following notations. For vectors x, E Rd we denote by 2; the 1th element. We denote the & -norm of as || oc the Euclidean norm as || as ||æ||1,, and the fg-norm as || lo- Recall ||æ||o {a 2; A OH}- Also, we let | min denote the minimal entry of æ in the sense of absolute value. The inner product between two vectors denoted as (a, v). Matrices such as X E Rdxd, are capitalized. For any matrix E, we denote its largest singular value by maxx (X) and its smallest singular value by min (E). Similarly, if 2 E Rddd diagonalizes, we use Amax(2) (2:) and Amin (2:) to denote the largest and smallest eigenvalues of Σ respectively. The Frobenius norm of X is given by ||.X] Ir. and the nuclear norm by ||
, the & f--orrm
2 Setup and background
We make the following assumptions regarding the objective function f.
Assumption (RSC Agarwal et al. 20112). The objective function f is (p, (), restricce strongly convex in Rdd i.e.
10, 0211 _ _100 02l1 < ((0)) s(02) (VI((),)) 02), VO,, O2 E ER..
(2)
Assumption 2 (RSSAgarwal et al. (2012). ) The objective function f is -rettricee smooth in Rd, i.e.
S(0,) - s(02) - (V ((2)), 02)
2 1100 - 0211 - 110 02.13 V0,, 02 E eRd.
Assumptions 1 |and extend the classical notions of strong convexity and L--Lisschtz smoothness.
These assumptions reduce to their classical counterparts when T is sufficiently small and the direction
01 02 is appropriately sparse. Observe that when 0, 02 is dense and f is convex, Assumption becomes vacuous and the upper bound in Assumption 2 scales linearly with the problem dimension d. This highlights that the RSC and RSS depend not only the magnitude of the direction 0 0₂ but also its structure.
In high-dimensional statistical learning settings where + 0c, standard strong convexity and smoothness assumptions fail to hold. However many important problems still satisfy variants of the RSC and RSS, with both y and L. remaining
and with exhibiting only moderate dependence on d, e.g., Jain al. (2014); :): Agarwal et al. (2012). We leverage this in Section high-dimensional statistical learning problems.
4 to establish fast computational and near optimal statistical guarantees for variety of
We further highlight that our results can be generalized by adopting weaker RSC condition, where we assume that (2) holds only for pairs 0, H" satisfying || - 6* < where 0 denotes the ground truth, rather than requiring it to hold globally. This relaxation broadens the applicability of our approach, allowing to accommodate a wider class of functions. Notably, for some generalized linear models (GLMs), the loss function does not necessarily satisfy (2) without imposing additional constraints. We provide detailed discussion of these results Appendix C.2
3
Main Result
In this section we present our main theoretical result, which establishes convergence guarantees for Sparse Polyak (c.f. Algorithm l). In this section, we focus the deterministic setting of f. while

---

the statistical case will be addressed in Section 4 Note that IHT corresponds to projected gradient descent when applied to (T). The projection onto the Eg given by the hard thresholding operator, already discussed in the introduction and formally defined as follows.
Definition (Hard Thresholding Operator), For any > O, and z E IRd, we define (T,(2) as the projection of on Bs(s) = {0 IRd
|10||o S $}. i.e.,
HT, (z) argmin |0 =1|2, BEBa(s)
where ties are broken lexicographicaly:
Algorithm 1 Iterative Hard Thresholding (IHT) with Polyak Step-Size
1: Input: Function, f. target function value 7. sparsity parameter 8, number of iterations T
2: Initialize: 0g E IR: with ||0ollo I
3: for t= to 7 - 1 do
4: Compute step-size 7/2 max max(s(0)) ((0,) f
(V ((B.))]
5: Update: 0++1 = HT, (0.- 7/1 V v(()))
6: end for
7: Output: 8T
The adaptive step-size rule in Algorithm differs from the classical Polyak rule by replacing || [0,))] with HT, ( f(e)))l.. This approach contrasts with the low-dimensional case and with the work in Li et (2024). Observe that even if the current iterate is sparse, sparsity of f(..)
can not be guaranteed. In fact, the worst-case relationship // f(+))l S V (ITT.(V((0))).. which holds with equality for vector in which all coordinates are identical, may hold. As result, using the full gradient norm can lead to using overly conservative step-sizes, slowing down convergence dramatically as d increases. Unless strong additional conditions are imposed on the problem, convergence may be even be jeopardized (see discussion Section for more details).
Before providing our main result we introduce some notation. Consider fixed values 2 s* > 0,
and let / A mine. [()). Assume the chosen target function value satisfies f2 l" Let 3TS, and L/ denote the restricted condition number.
== L + 3TS, Theorem 1. Let 10th121 denote the sequence of iterates generated by Algorithm Suppose the objective function I satisfies the RSC and RSS in Assumptions' and respectively. Let Ô be any s* s'-spase vector such that f [()) f, and assume µ > and (240R) s" Then, for any iterate O
such that 0,0c 0112
36||HT, (V (V11))11 we can guarantee
a²
10 012 <
110, - 0112.
80R
Moreover, let to 2 0 be the first iteration for which 110 - 0|1 <
36 HT (V/
Then, for all
t 2 to. 10, 0|12 < +
1
36 (VTHHH
80R
Theorem implies linear convergence at a rate scaling with R up to precision O (V This result is equivalent to that in Theorem 3 in Jain al. (014)), where the RSS constant assumed to be known, up to a constant factor. Thus, we successfully answer in the affirmative question (ii) posed in Section Further, observe that, if L, j and |[H,, (V.T, (V)( remain constant as the size the problem grows, the rate remains unchanged and the achievable precision does so too. This implies that for variety of statistical models the rate and final precision will remain invariant as the problem size increases, as long as the aforementioned quantities do not change. Thus, we also answer in the affirmative question (iii) posed in Section I We particularize the result to specific statistical models and provide further discussion Section and Appendix C
f2
Here, / is a user-defined target value that reflects the desired level optimization, which can be set above equal to the statistical accuracy of the problem. Such relaxation is natural in learning
5

---

problems, as it is often counterproductive to optimize to full precision. The continuity and restricted strong convexity of f imply that s" sS-spasee vector 0 such that f(0) always exists.
The additional factor 80R in the final precision (c.f. Theorem 1 when t to) stems from the expansiveness of the hard thresholding operator. The removal of this factor and support recovery guarantees can be achieved under the signal to noise ratio (SNR) condition (c.f. and are provided in the following Corollary. Such condition is widely used in hard thresholding and support recovery studies, as seen (Shen and Li 2017a Prop 2) Bouchot al. (2016); [hhnn and 2017ch: Yuan et al. (0166).
Corollary 1. Under the assumptions stated in Theorem 7 if. further. the SNR condition:
(3)
holds. for any 2 to, where to is defined in Theorem the support of 0, contains that of 0, and the sequence ||0 - 0 is non-increasing and upper bounded by 36||ETT,[ (V)
Algorithm [ offers an approach obtain a target accuracy that we assume known in advance to (1) without requiring precise knowledge of L, μ4, T. For scenarios in which we only have access to lower bound on the problem, Algorithm 2 serves as an alternative: in most learning problems, the bound can be simply set to 0, making the method broadly applicable. This adaptive variant of Algorithm builds on the framework of Hazan and Kakade (2019), which reviews gradient descent with Polyak's s step size and its double-loop counterpart. By updating the lower bound adaptively in an outer loop, the method ensures that either the accuracy O [V (()))) is attained or the updated lower bound remains valid.
Algorithm 2 IHT with Adaptive Polyak
1: Input: Function f a lower bound fi, and sparsity parameter s. E Rd
2: Initialize:
3: for k=1 to K do
4: for 0
do
5:
Compute step-size 7/8 = T0] (V/(..
[(0,)-
6:
Update: 0.+1 HT, (6. 7e9c)
7:
end for
8:
O2 argmingsT f(0.)
JA+1 f(
9:
2
10:
Og- = 0k
11:
end for
12:
Output: 0 argmin f f(x)
Let 2 8 O. and define f = minloollsss f(0), attained by some s" sparse vector 6*
Theorem 2. Consider the iterates {0x) generated by Algorithu
2 Assume that the func-
(V when 36
tion f fulfills Assumptions and Then for (1 + 1606
j > 0, 2 (480R)>s*, Algorithm requires at most T = log, (()(0.)-1(0))) T gradi-
ent evaluations to achieve f(0) f(0) I and ||@ 1- < (1
H !! Here
160R
8² 100-0" T |o<(1/(-1//600)) log 36(1+1/660 HT.
This theorem focuses on the distance of the iterates to 0". The quantity 7 can be interpreted as the number of iterations required reach the desired accuracy when applying Algorithm with f* The additional term O [((w)-f(0") in the definition of T corresponds to the number of outer iterations needed to obtain sufficiently tight lower bound for the targeted accuracy. Similar order guarantees are established in Theorems and 2 The proof of Theorem is provided in Appendix B

---

4 Statistical Guarantees
The results in Section 3 are deterministic in nature and consequently, do not depend on a data generation model. In, contrast, in this section we use Theorem I provide guarantees for specific statistical models. for sparse logistic regression and low-rank matrix regression respectively. We provide guarantees for additional statistical models, including sparse linear regression, in Appendix
and
the camputational-ssatistte
performance guarantees
4.1 Logistic Regression
We consider dataset consisting of observations 2; (zi, 14) for I,..., where æ; E Rd denote the feature vectors, and yi E {0, 1) denote the corresponding responses. The feature vectors organized into the design matrix
x x (21,...,
Raxd
We assume that the relationship between yi and z; follows the model
1
Pr(ys 2.) xx (')*
(4)
where 6* is an s -spase vector representing the underlying ground truth parameter. The objective function is defined
1 71
f())
log (1 + exp(x?0)) y. 0.
We assume that each covariate vector 2; is drawn independently from a multivariate normal distribu- tion (0, Σ), where Σ is non-singular covariance matrix. By invoking Corollary 1 in Yuan and Li (2021), ) we can establish that the objective function f(0) satisfies the RSS and RSC conditions, as formalized in the following lemma.
Lemma 1. Consider the sparse linear logistic regression problem described above. Suppose the covariates 2, are uniformly bounded such that ||>,||| < 1 for all [n]. Then f(0) L-smooth with T, = 1. Moreover, with probability at least - e con, the RSC condition holds with curvature parameter f := cxp(--4R) Quin (E) and tolerance == cxp(-4R) ((E) logd where R ==
|10"1, ((E)= max;=1....dd Eis, and Co, C1 > 0 are universal constants.
Corollary 2. Consider the sparse linear logistic regression problem described above. Under the assumptions of Lemma further suppose that the sample size is sufficiently large so that jı > O. Further, assume the design matrix X E Baxd is normalized such that \X,/ ///l C for all
j = 1,...,d. dd Let {0:)220 be the sequence of iterates produced by Algorithm 7 when applied to the sparse logistic regression problem. Assume the sparsity parameter satisfies (240R)² st, and
=((0). Then, with probability least
con 2. the following hold:
(ff 110, - 112 2 720 nat the iterates exhibit contraction toward 0", i.e., 110+11 - 0-12 S
(1 s)x) 11. (ii) Let to denote the first iteration at which ||0., -0"12 2|2 72C log d ap²
Then for
all t 2 to. the iterates remain confined in a neighborhood of 0": 1|0.-0 < (1 SUL) 72 elogd
ng
(ii) 0 satisfies the SNR condition B then the iterates remain confined in neighborhood of l ||0 1- 2 < 72 nng for all to and the support 0f 0 is exactly recovered and preserved for all subsequent iterations.
The proof of Corollary is provided in Appendix A.3
) to provide performance guarantees of HT with fixed step-size on Logistic Regression. However, we note that this assumption extremely restrictive in the higg-diimesional setting. We provide additional results that do not require ||æa|| < 1, Vi E [n]] in Appendix apply to additional GLMs that do not satisfy the RSC condition (2) globally, and are aligned with the results in Loh and Wainwright (2015) both in terms of sumple and asymptotic convergence rates. In these cases, convergence at rate of (1 - co/ k) can be guaranteed if the algorithm is suitably initialized. However, a stricter requirement the sparsity level, specifically O(R')s".
Remark The assumption || S Vi [n] is required in Yuan and Li (022))
Our results in Appendix C.2 further
7

---

is necessary otherwise, and rate of (1 1/2 can be guaranteed. Here co, and 01 are universal
constants.
The result above and the result in Appendix C.2 match the observed behavior of Sparse Polyak in Fig Namely, we observe that Sparse Polyak indeed achieves rate invariance, as d 00, the rate and final precision remain constant long as slog is left unchanged.
4.2 Matrix Regression
Consider the data generation model
y/ (X;, (*) + Ex. for ,,2...,7e
where -- E RRdd is a matrix of rank at most s". We assume, X; E Raxd vec(X.) ~ NN E) i.i.d, and E 0. Further, 6; o2) are i.i.d. and independent of X,. Define s(0) 2n
En_1 U/ O)).
By invoking (Agarwal et al. 2012 Lemma 7), we establish the RSS and RSC properties of f()), as formalized in the tollowing lemma.
Lemma 2. Consider the low-rank matrix regression problem described above. Then, with probability at least 1 --c e -on, f(O) satisfies the RSS and RSC conditions with respect to the Frobenius norm and the nuclear norm. The corresponding parameters are given by:
max (2)), p Fmin (2), and T T= - d
where ((Σ) == suP| =1/[|||=1 Var(u X1v), and co, c 0 are universal constunts.
To enforce the suitable low-rank structure the iterates define
PM.(W) = Σ auun?
where o,, = 1..., are the s largest singular values of W, and u, the corresponding singu- lar vectors. We substitute all instances of the HT, operator by PM, in Algorithm [ and [ (c.f.
Appendix B
Corollary 3. Consider the low rank matrix regression problem described above. Let (Othezo be the sequence of iterates generated by Algorithm [hen applied low rank matrix regression problem. Suppose that f(@*), n is sufficiently large such that p > O, and s 2 (240R)²**. Then, with
probability least 1- -aan -1d- the following holds: (i) O 0* 10 7200,2 then
1dd
nu
the iterates contract relative to è* as ||0++1 -0 1/ < (1 80R
10 1 (ii) Let to be the
first iteration for which ||Ot, θ || 22000><(E)es Then, for all 2 to. the iterates remain in na
stable neighborhood around θ* with ||0: - -
SOR)
7200o² ((E)ed
np
The proof of Corollary 3 can be found Appendix A.4
Remark 2. The result in Corollary also extends Algorithm 2 f satisfies the SNR condition
3) the iterates of Algorithm 2 recover the support after T gradient descent steps.
Sketch of the Proof for Theorem and Corollary 1
We provide sketch of the proof of the main results of the paper. We refer the reader to the appendix for the complete proof. The proof of Theorem [ follows the outline: (i) study the behavior of
||A,-1 0|12 given "t and ||e, Ø|| (ii) establish that under the assumption that % is sufficiently large the expansive effect of the Hard Thresholding operator can be offset by the contractive effect of the gradient update, (iii) show that 7/. is sufficiently large until we reach optimal statistical precision.
Both to finalize the proof of Theorem [and to establish Corollary [ (iv) the iterates remain confined within a neighborhood of 0² and, given that the SNR condition holds, the support can identified, providing further benefits to the algorithm's performance. We elaborate on points (i) through (iv).

---

(i) To understand the dynamics of ||0; 0 we analyze the combined effect of gradient descent under the RSC and RSS and the Hard Thresholding operator.
For this we suitably apply Lemma from Jain et al. 2014), exploit the RSC, and the properties of 0.41 in relation HT, yield:
10 0
12
1 +
((1 jime) 110, 1-1
(5a)
2% (f(0.) f(0))) + 10 |IT.((
(5b)
Setting 7/- ((0)) makes 5b nonpositive and yields a choice that is invariant with d. in ' and the RSS are critical to avoiding a step-size that scales with
The use of (V the Lipschitz smoothness constant (which scales as O(d) in the high-dimensional setting). Setting Tt = 1/4LL) in (5) recovers the result in Jain et al. (2014), which establishes that IHT requires at most O log([/e)) to be within e accuracy of near optimal statistical precision.
(ii) To achieve the optimal linear rate under our choice of step size, we require that
[++
8
(1 (1- 1m) S (1) I
(6)
Sufficient conditions for (6) are 8 2 22x28" and 7t 2 - for some universal constants cp: et and c2- However, if y/ were to scale with the Lipschitz smoothness constant, i.e. with d we would require 2 d's" to establish linear convergence. Observe that this requirement can not be fulfilled as S d.
(iii) We exploit the RSS and RSC to show that ?/2 2 407 when ||0 0| > 3 1 V
If, further, 66 holds with co 1/160, yielding ||0++ - 012 < (1 - co 1 0|
(240R)²s.,
The condition |00 - 21 (6/lHT,((//001 stems from I being a constrained problem. If the
|0,-
Polyak step-size had been left unaltered, additional regularity conditions are required to establish convergence in the constrained case. As established Theorem in Polyak (1969), one such condition is (())) 2 for some > 0 and any 0, Rd. which does not uniformly in high- ||0-0]| dimensional M-estimation. For some GLMs, we can show the condition holds locally and exploit this fact to provide more general results in Appendix C.2
(iv) From (i)-(ii)) it follows that there exists to at which 110, 012 < 36 HT. ((0) 2
Since 3 we show that
the potential expansion is at most (1 + SOR Theorem follows. If 0 satisfies
the inequality |00t, 0| 36 (V)(0))) guarantees that S c Stg. From here, we establish
that, this results in possible scenarios: (a) 72 JU. implying $ Stot ensuring that 11 0 S... 2 = |l@ 2 and eliminating the expansion term in (5) : or, (b) 2 40Z and
(6) holds. Based on (a) and (b), Corollary D is established by induction.
Numerical experiments
We first consider sparse linear regression and sparse logistic regression on synthetic data. This done to illustrate the algorithm's 's performance as the size of the problem grows while the problem condi- tioning is kept the same. In all scenarios that rely on synthetic data, we set (5000, 10000, 20000). and s* 300. The design matrix X Rdxa is generated to reflect a time-series structure with correlation parameter w = 0.5. We set the sample size n according to = (as log df, where > 0 a constant. In this section, we use 5. For each column index j E {l...., n), we generate sequence of i.i.d. standard normal variables eI, and construct T1,j = E1/ / subsequent entries are generated recursively Ir+1.j = wit. for / E {1...., d 1}. where &, (0, 1). The true parameter 0* created by sampling each entry from ((0 1), and assigning nonzero values to s* randomly chosen entries.
W The

---

Linear Regression
Logistic Regression
Logistic Regression
f(0)) f(
76
0.6
7
0.6
100
Algorithm
Al Algorithm
0.25
0.4
0.4
0.3
50
0.2
0
0.2
40
60
Iteration
Iteration É
Iteration
Figure Left and center: IHT with 2 (blue) vs. Algorithm (red) on linear and logistic regression
3Z
respectively. Right:
Choice of f on Algorithm In all scenarios 5, d 5000 and 700.
In the case of linear regression, each sample (where i (I.....n)) is generated according to the (0,00.20) 0.255...
model:
= - 0" + wis
For logistic regression, the relationship between ys and z; follows the model 4
(i) Comparison to fixed step-size: In Figure present comparison between IHT with a fixed step size (blue) and the adaptive step size used in Algorithm ! (red) when solving linear regression (left panel) and logistic regression (center panel) respectively. When working with constant step-size,
2014) for both linear and logistic regression. L for linear set the step size to 2 following Jain et al.
3Z
regression can be upper bounded as Amax (3 2(2*+=")) (Loh and Wainwright 2015 Appendix D.I), whereas L for logistic regression is one fourth of that of linear regression. In both settings, Amax(() (E) S [---)-1+"] Agarwal et al. 2012). Although both step size strategies share the same theoretical guarantees,
aa
observe that the adaptive step size speeds-up convergence. This advantage arises because an adaptive step size can adapt to the local curvature and therefore be significantly more aggressive, allowing for faster progress towards the solution. As expected, this effect is more pronounced when solving logistic regression where the functions' curvature will depend on the point. (ii) Rate invariance: In Figure (c.f. Section [ we compare the classical Polyak step size with the step size employed in Algorithm across different problem dimensions while maintaining The solid line represents the periormance of Algorithm [ whereas the dashed line corresponds to the classical Polyak step size. This demonstrates that when the condition number of Σ remains unchanged, the complexity of the method remains almost identical under our chosen step size. In contrast, Polyak's 's step size leads an increased number of iterations, even if the achievable statistical precision and (Amax (E)/ Amin (E)) remain the same.
5.
(iii) Choice of f Finally, Figure 6 (right) highlights the impact of the choice of F. The results confirm that f determines the best achievable accuracy. Additionally, from the formulation of 7s, we observe that f directly influences the step size magnitude, thus is impacting the convergence rate.
Experiments real world data are provided in Appendix D All experiments were conducted laptop equipped with 16 GB of RAM and 12th Gen Intel Core i5-12500H 3.10 GHz CPU.
10

---

References
Alekh Agarwal, Sahand Negahban, and Martin J. Wainwright. Fast global convergence of gradient methods for high-dimensional statistical recovery. The Annals of Statistics, 40(5):2452 (085)2245-22482 2012.
Kyriakos Axiotis and Maxim Sviridenko. Iterative hard thresholding with adaptive regularization: Sparser solutions without sacrificing runtime. In International Conference on Machine Learning, pages 177--197. PMLR, 2022.
Jonathan Barzilai and Jonathan M Borwein. Two-point step size gradient methods. IMA journal of numerical analysis, 8(1): (1)::441-44, 1988.
Thomas Blumensath and Mike E Davies. Iterative hard thresholding for compressed sensing, Applied and computational harmonic analysis, 27(():265--77 2009.
Jean-Luc Bouchot, Simon Foucart, and Pawel Hitczenko. Hard thresholding pursuit algorithms:
Number of iterations. Applied and Computational Harmonic Analysis, 41(2):412-435, 2016. doi:
https://doi.org/100 10.111000111661.
David Chapman and Ajay Jain. Musk (Version 2) data set. UCI Machine Learning Repository, 1994. 142244/3/ 1608.
DOL: https://doi.orr//100 24432/C5
Wanyou Cheng and Donghui Li. An active set modified Molak-Ribiére-Polyyee method for large-scale nonlinear bound constrained optimization. Journal of Optimization Theory and Applications, 155: 1084-1094, 2012.
Nikhil Devanathan and Stephen Boyd. Polyak minorant method for convex optimization. Journal of Optimization Theory and Applications, pages 1-20, 2024.
Elad Hazan and Sham Kakade. Revisiting the Polyak step size. arXiv preprint 1005:319,. 2019.
Prateek Jain, Ambuj Tewari, and Purushottam Kar. On iterative hard thresholding methods for high-dimensional M-estimation. Advances in neural information processing systems, 27:685-693, 2014.
Rajiv Khanna and Anastasios Kyrillidis. IHT dies hard: Provable accelerated iterative hard threshold- ing. In International Conference on Artificial Intelligence and Statistics, pages 188-198. PMLR, 2018.
Puya Latafat, Andreas Themelis, Lorenzo Stella, and Panagiotis Patrinos. Adaptive proximal algorithms for convex optimization under local Lipschitz continuity of the gradient. Mathematical Programming, pages 1-39, 2024.
Changhao Li, Zhixin Ma, Dazhi Sun, Guoming Zhang. and Jinming Wen. Stochastic IHT with stochastic Polyak step-size for sparse signal recovery. IEEE Signal Processing Letters, 31:2035- 2039, 2024. doi: 10.1 19190999..222.3422555.
Xingguo Li, Tuo Zhao, Raman Arora, Han Liu, and Jarvis Haupt. Stochastic variance reduced optimization for nonconvex sparse learning. In International Conference on Machine Learning, pages 917-925. PMLR, 2016.
Po-Ling Loh and Martin J Wainwright. Regularized M-estimators with nonconvexity: Statistical and algorithmic theory for local optima. The Journal of Machine Learning Research, 16(1):559-615 2015.
Nicolas Loizou, Sharan Vaswani, Issam Hadj Laradji, and Simon Lacoste-Julien. Stochastic Polyak step-size for SGD: An adaptive learning rate for fast convergence. In International Conference on Artificial Intelligence and Statistics, pages 1306-1314, PMLR, 2021.
Yura Malitsky and Konstantin Mishchenko. Adaptive gradient descent without descent. In Proceed- ings of the 37th International Conference Machine Learning, pages 17, 2020.

---

Yura Malitsky and Konstantin Mishchenko. Adaptive proximal gradient method for convex optimiza- tion, Advances in Neural Information Processing Systems, 37:100670-100697, 2024.
Stéphane G Mallat and Zhifeng Zhang. Matching pursuits with time-frequency dictionaries. /EEE
Transactions on signal processing, 41(12)3397-3415, 1993.
Aaron Mishkin, Ahmed Khaled, Yuanhao Wang, Aaron Defazio, and Robert Gower. Directional smoothness and gradient methods: Convergence and adaptivity. Advances in Neural Information Processing Systems, 37:148 O--4848, 2024.
Balas Kausik Natarajan. Sparse approximate solutions to linear systems. SIAM journal computing, 24(2):22--233 1995.
Deanna Needell and Joel A Tropp. CoSaMP: Iterative signal recovery from incomplete and inaccurate samples. Applied and computational harmonic analysis, 26(3):301-32,, 2009.
Sahand N. Negahban, Pradeep Ravikumar, Martin J. Wainwright, and Bin Yu. A unified framework for digh-dimensional analysis of M-estimators with decomposable regularizers. Statistical Science, 27(4):538-555 557, 2012.
Mchdi Neshat, Bradley Alexander, Nataliia Sergiienko, and Markus Wagner. Large-scale Wave Energy Farm data set. UCI Machine Learning Repository, 2020. DOI: https//dbbi.rg/10024443///
Yagyensh Chandra Pati, Ramin Rezaiifar, and Perinkulam Sambamurthy Krishnaprasad. Orthogonal matching pursuit: Recursive function approximation with applications to wavelet decomposition. In Proceedings of 27th Asilomar conference on signals, systems and computers, pages 40 44. IEEE, 1993.
B.T. Polyak. Minimization of unsmooth functionals. USSR Computational Mathematics and http:////....0010011511000011010
Mathematical Physics, ((3):1429,, 1969. doi:
Tongzheng Ren, Fuheng Cui, Alexia Atsidakou, Sujay Sanghavi, and Nhat Ho. Towards statistical and computational complexities Polyak step size gradient descent. In International Conference on Artificial Intelligence and Statistics, pages 3930-3961. PMLR, 2022.
Jie Shen and Ping Li. Partial hard thresholding: Towards principled analysis of support recovery. In Advances in Neural Information Processing Systems, volume 30, 2017a
Jie Shen and Ping Li. Partial hard thresholding: Towards a principled analysis of support recovery. Advances in Neural Information Processing Systems, 30, 2017b.
Jie Shen and Ping Li. On the iteration complexity of support recovery via hard thresholding pursuit. In Proceedings of the 34th International Conference on Machine Learning, volume 70 of Proceedings of Machine Learning Research, pages 31 31-2-1124, 2017c
Jie Shen and Ping Li. tight bound of hard thresholding. Journal of Machine Learning Research, 18 1208):1-42, 2018.
Martin J. Wainwright. High- Bigh--imensional Statistics: Non-Asymptotic Viewpoint. Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 2019.
Xiaoyu Wang, Mikael Johansson, and Tong Zhang, Generalized Polyak step size for first order optimization with momentum. In International Conference on Machine Learning, pages 35836- 35863. PMLR, 2023.
Zhaoran Wang, Han Liu, and Tong Zhang. Optimal computational and statistical rates of convergence for sparse nonconvex learning problems. Annals of statistics, 42(6):2164, 2014.
Lin Xiao and Tong Zhang. A groximal-gradices homotopy method for the sparse least-squares problem. SIAM Journal on Optimization, 23(2):1062-1099, 2013.
Tao-Tong Yuan, Ping Li, and Tong Zhang. Exact recovery of hard thresholding pursuit. In Advances in Neural Information Processing Systems, volume 29, pages 3558-3566, 2016.
12

---

Xiao-Tong Yuan, Ping Li, and Tong Zhang, Gradient hard thresholding pursuit. Journal of Machine Learning Research, 18(166):1-4, 2018.
Xiaotong Yuan and Ping Li. Stability and risk bounds iterative hard thresholding. Proceedings of The 24th International Conference Artificial Intelligence and Statistics, volume 130 of Proceedings of Machine Learning Research, pages 1702-1710, 2021.
Moslem Zamani and François Glineur. Exact convergence rate of the subgradient method by using Polyak step size. arXiv preprint arXiv:2407.15195, 2024.
Yanhang Zhang. Zhifan Li, Shixiang Liu, Xueqin Wang, and Jianxin Yin. Rethinking hard threshold- ing pursuit: Full adaptation and sharp estimation. arXiv preprint arXiv:2501.02554, 2025.
Danqing Zhou, Shiqian Ma. and Junfeng Yang. AdaBB:
Adaptive MiriliiBooowein method for convex optimization. Mathematics of Operations Research, 2025.
Pan Zhou, Xiaotong Yuan, and Jiashi Feng. Efficient stochastic gradient hard thresholding. Advances in Neural Information Processing Systems, 31:1988-11999 1997, 2018.
13

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: Yes. The main contributions, including the novel theoretical analysis and as- sumptions, explicitly stated in the abstract and introduction, and are rigorously supported in the main sections.
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 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.
is fine to include aspirational goals as motivation as long as is clear that these goals not attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: The paper discusses the limitations of the proposed method, particularly its reliance on prior knowledge of the optimal objective value. This limitation clearly stated in the main text. To address this issue, a double-loop variant of the algorithm proposed in the appendix, which can be applied when the required information is unavailable.
Guidelines:
The answer NA means that the paper has 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 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 few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or 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 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, 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 impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed not penalize honesty concerning limitations.
3. Theory assumptions and proofs
14

---

Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
Justification: For each theoretical result, the paper explicitly states all required assumptions.
A proof sketch of the main theorem is included in the main text, while complete and rigorous proofs of all results are provided in the appendix.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Answer: [Yes]
Justification; The paper provides sufficient details to reproduce the main experimental results, including the data generation process, parameter choices, and implementation details. These are described in the numerical experiments section.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is a dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution a specific model and empirical evaluation, may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, other means that appropriate to the research performed.
also be provided via detailed
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how to reproduce that algorithm.
(b) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution a new model (e.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset instructions for how to construct the dataset).
15

---

(d)
The authors should provide scripts to reproduce all experimental results for the new reproducible, they We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [No]
Justification: The code is not provided, as the implementations involve standard optimization routines with no algorithmic novelty. All data used the experiments are either synthetically generated using clearly described procedures or drawn from publicly available sources. The paper includes sufficient details to allow reproduction of the results without requiring access to the code.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips.... cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (c.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https
/nips..ccppubic/guuu
for more details.
The authors should provide instructions on data access and preparation, including how access the raw data, preprocessed data, intermediate data, and generated data, etc.
proposed method and baselines. If only subset of experiments should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: The paper specifies all experimental settings. Full details of the synthetic data and standard optimization benchmarks are provided the numerical experiments section of the main body, while the settings for the real-data experiments are presented in the appendix.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core the paper to 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, as supplemental material.
7. 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?
16

---

Answer: [No]
Justification: This paper focuses on deterministic optimization problems, where the goal is evaluate algorithmic performance rather than statistical model accuracy. As such, the use of error bars or statistical significance tests is not applicable or necessary in this context.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results accompanied by emror bars, confi- dence 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 a library function, bootstrap, etc.)
The assumptions made should 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.
is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2- 2-sigma error bar than state that they have 96% CL, 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 range (e.g. negative eтror rates).
If error bars reported in tables plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Question:
Answer: [Yes]
Justification: This paper does not involve large-scale computation; the experiments are limited evaluating optimization performance moderate-size synthetic problems and standard benchmarks. All relevant information is provided in the Numerical Experiments section.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (c.g., preliminary or failed experiments that didn't 't make into the paper).
9. Code of ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics
Answer: [Yes]
Justification: Yes, the research conducted in the paper conforms to the NeurIPS Code of Ethics in every respect.
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
17

---

If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g.. there a special consid- eration due to laws or regulations 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: The focus of the paper is purely on theoretical advancements in optimization methods. The primary objective is to improve algorithmic performance and theoretical understanding, which does not have direct societal impacts.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative socictal impacts include potential malicious unintended uses (c.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 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 system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification:
The paper poses no such risks.
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have 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 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 best faith effort.
12. Licenses for existing assets
18

---

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: [Yes]
Justification: The paper uses open-source datasets, and proper credit is given the creators original owners of these datasets. The license and terms of use for these datasets are respected, as specified by the respective sources.
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 dataset.
The authors should state which version of the asset used and, if possible, include 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, paperswi thcode. com/ datasete has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification:
The paper does not release new assets.
Guidelines:
The answer NA means that the paper does release new assets.
Researchers should communicate the details the dataset/coddemmodll as part their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether how consent was obtained from people whose asset used.
At submission time, remember anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does involve crowdsourcing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According to the NeurlPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
19

---

15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/teview based the requirements of your country institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with
human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent)
may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
Declaration LLM usage
Does the paper describe the usage of LLMs if it is important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, formatting purposes and does not impact the core methodology, scientific rigorousness, or originality of the research, declaration is not required.
Question:
Answer: [NA]
Justification: The core method development this research does not involve LLMs as any important, original, or non-standard components.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, non-standard components.
Please refer to our LLM for what should or should not be described.
20

---

Appendix
Main Theorems
A.I
Proof of Theorem
A.2 Proof of Corollary
24
A.3 Proof of Corollary
25
A.4 Proof of Corollary
26
B
Adaptive Lower Bound
26
B.1 Proof of Theorem 2
26
C Other Statistical Guarantees
C.I Sparse Linear Regression..
27
C.2 Generalized Linear Models
29
Experiments on real data
D.I Linear Regression
D.2 Logistic Regression...
A Main Theorems
In this section provide the formal proof the statements included in the main body of the paper. To formally establish Theorem [ we build on (Jain et al. 2014, Lemma 1) and Lemmas 4 6 (Jain control the expansive properties of the Hard Thresholding operator, 4 allows us to establish (5) (c.f. Section (5) Further, in Lemma we establish instrumental in establishing lower bound the step-size "/s. We lower bound "r in Lemma 6 and establish conditions under which this lower bound holds. We combine these results in the proof of Theorem iin Appendix A.] followed by formal proof of the corollaries included in the main body of the paper in the remaining sections of Appendix
al. 2014 Lemma allows whereas Lemma
consequences of the RSS (Assumption 2 that
For simplicity, let S; U. Also, we let 9t f(0)), and g ( When discussing specific statistical models we denote by 0* the ground truth, g" = (f(0)) and S supp(0"). For any index set S and vector f E Rd we define [0]s as the vector that retains the entries indexed by S, while setting all other entries zero.
We include the following fundamental lemma, which plays a key role in our analysis. corresponds to Lemma in Jain et al. (2014), and is presented here for completeness.
Lemma For any index set I, any 2 RR, let HT,(2). Then for any Ø E Rd such that |10% I 8", we have
||0 2112
[ s² 110 2112.
Lemma 4. Let 0 be any -sparse vector, and 02 be any s-sparse vector. Assume that the function fulfills Assumption / with = μ 3rs Let 0.11 HT, (0, 7xtt), and 0.+1 0 -ie9t-
For any " < ((6.)-f(8) (T)l
have
2
2
|10
0|12 <
1
00+1
0
S.+1
1 +
g*
(1
110:
-012.
Proof. For the first inequality,
(it) Ø
10+t1 01 11+++1 (0+1ls...
||
110
[++
1
11
- Ø11
where in (i) the triangle inequality and in (ii) used Lemma
21

---

For the second inequality, we consider the expansion
[10++1]s... 01 [[0r+1 - 50)s... 112
=11 - -Os.... 11? 200 0.,91l.... + 111
where
-7 - -0,
s.t1 (..101...11
-2 - 00) + 2 (0. 0, 1 s. /,,11
and consequently
00+1 03 -| S 10, - 0112 -2 - 0,91 1911s,18. 7 [ 12
2
Using the RSC yields
|l0
8 112 <(1 - fin))10: - 2 ( + + (21HHT++0-(991))
+ 2% - 0, 190ls,18.11
To obtain (5) we must upper bound the inner product in 7). for which we have:
((0e)s.)S..11 i m ([0o)s.)8... 71 n1 ( (
([0rss//Be1 - 11 (.. + nlaes, ( (1lors,,55.1)
I |10:|s -ll..s. 111. ll I1al 13
Then,
(4)
(H) S |10, - no9l. /s II
|10, [10:->79iss,8,,,,
S |1 Vedds.
where in (i) that S. \ S4+1 S Se \ S1+1, and in (ii) we exploit that |S, \ S SS \ S. and that S.+1 contains the indexes of the s largest elements of 02 7hh Thus, we obtain the overall upper bound
2% (0, - 0, 1,,.,,1 (9ols,)....)) < 27 271llllsss I1llblls. Ill + 27
which together with (7) yields
110:+1-0s, - 0واg 21 < (1 mn) - - -
+ (9T))
(8)
Given the upperbound on ns the two right most terms together are negative, and we thus the proof is complete.
Observe that in the proof of Lemma we use the iterate 0.+1 to treat the effect of the hard thresholding operator and gradient descent separately. We then restrict to the support S.41 to avoid the scaling of any bound with the ambient dimension d.
Lemma 5. (RSS- gras-gradient bound) Assume that f fulfills Assumptions 7 Where jt = 11- - -33 > 0, and
Assumption 2
Then, for any pair x, y of s-sparse vectors there holds
HT, (Vf(x)
vf(x))
- (V)))() x))
2(L + 378)
Proof. We define
(()-
From its formulation, we know o(l) inherits the RSS and RSC property of f-
As result, for any 28-sparse vector 2, have
o(x)= S 6(z) (0, x)
2 r 2|
=(2)
マ || 2/12
á
< 6(2)
L + 3TS
s ((v) + (Vo(y), 3) +
2|
2
22

---

Set 2 = y L43ra
1 H (Vo(y)), the inequality above indicates
((x) S Φ 2(L + 3T8)
HT, (Vo()))
which is equivalent to
f(x) (VI(t)) x))
[HT,(V V(())) f(w).
2(I. 3T8)
Lemma 6. Consider the iterates {8: generated by Algorithm 7 to solve ① Assume that fulfills Assumptions where ā = 3rs > O, and Further, denote by an arbitrary s -ppare vector for which f() is known and desirable. Then, the step-size
72 f(0.) 5(ō) 1
40L
for each 2 for which
18 HT s+s (3)
10, 011 N
112
f(0.) f(0).
Proof. Assume that f(0.) - f()) > 0, then
f(0.) - f(0) 12 + HT. O
" HT (
Given that f fulfills Assumption we may invoke Lemma 5 yielding the bound
f(0.) f(0)
7t 2 20I(S(0.) ((0) 0. 02 0))
10 (T. Ø 0.) + ||HT,(0) < -
we can guarantee that
2 40
Rearranging
( we have that the condition can be equivalently written as
HT, 11 < 2L (S(O) - (()) + (.. - -))
Invoking the RSC, sufficient condition for the above is
(T) < 21 - 10,
012 + 2( θ))
which can be guaranteed long
|HT, (+)
2L
f 10:
0
2|| .+. (0) Ø1)
2ī 10:
1-
1
Ø| - 2 HT ( |10 Ø11
2
To guarantee that the above holds it is sufficient to request that ||0 0| 2 thus the result follows.
106
H [T, ( ||.
and
25p
23

---

A.1 Proof of Theorem
As a consequence of Lemma we distinguish three cases for any L: (i) 10, - 01 > 36
and
f(O:) f > 0, (ii) 10: 0
[ and f(()) f 0, (iii) f(0.) f s
r²
In case (iii) no progress is made, i.e. B441 and by the RSC there holds
#
00
01 110 - || s+s* (i)
2
and thus
00+1 012 = 100 012 36 ,
Further, from the above it follows that if ||0, 01
36 HOI we can guarantee that > 0.
For case (i) we begin by invoking Lemma
4 which guarantees that
||0+11 0 < 1 + 3
s"
(1 - (1--m)ll 10: - 01(.
Using Lemma 6 we guarantee a lower bound on the step-size Y1 40z Further, under assumption s, namely, 2 (240R)2s*, we can bound the contraction factor:
s* (1 (m)) S (1+ 807 (
1 80R
1+
(1 - mm)) -
1 +
S
40R
and therefore,
11
0||² 1
110, - -112.
S0R
Thus, the first part of the theorem's 's statement follows, i.e. when 10, - 0| is sufficiently large, we can guarantee that 0.+1 will approach 0. We now left with establishing the veracity of the second statement. For this, let ty be the time defined in the theorem's 's statement. Then, we are under case (ii) or case (iii). we are under case (iii) there is nothing left to proof. If we are in case (ii), there are two further cases: (a) the iterates remain confined within a ball of radius 6 HT 2 //p i.e.
Vt N to, |0 0 36 HT
J12
and the theorem's 's second statement is therefore true, (b) there exists time > at which for the first time
36 HT, 1
S 10, - -
By Lemma 4 and and by definition of t1 have
36 HT 1
36 36|HT HT (HTHT())) 20)
|00, - 01 <
1 + S0R
(1
(im))
< (1+
80R
This implies, we find ourselves again case (i). Observe that going forward, no iterate can escape the above ball and therefore the second statement of the theorem holds.
A.2 Proof of Corollary 1
If holds for 0, for any fulfilling
12
36 HT. H n
there holds
7 H
HT
min e
|00cal 2
> 0.
(10)

---

This can be established by contradiction, i.c. if the condition above is violated, either ||0, 0|1
36 [HT, (7)l
or 33 can not hold. Consequently, S C St. Now we consider the iterate 8(+1- Note that
122
Sc SCS++_1 if
|HT,
Th 10
(11)
2/1
To see this note that for any $ St we have
|TT
[[0,]: vloraal = >
2/
and for any E S, given that |10 01 ( have
7 7|H
(
6
(
(TIl
|HT,
[0.] [[
2
2j
2/z
Because the Hard Thresholding operator selects the largest components of 0t+1, in the selection of the elements that should go into S the operator will not deselect elements from S in benefit of any outside of St. Thus, S C St-1
We now find conditions on " under which (TT) can be guaranteed. Observe that
2e S 7/ (lğ) + 11 - ol
7h | + T1102 Ø1)
(0) - 7t
6L |HT, [IHT, (71) +
where (i) follows from Lemma
5 and (ii) follows from the assumption that |0, -0 0| <
36 (T))"
A sufficient condition for to hold is thus given by
N 2 127
As a result, when % < 12T we have S c st+1, and thus [|0,+1 - 0 < 100 0| by
|10-11
0| ||00+1 0,g,_, S mw)ll 0, 0| Otherwise, V/ 2 2p+
1
still
405.
have ||0++1 0 2 00 0| by Lemma Thus, the proof is completed by induction.
A.3 Proof of Corollary 2
By invoking Lemma I we can guarantee that fulfills the RSC and RSM with probability least 1 - The function is convex by construction, and assuming that Algorithm [ is provided suitable parameters f and stipulated by the Corollary, may invoke Theorem Thus, with probability at least e- -can, the iterates satisfy a contractive relation for ||0, 0*| until the point where 10: 8* 12 (T
p
To complete the proof, it remains to establish that, with probability at least 1 - 3. it holds that
36 (T" S 720 s log d
(12)
np²
The proof of this claim essentially follows the arguments in (Wainwright) 2019 Example 7.14); for completeness, we include the full details here.
Define (xx) TTexp(-z) The gradient of f evaluated at 0* can be expressed as
g" =
(a 0*) y1) Ti.
Recall the relation between 24 and ys is governed by the model
1
Pr(y:
x;)) 1 + exp - 0*)
25

---

Under this model, it follows that each term oz (t) is a zero-mean sub-Gaussian random variable with sub-Gaussian parameter o² =
Thus, 119* oo is the maximum of d independent zero-mean sub-Gaussian random variables, each with variance proxy most o²
C² In By standard sub-Gaussian maximal inequalities, we have
Pr l 2 2C
log d 2n
+
C6 2
- 2e
nd2
for any 8 > Setting
V 2 d n
completes the proof of the bound in (2)
Applying the union bound allows us claim that the guarantees provided by Theorem hold with probability at least oye Guarantees on support recovery follow then by direct application of Corollary [
A.4 Proof to Corollary
By (Jain et al. 2014 Lemma 2), we can verify the result of Theorem translates to the matrix case, where the vector 22-norm is replaced by the Frobenius norm, the -noorm replaced by the nuclear norm, and the HT, operator is substituted by PM,.
By applying Theorem I and Lemma
2 we conclude that, with probability least 1 the iterates
satisfy contractive relation for || ( "1 until the point where ||O, O* 36 Ill
/
To complete the proof, it remains establish that, with probability least 1 3.
36 36|IPM
7200o
(13)
n1
npp2
Note that g" n e,, by Wainwright 2019
Corollary 10.10), we have
[1g 2 An
S 2c 2no²
Pr
where A, = 100 (x))
(V
2d + +0)
By setting 8
2d
we have
2d
l1g 112 S 10 / (E)
with probability at least 2e 4dd We thus apply the union bound to complete the proof of claim
B Adaptive Lower Bound
B.1 Proof of Theorem
Let 0r = f(0,)- f(0") (T.)T Suppose the step size % used in Algorithm
satisfies Yr = bag for some
scalar E 1 1|. By invoking Lemmas and 6 as long as |10: - 0-1 > HT (o"") we have:
b 110. 0-1
1
110 1- <
110, 0-1
160%
40
160
This establishes that the iterates exhibit contractive behavior until they enter small neighborhood of the optimum.
We now consider two possible cases depending on whether the lower bound surrogate fs is valid and how the step size compares to a, during epoch k.
26

---

Case (i): Suppose fRr is valid lower bound for f(d"), and that " I holds for all iterations
t -....,, (()) within epoch k. Then, the contractive relation applies repeatedly, and we obtain
P
T(a)
100 0-1 - (1 a)) (1
36 (9"(1")11
(0r(a)-00)) 0-
160L
160%
J12
By the restricted smoothness and strong convexity assumptions, this implies that the function suboptimality satisfies f(Or(a)) I" s E(x),, thus completing the proof for this case.
Case (ii): Alternatiely, suppose that f& a valid lower bound, but 7 > ag for some t epoch
This condition implies that
f(0,) - Fi > 2(f(0.) f(0))),
which in turn yields
J&41 f(0x) - Jn ((0)) + < f(0").
2
2
Hence, fa+1 is also valid lower bound. By induction, we conclude that if Case I never occurs, then all Jú. for k 1,...y K, remain valid lower bounds.
Moreover, under this scenario, the sequence f decreases geometrically. In particular,
FRtt S fr
fk
- fk
2
The geometric decrease of f" - Fu. in case (ii) ensures that if case (i) never occurs, then there exists some kn such that f(0*) - fh., < =(a). In that case, 8R, is either an output corresponding to case (i) (which completes the proof) or an output under case (ii), i.e., satisfies
1(0%)) f(00) < 5(0)) Jro < <((a).
C Other Statistical Guarantees
In this section we provide guarantees for additional statistical models not provided in the body of the paper. The results in Appendix Jopendix[C.Ihholll Chll for sparse linear regression. The results in Appendix C.2 apply to some GLMs and require the analysis of the behavior of Sparse Polyak under different regularity conditions. These
entirely captured in Theorem T
C.1 Sparse Linear Regression
In this section, we assume the dataset consists of data points = (I,, 31) for i = I,..., n, where I; E Rd denote the feature vectors, and ys E IR denote the responses. The feature vectors are aggregated into the design matrix
Rad X
and the responses are aggregated in IR" 3 A (yn,..., .ya) Let 0" denote the ground truth of the statistical model, with ||0||o S s", and f denote the corresponding objective value. Specifically, we assume that the responses and feature vectors x; are related by y: = x f²* + where x; are drawn from (0, 2) distribution, Σ is non singular, E, N(0, d²), and x; and s; are i.i.d and independent of one another. Additionally, the objective function f(0) 2n X - y/l².
From Agarwal al. (Lemaa 6), it follows that the RSS and RSC conditions hold with f0 with coefficients L = 20 max (2), where ((E) = max;=1..... Eg. Here op and C₁ are universal constants.
probability at least
J (2), and ci((E) logd
Corollary 4. Consider the sparse linear regression problem described above. Let (0) be the sequence of iterates generated by Algorithm 1 or Algorithm 2 when employed to solve sparse linear
27

---

regression problem. Suppose that = for Algorithm Assume we have sufficient samples for T O, with 2 (240R)ss* for Algorithm [ and 2 (480R) s* for Algorithm 2 Further, assume that each column of X is C-normalized, i.e., column of X. Then, with probability at least 1 e con 3. for any a Z 80 after T(a) iterations, mingsT (a) 10, 0* 1 is upper bounded by
X p I C for j = 1,...,d. d. Here X, denotes the j-th
E(a) (1+aa 28
n/722
The required number of iterations ([a)ffillll
1 log (1/(1 - 1/80R))
|100 110-00112
T(a) <
log
.nd
=(()
T(a)= < 1og2 4(f(0)) (E(aa
f(0)))
1
|00 1 112 E(()
log2
log
log (1/(1 (/1100))))
for Algorithm Hadd Algorithm
respectively.
Moreover, if satisfies the SNR condition 33, we can ensure that after T(0) iterations, the error is upper bounded by e(0), and the support of 60 has been identified, i.e. S C St 4I 2 T(0).
Corollary 4 establishes the convergence properties of Algorithm D and Algorithm 2
The error term is
of order O R²*" log'' df
which of the same order as that in (Jain et al. 2014 Theorem 3), where fixed step size considered under the assumption that L is known.
nift
Proof. For Algorithm [ the proof follows the same steps the proof of Corollary difference is the upper bound for (0"|HT.("")Il" which we provide next.
The only
When the columns of X are C-normalized, by (Wainwright 2019 Example 7.14), with probability 1 2-
log d
l|g 113. | I 8C
Using the union bound together with (Agarwal et al. 2012 Lemma 6), and the assumptions stated in the Corollary, yield the required assumptions for Theorem to hold. Thus, this completes the proof for Algorithm T
As we can see from the proof of Theorem I and Theorem 2 the accuracy level =(0) is determined by the point at which lower bound the step size can be established. According to the formulation of ?,. in Algorithm 2 it is guaranteed to at least half the step size used Algorithm [ This observation implies that the accuracy level e(O) also be achieved by Algorithm conditions of Case (i) in the proof of Theorem 2 are satisfied.
2 when the
To complete the proof for Algorithm
we only need to show that f(0)) f
18 H8 II implies
P
||0, 0-1|² < s0)).
By RSC,
[(0.) l 2 pil00 0-
- 0-1|-
When f(0.) f <
18||HT, ( it implies
HT (g))( Il.
(3) HTTT (T'('"')ll
+ 110 +100--001)
2 110 11.
J
A necessary condition for the inequality above is
36 H6 11 12
||: -
E00).
2
28

---

C.2 Generalized Linear Models
In this section, we consider dataset consisting observations 2; = (I ys) for = 1,.... n, where Jy R&d denote the feature vectors, and yi E IR denotes the corresponding responses. The feature vectors are organized into the design matrix
#1
X
Raxd
and the responses collected in the vector y (y) Un) T E RR.
For notational convenience, we denote the true underlying parameter of the statistical model and f the corresponding objective function value. We assume the relationship between 2; and yi is characterized by the conditional distribution
y1u 0- [1 0*)
Pr(s I4, "7, =exp
c(o)
where a > 0 a scale parameter, and &i is the cumulant function. Given this data generation model, we define the objective function
((0) = (w(a you
We assume that & is infinitely differentiable with (t) > 0 and uniformly bounded for all R. These assumptions are satisfied in a variety of settings, including logistic regression and multinomial regression Loh and Wainwright (2015). We assume the feature vectors I; are i.i.d. and drawn from multivariate normal distribution N70, Σ), where Σ is non-singular. Under the setting described above, the RSC condition does hold everywhere. However, the described setting it can be shown that the following milder RSC condition holds Loh and Wainwright (2015)
" 2
yy
x
x12
T bl xm|
if fly - x112 S 1, (14a)
2
f(v) - f()) (Vf()) (V))))
T |l - x
|1y 2112
if f - 2112 > 1. (14b)
2 2 ly x0133
We present lemma grouping the results that need to proceed: the RSS condition, the condition (14), and the order of achievable statistical accuracy for GLMs in our setting. This lemma aggregates results from (Loh and Wainwright, 2015 Proof of Corollary 2, Appendix D.1) and (Negahban et al.) 2012, Proposition 2). Notably, while (Negahban et al. 2012 Proposition 2) is originally stated centered only around the ground truth 8" (c.f. (4) its proof extends any given x. This implies that while our results are currently stated for = 0 to achieve optimal
statistical precision, we can instead state equivalent results to those in Theorem where 0 is such that ʃ(0) f.
Lemma 7. For the statistical models described above, with probability at least cid -1 o2e Assumption and 14) hold, and
log
IV Sco
where C2, > are universal constants. The constants μ. and L in (14) and Assumption 2
respectively, depend on v and Σ. Further = C3 log where C3 > a universal constans.
In the following, we keep the definition T= L + 3T8, and j 378, where is that of 14
Observe that Corollary stated for Algorithm for simplicity but an analogous statement Algorithm
f follows.
29

---

Let (22 denote the iterates generated by Algorithm 7 when applied to the general- ized linear models described above, Set the step size rule according
Corollary 5.
mmax [ [.0]
5 (g:
Define == 1|0-11, and Ro == 4R + 1. Assume we set f = f 0o 0, and suppose 2
(4s'. Further, assume the sample size is large enough to ensure PL > 5co V logg
Then, with probability at least
a oze " we guarantee that for all 2 T where
n s slog(d)) ()log()))
))
T O klog
+ +0 (o=2
there holds
110 1-2 < 1 1
3600s log
160
Moreover, if satisfies the SNR condition (3) we can guarantee that for all t 2 T, 100 - 00112 s 36cps logd and that the support of 6* has been recovered, i.e. S C..
Proof. To establish this result we leverage that the condition (14 combined with Sparse Polyak allow us to establish convergence despite the lack of RSC. More specifically, we will establish that Algorithm exhibits, under the conditions described the Corollary, most three modes of convergence:
110: - 11,,
if 110 6 2 1,
160
||0,+1 0-12 S
) 110
8- 11?
if 11: - 0*1 < 1, and 110 - 0-112 36c s log
1ooR
np2
360 np²
otherwise.
R60
(15)
Observe that under our current assumptions 36c log
18 and therefore, the list above is exhaustive nAa2
25
and the conditions on the second case are compatible.
We begin assuming that R 2 1. We will establish that if the modes of convergence provided above hold for R 2 1, when R < 1 we only observe the two last cases.
(i) We start off by exploiting [1). Assuming that ||0, 0* 2 may follow the strategy in Lemma 4 and exploit that f is convex (but not RSC) to obtain
||0+11 0-112 2
1
s*
(110, 1- - -2)((f(0)) - f(())) + 5 (16)
s
With the choice of step size
max(f((.) f(0)),0}
7t
5 26
16
simplifies to:
2 (100--0112 0-111
-(00))
1 - 00- 5 1
V s
2 [1
(17)
From Lemma and under the assumption that µt 50u
2a logd it follows that with probability
least 1 - d
c2e -n there holds
f(0,) f(0)) > 10 - 8" |TT2 - 0- 2 10: 11
30

---

thus by applying the above bound onto (17 we have
- 0010 - 1 +
8*
10: - 1
10: - 0²(3)
(18)
80 HT (gt))) (gt)l
Further, we may upper bound the norm of the gradient as
(H:))(9))11
< 2 (| (g ) + HT 2x (9:: 9" )²)
(19)
x 2 (I + +223110. 1-
where in (i) we invoke Lemma
2 5
Combining
with 19 yields
n12
|0 11 v
100 1=1 (20)
160 (1 (
212 13)
To guarantee that the above implies our first regime of convergence we need to establish that ||0 112 < R for all t 2 O. We return to this point after exploring the second and third regime which will be useful establishing the first.
(ii) On the other hand, assume instead that ||0, 112 < 1 and |0 - 101| 2 3Ge logd then, from ap²
14 the RSC holds and therefore we can invoke the result of Theorem in which
1
|0c+1 0-
110, 8-11
(21)
160R
(iil) If, instead |00 -0 0-112=
10-12 <
38c log ni
we have from Theorem that
36c log
||0++1 - 11 1+
< 1.
(22)
160R
nip2
Observe that from the three cases we consider ((0))
(2
21 and 22 already correspond to one
of our stated modes of convergence, and thus are to establish the first.
Clearly, when < 1, Oo is in a region which the RSC holds, and therefore, we will only observe the behavior in 21 and (22). Thus, we only to prove the first regime of convergence for R > 1. To establish that the behavior in the first regime holds, we need to establish that |6, - 0 < R holds for all 2 O. We proceed to establish this and consequently the behavior in the first regime by induction. Note that by our initial condition ||0, - 1² = =0011 R and thus the condition holds for t 0. Suppose ||0: 1 R for some t. If [|6, 11
then either (2T) and (22) hold
and the proof by induction is complete. If instead, 00 - 112 2 1, (20) holds. Then, by induction
hypothesis and under the Corollary's 's assumption the sample size there holds
HT (g))) + 2T, 1-112 < 16 + 2R 1 S RgL2
thus, combining with 20 and using that by assumption
> RA8 obtain
(1
110-1-00 1-113
100 0" 112 R.
R60
We have thus established the veracity imply that there exists to fulfilling
15 Observe that 15 together with 0, = 0 and 101122 =
to ----log (R) log(1 1/160RpR²)]
such that for all t to
||0: - 0012 1.
Further, this implies that in most
log 36cfs(1 1/(160R) log(d)
log(1 1/160R)
2
additional iterations optimal statistical precision is reached. Finally, the SNR condition holds, the term (1 + 160% in the third regime is replaced by I as a consequence of Corollary D
31

---

Theorem 3) with the following similarities and differences. To achieve optimal statistical precision, as a log 0 both Sparse Polyak and (Loh and Wainwright 2015, Theorem 3) require ( ( )) iterations. However, we require additional iterations O log( when under more general conditions, we do not make the assumption that ||0*| < 1 which is necessary in (Loh and Wainwright/ 2013 Theorem 3), where this condition can relaxed the expense of requiring the RSC to hold within larger radius.
Corollary 5 recovers the result in (Loh and Wainwright 2015
1. We observe however, that our result holds
We conclude this Appendix by highlighting the fact that for both sparse linear regression and sparse GLMs be verify the rate invariance of Sparse Polyak theoretically. When the problem size increases much faster than the sample size
00 but > log and Σ remain constant, IHT with Sparse Polyak will reach E neighborhood of the optimal statistical precision within at most O log(1/s)) for + O - log(R)) when R > in the case of GLMs. In both cases, this number does not change with increasing and n. Observe that these results allow us to answer in the affirmative questions (ii) and (iii) posed in Section
linear regression and at most ( log(=/log//////
D Experiments on real data
Linear Regression We consider a linear regression task using the Large-scale Wave Energy Farm dataset from the UCI Machine Learning Repository Neshat et al. (2020), which is publicly available under the CC BY 4.0 license. The terms of use are described #terms The goal is to predict the total power output of the wave farm based on sparse linear model. We randomly select 120 samples from the dataset, each containing 149 features. In our experiment, we set the sparsity level to 20. For the IHT method with fixed step size, we choose the step size as 10-12 This value determined via grid search over the range [10-13, 10 10-12), as step sizes outside this interval result in poor convergence or divergence. The results are presented in Figure (left).
ics. uci edu/
Logistic Regression We evaluate sparse logistic regression using the Molecule Musk dataset Chap man and Jain 1994) from the UCT Machine Learning Repository, which is publicly available under the BC 4.0 license. The terms of use described at pattts://archhn
ics. uci edu/ #terms The task is to classify molccules as musks or non-musks. We randomly select 120 samples from the dataset, each with 166 features. In our experiment, we set the sparsity level to s 20. For the IHT method with a fixed step size, we select a step size of 1.9 x 10-5, chosen via grid search over the interval [3 10-6 4 x 10 Step sizes outside this range lead to poor convergence or divergence. The results are shown in Figure (right).
We observe that in Fig 3 (right) Sparse Polyak performs better than both classic Polyak and IHT with the fixed siep-size, even if the step-size is optimized by grid search. This is expected, an adaptive step-size can adapt to the curvature any point in the algorithm's 's trajectory, whereas a fixed step-size cannot. We observe that Sparse Polyak performs better than Classic Polyak. On the other hand, in Fig 3 (left) observe that the best fixed step-size performs better than adaptive step-size. We conjecture that this is due to the factor 1/5 in Sparse Polyak which we presume is an artifact of our analysis and is not in fact necessary. Despite this fact, sparse Polyak consistently outperforms the classical Polyak rule in the sigh--imensional setting.
32

---

Lincar regression
Logistic regression
-1012
Tt fixed
fixed
0.6
Algorithm
Algorithm
Polyak
Polyak
0.5
(1ea) - (o8)
f(0e) 1(0
0.4
0.3
2
0.2
20
30
40
0
50
100
150
200
Iteration
Iteration
Figure 3: Performance comparison of IHT with optimal constant step size, Sparse Polyak and classical Polyak when performing: (left) linear regression on the Wave Energy Farm data set, and (right) logistic regression the Molecule Musk data set.
33


Reviews:
[{"review_id": "3BhPCqb7Ob", "reviewer": "Reviewer_LwtJ", "summary": "The problem studies minimizing the average of convex losses under the constraint of having at most $s$ non-zero element. The Iterative Hard Thresholding algorithm has been found to a some desirable property called rate invariance: under commonly used statistical models, and as long as the number of samples scales as $s \\log d$ (instead of $d$ which is the usual scaling in the dense setting), the number of steps required for achieving optimal statistical precision is independent of $d$ and $n$. This however requires a step size that depends on the restricted Lipschitz smoothness (RSS) constant which is in general unknown to the algorithm. The question then is if existing methods with adaptive step size already achieve this kind of guarantee.\n\nThe paper shows empirically that using Polyak’s step size leads to an error that grows with the dimension, thereby ruling out the possibility that the rate invariance property holds in this setting.\n\nInstead, the proposed method  tries to estimate the restricted Lipschitz smoothness. The resulting algorithm is termed sparse Polyak, as it can be viewed as projected gradient descent with a modified version of Polyak’s step size, where the gradient in the denominator is replaced by its hard-thresholded version. The main theorem shows that under the standard assumptions on $f$ used in prior work, convergence is similar to when the RSS constant is known and that the rate invariance property holds. The corollaries specialize the results for specific statistical learning models. The corollaries include logistic regression and matrix regression. The authors also relax some of the conditions in order to get results for generalized linear models for which the original conditions do not hold. The paper also contains numerical experiments that show an advantage of the sparse Polyak method.", "questions": "- Does the algorithm work for any possible initialization $\\theta_0$ being an $s$-sparse vector?\n- If one wants to achieve $f(\\theta_t) - \\hat{f}$ being upper bounded by a desired value, how do they choose $t$ (this is related to the question about the main theorem in the \"Weaknesses\" paragraph)? Similar question for the Corollaries. \n- For Lemma 1, $\\|x_i\\| \\leq 1$ implicitly means that $\\Sigma$ must be upper bounded. Is there a way to relax this condition and possibly have the norm of $\\Sigma$ appear in more places in the lemma guarantee?", "limitations": "yes", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\n- The paper makes progress on an important problem. Taking advantage of sparsity to avoid performance degrading by the dimension is highly desirable.\n-  The main result is general and applicable to many specific problems, as shown in Section 4.\n- I liked the evidence shown in Figure 6 of the experiments showing that the adaptive step-size makes more - aggressive progress because it uses local information.\n\nWeaknesses (mostly about writing and presentation):\n-  Is question (i) from the beginning of page 2 intended to be a more general question that subsumes (ii) and (iii)? If not, what does (i) ask that is not already covered by (ii) and (iii)?\n-  On my first read, I was confused by the description of problem (1) as “statistical,” given that no references to the underlying data distribution were made early on. It would help to clarify from the beginning which quantities (e.g., the RSS parameter) depend on the distribution of the $z_i$ and in what way (e.g., in expectation vs. with high probability). In the current version of the paper the first time these aspects are talked about is in the paper is in Section 4.\n- It is stated that the main theorem is almost equivalent to that of Jain et al. 2014. However the theorem there is stated in a different way (it talks about closeness in $f$ values rather than closeness of the parameters). It would be helpful to explain how to translate between the two and to clarify what the additive term $1/(80\\kappa)$ corresponds to in the guarantee from Jain et al.\n- In Lemma 1, $\\hat{f} = f(\\theta*)$ is input to the algorithm thus it is implicitly assumed that $f(\\theta*)$ is known.\n\n\nOverall I found some aspects of the presentation confusing but I am overall satisfied with the results of this paper. That said, I am not familiar with the literature, so my opinion is based mostly on the main body of the paper and a quick look at some of the cited works.", "was_revised": true, "final_justification": "The authors have provided clear explanations to my points. After reviewing the other rebuttals and discussions my view of the paper remains positive, I thus keep my original score.", "rebuttal": "We thank the reviewer for their insightful comments and their overall positive assessment of our work. We are glad the reviewer liked and it and answer the reviewer's concerns in order:\n\n$\\textbf{Weakness 1:}$  We agree with the reviewer in that the questions can be misleading. Questions (ii) and (iii) are subsumed by question (i) in the sense that they concretize in which sense Question (i) is answered in the negative. We will address this in the final version.\n\n$\\textbf{Weakness 2:}$ We agree with the reviewer that the distinction between deterministic results and statistical results can be made clearer. To clarify this, we will highlight prior to discussing our results in Section 3 that they are deterministic in nature. We will direct the reader to Section 4 for the statistical setting.\n\n$\\textbf{Weakness 3:}$ As the reviewer correctly points out, the deterministic result in Jain et. al. 2014 (Theorem 1) is indeed provided in terms of the objective function value. We highlight that in general, results in terms of estimation error are preferred. To obtain, based on Theorem 1 in Jain et. al. 2014, a result regarding the estimation error (which is preferred in the statistical context) one can exploit the RSC (where here $\\bar{\\mu}$ includes the terms $\\mu,\\,\\tau,\\,s$):\n$$\\frac{\\bar{\\mu}}{2}||\\theta_t - \\theta^{\\star}||^2 + f(\\theta^{\\star}) + \\langle \\nabla f(\\theta^{\\star}), \\theta_t - \\theta^{\\star} \\rangle \\leq f(\\theta_t)\\quad (1)$$  \nyielding by exploitation of sparsity of $\\theta_t$ and $\\theta_{\\star}$ and Cauchy-Schwarz followed by Young's inequality on the inner product in $(1)$, a result of the form in Theorem 3 in Jain et. al. 2014. Conversely, we can use the RSM to connect the estimation error with the objective function gap\n$$f(\\theta_t) - f(\\theta^{\\star}) \\leq \\langle \\nabla f(\\theta^{\\star}),\\theta_t - \\theta^{\\star}\\rangle + \\frac{\\bar{L}}{2}||\\theta^t - \\theta^{\\star}||^2$$\nfollowing an analogous approach as for $(1).$ \n\nWe agree that at face value the addition of the $1/(80\\bar{\\kappa})$ factor multiplying $\\frac{6||HT_s(\\nabla f(\\hat{\\theta}))||^2}{\\bar{\\mu}^2}$ may appear odd (final precision, Theorem 1), especially given that we can guarantee that the estimation error is contractive as long as its larger than $\\frac{6||HT_s(\\nabla f (\\hat{\\theta}))||^2}{\\bar{\\mu}^2},$ and given that an analogous to this does not appear in Theorem 3 of Jain et. al. 2014. This apparent oddity stems from the adaptive step-size rule. We can guarantee that the adaptive step-size is sufficiently large when we are far away from statistical optimality. This allows to guarantee that the estimation error remains contractive even in the presence of the hard thresholding operator (which can be expansive if the support has not been identified). However, as we get close to statistical precision, we can not guarantee that the step-size remains sufficiently large. This implies that the potential expansive effects of the projection onto the non-convex set can not be controlled by the gradient step, and we may exit the ball of radius $\\frac{6||HT_s(\\nabla f(\\hat{\\theta}))||^2}{\\bar{\\mu}^2}$ by a factor at most $1/(80\\bar{\\kappa}).$ If this happens, we contract back and thus the guarantee provided holds. This discussion is partially reflected in (iv) in the proof sketch (Section 5) and its formal counterpart can be found in the proof of Theorem 1 in Appendix A.1. Observe also that as a consequence, this term does not appear when support recovery is possible due to sufficiently high SNR. This is because if the support does not change, the HT operator does not introduce expansion.\n\n$\\textbf{Weakness 4:}$ We specify the results for the statistical models under the assumption that $f(\\theta^{\\star})$ is known. We provide the analogous result to Theorem 1 when $f(\\theta^{\\star})$ is unknown in Theorem 2 (Appendix B). Here, we do not require knowledge of $f(\\theta^{\\star})$ but a lower bound to it. For statistical learning problems, 0 is generally a valid lower bound, i.e. $0 \\leq f(\\theta^{\\star}).$ In this case, optimal statistical precision is achievable at the price of computing additional $\\log(1/\\varepsilon)$ factor in terms of number of iterations. The analogous results to those in Section 4 in the case where $f(\\theta^{\\star})$ is unknown can then be derived as Corollaries to Theorem 2 in the same way as the results in Section 4 are corollaries of Theorem 1.\n\n$\\textbf{Question 1:}$ The reviewer is correct in that the algorithm will converge as long as $\\theta_0$ is s-sparse. One may select $\\theta_{-1}$ to be an arbitrary vector and set $\\theta_0$ to include only its $s$ largest elements (in absolute value) and zero elsewhere,  or alternatively, and more conveniently (unless prior knowledge on where the solution may lie is) one may set $\\theta_0 = 0.$\n\n$\\textbf{Question 2:}$  This is a good point and we acknowledge that setting a practical termination criteria is non trivial. If an upper bound to the quantity $||\\nabla f(\\hat{\\theta})||^2_{\\infty}$ is known (in terms of order this is not unreasonable, i.e. for $\\hat{\\theta} = \\theta^{\\star}$ $||\\nabla f(\\hat{\\theta})||^2_{\\infty} \\leq \\mathcal{O}(\\frac{\\log d}{N})$) and therefore from $(1)$ it follows that\n\n$$||\\theta^t - \\theta^{\\star}||^2 \\leq \\frac{1}{2\\bar{\\mu}}\\left(f(\\theta^t) - f(\\theta^{\\star})\\right) + \\frac{s}{2\\bar{\\mu}}||\\nabla f(\\theta^{\\star})||^2_{\\infty}.$$\n\nThus, setting $\\varepsilon = c s ||\\nabla f(\\theta^{\\star}) ||^2_{\\infty}$ (which is known up to a constant factor $c$) and requesting $f(\\theta_t)-f(\\hat{\\theta}) \\leq \\varepsilon $ guarantees an upper bound on the estimation error of the order of optimal statistical precision. We note that the total number of iterations this will require is given by Theorem 1. Observe that this explanation holds for all statistical models as they are all based on the same deterministic result.\n\nWe want to highlight that one can not establish that the objective function value decreases monotonically as in standard gradient descent when the step-size is chosen to be sufficiently small. This is because Polyak's step size can be larger than that yielding faster convergence in practice but a loss of monotonicity. Further, in the standard (low-dimensional) scenario, with strongly convex functions and Lipschitz smooth gradients, if we were to employ Polyak's step-size, the termination criteria $f(\\theta_t) - f^{\\star} \\leq \\varepsilon$ guarantees $\\frac{\\mu}{2}||\\theta_t - \\theta^{\\star}||^2 \\leq \\varepsilon,$ where $\\mu$ is unknown. This implies that even in the standard case, we have a sense of the precision at the level of the iterates except for a constant factor.\n\n\n$\\textbf{Question 3: }$ This an important point. We agree that the assumption $\\|x_i\\|\\leq 1 $ is extremely restrictive. In fact, this will in general not hold in the high dimensional setting.  We included this result as is the one present for logistic regression among works considering the Hard Thresholding algorithm and thus facilitates comparison to prior works. This restriction originates from Corollary 1 of Yuan and Li 2021, where it was used to ensure that the RSC condition holds. This is necessary because GLMs generally do not satisfy the RSC without additional structural constraints. In contrast, our analysis relies on a milder version of the RSC condition (Loh and Wainwright 2015) allowing us to avoid such restrictive assumptions. This and guarantees for our algorithm's convergence are stated and proven in Appendix C.2, as a function of the covariance matrix $\\Sigma$. We highlight that this is the first IHT (regardless of step-size) result to work for logistic regression without artificial assumptions on the statistical model. This is because Polyak's step-size can guarantee linear convergence under milder conditions on the lower bound on the function's curvature.\n\nPlease let us know if further clarifications are needed."}, {"review_id": "L56TpY4F6H", "reviewer": "Reviewer_XXZH", "summary": "This paper proposes a sparse adaptation of the classical Polyak stepsize rule for high-dimensional sparse estimation problems. Specifically, the classical Polyak stepsize is given by $\\gamma_t=\\frac{f(\\theta_t)-f(\\theta^\\ast)}{\\|\\nabla f(\\theta_t)\\|^2}$, and the proposed modification replaces the denominator with the squared norm of the top-$s$ hard-thresholded gradient, i.e., $\\gamma_t=\\frac{f(\\theta_t)-f(\\theta^\\ast)}{\\|\\mathrm{HT}_s(\\nabla f(\\theta_t)) \\|^2}$. This modification is intended to better align the stepsize with the underlying sparsity structure of the target solution. The paper provides a convergence theorem showing that, under this modified stepsize rule, gradient descent converges linearly to a neighborhood of the optimum. The authors also present numerical experiments to examine empirically the proposed method.", "questions": "- Could the authors clarify where the novelty in the convergence proof lies?\n- Have the authors considered comparing the proposed method to other sparse optimization methods, such as top-$k$ SGD?\n- In practice, how should $s$ be chosen, especially when $s^*$ is unknown? Is there a data-driven or adaptive way to tune this parameter?\n- On line 39, “the same” appears to be repeated twice.\n- I found the sentence on lines 67–69 a bit unclear. Could the authors consider rephrasing or elaborating on its intended meaning?", "limitations": "The work is theoretical and methodological in nature and does not appear to raise direct societal concerns.", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\n- The motivation for adapting the Polyak stepsize to sparse settings is clearly explained, and the method is easy to implement.\n- Using the top-$s$ hard-thresholded gradient in the stepsize rule is a smart and intuitive adaptation that aligns well with the underlying sparse structure of the problem. It allows the stepsize to focus on informative coordinates and avoids distortion from noise in irrelevant dimensions.\n- The authors prove linear convergence for their algorithm. This provides theoretical support for the proposed method in high-dimensional settings.\n\nWeaknesses:\n- The paper does not clearly define the original Polyak stepsize rule. For readers who may not be familiar with the literature, it would be helpful to include the classical definition early on to better contextualize the proposed modification.\n- While the sparse adaptation is well-motivated, the linear convergence result is somewhat expected from the form of the update. The theoretical result, though important, may not be surprising to readers familiar with sparse optimization theory. While the authors mention that the proof technique may itself be of interest, it is not clear where the novelty lies. It would strengthen the paper if the authors\ncould more clearly highlight any innovative aspects of the analysis.\n- It would be helpful to include comparisons with other sparse-aware baselines, such as top-$k$ SGD, to better contextualize the proposed method's benefits.", "was_revised": true, "final_justification": "The proposal to use an adaptive stepsize, while intuitive, is nonetheless innovative and well-motivated for the sparse optimization setting. The authors provide theoretical convergence guarantees that support their proposed algorithm. Notably, the paper is among the first to emphasize the need for modifying adaptive stepsize rules to address the challenges posed by high-dimensional statistical learning problems, where naive adaptations may degrade as the dimension increases. The authors have also adequately addressed my concerns regarding the lack of comparison to other sparsity-aware methods. Overall, I believe the merits of this paper outweigh the reasons to reject, and I have thus accordingly raised my score.", "rebuttal": "We thank the reviewer for their insightful comments. We will address the reviewer's concerns in the order they were brought up. However, for both clarity and convenience, we merge the responses to Weakness 2 with Question 1 and Weakness 3 with Question 2.\n\n$\\textbf{Weakness 1:}$ Thank you for pointing this out. We agree that including the classical definition of Polyak's stepsize would improve clarity. We will include this in the revision.\n\n$\\textbf{Weakness 2 and Question 1:}$ While the sparse adaption we propose is relatively simple in form, and the resulting convergence may appear intuitive, to the best of our knowledge, this is the first adaptive step-size for IHT  that achieves such a result. Not only this, but to our knowledge we are the first to highlight that when dealing with high-dimensional statistical learning problems and employing projected gradient descent, adaptive step-size rules need to be modified to avoid degradation as the ambient dimension $d$ grow. In terms of proof technique, the novelty of our contribution lies in the following:\n\n(i) Our analysis decomposes the dynamics of IHT iterations into the expansiveness of the HT operator and the contractive effect of the gradient descent step. Compared to the original analysis in Jain et. al. 2014 our approach improves clarity and provides a clear understanding of how convergence of IHT is achieved despite not having a fixed step-size.\n\n(ii) We establish that the use of Polyak's step-size is sufficiently large throughout most of the objective function's landscape to guarantee fast convergence, after this stage is completed, we establish that the expansiveness introduced by projection onto the non-convex $\\ell_0$ does not harm the overall performance. This aspect is made obvious in the proof of Theorem 1 (Appendix A.1) and discussed briefly in point (iv) in the proof sketch (Section 5). We want to highlight here that Polyak's step-size is typically not used in combination with projected gradient descent over non-convex sets and find that the use of it here yielding optimal statistical guarantees was not necessarily obvious a priori. \n\n$\\textbf{Weakness 3 and Question 2:}$ We thank the reviewer for the suggestion. We would like to clarify that a comparison to the classical Polyak stepsize and our sparse variant is already provided in Figure 1, where we demonstrate the invariance of the proposed method.\n\nWe agree that comparing to other methods could further contextualize the method's performance. However, our focus in this work is to highlight the need for the design of adaptive step-size rules for the sparse high-dimensional setting. For this, we focused on analyzing a step-size strategy for IHT, which targets an $\\ell_0$ constrained optimization problem. In contrast, methods such as top-$k$ SGD address a different class of problems, where the setup, data availability/gradient oracle, and problem geometry differ from those considered in our work. Most importantly, our focus is on problems with explicit sparsity constraints on the solution, which top$-k$ SGD does not enforce. Instead top$-k$ SGD selects the larges $k$ coordinates (in absolute value) with the purpose of approximating the negative gradient direction (this is interesting in the distributed context, as can allow to communicate less). Top-$k$ SGD will therefore be well suited for problems in which the gradient is well approximated by its truncated counterpart, and is thus in a sense, approximately sparse. Observe that in our case, the gradients can be dense (and consequently not well approximated by a sparse quantity) and it's the iterates that are forced to be sparse. We further want to highlight the fact that not controlling for the sparsity of the iterates (as done by top$-k$ SGD) will force the algorithm to traverse the landscape in regions with large Lipschitz smoothness constants potentially implying extremely slow convergence rates. For these reasons, we believe that top$-k$ SGD and related sparsified/quantized gradient methods are not appropriate baselines for direct comparison in our setting. Accordingly, we restricted our empirical evaluations to algorithms using IHT. That said, we acknowledge the importance of comparisons to a broader class of algorithms/formulations and will consider alternative step-size rules (that should be tailed for the high-dimensional regime) in future work. \n\n$\\textbf{Question 3:}$ In practice the parameter $s$ is to be chosen by for instance cross validation. We acknowledge that this is not a desirable feature, but it is almost universal in the sparse high-dimensional setting. Either one needs to have a lower bound for $s,$ or knowledge of the noise variance $\\sigma^2$ (penalty case). We relax the requirement of IHT of having access to both $s$ and $\\bar{L}$ to only $s.$ We want to highlight that this implies that we will be searching for a single parameter as opposed to two as when performing standard IHT.\n\n$\\textbf{Question 4:}$ Thank you for pointing this out. We will correct the repetition on line 39 in the revised version.\n\n$\\textbf{Question 5:}$ Thank you for mentioning this. The intended meaning of these lines was to convey that the results in Jain et. al. 2014 require the iterate sparsity level to satisfy $s \\geq \\mathcal{O}(s^{\\star}\\bar{\\kappa}^2).$ This is highlighted here as it plays an important role in highlighting why an adaptive step-size rule is desirable. The scaling $\\mathcal{O}(s^{*}\\bar{\\kappa}^2)$ is only guaranteed in Jain et. al. 2014 if the step size is chosen $\\mathcal{O}(\\frac{1}{\\bar{L}}).$ If the step-size is smaller, the sparsity parameter degrades, implying that the final solution will be unnecessarily dense. Our scheme avoids this issue altogether as the step-size will always be sufficiently large to guarantee that the final solution is sparse in the same sense as in Jain et. al. 2014. We will revise lines 67-80 to clarify this point and enhance readability.\n\nPlease let us know if further clarifications are needed."}, {"review_id": "oS4cXuKrJZ", "reviewer": "Reviewer_drVv", "summary": "The paper considers a high-dimensional statistical estimation problem $\\min_{\\theta} f(\\theta) = \\frac{1}{n} \\sum_{i=1}^n \\ell(z_i,\\theta)$ where sparsity of $ \\theta $ is at most $s$. They focus on the case when the dimensionality grows much faster than the sample size, i.e. $d/n \\rightarrow \\infty$. The problem is non convex even if $\\ell$ is strongly convex/smooth due to the sparsity constraint. IHT is an algorithm that has been known to efficiently find solutions to this. The algorithm repeatedly applies projected gradient descent where the projection is via a hard thresholding operator. The stepsize used is $1/\\overline{L}$ where $\\overline{L}$ is the restricted Lipschitz constant, i.e. smoothness constant corresponding to certain important directions.\n\nThe paper looks at 3 questions: 1. Do already existing approaches to tune the step size work? 2. Can they achieve the same of better guarantess by choosing an optimal step size? 3. Can they guarantee rate invariance as $d/n \\rightarrow \\infty$?\n\nThe main contribution is an adaptive step size rule that performs well in high-dimensions and maintaines rate invariance. They show simulations that prove existing approaches do not work well in high-dimensions.\n\nKey idea & algorithm\n1. Measure smoothness only along sparse directions -- The authors replace the full gradient norm in Polyak’s rule with the norm of the gradient hard-thresholded to the allowed sparsity level s, yielding an estimate of the restricted Lipschitz constant.\n2. Sparse Polyak step size  followed by one projected-gradient (Iterative Hard Thresholding) step .\n\nTheoretical guarantees\n1. Linear convergence to a neighbourhood whose radius matches the optimal statistical precision, at a rate $O(1/\\kappa)$ where $\\kappa $is the restricted condition number (Theorem 1) .\n2. Rate invariance: if $\\kappa$ and the achievable precision stay constant while $ d/n \\rightarrow \\infty$, the required number of iterations remains constant—matching fixed-step IHT with oracle knowledge of the restricted smoothness constant .", "questions": "1. Your theory assumes access to a lower bound on $f(\\theta^*)$. In practice you estimate it once per outer loop. How sensitive are convergence and statistical accuracy to misspecifying that bound?\n2. Convergence is linear with rate $1/\\kappa$. the restricted condition number. Have you observed cases where $\\overline{\\kappa}$ is large even though sparsity holds, and does Sparse Polyak still help relative to classical Polyak there?", "limitations": "Yes", "rating": 4, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths: \n\n 1. First adaptive step-size rule that retains IHT’s constant-iteration property in the high-dimensional regime.\n2. General proof technique that may extend to other adaptive rules.", "was_revised": true, "final_justification": "Thanks for the rebuttal and for answering my questions. Based on other reviews, responses, and prior NeurIPS theory papers, I would like to retain my score.", "rebuttal": "We thank the reviewer for their insightful comments and their overall positive assessment of our work. In what follows we do our best to address the reviewer's questions.\n\n$\\textbf{Question 1:}$ Thank you for the question. As noted, our theory assumes access to a lower bound on $f(\\theta^{\\star}).$ In practice, this may not be known, but 0 is typically a lower bound for $f(\\theta^{\\star})$ for many learning problems. In Appendix B, we describe a double-loop scheme that relies on such a lower bound, and provide theoretical guarantees (Theorem 2) on both computational cost and accuracy, together with its corresponding proof. Corollaries for specific statistical models can be derived based on Theorem 2 in the same way as the Corollaries in Section 4 are derived of Theorem 1. \n\nImportantly, the final accuracy achieved is not sensitive to the tightness of the initial lower bound. Regardless of how loose the estimate is, the same level of accuracy can still be attained, at the cost of an additional $\\mathcal{O}\\left( \\log ((f(\\hat{\\theta}) - \\tilde{f}_0)/\\varepsilon) \\right).$ Here $\\hat{\\theta}$ is a solution to the $s^{\\star}$ constrained minimization problem and $\\tilde{f}_0$ denotes the lower bound (in general, this can be set to zero).\n\n$\\textbf{Question 2:}$ Yes, the reviewer is correct in that convergence is linear and the rate scales with $\\bar{\\kappa}.$ We agree that when $\\bar{\\kappa}$ is larger the convergence rate will be slower. In such cases, Sparse Polyak still offers advantages over classical Polyak, particularly due to its sparsity-aware scaling. Specifically, within our set-up the objective function is not strongly convex, and the Lipschitz smoothness parameter (non-restricted) scales with $\\mathcal{O}(d/N)$ implying, according to existing theory, potentially sub-linear rates, that will degrade as $d/N$ increases. In short, Sparse Polyak will offer an advantage whenever $\\bar{\\kappa} < \\bar{\\kappa}'$ (observe that this always holds in our set-up)  where $\\bar{\\kappa}'$ denotes the non-restricted condition number, will offer an advantage when the objective function is not strongly convex in all directions, and will yield near (due to the constant factor in our step-size) identical performance when $\\bar{\\kappa} = \\kappa'.$\n\nPlease let us know if further clarifications are needed."}, {"review_id": "9TpKc1uGQO", "reviewer": "Reviewer_yL2Y", "summary": "This paper proposes Sparse Polyak, a novel adaptive step-size rule for the Iterative Hard Thresholding algorithm, aimed at solving sparse high-dimensional M-estimation problems, particularly where the dimension $d$ significantly exceeds the sample size $n$. Traditional Polyak step sizes degrade in such settings due to poor estimation of the global Lipschitz constant. The proposed method instead estimates a restricted Lipschitz smoothness constant by focusing on significant descent directions relevant to the sparsity constraint. The authors present a thorough theoretical analysis that establishes the linear convergence of the proposed algorithm under restricted strong convexity and smoothness assumptions, a result comparable to that of IHT with an optimal fixed step size. Notably, the convergence rate only depends on the restricted condition number, which does not grow with dimension in many specific statistical problems, enabling rate invariance under high-dimensional settings. Exact recovery of the support set is established under additional SNR conditions. Experiments are conducted on both synthetic and real datasets to support the theoretical findings.", "questions": "1.  Can the authors discuss whether partial or approximate support recovery guarantees could be obtained under relaxed assumptions (e.g., $\\Vert \\theta^\\star_{\\min} \\Vert$ below the SNR threshold)?\n2. Can the authors comment on whether Sparse Polyak could be effective in dense or moderately sparse settings? Will it degrade significantly?\n3. How to choose an appropriate $\\hat{f}$ in practice? When the loss function is always nonnegative, can we always set it to zero and control the iteration complexity only by stopping criteria?", "limitations": "Yes", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 4, "significance": 4, "originality": 3, "strengths_and_weaknesses": "Strengths:\n1. The paper provides nontrivial extensions to the convergence theory of IHT with Polyak step-size in the high-dimensional regime. The main theorem and its corollaries are rigorously proved.\n2. The experiments are carefully designed and demonstrate the algorithm's advantage over the original Polyak step size, consistency with the optimal fixed step size, and rate invariance under high dimensions, thereby providing sufficient support for the theory.\n3. The paper connects its results to classical models, such as sparse logistic regression and low-rank matrix regression, demonstrating broad applicability.\n\nWeaknesses:\n1. Corollary 1 requires a relatively strong minimum signal condition for exact support recovery. In many practical applications, signal amplitudes are not well-separated.\n2. The paper assumes that each per-sample loss function is convex, which ensures that the empirical risk is convex. However, this assumption excludes a large class of modern M-estimation problems involving non-convex losses (e.g., Tukey’s biweight, logistic loss with non-convex regularizers, or robust statistics). In such settings, the optimization problem becomes non-convex both due to the objective and the constraint, and it is unclear whether the proposed Sparse Polyak step-size remains effective or stable. This limitation restricts the scope of the method to convex loss settings, and the authors do not discuss potential extensions or failure modes under non-convex losses, which are increasingly relevant in robust learning and signal processing applications.", "was_revised": true, "final_justification": "I have carefully reviewed the authors’ replies, and I find that they have adequately addressed all of my comments and questions. I do not have any further concerns. I would like to retain my original score.", "rebuttal": "We thank the reviewer for their insightful comments and overall positive evaluation of our work. We address the reviewer's concerns in order but merge Weakness 1 and Question 1 for clarity.\n\n$\\textbf{Weakness 1 and Question 1:}$ We thank the reviewer for their feedback and question. The minimum signal condition in Corollary 1 is a standard assumption in the literature and has been used prior to claim support identification for IHT-like algorithms. We agree that in some practical scenarios where signal amplitudes are not well-separated, exact support recovery may be more challenging. Upon revision of our proof, we observe that partial support recovery guarantees naturally follow from our analysis. Specifically, any coordinate whose magnitude exceeds the threshold on the right hand side of equation (3) can be correctly recovered, even if the minimum signal condition for exact support recovery is not satisfied. We will clarify this observation in the revision and believe it offers a useful perspective on how our method performs in practice when ideal conditions are not met.\n\n$\\textbf{Weakness 2:}$ The reviewer correctly notes that our current analysis assumes convexity of the empirical risk. We agree that extending the proposed Sparse Polyak stepsize to non-convex loss settings--such as robust M-estimation or non-convex regularization (other than $\\ell_0$)--is an important and challenging direction. A natural extension would be to consider non-convex problems that satisfy the PL condition. Results for the (standard) Polayk step-size for the PL exist, and are sufficiently fast to guarantee that the expansiveness of the HT operator can be absorbed. For non-convex problems with less structure, sublinear rates are the best one will potentially be able to expect, complicating the absorption of the effect of the HT operator. This will require a more careful analysis of the amount of coordinates that enter/exit the support per iteration. These are two promising directions that we leave as future work.\n\n$\\textbf{Question 2:}$ The performance of Sparse Polyak will improve upon that of classical Polyak whenever the restricted condition number is smaller than the condition number (observe that smaller or equal is always true). When restricted strong convexity does not hold, we expect sublinear convergence from both schemes with Sparse Polyak benefiting as the Restricted Lipschitz smoothness constant will be smaller than the Lipschitz smoothness constant. In the case where no sparsity is to be included in the problem (namely gradient descent) both step-size rules collapse and thus yielding the same guarantees but with a slight degradation due to the constant factor in our step size rule.\n\n$\\textbf{Question 3:}$ Thank you for the question. While it may seem appealing to use 0 to substitute the value $\\widehat{f},$ this approach generally does not guarantee optimal accuracy. The reason is that doing so leads to a more aggressive stepsize which may disrupt convergence altogether. To alleviate the downside of requiring knowledge of $\\widehat{f},$ in Appendix B, we desribe a double-loop scheme that adaptively updates a lower bound for $f(\\theta^{\\star}),$ which substitutes the unknown $\\widehat{f}.$ This approach guarantees the same level of accuracy, with an additional computation overhead of at most $\\mathcal{O}\\left(  \\log((\\widehat{f} - \\widetilde{f}_0)/\\varepsilon)\\right)$ outer iterations, where $\\widetilde{f}_0$ denotes the lower bound, typically set to zero. \n\nPlease let us know if further clarifications are needed."}]

Author rebuttal:

Meta-review: **Summary:** The paper studies a variant of the Polyak step-size for high-dimensional statistical estimation problems, where $d >> n$. The paper shows that using the Polyak step-size leads to an error that grows with the dimension. Instead, the proposed method tries to estimate the restricted Lipschitz smoothness. The resulting algorithm is termed sparse Polyak, as it can be viewed as projected gradient descent with a modified version of Polyak’s step size, where the gradient in the denominator is replaced by its hard-thresholded version. The main theorem shows that under the standard assumptions used in prior work, convergence is similar to when the RSS constant is known. The corollaries specialize the results for specific statistical learning models, and include logistic regression and matrix regression. The paper also contains numerical experiments that show an advantage of the sparse Polyak method.

**Strengths:**
- Well-written paper, with a clear contribution and a novel modification to a classical method. 
- Strong theoretical results showing linear convergence in high-dimensional settings. 
- Carefully designed experimental evaluation to validate the theory 

**Weaknesses:**
- Needing to choose the sparsity level $s$ via cross-validation is not ideal. 

**Decision and Suggested Changes:** The reviewers unanimously agree that the paper is well-written, and its contributions merit acceptance. After carefully reading the rebuttal/discussion, I tend to agree. Please incorporate the reviewer comments in the final version of the paper. In particular, addressing the following concerns will help strengthen the current version of the paper:
- Explanation of the comparison to Jain et al, 2024 and discussion about the termination criterion (response to Rev. LwtJ)
- Clarification about the comparison to classical Polyak (response to Rev. drVv, Rev. yL2Y)
Decision: Accept (poster)