Bounce: Reliable High-Dimensional Bayesian Optimization for Combinatorial and Mixed Spaces
Authors: Leonard Papenmeier
Luigi Nardi
Matthias Poloczek
Venue: NeurIPS 2023 poster (2023)
Keywords: Bayesian optimization
global optimization
Gaussian process
combinatorial optimization
high-dimensional
TL;DR: We propose a Bayesian optimization algorithm for combinatorial, mixed, and continuous spaces that is robust and gives state-of-the-art performance on a wide set of benchmarks.

Abstract: Impactful applications such as materials discovery, hardware design, neural architecture search, or portfolio optimization require optimizing high-dimensional black-box functions with mixed and combinatorial input spaces.
While Bayesian optimization has recently made significant progress in solving such problems, an in-depth analysis reveals that the current state-of-the-art methods are not reliable. 
Their performances degrade substantially when the unknown optima of the function do not have a certain structure. 
To fill the need for a reliable algorithm for combinatorial and mixed spaces, this paper proposes Bounce that relies on a novel map of various variable types into nested embeddings of increasing dimensionality.
Comprehensive experiments show that Bounce reliably achieves and often even improves upon state-of-the-art performance on a variety of high-dimensional problems.

Bounce: Reliable High-Dimensional Bayesian Optimization for Combinatorial and Mixed Spaces
Leonard Papenmeier
Luigi Nardi Lund University, Stanford University, DBtune juigi.nardia 1tth.se
Lund University
leonard. pap.ardee.peeeemeeeeeeees. se
Matthias Poloczek
Amazon
San Francisco, CA 94105, USA
matpol@amazon.com com
Abstract
Impactful applications such as materials discovery, hardware design, neural ar- chitecture search, portfolio optimization require optimizing high-dimensional black-box functions with mixed and combinatorial input spaces. While Bayesian optimization has recently made significant progress in solving such problems, an in-depth analysis reveals that the current state-of-the-art methods are not reliable. Their performances degrade substantially when the unknown optima of the func- tion do not have a certain structure. To fill the need for a reliable algorithm for combinatorial and mixed spaces, this paper proposes Bounce that relies on novel map of various variable types into nested embeddings of increasing dimensional- ity. Comprehensive experiments show that Bounce reliably achieves and often even improves upon state-of-the-an performance on a variety of high-dimensional problems.
Introduction
Bayesian optimization (BO) has become 'go-to' method for optimizing expensive-to-evaluute black-box functions [27, 83] that have numerous important applications, including hyperparame- ter optimization for machine learning models [9, 27], portfolio optimization finance [7]. chem- ical engineering and materials discovery [13. 15, 26, 33, 37, 38, 40, 57, 68, 71, 81]. hardware de- sign [22, 36, 51], or scheduling problems [42]. These problems are challenging for variety of reasons. Most importantly, they may expose hundreds of tunable parameters that allow for granular optimization of the underlying design but also lead to high-dimensional optimization tasks and the 'curses of dimensionality" [10, 63]. Typical examples are drug design [53, 76] and combinatorial testing [55]. Moreover, real-world applications often have categorical or ordinal tunable parameters, in addition to the real-valued parameters that BO has traditionally focused [10, 27, 70J. Recent efforts have thus extended BO to combinatorial and mixed spaces. Casmopol itan of Wan al. [79] uses trust regions (TRs) to accommodate high dimensionality, building upon prior work of Eriksson et al. [25] for continuous spaces. COMBO of Oh et al. [56] constructs a surrogate model based on a combinatorial graph representation of the function. Recently, Deshwal al. [21] presented BODi embedding and showed that it outperforms the prior
that employs a novel type of dictionary-bast work. However, the causes for BODi 's excellent performance not yet well-understood and require a closer examination. Moreover, the ability of methods for mixed spaces to scale to higher dimen- sionalities trails behind BO for continuous domains. In particular, Papenmeier et al. [60] showed that nested embeddings allow BO to handle a thousand input dimensions, thus outperforming vanilla TR-based approaches and raising the question of whether similar performance gains are feasible for combinatorial domains.
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

