We sincerely thank the reviewer for the thoughtful feedback. We especially appreciate your recognition of the real-world in-vitro results. Below, we respond to each of your questions:

Q1: Batch Extension and Equation 6

You’re absolutely right to flag this subtlety. The form in Equation 6 approximates the marginal gain in coverage under the simplifying assumption that only the best point in the batch contributes to coverage, rather than the full set. This is computationally much simpler to optimize in practice, but as you point out, it does not fully reflect the potential joint contribution of all q points. Your suggested form would indeed capture the “true” expected joint gain more faithfully. We chose the current formulation as a practical approximation. That said, exploring more accurate batch ECI objectives, like the one you suggest, is a very interesting direction for future work.

Q2: Redundancy in the ECI Expression

Thank you for this insightful question. You’re correct that, in theory, the “max(0, …)” term in the ECI expression may not be strictly necessary, since the coverage score can only be increased or remain the same after a new data point is added. Similarly, since the baseline coverage value (before adding the new point) is fixed, one could factor it out of the expectation and simplify the expression. We kept the original form primarily for consistency with standard Expected Improvement-style notation and to make the marginal gain in coverage explicit. In practice, the predictive distribution over objective values still introduces enough variability that the expectation continues to provide a meaningful acquisition signal, so we haven’t observed any under-exploration behavior due to this formulation. That said, your observation is well-taken, and this is something we plan to explore further in future work.

Q3: Single-Objective Baselines and Figure 2

Thank you for asking for clarification here. For each baseline (e.g., TuRBO), we perform T independent runs, one per objective, each with its own budget. To ensure a fair comparison on sample efficiency, we combine the results from these T runs by selecting the best K solutions (among the T collected ones) that maximize coverage. We then report the cumulative number of queries across all objectives (i.e., total function evaluations across the T runs) on the x-axis of Figure 2. So, there is no meta-algorithm cycling through objectives in a round-robin fashion. Instead, we assume access to T independent runs and aggregate the solutions post hoc to compute coverage, essentially giving the baselines every possible advantage while still showing that MOCOBO outperforms them.

Formatting and SI Units

Thank you for pointing out the SI unit formatting issues and appendix structure. We’ll take your suggestions into account when preparing the camera-ready version, including exploring the use of the siunitx package and reorganizing the appendix for improved readability.

We thank the reviewer for their thoughtful feedback and for highlighting the practical relevance and strong empirical performance of MOCOBO. We respond to each point below:

1. Coverage Score

We appreciate the reviewer’s suggestion to examine alternative definitions of coverage. However, we respectfully clarify that our coverage score is not a design choice, but rather the problem definition itself, grounded in prior work. Specifically, the continuous coverage score we use was originally formalized in the CluSO paper and is a principled way to pose the goal of finding a small set of K solutions such that each of the T objectives is optimized as much as possible by at least one solution in the set. This captures the setting we care about, for example, discovering K antibiotics such that each pathogen is targeted as effectively as possible, not just adequately. By contrast, threshold-based coverage is indeed a meaningful formulation for some domains but it constitutes a fundamentally different problem: it assumes that we know, a priori, a threshold value for each objective beyond which we do not care to improve performance further. This makes sense in some safety-critical settings (e.g., drug toxicity screening), but is not aligned with the domains we focus on where better objective values are always desirable and thresholds are typically unknown or unhelpful. For these reasons, we adopted the coverage formulation from CluSO and did not include threshold-based comparisons. That said, we appreciate this point and now explicitly clarify in the paper that the problem definition we adopt is the same as has been used in the black-box coverage literature (CluSO), and we better motivate its relevance to our domains of interest. We also agree with the reviewer that methods like MORBO are not designed for coverage, and we will further emphasize in the main text that these are adapted baselines used to highlight that off-the-shelf BO methods do not trivially solve the coverage problem. We include them not to criticize, but to underscore the importance of explicitly addressing this problem class. We do directly compare to CluSO, which is the only existing method we are aware of that does try to solve the black-box coverage optimization problem directly.

2. Trust Regions (TRs)

Thank you for raising this point. We have added an ablation study comparing MOCOBO with and without trust regions (TRs) on the rover task with 4 obstacle courses, and on the template free peptide design task. As expected, MOCOBO performs substantially better with TRs, confirming that TRs significantly improve performance in the high-dimensional settings we consider. We have added plots and discussion of this to the appendix. We report below the average coverage score achieved by MOCOBO with and without TRs after N function evaluations (higher is better).

