1 Introduction
Historically, the fields of statistical inference and stochastic optimization have often developed their own specific methods and approaches. Recently, however, there has been a growing interest in applying inferencebased methods to optimization problems and vice versa [1, 2, 3, 4]. Here we consider stochastic optimization problems where we observe noisecontaminated values from an unknown nonlinear function and we want to find the input that maximizes the expected value of this function.
The problem statement is as follows. Let be a metric space. Consider a stochastic function mapping a test point to real values characterized by the conditional pdf . Consider the mean function
(1) 
The goal consists in modeling the optimal test point
(2) 
Classic approaches to solve this problem are often based on stochastic approximation methods [5]. Within the context of statistical inference, Bayesian optimization methods have been developed where a prior distribution over the space of functions is assumed and uncertainty is tracked during the entire optimization process [6, 7]. In particular, nonparametric Bayesian approaches such as Gaussian Processes have been applied for derivativefree optimization [8, 9], also within the context of the continuumarmed bandit problem [10]. Typically, these Bayesian approaches aim to explicitly represent the unknown objective function of (1) by entertaining a posterior distribution over the space of objective functions. In contrast, we aim to model directly the distribution of the maximum of (2) conditioned on observations.
2 Description of the Model
Our model is intuitively straightforward and easy to implement^{1}^{1}1Implementations can be downloaded from http://www.adaptiveagents.org/argmaxprior. Let be an estimate of the mean constructed from data (Figure 1a, left). This estimate can easily be converted into a posterior pdf over the location of the maximum by first multiplying it with a precision parameter and then taking the normalized exponential (Figure 1a, right)
In this transformation, the precision parameter controls the certainty we have over our estimate of the maximizing argument: expresses almost no certainty, while expresses certainty. The rationale for the precision is: the more distinct inputs we test, the higher the precision—testing the same (or similar) inputs only provides local information and therefore should not increase our knowledge about the global maximum. A simple and effective way of implementing this idea is given by
(3) 
where , , , and are parameters of the estimator: is the precision we gain for each new distinct observation; is the number of prior points; is a finite, symmetric kernel function; is a prior precision function; and is a prior estimate of .
In (3), the mean function is estimated with a kernel regressor [11], and the total effective number of locations is calculated as the sum of the prior locations and the number of distinct locations in the data . The latter is estimated by multiplying the number of data points with the coefficient
i.e. the ratio between the trace of the Gramian matrix and the sum of its entries. Inputs that are very close to each other will have overlapping kernels, resulting in large offdiagonal entries of the Gramian matrix—hence decreasing the number of distinct locations (Figure 1b).
The expression for the posterior can be calculated, up to a constant factor, in quadratic time in the number of observations. It can therefore be easily combined with Markov chain Monte Carlo methods (MCMC) to implement stochastic optimizers as illustrated in Section
4.3 Derivation of the Model
3.1 FunctionBased, Indirect Model
Our first task is to derive an indirect Bayesian model for the optimal test point that builds its estimate via the underlying function space. Let be the set of hypotheses, and assume that each hypothesis corresponds to a stochastic mapping . Let be the prior^{2}^{2}2For the sake of simplicity, we neglect issues of measurability of . over and let the likelihood be . Then, the posterior of is given by
(4) 
For each , let be the subset of functions such that for all , ^{3}^{3}3Note that we assume that the mean function is bounded and that it has a unique maximizing test point.. Then, the posterior over the optimal test point is given by
(5) 
This model has two important drawbacks: (a) it relies on modeling the entire function space , which is potentially much more complex than necessary; (b) it requires calculating the integral (5), which is intractable for virtually all realworld problems.
3.2 DomainBased, Direct Model
We want to arrive at a Bayesian model that bypasses the integration step suggested by (5) and directly models the location of optimal test point . The following theorem explains how this direct model relates to the previous model.
Theorem 1.
The Bayesian model for the optimal test point is given by
(prior)  
(likelihood) 
where is the set of past tests.
Proof.
Using Bayes’ rule, the posterior distribution can be rewritten as
(6) 
Since this posterior is equal to (5), one concludes (using (4)) that
Note that this expression corresponds to the joint . The prior is obtained by setting . The likelihood is obtained as the fraction
where it shall be noted that the denominator doesn’t change if we add the condition . ∎
From Theorem 1 it is seen that although the likelihood model for the indirect model is i.i.d. at each test point, the likelihood model for the direct model depends on the past tests , that is, it is adaptive. More critically though, the likelihood function’s internal structure of the direct model corresponds to an integration over function space as well—thus inheriting all the difficulties of the indirect model.
3.3 Abstract Properties of the Likelihood Function
There is a way to bypass modeling the function space explicitly if we make a few additional assumptions. We assume that for any , the mean function is continuous and has a unique maximum. Then, the crucial insight consists in realizing that the value of the mean function inside a sufficiently small neighborhood of is larger than the value outside of it (see Figure 2a).
We assume that, for any and any , let denote the open ball centered on . The functions in fulfill the following properties:

