【机器学习】线性回归实验

任何一门机器学习的课，第一个接触的算法多半都是线性回归。这篇文章简单总结了相关的两个算法——最小二乘法和梯度下降法，以及它们的步骤，并用Python实现了相关算法。

步骤如下：

1. 首先用随机函数构造数据，在这里我构造了一百组数据。

2. 建立线性模型$y=wx+b$。

3. 分别用最小二乘法和梯度下降法训练模型，并与流行的机器学习第三方库得出的结果做对比。

4. 绘制样本点和预测的直线。

数据准备

方便起见，这里直接用numpy.random生成随机数作为数据，让y和x保持一定关系并添加噪声。

arr = []
for i in range(100):
    x = np.random.rand()
    y = 1.5 * x + 0.1 + np.random.rand() - 0.5
    arr.append((x, y, x * y, x * x))
    df = pd.DataFrame(arr, columns=['x', 'y', 'x*y', 'x*x'])

在这里，令$y=1.5x+0.1$并加上扰动项。

最小二乘法

相信大家对最小二乘法并不陌生，在高中数学里就已经多次接触过。在这里不过多说明原理，实现很简单，根据公式即可。
$$\begin{gathered} w=\frac{\sum xy-n\bar{x}\bar{y}}{\sum x^2-n{\bar{x}}^2} \\ b=\bar{y}-w\bar{x} \end{gathered}$$
Python实现：

# 最小二乘法
# 计算均值及其它变量
x_mean = np.mean(df['x'])
y_mean = np.mean(df['y'])
xy_sum = np.sum(df['x*y'])
xx_sum = np.sum(df['x*x'])
# 计算斜率与截距
w1 = (xy_sum - 100 * x_mean * y_mean) / (xx_sum - 100 * x_mean * x_mean)
b1 = y_mean - w1 * x_mean

梯度下降法

最小二乘法可以一步到位，但梯度下降法就不是这么回事了。下面简单介绍一下梯度下降法，这里我们只讨论简单的线性回归形式，输入x只有一个属性输出为y，共有m个样本$(x_i,y_i),i=1,\mathrm{2...}m$，我们可以使用均方误差确定w和b，均方误差是回归任务中最常用的性能度量，我们试图通过均方误差最小来求解w和b即
$$\begin{array}{rl}(w,b)& =\underset{(w,b)}{\mathrm{argmin}}\frac{1}{2m}\sum _{i=1}^m(f(x_i)-y_i{)}^2\\ & =\underset{(w,b)}{\mathrm{argmin}}\frac{1}{2m}\sum _{i=1}^m(w{x}_i+b-y_i{)}^2\end{array}$$
我们称其为代价函数或损失函数,我们只需要使其最小即可。函数沿着导数方向是变化最快的，为了更快的达到优化目标，沿着负梯度方向搜寻w,b使其最小，即使用梯度下降法更新权重即可求出w和b。

求导以后具体的转移方程为：
$$\begin{gathered} w^{\ast }=w-\eta \frac{1}{m}\sum _{i=1}^m(w{x}_i+b-y_i)x_i \\ {b}^{\ast }=b-\eta \frac{1}{m}\sum _{i=1}^m(w{x}_i+b-y_i) \end{gathered}$$
$\eta$为学习率（步长），随机初始化w,b,通过不断的迭代这两个式子计算w,b的最优值。

Python实现：

# 梯度下降法
e = []  # 残差
w2 = 100
b2 = 100
alpha = 0.1
for i in range(1000):
    temp = df['x'] * w2 + b2 - df['y']
    e.append((i, np.sum(temp)))
    w2 -= alpha * np.sum(temp * df['x']) / 100
    b2 -= alpha * np.sum(temp) / 100

与库函数的综合比较

说了这么多，怎么知道自己的算法实现得如何呢？究竟准不准确呢？这就需要我们画图来观察了。

准备流行库scipy作为对照，并用matplotlib.pyplot画图检验结果：

可以看到，三者得出的曲线相当接近，以至于在最后的综合对比图中三线合一，互相重叠而导致看不出区别了。

但单单从曲线上看是看不到数值的最终差异的，于是最后在控制台输出三种方法得出的斜率与截距。

可以看到，在最后得到的结果中，调库与最小二乘法得出的结果几乎一致，而梯度下降法得出的结果略有差异。由此可以得出结论，最小二乘法比梯度下降法更加准确。

源代码展示

注：有些代码只是我为了方便处理或表示，并没有什么重要的含义。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

plt.rcParams['font.sans-serif'] = ['SimHei']  # 修改字体配置
plt.rcParams['axes.unicode_minus'] = False  # 解决负号编码问题

if __name__ == '__main__':
    # 构造数据
    arr = []
    for i in range(100):
        x = np.random.rand()
        y = 1.5 * x + 0.1 + np.random.rand() - 0.5
        arr.append((x, y, x * y, x * x))
    df = pd.DataFrame(arr, columns=['x', 'y', 'x*y', 'x*x'])

    # 最小二乘法
    # 计算均值及其它变量
    x_mean = np.mean(df['x'])
    y_mean = np.mean(df['y'])
    xy_sum = np.sum(df['x*y'])
    xx_sum = np.sum(df['x*x'])
    # 计算斜率与截距
    w1 = (xy_sum - 100 * x_mean * y_mean) / (xx_sum - 100 * x_mean * x_mean)
    b1 = y_mean - w1 * x_mean

    # 梯度下降法
    w2 = 100
    b2 = 100
    alpha = 0.1
    for i in range(1000):
        temp = df['x'] * w2 + b2 - df['y']
        w2 -= alpha * np.sum(temp * df['x']) / 100
        b2 -= alpha * np.sum(temp) / 100

    # 调库
    popt, pcov = curve_fit(lambda X, w, b: w * X + b, df['x'], df['y'])
    w3 = popt[0]
    b3 = popt[1]

    # 效果展示
    w = [w1, w2, w3]
    b = [b1, b2, b3]
    color = ['red', 'green', 'blue']
    text = ['最小二乘法', '梯度下降法', '库函数', '综合对比']
    for i in range(4):
        plt.subplot(2, 2, i + 1)
        plt.scatter(df['x'], df['y'], 15)
        if i < 3:
            plt.plot(df['x'], df['x'] * w[i] + b[i], color=color[i])
            plt.title(text[i])
        else:
            plt.plot(df['x'], df['x'] * w1 + b1, color='red')  # 最小二乘法
            plt.plot(df['x'], df['x'] * w2 + b2, color='green')  # 梯度下降法
            plt.plot(df['x'], df['x'] * w3 + b3, color='blue')  # 调库
            plt.title(text[3])
    print(f"最小二乘法斜率为{w1:.4f}，截距为{b1:.4f}")
    print(f"梯度下降法斜率为{w2:.4f}，截距为{b2:.4f}")
    print(f"调用三方库斜率为{w3:.4f}，截距为{b3:.4f}")
    plt.show()