Rover Task

• N=40K, with TRs: $17.52 \pm 0.32$, without TRs: $-6.95 \pm 0.92$.
• N=80K, with TRs: $17.76 \pm 0.31$, without TRs: $-4.64 \pm 0.69$.

Peptide Task

• N=300K, with TRs: $-45.30 \pm 6.38$, without TRs: $-242.53 \pm 11.56$.
• N=600K, with TRs: $-39.50 \pm 5.74$, without TRs: $-225.84 \pm 13.30$.

To clarify the motivation for using TRs: it is well-established in the BO literature that standard BO (without TRs or other high-dim adaptations) performs poorly in high dims. [1] originally introduced TRs as a principled way to improve performance of high-dim BO. Since then, TRs have become a standard tool for any high-dim BO task. Since all of the tasks we consider are high-dimensional, this precisely why we adopt TRs. We note that recent work has proposed alternative approaches to improve high-dim BO without TRs [2]. However, there remain concerns about the applicability of these approaches to structured domains (e.g., see Table 2 in [3], results show that [2] performs surprisingly poorly on structured tasks).

3. Additional Ablations

We added an ablation comparing to “MOCOBO-EIT": MOCOBO with the reviewer’s suggested "maximum of standard EI per objective" acquisition function. We ran this ablation on the rover task with 4 obstacle courses and on the template free peptide design task. We have added plots and discussion of this to the appendix. We report below the average coverage score achieved by MOCOBO (as proposed with ECI), and by MOCOBO-EIT, after N function evaluations (higher is better).

Rover Task

• N=40K, MOCOBO: $17.52 \pm 0.32$, MOCOBO-EIT: $12.33 \pm 0.92$.
• N=80K, MOCOBO: $17.76 \pm 0.31$, MOCOBO-EIT: $16.00 \pm 1.07$.

Peptide Task

• N=300K, MOCOBO: $-45.30 \pm 6.38$, MOCOBO-EIT: $-180.32 \pm 5.44$.
• N=600K, MOCOBO: $-39.50 \pm 5.74$, MOCOBO-EIT: $-175.71 \pm 4.88$.

MOCOBO with our proposed ECI acquisition function performed significantly better on both tasks, demonstrating the importance of ECI to the performance of MOCOBO. Additionally, MOCOBO-EIT still performed fairly well, outperforming ALL other baseline methods considered in the paper after the full budget of 80K function evaluations on the rover task. This highlights the importance of other aspects MOCOBO.

Additionally, we agree that ablations of q (batch size), initialization dataset size, and K, will further strengthen the paper. We plan to run each of these and add results to the camera ready version.

4. Novelty

We appreciate the reviewer’s concerns regarding novelty and the opportunity to better clarify our contributions. We agree that Algorithm 1 is based on classic submodular maximization (e.g., Nemhauser et al. 1978), and will certainly clarify in the paper that we did not invent submodular optimization! Our novelty lies not in the greedy algorithm itself, but in how it is integrated with other components of MOCOBO to form a practical, scalable solution to black-box coverage optimization. Additionally, we have revised the intro and related work sections to better position MOCOBO as the first method for coverage optimization that is scalable to high-dimensional and structured domains. While CluSO addresses the same problem formulation, it was not designed for high-dimensional or structured input spaces. In fact, as our experiments show, CluSO performs poorly, consistently worse than baselines not designed for coverage, on all of the domains we consider. This is not a criticism of CluSO, but rather a reflection of its limitations in structured or high-dimensional domains. In contrast, MOCOBO is designed with scalability in mind and supports structured and high-dimensional inputs, making it, to our knowledge, the first black-box coverage optimization method applicable to real-world domains like protein and molecular design. We have revised the manuscript to better emphasize this contribution and to clearly differentiate our method from prior work.

5. Convergence Guarantees

We agree that regret bounds for MOCOBO would be interesting and we plan on investigating this in future work. We do kindly note that a lot of recent popular BO papers do not include regret bounds (e.g., [6-15]) or include ones that reduce to proving convergence via running random search as a subsequence (e.g., [4,5]).

6. Other

We have corrected the typo, revised writing to reduce redundancy in coverage definition, fixed caption formatting, and added discussion of limitations regarding computation required to pre-train domain specific VAEs. We thank the reviewer for pointing out that a mixed search space domain could be interesting to add, this is something we will consider for future work.

