贝叶斯优化
第三部分：贝叶斯优化
EN

采集函数

第 11 章构建了贝叶斯优化循环，并用一条规则选择下一次评估：后验均值加上后验标准差的某个倍数。这条规则有效，但它的权重相当于一个旋钮，而图 11.3 表明，旋钮向任一方向调得过头，代价都很大。读者自然会问：能否以有原则的方式确定不确定性的价值？

本章给出几种回答。每个采集函数都从评估的目的出发：超过迄今见到的最优值、改进最终的推荐，或者确定最大值的位置。由这一目的和高斯后验即可推出公式，探索与利用之间的平衡也由推导得出，无须手工设定。本章逐个推导经典的采集函数，先后在一维和二维中基于同一后验比较它们，最后讨论它们共同遗留的问题：求采集函数本身的最大值。维度很高时，这本身就是一个困难的优化问题。

12.1 什么是采集函数 #

先回顾设定。经过 nn 次评估后，数据为 Dn={(xi,yi)}i=1n\D_n = \{(\vx_i, y_i)\}_{i=1}^n；高斯过程后验对每个输入 x\vx 给出关于 f(x)f(\vx) 的高斯信念，均值为 μn(x)\mu_n(\vx)，标准差为 σn(x)\sigma_n(\vx)（第 8 章）。采集函数 an(x)a_n(\vx) 为每个输入打分，循环在得分最高处评估目标函数（定义 11.1）。

构造这种得分最清晰的办法，是先说明最终希望得到什么，再问多做一次评估预期能增加多少。记 u(D)u(\D) 为数据集的效用（utility），它是一个数，衡量带着数据 D\D 停止时处境有多好。如果在 x\vx 处评估并观测到 yy，效用就从 u(Dn)u(\D_n) 变为 u(Dn∪{(x,y)})u(\D_n \cup \{(\vx, y)\})。评估之前 yy 是未知的，但后验给出了它可能的取值，因此可以对它取平均。

定义 12.1 一步前瞻采集函数

给定效用 uu，一步前瞻（one-step lookahead）采集函数是在 x\vx 处再做一次评估所带来的效用期望增益：

an(x)=En ⁣[ u(Dn∪{(x,y)})−u(Dn) ],a_n(\vx) = \E_n\!\left[\,u(\D_n \cup \{(\vx, y)\}) - u(\D_n)\,\right],
(12.1)

其中 En\E_n 表示给定 Dn\D_n 时，按后验预测分布对结果 yy 取平均。

不同的效用给出不同的采集函数。表 12.1 预先列出本章讨论的几种。其中上置信界与 Thompson 采样并非来自效用，而是来自第 13.2 节的赌博机问题，在那里它们具有另一类理论保证（第 13 章）。

表 12.1 各采集函数及其回答的问题。
采集函数 一次评估的价值 所在节
改进概率（PI） 超过迄今最优值的概率 第 12.2 节
期望改进（EI） 迄今最优值预期上升的幅度 第 12.3 节
上置信界（UCB） 对该输入处函数值的乐观估计 第 12.4 节
Thompson 采样（TS） 该输入是最大值点的概率 第 12.5 节
知识梯度（KG） 最终推荐的价值预期上升的幅度 第 12.6 节
熵搜索（ES、PES、MES） 关于最大值点或最大值的期望信息量 第 12.7 节

12.1.1 为什么只向前看一步 #

式（12.1）只向前看一次评估，仿佛下一次评估就是最后一次。这是一种简化。如果还剩 N−nN - n 次评估，当前的最优选择就取决于此后还能做什么；下一次评估若能为后续的评估创造条件，就更有价值。把这一点写成式子，便得到一个动态规划，其中每个未来的结果都分叉出每个未来的选择。它给出最优策略，但用 Kushner 的话说，也“几乎不可能计算”（virtually impossible to compute）（Kushner，1964；转引自 Garnett，2023，第 12.3 节）。

因此，实际使用的采集函数几乎都是短视（myopic）的：若下一次评估就是最后一次，它们是最优的，否则只是近似最优（Frazier，2018）。这种近似可能比听上去更好。在少数能够算出最优多步策略的特殊问题中，短视规则的表现接近最优；Frazier（2018）引用的一项此类研究中，知识梯度达到了最优值的 98% 以上。不过，短视也有显而易见的代价：忽视未来的规则会低估探索，因为探索的回报来得较晚。下文有几条规则以不同方式把探索重新纳入考虑。

第 12.1 节引用的文献 3
  1. Kushner（1964）A New Method of Locating the Maximum Point of an Arbitrary Multipeak Curve in the Presence of Noise
  2. Garnett（2023）Bayesian Optimization
  3. Frazier（2018）A Tutorial on Bayesian Optimization

12.2 改进概率 #

最古老的规则提出的问题最简单：x\vx 超过迄今最优值的可能性有多大（Kushner，1964）？记这个值为 fn∗=max⁡i≤nyif^*_n = \max_{i \le n} y_i，称为当前最优值（incumbent）。评估精确时，它就是此刻所能推荐的最佳输入的值。

按照后验，f(x)f(\vx) 服从均值为 μn(x)\mu_n(\vx)、标准差为 σn(x)\sigma_n(\vx) 的高斯分布。它超过当前最优值至少一个裕量 ξ≥0\xi \ge 0 的概率，可分三步求得。

推导改进概率
  1. 标准化：Z=(f(x)−μn(x))/σn(x)Z = (f(\vx) - \mu_n(\vx)) / \sigma_n(\vx) 是标准正态变量（第 4.3 节）。
  2. 改写事件：f(x)>fn∗+ξf(\vx) > f^*_n + \xi 等价于 Z>(fn∗+ξ−μn(x))/σn(x)Z > (f^*_n + \xi - \mu_n(\vx)) / \sigma_n(\vx)。
  3. 由标准正态密度的对称性，P(Z>a)=1−Φ(a)=Φ(−a)\Prob(Z > a) = 1 - \Phi(a) = \Phi(-a)。
PI⁡n(x)=Φ ⁣(μn(x)−fn∗−ξσn(x)).\PI_n(\vx) = \Phi\!\left(\frac{\mu_n(\vx) - f^*_n - \xi}{\sigma_n(\vx)}\right).
(12.2)

用定义 12.1 的语言表述，改进概率是如下效用的期望增益：最优值提升至少 ξ\xi 时效用为 1，否则为 0。

改进概率有一个众所周知的缺陷：只计较是否改进，而不问改进多少。当 ξ=0\xi = 0 时，若某个输入的均值略高于当前最优值、标准差又极小，其改进概率接近 1，规则会优先选它，而不是选一个可能大幅超过当前最优值的不确定输入。结果是搜索在最优点附近以极小的步长缓慢移动。裕量 ξ\xi 正是为此而设。增大裕量，就是要求值得一取的改进，这会使大多数输入在式（12.2）中的均值项小于零；此时只有较大的 σn(x)\sigma_n(\vx) 才能带来可观的概率，规则便转向探索。

ξ\xi 应当取多大？Kushner 建议开始时取较大的值，随搜索进行逐渐减小（Brochu 等，2010），他的论文也讨论过在搜索过程中手动调整它（Garnett，2023，第 12.3 节）。Jones 发现这一方法“对目标值的选择极其敏感”（extremely sensitive to the choice of the target）：太小则搜索局限于局部，太大则永远不会细化有希望的解（Jones，2001；转引自 Brochu 等，2010）。下一节的规则直接考虑改进的大小，因而裕量的作用要小得多。

第 12.2 节引用的文献 4
  1. Kushner（1964）A New Method of Locating the Maximum Point of an Arbitrary Multipeak Curve in the Presence of Noise
  2. Brochu 等人（2010）A Tutorial on Bayesian Optimization of Expensive Cost Functions, with Application to Active User Modeling and Hierarchical Reinforcement Learning
  3. Garnett（2023）Bayesian Optimization
  4. Jones（2001）A Taxonomy of Global Optimization Methods Based on Response Surfaces

12.3 期望改进 #

与其问当前最优值是否会改进，不如问平均改进多少。在 x\vx 处评估，会把见到的最优值从 fn∗f^*_n 提升到 max⁡(fn∗,f(x))\max(f^*_n, f(\vx))。增益即改进（improvement）max⁡(f(x)−fn∗,0)\max(f(\vx) - f^*_n, 0)，当 f(x)f(\vx) 不够高时为零；它在后验下的期望就是期望改进：

EI⁡n(x)=En ⁣[max⁡ ⁣(f(x)−fn∗,  0)].\EI_n(\vx) = \E_n\!\left[\max\!\left(f(\vx) - f^*_n,\; 0\right)\right].
(12.3)

这正是取效用 u(D)=max⁡iyiu(\D) = \max_i y_i（观测到的最优值）时的式（12.1）。期望改进通常归功于 Močkus 及其同事（Močkus，1975；Frazier，2018），更早的显式公式可追溯到 Šaltenis 1971 年的工作（Garnett，2023，第 12.3 节）；它借由 Jones 等人（1998）的 EGO 算法成为标准方法。

