基于Python的机器学习系列（2）：梯度下降-CSDN博客

本文链接：https://blog.csdn.net/ljd939952281/article/details/141022625

在我们完成了机器学习系列的第一个部分——交叉验证之后，现在让我们继续深入探讨机器学习中的另一个关键主题：梯度下降。在本篇文章中，我们将详细介绍梯度下降及其在线性回归中的应用。

线性回归

线性回归是一种监督学习算法，用于根据连续的y值进行预测。

例如，给定以下数据：

样本	鸡蛋价格	黄金价格	石油价格	GDP
1	3	100	4	21
2	4	500	7	43

我们希望使用鸡蛋价格、黄金价格和石油价格来预测GDP。我们将鸡蛋价格、黄金价格和石油价格称为特征或X。我们将希望预测的值称为标签或目标或y。每一行称为样本或观测。

术语和符号

在本系列博文中，我们将使用以下符号：

x_j^(i) 表示第i个样本中的第j个特征。例如，x_1^(1) 表示第一个样本的鸡蛋价格（即3），x_2^(1) 表示第一个样本的黄金价格（即100），x_3^(1) 表示第一个样本的石油价格（即4）。
粗体大写字母X表示特征的整个矩阵，其中包含m行样本和n列特征。X的形状为(m, n)。
粗体小写字母x表示单个向量（即列）特征。任何向量默认都是列向量，形状为(?, )，其中？表示一个整数。表示行向量时，我们用x^T。
y^(i) 表示第i个样本的标签/目标，y表示所有标签的向量，形状为(m, )。

假设函数

假设函数将给定的输入X映射到预测的y。我们必须学习/训练这个函数。为了区分实际的y和预测的y，我们通常将预测的y称为yhat。

对于线性回归，假设函数（记作h_theta(x)，表示h取决于参数theta的x）定义如下：

h_theta(x^(i)) = theta_0 + theta_1*x^(i)_1 + theta_2*x^(i)_2 + ... + theta_n*x^(i)_n
               = sum(theta_j * x^(i)_j) for j from 0 to n
               = theta^T * x
               = X * theta

这里theta称为参数或权重或系数，用于将线性映射从X到y。我们通常不写成第一种形式，而是喜欢用第二种求和形式、第三种向量化形式（theta的形状为(n, 1)，x的形状也为(n, 1)），以及第四种矩阵形式（将所有样本组合在一起，其中X的形状为(m, n)，theta的形状为(n, 1)）。我们喜欢写成第四种形式，因为它容易实现。

注意：当theta_0存在时，X的形状为(m, n+1)，theta的形状为(n+1, 1)。

假设函数的结果称为模型，这个过程称为训练模型。测试模型的过程称为推断。

如何找到最佳theta ==> 梯度下降

如何找到最佳参数theta？

首先，我们必须定义一个损失/目标函数（记作J），用于衡量在给定的theta下，h(x^(i))与对应的y^(i)的接近程度：

J(theta) = (1/2) * sum((h_theta(x^(i)) - y^(i))^2) for i from 1 to m

注意，这里的1/2只是为了后续数学推导的方便。

我们希望选择theta以最小化J(theta)。这个过程称为优化。

theta* = argmin_theta J(theta)

一种非常成功的优化算法称为梯度下降算法，该算法基于根据导数（也称为梯度）更新theta。

这个图表说明了为什么梯度下降有效：

正式定义

我们从随机的theta开始，并反复执行以下更新：

theta_j := theta_j - alpha * (partial derivative of J with respect to theta_j)

这个更新同时对所有j=0, 1, ..., n进行（如果你忘了，j表示每个特征）。这里，alpha称为学习率（lr），范围从0到1。通常，我们尝试0.001作为默认值，这个值必须手动选择。任何参数，例如学习率，称为超参数。

为了实现这个算法，我们必须找到损失函数关于每个theta_j的偏导数。让我们首先对单个训练样本进行推导：

partial derivative of J with respect to theta_j = (h_theta(x) - y) * x_j

要修改更新规则以包含所有训练样本，我们将更新规则修改为包括求和：

