Hutchinson 迹估计器的现代分析.

摘要

本文建立了 Hutchinson 迹估计器精度分析的最新技术水平。 利用以前在此背景下未曾使用过的工具，特别是超收缩不等式和次伽马分布的集中性质，我们提供了一种优雅且模块化的分析，以及数值上更优越的界限。 除了这些改进之外，这项工作旨在更好地普及计算机科学界中提到的技术。

I 引言

I-A 背景

估计矩阵迹是一个具有基本意义的问题 [1, 2, 3] 并出现在许多问题中，例如逼近矩阵的谱特性 [4, 5, 6]， 求解偏微分方程 [7, 8, 9, 10, 11]，机器学习中的误差评估 [12]，以及组合计数 [13]，仅举几例。 对于特别关注数据科学或优化的读者来说，这至关重要，因为 Hessian 迹存储了关于曲率的有价值的信息。 为了提供一个有意义的例子，考虑在 MNIST 数据集 [14] 上拟合一个线性分类器，该数据集是机器学习中广泛用作基准和玩具问题的。 由于图像是 $28\times 28$ 分辨率，并且被分组到 $10$ 类中，参数数量为 $m=28\cdot 28\cdot 10=7840$，Hessian 矩阵的大小为 $m^2\approx 6\cdot {10}^6$。 与迹成正比的对角线平均值等于平均特征值（有价值的信息），而单个行或列几乎为零（直到随机噪声）。 这在 脚注 2 中进行了说明。

图 1： 在 MNIST 数据集上训练的线性分类器的 Hessian 矩阵。 与行和列相比，对角线存储了有价值的信息：行和列中的条目几乎为零，而对角线的平均值为 $\approx 0.25$。 Tensorflow 代码以 Python Notebook 的形式在 GitHub 上提供 。

从上面的例子可以看出，即使是玩具例子也会导致巨大的矩阵；对于更大的问题，存储和检查这样的矩阵是不可能的。 幸运的是，可以估计迹（等于对角线之和，或等效地等于特征值之和）而无需知道整个矩阵。 要求是能够有效地计算 矩阵向量 乘积；在 Hessian 的情况下，此类乘积可以通过自动微分 [15] 在所有流行的机器学习框架中计算，例如 Tensorflow [16]、PyTorch [17] 和 JAX [18]。 在假设一个 $m\times m$ 矩阵 $A$ 是 对称半正定（这适用于所有 Hessian）的情况下，Hutchinson [1] 提出的流行估计器为

${\mathrm{tr}}_H(A)\triangleq z^{T}Az,z{\sim }^{iid}{\{-1,1\}}^m,$ (1)

其中向量 $z$ 的分量是 Rademacher 随机变量，等于 $+1$ 或 $-1$，概率相等。 上述估计器仅使用一个采样的 $v$，因此是噪声的（方差很大），因此通常通过平均 $n$ 个独立试验来增强（对于适当大的 $n$）。 正式地：

${\mathrm{tr}}_{H^{(n)}}(A)\triangleq {\displaystyle \frac{1}{n}}{\displaystyle \sum _{i=1}^n}{z}_i^{T}A{z}_i,z_i{\sim }^{iid}{\{-1,1\}}^m.$ (2)

估计器是无偏的，即平均而言是正确的：

$𝔼{\mathrm{tr}}_H(A)=𝔼{\mathrm{tr}}_{H^{(n)}}(A)=\mathrm{tr}(A),$ (3)

这很容易证明。 本文关注的是一个更具挑战性的问题，即理解它们的 集中性质。 在这里，我们问以下问题：

什么样本量 $n$ 保证 在置信水平为 $1-\delta$ 时达到所需的相对误差 $ϵ$？

形式上，对于误差 ${\mathrm{err}}_{H^{(n)}}(A)\triangleq \frac{{\mathrm{tr}}_{H^{(n)}}(A)}{\mathrm{tr}(A)}-1$ 和固定的 $ϵ,\delta$，我们希望找到一个可能很小的 $n=n(ϵ,\delta )$ 使得