In this work, we assess and improve upon the state-of-the-arr in combinatorial BO. In particular, we make the following contributions:
1. We conduct an in-depth analysis of two state-of-the-art algorithms for combinatorial BO, COMBO [56] and BODi [21]. The analysis reveals that their performances often degrade consid- erably when the optimum of the optimization problem does not exhibit a particular structure common for synthetic test problems.
2. We propose Bounce (Bayesian optimization using increasingly high-dimensional combinatorial and continuous embeddings), a novel sigh-dimensional Bayesian optimization (HDBO) method that effectively optimizes over combinatorial, continuous, and mixed spaces. Bounce leverages parallel function evaluations efficiently and uses nested random embeddings to scale to high- dimensional problems.
We provide a comprehensive evaluation on a representative collection of combinatorial, contin- uous, and mixed-space benchmarks, demonstrating that Bounce is on par with or outperforms state-of-the thr methods.
2 Background and related work
Bayesian optimization. Bayesian optimization aims to find the global optimum x² E X of a black-box function f : X' + R, where X' is the D-dimensional search space imput space. Throughout this paper, we consider minimization problems, i.e., we aim to find x² E X such that f(x") S f(x) for all æ E X. The search space X' may contain variables of different types: continuous, categorical, and ordinal. We denote the number of continuous variables in X' by Noom and the number of combinatorial variables by Nommb = Teat + Foord D Acontt where we denote the number of categorical variables by Tlaa and the number of ordinal variables by Yoord-
Combinatorial domains. Extending BO to combinatorial spaces is challenging, for example, be- cause the acquisition function is only defined at discrete locations the dimensionality of the space grows drastically when using one-hot encoding for categorical variables. Due to its numerous ap- plications, combinatorial BO has received increased attention in recent years. BOCS [6] handles the exponential explosion of combinations by only modeling lower-order interactions of combinatorial variables and imposing sparse prior the interaction terms. COMBO [56] models each variable as a graph and uses the graph-Cartese performance on categorical problems in Appendix H. CoCaBo [67] combines multi-armed bandits and BO to allow for optimization in mixed spaces. It uses two separate kernels for continuous and combinatorial variables and proposes weighted average of product and sum kernel to model mixed spaces. Liu and Wang [47] show that certain, possibly combinatorial functions can be mod- eled by parametric function approximators such as neural networks or random forests. As general method to optimize the acquisition function with gradient-based methods, Daulton et al. [19] pro- pose probabilistic reparametrization.
product to represent the search space. We revisit COMBO's s
High-dimensional continuous spaces. Subspace-based methods primarily used for continu- ous spaces. Wang et al. [84] propose REMBO for HDBO in continuous spaces using Gaussian random projection matrices. REMBO suffers from distortions and projections outside the search domains that the corrections of Binois et al. [11, 12] address. The HeSBO algorithm of Nayebi et al. [52] avoids the need for corrections by using the CountSketch embedding [87]. Alebo of Letham et al. [46] builds upon REMBO, learning suitable corrections of distortions. TuRBo [25] is a method that operates in the full-dimensional input space X', relying on trust region (TR) focus the search on promising regions of the search space. BAXUS of Papenmeier et al. [60] combines the trust region approach of TuRBO with the random subspace idea of HeSBO. BAXUS uses novel family of nested random subspaces that exhibit better theoretical guarantees than the CountSketch embedding. While BAxUS handled a 1000D problem, it only considers continuous problems and cannot leverage parallel function eval- uations. GTBO [35] assumes the existence of an axis-aligned active subspace. The algorithm first identifies "active" variables and optimizes in the full-dimensional space by placing separate strong length scale priors onto active and inactive variables. Another line of recent approaches employs Monte-Carlo tree search (MCTS) to reduce the complexity of the problem. Wang et al. [82] use MCTS to learn a partitioning of the continuous search space to focus the search on promising re- gions the search space. Song et al. [74] use a similar approach, but instead of learning promising regions in the search space, they assume an axis-aligned active subspace and use MCTS to select important variables. Linear embeddings and random linear embeddings [14, 46, 52, 60, 84] require little no training data to construct the embedding but assume lincar subspace.

---

Algorithm The Bounce algorithm
Input: initial target dimensionality duni. evaluation budget batch size B, evaluation budget to input dimensionality mp. # new bins added per dimension b, number design experiment
(DOE) points Ninit
Output: optimizer x" E arg min (()).
I: 0, d dinit
2: < ADIUSTBINS(), dinit, mp)
D Section
homp dua
3: m,
4: S < INITIALEMEEDDINGGOO Teonts Theats Thcombs Thbin) 1(5 6: for -1,,..,, do
Section 3.1
5: D < ((=))))(--(())))))
P Sample and evaluate initial points.
7:
Lcont 6 6.88, Leomb f- min Reomb)
8:
while Leout > Lount Leomb L Leamb do
9: Find B candidates according to Sec. 3.2
10:
Evaluate B and update D: 'D D U {(z) (2))
11:
Update Teemu and Teomb
D Section
12:
d < D then
i i+1
Increase index for target dimensionality. Section 3.1
14:
S < INCREASHEMBEDDINISS b)
< # target variables in m.
bm
16:
dui (1- (1-(6+1)
17:
else
Reset D by resampling and evaluating new initial points
19:
Resample S, reset Leont and Leomb, j < + Minit
Combinatorial hign-dimensional domains. These works optimize black-box functions defined combinatorial mixed space with dozens of input variables. Thebelt et al. [77] ensemble kernel to model the Gaussian process (GP) prior and derive a formulation of the upper- confidence bound (UCB) that allows it to be optimized globally and to incorporate constraints. RDUCB [88] relies on random additive decompositions of the GP kernel to model the correlation between variables. Casmopolitan itaannttnn [79] follows TuRBO in using TRs to focus the search promis- ing regions of the search space and uses the Hamming distance to model TRs for combinatorial variables. For mixed spaces, Casmopolitan uses interleaved search and models continuous and categorical variables with two separate TRs. Kim et al. [44] use random projection matrix to ap- proach combinatorial problems continuous embedded subspace. When evaluating point, their approach projects the continuous candidate point to the high-dimensional search space and then rounds the next feasible combinatorial solution. Deshwal al. [20] propose two algorithms for permutation spaces, which occur in problems such as compiler optimization [36] and pose special challenges due to the superexponential explosion of solutions. BODi [21] proposes an embedding type based on a dictionary of anchor points in the search space. The pairwise Hamming distances between the point and each anchor point a, in the dictionary represent a point in the search space. The anchor points in the dictionary change at each iteration of the algorithm. They are sampled from the search space to cover a wide range of 'sequencies', i.c., the number of changes from 0 to 1 (and vice versa) in the binary vector. The authors hypothesize that the diverse sampling procedure leads to BODi's remarkable performance in combinatorial spaces with up to 60 dimensions. To our knowledge, BODi is the only other method combining embeddings and combinatorial spaces. We show Section 4.6 that BODi's reported good performance relies on an artificial structure of the optimizer x² and that its performance degrades considerably when this structure is violated.
tree-
The Bounce algorithm
To overcome the aforementioned challenges in HDBO for real-world applications. we propose Bounce, a new algorithm for continuous, combinatorial, and mixed spaces. Bounce uses GP [85] surrogate in lower-dimensional subspace, the target space, that is realized by partitioning input variables into "bins', the so-called target dimensions. Bounce only bins variables of the same type (categorical, binary, ordinal, and continuous). When selecting points to evaluate, Bounce sets

---

all input variables within the same bin to a single value. It thus operates in a subspace of lower dimensionality than the imput space and, in particular, maximizes the acquisition function in a sub- space of lower dimensionality. Bounce iteratively refines its subspace embedding by splitting bins into smaller bins, allowing for more granular optimization at the expense of higher dimensional- ity. Note that by splitting up bins, Bounce asserts that observations taken in earlier subspaces are contained in the current subspace; see Papenmeier et al. [60] for details. Thus, Bounce operates in a series of nested subspaces. It uses novel TR management to leverage batch parallelism efficiently, improving over the single point acquisition of BAxUS [60].
The nested subspaces. To model the GP in low-dimensional subspaces, Bounce leverages BAxUS' family of nested random embeddings [60]. In particular, Bounce employs the sparse count-sketch embedding [87] in which each input dimension is assigned to exactly one target dimension. When increasing the target dimensionality, Bounce creates b new bins for every existing bin and re- distributes the input dimensions that had previously been assigned to that bin across the now + 1 bins. Bounce allocates an individual evaluation budget m, to the current target space X; that is pro- portional to the dimensionality of X, When the budget for the current target space is depleted, and Bounce has not found a better solution, Bounce will increase the dimension of the target space until it reaches the input space of dimensionality D. Let dinit denote the dimensionality of the first target space, i.e., the random embedding that Bounce starts with. Then Bounce has to increase the target dimension [log++1 didiit =: k-times to reach the input dimensionality D. After calculating k, Bounce re-sets the split factor b such that the distance between the predicted final target dimension- ality d2 = dinit (b + 1)x and the input dimensionality D is minimized: b = [loga (D/dinit) - 11, where [x]] denotes the integer closest to x. This ensures that the predetermined evaluation budget for each subspace will be approximately proportional to its dimensionality. This contrasts to BAXUS [60] that uses a constant split factor and adjusts the initial target dimensionality dinit- The evaluation
budget m, for the i-th subspace X; is mi dl-ic-(- 16+1)
hmp:d:
where mp is the budget until D
is reached and is the maximum number of bins added per split.
Trust region management. Bounce follows TuREO [25] and Casmopolitan [79] in using trust regions (TRs) to efficiently optimize over target spaces of high dimensionality. TRs allow focus- ing on promising regions of the search space by restricting the next points to evaluate
region centered at the current best function value [25]. TR-based methods usually expand their TR if they find better points and conversely shrink it if they fail to make progress. If the TR falls below the threshold given by the base length, TuRBO and Casmopolitan restart with a fresh TR elsewhere. Casmopolitan [79] uses different base lengths for combinatorial and continuous variables. For combinatorial variables, the distance to the currently best function value is defined in terms of the Hamming distance, and the base length is an integer. For continuous variables, Casmopol i de- fines the base length in terms of the Euclidean distance, i.c., a real number. Similarly, Bounce has separate base lengths contt and commb for continuous and combinatorial variables but does not fix the factor by which the TR volume is increased or decreased upon successes or failures. Instead, the factor is adjusted dynamically so that the evaluation budget m, for the current target space X, is adhered to. This design is crucial to enable batch parallelism, as we describe next.
Batch parallelism. We allow Bounce to efficiently evaluate batches of points in parallel by using a scalable TR management strategy and g-expected improvement (qEI) [5, 80, 86] as the acquisition function for batches of size B 1. When Bounce starts with fresh TR, we sample Aimit initial points to initialize the GP, using a Sobol sequence for continuous variables and uniformly random values for combinatorial variables.
The TR management strategy of Bounce differs from previous strategies [25, 60, 65, 79] in that it uses dynamic factor determine the TR base length. Recall that Bounce shrinks the TR if it fails to make progress and starts a fresh TR if the TR falls below the threshold given by the base length. Bounce's 's rule based on the idea that the minimum admissible TR base length should be reached when the current evaluation budget is exhausted. If one employed the strategies of TuRBo [25], Casmopolitan [79], or BAXUS [60] for larger batch sizes B and Bounce's nested subspaces, then one would spend a large part of the evaluation budget in early target spaces. For example, suppose a continuous problem, the common values for the initial, minimum, and maximum TR base length, and the constant shrinkage factor of [65]. Then, such a method has to shrink the TR base length at least seven times (i.e., evaluate ʃ 7B-times) before it would increase the dimensionality of the target space. Thus, the method would risk depleting its budget before reaching target space suitable

---

for the problem. On the other hand, we will see that Bounce chooses an evaluation budget that is smaller in wiwdimmnsional target spaces and higher for later target spaces of higher dimensionality, Considering a 1000-dimeniional problem with an evaluation budget of 1000, it uses only 3, 12, and 47 samples for the first three target spaces of dimensionalities 2, 8, and 32.
Bounce's s strategy permits flexible TR shrinkage factors and base lengths, i.c., TR base lengths to vary within the range [Lmnnn Lmax] This allows Bounce comply with the evaluation budget my for the current target space XG Suppose that Bounce has evaluated j batches of B points each since it last increased the dimensionality of the target space, and let L, denote the current TR base length. Observe that hence m - j evaluations remain for X, Then Bounce sets the TR base length to L++1 == X Ly with X, < 1. if it found new best point whose objective value improves upon the incumbent by least e. We call this 'success'. Otherwise, Bounce observes 'failure' and sets Li+1 = X3 Lij. The rationale of this rule is that if the algorithm is in iteration j and only observes failures subsequently. then we apply this factor (m, - (i--jBB)timeee which the remaining number of function evaluations in the current subspace Xz Hence, the last batch of the i-th target space X, will have the minimum TR base length, and Bounce will increase the target dimensionality afterward. If the TR is expanded upon a 'success', we need to adjust ' not to use more than the allocated number of function evaluations in a target space. At each iteration, we therefore set adjustment factor Xj = (Lmin/ Lgg///mm 1/1m Note that Xj remains unchanged under this rule unless the TR expanded in the previous iteration.
The kernel choice. To harvest the sample efficiency of a low-dimensional target space, we would like to combine categorical variables into a single bin, even if they vary in the number of categories. This is not straightforward. For example, note that the popular one-hot encoding of categorical variables would give rise to multiple binary input dimensions, which would not be compatible with the above strategy binning variables to form nested subspaces. Bounce overcomes these obstacles and allows variables of the same type to share a representation in the target space. We provide the details in Sect. 3.1.
For the GP model, we use the CoCaBo kernel [67]. In particular, we model the continuous and combinatorial variables with two separate 6/2 1/2-Matérn kernels where we use automatic relevance determination (ARD) for the continuous variables and share one length scale for all combinatorial variables. Following Ru al. [67], we use mixture of the sum and the product kernel:
k( x²) = of (pkcabb (Famb: <mm ) Kemi (Zeny a'mm) + (1- (1-p) (Kemb (Zembr F(emb) + (Fenss F'mm))),
where Fen and Temb are the continuous and combinatorial variables in æx, respectively, and of the signal variance. The trade-off parameter & E [0, 1] is learned jointly with the other hyperparameters via likelihood maximization. See Appendix G.2 for additional details.
Algorithm gives a high-level overview of Bounce. In Appendix A, we prove that Bounce con- verges to the global optimum under mild assumptions. We now explain the different components of Bounce in detail.
3.1 The subspace embedding of mixed spaces
Bounce supports mixed spaces of four types of input variables: categorical. ordinal, binary, and continuous variables. We discuss binary and categorical variables separately because we model them differently. The proposed embedding maps only variables of single type to each "bin', i.e., no target dimension of the embedding has variables of different types. Thus, target dimensions are homogeneous in this regard. Note that the number of target dimensions of each type is implied by the current bin size of the embedding that may grow during the execution. The proposed embedding can handle categorical or ordinal input variables that differ in the number of discrete values they can take.
Continuous variables. As common in BO, we suppose that each continuous variable takes values in bounded interval. Thus, we may normalize each interval to [-,, The embedding of contin- uous variables, i.e., input dimensions, follows BAxUS [60]: each input dimension D, is associated with a random sign E [-1,+1) and one or multiple input dimensions can be mapped to the same target dimension of the low-dimensional embedded subspace. Recall that Bounce works on the low-dimensional subspace and thus decides an assignment v for every target dimension d, of the embedding. Then, all input variables mapped to this particular target dimension are set to this value tiy-
5

---

d d.
di
d di
d
Uy
U3
0
0
0 0
f// = 2
0
0
0
0
0g 3
0
Figure 1: The mapping (or binning) of categorical and ordinal variables. Suppose that variable uk has two categories and that 2y has three categories. Both are mapped to the target dimension d, that has cardinality 3 = max [2, 3}. While the mapping of to di is straightforward bijection, has fewer categories than dj. Thus, the mapping of uk to d repeats label Ordinal variables are mapped similarly.
Binary variables. Binary dimensions are represented by the values 1 and +1. Each input dimen- sion D; is associated with random sign 8; f-1,+11 +1), and the subspace embedding may map one or more binary input dimensions to the same binary target dimension. While the embedding for binary and continuous dimensions is similar, Bounce handles binary dimensions differently when optimizing the acquisition function.
Categorical and ordinal variables. Categorical variables that differ in the number of categories may be mapped to the same target dimension (bin). Suppose that the categorical variables u1.... Uy with cardinalities c1...., ..... Ct mapped to single bin that is associated with the target dimension dj of the subspace embedding. Then dj is of categorical type and has max(c; Cmax distinet categories, that is, its cardinality is given by the maximum cardinality of any variable mapped to it. Thus, the bin dj can represent every category of these input variables.
S S 4
Suppose that Bounce assigns the category l {L..... Cmax to the categorical bin (target dimen- sion) dy. We transform this label to a categorical assignment to each input variable U1,.... Ue, setting oy [k cm.// Recall that Bounce may split up bins, i.c., target dimensions, to increase the dimensionality of its subspace embedding. nality of the parent bin. This allows this point. Analogously to the random sign for binary variables, we randomly shuffle the categories before the embedding. This reduces the risk of Bounce being biased towards specific structure of the optimizer (see Appendix E).
such an event, every derived bin inherits the cardi- retain any observations the algorithm has taken up to
We treat ordinal variables as categorical variables whose categories correspond to the discrete values the ordinal variable can take. For the sake of simplicity, we suppose here that an ordinal variable is has range {1, .. ) and G 2 2 for all E {I....,E}. Figure I shows examples of the binning of categorical and ordinal variables.
3.2 Maximization of the acquisition function
We use expected improvement (EI) [43] for batches of size B 1 and qEI [5, 80, 86] for larger batches. We optimize the El using gradient-based methods for continuous problems and local search for combinatorial problems. We interleave gradient-based optimization and local search for func- tions defined over a mixed space; see Appendix G. for details.
4 Experimental evaluation
We evaluate Bounce empirically various benchmarks whose inputs combinatorial, con- tinuous, or mixed spaces. The evaluation comprises the ofate-of-the-art algorithms BODi [21], Casmopolitan [79), COMBO [56], SMAC [41], and RDUCE [88], using code provided by the authors. We also report Random Search [8] as a baseline. For categorical problems, COMBO's 's implementa- tion suffers from bug explained in Appendix H.2. We report the results for COMBO with the correct benchmark implementation as "COMBO (fixed)".
The experimental setup. We initialize every algorithm with five initial points. The plots show the performances of all algorithms averaged over 50 repetitions except BODi, which has 20 repetitions due to resource constraints caused by its high memory demand. The shaded regions give the stan- dard error of the mean. We use common random seeds for all algorithms and for randomizing the

---

50D-LABS
50D-LABS (publ.)
Bounce
Objective valuee
Casmopolitan
-2
BODI
Search
-3
COMBO
ROUCB
SMAC
200
400
200
400
Iteration
Iteration
Figure 2: The 50D low-autocoorreaatio binary sequence problem. Bounce finds the best solutions,
followed by COMBO.
125D-Cluse ClusterExpaanioonn
1c14
DDDDClusterEExansion (publ.)
-2.7
Bounce
-2.7
Casmopalitan
BODi
Objective
-2.76
Random Search
-..77
COMBO
RDUCB
-2.7
SMAC
200
400
200
400
Iteration
Iteration
Figure 3: The 125D weighted ClusterExpansion maximum satisfiability problem. We plot the total negative weight of clauses. Bounce produces the best assignments.
benchmark functions. We run all methods for 200 function evaluations unless stated otherwise. The Labs (Section 4.1) and MaxSat125 (Section 4.2) benchmarks are rum for 500 evaluations.
The benchmarks. 125D ClusterExpansi on [3, 4], 60D MaxSAT60 [21, 56], 25D PestControl, 53D Ackley53. and 25D Contamination [6, 39, 56]. Due to space constraints, we moved the results for the MaxSAT60, Contamination, and Ackley53 benchmarks to Appendix B.1. For each benchmark, results for the originally published formulation and for a modification where we move the opti- mal point to random location. The randomization procedure is fixed for each benchmark for all algorithms and repetitions. For binary problems, we flip each input variable independently with probability 0.5. For categorical problems, randomly permute the order of the categories. We motivate this randomization in Section 4.6.
The evaluation uses seven established benchmarks [21]: 53D SVM, 50D LABS,
report
4.1 50D Low-Autoccrrrelttion Binary Sequences (LABS)
LABS has n 50 binary dimensions. It has important applications in communications engineering and mathematics; see [58] for details. LABS is a hard combinatorial problem and currently solved exhaustive search. The goal is to find a sequence æ E (-1, with a maximum merit factor >²
F(x) = > where E(x) = Σx=1 k= C (æ) and Cx(x) = Σxi Z,T4+ for k = O,..., n - 1 are the autocorrelations of æ [58]. Figure 2 summarizes the performances. We observe that Bounce outperforms other algorithms on the benchmark's original and randomized versions.
4.2 Industrial Maximum Satisfibbility: 125D ClusterExpansior benchmark
MaxSat is notoriously hard problem for which various approximations and exact (exponential time) algorithms have been developed; see [32, 61] for an overview. We evaluate Bounce and the other algorithms the 125 -25-dimensional ClusterExpansion benchmark, real-world MaxSAT instance with many applications in materials science [2]. Unlike the MaxSAT60 benchmark (see Appendix 1.1.). ClusterExpansion is not crafted benchmark, and its optimum has synthetic structure [1, 3). We treat the MaxSat problems as black-box problems; hence, algorithms do not have access the clauses, and cannot use the usual algorithms.
Figure 3 shows the total negative weight of the satisfied clauses as function of evaluations. We cannot plot regret curves since the optimum is unknown [4]. We observe that Bounce finds better solutions than all other algorithms. BODi is the only algorithm for which we observe sensitivity to
7

---

2SD-PestControl
25D-PestControl (publ.)
Bounce
Objectivve value
BODI
166
Search
Casmopelltan
COMBO
14-
COMBO (fixed)
RDUCB
SMAC
100 Iteration
200
1.00
Iteration
Figure 4: The 25D categorical pest control problem. Bounce obtains the best solutions, followed by Casmopolitan. BODi's 's performance degrades significantly when shuffling the order of categories.
53D-SVM
53D-SVM (publ.)
x 10
Bounce
3 10-1
BODi
Objective
Random Search
Casmopolitan
x 10-1
RDUC8
SMAC
100 Iteration
200
100
200
Iteration
Figure 5: The 53-dimensional SVM benchmark. Bounce, BODi, and Casmopolitan achieve compa- rable solutions.
the location of the optimal assignment: for the published version of the benchmark, BODi quickly jumps to a moderately good solution but fails to make further progress.
4.3 25D Categorical Pest Control
PestControl is more complex version the Contaminat intaaminiiinn on benchmark and has 25 categorical
variables with five categories each [56]. The task is to select out of five actions {1, 2.... 5) each of 25 stations to minimize the objective function that combines total cost and a measure of the spread of the pest. We note the setting æ (5, 5, ......,) 5) achieves a good value of 12.57, while the best value found in our evaluation is 12.07 is x (5, 5,..,,.. 5, 1) and thus has a Hamming distance Figure 4 summarizes the performances of
of one. The random seed used in experiments is the algorithms. Bounce is robust to the location of the global optimum and consistently obtains the best solutions. In particular, the performances of COMBO and BODi depend on whether the optimum has a certain structure. We discuss this issue in detail in Appendix H.l.
4.4 SVM 53D AutoML task
In the SVM benchmark, we optimize over a mixed space with 50 binary and continuous parame- ters to tune an s-support vector regression (SVR) model [72]. The 50 binary parameters determine whether to include or exclude an input feature from the dataset. The 3 continuous parameters cor- respond to the regularization parameter C, the kernel width 7, and the = parameter of the =-SVR model [72]. Its root mean squared eтror on held-out dataset gives the function value. Figure 5 sum marizes the performances of the algorithms. We observe that Bounce, BODi, and Casmopolita tan achieve comparable solutions. BODi performs slightly worse if the ordering of the categories is shuffled and slightly better if the optimal assignment to binary variables one. COMBO does not support continuous variables and thus was omitted.
4.5 Bounce's efficacy for batch acquisition
We study the sample efficiency of Bounce when it selects a batch of B points in each iteration evaluate parallel. Figure shows the results for B 1,3, 10, and 20, where Bounce was run for min 200 B) function evaluations. We configure Bounce to reach the input dimensionality after 100 evaluations for B 1,3, and after 25B for = 10, 20. We observe that Bounce leverages parallel function evaluations effectively: it obtains a comparable function value considerably smaller number of iterations, thus saving wall-clock time for applications with time-consuming function evaluations. We also studied batch acquisition for continuous problems

---

250-PestControl
6OD-MaxSAT60
2SD-Contamination
23.0
Cojective
166
22.5
14
-150
22.0
12
-200
21.5
bel
Bounce, bm3
53D-SVM
S3DDAckeey53
S0D-LABS
b-5
Bounce, b-20
velue veluv
-2
colctive
0.3
-3
0
100
2000
0
100
200
1.00
200
batch wve Naliaations
batch evaluations
batch evaluations
Figure 6: Bounce benefits from the batch acquisition that allows parallelizing function evaluations.
We show the best function values obtained after each batch for batch sizes 1,3, 5, 10, and 20.
and found that Bounce also provides significant speed-up. Due to space constraints, we deferred the discussion to Appendix C.
4.6
The sensitivity of BODi and COMBO to the location of the optima
The empirical evaluation reveals that the performances of BODi [21] and COMBO [56] are sensitive to the location of the optima. Both methods degrade on least one benchmark when the optimum is moved to a randomly chosen point. This is particularly unexpected for categorical variables where moving the optimum to a random location is equivalent to shuffling the labels of the categories of each variable. Such a change of representation should not affect the performance of an algorithm.
BODi is more susceptible the location of the optimizer than COMBO. The performance of COMBO degrades only on the categorical PestControl benchmark, whereas BODi degrades on five out of seven benchmarks. Looking closer, we observe that BODi's 's performance degradation is particularly large for synthetic benchmarks like Ackley53 and MaxSAT60, where setting all variables to the same value is optimal. Figure summarizes the effects of moving the optimum on BODi. Due to space constraints, we moved the details and a discussion of categorical variables to the appendix. Similarly, setting all binary variables of the SVM benchmark to one produces a good objective value. It is not surprising, given that the all-one assignment corresponds to including all features previously selected for the benchmark because of their high importance.
We show in Appendix H.1 that BODi adds a point to its dictionary with zero Hamming distance to an all-zero all-one solution, with probability that increases with the dictionary size. Deshwal et al. [21, p. 7] reported that BODi's performance 'tends to improve' with the size of the dictionary. Moreover, BODi samples new dictionary in each iteration, eventually increasing the chance of having such a point in its dictionary. Thus, we hypothesize that BODi benefits from having a near- optimal solution in its dictionary, likely for all-zero all-one solutions. For COMBO, Figure 4 shows that the performance on PestControl degrades substantially if the labels of the categories are shuffled. Then COMBO's s sample-efficiency becomes comparable to Random Search.
Discussion
BO in combinatorial spaces has many exciting and impactful applications. Its applicability to real- world problems, such LABS that defy closed-form solution, makes it valuable tool for practi- tioners. Our empirical evaluation reveals that state-of-the-an methods fail to provide good solutions reliably. In particular, it finds that BODi and COMBO, which performed best in recent publications, are sensitive to the location of the optimizer. We identified design flaws in BODi and an implementation bug in COMBO as the root causes of the performance degradations.
The proposed Bounce algorithm reliable for high-dimensional black-box optimization combi- natorial, continuous, and mixed spaces. The empirical evaluation demonstrates that Bounce reliably outperforms the state-of-the-an on a diverse set of problems. Using a novel TR management strat- egy, Bounce leverages parallel evaluations of the objective function to improve its performance. We

---

25D-PestControl
SOD-LABS
53D-SVM
53D-AckleyS3
18-
Osgecde vakbe
16-
0.4
0.3
-3
0.2
100
200
100
200
100
200
0
100
200
125D- ClusterExpaasion
60D-MaxxAAT6
Z5DDConntminnttoo
Otpertiv value
23.0
-2.72
BOS (hand))
22.5
-2.74
- BOOK
-2.75
-200
100
200
itoration
Iterattion
Figure 7: BODi's 's performance degrades on five out of seven benchmarks when randomizing the location of the optimal solution: BODi the default version (orange) of the benchmark and on the modified version (green, dashed) where the optimum was moved.
anticipate headroom by tailoring the modeling of combinatorial objects, e.g., arising in the search for peptides or materials discovery [28, 34, 37, 55, 75, 78). Here seems particularly interesting to incorporate prior belief on the importance of decision variables while maintaining the overall scala- bility. Moreover, extending the present work to black-box constraints [24, 30), multiple objectives, and multiple information sources [18, 31, 62] will considerably expand the applicable use cases.
Limitations. Bounce is not designed to handle noisy evaluations of the objective function. While seems straightforward extend Bounce to handle noisy evaluations, c.g., by using a Gaussian process with a noise term and acquisition functions that account for noise [70], we leave this for future work. Moreover, in applications where the categorical or ordinal variables vary substantially in the number of values they can take, there may be better ways to 'bin" them.
Societal impact. Bayesian optimization has recently gained wide-spread popularity for tasks in drug discovery [53], chemical engineering [15, 33, 38, 68, 71], materials discovery [28, 34, 37, 59, 75, 78]. aerospace engineering [6, 45, 49], robotics [16, 17, 48, 50, 64), and many more. This highlights the Bayesian optimization communitys progress toward providing reliable 'off-the- shelf optimizer. However, this promise is not yet fulfilled for the newer domain of mixed-variable Bayesian optimization that allows optimization over hundreds of 'tunable levers', some of which are discrete, while others are continuous. This domain is of particular relevance for the tasks above. Bounce's 's ability incorporate more such levers in the optimization significantly impacts the above practical applications, allowing for more granular control of a chemical reaction or a processing path, to give some examples. The empirical evaluation shows that the performance of state-of-the- art methods is highly sensitive to the location of the unknown global optima and often degenerates drastically, thus putting practitioners risk. The proposed algorithm Bounce, however, achieves robust performance over broad collection of tasks and thus will become 'goto" optimizer for practitioners other fields. Therefore, we open-source the Bounce code.
Acknowledgments and Disclosure of Funding
Leonard Papenmeier and Luigi Nardi were partially supported by the Wallenberg AI, Autonomous Systems and Software Program (WASP) funded by the Knut and Alice Wallenberg Foundation. Luigi Nardi was partially supported by the Wallenberg Launch Pad (WALP) grant Dnr 2021.0348 and by affiliate members and other supporters of the Stanford DAWN project Ant Financial, Face- book, Google, Intel, Microsoft, NEC, SAP, Teradata, and VMware. The computations were en- abled by resources provided by the National Academic Infrastructure for Supercomputing in Swe- den (NAISS) at the Chalmers Centre for Computational Science and Engineering (C3SE) and the National Supercomputer Centre Linköping University, partially funded by the Swedish Research Council through grant agreement no. 2022-06725.
10

---

References
[1] Wehpage of the Fourth Max-SAT Evaluation. http://wwx. www.:w.www................. www.:aawwtt.d.ccc..mmm wuw.:tt/lllclaa/////// /iddeeex.nnn////nee/ herprabbbbtttte--bebcnaaararss php?disp=submitten benchmarks. Last access: 2023-05-02.
[2] MaxSAT Evaluation 2018 Solver and Benchmark Descriptions. 2018. URL http://hdl.
handle met/110338/277139 Publisher: Department of Computer Science, University of Helsinki.
[3] André Abramé and Djamal Habet. AHMAXSAT: Description and evaluation of branch and bound Max-SAT solver. Journal on Satisfibiliiy, Boolean Modeling and Computation, 9(1) 89-128, 2014.
[4] Fahiem Bacchus, Matti Järvisalo, and Ruben Martins. MaxSAT Evaluation 2018: New devel- opments and detailed results. Journal on Satisfiability, Boolean Modeling and Computation, 11(1):99-131 2019.
[5] Maximilian Balandat, Brian Karrer, Daniel Jiang, Samuel Daulton, Benjamin Letham, An- drew Gordon Wilson, and Eytan Bakshy. BoTorch: A Framework for Efficient Monte-Carlo Bayesian Optimization. In Advances in Neural Information Processing Systems (NeurIPS), volume 33, 2020.
[6] Ricardo Baptista and Matthias Poloczek. Bayesian optimization of combinatorial structures.
In International Conference on Machine Learning, pages 462-471. PMLR, 2018.
[7] Petr Baudiš and Petr Pošík. Online Black-Box Algorithm Portfolios for Continuous Opti- mization. In Parallel Problem Solving from Nature PPSN XIII, pages 40-49, Cham, 2014. Springer International Publishing.
[8] James Bergstra and Yoshua Bengio. Random Search for nal of Machine Learning Research, 13(2), 2012.
Optimization. Jour-
[9] James Bergstra, Rémi Bardenet, Yoshua Bengio, and Balázs Kégl. Algorithms for Hyper- Parameter Optimization. In Advances in Neural Information Processing Systems (NeurIPS), volume 24, 2011.
[10] Mickael Binois and Nathan Wycoff. A survey on high-dimensional Gaussian process modeling with application to Bayesian optimization. ACM Transactions Evolutionary Learning and Optimization, 2(2):1-26, 2022.
[11] Mickaël Binois, David Ginsbourger, and Olivier Roustant. warped kernel improving ro- bustness in Bayesian optimization via random embeddings. In International Conference Learning and Intelligent Optimization (LION), pages 281-286. Springer, 2015.
[12] Mickaël Binois, David Ginsbourger, and Olivier Roustant. On the choice of the low- dimensional domain for global optimization via random embeddings. Journal of global op- timization, 76(1) 76(1))69-99)) 90, 2020.
[13] Jim Boelrijk, Bernd Ensing, Patrick Forré, and Bob WJ Pirok. Closed-loop automatic gradient design for liquid chromatography using Bayesian optimization. Analytica Chimica Acta, 1242, 2023.
[14] Mohamed Amine Bouhlel, Nathalic Bartoli, Rommel G. Regis, Abdelkader Otsmane, and Joseph Morlier. Efficient global optimization for high-dimensional constrained problems by using the Kriging models combined with the partial least squares method. Engineering Opti- mization, 50( (001))2038-222 20532 2018.
[15] Benjamin Burger, Phillip M Maffettone, Vladimir V Gusev, Catherine M Aitchison, Yang Bai, Xiaoyan Wang, Xiaobo Li. Ben M Alston, Buyi Li, Rob Clowes, A mobile robotic chemist. Nature, 583(7815):237-24,, 2020.
[16] Roberto Calandra, Nakul Gopalan, André Seyfarth, Jan Peters, and Mare Peter Deisenroth. Bayesian Gait Optimization for Bipedal Locomotion. Learning and Intelligent Optimization, pages 274-290. Springer International Publishing, 2014.