theta_j := theta_j - alpha * sum((h_theta(x^(i)) - y^(i)) * x_j^(i)) for every j

或者

theta := theta - alpha * X^T * (h - y)

由于这种梯度下降算法使用整个训练集中的每个样本来计算梯度，我们称之为批量梯度下降。

实现

梯度下降有以下步骤：

准备数据
- 添加截距
- 将X、y和w调整到正确的形状
  - X -> (m, n)
  - y -> (m, )
  - w -> (n, )
- 训练-测试集拆分
- 特征缩放
- 清除任何缺失数据
- （可选）特征工程
预测并计算损失
- 损失函数是均方误差
```
J = ((h - y)^2) / 2
```
其中h为
```
h = X * theta
```

计算基于损失的梯度

损失函数的梯度为

partial derivative of J with respect to theta_j = X^T * (h - y)

使用更新规则更新theta
```
theta = theta - alpha * X^T * (h - y)
```
其中alpha是典型的学习率，范围在0到1之间
循环步骤2-4直到达到最大迭代次数或旧损失与新损失之间的差异小于预定义的阈值

让我们动手实践吧！

1. 准备数据

1.1 获取正确形状的X和y

# 1. 加载糖尿病数据集 
# 作为我们的回归案例研究
from sklearn.datasets import load_diabetes
# 类型 - Bunch
# Bunch - numpy数据的字典
# diabetes.feature_names
# print(diabetes)
diabetes = load_diabetes()

print("特征: ", diabetes.feature_names)

X = diabetes.data
X.shape # 样本数量，特征数量

m = X.shape[0]  # 样本数量
n = X.shape[1]  # 特征数量

print(m, n)

y = diabetes.target

# X中的行数与y中的行数相同
# 因为我们对所有y都有预测值
assert m == y.shape[0]

1.2 训练-测试集拆分

# 测试数据的适当大小是多少
# 70/30（小数据集）；80/20（中等数据集）；90/10（大数据集）；
# 为什么大数据集可以设置测试集大小为10，
# 因为大数据集的10％已经足够测试准确性
from sklearn.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size = 0.3)

assert len(X_train)  == len(y_train)
assert len(X_test)   == len(y_test)

1.3 特征缩放

# 我想对我的数据进行标准化，使均值为0，方差为1
# 在每个特征上平均，而不是在每个样本上
# 为什么我们需要标准化
# 因为标准化通常可以让我们更快地收敛
# 为什么 -> 因为数值在较小范围内
# 因此，梯度也在有限范围内，不会发散

from sklearn.preprocessing import StandardScaler

# 1. StandardScaler.fit(X)  # 这个scaler（或self）知道均值和标准差，因此现在
# 它知道如何转换数据
# 2  X = StandardScaler.transform(X)  # 非原地；将返回一些东西

# 1. StandardScaler.fit_transform(X) -> 依次执行1和2

# 创建一个StandardScaler对象
# StandardScaler是一个类
# scaler称为实例/对象

# 几乎总是，使用归一化或标准化来特征缩放数据
# 如果你假设你的数据是高斯分布的，使用标准化，否则，使用归一化

scaler = StandardScaler()

X_train = scaler.fit_transform(X_train)
X_test = scaler.transform(X_test)

#注意：你必须先拆分，然后再缩放....否则，你会得到数据泄漏

1.4 添加截距

# 他们想要的X的形状是什么
# (样本数量，特征数量) --> 正确的形状
# 截距呢
# w0是我们的截距
# w的形状是什么 -->(n+1, )
# 截距的形状是什么 -->(m, 1)
#X = [1 2 3     @  [w0
#     1 4 6         w1
#     1 9 1         w2 
#     1 10 2 ] 

# np.ones((形状))
intercept = np.ones((X_train.shape[0], 1))

# 根据轴=1连接截距
X_train = np.concatenate((intercept, X_train), axis=1)

# np.ones((形状))
intercept = np.ones((X_test.shape[0], 1))

# 根据轴=1连接截距
X_test = np.concatenate((intercept, X_test), axis=1)

