不确定性的线性代数

第 2.6.3 节介绍了一张表，记录每一对不确定量之间的协方差，即协方差矩阵。目标函数在已评估的每个输入、以及正在考虑的每个输入处都有一个值，贝叶斯优化器为这些值保存这样一张表，因此 50 次评估之后，这张表至少有 50 行。信念的每一次更新、每一次预测、对模型自身设置的每一次拟合，都是在这张表上的计算。本章讨论的就是这些计算。

线性代数是一门很大的学科，本章只选取本书后文用到的部分。本章把矩阵看作移动空间中各点的映射，因为这一图景能够解释几个关键性质：哪些表可以是协方差矩阵，协方差赋予一团不确定性什么样的形状，如何完全不求逆矩阵而解出高斯过程的方程，以及行列式衡量的是什么。大部分直觉可以由二维的图建立，但有几个事实要到三维才会出现，因此本章还将让读者亲手转动一个三维协方差。

这些内容的作用在下一章第 4 章中显现：这里的每个结果都会成为关于高斯分布的结论；在第 8 章中，这些结论又转化为一个算法。

3.1 向量 #

向量（vector）是一组有序的数。向量用粗体小写字母表示，按惯例写成列的形式：

符号 $\R^d$ 表示由 $d$ 个实数组成的所有这类列表的集合，$d$ 称为向量的维度（dimension）。写成 $\vx^\T = (x_1, \dots, x_d)$（读作“$\vx$ 转置”），列就变成了行。

向量在本书中有两种用法，区分二者可以避免一个常见的混淆。第一种用法中，向量是搜索空间中的一个点：一组超参数配置、一套外骨骼的设置、一个照片滤镜的参数。例如，Snoek 等人（2012）同时调节图像分类器训练过程的九个设置，因此他们尝试的每个配置都是 $\R^9$ 中的一个向量（第 9.2 节会再讨论他们的发现）。第二种用法中，向量保存目标函数在一组输入处的值，$\vf = (f(\vx_1), \dots, f(\vx_n))^\T$。它的维度是输入的个数，而不是超参数的个数，并且每评估一次就增加一维。高斯过程的协方差矩阵就位于这第二种空间中。

向量按分量相加；向量乘以一个数，即每个分量都乘以这个数。从几何上看，$\vx + \vy$ 是把 $\vy$ 的箭头接在 $\vx$ 的箭头末端，$2\vx$ 则是方向相同、长度加倍的箭头。

3.1.1 内积、长度与角度 #

两个维度相同的向量逐个分量相乘再求和，得到它们的内积（inner product），又称点积：

这一记号读作一行乘以一列，这正是第 3.2 节将介绍的矩阵乘法规则。长度和距离都由内积而来。向量的范数（norm）即其长度，$\lVert\vx\rVert = \sqrt{\vx^\T\vx}$，也就是 $d$ 维空间中的勾股定理；两点之间的距离是二者之差的范数，$\lVert\vx - \vy\rVert$。

内积还能衡量角度。对任意两个向量，

其中 $\theta$ 是两者的夹角。内积为零的两个向量是正交的（orthogonal），即互相垂直。除以两个长度就得到 $\cos\theta$：方向相同时为 1，互相垂直时为 0，方向相反时为 $-1$。第 2.6.3 节的相关系数就是这个余弦，只是计算对象换成了中心化的随机变量，而不是数的列表。

距离对贝叶斯优化很重要。核函数决定目标函数在两个输入处的值相关程度有多强（第 9 章），而最常用的核函数只通过两个输入之间的距离依赖于它们。按照径向基函数（RBF）核 $k(\vx, \vx') = \exp\!\left(-\lVert\vx - \vx'\rVert^2 / 2\ell^2\right)$，相近的输入取值相似，相距遥远的输入互不相关，“相近”的含义则由长度尺度 $\ell$ 决定。

3.1.2 高维空间中的距离 #

我们对距离的直觉来自二维和三维，到了更高的维度，这种直觉会产生误导。在单位立方体 $[0, 1]^d$ 中均匀随机地撒点。每个坐标对两点间距离的平方平均贡献 $1/6$，因此典型距离按 $\sqrt{d/6}$ 增长：二维中约为 0.58，六维中为 1.0，五十维中为 2.9。距离的离散程度却不随之增长。在一次用 200 个随机点所做的模拟中，从二维到一百维，两两距离的标准差始终保持在 0.24 左右。结果，对一个典型的点，它到最近邻的距离与到最远邻的距离之比，在二维中是 0.04，在六维中升至 0.22，在五十维中升至 0.65，在一百维中升至 0.74。

