AI 工程地基 17｜线性系统：七种解法，一条警戒线

AI 工程地基 17｜线性系统：七种解法，一条警戒线

「AI 工程地基」系列第十八篇。把 GitHub 上 58k star 的开源课 AI Engineering from Scratch（MIT 协议）一课一课啃完、核完、讲完。今天第 17 课：线性系统，Build 型，120 分钟，前置是线性代数前三课。

你每次 fit，都在做一件两千年前的事

每次写下 model.fit(X, y)，底下都在干同一件事：解一个线性方程组 Ax = b。A 是你的数据矩阵，b 是你的目标值，x 是要求的权重。线性回归是这个方程，最小二乘是这个方程，神经网络每层的 y = Wx + b 是这个方程的一边，高斯过程要对核矩阵做分解，马氏距离要对协方差求逆——全是这一个方程的变体。

这课开场白说，Ax = b 是数学里最古老、却还在驱动你神经网络的问题。多老？中国《九章算术》的方程章就在系统地做消元，高斯还要再等一千八百年才出生。

这一课把解这个方程的所有主流方法从零写一遍，最后收在七个方法的一张表上。看完你会发现它们其实在回答两个问题：怎么解，以及解出来的东西能不能信。

先换个读法：方程是配方

Ax = b 有两种读法。按行读，每个方程是一张超平面，解是所有平面的交点。按列读，A 的每一列是一种原料，解方程是找配方：几份原料一、几份原料二，能调出 b 这个味道。

两行两列的小例子：A 的第一列是 [2, 1]，第二列是 [1, -1]，要调出 b = [5, 1]。答案恰好是 2 份第一列加 1 份第二列。列的读法更根本——b 落在列能调出的范围里，方程有解；不落在里面，就找范围内离它最近的那个点。这个「最近」，就是后面最小二乘的全部。

消元法，和它容易翻车的地方

高斯消元你大概见过：一列一列往下消，把矩阵变成上三角，再从最后一个方程往回代。一百乘一百的矩阵，大约十亿次浮点运算，今天的电脑眨眼就完。

真正的坑在主元。消元要除以对角线上的元素，这个数叫主元。它是零，除以零直接崩。它很小，你就在放大舍入误差。

这课给了一个经典对比：左边不选主元直接消，乘数是 1000；右边先换行再消，乘数是 0.001。左边看着凶险。但有意思的是，配图两边最后都算出了 x = 1.000，各自标着「正确」。想演示的事故，在图里没有发生——在 float64 双精度下，这道题的误差藏在第八位有效数字之后，图上那点位数根本看不出来。

我把同样的计算换成 float32 单精度重跑：不选主元的 x₁ = 0.99993，误差 7.3 × 10⁻⁵。float32 满打满算七位有效数字，这一下丢了两位半。坏是真的，只是这张图没把它演出来。想亲眼看见，得自己动手把精度砍下来。

所以部分主元法（每次挑该列绝对值最大的行换上来）不是仪式，是保险。它把乘数压到最小，把误差放大压到最轻。

分解一次，用一千次

消元有个浪费：换一个 b，从头再来，又是 O(n³)。LU 分解把消元的中间产物存下来：A = LU，L 存各行被消时的乘数，U 是消完的上三角。之后每来一个新 b，先解 Ly = b（往下代），再解 Ux = y（往回代），两趟三角求解都是 O(n²)，便宜得多。

课里说，对同一个 A 解 1000 个系统，LU 在总工作量上「省下 1000/3 倍」。这个数我推不出来。实际省多少取决于矩阵大小：n = 100 时省 16 倍，n = 1000 时 143 倍，n = 10000 时 625 倍。矩阵越大，分解的固定成本摊得越薄，省的倍数越接近解的个数。方向没错：解得越多越划算。但「1000/3」这个具体数字，没有来历。

LU 的近亲还有两个。QR 分解把 A 拆成「正交矩阵 × 上三角」，正交矩阵乘上去不改变长度和角度，数值上格外稳，是最小二乘的标配。Cholesky 分解只伺候一种矩阵——对称正定（协方差矩阵、核矩阵、凸函数的 Hessian 全是这类），换来一半的运算量和一半的存储，高斯过程和 ridge 回归天天用它。

Cholesky「比 LU 快一倍」，这个说法在运算量层面是真的：n³/3 对 2n³/3。但课后练习 5 让你拿课程自带的代码验证这 2 倍，你会验不出来。我实测了：n = 50 时 Cholesky 快 1.28 倍，n = 200 时反而慢了四分之一。手写 Python 的循环开销和内存读写，把运算量省下的钱全吃回去了。同样的对比换到 NumPy 背后的专家级实现（LAPACK）上：2.35 倍、2.47 倍，「快一倍」这才兑现。一个只在专家实现里成立的数字，让初学者用自己的代码去验证，这道练习题出得欠考虑。想看真实的 2 倍，直接计 LAPACK 的时。

没有精确解的时候

机器学习里 A 通常瘦高——一百个数据点、三个参数，方程比未知数多。没有 x 能让所有方程同时成立，只能找最接近的：最小化 ‖Ax − b‖²，残差的平方和。

