Optimistic Optimization of Gaussian Process Samples

We consider maximizing a function $f:X\to R$. Throughout, $f^{\ast }$ denotes the maximum of the function, and $x^{\ast }$, the point where it is attained. $k:X\times X\to R$ is a kernel function. The function $f$ is assumed to be a sample from a GP $GP(\mu ,k)$. We assume that the GP is centered, i.e. $\mu =0$. One has access to observations $y_i$, where $y_i\sim f\left(x_i\right)+N(0,\lambda )$. The noiseless setting amounts to the assumption $\lambda \to 0$. After each new observation the posterior over the function is updated:

where $y_n=[y_1,...,y_n{]}^T,k_n(x)=[k(x_1,),...,k(x_n,){]}^T$, and $K_n$ is the Gram matrix with $K_n{}_{ij}=k_n(x_i,x_j)$. The posterior is used to define an acquisition function $a_n\left(x\right)$ at whose maximum the function is evaluated next. For an overview of different options of acquisition functions see Frazier [2018]. We will focus our analysis on the choice of GP-UCB Srinivas et al. [2009], where

An advantage of GP-UCB over other choices is its comparatively simple structure, which facilitates not just implementation but also analysis. While other aquisition functions may well be more sample-efficient in concrete settings, this will not be important for what is to follow, because we are primarily concerned with run-time complexity, where GP-UCB is arguably indeed the fastest possible choice.

The optimistic optimization principle is used to optimize functions that are known to fulfill a local Lipschitz assumption with respect to a dissimilarity $d:X\times X\to R$:

The method revolves around a hierarchical partitioning of the search space $X$ that can be described by an infinite binary tree. Search can thus be implemented as a fast descent through the tree, thus only involves evaluations on a countable set of mesh points, in contrast to the local numerical optimization done within Bayesian optimization. This is one of the two reasons for the drastically shorter run-time costs of OO over BO (the other being that it does not require computing a GP posterior).

The root node corresponds to the entire search space $X_{(0,1)}=X$ and is named $(0,1)$. Consider a node $(t,i)$ at depth $t$. The left child $(t+1,2i-1)$ and right child $(t+1,2i)$ represent two subregions $X_{(t+1,2i-1)}$ and $X_{(t+1,2i)}$ such that $X_{(t,i)}=X_{(t+1,2i-1)}\cup X_{(t+1,2i)}$, i.e. the tree covers the entire space. To indicate that a cell was explored at iteration $n$, we refer to it as $X_n$. Intuitively, it makes sense to select the cells in such a way that all points in a cell are similar to each other and all similar points are in the same cell. Formally, this can be expressed as

• $\forall x,y\in X_{(t,i)}:\left\{d\right(x,y)<\delta (t\left)\right\}$

• $\exists x_{(t,i)}\in X_{(t,i)}:\{y\in X:d(x_{(t,i)},y)<c\delta (t\left)\right\}\subset X_{(t,i)}$,

where $\delta \left(t\right)$ is a decreasing sequence of diameters and $c$ is a global constant. During search, the tree is build incrementally by adding the two children of a selected node. When a new node ($t,i$) is added, an observation is made at the center $x_{(t,i)}$ of the region $X_{(t,i)}$. In each round, the leaf with the highest upper bound $U_{(t,i)}=f\left(x_{(t,i)}\right)+\delta \left(t\right)$ is selected for expansion. Since $f$ is assumed to be locally Lipschitz with respect to $d$, selecting nodes by $U_{(t,i)}$ is a valid upper bound strategy (assuming noiseless observations). The sequence of $\delta \left(t\right)$ controls exploration. Therefore, the partitioning should be chosen in such a way, that the $\delta \left(t\right)$ can be as small as possible. The optimistic optimization principle is summarized in Algorithm 1 and Figure 1. Using a binary heap, the priority queue for the leaf nodes can be realized in $O\left(N\log N\right)$ if the budget $N$ is known in advance.

This algorithm is a batch method in the sense that all children of a newly explored node are evaluated in each step. We assume batch size two, but other choices are possible too and might be preferable in some applications, e.g. when parallel evaluation is possible. A conceptual difference to the upper bound in GP-UCB is that one here works with upper bounds for entire regions of the search space, not for individual points in it. Also, the exploration term $\delta \left(t\right)$ in the upper bound is not updated based on new observations, but derived a priori. Nevertheless, the exploration term, aka. the uncertainty decreases during the search as the tree grows.

The hierarchical optimistic optimization principle has its origins in the Bandit setting. Bubeck et al. [2011] apply it in the noisy setting, Munos [2011] apply it in the noiseless setting and Kleinberg et al. [2008] uses it with a slightly different Lipschitz assumption. The assumptions on the dissimilarity $d$ vary, e.g. some theoretical analyses require that it additionally is a semi-metric, or metric. Munos [2011] introduces a variant of the principle, simultaneous optimistic optimization (SOO), for situations in which the dissimilarity function is unknown. Valko et al. [2013] extend this idea to the noisy setting.