Continuous: Every function is such that its mean is continuous and bounded.

Maximum: For any , the functions are such that for all and all , .
Furthermore, we impose a symmetry condition on the likelihood function. Let and be in , and consider their associated equivalence classes and . There is no reason for them to be very different: in fact, they should virtually be indistinguishable outside of the neighborhoods of and . It is only inside of the neighborhood of when becomes distinguishable from the other equivalence classes because the functions in systematically predict higher values than the rest. This assumption is illustrated in Figure 2b. In fact, taking the loglikelihood ratio of two competing hypotheses
for a given test location should give a value equal to zero unless is inside of the vicinity of or (see Figure 2c). In other words, the amount of evidence a hypothesis gets when the test point is outside of its neighborhood is essentially zero (i.e. it is the same as the amount of evidence that most of the other hypotheses get).
3.4 Likelihood and Conjugate Prior
Following our previous discussion, we propose the following likelihood model. Given the previous data and a test point , the likelihood of the observation is
(7) 
where: is a normalizing constant;
is a posterior probability over
given and the data ; is a precision measuring the knowledge we have about the whole function given bywhere is a precision scaling parameter; is a parameter representing the number prior locations tested; and is an estimate of the mean function given by
In the last expression, corresponds to a prior estimate of with prior precision . Inspecting (7), we see that the likelihood model favours positive changes to the estimated mean function from new, unseen test locations. The pdf does not need to be explicitly defined, as it will later drop out when computing the posterior. The only formal requirement is that it should be independent of the hypothesis .
We propose the conjugate prior
(8) 
The conjugate prior just encodes a prior estimate of the mean function. In a practical optimization application, it serves the purpose of guiding the exploration of the domain, as locations with high prior value are more likely to contain the maximizing argument.
Given a set of data points , the prior (8) and the likelihood (7) lead to a posterior given by
(9) 
Thus, the particular choice of the likelihood function guarantees an analytically compact posterior expression. In general, the normalizing constant in (9) is intractable, which is why the expression is only practical for relative comparisons of test locations. Substituting the precision and the mean function estimate yields
4 Experimental Results
4.1 Parameters.
We have investigated the influence of the parameters on the resulting posterior probability distribution. Figure
3 shows how the choice of the precision and the kernel width affect the shape of the posterior probability density. We have used the Gaussian kernel(10) 
In this figure, 7 data points are shown, which were drawn as , where the mean function is
(11) 
The functions and were chosen as
(12) 
where the latter corresponds to the logarithm of a Gaussian with mean
and variance
. Choosing a higher value for leads to sharper updates, while higher values for the kernel width produce smoother posterior densities.4.2 Application to Optimization.
Comparison to Gaussian Process UCB.
We have used the model to optimize the same function (11) as in our preliminary tests but with higher additive noise equal to one. This is done by sampling the next test point directly from the posterior density over the optimum location , and then using the resulting pair
to recursively update the model. Essentially, this procedure corresponds to Bayesian control rule/Thompson sampling
[12, 13].We compared our method against a Gaussian Process optimization method using an upper confidence bound (UCB) criterion [10]. The parameters for the GPUCB were set to the following values: observation noise and length scale . For the constant that trades off exploration and exploitation we followed Theorem in [10] which states with . We have implemented our proposed method with a Gaussian kernel as in (10) with width . The prior sufficient statistics are exactly as in (12). The precision parameter was set to .
Simulation results over ten independent runs are summarized in Figure 4. We show the timeaveraged observation values of the noisy function evaluated at test locations sampled from the posterior. Qualitatively, both methods show very similar convergence (on average), however our method converges faster and with a slightly higher variance.
HighDimensional Problem.
To test our proposed method on a challenging problem, we have designed a nonconvex, highdimensional noisy function with multiple local optima. This Noisy Ripples function is defined as
where is the location of the global maximum, and where observations have additive Gaussian noise with zero mean and variance . The advantage of this function is that it generalizes well to any number of dimensions of the domain. Figure 5a illustrates the function for the 2dimensional input domain. This function is difficult to optimize because it requires averaging the noisy observations and smoothing the ridged landscape in order to detect the underlying quadratic form.
We optimized the 50dimensional version of this function using a MetropolisHastings scheme to sample the next test locations from the posterior over the maximizing argument. The Markov chain was started at , executing 120 isotropic Gaussian steps of variance before the point was used as an actual test location. For the argmax prior, we used a Gaussian kernel with lengthscale , precision factor , prior precision and prior mean estimate . The goal was located at the origin.
The result of one run is presented in Figure 5b. It can be seen that the optimizer manages to quickly ( samples) reach nearoptimal performance, overcoming the difficulties associated with the highdimensionality of the input space and the numerous local optima. Crucial for this success was the choice of a kernel that is wide enough to accurately estimate the mean function. The authors are not aware of any method capable of solving a problem is similar characteristics.
5 Discussion & Conclusions
We have proposed a novel Bayesian approach to model the location of the maximizing test point of a noisy, nonlinear function. This has been achieved by directly constructing a probabilistic model over the input space, thereby bypassing having to model the underlying function space—a much harder problem. In particular, we derived a likelihood function that belongs to the exponential family by assuming a form of symmetry in function space. This in turn, enabled us to state a conjugate prior distribution over the optimal test point.
Our proposed model is computationally very efficient when compared to Gaussian processbased (cubic) or UCBbased models (expensive computation of ). The evaluation time of the posterior density scales quadratically in the size of the data. This is due to the calculation of the effective number of previously seen test locations—the kernel regressor requires linear compuation time. However, during MCMC steps, the effective number of test locations does not need to be updated as long as no new observations arrive.
In practice, one of the main difficulties associated with our proposed method is the choice of the parameters. As in any kernelbased estimation method, choosing the appropriate kernel bandwidth can significantly change the estimate and affect the performance of optimizers that rely on the model. There is no clear rule on how to choose a good bandwidth.
In a future research, it will be interesting to investigate the theoretical properties of the proposed nonparametric model, such as the convergence speed of the estimator and its relation to the extensive literature on active learning and bandits.
References