由于 fn∗f^*_n 是已知的数，差值 f(x)−fn∗f(\vx) - f^*_n 是高斯变量，问题于是归结为高斯分布的一个性质。下面对任意高斯变量 DD 陈述这一性质，因为第 19.4 节会把它用于另一个高斯变量。

推导高斯变量正部的期望

设 DD 服从均值为 δ\delta、标准差为 s>0s > 0 的高斯分布，求 E[max⁡(D,0)]\E[\max(D, 0)]。

  1. 写成 D=δ+sZD = \delta + s Z，其中 ZZ 为标准正态变量（第 4.3 节）。
  2. 仅当 D>0D > 0，即 Z>−δ/sZ > -\delta/s 时，max⁡(D,0)\max(D, 0) 才不为零。因此 E[max⁡(D,0)]=∫−δ/s∞(δ+sz) ϕ(z) dz\E[\max(D, 0)] = \int_{-\delta/s}^{\infty} (\delta + s z)\, \phi(z)\, \dd z，其中 ϕ\phi 为标准正态密度。
  3. 把积分拆成两部分：δ∫−δ/s∞ϕ(z) dz+s∫−δ/s∞z ϕ(z) dz\delta \int_{-\delta/s}^{\infty} \phi(z)\, \dd z + s \int_{-\delta/s}^{\infty} z\, \phi(z)\, \dd z。
  4. 由 ϕ\phi 的对称性，第一个积分为 1−Φ(−δ/s)=Φ(δ/s)1 - \Phi(-\delta/s) = \Phi(\delta/s)。
  5. 对于第二个积分，密度满足 ϕ′(z)=−z ϕ(z)\phi'(z) = -z\,\phi(z)，故 z ϕ(z)z\,\phi(z) 的原函数为 −ϕ(z)-\phi(z)，从而 ∫a∞z ϕ(z) dz=ϕ(a)\int_{a}^{\infty} z\,\phi(z)\, \dd z = \phi(a)。取 a=−δ/sa = -\delta/s，并利用 ϕ(−a)=ϕ(a)\phi(-a) = \phi(a)，得第二个积分为 ϕ(δ/s)\phi(\delta/s)。
  6. 合并两项。
E[max⁡(D,0)]=δ Φ ⁣(δs)+s ϕ ⁣(δs).\E[\max(D, 0)] = \delta\, \Phi\!\left(\frac{\delta}{s}\right) + s\, \phi\!\left(\frac{\delta}{s}\right).
(12.4)

当 s=0s = 0 时，DD 为常数 δ\delta，期望为 max⁡(δ,0)\max(\delta, 0)，这也是 s→0s \to 0 时式（12.4）的极限。

对于期望改进，取 D=f(x)−fn∗−ξD = f(\vx) - f^*_n - \xi，其中可选的裕量 ξ\xi 与改进概率中的相同。其均值为 δ=μn(x)−fn∗−ξ\delta = \mu_n(\vx) - f^*_n - \xi，标准差为 s=σn(x)s = \sigma_n(\vx)：

EI⁡n(x)=(μn(x)−fn∗−ξ)Φ(z)+σn(x) ϕ(z),z=μn(x)−fn∗−ξσn(x).\EI_n(\vx) = \left(\mu_n(\vx) - f^*_n - \xi\right) \Phi(z) + \sigma_n(\vx)\, \phi(z), \qquad z = \frac{\mu_n(\vx) - f^*_n - \xi}{\sigma_n(\vx)}.
(12.5)

这就是经 Jones 等人（1998）而成为标准的闭式解，裕量 ξ\xi 是后来加入的。Brochu 等人（2010）报告了 Lizotte 的实验，结果提示 ξ=0.01\xi = 0.01（必要时按信号方差缩放）在几乎所有情形下都表现良好。

12.3.1 解读公式 #

式（12.5）的两项在同一个表达式中分别体现利用与探索。第一项在均值高于当前最优值处较大，第二项在标准差较大处较大。与上置信界不同，两者之间的兑换率并非人为选定，而是由积分得出。

以下三条性质使这一点更加精确。记 EI⁡(δ,s)\EI(\delta, s) 为式（12.4）的右端。

  • 期望改进不低于均值的改进。max⁡(⋅,0)\max(\cdot, 0) 是凸函数，由 Jensen 不等式得 E[max⁡(D,0)]≥max⁡(E[D],0)=max⁡(δ,0)\E[\max(D, 0)] \ge \max(\E[D], 0) = \max(\delta, 0)。不确定性只会增加价值。
  • 均值等于当前最优值时，期望改进约为 0.4 s0.4\,s。取 δ=0\delta = 0，由式（12.4）得 EI⁡=s ϕ(0)=s/2π≈0.399 s\EI = s\,\phi(0) = s / \sqrt{2\pi} \approx 0.399\,s。这正是第 2.6 节中的例子：均值的改进为零，期望改进却不为零。
  • 期望改进随均值和不确定性的增大而增大。偏导数为 ∂EI⁡/∂δ=Φ(δ/s)\partial \EI / \partial \delta = \Phi(\delta/s) 和 ∂EI⁡/∂s=ϕ(δ/s)\partial \EI / \partial s = \phi(\delta/s)（习题 12.1），两者均为正。无论均值处于何处，不确定性越大，期望改进越大。

期望改进与改进概率正是在最后一条性质上分道扬镳。改进概率对 ss 的导数为 −(δ/s2) ϕ(δ/s)-(\delta/s^2)\,\phi(\delta/s)，只要 δ>0\delta > 0 就为负：一旦均值高于当前最优值，改进概率便偏好确定性。图 12.1 展示了单个输入处的这两个量。

D 的密度D 大于 0 的概率（面积 = PI）max(D, 0) × 密度（面积 = EI）δ = −0.30，s = 0.50：PI = 0.274，EI = 0.0840.00.51.0−3−2−10123相对当前最优值的改进 D = f(x) − f*ₙf*ₙδPI 随 s 的变化（δ 固定）0.00.20.40.60.81.00.00.51.01.5标准差 sEI 随 s 的变化（δ 固定）0.00.20.40.00.51.01.5标准差 s
D 的密度D 大于 0 的概率（面积 = PI）max(D, 0) × 密度（面积 = EI）δ = −0.30，s = 0.50：PI = 0.274，EI = 0.0840.00.51.0−3−2−10123D = f(x) − f*ₙf*ₙδPI 随 s 的变化（δ 固定）0.00.20.40.60.81.00.00.51.01.5标准差 sEI 随 s 的变化（δ 固定）0.00.20.40.00.51.01.5标准差 s
图 12.1 单个输入处的改进。上：D=f(x)−fn∗D = f(x) - f^*_n 的密度，该变量服从高斯分布，均值 δ\delta 与标准差 ss 由滑块设定。零右侧的阴影面积为改进概率；虚线为改进 max⁡(D,0)\max(D, 0) 与密度的乘积，其下面积为期望改进，即式（12.4）。下：δ\delta 固定时，改进概率（PI）与期望改进（EI）随 ss 的变化，圆点为当前的 ss。当 δ>0\delta > 0 时，EI 面板中的虚线标出 δ\delta，即 ss 缩小时期望改进趋近的值。

以下三种设置显示了两者的差别。

均值低于当前最优值。图的初始设置为 δ=−0.3\delta = -0.3。改进概率与期望改进都随 ss 增大而增大：要超过当前最优值，函数值必须高于模型的预期，而不确定性越大，这种情况越可能出现。

均值高于当前最优值。设 δ=0.5\delta = 0.5。此时改进概率反而随 ss 增大而下降，因为更宽的信念把更多权重放在了达不到当前最优值的情形上。期望改进仍在上升：大幅改进所获得的额外权重，超过了达不到的情形所获得的额外权重，而后者没有任何代价，因为改进从不为负。

微小而确定的改进。设 δ=0.05\delta = 0.05，s=0.05s = 0.05。改进概率约为 0.84，而期望改进只有约 0.05。另一个输入取 δ=−0.3\delta = -0.3、s=1s = 1，改进概率只有 0.38，但期望改进约为 0.27，约为前一个输入的五倍。期望改进把第二个输入排在前面，改进概率则把第一个排在前面。

易错点期望改进会下溢

在大定义域中远离数据处，或在运行后期当前最优值已经很高时，δ/s\delta / s 是绝对值很大的负数，式（12.5）成为一些极小的数之差，浮点运算会把它们舍入为零。于是采集函数在定义域的大部分区域值和梯度都恰好为零，从那里出发的基于梯度的优化器无法移动。Ament 等人（2023）认为，文献中期望改进表现不稳定且往往偏弱，原因在于这一数值问题，而不在于期望改进这一思想本身。他们提出了 LogEI，用在尾部仍保持精确的公式计算 log⁡EI⁡\log \EI。LogEI 的最大值点与期望改进的相同或几乎相同；在他们的实验中，它与较新的采集函数持平或更优。BoTorch 以 LogExpectedImprovement 提供这一方法。

