说个真实经历。有阵子我在调一个高斯过程回归的代码深夜弹出一个报错matrix is not positive definite。第一反应是数据预处理出了问题查了半天没查到最后才发现是我构造的核矩阵在数值上差了一丁点最小特征值到了-3e-15这个量级。正是从那天起我认真把半正定矩阵PSDPositive Semi-Definite Matrix从头梳理了一遍。这个东西看似只是矩阵论里的小概念实际上统计、机器学习、最优化、数值计算全都在用它你会在各种论文和开源代码里反复撞见它的身影。这篇文章我打算用一篇完整的篇幅从直觉、定义、判定方法、性质到工程里的实际坑位把半正定矩阵一次讲透适合正在学矩阵论、机器学习算法或优化理论的朋友参考。1. 从二次型开始PSD矩阵到底在描述什么1.1 二次型的几何直觉半正定矩阵的严格定义绕不开二次型。给定一个对称矩阵 A二次型就是这样的标量表达式x^T A x其中 x 是一个列向量。展开到二维就是a11 * x1^2 2 * a12 * x1 * x2 a22 * x2^2这个式子背后有非常直观的物理含义可以把它看成在方向 x 上测量的能量或长度。如果它恒不小于零说明这个矩阵在任意方向上给出的测量结果都不会出现负数。我自己的理解方式是用椭圆。假设 A 是正定矩阵集合{x : x^T A x ≤ 1}在二维平面上是一个椭圆在三维空间里是一个椭球。这个椭球的每一个轴的方向由 A 的特征向量决定轴的长度和对应特征值相关。所以 A 的特征值大椭球在该方向就拉长特征值小就压扁如果某个特征值正好是零椭球在该方向就直接塌掉变成了一个低维薄片。半正定矩阵的含义就是允许出现这种塌掉的方向但绝不允许任何一个方向出现负数能量。这是它和正定矩阵唯一的区别也是半字的由来。1.2 为什么非负是底线很多刚接触这个概念的人会疑惑为什么我们这么在意一个二次型是否非负原因是在实际应用里这个量通常对应方差、距离平方、能量、KL散度或概率密度的指数项这些东西在最基本的物理和数学直觉上就不能为负。举个最直白的例子马氏距离的平方是d^2 (x - μ)^T Σ^{-1} (x - μ)如果 Σ 不是正定矩阵这个值就可能为负距离变成负数就没有任何意义了。再比如多元高斯的概率密度指数部分同样是二次型分母里还有det(Σ)^(1/2)如果 Σ 不正定整个密度函数就失去了合法性。在优化问题中Hessian 矩阵必须是 PSD 才能保证函数局部是凸的。Hessian 如果在某个方向给出负的曲率那就说明你站的不是一个局部极小点而是一个鞍点。所以非负不是一个多余的数学洁癖而是现实世界对合法测量的基本要求。2. 定义与特征值为什么对称和非负是这家的门牌号2.1 对称性是前提教科书上写的定义是一个对称矩阵 A如果对任意非零向量 x 都有x^T A x ≥ 0就称 A 是半正定矩阵记作A ⪰ 0。之所以必须强调对称是因为非对称矩阵的二次型只取决于它的对称部分。任何一个方阵 M 都可以拆成M (M M^T) / 2 (M - M^T) / 2右边第一项叫对称部分第二项叫反对称部分。代入二次型x^T M x反对称部分贡献恰好为零。这个可以用一维类比来理解奇数函数在整个实数轴上的积分是零反对称矩阵和向量的二次型也类似绕一圈就抵消掉了。所以工程实践里有一个很常见的低级错误传入一个几乎对称但又不完全对称的矩阵比如数据对齐时产生了浮点误差然后直接用特征值分解去判断半正定性。这时候不仅特征值会变成复数后续的 Cholesky 分解也可能直接报错。正确做法是先用(A A.T) / 2做一个显式的对称化再往下走。2.2 谱分解把定义翻译成特征值对称矩阵有个极其重要的性质就是一定能正交对角化A Q Λ Q^T其中 Q 是正交矩阵Q^T Q IΛ 是对角矩阵。把它代入二次型x^T A x x^T Q Λ Q^T x令y Q^T x由于 Q 是正交变换y 的长度和 x 一样但方向被旋转到 A 的特征向量坐标系里。于是x^T A x Σ λ_i * y_i^2注意y_i^2永远是非负的。所以一个对称矩阵是半正定的当且仅当它的所有特征值都是非负的是正定的当且仅当所有特征值都严格大于零。这个结论是理解 PSD 矩阵最核心的桥梁后面的主子式法、Cholesky 法等全部都是这个结论的旁支。2.3 实际算一个例子我常用的例子是A [[2, -1], [-1, 2]]特征多项式是(2 - λ)^2 - 1 0展开得到λ^2 - 4λ 3 0所以特征值是 1 和 3都大于零这是一个正定矩阵。几何上它对应的椭圆长轴短轴比例恰好是√3 : 1。再来看一个半正定但奇异的例子B [[1, 1], [1, 1]]特征值是 2 和 0。有一个方向的能量是零——具体来说是沿着(1, -1)方向x^T B x 0。这种情况正是半正定和正定之间最本质的差别半正定允许零特征值也就是允许信息在某些方向上完全消失。3. 四种实操判定法别再只背定义了3.1 特征值检验法最通用、最稳妥的判断方法就是直接求特征值。Python 里要用eigvalsh而不是eigvals前者专门针对对称矩阵做了优化输出严格实数速度也更快。import numpy as np A np.array([[2., -1.], [-1., 2.]]) eigvals np.linalg.eigvalsh(A) tol 1e-8 is_psd np.all(eigvals -tol) print(eigenvalues:, eigvals) print(PSD:, is_psd)容差tol的选取很有讲究。严格数学意义上只要特征值不是负数就是 PSD但浮点计算误差会导致本来为 0 的特征值变成-1e-15这种小负数。如果直接用 0判断一个数值上完美的半正定矩阵都可能被判负。所以在实际工程里我习惯把容差设成1e-8或者max(1e-8, n * eps * max_eigval)这种尺度。3.2 主子式判定法与一个经典陷阱正定矩阵的判断有一个常用的简化办法顺序主子式全为正。但很多人把它直接搬到半正定上写成顺序主子式全为非负这其实是不对的。来看一个反例C [[0, 0], [0, -1]]它的第一个顺序主子式是 0第二个顺序主子式是行列式 0都满足非负但这个矩阵显然不是半正定因为特征值是 0 和 -1。问题出在哪顺序主子式只检查了从左上角开始逐层扩大的子矩阵而没有检查删除任意行列后留下的子矩阵。半正定的正确判定条件是所有主子式包括去掉某些行和对应列后剩下的子式都必须非负。对对角矩阵来说主子式里包含了每一个对角线元素本身所以一旦某个对角线元素为负检查就能发现。但如果你的代码只写了顺序主子式这个 -1 就会被漏掉。3.3 Cholesky 分解与 Gram 矩阵稍微懂一点数值线性代数的朋友会喜欢用 Cholesky 分解。一个正定矩阵存在唯一的分解A L L^T其中 L 是下三角矩阵且对角元素为正。如果矩阵只是半正定但奇异标准的 Cholesky 分解会失败因为分解过程中会遇到零主元。判断逻辑很简单from scipy.linalg import cholesky try: L cholesky(A, lowerTrue) print(PD) except Exception: print(not PD, could be PSD or indefinite)注意这个方法检测的是正定性不是半正定性。对于半正定矩阵更稳妥的做法是 LDLT 分解scipy.linalg.ldl可以把A分解成L D L^T然后检查 D 的对角元素是否全部非负。另外还有一个很实用的等价条件任意矩阵 MM^T M一定是半正定矩阵反过来任何半正定矩阵 A 都可以被分解成B^T B的形式。这条性质在机器学习里几乎随处可见协方差矩阵是X^T X核矩阵是特征映射的内积矩阵它们天生就是 PSD 结构。3.4 快速判定流程建议如果只想快速知道这个矩阵能不能用我建议按这个顺序来检查矩阵是否对称必要时做(A A.T) / 2对称化。用numpy.linalg.eigvalsh看最小特征值。如果最小特征值大于1e-8直接当正定用。如果最小特征值落在[-1e-8, 1e-8]区间按半正定处理注意是否需要加扰动来解决奇异性问题。如果最小特征值明显小于零先去查数据生成流程而不是急着修复矩阵。4. 一张表理清PSD矩阵的家底性质与不等式4.1 性质对照表半正定矩阵的性质多而杂我把日常最常用到的整理成一张表格。建议收藏用到哪条翻哪条。性质说明特征值全非负定义的等价表述也是判断的黄金标准迹非负行列式非负必要不充分条件只能用来排除不能用来确认所有主子式非负等价条件比顺序主子式严格可写成B^T B形式构造性理解机器学习中大量使用PSD 矩阵之和仍是 PSD对加法封闭加正则项后仍合法M^T A M保持 PSD合同变换不破坏半正定性数据变换后协方差矩阵依然合法存在唯一的 PSD 平方根A A^{1/2} A^{1/2}且A^{1/2}也半正定奇异值等于特征值对 PSD 矩阵来说SVD 和特征分解结果一致秩等于正特征值个数零特征值的数量就是秩亏缺的维度第 6 条M^T A M保持 PSD 特别重要。比如原始协方差矩阵是 PSD你对数据做了线性变换乘以矩阵 M新的协方差矩阵M^T A M依然是 PSD。这保证了经过 PCA、白化、旋转等一系列操作后协方差矩阵不会突然变得不合规。4.2 Rayleigh 商与特征值区间在做 PCA、谱聚类、特征值问题的时候Rayleigh 商是个绕不开的工具。它的定义是R(x) x^T A x / x^T x对于对称矩阵Rayleigh 商有一个很漂亮的性质它总是落在最小特征值和最大特征值之间λ_min ≤ R(x) ≤ λ_max证明非常简单把A Q Λ Q^T代入R(x) 实际上就是各个特征值按y_i^2加权平均的结果权重之和为 1所以不可能跳出特征值的区间范围。这也是 PCA 里第一个主成分方向对应最大特征值特征向量的理论来源。顺便提一个不太常见但很有用的不等式如果 A 和 B 都是 PSD 矩阵那么Tr(AB) ≥ 0。注意 AB 本身不一定是对称矩阵但它的迹一定非负。这个结论在统计里面的 Stein 引理、信息几何和协方差估计分析中经常出现。证明思路是把 B 写成B^{1/2} B^{1/2}然后利用Tr(AB) Tr(B^{1/2} A B^{1/2}) ≥ 0被作用后的矩阵仍然是 PSD迹自然非负。5. PSD不是摆设统计、优化、机器学习里它都在场5.1 协方差矩阵与高维统计统计学里最常见的 PSD 矩阵就是样本协方差矩阵。假设数据矩阵 X 已经中心化维度是 n×p那么样本协方差S X^T X / (n-1)天然就是半正定矩阵因为任意向量 v 都有v^T S v ||X v||^2 / (n-1) ≥ 0那句任意向量的二次型等于某个范数的平方直接保证了非负性。当样本量 n 大于特征数 p 的时候X 通常列满秩S 是正定的。但在高维统计里经常出现 p 比 n 还大的情况此时 X 一定是秩亏的S 必然奇异。这意味着数据在某些方向上完全没有波动协方差矩阵有零特征值。很多初学者遇到numpy.linalg.LinAlgError: singular matrix时以为是 bug其实这是高维数据的数学本质不是程序写错了。常规修法是加一个小的对角扰动比如用S λI也就是统计里的岭正则化。也可以使用 Ledoit-Wolf 收缩估计把协方差矩阵往单位阵方向收缩。我自己在基因表达数据、用户行为矩阵这类场景里几乎每次都要处理这种问题。5.2 马氏距离与核方法马氏距离是协方差矩阵 PSD 性质的直接受益者。它用协方差矩阵的逆来重新定义距离可以理解成在数据真实的分布形状下测量距离。数据如果在一个方向上方差很小那么沿这个方向的微小偏差都应该被视为很大的距离变化。这在异常检测、聚类和图像检索里都非常实用。核方法则是把 PSD 用在了 Gram 矩阵上。给定一组样本和核函数k(x_i, x_j)Gram 矩阵K_ij k(x_i, x_j)半正定是 Mercer 条件的要求。只有 Gram 矩阵是 PSDSVM、高斯过程回归等算法里的二次规划才有解后验分布才合法。高斯核k(x, y) exp(-γ ||x-y||^2)在理论上是正定核但在数值计算里有个老毛病当数据点彼此距离极近时Gram 矩阵会接近一个全 1 矩阵特征值出现一个大的和一堆接近零的数值上很容易变成半正定但不严格正定。高斯过程社区管这个叫jitter解决办法是在对角线上加一个很小的数比如1e-8 * I。5.3 Hessian 矩阵与凸优化在优化理论中函数 f 的可微性和凸性之间有一个非常经典的等价关系f 是凸函数 ⟺ 对任意自变量 xHessian 矩阵 ∇²f(x) ⪰ 0如果 Hessian 矩阵的所有特征值都存在一个正的下界那么函数是强凸的。强凸的性质比凸更强它保证目标函数有一个唯一的极小点而且梯度下降能收敛得很快。以最小二乘问题为例目标函数是f(x) 1/2 * ||Ax - b||^2它的梯度是A^T A x - A^T bHessian 是A^T A。由于A^T A天然是 PSD所以最小二乘问题无条件是一个凸问题。这个事实解释了为什么我们从来不用担心最小二乘会收敛到局部极小点——它的所有局部极小点都是全局极小点。但是当 A 列秩不满时A^T A奇异最优解不唯一最小二乘的解有无穷多个。加上 L2 正则之后目标 Hessian 变成A^T A λI因为λI是正定的所以整体变成正定解就唯一化了。这就是岭回归能从数学上强制解唯一的原因。5.4 半正定规划SDP的简要面貌到了优化领域再往前走一步就遇到了半正定规划其中的变量本身就是一个 PSD 矩阵。这类问题的典型形式是最小化Tr(CX)条件包括Tr(A_i X) b_i和X ⪰ 0注意这里的未知量 X 是一个矩阵不是向量而约束X ⪰ 0是一个矩阵不等式。这个条件看似抽象其实等价于要求 X 的所有特征值非负。SDP 之所以重要是因为它把线性规划从非负的数轴锥推广到了非负的矩阵锥从而可以优雅地处理矩阵结构问题。组合优化里的 MaxCut 松弛、控制系统里的 LMI 分析、量子信息里的密度矩阵约束最后都可以转化成 SDP。如果你以后读优化方向的论文遇到把问题松弛成 SDP不必恐慌它背后的核心约束就是这一章讲的东西找到一个 PSD 矩阵同时满足若干线性等式。6. 分块矩阵视角Schur补与半正定的等价关系6.1 分块 PSD 判定Schur 补实际项目中遇到的大型 PSD 矩阵几乎都是分块结构的。协方差矩阵天然按变量分组分块高斯过程里的核矩阵按数据批次分块控制论里的李雅普诺夫矩阵也是分块形式。所以掌握分块视角非常实用。设一个对称分块矩阵A [[A11, A12], [A12^T, A22]]其中 A11 是 k×k 的可逆子块。如果 A 是 PSD那么 A11 和 A22 必然都是 PSD。这个可以直接验证取x (v, 0)二次型就退化成v^T A11 v非负性自动保留。但如果 A11 正定A 整体是否 PSD取决于一个关键对象——Schur 补S A22 - A12^T A11^{-1} A12结论是当 A11 正定时A 半正定当且仅当 S 半正定。证明思路是做一个合同变换用一个下三角消元矩阵 L 对 A 做变换L [[I, 0], [-A12^T A11^{-1}, I]]可以得到L A L^T [[A11, 0], [0, S]]因为合同变换不改变矩阵的正负惯性指数所以 A 的 PSD 性质完全由对角块 A11 和 S 决定。这个等价关系非常干净实际用起来也很顺手。6.2 分块求逆公式中的老熟人分块求逆公式大家可能都见过但未必注意到 Schur 补在里面反复出现。对于一个分块矩阵A^{-1} [[A11^{-1} A11^{-1} A12 S^{-1} A12^T A11^{-1}, -A11^{-1} A12 S^{-1}], [-S^{-1} A12^T A11^{-1}, S^{-1}]]右下角正好是 Schur 补的逆。这个结构在多元高斯条件分布中大有用处已知部分变量后剩余变量的条件协方差矩阵恰好就是 Schur 补。通俗地说Schur 补量化了把 A11 的信息消掉之后A22 还剩下多少自己的波动。如果 S 变成零矩阵说明 A22 的信息完全被 A11 解释了这对应数据里的确定性关系。顺手提一个热搜词里反复出现的分块矩阵的 n 次方。如果分块矩阵是对角的幂次直接对每个块分别求即可。一般分块矩阵求幂没有简单普适公式但可以对角化后分别求幂值。真正常用的场景是分块对角矩阵加上低秩扰动用 Woodbury 公式来做迭代更新。6.3 Sylvester 定律乘积特征值的非零部分分块和乘积常常绑定在一起这里有一个少为人知但极其有用的结论对任意两个方阵 A 和 BAB 和 BA 的非零特征值完全相同这就是 Sylvester 定律。如果 A、B 都是 PSD 矩阵AB 虽然不一定对称但它的特征值全是实数且非负。原因可以构造AB 与 A^{1/2} B A^{1/2}有相同的非零特征值而后者是 PSD 矩阵特征值自然全非负。这条性质在控制论判稳、卡尔曼滤波的可观测性分析以及很多统计检验中都有应用。如果你以后遇到某个算法的稳定性分析中出现AB 的特征值必须落在右半平面之类的说法它背后的逻辑很可能就是这条定律。7. 数值实践中的PSD陷阱与应对策略7.1 那些幽灵般的负特征值从哪来工程里最常见的症状是理论上应该 PSD 的矩阵一算特征值却出现几个小的负数。我根据经验把来源分成三类。第一类是浮点舍入误差。对称矩阵经过复杂的浮点运算后本来精确为零的特征值可能变成-1e-15这是完全正常的不需要过度处理。第二类是秩亏。样本量不足、特征重复、中心化后数据之间存在线性关系都会让矩阵产生精确的零特征值浮点运算下表现为微负值。第三类是超参数选择不当。比如高斯核的带宽取得过小Gram 矩阵趋近于全 1 矩阵数值上出现大量接近零的特征值再加上舍入误差就会报not positive definite。判断一个小负数到底是数值误差还是真正的负特征值我的经验阈值是负值绝对量级小于n * eps * λ_max就可以当作零来处理如果负值达到1e-3这个量级那一定是数据或模型结构出了问题加正则项只是掩耳盗铃。7.2 处理前一定要先对称化很多从计算流程里拿出来的矩阵细微处已经不对称了。比如两个矩阵乘积M^T N如果不相等再做一些浮点运算结果会和严格预期差出1e-16左右。此时不管用eigvalsh还是cholesky行为都不可预测因为eigvalsh默认只读下三角部分矩阵不对称时它会忽略掉上三角的信息。所以我的处理习惯是任何矩阵在做 PSD 判断之前先执行A_sym (A A.T) / 2这不仅让后续分解稳定也避免了因为A和A.T细微不一致而导致的诡异行为。7.3 三个常用修复方案如果确认矩阵需要修复我按优先级推荐三种方案。第一种是对角扰动也就是 jitterA_fixed A_sym 1e-8 * np.eye(n)这是最简单、最常用的方案适用于核矩阵、协方差矩阵这类加一点噪声不影响结果的场景。高斯过程里默认这么做原因就是高斯核在数据点靠得很近时协方差矩阵数值上几乎奇异。第二种是特征值截断。把负特征值强行置零再重构矩阵eigvals, eigvecs np.linalg.eigh(A_sym) eigvals_fixed np.maximum(eigvals, 0) A_fixed eigvecs np.diag(eigvals_fixed) eigvecs.T这种做法的代价是改变了原始矩阵的谱适合你不在乎特征值结构只在乎矩阵合法且距离原矩阵尽量近的场景。第三种是更精细的收缩估计类似 Ledoit-Wolf。它不是对每个矩阵都加同一个常数而是按特征值大小做自适应收缩A_shrunk (1 - α) A α * (Tr(A) / n) I其中 α 可以通过交叉验证选择。这个方法在协方差估计中表现远好于固定对角扰动因为它在保留信号和消除噪声之间自动找平衡。7.4 奇异半正定矩阵求解线性系统当 A 只是半正定但奇异时np.linalg.solve会报错但问题本身往往有解只是解不唯一。最常用的替代方案是伪逆x np.linalg.pinv(A) b它给出的是最小范数意义下的解也就是所有可行解中||x||最小的那个。在统计中这个解恰好对应线性回归里当X^T X奇异时的岭估计极限形式。另一个选择是np.linalg.lstsq(A, b)它在最小二乘意义下自动处理秩亏问题。如果矩阵规模很大不能用直接法可以考虑迭代法。共轭梯度法要求矩阵对称正定半正定加上对角扰动之后就可以安全使用了。注意扰动不要太大否则解会偏离原问题太多。7.5 一次真实的排查经历回到开头那个高斯过程回归的例子。我当时的排查链路是这样的先打印最小特征值发现是-3e-15。检查数据确认没有缺失值量纲差异不大。检查核函数带宽发现很小导致 Gram 矩阵几乎全 1。于是判定这是数值秩亏而不是真实负特征值。在核矩阵对角线上加了1e-8的 jitter问题解决。但如果当时看到的最小特征值是-0.1我会回头检查原始数据是否存在重复列、是否中心化之后仍有线性相关变量、核带宽是否小到令两点间相似度几乎为 1。加 jitter 能救数值问题但救不了真负特征值因为后者说明你的模型或数据在某个方向上能量为负这在物理上是不应该存在的。最后分享一个我在工程里的习惯任何矩阵到我手里先做三步检查——对称化、看最小特征值、尝试 Cholesky。如果看到一个-1e-15量级的负特征值我直接当它是零不需要大动干戈但如果负特征值到了-1e-3甚至更大先别急着加 jitter回头检查数据生成流程十有八九是秩亏或数据没对齐。理解半正定矩阵的真正价值不是让你背一堆定理而是让你在碰到matrix is not positive definite时能准确判断这到底是数学问题、数值问题还是数据问题并且知道下一步该怎么做。 半正定矩阵（PSD）从直觉到工程：判定方法、性质与数值陷阱
说个真实经历。有阵子我在调一个高斯过程回归的代码深夜弹出一个报错matrix is not positive definite。第一反应是数据预处理出了问题查了半天没查到最后才发现是我构造的核矩阵在数值上差了一丁点最小特征值到了-3e-15这个量级。正是从那天起我认真把半正定矩阵PSDPositive Semi-Definite Matrix从头梳理了一遍。这个东西看似只是矩阵论里的小概念实际上统计、机器学习、最优化、数值计算全都在用它你会在各种论文和开源代码里反复撞见它的身影。这篇文章我打算用一篇完整的篇幅从直觉、定义、判定方法、性质到工程里的实际坑位把半正定矩阵一次讲透适合正在学矩阵论、机器学习算法或优化理论的朋友参考。1. 从二次型开始PSD矩阵到底在描述什么1.1 二次型的几何直觉半正定矩阵的严格定义绕不开二次型。给定一个对称矩阵 A二次型就是这样的标量表达式x^T A x其中 x 是一个列向量。展开到二维就是a11 * x1^2 2 * a12 * x1 * x2 a22 * x2^2这个式子背后有非常直观的物理含义可以把它看成在方向 x 上测量的能量或长度。如果它恒不小于零说明这个矩阵在任意方向上给出的测量结果都不会出现负数。我自己的理解方式是用椭圆。假设 A 是正定矩阵集合{x : x^T A x ≤ 1}在二维平面上是一个椭圆在三维空间里是一个椭球。这个椭球的每一个轴的方向由 A 的特征向量决定轴的长度和对应特征值相关。所以 A 的特征值大椭球在该方向就拉长特征值小就压扁如果某个特征值正好是零椭球在该方向就直接塌掉变成了一个低维薄片。半正定矩阵的含义就是允许出现这种塌掉的方向但绝不允许任何一个方向出现负数能量。这是它和正定矩阵唯一的区别也是半字的由来。1.2 为什么非负是底线很多刚接触这个概念的人会疑惑为什么我们这么在意一个二次型是否非负原因是在实际应用里这个量通常对应方差、距离平方、能量、KL散度或概率密度的指数项这些东西在最基本的物理和数学直觉上就不能为负。举个最直白的例子马氏距离的平方是d^2 (x - μ)^T Σ^{-1} (x - μ)如果 Σ 不是正定矩阵这个值就可能为负距离变成负数就没有任何意义了。再比如多元高斯的概率密度指数部分同样是二次型分母里还有det(Σ)^(1/2)如果 Σ 不正定整个密度函数就失去了合法性。在优化问题中Hessian 矩阵必须是 PSD 才能保证函数局部是凸的。Hessian 如果在某个方向给出负的曲率那就说明你站的不是一个局部极小点而是一个鞍点。所以非负不是一个多余的数学洁癖而是现实世界对合法测量的基本要求。2. 定义与特征值为什么对称和非负是这家的门牌号2.1 对称性是前提教科书上写的定义是一个对称矩阵 A如果对任意非零向量 x 都有x^T A x ≥ 0就称 A 是半正定矩阵记作A ⪰ 0。之所以必须强调对称是因为非对称矩阵的二次型只取决于它的对称部分。任何一个方阵 M 都可以拆成M (M M^T) / 2 (M - M^T) / 2右边第一项叫对称部分第二项叫反对称部分。代入二次型x^T M x反对称部分贡献恰好为零。这个可以用一维类比来理解奇数函数在整个实数轴上的积分是零反对称矩阵和向量的二次型也类似绕一圈就抵消掉了。所以工程实践里有一个很常见的低级错误传入一个几乎对称但又不完全对称的矩阵比如数据对齐时产生了浮点误差然后直接用特征值分解去判断半正定性。这时候不仅特征值会变成复数后续的 Cholesky 分解也可能直接报错。正确做法是先用(A A.T) / 2做一个显式的对称化再往下走。2.2 谱分解把定义翻译成特征值对称矩阵有个极其重要的性质就是一定能正交对角化A Q Λ Q^T其中 Q 是正交矩阵Q^T Q IΛ 是对角矩阵。把它代入二次型x^T A x x^T Q Λ Q^T x令y Q^T x由于 Q 是正交变换y 的长度和 x 一样但方向被旋转到 A 的特征向量坐标系里。于是x^T A x Σ λ_i * y_i^2注意y_i^2永远是非负的。所以一个对称矩阵是半正定的当且仅当它的所有特征值都是非负的是正定的当且仅当所有特征值都严格大于零。这个结论是理解 PSD 矩阵最核心的桥梁后面的主子式法、Cholesky 法等全部都是这个结论的旁支。2.3 实际算一个例子我常用的例子是A [[2, -1], [-1, 2]]特征多项式是(2 - λ)^2 - 1 0展开得到λ^2 - 4λ 3 0所以特征值是 1 和 3都大于零这是一个正定矩阵。几何上它对应的椭圆长轴短轴比例恰好是√3 : 1。再来看一个半正定但奇异的例子B [[1, 1], [1, 1]]特征值是 2 和 0。有一个方向的能量是零——具体来说是沿着(1, -1)方向x^T B x 0。这种情况正是半正定和正定之间最本质的差别半正定允许零特征值也就是允许信息在某些方向上完全消失。3. 四种实操判定法别再只背定义了3.1 特征值检验法最通用、最稳妥的判断方法就是直接求特征值。Python 里要用eigvalsh而不是eigvals前者专门针对对称矩阵做了优化输出严格实数速度也更快。import numpy as np A np.array([[2., -1.], [-1., 2.]]) eigvals np.linalg.eigvalsh(A) tol 1e-8 is_psd np.all(eigvals -tol) print(eigenvalues:, eigvals) print(PSD:, is_psd)容差tol的选取很有讲究。严格数学意义上只要特征值不是负数就是 PSD但浮点计算误差会导致本来为 0 的特征值变成-1e-15这种小负数。如果直接用 0判断一个数值上完美的半正定矩阵都可能被判负。所以在实际工程里我习惯把容差设成1e-8或者max(1e-8, n * eps * max_eigval)这种尺度。3.2 主子式判定法与一个经典陷阱正定矩阵的判断有一个常用的简化办法顺序主子式全为正。但很多人把它直接搬到半正定上写成顺序主子式全为非负这其实是不对的。来看一个反例C [[0, 0], [0, -1]]它的第一个顺序主子式是 0第二个顺序主子式是行列式 0都满足非负但这个矩阵显然不是半正定因为特征值是 0 和 -1。问题出在哪顺序主子式只检查了从左上角开始逐层扩大的子矩阵而没有检查删除任意行列后留下的子矩阵。半正定的正确判定条件是所有主子式包括去掉某些行和对应列后剩下的子式都必须非负。对对角矩阵来说主子式里包含了每一个对角线元素本身所以一旦某个对角线元素为负检查就能发现。但如果你的代码只写了顺序主子式这个 -1 就会被漏掉。3.3 Cholesky 分解与 Gram 矩阵稍微懂一点数值线性代数的朋友会喜欢用 Cholesky 分解。一个正定矩阵存在唯一的分解A L L^T其中 L 是下三角矩阵且对角元素为正。如果矩阵只是半正定但奇异标准的 Cholesky 分解会失败因为分解过程中会遇到零主元。判断逻辑很简单from scipy.linalg import cholesky try: L cholesky(A, lowerTrue) print(PD) except Exception: print(not PD, could be PSD or indefinite)注意这个方法检测的是正定性不是半正定性。对于半正定矩阵更稳妥的做法是 LDLT 分解scipy.linalg.ldl可以把A分解成L D L^T然后检查 D 的对角元素是否全部非负。另外还有一个很实用的等价条件任意矩阵 MM^T M一定是半正定矩阵反过来任何半正定矩阵 A 都可以被分解成B^T B的形式。这条性质在机器学习里几乎随处可见协方差矩阵是X^T X核矩阵是特征映射的内积矩阵它们天生就是 PSD 结构。3.4 快速判定流程建议如果只想快速知道这个矩阵能不能用我建议按这个顺序来检查矩阵是否对称必要时做(A A.T) / 2对称化。用numpy.linalg.eigvalsh看最小特征值。如果最小特征值大于1e-8直接当正定用。如果最小特征值落在[-1e-8, 1e-8]区间按半正定处理注意是否需要加扰动来解决奇异性问题。如果最小特征值明显小于零先去查数据生成流程而不是急着修复矩阵。4. 一张表理清PSD矩阵的家底性质与不等式4.1 性质对照表半正定矩阵的性质多而杂我把日常最常用到的整理成一张表格。建议收藏用到哪条翻哪条。性质说明特征值全非负定义的等价表述也是判断的黄金标准迹非负行列式非负必要不充分条件只能用来排除不能用来确认所有主子式非负等价条件比顺序主子式严格可写成B^T B形式构造性理解机器学习中大量使用PSD 矩阵之和仍是 PSD对加法封闭加正则项后仍合法M^T A M保持 PSD合同变换不破坏半正定性数据变换后协方差矩阵依然合法存在唯一的 PSD 平方根A A^{1/2} A^{1/2}且A^{1/2}也半正定奇异值等于特征值对 PSD 矩阵来说SVD 和特征分解结果一致秩等于正特征值个数零特征值的数量就是秩亏缺的维度第 6 条M^T A M保持 PSD 特别重要。比如原始协方差矩阵是 PSD你对数据做了线性变换乘以矩阵 M新的协方差矩阵M^T A M依然是 PSD。这保证了经过 PCA、白化、旋转等一系列操作后协方差矩阵不会突然变得不合规。4.2 Rayleigh 商与特征值区间在做 PCA、谱聚类、特征值问题的时候Rayleigh 商是个绕不开的工具。它的定义是R(x) x^T A x / x^T x对于对称矩阵Rayleigh 商有一个很漂亮的性质它总是落在最小特征值和最大特征值之间λ_min ≤ R(x) ≤ λ_max证明非常简单把A Q Λ Q^T代入R(x) 实际上就是各个特征值按y_i^2加权平均的结果权重之和为 1所以不可能跳出特征值的区间范围。这也是 PCA 里第一个主成分方向对应最大特征值特征向量的理论来源。顺便提一个不太常见但很有用的不等式如果 A 和 B 都是 PSD 矩阵那么Tr(AB) ≥ 0。注意 AB 本身不一定是对称矩阵但它的迹一定非负。这个结论在统计里面的 Stein 引理、信息几何和协方差估计分析中经常出现。证明思路是把 B 写成B^{1/2} B^{1/2}然后利用Tr(AB) Tr(B^{1/2} A B^{1/2}) ≥ 0被作用后的矩阵仍然是 PSD迹自然非负。5. PSD不是摆设统计、优化、机器学习里它都在场5.1 协方差矩阵与高维统计统计学里最常见的 PSD 矩阵就是样本协方差矩阵。假设数据矩阵 X 已经中心化维度是 n×p那么样本协方差S X^T X / (n-1)天然就是半正定矩阵因为任意向量 v 都有v^T S v ||X v||^2 / (n-1) ≥ 0那句任意向量的二次型等于某个范数的平方直接保证了非负性。当样本量 n 大于特征数 p 的时候X 通常列满秩S 是正定的。但在高维统计里经常出现 p 比 n 还大的情况此时 X 一定是秩亏的S 必然奇异。这意味着数据在某些方向上完全没有波动协方差矩阵有零特征值。很多初学者遇到numpy.linalg.LinAlgError: singular matrix时以为是 bug其实这是高维数据的数学本质不是程序写错了。常规修法是加一个小的对角扰动比如用S λI也就是统计里的岭正则化。也可以使用 Ledoit-Wolf 收缩估计把协方差矩阵往单位阵方向收缩。我自己在基因表达数据、用户行为矩阵这类场景里几乎每次都要处理这种问题。5.2 马氏距离与核方法马氏距离是协方差矩阵 PSD 性质的直接受益者。它用协方差矩阵的逆来重新定义距离可以理解成在数据真实的分布形状下测量距离。数据如果在一个方向上方差很小那么沿这个方向的微小偏差都应该被视为很大的距离变化。这在异常检测、聚类和图像检索里都非常实用。核方法则是把 PSD 用在了 Gram 矩阵上。给定一组样本和核函数k(x_i, x_j)Gram 矩阵K_ij k(x_i, x_j)半正定是 Mercer 条件的要求。只有 Gram 矩阵是 PSDSVM、高斯过程回归等算法里的二次规划才有解后验分布才合法。高斯核k(x, y) exp(-γ ||x-y||^2)在理论上是正定核但在数值计算里有个老毛病当数据点彼此距离极近时Gram 矩阵会接近一个全 1 矩阵特征值出现一个大的和一堆接近零的数值上很容易变成半正定但不严格正定。高斯过程社区管这个叫jitter解决办法是在对角线上加一个很小的数比如1e-8 * I。5.3 Hessian 矩阵与凸优化在优化理论中函数 f 的可微性和凸性之间有一个非常经典的等价关系f 是凸函数 ⟺ 对任意自变量 xHessian 矩阵 ∇²f(x) ⪰ 0如果 Hessian 矩阵的所有特征值都存在一个正的下界那么函数是强凸的。强凸的性质比凸更强它保证目标函数有一个唯一的极小点而且梯度下降能收敛得很快。以最小二乘问题为例目标函数是f(x) 1/2 * ||Ax - b||^2它的梯度是A^T A x - A^T bHessian 是A^T A。由于A^T A天然是 PSD所以最小二乘问题无条件是一个凸问题。这个事实解释了为什么我们从来不用担心最小二乘会收敛到局部极小点——它的所有局部极小点都是全局极小点。但是当 A 列秩不满时A^T A奇异最优解不唯一最小二乘的解有无穷多个。加上 L2 正则之后目标 Hessian 变成A^T A λI因为λI是正定的所以整体变成正定解就唯一化了。这就是岭回归能从数学上强制解唯一的原因。5.4 半正定规划SDP的简要面貌到了优化领域再往前走一步就遇到了半正定规划其中的变量本身就是一个 PSD 矩阵。这类问题的典型形式是最小化Tr(CX)条件包括Tr(A_i X) b_i和X ⪰ 0注意这里的未知量 X 是一个矩阵不是向量而约束X ⪰ 0是一个矩阵不等式。这个条件看似抽象其实等价于要求 X 的所有特征值非负。SDP 之所以重要是因为它把线性规划从非负的数轴锥推广到了非负的矩阵锥从而可以优雅地处理矩阵结构问题。组合优化里的 MaxCut 松弛、控制系统里的 LMI 分析、量子信息里的密度矩阵约束最后都可以转化成 SDP。如果你以后读优化方向的论文遇到把问题松弛成 SDP不必恐慌它背后的核心约束就是这一章讲的东西找到一个 PSD 矩阵同时满足若干线性等式。6. 分块矩阵视角Schur补与半正定的等价关系6.1 分块 PSD 判定Schur 补实际项目中遇到的大型 PSD 矩阵几乎都是分块结构的。协方差矩阵天然按变量分组分块高斯过程里的核矩阵按数据批次分块控制论里的李雅普诺夫矩阵也是分块形式。所以掌握分块视角非常实用。设一个对称分块矩阵A [[A11, A12], [A12^T, A22]]其中 A11 是 k×k 的可逆子块。如果 A 是 PSD那么 A11 和 A22 必然都是 PSD。这个可以直接验证取x (v, 0)二次型就退化成v^T A11 v非负性自动保留。但如果 A11 正定A 整体是否 PSD取决于一个关键对象——Schur 补S A22 - A12^T A11^{-1} A12结论是当 A11 正定时A 半正定当且仅当 S 半正定。证明思路是做一个合同变换用一个下三角消元矩阵 L 对 A 做变换L [[I, 0], [-A12^T A11^{-1}, I]]可以得到L A L^T [[A11, 0], [0, S]]因为合同变换不改变矩阵的正负惯性指数所以 A 的 PSD 性质完全由对角块 A11 和 S 决定。这个等价关系非常干净实际用起来也很顺手。6.2 分块求逆公式中的老熟人分块求逆公式大家可能都见过但未必注意到 Schur 补在里面反复出现。对于一个分块矩阵A^{-1} [[A11^{-1} A11^{-1} A12 S^{-1} A12^T A11^{-1}, -A11^{-1} A12 S^{-1}], [-S^{-1} A12^T A11^{-1}, S^{-1}]]右下角正好是 Schur 补的逆。这个结构在多元高斯条件分布中大有用处已知部分变量后剩余变量的条件协方差矩阵恰好就是 Schur 补。通俗地说Schur 补量化了把 A11 的信息消掉之后A22 还剩下多少自己的波动。如果 S 变成零矩阵说明 A22 的信息完全被 A11 解释了这对应数据里的确定性关系。顺手提一个热搜词里反复出现的分块矩阵的 n 次方。如果分块矩阵是对角的幂次直接对每个块分别求即可。一般分块矩阵求幂没有简单普适公式但可以对角化后分别求幂值。真正常用的场景是分块对角矩阵加上低秩扰动用 Woodbury 公式来做迭代更新。6.3 Sylvester 定律乘积特征值的非零部分分块和乘积常常绑定在一起这里有一个少为人知但极其有用的结论对任意两个方阵 A 和 BAB 和 BA 的非零特征值完全相同这就是 Sylvester 定律。如果 A、B 都是 PSD 矩阵AB 虽然不一定对称但它的特征值全是实数且非负。原因可以构造AB 与 A^{1/2} B A^{1/2}有相同的非零特征值而后者是 PSD 矩阵特征值自然全非负。这条性质在控制论判稳、卡尔曼滤波的可观测性分析以及很多统计检验中都有应用。如果你以后遇到某个算法的稳定性分析中出现AB 的特征值必须落在右半平面之类的说法它背后的逻辑很可能就是这条定律。7. 数值实践中的PSD陷阱与应对策略7.1 那些幽灵般的负特征值从哪来工程里最常见的症状是理论上应该 PSD 的矩阵一算特征值却出现几个小的负数。我根据经验把来源分成三类。第一类是浮点舍入误差。对称矩阵经过复杂的浮点运算后本来精确为零的特征值可能变成-1e-15这是完全正常的不需要过度处理。第二类是秩亏。样本量不足、特征重复、中心化后数据之间存在线性关系都会让矩阵产生精确的零特征值浮点运算下表现为微负值。第三类是超参数选择不当。比如高斯核的带宽取得过小Gram 矩阵趋近于全 1 矩阵数值上出现大量接近零的特征值再加上舍入误差就会报not positive definite。判断一个小负数到底是数值误差还是真正的负特征值我的经验阈值是负值绝对量级小于n * eps * λ_max就可以当作零来处理如果负值达到1e-3这个量级那一定是数据或模型结构出了问题加正则项只是掩耳盗铃。7.2 处理前一定要先对称化很多从计算流程里拿出来的矩阵细微处已经不对称了。比如两个矩阵乘积M^T N如果不相等再做一些浮点运算结果会和严格预期差出1e-16左右。此时不管用eigvalsh还是cholesky行为都不可预测因为eigvalsh默认只读下三角部分矩阵不对称时它会忽略掉上三角的信息。所以我的处理习惯是任何矩阵在做 PSD 判断之前先执行A_sym (A A.T) / 2这不仅让后续分解稳定也避免了因为A和A.T细微不一致而导致的诡异行为。7.3 三个常用修复方案如果确认矩阵需要修复我按优先级推荐三种方案。第一种是对角扰动也就是 jitterA_fixed A_sym 1e-8 * np.eye(n)这是最简单、最常用的方案适用于核矩阵、协方差矩阵这类加一点噪声不影响结果的场景。高斯过程里默认这么做原因就是高斯核在数据点靠得很近时协方差矩阵数值上几乎奇异。第二种是特征值截断。把负特征值强行置零再重构矩阵eigvals, eigvecs np.linalg.eigh(A_sym) eigvals_fixed np.maximum(eigvals, 0) A_fixed eigvecs np.diag(eigvals_fixed) eigvecs.T这种做法的代价是改变了原始矩阵的谱适合你不在乎特征值结构只在乎矩阵合法且距离原矩阵尽量近的场景。第三种是更精细的收缩估计类似 Ledoit-Wolf。它不是对每个矩阵都加同一个常数而是按特征值大小做自适应收缩A_shrunk (1 - α) A α * (Tr(A) / n) I其中 α 可以通过交叉验证选择。这个方法在协方差估计中表现远好于固定对角扰动因为它在保留信号和消除噪声之间自动找平衡。7.4 奇异半正定矩阵求解线性系统当 A 只是半正定但奇异时np.linalg.solve会报错但问题本身往往有解只是解不唯一。最常用的替代方案是伪逆x np.linalg.pinv(A) b它给出的是最小范数意义下的解也就是所有可行解中||x||最小的那个。在统计中这个解恰好对应线性回归里当X^T X奇异时的岭估计极限形式。另一个选择是np.linalg.lstsq(A, b)它在最小二乘意义下自动处理秩亏问题。如果矩阵规模很大不能用直接法可以考虑迭代法。共轭梯度法要求矩阵对称正定半正定加上对角扰动之后就可以安全使用了。注意扰动不要太大否则解会偏离原问题太多。7.5 一次真实的排查经历回到开头那个高斯过程回归的例子。我当时的排查链路是这样的先打印最小特征值发现是-3e-15。检查数据确认没有缺失值量纲差异不大。检查核函数带宽发现很小导致 Gram 矩阵几乎全 1。于是判定这是数值秩亏而不是真实负特征值。在核矩阵对角线上加了1e-8的 jitter问题解决。但如果当时看到的最小特征值是-0.1我会回头检查原始数据是否存在重复列、是否中心化之后仍有线性相关变量、核带宽是否小到令两点间相似度几乎为 1。加 jitter 能救数值问题但救不了真负特征值因为后者说明你的模型或数据在某个方向上能量为负这在物理上是不应该存在的。最后分享一个我在工程里的习惯任何矩阵到我手里先做三步检查——对称化、看最小特征值、尝试 Cholesky。如果看到一个-1e-15量级的负特征值我直接当它是零不需要大动干戈但如果负特征值到了-1e-3甚至更大先别急着加 jitter回头检查数据生成流程十有八九是秩亏或数据没对齐。理解半正定矩阵的真正价值不是让你背一堆定理而是让你在碰到matrix is not positive definite时能准确判断这到底是数学问题、数值问题还是数据问题并且知道下一步该怎么做。