Tree ensemble kernels for Bayesian optimization with known constraints over  mixed-feature spaces
Authors: Alexander Thebelt
Calvin Tsay
Robert Matthew Lee
Nathan Sudermann-Merx
David Walz
Behrang Shafei
Ruth Misener
Venue: NeurIPS 2022 Accept (2022)
Keywords: Bayesian Optimization
Tree Ensembles
Global Optimization
Known Constraints
Black-box Optimization
Mixed-Variable Spaces
Hybrid Spaces

Abstract: Tree ensembles can be well-suited for black-box optimization tasks such as algorithm tuning and neural architecture search, as they achieve good predictive performance with little or no manual tuning, naturally handle discrete feature spaces, and are relatively insensitive to outliers in the training data. Two well-known challenges in using tree ensembles for black-box optimization are (i) effectively quantifying model uncertainty for exploration and (ii) optimizing over the piece-wise constant acquisition function. To address both points simultaneously, we propose using the kernel interpretation of tree ensembles as a Gaussian Process prior to obtain model variance estimates, and we develop a compatible optimization formulation for the acquisition function. The latter further allows us to seamlessly integrate known constraints to improve sampling efficiency by considering domain-knowledge in engineering settings and modeling search space symmetries, e.g., hierarchical relationships in neural architecture search. Our framework performs as well as state-of-the-art methods for unconstrained black-box optimization over continuous/discrete features and outperforms competing methods for problems combining mixed-variable feature spaces and known input constraints.

Tree ensemble kernels for Bayesian optimization with known constraints over mixed-feature spaces
Alexander Thebelt Imperial College London London, UK
Calvin Tsay Imperial College London London, UK
Robert M. Lee BASF SE Ludwigshafen, Germany
Nathan Sudermann-Mers Cooperative State University Mannheim, Germany
David Walz
Behrang Shafei
BASF SE
BASF SE
Ludwigshafen, Germany
Ludwigshafen, Germany
Ruth Misener
Imperial College London London, UK
Abstract
Tree ensembles can be well-suited for black-box optimization tasks such algo- rithm tuning and neural architecture search, as they achieve good predictive perfor- mance with little no manual tuning, naturally handle discrete feature spaces, and are relatively insensitive to outliers in the training data. Two well-known challenges in using tree ensembles for black-box optimization are (i) effectively quantifying model uncertainty for exploration and (ii) optimizing over the piece-wise wise constant acquisition function. To address both points simultaneously, we propose using the kernel interpretation of tree ensembles as Gaussian Process prior to obtain model variance estimates, and we develop a compatible optimization formulation for the acquisition function. The latter further allows constraints to improve sampling efficiency by considering domain-know knowleeg engineering settings and modeling search space symmetries, e.g., hierarchical relationships in neural architecture search. Our framework performs as well as state-of- offtth-artttt methods for unconstrained black-box optimization over continu- ous/discrete features and outperforms competing methods for problems combining mixed-variable feature spaces and known input constraints.
seamlessly integrate known
Introduction
Many black-box optimization problems contain feature relationships known priori based on domain knowledge, such as hierarchies or constraints [5] 68 For example, hierarchical structures arise in neural architecture search 72) where hyperparameters such as kernel size are only relevant a convolutional layer is selected. Explicit constraints can also arise, e.g., matching kernel size and stride to input channel size and padding. In many cases, Bayesian optimization can incorporate known hierarchies and/or constraints, given that a suitable surrogate model selected. To this end, Fromont et al. 2 and Nijssen and Fromont [33] impose variable hierarchy by constraining the splitting order of decision trees, i.e., certain attributes must selected before others.
Tree-based models, such as random forests or cradieet-boostee trees, remain popular in many appli- cations, as they inherit the innate ability of simple decision trees seamlessly handle categorical and
"Coresponding author: alexander.theeell11888. imperiall........kkkkkkkk
36th Conference on Neural Information Processing Systems (NcurIPS 2022).

---