因此，在高维空间中，每个点都远离其他所有点，而且所有点之间的距离都相差无几。对于五十维问题，长度尺度适合二维问题的核函数会认为每一对点都互不相关，模型从一次评估中得不到关于其他任何评估的信息。这是贝叶斯优化在高维中困难的根源之一，种种补救办法，例如随维度增长的长度尺度、为每个输入分别设定长度尺度（第 9.2 节），也都从这里出发（第 30 章）。

3.2 作为变换的矩阵 #

矩阵（matrix）是一张矩形的数表，用粗体大写字母表示。$m \times n$ 矩阵 $\mA$ 有 $m$ 行、$n$ 列，$A_{ij}$ 是第 $i$ 行第 $j$ 列的元素。这张表是矩阵的存储格式，而矩阵的含义在于它所定义的映射。

3.2.1 矩阵与向量的乘积 #

$m \times n$ 矩阵乘以 $n$ 维向量，得到 $m$ 维向量，其第 $i$ 个分量是第 $i$ 行与该向量的内积：

同一个乘积还有第二种读法，也是应当在头脑中想象的那一种。把 $\mA$ 的各列记作 $\mathbf{a}_1, \dots, \mathbf{a}_n$，则

向量 $\mathbf{e}_1 = (1, 0, \dots, 0)^\T$ 选出第一列：$\mA\mathbf{e}_1 = \mathbf{a}_1$。可见矩阵的各列记录了各坐标方向被映射到何处；按照式（3.2），其他每个向量都按其坐标移动：$\vx$ 含有多少 $\mathbf{e}_1$，就取多少 $\mathbf{a}_1$，依此类推。

这类映射是线性的（linear）：它把和映射为和，把倍数映射为倍数，$\mA(\vx + \vy) = \mA\vx + \mA\vy$，$\mA(c\vx) = c\,\mA\vx$。在平面上，线性映射把直线映射为直线，保持原点不动，并把正方形网格变为平行四边形网格。单位圆，即所有长度为一的向量构成的集合，变为一个椭圆。

可以尝试以下操作：

• 拖动 $\mA\mathbf{e}_1$ 的末端。只有第一列改变，整个网格随之移动：每个点的移动量，是它的第一个坐标乘以这一改变量。
• 依次选择预设：旋转、剪切、拉伸。旋转使圆转动而形状不变；拉伸以不同的倍数缩放两个轴；剪切使一个轴沿另一个轴滑动，把正方形变为面积不变的斜平行四边形。
• 选择预设：奇异。两列指向同一方向，椭圆坍缩成一条线段，整个平面被压扁到一条直线上。这条直线上的每个点都来自无穷多个点，因此这个映射无法撤销。

3.2.2 复合、单位矩阵与逆矩阵 #

先作用 $\mathbf{B}$、再作用 $\mA$，仍是一个线性映射，其矩阵是乘积 $\mA\mathbf{B}$，各列是 $\mA$ 作用于 $\mathbf{B}$ 各列的结果。逐个元素来看，$(\mA\mathbf{B})_{ij} = \sum_k A_{ik} B_{kj}$，即 $\mA$ 的第 $i$ 行与 $\mathbf{B}$ 的第 $j$ 列的内积。顺序很重要：先拉伸再旋转，与先旋转再拉伸并不相同，因此一般而言 $\mA\mathbf{B} \ne \mathbf{B}\mA$。按这种方式相乘两个 $n \times n$ 矩阵需要 $n^3$ 次乘法，$i$、$j$、$k$ 的每种组合各一次。

单位矩阵（identity matrix）$\mI$ 对角线上为一、其余为零，它使每个向量保持原位。如果另一个矩阵 $\mA^{-1}$ 能撤销方阵 $\mA$ 的作用，即 $\mA^{-1}\mA = \mA\mA^{-1} = \mI$，则称这个方阵可逆（invertible）。矩阵可逆，当且仅当它不把空间压扁：没有非零向量被映射为零，因而没有两个点映射到同一个像。图 3.1 中的奇异预设正是相反的情形。乘积求逆要颠倒顺序，$(\mA\mathbf{B})^{-1} = \mathbf{B}^{-1}\mA^{-1}$，因为最后施加的映射必须最先撤销。