Q1

Addressed above via new ablation study.

Q2

MOCOBO’s greedy covering set construction does lead to longer overhead time compared to baseline approaches. MOCOBO overhead takes approximately 1.4 times longer than TuRBO and MORBO and approximately 1.2 times longer to run than ROBOT (ROBOT is slower than other baselines due to additional time needed to compute diversity constraints between observations on each iteration). However, the overall additional computation time required by a run of MOCOBO is negligible when optimizing moderately expensive black-box functions, where black-box function evaluation time dominates overall computational cost.

Q3

As we mentioned above, CluSO performs poorly on the tasks considered in this paper because it was not designed to scale to high-dimensional or structured input spaces. We have added discussion of this to the paper.

Q4

See section 2 above.

Cites

[1] “Scalable Global Optimization via Local Bayesian Optimization” Eriksson et al. NeurIPS 2019. [2] “Vanilla Bayesian Optimization Performs Great in High Dimensions” Hvarfner et al. ICML 2024. [3] “A survey and benchmark of high-dimensional Bayesian optimization of discrete sequences” González-Duque et al. NeurIPS 2024. [4] “Scalable Constrained Bayesian Optimization” Eriksson et al. AISTATS 2021. [5] “Discovering Many Diverse Solutions with Bayesian Optimization” Maus et al. AISTATS 2023. [6] “Bayesian optimization and attribute adjustment” Eismann et al. UAI 2018. [7] “Sample-Efficient Optimization in the Latent Space of Deep Generative Models via Weighted Retraining” Tripp et al. NeurIPS 2020. [8] “Max-value Entropy Search for Multi-Objective Bayesian Optimization with Constraints” Belakaria et al. NeurIPS 2020. [9] “High-dimensional Bayesian optimization with sparse axis-aligned subspaces” Eriksson et al. UAI 2021. [10] “Combining Latent Space and Structured Kernels for Bayesian Optimization over Combinatorial Spaces” Deshwal et al. NeurIPS 2021. [11] “Local Latent Space Bayesian Optimization over Structured Inputs” Maus et al. NeurIPS 2022. [12] “Advancing Bayesian Optimization via Learning Correlated Latent Space” Lee et al. NeurIPS 2023. [13] “Joint Composite Latent Space Bayesian Optimization” Maus et al. ICML 2024. [14] “Approximation-Aware Bayesian Optimization” Maus et al. NeurIPS 2024. [15] “Latent Bayesian Optimization via Autoregressive Normalizing Flows” Lee et al. ICLR 2025.

We thank the reviewer for the thoughtful and encouraging feedback. We appreciate the recognition of MOCOBO’s practical relevance, principled design, strong empirical performance, and in vitro validation. We respond to each point below:

MAJOR POINT: Code Release and Reproducibility

We fully agree that releasing code is essential for reproducibility and for supporting the community. While the code was not public at submission time, we want to clarify that we have a complete GitHub repository that includes all code to reproduce all MOCOBO results on all tasks in the paper. We also recognize the importance of accessibility and usability. To that end, we’ve included a detailed README and environment setup to make it easy for others to run MOCOBO and apply it to new problems. We have invested significant effort in this because we care deeply about making our method practical and usable by the broader research and practitioner community.

Our plan is to make this repository fully public upon paper acceptance, and we will include the link in the camera-ready version of the paper.

Ordinarily, we would have liked to share an anonymous version of the repository during the review process so that reviewers could verify results without compromising anonymity. However, due to new NeurIPS rules this year, we are unfortunately not allowed to include links of any kind (even anonymized) in the author rebuttal. We appreciate your understanding and assure you that full code release is planned and prioritized.

MINOR POINTS

Thank you for your detailed proofreading suggestions, and in particular, for catching numerous typos, missing capitalizations, and formatting inconsistencies. These comments will be very helpful in polishing the camera-ready version.

QUESTIONS

Q1: How is the reward computed in Figure 1?