An interim summary: OO descends along a search tree, i.e. in a discrete sequence of steps, without numerical optimization, and only using local summary statistics, rather than updating a global posterior. This makes OO very fast, at least compared to Bayesian optimization. But the approach also has a downside: At least at first sight, it is not clear how to encode salient prior information about the global optimization problem into the algorithm. By contrast, many Bayesian optimization experts see the rich language of GP prior models as a key strength of their framework. The following section thus investigates formal connections between OO and BO. The goal is to understand to which degree the structural language of a GP prior can be translated into the algorithmic efficiency of OO.

The canonical pseudo-metric $d:X\times X\to R$ for a centered GP is defined as

The fundamental relevancy of the canonical pseudo-metric arises from the fact that it controls the increments of the GP, in the sense of the following deviation inequality (Pisier [1999], Theorem 4.7):

An important class of kernel-induced metrics is formed by monotonic transformations of the Euclidean metric, i.e. $d(x,y)=g(\parallel x-y{\parallel }_2)$ where $g:R\mapsto R$ is monotonically increasing. Many kernels used in practice are in this class, e.g. the square-exponential kernel, the Matérn class of kernels, the rational-quadratic kernels and the Wiener kernel as well as sums and products thereof. We refer to this class of kernels as $K$. Kernels violating this assumptions are e.g. polynomial or periodic kernels.

The first challenge in applying the optimistic optimization principle on samples of a GP consists in the probabilistic nature of the deviation inequality. To obtain a valid upper bound for the maximal deviation in a cell ${sup}_{x\in X_n}|f\left(x\right)-f\left(x_n\right)|$, the deviation for all points in the cell has to be bounded. We approximate such an upper bound by introducing two simplifications: We discretize the search space $X$ into a finite number of points $\hat{X}$ and we introduce an independence assumption between the $|f\left(x\right)-f\left(x_n\right)|$. Then we take a union bound approach:

The bounds have to hold at each step $n$, so we additionally take a union bound over the number of steps. This implies the following statement, that holds with high probability :

where ${\beta }_n$ are appropriate constants specified in the Supplements. The union bound approximation will be good if the $|f\left(x\right)-f\left(x_n\right)|$ are (nearly) uncorrelated, or the size of the discretization $|\hat{X}|$ is chosen sufficiently small. Otherwise the bounds are loose, which leads to over-exploration. Though approximate, this idea of neglecting correlations to simplify the calculation of the expected supremum of dependent Gaussian variables is, for example, also done in Maximum Value Entropy Search [Wang and Jegelka, 2017] for BO or in related settings [Grosse et al., 2021]. The other extreme, a greedy approach with ${\beta }_n=1$ has also been taken in recent work [Rando et al., 2022]. With generic chaining [Talagrand, 1996] it is possible to improve over the union bound approach. However, to the best of our knowledge state-of-the-art algorithms [Borst et al., 2020] to optimize for tighter bounds require polynomial run-time in the number of points per cell for arbitrary kernel functions. For special cases, like a Wiener kernel [Talagrand, 2021] or Matérn kernel functions [Shekhar and Javidi, 2018] analytical attempts to derive chaining based upper bounds exist.

The second main step in applying the OO principle is to choose the cells and location of the centers in such a way that the diameters of the cells are as small as possible. For $k$ children nodes, this is a metric $k$-center problem – one of the classical NP-hard problems [Gonzalez, 1985]. A greedy approximation consists in iteratively picking the $k$ centers with the largest distance to the previously picked centers, and requires $O\left(|\hat{X}|k\right)$ time. The greedy procedure is guaranteed to result in a $2$-approximation, and there is no polynomial time algorithm doing better (unless P$=$NP). Working with a greedy instead of the optimal partitioning scheme thereby leads to an additional factor of $2$ in the below regret bound, but is not harmful in the sense that the search gets stuck in a local optimum. NP-hardness also appears in the context of BO, e.g. the exploration term used in GP-UCB. And the acquisition functions of information-theoretic BO methods [Hennig and Schuler, 2012], [Wang and Jegelka, 2017] are related to the maximization of information gain, which is also a NP-hard problem.
From an implementation perspective, it is desirable to constrain the partitioning to axis-parallel boxes. For some kernel functions, e.g. the polynomial kernel, requirement (b) of the OO principle cannot be fulfilled with axis parallel boxes. One can nevertheless run the algorithm, but it will clearly be less information-efficient. One may even consider a randomized choice of centers, e.g. as done in Monte Carlo Tree Search [Chaslot et al., 2008]. However, it still remains to calculate the maximal distance from a point in the cell $X$ to the center $x_n$. The computational complexity of this part is comparable to the numerical optimization of the acquisition function in BO. An advantage is that the domain over which one optimizes shrinks in each step. A disadvantage, though, is that if the numerical optimization is suboptimal and the maximal distance within a cell is underestimated, cells containing the optimum might get irreversibly pruned.
An important observation is that the partitioning scheme itself does not depend on the objective $f$, but only the kernel/distance function (how the search tree grows, however, does depend on $f$). This opens up the possibility of finding a good partitioning and the corresponding maximal distances analytically. For distances derived from kernels in $K$, one can apply the following regular partitioning scheme: At each step, cut along the longest dimension in order to obtain the two children cells. Use the euclidean centers as centers. A point that maximizes the distance to the center will always be one of the $2^m$ corner points in $m$ dimensions. Thus, for this type of kernels, the costs reduce to $O\left(N{2}^m\right)$ for the partitioning, or $O(N\log N+N{2}^m)=O\left(N\log N\right)$ in total.

