Learning Theory for Kernel Bilevel Optimization

Dear Reviewer,

Thank you very much for your review. We are glad to hear that you found the paper clear and that it addresses an important problem. We appreciate your feedback and understand your concerns, which we address in detail below. Please feel free to reach out with any further questions.

1- Main challenge in KBO: Beyond empirical processes

Our analysis required considering objects that are not standard empirical processes. The functional setting involves processes that take values in an infinite-dimensional space (as opposed to real-valued processes). Hence, we cannot apply classical empirical process theory to those: one could apply the classical result to get an error per dimension, but then summing the errors yields a divergent sum. Instead, we had to use results on $U$-processes (Step 3 in Section 4.2), which, to our knowledge, have not been used in the analysis of finite-dimensional problems. Considering these $U$-processes was key for obtaining the $1/\sqrt{n}$ and $1/\sqrt{m}$ rates. As discussed in the last paragraph of Section 4.2, using other "simpler" techniques results in slower rates. We will clarify this technical novelty in the introduction and throughout the text.

2- Convergence rates of (projected) gradient algorithms

While our convergence rates resemble those of classical gradient algorithms, their derivation in our setting is far from standard due to the functional nature and statistical dependencies specific to the bilevel kernel framework. In particular, the parameter updates do not rely on the exact gradient, but rather on a data-dependent gradient estimator arising from the empirical bilevel problem $\widehat{\text{KBO}}$. This introduces two main sources of errors: (1) a statistical error, which stems from approximating the population gradient and depends on the sample sizes, (2) an optimization error, which arises from performing a finite number of (projected) gradient descent steps.

3- Optimality of the convergence rates

We did not investigate the optimality of the convergence rates, as the primary goal of our paper is to provide the first finite-sample generalization bounds for nonparametric bilevel optimization. We agree that studying the optimality of these rates is an interesting future research direction and plan to include this point in the conclusion of the revised version.

Best regards,

Authors

Dear Reviewer,

Thank you for your detailed review and thoughtful comments. We are pleased that you found our paper strong and appreciate the constructive suggestions to improve it. Below, we provide detailed responses to your remarks, and we would be happy to address any further questions you may have, should any of our responses require clarification.

1- Discussion of the constants

We agree that a more explicit expression would be beneficial. In our case, the most interesting dependence is w.r.t. $\lambda$ as it (1) controls the strong convexity of the inner objective and the smoothness of the outer objective (see Proposition C.4 of Appendix C.2), crucial for optimization, and (2) allows providing generalization for the un-regularized problem ($\lambda= 0$) by controlling the behavior of such constant as $\lambda$ approaches 0, provided additional source assumptions are made (e.g., Smale & Zhou, 2007; Sriperumbudur et al., 2017). Under our current relatively weak assumptions, we are not able to characterize the exact dependence of the constants on $\lambda$, especially to control the behavior of the constant as $\lambda$ decreases to 0: parts of the terms in the constant have only a qualitative dependence on $\lambda$. Providing more quantitative estimates would require stronger regularity assumptions such as Lipschitzness of some of the derivatives. Our choice of assumptions is partly motivated by the simplicity of the resulting proofs; however, it is possible to adapt the proofs to get more quantitative dependence on $\lambda$ and provide generalization errors for un-regularized problems, which would be a great future work direction. We will still add a discussion on the constants and, in particular, discuss how the regularization parameter affects both generalization and optimization: large $\lambda$ yields smoother problems that are faster to optimize using large step sizes, but introduces a large bias in the statistical estimation. When $\lambda$ is a parameter to select, a trade-off that accounts for the bias (un-regularized vs regularized problems), variance (finite sample vs population as done in this study) and optimization (step size) must be found.

2- Remarks on the gradient scheme

Corollary 4.2 informs practical algorithm design choices. Indeed, it underscores that the convergence of the bilevel algorithm requires a careful balance between data availability (sample sizes) and computational budget (number of gradient steps). We will make sure to discuss that right after Corollary 4.2.

3- Connection with existing results in IV

Thank you for raising this point, we will include a discussion on that based on the following. Indeed, our problem formulation is related to IV, which is a particular case, with the important difference that our setting currently only considers regularized inner problems (fixed $\lambda$, independent of sample size), while typical studies on IV consider un-regularized population problems. It is more challenging to study the latter, and existing literature on the topic made important effort in that direction, obtaining minimax optimal rates depending on source assumptions, as mentioned in your questions. The challenge in our case, lies in an orthogonal direction: by considering more general inner objectives, even the study of generalization for regularized problems presents additional obstacles that have not been addressed before (see also the answer to the general question). A next natural step would be to study un-regularized problems, after which a more direct comparison with existing results on IV could be made. This is unfortunately beyond the scope of this work.