Yes, each curve represents an objective function over the same domain. The “reward” in this case is actually just the coverage score as defined in equation 1. For the single pareto optimal point (the star), the coverage score is just the sum over the three objective function values for that single point. For the covering set of K=2 points (the two black squares), this coverage score is the sum of the higher objective value obtained between the two points for each one of the three objective functions. We decided to use “reward” in Figure 1 rather than “coverage score” just because Figure 1 is discussed at the start of the paper before the “coverage score” is formally introduced. We do think that the “pareto frontier” given by the black dotted line is accurate in the sense that if we move along the x-axis in any direction away from this dotted line, all objective function values decrease, and if we move towards the black dotted line, all objective function values increase. When moving between any pair of x-axis points that are both along the dotted line, some objective functions increase while others decrease. To our knowledge, this defines the black dotted line as the pareto frontier for the three objective functions. Please let us know if that answers your question, we are happy to discuss this with you further.

Q2: Could a logEI version of ECI improve performance?

Indeed using logEI has been shown to lead to substantial improvements in BO over using standard EI. It therefore worth investigating whether some version of logECI could be used instead of ECI to improve MOCOBO similarly. Thank you for the suggestion, this is an exciting direction for future work!

Q3: Clarify mapping between design points and objectives in relation to the MCP

Thank you for this insightful question. We agree that this portion of Appendix A9 would benefit from clearer explanation.

To clarify: in our reduction from the Maximum Coverage Problem (MCP), the connection between a subset $A_i$ and a design point $x_i \in D_s$ is encoded via the binary function values $f_t(x_i) \in$ {$0,1$}. That is, each objective $f_t$ plays the role of an element $e_t \in U$, and each design point $x_i$ corresponds to a subset $A_i \subseteq U$, where an element $e_t$ is “covered” by $A_i$ if and only if $f_t(x_i) = 1$.

This construction ensures that while the design points are not literally subsets, the induced binary coverage behavior matches the structure of the MCP instance. The design point $x_i$ “covers” objective $f_t$ if it achieves a value of 1 under that objective. Thus, the coverage score $c(S)$ counts how many objectives are covered (i.e., attain a value of 1 for at least one point in S), just like MCP counts how many elements are covered by a set of subsets.

We have revised Appendix A9 to more clearly explain this mapping, especially the fact that the “subsets” in our reduction are represented functionally via binary-valued objectives over points in $D_s$, rather than via explicit set membership.

Q4: Why does MOCOBO outperform baselines even when they solve easier sub-problems?

Indeed, the baselines (TuRBO, LOL-BO, ROBOT) are run independently per objective, so they produce T solutions total (one solution specifically designed to optimize each of the T objectives). The best K-covering subset is then selected from those. While this does provide an advantage in terms of raw objective optimization, it lacks joint optimization: the T solutions are optimized independently, so coverage suffers when we select only the best K<T from among those T solutions. In contrast, MOCOBO jointly selects K diverse solutions to maximize total objective coverage. Indeed your comment would be correct if K=T, but when K<T, there is a need to design the small set of K solutions specifically such that some solutions in the set of K are good optimizers for more than one of the T objectives, in order to obtain a set of K<T solutions that together provide at least one good optimizer for each of the T objectives. Please let us know if this answers your question, we are happy to discuss this with you further.

Q5: Practical use case for the rover task?

The rover task is a well-established high-dimensional benchmark in BO. While our adaptation to multi-objective coverage is synthetic, it is intended as a proof-of-concept for settings where agents must be robust to diverse terrain types, e.g., autonomous vehicles navigating multiple environments. Thank you! We thank the reviewer again for their detailed, constructive feedback. We are excited to make our code public and are grateful for your thoughtful review.

We thank the reviewer for their thoughtful and constructive feedback. Below we address all points raised:

Justification for K < T

Thank you for highlighting this point. We agree that the constraint K < T deserves more discussion. This constraint reflects many real-world scenarios with strict evaluation or deployment budgets, such as: 1) In drug discovery, where only a small number of compounds can be synthesized or tested, 2) In materials design, where the cost of fabricating each candidate is high, 3) In robotics or system design, where one might only be able to deploy a handful of control policies, etc. In settings such as these, K < T encodes a practical resource constraint: the user cannot afford to optimize each objective independently and must instead find a small set of solutions with broad utility. We will expand our discussion of the problem setting to clearly articulate this point and further justify the K<T constraint with concrete examples.

Problem Motivation and Presentation

We appreciate the reviewer’s request for greater clarity in the problem framing and example. We note that our problem definition is defined in precise, mathematical terms in the Method’s section (see equation 2). However, we agree with the reviewer that the explanation of the problem setting in the abstract and introduction (text preceding this precise mathematical definition) can be improved to more clearly define the problem in words. For the camera ready version, we will revise the abstract and introduction to more clearly define the goal: namely, to identify a set of K < T solutions that collectively provide strong coverage over all T objectives, such that each objective is well-optimized by at least one solution in the set. This differs from the standard Pareto setting, which typically emphasizes trade-offs at the level of individual solutions, rather than coverage across objectives via a set.

