BayeSQP: Bayesian Optimization through Sequential Quadratic Programming

We thank the reviewer for the comments and time taken for the review. We address their questions and stated weaknesses below.

Q1 (and W2): Can you provide runtime breakdowns (e.g., % time in SOCP solve vs. GP inference vs. sampling)?

We provide runtime breakdowns for BayeSQP on within-model comparisons across different dimensions below. Note that we identified a more efficient Hessian evaluation using Einsum notation, which significantly reduces runtime by avoiding nested loops. Accordingly, an updated Figure 5b reflects these improved runtimes.

From the results above, we can observe that the active sampling-based line search (TS) is the most compute intensive part of BayeSQP, followed by the building the Hessian, and lastly solving the SOCP subproblem. The remaining runtime is attributed to logging overhead (note that with only 3sec of total runtime already 1.5sec of logging are quite significant) and initializing the GP models which we reinitialize in each iteration following various BoTorch tutorials. Maintaining the same model over iterations and appending the data should further decrease runtime. The time for initialization is however similar for all algorithms as it all is done through the same constructors from BoTorch. Also the time for logging is similar across methods. We will include this runtime breakdown in a similar form to Figure 5b in the paper. However, for this figure, we will report the time per BayeSQP steps. Note that above, for dimension 4, we have a longer runtime than for dimension 8 because the number of substeps is larger. The table below shows the runtime breakdown per substep. The percentages are based on the sum of the core components.

The per-step analysis reveals expected dimensional scaling patterns that provide insights into algorithmic behavior. While Hessian evaluation time increases with dimension (from 8.1% to 20.1%), absolute computation times remain efficient even at higher dimensions (0.05 seconds at 96 dimensions). We appreciate this suggestion, as the breakdown provides valuable insights into our algorithm's computational characteristics of our core components.

Note that learning hyperparameters adds the same runtime for all baselines as the log marginal likelihood objective only operates on the marginal standard GP. Upon request of reviewer Br7S we now also have added a discussion on the complexity of the different parts of the algorithm.

Q2: Would first-order-only BayeSQP (without Hessians) perform similarly? An ablation would clarify the benefit of the full model.

This is an interesting question and follows along the lines of the request by reviewer DWmx in W2. In the revised version, we have now added an additional first-order BO method, namely MPD [1], which has the substeps of learning the current gradient and then performing gradient descent using the predicted gradients. To an extent, BayeSQP is a second-order extension of this method. We can observe that after 8 dimensions, BayeSQP also outperforms this baseline leveraging the second-order approximation in a SQP style.

Q3: Have you considered using initialization heuristics or multi-start strategies to mitigate local-only behavior?

We have not explored initialization heuristics or multi-starts but see it as an exciting future direction of research. In general, we believe that BayeSQP can be also be applied in a batched setting using strategies similar to TuRBO [2], which uses for each entry in a batch a separate trust region.

Q4 (and W3): Why not include at least one CEC-style constrained benchmark, even if surrogate mismatch is expected, to demonstrate robustness?

Thank you for the comment. By including Ackley, we aimed to demonstrate our method's limitations on such “adversarial” objective functions for purely local approaches. For this revision, we prioritized implementing two additional baselines (SAASBO as requested by reviewer Br7S, and MPD discussed above) as well as performing ablations to further provide insights into the behavior of BayeSQP. As you note, the CEC-style functions are often stated in a way in which these local approaches tend to yield not meaningful results. We will follow your advice and further expand the discussion of BayeSQP's limitations in such settings within the main paper.

Q5: Could you expand on when the uncertainty-aware formulation (δ < 0.5) is most beneficial — e.g., high noise, narrow constraints, ill-conditioning?

We see the benefits of our uncertainty-aware formulation in two different mechanisms. First, incorporating uncertainty typically improves optimization behavior by providing more reliable search directions when gradient estimates are uncertain. Second, the choice of $\delta$ will have a direct influence on the robustness of the final solution. The magnitude of $\delta$ will likely depend on the application at hand and will mostly influence the final optimal solution. If, e.g., one aims to find a robust constraint satisfying solution one can directly encode this through the choice of very small $\delta$. By design, it will yield a solution that has some safety margin. This touches also on the high-noise setting, where we naturally want some margin to the constraint. With Figure 3, we aimed to provide some qualitative insights on the influence of delta.

References

[1] Nguyen, Quan, et al. "Local Bayesian optimization via maximizing probability of descent." Advances in neural information processing systems 35 (2022): 13190-13202.

[2] Eriksson, David, et al. "Scalable global optimization via local Bayesian optimization." Advances in neural information processing systems 32 (2019).

We thank the reviewer for the comments and time taken for the review. We address their questions and stated weaknesses below.