最小化点满足一个漂亮的方程：AᵀA x = Aᵀb，叫正规方程。对 ‖Ax − b‖² 求梯度令它为零，三行推导就到。而它就是线性回归本身——你的数据矩阵代进去，解出来的就是回归权重。

这课的数值例子是四个点 (1,3)、(2,5)、(3,6)、(4,8) 拟合直线。算出 AᵀA = [[4, 10], [10, 30]]，Aᵀb = [22, 63]，书上给的解是 x = [1.5, 1.7]，截距 1.5、斜率 1.7。

我把这四个方程亲手解了一遍。斜率不是 1.7，是 1.6。克莱姆法则算这一项，分子 32、分母 20，32 除以 20 等于 1.6。代回去检验：斜率 1.6 的残差平方和是 0.2，换成 1.7 是 0.5；正规方程解的标志是残差和两个列向量内积为零，1.6 过检，1.7 不过。数据没问题，推导没问题，最后一步印错了。

这种错误不致命，但值得亲手算一遍。机器学习里危险的从来不是没有答案，是拿着次优解当最优解用。

警戒线：条件数

解算得快不算完，还得问一句：这个解能信吗？衡量的是条件数 κ，矩阵最大奇异值除以最小奇异值。经验法则直白：κ ≈ 10^k，你就丢大约 k 位精度。float64 有 16 位精度，所以 κ 到 10¹⁶，答案等于全丢光，矩阵在浮点世界里就是奇异的。

教科书都说解最小二乘用 QR 比走正规方程稳，原因通常一句带过。机制其实就一行字：正规方程把条件数平方。κ(AᵀA) 恰好等于 κ(A) 的平方。我造了个 50 行 3 列、其中两列近乎重复的数据实测：κ(X) = 37.0，κ(XᵀX) = 1371——正好 37 的平方。条件数一平方，丢的位数跟着翻倍：X 本身丢一两位，进了正规方程就丢三四位。QR 直接作用在 X 上，只承受那个 37。这就是「稳」的全部来历。

还有个补丁式的救法：伪逆，A⁺ = VΣ⁺Uᵀ，经 SVD 算出来。系统有唯一解它给唯一解，无解它给最小二乘解，无穷多解它给范数最小的那个。任何形状、任何秩的矩阵都伺候，NumPy 的 lstsq 和 pinv 内部都走它。全能，但慢——所以它是兜底，不是首选。

加 0.01，把三亿压回零点一

课后练习 3 是全课最好的一道题，数字我从头跑了一遍。

造一个 50 行 2 列的矩阵，第二列几乎是第一列的复制品，只差 10⁻¹⁰ 量级的噪声。此时 XᵀX 的条件数是 2.7 × 10¹⁶——不偏不倚踩在 float64 的死刑线上，16 位精度理论上全灭。

不加正则直接解：权重是 −3.39 亿和 +3.39 亿。发生了什么？两列几乎一样，模型学会了「用这列减那列」，把 10⁻¹⁰ 的微小差异放大三亿倍去拟合噪声，两个巨数再互相抵消。数学上完全成立，实际上是垃圾。

往矩阵加 0.01 × I 再解——这就是 ridge 回归：条件数从 2.7 × 10¹⁶ 落回 7778，权重变成 −0.082 和 −0.082。两列几乎一样，就给几乎一样的权重，说得通。关键看代价：残差从 42.388 涨到 42.428，涨了千分之一。

正则化买的从来不是更好的拟合，是一个值得相信的解。

迭代法：矩阵大到分解不动时

上百万未知数的稀疏矩阵（比如大规模 PDE 或超大核矩阵），LU 和 Cholesky 都太贵。共轭梯度（CG）换思路：不分解，从零向量出发，每轮沿着一个精心挑选的方向挪一步，残差不够小就再来。

课里说 CG 最多 n 步找到精确解，括号里注明 in exact arithmetic——精确算术下。这个括号是活的。我用 100 × 100 的对称正定矩阵实测：条件数约 5 时，21 步收敛，远低于 100；条件数 380 时，138 步——超过了 100。浮点误差把「最多 n 步」的保证戳破了。条件数越好收敛越快，所以正则化在这又立了一功。这个系列读到这里，见过太多不写的括号了；这次写上了，还写得这么必要。

收摊之前，对一遍数

课尾的实战代码值得一提。同一份数据，课程自写的路线（手搓消元解正规方程、手搓 Cholesky 解 ridge）跟 NumPy 的 SVD 路线、sklearn 的路线对到同一个数上：普通最小二乘两边差 8.9 × 10⁻¹⁶，ridge 差 2.2 × 10⁻¹⁶，都是机器精度。两条完全独立的算法路径落到同一个数，这是这个系列里少见的一手好验证。

课后练习 1 的系统 [[1,2,3],[4,5,6],[7,8,10]] x = [6, 15, 27]，用课里自写的消元、自写的 LU、NumPy 官方求解器分别跑，三个都是 [3, −3, 3]。你可以再补一手：代回原方程逐行验，6、15、27 一个不差。第 18 课是凸优化，从「解方程」走向「带约束地解方程」。下一课见。