b

DiscoverSearch
About
Gaussian Process Optimization in the Bandit Setting: No Regret and Experimental Design
2009·arXiv
Abstract
Abstract

Many applications require optimizing an unknown, noisy function that is expensive to evaluate. We formalize this task as a multi-armed bandit problem, where the payoff function is either sampled from a Gaussian process (GP) or has low RKHS norm. We resolve the important open problem of deriving regret bounds for this setting, which imply novel convergence rates for GP optimization. We analyze GP-UCB, an intuitive upper-confidence based algorithm, and bound its cumulative regret in terms of maximal information gain, establishing a novel connection between GP optimization and experimental design. Moreover, by bounding the latter in terms of operator spectra, we obtain explicit sublinear regret bounds for many commonly used covariance functions. In some important cases, our bounds have surprisingly weak dependence on the dimensionality. In our experiments on real sensor data, GP-UCB compares favorably with other heuristical GP optimization approaches.

In most stochastic optimization settings, evaluating the unknown function is expensive, and sampling is to be minimized. Examples include choosing advertisements in sponsored search to maximize profit in a click-through model (Pandey & Olston, 2007) or learning optimal control strategies for robots (Lizotte et al., 2007). Predominant approaches to this problem include the multi-armed bandit paradigm (Robbins, 1952), where the goal is to maximize cumulative reward by optimally balancing exploration and exploitation, and experimental design (Chaloner & Verdinelli, 1995), where the function is to be explored globally with as few evaluations as possible, for example by maximizing information gain. The challenge in both approaches is twofold: we have to estimate an unknown function f from noisy samples, and we must optimize our estimate over some high-dimensional input space. For the former, much progress has been made in machine learning through kernel methods and Gaussian process (GP) models (Rasmussen & Williams, 2006), where smoothness assumptions about f are encoded through the choice of kernel in a flexible nonparametric fashion. Beyond Euclidean spaces, kernels can be defined on diverse domains such as spaces of graphs, sets, or lists.

We are concerned with GP optimization in the multi-armed bandit setting, where f is sampled from a GP distribution or has low “complexity” measured in terms of its RKHS norm under some kernel. We provide the first sublinear regret bounds in this nonparametric setting, which imply convergence rates for GP optimization. In particular, we analyze the Gaussian Process Upper Confidence Bound (GP-UCB) algorithm, a simple and intuitive Bayesian method (Auer et al., 2002; Auer, 2002; Dani et al., 2008). While objectives are different in the multi-armed bandit and experimental design paradigm, our results draw a close technical connection between them: our regret bounds come in terms of an information gain quantity, measuring how fast f can be learned in an information theoretic sense. The submodularity of this function allows us to prove sharp regret bounds for particular covariance functions, which we demonstrate for commonly used Squared Exponential and Mat´ern kernels.

Related Work. Our work generalizes stochastic linear optimization in a bandit setting, where the unknown function comes from a finite-dimensional linear space. GPs are nonlinear random functions, which can be represented in an infinite-dimensional linear space. For the standard linear setting, Dani et al. (2008) provide a near-complete characterization1

(also see Auer 2002; Dani et al. 2007; Abernethy et al. 2008; Rusmevichientong & Tsitsiklis 2008), explicitly dependent on the dimensionality. In the GP setting, the challenge is to characterize complexity in a differ-ent manner, through properties of the kernel function. Our technical contributions are twofold: first, we show how to analyze the nonlinear setting by focusing on the concept of information gain, and second, we explicitly bound this information gain measure using the concept of submodularity (Nemhauser et al., 1978) and knowledge about kernel operator spectra.