3.2.3 转置 #

转置（transpose）$\mA^\T$ 把矩阵沿对角线翻转，使 $(\mA^\T)_{ij} = A_{ji}$，$m \times n$ 矩阵变为 $n \times m$ 矩阵。列向量是 $n \times 1$ 矩阵，其转置是 $1 \times n$ 的行，内积 $\vx^\T\vy$ 就是一行与一列的乘积。反过来的乘积 $\vx\vy^\T$，即一列乘以一行，是一个 $n \times n$ 矩阵，称为外积（outer product）。乘积的转置同样要颠倒顺序，$(\mA\mathbf{B})^\T = \mathbf{B}^\T\mA^\T$。

第 8 章中各个量的形状都遵循这些规则。设有 $n$ 个已观测的输入和 $m$ 个待预测的输入，已观测输入之间的核矩阵 $\mK$ 是 $n \times n$ 的，已观测输入与新输入之间的矩阵 $\mK_*$ 是 $n \times m$ 的；后验均值 $\mK_*^\T\mK^{-1}\vy$ 把一个 $m \times n$ 矩阵、一个 $n \times n$ 矩阵与一个 $n$ 维向量相乘，为 $m$ 个新输入各给出一个预测。

3.3 对称矩阵与正定矩阵 #

并非每张方形数表都能作为协方差矩阵。本节找出其中的条件，结果表明它是几何性的：协方差矩阵不能让任何方向的方差为负。

3.3.1 二次型 #

矩阵等于其转置，即 $A_{ij} = A_{ji}$ 时，称为对称的（symmetric）。协方差矩阵是对称的，因为 $\Cov[x_i, x_j] = \Cov[x_j, x_i]$；核矩阵也是，因为 $k(\vx, \vx') = k(\vx', \vx)$。

对于对称矩阵 $\mA$，表达式 $\vw^\T\mA\vw$ 对每个向量 $\vw$ 给出一个数，称为二次型（quadratic form）。在二维中，取 $\mA = \begin{bmatrix} a & b \\ b & d \end{bmatrix}$，

这是 $\vw$ 各分量的二次多项式。协方差矩阵与加权和相遇时，就会出现二次型，原因在于下面的恒等式。

3.3.2 每个方向的方差都非负 #

方差不可能为负。因此由上面的推导，协方差矩阵必须对每个 $\vw$ 都满足 $\vw^\T\mSigma\vw \ge 0$。这一性质有专门的名称。

取 $\mA = \mSigma$，正定性的含义是：空间中每个方向的方差都为正，变量的任何加权组合都不是确切已知的。仅为半正定的协方差矩阵，存在某个方差为零的组合，即一个固定不变的组合。在高斯过程中，两个观测输入重合时就会出现这种情况：此时 $f(\vx)$ 与 $f(\vx')$ 是同一个随机变量，$f(\vx) - f(\vx')$ 恰好为零。下面的每个结果都假定严格正定，因此各种实现都会在对角线上加一点噪声或“抖动项”（第 8.4 节）。

二维中这一条件很容易检验。标准差为 $\sigma_1, \sigma_2$、相关系数为 $\rho$ 的协方差矩阵是正定的，当且仅当 $\sigma_1, \sigma_2 > 0$ 且 $-1 < \rho < 1$，而且这一范围内的每个相关系数都是可能的。从三维开始，这一条件就不够了。

3.3.3 无法共存的三个相关系数 #

假设三个量的方差都是 1，第一个与第二个的相关系数为 $0.8$，与第三个的相关系数也为 $0.8$，而第二个与第三个不相关：

每个相关系数都在 $-1$ 与 $1$ 之间，每一对单独来看都是有效的二维协方差。然而，取权重 $\vw = (-2, 1, 1)^\T$，得

即量 $x_2 + x_3 - 2x_1$ 的方差为负。任何三个随机变量都不可能具有这样的相关系数。直观地说，如果 $x_1$ 与 $x_2$ 紧密地共同变化，又与 $x_3$ 紧密地共同变化，那么 $x_2$ 与 $x_3$ 必然在一定程度上共同变化，不可能互不相关。二者必须相关到什么程度，由习题 3.2 求出。

