Accelerating first-order methods for nonconvex-strongly-convex bilevel optimization under general smoothness

We thank Reviewer aCYi for the detailed and thoughtful review. We are grateful for the recognition that:

• Our paper is clearly written and easy to follow;
• We provide a fully first-order accelerated method with provable complexity guarantees under general Hölder smoothness, extending recent nonconvex optimization techniques to the bilevel setting;
• Our experiments validate the theoretical results.

We respectfully acknowledge the concerns raised regarding novelty and literature coverage. Below, we provide point-by-point clarifications and improvements to address these issues.

W1. On related works in bilevel optimization not covered (e.g., [1], [2])

We appreciate the reviewer’s suggestion and fully agree that [1] and [2] are highly relevant.

[1] proposes the Moving-Average Stochastic Bilevel Optimization (MA-SOBA) framework, which achieves an optimal sample complexity of $\mathcal{O}(\epsilon^{-2})$ under relaxed smoothness assumptions. [2] studies the complexity of first-order methods in stochastic bilevel optimization, providing tight upper and lower bounds under different assumptions.

These works address stochastic bilevel optimization, whereas our method focuses on deterministic problems under general Hölder continuity, using fully first-order updates. We will revise the related work section to incorporate a detailed comparison with [1, 2] in terms of oracle complexity, smoothness assumptions, and algorithmic distinctions.

W2. Novelty beyond [3, "PRAF2BA"]

We thank the reviewer for highlighting this. Our approach is distinguished from [3] in several key aspects:

1. Restart scheme: Our scheme (Eq. (5)) is tailored to restricted Hölder continuity, and becomes independent of $\epsilon$ when $\nu_f = \nu_g = 1$. In contrast, [3] uses an $\epsilon$-dependent threshold.

2. Descent analysis: We introduce a novel potential function to handle errors in $\nabla L_\lambda^\star(x)$ (see Lemma 7), which is not present in [3].

3. Weaker precision requirement: Our method only requires $\sigma \propto \epsilon$, while [3] requires $\sigma \propto \epsilon^2$.

4. Same convergence rate in Lipschitz case: Under $\nu_f = \nu_g = 1$, our rate matches that of PRAF2BA: $\widetilde{\mathcal{O}}(\epsilon^{-7/4})$.

5. Hölder-smooth complexity: We are, to our knowledge, the first to achieve:

   $$\mathcal{O}\left(\Delta \ell^{\frac{2+2\nu_f-\nu_f\nu_g}{2+2\nu_f}} \kappa^{\frac{6+7\nu_f-2\nu_f\nu_g}{2+2\nu_f}} \epsilon^{-\frac{4+4\nu_f-\nu_f\nu_g}{2+2\nu_f}}\right)$$

   under general Hölder smoothness using only first-order information.

We will highlight these contributions more clearly in the revised introduction and discussion sections.

W3. Assumption

Lipschitz assumptions on $f$ and $g$ are common in the bilevel literature (e.g., [3–5]). Our generalization to Hölder continuity enables the inclusion of more realistic and structured objectives without compromising tractability. See the response to Q3 for details. We will clarify this motivation in the revised assumptions section.

W5 & Q1. $\kappa$ clarification

Thank you for the opportunity to clarify. Our definition of $\kappa$ includes all relevant smoothness constants up to second order:

This generalization is aligned with [4]. In contrast, [3] uses a simplified $\kappa = \max\{L_f, L_g\}/\mu$. In practice, [3]'s convergence bounds depend on higher-order constants even though they are omitted.

Accordingly, we refine Table 1 to use our definition of $\kappa$:

This ensures consistency and transparency in how $\kappa$ is used across methods.

W4 & Q2. On including comparisons with PRAF2BA and other baselines

We agree that including PRAF2BA is essential. Due to space limits, we included results on the Data Cleaning task ($p=0.2$). Other baselines such as AID-BIO, BA-CG, RAHGD, and ITD-BIO were also aligned with those used in [3].

Q3. On motivation for Hölder continuity generalization

We appreciate this thoughtful question. Our motivation for Hölder continuity generalization is threefold:

1. Theoretical Generalization: Hölder continuity is a strict generalization of Lipschitz continuity, allowing us to cover a broader class of functions. For example, $f(x) = \sqrt{x}$ is (0.5,1)-Hölder continuous but not Lipschitz continuous. Moreover, the optimization community has long studied such settings from a theoretical perspective [6–11].
2. Emerging Applications: Recent studies in related areas, such as positional encodings in graph learning, have demonstrated Hölder continuity properties [11]. These suggest that such regularity conditions may naturally arise in structured models or neural network components. This opens the door to identifying or designing architectures where Hölder-type conditions are more appropriate.
3. Timeliness and Potential Impact: As evidenced in works like [10], there is growing interest in developing optimization algorithms under relaxed smoothness assumptions. By establishing provable guarantees for bilevel optimization under Hölder continuity, we aim to contribute foundational insights that may inspire further investigation into both theoretical properties and practical relevance. Our decision to generalize only Assumption 1(v) from [1](i.e., the Hölder continuity of $\nabla^2\_{xx} f$) was motivated by both analytical tractability and practical relevance. Specifically, this term governs the curvature of the upper-level objective with respect to the decision variable $x$, which plays a central role in our the property of $\nabla^2 L\_\lambda^*(x)$. On the other hand, extending the Hölder continuity to cross-derivatives such as $\nabla^2\_{xy} f$ introduces significant technical challenges. In particular, under Hölder continuity of $\nabla^2\_{xy} f$, it becomes difficult to tightly control terms like $||\nabla^2\_{xy} f||$ in Lemma 2, which would be necessary to ensure that the gradient approximation error remains bounded. Therefore, to maintain clarity and ensure provable convergence under the general Hölder setting, we chose to relax only the most influential smoothness assumption—while keeping others Lipschitz continuous to avoid overly technical complications. We leave further generalizations as future work.

Summary

We sincerely thank the reviewer for their thoughtful and detailed feedback. In the final version, we will (i) expand the related work to include [1,2], (ii) clarify our technical contributions beyond [3], and (iii) include PRAF2BA baselines and revise our complexity table accordingly. We hope these clarifications help demonstrate the novelty and rigor of our work.

References

[1] Chen et al, "Optimal algorithms for stochastic bilevel optimization under relaxed smoothness conditions." Journal of Machine Learning Research 25, no. 151 (2024): 1-51. [2] Kwon et al, "On the complexity of first-order methods in stochastic bilevel optimization." arXiv preprint arXiv:2402.07101 (2024). [3] Li, "Accelerated Fully First-Order Methods for Bilevel and Minimax Optimization." arXiv preprint arXiv:2405.00914 (2024). [4] Chen et al, "Near-optimal nonconvex-strongly-convex bilevel optimization with fully first-order oracles." arXiv preprint arXiv:2306.14853 (2023). [5] Yang et al, "Accelerating inexact hypergradient descent for bilevel optimization." arXiv preprint arXiv:2307.00126 (2023). [6] Berger et al, "On the quality of first-order approximation of functions with Hölder continuous gradient", JOTA 2020. [7] Nesterov, "Universal gradient methods for convex optimization problems", Mathematical Programming 2015. [8] Lan, Bundle-level type methods uniformly optimal for smooth and nonsmooth convex optimization, Mathematical Programming 2015. [9] Chien et al, Convergent Privacy Loss of Noisy-SGD without Convexity and Smoothness, ICLR 2025. [10] Marumo et al, "Universal heavy-ball method for nonconvex optimization under hölder continuous hessians". Mathematical Programming, pages 1–29(2024). [11] Huang et al., "On the Stability of Expressive Positional Encodings for Graphs", ICLR (2024).

We thank reviewer Kp5u for their constructive comments and for recognizing key strengths of our work, including the fully first-order nature of the proposed algorithm and its faster convergence rate compared to prior first-order methods such as [7] and [10]. Below, we address the reviewer’s concerns in detail.

Q1. On the Novelty and Contribution Beyond Prior Works

We appreciate the opportunity to clarify the novelty of our approach. While our work is indeed inspired by ideas from [8], [7], and [10], it is not a mere combination of existing techniques. Our contributions introduce several key innovations, both algorithmic and analytical, which we summarize below:

1. Extension of acceleration to bilevel problems under general smoothness assumptions:

   • In contrast to [9], which addresses single-level nonconvex optimization, our work targets the more challenging bilevel setting.
   • We generalize smoothness assumptions and accelerate first-order methods beyond the frameworks in [7] and [10].

2. New algorithmic design for bilevel acceleration:

   • Our restart scheme (see Eq. (5)) differs from that in [9], as it is designed to operate under restricted Hölder continuity.
     Notably, when $\nu_f = \nu_g = 1$, our restart threshold becomes independent of $\epsilon$, whereas [9] uses a threshold that explicitly depends on $\epsilon$.
   • We introduce a potential function analysis that carefully accounts for gradient estimation errors in $\nabla L\_\lambda^\star(x)$ (see Lemma 7).
   • Our analysis allows the inner-loop solvers to compute approximate gradients with precision $\sigma \propto \epsilon$, in contrast to [8], which requires $\sigma \propto \epsilon^2$.

3. Complexity improvements using only first-order information:

   • Compared to [7], our method achieves a faster convergence rate of $\widetilde{O}(\epsilon^{-7/4})$.
   • While [8] achieves the same $\epsilon$-dependence, it relies on second-order information (Hessian-vector products), whereas our method uses only gradients.

We will explicitly summarize these technical distinctions in the revised version, particularly in the introduction and discussion sections.

Q2. Motivation for Hölder Continuity in Bilevel Optimization

We appreciate this insightful question. While Lipschitz smoothness is a common assumption, it can be overly restrictive in practice, especially in bilevel settings with complex or structured objectives.

Our motivation for considering Hölder continuity is threefold:

1. Theoretical generalization:
   Hölder continuity strictly generalizes Lipschitz smoothness. For example, $f(x) = \sqrt{x}$ is $(0.5, 1)$-Hölder continuous but not Lipschitz continuous. Theoretical works have long studied optimization under such assumptions [1–5].

2. Emerging applications:
   Recent studies in graph learning have identified Hölder-type behavior in positional encodings [6]. While these do not directly correspond to gradient Hölder continuity, they suggest that such regularity may naturally arise in structured models or neural networks.

3. Timeliness and future potential:
   The optimization community is actively exploring relaxed smoothness assumptions (e.g., [5]), indicating a trend toward broadening the theoretical foundation of optimization algorithms. Our work provides provable guarantees under Hölder continuity for bilevel problems, paving the way for future exploration of this setting.

In summary, Hölder continuity allows us to study a broader class of bilevel problems beyond standard smoothness. We view our work as an important step toward a more general theory of bilevel optimization.

Q3. On Experiments Showing Gradient Norm Convergence

We agree that convergence in $\|\|\nabla \varphi(x)\|\|$ (or the surrogate $\|\|\nabla L^\star_\lambda(x)\|\|$) is an important empirical indicator to support the theory.
We have already computed these gradient norms during training in the Data Cleaning experiment. The table below summarizes the results for two noise levels ($p = 0.2$ and $p = 0.4$), showing clear convergence toward stationarity:

We will include these results and their visualizations in the final version to better support the theoretical claims.

Summary

We thank the reviewer for highlighting important points around novelty, smoothness assumptions, and empirical validation. In response, we will:

• Clearly articulate our technical contributions beyond prior work;
• Expand the discussion of Hölder continuity and its relevance;
• Include convergence plots for gradient norm to support theoretical guarantees.

We hope these clarifications demonstrate the significance, rigor, and originality of our work.

References
[1] Berger et al. On the quality of first-order approximation of functions with Hölder continuous gradient. JOTA, 2020.
[2] Nesterov. Universal gradient methods for convex optimization problems. Mathematical Programming, 2015.
[3] Lan. Bundle-level methods for smooth and nonsmooth convex optimization. Math Programming, 2015.
[4] Chien et al. Convergent Privacy Loss of Noisy-SGD without Convexity and Smoothness. ICLR, 2025.
[5] Marumo & Takeda. Universal heavy-ball method under Hölder continuous Hessians. Mathematical Programming, 2024.
[6] Huang et al. On the Stability of Expressive Positional Encodings for Graphs. ICLR, 2024.
[7] Chen et al. Near-optimal nonconvex-strongly-convex bilevel optimization with fully first-order oracles. arXiv, 2023.
[8] Yang et al. Accelerating inexact hypergradient descent for bilevel optimization. arXiv, 2023.
[9] Li et al. Restarted Nonconvex Accelerated Gradient Descent: No More Polylogarithmic Factor in the $O(\epsilon^{-7/4})$ Complexity. PMLR, 2022.
[10] Kwon et al. A fully first-order method for stochastic bilevel optimization. ICML, 2023.

