基于高斯过程的贝叶斯优化(三)GP超参数的估计
前面的文章大致描述了基于高斯过程(GP)贝叶斯优化的原理框架，该框架中也存在了几个参数，本篇文章简单介绍如何对他们进行估计。
首先介绍一下贝叶斯优化框架的超参数有哪些：
回忆我们将高斯过程表述为以下形式:
\[f ( x ) \sim G P \left( m ( x ) , k \left( x , x ^ { \prime } \right) \right)\]
其中$m(x)$表示均值函数，一般都设为0，不需要更新，我们更关心的是核函数k，核函数的选取主要有两种：squared exponential kernel以及Matern kernel
下面给出两种核函数的具体形式：
squared exponential kernel:
\[k \left( \mathbf { x } _ { i } , \mathbf { x } _ { j } \right) = \exp \left( - \frac { 1 } { 2 } \left( \mathbf { x } _ { i } - \mathbf { x } _ { j } \right) ^ { T } \operatorname { diag } ( \boldsymbol { \theta } ) ^ { - 2 } \left( \mathbf { x } - \mathbf { x } ^ { \prime } \right) \right)\]
其中$\operatorname { diag }\boldsymbol( { \theta })$表示对角阵
Matern kernel:
\[k \left( \mathbf { x } _ { i } , \mathbf { x } _ { j } \right) = \frac { 1 } { 2 ^ { < - 1 } \Gamma ( \zeta ) } \left( 2 \sqrt { \zeta } \left\| \mathbf { x } _ { i } - \mathbf { x } _ { j } \right\| \right) ^ { \zeta } H _ { \zeta } \left( 2 \sqrt { \zeta } \left\| \mathbf { x } _ { i } - \mathbf { x } _ { j } \right\| \right)\]
实践中，Jasper Snoek[1]推荐采用Matern kernel。
问题在于如何对核函数中存在的超参数进行估计呢？
实际中一般有两种方法，极大似然与MCMC估计法，先介绍采用极大似然方法更新超参数。
极大似然方法
极大似然通过直接使得GP超参数$\theta$的后验分布最大化进行计算。
\[L = \log p ( y | \theta ) = - \frac { 1 } { 2 } y ^ { T } \left( K + \sigma _ { n } ^ { 2 } I \right) ^ { - 1 } y - \frac { 1 } { 2 } \log \left| K + \sigma _ { n } ^ { 2 } I \right| - \frac { n } { 2 } \log 2 \pi\]
注意这里的形式是加了噪声项的高斯模型。其中$y = f \left( \mathbf { x } \right) + \epsilon,$，并且$\epsilon \sim \mathcal { N } \left( 0 , \sigma _ { \text { noise } } ^ { 2 } \right)$，此时的协方差矩阵K为
\[\mathbf { K } = \left[ \begin{array} { c c c } { k \left( \mathbf { x } _ { 1 } , \mathbf { x } _ { 1 } \right) } & { \dots } & { k \left( \mathbf { x } _ { 1 } , \mathbf { x } _ { t } \right) } \\ { \vdots } & { \ddots } & { \vdots } \\ { k \left( \mathbf { x } _ { t } , \mathbf { x } _ { 1 } \right) } & { \dots } & { k \left( \mathbf { x } _ { t } , \mathbf { x } _ { t } \right) } \end{array} \right] + \sigma _ { \text { nois } } ^ { 2 } I\]
此时：
\[P \left( y _ { t + 1 } | \mathcal { D } _ { 1 : t } , \mathbf { x } _ { t + 1 } \right) = \mathcal { N } \left( \mu _ { t } \left( \mathbf { x } _ { t + 1 } \right) , \sigma _ { t } ^ { 2 } \left( \mathbf { x } _ { t + 1 } \right) + \sigma _ { \text { noise } } ^ { 2 } \right)\]
其中
\[\mu _ { t } \left( \mathbf { x } _ { t + 1 } \right) = \mathbf { k } ^ { T } \left[ \mathbf { K } + \sigma _ { \text { noise } } ^ { 2 } I \right] ^ { - 1 } \mathbf { y } _ { 1 : t }\]
\[\sigma _ { t } ^ { 2 } \left( \mathbf { x } _ { t + 1 } \right) = k \left( \mathbf { x } _ { t + 1 } , \mathbf { x } _ { t + 1 } \right) - \mathbf { k } ^ { T } \left[ \mathbf { K } + \sigma _ { \text { noise } } ^ { 2 } I \right] ^ { - 1 } \mathbf { k }\]
MCMC方法
MCMC方法是通过马尔科夫链实现对一个复杂分布p(x)进行采样的技术。这里直接略过对MCMC方法的基本介绍。
注意到我们的最终目的仍然是求模型M中某个特定超参数下acquisition function值，因此，可以通过条件概率方法进行展开
\[\hat { a } \left( \mathbf { x } ; \left\{ \mathbf { x } _ { n } , y _ { n } \right\} \right) = \int a \left( \mathbf { x } ; \left\{ \mathbf { x } _ { n } , y _ { n } \right\} , \theta \right) p ( \theta | \left\{ \mathbf { x } _ { n } , y _ { n } \right\} _ { n = 1 } ^ { N } ) \mathrm { d } \theta\]
这里的$\theta$表示GP中的超参数，而 $\left\{ \mathbf { x } _ { n } , y _ { n } \right\} _ { n = 1 } ^ { N }$表示的针对模型M选择超参数(x)并进行评估的结果(y)。为了估计这个积分的值，我们需要实现对分布$ p ( \theta | \left\{ \mathbf { x } _ { n } , y _ { n } \right\} _ { n = 1 } ^ { N } )$的采样。注意到这里
\[p ( \theta | \left\{ \mathbf { x } _ { n } , y _ { n } \right\} _ { n = 1 } ^ { N } ) \propto p(\theta) p(\left\{ \mathbf { x } _ { n } , y _ { n } \right\} _ { n = 1 } ^ { N } | \theta)\]
也就是说只要再给GP超参数$\theta$提供一个先验让上式右边可以具体计算，分布$ p ( \theta | \left\{ \mathbf { x } _ { n } , y _ { n } \right\} _ { n = 1 } ^ { N } )$就能通过MCMC方法就行采样。在这里还存在一个问题，当对$\theta$进行采样时，如果采用Metropolis-Hastings采样方法，由于接受率可能很低，而在计算接受概率时势必需要根据新采样的$\theta$重新计算其后验概率，这可能造成很大的计算开销，因此Iain Murray[2]等人改进了Elliptical Slice sampling方法进行采样，并取得了较好的效果。
References
[1] Snoek J , Larochelle H , Adams R P . Practical Bayesian Optimization of Machine Learning Algorithms[J]. Advances in neural information processing systems, 2012.
[2] Murray I , Adams R P . Slice sampling covariance hyperparameters of latent Gaussian models[C]// International Conference on Neural Information Processing Systems. Curran Associates Inc. 2010.
[3] Brochu E , Cora V M , De Freitas N . A Tutorial on Bayesian Optimization of Expensive Cost Functions, with Application to Active User Modeling and Hierarchical Reinforcement Learning[J]. Computer Science, 2010.