这并非罕见的特例。在金融中，许多资产之间的相关系数表可能不是半正定的；如何把这样的表修正为最接近的有效相关矩阵，已有专门的文献（Higham，2002）。对高斯过程而言，这正是核函数不能是任意相似度分数的原因。一条看似自然的规则，例如“两个输入距离小于 1 时完全相关，否则不相关”，在输入 $0$、$0.6$ 与 $1.2$ 上就会失效：前两个完全相关，后两个完全相关，而第一个与最后一个不相关，这是同一类不可能的三元组。第 9 章的核函数都经过构造，保证这种情况永远不会发生：一个函数是有效的核函数，当且仅当对任意选取的输入，其矩阵都是半正定的（第 7.2 节）。

在下一节的图 3.2 中，读者可以尝试构造这样的矩阵，观察问题出在何处。

3.4 特征向量：矩阵拉伸的轴 #

仅从矩阵的元素，很难一眼看出它的作用。在图 3.1 中，大多数向量 $\vx$ 被 $\mA$ 既转动又拉伸。不过，有些方向只被拉伸：像 $\mA\vx$ 沿着 $\vx$ 本身的方向。找到这些方向，就能把这个映射描述为一组互相独立的拉伸。

特征向量是矩阵不转动的方向；特征值是这一方向被拉伸的倍数，负的特征值使方向反转。在图 3.1 中，沿圆周拖动 $\vx$，直到 $\vx$ 与 $\mA\vx$ 对齐，读数就会给出特征值。默认矩阵在 $18°$ 附近（特征值为 1.4）和 $135°$ 附近（特征值为 0.6）出现对齐。打开显示特征向量，可以看到这两个方向。旋转预设根本没有实特征向量，因为它转动了每一个方向。

3.4.1 对称矩阵 #

一般矩阵的特征向量不一定互相垂直，默认矩阵就是一例。协方差矩阵只可能是对称矩阵，而对称矩阵的性质要好得多。

证明特征值为实数所需的篇幅，超出了它在这里应占的分量（Strang，2016），但垂直性的证明很简短。

这个分解从右往左读，就是一组操作步骤。$\mathbf{U}^\T$ 转动空间，使特征向量与坐标轴对齐；$\bm{\Lambda}$ 把每个轴按对应的特征值拉伸；$\mathbf{U}$ 再把空间转回去。对称矩阵就是一组互相垂直的拉伸，仅此而已。在图 3.1 中选择对称预设，并打开特征向量：它们互相垂直，而且正是椭圆的轴。

3.4.2 特征值与正定性 #

谱定理把正定性的定义转化为关于特征值的条件。

对于协方差矩阵，这个推导还能说明更多。取长度为一的 $\vw$，即一个方向。它沿各特征向量的坐标满足 $\sum_i c_i^2 = 1$，所以方向 $\vw$ 上的方差是各特征值以 $c_i^2$ 为权重的加权平均。特征值就是沿各特征向量的方差，其他任何方向的方差都介于最小与最大的特征值之间。特征值最大的特征向量，就是不确定性最宽的方向。不确定性的形状是一个椭球，其轴沿各特征向量，半轴长为 $\sqrt{\lambda_i}$（即沿各轴的标准差）；第 4.2.1 节将针对高斯分布推导这一点。

3.4.3 三维中的协方差 #

在二维中，椭圆概括了一切。在三维中，会出现两种平面无法展示的现象。图 3.2 画出了一个三维协方差的椭球 $\{\vx : \vx^\T\mSigma^{-1}\vx = 1\}$，即距中心一个标准差的点的集合（这里的标准差在什么意义下成立，由第 4.2.1 节精确说明），同时画出它的三根轴，以及它在立方体各面墙上的影子。旁边是三个成对视图：每次取三个坐标中的两个，画出 $\mSigma$ 中相应 $2 \times 2$ 子块的椭圆。

可以尝试以下操作：

• 转动立方体。椭球有三根互相垂直的轴，即特征向量；从大多数角度看，它们都不与任何坐标轴对齐。读数列出各特征值，其平方根就是各轴的半长。
• 比较影子与成对视图。二者是相同的椭圆。协方差椭球在两个坐标所在平面上的影子，就是相应 $2 \times 2$ 子块的椭圆；因此，忽略高斯分布的一个变量，相当于删去它所在的行和列（第 4.4 节）。
• 按下不可能的三元组。成对视图显示三个有效的椭圆，相关系数分别为 0.8、0.8 和 0。椭球消失，红色虚线标出方差会为负的方向，即特征值 $-0.13$ 对应的特征向量。现在慢慢调高 $\rho_{23}$：一旦 $\rho_{23}$ 超过 0.28（滑块上为 0.29），椭球就会重新出现，起初是一张扁平的圆盘，随 $\rho_{23}$ 增大而变厚。
• 按下链式。这里 $x_1$ 与 $x_2$ 相连，$x_2$ 与 $x_3$ 相连，相关系数都是 0.8，而 $\rho_{13} = 0.64$ 是两者之积。椭球变成沿对角线的一根雪茄，有一个特征值远大于其他特征值：大部分不确定性都集中在单一方向上。