最后补充两点实际说明。第一，式（12.5）假定评估精确，因而 fn∗f^*_n 是已知值。有噪声时，观测到的最优值本身不确定且偏高，有几种带噪声的变体用其他量取代它，第 14.2 节对这些变体做了比较。第二，一个真实应用展示了这一公式的实际效果。Snoek 等人（2012）用期望改进为机器学习算法调参，他们没有固定高斯过程的超参数，而是在超参数的样本上对式（12.5）取平均。他们还用第二个高斯过程为训练时间建模，并最大化每秒期望改进（expected improvement per second），这一准则偏好既有希望、又能快速评估的输入。第 22 章重现了一个这类问题。

第 12.3 节引用的文献 7
  1. Močkus（1975）On Bayesian Methods for Seeking the Extremum
  2. Frazier（2018）A Tutorial on Bayesian Optimization
  3. Garnett（2023）Bayesian Optimization
  4. Jones 等人（1998）Efficient Global Optimization of Expensive Black-Box Functions
  5. Brochu 等人（2010）A Tutorial on Bayesian Optimization of Expensive Cost Functions, with Application to Active User Modeling and Hierarchical Reinforcement Learning
  6. Ament 等人（2023）Unexpected Improvements to Expected Improvement for Bayesian Optimization
  7. Snoek 等人（2012）Practical Bayesian Optimization of Machine Learning Algorithms

12.4 上置信界 #

上置信界就是第 11 章已经用过的规则：

UCB⁡n(x)=μn(x)+β1/2 σn(x).\UCB_n(\vx) = \mu_n(\vx) + \beta^{1/2}\, \sigma_n(\vx).
(12.6)

它的原则来自赌博机文献（第 13.2 节），称为面对不确定性时的乐观原则（optimism in the face of uncertainty）：假定世界在合理范围内尽可能好，再由评估来纠正。如果乐观有道理，评估就找到了好的输入；否则，评估会缩小该处的 σn(x)\sigma_n(\vx)，被抬高的估计随之回落，规则便转向别处。

这一权重有概率上的解释。在后验下，P(f(x)≤μn(x)+β1/2σn(x))=Φ(β1/2)\Prob(f(\vx) \le \mu_n(\vx) + \beta^{1/2} \sigma_n(\vx)) = \Phi(\beta^{1/2})，因此式（12.6）是关于 f(x)f(\vx) 的信念的分位数。当 β1/2=2\beta^{1/2} = 2 时，它是 97.7% 分位数；最大化它，就是选出合理范围内最好情形最高的输入。

β\beta 应当取多少？Srinivas 等人（2010）给出的回答是一个随评估次数 tt 缓慢增长的取值方案。对于由候选输入构成的有限定义域 X\X 和置信参数 δ∈(0,1)\delta \in (0, 1)，他们的 GP-UCB 规则取

βt=2log⁡ ⁣(∣X∣ t2π26δ),\beta_t = 2 \log\!\left(\frac{|\X|\, t^2 \pi^2}{6 \delta}\right),
(12.7)

并证明了在这一方案下遗憾次线性增长，因而与最大值的平均差距趋于零（定理及其证明见第 13.4 节）。该方案的取值恰好大到足以保证：以至少 1−δ1 - \delta 的概率，区间 μ±βt1/2σ\mu \pm \beta_t^{1/2} \sigma 在每个候选点、每一步同时包含真实函数。对于连续定义域，方案中会多出一个与 dlog⁡td \log t 成正比的项（Srinivas 等，2010）。

具体数字显示出理论有多保守。取 ∣X∣=1000|\X| = 1000 个候选点、δ=0.1\delta = 0.1，由式（12.7）得 t=10t = 10 时 βt1/2≈5.4\beta_t^{1/2} \approx 5.4，t=100t = 100 时 ≈6.2\approx 6.2，相当于五到六个标准差的乐观。由图 11.3 可见，在贯穿全书的示例目标函数上，权重取 1 至 2 最好，取 6 左右则明显更差。原作者也发现，把 βt\beta_t 缩小为五分之一后，算法表现更好（Srinivas 等，2010）。实践中，软件库由用户指定常数 β\beta；BoTorch 的 UpperConfidenceBound 计算均值加上 β\sqrt{\beta} 倍的标准差。取值方案对证明很重要，因为证明需要永不停止的探索。在固定预算下什么常数表现良好，则是另一个实际问题。

第 12.4 节引用的文献 1
  1. Srinivas 等人（2010）Gaussian Process Optimization in the Bandit Setting: No Regret and Experimental Design

12.5 Thompson 采样 #

本章最古老的思想出自 1933 年。Thompson（1933）考虑两种疗效未知的治疗方案，提出为每位新病人分配某种治疗的概率，应等于根据已有证据该治疗较优的概率。对于高斯过程，这条规则为：

  1. 从给定 Dn\D_n 的后验中抽取一个函数 gg（第 8.5 节）。
  2. 在该样本取最大值处评估目标函数：xn+1∈arg max⁡xg(x)\vx_{n+1} \in \argmax_{\vx} g(\vx)。

第一步是唯一带随机性的一步，这条规则之所以有效，正在于此。给定数据，样本 gg 与未知的 ff 同分布，因此 gg 的最大值点与 ff 的最大值点也同分布：

P ⁣(xn+1∈A∣Dn)=P ⁣(x⋆∈A∣Dn)for every region A.\Prob\!\left(\vx_{n+1} \in A \mid \D_n\right) = \Prob\!\left(\vx^\star \in A \mid \D_n\right) \quad \text{for every region } A.
(12.8)

这称为概率匹配（probability matching）：Thompson 采样评估每个区域的频率，恰好等于模型认为最大值位于该区域的概率。模型确信较差的区域几乎从不被选中。可能藏有最大值的区域，被选中的频率与这种可能性成正比，无论模型在其余部分有多不确定。

Thompson 采样没有需要调节的权重或裕量，而且天然适合并行：要一次选出十次评估，只需抽取十个样本，分别取其最大值点。它的理论保证与上置信界接近。Russo 与 Van Roy（2014）建立了后验采样与上置信界算法之间的联系，借此可把为上置信界算法证明的遗憾界转化为后验采样的贝叶斯遗憾界，其中包括高斯过程模型的一个界；Chowdhury 与 Gopalan（2017）则在未知函数固定、而非从先验中抽取的情形下，为高斯过程版本证明了一个遗憾界。

代价在第 1 步。在 mm 个输入构成的网格上抽取一个样本，需要对 m×mm \times m 的后验协方差做 Cholesky 分解，计算量为 O(m3)O(m^3)，因此精确采样只限于几千个输入以内。维度更高时，实现上有两种做法：一是在候选集上抽样，候选集由几千个覆盖有希望区域的点组成；二是抽取可在任意位置求值的近似样本函数，这类函数由有限个随机基函数构建（随机特征），或由经数据修正的先验样本构建（路径式更新），再用梯度方法求其最大值（Rahimi 与 Recht，2007；Wilson 等，2020）。为什么在维度很多时候选集的选择会成为难点，见第 12.9 节。

第 12.5 节引用的文献 5
  1. Thompson（1933）On the Likelihood that One Unknown Probability Exceeds Another in View of the Evidence of Two Samples
  2. Russo 与 Van Roy（2014）Learning to Optimize via Posterior Sampling
  3. Chowdhury 与 Gopalan（2017）On Kernelized Multi-armed Bandits
  4. Rahimi 与 Recht（2007）Random Features for Large-Scale Kernel Machines
  5. Wilson 等人（2020）Efficiently Sampling Functions from Gaussian Process Posteriors

12.6 知识梯度 #

期望改进暗含一个假定：最终推荐的是已评估输入中观测值最好的那一个（Frazier，2018）。但很多时候，只要模型确信某个从未评估过的输入很好，也乐于推荐它。而在评估有噪声时，任何观测值本来就不能照单全收，因此自然的推荐是后验均值最高的输入（第 11.2.3 节）。

知识梯度认真对待这种推荐。它的效用是最终推荐的价值：如果现在停止，推荐的将是后验均值的最大值点，其在后验下的期望值为

u(Dn)=max⁡x′μn(x′)=:μn∗.u(\D_n) = \max_{\vx'} \mu_n(\vx') =: \mu^*_n.

在 x\vx 处多做一次评估，改变的是各处的后验均值，而不只是 x\vx 处的，因此即使 f(x)f(\vx) 本身平平无奇，也可能抬高 μ∗\mu^*。把这一效用代入式（12.1），即得知识梯度（knowledge gradient）：

KGn(x)=En ⁣[ max⁡x′μn+1(x′)  |  xn+1=x]−max⁡x′μn(x′).\mathrm{KG}_n(\vx) = \E_n\!\left[\,\max_{\vx'} \mu_{n+1}(\vx') \;\middle|\; \vx_{n+1} = \vx\right] - \max_{\vx'} \mu_n(\vx').
(12.9)