Dear Reviewer Kp5u,

We hope this message finds you well.

We truly appreciate your earlier indication that you might consider increasing your score, and we would be very grateful if you feel our responses have resolved the issues and could kindly reflect this in the system (e.g. via the acknowledgement / final justification).

If there remain any points that are unclear or that you would like us to expand on, we would be glad to provide further clarification at your convenience.

Thank you again for your time and constructive feedback! It has been very helpful to us.

Best regards,

The Authors

We thank reviewer bRp6 for the thoughtful review and constructive suggestions. We are encouraged that the reviewer acknowledges the importance of accelerating bilevel optimization using only first-order information and the generality of our Hölder-type smoothness assumptions. Below, we address each concern in detail.

W1. Clarity and Interpretation of Lemmas

We appreciate the reviewer’s feedback regarding the lack of clarity in the presentation. We will substantially revise the narrative in the final version. Here, we clarify the reviewer’s questions point by point:

(Q1) Lemma 2 – Interpretation of the trade-off:

Lemma 2 provides a Hölder-type bound that consists of two terms:

• The first term, $O(\ell \kappa^{\nu_f})$, increases with $\nu_f$, and reflects the contribution from the upper-level Hölder continuity;
• The second term, $O(\lambda^{1-\nu_g} \ell \kappa^{4+\nu_g}) R^{1-\nu_f}$, decreases with $\nu_f$, and stems from the coupling between the upper and lower objectives via the penalty formulation (see Equation (16)).

To provide more intuition: in the proof of Lemma 2, we decompose $\nabla^2 L_\lambda^\star(x)$ as:

where $A(x)$ corresponds to the Hessian of the penalized upper-level objective, and $B(x)$ arises from the perturbation of the lower-level optimal solution $y_\lambda^\star(x)$. The two terms yield:

• $O(\ell \kappa^{\nu_f})$ (independent of $\lambda$),
• $O(\lambda^{1-\nu_g} \ell \kappa^{4+\nu_g}) R^{1-\nu_f}$ (dependent on $\lambda$ unless $\nu_g = 1$).

When $\nu_f$ decreases, the first term becomes smaller, as it relates to the Hölder constant of the composition of $\nabla_{xx}^2 f(x, y)$ with $y_\lambda^\star(x)$. However, the second term increases due to sensitivity to the radius $R$:

In Algorithm 1, we set $\mathcal{R}$ (Eq. (6)) and $\lambda$ (Theorem 1) such that the second term dominates. We will emphasize this intuition more clearly in the revised version.

(Q2) Role of Lemmas 5, 6, and 8:

• Lemma 5 bounds the error between the average gradient over an epoch and the true gradient at the averaged iterate $\bar{w}_k$. This is essential because our algorithm outputs $\bar{w}_k$ and we need this estimate to control $\|\|\nabla L\_\lambda^\star(\bar{w}_k)\|\|$.
• Lemma 6 provides a Hessian-free quadratic surrogate for the function difference, used to show descent in the potential function (Eq. (9)).
• Lemma 8 gives the aggregate descent over one epoch. Together with Lemma 9, it ensures a convergence rate on $\min_i \|\|\nabla L^\star_\lambda(\bar{w}_i)\|\|$.

We agree that the current version lacks intuitive discussion. In the final version, we will add brief commentary before each lemma to clarify its role in the convergence proof.

W1. When does Lemma 8 ensure descent?

As shown in the proof of Proposition 1, whenever $\bar{w}\_{t,k}$ is not an $\epsilon$-stationary point, Lemma 8 ensures that $L\_\lambda^\star(x)$ decreases in the t-th epoch.

