1.1 线性回归的含义

线性回归试图找到一个线性方程，来描述自变量 $X$ 与因变量 $y$ 之间的关系。对于单变量情况，就是一条直线：

$y = w_1 x + w_0$

其中 $w_1$ 是斜率，$w_0$ 是截距。在机器学习中，我们统称为模型的权重（weights） 。

对于多元情况，可以写成矩阵形式：

$y = Xw$

其中 $X$ 是 $n \times p$ 的设计矩阵（每行一个样本，每列一个特征），$w$ 是 $p \times 1$ 的权重向量。

1.2 最小二乘法

给定一条直线，每个数据点都有一个残差（residual）：真实值与预测值之差。我们需要一个指标来衡量“总体残差”的大小。

最经典的做法是：最小化残差平方和。这就是最小二乘法（Ordinary Least Squares, OLS）。

损失函数为：

$L(w) = \sum_{i=1}^{n} (y_i - X_i w)^2 = \|y - Xw\|_2^2$

1.3 矩阵求导推导闭式解

现在我们来手推这个最小化问题的解析解。

首先展开损失函数：

$L(w) = (y - Xw)^T(y - Xw) = y^T y - y^T X w - w^T X^T y + w^T X^T X w$

由于 $y^T X w$ 和 $w^T X^T y$ 都是标量且互为转置（相等），所以：

$L(w) = y^T y - 2w^T X^T y + w^T X^T X w$

对 $w$ 求梯度：

$\nabla_w L(w) = -2X^T y + 2X^T X w$

关于矩阵的梯度（求导）：

对于整理后的损失函数 $L(w) = y^T y - 2w^T X^T y + w^T X^T X w$

• 第1项 $y^T y$：与 $w$ 无关，梯度为 0。
• 第2项 $-2w^T X^T y$：这可以看成 $-2 \cdot w^T a$ 的形式（其中 $a = X^T y$）。使用公式 $\nabla_w (w^T a) = a$，所以这一项梯度为 $-2a = -2X^T y$。
• 第3项 $w^T X^T X w$：这是典型的二次型 $w^T A w$，其中 $A = X^T X$，且 $A$ 是对称矩阵。使用公式 $\nabla_w (w^T A w) = (A + A^T)w = 2A w$，所以这一项梯度为 $2X^T X w$。

三项相加，立即得到：$\nabla_w L(w) = -2X^T y + 2X^T X w$

令梯度为零：

$-2X^T y + 2X^T X w = 0$

$X^T X w = X^T y$

这就是著名的正规方程（Normal Equation） 。当 $X^T X$ 可逆时：

$\boxed{w = (X^T X)^{-1} X^T y}$

关于闭式解的存在前提

这个闭式解要求 $X^T X$ 是满秩的（即特征之间不存在完全共线性）。如果特征数量大于样本数量（$p > n$），或者存在多重共线性，$X^T X$ 就是奇异矩阵，无法求逆。

1.4 Python 手动实现

1
import numpy as np
2
import matplotlib.pyplot as plt
3

4
class LinearRegression:
5
    """手动实现线性回归（最小二乘法）"""
6

7
    def __init__(self):
8
        self.w = None
9
        self.intercept_ = None
10

11
    def fit(self, X, y):
12
        """
13
        X: shape (n_samples, n_features)
14
        y: shape (n_samples,)
15
        """
16
        # 添加一列全1，用于处理截距项
17
        X_b = np.c_[np.ones((X.shape[0], 1)), X]  # (n, p+1)
18

19
        # 正规方程求解
20
        # w = (X^T X)^(-1) X^T y
21
        self.w = np.linalg.inv(X_b.T @ X_b) @ X_b.T @ y
22

23
        self.intercept_ = self.w[0]
24
        self.coef_ = self.w[1:]
25
        return self
26

27
    def predict(self, X):
28
        X_b = np.c_[np.ones((X.shape[0], 1)), X]
29
        return X_b @ self.w
30

31
# 生成示例数据
32
np.random.seed(42)
33
X = 2 * np.random.rand(100, 1)
34
y = 4 + 3 * X + np.random.randn(100, 1)
35

36
# 训练模型
37
model = LinearRegression()
38
model.fit(X, y.flatten())
39
print(f"截距: {model.intercept_:.4f}")
40
print(f"系数: {model.coef_[0]:.4f}")
41
# 输出: 截距: 4.2221, 系数: 2.7709
42