$|{\mathrm{err}}_{H^{(n)}}(A)|⩽ϵ\text{with prob. at least}1-\delta .$ (4)

I-B 相关工作

该估计器已有 30 多年的历史 [1]，尽管存在替代方案（例如基于 Lanczos 积分的方法 [19]），但它在方差方面已被证明是最好的 [20]，并且可以说是最简单的，受到统计/学习软件开发者的青睐 [21]。 令人惊讶的是，它在 Avron 和 Toledo 的工作 [3] 之前，一直被使用而没有进行严格的准确性评估，他们的工作确定了有限样本大小的保证。 他们的结果后来被 Roosta-Khorasani 和 Ascher [20] 改善，他们基本上摆脱了涉及粗略并集界定的损失证明步骤。 他们的方法基于 Chernoff 式的直接计算，并提供了边界 $n=O({ϵ}^{-2}\log (1/\delta ))$。

II 贡献

II-A 摘要

在这项工作中，我们对 Hutchinson 的估计器进行了 模块化、现代化和更准确 的分析。 与之前工作的改进可以概括如下：

• 技术的模块化和新颖性。 在第一步中，我们将估计器精度明确地与 Rademacher 二次混沌 的色散联系起来，该混沌具有酉矩阵核。 然后，我们通过 超收缩不等式 获得了对这种混沌的界限，这是一个布尔傅立叶分析的重要结果 [22]，最初由 Aline Bonami 提出 [23]。 作为最后一步，我们将此界限表示为 次伽马性质（见 Boucheron 的著作 [24]），这使得得出最终的集中结果变得方便（且准确）。 此外，我们给出了多样本 和 单样本估计器的结果；这很有趣，因为基本构建块可能比通过平均（例如，参见中值算法 [25]）获得的提升有所不同。

  因此，与之前工作中的特设计算相反，我们能够使用来自高维概率和傅里叶分析领域的既定工具。 我们透明的方法为进一步改进开辟了道路。

• 优越的精度。 通过专用工具，我们获得了数值上好一个数量级的界限。 它们被表述为基本公式，便于在迹估计的实际应用中使用。

II-B 预备知识

下面我们解释在下一节中提出结果时使用的符号和定义。

r.v. 的 $d$-范数 $X$ 为 ${\Vert X\Vert }_d\triangleq {(𝔼{|X|}^d)}^{1/d}$；当 $d⩾1$ [25] 时，它是一个有效的范数（例如，它服从三角不等式）。 研究集中性的一个有用工具是 矩生成函数，定义为 ${\mathrm{MGF}}_X(t)=𝔼\exp (tX)$（实参数 $t$ 的函数）。 在现代高维概率中，人们根据 $\mathrm{MGF}$ 的行为对分布进行分类；特别是当 $\log {\mathrm{MGF}}_{|X|}(t)⩽\frac{v{t}^2}{2(1-ct)}$ 时，我们称 $X$ 为 次伽马，方差因子为 $v$，尺度为 $c$，并用 $X\in \mathrm{\Gamma }(v,c)$ [24] 表示； 相同的公式适用于中心伽马分布，因此得名。

II-C 结果

以下我们假设 $A$ 是非零对称半正定矩阵。 然后，必然地，$\mathrm{tr}(A)>0$。

II-C1 单样本估计器

下面我们给出估计器在 (1) 中的精度的以下界限：

Theorem 1 (哈奇森估计器 1 样本界限)。

对于任何整数 $d⩾2$，哈奇森估计器 (1) 的相对误差服从以下界限

${\Vert {\mathrm{err}}_H(A)\Vert }_d⩽d-1.$ (5)

特别地，这意味着次伽玛行为

${\mathrm{err}}_H(A)\in \mathrm{\Gamma }(1,{\displaystyle \frac{8}{3}}),$ (6)