4- Can existing mimimax results on IV be used for KBO?

It would be great to obtain minimax rates for KBO. However, it is unclear to us how the techniques used in the context of IV can be used for the KBO. Results for IV heavily rely on the closed-form expression of the first stage, while KBO assumes general loss structure for the inner problem (no closed-form of the minimizer in general). We believe the non-linear dependence of the inner solution on the outer parameter, in the general case, is precisely what makes it hard to apply existing results on IV and two stage regression (see also the response to the general question).

5- Measurability of the kernel

We agree that explicitly including this assumption would make the exposition more rigorous and self-contained. We will add a corresponding measurability condition to our list of formal assumptions in the revised version.

6- General question

Thank you for this insightful question. We will add a discussion based on the following points.

• Spectral techniques. The idea of spectral filtering is powerful as it allows treating several algorithms in a unified manner as illustrated in (Meunier et al., 2024a). To our knowledge, using the results on spectral filtering for the outer loss suggests that such loss must be a quadratic loss in $\omega$. This is true when using a regression loss and if the inner solution has a linear dependence in the outer parameter. Such linearity arises when the inner objective is also a regression loss, as in KIV, but is no longer true in general: using general losses in the inner level, yields a non-linear dependence on the outer parameter, that is even implicit. Because of this, it is unclear to us at this stage how to directly transfer the minimax results in (Meunier et al., 2024b) to our setting.

• Surrogate loss and calibration. The approach based on surrogate loss and calibration, suggests two directions:

  • Using a regression loss as surrogate for either both or one of the losses. However, this would imply using the estimators arising from the use of such calibrated regression losses. It is unclear to us how to construct computable calibrated losses in the form of a quadratic loss from the general losses we considered. Additional assumptions on the structure of these losses would be needed.
  • The other direction is to treat the general losses as surrogate for the quadratic losses, but this means one has to derive statistical bounds for these general losses (which is what the current work does), and then transfer those to control quadratic errors.

Overall, it remains unclear to us how to directly leverage the results in the case of KIV on KBO without restriction on the form of the loss to be essentially quadratic.

Best regards,

Authors

---

References

(Smale & Zhou, 2007) S. Smale and D.-X. Zhou. Learning theory estimates via integral operators and their approximations. Constructive Approximation 2007.

(Sriperumbudur et al., 2017) B. Sriperumbudur, K. Fukumizu, A. Gretton, A. Hyvärinen, and R. Kumar. Density estimation in infinite dimensional exponential families. JMLR 2017.

(Meunier et al., 2024a) D. Meunier, Z. Shen, M. Mollenhauer, A. Gretton, and Z. Li. Optimal rates for vector-valued spectral regularization learning algorithms. NeurIPS 2024.

(Meunier et al., 2024b) D. Meunier, Z. Li, T. Christensen, and A. Gretton. Nonparametric Instrumental Regression via Kernel Methods is Minimax Optimal. arXiv:2411.19653 2024.

Dear Reviewer,

Thank you for your thoughtful review of our paper and for your valuable suggestions for improvement. We answer your concerns below and are happy to clarify any further concerns you may have.

1- Clarifying the purpose of Section 2.3

Thank you for this helpful comment. As suggested, we will explain at the beginning of the section the following idea: Section 2.3 adresses the challenge of performing implicit differentiation for KBO. This first requires deriving an abstract, a priori intractable, expression for the implicit gradient, from which we derive a more concrete expression in Proposition 2.1. By concrete we mean that it requires solving a regression problem in the RKHS (the adjoint problem), which can be approximately solved using finite samples, and plug it into an expectation (which can be also approximated using finite samples). Hence, Proposition 2.1 lays the ground for defining a plug-in estimate, later used in our statistical analysis, and which allows getting around the intractability of the implicit gradient.

2- Abstract operators in L156

Thank you for noticing this point, the phrasing is indeed vague. The operators we meant include $\partial_h$, but also $\partial_h^2$ and $\partial_{\omega, h}^2$, which appear when replacing the Jacobian $\partial_\omega h^\star_\omega$ by its expression in Equation 1. We will make sure to be precise about which operators we had in mind in the revised version.

3- Literature review

Thank you for the suggestion. We will add the following in the appendix of the revised version of our paper, which also answers your question about L213.