43
# 可视化
44
X_test = np.array([[0], [2]])
45
y_pred = model.predict(X_test)
46
plt.scatter(X, y, alpha=0.6)
47
plt.plot(X_test, y_pred, 'r-', linewidth=2)
48
plt.xlabel("X")
49
plt.ylabel("y")
50
plt.show()

2.1 正则化

当 $X^T X$ 不可逆时（特征数 > 样本数，或存在多重共线性），OLS 无法求解。即使可逆，模型也可能过拟合 —— 在训练集上表现很好，但在新数据上很差。

正则化（Regularization） 的核心思想是：在损失函数中加入一个惩罚项，限制模型参数的大小，从而防止过拟合。

2.2 Ridge 的数学形式

Ridge 回归（又称 L2 正则化）在 OLS 损失函数上增加了权重平方和作为惩罚项：

$L_{\text{ridge}}(w) = \|y - Xw\|_2^2 + \lambda \|w\|_2^2$

其中 $\lambda \geq 0$ 是正则化强度超参数。

写成展开形式：

$L_{\text{ridge}}(w) = y^T y - 2y^T X w + w^T X^T X w + \lambda w^T w$

对 $w$ 求梯度：

$\nabla_w L_{\text{ridge}}(w) = -2X^T y + 2X^T X w + 2\lambda w$

令梯度为零：

$(X^T X + \lambda I) w = X^T y$

$\boxed{w = (X^T X + \lambda I)^{-1} X^T y}$

Ridge 能处理不可逆矩阵的根本原因

即使 $X^T X$ 不可逆，加上 $\lambda I$ 后矩阵一定可逆（因为 $\lambda > 0$ 时，$X^T X + \lambda I$ 是正定矩阵）。这就是 Ridge 能处理不可逆矩阵的根本原因。

2.3 Ridge 的几何解释

Ridge 等价于在以下约束下最小化残差平方和：

$\sum_{j=1}^{p} w_j^2 \leq t$

在二维空间中，L2 约束区域是一个圆形。圆形的边界是光滑的，最优解可能落在圆周的任何位置，但很难正好落在坐标轴上，因此 Ridge 通常不会将系数压缩到 0 —— 它做的是收缩（shrinkage），而不是选择（selection）。

2.4 Python 实现

1
class RidgeRegression:
2
    """手动实现Ridge回归（L2正则化）"""
3

4
    def __init__(self, alpha=1.0):
5
        self.alpha = alpha  # 即 lambda
6
        self.w = None
7

8
    def fit(self, X, y):
9
        X_b = np.c_[np.ones((X.shape[0], 1)), X]
10
        n_features = X_b.shape[1]
11

12
        # Ridge闭式解: (X^T X + alpha * I)^(-1) X^T y
13
        # 注意：截距项通常不参与正则化，这里为了简便对所有参数都加了惩罚
14
        I = np.eye(n_features)
15
        self.w = np.linalg.inv(X_b.T @ X_b + self.alpha * I) @ X_b.T @ y
16
        return self
17

18
    def predict(self, X):
19
        X_b = np.c_[np.ones((X.shape[0], 1)), X]
20
        return X_b @ self.w
21

22
# 对比不同alpha的效果
23
from sklearn.datasets import make_regression
24
X, y = make_regression(n_samples=50, n_features=10, noise=0.1, random_state=42)
25

26
for alpha in [0, 0.1, 1, 10]:
27
    model = RidgeRegression(alpha=alpha)
28
    model.fit(X, y)
29
    coef_norm = np.linalg.norm(model.w[1:])  # 忽略截距
30
    print(f"alpha={alpha:4.1f}, 系数L2范数={coef_norm:.4f}")
31
# 输出: alpha越大，系数范数越小（收缩越强）

3.1 Lasso 的数学形式

Lasso（Least Absolute Shrinkage and Selection Operator）使用 L1 正则化，惩罚项是权重的绝对值之和：

$L_{\text{lasso}}(w) = \|y - Xw\|_2^2 + \lambda \|w\|_1 = \|y - Xw\|_2^2 + \lambda \sum_{j=1}^{p} |w_j|$

与 Ridge 不同，Lasso 的损失函数在 $w_j=0$ 处不可微，因此没有简单的闭式解。

3.2 Lasso 产生稀疏解的原因

1. 几何解释

Lasso 等价于在约束 $\sum_{j=1}^{p} |w_j| \leq t$ 下最小化残差平方和。