Q1: The method uses constrained Thompson sampling with a fixed number of samples (M=3) for line search. Could the authors elaborate on how sensitive BayeSQP is to this choice? Would a larger budget significantly improve performance?

Thank you for the interesting question. We performed a sweep over different values with $M \in \\{ 1, …, 10\\}$ to test the influence on the final performance on the 64D constrained within model objective functions. We find that for a within-model setting the median performance does not significantly change but by increasing $M$ the 90% upper quantile of the performance can be reduced. In general, we hypothesize that there is a deep connection of $M$ and the number of subsampling steps $K$. We therefore conducted a second sweep over combinations (similar to Figure 5c) with $K \in \\{4, 8, 16, 32, 64, 128, 256\\}$ and found that one can reduce the number of subsampling steps when slightly increasing $M$ for a “pure” within-model (WM) setting, i.e., with perfect model knowledge. We performed the same sweep also with online hyperparameter optimization (OHO). Here, we find that increasing $M$ can improve performance but decreasing $K$ comes at a higher cost as a higher number of subsamples facilitates better learning of hyperparameters. Increasing $K$ beyond D also does not improve performance as it reduces the number of line searches given a fixed budget. That said, using D+1 subsampling points is a compelling heuristic as it offers a clean theoretical justification, as discussed in the paper, however, if more inductive biases are available one can probably reduce this number. We will include these additional sensitivity ablations in the main body of the paper.

Q2: However, one could imagine alternative approaches, such as using Conditional Value-at-Risk (CVaR) or distributionally robust optimization formulations. Have the authors considered or tested such alternatives? We have not yet experimented with CVaR or distributionally robust formulations, but we agree that this is a promising direction. We see BayeSQP as a flexible framework that could benefit from various extensions and refinements, including improvements to the subproblem formulation. In particular, incorporating risk-sensitive objectives like CVaR or potentially leveraging approximate Hessian uncertainties for more informed subsampling are exciting directions which seem promising for future research. We hope that this initial version of BayeSQP can serve as a foundation for further research in these directions.

Q3: Although the authors include an illustrative analysis of the role of and , a more systematic sensitivity study would be helpful. Providing practical recommendations for setting these parameters could make the method more accessible to non-expert users.

Thank you for that comment. With the illustrative example, we wanted to provide intuition on the behavior of BayeSQP for different confidence parameters and then highlight with Figure 5 c) that the exact choice of the delta does not seem to be critical as long as some delta_c < 0.5 is specified. Still, we think that we can further improve in explicit recommendations for practitioners as well as researchers. We believe that the above mentioned results will help with this and we will make sure to state suggestions more explicitly in a final version of the manuscript.

Q4: Could the authors clarify the meaning of this expectation, or consider reformulating the objective if the method is not intended to handle true stochastic optimization problems?

This was a mistake from our side and thank you for spotting this! We have fixed this in the manuscript. To clarify: we assume a fixed f from which we can only obtain noisy zero-order observations.

W1: The experiments are well-executed, but they’re limited to synthetic functions. Including one or two real-world applications would help demonstrate the practical impact and usability of the method.

We agree that such an application would further strengthen the contributions of our paper. Unfortunately, we were not able to set up a new experiment in the limited time of the rebuttal and rather focused on adding baselines (MPD as requested in W2 and SAASBO requested by reviewer Br7S to which we also refer to for a longer discussion on this baseline) as well as conducting experiments that further provide insights into the behavior of BayeSQP as discussed in Q1. We also added additional out-of-model comparisons (similar to the mentioned experiments in Q1) to highlight BayeSQPs performance also in the OHO setting. In these experiments, we observe a slight drop in performance, which is to be expected, but the rankings among the algorithms remains consistent with BayeSQP showing the best performance in the high-dimensional setting.

W2: The comparison focuses on established baselines like TuRBO and logEI, but doesn’t include some newer methods for local BO, such as the algorithm proposed by Nguyen et al. (2022). A broader comparison would help position BayeSQP more clearly in the current landscape.

Thank you for the comment and suggestion. We agree that including more recent local BO methods would better position BayeSQP within the current landscape and can improve the empirical section of our paper. Adding the method from the paper above is especially interesting as, in a sense, BayeSQP is an extension of this class of first-order BO methods to a second order method. We therefore added MPD (method from the mentioned paper). We find that BayeSQP outperforms this baseline in all experiments with dimension greater than 8 and we furthermore observe that in the high-dimensional experiments, also TuRBO outperformed this baseline. We should also note that extending MPD (and similar methods like GIBO) to constrained optimization settings is not straightforward.

Dear Reviewer DWmx, given the authors' response, please raise any remaining questions and/or concerns in a timely fashion so the authors have a chance to reply. Thanks!