• Computational aspects of KBO. The computational complexity of KBO algorithms typically scales cubically in sample size as for many kernel methods. For example, in the case of kernel instrumental variable regression (Singh et al., 2019), solving the inner level system costs $\mathcal{O}(n^3)$, and performing vector matrix products operations costs $\mathcal{O}(mnd)$, leading to a total complexity of $\mathcal{O}(n^3+mnd)$. To improve scalability, potential strategies rely on classical techniques in kernel methods such as Random Fourier Features, which approximate kernel functions in a finite-dimensional feature space and enable more efficient gradient computations (Rahimi & Recht, 2007), and Nyström approximations, which mitigate the computational burden of full kernel matrices by using a low-rank approximation of the kernel (Williams & Seeger, 2001).

• Existing bounds for the non-kernel case. Bao et al. (2021) laid foundational work in this direction by analyzing uniform stability in full-batch bilevel optimization. Their criterion for generalization compares the population outer loss evaluated at the output of a randomized algorithm to the empirical outer loss evaluated using the same algorithm. Given a number $\kappa$ between 0 and 1, they obtain a decay at a rate of $\mathcal{O}(T^\kappa / m)$ for unrolled optimization, which decreases as $1/m$ in outer sample size, but increases with the number of outer iterations made. The generalization criterion is unlike the one we consider and which compares the population outer objective at the theoretically optimal inner solution $h^\star_\omega$ to the empirical loss evaluated at the empirical solution $\hat{h}_\omega$. Complementing this upper bound, Wang et al. (2024) established lower bounds on the uniform stability of gradient-based bilevel algorithms, demonstrating a rate of $\Omega(1/m)$. Building on (Bao et al., 2021), Zhang et al. (2024) extended the analysis to stochastic bilevel optimization, establishing on-average stability bounds and deriving a generalization rate of $\mathcal{O}(1/\sqrt{m})$ due to the presence of stochastic gradients. In a related context, Oymak et al. (2021) studied generalization in neural architecture search using a bilevel formulation, showing that approximate inner solutions and Lipschitz continuity of the outer loss yield a generalization bound of $\mathcal{O}(1/\sqrt{n}+1/\sqrt{m})$. Arora et al. (2020) investigated representation learning for imitation learning via bilevel optimization, offering generalization bounds of order $\mathcal{O}(1/\sqrt{m})$ that depend both on the size of the dataset and the stability of learned representations.

4- Clarification on $\nabla \widehat{\mathcal{F}}(\omega)$

We apologize for the confusion at L224 and thank you for pointing this out. We considered two estimators: $\nabla \widehat{\mathcal{F}}(\omega)$ (first defined in L199 in the text as the gradient of $\widehat{\mathcal{F}}(\omega)$) and $\widehat{\nabla\mathcal{F}}(\omega)$ which is indeed defined after L224 in Equation 4. In L224, we were referring to the first estimator $\nabla \widehat{\mathcal{F}}(\omega)$. To make this clear in the text, we plan to add in Section 3 a reference to Proposition B.1 in Appendix B, where the detailed expression of the first gradient estimator $\nabla \widehat{\mathcal{F}}(\omega)$ is provided.

5- Inner and outer distributions in Section 2.2

Thank you for the suggestion. You are correct: the expectations correspond to different distributions (training for the inner and test for the outer). We will update the equation in the revised version of our paper to explicitly highlight this distinction for clarity.

6- Experiments

Thank you for suggesting the addition of experimental results in the case where $m \neq n$: we find this to be a valuable and interesting experiment. To provide these results, we retain the same experimental setup as in the main paper and vary both $m$ and $n$ simultaneously between 100 and 5000. For each distribution (standard Gaussian and Student’s $t$-distribution with degrees of freedom 2.1, 2.5, and 2.9), we plot heatmaps of the quantities of interest, with $n$ on the $x$-axis and $m$ on the $y$-axis. These quantities include $|\mathcal{F}(\omega_0)-\widehat{\mathcal{F}}(\omega_0)|$, $\lVert\nabla\mathcal{F}(\omega_0)-\widehat{\nabla\mathcal{F}}(\omega_0)\rVert$, $\lVert\nabla\mathcal{F}(\omega_T)\rVert$, and $\min_{i=0,\ldots, T}\lVert \nabla\mathcal{F}(\omega_i)\rVert$. We plan to include them in the appendix of the revised version of our paper. Here, since all cases show similar trends, we report only the results for the Gaussian case and focus on the quantity $\lVert\nabla\mathcal{F}(\omega_0)-\widehat{\nabla\mathcal{F}}(\omega_0)\rVert$. The following table summarizes the results for $m, n \in$ {$100, 500, 1000, 5000$}. As observed, the quantity decreases as both $m$ and $n$ increase, and so is the upper-bound. However, we are unable to experimentally assess the tightness of the upper-bound at this stage.