在二维空间中，L1 约束区域是一个菱形。菱形的顶点在坐标轴上。当目标函数的等高线（椭圆）与约束区域相切时，切点很可能落在菱形的顶点上，此时对应的系数就被压缩为 0。

随着维度增加，多维菱形的“角”越来越多，系数被置零的概率就越大。这就是 Lasso 能够自动进行特征选择的原因——它会产生稀疏解。

2. L1 产生稀疏解而 L2 不产生的原因

• L2（圆形）：边界光滑，切点落在坐标轴上的概率为 0
• L1（菱形）：边界有尖角，切点很容易落在尖角上（即坐标轴）

总结：L1 的几何形状有“角”，所以容易产生零系数。

在3.3节现有推导基础上，我从专业角度为你补充了三块核心内容，旨在将“跳跃的数学步骤”彻底展开，并赋予其直观的统计含义。请查收修订后的完整 3.3节：

3.3 软阈值算子（Soft-Thresholding Operator）

虽然 Lasso 没有闭式解，但在单变量情况下，我们可以推导出解析解。

考虑单变量 Lasso 问题（没有截距）：

$\min_{w} \frac{1}{2}(y - Xw)^2 + \lambda |w|$

其中 $X$ 和 $y$ 都是标量（为简单起见，假设 $X^T X = 1$，即特征已标准化）。

目标是：

$\min_{w} \frac{1}{2}(y - Xw)^2 + \lambda |w|$

次梯度（Subgradient）条件的完整展开

由于绝对值函数在 $w=0$ 处不可导，我们需要引入次梯度 $\partial |w|$。对于绝对值函数，其次导数为：

令 $z = X^T y$ 且利用标准化条件 $X^T X = 1$，损失函数的次梯度为：

$\partial_w L_{\text{lasso}} = -X(y - Xw) + \lambda \partial |w| = (w - z) + \lambda \partial |w|$

根据凸优化的最优性条件（$0 \in \partial f(w)$），分三种情况讨论：

1. 当 $w > 0$ 时：$\partial |w| = 1$，令 $w - z + \lambda = 0$，得 $w = z - \lambda$。为保证 $w>0$，须有 $z > \lambda$。
2. 当 $w < 0$ 时：$\partial |w| = -1$，令 $w - z - \lambda = 0$，得 $w = z + \lambda$。为保证 $w<0$，须有 $z < -\lambda$。
3. 当 $w = 0$ 时：$\partial |w| = [-1, 1]$，令 $0 \in (0 - z) + \lambda[-1, 1]$，即 $z \in [-\lambda, \lambda]$ 时成立，等价于 $|z| \leq \lambda$。

综合三种情况，得到软阈值算子（Soft-Thresholding Operator）：

$\boxed{S(z, \lambda) = \text{sign}(z) \cdot \max(|z| - \lambda, 0)}$

其中 $z = X^T y$（在标准化特征下，$z$ 也等于普通最小二乘的单变量解 $w_{OLS}$）。

z 的统计直觉与“软阈值”的物理意义

$z$ 的本质：在标准化特征下，$z = X^T y$ 实际上就是特征 $X$ 与响应变量 $y$ 的协方差（相关性）。它衡量了当前特征对预测结果的边际贡献强度。

阈值的逻辑：

• 如果相关性 $|z| \leq \lambda$，说明该特征的贡献强度被正则化惩罚完全“淹没”，模型判定其为噪声，系数直接置零（实现特征筛选）。
• 如果相关性 $|z| > \lambda$，说明该特征有效，但系数不是直接取 $z$，而是向零收缩 $\lambda$ 的量（减去 $\lambda$）。这种“收缩”不同于 Ridge 的等比缩放，它是绝对值上的硬性削减，因此得名“软阈值”（Soft-Thresholding）。

代码实现：

1
def soft_threshold(z, lam):
2
    """软阈值算子"""
3
    return np.sign(z) * np.maximum(np.abs(z) - lam, 0)
4

5
# 演示
6
z_values = np.linspace(-3, 3, 100)
7
for lam in [0, 0.5, 1, 2]:
8
    s_values = [soft_threshold(z, lam) for z in z_values]
9
    plt.plot(z_values, s_values, label=f'λ={lam}')