We will revise all relevant text to provide a clearer explanation of the problem definition and to clarify that we are not assuming a single solution can cover all objectives equally well. Rather, the goal is to find a small set of solutions such that each objective is well-optimized by at least one member of the set. For example, as you rightly point out, the original wording of the statement “we may have T pathogens and aim to identify a set of K < T antibiotics such that at least one antibiotic can be used to treat each pathogen” is potentially misleading. To more accurately convey our intent, we plan to replace this with a revised version such as:

“In drug discovery, we may be faced with the challenge of designing antibiotics for T different pathogens. Developing T individual antibiotics (a separate antibiotic specifically designed to treat each pathogen) is not ideal, as synthesizing and testing each antibiotic is highly expensive (a smaller number of antibiotics is ideal). However, we still want to design more than one antibiotic, as a single antibiotic is highly unlikely to be able to treat all T pathogens effectively. We therefore aim to identify a smaller set of K < T antibiotics, such that each pathogen can be effectively treated by at least one of the K selected antibiotics. The goal is not for every antibiotic to treat every pathogen, but for the set of K antibiotics as a whole to provide ‘coverage’ of the T pathogens, such that if we encounter any one of the T pathogens, we will be able to treat it with at least one of the antibiotics in the optimized set.”

Clarification of Figure 1

We appreciate the feedback that the illustrative example in Figure 1 could be clearer. You are absolutely right that in this simple 2D example with K = 2 and T = 3, the two square points that yield good coverage could indeed be recovered by searching over the Pareto frontier for a pair that maximizes the coverage score. However, in the real-world tasks we consider, where the objectives are high-dimensional, black-box, and highly non-convex, identifying such a set of K coverage-maximizing points is non-trivial and likely cannot be accomplished by naïvely sampling the Pareto frontier.

Additionally, as you note, many standard multi-objective methods do output multiple Pareto-optimal points. However, these points are typically meant to serve as alternative trade-off options, from which a practitioner selects one solution to apply, depending on their preferences. The methods do not aim to maximize the combined utility of the full set; instead, they ensure that each individual point lies on the Pareto frontier. In contrast, our formulation explicitly seeks a set of solutions that work together: the goal is to maximize aggregate coverage across all objectives, not to offer mutually exclusive alternatives.

The star point in Figure 1 is meant to illustrate the limitations of single-solution approaches, namely, that no individual solution can provide high performance across all objectives when trade-offs are extreme. This motivates the need for multiple jointly selected solutions, each covering different subsets of objectives.

We will revise the text to clarify this distinction, explicitly explain how coverage optimization differs from Pareto-based approaches (even those that return multiple points), and emphasize that Figure 1 is a simplified schematic meant to illustrate intuition, not a literal depiction of how Pareto methods operate.

Related Work

Thank you for pointing us to Malkomes et al. (ICML 2021) [1] and the related diversity-in-objective-space paper. We agree that while their definition of “coverage” differs from ours (input vs. objective space), their use of an Expected Coverage Improvement (ECI) acquisition function is relevant and valuable to acknowledge. We will update the Related Work section to include discussion of these papers.

Minor Comments

Thank you for the detailed comments on writing and presentation. These comments will be very helpful in polishing the camera-ready version.

Questions

Q1: Have you considered using a lazy greedy approach to accelerate selection?

While we did not use a lazy greedy approach in this version of MOCOBO, we plan to explore this in future work. Thank you for the suggestion!

Q2: Are you really doing 20k function evaluations? Are you training GPs on that scale?

To clarify: while some of our large-scale tasks involve collecting up to tens or even hundreds of thousands of data points (e.g., the peptide design task), we do not train exact Gaussian processes on the full dataset, as this would be computationally infeasible.

Instead, as described in Appendix A.7.2, we use PPGPR, a scalable approximate GP surrogate model. This setup allows us to handle large datasets efficiently while maintaining strong predictive performance. Additionally, at each step of optimization, we update the surrogate using only a small subset of 1000 data points: this includes the most recent batch of observations and the highest-performing points collected so far. This ensures that the model remains up to date without retraining on the full history on each step. Full details can be found in Appendix A.7.2.