指数均值方差爆炸：最小二乘法的救援

# 指数均值方差的爆炸性：最小二乘法的救援 – 机器学习研究博客 来源：https://francisbach.com/spectral_log_density_estimation/ 机器学习中一项常见任务是估计或优化具有（可能连续）多项的“log-sum-exp”（https://fr.wikipedia.org/wiki/LogSumExp）函数，例如： $$ \\log \\Big\( \\int\_{\\mathcal\{X\}\} e^\{v\(x\)\} dq\(x\) \\Big\)，$$ 其中 \\\(v: \\mathcal\{X\} \\to \\mathbb\{R\}\\\) 是某个势函数，\\\q\\\) 是集合 \\\(\\mathcal\{X\}\\\)\ 上的概率分布。这在数据科学中有**众多**应用，通常通过概率模型（https://en.wikipedia.org/wiki/Exponential_family#log-partition_function）的归一化来实现，但也可作为最大值的平滑近似（https://en.wikipedia.org/wiki/Smooth_maximum），在变换器（https://en.wikipedia.org/wiki/Transformer_(deep_learning)）中通过其导数使用，或在强化学习中使用熵正则化 \[19 (https://cdn.aaai.org/AAAI/2008/AAAI08-227.pdf)\]。有时集合 \\\(\\mathcal\{X\}\\\)\ 是有限的（可能很大），可以通过显式求和完成积分，但通常精确计算不可行，因此改用从概率分布 \\\(q\\\) 中采样。 关键难点来自这类估计的方差，特别是当 \\\(v\\\) 取值较大时。在最简单的例子中，对于独立且服从均值为 \\\(\\mu\\\)、方差为 \\\(\\sigma^2\\\)\ 的正态分布的 \\\(z\_1,\\dots,z\_n \\in \\mathbb\{R\}\\\)，估计 \\\(\\mathbb\{E\}\[e^z\]\)\ 的相对平方误差为： $$\\frac\{ {\\rm var}\big\( \\frac\{1\}\{n\} \\sum\_\{i=1\}^n e^\{z\_i\} \\big\) \}\{\( \\mathbb\{E\}\[ e^\{z\} \]\)^2\} = \\frac\{1\}\{n\} \\frac\{ {\\rm var}\(e^z\) \}\{\( \\mathbb\{E\}\[ e^\{z\} \]\)^2\} = \\frac\{ e^\{\\sigma^2\}\-1\}\{n\}。$$ 当 \\\(n\\\) 增长时它趋近于零（正如大数定律（https://en.wikipedia.org/wiki/Law_of_large_numbers）所预期），但当 \\\(\\sigma\\\) 增长时呈指数级爆炸。即使取对数也无法改变这种方差爆炸，即 \\\({\\rm var}\big\( \\log \\big\( \\frac\{1\}\{n\} \\sum\_\{i=1\}^n e^\{z\_i\} \\big\)\\big\)\)\ 的渐近增长也大致类似 \\\(\\frac\{ e^\{\\sigma^2\}\-1\}\{n\}\\\)\（当 \\\(n\\\) 很大时，可通过delta方法（https://en.wikipedia.org/wiki/Delta_method)得到）。 尽管难以估计，但 log-sum-exp 函数具有许多优良性质（这也是人们喜爱它的原因）；我特别喜欢以下两点：(1) 它是最大值的平滑近似（参见例如这篇早期文章（https://francisbach.com/integration-by-parts-randomized-smoothing-score-functions/）），(2) 它是适合最大似然估计（https://en.wikipedia.org/wiki/Maximum_likelihood_estimation）的概率模型归一化方式，特别是在分层概率模型中，（条件）独立性假设使得相关损失函数可分离（正如概率图模型（https://en.wikipedia.org/wiki/Graphical_model) 中充分应用的那样）。 本文试图回答的主要问题是： *我们能否在保留优化 log-sum-exp 函数优势的同时，减少其计算/统计劣势的影响？* 在谱系的另一端是**最小二乘回归**（https://en.wikipedia.org/wiki/Least_squares），其特性本质上与上述完全相反： - 在*正面*方面，我们在多种形式下获得计算和统计的简便性，例如，它通过线性代数为线性模型提供闭式估计，它基于计算具有固定受控方差的矩，并且在各种设置（加速、随机梯度下降等）中提供精确分析。参见例如这篇关于加速（https://francisbach.com/jacobi-polynomials/) 和这篇关于平均（https://francisbach.com/the-sum-of-a-geometric-series-is-all-you-need/) 的文章。 - 在*负面*方面，将最小二乘回归用于所有预测问题，特别是具有离散输出的问题，会产生一些伪影。传统的例子是具有高斯类条件数据（具有相同协方差矩阵）的分类，其中对独热编码输出使用最小二乘法会有问题，例如“掩蔽”效应（参见 \[13, Section 2\.4\] 和下面的例子），或者与使用多项逻辑回归（https://en.wikipedia.org/wiki/Multinomial_logistic_regression)（又名 softmax 回归）相比具有较高的近似误差，因为后者的对数条件概率是仿射的（https://en.wikipedia.org/wiki/Linear_discriminant_analysis)。 我们能否调和它们？换句话说，最小二乘法真的是我*唯一*需要的吗？（我的同事们有时会嘲笑我对最小二乘法的热爱）。 请注意，还有另一种（经典的）通过最小二乘法看世界的方式：通过牛顿法（https://en.wikipedia.org/wiki/Newton%27s_method) *串行*进行，在此上下文中导致迭代加权最小二乘法（https://en.wikipedia.org/wiki/Generalized_linear_model)，但这仅用于计算，没有统计改进。我们追求的是更强的东西：我们能否为通常需要优化凸函数（如逻辑或 softmax 回归）的最大似然问题，获得基于最小二乘法的**闭式**估计器？ 有趣的是，我的新尝试可以总结为一个积分方程： $$ t \\log t\\, – t \+ 1 = \\int\_0^1 \\\!\\\! \\frac\{ \(t\-1\)^2\}\{\\rho t \+ 1\-\\rho\} \(1\-\\rho\) d\\rho，$$ 这可以通过通常的积分技巧验证。让我们看看为什么以及如何做到！ ## 相对密度估计作为测试平台 在本文中，我研究一个简单的基本问题，我们可以在此研究并比较各种估计框架，注意到它可以通过多种方式扩展（特别是通过互信息，见下文）。 我们考虑集合 \\\(\\mathcal\{X\}\\\)\ 上的两个概率分布 \\\(p\\\) 和 \\\(q\\\)\；我们的目标是估计相对密度的对数 \\\(\\log \\big\(\\frac\{dp\}\{dq\}\(x\)\\big\)\\\)\。由于变分公式 \[1 (https://arxiv.org/pdf/0809.0853)\]： $$\{\\rm KL\}\(p\\\|q\) = \\int\_\{\\mathcal\{X\}\} \\log \\big\(\\frac\{dp\}\{dq\}\(x\)\\big\) dp\(x\) = \\sup\_\{v: \\mathcal\{X\} \\to \\mathbb\{R\}\} \\int\_\{\\mathcal\{X\}\} v\(x\) dp\(x\) \+ 1\\, – \\int\_\{\\mathcal\{X\}\} e^\{v\(x\)\} dq\(x\)。 \\tag\{1\}$$ 这被证明等价于估计 Kullback-Leibler（KL）散度（https://en.wikipedia.org/wiki/Kullback%E2%80%93Leibler_divergence)。 这是 \\\(f\\\)-散度（https://en.wikipedia.org/wiki/F-divergence)（参见例如 \[2 (https://people.lids.mit.edu/yp/homepage/data/itbook-export.pdf)）一个特别重要的实例，其定义和基于 \\\(f\\\) 的 Fenchel 共轭（https://en.wikipedia.org/wiki/Convex_conjugate) \\\(f^\\ast\\\) 的变分公式如下： $$D\(p\\\|q\) = \\int\_\{\\mathcal\{X\}\} f \\big\( \\frac\{dp\}\{dq\}\(x\) \\big\) dq\(x\)= \\sup\_\{v: \\mathcal\{X\} \\to \\mathbb\{R\}\} \\int\_\{\\mathcal\{X\}\} v\(x\) dp\(x\) \\, – \\int\_\{\\mathcal\{X\}\} f^\\ast\(v\(x\)\) dq\(x\)，$$ 该表示源于对每个 \\\(t = \\frac\{dp\}\{dq\}\(x\)\\\)\ 应用 \\\(f\(t\) = \\sup\_\{ u \\in \\mathbb\{R\}\} ut\-f^\\ast\(u\)\\\)。KL 散度对应于 \\\(f\(t\) = t \\log t \, – t \+ 1\\\) 和 \\\(f^\\ast\(u\) = e^u \, – 1\\\)。 注意，对于 KL 散度的特殊情况，当对 \\\(v\\\) 上的常数进行优化时，我们得到 Donsker-Varadhan 表示 \[3 (https://onlinelibrary.wiley.com/doi/epdf/10.1002/cpa.3160290405)\]： $$\{\\rm KL\}\(p\\\|q\) = \\sup\_\{v: \\mathcal\{X\} \\to \\mathbb\{R\}\} \\int\_\{\\mathcal\{X\}\} v\(x\) dp\(x\)\\, – \\log \\Big\( \\int\_\{\\mathcal\{X\}\} e^\{v\(x\)\} dq\(x\) \\Big\)。 \\tag\{2\}$$ 我们看到 log-sum-exp 函数显式出现。为了从 \\\(p\\\) 和 \\\(q\\\) 的 i.i.d. 样本估计势 \\\(v\\\)，传统的变分方法对应于用经验平均替代积分。对于 \\\(q\\\)，当只有样本可用时，这会导致一个可能不稳定的经验平均。本文的目标是探索另一种方式（我介绍这一新框架背后的主要原理；更多细节参见 \[7 (https://arxiv.org/abs/2605.10668)）。 通过 \\\(f\\\)-散度的框架对于新方法确实至关重要，因为其他散度将在定义中发挥工具性作用。在继续之前，我们陈述另一个（等价的）具有*两个*势 \\\(v\\\) 和 \\\(w\\\) 的变分公式，我们稍后会用到：\\\（D\(p\\\|q\)\\）等于： $$\\sup\_\{v,w: \\mathcal\{X\} \\to \\mathbb\{R\}\} \\int\_\{\\mathcal\{X\}\} v\(x\) dp\(x\) \+ \\int\_\{\\mathcal\{X\}\} w\(x\) dq\(x\) \\mbox{ such that } \\forall x \\in \\mathcal\{X\}, w\(x\) \\leqslant \-f^\\ast\(v\(x\)\)。 \\quad \\tag\{3\}$$ 在最优点，我们得到 \\\(w\(x\) = \-f^\\ast\(v\(x\)\)\\\)，并在 KL 情况下恢复方程 (1)。约束是凸的，但对于 \\\(f\(t\) = t \log t – t \+ 1\\\)，它远非传统凸优化方法通常允许的。该公式出现在 \[25, Theorem 4\.4\] 中，并具有保持问题对称性的良好性质（即，如果交换 \\\(p\\\) 和 \\\(q\\\)，这相当于将 \\\(f\\\) 替换为 \\\(t \\mapsto t f(1/t\)\\\)，并且这对应于交换 \\\(v\\\) 和 \\\(w\\\)）。在下文中，我们将获得满足约束 \\\(\\forall x \\in \\mathcal\{X\}, \\ w\(x\) \\leqslant \-f^\\ast\(v\(x\)\)\\\) 的函数 \\\(v\\\) 和 \\\(w\\\) 的候选者，通常不取等号。 ## 加权卡方散度 另一个相关的 \\\(f\\\)-散度函数是，对于参数 \\\(\\rho \\in \[0,1\]\\\)， $$ f\(t\) = \\frac\{1\}\{2\} \\frac\{ \(t\-1\)^2\}\{ \\rho t \+ 1\-\\rho\}。$$ 它导致加权卡方（https://en.wikipedia.org/wiki/Chi-squared_test) 散度： $$ D\(p\\\|q\) = \\frac\{1\}\{2\} \\int\_\{\\mathcal\{X\}\} \\frac\{ \\big\(\\frac\{dp\}\{dq\}\(x\)\-1 \\big\)^2\}\{ \\rho \\frac\{dp\}\{dq\}\(x\) \+ 1\-\\rho\} dq\(x\)。$$ 它已在应用数学的各个领域 \[4 (https://doi.org/10.1524/strm.2001.19.3.239),5 (https://link.springer.com/book/10.1007/978-1-4612-4946-7)\] 中使用，并且对于特殊情况有多个名称，例如当 \\\(\\rho =0\\\) 时为 Pearson 卡方散度，或当 \\\(\\rho=1\\\) 时为 Neyman 卡方（或反向 Pearson），或当 \\\(\\rho=1/2\\\) 时为 Le Cam 散度。 上述函数 \\\(f\\\) 的形式为“二次除以仿射”，它通过“二次加仿射”函数具有变分表示： $$ \\frac\{1\}\{2\} \\frac\{ \(t\-1\)^2\}\{ \\rho t \+ 1\-\\rho\} = \\sup\_\{u \\in \\mathbb\{R\}\} \\ \(t\-1\) u \\, – \\frac\{1\}\{2\} \( \\rho t \+ 1 – \rho\) u^2，$$ 其中最优 \\\(u = \\frac\{t\-1\}\{\\rho t \+ 1 \\, – \rho\}\\\)（顺便说一句，这*不是* Fenchel 表示）。 因此，对每个 \\\(x \\in \\mathcal\{X\}\\\)\ 应用此式，令 \\\(t = \\frac\{dp\}\{dq\}\(x\)\\\)，对于这个函数 \\\(f\\\)，我们有： $$ D\(p\\\|q\) = \\\!\\\! \\sup\_\{u\(\\rho,\\cdot\):\\mathcal\{X\} \\to \\mathbb\{R\}\} \\int\_\{\\mathcal\{X\}\} \\Big\\\{ \\big\(\\frac\{dp\}\{dq\}\(x\) \-1\\big\) u\(\\rho,x\) \\, – \\frac\{1\}\{2\} \\big\( \\rho \\frac\{dp\}\{dq\}\(x\) \+ 1 – \\rho\\big\) u\(\\rho,x\)^2 \\Big\\\} dq\(x\)。$$ 这恰好是一个二次变分问题，因为函数 \\\(u\(\\rho,\\cdot\): \\mathcal\{X\} \\to \\mathbb\{R\}\\\)\ 仅以二次形式出现。最优变分函数为 \\\(\\displaystyle u\(\\rho,x\) = \\frac\{ \\frac\{dp\}\{dq\}\(x\) \\, – 1\}\{ \\rho \\frac\{dp\}\{dq\}\(x\) \+ 1\-\\rho\}\)。\\\) 请注意，我们刚刚定义的二次代价函数显式地是一个最小二乘预测问题，即在给定 \\\(x\\\) 的情况下预测 \\\(y\\\)，其中以概率 \\\(\\rho\\\)，\\\ (y\\\) 取值 \\\(1/\\rho\\\) 且 \\\(x\\\) 从 \\\(p\\\) 中采样，以概率 \\\(1\-\\rho\\\)，\\\ (y\\\) 取值 \\\(\-1/\(1\-\\rho\)\\\) 且 \\\(x\\\) 从 \\\(q\\\) 中采样（因此这让人想起噪声对比估计 \[26\]）。因此我们看到在 \\\(\\rho=0\\\) 或 \\\(\\rho=1\\\) 处存在潜在的不稳定性（我们稍后会看到并非如此）。 \\\（D\(p\\\|q\)\\）的这种变分公式可以通过将对 \\\(p\\\) 或 \\\(q\\\) 的期望组合在一起重写为： $$ \\sup\_\{u\(\\rho,\\cdot\):\\mathcal\{X\} \\to \\mathbb\{R\}\} \\int\_\{\\mathcal\{X\}\} \\big\( u\(\\rho,x\) \\, – \\frac\{\\rho\}\{2\} u\(\\rho,x\)^2 \\big\) dp\(x\) \+ \\int\_\{\\mathcal\{X\}\} \\big\( \-u\(\\rho,x\)\\, – \\frac\{1\-\\rho\}\{2\} u\(\\rho,x\)^2 \\big\) dq\(x\)。$$ 那么我们恰好得到一个双势公式，其中： $$v\(x\) = u\(\\rho,x\) \\, – \\frac\{\\rho\}\{2\} u\(\\rho,x\)^2 \\ \\ \\mbox{ and } \\ \\ w\(x\) = \-\\, u\(\\rho,x\) \\, – \\frac\{1\-\\rho\}\{2\} u\(\\rho,x\)^2$$ 作为两个势（严格的读者可以验证 \\\(\\forall x \\in \\mathcal\{X\}, w\(x\) \\leqslant \-f^\\ast\(v\(x\)\)\\\)，这并不明显，因此它们是方程 (3) 的优化器）。 总结一下，对于函数 \\\(f: t \\mapsto \\frac\{1\}\{2\} \\frac\{ \(t\-1\)^2\}\{ \\rho t \+ 1\-\\rho\}\\\)，我们恰好得到了我们想要的东西，即通过最小二乘法的变分公式。如何将其扩展到更一般的 \\\(f\\\)-散度？ ## 通过积分扩展 如果我们可以写成 \\\(\\displaystyle f\(t\) = \\frac\{1\}\{2\} \\int\_0^1 \\frac\{ \(t\-1\)^2\}\{\\rho t \+ 1\-\\rho\} d\\nu\(\\rho\)\\\)，其中 \\\(\\nu\\\) 是区间 \\\(\[0,1\]\\\)\ 上的某个非负测度，那么我们可以直接利用上述推导，将 \\\(D\(p\\\|q\) \\\) 表示为： $$\\\!\\\!\\sup\_\{u:\[0,1\] \\times \\mathcal\{X\} \\to \\mathbb\{R\}\} \\int\_\{0\}^1\\\!\\\!\\\!\\\! \\int\_\{\\mathcal\{X\}\} \\\!\\\! \\Big\[\\\! \\, \(\\frac\{dp\}\{dq\}\(x\) \-1\) u\(\\rho\)\(x\) \\, – \\frac\{1\}\{2\} \( \\rho \\frac\{dp\}\{dq\}\(x\) \+ 1 \\, – \\rho\) u\(\\rho\)\(x\)^2 dq\(x\) \\\!\\Big\] d\\nu\(\\rho\)， \\tag\{4\}$$ 现在 \\\(u\\\) 是从 \\\(\[0,1\] \\times \\mathcal\{X\}\\\)\ 到 \\\(\\mathbb\{R\}\\\)\ 的函数。 那么我们恰好得到一个双势公式，如方程 (3) 所述，其中： $$v\(x\) = \\int\_0^1 \\big\[ u\(\\rho,x\) \\, – \\frac\{\\rho\}\{2\} u\(\\rho,x\)^2 \\big\] d\\nu\(\\rho\) \\tag\{5\}$$ 和 $$w\(x\) = \\int\_0^1 \\big\[ \-u\(\\rho,x\) \\, – \\frac\{1\-\\rho\}\{2\} u\(\\rho,x\)^2 \\big\] d\\nu\(\\rho\)， \\tag\{6\}$$ （它们满足约束 \\\(w\(x\) \+ f^\\ast\(v\(x\)\) \\leqslant 0\\\)，这同样不简单），其中 \\\(u\(\\rho,\\cdot\)\\)\ 是方程 (4) 的最大化器。这对应于并行执行**连续的最小二乘问题**。 这些推导对所有具有积分表示的 \\\(f\\\)-散度有效，特别是 KL 散度，因为我们有： $$t \log t\\, – t \+ 1 = \\int\_0^1 \\\!\\\! \\frac\{ \(t\-1\)^2\}\{\\rho t \+ 1\-\\rho\} \(1\-\\rho\) d\\rho，$$ 即 \\\(d\\nu\(\\rho\) = 2 \(1\-\\rho\) d\\rho\\\)。这个积分表示并非凭空而来；事实上，它来自于我们在早期文章中探讨的算子凸和算子单调（https://francisbach.com/matrix-monotony-and-convexity/) 函数理论。它包括 KL 散度，显然还有所有加权卡方散度（\\\（\\nu\\\) 是狄拉克测度），以及所有 \\\(\\alpha\\\)-散度 \[14 (https...