最小二乘法是回归分析中最基础、最常见的参数估计方法。它的目标可以概括为一句话:给定一组观测数据,找到一条曲线(在线性回归中通常是直线),使模型预测值与真实值之间的误差平方和达到最小。误差平方和写成公式就是 S = Σ(y_i - ŷ_i)²,其中 y_i 是第 i 个真实值,ŷ_i 是模型给出的预测值。之所以用平方而不是绝对误差,一个重要原因是平方函数连续可导,极值点可通过导数或矩阵求导直接求解;同时平方对大误差更敏感,惩罚力度随偏离程度平方增长,这让拟合结果在大多数场景下更合理。

如果接受误差平方和这一标准,最小二乘问题就变成了一个凸优化问题。对一元线性回归 y = ax + b,误差平方和 S 对参数 a、b 分别求偏导并令偏导数为零,可以解出唯一最优解。多元情况则更适合用矩阵描述:设 X 是设计矩阵,每一行是一个样本的特征,y 是目标向量,w 是待求系数。模型预测为 Xw,误差向量为 Xw - y。最小化 L(w) = ||Xw - y||²,展开后对 w 求梯度并令其为零,可得到正规方程 X^T X w = X^T y。当 X^T X 可逆时,解为 w = (X^T X)^(-1) X^T y。接下来我们从不同角度拆解这个公式。
一、误差函数构造与几何意义
选择平方误差并非唯一选择,但它在概率论和几何上都有良好性质。从概率角度看,如果误差服从均值为零的正态分布,最大化似然函数等价于最小化平方误差,因此最小二乘是高斯噪声假设下的最大似然估计。从几何角度看,线性回归可以理解为在 X 的列空间里寻找一个向量 Xw,使其与观测向量 y 的欧氏距离最短。这个距离最短的点正是 y 在列空间上的正交投影,最小二乘解就是在求这个投影坐标。
正交投影的视角解释了为什么正规方程是 X^T X w = X^T y。残差向量 r = y - Xw 应当垂直于 X 的每一个列向量,即 X^T r = 0。代入 r 后立即得到 X^T (y - Xw) = 0,整理后就是正规方程。这说明最小二乘不是简单地把点往线附近推,而是在一个由特征张成的子空间里做垂直投影。因此当特征之间存在强线性相关时,X 的列向量近似共线,列空间维度下降,投影变得不稳定,这也是共线性问题会影响最小二乘解的几何原因。
还需要区分竖直误差和垂直距离。教科书上的散点图中,最小二乘通常最小化的是竖直方向误差,也就是真实 y 与预测 y 之间的差,而不是点到直线的垂直距离。竖直误差假设 x 是相对确定的输入,y 是含有噪声的输出。如果 x 和 y 都含有误差,就需要使用总体最小二乘,它最小化点到直线的正交距离。理解这一点可以避免在特定测量误差场景下误用普通最小二乘。
二、Python 正规方程实现与验证
在 Python 中实现最小二乘最直接的方式是使用 numpy 构造矩阵并求解正规方程。下面生成一组带有线性关系的数据,并加入正态噪声,然后通过公式 w = (X^T X)^(-1) X^T y 计算系数。
import numpy as np
rng = np.random.default_rng(42)
x = rng.uniform(0, 10, 100)
noise = rng.normal(0, 1.5, 100)
y = 3.2 * x + 1.7 + noise
# 构建设计矩阵:第一列全为 1,对应截距项
X = np.column_stack([np.ones(x.shape[0]), x])
# 正规方程求解
theta = np.linalg.inv(X.T @ X) @ X.T @ y
print('正规方程系数:', theta)
print('理论系数约等于: [1.7 3.2]')
运行后得到的两个系数分别接近截距 1.7 和斜率 3.2。这里 X.T @ X 计算了 X^T X,np.linalg.inv 用于求逆,@ 表示矩阵乘法。由于数据中加入了噪声,估计值不会完全等于理论值,但会围绕真值小幅波动。如果只关心数值稳定性和代码简洁性,numpy 提供了 np.linalg.lstsq,它内部使用奇异值分解,可以处理 X^T X 不可逆或接近奇异的情况。
theta_lstsq, residuals, rank, singular = np.linalg.lstsq(X, y, rcond=None)
print('lstsq 系数:', theta_lstsq)
print('残差平方和:', residuals[0] if residuals.size else '未返回')
正规方程的计算量主要集中在 X^T X 的乘法与矩阵求逆。设特征数量为 p,求解线性方程组的时间复杂度约为 O(p³)。当 p 较小时,例如几百维以内,正规方程非常高效、实现简单,适合快速验证和教学演示。但当 p 很大,比如文本分类模型中的数万维特征,矩阵求逆会变得极其昂贵,而且 X^T X 可能因为样本数少于特征数而不可逆。这时梯度下降等迭代方法会成为更实际的选择。
三、梯度下降与最小二乘的迭代视角
梯度下降不直接求解闭式解,而是从一组初始系数出发,沿着损失函数负梯度方向逐步更新参数。对于平方误差损失 L(w) = (1/n) ||Xw - y||²,其梯度为 ∇L(w) = (2/n) X^T (Xw - y)。每次迭代让 w 减去学习率乘以梯度即可。虽然梯度下降需要设置学习率、迭代次数等超参数,但它把计算从 O(p³) 降到每次迭代 O(np),更适合高维稀疏数据。
# 梯度下降求解最小二乘
def gradient_loss(theta, X, y):
n = X.shape[0]
return (2.0 / n) * X.T @ (X @ theta - y)
theta_gd = np.zeros(X.shape[1])
learning_rate = 0.01
for epoch in range(5000):
grad = gradient_loss(theta_gd, X, y)
theta_gd -= learning_rate * grad
if epoch % 1000 == 0:
loss = np.mean((X @ theta_gd - y) ** 2)
print(f'epoch={epoch}, loss={loss:.4f}')
print('梯度下降系数:', theta_gd)
在数据量不大、特征维度不高时,正规方程和梯度下降都会收敛到相近的参数。二者的主要差别不在结果精度,而在计算资源与工程场景。正规方程一步到位但依赖矩阵可逆;梯度下降需要调参和观察损失曲线,却可以在线更新、支持海量数据。实际项目中,如果特征数量在几千以内,通常直接用 np.linalg.lstsq 或 sklearn.linear_model.LinearRegression;如果特征维度极高且数据是稀疏矩阵,则优先考虑 SGDRegressor 这类随机梯度下降实现。
此外,最小二乘还可以自然扩展为正则化回归。在损失函数后追加 L2 惩罚项,得到岭回归,其正规方程变为 w = (X^T X + λI)^(-1) X^T y,其中 λ 是正则化强度,I 是单位矩阵。这样即使 X^T X 不可逆,加上 λI 后通常也可逆,且能收缩系数、抑制共线性带来的方差放大。Python 中可以用 sklearn.linear_model.Ridge 或手动修改正规方程来验证这一效果。
四、常见误区与边界条件
第一个常见误区是认为最小二乘只能拟合直线。实际上最小二乘是一种损失函数和参数估计原则,可以用于多项式、指数、对数等多种曲线。只要模型对参数是线性的,例如 y = a + b x + c x²,仍可通过设计矩阵把 x、x² 作为新特征,用同样的正规方程求解。真正的限制在于模型必须对参数线性,而不是对变量线性。也正因如此,最小二乘可以配合基函数扩展处理很多非线性趋势。
第二个误区是只看 R² 高就认为模型好。最小二乘总是最小化训练集上的平方误差,但在特征过多或存在共线性时,模型可能在训练集上表现很好,却把噪声也拟合进去。解决方法包括使用验证集、交叉验证、正则化,以及检查特征相关性。共线性会让 X^T X 的条件数变得很大,系数符号可能与业务直觉相反,此时删除冗余特征或使用岭回归通常比强行解释普通最小二乘系数更可靠。
第三个容易忽略的边界条件是异常值。平方误差对较大的偏离非常敏感,一个远离主体的异常点就可能把拟合线明显拉向自己。若数据中存在明显异常值,可以考虑先用稳健回归方法,如 Huber 回归或 RANSAC,或者检查数据收集过程。最小二乘并非万能,但它为理解更复杂的模型提供了一个清晰基准:先比较你的模型是否优于一个简单线性回归,再判断复杂结构是否值得引入。