这反过来又给出了以下概率尾界

$\mathrm{Pr}[|{\mathrm{err}}_H(A)|⩾ϵ]⩽{\mathrm{e}}^{-\frac{{ϵ}^2}{2\left(1-\frac{8}{3}ϵ\right)}}.$ (7)

证明过程如下：a) 我们用 Rademacher 混沌的形式表示误差 b) 我们使用超收缩不等式来限制它的矩；这建立了部分 (5) c) 我们将矩界限代入矩生成函数的泰勒展开式，然后估计级数，得到次伽玛条件 (6)。 根据次伽玛分布的尾部性质，我们得出尾部界限 (7)。

II-C2 多样本估计器

接下来，我们转向 (2) 中的多样本估计器。 我们建立以下

Theorem 2 (哈奇森估计器 $n$样本界限)。

$n$样本哈奇森估计器的相对误差具有以下次伽玛行为

${\mathrm{err}}_H(A)\in \mathrm{\Gamma }({\displaystyle \frac{1}{n}},{\displaystyle \frac{8}{3}}).$ (8)

特别地，这意味着误差尾部概率的以下界限，对于任何 均有效：

$\mathrm{Pr}[|{\mathrm{err}}_H^{(n)}(A)|⩾ϵ]⩽{\mathrm{e}}^{-\frac{n{ϵ}^2}{2\left(1-\frac{8}{3}ϵ\right)}}.$ (9)

此结果通过使用关于子伽马随机变量之和的集中性的额外事实，从 定理 1 中获得：实质上，方差因子从 $1$ 减小到 $\frac{1}{n}$，用于长度为 $n$ 的总和（就像 $n$ 个 iid 项的方差一样）。

Corollary 1 (Hutchinson 估计量的样本大小).

相对误差绝对以 $ϵ$ 为界限，概率为 $1-\delta$，前提是样本大小 $n$ 大于或等于

$n(ϵ,\delta )={\displaystyle \frac{2\left(1-\frac{8}{3}ϵ\right)\log \left(\frac{1}{\delta }\right)}{{ϵ}^2}}.$ (10)

这对任何 和 都成立。

II-D 与相关工作的数值比较

下面我们将展示相对于先前工作  [3, 20] 的数值改进。 具有精确常数的界限在下面的 表 I 中进行了总结。

Author	Sample Size $n=n(ϵ,\delta )$
this work	$\frac{2\left(1-\frac{8}{3}ϵ\right)\log \left(\frac{1}{\delta }\right)}{{ϵ}^2}$
[20]	$\frac{6\log \left(\frac{2}{\delta }\right)}{{ϵ}^2}$
[3]	$\frac{6\log \left(\frac{2}{\delta }\right)+6\mathrm{r}\mathrm{a}\mathrm{n}\mathrm{k}(A)}{{ϵ}^2}$

表 I： Hutchinson 估计器精度的界限。

此外，在 图 2 中，我们详细评估了样本大小 $n$ 如何取决于误差 $ϵ$，当置信参数 $\delta$ 固定时。 此实验的设置是讨论的 MNIST 上玩具分类器的 Hessian。

图 2： Hutchinson 估计量的样本大小。 我们假设 $\delta ={10}^{-3}$，这导致 $0.999$ 的置信度。 矩阵大小为 $7840$，对应于线性 MNIST 分类模型的海森矩阵（输入数据是 $28\times 28$ 分辨率的黑白图像，分组为 $10$ 类）。

III 辅助结果

III-1 矩阵理论

我们将需要以下关于对称矩阵分解的事实 [26]。 回想一下，当 $B^T=B^{-1}$ 时，矩阵 $B$ 是正交的 (${\cdot }^T$ 表示转置)；特别地，$B$ 的列和行是单位长度的。

Lemma 1 (对称矩阵分解).

