自求导方法实现线性回归算法-优快云博客

线性回归是统计学中用于建立因变量（目标值）与一个或多个自变量（特征）之间关系的预测模型。在众多预测模型中，线性回归因其简单、直观、易于理解而广受欢迎。它的核心思想是寻找一个线性方程，使得预测值与实际值之间的差异最小。在机器学习中，这种差异通常通过损失函数来量化，而自求导方法，尤其是梯度下降算法，是寻找最优解的有效手段。本文将详细探讨如何使用自求导方法来实现线性回归算法，并提供代码示例。

线性回归模型

线性回归模型基于一个简单的假设：因变量 (Y) 与自变量 (X) 之间存在线性关系。这种关系可以用数学公式表示为：

[ $Y = \beta_0 + \beta_1 X + \epsilon$ ]

其中，( $\beta_0$ ) 是模型的截距项，表示当所有自变量为零时因变量的期望值；( $\beta_1$ ) 是模型的斜率，表示自变量 ( $X$ ) 每变化一个单位，因变量 ( $Y$ ) 预期将如何变化；( $\epsilon$ ) 是误差项，表示模型无法解释的随机变异。

损失函数

损失函数是衡量模型预测值与实际值之间差异的函数。在线性回归中，最常用的损失函数是均方误差（MSE），它计算所有样本的预测误差平方的平均值。MSE 的定义如下：

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

其中，( $y_i$ ) 是第 ( $i$ ) 个样本的实际观测值，( $\hat{y}_i$ ) 是模型对第 ( $i$ ) 个样本的预测值，( $n$ ) 是样本总数。MSE 通过惩罚预测误差的平方，鼓励模型寻找能够最小化预测误差的参数。

参数估计

在确定了损失函数之后，下一步是估计模型参数，即找到能够最小化损失函数的 ( $\beta_0$ ) 和 ( $\beta_1$ ) 的值。这通常通过计算损失函数关于参数的梯度并更新参数来实现。梯度是一个向量，指向损失函数增长最快的方向，因此，梯度的反方向就是损失函数下降最快的方向。

梯度计算

对于MSE损失函数，梯度计算如下：

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

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

这些梯度表达式告诉我们，为了最小化损失函数，我们需要如何调整每个参数。具体来说，我们需要对每个参数的梯度乘以一个负的学习率，然后更新参数。

参数更新

使用梯度和学习率 ( $\alpha$ ) 更新参数的过程如下：

[ $\beta_0 = \beta_0 - \alpha \frac{\partial \text{MSE}}{\partial \beta_0}$ ]

[ $\beta_1 = \beta_1 - \alpha \frac{\partial \text{MSE}}{\partial \beta_1}$ ]

这里，( $\alpha$ ) 是学习率，它控制着我们在梯度下降方向上迈出的步长。如果学习率太大，我们可能会跳过最小值；如果学习率太小，收敛可能会非常缓慢。

迭代优化

参数更新的过程需要重复进行，直到满足某些停止条件。这些条件可能包括梯度的大小低于某个阈值，或者达到预设的最大迭代次数。每次迭代都包括计算当前参数下的梯度，并更新参数，目的是逐步逼近损失函数的最小值。

模型评估

一旦模型训练完成，我们使用测试数据集来评估模型的性能。常用的评估指标包括均方误差（MSE）、平均绝对误差（MAE）和决定系数（R²）。这些指标提供了模型预测准确性的量化度量，帮助我们理解模型在未见数据上的表现。

代码实现

以下是使用Python实现简单线性回归的代码示例：

import numpy as np

# 生成模拟数据
np.random.seed(0)
X = 2 * np.random.rand(100, 1)
y = 4 + 3 * X + np.random.randn(100, 1)

# 添加截距项
X_b = np.c_[np.ones((100, 1)), X]

# 定义学习率和迭代次数
alpha = 0.1
iterations = 1000

# 初始化参数
theta = np.random.randn(2, 1)

# 梯度下降
for iteration in range(iterations):
    predictions = X_b.dot(theta)
    errors = predictions - y
    gradients = X_b.T.dot(errors) / len(X)
    theta = theta - alpha * gradients

# 打印结果
print(f"截距项（beta_0）: {theta[0][0]}")
print(f"斜率（beta_1）: {theta[1][0]}")

# 预测新数据
X_new = np.array([[0], [2]])
X_new_b = np.c_[np.ones((2, 1)), X_new]
y_new_pred = X_new_b.dot(theta)
print(f"预测值: {y_new_pred}")

代码解释

数据生成：首先，我们使用 np.random.rand 生成100个随机样本作为自变量 ( $X$ )，并根据线性关系 ( $y = 4 + 3x + \epsilon$ ) 生成对应的因变量 ( $y$ )，其中 ( $\epsilon$ ) 是正态分布的随机噪声。
添加截距项：在线性回归中，我们通常需要在模型中添加一个截距项。这里，我们通过 np.c_[np.ones((100, 1)), X] 将截距项（全为1的列向量）与自变量 ( $X$ ) 合并。
参数初始化：我们随机初始化模型参数 ( $\theta$ )，它包含了截距项 ( $\beta_0$ ) 和斜率 ( $\beta_1$ )。
梯度下降：在迭代过程中，我们首先计算模型的预测值，然后计算预测值与实际值之间的误差。接着，我们计算梯度，即损失函数关于参数的导数，并使用梯度下降法更新参数。
结果输出：最后，我们输出训练后的截距项和斜率，并使用这些参数预测新的数据点。

损失函数的计算

在代码中，我们没有直接实现损失函数的计算，但是可以通过以下代码添加对MSE损失函数的计算，以便在每次迭代中监控损失函数的变化：

# 损失函数计算
def compute_mse(y_true, y_pred):
    return ((y_true - y_pred) ** 2).mean()

# 监控损失函数的变化
mse_history = []

for iteration in range(iterations):
    predictions = X_b.dot(theta)
    errors = predictions - y
    gradients = X_b.T.dot(errors) / len(X)
    theta = theta - alpha * gradients
    mse = compute_mse(y, predictions)
    mse_history.append(mse)

# 打印最终的MSE
print(f"最终的MSE: {mse}")

这段代码定义了一个 compute_mse 函数来计算MSE，并在每次迭代中记录MSE的值。这样，我们不仅可以观察到参数的变化，还可以监控模型性能的改进。

参数更新的可视化

为了更直观地理解参数是如何随着迭代次数变化的，我们可以将参数的变化绘制成图表：

import matplotlib.pyplot as plt

# 绘制参数变化
plt.plot(theta[0], theta[1], 'r--')
plt.xlabel('Beta 0 (Intercept)')
plt.ylabel('Beta 1 (Slope)')
plt.title('Parameter Convergence')
plt.show()

这段代码使用 matplotlib 库来绘制参数 ( $\beta_0$ ) 和 ( $\beta_1$ ) 的变化。通过观察图表，我们可以直观地看到参数是如何随着迭代次数的增加而收敛的。