---

[17] Roberto Calandra, André Seyfarth, Jan Peters, and Mare Peter Deisenroth. Bayesian optimiza- tion for learning gaits under uncertainty. Annals of Mathematics and Artificial Intelligence, 76 (1):5-23, 2016.
[18] Zachary Cosenza, Raul Astudillo, Peter Frazier, Keith Baar, and David E Block. Multi- Information Source Bayesian Optimization of Culture Media for Cellular Agriculture. Biotech- nology and Bioengineering. 2022.
[19] Samuel Daulton, Xingchen Wan, David Eriksson, Maximilian Balandat, Michael A. Osborne, and Eytan Bakshy. Bayesian Optimization over Discrete and Mixed Spaces via Probabilistic Reparameterization. In Advances in Neural Information Processing Systems (NeurIPS), vol- ume 35, 2022.
[20] Aryan Deshwal, Syrine Belakaria, Janardhan Rao Doppa, and Dae Hyun Kim. Bayesian opti- mization over permutation spaces. In Proceedings of the AAAI Conference on Artificial Intelli- gence, volume 36, pages 6515-6523. 2022.
[21] Aryan Deshwal, Sebastian Ament, Maximilian Balandat, Eytan Bakshy, Janardhan Rao Doppa, and David Eriksson. Bayesian Optimization over High-Dimensional Combinatorial Spaces via Dictionary-based Embeddings. In International Conference on Artificial Intelligence and Statistics, pages 7021-7039, PMLR, 2023.
[22] Adel Ejjeh, Leon Medvinsky, Aaron Councilman, Hemang Nehra, Suraj Sharma, Vikram Adve, Luigi Nardi, Eriko Nurvitadhi, and Rob A. Rutenbar. HPVM2FPGA: Enabling True Hardware-Agnostic FPGA Programming. Proceedings of the 33rd IEEE International Con- ference on Application-speccifc Systems, Architectures, and Processors, 2022.
[23] David Eriksson and Martin Jankowiak. High-dimensional Bayesian optimization with sparse axis-aligned subspaces. In Cassio de Campos and Marloes H. Maathuis, editors, Proceedings of the Thirty-Seventh Conference on Uncertainty in Artificial Intelligence, volume 161 of Pro- ceedings of Machine Learning Research, pages 493-503. PMLR, 27-30 Jul 2021.
[24] David Eriksson and Matthias Poloczek. Scalable Constrained Bayesian Optimization. In Pro- ceedings of The 24th International Conference Artificial Intelligence and Statistics, volume 130 of Proceedings of Machine Learning Research, pages 730-738. PMLR, 13-15 Apr 2021.
[25] David Eriksson, Michael Pearce, Jacob Gardner, Ryan D Turner, and Matthias Poloczek. Scal- able Global Optimization via Local Bayesian Optimization. In Advances in Neural Information Processing Systems (NeurIPS), volume 32, pages 1499-5507, 2019.
[26] CA Floudas, HK Fung, SR McAllister, M Mönnigmann, and R Rajgaria. Advances in protein structure prediction and de novo protein design: A review. Chemical Engineering Science, 61 (3):96-9988, 2006.
[27] Peter Frazier. A tutorial Bayesian optimization. arXiv preprint arXiv:1807.028114 2018.
[28] Peter I. Frazier and Jialei Wang. Bayesian Optimization for Materials Design, pages 45-75.
Springer International Publishing, Cham, 2016. ISBN 978-3-319-23871-5.
[29] Jacob R Gardner, Geoff Pleiss, David Bindel, Kilian Q Weinberger, and Andrew Gordon Wil- son. GPyTorch: Blackbox Matrix-Matrix Gaussian Process Inference with GPU Acceleration. In Advances in Neural Information Processing Systems (NeurIPS), volume 31, 2018
[30] Seyede Fatemeh Ghoreishi and Douglas Allaire. Multi-information source constrained Bayesian optimization. Structural and Multidisciplinary Optimization, 59:977-991, 2019.
[31] Seyede Fatemeh Ghoreishi and Douglas I. Allaire. A fusion-based multi-information source op- timization approach using knowledge gradient policies. In 2018 AIAA/ASCF/AFSSSSC Struc- tures. Structural Dynamics, and Materials Conference, page 1159, 2018.
[32] Pierre Hansen and Brigitte Jaumard. Algorithms for the maximum satisfiability problem. Com- puting, 44(4):279-33) -03,, 1990.
12

---

[33] Florian Hase, Loic M Roch, Christoph Kreisbeck, and Alán Sspuru-Guzik Phoenics: Bayesian optimizer for chemistry. ACS central science, (19):1134-111 2018.
[34] Florian Häse. Matteo Aldeghi, Riley Hickman, Loic M Roch, and Alán Aspuru- Guzik. Gryf- fin: An algorithm for Bayesian optimization of categorical variables informed by expert knowl- edge. Applied Physics Reviews, 8(3):011406, 2021.
[35] Erik Hellsten, Carl Hvarfner, Leonard Papenmeier, and Luigi Nardi. High-dimensions Bayesian Optimization with Group Testing, arXiv preprint arXiv:2310.03515, 2023.
[36] Erik Hellsten, Artur Souza, Johannes Lenfers. Rubens Lacouture, Olivia Hsu, Adel Ejjeh. Fredrik Kjolstad, Michel Steuwer, Kunle Olukotun, and Luigi Nardi. BaCO: A Fast and Portable Bayesian Compiler Optimization Framework. ACM International Conference on Architectural Support for Programming Languages and Operating Systems, 2023.
[37] Henry C Herbol, Weici Hu, Peter Frazier, Paulette Clancy, and Matthias Poloczek. Efficient perovskites via Bayesian optimiza-
search of compositional space for hybrid tion. npj Computational Materials, 4(1):1-7, 2018.
(38) José Miguel Hernández-LLohato. James Requeima, Edward O. Pyzer-Knaphe and Alán Aspuru- Guzik. Parallel and Distributed Thompson Sampling for Large-scale Accelerated Exploration of Chemical Space. In Proceedings of the 34th International Conference on Machine Learning. volume 70, pages 1470-1479, PMLR, 06-11 Aug 2017.
[39] Yingjie Hu, Jian- MaanQQiag Hu, Yifan Xu, Fengchun Wang, and Rong Zeng Cao. Contamination control food supply chain. In Proceedings of the 2010 Winter Simulation Conference, pages 2678-2681 IEEE, 2010.
[40] Zak E Hughes, Michelle A Nguyen, Jialei Wang. Yang Liu, Mark T Swihart, Matthias Poloczek, Peter Frazier, Mare R Knecht, and Tiffany R Walsh. Tuning Materials-Binding Peptide Sequences toward Gold-and Silver-Binding Selectivity with Bayesian Optimization. 2021.
ACS nano, 15(11):18200-188
[41] Frank Hutter, Holger H. Hoos, and Kevin Leyton- Sequential Model-Based Optimiza- tion for General Algorithm Configuration. In Carlos A. Coello Coello, editor, Learning and Intelligent Optimization, pages 507-523, Berlin, Heidelberg, 2011. Springer Berlin Heidelberg. ISBN 978-3-642-25566-3
[42] Anant Singh Jain and Sheik Meeran. Deterministic job-shop scheduling: Past, present and future. European journal of operational research, 1 133())390-4434 1999.
[43] Donald R Jones, Matthias Schonlau, and William Welch. Efficient global optimization of expensive black-box functions. Journal of Global optimization, 13(4):455, 1998.
[44] Jungtaek Kim, Seungjin Choi, and Minsu Cho. Combinatorial Bayesian optimization with random mapping functions to convex polytopes. In Uncertainty in Artificial Intelligence, pages 1001-1011. PMLR, 2022.
[45] Rémi Lam, Matthias Poloczek, Peter Frazier, and Karen E Willcox. Advances in Bayesian optimization with applications in aerospace engineering. In 2018 AIAA Non-Deterministic Approaches Conference, page 1656, 2018.
[46] Ben Letham, Roberto Calandra, Akshara Rai, and Eytan Bakshy. Re-Examining Linear Em- beddings for High-Dimensional Bayesian Optimization. In Advances in Neural Information Processing Systems (NeurIPS), volume 33, pages 1546-1558, 2020.
[47] Chong Liu and Yu-Xiang Wang. Global optimization with parametric function approximation. In Proceedings of the 40th International Conference on Machine Learning, volume 202, pages 22113-22136, 2023,
[48] Daniel Lizotte, Tao Wang, Michael H Bowling, Dale Schuurmans, al. Automatic Gait Optimization With Gaussian Process Regression. In IJCAI, volume 7, pages 944-949, 2007.
13