任何对称实矩阵 $A$ 可以写成 $A=B^T\mathrm{\Lambda }B$，其中 $B$ 是一个正交实矩阵，而 $\mathrm{\Lambda }$ 是对角的，由 $A$ 的（必然是实数的）特征值组成。

III-2 二次混沌

我们将需要以下关于著名的超对比不等式 [23, 22] 的特例：

Lemma 2 (($2,d)$-超对比).

令 $F$ 是 Rademacher 变量 $Z_i$ 中的一个 $2$ 次多项式（例如，对于某些权重 $a_{i,j}$，$F={\sum }_{i\ne j}{a}_{i,j}{Z}_i{Z}_j$）。 那么

${\Vert F\Vert }_d⩽(d-1){\Vert F\Vert }_2.$ (11)

III-3 次伽马随机变量

次伽马分布的尾部概率可以通过 Crammer-Chernoff 方法从 MGF 界限估计，如下所示：

Lemma 3 (次伽马尾部).

令 $X\in \mathrm{\Gamma }(v,c)$，则

$\forall ϵ>0:\mathrm{Pr}[|X|⩾ϵ]⩽{\mathrm{e}}^{-\frac{{ϵ}^2}{v+cϵ}}.$ (12)

此外，对于次伽马分布，以下合成性质对于独立项的总和成立：

Lemma 4 (次伽马和).

令随机变量 $X_i\in \mathrm{\Gamma }(v_i,c_i)$ 独立，则

${\displaystyle \sum _i}{X}_i\in \mathrm{\Gamma }(v^{\prime },c^{\prime }),v^{\prime }\triangleq {\displaystyle \sum _i}{v}_i,c^{\prime }\triangleq \underset{i}{\max }{c}_i.$ (13)

这些引理的证明可以在 [24] 中找到。

IV 证明

IV-A 定理 1 的证明

定义以下随机变量

$Y=z^{T}Az-𝔼[z^{T}Az]=z^{T}Az-\mathrm{tr}(A).$ (14)

然后我们的任务是限制 $Y$ 的集中度。

IV-A1 归约到非对角二次混沌

由于 $A$ 是对称的，根据 引理 1，我们有 $A=B^T\mathrm{\Lambda }B$，其中 $B$ 是正交的，$\mathrm{\Lambda }$ 是由 $A$ 的特征值组成的对角矩阵。 此外，由于矩阵 $A$ 是半正定的，因此它的特征值是非负的。 因此，我们可以写成

$z^{T}Az={(Bz)}^T\mathrm{\Lambda }(Bz)={\Vert {\mathrm{\Lambda }}^{1/2}Bz\Vert }_2^2.$ (15)

这可以分解为

$z^{T}Az={\displaystyle \sum _i}{Y}_i,Y_i$ $\triangleq {({\mathrm{\Lambda }}^{1/2}Bz)}_i^2.$ (16)

为了进一步简化，用 ${\lambda }_1,\mathrm{\dots },{\lambda }_m$ 表示 $A$ 的所有特征值；这些正是 $\mathrm{\Lambda }$ 的对角线元素。 我们注意到 $Y_i={\lambda }_i{({\sum }_j{B}_{i,j}{z}_j)}^2$， 以及 $𝔼Y_i={\lambda }_i{\sum }_j{B}_{i,j}^2$（因为 $𝔼z_j^2=1$）。 因此我们得到

$Y_i-𝔼Y_i={\lambda }_i{\displaystyle \sum _{j\ne j^{\prime }}}{B}_{i,j}{B}_{i,j^{\prime }}{z}_j{z}_{j^{\prime }}.$ (17)

因此，每个 $Y_i$ 都是 $z_j$ 中具有零对角线的二次型。

IV-A2 限制二次混沌

对于每个固定的 $i$，我们将 引理 2 应用于 $Y_i-𝔼Y_i$，它是一个关于 Rademacher 随机变量 $z_i$ 的 $2$ 次多项式（注意 非对角线属性保证次数正好是 2）。