We thank the reviewer for the comments and time taken for the review. We address their questions and stated weaknesses below.

Q1: Is TuRBO incorrectly labeled as SCBO on line 96?

Yes, thank you very much for catching this! We have fixed this typo in the document.

Q2 (and Weakness 1): There is substantial work on scalable BO which is not discussed and never compared with. For example [1] and [3]. While the focus is on constrained optimization, there is a section with unconstrained problems, so more baselines need to be included.

Thank you for pointing out these additional references (including [2]). We agree that there is a large body of work on making BO more scalable. In our submitted version, we focused primarily on local methods—a subset of scalable BO approaches—since BayeSQP extends prior work that leverages first-order optimization techniques. However, your comment highlights that we should expand our related work discussion. We therefore include [1-3] and subsequent approaches in the related work section, adding a paragraph on scalable BO methods before discussing the local approaches. We also provide explicit experimental comparisons by incorporating an additional baseline from this set of methods to our experiments.

Here, we decided to include SAASBO in our unconstrained experiments as a popular baseline in the sparse subspace setting. We find that it outperforms BayeSQP up to 8D, performs similarly on 16D, and afterwards BayeSQP also outperforms this baseline. Thank you for suggesting this baseline in Q4.

Q3 (and Weakness 2): I didn't see the assumptions on the objective and constraints in the main text. Is this included somewhere and if not can the authors please include this information?

No, you are correct that we did not explicitly formulate assumptions on the objective and constraints in our submitted version. While BayeSQP does not impose explicit structural assumptions (such as sparse subspaces) on the objective and constraints, our modeling approach does embed implicit assumptions about these functions. The sub-problem is based on the assumption that the modeling approach of the objective with a GP is reasonable. Since we rely also on derivatives of the GP, severe model misspecification could potentially harm performance more than other approaches which we in part see on the Ackley experiments. We will make sure to make this more explicit in the final version. Specifically, we will add the following to Section 2: “Here and in the following, we will assume that $f(x)$ and $c(x)$ are samples from Gaussian processes and we will base our algorithm as well as modelling approach on this assumption”. We will further add a small discussion to the existing one on limitations stating that a severe model mismatch could harm the performance of BayeSQP.

Q4: It would be nice to see more baselines included to know how much benefit one gets from incorporating this over existing, potentially simpler methods. For the unconstrained problems, why not compare with methods assuming an additive structure to the function [1], SAASBO [2], or random projections [3]?

As discussed above, we agree that a discussion on these methods was missing and we now have added it to the manuscript. We furthermore agree that a comparison to one of the methods mentioned above would benefit the current presentation. Therefore, we included SAASBO as a representative recent method from this category. SAASBO matched or exceeded BayeSQP's performance up to 16 dimensions, but BayeSQP significantly outperformed SAASBO at 32 dimensions. For dimensions 64 and 96, SAASBO could not complete runs within our 24-hour computational budget. This computational limitation has also been noted in prior evaluations of SAASBO [4].

Q5 (and Weakness 3): What is the computational complexity of the algorithm? This is not explicitly stated in the text as far as I can see.

The computational complexity with respect to the number of data points follows standard BO at $O(N^3)$ for Gram matrix inversion. Since this inverse is cached, it can be reused for evaluating both the gradient GP and the Hessian mean. The training time is identical to a standard GP as the log marginal likelihood objective only operates on the marginal standard GP. Evaluating the mean of the Hessian at a test point scales as $O(N^2 d^2)$. The $d^2$ factor arises from computing all d×d Hessian elements, while $N^2$ comes from the matrix operations with the cached Gram matrix inverse. The optimization subproblem involves solving a second-order cone program (SOCP) with $d+1+m$ decision variables (recall that m is the number of constraints), which has polynomial complexity independent of the number of data points. For our problem scale, CVXOPT provides adequate performance for solving the optimization subproblem, as demonstrated in the runtime results below.

On that note, we want to highlight that we found a more efficient implementation to evaluate the Hessian leveraging Einsum notation vs. iteratively building up the Hessian. With this, we were able to significantly reduce runtime by eliminating redundant computations and leveraging vectorized tensor operations. The optimized implementation achieves almost a $d^2$ speed up as it avoids nested loops and we can now report runtimes on, e.g., the 96D within-model functions setting of TuRBO (WM: 310sec, OHO: 12 206sec), logEI (WM: 186sec, OHO: 1743 sec), BayeSQP (WM: 7sec, OHO: 110sec) where WM is pure within-model, i.e., perfect model knowledge and OHO is with online hyperparameter optimization. As stated above, hyperparameter optimization is similar for all methods.

References