---

[49] Trent W Lukaczyk, Paul Constantine, Francisco Palacios, and Juan J Alonso. Active subspaces for shape optimization. In 10th AlAA multidisciplinary design optimization conference, page
11 1171, 2014.
[50] Matthias Mayr, Faseeh Ahmad, Konstantinos I. Chatzilygeroudis, Luigi Nardi, and Volker Krüger. Skill-based Multi-objective Reinforcement Learning of Industrial Robot Tasks with Planning and Knowledge Integration. CoRR, abs/2203. 10033, 2022.
[51] Luigi Nardi, David Koeplinger, and Kunle Olukotun. Practical design space exploration. In 2019 IEEE 27th International Symposium on Modeling, Analysis, and Simulation of Computer and Telecommunication Systems (MASCOTS), pages 347-358, 2019.
[52] Amin Nayebi, Alexander Munteanu, and Matthias Poloczek. A framework for Bayesian Op- timization in Embedded Subspaces. In Proceedings of the 36th International Conference Machine Learning, volume 97 of Proceedings of Machine Learning Research (PMLR), pages 4752-4761, 09-15 Jun 2019.
[53] Diana M. Negoescu, Peter Frazier, and Warren B. Powell. The Knowledge-Gradie Algo- rithm for Sequencing Experiments in Drug Discovery. INFORMS Journal Computing, 23 (3):346366) 363, 2011.
[54] Kai Wang Ng, Guo-Liang Tian, and Man-Lai Tang. Dirichlet and related distributions: Theory, methods and applications. 2011.
[55] Changhai Nie and Hareton Leung. A survey of combinatorial testing. ACM Computing Surveys (CSUR), 43(2):1-22 2011.
[56] Changyong Oh, Jakub Tomezak, Efstratios Gavves, and Max Welling. Combinatorial Bayesian Optimization using the Graph Cartesian Product. Advances in Neural Information Processing Systems (NeurIPS), 32, 2019.
[57] Steve O'Hagan, Warwick B Dunn, Marie Brown, Joshua Knowles, and Douglas Kell. Closed-loop, multiobjective optimization of analytical instrumentation: gas chromatoraaplytiim-o-ffiiiee mass spectrometry of the metabolomes of human serum and of yeast fermentations. Analytical Chemistry, 77()):90-330, 2005.
[58] Tom Packebusch and Stephan Mertens. Low autocorrelation binary sequences. Journal of Physics A: Mathematical and Theoretical, 49(16):165001, 2016.
[59] Daniel Packwood. Bayesian Optimization for Materials Science. Springer, 2017.
[60] Leonard Papenmeier, Luigi Nardi, and Matthias Poloczek. Increasing the Scope You Learn: Adaptive Bayesian Optimization in Nested Subspaces. In Advances in Neural Information Processing Systems (NeurIPS), volume 35, 2022.
[61] Matthias Poloczek and David P Williamson. An experimental evaluation of fast approximation algorithms for the maximum satisfiability problem. Journal of Experimental Algorithmics (JEA), 22:1-18, 2017.
[62] Matthias Poloczek, Jialei Wang, and Peter Frazier. Multi-information source optimization. Advances in neural information processing systems, 30, 2017.
[63] Warren B Powell. A unified framework for stochastic optimization. European Journal of Operational Research, 275(3):755-821, 2019,,
[64] Akshara Rai, Rika Antonova, Seungmoon Song, William Martin, Hartmut Geyer, and Christo- pher Atkeson. Bayesian optimization using domain knowledge the ATRIAS biped. In 2018 IEFE International Conference on Robotics and Automation (ICRA), pages 1771-1778, 2018.
[65] Rommel G Regis and Christine A Shoemaker. Combining radial basis function surrogates and dynamic coordinate search high- high-dimensional expensive black-box optimization. Engineer- ing Optimization, 45(5):529-555, 2013.
14

---

[66] Danilo Jimenez Rezende, Shakir Mohamed, and Daan Wierstra. Stochastic backpropagation and upproximate inference deep generative models. In International conference on machine learning, pages 1278-1 186. PMIR, 2014.
[67] Binxin Ru, Ahsan Alvi, Vu Nguyen, Michael A Osborne, and Stephen Roberts. Bayesian optimisation over multiple continuous and categorical inputs. In International Conference on Machine Learning, pages 8276-8285. PMLR, 2020.
[68] Artur M Schweidtmann, Adam D Clayton, Nicholas Holmes, Eric Bradford, Richard A Bourne, and Alexei A Lapkin. Machine learning meets continuous flow chemistry: Automated opti- mization towards the Pareto front of multiple objectives. Chemical Engineering Journal, 352: 277-282, 2018.
[69] Kenan Sehic, Alexandre Gramfort, Joseph Salmon, and Luigi Nardi. LassoBench: A High- Dimensional Hyperparameter Optimization Benchmark Suite for Lasso. In First Conference on Automaled Machine Learning (Main Track), 2022.
[70] Bobak Shahriari, Kevin Swersky, Ziyu Wang, Ryan Adams, and Nando De Freitas. Taking the human out of the loop: A review of Bayesian optimization. Proceedings of the IEEE, 104 (1):148-175. 2015.
[71] Benjamin J Shields, Jason Stevens, Jun Li, Marvin Parasram, Farhan Damani, Jesus Martinez Alvarado, Jacob M Janey, Ryan Adams, and Abigail G Doyle. Bayesian reaction optimiza- tion as tool for chemical synthesis. Nature, 590(7844):89
2021.
[72] Alex J Smola and Bernhard Schölkopf. A tutorial on support vector regression. Statistics and computing, 14:199-222, 2004.
[73] Francisco J. Solis and Roger J-B. Wets. Minimization by Random Search Techniques. Mathe- matics of Operations Research, 6(1):19-30, 1981.
[74] Lei Song, Ke Xue, Xiaobin Huang. and Chao Qian. Monte Carlo Tree Search based Vari- able Selection for High Dimensional Bayesian Optimization. Advances in Neural Information Processing Systems (NeurIPS), 35, 2022.
[75] Artur Souza, Leonardo B Oliveira, Sabine Hollatz, Matt Feldman, Kunle Olukotun, James M Holton, Aina E Cohen, and Luigi Nardi. DeepFreak: Learning crystallography diffraction patterns with automated machine learning. arXiv preprint arXiv: /904. 1834, 2019.
[76] Hampus Gummesson Svensson, Esben Jannik Bjerrum, Christian Tyrchan, Ola Engkvist, and Morteza Haghir Chehreghani. Autonomous drug design with multi-armed bandits. 2022 IEEE International Conference on Big Data (Big Data), pages
IEEE, 2022.
[77] Alexander Thebelt, Calvin Tsay, Robert Lee, Nathan Sudermann-Mert David Walz, Behrang Shafei, and Ruth Misener. Tree ensemble kernels for Bayesian optimization with known con- straints over fixtu-feature spaces. Advances in Neural Information Processing Systems, 35: 37401-37415, 2022.
[78] Tsuyoshi Ueno, Trevor David Rhone, Zhufeng Hou, Teruyasu Mizoguchi, and Koji Tsuda. COMBO: An efficient Bayesian optimization library for materials science. Materials Discov- ery, 21, 2016.
[79] Xingchen Wan, Vu Nguyen, Huong Ha, Binxin Ru, Cong Lu, and Michael A. Osborne. Think Global and Act Local: Bayesian Optimisation over High-Dimensional Categorical and Mixed Search Spaces. In Proceedings of the 38th International Conference on Machine Learning. volume 139 of Proceedings of Machine Learning Research, pages 10663-10674, PMLR, 18- 24 Jul 2021.
[80] Jialei Wang, Scott C. Clark, Eric Liu, and Peter Frazier. Parallel Bayesian global optimization of expensive functions. Operations Research, 68(6): (666)))850-8666 2020.
[81] Ke Wang and Alexander Dowling. Bayesian optimization for chemical products and func- tional materials. Current Opinion in Chemical Engineering, 36: 100728, 2022.
15

---

[82] Linnan Wang, Rodrigo Fonseca, and Yuandong Tian. Learning Search Space Partition for Black-box Optimization using Monte Carlo Tree Search. Advances Neural Information Processing Systems (NeurIPS), 33: 13119111-1 2020.
[83] Xilu Wang, Yaochu Jin, Sebastian Schmitt, and Markus Olhofer. Recent Advances in Bayesian Optimization. ACM Comput. Surv., 2023.
[84] Ziyu Wang, Frank Hutter, Masrour Zoghi, David Matheson, and Nando de Feitas. Bayesian Optimization in Billion Dimensions via Random Embeddings. Journal of Artificial Intelli- gence Research (JAIR), 55:361-387, 2016.
[85] Christopher K Williams and Carl Edward Rasmussen. Gaussian processes for machine learn- ing, volume MIT press Cambridge, MA, 2006.
[86] James T Wilson, Riccardo Moriconi, Frank Hutter, and Marc Peter Deisenroth. The reparame- terization trick for acquisition functions. NeurlPS Workshop Bayesian Optimization, 2017.
[87] David P Woodruff. Sketching tool for numerical linear algebra. Foundations and Trends@ in Theoretical Computer Science, 10(1-2)11-17, 2014.
[88] Juliusz Krzysztof Ziomek and Haitham Bou Ammar. Are Random Decompositions all we need in High Dimensional Bayesian Optimisation? In International Conference Machine Learning, pages 43347-43368. PMLR, 2023.
16

---

A Consistency of Bounce
In this section, we prove the consistency of the Bounce algorithm. The proof is based on Papenmeier al. [60] and Eriksson and Poloczek [24].
Theorem (Bounce consistency). With the following definitions
Def. 1. (an))== is sequence of points of decreasing function values;
Def. 2. x² arg mingex is a minimizer of in X;
and under the following ussumptions:
Ass. 1. D is finit;;
Ass. 2. is observed without noise;
Ass. 3. The range of f is bounded in X, i.e., 30 R.4 [f()) < C Vo E Xx;
Ass. For at least one of the minimizers æ; the (partial) assignment corresponding the continu- ous variables lies in a (continuous) region with positive measure;
Ass. 5. One Bounce reached the imput dimensionality D. the continuous elements of the initial points { {Icont, 1 after each TR restart are chosen
(a) uniformly at random for continuous variables; and
(b) such that every realization of the combinatorial variables has positive probability;
then the Bounce algorithm finds a global optimum with probability as the number of samples N goes
Proof. The range of f is bounded per Assumption 3, and Bounce only considers function evalua- tion a 'success' the improvement over the current best solution exceeds a certain constant thresh- old. Bounce can only have finite number of 'successful" evaluations because the range of f is bounded per Assumption For the sake of a contradiction, we suppose that Bounce does not obtain an optimal solution its number of function evaluations N oc. Thus, there must be a sequence of failures, such that the TRs in the current target space, i.c., the current subspace, will eventually reach its minimum base length. Recall that in such an event, Bounce increases the target dimension by splitting up the 'bins', thus creating a subspace of (b + 1)-times higher dimensionality. Then TR that again experiences sequence of failures that lead to another split, and so on. This series of events repeats until the embedded subspace eventually equals the input space and thus has dimensionality D. See lines 12 16 in Algorithm in Sect. 3.
Bounce creates
Still supposing that Bounce does not find an optimum in the input space, there must be a sequence of failures such that the side length of the TR again falls below the set minimum base length, now forcing restart of Bounce. Recall that at every restart, Bounce samples fresh set of initial points uniformly at random from the input space; see line in Algorithm Therefore, with probability 1, random sample will eventually be drawn from any subset ) C X with positive Lebesgue measure (v())= > 0):
lim (1 - (())* 11 k=000
where y is the uniform probability measure of the sampling distribution that Bounce employs for initial data points upon restart [73].
Let
a inf (t: E f(x) < <<>0> 0}
denote the essential infimum of ff on X' with being the Lebesgue measure [73].
Following Solis and Wets [73], we define the optimality region, i.c., the set of points whose function value is larger by most than the essential infimum:
R (xr X f(x) <a+e}
with e > and M < Because Ass. 4, at least one optimal point lies in a region of positive measure that is continuous for the continuous variables. Therefore, we have that = f(x"). Note that this is also the case if the domain of f only consists of combinatorial variables (Ass. 5). Then, R..M fx X f()) < f(x") +e).
17

---

