最小二乘问题详解7：正则化最小二乘

1. 引言

在之前的文章《最小二乘问题详解4：非线性最小二乘》、《最小二乘问题详解5：非线性最小二乘求解实例》和《最小二乘问题详解6：梯度下降法》中分别介绍了使用Gauss-Newton方法（简称GN方法）和梯度下降法求解最小二乘问题之后，让我们插入另一个基础知识：正则化最小二乘（Regularized Least Squares），也就是大家常说的岭估计（Ridge Estimator），因为接下来要介绍的 Levenberg-Marquardt方法会用到这个思想。

本文讨论的岭估计在机器学习中也常称为岭回归（Ridge Regression）

2. 问题

复习《最小二乘问题详解2：线性最小二乘求解》中讨论的标准线性最小二乘问题：

其解为正规方程

的解。当

列满秩时，解唯一且为：

然而，在实际应用中，直接使用这个解可能会遇到一些问题。

1. 矩阵病态（Ill-conditioning）：
   • 矩阵

   A^T A

   的条件数可能非常大。

   • 这意味着即使数据

   b

   中有微小的扰动，也会导致解

   \theta^*

   发生剧烈变化，数值不稳定。

   • 条件数大的一个常见原因是特征之间高度相关（多重共线性），此时

   A^T A

   接近奇异，逆矩阵难以精确计算。

2. 过拟合（Overfitting）：
   • 当模型参数过多或特征维度很高时，标准最小二乘倾向于拟合训练数据中的噪声，导致泛化能力差。
   • 解

   \theta^*

   的模长（范数）可能非常大，表示模型对某些特征赋予了不合理的高权重。

3. 秩亏（Rank Deficiency）：
   • 若

   A

   不是列满秩（即

   \text{rank}(A) < n

   ），则

   A^T A

   是奇异的，无法求逆，正规方程有无穷多解。

矩阵的病态问题详见6.1节。

为了解决这些问题，我们引入正则化（Regularization）——通过在目标函数中加入一个额外的惩罚项，来“约束”解的行为，使其更稳定、更平滑，并避免过拟合。其实，正则化是一个通用的数学和机器学习思想：在优化问题中，通过向目标函数添加一个额外的“惩罚项”（penalty term），来引导解具有某种期望的性质，比如平滑性、稀疏性、小范数、结构简单等，从而提高模型的稳定性、泛化能力或可解释性。

3. 定义

最常用的正则化方法之一是Tikhonov 正则化，其核心思想是在原始的残差平方和基础上，加上一个关于参数

的L2范数平方的惩罚项：

其中：

是数据拟合项（data fidelity term），衡量模型对数据的拟合程度。

是正则化项（regularization term），也叫 Tikhonov 项，它惩罚较大的参数值。

是正则化参数（regularization parameter），控制正则化的强度：

：正则化作用消失，退化为标准最小二乘。

：正则化主导，迫使

，模型趋于常数零。

直观理解就是：通过加入正则项，不仅要保证预测误差小，还要保证模型参数不要太大。这相当于在“拟合数据”和“保持模型简单”之间做权衡。

4. 求解

将目标函数展开：

对

求梯度并令其为零：

整理得：

这就是岭估计的正规方程。

由于

，矩阵

是半正定的，而

是正定的，因此

是严格正定的，从而总是可逆的，无论

是否列满秩！于是，岭估计的闭式解为：

可见，岭估计只是在

的对角线上加了一个小量

，这相当于“抬升”了所有特征值，使得原本接近零的特征值变得远离零，从而显著改善了矩阵的条件数，提高了数值稳定性。

关于正定矩阵的问题详见6.2节。

从几何角度看，标准最小二乘是寻找使

投影最接近

的点；而岭估计则在此基础上“收缩”参数空间，使

向零靠近。这可以理解为给

的空间加上一个“弹簧”，防止参数跑得太远。

5. 分解

虽然有了闭式解，但在实际数值计算中，仍然不推荐直接计算

，因为

可能病态，显式构造

会损失精度。

5.1 QR分解

将正则化最小二乘问题转化为一个更大的最小二乘问题：

这是一个标准的最小二乘问题，可以用 QR 分解高效求解。

令：

对

做QR分解：

则解为：

其中

是

的前

行（对应

部分）。

5.2 SVD分解

设

是

的 SVD。

代入岭估计的闭式解：

利用

，得：

记

，则：

即：

对比标准最小二乘解（SVD 形式）：

可见，岭估计通过因子

对小奇异值方向进行了压制。当

时，标准解会爆炸（

），但岭估计中该项趋于 0，从而避免了对噪声方向的过度放大。

6. 补充

有一些《线性代数》、《矩阵论》相关的知识笔者也忘记了，这里就总结一下。如果有的读者很熟悉，可以直接略过。

6.1 病态矩阵

病态矩阵指的是矩阵条件数（Condition Number）很大的矩阵。矩阵病态会导致输入数据的微小扰动（如

或

的舍入误差），线性方程组

解就会剧烈变化。换句话说，系统对噪声极度敏感，数值计算中结果不可靠。

对于一个可逆矩阵

，其谱条件数（Spectral Condition Number）定义为：

其中：

是

的最大奇异值，

是

的最小奇异值。

注意这个定义也适用于非方阵

（

），只要

列满秩。条件数

的含义如下：