1.5 特征工程（可选）

特征工程有时是有用的（例如，多项式，核函数），以便创建与目标之间的非线性关系。

这里我们将跳过这个步骤。

2. 训练你的算法

2.1 定义你的算法

from time import time

# 步骤1：准备数据
# X_train, X_test已经连接了截距
# [1, 特征
#  1, 特征....]

# 确保我们的X_train与y_train具有相同的样本大小
assert X_train.shape[0] == y_train.shape[0]

# 初始化我们的w
# 我们不必执行X.shape[1] + 1，因为我们的X_train已经有了
# 截距
# w = theta/系数
theta = np.zeros(X_train.shape[1])

# 定义学习率
# 以后，你会知道最好让它慢慢减少
# 一旦我们执行了很多迭代，我们希望更新速度减慢，这样可以更好地收敛
lr = 0.1

# 定义我们的最大迭代次数
# 通常称之为epochs <---机器学习人喜欢这样称呼
max_iter = 50

def h_theta(X, theta):
    return X @ theta

def mse(yhat, y):
    return ((yhat - y)**2).sum() / yhat.shape[0]

def gradient(X, error):
    m = X.shape[0]
    return (X.T @ error) / X.shape[0]

start = time()

# 定义你的for循环
for i in range(max_iter):
    
    # 1. yhat = X @ w
    # 预测
    # yhat (m, ) = (m, n) @ (n, )
    yhat = h_theta(X_train, theta)

    # 2. error = yhat - y_train
    # 用于计算梯度的误差
    # error (m, ) = (m, ) - (m, )
    error = yhat - y_train

    # 3. grad = X.T @ error
    # grad (n, ) = (n, m) @ (m, )
    # 每个特征j的梯度
    grad = gradient(X_train, error)

    # 4. w = w - lr * grad
    # 更新w
    # w (n, ) = (n, ) - 标量 * (n, )
    theta = theta - lr * grad

time_taken = time() - start

2.2 计算准确率/损失

# 我们得到了我们可爱的w
# 现在是时候检查我们的准确性了
# 1. 进行预测
yhat = h_theta(X_test, theta)

# 2. 计算均方误差
mse = mse(yhat, y_test)

# 打印均方误差
print("MSE: ", mse)
print("Time used: ", time_taken)

绘制图表

import matplotlib.pyplot as plt

plt.figure(figsize=(5,5))
plt.scatter(y_test, yhat, c='crimson')

p1 = max(max(yhat), max(y_test))
p2 = min(min(yhat), min(y_test))
plt.plot([p1, p2], [p1, p2], 'b-')
plt.xlabel('True Values', fontsize=15)
plt.ylabel('Predictions', fontsize=15)
plt.axis('equal')
plt.show()

检查你的理解

使用鸡蛋价格、黄金价格等表格，请完成以下问题：

画出向量 x_2
画出向量 x^(2)
x^(2)_1 是多少
y^(1) 是多少
画出向量 y
画出整个矩阵 X
给定第一个样本 x^(1)，以及权重 θ 为 [2, 0.1, 3]，求 y_hat。请手动计算。
给定整个矩阵 X，并使用上述相同的权重，求 y_hat。请手动计算。
给定整个矩阵 X，给定权重 θ 为 [2, 0.1, 3]，以及 y，求 J(θ)。请手动计算。
对于训练样本1，即 x^(1)，求偏导数 ∂J/∂θ_2。请手动计算。
对于所有样本，即 X，求偏导数 ∂J/∂θ_2。请手动计算。
对于所有样本，求偏导数 ∂J/∂θ。请手动计算。
给定 α = 0.01，新 θ 是多少？
证明 θ_j := θ_j - α * ∑((h_θ(x^(i)) - y^(i)) * x_j^(i)) 与 θ := θ - α * X^T * (h - y) 相同。
在梯度下降的开始，初始权重是多少？
复习如何求导数/梯度，特别是乘法和链式法则。
用你自己的话解释为什么求导数可以帮助你找到更好的权重。