[1] Kandasamy, Kirthevasan, Jeff Schneider, and Barnabás Póczos. "High dimensional Bayesian optimisation and bandits via additive models." International conference on machine learning. PMLR, 2015.

[2] Eriksson, David, and Martin Jankowiak. "High-dimensional Bayesian optimization with sparse axis-aligned subspaces." Uncertainty in Artificial Intelligence. PMLR, 2021.

[3] Wang, Ziyu, et al. "Bayesian optimization in a billion dimensions via random embeddings." Journal of Artificial Intelligence Research 55 (2016): 361-387.

[4] Papenmeier, Leonard, Luigi Nardi, and Matthias Poloczek. "Increasing the scope as you learn: Adaptive Bayesian optimization in nested subspaces." Advances in Neural Information Processing Systems 35 (2022): 11586-11601.

Dear Reviewer Br7S, given the authors' response, please raise any remaining questions and/or concerns in a timely fashion so the authors have a chance to reply. Thanks!

We thank the reviewer for the comments and time taken for the review. We address their questions and stated weaknesses below.

Q1 (and W1): Could BayeSQP's strong performance simply come from its initial feasibility guarantee rather than the core algorithm itself?

Thank you for the question but we think there may be a misunderstanding. We think our explanation of feasibility may not have been sufficiently clear, as there are actually two distinct types of feasibility at play in BayeSQP that we should distinguish better. First, there is the feasibility of a sample point x_t with respect to the original problem constraints. Second, there is the feasibility of the SOCP subproblem that BayeSQP solves at each iteration. Upon reviewing our manuscript, we recognize that our terminology was not consistently precise in distinguishing these concepts, which we will address in a final version. To directly address the question: BayeSQP implements a slacked version of the SOCP subproblem to ensure computational feasibility of the subproblem itself, not to guarantee that the resulting sample points satisfy the original problem constraints. Therefore, BayeSQP's performance cannot be attributed to an initial feasibility guarantee for the sample points. In fact, we often have multiple iterations at the start, e.g., on the constrained Hartmann, where the sample points remain infeasible similar to the other algorithms.

Q2 (and W2): Why does logEI converge quickly and perform poorly in high-dimensional problems?

LogEI's poor performance in high-dimensional settings stems primarily from its global search approach. Essentially, as nicely stated in [1], when the true objective function has small lengthscales relative to domain, global search becomes increasingly ineffective as dimensionality grows. Given a small length scale for the unit hypercube and a large number of dimensions, a global search is bound to fail to a certain extent. Therefore, as problem complexity increases, logEI's behavior increasingly resembles random search which becomes particularly evident in the "pure" within-model setting as also highlighted in prior work (e.g., [2]).

To investigate whether enforcing larger lengthscales could mitigate this issue, we conducted an experiment where we initialized the lengthscales as sqrt(d) [3] on the 64D within model objective. We can observe that this does improve performance slightly however, both TuRBO and BayeSQP are able to also outperform this baseline with BayeSQP showing the strongest performance.

Q3: How many initial points were used for logEI in these experiments?

For all experiments and all baselines, we use D initial points. Such a scaling of the number of initial points with the number of dimensions is quite common.

W1: The performance of BayeSQP may come from guaranteeing feasible points at the start rather than from the rest of the pieces of its algorithm. It would be useful to see an ablation study without the fallback to isolate its effect.

As discussed above, we believe there might have been misunderstanding about which feasibility was meant in the submitted version when discussing the slacked version of the SOCP. We adjust the text to make the distinction more precise. If we however misunderstood the point you were making, please let us know so that we can further clarify.

W2: I think the number of repetitions, 20, is insufficient. They should use at least 30, and preferably 50 or 100, as is common in the literature.

We appreciate this feedback and have increased the number of repetitions from 20 to 32, following the logEI paper [3]. We are happy to increase this number even further for the final version but opted for 32 due to the limited time frame for the rebuttal. Importantly, our results remain consistent with the results reported in the submitted version of the paper.

Also thank you for spotting these typos and the recommendation for Table 1! We have fixed the typos and included bold highlighting in the manuscript!

References

[1] Wüthrich, Manuel, Bernhard Schölkopf, and Andreas Krause. "Regret bounds for Gaussian-process optimization in large domains." Advances in Neural Information Processing Systems 34 (2021): 7385-7396.

[2] Müller, Sarah, Alexander von Rohr, and Sebastian Trimpe. "Local policy search with Bayesian optimization." Advances in Neural Information Processing Systems 34 (2021): 20708-20720.

[3] Ament, Sebastian, et al. "Unexpected improvements to expected improvement for bayesian optimization." Advances in Neural Information Processing Systems 36 (2023): 20577-20612.

The authors have satisfactorily addressed my concerns and those of the other reviewers I wanted to understand.