discrete input spaces, Moreover, they are highly parallelizable and scalable to high- gigh-diennsional data. Despite these modelling advantages, the deployment of tree-based models in Bayesian optimization has been limited by challenges pertaining to (i) quantifying prediction uncertainty and (ii) optimizing acquisition functions defined by their discontinuous response surfaces [61]. Early works, e.g., the popular SMAC algorithm [34], addressed (i) using empirical variance within tree ensemble and (ii) via local and/or random search methods. Moreover, recent works 44 50 66 propose mixed-integer formulations for tree ensembles, enabling optimization over their mean functions.
Contributions. Sections and B present related work and methods used to derive the approach proposed in this paper. We present mixed-integer second cone optimization formulation for tree kernel Gaussian processes in Section 4 and show that that the free agreement ratio, i.e., the hyperparameter introduced by the tree ensemble kernel. sufficiently represents the model uncertainty in Section 5.1] Section 5.2 shows that solving the mixed-integer second-order cone optimization problem considerably outperforms sampling-bane -baned strategies in Bayesian optimization. Our approach of using tree ensemble kernels as Gaussian process prior particularly useful for applications combining mixed-varie variable spaces and known input constraints. A Python implementation of the proposed algorithm is available at: www.gi com/ cog- wipemimalttreeeeeee ceree_eeeeeel_el gp
2 Related work
x E arg max fxx), where f is an expensive-to- evaluate black-box function that can be queried at inputs X' derive the optimal solution x BO iteratively updates a surrogate model of f and optimizes corresponding acquisition function that balances exploitation and exploration. Maximizing the acquisition function produces a new query x" common
Bayesian optimization (BO) solves 22
61
which is evaluated and added to the set of observations. Gaussian processes (GPs) [57]] choice of BO surrogate due to their flexibility, e.g., domain-spectin knowledge can be built into the GP prior via mean and kernel functions, and reliable uncertainty quantification identify unexplored implement BO with GP surrogates and offer a wide selection of kernels mainly suited for continuous search spaces. Prior works developed GPs with modified kernels to integrate discrete features IC
search areas. Open-source tools such as BoTorch [2]
or considered conditional feature Nguyen et al. 152 integrate catgorical and category-specific continuous inputs by formulating the black-box optimization problem as multi-arm bandit problem for which each category corresponds an arm. Similarly, Gopakumar al. 26 handle mixed-type inputs using multi-armed bandits. Besides GPs, tree ensemble- surrogates show excellent performance for black-box optimization with vixid-variable settings, i.e., with continuous, integer and categorical variables, and for structured search spaces, e.g., hierarchical and conditional feature spaces 64]
25
30 59
spaces 2
33
f
4
Black-box optimization tools using tree ensembles, e.g.. SMAC 34 and Scikit-Optimize (SKOPT) are useful for applications such as neural architecture search (NAS) and algorithm tuning. Shahriari al. [6I) mention challenges in deploying tree ensembles for BO: (i) quantifying uncertainty for exploration purposes, and (ii) optimizing over the non non-differentiiale discrete acquisition function to determine the next query point. SMAC identifies uncertain search space regions using empirical variables across tree predictions of the random forest and optimizes the acquisition function combining local and random search. Bergstra et al. 5 proposes the Tree Parzen Estimator (TPE) to handle categorical variables and conditional structures by modeling individual input dimensions by a kernel density estimator. However, the TPE approach ignores dependencies between dimensions. For gradient-boosted tree ensembles, SKOPT derives uncertainty with quantile regression to fit two models for the 16th and 84th percentile and averages the predictions to estimate the standard deviation. For random forests, SKOPT uses an uncertainty strategy similar to SMAC. In general, SKOPT relies on random sampling to optimize the acquisition function.
Misié 49 proposed a mixed-integer optimization formulation for tree ensemble mean functions that has been used in several applications [|]] 50 65-667 Besides improving the solution to an acquisition function, mixed-integer formulations also allow explicit consideration of input constraints to incorporate domain knowledge. Some software tools, e.g.. BoTorrh, support linear equality and inequality constraints of continuous variables at the acquisition function optimization step, while tree ensemble-base algorithms do not support input constraints. Papalexopoulos et al. [54] use ReLU neural networks as surrogate models and deploy mixed-integer linear formulation to optimize the
acquisition function. The approach relies on random initialization and stochasticity in the model
training allow for exploration. Daxberger et al. [6 handle -axiiabee search spaces by using

---

a Bayesian linear regressor that uses an integer solver to search the discrete subspace. The authors introduce features capturing the discrete parts of the search space by using BOCS model B
while continuous parts are handled with random Fourier features [56]. Genetic Algorithms (GA) are another class of algorithms, which deploy evolution-based selection heuristics to maximize black-box
functions KG While there is no feasibility guarantee for input constraints, GA implementations like pymoo i support constraint optimization by minimizing constraint violation.
We compare our BO approach, which uses the kernel interpretation of tree ensembles as a Gaussian process prior, to: SMAC, the random forest and madient-bboster tree versions of SKOPT, (SXOPT-RF, and SKOPT-GBRT, respectively) in the Section S numerical studies as baseline for other tree ensemble-based algorithms. We also compare against the default upper-confident pected improvement BO implementations of BoTorch (UCB- CEE-MATERN and EI SI-MATERN, respectively) and the default algorithm of pymoo to include black-box algorithms that partially support con- strained optimization.
bound and ex-
3 Technical background on prior work
3.1 Tree ensemble kernel as a Gaussian process prior
Exploring the search space requires quantifying the uncertainty of the underlying surrogate model. We use the kernel interpretation of tree ensembles based on random partitions [5] 73 The tree kernel captures correlation between two input data points (x, x') R":
Teeee (x, x') = o [T|- Txx'')
(1)
To derive the tree kernel, we first train a gradient boosted tree ensemble T on data set X E R=** with denoting the dimensionality of the search space and m the size of the data set. Every tree in the tree ensemble T maps inputs x E IR" onto a leaf by sequentially evaluating splitting conditions. Each leaf 1 defines a subspace x; C R" restricted by active splits E splits(t). Two inputs (x, xx, are fully correlated in tree t if both end up in the leaf subspace Xt,1 and uncorrelated otherwise. The Eq. vector zxx) consists of binary elements 2 indicating if leaf [ Le is active for input x, with C, denoting the set of all leaves in tree t. The inner product 2 (x) ( z(x'), normalized by the total number of trees (7). gives the ratio of trees in the ensemble for which and x' fall into the same leaf. We modify the kernel by adding a trainable signal variance of 42 Note that the tree ensemble trained separately. Davies and Ghahramani 15 prove that the tree kernel is suitable GP prior. The resulting (non-stationary and supervised) tree kernel describes prior over piece-wise constant functions when used in a GP.
defining the kernel and the kernel hyperparameters
3.2 Posterior distribution
We approximate ff as a Gaussian process with zero mean and kernel kive: f(.) ~ GP((, Knee). Since Treee is a valid Mercer kernel, the mean (x) and variance (xx) of the GP at R" is 55
M(x) Kx.x (Kx,x)=1 y
(2a)
V V(x) Kx,x Kx.x (Kx,x) KI
(2b)
The Gram matrix Kx.x x E x has entries describing pairwise correlations computed based on the kernel function in Eq. ①. The entries of vector Kx,x E IR contain correlations between the input x and sampled data points, defined as [haree (x, x1), Kivee (X,, xx)).... &ive (x, x,,)] with vectors x; referring to in data set X. Target values y I the corresponding observations for X. Eq. 2 describes the noise-free case of the GP mean and variance. To fit GP based on the tree kernel function, we usually require noise term, i.c.. diagonal matrix o that is added to Kx,x- We set hyperparameters oū and oi by maximizing the log marginal likelihood. Lee et al. 4p compute the Eq. (2a) inverse the Gram matrix Kx.x efficiently by exploiting the property that the rank of kx.x is at most the number of leaves over trees. Fig. [ visualizes the Eq. 2 M (x), V(x) and the upper-confidence bound (UCB) [4] response surface. Fig. I(a) shows that the tree kernel-based GP mean gives a good piecewise-constart
approximation of negative one times the Branin function. Fig. shows that variance peaks reveal areas where data are sparse, reliably identifying uncertainty in the underlying surrogate model. Finally, Fig. shows how an acquisition function such as the UCB can effectively manage the iegloitationeexplratite trade-off.
3

---

(a) (x)
(b) V(x)
(c) UCB(x)
Figure I: Tree kernel GP trained on 40 random points of negative one times the Branin function for the intervals (I-5.0, 10.0], [0.0, 15.0))T. Function values increase with the colour brightness and white contour lines indicate the true functional shape of the negated Branin function.
3.3 Global optimization of tree ensembles
proposes mixed-integer linear optimization formulation that ensures that binary variables Zr, follow the logic of the tree ensemble:
Mišić 49
Σ 2t,1 =
11 E T, (3a)
JEL,
Σ
24,1 S
Σ W(()
At E T, splits(t),, (3b)
(Eleft(s)
jEC(s)
Σ
2, 1
Σ
W(s)dh
1A T, splits((), (3c)
(Eright(s)
jEC(s)
Σ Vi.i 1,
Ei (3d)
Vi.j S Miit+1
WE [K; - 1], (3c) vi [n] Vj [K]], (3f)
M/ E {0,1},
2x 2 0,
VIET, C. (3g)
Eq. Bap guarantees exactly one active leaf I in leaf set L, for tree t. Eqs. 3b
3c
ensure that
binary variables žτ, are only active all previous split binaries Pv(e), corresponding to continuous splitting thresholds are active. At any node in given tree, left and right(s) contain all leaves following the left and right branches, respectively. The mapping V(s) gives the splitting feature at node & E splits(t) in tree f, with splits(t) defining the set of all splits in tree f. Tree ensembles can handle continuous, integer and categorical data. Continuous splits conditions, where (v(s). is the learned splitting threshold. Therefore, C(s) only contains a single index j representing the threshold (V(()).- Categorical splits s are characterized by (subsets of) categories available for feature V(s) and define a splitting condition based on the inclusion of (v(()) in the category subset at split s. For categorical splits, C(s) includes the category indices comprising the category subset at split s. Eq. (3d) ensures that only one category is active per categorical variable i E C. Continuous splitting thresholds of all trees in the ensemble are ordered according to U.KK with K; denoting the index for the last split of continuous feature E N. To enforce this order, Eq. Beh ensures that binary variables v: Vijj corresponding to the split thresholds Ua.j. are activated sequentially. The model comprising Eqs. 3 has no direct dependency
defined by JV((s) S (V(s)
U111 U22
x and is fully defined by binary variables indicating which splits and leaves of the tree model active. However,
to allow the user to include extra equality and inequality constraints on the input vector x, we bound

---

the continuous variables based on the active splits by adding linking constraints 50
K (ii - 115-1) 1) (1 -vi)),
32 2
WWE,, (4a)
K
34 S o
(Ri,) vii+t)
Vi EiEW, (4b)
==1
32 [ot vi T; Lj E [Ki] Pb. 1}.
Vi E.,, (4c)
e C, (4d)
For categorical variables, Eq. (4d maps the indices of active categories onto 2; for C. These
optimization formulations are implemented in open-source software ENTMOOT 666 and OMLT
Tree ensemble kernels for Bayesian optimization
The technical details in Section 3 are insufficient to use the tree ensemble kernel in Bayesian optimization framework. While Eqs. S and allow optimization over the mean an associated acquisition function 49), Bayesian optimization also requires quantifying model uncertainty for exploration. This section proposes a mixed-integer second-order cone optimization formulation to capture the standard deviation of GP with tree kernel prior. Combining this optimization formulation with the already-develope the upper-confidence bound (UCB) of the tree kernel-based GP. The advantage of deriving mixed- integer second-order cone optimization formulation of the UCB acquisition function is that we can globally optimize the acquisition function. Additionally, we can easily incorporate explicit input constraints that capture domain knowledge and/or known search space relationships.
mixed-integer formulation of the mean function, we derive
Our optimization problem, which includes the UCB acquisition function, is:
Xii X X E arg max puxx) + kx()).
(5a)
h(x) 0, g(x) S O,
(5b)
(5c)
with μ(x) and (xx) denoting the surrogate model's mean prediction and standard deviation, re- spectively. Functions h(x) and gxx) are known constraints and handled similarly to E (9 In our implementation, the constraints can be linear, quadratic, or polynomial. Hyperparameter k 2 0 controls the peploitation-exploaaaaee
trade-off to determine the next black-box function query area. The solution to Eq. 5a) is defined by X and x for non-categorical variables i E ie.. continuous and integer features, and X a set valid category subsets for categorical variables i C. Two vectors (xj, and x define the non-categorical variables because the trees are piecewise-constant over intervals and the vectors define the lower and upper bounds these intervals.
Next, we formalize the mean and variance of the tree kernel- based GP. Eqs. 6 and are equivalent.
m(x)= M (x) kx.x (Kx.))-1 y
(6a)
x²(xx Vxx) = Kxxx Kx,x (Kx,x) KT x.X
(6b)
The Gram matrix Kx.x and target vectory are constants in the optimization model, as these quantities only depend on the data set X. The value Kx,x is directly related to signal variance hyperparameter og since there is full leaf overlap for two identical inputs:
Kx.x = o
The vector Xx.x contains the kernel output of x with individual data points x,- We compute tree is trained but before solving the
the constant matrix A E m= after the gradient-booste acquisition function. The entries of A are equal to 1 for all active leaves of data point and 0 otherwise. Given matrix A, Eqs. (8) capture the kernel output by summing over binary variables 21,1 that are active for data point i.
Xxxx [ x1), KIcc (x, x2), ... A Xm)]
(8a)
5

---

( oū |71-
'ET'EE- SET IEC,
A,,22
vi E [ml]
(8b)
To get an intuition for Eqs. CB. note that the largest possible value for each Teee (x. x;) oe (if the corresponding entries of matrix A are all equal to I, and there is full leaf overlap) and the smallest possible value for each horre (x, xi) O (if there is no leaf overlap). In general, values of each element of vector xx,x will range between O and oß: higher values in the elements of Xx.x indicate a higher degree of overlap between the next query location x and data point x; in the set X.
Without considering additional tree model constraints h and g. the resulting acquisition function is a mixed-integer quadratic optimization problem. The optimization problem mixed-intess because of binary variables v and quadratic because of Eq. (66). The quadratic Eq. (bb) components are o² and the terms frre [x X;) arising from the inner product of xx,x with itself. We only require one direction of the Eq. (6b) equality (<) and re- re-write Eq. (bb) as second-order cone constraint. Second-order cone programming 41] optimizes over a linear objective subject to both linear and second-order cone constraints (here, convex quadratic constraints, but the theory is more general), More recently, solvers integrate advanced methods solving second-order cone problems in the mixed- integer setting [4] 20
[6 Solver Gurobi 9 [2X] automatically finds that Eq. 6b (with S rather than =) can be represented as a second-order cone and makes the appropriate algorithm modifications, e.g,, as described by E3 [T). The proposed formulation is also compatible with open-source solver alternatives including Bonmin
MindtPy 6 Pajarito [3. and SHOT
47
The acquisition function, with Objective Sa and Constraints 3
8b is a mrixed-integer second-order cone program which can be solved with optimization solvers. valid solution the proposed model is a set of active leaves L and the intersection of all leaf subspaces [xi=
x
Xeat-
To derive the next black-box function query point x", propose some heuristics. The acquisition function value is constant for x E [ x X i.e., all contained points are equivalent from the perspective of the tree kernel-GP. For the continuous and integer features, we observe that tree models tend to learn split thresholds close to training data points. We propose the center of [xip Xinl as the next query point x* for continuous and integer features:
m = 2 (dib, + abb.)) V.
For integer features with fractional mid- mid-point between the upper and lower bound, we randomly select its floor or ceiling. For categorical features, we sample from the subset of available categories:
mid, uniform(x(a)) Vi E C.
(10)
For the unconstrained case, we use ii as the next query point x". When additional input constraints h(x) and/or g(x) are given, we know that at least one point in [xj,, xtpl x is feasible even if the heuristic X is infeasible. To repair the solution X we project it onto the feasible space:
x" E arg min
Σ (#ind, - xz) Σ
Σ
ViJ
(1la)
xc xE(Mi>>i]>>< iEN
iEC jE
h(x) 0,
(11b)
g(x) S O.
(11c)
Eq. (11) projects X onto the feasible set defined by h(x) and g(x). The time complexity of
solving Eq. 11 and that of solving Eqs. 3 - for the search space [x x x are both NP-hard.
However, in preliminary evaluations we found that Gurobi 9 often solves both problems to e-global optimality for moderately-sized tree models in less time compared to random sampling.
Hyperparameters. We train the kernel hyperparameters signal variance of and noise of by maximiz- ing the log marginal likelihood. Additional hyperparameters are introduced by the gradient-booste tree ensemble trained at every iteration, i.e., maximum tree depth and number of trees, and through e in the UCB. For the Section numerical studies, we leave k and gradient-boosting hyperparame- ters constant and only increase maximum tree depth and number of trees for the high-dimensional CIFAR-NAS benchmark to capture more complicated interactions. The appendix reports specific values for all hyperparameters.
Limitations. The method suffers from standard BO challenges where the tree kernel may not be good prior for the underlying black-box function, e.g., purely continuous feature spaces. More

---

limitations arise from solving an NP-hard problem to global optimality when working with large data sets in nigh-dimensional search spaces. For cases where solving the NP-hard problem difficult but BO is still applicable, Gurobi 9 can typically develop good feasible solutions as a heuristic.
5 Numerical studies
This section empirically evaluates the performance of tree kernel-based GPs using a wide variety of synthetic and real-world benchmark problems. We show (i) the ability of tree kernels to capture uncertainty of the underlying tree ensemble, (ii) the advantage of using global vs. local strategies for optimizing the acquisition function and (iii) the proposed algorithm's superior performance cases with constrained search spaces and mixed variable types. LEAF- HEAF-GP denotes the Section A proposed ootlinee the baseline of methods we compare against. For every run, visualize the median and confidence intervals of the first and third quartile based on 20 individual runs with varying random seeds. Further technical details be found in the appendix.
algorithm, and Section 2
5.1
Uncertainty metric
1200
2.0
Relative Prediction Error
1.5
Prediction Mean
800
1.0
600
0.5
0.0
200
-0.5
0.4
0.5
0.7
0.8
0.9
0.4
0.5
0.7
0.8
0.9
Tree Agreement Ratio
Tree Agroement Ratio
(a) Relative Model Error
(b) Prediction Mean
Figure 2: The relative prediction error (Eq. 13 and model prediction mean over the maximum tree agreement ratio R for benchmark problem Rastrigin (10D). Changing R is equivalent changing the maximum kernel covariance. Plot shows the median line and confidence intervals (first and third quartile) from 20 random seeds. Section 5.1
provides more details.
The tree kernel-based GP uses the leaf overlap measure to quantify correlation between two inputs. In BO, such measures help identify unexplored areas where correlation existing training data is low and we expect inaccurate model predictions. To empirically test the tree kernel's capability of identifying uncertainty in the underlying tree ensemble, we change the optimization formulation:
Xio Xbb X
Earg max (xx), X.8.M
(12a)
Eq.
Eq.
Eq.
Eq.
Eq.
(12b)
Σ
Σ
S
R. Vi E [m]
(12c)
BETIEC.
Eq. (2a) maximizes over the mean prediction of the tree kernel-based GP with Eq. 12ch restricting the ratio of tree agreement for the black-box query point x* defined as the Eq. center of the optimal area x X x A solution problem Eq. 112 guarantees maximum leaf overlap of R between x" and all available data points X. Limiting the tree agreement ratio corresponds to constraining the Eq. Q kernel correlation between the existing dataset and the proposed optimal area Xio, X X Here, we evaluate values of R E [0.35, 1.0] with increments of 0,05 and compute the model error according to:
[(x") firue
Cervor
(13)
μ(x*)
Fig. shows results for the Rastrigin optimization problem becomes more restricted, leading smaller objective values for solution
63
benchmark. With smaller values of R. the Eq. 12
Fig. a shows that growing leaf overlap (increasing R) reduces the model error defined in Eq. [3

---

This suggests that the kernel works intended and that low kernel correlation can reveal search space areas with high model uncertainty. The appendix provides results for additional benchmarks.
5.2 Local vs. global acquisition optimization
3.0
350
9 Back Box Objective
Back-Box Objective
300
2.0
250
E-MATERN
1.5
LEA--PP
200
LEAF-GN-RNC
1.0
150
SMAC
0.5
UCB-MATEIN
0
50
100
150
200
100
200
300
4000
500
teeratoon
Iterations
(a) Hartmann (6D)
(b) Styblinski-Tang (10D)
Figure 3: Black-box optimization progress of LEAF-GP baseline. Plot shows the median line and confidence intervals (first and third quartile) from 20 random seeds. Section 5.2 provides more details
This section (i) compares LEAF-GE to other state-of-the-art algorithms for common benchmarks and (ii) shows the advantage of global vs. local strategies for optimizing the acquisition function. We introduce a variation the proposed algorithm, LEAF-GP-RIE which optimizes the same acquisition function LEAF- GP, but
random sampling instead of Gurobi 9. Fig, shows results for Hartmann (6D) 6 and Styblinski Tang (10D) 55.. We do not expect tree model-based approaches to perform well since both benchmark functions are continuous. However, the benchmarks show that all approaches perform similarly for Hartmann (6D) with MOB-MATERN and EI MITMRTERN improving the black-box objective at the fastest rate. For Styliiski-Taag (10D), LEAF-GP significantly outperforms other algorithms. Although LEAF BEAF-GP is not specialized to this setting, observe that performs roughly equivalently to the state of the art these continuous, unconstrained optimization problems.
Although sampling-based LEAF-GP-RND and optimization LEAF-GP perform similarly on smaller Hartmann (6D), LEAF-GP is particularly strong the higher-dimensional Styblinski-Tang (10D) benchmark function, and two additional benchmarks the appendix. Moreover, LEAF-GP-RMN does not support explicit input constraints and performs particularly bad for BO with known con- straints. For the following constrained benchmarks, we remove LEAF-GP RRD from the comparison. The appendix has additional details and results.
5.3 Constrained spaces
This section presents numerical benchmarks with known input constraints. The acquisition opti- mization strategy of LEAF-GP allows explicit consideration of input constraints, i.c., logical and convex/non-convex n-th degree polynomial equality and inequality constraints for mixed variable spaces. The UCB-MATERN and EI -MAMRR implementations BoTorch supports only linear equality and inequality constraints at the acquisition optimization step. The GA pymoo algorithm has an interface for callable constraint functions which are considered when generating new generations of candidate points. When method does not support the specific input constraints, we penalize the objective:
(14) where (xx < 0 and h(x) = 0 are inequality and equality constraints, respectively. This penalty strategy allows methods that do not support explicit constraints to still produce feasible points given enough iterations. Eq. (14 introduces the hyperparameter À which weights the penalty, we test values E (1, 10, 100} and only plot the best run for all methods that rely on constraint penalization. To initialize every method with feasible points, we draw random samples from uniform distribution within the bounds that define the search space and compute the closest feasible point similar to For constrained benchmarks, introduce -AANDRR which simply projects random
fpenaity A (max(g(xx, h)x)+(())))))
Eq. II

---

Black-Box Objective
14.5
0.8
14.0
0.7
D-MATERN
E-MATERN
13.5
LEAF
13.0
SMAC
UCSMMATER
0.4
UCB-ATTER
20
40
60
80
100
20
40
60
80
100
120
Iterations
lterations
(a) (13D, 9IC)
(b) G3 (5D, IEC)
1
30.6K
Back-Box Objective
30.4K
Black Box Objective
1 1.6K
30.2K
30.0K
LI-AATE FEAS-RND
MATTEA
FEAS-AND
GA
29.6K
LEAF-G"
2K
LEAP-G"
SKOT1-GBBT
NOFF-C8M1
SKOF1-R2
SCOP1-9
29.4K
SMAC
SMAC
UCB MATERN
1K
UESMATERN
20
40
60
80
100
40
100
Iterations
Iterations
(c) G4 (5D, 6IC)
(d) Alkylation (7D, 14IC)
Figure 4: Feasible black-box optimization progress of LEAF-GP vs. baseline. Plot shows the median line and confidence intervals (first and third quartile) from 20 random seeds. Confidence intervals are neglected for methods that cannot improve the initial training data. Figure subtitles give the function name and number of: dimensions (D), equality constraints (EC), and inequality constraints (IC). Section 5.3 provides more details.
samples onto the set of feasible points that satisfy h(x) and g(x). Fig. plots feasible solutions
to four different continuous benchmark problems, G1, G3, G4 B2 and Alkylation 60) The G1 benchmark has linear inequality constraints only, which are supported by NCB-MATERN, MATERN and LEAF-GP Fig. 4(a) shows LEAF-GP making quick progress at the beginning with UCB-MATERN catching up towards the end. G3 and G4 have different combinations of nonlinear constraints, and LEAF-GP significantly outperforms competing methods and random feasible sampling. The Alkylation benchmark determines optimal operating conditions for simplified alkylation process considering complicated nonlinear constraints representing economic, physical and performance limits. Again, LEAF-GP outperforms other methods, which often struggle find feasible solutions. More details regarding the presented benchmark problems are given in the appendix.
5.4 Mixed-variable spaces
We now consider spaces that exhibit mixtures of continuous, integer and categorical variables for constrained problems. Tree model-basee algorithms naturally support categorical variables by replacing continuous splits with categorical splits. For all methods that do not support categorical features, we use one-hot encoding. The category with the highest corresponding auxiliary variable value is chosen for the subsequent query point. On the Pressure Vessel (4D) 12] benchmark, which comprises two continuous features, two integer features and three inequality constraints, Fig. ((a)
shows that LEAF-GP outperforms other black-box optimization methods.
The CIFAR-NAS (29D) a high-dimensional benchmark problem with one continuous, 23 integer and five categorical variables. It describes properties of different layers and training hyperparameters for convolutional neural network (CNN) trained on CIFAR-10 40 The problem's 's search space is hierarchical, as different layers, ie.. convolutional and fully fully-connected layers, can be activated.

---

0.650
Black-Bon Objective
-50.0K
Blac-Boxx Objective
0.600
100.0K
0.575
EIHATERN
WATERN
0.550
REASAND
-150.KK
LEAF
500PT-GBR
0.525
SHOPT GORT
LEAF
-200.KK
SWAC
UCB-MATERN
MON-MATERN
20
40
60
80
100
25
50
75 Iterations
100
125
150
Iterations
(a) Pressure Vessel (4D, 3IC)
(b) CIFAR -ASS (29D)
Figure 5: Feasible black-box optimization comparing LEAF-GP vs. baseline. Plot shows median line and confidence intervals (first and third quartile) from 20 random seeds. Figure subtitles give the number of dimensions (D) and inequality constraints (IC). Section 5.4 provides details.
Properties describing layer are only relevant if the layer is active. For this problem, we introduce two types of constraints guide the LEAF-GP: (i) constraints that allow for feasible network designs, i.e. different values for padding, stride, kernel size and max pooling in earlier layers affects what kernel size is feasible in following layers and (ii) constraints that ensure that hierarchies are respected. Towards the former, we include constraints to ensure the output size of every convolutional layer is positive. general, tree models are particularly good capturing hierarchical search space structures. To enforce hierarchical relationships during the acquisition optimization step, we introduce indicator constraints that force layer properties to take predefined default values if the layer inactive. This ensures that the optimizer avoids leaves that infer splitting conditions from feature properties of inactive layers, thus effectively reducing the search space. To expedite tests, we train the CNN on half of the CIFAR-10 training data and optimize for test accuracy. Fig. 5(b) summarizes the results of this study and shows that tree model-based algorithms generally outperform -ATTR and EI EMMMATEN. Utilizing search space constraints, SMAC and LEAF-GP find the neural architectures with the highest test accuracy with SMAC slightly outperforming LEAF-GP. Since finding feasible architectures for CIFAR-10 10 is not particularly challenging. Fig. in the appendíx shows benchmark for tuning variational autoencoders (VAE) adapted from Daxberger et al. [16], where LEAF significantly outperforms other algorithms. Finding feasible VAE architectures is more difficult given the requirement that the latent encoding must be decoded back to original input image size.
Conclusion
We present framework for black-box optimization based tree kernel Gaussian processes that simultaneously allows (i) reliable uncertainty quantification mixed feature spaces and (ii) incorpo- ration of explicit input constraints. Although these two needs have been considered separately, are able to address both simultancously through the sixeedinttger second-orde cone formulation of the acquisition function. The numerical studics show that the proposed strategy performs competi- tively with state-of-the-arr algorithms for unconstrained problems and may significantly outperform existing methods for constrained benchmarks, especially those with mixed feature spaces. We optimization constraints together with the acquisition functions to incorporate domain knowledge and leverage hierarchical search space structures, e.g., for neural architecture search.
7 Acknowledgements
This work was supported by BASF SE, Ludwigshafen am Rhein, EPSRC Research Fellowships to R.M. (EP/PO16871/I) and CT (EP/T001577/I), and an Imperial College Research Fellowship to CT.
10

---

References
[1] Alizadeh and D. Goldfarb. Second-order cone programming. Mathematical programming, 95(1):3-51, 2003.
[2] M. Balandat, B. Karrer, D. Jiang. S. Daulton, B. Letham, A. G. Wilson, and E. Bakshy. Botorch: a framework for efficient Monte-Carlo Bayesian optimization. Advances in Neural Information Processing Systems, 33:21524-21555 2020.
[3] R. Baptista and M. Poloczek. Bayesian optimization of combinatorial structures. In International Conference on Machine Learning, pages 462-471. PMLR, 2018.
[4] H. Y. Benson and U. Saglam. Mixed-integer second-order cone programming: A survey, In Theory Driven by Inftuential Applications, pages 13-36. INFORMS, 2013.
[5] Bergstra, R. Bardenet, Y. Bengio, and B. Kégl. Algorithms for hyper-parameter optimization. Advances in Neural Information Processing Systems, 24, 2011.
[6] D. E. Bernal, Q. Chen, F. Gong, and I. E. Grossmann. Mixed-integer nonlinear decomposition toolbox for Pyomo (MindtPy). In Computer Aided Chemical Engineering, volume 44, pages 895-900. Elsevier, 2018.
[7] J. Blank and K. Deb. pymoo: Multi-objective optimization in python. IEEE Access, 8:89497- 89509, 2020.
[8] P. Bonami, L. T. Biegler, A. R. Conn, G. Coruéjols, I. E. Grossmann, C. D. Laird, J. Lee, A. Lodi, F. Margot, N. Sawaya, and A. Wächter. An algorithmic framework for convex mixed integer nonlinear programs. Discrete Optimizationn, 5(2): 1022:18620444 204, 2008.
[9] Boukouvala and M. G. Ierapetritou. Derivative-free optimization for expensive constrained problems using novel expected improvement objective function. AIChE Journal, 60(7): 2462-2474, 2014.
[10] P. Buathong, D. Ginsbourger, and T. Krityakierne. Kernels over sets of finite sets using RKHS embeddings, with application Bayesian (combinatorial) optimization. International Conference on Artificial Intelligence and Statistics, pages 2731 PMLR, 2020.
[11] F. Ceccon, J. Jalving, J. Haddad, A. Thebelt, C. Tsay, C. D. Laird, and R. Misener. OMLT: Optimization machine learning toolkit. arXiv, 2202.02414, 2022.
[12] C. A. Coello and M. Montes. Constraint-handling in genetic algorithms through the use of dominance-based tournament selection. Advanced Engineering Informatics, 16(3):193-203 2002.
[13] C. Coey, M. Lubin, and J. P. Vielma. Outer approximation with conic certificates for mixed- integer convex problems. Mathematical Programming Computation, 12(2):249-295, 2020.
[14] D. D. Cox and S. John. A statistical method for global optimization. In Proceedings 1992 IEEE International Conference on Systems, Man, and Cybernetics, pages 1241-1246. IEEE, 1992.
[15] A. Davies and 7. Ghahramani. The random forest kernel and other kernels for big data from random partitions. arXivv, 1402.4293, 2014.
[16] E. Daxberger, A. Makarova, M. Turchetta, and Krause. Mixed-variable Bayesian optimization. In Proceedings of the Twenty-Ninth International Conference on International Joint Conferences on Artificial Intelligence, pages
2021.
[17] L. Deng. The MNIST database of handwritten digit images for machine learning research. IEEE Signal Processing Magazine, 29(6):141-142, 2012.
[18] A. Deshwal, S. Belakaria, and J. R. Doppa. Scalable combinatorial Bayesian optimization with tractable statistical models. arXiv, 2008.08177 2020.
[19] A. Deshwal, S. Belakaria, and J. R. Doppa. Bayesian optimization over hybrid spaces. In
International Conference Machine Learning, pages 2632-2643. PMLR, 2021.

---

[20] S. Drewes and S. Ulbrich. Mixed integer second order cone programming. Verlag Dr. Hut Germany, 2009.
[21] Elsken, H. Metzen, and Hutter. Neural architecture search: survey. The Journal of Machine Learning Research, 20(1):1997-2017, 2019.
[22] P. I. Frazier. A tutorial on Bayesian optimization. arXiv, 1807.02811 2018.
[23] É. Fromont, H. Blockeel, and J. Struyf. Integrating decision tree learning into inductive databases. In International Workshop on Knowledge Discovery in Inductive Databases, pages 81-96. Springer, 2006.
[24] R. Gardner, G. Pleiss, Bindel, Q. Weinberger, and A. G. Wilson. GPyTorch: Blackbox matrix-matrix Gaussian process inference with GPU acceleration. In Advances in Neural Information Processing Systems, 2018.
[25] E. C. Garrido-Merchán and D. Hernández-Lobato. Dealing with categorical and integer- variables in Bayesian optimization with Gaussian processes. Neurocompuling, 380:20-35, 2020.
[26] S. Gopakumar, S. Gupta, S. Rana, Nguyen, and S. Venkatesh. Algorithmic assurance: An active approach to algorithmic testing using Bayesian optimisation. Advances in Neural Information Processing Systems, 31, 2018.
[27] S. Greenhill, S. Rana, S. Gupta, P. Vellanki, and S. Venkatesh. Bayesian optimization for adaptive experimental design: A review. IEEE access, 8:13937-13948, 2020.
[28] Gurobi Optimization, LLC. Gurobi Optimizer Reference Manual, 2022.
[29] E. Han, I. Arora, and J. Scarlett. High-dimensional Bayesian optimization via tree-structured additive models. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 35, pages 7630-7638. 2021.
[30] F. Häse, M. Aldeghi, R. Hickman, L. M. Roch, and A. Aspuru-Guzik. Gryffin: An algorithm for Bayesian optimization of categorical variables informed by expert knowledge. Applied Physics Reviews, 8(3):031406, 2021.
[31] T. Head, MechCoder, G. Louppe, I. Shcherbatyi, fcharras, Z. Vinícius, emmalone, C. Schröder, nel215, and N. C. et al. sckiiroptiiizehhcckicccppiiie v0.5.2. Zenodo, Mar. 2018. doi: 10.5281/2emonl 1207017
1207017.
[32] A. Hedar. Test function web pages. Retrieved May 18, 2022, from httmm//tt::///wwttiiiiiii amp www.sfffeeesllll htm
i. ..kotoou.uaa httpa///wweebbeessttttttttaaaaaa
[33] H. Hijazi, P. Bonami, and A. Ouorou. An ape for separable mixed-integer nonlinear programs. INFORMS Journal on Computing, 26(1):31 2014.
[34] Hutter, H. H. Hoos, and K. Leyton-Brown. Sequential Model-Based Optimization for General Algorithm Configuration. In Proceedings of the 5th International Conference on Learning and Intelligent Optimization, pages 403. Springer- Verlag, 2011.
[35] R. Jenatton, C. Archambeau, González, and M. Seeger. Bayesian optimization with tree- structured dependencies. In International Conference on Machine Learning, pages 1655-1664. PMLR, 2017.
[36] Y. Jin and Branke. Evolutionary optimization uncertain environments-a survey. IEEE Transactions on Evolutionary Computation, 93)):003-317 2005.
[37] G. Ke, Q Meng, T. Finley, T. Wang, W. Chen, W. Ma, Q. Ye, and T. Liu. LightGBM: A Highly Efficient Gradient Boosting Decision Tree. In Proceedings of the 31st International Conference on Neural Information Processing Systems, pages 31 399-3157. Curran Associates Inc., 2017.
[38] D. P. Kingma and J. Ba. Adam: A method for stochastic optimization. arXiv, 1412.6988 2014.
[39] D. P. Kingma and M. Welling. Auto-encoding variational Bayes. arXiv, 1312 2013.
12

---

[40] A. Krizhevsky, V. Nair, and G. Hinton. CIFAR-I0 (Canadian Institute for Advanced Research). URL BRLLhtp::/ http://wwwr cs. toronto. wtt.oto...e.uutrrii//////// cifar .html
[41] Y.J. Kuo and H. D. Mittelmann. Interior point methods for second-order cone programming and or applications. Computational Optimization and Applications, 28(3):255-285, 2004.
[42] D. Lee, H. Park, and C. D. Yoo. Face alignment using cascade Gaussian process regression trees. In Proceedings of the IEEE Conference Computer Vision and Pattern Recognition, pages 2204-4212, 2015.
[43] J.-C. Lévesque, A. Durand, C. Gagné, and R. Sabourin. Bayesian optimization for conditional hyperparameter spaces. In 2017 International Joint Conference on Neural Networks (IJCNN), pages 286-293. IEEE, 2017.
[44] M. Lindauer, K. Eggensperger, M. Feurer, A. Biedenkapp. D. Deng, C. Benjamins, T. Ruhkopf, R. Sass, and F. Hutter. SMAC3: A Versatile Bayesian Optimization Package for Hyperparameter Optimization. The Journal of Machine Learning Research. 23:54-1, 2022.
[45] M. S. Lobo, L. Vandenberghe, S. Boyd, and H. Lebret. Applications of second-order cone programming. Linear algebra and its applications, 284(1-3):193-228, 1998.
[46] M. Lubin, E. Yamangil, R. Bent, and P. Vielma. Extended formulations in mixed-integer convex programming. In International Conference on Integer Programming and Combinatorial Optimization, pages 102-113. Springer, 2016.
[47] A. Lundell and J. Kronqvist. Polyhedral approximation strategies for nonconvex mixed-integer nonlinear programming in SHOT. Journal of Global Optimization, 82(4):863- 896, 2022.
[48] X. Ma and M. Blaschko, Additive tree-structured covariance function for conditional parameter spaces in Bayesian optimization. In International Conference on Artificial Intelligence and Statistics, pages 1015-1025. PMLR, 2020.
[49] V. V. Mišić. Optimization of tree ensembles. Operations Research, 68(5):1605-1624, 2020.
[50] M. Mistry, D. Letsios, G. Krennrich, R. M. Lee, and Misener, Mixed-integer convex nonlinear optimization with cradient-booster trees embedded. INFORMS Journal on Computing, 33(3): 1003-1119, 2021.
[51] G. Nanfack, P. Temple, and B. Frénay, Constraint enforcement on decision trees; A survey, ACM Computing Surveys (CSUR), 2022.
[52] D. Nguyen, S. Gupta, S. Rana, A. Shilton, and S. Venkatesh. Bayesian optimization for categorical and -spegory-specific continuous inputs. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 34, pages 2266, 2020.
[53] S. Nijssen and E. Fromont. Mining optimal decision trees from itemset lattices. In Proceedings of the 13th ACM SIGKDD international conference on Knowledge discovery and data mining, pages 530-539, 2007.
[54] T. P. Papalexopoulos, C. Tjandraatmadja, R. Anderson, J. Vielma, and D. Belanger. Con- strained discrete black-box optimization using mixed-integer programming. In International Conference on Machine Learning, pages 17295-17322. PMLR, 2022.
[55] A. Paszke, S. Gross, F. Massa, A. Lerer, J. Bradbury, G. Chanan, T. Killeen, Z. Lin, N. Gimelshein, and L. e. a. Antiga. Pytorch: An imperative style, high-performance deep learning library. In Advances in Neural Information Processing Systems, volume 32, pages 8024-8035. Curran Associates, Inc., 2019.
[56] A. Rahimi and B. Recht. Random features for large-scale kernel machines. Advances in Neural Information Processing Systems, 20, 2007.
[57] C. Rasmussen and C. K. I. Williams. Gaussian Processes for Machine Learning. The MIT Press, 2006.
13

---

[58] S. R. Rath. Convolutional
Variational Autoencoder in PyTorch on MNIST May 18. 2022, from https://aeuggercaafe com/
Dataset.
Retrieved
autoencoder- www.rder------poooro-ooo-------- dataset/
[59] B. Ru, A. Alvi, V. Nguyen, M. A. Osboree, and S. Roberts, Bayesian optimisation over multiple continuous and categorical inputs. In International Conference on Machine Learning, pages PMLR, 2020.
[60] R. Sauer, A. Colville, and C. Burwick. Computer points way more profits. Hydrocarbon Processing, 84(2), 1964.
[61] B. Shahriari, K. Swersky, Z. Wang, R. P. Adams, and N. de Freitas. Taking the human out of the loop: A review of Bayesian optimization. Proceedings of the IEEE, 104(1):148-175, 2016.
[62] solab ntu. opt-prob-coollectoon Retrieved May 18, 2022, from http:/ http:/www.gttttt. com/ solab-ntu bigabbattt/ppp-pro---olllleeeeed
[63] S. Surjanovic and D. Bingham. Virtual library of simulation experiments: Test functions and datasets. Retrieved May 18, 2022, from http::/ww.. sfu. ca/ "ssurjand
[64] D. Thain, T. Tannenbaum, and M. Livny. Distributed computing in practice: the Condor experience. Concurrency- Practice and Experience, 17(2-4):323-356 2005.
[65] A. Thebelt, Kronqvist, R. M. Lee, N. Sudermann-Merx, and Misener. Global optimization with ensemble machine learning models. In Computer Aided Chemical Engineering, volume 48, pages 1981-1986. Elsevier, 2020.
[66] A. Thebelt, J. Kronqvist, M. Mistry, R. M. Lee, Sudermann- MMrxx, and R. Misener. ENT- MOOT: A framework for optimization over ensemble tree models. Computers & Chemical Engineering, 151:107343, 2021.
[67] A. Thebelt, C. Tsay, R. M. Lee, N. Sudermann-Mers D. Walz, T. Tranter, and R. Misener. Multi-objective constrained optimization for energy applications via tree ensembles. Applied Energy, 306: 18061, 2022.
[68] A. Thebelt, J. Wiebe, J. Kronqvist, C. Tsay, and R. Misener. Maximizing information from chemical engineering data sets: Applications to machine learning. Chemical Engineering Science, 152:1 2022.
[69] A. Tran, J. Sun, J. M. Furlan, K. V. Pagalthivarthi. R. J. Visintainer, and Y. Wang. pBO- 2GP-3B: A batch parallel known/unknown constrained Bayesian optimization with feasibility classification and its applications in computational fluid dynamics. Computer Methods in Applied Mechanics and Engineering. 347:827-852 2019.
[70] M. Trencseni. Solving CIFAR-10 with Pytorch and SKL. Retrieved May 18, 2022, from https://bytepane wttp:/ppwwa.ccm.ooooooaeam wwwm:oolvngg-ciee.eeeao cifar 10- www.httthhpptthh.ooooooomoo and- sk1 html
[71] J. P.. Vielma, I. Dunning, J. Huchette, and M. Lubin. Extended formulations in mixed integer conic quadratic programming. Mathematical Programming Computation, 933):369-418, 2017.
[72] C. Ying, A. Klein, E. Christiansen, E. Real, K. Murphy, and F. Hutter. Nas-bench- 101: Towards reproducible neural architecture search. In International Conference on Machine Learning, pages 7105-7114 PMLR, 2019.
[73] A. Zafari, R. Zurita- Milla, and E. Izquierdo-Vertee Verdiguier. Exaluating the performance of random forest kernel for land cover classification. Remote sensing, 11(5):575, 2019.
14

---

Checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's
contributions and scope? [Yes] See Abstract and Section !
(b) Did you descrihe the limitations of your work? [Yes] See Section
(c) Did you discuss any potential negative societal impacts of your work? [No] There is no obvious societal impact of the method proposed.
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]
2. If you are including theoretical results...
(a) Did you state the full set of assumptions of all theoretical results? [N/A]
(b) Did you include complete proofs of all theoretical results? [N/A]
3. If you ran experiments..
(a) Did you include the code, data, and instructions needed to reproduce the main experi- mental results (either in the supplemental material or URL)? [Yes] We submit the code to run all experiments for the proposed method.
(b) Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Yes] See end of Section 4 and supplemental material.
(c) Did you report error bars (e-g., with respect the random seed after running experi- ments multiple times)? [Yes] We run 20 random seeds for every experiment and report median and confidence intervals, i.e., first and third quartile. See Section
(d) Did you include the total amount of compute and the type of resources used (e.g., type of GPUs, internal cluster, cloud provider)? [Yes] We report the computational resources used to conduct the experiments in the supplementary material.
4. If you are using existing assets (e.g., code, data, models) or curating/releasing new assets...
(a) your work uses existing assets, did you cite the creators? [Yes] Yes, cite all libraries that were used when preparing the numerical studies. For more information, we refer to the supplemental material.
(b) Did you mention the license of the assets? [No] License information for all assets can be found in the documentation of each of the assets.
(c) Did you include any new assets either in the supplemental material or URL? [No] The code to reproduce the paper's main results will be published open-source after the peer review process.
(d) Did you discuss whether and how consent was obtained from people whose data you're using/curating? [N/A]
(e) Did you discuss whether the data you using/curating contains personally identifiable information or offensive content? [N/A]
5. If you used crowdsourcing or conducted research with human subjects...
(a) Did you include the full text of instructions given participants and screenshots, if applicable? [N/A]
(b) Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, applicable? [N/A]
(c) Did you include the estimated hourly wage paid to participants and the total amount spent participant compensation? [N/A]
15


Reviews:
[{"review_id": "phWwRBnvuir", "reviewer": "Reviewer_94qs", "summary": "The paper develops a new framework for Bayesian optimization with mixed-type features and flexible feature constraints. The idea is to utilize a GP model with a kernel function based on tree ensembles that can naturally handle mixed-type features. A mixed-integer optimization technique is used to minimize the piecewise-constant LCB acquisition function. The proposed method shows promising performance in various benchmarks.", "questions": "Major comments:\n\n1. As far as I understand, there are two possible ways to train the proposed surrogate model with a tree ensemble kernel: (i) treat the tree structure as a kernel hyperparameter and learn it by maximizing marginal likelihood (which is an end-to-end approach and provides a coherent statistical model), and (ii) fit the tree structure first by gradient boosting, and then fit a GP model with the given tree structure (which does not seem to be a well-defined statistical model to me). It is unclear to me which approach is used. If (ii) is used, what are the responses used to train boosting trees?\n\n2. The experiment in Section 5.1 aims to verify the uncertainty metrics. However, in Eq. (12a) only the posterior mean is included while the posterior variance (which accounts for uncertainty) is missing. Moreover, the plots in Figure 2 only show the error for the mean prediction, which again does not convey any information on uncertainty. It is unclear to me how the conclusion in Section 5.1 is reached and why the constraint in Eq. (12c) is needed.\n\n3. Since the posterior mean and variance from the tree-kernel GP are piecewise constant functions corresponding to the tree partitions, I wonder how sensitive the performance is with respect to the tree ensemble hyperparameters such as the number of boosting iterations.\n\n4. How is the proposed method compared to its competitors in terms of computation time? I would expect the sparsity in the tree-based covariance matrix can reduce computation for GP.\n\nMinor comments:\n\n* I think the notations in Section 3.3 are quite confusing and should be clearly defined. For example, what do $v_{i, j}$, $\\nu_{i, j}$ and $\\mathrm{splits}(t)$ mean?\n* In Eq. (5a), I understand that the optimization is done over $\\mathbf{x}$. But since $\\mathbf{z}$ is a function of $\\mathbf{x}$ given the tree structures (as in Section 3.1), it should not be included under $\\mathrm{argmin}$. It is also not clear what $\\nu$ is.\n* Line 98: Leaf subspace $\\mathbf{x}_l \\in \\mathbb{R}^n$ -> $\\mathbf{x}_l \\subset \\mathbb{R}^n$.\n", "limitations": "Please see the major comments above.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strength:\n\n* A novel framework for Bayesian optimization to handle mixed-type features and flexible feature constraints.\n* Promising performance in various benchmarks.\n\nWeakness (see major comments below for details):\n\n* Training of the surrogate model is unclear.\n* The experiment in Section 5.1 does not seem solid.\n* The paper could benefit from additional sensitivity and timing results.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "J6aiqAxoGr", "reviewer": "Reviewer_NYXt", "summary": "This paper focuses on Bayesian Optimization and utilizes Gaussian Processes (GPs) with ensemble tree kernels. The proposed approach advocates a mixed-integer second-order cone optimization formulation, relying on the lower confidence bound (LCB) acquisition criterion, that allows for uncertainty quantification and readily accommodates the incorporation of input constraints. Numerical tests show the benefits of the proposed framework compared to existing methods in the literature. ", "questions": "Some minor comments are listed below:\n- In Eq. 6b, instead of $\\sigma^2$ it should be $\\sigma^2(\\mathbf{x})$.\n- In line 51, instead of “prior work” it should be “prior works”.\n- In line 118, instead of “where data is sparse” it should be “where data are sparse”.\n- In line 167, instead of “We compute constant matrix”, it should be “We compute the constant matrix”. \n- In line 199, in the phrase “we found that Gurobi 9 often solves” a citation is needed.\n", "limitations": "The limitations subsection is well written. Some additional concerns of mine are mentioned in the “Strengths and Weaknesses” section of the review. \n", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "STRENGTHS \n\n- The ensemble tree kernel GP - based approach is intuitive and seems to work well in practice. \n- The paper is well written and the ideas of the paper are well presented.\n\nWEAKNESSES \n\n- The novelty of the paper is somewhat limited. The ensemble tree kernel is known and the same holds true for Eqs. 3 and 4. In my mind, the advocated approach is just a combination of some intuitive well-known components. \n-  Additional experiments are needed to corroborate the benefits of the proposed approach compared to existing methods in the literature. Except the GP-UCB approach which is indeed a good baseline, there exist many other well-known acquisition criteria combined with GP surrogate model such as Expected Improvement (EI), Thompson Sampling (TS), Entropy Search (ES) and Predictive Entropy Search (PES) to list a few. So, in my opinion it is important to compare with some of these approaches to better demonstrate the benefits of the advocated method, and especially with the EI which seems to be very effective in several BO problems. In addition, there exist some very interesting ensemble approaches that can effectively handle continuous-categorical data; see e.g [1], [2]. A numerical comparison (or even discussion) with these approaches would also be nice. \n\n[1] Gopakumar, Shivapratap, Sunil Gupta, Santu Rana, Vu Nguyen, and Svetha Venkatesh. \"Algorithmic assurance: An active approach to algorithmic testing using bayesian optimisation.\" Advances in Neural Information Processing Systems 31 (2018).\n \n[2] Nguyen, Dang, Sunil Gupta, Santu Rana, Alistair Shilton, and Svetha Venkatesh. \"Bayesian optimization for categorical and category-specific continuous inputs.\" In Proceedings of the AAAI Conference on Artificial Intelligence, vol. 34, no. 04, pp. 5256-5263  (2020).\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "ARDHj4lTfhO", "reviewer": "Reviewer_PeTE", "summary": "This paper proposes a feasible way to perform Bayesian optimization using ensembles of decision trees. Although it is easy to use decision trees ensembles to predict a mean and variance (and therefore formulate acquisition functions like upper confidence bounds), _optimizing_ these acquisition functions is challenging because they are piecewise constant. The key insight of this paper is that if the decision tree predictions are used to form a GP kernel, upper confidence bound acquisition function can be formulated as a conic problem and be inputted into a generic solver to find the global optimum. A secondary advantage of this approach is that the solver can also tolerate a large number of constraints from a variety of families, which is practically useful.\n\nThe authors provide the technical details of this formulation and show that it works well for a variety of synthetic optimization problem and one real-world problem over mixed-variable, hierarchical spaces.", "questions": "- Trees are effectively parametric methods: once trained on data, updating the trees generally requires retraining. This is in contrast to the GP with the tree kernel, which will make different predictions conditioned on new data even if the underlying trees are not retrained. _In your baseline experiments, do you keep the trees fixed or do you retrain them?_ If you are not retraining them, this could explain the roughly constant performance of the baselines in e.g. Figure 4. I think that you should be retraining them.\n- _Do all tree-based methods use the same underlying trees?_ To minimize external sources of variation, I think that the trees in SMAC/skopt should be the exact same trees used in the GP kernel. If this is not possible/extremely hard to do then don't worry about it.\n- _What things are randomized in the random seeds?_ In general I think the main 2 sources of randomness are 1) randomness in the BO procedure (e.g. initialization of interior point methods) and 2) randomness in training of the tree ensembles. Do your error bars in Figures 2-4 cover both of these sources of variation, or just the randomness of BO?\n- _How necessary is the gurobi solver specifically?_ Could you use another similar solver? Would other solvers be able to handle the same constraints?\n- The solver seems to be able to handle a large family of constraints. _Is the solution exact/nearly exact for all these constraint types, or does the solver lose its guarantees for certain constraint families (e.g. non-convex constraints)?_ (I'm not very familiar with cone solvers)", "limitations": "As far as I can tell, the limitations of this work are discussed adequately.", "rating": 8, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Overall I thought that this paper was very good.\n\n### Strengths\n\n- **Very relevant problem choice**: trees are one of the few models that can easily handle mixed-variable/hierarchical data, _and_ work well on small datasets. To my knowledge, optimizing the acquisition function is a significant hurdle preventing the widespread use of tree ensembles in Bayesopt.\n- **Elegant solution**: formulating the acquisition function as a cone problem allows for decades of research on convex solvers to be applied. This is a clever insight. It is nice to see a (semi-)exact solution to an optimization problem presented at NeurIPS, where most papers tend to use approximations which are hard to verify.\n- **Extremely well-written**: this paper is definitely in the top 10% of writing quality. Everything was extremely clear. I especially appreciated the discussion of the related work and the technical background (section 3). The technical details in section 4 were very well presented. I don't know if I would 100% be able to implement the method based on the paper (the explanation of $\\nu$ and $v$ in equations 3-4 could be slightly clearer), but it is closer to this ideal than perhaps 90% of papers at NeurIPS. Also, unlike many ML conference papers, it did not feel like the authors crammed a 12 page paper into 9 pages by moving important explanations/equations into the appendix, which I really appreciated.\n- **Highly Usable**: this method seems like it should work \"out of the box\" (assuming access to a conic solver) because it is simple and could be plugged into any decision tree package. I highlight this because I feel like most methods presented at generally do not work \"out of the box\" because they are brittle, require a lot of tuning, or have undisclosed limitations.\n\n### Weaknesses\n\n- **Several technical details unclear**: see \"Questions\" section. These are not absolutely critical to the paper but I would appreciate clarification.\n- **Unclear effect of \"different objective\" vs \"exact solution\"**: the proposed method has 2 deviations from using the empirical mean/variance of a tree ensemble for Bayesopt: 1) making predictions with a tree kernel GP instead of empirical mean/variance 2) optimizing UCB exactly using a cone solver. In section 5 it was generally unclear to me how much of the difference in final performance came from using the trees to make predictions in a different way (GP instead of empirical mean/variance), and how much came from performing [nearly] exact maximization of the acquisition function. A good experiment to do would be use the GP formulation to find a promising set of candidate points (perhaps using multiple seeds), then use the original empirical mean/variance UCB objective to select a point to query in the end. This would be combining \"exact optimization\" with the standard way of making predictions.\n- **Code requires a non-public library** gurobi. This is not a deal-breaker, and I realize there are academic licenses available, but it would be really nice to have a fully open-source version for the camera-ready paper. I ask some additional questions about the choice of solver below.\n\n### Suggestions\n\n- Given that the paper focuses specifically on non-continuous spaces, why is $x$ defined to be a real number when introducing Bayesopt on line 44? I think $x$ should just be an element of a general space $\\mathcal{X}$\n- Many readers at NeurIPS may not be very familiar with cone problems / conic solvers. It would be nice to put some information on these in the background, even if it is just a short section saying that \"it is a very general way of formulating certain optimization problems, with a variety of established methods for solving them\"\n- \"LCB\" is a fairly uncommon formulation of Bayesopt; I think many readers are accustomed to seeing problems posed as UCB _maximization_. Perhaps it would be slightly helpful to the reader to frame the problem as maximization and multiple everything in the paper by $-1$?\n- I saw in the appendix that you use Adam to optimize GP kernel hyperparameters. One problem with Adam is that it may converge very slowly if the learning rate is too low, or not fully converge if it is too large. Especially since you are only optimizing 2 hyperparameters, I recommend using L-BFGS instead. In my personal experience this method is much more robust.\n- The objective in equation 12 would also be another interesting baseline to see: it would show how much the GP variance actually contributes, as opposed to just finding a point with a sufficient amount of disagreement among ensemble members.\n\n\n### Overall\n\nOriginality: seems high, but I am not very familiar with related work in this area\n\nQuality: high, good paper, good method, good evaluation\n\nClarity: very high!\n\nSignificance: high. The proposed method solves a real problem and is likely to be used in practice, especially if a good open-source library is published.\n\nI definitely recommend accepting this paper and my score reflects that. I would consider increasing my score slightly if my suggestions are addressed. I could also imagine decreasing my score if:\n- The novelty is not as high as I thought (I am not very familiar with related work in this area, so perhaps there is a paper that does something very similar which was not cited)\n- There were technical flaws that I did not catch (I know only the basics of convex optimization/cone problems, so this is possible)", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}]

Author rebuttal:

Decision: Accept