10
plt.axhline(0, color='black', linewidth=0.5)
11
plt.axvline(0, color='black', linewidth=0.5)
12
plt.xlabel('z')
13
plt.ylabel('S(z, λ)')
14
plt.legend()
15
plt.title('软阈值算子')
16
plt.show()

3.4 坐标下降法求解 Lasso

对于多维 Lasso，通常使用坐标下降法（Coordinate Descent） 求解。核心思想是每次只优化一个参数，固定其他所有参数。这样一来，原本复杂的 $p$ 维优化问题，就被拆解成了无数个一维的单变量 Lasso 问题。

在第 $t$ 轮迭代中，当我们准备更新第 $j$ 个特征对应的权重 $w_j$ 时：

• 第一步：计算当前残差与第 $j$ 个特征的“相关系数”（专业叫法叫偏残差）。这个值就相当于 3.3 节里的 $z$，它衡量了剔除了其他特征影响后，当前这个特征还有多大的贡献潜力。
• 第二步：把贡献潜力套入 3.3 节的软阈值公式 $S(z, \lambda)$ 里。
  • 如果潜力值很小（$|z| \le \lambda$），软阈值公式直接输出 0 —— 把这个特征踢出模型；
  • 如果潜力值足够大（$|z| > \lambda$），软阈值公式就输出 $z - \lambda$ 或 $z + \lambda$ —— 保留特征，但强制砍掉 $\lambda$ 的冗余分量。

即在坐标下降法中，软阈值算子会逐坐标地对特征进行“相关性截断”，从而在多维场景中高效地产生稀疏解。这使得 Lasso 在超高维数据（$p \gg n$）中依然具备强大的特征选择能力。

对于 Lasso 的第 $j$ 个坐标更新：

$w_j^{(new)} = S\left(\frac{1}{n} \sum_{i=1}^{n} X_{ij} (y_i - \sum_{k \neq j} X_{ik} w_k), \lambda\right)$

即计算当前残差与第 $j$ 个特征的相关性，然后应用软阈值。

和 Ridge 回归的直观对比：

• Ridge 的坐标更新是给系数乘以一个小于 1 的固定比例（如 (0.9 \times w_j)），系数只会无限逼近 0，但永远到不了 0。
• Lasso 的坐标更新是硬生生的“阈值截断”，只要贡献潜力不够大，系数直接精确等于 0。 这个“绝对零”的特性，就是坐标下降法让 Lasso 产生稀疏解的工程实现。

scikit-learn 的 Lasso 实现默认就使用坐标下降法。

1
from sklearn.linear_model import Lasso
2
from sklearn.datasets import make_regression
3

4
X, y = make_regression(n_samples=100, n_features=20, n_informative=5, noise=0.1, random_state=42)
5

6
# 注意：alpha越大，稀疏性越强
7
for alpha in [0.01, 0.1, 1]:
8
    lasso = Lasso(alpha=alpha, max_iter=10000)
9
    lasso.fit(X, y)
10
    n_nonzero = np.sum(np.abs(lasso.coef_) > 1e-6)
11
    print(f"alpha={alpha:.2f}, 非零系数个数={n_nonzero}")
12
# alpha越大，非零系数越少（越稀疏）

4.1 结合 L1 和 L2 的原因

• Lasso（L1）：能做特征选择，产生稀疏解。但在特征高度相关时，Lasso 会随机选择其中一个，导致不稳定。
• Ridge（L2）：系数收缩但不置零，在相关特征上表现稳定。

ElasticNet 结合了两者的优点。

4.2 ElasticNet 的数学形式

ElasticNet 的损失函数为：

$L_{\text{elastic}}(w) = \frac{1}{2n}\|y - Xw\|_2^2 + \alpha \cdot \rho \|w\|_1 + \frac{\alpha(1-\rho)}{2}\|w\|_2^2$

其中：

• $\alpha$：整体正则化强度
• $\rho$（即 l1_ratio）：L1 和 L2 的混合比例
  • $\rho = 0$：纯 Ridge（L2）
  • $\rho = 1$：纯 Lasso（L1）
  • $0 < \rho < 1$：两者的组合

ElasticNet 也可以用坐标下降法求解。scikit-learn 的 ElasticNet 默认使用坐标下降。

1
from sklearn.linear_model import ElasticNet
2

3
X, y = make_regression(n_samples=100, n_features=20, n_informative=5, noise=0.1, random_state=42)
4

5
for l1_ratio in [0, 0.5, 1]:
6
    enet = ElasticNet(alpha=0.1, l1_ratio=l1_ratio)
