线性回归中的正规方程教程

我们遇到的大多数问题都有多种解法。比如，从房间一侧走到另一侧，我们可以沿着房间边缘绕过去，也可以直接穿过去。

正规方程正是对此概念的一种体现。它只是解决问题的另一种方式。您问是解决什么问题？我们会在本文余下部分介绍。现在，您只需要知道，在某些条件下实现线性回归时，它是一种高效的方法，能帮您节省大量时间。

让我们深入了解一下……

什么是正规方程？

正规方程是一种闭式解，用于在普通最小二乘线性回归中找到使代价函数最小化的 θ 值。另一种表述是：正规方程是一种解析方法，用来寻找使损失函数最小化的系数。两种说法都对，但它们具体意味着什么？ 

线性回归通过计算输入特征的加权和再加上一个偏置项来进行预测 y_hat。用数学表达如下：

其中 θ 表示参数，n 为特征数量。

本质上，上式所做的是对 θ 与 x 的点积求和。因此，更简洁的表示方式是使用向量化形式：

h(θ) 为假设函数。

给定这一近似的目标函数，我们可以用模型进行预测。为判断模型是否学得良好，衡量其在训练数据上的表现非常重要。为此我们计算损失函数。训练过程的目标是找到使损失函数最小化的 theta（θ）取值。

下面是损失函数的数学表示：

在上述方程中，theta（θ）是一个 n + 1 维向量，损失函数是该向量的函数。因此，需要对损失函数 J 分别对每个参数 θ_j 求偏导，并令其等于零。按照这一过程求解从 θ_0 到 θ_n 的所有 θ 值，便可得到使损失函数最小化的参数。

按照上述过程从 θ_0 到 θ_n 推导参数的解会非常繁琐。确实存在更快的解法。

来看一下正规方程的公式：

其中：

θ → 使损失函数最小化的参数 X → 每个样本的输入特征值 y → 每个样本的输出值向量

正规方程 vs 梯度下降

虽然两种方法都旨在找到使损失函数最小化的 theta（θ），但它们的求解思路大相径庭。

我们已在上一节介绍了正规方程的工作原理，下面简要回顾一下梯度下降，并给出两种技术的差异点。

梯度下降

梯度下降是机器学习中最常用的算法之一。它通过迭代寻找使损失函数最小化的 theta（θ）。

该过程首先评估模型表现。接着，从损失函数计算偏导数，用以指示当前点的斜率。最后，沿负梯度方向按比例迈步，通过更新当前参数集来向损失函数的最小值下降——见下式。

重复这一过程，直至在损失函数的最小值处收敛。

它们有何不同？

正规方程与梯度下降最明显的区别是前者是解析法。梯度下降采用迭代方式，即逐步更新参数直至收敛。由此带来的另一个细微差别是，梯度下降需要我们设定学习率，以控制朝向损失函数最小值迈步的步长。正规方程无需设定学习率，因为我们不进行迭代——结果可直接求得。

此外，采用正规方程时通常不需要进行特征缩放；我们常进行特征缩放是为了让特征具有相近的取值范围。这是因为梯度下降对数据点的取值范围较为敏感。使用梯度下降而未对特征归一化，可能会使损失函数的等高线图出现偏斜，而正规方程不存在这一问题。

何时使用正规方程

判断是否应优先使用正规方程而非梯度下降的最佳方式，是了解其劣势。

当数据集中的特征数量很大时，计算正规方程在计算上会变得困难。原因在于，为求解参数 θ，必须计算项 (X’ X)^-1。计算 X’ X 会得到一个 n×n 的矩阵，而在大多数计算机实现中，对矩阵求逆的开销大致随矩阵维度的立方增长。这意味着逆运算的时间复杂度为 O(n^3)，当 n 非常大时，会使正规方程变得极其缓慢——了解更多关于时间复杂度的信息。

因此，当数据集的特征数较多时，最好使用梯度下降。知名机器学习与 AI 专家 Andrew Ng 建议，当特征数 n 大于 10,000 时应考虑使用梯度下降。若特征数不超过 10,000，您可能更适合使用正规方程等解析方法，因为无需选择学习率这一超参数，从而少了一个需要调优的超参数。不过需要说明的是，像 R 或 Python 这样的编程语言在求解普通最小二乘线性回归的系数时，实际上多采用QR 分解，因为 QR 分解在数值上比正规方程更稳定、更高效。 

用 Python 从零实现正规方程

让我们生成一个回归问题来测试该方程：

import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import make_regression

# Generate a regression problem
X, y = make_regression(
    n_samples=100,
    n_features=2,
    n_informative=2,
    noise = 10,
    random_state=25
    )

# Visualize feature at index 1 vs target
plt.subplots(figsize=(8, 5))
plt.scatter(X[:, 1], y, marker='o')
plt.xlabel("Feature at Index 1")
plt.ylabel("Target")
plt.show()
技术支持

下面我们来实现正规方程：

# adds x0 = 1 to each instance
X_b = np.concatenate([np.ones((len(X), 1)), X], axis=1)
# calculate normal equation
theta_best = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)
# best values for theta
intercept, *coef = theta_best
print(f"Intercept: {intercept}\n\
Coefficients: {coef}")

Intercept: 0.35921242677977794
Coefficients: [6.129199175400593, 96.44309685893134]
技术支持

让我们通过做一次预测来检验模型：

# making a new sample
new_sample = np.array([[-2, 0.25]])
# adding a bias term to the instance
new_sample_b = np.concatenate([np.ones((len(new_sample), 1)), new_sample], axis=1)
# predicting the value of our new sample
new_sample_pred = new_sample_b.dot(theta_best)
print(f"Prediction: {new_sample_pred}")
Prediction: [12.21158829]
技术支持

每当您从零实现一个机器学习算法，拥有一种验证结果的方法总是有帮助的；Scikit-learn 是 Python 中最流行的机器学习库之一。它包含多种算法的实现，包括线性回归，我们将用它来验证我们的正规方程。

from sklearn.linear_model import LinearRegression

lr = LinearRegression()
lr.fit(X, y)

print(f"Intercept: {lr.intercept_}\n\
Coefficients: {lr.coef_}")
print(f"Prediction: {lr.predict(new_sample)}")

Intercept: 0.3592124267797807
Coefficients: [ 6.12919918 96.44309686]
Prediction: [12.21158829]
技术支持

两者的结果近似相同，因此可以确认我们的解是正确的。