ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

Moore-Penrose广义逆:从定义到SVD计算与工程实践

Moore-Penrose广义逆:从定义到SVD计算与工程实践 矩阵这种东西工作时最烦的就是遇到“不可逆”。刚开始干活那几年我一看到pinv这个函数就觉得它是个“凑合用的逆”后来被各种病态问题折磨过几轮才真正意识到 Moore-Penrose 广义逆也常称为加号广义逆、伪逆记作 $A^{}$才是线性代数里最被低估的工具之一。这篇文章我不打算念念教材而是把 $A^{}$ 的定义、四条核心条件、唯一性直觉、怎么算、以及实际项目中容易踩的坑完整拆一遍。不管你是做数据分析、控制、图像、机器学习还是单纯被矩阵论折磨理解 $A^{}$ 都能让你在面对“不可逆矩阵”时多一张底牌。它不是万能的但在绝大多数需要“最小二乘解”或“最小范数解”的场合它都是那个正确的默认选择。1. 背景为什么我们反复需要“一个并不存在的逆”先从一个几乎人人都遇到过的场景说起。解线性方程组 $Ax b$ 的时候如果 $A$ 是方阵且可逆直接写 $x A^{-1}b$ 就结束了。可是现实里哪有这么便宜的事要么 $A$ 是 $m\times n$、$m\ne n$ 的长条矩阵要么 $A$ 方方正正但行列式为零——缺秩。很多教材在这里给出的答案是“此路不通”但工程问题不会因为矩阵“不可逆”就停止向你要结果。你必须有一个“尽量像逆、但又在退化方向上让步”的东西这正是 Moore-Penrose 广义逆存在的理由。1.1 三种让你写不出“$A^{-1}$”的情况先给矩阵“把把脉”。设 $A$ 是 $m\times n$ 矩阵左侧代表方程个数右侧代表未知量个数。$mn$ 且 $\det A\ne 0$天下太平唯一解 $xA^{-1}b$。$mn$超定方程组方程的数量比未知数多。系统往往没有严格解只能找一个让残差 $|Ax-b|$ 尽量小的“妥协解”。$mn$欠定方程组未知数的数量比方程多。解有无穷多个精确的“那一个解”根本不存在我们还要额外加一条标准比如找范数最小的那个。第二和第三种情况都有一个共同点它们都不是“真正的可逆系统”但你仍然希望用一个统一的框架去处理。Moore-Penrose 广义逆恰好就长在这个框架上。它不要求 $A$ 是方阵不要求满秩任何实矩阵或复矩阵都能算出 $A^{}$。这个“什么矩阵都定义得出来”的性质是它成为诸多算法底层标配的根本原因。1.2 核心思想在“有用的方向上做逆在退化方向上清零”我在给新人讲 $A^{}$ 时最喜欢用投影来比喻。方块矩阵可逆时$A^{-1}$ 把“像空间”的每一个向量毫无损失地映射回“定义域”。缺秩或非方阵时$A$ 不再是一对一的满射这时候强行求逆就会在没定义清楚的方向上“瞎猜”。$A^{}$ 的做法非常朴素把整个空间拆成两部分——一部分是“$A$ 能好好作用的方向”另一部分是“$A$ 直接压扁掉的方向”。在能好好作用的方向上$A^{}$ 老老实实做逆在被压扁的方向上$A^{}$ 干脆置零。用数学家的话说$A^{}$ 在 $A$ 的值域 ${\rm Im}(A)$ 上与 $A$ 互为逆在 ${\rm Im}(A)$ 的正交补上为零。这一条简单的准则就是后面所有性质的源头。顺便说一句历史。这个概念最早可以追溯到 Moore 在 1920 年前后的工作后来 Penrose 在 1955 年独立给出了我们现在常用的四条公理式定义。正因为 Penrose 把这四条条件写得这么干净后世才把它命名为 Moore-Penrose 广义逆。你写代码时调用pinv或pinv的实现背后基本都能溯源到 Penrose 的这套框架。2. 四条件定义Moore-Penrose 逆的“身份证”2.1 Penrose 四条等式给定实矩阵 $A$复矩阵只要把转置换成共轭转置 $A^{*}$其余逻辑完全一样满足下面四个条件的矩阵 $X$ 称为 $A$ 的 Moore-Penrose 广义逆$AXA A$$XAX X$$(AX)^{T} AX$$(XA)^{T} XA$。这四条叫 Penrose 条件。满足这四个条件的矩阵 $X$ 存在且唯一记作 $A^{}$。注意如果只满足其中部分条件得到的矩阵就不唯一统称为广义逆比如满足 (1) 的叫内逆或“${1}$-逆”满足 (1)(2) 的叫自反广义逆。我们在实际中用到的几乎都是完整满足四条的那个唯一对象。提示这里默认矩阵是实的所以转置用 $(\cdot)^{T}$如果是复数矩阵所有转置都要换成共轭转置 $(\cdot)^{H}$ 或 $(\cdot)^{*}$后面不再重复。2.2 四条条件到底在说什么刚接触这四条的人往往很晕为什么偏偏是这四个少一个行不行我逐条说。条件 (1) $AXA A$ 是“逆”这个身份的最低要求。把它改写成 $A X A A$意思是在某些方向上 $X$“还原”了 $A$ 的行为若 $y Ax$那么 $A X y A X A x A x y$也就是说 $AX$ 作用在 $A$ 的值域上时等于恒等。条件 (2) $XAX X$ 是对称地约束 $X$ 本身保证 $XA$ 在 $X$ 的值域上也等于恒等。这两个条件合起来已经让 $X$ 成为一个“代数上说得过去”的广义逆但此时 $X$ 依然不唯一。真正把唯一性锁死的是条件 (3) 和 (4)它们要求 $AX$ 与 $XA$ 都是对称矩阵。于是 $AX$ 和 $XA$ 都成了正交投影矩阵而不是乱七八糟的斜投影。几何上$AA^{}$ 是把空间投影到 ${\rm Im}(A)$ 的正交投影$A^{}A$ 是把空间投影到行空间 ${\rm Im}(A^{T})$ 的正交投影。正交投影意味着天然带有“垂线段最短”的优化属性这正是后面解决最小二乘问题的关键。2.3 存在性与唯一性的块状证明回到为什么非用四条不可。我用一个技巧性不强但非常直接的证明来说明。假设 $A$ 的奇异值分解为 $AU\Sigma V^{T}$其中 $\Sigma$ 是一个“对角块”矩阵$$ \Sigma\begin{pmatrix} \Sigma_r 0\ 0 0 \end{pmatrix}_{m\times n},\qquad \Sigma_r\mathrm{diag}(\sigma_1,\dots,\sigma_r),\ \sigma_i0. $$令 $XV^{T} X U$把 $X$ 按同样方式分块成 $\begin{pmatrix}X_{11} X_{12}\ X_{21} X_{22}\end{pmatrix}$。此时四个条件在正交基下等号不变我们把它们全部翻译到分块上。先看条件 (1) $\Sigma X\Sigma\Sigma$左边只剩左上角一块 $\Sigma_r X_{11}\Sigma_r$于是立刻得到 $X_{11}\Sigma_r^{-1}$条件 (3) 要求 $\Sigma X$ 对称。$\Sigma X\begin{pmatrix}\Sigma_r X_{11} \Sigma_r X_{12}\ 0 0\end{pmatrix}$对称性让右上角必须等于左下角的转置也就是 $\Sigma_rX_{12}0$所以 $X_{12}0$。条件 (4) 类似地推出 $X_{21}0$。条件 (2) $X\Sigma XX$ 再补上最后一块 $X_{22}$最后只剩 $X_{22}0$。也就是说在奇异值分解的坐标下满足条件的矩阵必须长成这样$$ X\begin{pmatrix} \Sigma_r^{-1} 0\ 0 0 \end{pmatrix}_{n\times m}. $$于是 $XV\begin{pmatrix}\Sigma_r^{-1} 0\0 0\end{pmatrix}U^{T}$ 唯一确定。同时我们反向把 $V\Sigma^{}U^{T}$ 代入四个条件也都能满足所以存在性和唯一性同时证毕。这个证明看起来是在走 SVD 的捷径但它很能说明一个道理$A^{}$ 的所有性质本质上都被奇异值分解控制住了。3. 性质全景从“正统逆”到“最优化解”3.1 与普通逆、左逆、右逆的关系如果 $A$ 本身可逆$A^{}$ 就是 $A^{-1}$。把 $A^{-1}$ 代入四个条件一个不漏地成立由于唯一性$A^{}A^{-1}$。这说明 $A^{}$ 是普通逆的推广普通逆只是它的特例。当矩阵不满秩但“有一侧信息是完整的”时$A^{}$ 还有更简洁的闭式$A$ 列满秩$m\ge n$ 且 ${\rm rank}(A)n$时$A^{T}A$ 可逆$A^{}(A^{T}A)^{-1}A^{T}$也就是通常说的左逆。它在超定最小二乘中非常常见。$A$ 行满秩$m\le n$ 且 ${\rm rank}(A)m$时$AA^{T}$ 可逆$A^{}A^{T}(AA^{T})^{-1}$也就是右逆常用于欠定问题。这两个公式不要死记记住它们来自“在对应方向上可逆”即可。如果左右都不满秩就只能老老实实走 SVD 或迭代路线。3.2 值域、零空间与投影关系$A^{}$ 对值域和零空间的刻画几乎是“教科书级”的$$ {\rm Im}(A^{}){\rm Im}(A^{T}),\qquad {\rm Ker}(A^{}){\rm Ker}(A^{T}). $$更有用的是两个投影恒等式$$ AA^{}P_{{\rm Im}(A)},\qquad A^{}AP_{{\rm Im}(A^{T})}. $$左边是到列空间的正交投影右边是到行空间的正交投影。这两个等式意味着什么意味着对任意 $b$$AA^{}b$ 就是 $b$ 在列空间上的正交投影。判断线性方程组 $Axb$ 是否有解可以直接看 $AA^{}b$ 是否等于 $b$等于则有解不等于则没有严格解只能找最小二乘意义下的近似解。这种“用投影看解的存在性”的视角在刚接触时觉得抽象但用熟之后能省下大量纠结。3.3 为什么 $A^{}b$ 是最小二乘和最小范数的“双料冠军”这是 $A^{}$ 最吸引人的性质也是它活跃在回归、控制、机器学习中的原因。考虑任意 $b$定义 $x^{}A^{}b$那么$x^{}$ 是最小二乘解$|Ax^{}-b|_2 \min_x |Ax-b|_2$在所有最小二乘解里$x^{}$ 的欧氏范数最小。原因其实很直白。把方程组 $Axb$ 的解集写成 $x_{\rm particular} {\rm Ker}(A)$其中 $x_{\rm particular}$ 是任意一个特解。由于 $A^{}b\in {\rm Im}(A^{}){\rm Im}(A^{T})$而 ${\rm Ker}(A)$ 与 ${\rm Im}(A^{T})$ 正交所以 $x^{}$ 自动与零空间正交再加上它确实落在正确的仿射空间里它就成了距离原点最近的解。这类解还被冠以“最小范数最小二乘解”minimum-norm least-squares solution的名字。很多书把这个性质写在很后面但我建议你最先记住它因为后面一堆应用都是围绕它展开的。4. 怎么算从 SVD 到手写实现4.1 通用的 SVD 配方$A^{}$ 的通用算法几乎总是走奇异值分解。设 $A$ 的瘦 SVD或者完整 SVD为$$ AU\Sigma V^{T}, $$其中 $U$、$V$ 是正交矩阵$\Sigma$ 是 $m\times n$ 的“对角矩阵”对角线上的 $\sigma_i$ 按从大到小排列。那么$$ A^{}V\Sigma^{}U^{T}, $$这里的 $\Sigma^{}$ 是一个 $n\times m$ 矩阵做法只有一句话把 $\Sigma$ 转置然后把每个非零奇异值换成倒数零保持为零。写成块状就是$$ \Sigma\begin{pmatrix} \Sigma_r 0\ 00 \end{pmatrix} \Longrightarrow \Sigma^{}\begin{pmatrix} \Sigma_r^{-1} 0\ 00 \end{pmatrix}. $$为什么非零的要取倒数在那些奇异值非零的方向上$A$ 是可逆的可是奇异值等于零的方向被压缩没了你没法恢复那些信息只能给零。这就呼应了第一节的几何思想有用的方向做逆退化方向清零。4.2 实际写代码时的三个细节我自己用 Python 给这类算法做实现时有几个细节值得拿出来说。第一不要为满秩情况重复设计“自己的伪逆”。如果你的矩阵满足列满秩直接解 $(A^{T}A \lambda I)x A^{T}b$ 或调用np.linalg.solve不要算显式的 $(A^{T}A)^{-1}$。显式求逆会引入不必要的数值误差。只有当你明确需要得到 $A^{}$ 这个矩阵本身才用np.linalg.pinv(A)。第二np.linalg.pinv默认会做奇异值截断。它内部把小于max(shape) * eps * max_σ的奇异值直接当零处理。这个默认阈值在新手期够用但在病态问题里建议自己传入rtol。遇到矩阵条件数特别大但并非真正缺秩时可能你想要保留一些特别小的奇异值也可能你想截掉更多——这取决于你面对的应用别一直沿用默认值。第三如果你基于 SVD 手动实现 $A^{}$绕不开的一步是判断哪些奇异值“足够大”。我自己习惯用比例而非绝对值算σ_i / σ_max小于某个阈值比如 $10^{-12}$ 或 $10^{-14}$就归零。绝对阈值在不同量纲的数据下会害死人。4.3 欠定和超定问题的快速公式如果不想调用通用伪逆两个特殊场景有专门公式。对于列满秩的超定问题$\hat\beta(X^{T}X)^{-1}X^{T}Y$ 是最小二乘解的核心在实际回归中我们通常不显式求逆而是用 QR 分解解正规方程。对于行满秩的欠定问题$A^{}A^{T}(AA^{T})^{-1}$它在机器人逆运动学里非常常见因为机器人雅可比矩阵 $J$ 是 $3\times n$ 或 $6\times n$想要最小范数关节速度用的就是 $J^{T}(JJ^{T})^{-1}$ 这类右逆。到这里你会发现这些特殊公式都只是 $A^{}$ 在满秩退化情况下的特例。真正到了缺秩场景公式没法用了通用 SVD 才顶上。5. 应用场景从最小二乘回归到逆运动学5.1 线性回归的“标准答案”其实是伪逆先看最常规的例子。线性回归写成矩阵形式是 $YX\beta$其中 $Y$ 是 $n$ 个观测值$X$ 是 $n\times p$ 的设计矩阵。当 $p\le n$ 且 $X$ 满列秩时教科书给你 $\hat\beta(X^{T}X)^{-1}X^{T}Y$。而这个表达式正是 $X^{}Y$。有些教材把伪逆藏在最小二乘解背后你完全没有意识到自己其实在用 $A^{}$。更关键的是当 $X$ 缺秩比如特征之间有完美共线性时$(X^{T}X)^{-1}$ 已经不存在但 $X^{}Y$ 依然稳定地给出一个解——这个解不仅是残差最小的解还是所有最小二乘解里范数最小的那个。所以很多现代机器学习库在做线性模型时底层用的其实是伪逆而不是单纯记正规方程。5.2 欠定系统越少的信息越需要最小范数欠定问题 $mn$ 在现实中比想象中常见。比如机器人手臂要到达目标位置关节数经常远大于末端自由度雅可比矩阵是“宽”的再比如多天线系统中的波束成形测量维度小于待估参数。这些场景里解不唯一额外要求“动作尽量小”或“能量尽量小”是很自然的工程诉求。用 $A^{}$ 就天然得到范数最小的那个解所以它成了这类问题的默认工具。我做控制系统时经常遇到类似的取舍。给定一个 $Axb$如果矩阵欠定直接解一个带约束的最小化问题当然也行但每次求解可能要多写几行优化代码当问题规模很大、需要在循环里反复调用时直接使用 $A^{}$ 的闭式解往往比内点法快一个数量级。代价是它给出的只是最小范数这一特定角度的解如果你需要的是“稀疏解”或“只有少数变量非零”那 $A^{}$ 就不是最优了。5.3 其他领域的典型身影$A^{}$ 在信号处理、图像去模糊、控制理论里都有经典身影。图像去模糊中的伪逆滤波本质上是把退化算子写成矩阵再用伪逆去反演系统辨识里多输入多输出系统的“最速下降”型控制律也会用到广义逆。在很多被叫做“Inverse Problem”的领域里第一步就是把算子离散化第二步就是求伪逆或带正则化的伪逆第三步才谈后续的约束条件。下面给一个极小的实现片段用来验证思想也方便你复制到本地跑import numpy as np A np.array([[1., 2.], [2., 4.], [3., 6.]]) Ap np.linalg.pinv(A) print(Ap) # 验证四个条件 X Ap print(np.allclose(A X A, A)) # 条件1 print(np.allclose(X A X, X)) # 条件2 print(np.allclose(A X, (A X).T)) # 条件3 print(np.allclose(X A, (X A).T)) # 条件4这段代码的输出里会有四个True说明np.linalg.pinv给出的确实是 Moore-Penrose 广义逆。自己手写 SVD 时最常出错的地方就是 $\Sigma^{}$ 的尺寸好好核对是 $m\times n$ 还是 $n\times m$。6. 一个完整的手算小案例6.1 问题与 SVD 路径我挑一个非常简单但能说明问题的矩阵让你看清步骤也方便你自己用手算一遍验证。设$$ A\begin{pmatrix} 1 2\ 2 4\ 3 6 \end{pmatrix}. $$这个矩阵第二列是第一列的 2 倍所以秩为 1既不是列满秩也不是行满秩普通左逆右逆公式全部失效只能用通用定义。先写成秩一分解$Au,v^{T}$其中 $u(1,2,3)^{T}$$v(1,2)^{T}$。这里注意 $u$ 是三维列向量$v$ 是二维列向量$uv^{T}$ 正好是一个 $3\times2$ 的秩一矩阵。$v$ 的范数是 $\sqrt{5}$$u$ 的范数是 $\sqrt{14}$$A$ 的奇异值 $\sigma_1 |u||v| \sqrt{70}$。于是$$ A^{}\frac{1}{\sigma_1}\frac{v}{|v|}\frac{u^{T}}{|u|} \frac{1}{\sqrt{70}}\frac{(1,2)}{\sqrt5}\frac{(1,2,3)}{\sqrt{14}} \frac{1}{70}\begin{pmatrix} 123\ 246 \end{pmatrix}. $$6.2 用 Penrose 条件验证拿到这个结果最稳妥的事就是代入四个条件验算。先看 $A^{}A$$$ A^{}A\frac{1}{70}\begin{pmatrix} 123\ 246 \end{pmatrix} \begin{pmatrix} 12\ 24\ 36 \end{pmatrix} \frac{1}{70}\begin{pmatrix} 1428\ 2856 \end{pmatrix} \begin{pmatrix} 1/52/5\ 2/54/5 \end{pmatrix}. $$它显然是对称的而且满足 $(A^{}A)^{2}A^{}A$所以是正交投影。再看 $AA^{}$$$ AA^{} \frac{1}{70}\begin{pmatrix} 12\ 24\ 36 \end{pmatrix} \begin{pmatrix} 123\ 246 \end{pmatrix} \frac{1}{70}\begin{pmatrix} 51015\ 102030\ 153045 \end{pmatrix}. $$也是一个对称幂等矩阵。条件 (1)(2) 也不难代入验证。四个条件成立所以这个 $A^{}$ 确实是唯一定义中的那个。6.3 用这个结果解一个最小二乘问题取 $b(1,1,1)^{T}$令$$ x^{}A^{}b\frac{1}{70}\begin{pmatrix} 123\ 246 \end{pmatrix}\begin{pmatrix} 1\1\1 \end{pmatrix} \begin{pmatrix} 6/70\12/70 \end{pmatrix} \begin{pmatrix} 3/35\6/35 \end{pmatrix}. $$算一下残差 $Ax^{}-b$。$Ax^{}\begin{pmatrix}3/35\cdot16/35\cdot2,\ 3/35\cdot26/35\cdot4,\ 3/35\cdot36/35\cdot6\end{pmatrix}\begin{pmatrix}15/35,\ 30/35,\ 45/35\end{pmatrix}\begin{pmatrix}3/7,\ 6/7,\ 9/7\end{pmatrix}$。残差为 $\begin{pmatrix}4/7,\ 1/7,\ -2/7\end{pmatrix}$它与列空间 $\mathrm{span}{(1,2,3)^{T}}$ 做内积$$ \frac{4}{7}\cdot1\frac{1}{7}\cdot2-\frac{2}{7}\cdot30. $$残差垂直于列空间这就是最小二乘解的正交性条件。同时如果往 $x^{}$ 上加上 ${\rm Ker}(A)$ 的任意非零向量例如 $(-2,1)^{T}$范数只会变大所以 $x^{}$ 确实是最小范数解。整个流程和理论完全对上了。7. 常见问题与避坑指南7.1 “为什么我的伪逆算出来巨大”这是我被问得最多的问题。多数时候不是实现错了而是矩阵太接近缺秩。当矩阵最小奇异值很小但不为零伪逆会把它放大成非常大的倒数同时只要奇异值发生极小扰动伪逆会发生剧烈变化。从数学上讲矩阵到伪逆的映射在缺秩边界上不连续这是数值上要不断小心的地方。对策很直白如果你的问题本质上病态用截断奇异值TSVD或岭正则化。岭估计 $(\lambda I A^{T}A)^{-1}A^{T}b$ 在 $\lambda$ 适当大时数值表现会比直接使用 $A^{}$ 稳得多。不要迷信“闭式解一定稳定”在病态问题上闭式解就摔得越狠。7.2 五个容易踩的坑把 $(AB)^{}B^{}A^{}$ 当成恒等式。这个等式并不恒成立只有当 $A$ 有正交列而 $B$ 有正交行等特殊条件下才成立。计算时别贪这条捷径。复数矩阵忘掉共轭。Penrose 条件里的转置在复数域必须换成共轭转置 $A^{*}$否则一定算出错。拿伪逆硬套稀疏解。$A^{}$ 给出的是最小范数意义上的解如果你想找“只有少数分量非零”的解请用 L1 正则化或匹配追踪。在列满秩问题里显式求 $(A^{T}A)^{-1}$。先求逆再乘向量数值误差通常比直接解线性方程组大得多。忽略了伪逆对统计推断的影响。如果设计矩阵缺秩伪逆给出的回归系数只是无穷多个最优解中的一个它的可解释性不一定好不要直接拿它解读系数。7.3 快速自检验证你的伪逆实现每当我实现完一个pinv功能都会用随机矩阵做一次批量化自检。流程是生成随机矩阵 $A$计算 $XA^{}$然后用四个条件逐一验证np.allclose。这一步能过滤掉绝大多数转置写反、实复域搞错、尺寸算错的基础 bug。值得养成习惯。提示随机生成时别只用满秩矩阵特意加上缺秩矩阵、长条矩阵、瘦高矩阵覆盖全类型。这样自检才有效。我在实际项目中还有个小偏好当矩阵不是特别庞大时我先用完全 SVD 做一遍再用瘦 SVD 做一遍对比结果。如果两种写法结果不一致八成是 $\Sigma$ 的尺寸或零行零列的裁剪出了问题。这类错误在纸面上看不出来一到真实数据就会变成莫名其妙的 NaN。最后再分享一个很实用的习惯遇到任何线性系统先别急着直接解花十秒钟看它的尺寸和秩。能分清楚“超定”“欠定”“方阵可逆”这三种情况再决定用伪逆、QR 还是直接求逆你会少走很多弯路。$A^{}$ 是线性代数工具箱里最通用的万能钥匙但用之前搞清楚锁的类型才不会把钥匙拧断。
返回列表