53D-AckleyS3
53D-Ackley53 (publ.)
Bounce
10a
BODI
Objective
Random Search
Casmopolitan
100
100
Iteration
Iteration
Figure 8: Bounce the other algorithms on the synthetic Ackley53 benchmark function. Bounce outperforms all other algorithms and quickly finds excellent solutions, BODi performance degrades upon randomization.
Let (#t))== denote the sequence of best points that Bounce discovers with ot being the best point up to iteration k. This sequence satisfies Def. 1 by construction. Note that o E Re, implies that oi E R. for all K + [73] because observations are noise-free. Then,
P P E R
1 - P c E X Re
21 -1-11 (1-1(Re,.M)) (()e,M)
and,
2 lim PP (art
E
Re,M] 2
lim (1 - je( ()) =1, k-+00
1000
(1)
i... 2 eventually falls into the optimality region [73]. By letting + 0, æ converges to the global optimum with probability 1 as ke
D0.
B Additional experiments
We compare Bounce to the other algorithms on three additional benchmark problems: Ackley53 and MaxSAT60 [21]. Moreover,
two additional studies investigate the performance of Bounce further. First, we run Bounce on a set of continuous problems from Papenmeier et al. [60] showcase the performance and scalability of Bounce on purely continuous problems. We then present "low-sequency" version of Bounce to showcase how such a version can outperform its competitors on the original benchmarks by introducing bias towards low-sequency solutions.
B.1 Bounce and other algorithms on additional benchmarks
B.1.1 The synthetic Ackley53 benchmark function
Ackley53 is 33-dimensional function with 50 binary and three continuous variables. Wan al. [79] discretized 50 continuous variables of the original Ackley function, requiring these variables to be either zero or one. This benchmark was designed such that optimal value of 0.0 is the ori- gin (0..... 0). Here, we perturb the optimal assignment of combinatorial variables by flipping each binary variable with probability 1/2. Figure 8 summarizes the performances the algorithms. Bounce outperforms all other algorithms and proves to be robust to the location of the optimum point. Casmopolitan taan a distanced runner-up. BODi initially outperforms Casmopolitan on the published benchmark version but falls behind later.
B.1.2 Contamination control
The Contaminati benchmark models supply chain with 25 stages [39]. At each stage, binary decision is made whether to quarantine food that has not yet been contaminated. Each such inter- vention is costly, and the goal is minimize the number of contaminated products and prevention cost [6, 56). Figure 9 shows the performances the algorithms.
Bounce, Casmopol tan, and BODi all produce solutions of comparable objective value. Bounce and Casmopoli i tan find better solutions than BOD initially, after about 100 function evaluations, the solutions obtained by the three algorithms are typically par.
18

---

ZSD-Contaminatton
ZSD-Contariination (publ.)
22.50
Bounce
22.25
BODI
22.00
Random Search
aammopalitan
21.75
RDUCE
SMAC
21.50
100
200
0
100
200
tenaaion
Iteration
Figure 9: Bounce and the other algorithms on the 25-dimensional contamination problem. Bounce performs par with Casmopolite itan and BODi on both versions of the benchmark.
B0D-MaxSAT60
BOD-MaxSATéO (publ.)
-50
Baunce
BODI
-100
Objective
Cosmopolitan
150
COMBO
RDUCB
SMAC
-200
100
200
0
100
200
Iteration
Iteration
Figure 10: Bounce and other algorithms on the 60-dimensional weighted maximum satisfiability problem. Bounce is the first to find an optimal solution (left). On the published version (right), Bounce comes in second after BOD1.
B.1.3 The MaxSAT60 benchmark
MaxSAT60 is a 60-dimensional, weighted instance of the Maximum Satisfiability (MaxSAT) prob- lem. MaxSAT AT is notoriously hard combinatorial problem that cannot be solved in polynomial time (unless PP NP). The goal find binary assignment to the variables that satisfies clauses of max- imum total weight. For every in {1, 2....... ... d}. this benchmark has one clause of the form 2; with weight of and 638 clauses the form -, V -xy with a weight of 61. Following [21, 56, 79], we normalize these weights to have zero mean and unit standard deviation. This normalization causes the one-variable clauses to have a negative weight, i.e., the function value improves such a clause not satisfied, which is atypical behavior for MaxSAT problem. Since the clauses with two variables are satisfied for 2; =
and the clauses with one variable negative weights are never satis- fied for 2; 0, the normalized benchmark version has global optimum at x (0, (O..... 0) by con- struction. The problem's 's difficulty is finding an assignment for variables such that all two-variable clauses satisfied and many one-variable clauses possible captured by normalized weights.
Figure 10 summarizes the performances the algorithms. The general version that attains the global optimum for randomly selected binary assignment is shown on the left. The special case where the global optimum to the all-zero assignment is shown the right.
We observe that Bounce requires the smallest number of samples to find an optimal assignment in general, followed by BODi and Casmopoli tann. Only in the special case where the optimum the all-zero assignment, BODi ranks first, confirming the corresponding result in Deshwal al. [21].
19

---

C An evaluation of Bounce on continuous problems
5000-Branin2
124D-Mopta08
500D-Hartmann6
Dbjective value Regret
300
20"
260
10a
20*
240
Bounce. b=1
188D-SVM
1B0D-LassoDNA
1000D Lassolaad
Bounce, b-3
Objedtre value valut Regret
b-10
0.20
0.305
Bounce. b=
0.300
0.15
0.255
0.10
0.290
0.05
200
400
200
400
400
batch evaluations
batch eva evaliattons
batch evaluations
Figure 11: Bounce continuous problems with different batch sizes (plotted in terms of batch evaluations).
500D-Branin2
124D-Mopta08
500DD-arrmamm6
Objective Regret
280
20²
260
220
1000
2000
1000
2000
1000
2000
Bounce. b-1
388D-SVM
1BB0D-LaSSODNA
1000D-LassoHard
Obledive vakue Regret
0.31
10
Bounce, b-
0.20
0.30
0.10
0.05
1000
2000
1000
2000
1000
2000
evaluations
function
function
Figure 12: Bounce on continuous problems with different batch sizes (plotted in terms of function evaluations).
To showcase the performance and scalability of Bounce, we run it on a set of continuous problems from Papenmeier al. [60]. The 124-dimensional Mopta08 benchmark is constrained vehicle optimization problem. We adopt the soft-constrained version from Eriksson and Jankowiak [23]. The 988-diimnsiionl SVM problem [23] concerns the classification performance with an SVR on the slice localization dataset. The 180-dimensional LassoDNA benchmark [69] is a sparse regression problem real-world dataset, and the 1000-dimensional LassoHard benchmark optimizes over a synthetic dataset. The 500-dimensional Branin2 and Hartmann6 problems are versions of the and 6-dimensional benchmark problems where additional dimensions with effect the function value were added.
We set the number function evaluations max(2000, 500B) for batch size of B and configure Bounce such that it reaches the input dimensionality after 500 function evaluations. Figures 11 and 12 show the simple regret for the synthetic Branin2 and Hartmann6 problems, and the best function value obtained after given number of batch (Figure 11) or function (Figure 12) evaluations for the remaining problems: Mopta08, SVM, LassoDNA, and LassoHard.
We observe that Bounce always benefits from more parallel function evaluations. The difference between smaller batch sizes, such as B and B or B and B 10, more remarkable than between larger batch sizes, like B 10 and B 20. Parallel function evaluations prove
20

---

especially effective on SVM and LassoDNA. Here, the optimization performance improves drastically. We conclude that small number of parallel function evaluations already helps considerably increase the optimization performance.
On the synthetic Branin2 and Hartmann6 problems, Bounce quickly converges to the global opti- mum. Here, we see that larger number of parallel function evaluations also helps in converging to a better solution.
In Figure 11, we also compare to BAxUS by Papenmeier et al. [60], which does not support parallel function evaluations and is therefore run with batch size of 1. We observe that Bounce with a batch size of 1 outperforms BAxUS on all benchmarks except for Hartmann6, showcasing Bounce's overall performance improvements even for continuous problems and single-element batches.
D Batched evaluations on mixed and combinatorial problems
In addition to Figure which shows the best objective value in relation to the number of batches, Figure 13 shows the objective value in terms of function evaluations.
25D-PestControl
60D-Max5AT60
25D- 3D-Contaminaatio
18
23.0
Obyccive
26
-100
22.5
344
-150
22.0
12
200
Bounce, b-1
1000
2000
Bounce. b=3
S3D-SVM
S3D-AckleyS:
50D-LABS
Bounce, b-5
Bounce, b-10
objecivve value
Bounce, b-20
4 10
3 30-1
2x
0
1000
2000
2000
2000
2000
2000
function evoluations
function evalustions
function Ieluaatoons
Figure 13: Bounce on mixed and combinatorial problems with different batch sizes (plotted in terms of function evaluations).
The figure confirms that Bounce leverages parallel function values efficiently.
E Low-sequency version of Bounce
We show how we can bias Bounce towards low-sequency solutions. A binary vector has low se- quency if there are few changes from 0 to or vice versa. Similarly, a solution to a categorical ordinal problem has low sequency if there are few changes from category level to other. We remove the random signs (for binary and continuous variables) and the random offsets (for categorical and ordinal variables) from the Bounce embedding. We conduct this study to show a) that Bounce can outperform BODi on the unmodified versions of the benchmark problems if we introduce similar bias towards low-sequency solutions and b) that the random signs empirically show to remove biases towards lew-sequency solutions. However, we want to emphasize that the results of this section are not representative of the performance of Bounce on arbitrary real-world problems. Nevertheless, if one knows that the problem has low-sequency structure, then Bounce be configured exploit this structure and outperform BODi.
Figure 14 shows the results of the low-sequency version of Bounce on the original benchmarks from Section 4. We observe that Bounce outperforms BODi and the other algorithms on the unmodified versions of the benchmark problems. This shows that Bounce can outperform B0Di on the unmodi- fied version of the benchmarks if we introduce a similar bias towards lewseequnccy solutions.
Figure 15 shows the results of the low-sequency version of Bounce on the flipped benchmarks from
Section 4. The sequency version of Bounce robust towards the randomization of the optimal point.
21

---

53D-5VM
53D-Ackley55
60D-MaxSAT50
CR value
40
10a
-100
2x10-
10
100
Bounce
8DDi
25D-PestControl
25D-Contamination
50D-LABS
Casmopolitan
23.5
COMBO
18
23.0
objedtve
26-
22.5
22.0
12 100
21.5
0
200
100
200
100
200
Iteration
Iteratioon
Iteration
Figure 14: 'Low-sequency" version Bounce on the original benchmarks from Section 4: with bias towards low-sequency solutions, Bounce outperforms BODi the original versions of the benchmark problems.
53D-SVM
53D Ackley53
600 MaxSAT60
-
Oboccive volue
109
4 4x
-100
30-
-150
2x
10
-200
100
200
100
200
100
200
Bounce
25D-PastControl
25D-Contamination
50D-LABS
Casmopolitan
23.5
CONBO
3.8
26
Objective
22.5
24
22.0
12
200
100
200
100
Figure 15: "Low-sequency" version of Bounce on the modified benchmarks from Section 4.
F Effect of trust-region management
Bounce uses novel trust-region management that differs from previous approaches in that it al- lows arbitrary trust region base lengths in [Lmin Lmma]. This strategy allows Bounce to efficiently leverage parallel function evaluations, which we refer to as batch acquisition.
We compare Bounce's trust-region management strategy with the strategy employed by BAXUS [60] for purely continuous benchmarks and an adapted version of Casmopolitan [79] strategy for mixed or discrete-space benchmarks. For the latter, the difference is that in Bounce we reduce the failure tolerance as described by Papenmeier al. [60]: We first calculate the number times the
TR base length needs be reduced to reach the minimum TR base length
log1.5-1
Linia
Limh
and then find the failure tolerance for the i-th target space by T/a max min ([m=
where
m is the budget for the i-th target space [60]. This change in regard to Casmopol taa is necessary because Casmopoli itan's failure tolerance of 40 [79] would cause the algorithm to spend a large part of the evaluation budget in initial target spaces of low dimensionality.
We show the effect of this strategy on discrete and mixed-space benchmarks in Figure 16. Figure 17 shows the effect on continuous benchmarks. All experiments are replicated times. observe that the proposed trust-region management method not only enables efficient batch parallelism but also improves the performance for single-function evaluations when compared to the respective baseline for the TR management that stated above. The TR management proposed here usually
22

---

Contaminaaionn
Dijectio calow heggret
-0
50
-6.
-100
22.5
14
-150
Bounce,
-200
Bounce, b-3
Bounce, bm5
53D-SVM
53D-AckeeS3
50D-LABS
Objedive value FRegre
20"
3xx10
b-
2.4
prav. TR method
22
2x
18x10
16x
100
200
100 evaiuations
Figure 16: Effect of the proposed trust region management the discrete and mixed benchmarks. Here, the baseline is the TR management of Casmopoli tan. We see that often, even large batch sizes for the old TR management strategy do not outperform the proposed new strategy with batch size of 1.
500D-Branin2
124D-Mopta08
500D-Hartmann6
objozive valuevevall////o Rogres
101
300
10
250
Bounce, bm1
Bounce, bm3
388D-SVM
180D LassaDNA
10s0D-Lssaaard
Odeccbe value RRgges
0.25
0.32
15.0
Bounce, bm
0.31
12.5
prev. TR method
10.0
0.15
0.30
0.10
0.29
0.0
400
evalustions
Figure 17: Effect of the proposed trust-region management the continuous benchmarks. Here, the baseline is the TR management of BAXUS. We observe that larger batches help to improve per- formance for both methods. The new strategy, however, outperforms the baseline in almost all cases with the exception of LassoDNA with a batch size of 3.
provides better solutions the same number of batches than the respective baseline TR management. There are a few exceptions. On Labs, the previous TR management strategy outperforms Bounce's 's strategy by small margin for batch sizes 5, 10, and 20, and on Ackley53, the previous strategy converges to a better solution for batch size of three, five, and ten.
On the continuous benchmarks (Figure 17), we observe that the new TR management strategy also improves the performance for single function evaluations with the exception of Branin2. Similar to the mixed and combinatorial-space benchmarks, large batch sizes (10 and 20) bring little advantage compared to batch size of 5, For Hartmann6, the previous strategy with batch size of 20 performs worse than 3.
G Implementation details
We implement Bounce in Python using the BoTorch [5] and GPyTorch [29] libraries.
We employ a T(1.5, (.1. prior on the lengthscales of both kernels and F(1.5, 0.5) prior on the signal variance. We further use a I'(1.1, 0.1) prior on the noise variance.
23

---

Motivated by Wan et al. [79] and Eriksson et al. [25], we use an initial trust region baselength of 40 for the combinatorial variables, and 0.8 for the continuous variables. We maintain two separate TR shrinkage and expansion parameters (Temb and *emt) for the combinatorial and continuous variables, respectively such that each TR base length reaches its respective minimum of 1 and 2- after given number of function evaluations. When Bounce finds a better or worse solution, we increase or decrease both TR base lengths.
We use the author's s implementations for COMBO², BODi 3, RDUCB*, SMAC", and Casmopolitan", We use the same settings as the authors for COMBO, RDUCB, SMAC, and BODi. For Casmopol tan, we use the same settings as the authors for benchmarks reported in Wan et al. [79] and set the initial trust region base length to 40 otherwise.
Due to its high-mun memmor footprint, we ran BODi on NVidia A100 80GB GPUs for 300 GPU/h. We ran Bounce on NVidia A40 GPUs for 2,000 GPU/h. We ran the remaining methods for 20,000 GPU/h on one core of Intel Xeon Gold 6130 CPUs with 60GB of memory.
G.1 Optimization of the acquisition function
We use different strategies to optimize the acquisition function depending on the type of variables present in problem.
Continuous problems. For purely continuous problems, we follow a similar approach as TuRBO [25]. In particular, we use the lengthscales of the GP posterior to shape the TR. We use gradient descent to optimize the acquisition function within the TR bounds with 10 random restarts analytical EI. For larger batch sizes,
and 512 raw samples. For batch size of 1,
the
BoTorch implementation of qEI [5, 66, 86].
Binary problems. Similar Wan et al. [79], we use discrete TRs centered the current best solution. A discrete TR describes all solutions with a certain Hamming distance to the current best solution.
We use a local search approach to optimize the acquisition function for all problems with combina- torial search space of only binary variables: When starting the optimization, we first create a set of min max (2000, 200 dz)) random solutions. The choice of the number of random solutions based on TuRBO [25]. For each candidate, we first draw L, indices uniformly at random from (1...., d without replacement, where L, is the TR length the i-th iteration. We then sample d; values in {0, }} and set the candidate at the sampled indices to the sampled values. All other values set to the values of the current best solution. Note that this construction ensures that each can- didate solution lies in the TR bounds of the current best solution. We add all neighbors (i.e., points with Hamming distance of 1) of the current best solution to the set of candidates. This is inspired by BODi [21]. We find the 20 candidates with the highest acquisition function value and use local search to optimize the acquisition function within the TR bounds: At each local search step, we create all direct neighbors that do not coincide with the current best solution or would violate the TR bounds. We then move the current best solution to the neighbor with the highest acquisition function value. We repeat this process until the acquisition function value does not increase anymore. Finally, we return the best solution found during local search.
Categorical problems. We adopt the approach for binary problems, i.c., we first create set of random solutions with the same size as for purely binary problems and start the local search on the 20 best initial candidates.
Suppose the number of categorical variables of the problem is smaller or equal to the current TR length. In that case, we sample, for each candidate and each categorical variable, index uniformly at random from 21.2.... [;]] where |ex] is the number of values of the i-th categorical variable. We then set the candidate at the sampled index and all other values to 0.
http:::/ggtaaa..cco//////////// unspecified license, last access: 2023-05-04
3 httpe:///gthbbboocmmaaaaaaaaaaaa license provided, last access: 2023-05-04
http:://gthabl.ccc//uueeeeoooa/r httpeeeeeeentttee//eeeee MIT license, last access: 2023-10-
20 htttpp//ggtaaal.cc///ingaaaaaaaa MIT license, last access: 2023-05-04
httpp://ggihaa.cc////ttaaa/aagea license, last access 2023-10-20
24

---

If the number of categorical variables of the problem is larger than the current TR length L, we first sample L, categorical variables uniformly at random from [d,] without replacement. For each initial candidate and each sampled categorical variable, we sample an index uniformly at random, for which we set the categorical variable to and all other values to O. The values for the variables that were not sampled are set to the values of the current best solution.
As for the binary case, we add all neighbors of the current best solution to the of candidates, and we sample the 20 candidates with the highest acquisition function value.
We then multi-start local search to optimize the acquisition function within the TR bounds while neighbors are created by changing the index of one categorical variable. Again, we repeat until convergence and return the best solution found during the local search.
Ordinal problems. The construction for ordinal problems is similar to the one for categorical problems.
Suppose the number of ordinal variables the problem is smaller or equal to the current TR length. In that case, sample an ordinal value uniformly at random to set the ordinal variable for each candidate and each ordinal variable. Otherwise, we choose as many ordinal variables as each candi- date's 's current TR length and sample an ordinal value uniformly random to set the ordinal variable. We add all neighbors of the current best solution, all solutions where the distance to the current best solution is for ordinal variable, to the set of candidates. We then sample the 20 candidates with the highest acquisition function value and use local search to optimize the acquisition function within the TR bounds. In the local search, we increment decrement the value of single ordinal variable.
Mixed problems. Mixed problems are effectively handled by treating every variable type sepa- rately. Again, we create a set of initial random solutions where the values for the different variable types are sampled according to the abovementioned approaches. This can lead solutions outside the TR bounds. We remove these solutions and find the 20 best candidates only across the solutions within the TR bounds.
When optimizing the acquisition function, we differentiate between continuous and combinatorial variables. We optimize the continuous variables by gradient descent with the same settings as purely continuous problems. When optimizing, we fix the values for the combinatorial values.
We use local search optimize the acquisition function for the combinatorial variables. In this step, we fix the values for the continuous variables and only optimize the combinatorial variables. We create the neighbors by creating neighbors within Hamming distance of 1 for each combinatorial variable type and then combining these neighbors. Again, we run a local search until convergence.
We do five interleaved steps, starting with the continuous variables and ending with the combinatorial variables.
G.2 Kernel choice
We the CoCaBo kernel [67] with one global lengthscale for the combinatorial and ARD for the continuous variables:
k(x,a²)= of (phicmb (DremmbbDeenn a(mb)) Himm (Fems alcm))
(1 p) (Kimb (Femb; amb)) + Kent (Featt al(o)))
V5remb
5r
V5ramb
Kembl (Femb; (mmb))
1
lenb
3 36
exp
Lemb
Temb
Lemb
x'mblll ll
Rem ( Fremt (mm)
Vorca + 5r
exp (5rrcm))
(I, x2)2
Fent
xE[d] d, continuous
ent.
where Dent and Femb are the continuous and combinatorial variables in æ, respectively, i.e.
Lemb (xx d, is combinatorial)
25

---

ent = (zz: d; is continuous)
and p [0, is tradeoff parameter learned jointly with the other hyperparameters during the likelihood maximization. Here, femb and kents are the lengthscale hyperparameters.
H Additional analysis of BODi and COMBO
H.1 Analysis of BODi
Binary problems. BODi by Deshwal et al. [21] uses a dictionary of anchor points A am) to encode a candidate point In particular, the i-th entry of the m-dimensional embedding @ (2(2) is obtained by computing the Hanming-distance between 2 and a, Notably, Deshwal et al. [21] choose the dimensionality of the embedding as 128, which is larger than the dimensionality of the objective functions themselves.
The sampling procedure for dictionary elements a, is chosen to yield wide range of sequencies.
The sequency of a binary string defined the number of times the string changes from 0 to 1 and vice versa. Deshwal et al. [21] propose two approaches to generate the dictionary elements: (i) using binary wavelets, and (ii) by first drawing a Bernoulli parameter ~ (4( 1) for each E [m] and then drawing binary string a, from the distribution B(0;). The latter approach is their default method.
We prove that BODi generates an all-zeros (or, by a similar symmetry argument, all-ones) representer point with a probability that is significantly higher than 2 2- Moreover, we claim without proof that a similarly increased probability holds for points with low Hamming distance to all-zeros or all-ones. This is consistent with the intention of BODi to sample points of diverse sequency and necessarily an issue.
However, claim without proof that having such points in the dictionary substantially increases the probability that BODi evaluates the all-zeros (and, by symmetry, the all-ones) points. That hypothesis consistent with observation in Section 4.6 that BODi has much higher chance of finding good or optimal solutions when they are near the all-zero point.
Deshwal et al. [21] choose the dictionary to have m 128 dictionary elements. Given Bernoulli parameter 0,, the probability that the i-th dictionary point a; a point of sequency zero given by 0," + (1 - P:
IP("zero sequency" |0,)
II op +1 (1
Then, since 0g follows uniform distribution, the overall probability for point of zero sequency is given by
m
Pe"zeero sequency") =
(1-0P - (1 0)) p(.)) d0, i=1
prob. of times not sequency
DD+1
(1 0,) D+1
0 D+1
D -
i=1
00 m1
=1 D +
i.e., the probability of at least one dictionary element being of sequency zero is (1 The probability of BODi's 's dictionary to contain zero-sequency point increases with the number of dictionary elements m and decreases with the function dimensionality D (see Figure 18).
D+I m
For instance, for the 60-dimensional MaxSAT60 benchmark, the probability that at least one dictio- 128
nary element is of sequency zero is 1
2
RJ 0.986 (see Figure 18).
Note that at least one point 2" has probability of S 1/2d to be drawn. The probability of the dictionary containing that x is less than equal to 1 (1 1/2d) m which is already less than
26

---

1000
1.000
problem dimensions
800
0.750
600
0.625
0.500
400
0.375
0.250
200
0.125
MaxSAT60
0.000
200
400
600
800
1000
# dictionary elements
Figure 18: Probabilities BODi to contain a vero-sequency solution for different choices of the dictionary size and the function dimensionality D).
0.01 for = 14 and m = 128. In Section 4, we show that randomizing the optimal point structure leads to performance degradation for BODi. We hypothesize this due to the reduced probability of the dictionary containing the optimal point after randomization.
Categorical problems. We calculate the probability that BODi contains vector in its dictionary where all elements are the same. For categorical problems, BODi first samples a vector O from the Tmax-eimplex A Tiaa for each vector a, in the dictionary, with Tmax being the maximum number of categories across all categorical variables of problem. We assume that all variables have the same number of categories as is the case for the benchmarks in Deshwal et al. [21]. of categories of the variables. For each element in a, BODi draws value from the categorical distribution with probabilities 0. While line 7 in Algorithm that the elements in 0 are shuffled for every element in a, we observe that 0 remains fixed based on the implementation provided by the authors?. The random resampling of elements from @ is probably only used for benchmarks where the number of realizations differs between categorical variables.
7 the number
Deshwal et al. [21] might suggest
Then, for fixed 0, the probability that all D elements in for any are equal some fixed value I E (1,2,...,,7) 7 is given by 0f. The probability that, for any the m dictionary elements, all D elements in are equal to some fixed value i. E 21,2,..., is given by
/
PP(alll one specific category") = II
0 (()))dd)
(2)
1=1
We note that 0 follows a Dirichlet distribution with = [54]. Then, O, is marginally Beta(1, T- distributed [54]. With that, Eq. (2) becomes
P"aall one specific category")) =1 II
Es Bct [1 0.]
1m1
Bot [0P]
now, by using the equality E[x"1 bution [54]
D-1 atr a+
for the D-th raw moment of a Beta(a, 3) distri-
D
m
1
++r
T +1
T ++ -1 -
"See https://githut: wttpaau...oo///aaaaaaaaaaa bodi/categorical iitirgoorrca_ddictinaaaryy kernel peernel.pyHL1..
27

---

T-2
T-3
T=4
T+5
T-0
prod probicon dimorsioos D
0.573
20
20
20
20
20
0.753
0.590
10
10
10
30
10
0.375
0 252
0.123
0.000
50 dictinary
102 dictionary omonss
50
dictionary elements
dictioaay
dictonary
Figure 19: Probabilities of BODi's dictionary to contain at least categorical point where each category has the same value. The probability increases with the number of dictionary elements m but decreases with the number of categories and the number of problem dimensions D.
D!
m
(T+D-
We discussed in Section 4.3 that the PestControl benchmark obtains a good solution at a 5. One could assume that BODi performs well on this benchmark because its dictionary will likely contain this point. However, we observe that the probability is effectively zero for 5, = 128, and 25 (see Figure 19), which are the choices for the PestControl benchmark in Deshwal et al. [21]. This raises the question of (i) whether our hypothesis is wrong and (ii) what the reason for BODi 's performance degradation on the PestControl benchmark is.
We show that BODi's reference implementation differs from the algorithmic description in an impor- tant detail, causing BODi to be considerably more likely sample category five on PestControl (or the "last" category for arbitrary benchmarks) than any other category.
12000
Category 1
Category 2
10000
Category 3
Category 4
Category 5
8000
Count
6000
4000
2000
10 11
12 13
14
16 17
18
19
20
21
22
24
number elements set to value
Figure 20: Histograms over the number of dictionary element entries set each category for 20,000 repetitions of the sampling of dictionary elements for the PestControl benchmark. For each of the five categories and each value the r-axis, the figure shows how often the number of entries in dictionary element equals the value on the z-axis for the given category. For example, the count for and category zero, indicating that all 20,000 dictionary points had least one entry '5', There is considerably higher chance for a dictionary element entry to be set to category five than any of the other categories.
Figure 20 shows five histograms over the number of dictionary elements each category. The values on the z-axis give the number of elements in a 25-dimensional categorical vector being set specific category. One would expect that the histograms have similar shape regardless of the category. However, for category 5, we see that more elements are set to this category than for the other categories: The probability of elements being set to category is almost twice as high as probability of being set to another category for k 2 3. In contrast, the probability that element in the vector belongs to category 5 is virtually zero. This behavior is beneficial for the PestControl benchmark, which obtained the best value found during our experiments for xx 55.5.... ...,
28

---

(see Section 4). While we see that the probability of each dictionary entry being set to category 5 is very low, we assume that we sample sufficiently many dictionary elements within small Hamming distance to the optimizer such that BODi's GP can use this information to find the optimizer.
The reason for oversampling of the last category lies in rounding issue in sampling dictionary elements. In particular, for a given dictionary element a; and corresponding vector O with |0] т, for each {1,2, ......" T remaining D E Do: elements are then set to category т. This causes the last category to be overrepresented in the dictionary elements. For the choices of the PestControl benchmark, D = 25 and r = 5, the first four categories had a probability of z 0.1805. In contrast, the last one had a probability of = 0.278 for 10 simulations". We assume that the higher probability of the last category is the reason for the performance difference between the modified and the unmodified version of the PestControl benchmark.
1}, Deshwal et al. [21] set [D0,] elements to category i. The
H.2 COMBO on categorical problems
On the categorical PestControl benchmark, COMBO [56] behaved similar to BODi [21]
COMBO
COMBO (flip)
140000
# element entries set to category
120000
100000
80000
60000
40000
20000
0
category
category
The histogram shows how many element entries of a 25-dimensional dictionary element are set to each of the five categories. There a considerably higher chance for a dictionary element entry to be set to category or 5 than one of the other categories.
Figure 21:
The histogram in Figure 21 shows how many element entries of 25-dimensional dictionary element are set to each of the five categories. We see that the first and the last categories are over-represented on both benchmark versions. As discussed in Section 4, this benchmark attains the best observed value for x = (5, 55.5.....5 5, 1). We observe that on the original and modified version of the bench- mark, where the categories are shuffled, COMBO sets disproportionately many elements categories one and five. Note that for the modified version of the benchmark, all categories are equally likely in the optimal solution.
We argue that this behavior is at least partially caused by implementation error in the construction of the adjacency matrix and the Laplacian for categorical problems This error causes categorical variables to be modeled like ordinal variables. According to Oh et al. [56], categorical variables are modeled complete graph (see Figure 22).
"The 95% confidence intervals for categories 1-5 are (0.1799, 0817), 18.1802, 0.1810), 18.1, 0.1811). (0.1801, 0.1809), (0.2775, 0.2783). Pairwise Wilcoxon signed-rant tests between categories 1-4 and category gives p values of (W z 4.7 1o each).
con/ QUV http://.....ooooooooooooooooooaa
hotpD//errpriiinneeeeeeeeeeennns fwwct cttp:eipppppeeatttrrccccc..cclca PyW1137, last access: 2023-04-25
29

---

al
a2
a0
43
a4
solloringg the
riceipttoon in the paper of Oh et al. [56], a carfggriali
Tigune 22.
variabee
five with
cureesortes is modella as compller a srappit.
of CONED,, we find the aellaceeay matth: for the first of carlggyy
mno Looking
code source the Mutzable piggorccii with five calteeutte is consttuucr as
0
0
0
1
0
0
0
0
0
0
0
0
0
the adijuceney matrix for a Sive-vettee patti greppi.
is which
30


Reviews:
[{"review_id": "Dp42J6RrJl", "reviewer": "Reviewer_ACyR", "summary": "This paper proposes an algorithm called Bounce, using nested embedding for mixed and combinatorial search space. Bounce partitions input variables into ‘bins’ and set all variables within the same bin to a single value to reduce the dimension. During the optimization, Bounce splits up bins into smaller bins for more refined optimization. ", "questions": "I am puzzled about why the performance of BODi and COMBO can be sensitive to the location of the optima, which is equivalent to shuffling the labels of the categories of each variable. BODi uses Hamming distance and COMBO uses the combinatorial graph to model the relation of different labels for each variable. Different labels for the same variables are the same, so shuffling the labels will not impact the performance. Can you provide more explanation?\n\nI think the variable selection methods [4-5] maybe also can be applied in the mixed spaces tasks. I want to see some discussions. \n\nI'm glad to increase my score if you address the proposed isssues.", "limitations": "Yes.", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1/ The combinatorial and mixed optimization is significant and has wide range of applications in the real world. \n\n2/ The algorithm is reasonable. Bounce designs a reasonable mapping for binary, categorical and ordinal variables, while a similar algorithm, BAXUS, which also uses nested random subspaces, has shown good performance on continuous optimization. \n", "weaknesses": "I don't have any major complaints. The proposed algorithm is a natural extension of BAXUS for categorical and mixed space. The experiments show that Bounce is reliable and can achieve good performance.   \n\nHere are some minor comments: \n\n1/ One minor weakness is the the high-dimensional continuous spaces part in Section 2. A recent survey paper [1] categorizes these works into several categories, i.e., low-dimensional embeddings, decomposition, and variable selection. This paper mainly discuss the work in low-dimensional embeddings I think discuss the recent works in decompostion [3-4] and variable selection [5-6] would further strengthen this paper. \n\n[1] A Survey on High-dimensional Gaussian Process Modeling with Application to Bayesian Optimization. ACM TELO, 2022. \n\n[2] High-Dimensional Bayesian Optimization via Tree-Structured Additive Models. AAAI, 2021. \n\n[3] Are Random Decompositions all we need in High Dimensional Bayesian Optimisation? ICML, 2023.\n\n[4] Monte Carlo Tree Search based Variable Selection for High Dimensional Bayesian Optimization. NeurIPS, 2022.\n\n[5] Fast and Scalable Spike and Slab Variable Selection in High-Dimensional Gaussian Processes. AISTATS, 2022.\n\n2/ It would benificial to see the the disscussion (even experimental comparison) with a work with similar topic [6]. \n\n[6] Tree ensemble kernels for Bayesian optimization with known constraints over mixed-feature spaces", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their valuable remarks and are happy to address them.\n\n> 1/ One minor weakness is the the high-dimensional continuous spaces part in Section 2. A recent survey paper [1] categorizes these works into several categories, i.e., low-dimensional embeddings, decomposition, and variable selection. This paper mainly discuss the work in low-dimensional embeddings I think discuss the recent works in decompostion [3-4] and variable selection [5-6] would further strengthen this paper.\n\nWe agree that comparing against decomposition-based techniques is a valuable addition to the paper. We, therefore, compare against the recent RDUCB [1] algorithm that relies on additive decompositions. Furthermore, we compare against the SMAC [2]  Random Forests surrogate model. We will further discuss the suggested papers in the related work section. The appendix submitted in April in the supplementary material contains an extended related work section that discusses Monte Carlo Tree Search (MCTS) based techniques [3, 4]. We plan to integrate this section into the main text if there is sufficient space and will expand the discussion on the recent work on decomposition and variable selection as suggested by the reviewer.\n\n> I am puzzled about why the performance of BODi and COMBO can be sensitive to the location of the optima, which is equivalent to shuffling the labels of the categories of each variable. BODi uses Hamming distance and COMBO uses the combinatorial graph to model the relation of different labels for each variable. Different labels for the same variables are the same, so shuffling the labels will not impact the performance. Can you provide more explanation?\n\nWe hope Appendix D in the supplementary material answers the open questions about BODi and COMBO’s behavior. In the camera-ready, we will refer to this appendix more prominently in the main text since it contains important insights into these methods. Due to space constraints, we cannot move it to the main text. \n\n[1] Are Random Decompositions all we need in High Dimensional Bayesian Optimisation? ICML, 2023 \\\n[2] Sequential Model-Based Optimization for General Algorithm Configuration, LION, 2011 \\\n[3] Learning Search Space Partition for Black-box Optimization using Monte Carlo Tree Search, NeurIPS, 2020 \\\n[4] Monte Carlo Tree Search based Variable Selection for High Dimensional Bayesian Optimization, NeurIPS, 2022"}, {"review_id": "L416EaPUW4", "reviewer": "Reviewer_nHAf", "summary": "The paper considers the problem of optimizing black-box functions defined over high-dimensional combinatorial and mixed continuous-combinatorial spaces. The key idea is to use the Bayesian optimization framework specialized to increasingly large nested embeddings of input dimensions in order to tackle the high-dimensionality challenge. The paper proposes separate embeddings for different class of input variables i.e. continuous, binary, categorical and ordinal. A trust region based approach is employed to enable parallel candidate evaluations in the proposed approach. Experiments are performed on multiple benchmarks to demonstrate the efficacy of the approach.", "questions": "The intuition or principle behind the embedding specific to categorical variable is not entirely clear to me. While for continuous variables, the embedding comes from the idea of count-sketches, it is not clear immediately what is captured by assigning categorical variables of eve different cardinalities to the same bin. For example, if there is an outlier variable with very large cardinality or if the range of cardinalities is quite large, the embedding would be very sensitive to the outlier? Please expand the description about the categorical embedding more if possible.", "limitations": "N/A", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "- The problem space considered in the paper is quite relevant and arises in multiple real-world applications. The paper is written really well describing the scope of the problem and associated applications. \n\n- The proposed idea (although building on BAxUS and TurBO) is principled and explained well. The method clearly outperforms all state-of-the-art baselines on several benchmarks demonstrating its efficacy. \n\n- The related work discussion and the amount of effort to compare with all the baselines in proper way is commendable and deserves credit. \n\nOverall, I think the paper will be useful and practical contribution to this problem space of high dimensional combinatorial spaces. ", "weaknesses": "Please see my question below.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We're grateful for the reviewer's input and will be addressing their comments.\n\nWe expanded the empirical section further and added RDUCB [1] and SMAC [2].\n\n> The intuition or principle behind the embedding specific to categorical variable is not entirely clear to me. While for continuous variables, the embedding comes from the idea of count-sketches, it is not clear immediately what is captured by assigning categorical variables of eve different cardinalities to the same bin. For example, if there is an outlier variable with very large cardinality or if the range of cardinalities is quite large, the embedding would be very sensitive to the outlier? Please expand the description about the categorical embedding more if possible.\n\nIn the camera-ready version, we will add a visualization of the proposed binning procedure. We agree with the reviewer that binning variables of very different cardinalities can lead to undesired effects. Consider the case of ordinal variables. Representing a variable of low cardinality with a high number of ‘levels’ may cause the surrogate model to exhibit variability where there is none in the unknown black-box function: many levels correspond to the same function value since they are mapped to the same x-value of the low-cardinality variable. However, there is a sudden change between two levels if the next level corresponds to another x-value of the low-cardinality variable. Therefore, the surrogate model has wide flat regions with jumps between two levels. We have ideas to tackle this problem but left them for future work. \n\nWe will discuss this limitation in the paper.\n\n[1] Are Random Decompositions all we need in High Dimensional Bayesian Optimisation?, ICML, 2023 \\\n[2] Sequential Model-Based Optimization for General Algorithm Configuration, LION, 2011"}, {"review_id": "KqRZptPEng", "reviewer": "Reviewer_HPPy", "summary": "The paper proposes a new Bayesian optimization algorithm for combinatorial and mixed search spaces containing input variables of different types (continuous, binary, categorical, ordinal), which promises to be relevant for problems as varied as materials discovery, hardware design, neural architecture search, and portfolio optimization.\nBounce relies on a Gaussian process (GP) model of the objective function, which is based on lower-dimensional subspaces of the original search space, generated by partitioning input variables into “bins”.  Notably, the bins only contain input variables of the same type (e.g. continuous) and all variables in a bin are forced to take the same value during the optimization of the acquisition function, in effect operating in a lower-dimensional subspace. During the course of an optimization run, Bounce splits up bins into smaller ones, enabling the algorithm to propose candidates with increasingly finer structure. Bounce also leverages existing techniques for high-dimensional BO, including the trust-region approach.", "questions": "- In the first 60 or so iterations, Bounce is outperformed by BODi on the non-modified 125d MaxSAT problem, as well as the  non-modified PestControl problem. Since sample-efficiency is a chief concern for BO methods, is there a way to leverage  BODi’s inductive biases for Bounce? On a high level, something like this might be possible since both methods leverage embeddings of the input space.\n- Is there something more you can say about any theoretical aspects of the method? Does a global convergence guarantee hold? What can you say about the convergence rate, either local or global? How does the choice of embedding influence the convergence?\n- Are there other embeddings approaches that could be combined with Bounce?\n\nSuggestions:\n\n- The description of the cardinal and ordinal embeddings are a bit repetitive and verbose, especially since they formally use the same approach. Can you unify the presentation, possibly include a single math-mode version of the embedding formula? That would aid the readability of the section.\n- Line 191: “The embedding of ordinal variables follows [that of the] categorical variables”", "limitations": "Yes", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "- Strong empirical results, comparing against published baselines, in addition to ablations on the batch optimization performance.\n- Performance is more robust to shifts of the solution than previously published methods.\n- Sheds light on a non-trivial structural assumption of prior methods.\n- Generally a well written paper.", "weaknesses": "- No theoretical analysis. \n- Without a “modified” or “shifted” label for the left subplots of Figures 1, 2, 3, 4, the presentation is confusing and can easily be misread to suggest that the left subplots report results on the problems as defined in the literature, but they do not. Please add a “modified” or “shifted” label for the left subplots.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate the reviewer's insights and will now discuss their remarks.\n\n> Without a “modified” or “shifted” label for the left subplots of Figures 1, 2, 3, 4, the presentation is confusing and can easily be misread to suggest that the left subplots report results on the problems as defined in the literature, but they do not. Please add a “modified” or “shifted” label for the left subplots.\n\nWe agree with their assessment that the Figure captions can be misleading when scanning the paper. We will add a label as you suggested. We submitted a PDF file with the rebuttal that shows the updated figures.\n\n> In the first 60 or so iterations, Bounce is outperformed by BODi on the non-modified 125d MaxSAT problem, as well as the non-modified PestControl problem. Since sample-efficiency is a chief concern for BO methods, is there a way to leverage BODi’s inductive biases for Bounce? On a high level, something like this might be possible since both methods leverage embeddings of the input space.\n\nRegarding the question about the inductive bias, we would like to stress that BODi is biased towards problems where the optimum lies in the origin or, for categorical problems, where the optimum is realized by setting all categories to the same ‘value’. In Appendix D, in the supplementary material submitted in April 2023, we discuss that most benchmarks in the BODi paper have this property due to their synthetic nature. We do not believe that this bias is reasonable for most practical applications. \n\nHowever, suppose the user has a prior belief that the optimal solution is located in or close to the origin or realized by setting all categorical variables to the same ‘value’. In that case, the user can similarly bias Bounce. We discuss this in Appendix B.3 in the supplementary material and show that this ‘low-sequency version’ of Bounce outperforms or is on par with BODi even on the version of the benchmark where the optimum is at the origin.\n\n> Is there something more you can say about any theoretical aspects of the method? Does a global convergence guarantee hold? What can you say about the convergence rate, either local or global? How does the choice of embedding influence the convergence?\n\nWe refer to Appendix A in the supplementary material, where we prove that Bounce is consistent, i.e., converges to the global optimum in the limit. We did not prove regret bounds and would like to point out that no regret bounds are known for this line of work.\n\n> Are there other embeddings approaches that could be combined with Bounce?\n\nFor the general case, we are not aware of any embedding that can generally be combined with Bounce. The Bounce embedding construction that conveys observations from a lower-dimensional to a higher-dimensional subspace is a key contribution of the paper. However, for the special case of purely continuous problems, Bounce could also use the HeSBO [1] embedding since it also uses a many-to-one mapping of input dimensions to target dimensions. Note that [2] showed that the HeSBO embedding has a lower worst-case probability of containing the optimum, so we opted for the BAxUS [2] embedding construction. \n\n> The description of the cardinal and ordinal embeddings are a bit repetitive and verbose, especially since they formally use the same approach. Can you unify the presentation, possibly include a single math-mode version of the embedding formula? That would aid the readability of the section.\n\nThank you for the suggestion! We agree and will revise the section accordingly for the camera-ready version.\n\n[1] A Framework for Bayesian Optimization in Embedded Subspaces, ICML, 2019 \\\n[2] Increasing the Scope as You Learn: Adaptive Bayesian Optimization in Nested Subspaces, NeurIPS, 2022"}, {"review_id": "OBTn13jm98", "reviewer": "Reviewer_RWfy", "summary": "The paper proposes a new BO method, namely Bounce, to tackle the problem of BO with combinatorial and mixed space. The key idea is based on the trust region approach (as with TurBO) and adaptive space (as with BAxUS), but is applied to the combinatorial and mixed variables. The proposed method is also enabled to work with batch. Experiments are conducted on various synthetic and real-world problems to evaluate the efficacy of the proposed method.", "questions": "Please address my comments in the previous section (section Weakness)", "limitations": "I don’t find any dedicated section describing the limitations of the work.", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "+ The paper’s writing is generally clear and easy to understand\n+ The paper tackles an interesting problem, which is to solve BO problems with categorical and mixed variables\n+ The methodology developed in the paper seems to be sound to me. \n+ The experiments are conducted with various benchmark optimization problems (synthetic and real-world)\n", "weaknesses": "+ The proposed method seems to be largely an extension from existing methods: TurBO and BAxUS – normally this is alright if the results are impressive but there are some issues with the experiments which I will describe in more detail in the below bullets.\n+ There is not much insight on why the proposed method, Bounce, works. There is a theoretical analysis in the appendix (Theorem 1) that shows the consistency of the Bounce. However, I think this theorem doesn’t have much meaning. It assumes that the search domain is finite, and the objective function is noiseless, and it states that Bounce will find a global optimum with probability of 1 when the number of samples N goes to infinity. But for a finite domain and noiseless observations, any algorithm will be able to find the global optimum?\n+ One of the contributions of the paper is to conduct an in-depth analysis of two state-of-the-art algorithms for combinatorial BO, COMBO [45] and BODi [17], however, I don’t really find this to be interesting. Unless I missed it, it seems to be discussed in Section 4.6? But this is largely empirical and not much insight.\n+ The experiments seem to lack of some well-known baselines for categorical and mixed search space, for example, SMAC, TPE, and some new BO methods like the work “Bayesian Optimization over Hybrid Spaces” by Deshwal et al (ICML 2021), and the work “Bayesian Optimization over Discrete and Mixed Spaces via Probabilistic Reparameterization” by Daulton et al (NeurIPS 2022).\n+ Also, in the experiments, the number of iterations is quite small very high-dimensional problems. The number of iterations evaluated in all problems is just 200, which is too small.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We are thankful for the reviewer's valuable input and will ensure their remarks are duly considered.\n\n> There is not much insight on why the proposed method, Bounce, works. There is a theoretical analysis in the appendix (Theorem 1) that shows the consistency of the Bounce. However, I think this theorem doesn’t have much meaning. It assumes that the search domain is finite, and the objective function is noiseless, and it states that Bounce will find a global optimum with probability of 1 when the number of samples N goes to infinity. But for a finite domain and noiseless observations, any algorithm will be able to find the global optimum?\n\nIt is a crucial property of an algorithm to be consistent. We politely disagree with the assessment that any algorithm can find the optimum with an unlimited evaluation budget and a finite domain. In particular, most algorithms in the subspace BO literature (REMBO [1], HeSBO [2], Alebo [3], …) lack this property because they make a random bet on the embedding, and they are unable to recover from choosing a wrong embedding. See, for example, the discussion in [3].\n\nWhile we agree that regret bounds would be an interesting addition, no regret bounds have been proven for this line of work. Nevertheless, our design choices are well motivated, and the fact that Bounce first optimizes over a low-dimensional subspace and eventually reverts to optimization in the input space is crucial to its performance.\n\n> One of the contributions of the paper is to conduct an in-depth analysis of two state-of-the-art algorithms for combinatorial BO, COMBO [45] and BODi [17], however, I don’t really find this to be interesting. Unless I missed it, it seems to be discussed in Section 4.6? But this is largely empirical and not much insight.\n\nRegarding the deeper analysis of BODi and COMBO, we would like to refer the reviewer to Appendix D, which discusses the causes of the performance degradation of those algorithms in detail. In the camera-ready version of the manuscript, we will refer to this appendix more prominently in the main text —  We had to move this detailed analysis to the appendix due to space constraints. \n\n> The experiments seem to lack of some well-known baselines for categorical and mixed search space, for example, SMAC, TPE, and some new BO methods like the work “Bayesian Optimization over Hybrid Spaces” by Deshwal et al (ICML 2021), and the work “Bayesian Optimization over Discrete and Mixed Spaces via Probabilistic Reparameterization” by Daulton et al (NeurIPS 2022).\n\nWe evaluated RDUCB [4], a recent additive method proposed at this year’s ICML, as an additional algorithm. We further added SMAC [5] to the comparison. Similarly to the CASMOPOLITAN [6] paper, we see SMAC performing poorly. We did not compare against probabilistic reparametrization (PR) [7] as we see PR falling in a different category, i.e., a meta-algorithm to optimize the acquisition function. However, we see the potential for improvement by combining Bounce with PR and would like to explore this in the future. \n\n> Also, in the experiments, the number of iterations is quite small very high-dimensional problems. The number of iterations evaluated in all problems is just 200, which is too small.\n\nWe use the same evaluation budget as BODi [8], which also tackles high-dimensional problems, and Bounce converges after 200 iterations on most benchmarks. We want to emphasize that sample efficiency and the ability to find good solutions quickly are important properties of Bounce. To study the behavior for larger sample budgets, we increased the number of function evaluations to 500 for the Labs and ClusterExpansion benchmarks in the figures we submit with the rebuttal. We observe that Bounce continues to outperform the other algorithms and do not find a qualitative difference.\n\n> I don’t find any dedicated section describing the limitations of the work.\n\nFor the camera-ready version, we will discuss limitations more prominently in the discussion section.\n\n[1] Bayesian optimization in a billion dimensions via random embeddings, JAIR, 2016 \\\n[2] A Framework for Bayesian Optimization in Embedded Subspaces, ICML, 2019 \\\n[3] Re-Examining Linear Embeddings for High-Dimensional Bayesian Optimization, NeurIPS, 2020 \\\n[4] Are Random Decompositions all we need in High Dimensional Bayesian Optimisation? ICML, 2023 \\\n[5] Sequential Model-Based Optimization for General Algorithm Configuration, LION, 2011 \\\n[6] Think Global and Act Local: Bayesian Optimisation for Categorical and Mixed Search Spaces, ICML, 2021 \\\n[7] Bayesian Optimization over Discrete and Mixed Spaces via Probabilistic Reparameterization, NeurIPS, 2022 \\\n[8] Bayesian Optimization over High-Dimensional Combinatorial Spaces via Dictionary-based Embeddings, AISTATS, 2023"}, {"review_id": "dioRpW0vAN", "reviewer": "Reviewer_FW4t", "summary": "Most real-world optimization problems contain a mix of continuous, categorical, binary, and ordinal variables and might be high-dimensional. Although there are some methods that have investigated this problem, they are not always reliable with regard to finding a satisfactory optimum. In order to tackle this problem, the authors propose Bounce (Bayesian Optimization Using iNcreasingly high-dimensional Combinatorial and continuous Embeddings). This work seems to be a natural extension of BAxUS, which used a similar embedding strategy as for Bounce but only for continuous variables. It uses a novel trust region management system to grow or shrink the trust regions. The proposed method is tested on a representative range of test problems and compared with state-of-the-art benchmarks and shows convincing results. ", "questions": "- Although I do understand that this is a typical question in BO paper reviews, for most problems, benchmarks are performed for 200 iterations. Would you say that this is enough for the dimensionalities of the benchmark problems? For instance, for the 125-D MaxSAT, you would say you probably need a lot more data points for the surrogate model to model the objective function effectively. What are your thoughts on this and do you think some of the benchmarks can be better than Bounce at a higher budget?\n\n- It is mentioned in Appendix C that the original implementations for COMBO, BODi and Casmopolitan are used. And that you use the same setting as in those works, what kind of settings are these exactly? Do these approaches use the same kernels for the GP models for example? Could this influence the performance on the benchmarks?\n\n- Have the authors contacted Oh et al. regarding the bug in COMBO? If it is possible, it would of course be ideal if results using a bug-less version can be provided.\n\n- I really like the results shown in Figure 5 regarding the efficacy of batch acquisition. In general, I like that the # batch evaluations are used on the x-axis as this directly shows what is good to use for users who are optimizing parallel processes. However, I'm also interested in how these plots look as a function of function evaluations. My intuition would say that b=1 is always better, as in this way the model is as informed as it can for the next iteration. Just out of interest, could you elaborate on this here or could you add a plot to the appendix that shows these results as a function of iterations?\n\n- As a suggestion to your citations regarding chemical engineering and materials discovery. There also is a range of works regarding the optimization of lab equipment for sample analysis. See for instance [Hagan et al.](https://pubs.acs.org/doi/full/10.1021/ac049146x?casa_token=dnIitFW7lO4AAAAA%3Ado8InveDdcPw3TwtOVSQuRvM5NQhSFEo3M1jmpdpEHvRbzA0f1jYTJ2_bloYb7he8-Ofb8u96oMXvCvO), [Boelrijk et al.](https://scholar.google.nl/citations?view_op=view_citation&hl=nl&user=1z-BBwkAAAAJ&citation_for_view=1z-BBwkAAAAJ:9yKSN-GCB0IC). These problems can typically contain mixed spaces and many variables.\n\nSome small textual remarks:\n- Could it be that d_0 and d_{init} are interchangeably used? For instance lines 2 and 3 of Algorithm 1. \n- Line 90, 'sequencies' should be 'sequences' or was this the originally proposed name by the BODi paper?\n- Lines 277-278, as well should be as well as?\n- Typo in Ass. 5. in Appendix Section a. One should be once?", "limitations": "The authors describe the societal impact of their work.\nThe authors do not describe the limitations of their work. Perhaps the authors could dedicate some sentences to this, for example, would the algorithm handle a noisy setting? Or does this violate the binning procedure to some extent?", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- Most real-world optimization problems have mixed spaces and are high-dimensional for which vanilla Bayesian optimization methods are not suited, this work thus addresses an important problem.\n- The proposed algorithm is tested on a broad range of problems and compared with state-of-the-art benchmarks and has a strong performance.\n- The authors provide proof that the algorithm converges in the limit of infinite iterations, making it a reliable method.\n- The paper is well-written and well-structured.", "weaknesses": "- In general, the paper is very complete and well-written and contains a very complete related work section. However, the paper does expect the reader to have a significant amount of knowledge of previous work. Especially about BAxUS and papers using TR management strategies. Perhaps the authors could provide a little bit more information and illustrations in this work to make it more easily accessible. For instance, provide an illustration of the binning procedures and/or about the TR management.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate the reviewer's remarks and will make sure they are appropriately addressed.\n\n> Perhaps the authors could provide a little bit more information and illustrations in this work to make it more easily accessible. For instance, provide an illustration of the binning procedures and/or about the TR management.\n\nThank you for the suggestion! We will add a figure that explains the binning procedure and embedding to the camera-ready to facilitate the understanding of this crucial part of the algorithm.\n\n> Although I do understand that this is a typical question in BO paper reviews, for most problems, benchmarks are performed for 200 iterations. Would you say that this is enough for the dimensionalities of the benchmark problems?\n\nWhen running our experiments, we saw Bounce converging on most benchmarks after 200 function evaluations. A lower evaluation budget is common for discrete BO papers due to the increased cost of optimizing the acquisition function. For instance, the authors of BODi [1] worked with the same evaluation budget. However, we agree that running a higher evaluation budget would be interesting for the Labs and MaxSAT125 benchmarks where Bounce has not yet converged. Due to the high cost of these experiments, we will add these analyses to the camera-ready version.\n\n> It is mentioned in Appendix C that the original implementations for COMBO, BODi and Casmopolitan are used. [...] Do these approaches use the same kernels for the GP models for example? Could this influence the performance on the benchmarks?\n\nWe agree that the choice of the kernel is a key aspect of these methods and is instrumental to their performance. Thus we use the exact settings and implementations of the respective authors. The core contributions of BODi and COMBO are their kernel constructions based on dictionary and diffusion kernels, respectively. We do not think replacing the kernel for these methods is reasonable. \n\nNote that Bounce uses the CoCaBO [2] kernel construction, also used in CASMOPOLITAN [3]. We will clarify this in the appendix.\n\n> Have the authors contacted Oh et al. regarding the bug in COMBO? If it is possible, it would of course be ideal if results using a bug-less version can be provided.\n\nWe recently contacted the COMBO authors and are currently waiting for a response. We fixed the bug in COMBO independently, and the PDF uploaded with this rebuttal contains an updated figure for the PestControl benchmark, showing the original version of COMBO and a fixed version (marked with “(fixed)”). We observe that fixing this bug improves COMBO’s performance on this benchmark considerably so that COMBO outperforms BODi on the modified benchmark. As expected, fixing the bug makes COMBO agnostic towards the modification of the benchmark. Note that Bounce still outperforms all other algorithms. PestControl is the only benchmark in our experiments where COMBO is only affected by this bug (see Appendix D.2). We will update our discussion to reflect these new findings.\n\n> I really like the results shown in Figure 5 regarding the efficacy of batch acquisition. [...] However, I'm also interested in how these plots look as a function of function evaluations. [...] Just out of interest, could you elaborate on this here or could you add a plot to the appendix that shows these results as a function of iterations?\n\nA figure with the number of function evaluations on the x-axis for the batched version of Bounce is an interesting addition. As expected, large batches perform ‘worse’ when plotting against the number of function evaluations. There is almost no difference between small batches, highlighting the efficacy of our batching strategy. We added this figure to the PDF in the global response and will add it to the appendix in the paper.\n\n> As a suggestion to your citations regarding chemical engineering and materials discovery. [...]See for instance Hagan et al., Boelrijk et al.\n\nWe agree on the relevance of this application area and will gladly add Hagan et al. and Boelrijk et al. as references for the practical applications.\n\n> Could it be that $d_0$ and $d_{\\textrm{init}}$ are interchangeably used? For instance lines 2 and 3 of Algorithm 1.\n\n$d_0$ and $d_{\\textrm{init}}$ are the same variable. We will revise this and keep only one of them.\n\n> Line 90, 'sequencies' should be 'sequences' or was this the originally proposed name by the BODi paper?\n\nThe term ‘sequencies’ was introduced in the BODi paper [1] and describes the number of changes from 0 to 1 and vice versa in a bit vector — This is equivalent to ‘frequency’ over a continuous domain. We will clarify this.\n\n> The authors describe the societal impact of their work. The authors do not describe the limitations of their work. Perhaps the authors could dedicate some sentences to this, for example, would the algorithm handle a noisy setting? Or does this violate the binning procedure to some extent?\n\nWhile we discuss limitations in the discussion, we will expand this section, as suggested. We did not consider noisy function evaluations. In future work, we would like to explore the performance of a variant of Bounce on noisy problems where we replace EI or qEI with noisy EI or noisy qEI [4]. With these changes, we expect Bounce to also perform well on noisy problems.\n\n[1] Bayesian Optimization over High-Dimensional Combinatorial Spaces via Dictionary-based Embeddings, AISTATS, 2023 \\\n[2] Bayesian Optimisation over Multiple Continuous and Categorical Inputs, ICML, 2020 \\\n[3] Think Global and Act Local: Bayesian Optimisation for Categorical and Mixed Search Spaces, ICML, 2021 \\\n[4] BoTorch: A Framework for Efficient Monte-Carlo Bayesian Optimization, NeurIPS, 2020"}]

Author rebuttal:
We thank all the reviewers for their valuable feedback. We are pleased that the reviewers appreciated the problem's relevance, as acknowledged by reviewers ACyR, nHAf, and HPPy, along with their positive assessment of Bounce's performance across various benchmarks, as highlighted by reviewers ACyR, nHAf, HPPy, and FW4t.  We are happy to read that the reviewers, including reviewers nHAf, HPPy, and FW4t, found the paper to be well-written. 

Several reviewers asked for a deeper analysis of the causes of the performance degradation we observed for BODi [1] and COMBO [2]. We want to point the reviewers to the in-depth analysis in Appendix D of the supplementary material that provides such an analysis — please notice that this year’s author guidelines do not allow the supplementary material to be submitted together with the main text. 

Some reviewers asked for a comparison with additional methods. We submit a PDF file with the rebuttal that shows the updated figures. We added an evaluation of RDUCB [3], a recent algorithm presented at ICML 2023 that belongs to the class of additive methods. We also added a comparison with the well-known SMAC [4]. Bounce outperforms both the recent RDUCB algorithm and SMAC on all benchmarks.

Furthermore, we increased the number of function evaluations to 500 for the 50D-Labs and 125D-ClusterExpansion benchmarks. Due to the high computational cost, we could only run ten repetitions for BODi but will increase the number of repetitions for the camera-ready version. We don’t expect the results to change significantly since BODi has shown low variability in our experiments. Furthermore, not all 500 iterations have finished yet for COMBO due to its large running time. The plots show the mean averaged over  fifty runs. We observe that Bounce continues to outperform the other algorithms.

Based on the comment by reviewer FW4t, we fixed COMBO’s bug that affected its performance on the PestControl benchmark. Fig. 1 (c) in the PDF shows both the performance of COMBO with (marked “COMBO” in the legend) and without the bug (marked “COMBO (fixed)” in the legend). Fixing the bug improves COMBO’s performance on this benchmark and, as expected, makes COMBO agnostic towards the modification of the benchmark. Note that Bounce still outperforms COMBO and all other algorithms.

We note that we contacted the authors of BODi, and they acknowledged the sensitivity of BODi toward the location of the optimal point. We further contacted the authors of COMBO and are currently awaiting a response.


[1] Bayesian Optimization over High-Dimensional Combinatorial Spaces via Dictionary-based Embeddings, AISTATS, 2023 \
[2] Combinatorial Bayesian Optimization using the Graph Cartesian Product, NeurIPS, 2019 \
[3] Are Random Decompositions all we need in High Dimensional Bayesian Optimisation?, ICML, 2023 \
[4] Sequential Model-Based Optimization for General Algorithm Configuration, LION, 2011

Meta-review: This paper considers the important problem of optimizing expensive black-box functions over high-dimensional combinatorial spaces with many real-world applications. The paper proposes an effective Bayesian optimization algorithm using a combination of ideas including BAxUS and trust-region methods to handle the high-dimensional challenge, and different embeddings for continuous, discrete, and ordinal variables. The paper showed good empirical performance on multiple benchmarks.

All the reviewers' were positive about this paper. They raised some clarification questions and the authors' rebuttal addressed many of them satisfactorily. 

I recommend accepting the paper and strongly encourage the authors' to incorporate reviewers' comments and rebuttal discussion into the final paper.
Decision: Accept (poster)