良态，解非常稳定

良态

中等病态

严重病态，浮点计算可能完全失效

双精度浮点数（约16位有效数字）完全失效

任何实矩阵

都可以进行奇异值分解（SVD）：

其中：

是左奇异向量矩阵（正交），

是右奇异向量矩阵（正交），

是对角矩阵，对角线元素为奇异值(Singular Value)

，

。

奇异值表示矩阵

在各个正交方向上的“拉伸”程度。如果

，说明

把某个方向“压扁”了，接近奇异，也就是矩阵病态。

6.2 正定矩阵

设

是一个实对称矩阵（即

），我们有如下定义：

如果对所有非零向量

，都有：

则称

是正定矩阵。

如果对所有非零向量

，都有：

则称

是半正定矩阵。

对于实对称矩阵

，我们可以进行谱分解（特征值分解）：

其中

是正交矩阵（列是特征向量），

是对角矩阵，对角元素是特征值。

那么：

正定

所有特征值

半正定

所有特征值

既然正定矩阵的所有特征值都 大于零，那么意味着正定矩阵是满秩，也就是正定矩阵一定可逆。

回到岭估计中的

：

是半正定的：

• 因为对任意

，有

• 所以它的所有特征值

1. 加上 \(\lambda I\)（其中 \(\lambda > 0\)）后： 新矩阵是 \(A^T A + \lambda I\) 它的特征值变为：原来的每个特征值

因此，

总是存在的，无论

是否列满秩。这正是岭估计的精髓所在：通过加一个小的正数

到对角线上，把原本可能奇异（不可逆）的

“修复”成一个严格正定、可逆的矩阵，从而保证了解的唯一性和数值稳定性。

7. 实例

如果线性最小二乘问题的设计矩阵

接近线性相关，那么普通方法求得的解不稳定，可以使用岭估计来给出稳定解。代码实现如下：

#include <Eigen/Dense>
#include <iostream>
#include <random>
#include <vector>

using namespace Eigen;
using namespace std;

int main() {
  // 设置随机数生成器
  default_random_engine gen;
  normal_distribution<double> noise(0.0, 0.5);  // 噪声 ~ N(0, 0.5)

  // 数据生成：y = x1 + x2 + 0.1*x3 + noise
  // 但 x3 ≈ x1 + x2 （高度相关，造成病态）
  int n_samples = 20;
  int n_features = 3;
  MatrixXd A(n_samples, n_features);
  VectorXd b(n_samples);

  for (int i = 0; i < n_samples; ++i) {
    double x1 = (double)rand() / RAND_MAX * 10;
    double x2 = (double)rand() / RAND_MAX * 10;
    double x3 = x1 + x2 + (double)rand() / RAND_MAX * 0.1;  // x3 ≈ x1 + x2

    A(i, 0) = x1;
    A(i, 1) = x2;
    A(i, 2) = x3;

    double true_y = 1.0 * x1 + 1.0 * x2 + 0.1 * x3;
    b(i) = true_y + noise(gen);  // 加噪声
  }

  //使用 SVD 计算条件数
  BDCSVD<MatrixXd> svd(A, ComputeThinU | ComputeThinV);
  cout << "Condition number of A: "
       << svd.singularValues()(0) / svd.singularValues()(2)
       << endl;  // 假设 3 个奇异值

  // 普通最小二乘解：theta = (A^T A)^{-1} A^T b
  VectorXd theta_ols = A.bdcSvd(ComputeThinU | ComputeThinV).solve(b);

  // 岭回归解：theta = (A^T A + λI)^{-1} A^T b
  MatrixXd AtA = A.transpose() * A;
  VectorXd Atb = A.transpose() * b;
  double lambda = 0.01;
  MatrixXd ridge_matrix =
      AtA + lambda * MatrixXd::Identity(n_features, n_features);
  VectorXd theta_ridge = ridge_matrix.ldlt().solve(Atb);

  // 输出结果
  cout << "True weights: [1.0, 1.0, 0.1]" << endl;
  cout << "OLS weights:  " << theta_ols.transpose() << endl;
  cout << "Ridge weights:" << theta_ridge.transpose() << endl;

  return 0;
}

这里可以看到，我们使第三个已知参数向量，可以用第一个已知参数向量和第二个已知参数向量接近线性表示（

），那么设计矩阵就接近线性相关，这个设计矩阵就是病态的。分别用普通最小二乘求解和岭估计来求解，最后结果如下：

Condition number of A: 715.286
True weights: [1.0, 1.0, 0.1]
OLS weights:    2.00569   2.03977 -0.938346
Ridge weights: 1.02399  1.05935 0.038262

理论上来说，使用QR/SVD方法可以一定程度上解决矩阵的病态问题。但是从这个实例结果可以看到，即使使用最稳健的 SVD 方法求解普通最小二乘，参数估计仍严重偏离真实值。这是因为问题本身的结构特征高度相关，导致参数不可识别。

相比之下，岭估计通过引入正则化项

，显著改善了矩阵的条件数，给出了稳定且接近真实的参数估计。岭估计通过牺牲一定的精度，提升模型的泛化能力，保证在噪声存在下仍能稳定预测。岭估计的另外一个问题是正则化参数

不太好进行给定，往往需要一些先验经验的辅助才能确定，有机会再进一步进行讨论了。