Brochu et al. [2009]
E. Brochu, V. Cora, and N. de Freitas.
A tutorial on bayesian optimization of expensive cost functions, with application to active user modeling and hierarchical reinforcement learning.
Technical Report TR2009023, University of British Columbia, Department of Computer Science, 2009.  Rawlik et al. [2010] K. Rawlik, M. Toussaint, and S. Vijayakumar. Approximate inference and stochastic optimal control. arXiv:1009.3958, 2010.
 Shapiro [2000] A. Shapiro. Probabilistic Constrained Optimization: Methodology and Applications, chapter Statistical Inference of Stochastic Optimization Problems, pages 282–304. Kluwer Academic Publishers, 2000.
 Kappen et al. [2012] H.J. Kappen, V. Gómez, and M. Opper. Optimal control as a graphical model inference problem. Machine Learning, 87(2):159–182, 2012.
 Kushner and Yin [1997] H.J. Kushner and G.G. Yin. Stochastic Approximation Algorithms and Applications. SpringerVerlag, 1997.
 Mockus [1994] J. Mockus. Application of bayesian approach to numerical methods of global and stochastic optimization. Journal of Global Optimization, 4(4):347–365, 1994.
 Lizotte [2008] D. Lizotte. Practical Bayesian Optimization. Phd thesis, University of Alberta, 2008.
 Jones et al. [1998] D.R. Jones, M. Schonlau, and W.J. Welch. Efficient global optimization of expensive blackbox functions. Journal of Global Optimization, 13(4):455–492, 1998.
 Osborne et al. [2009] M.A. Osborne, R. Garnett, and S.J. Roberts. Gaussian processes for global optimization. In 3rd International Conference on Learning and Intelligent Optimization (LION3), 2009.
 Srinivas et al. [2010] N. Srinivas, A. Krause, S. Kakade, and M. Seeger. Gaussian process optimization in the bandit setting: No regret and experimental design. In International Conference on Machine Learning, 2010.
 Hastie et al. [2009] T. Hastie, R. Tbshirani, and J. Friedman. The Elements of Statistical Learning. Springer, second edition, 2009.
 May and Leslie [2011] B.C. May and D.S. Leslie. Simulation studies in optimistic Bayesian sampling in contextualbandit problems. Technical Report 11:02, Statistics Group, Department of Mathematics, University of Bristol, 2011.