第一个教训是，逐对检查并不能证明协方差矩阵有效；是否有效，是整张表的性质。第二个教训是，协方差可能在某个方向上几乎是扁平的，而这个方向无法通过任何单个变量察觉：一个接近奇异的矩阵，其成对视图未必有哪一个看起来是退化的。

3.4.4 核矩阵的特征值 #

第二个教训正是高斯过程的常态。取 $[0, 1]$ 中等间距的 100 个输入，在长度尺度为 0.1 的径向基函数核下，核矩阵在精确算术中是一个 $100 \times 100$ 的正定矩阵。它最大的几个特征值是 23.9、21.2 和 17.5，但特征值衰减极快，100 个中只有 28 个超过 $10^{-10}$。在双精度算术中，最小的那些特征值淹没在舍入误差里；在我们用 NumPy 做的一次运行中，算出的最小特征值甚至略微为负，$-4 \times 10^{-15}$。相邻输入处的函数值相关性极强，以至于 100 维空间中的大多数方向几乎没有方差。在同一次运行中，NumPy 对这个矩阵的 Cholesky 分解失败了；在对角线上加 $10^{-6}$，即第 8.4 节的抖动项，分解即可成功。函数上的光滑先验实际只有远少于 100 个的独立不确定方向，这也是高斯过程能够从少数几次评估中学到光滑函数的部分原因。

最大特征值与最小特征值之比称为条件数（condition number），它衡量矩阵有多接近奇异。双精度数约有 16 位有效十进制数字，用条件数为 $10^k$ 的矩阵求解方程组，可能损失约 $k$ 位（Golub 与 Van Loan，2013）。光滑核函数在间隔很近的输入上，核矩阵的条件数常常接近 $10^{16}$，此时有效数字已一位不剩，因此抖动项或观测噪声是数值上的必需，而不是建模上的选择。

恰好为零的特征值同样携带信息。在第 18.5 节中，由一组成对比较构建的矩阵，其特征值可以统计出这样的方向有多少个：在这些方向上，这些比较对一个人的偏好没有提供任何信息。

3.5 不求逆而解方程组 #

高斯过程回归的公式中满是逆矩阵：式（8.3）中的后验均值 $\vk(\vx)^\T\mK^{-1}\vy$ 与方差 $k(\vx, \vx) - \vk(\vx)^\T\mK^{-1}\vk(\vx)$。按字面理解，就是先求 $\mK$ 的逆，再相乘。数值计算的惯例是绝不这样做。本节解释原因，并说明应当怎样做。

3.5.1 为什么不求逆 #

含有 $\mK^{-1}\vy$ 的公式，从来不需要 $\mK^{-1}$ 本身，需要的是线性方程组 $\mK\bm{\alpha} = \vy$ 的解向量 $\bm{\alpha}$。宁可求解也不求逆，理由有两个。其一，构造逆矩阵的代价更高，是下文所用分解的数倍。其二，算出的逆矩阵不够精确：乘以经过舍入的逆矩阵，比直接用分解求解损失更多的有效数字（Golub 与 Van Loan，2013）。矩阵病态时（核矩阵常常如此），这一差别决定了答案是否可用。

3.5.2 三角方程组 #

有些方程组很容易求解。对角线以上的元素全为零的矩阵是下三角的（lower triangular）。以这样的矩阵构成的方程组 $\mL\vz = \mathbf{b}$，可以从上到下逐个求出分量：第一个方程只含 $z_1$，第二个只含 $z_1$ 和 $z_2$，依此类推。

这种前代（forward substitution）约需 $n^2$ 次运算，而一般方法需要 $n^3$ 次。上三角（upper triangular）方程组的对角线以下为零，用同样的方式从下往上求解，即回代（back substitution）。求解对称正定方程组的整体策略，就是把它化为两个三角方程组。

