手撕岭回归

岭回归

在线性回归模型中，其参数估计公式为 β=(X^TX)^-1X^TY ，当 X^TX 不可逆时无法求出 β ，另外如果 |X^TX| 越趋近于0，会使得回归系数趋向于无穷大，此时得到的回归系数是无意义的。解决这类问题可以使用岭回归

推导过程

线性回归的目标函数

J(β)=∑(y-Xβ)²

为了保证回归系数 β 可求，岭回归模型在目标函数上加了一个L2范数的惩罚项

$$\begin{gathered} J(\beta )=\sum (y-X\beta {)}^2+\lambda ||\beta |{|}_2^2 \\ =\sum (y-X\beta {)}^2+\sum {\lambda \beta }^2 \end{gathered}$$

其中 λ 为非负数， λ 越大，则为了使 J(β) 最小，回归系数 β 就越小。

J(β)=(y-Xβ)^T(y-Xβ)+λβ^Tβ

=y^Ty-y^TXβ-β^TX^Ty+β^TX^TXβ+λβ^Tβ

使用微分方程

$$\begin{gathered} \frac{\partial J(\beta )}{\partial \beta }=0 \\ \Rightarrow 0-X^{T}y-X^{T}y+2X^{T}X\beta +2\lambda \beta =0 \\ \Rightarrow \beta =(X^{T}X+\lambda I{)}^{-1}{X}^{T}y \end{gathered}$$

L2范数惩罚项的加入使得 (X^TX+λI) 满秩，保证了可逆，但是也由于惩罚项的加入，使得回归系数 β 的估计不再是无偏估计。所以岭回归是以放弃无偏性、降低精度为代价解决病态矩阵问题的回归方法。

单位矩阵 I (eye) 的对角线上全是1，像一条山岭一样，这也是岭回归名称的由来。

λ的选择

• 模型的方差：回归系数的方差
• 模型的偏差：预测值和真实值的差异

随着 λ 的增大
模型的方差就越小；
β 的估计值更加偏离真实值，模型的偏差就越大；
岭回归的关键是找到一个合理的 λ 来平衡模型的方差和偏差。

后续换专门介绍方差和偏差

单位矩阵

在线性代数中，n阶单位矩阵，是一个n*n的方形矩阵，其主对角线元素为1，其余元素为0

记为 I (或者E)

代码实现

这次我们在之前的matrix类创建一个静态方法

各大语言里面基本都是定义为eye，估计是怕跟逆混淆

class Matrix {
  ……
  /**
   * 生成单位矩阵
   * @param n 尺寸
   */
  static eye(n) {
    const mat = new Matrix(n, n);
    for (let i = 0; i < n; i++) {
      mat.arr[i][i] = 1;
    }
    return mat;
  }
  ……
}

调用一下这个静态方法

// node --experimental-modules mat-advance.mjs
import Matrix from './matrix.mjs';
console.log('单位矩阵', Matrix.eye(3));

运行结果

单位矩阵 Matrix { _arr:
  [
    [ 1, 0, 0 ],
    [ 0, 1, 0 ],
    [ 0, 0, 1 ]
  ]
}

岭回归开发

代码实现

其实岭回归就是，在之前的 X^TX 之中多加了一个 λI

所以我们增加一个方法 ridge ，代表设置岭回归，然后之前的 getBeta 增加对 λ 的判断

import Matrix from '../matrix/matrix.mjs'
class Regression {
  ……
  /**
   * 岭回归
   * @param lambda 系数
   */
  ridge(lambda) {
    this.lambda = Number(lambda);
  }
  ……
  /**
   * 获取beta
   * @returns {Matrix}  参数矩阵
   */
  getBeta() {
    const xMat = new Matrix(this.x);
    const yMat = new Matrix(this.y);
    let denom = xMat.T.multiply(xMat);
    if (this.lambda) denom = denom.plus(Matrix.eye(denom.m).multiply(this.lambda));
    const beta = denom.I.multiply(xMat.T.multiply(yMat)); // beta矩阵公式
    return beta;
  }
}

还是之间的数据，我们增加一个设置

原始代码看之前的手撕多元线性回归

// node --experimental-modules ridge.mjs
import Regression from './regression.mjs'
const model = new Regression() // 构造回归模型
……
model.ridge(0.1); // 岭回归
model.fit(x, y); // 训练模型
const predictions = model.predict(testX); // 预测数据
console.log(predictions);

这时候再看看结果

[
  [ 360.7185157070388 ],
  [ 245.89242970494647 ],
  [ 441.9219487447597 ],
  [ 443.09254205964197 ],
  [ 217.6102287623099 ]
]

测试数据

我们需要选择最佳的 λ 所以需要测试

下一篇我会详细介绍一下数据判断

---

完全自己用代码写一遍，是不是更能有助于理解机器学习呢？

源码:

• https://gitee.com/thales-ucas/manuai.git
• https://github.com/thales-ucas/manuai.git