6.2-罗宾斯-门罗算法
随机近似 (stochastic approximation)指的是一大类求解方程或者优化问题的随机迭代算法24。与许多其他算法如基于梯度的算法(梯度上升或梯度下降)等优化方法相比，随机近似具有显著的优势，因为它不需要目标函数或导数的解析表达式。
罗宾斯-门罗(Robbins-Monro，RM)算法是随机近似领域的一项开创性工作[24,25,26,27]]。著名的随机梯度下降算法就是RM算法的一种特殊形式(\(6.4\)节)。
下面我们介绍\(RM\)算法，假设我们想要求解如下方程：
其中\(w\in\mathbb{R}\)是未知变量，\(g:\mathbb{R} \rightarrow \mathbb{R}\)是一个函数。许多问题都可以表述为上面方程的求根问题。例如，如果\(J(w)\)是一个需要优化的目标函数，这个优化问题可以转换为求解\(g(w)=\nabla_wJ(w)=0.\)。另外，如果方程右边有一个常数，例如\(g(w) = c\)，这样的方程仍然可以通过将\(g(w)-c\)转换为上述的形式。
Note
梯度等于\(0\)是\(J(w)\)达到最大或最小的一个必要条件。
如果\(g\)的表达式是已知的，那么可以利用很多数值算法来求解\(g(w)=0\)。然而，我们面临的问题是函数\(g\)的表达式可能是未知的。例如，函数可以用人工神经网络来表示，而人工神经网络的结构和参数都是未知的。具体来说，当输入\(w\)时，我们只能获得输出值，并且该输出值还可能包含观测噪声：
其中\(\eta\in\mathbb{R}\)是观测噪声，而且可能不是高斯分布的。这个函数是一个一个黑盒系统，只有输入\(w\)和噪声观测值\(\tilde{g}(w,\eta)\)是已知的(见图\(6.2\))。我们的目标是利用\(w\)和\(\tilde{g}\)求解\(g(w) = 0\)。
图\(6.2\)：由\(w\)和\(\tilde{g}\)求解 \(g(w) = 0\)问题的示意图。
求解\(g(w) = 0\)的\(RM\)算法是
其中，\(w_k\)是第\(k\)个方程解的估计值，\(\tilde{g}(w_k, \eta_k)\)是输入\(w_k\)后输出的观测值，\(a_k\)是一个正系数。可以看到，\(RM\)算法不需要关于函数表达式的任何信息。而仅需要输入和输出。
图\(6.3\)：RM 算法示例。
为了分析该算法的理论性质，我们先看一个\(g(w) = w^3 - 5\)的例子以直观理解\(RM\)算法。其真实的解为\(5^{1/3} ≈ 1.71\)。假设现在我们只能观察到输入\(w\)和输出 \(\tilde{g}(w) = g(w) + \eta\)，其中\(\eta\)是独立同分布的，并且服从均值为\(0\)且标准差为\(1\)的正态分布，初始值为\(w_1 = 0\)，系数是\(a_k = 1/k\)。\(w_k\)的变化过程如图\(6.3\)所示。尽管观测结果受到噪声\(\eta_k\)的干扰，但估计的\(w_k\)仍然可以收敛到真实的解。需要注意的是，如果初始值\(w_1\)选择的不合适，\(RM\)算法可能会发散。\(RM\)算法收敛性在下面将详细介绍。
6.2.1 收敛特性¶
为什么\((6.5)\)中的\(RM\)算法可以找到\(g(w) = 0\)的根呢？下面我们用一个例子来直观地解释，之后进行严格的收敛性分析。
请看图\(6.4\)所示的例子，\(g(w) = tanh(w-1)\)。\(g(w) = 0\)的真实解是\(w^∗ = 1\)。我们应用\(RM\)算法看能否成功得到这个解，选取\(w_1 = 3\)和\(a_k = 1/k\)。为了更好地说明收敛的原因，我们考虑最简单的没有观测噪声的情况：\(\eta_k ≡ 0\)，即\(\tilde{g} (w_k, \eta_k) = g(w_k)\)。此时，\(RM\)算法变为\(w_{k+1} = w_k -a_k g(w_k)\)。由\(RM\)算法生成的\({w_k}\)如图\(6.4\)所示。可以看出，\(w_k\)趋近于真正的根\(w^∗ = 1\)。这个简单的例子能从直观上说明\(RM\)算法为什么收敛。
图\(6.4\)：说明RM算法收敛性的示例。
这个简单的例子可以说明RM算法收敛的原因。
- 当\(w_k>w^*\)时，有\(g(w_k)>0\)，此时 \(w_{k+1}=w_k-a_kg(w_k)<w_k\)。如果\(|a_k g(w_k)|\)足够小，我们有\(w^*<w_{k+1}<w_k\)。因此\(w_{k+1}\)比\(w_k\)更接近于\(w^*\)。
- 当\(w_k<w^*\)时，我们有\(g(w_k)<0\)，此时\(w_{k+1}=w_k-a_kg(w_k)>w_k\)。如果\(|a_k g(w_k)|\)足够小，我们有\(w^*>w_{k+1}>w_k\)。因此\(w_{k+1}\)比\(w_k\)更接近于\(w^*\)。
无论哪种情况，\(w_{k+1}\)都比\(w_k\)更接近\(w^∗\)。所以从直观上可以知道\(w_k\)收敛于\(w^∗\)。
上面的例子很简单，并没有考虑到观测误差。若存在观测误差时，下面的定理给出了严格的收敛性结果。
Info
定理6.1. (罗宾斯-门罗定理)。在\((6.5)\)中的罗宾斯-门罗算法中，如果满足
-
\(0<c_{1}\leq\nabla_{w}g(w)\leq c_{2}\;for\;all\;w\);
-
\(\sum_{k=1}^{\infty}a_{k}=\infty \; and\; \sum_{k=1}^{\infty}a_{k}^{2}<\infty\);
-
\(\mathbb{E}[\eta_{k}|\mathcal{H}_{k}]=0\; and\; \mathbb{E}[\eta_{k}^{2}|\mathcal{H}_{k}]<\infty\);
其中\(\mathcal{H}_{k}=\{w_{k},w_{k-1},\ldots\}\),那么\(w_k\)几乎必然收敛到\(g(w^∗) = 0\)的根\(w^∗\)。
Note
几乎必然收敛。相关概念在附录\(B\)中有详细介绍。
我们将该定理的证明放到第\(6.3.3\)节。首先讨论定理\(6.1\)中的三个条件。
-
在第一个条件中，\(0<c_{1}\leq\nabla_{w}g(w)\)要求\(g(w)\)是一个单调递增函数。这个条件确保了\(g(w) = 0\)的根存在且唯一。如果\(g(w)\)是单调递减函数，我们只需将\(-g(w)\)视为单调递增的新函数即可。
如果 \(g(w)\) 是一个目标函数 \(J(w)\) 的梯度，即 \(g(w) = \nabla_w J(w)\)，那么 \(c_1 \le \nabla_w J(w)\) 等价于 \(c_1 \le \nabla_w J(w)\)，这意味着 \(J(w)\) 是凸的 (convex)，这是在优化问题中经常采用的假设，也说明该条件是比较宽松的。 此外，如果 \(g(w)\) 是单调递减的，我们可以简单地将 \(-g(w)\) 视为一个新的单调递增的函数。
不等式 \(\nabla_w g(w) \le c_2\) 表明 \(g(w)\) 的梯度不能超过无穷大. \(g(w) = \tanh(w - 1)\) 是满足这个条件的，但是目前面临到的另一个例子 \(g(w) = w^3 - 5\) 则不满足这个条件。因此，当我们用 RM 算法求解 \(g(w) = w^3 - 5 = 0\) 时，需要选择合适的初始值，否则可能发散。
-
关于\(\{a_k\}\)的第二个条件很有趣。我们在强化学习算法中经常看到类似的条件。其中 \(\sum_{k=1}^{\infty} a_k^2 < \infty\) 意味着 \(\lim_{n \to \infty} \sum_{k=1}^{n} a_k^2\) 是有上界的，这要求 \(a_k\) 随着 \(k \to \infty\) 收敛到 0。条件 \(\sum_{k=1}^{\infty} a_k = \infty\) 意味着 \(\lim_{n \to \infty} \sum_{k=1}^{n} a_k\) 是无限大的，这要求 \(a_k\) 不能太快地收敛到 0。 这些条件很有意思，稍后将做详细分析。
-
第三个条件是关于噪声 \(m_k\) 的，它并不要求观测误差 \(m_k\) 是高斯分布的。该条件成立的一个重要特例是 \(\{m_k\}\) 是一个独立同分布的随机序列，因此满足 \(E[m_k] = 0\) 且 \(E[m_k^2] < \infty\)。此时，因为 \(m_k\) 与 \(H_k\) 是独立的，所以 \(E[m_k | H_k] = E[m_k] = 0\) 且 \(E[m_k^2 | H_k] = E[m_k^2]\)。
接下来，我们将更仔细地研究关于系数\({ak}\)的第二个条件。
-
为什么第二个条件对RM算法的收敛是必要的？
当然，我们稍后给出的定理的严格证明可以很自然地回答这个问题。不过，我们先简单介绍一下其基本思想。
第一，条件 \(\sum_{k=1}^{\infty} a_k^2 < \infty\) 要求随着 \(k \to \infty\) 时 \(a_k \to 0\)。这个条件为什么重要呢？由于 [ w_{k+1} - w_k = -a_k \tilde{g}(w_k, \eta_k) ] 如果 \(a_k \to 0\)，那么 \(a_k \tilde{g}(w_k, \eta_k) \to 0\) (假设 \(\tilde{g}(w_k, \eta_k)\) 总是有界的)。因此，\(w_{k+1} - w_k \to 0\)，即当 \(k \to \infty\) 时 \(w_{k+1}\) 与 \(w_k\) 非常接近。反过来说，如果 \(a_k\) 不会收敛到 0，那么 \(w_{k+1}\) 可能不会收敛，因此 \(w_k\) 会一直波动而无法收敛。
第二，\(\sum_{k=1}^{\infty} a_k = \infty\)表示\(a_k\)不能太快收敛到0。为什么这个条件重要呢？总结方程两边\(w_2-w_1=a_1\tilde{g}(w_1,\eta_1)\),\(w_3-w_2=-a_2\tilde{g}(w_2,\eta_2)\)，\(w_4-w_3=-a_3\tilde{g}(w_3,\eta_3)\)得到
\[w_1 - w_\infty = \sum_{k=1}^{\infty} a_k g(w_k, \eta_k)。\]如果\(\sum_{k=1}^{\infty} a_k < \infty\)，则 \(|\sum_{k=1}^{\infty} a_k g(w_k, \eta_k)|\) 也是有界的。令\(b\)表示这个上界，因此可以得到
\[|w_1 - w_\infty| = \left|\sum_{k=1}^{\infty} a_k g(w_k, \eta_k)\right| \leq b.\tag{6.6}\]此时，如果我们选择的初始值 \(w_1\) 离方程的解 \(w^*\) 很远，以至于 \(|w_1 - w^*| > b\)，那么 (6.6) 告诉我们 \(w_\infty = w^*\) 是不可能的，所以在这种情况下 RM 算法无法找到 \(w^*\)。
-
什么样的序列满足\(\sum_{k=1}^{\infty} a_k = \infty\)和\(\sum_{k=1}^{\infty} a_k^2 < \infty\)呢？
一个典型的序列是
\[a_k = \frac{1}{k}.\]一方面，有
\[\lim_{n \to \infty} \left( \sum_{k=1}^{n} \frac{1}{k} - \ln n \right) = \kappa,\]其中\(\kappa \approx 0.577\)被称为欧拉-马歇罗尼常数(或欧拉常数)[28]。由于当\(\lim_{n \to \infty}\)时\(\ln n = \infty\)，因此
\[\sum_{k=1}^{\infty} \frac{1}{k} = \infty.\]实际上，\(H_n = \sum_{k=1}^{n} \frac{1}{k}\)被称为微积分中的调和数 (harmonic number)。另一方面，有
\[\sum_{k=1}^{\infty} \frac{1}{k^2} = \frac{\pi^2}{6} < \infty.\]找到\(\sum_{k=1}^{\infty} \frac{1}{k^2}\)的值被称为巴塞尔问题30。
综上所述，序列 \(\{a_k = 1/k\}\) 满足定理\(6.1\)中的第二个条件。值得注意的是，稍微修改一下(例如\(a_k = 1/(k+1)\)或\(a_k = c_k/k\))其中\(c_k\)是有界的，仍然可以满足条件2。
值得一提的是，\(a_k\) 在实际中常常会被选择为一个小的正数。此时，条件 \(\sum_{k=1}^{\infty} a_k = \infty\) 仍然成立，但是条件 \(\sum_{k=1}^{\infty} a_k^2 < \infty\) 不再成立。这样选择 \(a_k\) 的原因是它能够很好地利用后面(\(k\) 比较大时)得到的样本。否则，如果 \(a_k\) 逐渐收敛到 0，那么当 \(k\) 较大时，得到的数据对 \(w_k\) 的影响已经微乎其微了。实际上，当 \(a_k\) 恒等于一个正数时，算法仍然可以在某种意义上收敛，详情参见 [24, 第 1.5 节]]。实际上，我们之所以希望 \(k\) 较大时样本数据仍有效，基本原因是这可以应对时变系统（例如函数 \(g(w)\) 随时间缓变化）。
Note
如果\(a_k=1/k\)，当\(k\)比较大时，后面的数据所起到的作用就非常小了，所以希望未来的数据仍然可以有用，就让\(a_k\)趋近于一个非常小的数而非\(0\)。
6.2.2 在期望值估计问题中的应用¶
第 6.1 节已经对增量最小期望值计算法有过讨论。我们先简要回忆一下。式（6.4）给出了一个迭代算法：\(w_{k+1} = w_k + \alpha_k (x_k - w_k)\)。当 \(\alpha_k = 1/k\) 时，可以得到 \(w_{k+1}\) 的解析表达式为 \(w_{k+1} = 1/k \sum_{i=1}^{k} x_i\)。然而，当 \(\alpha_k\) 取一般值时，我们将无法得到 \(w_{k+1}\) 的解析表达式。此时这个算法仍然能解决期望值计算问题吗？
下面展示算法（6.4）实际上是一个特殊的 RM 算法，因此其收敛性也可以得到保证。具体来说，定义如下函数：
此时，求解期望值\(\mathbb{E}[X]\)的问题可以转换为求解方程\(g(w) = 0\)的问题。对任意给定的\(w\)值，我们能获得的带噪声的观测值为\(\tilde{g}(w, \eta) = w - x,\)其中\(x\)是\(X\)的一个随机样本。注意\(\tilde{g}\)可以写为
其中\(\eta \doteq \mathbb{E}[X] - x\)。这样我们就把期望值估计的问题转换为了一个可以用RM算法求解的问题。相应RM算法求解该问题的更新公式为
上式正是\((6.4)\)中的算法。根据定理\(6.1\)，如果 \(\sum_{k=1}^{\infty} \alpha_k = \infty\) 且 \(\sum_{k=1}^{\infty} \alpha_k^2 < \infty\)，并且 \(\{x_k\}\) 是独立同分布的，那么 \(w_k\) 将几乎肯定收敛到 \(\mathbb{E}[X]\)。值得一提的是，该收敛性不依赖于 \(X\) 的分布假设。