3.5.3 Cholesky 分解 #

这个因子是矩阵的一种平方根。$1 \times 1$ 矩阵 $[a]$ 的因子是 $[\sqrt{a}]$，当且仅当 $a > 0$ 时存在。$2 \times 2$ 的情形展示了一般算法如何运作，以及它为什么需要正定性。

第 4 步中的量 $d - b^2/a$，是第 3.7 节中 Schur 补的初次出现。对于 $a = \sigma_1^2$、$b = \rho\sigma_1\sigma_2$、$d = \sigma_2^2$ 的协方差矩阵，它等于 $\sigma_2^2(1 - \rho^2)$，即已知第一个变量之后，第二个变量剩余的方差，第 4.5 节将说明这一点。

一般算法以同样的方式，每次填写 $\mL$ 的一列。

第 2 步使这一算法既是分解，也是检验：它成功当且仅当矩阵正定，而且在实践中是代价最低的正定性检验。在图 3.2 中，不可能的三元组在第三列失败，此处 $s$ 为负。这一算法约需 $n^3/3$ 次浮点运算（Golub 与 Van Loan，2013）。

3.5.4 用因子求解 #

有了 $\mA = \mL\mL^\T$，方程组 $\mA\vx = \mathbf{b}$ 就拆成两个三角方程组。先用前代解 $\mL\vz = \mathbf{b}$，再用回代解 $\mL^\T\vx = \vz$。于是 $\mA\vx = \mL(\mL^\T\vx) = \mL\vz = \mathbf{b}$，正是所求。$\mL$ 一旦已知，每个新的右端项只需 $O(n^2)$。这就是算法 8.1：对核矩阵做一次分解，再通过三角求解得到权重 $\bm{\alpha}$ 和每个预测方差。

import numpy as np
from scipy.linalg import cho_factor, cho_solve, solve_triangular

A = np.array([[4.0, 2.0], [2.0, 3.0]])
b = np.array([2.0, 4.0])

L = np.linalg.cholesky(A)                  # 下三角，L @ L.T == A
z = solve_triangular(L, b, lower=True)     # 前代
x = solve_triangular(L.T, z, lower=False)  # 回代
print(x)                                   # [-0.25  1.5 ]

c = cho_factor(A)                          # 同样的计算，打包成现成的函数
print(cho_solve(c, b))                     # [-0.25  1.5 ]

3.5.5 实际的计算代价 #

精确高斯过程能够处理多少个观测，取决于这一立方级的代价。表 3.1 给出了用 NumPy 以双精度测得的耗时，测量所用的是编写本章的笔记本电脑（一台 Apple M5 Pro）；在其他机器上，数字会有所不同，但增长速率不变。

$n$ 翻倍，时间大约变为八倍，与 $n^3$ 的预测一致；求逆的代价约为分解的五倍。一次贝叶斯优化运行的观测很少超过几百个，因此单次分解的代价很低。真正使代价成倍增加的，是拟合核函数自身的设置（例如长度尺度），因为拟合的每一步都要重新分解矩阵（第 9.3 节）。对于更大的数据集，GPyTorch 库用共轭梯度法代替分解，并在图形处理器上运行；共轭梯度法是一种迭代方法，只需要核矩阵与向量的乘积。据其作者报告，这把精确推断的渐近代价从 $O(n^3)$ 降到了 $O(n^2)$（Gardner 等，2018）。

3.5.6 生成样本的平方根 #

Cholesky 因子还提供了一种由独立的随机性构造相关的随机性的方法。如果 $\vz$ 的各坐标相互独立、方差为 1，那么 $\mL\vz$ 的协方差是 $\mL\mL^\T = \mSigma$，这是第 4.3 节将推导的规则 $\Cov[\mA\vz] = \mA\Cov[\vz]\mA^\T$ 的推论。在图 3.2 中打开显示样本，可以看到 160 个这样的点：每个点都是 $\mL\vz$，其中 $\vz$ 是一次抽取的三个均值为 0、方差为 1 的独立随机数（标准正态随机数，即第 4.1.1 节的钟形曲线），这些点合在一起填满了椭球的形状。由于 $\mL$ 是下三角的，第一个坐标只用到 $z_1$，第二个用到 $z_1$ 和 $z_2$，第三个三者都用到，因此每个新坐标都由前面各坐标的随机性加上自身的一份新的随机性构成。第 8.5 节的后验样本就是这样抽取的。

