Practical Bayesian Optimization of Objectives with Conditioning Variables
Year: 2021

Abstract:
Bayesian optimization is a class of data efficient model based algorithms typically 1 focused on global optimization. We consider the more general case where a user is 2 faced with multiple tasks that each need to be optimized, find the global optima 3 within each task, all the task conditional optima. For example given a range of 4 cities with different patient distributions, we optimize the ambulance locations for 5 each and every city; given subclass partitions of CIFAR-10, we optimize CNN 6 hyperparameters for each partition. Similarity across tasks boosts optimization of 7 each task in two ways: in modelling by data sharing across objectives, and also in 8 acquisition by quantifying how a single point on one task can help learn the optima 9 of similar tasks. For this we propose a framework for conditional optimization: 10 ConBO. This can be built on top of a range of acquisition functions and we propose 11 a new Hybrid Knowledge Gradient acquisition function. The resulting method is 12 intuitive and theoretically grounded, either matches or significantly outperforms 13 recently published works on a range of problems, and thanks to the unique nature 14 of conditional optimization, is easily parallelized to collect a batch of points. 15

Bayesian optimization is a class of data efficient model based algorithms typically1 focused on global optimization. We consider the more general case where a user is2 faced with multiple tasks that each need to be optimized, find the global optima3 within each task, all the task conditional optima. For example given a range of4 cities with different patient distributions, we optimize the ambulance locations for5 each and every city; given subclass partitions of CIFAR-10, we optimize CNN6 hyperparameters for each partition. Similarity across tasks boosts optimization of7 each task in two ways: in modelling by data sharing across objectives, and also in8 acquisition by quantifying how a single point on one task can help learn the optima9 of similar tasks. For this we propose a framework for conditional optimization:10 ConBO. This can be built on top of a range of acquisition functions and we propose11 a new Hybrid Knowledge Gradient acquisition function. The resulting method is12 intuitive and theoretically grounded, either matches or significantly outperforms13 recently published works on a range of problems, and thanks to the unique nature14 of conditional optimization, is easily parallelized to collect a batch of points.15
1 Introduction16
Expensive stochastic black box functions arise in many fields such as fluid simulations [1], engineering17 wing design [2], and machine learning parameter tuning [3]. Bayesian optimization is a powerful set18 of tools to optimize such functions, finding the input with highest long term average performance19 x∗ = arg maxx E[f(x)], (the expectation represents averaging over the performance noise). In this20 work we consider an under-explored generalization of the standard setting previously referred to as21 “conditional optimization" [4] where a user has a collection of functions, or tasks, and simply seeks22 the peak of each function/task. Formally, we have an expensive black box f that takes as arguments23 both a task s ∈ S and an input x ∈ X , typically box-constrained continuous variables, and returns a24 noisy scalar performance25 f(s, x) : S ×X → R. (1) At each iteration an algorithm determines both task and input (s, x) then observes performance26 y = f(s, x) and the goal is to learn the input with highest average performance for each task27
x∗(s) = arg max x E[f(s, x)]. (2)
x∗(s) is referred to as the optima conditioned on s, see Figure 1 black line. In certain applications,28 one may want to give higher priority to particular tasks hence a task weighting function W (s) may29 also be specified. In this work we consider the following applications.30
CNN hyperparameters: the CIFAR-10 dataset contains 10 classes, this is split into five mutually31 exclusive binary classification datasets and for each we train a CNN, each CNN is a task in S =32
Submitted to 35th Conference on Neural Information Processing Systems (NeurIPS 2021). Do not distribute.
{1, ..., 5}. For each CNN, we optimize dropout rates, batch size and Adam parameters, so X ⊂ R733 and f(s, x) is the validation accuracy. We assume all five CNNs have equal priority, W (s) = 1/5.34
1.0
0.5
0.0
0.5
1.0
1.5
2.0
2.5
3.0
Ac tio
n x
Negative Rosenbrock Function
data maxx n(s, x) (s, x)n+1 states
2.0 1.5 1.0 0.5 0.0 0.5 1.0 1.5 2.0 State s
0
2
4
St at
e Va
lu e
tasks
T sk
Figure 1: Top: GP model. Each task (vertical slice) is an objective function to maximise. A new sample (s, x)n+1 provides information about the optima of similar tasks. Bottom: the acquisition function value of each task using hybrid Knowledge Gradient.
Ambulances in a square: [5] given a range of35 30km×30km cities, each city is a task with a36 different population centre s ∈ [0, 30]2. For37 a given city, we optimize the 2D location of38 three ambulance bases, x ∈ [0, 30]3×2 ⊂ R6.39 Given a city s and ambulance bases x, a vir-40 tual environment randomly generates patients41 for a simulated day and the average ambulance42 journey time, f(s, x), is returned. Inland cities,43 like Paris or London, with population centres44 in the middle are more common than coastal45 cities, like Singapore or Dubai, with a popula-46 tion centre on a boundary. Thus, W (s) mirrors47 this distribution of city centres, a Gaussian cen-48 tred at s = (15, 15) with 7km standard deviation49 and truncated to [0, 30]2.50
Assemble to order: [6, 7] a company owns51 many stock warehouses and each one faces a52 different level of demand s ∈ [0.5, 1.5]. Equal53 priority is given to all warehouses and so W (s)54 is uniform. At each warehouse, stock levels are55 controlled by setting targets of a control pol-56 icy x = [0, 20]8. Given a demand level s and57 control parameters x, a simulator generates cus-58 tomer orders according to demand, the ware-59
house stock is sold and replenished according to policy x, and profit over a simulated month is60 returned as f(s, x).61
The closely related multi-task Bayesian optimisation or Bayesian Quadrature optimisation methods62 [8, 9, 10, 11] aim to find a single peak that maximises the weighted sum/integral over tasks, x∗ =63 arg maxx ∫ E[f(s, x)]W (s)ds. However, as they are fundamentally designed for a different purpose,64 we observe that they struggle to find the peak of each and every task, see experiments in Supplementary65 Material 6 (SM 6).66
Also closely related are contextual optimization algorithms [12, 13, 14], conditional optimization67 has also been referred to as “offline contextual" [15]. In contextual applications, at iteration n the68 next task (called context) sn is passed from the black box to a contextual algorithm that intelligently69 determines xn. Then the black box returns both performance yn and the next task/context sn+1. One70 could “hotfix" a contextual algorithm for conditional optimization problems by randomly choosing71 sn+1 ∼ P[s] ∝W (s) at each iteration. However, we empirically show that, unsurprisingly, randomly72 choosing the task in each iteration is significantly worse than intelligently choosing the task.73
Multi-fidelity and multi-information source methods [16, 17, 18, 19] assume there is a single known74 constant target task s∗ ∈ S corresponding to the highest fidelity level or most expensive information75 source that must be optimized, x∗ = arg maxx E[f(s∗, x)]. If the cost of evaluation varies across the76 function domain, c(s, x), cheaper regions of the domain can be evaluated improving sample efficiency.77 Again, such methods may be hotfixed for the setting we consider by artificially designating a target78 task, e.g. the highest priority task s∗ = arg maxsW (s). However, we observe that if cost is constant79 across the domain (as is commonly assumed in the BO literature), these methods greedily optimise80 the artificial target task s∗ and blindly neglect all other tasks, see SM4. This desirable behaviour in81 multi-fidelity optimization problems leads to failure in conditional optimization problems.82
To the best of our knowledge, there exist only a small number of works that propose algorithms83 specifically designed for the conditional setting.84
The Surrogate Collaborative Tuning (SCoT) algorithm [20] optimizes a finite set of tasks. It iteratively85 visits each task, s, in a round-robin fashion and determines the input x by expected improvement (EI).86 The authors apply this to optimize the hyperparameters of an ML model for multiple datasets. The87
Profile Expected Improvement (PEI) and Profile Expected Quantile Improvment (PEQI) algorithms88 [21, 4] consider continuous tasks and significantly improve upon SCoT by dynamically determining89 the task in each iteration. The acquisition value of a given point (s, x) is the expected improvement90 of a new output over the best predicted output within the same task.91
The REVI algorithm [22] at each iteration discretizes both task S and input X spaces. Given a new92 hypothetical point, (s, x)n+1, for each task in the discretization, si, the discrete Knowledge Gradient93 over the X discretization quantifies how (s, x)n+1 will benefit task si. Summing the benefits over the94 tasks yields the acquisition value of (s, x)n+1. REVI was designed to account for how all tasks can be95 optimized by each sample. However, this comes with exponential cost, increasing problem dimension96 leads to exponentially increasing discretization size with corresponding exponential space and time97 requirements. On the other hand, sparse discretizations can lead to poor and arbitrary measures98 of acquisition value. Consequently, REVI has only successfully been applied to low-dimensional99 synthetic problems.100
Most recently, the Multi-task Thompson sampling (MTS) method [15] uses a novel kernel with a101 length scale that varies across tasks. This is combined with a Thompson sampling method for col-102 lecting new data and showed performance improvements over REVI. However, Thompson sampling103 also does not account for how one sample provides benefit for optimizing all tasks, a fundamental104 structural property of conditional optimization.105
We propose a method that exploits the basic structure of conditional optimization while also being106 highly scalable, and therefore applicable to more challenging real-world problems. We make the107 following contributions:108
1. A framework for conditional Bayesian optimization with theoretical guarantees: ConBO.109 2. A new, fast global optimization method: Hybrid Knowledge Gradient.110 3. State-of-the-art performance on open source problems including CNNs and simulators.111
2 The Conditional Bayesian Optimization Algorithm112
We first discuss the fitting of the Gaussian process model and the predicted conditional optima. We113 then motivate the acquisition value for a single task and how this is integrated over tasks yielding the114 acquisition function. Because integrating over tasks multiplies the computational burden, we propose115 Hybrid Knowledge Gradient as a solution. See Algorithm 1 in the SM 2 for a high level summary.116
At a stage after having observed n data points, {(si, xi, yi)}ni=1 where yi = f(si, xi), a Gaussian117 process is fit over the joint space S × X to scalar outputs y. Let X̃n = ((s, x)1, ..., (s, x)n) and118 Y n = (y1, ..., yn). A Gaussian process is defined by a prior mean and prior covariance function,119 µ0(s, x), k0((s, x), (s′, x′)) which are chosen for each application, for more information see [23].120 Let the data covariance matrix be K = k0(X̃n, X̃n) ∈ Rn×n, the posterior mean is given by121
µn(s, x) = µ0(s, x) + k0((s, x), X̃n) ( K + σ20I )−1 (Y n − µ0(X̃n))
and the posterior covariance is given by122 kn((s, x),(s′, x′)) = k0((s, x), (s′, x′)) + k0((s, x), X̃n) ( K + σ20I )−1 k0(X̃n, (s′, x′)).
At the end of data collection, in standard BO, often the peak posterior mean is returned as the123 predicted best input. Generalizing to the conditional setting we simply condition on s,124
x∗N (s) = arg max x µN (s, x). (3)
During data collection, a new data point yn+1 at (s, x)n+1 will update the Gaussian process over the125 whole domain S ×X . To construct an acquisition function for conditional optimization, we start126 by looking for standard acquisition functions that account for how the model changes at unsampled127 locations (s′, x′) 6= (s, x)n+1. Specifically, the popular Expected improvement (EI) [24] and upper128 confidence bound (UCB) [25] methods are both functions of the mean and kernel at the sampled point129 only, i.e. of ( µn((s, x)n+1), kn((s, x)n+1, (s, x)n+1) ) , hence would require non-trivial modification130 to be able to account for how a sample affects similar tasks. Methods that do utilise the mean and131 kernel at unsampled points include Entropy search (ES [26, 27] and PES [28]) that measures the132
mutual information between the new output P[yn+1|xn+1] and the (uncertain) location of the peak133 P[x∗|X̃n, Y n], Max-value entropy search (MES) [29] that measures mutual information between134 the new output P[yn+1|xn+1] and the (uncertain) largest output P[max y|X̃n, Y n], and Knowledge135 Gradient (KG) [30] that measures the expected peak of the new posterior mean E[maxµn+1(x)]136 caused by a new yn+1 at xn+1.137
Each single task s defines a single global optimization problem over x ∈ X . Given a proposed138 sample (s, x)n+1 = (sn+1, xn+1), the acquisition benefit for task s may be computed using ES, PES,139 MES or KG. In this work we adopt KG for its Bayesian decision theoretic derivation that extends140 seamlessly to the conditional setting. For KG, the benefit for a given task is the expected increase in141 peak predicted performance within the task. We denote the task-conditioned KG as142
KGc(s; (s, x)n+1) = Eyn+1 [max x′ µn+1(s, x′)|(s, x)n+1]−max x′′ µn(s, x′′).
We discuss numerical evaluation of KGc(·) in Section 2.1. Similar expressions for conditioned143 entropy methods, ESc(·), PESc(·), MESc(·), are derived in SM 7. Integrating over all tasks s yields144 the total acquisition value145 ∫
S
KGc(s; (s, x)n+1)W (s)ds. (4)
For discrete S the integral is replaced by summation. For continuous S, the integral over tasks s146 cannot be computed analytically so we use Monte Carlo with importance sampling. When using147 a kernel that factorises k(s, x, s′, x′) = σ20kS(s, s ′)kX(x, x ′), like squared exponential or Matérn,148 similarity across tasks is encoded in kS(s, s′). This naturally leads to the proposal distribution149 q(s|sn+1) ∝ kS(s, sn+1). In our continuous task experiments, we use the Matérn kernel and a150 Gaussian proposal distribution with mean sn+1 and the task kernel length scales, ls, as standard151 deviations,152 q(s|sn+1) ∼ N (s|sn+1, diag(l2s)). (5) We generate ns = 20 samples SMC = {s1, ..., sns}, finally the acquisition function is153
ConBO(sn+1, xn+1) = ∑
si∈SMC
W (si)
q(si|sn+1) KGc(si; (s, x)n+1).
Figure 1 shows a set of sampled tasks and the KGc(·) for each one. Each KG term directly measures154 increase in predicted performance for one task, e.g. if y values are dollar amounts, ConBO with155 KGc(·) is the sum of dollar increases over all tasks. However, for entropy based methods, ConBO156 becomes a sum of Shannon information units thereby indirectly optimizing the dollar amounts.157 The randomly sampled tasks SMC may be resampled with each call to ConBO(s, x) and gradients158 estimated enabling the optimal (s, x)n+1 to be found with a stochastic gradient ascent optimizer such159 as Adam [31]. Selecting each point according to maximising ConBO is also myopically optimal in a160 value of information framework:161
Theorem 1 Let (s∗, x∗) ∈ arg max ConBO(s, x) be a point chosen for sampling. (s∗, x∗) is also162 the point that maximises the myopic Value of Information, the increase in predicted performance.163
Further, in finite search space, with an infinite sampling budget all points will be sampled infinitely:164
Algorithm 1 Computing ConBO(s, x). For each call with a candidate (task, input) point, similar tasks are sampled. For each sampled task a cheap acquisition function, hybrid KG, is evaluated. The output is the importance weighted average of hybrid KG values. Details and full numerical expressions are given in SM 3. Require: candidate point (s, x)n+1, parameters ns, nz
Sample ns tasks similar to sn+1, SMC ∼ q(s|sn+1) Initialize output value Q← 0 parfor si in SMC do
Initialize X̃∗i ← {} parfor j in 1, ..., nz do zj ← Φ−1 ( (2j − 1)/2nz ) yn+1j computed with (s, x)
n+1 and zj µn+1j (s, x) constructed by rank-1 update with (s, x, yj) n+1 x∗j ← arg maxx µ n+1 j (si, x) with Optimizer()
X̃∗i ← X̃∗i ∪ (si, x∗j ) end parfor αi ← KGd(X̃∗i ) Q← Q+ αiW (si) / q(si|sn+1)
end parfor return average over tasks Q/ns
Theorem 2 Let S and X be finite sets and N the budget to be sequentially allocated by ConBO.165 Let n(s, x,N) be the number of samples allocated to a point (s, x) within budget N . Then for all166 (s, x) ∈ S ×X we have that limN→∞ n(s, x,N) =∞.167
The law of large numbers ensures that the algorithm learns the true expected performance for all168 points. Proofs are given in the SM 1.169
2.1 Hybrid Knowledge Gradient170
By definition, KG is more expensive than EI and UCB. Further, the function KGc(si, (s, x)n+1)171 must be computed once for each sampled task si, the computational cost is therefore ns times the172 global acquisition function equivalent. To alleviate this cost, we propose a novel, efficient algorithm173 for computing KG. In the following section we assume constant s for brevity, reducing to the174 global optimization setting. Given a hypothetical location xn+1, KG quantifies the value of a new175 hypothetical observation yn+1 by the expected increase in the peak of the posterior mean176
KG(xn+1) = Eyn+1 [
max x′
µn+1(x′) ∣∣xn+1]−max
x′′ µn(x′′). (6)
However, maxx′ µn+1(x′) has no explicit formula and approximations are required which we describe177 next. At time n, the next posterior mean is unknown, however, it may be written as µn+1(x) =178 µn(x) + σ̃(x;xn+1)Z where σ̃(x;xn+1) is a deterministic function and the scalar Z ∼ N (0, 1)179 captures the randomness of yn+1, see SM 2. Previously, KG(x) has been computed in two ways.180
KG by discretization [30, 32]: in Equation 6, the maximizations may be performed over a discrete set of d points Xd ⊂ X . Denoting µ = µn(Xd) ∈ Rd and σ̃(xn+1) = σ̃(Xd;xn+1) ∈ Rd, then
KGd(xn+1) = EZ [ max{µ+ σ̃(xn+1)Z} ] −maxµ.
The max{µ + σ̃(xn+1)Z} is a piece-wise linear function of Z and the expectation is analytically181 tractable. The output is a lower bound of the true KG(x). REVI [22] and the MiSo algorithm182 [17] used KGd with 3000 uniformly random distributed points. This method suffers the curse of183 dimensionality as Xd must grow exponentially with dimension. Further, the discretization may likely184 contain many useless points in uneventful regions of X , see Figure 2 centre-left plot.185
Hybrid Knowledge Gradient
186
KG by Monte Carlo [33, 34]: given xn+1, the method samples up to nz = 1000 Gaussian values187 of Z. For each Zj , construct the posterior mean, µn+1j (x), find the maximum with a continuous188 numerical Optimizer() like L-BFGS or CG. Averaging the maxima from all Zj yields189
KGMC(xn+1) = 1
ns ∑ j Optimizer x′ ( µn+1j (x) ) − Optimizer x′′ ( µn(x′′) ) .
The result is an unbiased estimate of true KG(x) and scales better to higher dimensional X , the190 univariate Z is discretized by Monte Carlo samples instead of X . However, for a good estimate, nz191 must be large, many Optimizer() calls are required. See Figure 2 centre right.192
We instead propose a simple mixture of the two approaches above that both scales to higher dimen-193 sional X and drastically reduces the number of Optimizer() calls.194
Hybrid KG: given xn+1, following KGMC we use nz = 5 values of Zj , construct µn+1j (x) and use Optimizer() to find the peak location x∗j . The set of peak locations X ∗ d = {x∗1, ..., x∗nz} is used as a dynamic optimized discretization in KGd thus analytically computing an extremely tight lower bound of the true KG(x). Let µ∗ = µn(X∗d ) and similarly for σ̃ ∗(xn+1), Hybrid KG is given by
KGh(xn+1) = EZ [maxµ∗ + σ̃∗(xn+1)Z ] −maxµ∗.
Compared to KGd, the hybrid method removes redundant points in the discretization Xd, all X∗d195 points contribute to maxµn+1(X∗d ) and there are far fewer points. Compared to KGMC that samples196 many Zj and optimizes many x∗j , Hybrid KG uses far fewer Zj and optimizes far fewer x ∗ j . See197 Algorithm 1 for a summary of evaluating ConBO with Hybrid KG. ConBO can be reduced to the198 REVI algorithm by replacing the dynamic importance sampled tasks with a pre-frozen discretization199 of tasks and the (dynamic optimized) hybrid KG with discrete KG over a pre-frozen discretization of200 inputs. These changes drastically improve scalability enabling ConBO to be applied to a far broader201 range of applications.202
To ensure asymptotic convergence, in a discrete domain, we require that the acquisition function203 is non-negative, KGh(x) ≥ 0, and the acquisition function is zero where GP variance is zero,204 KGh(x) = 0 ⇐⇒ kn(x, x) = 0. Therefore, always choosing xn+1 = arg max KGh(x) ensures205 only points with GP variance will be revisited until all points have no variance i.e. the true function is206 known for all points. We can ensure these properties by setting nz ≥ 2 and at least one Zj is equal to207 zero.208
Theorem 3 Let nz ≥ 2 and let Z = {Zj |j = 1, ..., nz}. If 0 ∈ Z then KGh(x) ≥ 0 for all x ∈ X209 and if x is sampled infinitely often KGh(x) = 0.210
Proof is in the SM 1. The Zj values can be fixed, for nz = 5 we use equal Gaussian quantiles211 Z = { Φ−1(0.1),Φ−1(0.3), . . . ,Φ−1(0.9) } where Φ(·) is the Gaussian CDF. Using quantile spacing212 and odd nz ensures Zj = Φ−1(0.5) = 0 is included which satisfies the assumptions of asymptotic213 convergence. See Fig.2 (right).214
The computational complexity of a single call to ConBO requires the posterior variance (O(n2)) and215 nsnz runs of Optimizer(). Let ncalls be the number of times that Optimizer() calls the posterior216 mean costing O(n). Thus, ConBO total complexity is O(n2 + nncalls nsnz). Note this is linear in217 nsnz , the size of the (small) dynamic optimal discretization over S ×X . Thompson sampling with218
discretization uses one operation that scales cubically as O(n2(nsnz) + n(nsnz)2 + (nsnz)3) in219 the worst case and to reduce cost special techniques are required e.g. Fourier features, CG matrix220 inversion.221
3 Experiments222
We consider synthetic benchmarks and the three applications described in Section 1. In SM 4 we223 also present parallelization batch sampling results, and experiments with contextual, multi-task and224 multi-fidelity methods in SM 6. We also briefly test Hybrid KG for global optimization, SM 5, and225 observe that KGMC performs worse in the same computation time hence we exclude KGMC methods226 from the conditional experiments. For each benchmark, for evaluation, held out test tasks are sampled227 from P[s] ∝W (s), for each test task, stesti , the predicted optimal input is computed, x∗(stesti ) and228 we report the true black box output averaged over all test tasks, see SM 3 for details.229
3.1 How Accurate is Hybrid Knowledge Gradient?230
The theoretical Knowledge Gradient cannot be analytically computed and must be approximated231 (similarly, entropy based methods are not fully tractable). To illustrate the quality of the approxima-232 tions, we collected 20 points from the Rosenbrock test function and fit a Gaussian process. Another233 point was then randomly selected and true KG was estimated by the different methods, keeping the234 20 + 1 points fixed and only varying the nz parameter. Each estimate of KG was recalculated 50235 times and Table 1 shows the mean and two standard deviations. Discrete KG results in an extremely236 loose, volatile bound. The (state of the art) MC method with nz = 50 has a run-to-run variance of237 ±35%, variance reduction techniques used in MC++ (latin hypercube inverse sampling and control238 variates) can reduce this to ±6%. Hybrid KG with nz = 3 is extremely stable, run-to-run variance239 is too small to show, and is within the error margins of MC++ with nz = 50 which is ∼ 17× more240 expensive to compute. Increasing nz tightens the bound, hybrid KG with nz = 5 has 98.2% of the241 value of nz = 50. As a result, Hybrid KG is much easier to optimize, the Adam optimizer may be242 used with a much higher learning rate and lower momentum and thus converges much faster.243
As an aside, the recently proposed one-shot KG [34] enhances the optimization of KGMC . By244 freezing the Z values between calls to KGMC thus X̃∗ may be reused, this enables joint gradient245 ascent of (xn+1, X̃∗). In the global optimization use case (a constant single task) one-shot KG and246 Hybrid KG can be combined. In the conditional setting, the X̃∗i of past sampled tasks may be saved,247 however each call to ConBO must sample new tasks that are not in the saved history, one-shot KG248 style joint optimization is not possible. In our implementation we utilise caching of X̃∗i from old249 tasks to heuristically warm start finding X̃∗i for new tasks, see SM 3.250
3.2 Synthetic Functions251
We perform low-dimensional toy experiments in an ideal setting as a sanity check where we expect252 all conditional methods to perform similarly. We use the popular Branin-Hoo and Rosenbrock test253 functions in 2D defining the (task, input) domain as displayed in Figure 1.254
Synthetic Functions Using the Rosenbrock and Branin-Hoo functions, we consider a uniform and a255 triangular task weighting W (s). For baselines, we adopt two policy based methods. KNN: (dummy256 baseline) randomly collect data, x∗KNN (s) takes a task, s, and returns the best input from 10 nearest257 neighbor tasks. PG: policy gradient, a parametric quadratic policy x∗PG(s) = πθ(s) is learnt by258 maximising observed performance values, each iteration samples a task from P[s] ∝W (s) then x is259 sampled with an -greedy strategy. For a controlled ablation study, all the BO methods fit the same260 GP and for evaluation use the same definition of x∗N (s) (Equation 3), methods only differ by their261 data acquisition strategy. UNI: random data collection, the most recent conditional methods PEQI,262 REVI, MTS with all parameters given in SM 3. ConBO-nz: given (s, x), 20 tasks are importance263 sampled, KGh with nz = 3, 5 points is used. EI: expected improvement that treats (s, x) as inputs to264 be optimized.265
Results are shown in Figure 3. Policy based methods KNN and PG consistently perform worse than266 the Gaussian process methods. Surprisingly, the conditional BO algorithm PEQI performs similarly to267 UNI and much worse than EI. All other conditional methods outperform all non-conditional methods.268
3.3 CNN Training Hyperparameters269
We apply MTS, REVI, and ConBO variants, and we adopt the recently proposed kernel used for BO270 with Common Random Numbers [35],271
k((s, x), (s′, x′)) = σ20M(x, x ′; l) + δs′s(σ 2 1M(x, x ′; l) + σ23).
where M(x, x′; l) is a Matérn 52 kernel with length scales l. The first term models a common trend272 function across all tasks and the second term models how each task independently differs from the273 trend. The differences are composed of another Matérn and the constant kernel to model a global274 offset e.g. one dataset may have universally higher validation accuracy. This kernel has far fewer275 parameters than a full multi-task product kernel, it is easy to fit and scales to an arbitrary number of276 tasks (or datasets) without adding extra parameters.277
In this problem setting, learning hyperparameters over similar datasets, one may expect that the278 optimal hyperparameters would be the same for all datasets. Therefore, as a baseline we apply EI279 to learn the hyperparameters of the first dataset (task 1). We then evaluate the objective function280 (validation accuracy) on the rest of the datasets using the best observed hyperparameters from dataset281 1, we refer to this as EI + Transfer.282
Argument Optimization versus Data Optimization In continuous task settings, it is not possible to283 evaluate every task. In discrete task settings with large sampling budgets N |S|, a user may desire284 a single high (although stochastic) output value, max y, for every task. For example, in network285 hyperparameter optimization, the network with the best validation error, max y, will be deployed286 (and the hyperparameters x may or may not be reused). We refer to this is data optimization (DO).287 For simulated environments, the input (or “action") that generalizes providing the best long-term288 average performance, maxx E[f(s, x)], is deployed and we refer to this is argument optimization289 (AO). Past work [22] has shown that argument optimization finds inputs that generalize better but290
may not provide optimal max y for all tasks during the optimization run. The authors propose a “DO291 trick”, use an AO method for N − |S| iterations, then finally allocate one sample per task with input,292 xn+1, determined by EI within the task. We apply this trick to all algorithms in this experiment.293
Results are shown in Figure 4. For the medium budget of 50 samples, ConBO performs best of294 the standard algorithms yet it is still worse than the EI+Tr baseline. Applying the final round of295 DO improves all results to match the baseline. For the large budget of 100 samples, all methods296 outperform the baseline suggesting dataset specific fine-tuning of hyperparameters is required to297 achieve best results. Again, DO provides a significant boost to performance for all methods.298
Gaussian process kernel parameter learning required approximately 2–5 seconds using Tensorflow. In299 Figure 4 (right) we show the runtime of each algorithm excluding model fitting and network training,300 purely acquisition function optimization time. MTS and ConBO-3 are quickest while Conbo-5301 increases linearly over ConBO-3 and REVI takes much longer.302
3.4 Ambulances and Warehouses303
We apply all methods from Section 3.2 to two benchmarks from the www.SimOpt.org library for304 simulation optimization problems. The ambulance problem (AMB) is 8-dimensional and consists of305 a range of cities and one must optimize ambulance locations for each city. The Assemble-to-order306 problem (ATO) is 9-dimensional consisting of a range of warehouses and one must optimize target307 stock level for each warehouse. Results are shown in Figure 5. Of the policy based methods, PG308 performs poorly and does not show on the plots whilst KNN performs poorly on AMB and performs309 well on ATO suggesting that AMB is a more difficult problem. Of the GP based methods, EI performs310 well for smaller budgets. Although it is not a conditional algorithm we include it to highlight that311 sometimes the simplest idea can also work. Of the conditional methods, MTS, REVI, and ConBO-3312 all perform similarly, either slightly (AMB) or largely (ATO) outperforming UNI. These methods313 struggle in higher dimensions while ConBO-5 uses a more accurate acquisition function and is the314 only method that consistently performs well across all problems. We hypothesize that these problems315 are more difficult than the synthetics and CNN and truly stress test conditional algorithms.316
4 Conclusion317
Potential Limitations and Broader Impact there are multiple ways in which ConBO may fail, in318 this work we have not investigated how ConBO or Hybrid KG suffers with poorly learnt Gaussian319 process hyperparameters. In many applications, a poorly chosen kernel or unoptimized hyperparame-320 ters can lead to poor performance and our proposed methods may be more sensitive or more robust to321 these failures than alternative approaches. We propose a general purpose optimization algorithm and322 analysis, we do not currently see any immediate societal impact.323
We investigate Conditional Bayesian optimization and propose ConBO. ConBO is designed from324 the ground up to fully exploit the structure of conditional problems, namely that optimizing one task325 helps optimize similar tasks. Hence every point should be collected to maximise the benefit of all326 tasks. However, this can lead to excessive computational cost, particularly in higher dimensions.327 Thus we also propose Hybrid KG that mixes past methods to be both fast and scalable.328
References329 [1] Victor Picheny and David Ginsbourger. A nonstationary space-time gaussian process model for330 partially converged simulations. SIAM/ASA Journal on Uncertainty Quantification, 1(1):57–78,331 2013.332
[2] Shinkyu Jeong, Mitsuhiro Murayama, and Kazuomi Yamamoto. Efficient optimization design333 method using kriging model. Journal of aircraft, 42(2):413–420, 2005.334
[3] Jasper Snoek, Hugo Larochelle, and Ryan P Adams. Practical bayesian optimization of machine335 learning algorithms. In Advances in neural information processing systems, pages 2951–2959,336 2012.337
[4] Diariétou Sambakhé, Lauriane Rouan, Jean-Noël Bacro, and Eric Gozé. Conditional opti-338 mization of a noisy function using a kriging metamodel. Journal of Global Optimization,339 73(3):615–636, 2019.340
[5] Naijia Anna Dong, David J Eckman, Xueqi Zhao, Shane G Henderson, and Matthias Poloczek.341 Empirically comparing the finite-time performance of simulation-optimization algorithms. In342 2017 Winter Simulation Conference (WSC), pages 2206–2217. IEEE, 2017.343
[6] Jing Xie, Peter I Frazier, and Stephen E Chick. Bayesian optimization via simulation with344 pairwise sampling and correlated prior beliefs. Operations Research, 64(2):542–559, 2016.345
[7] Matthias Poloczek, Jialei Wang, and Peter I Frazier. Warm starting bayesian optimization. In346 2016 Winter Simulation Conference (WSC), pages 770–781. IEEE, 2016.347
[8] Kevin Swersky, Jasper Snoek, and Ryan Prescott Adams. Multi-task bayesian optimization. In348 Advances in Neural Information Processing Systems. Curran Associates, Inc., 2013.349
[9] Supratik Paul, Konstantinos Chatzilygeroudis, Kamil Ciosek, Jean-Baptiste Mouret, Michael350 Osborne, and Shimon Whiteson. Alternating optimisation and quadrature for robust control. In351 Proceedings of the AAAI Conference on Artificial Intelligence, volume 32, 2018.352
[10] Saul Toscano-Palmerin and Peter I Frazier. Bayesian optimization with expensive integrands.353 arXiv preprint arXiv:1803.08661, 2018.354
[11] Michael Pearce and Juergen Branke. Bayesian simulation optimization with input uncertainty.355 In 2017 Winter Simulation Conference (WSC), pages 2268–2278. IEEE, 2017.356
[12] Wei Chu, Lihong Li, Lev Reyzin, and Robert Schapire. Contextual bandits with linear payoff357 functions. In Geoffrey Gordon, David Dunson, and Miroslav Dudík, editors, Proceedings of358 the Fourteenth International Conference on Artificial Intelligence and Statistics, volume 15 of359 Proceedings of Machine Learning Research, pages 208–214, Fort Lauderdale, FL, USA, 11–13360 Apr 2011. JMLR Workshop and Conference Proceedings.361
[13] Wei Chu, Lihong Li, Lev Reyzin, and Robert Schapire. Contextual bandits with linear payoff362 functions. In Proceedings of the Fourteenth International Conference on Artificial Intelligence363 and Statistics, pages 208–214. JMLR Workshop and Conference Proceedings, 2011.364
[14] Supratik Paul, Michael A Osborne, and Shimon Whiteson. Contextual policy optimisation.365 arXiv preprint arXiv:1805.10662, 2018.366
[15] Ian Char, Youngseog Chung, Willie Neiswanger, Kirthevasan Kandasamy, Andrew Oakleigh367 Nelson, Mark Boyer, Egemen Kolemen, and Jeff Schneider. Offline contextual bayesian368 optimization. In Advances in Neural Information Processing Systems, pages 4629–4640, 2019.369
[16] Deng Huang, Theodore T Allen, William I Notz, and R Allen Miller. Sequential kriging370 optimization using multiple-fidelity evaluations. Structural and Multidisciplinary Optimization,371 32(5):369–382, 2006.372
[17] Matthias Poloczek, Jialei Wang, and Peter Frazier. Multi-information source optimization. In373 Advances in Neural Information Processing Systems, pages 4289–4299, 2017.374
[18] Kirthevasan Kandasamy, Gautam Dasarathy, Jeff Schneider, and Barnabás Póczos. Multi-fidelity375 bayesian optimisation with continuous approximations. In Proceedings of the 34th International376 Conference on Machine Learning-Volume 70, pages 1799–1808. JMLR. org, 2017.377
[19] Aaron Klein, Stefan Falkner, Simon Bartels, Philipp Hennig, and Frank Hutter. Fast Bayesian378 Optimization of Machine Learning Hyperparameters on Large Datasets. In Aarti Singh and379 Jerry Zhu, editors, Proceedings of the 20th International Conference on Artificial Intelligence380 and Statistics, volume 54 of Proceedings of Machine Learning Research, pages 528–536, Fort381 Lauderdale, FL, USA, 20–22 Apr 2017. PMLR.382
[20] R. Bardenet, Mátyás Brendel, Balázs Kégl, and Michele Sebag. Collaborative hyperparameter383 tuning. In International Conference on Machine Learning, pages 199–207, 2013.384
[21] David Ginsbourger, Jean Baccou, Clément Chevalier, Frédéric Perales, Nicolas Garland, and385 Yann Monerie. Bayesian adaptive reconstruction of profile optima and optimizers. SIAM/ASA386 Journal on Uncertainty Quantification, 2(1):490–510, 2014.387
[22] Michael Pearce and Juergen Branke. Continuous multi-task bayesian optimisation with correla-388 tion. European Journal of Operational Research, 270(3):1074–1085, 2018.389
[23] C. E. Rasmussen and C. K. I. Williams. Gaussian Processes for Machine Learning. MIT Press,390 2004.391
[24] D. R. Jones, M. Schonlau, and W. J. Welch. Efficient global optimization of expensive black-box392 functions. Journal of Global optimization, 13(4):455–492, 1998.393
[25] Andreas Krause and Cheng S Ong. Contextual gaussian process bandit optimization. In394 Advances in neural information processing systems, pages 2447–2455, 2011.395
[26] Philipp Hennig and Christian J Schuler. Entropy search for information-efficient global opti-396 mization. Journal of Machine Learning Research, 13(Jun):1809–1837, 2012.397
[27] J. Villemonteix, E. Vazquez, and E. Walter. An informational approach to the global optimization398 of expensive-to-evaluate functions. Journal of Global Optimization, 44(4):509–534, 2009.399
[28] J. M. HernÂŽandez-Lobato, M. W. Hoffman, and Z. Ghahramani. Predictive entropy search400 for efficient global optimization of black-box functions. In Advances in Neural Information401 Processing Systems, pages 918–926. Curran Associates, Inc., 2014.402
[29] Zi Wang and Stefanie Jegelka. Max-value entropy search for efficient bayesian optimization.403 In Proceedings of the 34th International Conference on Machine Learning-Volume 70, pages404 3627–3635. JMLR. org, 2017.405
[30] P. Frazier, W. Powell, and S. Dayanik. The knowledge-gradient policy for correlated normal406 beliefs. INFORMS Journal on Computing, 21(4):599–613, 2009.407
[31] Diederik P Kingma and Jimmy Ba. Adam: A method for stochastic optimization. arXiv preprint408 arXiv:1412.6980, 2014.409
[32] Warren Scott, Peter Frazier, and Warren Powell. The correlated knowledge gradient for410 simulation optimization of continuous parameters using gaussian process regression. SIAM411 Journal on Optimization, 21(3):996–1026, 2011.412
[33] Jian Wu and Peter I Frazier. Discretization-free knowledge gradient methods for bayesian413 optimization. arXiv preprint arXiv:1707.06541, 2017.414
[34] Maximilian Balandat, Brian Karrer, Daniel Jiang, Samuel Daulton, Ben Letham, Andrew G Wil-415 son, and Eytan Bakshy. Botorch: A framework for efficient monte-carlo bayesian optimization.416 Advances in Neural Information Processing Systems, 33, 2020.417
[35] Michael Pearce, Matthias Poloczek, and Juergen Branke. Bayesian optimization allowing for418 common random numbers. arXiv preprint arXiv:1910.09259, 2019.419
Checklist420
1. For all authors...421 (a) Do the main claims made in the abstract and introduction accurately reflect the paper’s422 contributions and scope? [Yes] ConBO is theoretically grounded, see Theorems 1-3.423 ConBO outperforms recently published works on a range of problems, see Section 3,424 ConBO is easily parallelizable, see SM Section 4.425
(b) Did you describe the limitations of your work? [Yes] see Section 4426 (c) Did you discuss any potential negative societal impacts of your work? [Yes] see Section427 4428 (d) Have you read the ethics review guidelines and ensured that your paper conforms to429 them? [Yes]430 2. If you are including theoretical results...431
(a) Did you state the full set of assumptions of all theoretical results? [Yes] see SM Section432 1433
(b) Did you include complete proofs of all theoretical results? [Yes] see SM Section 1434 3. If you ran experiments...435
(a) Did you include the code, data, and instructions needed to reproduce the main exper-436 imental results (either in the supplemental material or as a URL)? [Yes] All code is437 submitted with the SM438
(b) Did you specify all the training details (e.g., data splits, hyperparameters, how they439 were chosen)? [Yes] all parameters for all baselines are given in SM section 3.440
(c) Did you report error bars (e.g., with respect to the random seed after running experi-441 ments multiple times)? [Yes]442
(d) Did you include the total amount of compute and the type of resources used (e.g., type443 of GPUs, internal cluster, or cloud provider)? [Yes] SM section 3444
4. If you are using existing assets (e.g., code, data, models) or curating/releasing new assets...445 (a) If your work uses existing assets, did you cite the creators? [Yes] code from Simopt.org446 is cited.447 (b) Did you mention the license of the assets? [N/A]448 (c) Did you include any new assets either in the supplemental material or as a URL? [Yes]449 the benchmark test problems.450 (d) Did you discuss whether and how consent was obtained from people whose data you’re451 using/curating? [N/A]452 (e) Did you discuss whether the data you are using/curating contains personally identifiable453 information or offensive content? [N/A]454 5. If you used crowdsourcing or conducted research with human subjects...455
(a) Did you include the full text of instructions given to participants and screenshots, if456 applicable? [N/A]457
(b) Did you describe any potential participant risks, with links to Institutional Review458 Board (IRB) approvals, if applicable? [N/A]459
(c) Did you include the estimated hourly wage paid to participants and the total amount460 spent on participant compensation? [N/A]461

Review summaries:
This paper proposes a new algorithm for a problem entitled "conditional Bayesian optimization." The goal is to search for a *series of input locations* that maximizes the performance of a given set of tasks. Unlike multi-fidelity optimization, the objective is to find good values for all tasks rather than the high-fidelity task. The main contribution is a Hybrid Knowledge Gradient (KG) algorithm that mixes two ways of computing KG: one that discretizes the input space (x, t) and computes KG for a discrete set of points,  and another one that uses the reparametrization trick to perform stochastic gradient descent over Monte Carlo samples. Formal guarantees about this method are also given. Empirical results show that the proposed method effectively reduces computation time and improves the quality of the solution.

This paper presents a framework for conditional optimization, which is called ConBO. Authors also proposed a new acquisition function: hybird knowledge gradient, which is built on KG. Empirical experiments show their proposed method outperform many competing algorithms. 

The authors tackle an underexplored problem they call conditional BO, which is a variant of multi-task BO in which one seeks to identify the optimum of each task, instead of the traditional multi-task BO seeks to identify the optimum of the average over all tasks.   The authors introduce ConBO, which uses the knowledge gradient acquisition function to determine the next point to evaluate, and discuss its efficient computation through a combination of intelligent sampling and discretization.   ConBO is evaluated on a variety of problems, including three good applications: hyperparameter optimization, ambulance base placement, and assemble-to-order policy selection. Results are promising, and suggest that ConBO outperforms a well-chosen set of sensible baselines.

The authors propose a method for conditional optimization problems: finding the optimum of each function in set of objective functions. In this setting, the method jointly selects the parameters to evaluate and the task to evaluate. The authors propose a theoretically-grounded variant of the Knowledge Gradient to address this class of problems and propose a new method of approximating the KG acquisition function.