7
    enet.fit(X, y)
8
    n_nonzero = np.sum(np.abs(enet.coef_) > 1e-6)
9
    print(f"l1_ratio={l1_ratio:.1f}, 非零系数个数={n_nonzero}")
10
# l1_ratio=0 (Ridge): 全部非零
11
# l1_ratio=0.5: 部分稀疏
12
# l1_ratio=1 (Lasso): 最稀疏

5.1 MSE（均方误差）

$\text{MSE} = \frac{1}{n} \sum_{i=1}^{n} (y_i - \hat{y}_i)^2$

MSE 是线性回归最常用的损失函数。它对大误差给予更大的惩罚（平方放大）。

5.2 RMSE（均方根误差）

$\text{RMSE} = \sqrt{\text{MSE}} = \sqrt{\frac{1}{n} \sum_{i=1}^{n} (y_i - \hat{y}_i)^2}$

RMSE 只是 MSE 开方，量纲与原始数据一致，更易解释。

5.3 MAE（平均绝对误差）

$\text{MAE} = \frac{1}{n} \sum_{i=1}^{n} |y_i - \hat{y}_i|$

MAE 是线性的，所有误差被同等对待。

5.4 异常值敏感性对比

这是三者最核心的区别：

• 直观解释：最小化 MSE 等价于预测均值 —— 均值对异常值敏感；最小化 MAE 等价于预测中位数 —— 中位数对异常值鲁棒
• 实际建议：如果异常值代表重要的业务场景（如欺诈检测），用 MSE/RMSE；如果异常值只是噪声/脏数据，用 MAE；如果想兼顾两者，可以用 Huber 损失（在误差小时用平方，误差大时用绝对值）

代码实现：

1
from sklearn.metrics import mean_squared_error, mean_absolute_error
2

3
# 模拟正常数据和异常值
4
y_true = np.array([1, 2, 3, 4, 5])
5
y_pred_normal = np.array([1.1, 1.9, 3.2, 3.8, 5.1])
6
y_pred_outlier = np.array([1.1, 1.9, 3.2, 3.8, 50])  # 最后一个点是异常预测
7

8
print("=== 正常情况 ===")
9
print(f"MSE: {mean_squared_error(y_true, y_pred_normal):.4f}")
10
print(f"RMSE: {np.sqrt(mean_squared_error(y_true, y_pred_normal)):.4f}")
11
print(f"MAE:  {mean_absolute_error(y_true, y_pred_normal):.4f}")
12

13
print("\n=== 存在异常预测 ===")
14
print(f"MSE: {mean_squared_error(y_true, y_pred_outlier):.4f}")
15
print(f"RMSE: {np.sqrt(mean_squared_error(y_true, y_pred_outlier)):.4f}")
16
print(f"MAE:  {mean_absolute_error(y_true, y_pred_outlier):.4f}")
17
# 可以看到：MSE/RMSE被异常值大幅拉高，而MAE变化相对温和

总结

方法	正则项	闭式解	稀疏性	适用场景
OLS	无	✅ $(X^TX)^{-1}X^Ty$	❌	特征数 < 样本数，无共线性
Ridge	L2: $\lambda\|w\|_2^2$	✅ $(X^TX + \lambda I)^{-1}X^Ty$	❌	多重共线性，特征相关
Lasso	L1: $\lambda\|w\|_1$	❌（坐标下降）	✅	特征选择，高维数据
ElasticNet	L1 + L2	❌（坐标下降）	✅（可控）	相关特征分组，兼顾选择与稳定

核心要点回顾

1. 最小二乘法的闭式解 $w = (X^TX)^{-1}X^Ty$ 是线性回归的基石，但要求 $X^TX$ 可逆。
2. Ridge（L2） 通过在 $X^TX$ 对角线加 $\lambda$ 保证矩阵可逆，系数收缩但不归零。
3. Lasso（L1） 的稀疏性来自两个层面：

• 几何：L1 约束是菱形，更容易在坐标轴上相切
• 算法：软阈值算子 $S(z, \lambda) = \text{sign}(z)\max(|z|-\lambda, 0)$ 会将小系数直接置零

4. ElasticNet 结合 L1 和 L2，在特征选择的同时保持分组稳定性。
5. MSE/RMSE 对异常值敏感（平方放大误差），MAE 相对鲁棒（绝对值不放大）。