3.6 作为体积的行列式 #

多元高斯分布的密度（第 4.2 节）和用于拟合核函数的边际似然（第 9.3 节）中都含有行列式。行列式的含义是几何性的。

方阵的行列式（determinant）记作 $\det\mA$ 或 $\lvert\mA\rvert$，是映射 $\mA$ 对面积（二维）、体积（三维）以及一般 $d$ 维体积的放大倍数。它的符号记录了这一映射是否像镜子那样翻转了方向。对于 $2 \times 2$ 矩阵，

即两列张成的平行四边形的有向面积。这就是图 3.1 中的黄色平行四边形，即单位正方形的像；奇异预设把它压扁为零面积，反射预设则使行列式变为负。

从这幅图中可以得出三条性质。先后施加两个映射，它们的体积倍数相乘，所以 $\det(\mA\mathbf{B}) = \det\mA \cdot \det\mathbf{B}$，特别地，$\det(\mA^{-1}) = 1/\det\mA$。矩阵可逆，当且仅当其行列式不为零，即它不把空间压扁。对于对称矩阵 $\mA = \mathbf{U}\bm{\Lambda}\mathbf{U}^\T$，旋转 $\mathbf{U}$ 与 $\mathbf{U}^\T$ 保持体积不变，而 $\bm{\Lambda}$ 把第 $i$ 个轴拉伸 $\lambda_i$ 倍，所以

因此，对于协方差矩阵，行列式是衡量不确定性总体大小的单个数，有时称为广义方差（generalized variance）。图 3.2 的椭球半轴为 $\sqrt{\lambda_i}$，所以体积是 $\frac{4}{3}\pi\sqrt{\lambda_1\lambda_2\lambda_3} = \frac{4}{3}\pi \sqrt{\det\mSigma}$。相关性会使它缩小：对于方差为一、相关系数为 $\rho_{12}, \rho_{13}, \rho_{23}$ 的三个变量，

它对于独立的变量等于 1，随椭球变扁而趋向 0，对于那个不可能的三元组则变为负数（$1 - 0.64 - 0.64 = -0.28$）。在第 6.5 节中，这类行列式衡量一组评估能揭示多少关于目标函数的信息，使行列式最大的那组评估，信息量也最大。

3.6.1 由 Cholesky 因子求对数行列式 #

三角矩阵的行列式是其对角元之积。（在 $2 \times 2$ 的公式中，$b = 0$ 时 $ad - bc$ 就是 $ad$。一般地，三角矩阵把第 $i$ 个坐标轴按其第 $i$ 个对角元拉伸，然后做剪切；剪切使点沿平行于其他轴的方向移动，不改变任何体积，正如图 3.1 的剪切预设所展示的。）对 $\mA = \mL\mL^\T$，由行列式的乘积法则得

取对数并非多此一举。有许多极小特征值的核矩阵，其行列式远小于最小的正双精度数（约 $5 \times 10^{-324}$）。对于第 3.4.4 节中 100 个输入的例子，加上使之可以分解的抖动项之后，算出的特征值之积会下溢为零，而式（3.6）给出 $\log\det\mA \approx -1139$，这是一个大小适中的数，且一旦得到因子，计算它几乎没有额外代价。高斯对数密度和对数边际似然的每一种实现，都用这种方式计算行列式（第 4.7 节）。

3.7 分块矩阵与 Schur 补 #

高斯过程同时处理两组变量：已观测的函数值和待预测的函数值。因此它的矩阵自然地分成若干块，如式（8.1）所示。本节建立这类分块的代数，最终得到第 4.5 节所需的那个公式。

3.7.1 矩阵的分块 #

分块的对称矩阵形如

其中 $\mA$ 是 $p \times p$ 的，$\mathbf{D}$ 是 $q \times q$ 的，$\mathbf{B}$ 是 $p \times q$ 的。各块可以像数一样相乘，只是每个乘积都要保持顺序，因为矩阵乘法不满足交换律。对于协方差，$\mA$ 与 $\mathbf{D}$ 是各组内部的协方差，$\mathbf{B}$ 存放两组之间的协方差。

3.7.2 分块消元 #

手工求解二元一次方程组时，用第一个方程从第二个方程中消去第一个未知数。用分块做同样的一步，就得到本节的核心对象。