Kleinberg et al. (2008) provide regret bounds under weaker and less configurable assumptions (only Lipschitz-continuity w.r.t. a metric is assumed; Bubeck et al. 2008 consider arbitrary topological spaces), which however degrade rapidly with the dimensionality of the problem (Ω(T

linearity w.r.t. a fixed basis is often too stringent

an assumption, while Lipschitz-continuity can be too coarse-grained, leading to poor rate bounds. Adopting GP assumptions, we can model levels of smoothness in a fine-grained way. For example, our rates for the frequently used Squared Exponential kernel, enforcing a high degree of smoothness, have weak dependence on the dimensionality:  O(�T(log  T)d+1) (see Fig. 1).

There is a large literature on GP (response surface) optimization. Several heuristics for trading off exploration and exploitation in GP optimization have been proposed (such as Expected Improvement, Mockus et al. 1978, and Most Probable Improvement, Mockus 1989) and successfully applied in practice (c.f., Lizotte et al. 2007). Brochu et al. (2009) provide a comprehensive review of and motivation for Bayesian optimization using GPs. The Efficient Global Optimization (EGO) algorithm for optimizing expensive black-box functions is proposed by Jones et al. (1998) and extended to GPs by Huang et al. (2006). Little is known about theoretical performance of GP optimization. While convergence of EGO is established by Vazquez & Bect (2007), convergence rates have remained elusive. Gr¨unew¨alder et al. (2010) consider the pure exploration problem for GPs, where the goal is to find the optimal decision over T rounds, rather than maximize cumulative reward (with no exploration/exploitation dilemma). They provide sharp bounds for this exploration problem. Note that this methodology would not lead to bounds for minimizing the cumulative regret. Our cumulative regret bounds translate to the first performance guarantees (rates) for GP optimization.

Summary. Our main contributions are:

We analyze GP-UCB, an intuitive algorithm for GP optimization, when the function is either sam-

image

Figure 1. Our regret bounds (up to polylog factors) for linear, radial basis, and Mat´ern kernels — d is the dimension, T is the time horizon, and  νis a Mat´ern parameter.

image

We bound the cumulative regret for GP-UCB in terms of the information gain due to sampling, establishing a novel connection between experimental design and GP optimization.

By bounding the information gain for popular classes of kernels, we establish sublinear regret bounds for GP optimization for the first time. Our bounds depend on kernel choice and parameters in a fine-grained fashion.

We evaluate GP-UCB on sensor network data, demonstrating that it compares favorably to existing algorithms for GP optimization.

Consider the problem of sequentially optimizing an unknown reward function  f : D → R: in each round t, we choose a point  xt ∈ Dand get to see the function value there, perturbed by noise:  yt = f(xt)+ϵt. Our goal is to maximize the sum of rewards �Tt=1 f(xt), thus to perform essentially as well as  x∗ = argmaxx∈D f(x) (as rapidly as possible). For example, we might want to find locations of highest temperature in a building by sequentially activating sensors in a spatial network and regressing on their measurements. D consists of all sensor locations, f(x) is the temperature at x, and sensor accuracy is quantified by the noise variance. Each activation draws battery power, so we want to sample from as few sensors as possible.

Regret. A natural performance metric in this context is cumulative regret, the loss in reward due to not knowing f’s maximum points beforehand. Suppose the unknown function is f, its maximum point1x∗ = argmaxx∈D f(x). For our choice  xtin round t, we incur instantaneous regret  rt = f(x∗) − f(xt). The  cumulative regret RTafter T rounds is the sum of instantaneous regrets:  RT = �Tt=1 rt. A desirable asymptotic property of an algorithm is to be no-regret: limT →∞ RT /T = 0.Note that neither  rtnor  RTare ever revealed to the algorithm. Bounds on the average regret  RT /Ttranslate to convergence rates for GP optimization: the maximum maxt≤T f(xt) in the first T rounds is no further from  f(x∗) than the average.

2.1. Gaussian Processes and RKHS’s

Gaussian Processes. Some assumptions on f are required to guarantee no-regret. While rigid parametric assumptions such as linearity may not hold in practice, a certain degree of smoothness is often warranted. In our sensor network, temperature readings at closeby locations are highly correlated (see Figure 2(a)). We can enforce implicit properties like smoothness without relying on any parametric assumptions, modeling f as a sample from a Gaussian process (GP): a collection of dependent random variables, one for each x ∈ D, every finite subset of which is multivariate Gaussian distributed in an overall consistent way (Ras- mussen & Williams, 2006). A  GP(µ(x), k(x, x′)) is specified by its mean function  µ(x) = E[f(x)] and covariance (or kernel) function  k(x, x′) = E[(f(x) −µ(x))(f(x′) − µ(x′))]. For GPs not conditioned on data, we assume2that  µ ≡0. Moreover, we restrict k(x, x) ≤1,  x ∈ D, i.e., we assume bounded variance. By fixing the correlation behavior, the covariance function k encodes smoothness properties of sample functions f drawn from the GP. A range of commonly used kernel functions is given in Section 5.2.

In this work, GPs play multiple roles. First, some of our results hold when the unknown target function is a sample from a known GP distribution GP(0, k(x, x′)). Second, the Bayesian algorithm we analyze generally uses GP(0, k(x, x′)) as prior distribution over f. A major advantage of working with GPs is the existence of simple analytic formulae for mean and covariance of the posterior distribution, which allows easy implementation of algorithms. For a noisy sample  yT = [y1 . . . yT ]Tat points  AT = {x1, . . . , xT }, yt = f(xt)+ϵtwith  ϵt ∼ N(0, σ2) i.i.d. Gaussian noise, the posterior over f is a GP distribution again, with mean  µT (x), covariance  kT (x, x′) and variance  σ2T(x):

µT (x) = kT (x)T(KT + σ2I)−1yT ,(1) kT (x, x′) = k(x, x′) − kT (x)T(KT + σ2I)−1kT (x′),

image

where  kT (x) = [k(x1, x) . . . k(xT , x)]Tand  KTis the positive definite kernel matrix [k(x, x′)]x,x′∈AT.

RKHS. Instead of the Bayes case, where f is sampled from a GP prior, we also consider the more agnostic case where f has low “complexity” as measured under an RKHS norm (and distribution free assumptions on the noise process). The notion of reproducing kernel Hilbert spaces (RKHS, Wahba 1990) is intimately related to GPs and their covariance functions  k(x, x′). The RKHS  Hk(D) is a complete subspace of  L2(D) of nicely behaved functions, with an inner product  ⟨·, ·⟩kobeying the reproducing property: ⟨f, k(x, ·)⟩k = f(x) for all  f ∈ Hk(D). It is literally constructed by completing the set of mean functions µTfor all possible  T, {xt}, and  yT. The induced

RKHS norm  ∥f∥k =�⟨f, f⟩kmeasures smoothness of f w.r.t. k: in much the same way as  k1would generate smoother samples than  k2as GP covariance functions, ∥·∥k1assigns larger penalties than  ∥·∥k2. ⟨·, ·⟩kcan be extended to all of  L2(D), in which case  ∥f∥k < ∞iff f ∈ Hk(D). For most kernels discussed in Section 5.2, members of  Hk(D) can uniformly approximate any continuous function on any compact subset of D.

2.2. Information Gain & Experimental Design

One approach to maximizing f is to first choose points  xtso as to estimate the function globally well, then play the maximum point of our estimate. How can we learn about f as rapidly as possible? This question comes down to Bayesian Experimental Design (henceforth “ED”; see Chaloner & Verdinelli 1995), where the informativeness of a set of sampling points  A ⊂ Dabout f is measured by the information gain (c.f., Cover & Thomas 1991), which is the mutual information between f and observations  yA = f A+ϵAat these points:

image

quantifying the reduction in uncertainty about f from revealing  yA. Here,  f A = [f(x)]x∈Aand εA ∼ N(0, σ2I). For a Gaussian, H(N(µ, Σ)) =

2log  |2πeΣ|, so that in our setting I(yA; f) =I(yA; f A) = 12log  |I + σ−2KA|, where  KA =[k(x, x′)]x,x′∈A. While finding the information gain maximizer among  A ⊂ D, |A| ≤ Tis NP-hard (Ko et al., 1995), it can be approximated by an efficient greedy algorithm. If  F(A) = I(yA; f), this algorithm picks  xt = argmaxx∈D F(At−1∪{x}) in round t, which can be shown to be equivalent to

image

where  At−1 = {x1, . . . , xt−1}. Importantly, this simple algorithm is guaranteed to find a near-optimal solution: for the set  ATobtained after T rounds, we have that

image

at least a constant fraction of the optimal information gain value. This is because F(A) satisfies a diminishing returns property called submodularity (Krause & Guestrin, 2005), and the greedy approximation guarantee (5) holds for any submodular function (Nemhauser et al., 1978).

While sequentially optimizing Eq. 4 is a provably good way to explore f globally, it is not well suited for function optimization. For the latter, we only need to identify points x where f(x) is large, in order to concentrate sampling there as rapidly as possible, thus exploit our knowledge about maxima. In fact, the ED rule (4) does not even depend on observations  ytobtained along the way. Nevertheless, the maximum information gain after T rounds will play a prominent role in our regret bounds, forging an important connection between GP optimization and experimental design.

For sequential optimization, the ED rule (4) can be wasteful: it aims at decreasing uncertainty globally, not just where maxima might be. Another idea is to pick points as  xt = argmaxx∈D µt−1(x), maximizing the expected reward based on the posterior so far. However, this rule is too greedy too soon and tends to get stuck in shallow local optima. A combined strategy is to choose

image

where  βtare appropriate constants. This latter objective prefers both points x where f is uncertain (large σt−1(·)) and such where we expect to achieve high rewards (large  µt−1(·)): it implicitly negotiates the exploration–exploitation tradeoff. A natural interpretation of this sampling rule is that it greedily selects points x such that f(x) should be a reasonable upper bound on  f(x∗), since the argument in (6) is an upper quantile of the marginal posterior  P(f(x)|yt−1). We call this choice the Gaussian process upper confidence bound rule (GP-UCB), where  βtis specified depending on the context (see Section 4). Pseudocode for the GP-UCB algorithm is provided in Algorithm 1. Figure 2 illustrates two subsequent iterations, where GP-UCB both explores (Figure 2(b)) by sampling an input x with large  σ2t−1(x) and exploits (Figure 2(c)) by sampling x with large  µt−1(x).

The GP-UCB selection rule Eq. 6 is motivated by the UCB algorithm for the classical multi-armed bandit problem (Auer et al., 2002; Kocsis & Szepesv´ari, 2006). Among competing criteria for GP optimization (see Section 1), a variant of the GP-UCB rule has been demonstrated to be effective for this application (Dorard et al., 2009). To our knowledge, strong theoretical results of the kind provided for GP-UCB in this paper have not been given for any of these search heuristics. In Section 6, we show that in practice GP-UCB compares favorably with these alternatives.

If D is infinite, finding  xtin (6) may be hard: the upper confidence index is multimodal in general. However, global search heuristics are very effective in practice (Brochu et al., 2009). It is generally assumed

image

that evaluating f is more costly than maximizing the UCB index.

UCB algorithms (and GP optimization techniques in general) have been applied to a large number of problems in practice (Kocsis & Szepesv´ari, 2006; Pandey & Olston, 2007; Lizotte et al., 2007). Their performance is well characterized in both the finite arm setting and the linear optimization setting, but no convergence rates for GP optimization are known.

We now establish cumulative regret bounds for GP optimization, treating a number of different settings: f ∼GP(0, k(x, x′)) for finite  D, f ∼GP(0, k(x, x′)) for general compact D, and the agnostic case of arbitrary f with bounded RKHS norm.

GP optimization generalizes stochastic linear optimization, where a function f from a finite-dimensional linear space is optimized over. For the linear case, Dani et al. (2008) provide regret bounds that explicitly depend on the dimensionality3 d. GPs can be seen as random functions in some infinite-dimensional linear space, so their results do not apply in this case. This problem is circumvented in our regret bounds. The quantity governing them is the maximum information gain γTafter T rounds, defined as:

image

where I(yA; f A) = I(yA; f) is defined in (3). Recall that I(yA; f A) = 12log  |I + σ−2KA|, where  KA =[k(x, x′)]x,x′∈Ais the covariance matrix of  f A =[f(x)]x∈Aassociated with the samples A. Our regret bounds are of the form  O∗(√TβT γT), where  βTis the confidence parameter in Algorithm 1, while the bounds of Dani et al. (2008) are of the form  O∗(√TβT d) (d the dimensionality of the linear function space). Here and below, the  O∗notation is a variant of O, where log factors are suppressed. While our proofs – all provided in the Appendix – use techniques similar to those of Dani et al. (2008), we face a number of additional

image

Figure 2. (a) Example of temperature data collected by a network of 46 sensors at Intel Research Berkeley. (b,c) Two iterations of the GP-UCB algorithm. It samples points that are either uncertain (b) or have high posterior mean (c).

significant technical challenges. Besides avoiding the finite-dimensional analysis, we must handle confidence issues, which are more delicate for nonlinear random functions.

Importantly, note that the information gain is a prob- lem dependent quantity — properties of both the kernel and the input space will determine the growth of regret. In Section 5, we provide general methods for bounding  γT, either by efficient auxiliary computations or by direct expressions for specific kernels of interest. Our results match known lower bounds (up to log factors) in both the K-armed bandit and the d-dimensional linear optimization case.

Bounds for a GP Prior. For finite D, we obtain the following bound.

Theorem 1 Let δ ∈ (0,1) and βt =2 log(|D|t2π2/6δ). Running GP-UCB with  βtfor a sample f of a GP with mean function zero and covariance function  k(x, x′), we obtain a regret bound of  O∗(�TγTlog |D|) with high probability. Precisely,

image

where C1 = 8/log(1 +  σ−2).

The proof methodology follows Dani et al. (2007) in that we relate the regret to the growth of the log volume of the confidence ellipsoid — a novelty in our proof is showing how this growth is characterized by the information gain.

This theorem shows that, with high probability over samples from the GP, the cumulative regret is bounded in terms of the maximum information gain, forging a novel connection between GP optimization and experimental design. This link is of fundamental technical importance, allowing us to generalize Theorem 1 to infinite decision spaces. Moreover, the submodularity of I(yA; f A) allows us to derive sharp a priori bounds, depending on choice and parameterization of k (see Section 5). In the following theorem, we generalize our result to any compact and convex  D ⊂ Rdunder mild assumptions on the kernel function k.

Theorem 2 Let  D ⊂[0, r]dbe compact and convex, d ∈ N, r >0. Suppose that the kernel  k(x, x′) satisfies the following high probability bound on the derivatives of GP sample paths f: for some constants a, b > 0,

Pr {supx∈D |∂f/∂xj|> L} ≤ae−(L/b)2, j = 1, . . . , d.

image

βt = 2 log(t22π2/(3δ)) + 2d log�t2dbr�log(4da/δ)�.

Running the GP-UCB with  βtfor a sample f of a GP with mean function zero and covariance function k(x, x′), we obtain a regret bound of  O∗(√dTγT) with high probability. Precisely, with  C1 = 8/log(1 +  σ−2) we have

image

The main challenge in our proof (provided in the Appendix) is to lift the regret bound in terms of the confidence ellipsoid to general D. The smoothness assumption on  k(x, x′) disqualifies GPs with highly erratic sample paths. It holds for stationary kernels k(x, x′) = k(x − x′) which are four times differen-tiable (Theorem 5 of Ghosal & Roy (2006)), such as the Squared Exponential and Mat´ern kernels with  ν >2 (see Section 5.2), while it is violated for the OrnsteinUhlenbeck kernel (Mat´ern with  ν = 1/2; a stationary variant of the Wiener process). For the latter, sample paths f are nondifferentiable almost everywhere with probability one and come with independent increments. We conjecture that a result of the form of Theorem 2 does not hold in this case.

Bounds for Arbitrary f in the RKHS. Thus far, we have assumed that the target function f is sampled from a GP prior and that the noise is  N(0, σ2) with known variance  σ2. We now analyze GP-UCB in an agnostic setting, where f is an arbitrary function from the RKHS corresponding to kernel  k(x, x′). Moreover, we allow the noise variables  εtto be an arbitrary martingale difference sequence (meaning that E[εt | ε<t] = 0 for all t ∈ N), uniformly bounded by  σ. Note that we still run the same GP-UCB algorithm, whose prior and noise model are misspecified in this case. Our following result shows that GP-UCB attains sublinear regret even in the agnostic setting.

Theorem 3 Let  δ ∈(0, 1). Assume that the true underlying f lies in the RKHS  Hk(D) corresponding to the kernel  k(x, x′), and that the noise  εthas zero mean conditioned on the history and is bounded by  σalmost surely. In particular, assume  ∥f∥2k ≤ Band let  βt = 2B+ 300γtlog3(t/δ). Running GP-UCB with βt, prior  GP(0, k(x, x′)) and noise model  N(0, σ2), we obtain a regret bound of  O∗(√T(B√γT +γT)) with high probability (over the noise). Precisely,

image

where C1 = 8/log(1 +  σ−2).

Note that while our theorem implicitly assumes that GP-UCB has knowledge of an upper bound on  ∥f∥k, standard guess-and-doubling approaches suffice if no such bound is known a priori. Comparing Theorem 2 and Theorem 3, the latter holds uniformly over all functions f with  ∥f∥k < ∞, while the former is a probabilistic statement requiring knowledge of the GP that f is sampled from. In contrast, if  f ∼GP(0, k(x, x′)), then  ∥f∥k = ∞almost surely (Wahba, 1990): sample paths are rougher than RKHS functions. Neither Theorem 2 nor 3 encompasses the other.

Since the bounds developed in Section 4 depend on the information gain, the key remaining question is how to bound the quantity  γTfor practical classes of kernels.

5.1. Submodularity and Greedy Maximization

In order to bound  γT, we have to maximize the information gain  F(A) = I(yA; f) over all subsets  A ⊂ Dof size T: a combinatorial problem in general. However, as noted in Section 2, F(A) is a submodular function, which implies the performance guarantee (5) for maximizing F sequentially by the greedy ED rule (4). Dividing both sides of (5) by 1−1/e, we can upper-bound γTby (1  − 1/e)−1I(yAT ; f), where  ATis constructed by the greedy procedure. Thus, somewhat counterintuitively, instead of using submodularity to prove that F(AT) is near-optimal, we use it in order to show that γTis “near-greedy”. As noted in Section 2, the ED rule does not depend on observations  ytand can be run without evaluating f.

The importance of this greedy bound is twofold. First, it allows us to numerically compute highly problem-specific bounds on  γT, which can be plugged into our results in Section 4 to obtain high-probability bounds on  RT. This being a laborious procedure, one would prefer a priori bounds for  γTin practice which are simple analytical expressions of T and parameters of k. In this section, we sketch a general procedure for obtaining such expressions, instantiating them for a number of commonly used covariance functions, once more relying crucially on the greedy ED rule upper bound. Suppose that D is finite for now, and let  f = [f(x)]x∈D, KD = [k(x, x′)]x,x′∈D. Sampling f at  xt, we obtain  yt ∼ N(vTt f , σ2), where  vt ∈ R|D|is the indicator vector associated with  xt. We can upper-bound the greedy maximum once more, by relaxing this constraint to  ∥vt∥ = 1 in round tof the sequential method. For this relaxed greedy procedure, all  vtare leading eigenvectors of  KD, since successive covariance matrices of  P(f |yt−1) share their eigenbasis with  KD, while eigenvalues are damped according to how many times the corresponding eigenvector is selected. We can upper-bound the information gain by considering the worst-case allocation of T samples to the min{T, |D|} leading eigenvectors of  KD:

image

subject to �t mt = T, and spec(KD) = {ˆλ1 ≥ ˆλ2 ≥. . . }. We can split the sum into two parts in order to obtain a bound to leading order. The following Theorem captures this intuition:

Theorem 4 For any  T ∈ Nand any  T∗ = 1, . . . , T:

image

where  nT = �|D|t=1ˆλtand  B(T∗) = �|D|t=T∗+1ˆλt.

Therefore, if for some  T∗ = o(T) the first  T∗eigenvalues carry most of the total mass  nT, the information gain will be small. The more rapidly the spectrum of  KDdecays, the slower the growth of  γT. Figure 3 illustrates this intuition.

5.2. Bounds for Common Kernels

In this section we bound  γTfor a range of commonly used covariance functions: finite dimensional linear, Squared Exponential and Mat´ern kernels. Together with our results in Section 4, these imply sublinear regret bounds for GP-UCB in all cases.

image

Figure 3. Spectral decay (left) and information gain bound (right) for independent (diagonal), linear, squared exponential and Mat´ern kernels (ν = 2.5.) with equal trace.

Finite dimensional linear kernels have the form k(x, x′) = xT x′. GPs with this kernel correspond to random linear functions  f(x) = wT x, w ∼ N(0, I).

The Squared Exponential kernel is k(x, x′) =exp(−(2l2)−1∥x − x′∥2), la lengthscale parameter. Sample functions are differentiable to any order almost surely (Rasmussen & Williams, 2006).

The Mat´ern kernelis given by k(x, x′) =(21−ν/Γ(ν))rνBν(r), r = (√2ν/l)∥x − x′∥, where  νcontrols the smoothness of sample paths (the smaller, the rougher) and  Bνis a modified Bessel function. Note that as  ν → ∞, appropriately rescaled Mat´ern kernels converge to the Squared Exponential kernel.

Figure 4 shows random functions drawn from GP distributions with the above kernels.

Theorem 5 Let  D ⊂ Rdbe compact and convex,  d ∈N. Assume the kernel function satisfies  k(x, x′) ≤1.

1. Finite spectrum. For the d-dimensional Bayesian linear regression case:  γT = O�dlog  T�.

2. Exponential spectral decay. For the Squared Exponential kernel:  γT = O�(log  T)d+1�.

3. Power law spectral decay. For Mat´ern kernelswith ν > 1: γT = O�T d(d+1)/(2ν+d(d+1))(log  T)�.

A proof of Theorem 5 is given in the Appendix, , we only sketch the idea here. γTis bounded by Theorem 4 in terms the eigendecay of the kernel matrix KD. If D is infinite or very large, we can use the operator spectrum of  k(x, x′), which likewise decays rapidly. For the kernels of interest here, asymptotic expressions for the operator eigenvalues are given in Seeger et al. (2008), who derived bounds on the information gain for fixed and random designs (in contrast to the worst-case information gain considered here, which is substantially more challenging to bound). The main challenge in the proof is to ensure the existence of discretizations  DT ⊂ D, dense in the limit, for which tail sums  B(T∗)/nTin Theorem 4 are close to corresponding operator spectra tail sums.

Together with Theorems 2 and 3, this result guarantees sublinear regret of GP-UCB for any dimension (see Figure 1). For the Squared Exponential kernel, the dimension d appears as exponent of log T only, so that the regret grows at most as  O∗(√T(log  T)d+12) – the high degree of smoothness of the sample paths effectively combats the curse of dimensionality.

We compare GP-UCB with heuristics such as the Expected Improvement (EI) and Most Probable Improvement (MPI), and with naive methods which choose points of maximum mean or variance only, both on synthetic and real sensor network data.

For synthetic data, we sample random functions from a squared exponential kernel with lengthscale parameter 0.2. The sampling noise variance  σ2was set to 0.025 or 5% of the signal variance. Our decision set D = [0, 1] is uniformly discretized into 1000 points. We run each algorithm for  T = 1000 iterations with δ = 0.1, averaging over 30 trials (samples from the kernel). While the choice of  βtas recommended by Theorem 1 leads to competitive performance of GP-UCB, we find (using cross-validation) that the algorithm is improved by scaling  βtdown by a factor 5. Note that we did not optimize constants in our regret bounds.

Next, we use temperature data collected from 46 sensors deployed at Intel Research Berkeley over 5 days at 1 minute intervals, pertaining to the example in Section 2. We take the first two-thirds of the data set to compute the empirical covariance of the sensor readings, and use it as the kernel matrix. The functions f for optimization consist of one set of observations from all the sensors taken from the remaining third of the

image

Figure 4. Sample functions drawn from a GP with linear, squared exponential and Mat´ern kernels (ν = 2.5.)

image

Figure 5. Comparison of performance: GP-UCB and various heuristics on synthetic (a), and sensor network data (b, c).

data set, and the results (for  T = 46, σ2 = 0.5 or 5% noise,  δ = 0.1) were averaged over 2000 possible choices of the objective function.

Lastly, we take data from traffic sensors deployed along the highway I-880 South in California. The goal was to find the point of minimum speed in order to identify the most congested portion of the highway; we used traffic speed data for all working days from 6 AM to 11 AM for one month, from 357 sensors. We again use the covariance matrix from two-thirds of the data set as kernel matrix, and test on the other third. The results (for  T = 357, σ2 = 4.78 or 5% noise,  δ = 0.1) were averaged over 900 runs.

Figure 5 compares the mean average regret incurred by the different heuristics and the GP-UCB algorithm on synthetic and real data. For temperature data, the GP-UCB algorithm and EI heuristic clearly outperform the others, and do not exhibit significant difference between each other. On synthetic and traf-fic data MPI does equally well. In summary, GP-UCB performs at least on par with the existing approaches which are not equipped with regret bounds.

We prove the first sublinear regret bounds for GP optimization with commonly used kernels (see Figure 1), both for f sampled from a known GP and f of low RKHS norm. We analyze GP-UCB, an intuitive, Bayesian upper confidence bound based sampling rule. Our regret bounds crucially depend on the information gain due to sampling, establishing a novel connection between bandit optimization and experimental design. We bound the information gain in terms of the kernel spectrum, providing a general methodology for obtaining regret bounds with kernels of interest. Our experiments on real sensor network data indicate that GPUCB performs at least on par with competing criteria for GP optimization, for which no regret bounds are known at present. Our results provide an interesting step towards understanding exploration–exploitation tradeoffs with complex utility functions.

We thank Marcus Hutter for insightful comments on an earlier version of this paper. This research was partially supported by ONR grant N00014-09-1-1044, NSF grant CNS-0932392, a gift from Microsoft Corporation and the Excellence Initiative of the German research foundation (DFG).

Abernethy, J., Hazan, E., and Rakhlin, A. An efficient algorithm for linear bandit optimization, 2008. COLT.

Auer, P. Using confidence bounds for exploitationexploration trade-offs. JMLR, 3:397–422, 2002.

Auer, P., Cesa-Bianchi, N., and Fischer, P. Finite-time

analysis of the multiarmed bandit problem. Mach. Learn., 47(2-3):235–256, 2002.

Brochu, E., Cora, M., and de Freitas, N. A tutorial on Bayesian optimization of expensive cost functions, with application to active user modeling and hierarchical reinforcement learning. In TR-2009-23, UBC, 2009.

Bubeck, S., Munos, R., Stoltz, G., and Szepesv´ari, C. On- line optimization in X-armed bandits. In NIPS, 2008.

Chaloner, K. and Verdinelli, I. Bayesian experimental de- sign: A review. Stat. Sci., 10(3):273–304, 1995.

Cover, T. M. and Thomas, J. A. Elements of Information Theory. Wiley Interscience, 1991.

Dani, V., Hayes, T. P., and Kakade, S. The price of bandit information for online optimization. In NIPS, 2007.

Dani, V., Hayes, T. P., and Kakade, S. M. Stochastic linear optimization under bandit feedback. In COLT, 2008.

Dorard, L., Glowacka, D., and Shawe-Taylor, J. Gaussian process modelling of dependencies in multi-armed bandit problems. In Int. Symp. Op. Res., 2009.

Freedman, D. A. On tail probabilities for martingales. Ann. Prob., 3(1):100–118, 1975.

Ghosal, S. and Roy, A. Posterior consistency of Gaussian process prior for nonparametric binary regression. Ann. Stat., 34(5):2413–2429, 2006.

Gr¨unew¨alder, S., Audibert, J-Y., Opper, M., and Shawe- Taylor, J. Regret bounds for gaussian process bandit problems. In AISTATS, 2010.

Huang, D., Allen, T. T., Notz, W. I., and Zeng, N. Global optimization of stochastic black-box systems via sequential kriging meta-models. J Glob. Opt., 34:441–466, 2006.

Jones, D. R., Schonlau, M., and Welch, W. J. Efficient global optimization of expensive black-box functions. J Glob. Opti., 13:455–492, 1998.

Kleinberg, R., Slivkins, A., and Upfal, E. Multi-armed bandits in metric spaces. In STOC, pp. 681–690, 2008.

Ko, C., Lee, J., and Queyranne, M. An exact algorithm for maximum entropy sampling. Ops Res, 43(4):684–691, 1995.

Kocsis, L. and Szepesv´ari, C. Bandit based monte-carlo planning. In ECML, 2006.

Krause, A. and Guestrin, C. Near-optimal nonmyopic value of information in graphical models. In UAI, 2005.

Lizotte, D., Wang, T., Bowling, M., and Schuurmans, D. Automatic gait optimization with Gaussian process regression. In IJCAI, pp. 944–949, 2007.

McDiarmid, C. Concentration. In Probabilistiic Methods for Algorithmic Discrete Mathematics. Springer, 1998.

Mockus, J. Bayesian Approach to Global Optimization. Kluwer Academic Publishers, 1989.

Mockus, J., Tiesis, V., and Zilinskas, A. Toward Global Optimization, volume 2, chapter Bayesian Methods for Seeking the Extremum, pp. 117–128. 1978.

Nemhauser, G., Wolsey, L., and Fisher, M. An analysis of the approximations for maximizing submodular set functions. Math. Prog., 14:265–294, 1978.

Pandey, S. and Olston, C. Handling advertisements of un- known quality in search advertising. In NIPS. 2007.

Rasmussen, C. E. and Williams, C. K. I. Gaussian Processes for Machine Learning. MIT Press, 2006.

Robbins, H. Some aspects of the sequential design of ex- periments. Bul. Am. Math. Soc., 58:527–535, 1952.

Rusmevichientong, P. and Tsitsiklis, J. N. Linearly param- eterized bandits. abs/0812.3465, 2008.

Seeger, M. W., Kakade, S. M., and Foster, D. P. Infor- mation consistency of nonparametric Gaussian process methods. IEEE Tr. Inf. Theo., 54(5):2376–2382, 2008.

Shawe-Taylor, J., Williams, C., Cristianini, N., and Kan- dola, J. On the eigenspectrum of the Gram matrix and the generalization error of kernel-PCA. IEEE Trans. Inf. Theo., 51(7):2510–2522, 2005.

Srinivas, N., Krause, A., Kakade, S., and Seeger, M. Gaus- sian process optimization in the bandit setting: No regret and experimental design. In ICML, 2010.

Stein, M. Interpolation of Spatial Data: Some Theory for Kriging. Springer, 1999.

Vazquez, E. and Bect, J. Convergence properties of the expected improvement algorithm, 2007.

Wahba, G. Spline Models for Observational Data. SIAM, 1990.

In this section, we provide details for the proofs of Theorem 1 and Theorem 2. In both cases, the strategy is to show that  |f(x) − µt−1(x)| ≤ β1/2t σt−1(x) for all t ∈ Nand all  x ∈ D, or in the infinite case, all x in a discretization of D which becomes dense as t gets large.

A.1. Finite Decision Set

We begin with the finite case,  |D| < ∞.

Lemma 5.1 Pick δ ∈ (0,1) and set βt =2 log(|D|πt/δ), where �t≥1 π−1t = 1, πt >0. Then,

image

holds with probability  ≥ 1 − δ.

Proof Fix  t ≥1 and  x ∈ D. Conditioned on  yt−1 =(y1, . . . , yt−1), {x1, . . . , xt−1}are deterministic, and f(x) ∼ N(µt−1(x), σ2t−1(x)). Now, if  r ∼ N(0,1), then

image

for c > 0, since  e−c(r−c) ≤1 for  r ≥ c. Therefore, Pr{|f(x) − µt−1(x)| > β1/2t σt−1(x)} ≤ e−βt/2, using r = (f(x)−µt−1(x))/σt−1(x) and  c = β1/2t. Applying the union bound,

image

holds with probability  ≥ 1 − |D|e−βt/2. Choosing |D|e−βt/2 = δ/πtand using the union bound for t ∈ N, the statement holds. For example, we can use πt = π2t2/6.

Lemma 5.2 Fix  t ≥1. If  |f(x) − µt−1(x)| ≤β1/2t σt−1(x) for all  x ∈ D, then the regret  rtis bounded by 2β1/2t σt−1(xt).

Proof By definition of  xt: µt−1(xt)+β1/2t σt−1(xt) ≥µt−1(x∗) +  β1/2t σt−1(x∗) ≥ f(x∗). Therefore,

rt = f(x∗) − f(xt) ≤ β1/2t σt−1(xt) +  µt−1(xt) − f(xt)

image

Lemma 5.3 The information gain for the points selected can be expressed in terms of the predictive variances. If  f T = (f(xt)) ∈ RT:

image

Proof Recall that I(yT ; f T ) = H(yT ) −(1/2) log  |2πeσ2I|. Now, H(yT ) = H(yT −1) + H(yT |yT −1) = H(yT −1) + log(2πe(σ2+  σ2t−1(xT)))/2. Here, we use that  x1, . . . , xTare deterministic conditioned on  yT −1, and that the conditional variance σ2T −1(xT) does not depend on  yT −1. The result fol- lows by induction.

Lemma 5.4 Pick  δ ∈(0, 1) and let  βtbe defined as in Lemma 5.1. Then, the following holds with probability ≥ 1 − δ:

image

where C1 := 8/log(1 +  σ−2) ≥ 8σ2.

Proof By Lemma 5.1 and Lemma 5.2, we have that {r2t ≤ 4βtσ2t−1(xt) ∀t ≥ 1}with probability  ≥ 1 − δ. Now,  βtis nondecreasing, so that

image

with C2 = σ−2/log(1 +  σ−2) ≥1, since s2 ≤ C2log(1 +  s2) for s ∈ [0, σ−2], and σ−2σ2t−1(xt) ≤ σ−2k(xt, xt) ≤ σ−2. Noting that C1 = 8σ2C2, the result follows by plugging in the representation of Lemma 5.3.

Finally, Theorem 1 is a simple consequence of Lemma 5.4, since  R2T ≤ T �Tt=1 r2tby the Cauchy- Schwarz inequality.

A.2. General Decision Set

Theorem 2 extends the statement of Theorem 1 to the general case of  D ⊂ Rdcompact. We cannot expect this generalization to work without any assumptions on the kernel  k(x, x′). For example, if k(x, x′) = e−∥x−x′∥(Ornstein-Uhlenbeck), while sample paths f are a.s. continuous, they are still very erratic: f is a.s. nondifferentiable almost everywhere, and the process comes with independent increments, a stationary variant of Brownian motion. The additional assumption on k in Theorem 2 is rather mild and is satisfied by several common kernels, as discussed in Section 4.

Recall that the finite case proof is based on Lemma 5.1 paving the way for Lemma 5.2. However, Lemma 5.1 does not hold for infinite D. First, let us observe that we have confidence on all decisions actually chosen.

Lemma 5.5 Pick  δ ∈(0, 1) and set  βt = 2 log(πt/δ), where �t≥1 π−1t = 1, πt >0. Then,

image

holds with probability  ≥ 1 − δ.

Proof Fix  t ≥1 and  x ∈ D. Conditioned on yt−1 = (y1, . . . , yt−1), {x1, . . . , xt−1}are deterministic, and  f(x) ∼ N(µt−1(x), σ2t−1(x)). As before, Pr{|f(xt) − µt−1(xt)| > β1/2t σt−1(xt)} ≤ e−βt/2. Since  e−βt/2 = δ/πtand using the union bound for t ∈ N, the statement holds.

Purely for the sake of analysis, we use a set of discretizations  Dt ⊂ D, where  Dtwill be used at time t in the analysis. Essentially, we use this to obtain a valid confidence interval on  x∗. The following lemma provides a confidence bound for these subsets.

Lemma 5.6 Pick δ ∈ (0,1) and set βt =2 log(|Dt|πt/δ), where �t≥1 π−1t = 1, πt >0. Then,

image

holds with probability  ≥ 1 − δ.

Proof The proof is identical to that in Lemma 5.1, except now we use  Dtat each timestep.

Now by assumption and the union bound, we have that

image

which implies that, with probability greater than 1  −dae−L2/b2, we have that

image

This allows us to obtain confidence on  x⋆as follows.

Now let us choose a discretization  Dtof size (τt)dso that for all  x ∈ Dt

image

where [x]tdenotes the closest point in  Dtto x. A suf-ficient discretization has each coordinate with  τtuniformly spaced points.

image

�t≥1 π−1t = 1, πt >0. Let  τt = dt2br�log(2da/δ) Let [x∗]tdenotes the closest point in  Dtto  x∗. Hence, Then,

|f(x∗) − µt−1([x∗]t)| ≤ β1/2t σt−1([x∗]t) + 1t2 ∀t ≥1

holds with probability  ≥ 1 − δ.

Proof Using (9), we have that with probability greater than 1  − δ/2,

image

Hence,

image

Now by choosing  τt = dt2br�log(2da/δ), we have that

image

This implies that  |Dt| = (dt2br�log(2da/δ))d. Using δ/2 in Lemma 5.6, we can apply the confidence bound to [x∗]t(as this lives in  Dt) to obtain the result.

Now we are able to bound the regret.

Lemma 5.8 Pick δ ∈ (0,1) and set βt =2 log(4πt/δ) + 4d log(dtbr�log(4da/δ)), where �t≥1 π−1t = 1, πt >0. Then, with probability greater than 1  − δ, for all  t ∈ N, the regret is bounded as follows:

image

Proof We use  δ/2 in both Lemma 5.5 and Lemma 5.7, so that these events hold with probability greater than 1  − δ. Note that the specification of  βtin the above lemma is greater than the specification used in Lemma 5.5 (with  δ/2), so this choice is valid.

By definition of  xt: µt−1(xt) +  β1/2t σt−1(xt) ≥µt−1([x∗]t)+β1/2t σt−1([x∗]t). Also, by Lemma 5.7, we have that  µt−1([x∗]t)+β1/2t σt−1([x∗]t)+1/t2 ≥ f(x∗), which implies  µt−1(xt)+β1/2t σt−1(xt) ≥ f(x∗)−1/t2. Therefore,

image

which completes the proof.

Now we are ready to complete the proof of Theorem 2. As shown in the proof of Lemma 5.4, we have that with probability greater than 1  − δ,

image

so that by Cauchy-Schwarz:

image

Hence,

image

(since �1/t2 = π2/6). Theorem 2 now follows.

image

states that if derivatives up to fourth order exists for (x, x′) �→ k(x, x′), then f is almost surely continuously differentiable, with  ∂f/(∂xj) distributed as Gaussian processes again. Moreover, there are constants  a, bj >0 such that

image

Picking  L = [log(da2/δ)/minj bj]1/2, we have that ae−bjL2 ≤ δ/(2d) for all j = 1, . . . , d, so that for K1 = d1/2L, by the mean value theorem, we have Pr{|f(x)−f(x′)| ≤ K1∥x−x′∥ ∀ x, x′ ∈ D} ≥ 1−δ/2.

Also, note that  K1 = O((log  δ−1)1/2).

This statement is about the joint distribution of  f(·) and its partial derivatives w.r.t. each component. For a certain event in this sample space, all  ∂f/(∂xj) exist, are continuous, and the complement of (10) holds for all j. Theorem 5 of Ghosal & Roy (2006), together with the union bound, implies that this event has probability  ≥ 1 − δ/2. Derivatives up to fourth order exist for the Gaussian covariance function, and for Mat´ern kernels with  ν >2 (Stein, 1999).

In this section, we detail a proof of Theorem 3. Recall that in this setting, we do not know the generator of the target function f, but only a bound on its RKHS norm  ∥f∥k.

Recall the posterior mean function  µT (·) and posterior covariance function  kT (·, ·) from Section 2, conditioned on data (xt, yt), t = 1, . . . , T. It is easy to see that the RKHS norm corresponding to  kTis given by

image

This implies that  Hk(D) = HkT (D) for any T, while the RKHS inner products are different:  ∥f∥kT ≥ ∥f∥k. Since  ⟨f(·), kT (·, x)⟩kT = f(x) for any  f ∈ HkT (D) by the reproducing property, then

image

by the Cauchy-Schwarz inequality.

Compared to our other results, Theorem 3 is an agnostic statement, in that the assumptions the Bayesian UCB algorithm bases its predictions on differ from how f and data  ytare generated. First, f is not drawn from a GP, but can be an arbitrary function from  Hk(D). Second, while the UCB method assumes that the noise  εt = yt − f(xt) is drawn independently from  N(0, σ2), the true sequence of noise variables  εtcan be a uniformly bounded martingale difference sequence:  εt ≤ σfor all  t ∈ N. All we have to do in order to lift the proof of Theorem 1 to the agnostic setting is to establish an analogue to Lemma 5.1, by way of the following concentration result.

Theorem 6 Let  δ ∈(0, 1). Assume the noise variables  εtare uniformly bounded by  σ. Define:

image

Then

Pr�∀T, ∀x ∈ D, |µT (x) f(x)| ≤  β1/2T +1σT (x)�1δ.

B.1. Concentration of Martingales

In our analysis, we use the following Bernstein-type concentration inequality for martingale differences, due to Freedman (1975) (see also Theorem 3.15 of Mc- Diarmid 1998).

Theorem 7 (Freedman) Suppose  X1, . . . , XTis a martingale difference sequence, and b is an uniform upper bound on the steps  Xi. Let V denote the sum of conditional variances,

image

Then, for every a, v > 0,

image

B.2. Proof of Theorem 6

We will show that:

image

Theorem 6 then follows from (11). Recall that  εt =yt − f(xt). We will analyze the quantity  ZT =∥µT − f∥2kT, measuring the error of  µTas approxi- mation to f under the RKHS norm of  HkT (D). The following lemma provides the connection with the information gain. This lemma is important since our concentration argument is an inductive argument — roughly speaking, we condition on getting concentration in the past, in order to achieve good concentration in the future.

image

�Tt=1min{σ−2σ2t−1(xt), α} ≤ 2αlog(1 +  α)γT , α > 0.Proof We have that min{r, α} ≤ (α/log(1 + α)) log(1+r). The statement follows from Lemma 5.3.

image

The next lemma bounds the growth of  ZT. It is for- mulated in terms of normalized quantities:  �εt = εt/σ, �f = f/σ, �µt = µt/σ, �σt = σt/σ. Also, to ease notation, we will use  µt−1, σt−1as shorthand for  µt−1(xt), σt−1(xt).

Lemma 7.2  For all T ∈ N,

image

Proof If αt = (Kt + σ2I)−1yt, then  µt(x) =αTt kt(x). Then, ⟨µT , f⟩k = f TT αT, ∥µT ∥2k =yTT αT − σ2∥αT ∥2. Moreover, for  t ≤ T, µT (xt) =δTt KT(KT+  σ2I)−1yT = yt − σ2αt. Since  ZT =∥µT −f∥k +σ−2 �t≤T(µT (xt)−f(xt))2, we have that

image

Now,  −yTT(KT +σ2I)−1yT .= 2 log P(yT), where “  .=”means that we drop determinant terms, thus concentrate on quadratic functions. Since log  P(yT ) =�tlog  P(yt|y<t) = �tlog  N(yt|µt−1(xt), σ2t−1(xt) + σ2), we have that

image

with  R = �t(µt−1 − f(xt))2/(σ2 + σ2t−1) ≥0. Dropping  −Rand changing to normalized quantities concludes the proof.

We now define a useful martingale difference sequence. First, it is convenient to define an “escape event”  ETas:

image

where I{·}is the indicator function. Define the random variables  Mtby

image

Now, since  �εtis a martingale difference sequence with respect to the histories  H<tand  Mt/�εtis deterministic given  H<t, Mtis a martingale difference sequence as well. Next, we show that with high probability, the associated martingale �Tt=1 Mtdoes not grow too large.

Lemma 7.3 Given  δ ∈(0, 1) and  βtas defined in in Theorem 6, we have that

image

The proof is given below in Section B.3. Equipped with this lemma, we can prove Theorem 6.

Proof [of Theorem 6] It suffices to show that the high-probability event described in Lemma 7.3 is contained in the support of  ETfor every T. We prove the latter by induction on T.

By Lemma 7.2 and the definition of  β1, we know that Z0 ≤ ∥f∥k ≤ β1. Hence  E0 = 1 always. Now supposethe high-probability event of Lemma 7.3 holds, in particular �Tt=1 Mt ≤ βT +1/2. For the inductive hypoth- esis, assume  ET −1 = 1. Using this and Lemma7.2:

image

The equality in the second step uses the inductive hypothesis. Thus we have shown  ET = 1, completingthe induction.

B.3. Concentration

What remains to be shown is Lemma 7.3. While the step sizes  |Mt|are uniformly bounded, a standard application of the Hoeffding-Azuma inequality leads to a bound of  T 3/4, too large for our purpose. We use the more specific Theorem 7 instead, which requires to control the conditional variances rather than the marginal variances which can be much larger.

Proof [of Lemma 7.3] Let us first obtain upper bounds

on the step sizes of our martingale.

image

where the first inequality follows from the definition of  Et. Moreover, r/(1 +  r2) ≤min{r, 1/2} for  r ≥0. Therefore,  |Mt| ≤ β1/2T, since  |�εt| ≤1 and  βtin nondecreasing. Next, we bound the sum of the conditional variances of the martingale:

image

In the last line, we used Lemma 7.1 with  α = 1/4, noting that 8α/log(1 +  α) ≤9. Since we have established that the sum of conditional variances,  VT, is always bounded by 9βT γT, we can apply Theorem 7 with parameters  a = βT +1/2, b = β1/2T +1and  v = 9βT γTto get

image

Note that our choice of  βT +1satisfies:

image

Therefore, the previous probability is bounded by δ/T 2, whereas the last inequality follows from the definition of  βT +1. With a final application of the union

bound:

image

completing the proof of Lemma 7.3.

In this section, we show how to bound  γT, the maximum information gain after T rounds, for compact D ⊂ Rd(assumptions of Theorem 2) and several commonly used covariance functions. In this section, we assume4that  k(x, x) = 1 for all x ∈ D.

The plan of attack is as follows. First, we note that the argument of  γT, I(yA; f A) is a submodular function, so  γTcan be bounded by the value obtained by greedy maximization. Next, we use a discretization  DT ⊂ Dwith  nT = |DT | = T τwith nearest neighbour distance o(1), consider the kernel matrix  KDT ∈ RnT ×nT, and bound  γTby an expression involving the eigenvalues {ˆλt}of this matrix, which is done by a further relaxation of the greedy procedure. Finally, we bound this empirical expression in terms of the kernel operator eigenvalues of k w.r.t. the uniform distribution on D. Asymptotic expressions for the latter are reviewed in Seeger et al. (2008), which we plug in to obtain our results. A key step in this argument is to ensure the existence of a discretization  DT, for which tails of the empirical spectrum can be bounded by tails of the process spectrum. We will invoke the probabilistic method for that.

C.1. Greedy Maximization and Discretization

In this section, we fix  T ∈ Nand assume the existence of a discretization  DT ⊂ D, nT = |DT |on the order of  T τ, such that:

image

We come back to the choice of  DTbelow. We restrict the information gain to subsets  A ⊂ DT:

image

Of course, ˜γT ≤ γT, but we can bound the slack.

Lemma 7.4 Under the assumptions of Theorem 2, the information gain  FT ({xt}) = (1/2) log |I + σ−2K{xt}|is uniformly Lipschitz-continuous in each component  xt ∈ D.

Proof The assumptions of Theorem 2 imply that the kernel  K(x, x′) is continuously differentiable. The result follows from the fact that  FT ({xt}) is continuously differentiable in the kernel matrix  K{xt}.

image

Lemma 7.5 Let  DTbe a discretization of D such that (13) holds. Under the assumptions of Theorem 2, we have that

image

Proof Fix  T ∈ N, and let  A = {x1, . . . , xT }be a maximizer for  γT. Consider neighbours [xt]T ∈ DTaccording to (13), [A]T = {[xt]T }. Then,

0  ≤ γT −˜γT ≤ γT −I(y[A]T ; f [A]T ) = FT (A)−FT([A]T ),

where  FT ({xt}) = (1/2) log  |I + σ−2K{xt}|. By Lemma 7.4,  FTis uniformly Lipschitz-continuous in each component, so that  |γT −I(y[A]T ; f [A]T )| =O(T maxt ∥xt − [xt]T ∥) = O(T 1−τ/d) by (13) and the mean value theorem.

We concentrate on ˜γTin the sequel. Let  KDT =[k(x, x′)]x,x′∈DTbe the kernel matrix over the entire  DT, and  KDT = UˆΛU Tits eigendecomposition, with ˆλ1 ≥ˆλ2 ≥ · · · ≥0 and  U = [u1 u2 . . .] orthonormal. Here, if  T > nT, define ˆλt = 0 fort = nT+ 1, . . . , T. Information gain maximization over a finite  DTcan be described in terms of a simple linear-Gaussian model over the unknown  f ∈ RnT, with prior  P(f ) = N(0, KDT) and likelihood potentials  P(yt|f ) = N(vTt f , σ2) with unit-norm features, ∥vt∥ = 1. With the following lemma, we upper-bound˜γTby way of two relaxations.

Lemma 7.6 For any  T ≥1, we have that

image

subject to  mt ∈ N, �t mT = T, where ˆλ1 ≥ ˆλ2 ≥ . . .is the spectrum of the kernel matrix  KDT. Here, if T > nT, then  mt = 0for  t > nT.

Proof As shown by Krause & Guestrin (2005), the function  F(A) = I(yA; f) is submodular. In the particular case considered here, this can be seen as follows: F(A) = H(yA) −H(yA | f), where the entropy H(yA) is a (not-necessarily monotonic) submodular function in A, and since the noise is conditionally independent given f , H(yA | f) is an additive (modular) function in A. Subtracting a modular function preserves submodularity, thus F(A) is submodular. Furthermore, the information gain is monotonic in A (i.e.,  F(A) ≤ F(B) whenever A ⊆ B) (Cover & Thomas, 1991). Thus, we can apply the result of Nemhauser et al. (1978)5which guarantees that ˜γTis upper-bounded by 1/(1 − 1/e) times the value the greedy maximization algorithm attains. The latter chooses features of the form vt = δxt = [I{x=xt}] in each round,  xt ∈ DT. We upper-bound the greedy maximum once more by relaxing these constraints to  ∥vt∥ = 1 only.In the remainder of the proof, we concentrate on this relaxed greedy procedure. Suppose that up to round t, it chose v1, . . . , vt−1. The posterior  P(f |yt−1) has inverse covariance matrix  Σ−1t−1 = K−1DT+  σ−2V t−1V Tt−1, V t−1 = [v1 . . . vt−1], and the greedy procedure selects v so to maximize the variance  vT Σt−1v: the eigenvector corresponding to  Σt−1’s largest eigenvalue (by the Rayleigh-Ritz theorem). Since  Σ0 = KDT, then  v1 = u1. Moreover, if all  vt′, t′ < t, have been chosen among U ’s columns, then by the inverse covariance expression just given,  KDTand  Σt−1have the same eigenvectors, so that  vtis a column of U as well. For example, if  vt = uj, then comparing  Σt−1and  Σt, all eigenvalues other than the j-th remain the same, while the latter is shrunk. Therefore, after T rounds of the relaxed greedy procedure: vt ∈ {u1, . . . , umin{T,nT }}, t = 1, . . . , T: at most the leading T eigenvectors of  KDTcan have been selected (possibly multiple times). If  mtdenotes the number that the t-th column of U has been selected, we obtain the theorem statement by a final bounding step.

C.2. From Empirical to Process Eigenvalues

The final step will be to relate the empirical spec- trum  {ˆλt}to the kernel operator spectrum. Since log(1 +  σ−2mtˆλt) ≤ σ−2mtˆλtin Theorem 7.6, we will mainly be interested in relating the tail sums of the spectra. Let  µ(x) = V(D)−1I{x∈D}be the uniform distribution on  D, V(D) =�x∈D dx, and assume that k is continuous. Note that�k(x, x)µ(x) dx = 1 byour assumption k(x, x) = 1, so that k is HilbertSchmidt on  L2(µ). Then, Mercer’s theorem (Wahba, 1990) states that the corresponding kernel operator has a discrete eigenspectrum  {(λs, φs(·))}, and

image

where  λ1 ≥ λ2 ≥ · · · ≥0, and  Eµ[φs(x)φt(x)] =δs,t. Moreover, �s≥1 λ2s < ∞, and the expansion of k converges absolutely and uniformly on  D ×D. Note that �s≥1 λs = �s≥1 λs Eµ[φs(x)2] =�K(x, x)µ(x) dx = 1. In order to proceed from The-orem 7.6, we have to pick a discretization  DTfor which (13) holds, and for which �t>T∗ˆλtis not much larger than �t>T∗ λt. With the following lemma, we deter- mine sizes  nTfor which such discretizations exist.

Lemma 7.7 Fix  T ∈ N, δ >0 and  ε >0. There exists a discretization  DT ⊂ Dof size

nT = V(D)(ε/√d)−d[log(1/δ)+dlog(√d/ε)+log V(D)]

which fulfils the following requirements:

 ε-denseness: For any  x ∈ D, there exists [x]T ∈DTsuch that  ∥x − [x]T ∥ ≤ ε.

If spec(KDT ) = {ˆλ1 ≥ˆλ2 ≥. . . }, then for any T∗ = 1, . . . , nT :

image

Proof First, if we draw  nTsamples ˜xj ∼ µ(x) independently at random, then  DT = {˜xj}is  ε-dense with probability  ≥ 1 − δ. Namely, cover D with N = V(D)(ε/√d)−dhypercubes of sidelength  ε/√d, within which the maximum Euclidean distance is  ε. The probability of not hitting at least one cell is upper-bounded by  N(1 − 1/N)nT. Since log(1  − 1/N) ≤−1/N, this is upper-bounded by  δif  nT ≥ Nlog(N/δ).

Now, let  S = n−1T �T∗t=1ˆλt. Shawe-Taylor et al. (2005) show that  E[S] ≥ �T∗t=1 λt. If C is the event  {DTis  ε−dense }, then Pr(C) ≥ 1 − δ. Since S ≤ n−1TtrKDT = 1 in any case, we have thatE[S|C] ≥ E[S] −Pr(Cc) ≥ �T∗t=1 λt − δ. By the probabilistic method, there must exist some  DTfor which C and the latter inequality holds.

The following lemma, the equivalent of Theorem 4 in the context here, is a direct consequence of Lemma 7.6.

Lemma 7.8 Let  DTbe some discretization of D,

nT = |DT |. Then, for any T∗ = 1, . . . , min{T, nT }:

image

Proof We split the right hand side in Lemma 7.6 at  t = T∗. Let  r = �t≤T∗ mt. For  t ≤ T∗: log(1 +  mtˆλt/σ2) ≤log(rnT /σ2), since ˆλt ≤ nT. For t > T∗: log(1+mtˆλt/σ2) ≤ mtˆλt/σ2 ≤ (T−r)ˆλt/σ2.

The following theorem describes our “recipe” for obtaining bounds on  γTfor a particular kernel k, given that tail bounds on  Bk(T∗) = �s>T∗ λsare known.

Theorem 8 Suppose that  D ⊂ Rdis compact, and k(x, x′) is a covariance function for which the additional assumption of Theorem 2 holds. Moreover, let  Bk(T∗) = �s>T∗ λs, where  {λs}is the operator spectrum of k with respect to the uniform distribution over D. Pick  τ >0, and let  nT = C4T τ(log T) with C4 = 2V(D)(2τ+ 1). Then, the following bound holds true:

image

for any T∗ ∈ {1, . . . , nT }.

Proof Let  ε = d1/2T −τ/dand  δ = T −(τ+1). Lemma 7.7 provides the existence of a discretization DTof size nTwhich is ε-dense, and for which n−1T �T∗t=1ˆλt ≥ �T∗t=1 λt − δ. Since n−1T �nTt=1ˆλt = 1 = �t≥1 λt, then �t>T∗ˆλt ≤ Bk(T∗) +  δ. The statement follows by using Lemma 7.8 with these bounds, and finally employing Lemma 7.5.

C.3. Proof of Theorem 5

In this section, we instantiate Theorem 8 in order to obtain bounds on  γTfor Squared Exponential and Mat´ern kernels, results which are summarized in Theorem 5.

image

For the Squared Exponential kernel  k, Bk(T∗) is given by Seeger et al. (2008). While  µ(x) was Gaussian there, the same decay rate holds for  λsw.r.t. uniform µ(x), while constants might change. In hindsight, it turns out that  τ = dis the optimal choice for the discretization size, rendering the second term in Theorem 5 to be O(1), which is subdominant and will be neglected in the sequel. We have that  λs ≤ cBs1/dwith B < 1. Following their analysis,

image

where  α = −log  B, β = αT 1/d∗. Therefore,  Bk(T∗) =O(e−β βd−1), β = αT 1/d∗.

We have to pick  T∗such that  e−βis not much larger than (TnT )−1. Suppose that  T∗ = [log(TnT )/α]d, so that  e−β = (TnT )−1, β = log(TnT). The bound becomes

image

with  nT = C4T d(log T). The first part dominates, so that r = T and  γT = O([log(T d+1(log T))]d+1) =O((log  T)d+1). This should be compared with E[I(yT ; f T )] = O((log  T)d+1) given by Seeger et al. (2008), where the  xtare drawn independently from a Gaussian base distribution. At least restricted to a compact set D, we obtain the same expression to leading order for max{xt}I(yT ; f T).

image

For Mat´ern kernels k with roughness parameter  ν, Bk(T∗) is given by Seeger et al. (2008) for the uniform base distribution  µ(x) on D. Namely,  λs ≤cs−(2ν+d)/dfor almost all  s ∈ N, and  Bk(T∗) =O(T 1−(2ν+d)/d∗). To match terms in the ˜γTbound, we choose  T∗ = (TnT )d/(2ν+d)(log(TnT))κ(κchosen below), so that the bound becomes

image

with  nT = C4T τ(log T). For  κ = −d/(2ν + d), we obtain that the maximum over r is  O(T∗log(TnT )) =O(T (τ+1)d/(2ν+d)(log T)). Finally, we choose  τ =2νd/(2ν+d(d+1)) to match this term with  O(T 1−τ/d). Plugging this in, we have  γT = O(T 1−2η(log T)), η = ν2ν+d(d+1). Together with Theorem 2 (for  ν >2), we have that  RT = O∗(T 1−η) (suppressing log factors): for any  ν >2 and any dimension d, the GPUCB algorithm is guaranteed to be no-regret in this case with arbitrarily high probability.

How does this bound compare to the bound on E[I(yT ; f T)] given by Seeger et al. (2008)? Here,  γT =O(T d(d+1)/(2ν+d(d+1))(log T)), while E[I(yT ; f T )] =O(T d/(2ν+d)(log  T)2ν/(2ν+d)).

image

For linear kernels  k(x, x′) = xT x′, x ∈ Rdwith  ∥x∥ ≤1, we can bound  γTdirectly. Let  XT = [x1 . . . , xT ] ∈Rd×Twith all  ∥xt∥ ≤1. Now,

log  |I + σ−2XTT XT | = log |I + σ−2XT XTT |≤log  |I + σ−2D|

with  D = diag diag−1(XT XTT), by Hadamard’s in- equality. The largest eigenvalue ˆλ1of  XT XTTis O(T), so that

image

and  γT = O(dlog T).


Designed for Accessibility and to further Open Science