Best regards,

Authors

---

References

(Singh et al., 2019) R. Singh, M. Sahani, and A. Gretton. Kernel Instrumental Variable Regression. NeurIPS 2019.

(Rahimi & Recht, 2007) A. Rahimi and B. Recht. Random Features for Large-Scale Kernel Machines. NIPS 2007.

(Williams & Seeger, 2001) C. K. I. Williams and M. Seeger. Using the Nyström Method to Speed Up Kernel Machines. NIPS 2001.

(Bao et al., 2021) F. Bao, G. Wu, C. Li, J. Zhu, and B. Zhang. Stability and generalization of bilevel programming in hyperparameter optimization. NeurIPS 2021.

(Wang et al., 2024) R. Wang, C. Zheng, G. Wu, X. Min, X. Zhang, J. Zhou, and C. Li. Lower bounds of uniform stability in gradient-based bilevel algorithms for hyperparameter optimization. NeurIPS 2024.

(Zhang et al., 2024) X. Zhang, H. Chen, B. Gu, T. Gong, and F. Zheng. Fine-grained analysis of stability and generalization for stochastic bilevel optimization. In Proceedings of the Thirty-Third International Joint Conference on Artificial Intelligence, 2024.

(Oymak et al., 2021) S. Oymak, M. Li, and M. Soltanolkotabi. Generalization guarantees for neural architecture search with train-validation split. ICML 2021.

(Arora et al., 2020) S. Arora, S. Du, S. Kakade, Y. Luo, and N. Saunshi. Provable representation learning for imitation learning via bi-level optimization. ICML 2020.

Dear Reviewer,

Thank you for providing feedback and highlighting the strengths of our paper. Below, we address your concerns and welcome any further questions.

1- Assumptions

We will add the following discussion on the assumptions.

• Assumption (A) – Boundedness of the kernel $K$. We require $K$ to be bounded, which holds for the most commonly used kernels (Gaussian, Laplacian, Matérn, spline, etc.). It also holds if $K$ is a Mercer kernel, i.e., a continuous, positive‐definite kernel on a compact domain $\mathcal{X}$. For instance, a kernel built using neural network features (e.g., Neural Tangent Kernel) satisfies this assumption on the space of images, which is compact (pixel values have bounded range).
• Assumption (B) – Compactness of the set of the labels $\mathcal{Y}$. In most supervised learning applications, this assumption is immediately met. For example, in regression tasks one often has $y\in [a,b]$ (e.g., probabilities in [0,1] or physical measurements on a known scale), which is compact. Similarly, in classification problems $\mathcal{Y}=${$1,\dots,C$} is compact.
• Assumption (C) – Regularity of pointwise losses. This assumption holds for most loss functions used in practice, including the squared loss, logistic loss, cross-entropy loss, and KL divergence. This assumption no longer holds when considering SVM with the hinge loss. To handle such a case, one would need to extend the analysis to a non-smooth setting, which would be an interesting future work direction.
• Assumption (D) – Convexity of the inner loss. What makes this convexity assumption mild is the fact that it is made on the variable $v$, which corresponds to the value taken by a prediction model $h$. This is weaker than making a convexity assumption with respect to the parameters of such a model. For instance, as discussed in (Petrulionyte et al., 2024), in the case of a neural network model, convexity in the parameters would not hold in general, but it would still hold in the values taken by the model for most commonly used losses in practice, such as the squared loss, logistic loss, cross-entropy loss, and KL divergence.

2- Extension to stochastic bilevel algorithms

We thank the reviewer for this insightful question. While our theoretical analysis focuses on a full-batch bilevel setting with exact gradients, the proposed framework could be extended to stochastic variants. A promising direction would be to consider approximate kernel representations, such as random Fourier features or Neural Tangent Kernels, which enable scalable learning using kernel methods while preserving useful theoretical properties. These extensions lie beyond the current scope of our work and represent an exciting direction for future research. We will make sure to clearly highlight this point in the conclusion of the revised version.

Best regards,

Authors

---

Reference

(Petrulionyte et al., 2024) I. Petrulionyte, J. Mairal, and M. Arbel. Functional Bilevel Optimization for Machine Learning. NeurIPS 2024.