对于元素为 $a, b, d$ 的 $2 \times 2$ 矩阵，$a$ 的 Schur 补就是 Cholesky 推导中的那个数 $d - b^2/a$。一般地，它是 $\mathbf{D}$ 去掉由第一组解释的部分之后剩余的部分。

由此得出三个推论，本书后文都会用到。

行列式。$\mathbf{E}$ 是对角块为单位矩阵的分块三角矩阵，所以它的行列式为 1，并且 $\det\mathbf{M} = \det\mA \cdot \det\mathbf{S}$。

正定性。对任意向量 $\vw$，令 $\mathbf{v} = \mathbf{E}^{-\T}\vw$；则 $\vw^\T\mathbf{M}\vw = \mathbf{v}_1^\T\mA\mathbf{v}_1 + \mathbf{v}_2^\T\mathbf{S}\mathbf{v}_2$，其中 $\mathbf{v}_1$ 与 $\mathbf{v}_2$ 是 $\mathbf{v}$ 的两部分。所以 $\mathbf{M}$ 正定，当且仅当 $\mA$ 与它的 Schur 补 $\mathbf{S}$ 都正定。这就是不可能的三元组的精确表述：考虑 $x_1$ 之后，$x_2$ 与 $x_3$ 的剩余部分仍必须构成有效的协方差（习题 3.2）。

逆矩阵。对第 3 步逐块求逆，利用 $(\mathbf{E}^{-1})^{-1} = \mathbf{E}$，得到 $\mathbf{M}^{-1} = \mathbf{E}^\T \begin{bmatrix} \mA^{-1} & \mathbf{0} \\ \mathbf{0} & \mathbf{S}^{-1} \end{bmatrix} \mathbf{E}$，乘开得

逆矩阵的右下块是 Schur 补的逆。第 4.5 节将用到这一点：取 $\mA$ 为已观测值的协方差、$\mathbf{D}$ 为未观测值的协方差，证明条件化之后的协方差就是 Schur 补 $\mathbf{D} - \mathbf{B}^\T\mA^{-1}\mathbf{B}$，即先验协方差减去观测所解释的部分。第 8.1 节是同一个公式，只是各块换成了核矩阵；第 B.1 节与第 B.2 节把它与相关的恒等式汇集在一起（Petersen 与 Pedersen，2012）。

3.7.3 增加一个观测 #

从分块的角度，还可以看出如何低成本地更新 Cholesky 因子，而贝叶斯优化器每次评估之后都需要这样做。假设对目前已观测的 $n$ 个输入有 $\mK = \mL\mL^\T$，现在来了一个新输入，它与旧输入的协方差为 $\vk$，自身方差为 $\kappa$。新的核矩阵及其因子具有如下分块形式

比较各块，得 $\mL\mathbf{l} = \vk$，即一次代价为 $O(n^2)$ 的前代；以及 $l_{\ast}^2 = \kappa - \mathbf{l}^\T\mathbf{l}$。由于 $\mathbf{l}^\T\mathbf{l} = \vk^\T\mL^{-\T}\mL^{-1}\vk = \vk^\T\mK^{-1}\vk$，新对角元的平方是 $\kappa - \vk^\T\mK^{-1}\vk$，即 $\mK$ 的 Schur 补，它也是高斯过程在新输入处的后验方差（式（8.3））。只要核函数的设置保持不变，把因子扩大一行只需 $O(n^2)$，而不必从头花费 $O(n^3)$。如果没有观测噪声，而新输入又与某个旧输入重合，$l_{\ast}$ 就为零，分解随之失败：矩阵已经变为奇异，原因见第 3.3.2 节。

3.8 习题 #

延伸阅读 #

• Strang（2016）耐心地介绍了作为映射的矩阵、特征值与正定矩阵，与本章一样侧重几何。
• Golub 与 Van Loan（2013）是讲解这些计算如何在浮点运算中完成的标准参考书，内容包括三角方程组、Cholesky 分解、运算量与条件数。
• Rasmussen 与 Williams（2006）的附录 A 收录了高斯过程回归所用的矩阵恒等式和基于 Cholesky 分解的计算。
• Petersen 与 Pedersen（2012）是一本紧凑的恒等式手册，包括分块求逆与行列式的公式，可用于核对推导。
• Sanderson（2016）是一套动画视频，展示矩阵在平面上和空间中的作用，也就是图 3.1 背后的图景。
• Higham（2002）讨论如何修正无效的相关矩阵，即第 3.3.3 节在实际中的一面。