由于由元组 $(j,j^{\prime })$ 索引的随机变量 $z_j{z}_{j^{\prime }}$ 对于 $j\ne j^{\prime }$ 是不相关的， 我们有 ${\Vert Y_i-𝔼Y_i\Vert }_2^2={\lambda }_i^2{\sum }_{j\ne j^{\prime }}{B}_{i,j}^2{B}_{i,j^{\prime }}^2$； 但是 $B$ 是正交的，所以 ${\sum }_{i,j}{B}_{i,j}^2=1$ 和 ${\Vert Y_i-𝔼Y_i\Vert }_2^2⩽{\lambda }_i^2$。 因此，我们得到

${\Vert Y_i-𝔼Y_i\Vert }_d⩽{\lambda }_i^2(d-1).$ (18)

根据三角不等式

${\Vert {\displaystyle \sum _i}(Y-𝔼Y_i)\Vert }_d⩽(d-1){\displaystyle \sum _i}{\lambda }_i.$ (19)

最后，根据 $Y_i$ 的定义和线性代数中的标准恒等式 $\mathrm{tr}(A)={\sum }_i{\lambda }_i$ [26]

${\Vert z^{T}Az-\mathrm{tr}(A)\Vert }_d⩽(d-1)\mathrm{tr}(A).$ (20)

将这个不等式的两边都除以 $\mathrm{tr}(A)$，我们就完成了 定理 1 第一部分的证明。

IV-A3 限制矩生成函数

令 $E=z^{T}Az/\mathrm{tr}(A)-1$ 为估计量的相对误差。 然后根据先前的结果，泰勒展开式 ${\mathrm{e}}^{tx}={\sum }_{d⩾0}{(tx)}^d/d!$ 以及 $𝔼E=0$ 的事实，我们有

$MGF(E)⩽1+{\displaystyle \sum _{d⩾2}}{\displaystyle \frac{t^d{(d-1)}^d}{d!}}.$ (21)

令 $a_d={(d-1)}^d/d!$。 那么我们有

${\displaystyle \frac{a_{d+1}}{a_d}}={\displaystyle \frac{d^{d+1}}{{(d-1)}^d}}\cdot {\displaystyle \frac{1}{d+1}}.$ (22)

当 $d⩾2$ 时，该比率递减，因此

${\displaystyle \frac{a_{d+1}}{a_d}}⩽{\displaystyle \frac{8}{3}}.$ (23)

因此，可以得出

$\mathrm{MGF}(E)⩽1+{\displaystyle \frac{t^2}{2(1-\frac{8}{3}t)}},$ (24)

其中，根据标准不等式 $\log (1+u)⩽u$，这意味着

$\log \mathrm{MGF}(E)⩽{\displaystyle \frac{t^2}{2(1-\frac{8}{3}t)}},$ (25)

因此，根据定义，$E\in \mathrm{\Gamma }(1,8/3)$。 这完成了 定理 1 第二部分的证明。

IV-A4 总结集中特性

定理 1 的第三部分直接由 引理 3 得出。

IV-B 定理 2 的证明

我们首先应用 引理 4 中关于独立次伽马分布之和的结果，这证明了定理的第一部分。 第二部分接着介绍了关于子伽马尾部的结果，如 引理 3 中所述。

IV-C 推论 1 的证明

推论的证明是通过将 定理 2 中尾部界限的右侧等同于 $\delta$ 来实现的。 关于 $n$ 求解，得到样本大小的公式。

V 结论

本工作对哈奇森提出的经典迹估计器建立了一个新的最先进界限。 核心思想是巧妙地利用了对 Rademacher 混乱的界限和子伽马分布理论。 这些结果在数值上优越，并且对于任何涉及迹估计的问题都具有直接的意义。 此外，作者希望为计算机科学界更好地普及在问题背景下新颖的讨论技术。