Motivated by the analysis above, we propose a variant of OO, which we call GP-OO. It consists in running Algorithm 1 with the utility $U_{\left(n\right)}=f\left(x_{\left(n\right)}\right)+{\beta }_n^{1/2}\Delta \left(X_n\right)$ in line 4, where $\Delta$ is as defined in Eq. (5). For kernels that match the observation above, one can then still choose the evaluation node (in $O\left(2^m\right)$ time) without numerical optimization. Figure 2 shows the algorithm running on GP samples from a square-exponential and a Matérn kernel with $\nu =2/3$ on the domain $[0,5{]}^2$, as well as from a quadratic kernel on the domain $[-1,1{]}^2$. In regions with higher function values the partitions are finer. The Matérn kernel yields higher distances than the square exponential, leading to more exploration, reflecting that the samples are less smooth. By optimizing the partitions and centers of the cells with respect to the canonical pseudo-metric, larger parts of the search space can be covered while keeping the cell diameters constant as shown by the two examples with the quadratic kernel.

Regular partitions. We begin with on-model experiments for a common setting in Bayesian optimization, where GP-OO can showcase its speed advantages: 20 samples from a GP with squared exponential kernel (lengthscale $l=0.1$) on the unit cube $X=[0,1{]}^3$. For all experiments, the noise level for GP-UCB is set to a very small constant $\lambda =0.001$ since we assume noiseless observations. For GP-OO, we use ${\beta }_n=2\log \left(|\widehat{X_n}|N/2ϵ|\right)$ with fixed horizon $N=10000$, and for GP-UCB, we use ${\beta }_n=2\log \left(|\hat{X}|n^2{\pi }^2/6ϵ\right)$ according to theory. We did a grid search over the size of the discretization $|\hat{X}|$ with values in $\{1,10,100,1000\}$. The confidence level $ϵ$ is set to $0.01$ for both methods, and we set $|\hat{X}|=1$. Even though GP-UCB is more sample-efficient than GP-OO (Figure 3, left) GP-OO is orders of magnitude faster per step (Figure 3, right).

This suggests there is a “sweet-spot”: We artificially added $c\in \{0.01,0.1,1,10\}$ seconds to the per-evaluation costs to simulate different evaluation costs (Figure 4). For evaluation time around 0.1 seconds and below, GP-OO becomes faster than GP-UCB. When the objective becomes “expensive”, BO can leverage its sample efficiency (at the very extreme end of nearly instant function-evaluations, random search, which has negligible overhead, is always a trivially optimal asymptotic baseline).

Non-regular partitions. For demonstration, we also explore a setting with non-regular optimal partitions, where GP-OO can not perform as well. We sample 100 functions from a GP with quadratic kernel $k(x,y)=(x^{T}y{)}^2$ with bias 0 on the domain $[-1,1{]}^3$. We consider GP-UCB and GP-OO, once with euclidean partitions and once with partitions optimized with respect to the canonical metric. The exploration constant was set to $\beta =1$. The results shown (Figure 5) show that it is advantageous to optimize the partitioning scheme with respect to the canonical metric. However, we restrict ourselves to partitions with axis-aligned cells, i.e. the perfectly correlated corner points with $k(x,y)=1$ and $d(x,y)=0$ end up in different cells. At this point, GP-UCB has a clear advantage, as information can be shared across the search space between the corner points. Also, the run-time advantage decreases due to the numerical optimization of the partitions.

We consider the minimization of common optimization benchmark functions, on domains with dimensions ranging from 2 to 10. Since the functions are deterministic, and GP-OO is deterministic up to random tie breaks, we randomly sample 10 domains $X$ to obtain variation (see Supplements). We use a Matérn kernel with $\nu =3/2$ and length-scales as recommended by Rando et al. [2022]. For the exploration constant $\beta$ we did a grid search over $\{0.1,1,10,100\}$. We use confidence level $ϵ=0.05$. Figure 6 shows minimal attained function value over time. Here, GP-OO outperforms GP-UCB in many cases. This is not just due to runtime, but in fact also sometimes holds in function evaluations (see Supplements). Here, both methods did not always find the global optimum (which can not be guaranteed here as these functions are not necessarily within the RKHS).