在这一意义下，一次查询的价值是再得到一个回答之后最终推荐的期望值。使其最大的查询是一步贝叶斯最优（one-step Bayes optimal）的：如果会话在这次评估后结束，任何其他选择在期望意义上都不会留下更好的推荐；知识梯度正是选出这种查询的采集函数。查询为成对选项时，第 19.4.1 节所依赖的正是这一性质。

计算式（12.9），需要知道新增一个观测时后验均值如何变化。

推导再增加一个观测后的后验均值

设下一个观测为 y=f(x)+εy = f(\vx) + \varepsilon，噪声方差为 σε2\sigma_\varepsilon^2。（全书其他地方把这一噪声方差记作 σn2\sigma_n^2；本节改用新记号，因为每个公式中后验方差 σn2(x)\sigma_n^2(\vx) 都紧挨着它。）记 kn(x′,x)k_n(\vx', \vx) 为给定 Dn\D_n 时 f(x′)f(\vx') 与 f(x)f(\vx) 的后验协方差。

  1. 给定 Dn\D_n，(f(x′),y)(f(\vx'), y) 服从联合高斯分布，均值分别为 μn(x′)\mu_n(\vx') 和 μn(x)\mu_n(\vx)，yy 的方差为 σn2(x)+σε2\sigma_n^2(\vx) + \sigma_\varepsilon^2；由于噪声与 ff 独立，协方差为 kn(x′,x)k_n(\vx', \vx)。
  2. 以 yy 为条件（第 4.5 节），得到 μn+1(x′)=μn(x′)+kn(x′,x)σn2(x)+σε2 (y−μn(x))\mu_{n+1}(\vx') = \mu_n(\vx') + \dfrac{k_n(\vx', \vx)}{\sigma_n^2(\vx) + \sigma_\varepsilon^2}\,(y - \mu_n(\vx))。
  3. 观测之前，y−μn(x)y - \mu_n(\vx) 服从均值为零、标准差为 σn2(x)+σε2\sqrt{\sigma_n^2(\vx) + \sigma_\varepsilon^2} 的高斯分布，因此可写成该标准差与标准正态变量 ZZ 的乘积。
  4. 代入第 2 步。
μn+1(x′)=μn(x′)+σ~n(x′,x) Z,σ~n(x′,x)=kn(x′,x)σn2(x)+σε2.\mu_{n+1}(\vx') = \mu_n(\vx') + \tilde\sigma_n(\vx', \vx)\, Z, \qquad \tilde\sigma_n(\vx', \vx) = \frac{k_n(\vx', \vx)}{\sqrt{\sigma_n^2(\vx) + \sigma_\varepsilon^2}}.
(12.10)

新后验均值的每一点都按同一个标准正态变量 ZZ 的某个倍数移动，因为推动它们的是同一个数，即结果 yy。与 x\vx 强相关的输入移动较多，远离 x\vx 的输入几乎不动。于是知识梯度为

KGn(x)=E ⁣[max⁡x′(μn(x′)+σ~n(x′,x) Z)]−max⁡x′μn(x′).\mathrm{KG}_n(\vx) = \E\!\left[\max_{\vx'} \left(\mu_n(\vx') + \tilde\sigma_n(\vx', \vx)\, Z\right)\right] - \max_{\vx'} \mu_n(\vx').
(12.11)

由这一形式可以得出三个推论。

  • 知识梯度非负。若干函数的最大值是凸的，且 ZZ 的均值为零，故由 Jensen 不等式，最大值的期望至少等于 Z=0Z = 0 时的最大值，即 μn∗\mu^*_n。在期望意义上，信息不会带来损失。
  • 精确的重复评估毫无价值。若 x\vx 已无噪声地评估过，则 σn(x)=0\sigma_n(\vx) = 0，且对每个 x′\vx' 都有 kn(x′,x)=0k_n(\vx', \vx) = 0，各处均值都不会移动，KGn(x)=0\mathrm{KG}_n(\vx) = 0。有噪声时（σε>0\sigma_\varepsilon > 0），重复评估仍可能有价值，这正是知识梯度能妥善处理带噪声问题的原因（Frazier，2018）。
  • 只涉及两个输入时，知识梯度即 Clark 公式。若式（12.11）中的最大值只在两个输入上取，它就是两个联合高斯值 a1+b1Za_1 + b_1 Z 与 a2+b2Za_2 + b_2 Z 的最大值，其期望由 Clark（1961）的公式给出（第 B.5 节）。同一公式将在偏好贝叶斯优化中再次出现：那里一次查询是一对选项，采集函数是较好选项的期望效用（EUBO，式（19.3））；回答无噪声时，它具有知识梯度的一步最优性（第 19.4.1 节）。

知识梯度也有助于理解期望改进。如果推荐必须是已评估的输入，且评估精确，那么式（12.1）中的效用就成为观测到的最优值，一步增益恰好就是期望改进（Frazier，2018）。期望改进即只推荐已测输入的决策者所用的知识梯度。

图 12.2 直观地展示了这种前瞻。

当前后验均值在 x 处再评估一次后的均值（7 种结果）各自的最大值x = 0.62：当前最高均值 0.52，评估后最高均值的期望 0.61，KG = 0.091−1.0−0.50.00.51.01.5f(x)当前最高均值 μ*0.00.10.20.00.20.40.60.81.0输入 x知识梯度
当前后验均值在 x 处再评估一次后的均值（7 种结果）各自的最大值x = 0.62：当前 μ* 0.52，评估后 0.61，KG = 0.091−1.0−0.50.00.51.01.5f(x)当前最高均值 μ*0.00.10.20.00.20.40.60.81.0输入 x知识梯度
图 12.2 示例目标函数上作为一步前瞻的知识梯度。用滑块或点击选择候选点 xx（橙色）。紫色曲线是在该处评估之后模型将得到的后验均值（式（12.10）），对应七个等可能的结果 ZZ（预测分布的分位数）；圆点标出每条曲线的最高处，蓝色虚线为当前的最高均值 μn∗\mu^*_n。知识梯度等于新最大值高度的期望减去 μn∗\mu^*_n，期望对所有结果计算，而不仅是画出的七个。下方横条给出每个候选点的知识梯度（精确计算），黑色三角形标出其最大值。“在 x 处评估”把候选点加入数据。

观察紫色曲线在何处散开。默认候选点位于右侧尚未探索的空隙中，各虚拟观测下的均值曲线散得很开：结果高，均值的新最大值就出现在那里；结果低，原来的最大值保持不变。低的结果无法降低 μ∗\mu^*（原来的最大值仍然可用），高的结果则会抬高它，因此平均值上升。

把候选点移到已评估的点上。紫色曲线全部收拢到蓝色曲线上，读数显示 KG 为零，即上文的第二个推论。再把候选点向任一侧稍稍移开，KG 立即回升。两个靠得很近的精确评估能揭示 ff 在两者之间的斜率，而在峰顶附近，斜率能够移动均值的最大值。横条中那些窄缺口正是这一效应。

与横条对照。知识梯度在当前最好的区域附近最大，而不是在最宽空隙的中间。在最优点附近，一次评估能直接移动均值的最大值；在远处的空隙中，结果必须很大才会产生影响。知识梯度的探索少于 β1/2=2\beta^{1/2} = 2 的上置信界，图 12.3 将再次显示这一规律。

加入噪声。把“噪声标准差 σ_n”设为 0.2。带噪声的结果对均值的影响较小，因此虚拟观测下的均值曲线散开得较少，那些缺口也变宽成了谷。在宽峰上的已评估输入处，知识梯度现在略大于零：在那里再做一次带噪声的测量，仍可能移动均值的最大值。在 x=0.45x = 0.45 附近凹陷处的已评估输入上，知识梯度仍为零，因为那里的任何结果都不能使该区域成为最好的区域。

知识梯度最初是为在信念相互独立的有限个备选项中做选择而提出的（Frazier 等，2008），后来推广到相关信念（Frazier 等，2009），即网格上高斯过程的情形。在有限集合上，式（12.11）可以精确计算：直线 aj+bjza_j + b_j z 的最大值是 zz 的分段线性函数，把每一段对正态密度积分即得闭式解（Frazier 等，2009）；图 12.2 正是这样计算的。在连续定义域上，内层最大值没有闭式解。实现中要么模拟结果，对每个结果重新求解内层最大化，以此估计知识梯度（Frazier，2018）；要么像 BoTorch 那样，在单个“一次性”（one-shot）问题中，把候选点与每个模拟结果对应的最大值点一起优化（Balandat 等，2020）。（知识梯度与期望改进在起源上如何交织在一起，见第 11.5 节。）

第 12.6 节引用的文献 5
  1. Frazier（2018）A Tutorial on Bayesian Optimization
  2. Clark（1961）The Greatest of a Finite Set of Random Variables
  3. Frazier 等人（2008）A Knowledge-Gradient Policy for Sequential Information Collection
  4. Frazier 等人（2009）The Knowledge-Gradient Policy for Correlated Normal Beliefs
  5. Balandat 等人（2020）BoTorch: A Framework for Efficient Monte-Carlo Bayesian Optimization

前面的规则都按一次评估对最终答案的作用来衡量其价值。另一类规则则按评估所提供的关于最大值的信息来衡量。ff 上的后验诱导出最大值位置的分布 p(x⋆∣Dn)p(\vx^\star \mid \D_n)：从后验中抽取一个函数，找出其最大值所在，如此重复。这个分布越分散，最大值的位置就越不确定。一次评估若预期能使这个分布更加集中，就有价值。

第 6.1 节中的熵度量分布的分散程度，它只依赖于概率，不依赖于距离。熵搜索（entropy search）选择在期望意义上使 p(x⋆∣D)p(\vx^\star \mid \D) 的熵降低最多的评估：

ESn(x)=H ⁣[p(x⋆∣Dn)]−En ⁣[H ⁣[p(x⋆∣Dn∪{(x,y)})]].\mathrm{ES}_n(\vx) = H\!\left[p(\vx^\star \mid \D_n)\right] - \E_n\!\left[H\!\left[p(\vx^\star \mid \D_n \cup \{(\vx, y)\})\right]\right].
(12.12)

熵的期望减少量就是结果 yy 与 x⋆\vx^\star 之间的互信息（第 6.3 节），而式（12.12）就是以 x⋆\vx^\star 的负熵为效用的式（12.1）。这一思想由 Villemonteix 等人（2009）提出，Hennig 与 Schuler（2012）也独立地提出了它，这一名称即出自后者（Garnett，2023，第 12.4 节）。

计算式（12.12）很困难：x⋆\vx^\star 的分布没有闭式解，其熵必须近似计算，而且这种近似要对每个假想的结果 yy 重复进行。有两种重新表述使这一思想变得实用。

预测熵搜索（predictive entropy search）利用了互信息的对称性：yy 所含关于 x⋆\vx^\star 的信息，等于 x⋆\vx^\star 所含关于 yy 的信息，因此

PESn(x)=H ⁣[p(y∣Dn,x)]−Ex⋆ ⁣[H ⁣[p(y∣Dn,x,x⋆)]].\mathrm{PES}_n(\vx) = H\!\left[p(y \mid \D_n, \vx)\right] - \E_{\vx^\star}\!\left[H\!\left[p(y \mid \D_n, \vx, \vx^\star)\right]\right].

第一项是高斯分布的熵，有闭式解；第二项对采样得到的最大值点取平均，每个最大值点都需要近似计算已知 x⋆\vx^\star 会如何改变 x\vx 处的预测（Hernández-Lobato 等，2014）。精确的熵搜索与预测熵搜索是同一个函数，区别只在于近似方法（Frazier，2018）。

最大值熵搜索（max-value entropy search，MES）换了目标：它寻求的是关于最大值本身 f⋆=f(x⋆)f^\star = f(\vx^\star) 这一个数的信息，而不是最大值的位置（Wang 与 Jegelka，2017）。已知 f⋆f^\star，关于 f(x)f(\vx) 只能得到一个简单的结论：它不能超过 f⋆f^\star。因此，给定 f⋆f^\star 时关于 f(x)f(\vx) 的信念，就是在 f⋆f^\star 处从上方截断的高斯后验，即截断高斯分布，其熵有闭式解。对 KK 个采样得到的最大值 f1⋆,…,fK⋆f^\star_1, \dots, f^\star_K 取平均，得到

MESn(x)≈1K∑k=1K[γk(x) ϕ(γk(x))2 Φ(γk(x))−log⁡Φ(γk(x))],γk(x)=fk⋆−μn(x)σn(x),\mathrm{MES}_n(\vx) \approx \frac{1}{K} \sum_{k=1}^{K} \left[\frac{\gamma_k(\vx)\, \phi(\gamma_k(\vx))}{2\, \Phi(\gamma_k(\vx))} - \log \Phi(\gamma_k(\vx))\right], \qquad \gamma_k(\vx) = \frac{f^\star_k - \mu_n(\vx)}{\sigma_n(\vx)},
(12.13)

即 Wang 与 Jegelka（2017）的式（6）。样本 fk⋆f^\star_k 可以取后验样本函数的最大值，也可以从其分布的更廉价近似中抽取。当 γk\gamma_k 较小，即以标准差计，采样得到的最大值比 x\vx 处的均值高出不多时，每一项都较大：在该处评估，可以揭示最大值是否大致就这么高。作者报告，最大值熵搜索的效果与熵搜索和预测熵搜索相当或更好，成本却只是后两者的一小部分，而且对样本数远没有那么敏感（Wang 与 Jegelka，2017）。

第 12.7 节引用的文献 6
  1. Villemonteix 等人（2009）An Informational Approach to the Global Optimization of Expensive-to-Evaluate Functions
  2. Hennig 与 Schuler（2012）Entropy Search for Information-Efficient Global Optimization
  3. Garnett（2023）Bayesian Optimization
  4. Hernández-Lobato 等人（2014）Predictive Entropy Search for Efficient Global Optimization of Black-box Functions
  5. Frazier（2018）A Tutorial on Bayesian Optimization
  6. Wang 与 Jegelka（2017）Max-value Entropy Search for Efficient Bayesian Optimization

12.8 相互比较 #

至此，每条规则都已单独推导。为了看清它们的行为有多大差别，图 12.3 把全部六条规则放在示例目标函数的同一个后验上。读者可以充当优化器：点击图中任意位置，即在该处评估目标函数，然后观察每条规则下一步会选在何处。

隐藏的目标函数后验均值95% 区间p(x*)：最大值可能的位置EI：下一步评估 x = 0.19−1.0−0.50.00.51.01.5f(x)PIEIUCBThompsonKGMES0.00.20.40.60.81.0输入 x
隐藏的目标函数后验均值95% 区间p(x*)：最大值可能的位置EI：下一步评估 x = 0.19−1.0−0.50.00.51.01.5f(x)PIEIUCBTSKGMES0.00.20.40.60.81.0输入 x
图 12.3 同一后验上的六个采集函数。上：示例目标函数（虚线）、迄今评估之后的后验（圆点为评估点），以及紫色条形，表示 64 个后验样本中每个输入成为最大值点的次数，是 p(x⋆∣Dn)p(x^\star \mid \mathcal{D}_n) 的估计。点击图中某处即在该处评估目标函数，点击圆点可将其删除。下：改进概率 PI 与期望改进 EI（式（12.2）、式（12.5），裕量 ξ\xi）、上置信界 UCB（式（12.6），权重 β1/2\beta^{1/2}）、一个 Thompson 样本、知识梯度 KG（式（12.11），精确计算），以及最大值熵搜索 MES（式（12.13），由 64 个样本的最大值得到）。每条横条都经过缩放以填满其高度，因此只有其形状与最大值位置（圆点）有意义。所选规则以橙色显示，它的下一个查询是后验图上的橙色线。“评估其最大值处”用所选规则运行一步循环；“新样本”重新抽取 Thompson 采样与最大值熵搜索所用的随机样本。

默认的后验有四次评估：两次在宽峰上，一次在其后的凹陷中，一次在最右端。x=0.73x = 0.73 附近的高峰位于尚未探索的空隙中。可以做以下几个实验：

各规则意见不一。改进概率、期望改进、最大值熵搜索与知识梯度都选在宽峰附近，xx 约在 0.19 至 0.27 之间；那里的均值已接近当前最优值，很可能获得适度的改进。β1/2=2\beta^{1/2} = 2 的上置信界则进入空隙，选在 x≈0.68x \approx 0.68 处，那里可信带的上沿最高。Thompson 采样取决于所抽的样本：按几次“新样本”，它的选择会在宽峰、空隙与左端之间跳动，频率大致与紫色条形成比例。

增大裕量。把 ξ\xi 从 0.01 往上调。ξ=0.1\xi = 0.1 时，期望改进已经转向空隙：要求至少 0.1 的改进，宽峰附近那些幅度小而很可能出现的收益便失去了价值。改进概率则一直停留在宽峰，直到 ξ\xi 达到 0.39。两条规则都会在某个裕量处突然改变选择，这正是 Jones 所警告的敏感性。

运行循环。选一条规则，反复按“评估其最大值处”，然后按“重置”再试另一条。每条规则最终都会到达高峰，只是路线不同：上置信界最先到达；期望改进、知识梯度与最大值熵搜索在宽峰上再评估一次之后到达；改进概率再评估两次之后到达；Thompson 采样更晚。在这个后验上，改进概率碰巧最快进入距最大值 0.05 以内；知识梯度追求的是提高均值的最大值，而非观测到的最优值，因此按观测最优值衡量，它最慢。一个问题上的一次运行不能说明任何排名。

观察紫色条形。高峰一经评估，x⋆x^\star 的估计分布便收缩到高峰上。Thompson 采样与最大值熵搜索随后把评估集中在那里；而只要其他区域的可信带仍宽，上置信界就会继续访问这些区域。

没有哪条规则在所有问题上都最好。在基准函数上的比较中，不同问题偏向不同的规则；第 13.5.2 节中的遗憾界针对不同规则、在不同设定下证明，因此也无法为这些规则排序（推断）。表 12.2 列出了确实成立的实际差别。

表 12.2 本章各采集函数的实际性质。
规则 对高斯过程有闭式解？ 需设定的参数 评估的价值依据 每个候选点的成本
PI 有，式（12.2） 裕量 ξ\xi，敏感 超过当前最优值的概率 一次预测
EI 有，式（12.5） 裕量 ξ\xi，常取 0 或很小 相对当前最优值的期望增益 一次预测
UCB 有，式（12.6） 权重 β\beta，影响大 乐观值 一次预测
Thompson 采样 无；为随机样本 无 成为最大值点的概率 在所有候选点上的一次联合采样
KG 在有限集合上有，式（12.11） 无 推荐值的增益 每个结果一次最大化
ES、PES 无 无 关于 x⋆\vx^\star 的信息 代价高昂的近似
MES 给定采样得到的最大值时有，式（12.13） 样本数 关于 f⋆=f(x⋆)f^\star = f(\vx^\star) 的信息 每个样本一次预测

12.8.1 二维情形 #

一维的图像掩盖了一个重要事实：需要查看的位置随维度呈指数增长，而后验只在数据附近才有信息量。图 12.4 在一个二维问题上重复这一比较。该问题是 Branin 函数，一个有三个等高全局最大值的标准测试函数（这里已缩放到单位正方形并取负，使值越大越好）。

6 次评估；最优值 0.62（最大值 0.99）；EI：选择 (0.63, 0.37)隐藏的目标函数（3 个最大值点，+）+++后验均值后验标准差EI
6 次评估；最优值 0.62（最大值 0.99）；EI：选择 (0.63, 0.37)隐藏的目标函数（3 个最大值点，+）+++后验均值后验标准差EI
图 12.4 二维 Branin 函数上的采集函数，已有六次初始评估（圆点）。左上：隐藏的目标函数，三个全局最大值以 + 标出；低于 −1.5 的值（有一角降到约 −4.7）以最浅的颜色绘制。右上：同一色标下的后验均值。左下：后验标准差。右下：所选的采集函数（Thompson 采样时为一个后验样本）。各面板中颜色越深，值越高。每个面板上的橙色圆环都是所选规则的下一个查询点。点击任一面板即在该处评估目标函数；“评估最大值处”由规则选择。网格有 31 × 31 个点，同时作为候选点。

比较下方两个面板。标准差只在六次评估周围的小圆盘内较低，几乎整个正方形都不确定。期望改进在较高的均值与较大的标准差相遇处较大，这里是下半部分各评估点之间的一片宽阔区域；在评估点处以及均值较低的右上方，期望改进都较小。

在规则之间切换。改进概率与知识梯度选在最好的评估点附近，因为在那里评估很可能抬高当前最优值或均值的最大值。上置信界与最大值熵搜索探得更远。Thompson 采样的样本是一整张曲面，在未探索的角落里有其自身的峰；其最大值点可能落在模型允许出现最大值的任何地方。

运行循环。选择“PI”，按十次“评估最大值处”。评估点从最好的初始点出发，以小步沿斜坡向下移动，到达 (0.54,0.15)(0.54, 0.15) 附近的最大值：这正是第 12.2 节中那种谨慎的行为，在这里碰巧奏效。按“重置”后，用“EI”或“UCB”做同样的操作。它们有几次评估落在正方形的边缘和角落，因为那些位置之外没有评估点，标准差始终很大；其余评估落在三个最大值中的两个附近。三个最大值等高，一次运行先找到哪一个，取决于最初的几次评估。

设想六维的情形。图中不确定区域占正方形的大部分，31 × 31 的候选网格能把它细密地覆盖。同样的网格在六维中有 316≈9×10831^6 \approx 9 \times 10^8 个点。采集函数仍然只在数据附近有信息量，而这部分如今只占体积的极小一部分。如何找到它的最大值，是下一节的主题。

12.9 采集函数的优化 #

本章的每条规则最后都归结为“在采集函数最大处评估”。前面的图通过检查网格上的每个点来做到这一点，这在一维或二维中可行，在十维中则不可能。最大化采集函数本身就是一个全局优化问题，而循环每一步都要求解一个。

这一问题之所以可解，在于代价低。计算一次期望改进或上置信界只需一次后验预测：每步做一次 Cholesky 分解之后，均值的计算量为 O(n)O(n)，方差为 O(n2)O(n^2)（第 8.4 节）。观测有几百个时，这只需几微秒，因此优化器每一步可以承受数万次采集函数求值，而耗时数小时的目标函数只评估一次。采集函数也是光滑的，对于高斯过程还能以闭式求导，因此可以使用梯度方法（Frazier，2018）。

这一问题之所以困难，在于形状。采集函数是非凸的，有许多局部极大值，每个值得关注的区域附近都有一个或多个。更麻烦的是，它们在远离数据处几乎是平的：使用平稳核时，在远离所有观测的地方，后验回到先验，后验的均值和标准差不再变化，采集函数及其梯度也随之不再变化（Garnett，2023，第 9.2 节）。在高维中，几乎整个定义域都远离数据，因此从随机点出发的梯度方法通常只得到几乎为零的梯度，停在原地。

标准的解决办法是多起点局部优化，Garnett 的教科书推荐这一方法（Garnett，2023，第 9.2 节），BoTorch 也采用了它（Balandat 等，2020）。

算法 12.1 用多起点梯度上升最大化采集函数

输入：带梯度的采集函数 ana_n；定义域 [0,1]d[0, 1]^d；数目 R≪MR \ll M。

  1. 筛选。在 MM 个拟随机点上计算 ana_n，例如加扰 Sobol 序列的前 MM 个点（第 11.4 节）。如果采集函数在别处很可能是平的，再加入一些最好观测附近的点。
  2. 选择。从中选出 RR 个起点，偏向取值高的点，同时保留一定的多样性。
  3. 爬升。从每个起点出发，运行满足箱形约束的基于梯度的局部优化器，例如 L-BFGS-B。（这是一种拟 Newton 法：用相继的梯度估计 ana_n 的曲率，而不计算二阶导数。）
  4. 返回找到的最佳局部极大值。

BoTorch 的默认设置遵循这一框架：筛选点取自加扰 Sobol 序列；起点按与 exp⁡(ηZ)\exp(\eta Z) 成正比的概率随机抽取，其中 ZZ 为标准化后的采集函数值，η\eta 为温度；局部爬升使用 L-BFGS-B。各次爬升相互独立，因此易于并行（Garnett，2023，第 9.2 节）。

采集函数的蒙特卡洛版本不使用闭式解，而是对样本取平均来估计期望，这类版本还需要另一个想法。Wilson 等人（2018）表明，如果把样本写成固定随机数的固定变换，蒙特卡洛估计就是输入的光滑函数，同样可以用梯度上升求解。他们还找出了一族采集函数（包括期望改进和上置信界），其性质保证可以贪心地构建一批评估：一次选一个点，每个点都在已选点给定的条件下取最大值。批量选择是第 14.3 节的主题。

12.9.1 从二维到二十维 #

图 12.4 中的网格有 961 个点，没有遗漏任何地方。当维度增长到典型调参问题的 6 至 20 维时，有三件事会发生变化。

网格不再可行，均匀筛选也是如此。每个输入取十个值的网格，在六维中有 10610^6 个点，在二十维中有 102010^{20} 个。均匀随机或 Sobol 筛选点不需要网格，却以另一种形式继承了网格的弱点：它们大多落在远离数据、采集函数平坦的地方。可以这样理解：若一次评估在约 0.20.2 的半径内为模型提供信息，那么一次评估只影响六维单位立方体的 0.033%（即第 11.3 节背后的计算）。有 60 次评估时，均匀放置的筛选点中至多约 2% 落在某个评估点的这一半径之内，在其余各点处，采集函数基本是常数。

图 12.5 在 Hartmann-6 上直接测量了这一点。Hartmann-6 是一个标准的六维测试函数，d>6d > 6 时混在若干无关输入之中。

d = 20，50 次评估：最好的均匀候选点 EI 为 0.0015；最好的局部候选点为 0.2（是前者的 130 倍）到最近评估点的距离（以长度尺度计）在立方体中均匀分布最好的 5 次评估附近01234期望改进（对数刻度）在立方体中均匀分布最大值 0.0015最好的 5 次评估附近最大值 0.210⁻¹²10⁻⁹10⁻⁶10⁻³1
d = 20，50 次评估。最大 EI：均匀 0.0015，局部 0.2（是前者的 130 倍）到最近评估点的距离（以长度尺度计）在立方体中均匀分布最好的 5 次评估附近01234期望改进（对数刻度）在立方体中均匀分布最大值 0.0015最好的 5 次评估附近最大值 0.210⁻¹²10⁻⁹10⁻⁶10⁻³1
图 12.5 在 2 至 20 维中寻找期望改进最大值的两种方式。对 Hartmann-6 做 10+2d10 + 2d 次拉丁超立方评估之后（d=2d = 2 时取穿过其最大值的切片；d>6d > 6 时加入 d−6d - 6 个无关输入），拟合长度尺度为 0.2d0.2\sqrt{d} 的高斯过程，再在两组点上计算期望改进：从立方体中均匀抽取的 1,000 个点（灰色），以及对最好的五次评估做随机扰动（每个坐标的标准差为 0.05）得到的 1,000 个点（橙色）。左：每个候选点到最近评估点的距离，以长度尺度为单位。右：对数刻度下的期望改进值，并标出每组的最大值；低于 10−1210^{-12} 的值计入最左边的区间。期望改进以标准化观测为单位。“重新抽取”更换评估点与候选点；数值仅作示意。

从 d=2d = 2 切换到 d=20d = 20。在二维中，最好的均匀候选点与数据的最好扰动一样好：两组点找到的是期望改进的同一个峰。在二十维中，按默认抽取，1,000 个均匀候选点中最好的那个，其期望改进不到最好扰动点的百分之一。长度尺度随 d\sqrt{d} 增长，这正是 Hvarfner 等人（2024）在高维中缩放长度尺度先验的速率；然而在二十维中，典型的均匀候选点到最近评估点仍约有 1.4 个长度尺度，而在二维中约为 0.4（左面板）。在那样的位置，后验均值平平无奇，当前最优值远在几个标准差之外，期望改进极小且几乎不变，因此算法 12.1 的筛选步骤会从没有信息量的点开始爬升。

候选点取自数据附近。因此，实用的实现把候选点放在采集函数有结构的地方。一种做法是把整个立方体上的拟随机样本与对迄今最好输入的随机扰动混合起来，每次只扰动部分坐标。TuRBO 是一种针对高维问题的方法，它在信赖域（trust region）内搜索，即围绕迄今最优点的一个箱形区域（第 14.6.2 节）。它的 Thompson 样本就在按这种方式构建、含 min⁡(100d,5000)\min(100d, 5000) 个点的候选集上抽取：候选点的每个坐标以 min⁡(1,20/d)\min(1, 20/d) 的概率在信赖域内取拟随机值，否则保持区域中心的值（Eriksson 等，2019）。正是这种扰动，使候选点在 dd 增大时仍留在有信息量的区域内。同样的想法也体现在 BoTorch 的一个选项中：在最好的输入周围采样并添加点，作为算法 12.1 第 1 步的筛选点。

数值上的平坦成为常态。维度很多、观测也很多时，期望改进在浮点运算中能与零区分开的区域会缩小，普通的期望改进交给优化器的，是一个几乎处处恰好为零、梯度也为零的函数。正是在这种情形下，计算 log⁡EI⁡\log \EI 而非 EI⁡\EI 的收益最大（Ament 等，2023）。

这些都没有改变采集函数的统计性质，改变的是能否找到它的最大值，而这在实践中很重要。第 13 章中的理论保证假定最大化是精确的（Srinivas 等，2010），而 Ament 等人（2023）发现，仅仅改进最大化，就改变了期望改进与较新采集函数的比较结果。在采信一项已发表的比较之前，值得先检查它是否使用了好的内层优化器（推断）。高维问题中代理模型方面的讨论见第 14.6 节，在那里，默认的长度尺度先验与采集函数同样重要。

第 12.9 节引用的文献 8
  1. Frazier（2018）A Tutorial on Bayesian Optimization
  2. Garnett（2023）Bayesian Optimization
  3. Balandat 等人（2020）BoTorch: A Framework for Efficient Monte-Carlo Bayesian Optimization
  4. Wilson 等人（2018）Maximizing Acquisition Functions for Bayesian Optimization
  5. Hvarfner 等人（2024）Vanilla Bayesian Optimization Performs Great in High Dimensions
  6. Eriksson 等人（2019）Scalable Global Optimization via Local Bayesian Optimization
  7. Ament 等人（2023）Unexpected Improvements to Expected Improvement for Bayesian Optimization
  8. Srinivas 等人（2010）Gaussian Process Optimization in the Bandit Setting: No Regret and Experimental Design

12.10 习题 #

习题 12.1

设 EI⁡(δ,s)=δ Φ(δ/s)+s ϕ(δ/s)\EI(\delta, s) = \delta\,\Phi(\delta/s) + s\,\phi(\delta/s)，见式（12.4）。证明 ∂EI⁡/∂δ=Φ(δ/s)\partial \EI / \partial \delta = \Phi(\delta/s)，∂EI⁡/∂s=ϕ(δ/s)\partial \EI / \partial s = \phi(\delta/s)，并由此得出期望改进随均值和标准差的增大而增大。再计算改进概率 Φ(δ/s)\Phi(\delta/s) 对 ss 的导数，并说明它何时为负。

解答

利用 ϕ′(z)=−z ϕ(z)\phi'(z) = -z\,\phi(z)，并记 z=δ/sz = \delta/s。对 δ\delta：∂δ[δ Φ(z)]=Φ(z)+δ ϕ(z)/s=Φ(z)+z ϕ(z)\partial_\delta[\delta\,\Phi(z)] = \Phi(z) + \delta\,\phi(z)/s = \Phi(z) + z\,\phi(z)，且 ∂δ[s ϕ(z)]=s ϕ′(z)/s=−z ϕ(z)\partial_\delta[s\,\phi(z)] = s\,\phi'(z)/s = -z\,\phi(z)。两者之和为 Φ(z)\Phi(z)。对 ss，由 ∂z/∂s=−δ/s2=−z/s\partial z / \partial s = -\delta/s^2 = -z/s：∂s[δ Φ(z)]=−δ ϕ(z) z/s=−z2ϕ(z)\partial_s[\delta\,\Phi(z)] = -\delta\,\phi(z)\,z/s = -z^2\phi(z)，且 ∂s[s ϕ(z)]=ϕ(z)+s ϕ′(z)(−z/s)=ϕ(z)+z2ϕ(z)\partial_s[s\,\phi(z)] = \phi(z) + s\,\phi'(z)(-z/s) = \phi(z) + z^2\phi(z)。两者之和为 ϕ(z)\phi(z)。Φ\Phi 与 ϕ\phi 都为正，所以期望改进对两个自变量都是递增的。对于改进概率，∂sΦ(z)=ϕ(z) (−z/s)=−(δ/s2) ϕ(δ/s)\partial_s \Phi(z) = \phi(z)\,(-z/s) = -(\delta/s^2)\,\phi(\delta/s)，它恰好在 δ>0\delta > 0 时为负：一旦均值高于目标，改进概率就偏好更小的不确定性。

习题 12.2

两个输入的后验信念相互独立，分别为 f(a)∼N(0.5,0.12)f(a) \sim \N(0.5, 0.1^2) 和 f(b)∼N(0.3,0.42)f(b) \sim \N(0.3, 0.4^2)，且不能推荐其他输入。只允许做一次精确评估。用式（12.11）计算评估 bb 的知识梯度，并证明它等于 bb 相对于值 0.50.5 的期望改进。为什么评估 aa 的知识梯度如此之小？如果 aa 与 bb 相关，这一点是否仍然成立？

解答

精确评估 bb 即得知 f(b)f(b)，因此 σ~(b,b)=σ(b)=0.4\tilde\sigma(b, b) = \sigma(b) = 0.4；由独立性，σ~(a,b)=0\tilde\sigma(a, b) = 0。新的均值对 aa 为 0.50.5，对 bb 为 0.3+0.4Z0.3 + 0.4Z，所以 KG(b)=E[max⁡(0.5,0.3+0.4Z)]−0.5=E[max⁡(0.3+0.4Z−0.5,0)]\mathrm{KG}(b) = \E[\max(0.5, 0.3 + 0.4Z)] - 0.5 = \E[\max(0.3 + 0.4Z - 0.5, 0)]，即 δ=−0.2\delta = -0.2、s=0.4s = 0.4 的高斯变量正部的期望。由式（12.4），它等于 −0.2 Φ(−0.5)+0.4 ϕ(−0.5)≈−0.2×0.3085+0.4×0.3521≈0.079-0.2\,\Phi(-0.5) + 0.4\,\phi(-0.5) \approx -0.2 \times 0.3085 + 0.4 \times 0.3521 \approx 0.079，也就是 bb 相对于当前最优值 0.50.5 的期望改进。评估 aa 只移动 aa 的均值：max⁡(0.5+0.1Z,0.3)\max(0.5 + 0.1Z, 0.3) 的期望为 0.5+E[max⁡(0.3−0.5−0.1Z,0)]=0.5+E[max⁡(−0.2+0.1Z′,0)]0.5 + \E[\max(0.3 - 0.5 - 0.1Z, 0)] = 0.5 + \E[\max(-0.2 + 0.1Z', 0)]，其中 Z′=−ZZ' = -Z，它等于 0.5+(−0.2 Φ(−2)+0.1 ϕ(−2))0.5 + (-0.2\,\Phi(-2) + 0.1\,\phi(-2))，小于 0.5+0.0010.5 + 0.001。所以 KG(a)<0.001\mathrm{KG}(a) < 0.001：评估 aa 几乎从不改变推荐，因为 aa 必须跌到其均值以下两个标准差，才会输给 bb。如果存在正相关，评估 aa 也会移动 bb 的均值，知识梯度会把这一点计算在内；知识梯度正是这样按评估对整个后验的影响来衡量其价值的。

习题 12.3

两个输入的信念相互独立：f(a)∼N(μa,sa2)f(a) \sim \N(\mu_a, s_a^2)，f(b)∼N(μb,sb2)f(b) \sim \N(\mu_b, s_b^2)。证明 Thompson 采样以概率 Φ ⁣((μa−μb)/sa2+sb2)\Phi\!\left((\mu_a - \mu_b) / \sqrt{s_a^2 + s_b^2}\right) 选择 aa，并验证这就是 aa 为较好输入的后验概率。当两个标准差都缩小、而 μa>μb\mu_a > \mu_b 保持不变时，这一概率如何变化？

解答

Thompson 采样独立地抽取 g(a)∼N(μa,sa2)g(a) \sim \N(\mu_a, s_a^2) 和 g(b)∼N(μb,sb2)g(b) \sim \N(\mu_b, s_b^2)，若 g(a)>g(b)g(a) > g(b) 就选择 aa。差值 g(a)−g(b)g(a) - g(b) 服从均值为 μa−μb\mu_a - \mu_b、方差为 sa2+sb2s_a^2 + s_b^2 的高斯分布，所以 P(g(a)>g(b))=Φ((μa−μb)/sa2+sb2)\Prob(g(a) > g(b)) = \Phi((\mu_a - \mu_b)/\sqrt{s_a^2 + s_b^2})。把 gg 换成 ff 做同样的计算，就得到 f(a)>f(b)f(a) > f(b) 的后验概率，这就是式（12.8）在这种情形下的形式。当标准差缩小时，Φ\Phi 的自变量无界增大，概率趋于 1：Thompson 采样停止探索 bb 的速度，恰好与证据排除它的速度一样快。

习题 12.4

证明最大化 μn(x)+β1/2σn(x)\mu_n(\vx) + \beta^{1/2}\sigma_n(\vx) 等价于最大化关于 f(x)f(\vx) 的后验信念的 qq 分位数，并对 β1/2=1\beta^{1/2} = 1、22 以及第 12.4 节中的 βt1/2≈6.2\beta_t^{1/2} \approx 6.2 求出 qq。最后这个值说明，理论预期真实函数超出区间的频率有多高？

解答

对于均值为 μ\mu、标准差为 σ\sigma 的高斯信念，qq 分位数为 μ+Φ−1(q) σ\mu + \Phi^{-1}(q)\,\sigma。所以 μ+β1/2σ\mu + \beta^{1/2}\sigma 就是 q=Φ(β1/2)q = \Phi(\beta^{1/2}) 的分位数，在每个输入处的 qq 都相同，最大化其中一个就是最大化另一个。对 β1/2=1\beta^{1/2} = 1，q≈0.841q \approx 0.841；对 22，q≈0.977q \approx 0.977；对 6.26.2，qq 与 1 相差约 3×10−103 \times 10^{-10}。理论把区间设得如此之宽，是为了以高概率让函数在每一个候选点、每一步同时都留在区间内，这就要求每个候选点、每一步的失效概率都极小。这就是式（12.7）背后的联合界（多个事件中任何一个发生的概率，至多为它们各自概率之和；第 13.4.3 节用到了它），也是取值方案随 log⁡∣X∣\log |\X| 和 log⁡t\log t 增长的原因。

延伸阅读 #

参考文献

  1. Ament, S., Daulton, S., Eriksson, D., Balandat, M., and Bakshy, E. (2023). Unexpected Improvements to Expected Improvement for Bayesian Optimization. Advances in Neural Information Processing Systems 36 (NeurIPS 2023). 引用于 §12.3 §12.9
  2. Balandat, M., Karrer, B., Jiang, D. R., Daulton, S., Letham, B., Wilson, A. G., and Bakshy, E. (2020). BoTorch: A Framework for Efficient Monte-Carlo Bayesian Optimization. Advances in Neural Information Processing Systems 33 (NeurIPS 2020). 引用于 §12.6 §12.9
  3. Brochu, E., Cora, V. M., and de Freitas, N. (2010). A Tutorial on Bayesian Optimization of Expensive Cost Functions, with Application to Active User Modeling and Hierarchical Reinforcement Learning. arXiv preprint. 预印本引用于 §12.2 §12.3
  4. Chowdhury, S. R., and Gopalan, A. (2017). On Kernelized Multi-armed Bandits. International Conference on Machine Learning. 引用于 §12.5
  5. Clark, C. E. (1961). The Greatest of a Finite Set of Random Variables. Operations Research. 引用于 §12.6
  6. Eriksson, D., Pearce, M., Gardner, J., Turner, R. D., and Poloczek, M. (2019). Scalable Global Optimization via Local Bayesian Optimization. Advances in Neural Information Processing Systems 32 (NeurIPS 2019). 引用于 §12.9
  7. Frazier, P. I. (2018). A Tutorial on Bayesian Optimization. arXiv. 预印本引用于 §12.1 §12.3 §12.6 §12.7 §12.9
  8. Frazier, P. I., Powell, W. B., and Dayanik, S. (2008). A Knowledge-Gradient Policy for Sequential Information Collection. SIAM Journal on Control and Optimization. 引用于 §12.6
  9. Frazier, P., Powell, W., and Dayanik, S. (2009). The Knowledge-Gradient Policy for Correlated Normal Beliefs. INFORMS Journal on Computing. 引用于 §12.6
  10. Garnett, R. (2023). Bayesian Optimization. Cambridge University Press. 引用于 §12.1 §12.2 §12.3 §12.7 §12.9
  11. Hennig, P., and Schuler, C. J. (2012). Entropy Search for Information-Efficient Global Optimization. Journal of Machine Learning Research. 引用于 §12.7
  12. Hernández-Lobato, J. M., Hoffman, M. W., and Ghahramani, Z. (2014). Predictive Entropy Search for Efficient Global Optimization of Black-box Functions. Advances in Neural Information Processing Systems 27 (NeurIPS 2014). 引用于 §12.7
  13. Hvarfner, C., Hellsten, E. O., and Nardi, L. (2024). Vanilla Bayesian Optimization Performs Great in High Dimensions. International Conference on Machine Learning. 引用于 §12.9
  14. Jones, D. R. (2001). A Taxonomy of Global Optimization Methods Based on Response Surfaces. Journal of Global Optimization. 引用于 §12.2
  15. Jones, D. R., Schonlau, M., and Welch, W. J. (1998). Efficient Global Optimization of Expensive Black-Box Functions. Journal of Global Optimization. 引用于 §12.3
  16. Kushner, H. J. (1964). A New Method of Locating the Maximum Point of an Arbitrary Multipeak Curve in the Presence of Noise. Journal of Basic Engineering. 引用于 §12.1 §12.2
  17. Močkus, J. (1975). On Bayesian Methods for Seeking the Extremum. Optimization Techniques IFIP Technical Conference. 引用于 §12.3
  18. Rahimi, A., and Recht, B. (2007). Random Features for Large-Scale Kernel Machines. Advances in Neural Information Processing Systems 20 (NeurIPS 2007). 引用于 §12.5
  19. Russo, D., and Van Roy, B. (2014). Learning to Optimize via Posterior Sampling. Mathematics of Operations Research. 引用于 §12.5
  20. Snoek, J., Larochelle, H., and Adams, R. P. (2012). Practical Bayesian Optimization of Machine Learning Algorithms. Advances in Neural Information Processing Systems 25 (NeurIPS 2012). 引用于 §12.3
  21. Srinivas, N., Krause, A., Kakade, S. M., and Seeger, M. (2010). Gaussian Process Optimization in the Bandit Setting: No Regret and Experimental Design. ICML 2010. 引用于 §12.4 §12.9
  22. Thompson, W. R. (1933). On the Likelihood that One Unknown Probability Exceeds Another in View of the Evidence of Two Samples. Biometrika. 引用于 §12.5
  23. Villemonteix, J., Vazquez, E., and Walter, E. (2009). An Informational Approach to the Global Optimization of Expensive-to-Evaluate Functions. Journal of Global Optimization. 引用于 §12.7
  24. Wang, Z., and Jegelka, S. (2017). Max-value Entropy Search for Efficient Bayesian Optimization. Proceedings of the 34th International Conference on Machine Learning (ICML 2017). 引用于 §12.7
  25. Wilson, J. T., Hutter, F., and Deisenroth, M. P. (2018). Maximizing Acquisition Functions for Bayesian Optimization. Advances in Neural Information Processing Systems 31 (NeurIPS 2018). 引用于 §12.9
  26. Wilson, J. T., Borovitskiy, V., Terenin, A., Mostowski, P., and Deisenroth, M. P. (2020). Efficiently Sampling Functions from Gaussian Process Posteriors. Proceedings of the 37th International Conference on Machine Learning (ICML 2020). 引用于 §12.5