W1. Algorithm Clarity

Thank you for the helpful suggestion. In the final version, we will:

• Add line-by-line comments to Algorithm 1 to explain the purpose of each operation (e.g., hypergradient estimation, extrapolation, restart check);
• Explicitly define the estimator $\hat \nabla L^\star_\lambda$ in the main text (used in Line 6), where:

with $y_{t,k}$ and $z_{t,k}$ computed in Lines 4 and 5.

Q3. Clarification on line 192 after Lemma 8:

The sentence was poorly phrased. Lemma 8 holds for any $\sigma \geq 0$, where $\sigma$ quantifies the inaccuracy of the gradient estimator.

To aid intuition: when $\sigma = 0$, i.e., the estimator is exact, the function value $L_\lambda^\star(x)$ decreases monotonically. This implies that for sufficiently small $\sigma$, the descent property remains approximately valid.

Importantly, Algorithm 1 does not require $\sigma = 0$—convergence is guaranteed under bounded inexactness. We will revise the phrasing for clarity.

Q4. Clarification on Comparison to Yang et al. (2023)

The statement on line 214 that our method is "better than" [Yang et al., 2023] refers to three key distinctions:

• Under standard Lipschitz continuity, we achieve the same convergence rate of $O(\epsilon^{-7/4})$ using only first-order information, while [Yang et al.] requires Hessian-vector products;
• Under general Hölder continuity, we provide the first provable acceleration result for nonconvex–strongly convex bilevel problems;
• Our analysis assumes that inner-loop solvers achieve gradient estimates with precision $\sigma = O(\epsilon)$, while [Yang et al.] and [Li et al.] require stricter precision $\sigma = O(\epsilon^2)$.

We will clarify this sentence to better reflect our advantages in oracle efficiency and smoothness generality.

W1. Experimental Comparison

We agree that RAGD-GS and RAHGD exhibit similar performance on small-scale tasks, which is expected since both achieve the same asymptotic rate $\widetilde{O}(\epsilon^{-7/4})$. Our empirical contribution is to show that RAGD-GS reaches this performance using only first-order gradients, while RAHGD requires second-order information. We provide additional quantitative results of Data Hypercleaning over five different random seeds to strengthen this point.

Final Remarks

We thank the reviewer again for the detailed and insightful feedback. In response, we will:

• Improve clarity by reorganizing and annotating technical content;
• Provide intuitive explanations for key lemmas and parameters;
• Clarify and strengthen comparisons to related work;
• Expand experimental analysis and reporting.

We hope these clarifications demonstrate the rigor and value of our contributions.

Thank you very much for your thoughtful follow-up and for considering an increased score based on our clarifications.

We truly appreciate your recognition that our contribution lies in achieving competitive performance using only first-order information under relaxed smoothness assumptions, while employing a distinct restart condition to enable acceleration. We are also grateful for your constructive suggestion on improving the writing. We will take this seriously and are fully committed to thoroughly revising the manuscript to improve clarity and incorporate all the reviewers’ comments.

Thank you again for your time and valuable feedback.

We thank Reviewer 3FGu for their thoughtful review and for recognizing the significance of our contributions, notably:

• The development of an efficient accelerated first-order method for the challenging nonconvex–strongly-convex bilevel optimization setting;
• The use of Hölder continuity, which generalizes standard Lipschitz assumptions and broadens applicability;
• The improved iteration complexity of $O(\epsilon^{-7/4})$ under first-order oracle access, matching or surpassing the performance of prior second-order oracles.

We appreciate the constructive feedback, and we respond to the key concerns in detail below.

W1. On the use and interpretation of the condition number $\kappa$ Thank you for raising this clarification. Indeed, our definition of $\kappa$ encompasses all the relevant smoothness constants ($\ell := \max\{C_f, L_f, \rho_f, H_f, L_g, \rho_g, M_g\}$), including first- and second-order Hölder/Lipschitz parameters for both $f$ and $g$, similar definition as in [1]. While in [2], $\kappa$ is defined as $\max{L_g, L_f}/\mu$. In [2, Lemma 2.6] and [2, Lemma 2.7], the detailed form of Lipschitz constant is actually dependent on second order or third order Lipschitz bounds. The authors of [2] omit the high-order smoothness constants and only consider $\{L_f,L_g,\mu\}$. Thus under our definition of $\kappa$, the complexity bound of [2] in our paper's Table 1 should be modified as following:

W2. On applicability to only deterministic settings

We agree that extending the method to the stochastic setting is an important direction. However, our focus in this work is to first close the gap in deterministic nonconvex–strongly convex bilevel optimization using fully first order oracle under Hölder continuity. Our current analysis lays the necessary groundwork for future stochastic extensions. [5] proposes the Moving-Average Stochastic Bilevel Optimization (MA-SOBA) framework, which achieves an optimal sample complexity of under relaxed smoothness assumptions. [6] studies the complexity of first-order methods in stochastic bilevel optimization, providing tight upper and lower bounds under different assumptions. We will add a discussion on this in the limitations section to better reflect this scope.

W3. On unknown parameters (e.g., smoothness constants) in Algorithm 1

This is a valid concern. In our current theoretical analysis, we assume known smoothness constants to simplify exposition and achieve clean complexity bounds. In practice, we use a grid search over reasonable ranges (see Section 5). Our parameters setting of Algorithm 1 is common in work before [1,2].

W4. On empirical results and reproducibility (seeds and error bars) We acknowledge this oversight. We provide experiments on Data Cleaning to report the mean ± standard deviation over 5 random seeds. Due to space limitations, we only report results for "RAGD-GS" and "RAHGD".

Q1. On Complexity Accounting: AGD Steps and Overall Iteration Complexity

Thank you for the question. In our paper, we denote the total number of epochs as $T$, where each epoch corresponds to one outer iteration of Algorithm 1 before an $\mathcal{O}(\epsilon)$-stationary point is found.

In the proof of Proposition 1, we assume that the algorithm has not yet found an $\mathcal{O}(\epsilon)$-stationary point after epochs $1$ to $t$. Let us denote by $k$ the total number of outer iterations at the end of epoch $t$. As shown in line 505, we derive the following descent guarantee for the $t$-th epoch:

By summing inequality (1) over all $t$ epochs, we obtain:

where $K$ is the total number of outer iterations. Since $\Delta_\lambda$ is finite, it follows that $K$ must also be finite, and therefore $t$ must be bounded. This implies that Algorithm 1 terminates in a finite number of outer iterations.

The resulting outer-loop complexity of $\widetilde{\mathcal{O}}(\epsilon^{-7/4})$ matches the known optimal complexity for nonconvex optimization [4]. Moreover, as established in Theorem 2, we provide the total first-order oracle complexity, which includes the cost of inner AGD iterations. The inner AGD subroutine enjoys linear convergence, so its per-call complexity is logarithmic in $1/\epsilon$.

We will revise the main text and appendix to better clarify this complexity accounting and ensure full transparency.

Reference [1] Lesi Chen, Yaohua Ma, and Jingzhao Zhang. "Near-optimal nonconvex-strongly-convex bilevel optimization with fully first-order oracles." arXiv preprint arXiv:2306.14853 (2023). [2] Yang, Haikuo, et al. "Accelerating inexact hypergradient descent for bilevel optimization." arXiv preprint arXiv:2307.00126 (2023). [3] Jeongyeol Kwon, Dohyun Kwon, Stephen Wright, and Robert D Nowak. "A fully first-order method for stochastic bilevel optimization". In International Conference on Machine Learning, pages 18083–18113. PMLR (2023). [4] Huan Li and Zhouchen Lin. "Restarted nonconvex accelerated gradient descent: No more polylogarithmic factor in the in the o (epsilonˆ(-7/4)) complexity". Journal of Machine Learning Research, 24(157):1–37 (2023). [5] Chen, Xuxing, Tesi Xiao, and Krishnakumar Balasubramanian. "Optimal algorithms for stochastic bilevel optimization under relaxed smoothness conditions." Journal of Machine Learning Research 25, no. 151 (2024): 1-51. [6] Kwon, Jeongyeol, Dohyun Kwon, and Hanbaek Lyu. "On the complexity of first-order methods in stochastic bilevel optimization." arXiv preprint arXiv:2402.07101 (2024).