Ortega and Braun [2010]
P.A. Ortega and D.A. Braun.
A minimum relative entropy principle for learning and acting.
Journal of Artificial Intelligence Research
, 38:475–511, 2010.
References

Brochu et al. [2009]
E. Brochu, V. Cora, and N. de Freitas.
A tutorial on bayesian optimization of expensive cost functions, with application to active user modeling and hierarchical reinforcement learning.
Technical Report TR2009023, University of British Columbia, Department of Computer Science, 2009.  Rawlik et al. [2010] K. Rawlik, M. Toussaint, and S. Vijayakumar. Approximate inference and stochastic optimal control. arXiv:1009.3958, 2010.
 Shapiro [2000] A. Shapiro. Probabilistic Constrained Optimization: Methodology and Applications, chapter Statistical Inference of Stochastic Optimization Problems, pages 282–304. Kluwer Academic Publishers, 2000.
 Kappen et al. [2012] H.J. Kappen, V. Gómez, and M. Opper. Optimal control as a graphical model inference problem. Machine Learning, 87(2):159–182, 2012.
 Kushner and Yin [1997] H.J. Kushner and G.G. Yin. Stochastic Approximation Algorithms and Applications. SpringerVerlag, 1997.
 Mockus [1994] J. Mockus. Application of bayesian approach to numerical methods of global and stochastic optimization. Journal of Global Optimization, 4(4):347–365, 1994.
 Lizotte [2008] D. Lizotte. Practical Bayesian Optimization. Phd thesis, University of Alberta, 2008.
 Jones et al. [1998] D.R. Jones, M. Schonlau, and W.J. Welch. Efficient global optimization of expensive blackbox functions. Journal of Global Optimization, 13(4):455–492, 1998.
 Osborne et al. [2009] M.A. Osborne, R. Garnett, and S.J. Roberts. Gaussian processes for global optimization. In 3rd International Conference on Learning and Intelligent Optimization (LION3), 2009.
 Srinivas et al. [2010] N. Srinivas, A. Krause, S. Kakade, and M. Seeger. Gaussian process optimization in the bandit setting: No regret and experimental design. In International Conference on Machine Learning, 2010.
 Hastie et al. [2009] T. Hastie, R. Tbshirani, and J. Friedman. The Elements of Statistical Learning. Springer, second edition, 2009.
 May and Leslie [2011] B.C. May and D.S. Leslie. Simulation studies in optimistic Bayesian sampling in contextualbandit problems. Technical Report 11:02, Statistics Group, Department of Mathematics, University of Bristol, 2011.

Ortega and Braun [2010]
P.A. Ortega and D.A. Braun.
A minimum relative entropy principle for learning and acting.
Journal of Artificial Intelligence Research
, 38:475–511, 2010.
Comments
There are no comments yet.