(3)线性回归

325 阅读5分钟

  线性回归是利用数理统计中回归分析,来确定两种或两种以上变量间相互依赖的定量关系的一种统计分析方法,运用十分广泛。回归用于预测输入变量和输出变量之间的关系,特别是当输入变量的值发生变化时,输出变量的值也随之发生变化。回归模型正是表示从输入变量到输出变量之间映射的函数。

一元线性回归

一元线性模型非常简单,假设我们有变量 xix_i 和目标 yiy_i,每个 i 对应于一个数据点,希望建立一个模型

y^i=wxi+b\hat{y}_i = w x_i + b

y^i\hat{y}_i 是我们预测的结果,希望通过 y^i\hat{y}_i 来拟合目标 yiy_i,通俗来讲就是找到这个函数拟合 yiy_i 使得误差最小,即最小化值使用函数表示为:

L=1ni=1n(y^iyi)2L = \frac{1}{n} \sum_{i=1}^n(\hat{y}_i - y_i)^2

L 就是俗称为代价函数,其即为预测值与真实值之间的平均的平方距离,统计中一般称其为MAE(mean square error)均方误差。把之前的函数式代入代价函数,并且将需要求解的参数w和b看做是函数L的自变量,可得

L(w,b)=1ni=1n(wxi+byi)2L(w,b) = \frac{1}{n} \sum_{i=1}^n(w x_i + b - y_i)^2

现在的任务是求解最小化L时w和b的值,即核心目标优化式为

(w,b)=argmin(w,b)i=1n(wxi+byi)2(w^*, b^*) = \arg \underset{(w,b)}{min} \sum_{i=1}^n(w x_i + b - y_i)^2

求解方式有两种

1)最小二乘法(least square method)

求解 w 和 b 是使损失函数最小化的过程,在统计中,称为线性回归模型的最小二乘“参数估计”(parameter estimation)。我们可以将 L(w,b) 分别对 w 和 b 求导,得到(最小二乘法的本质是什么)

Lw=w(i=1n((wxi)2+2(wxi)(byi)+(byi)2)=2i=1n(wxi+byi)xi=2(wi=1nxi2i=1nxi(yib))=0\frac{∂L}{∂w} = \frac{∂}{∂w}(\sum_{i=1}^n ((wx_i)^2 + 2(wx_i)(b-y_i) + (b - y_i)^2) \\ = 2 \sum_{i=1}^n (wx_i + b - y_i)x_i = 2 (w\sum_{i=1}^n x_i^2 - \sum_{i=1}^n x_i(y_i - b)) = 0

Lb=w(i=1n((wxiyi)2+2(wxiyi)b+b2)=2i=1n(wxi+byi)=2(nbi=1n(yiwxi))=0\frac{\partial L}{\partial b} = \frac{∂}{∂w}(\sum_{i=1}^n ((wx_i-y_i)^2 + 2(wx_i-y_i)b + b^2) \\ = 2 \sum_{i=1}^n (wx_i + b - y_i) = 2 (nb - \sum_{i=1}^n (y_i - wx_i)) = 0

令上述两式为0,可得到 w 和 b 最优解的闭式(closed-form)解:

w=i=1nyi(xix)i=1nxi21n(i=1nxi)2w = \frac{\sum_{i=1}^n y_i(x_i - \overline{x} ) }{\sum_{i=1}^n x_i^2 - \frac{1}{n}(\sum_{i=1}^n x_i )^2}
b=1ni=1n(yiwxi)b = \frac{1}{n} \sum_{i=1}^n (y_i - wx_i)

2)梯度下降(gradient descent)

我们接触到的第一个优化算法,非常简单,但是却非常强大,在深度学习中被大量使用,所以让我们从简单的例子出发了解梯度下降法的原理。

梯度

梯度在数学上就是导数,如果是一个多元函数,那么梯度就是偏导数。比如一个函数f(x, y),那么 f 的梯度就是

(fx, fy)(\frac{\partial f}{\partial x},\ \frac{\partial f}{\partial y})

可以称为 grad f(x, y) 或者 ∇𝑓(𝑥,𝑦)∇f(x,y)。具体某一点 (𝑥0, 𝑦0)(x0, y0) 的梯度就是 ∇𝑓(𝑥0, 𝑦0)∇f(x0, y0)。

梯度有什么意义呢?从几何意义来讲,一个点的梯度值是这个函数变化最快的地方,具体来说,对于函数 f(x, y),在点 (𝑥0,𝑦0)(x0,y0) 处,沿着梯度 ∇𝑓(𝑥0, 𝑦0)∇f(x0, y0) 的方向,函数增加最快,也就是说沿着梯度的方向,我们能够更快地找到函数的极大值点,或者反过来沿着梯度的反方向,我们能够更快地找到函数的最小值点。

梯度下降法

有了对梯度的理解,我们就能了解梯度下降发的原理了。上面我们需要最小化这个误差,也就是需要找到这个误差的最小值点,那么沿着梯度的反方向我们就能够找到这个最小值点。

我们可以来看一个直观的解释。比如我们在一座大山上的某处位置,由于我们不知道怎么下山,于是决定走一步算一步,也就是在每走到一个位置的时候,求解当前位置的梯度,沿着梯度的负方向,也就是当前最陡峭的位置向下走一步,然后继续求解当前位置梯度,向这一步所在位置沿着最陡峭最易下山的位置走一步。这样一步步的走下去,一直走到觉得我们已经到了山脚。当然这样走下去,有可能我们不能走到山脚,而是到了某一个局部的山峰低处。

1.png

类比我们的问题,就是沿着梯度的反方向,我们不断改变 w 和 b 的值,最终找到一组最好的 w 和 b 使得误差最小。

v2-b722c2fca0ea2c1bc71975dd965d0c97_720w.webp

在更新的时候,我们需要决定每次更新的幅度,比如在下山的例子中,我们需要每次往下走的那一步的长度,这个长度称为学习率,用 η\eta 表示,这个学习率非常重要,不同的学习率都会导致不同的结果,学习率太小会导致下降非常缓慢,学习率太大又会导致跳动非常明显。

最后我们的更新公式就是

w:=wηf(w, b)wb:=bηf(w, b)bw := w - \eta \frac{\partial f(w,\ b)}{\partial w} \\ b := b - \eta \frac{\partial f(w,\ b)}{\partial b}

通过不断地迭代更新,最终我们能够找到一组最优的 w 和 b,这就是梯度下降法的原理。

代码

import torch
import numpy as np
from torch.autograd import Variable
import matplotlib.pyplot as plt


def fun1():
    torch.manual_seed(2017)
    x_train = np.array([[3.3], [4.4], [5.5], [6.71], [6.93], [4.168],
                        [9.779], [6.182], [7.59], [2.167], [7.042],
                        [10.791], [5.313], [7.997], [3.1]], dtype=np.float32)

    y_train = np.array([[1.7], [2.76], [2.09], [3.19], [1.694], [1.573],
                        [3.366], [2.596], [2.53], [1.221], [2.827],
                        [3.465], [1.65], [2.904], [1.3]], dtype=np.float32)

    plt.ion()
    plt.figure()
    plt.plot(x_train, y_train, 'bo')
    plt.show()

    # 转换成 Tensor
    x_train = torch.from_numpy(x_train)
    y_train = torch.from_numpy(y_train)
    # w = Variable(torch.randn(1), requires_grad=True)  # 随机初始化
    # b = Variable(torch.zeros(1), requires_grad=True)  # 使用 0 进行初始化
    w = torch.randn(1, requires_grad=True)
    b = torch.zeros(1, requires_grad=True)
    print('w:', w)
    print('b:', b)

    def linear_model(x):
        return x * w + b

    y_ = linear_model(x_train)
    plt.figure()
    plt.plot(x_train.data.numpy(), y_train.data.numpy(), 'bo', label='real')
    plt.plot(x_train.data.numpy(), y_.data.numpy(), 'ro', label='estimated')
    plt.legend()
    plt.show()

    # 计算误差
    def get_loss(y_, y):
        return torch.mean((y_ - y_train) ** 2)

    loss = get_loss(y_, y_train)
    print(loss)
    # 自动求导
    loss.backward()
    # 查看 w 和 b 的梯度
    print(w.grad)
    print(b.grad)
    # 更新一次参数
    w.data = w.data - 1e-2 * w.grad.data
    b.data = b.data - 1e-2 * b.grad.data
    y_ = linear_model(x_train)
    plt.figure()
    plt.plot(x_train.data.numpy(), y_train.data.numpy(), 'bo', label='real')
    plt.plot(x_train.data.numpy(), y_.data.numpy(), 'ro', label='estimated')
    plt.legend()

    for e in range(100):  # 进行 10 次更新
        y_ = linear_model(x_train)
        loss = get_loss(y_, y_train)

        w.grad.zero_()  # 记得归零梯度
        b.grad.zero_()  # 记得归零梯度
        loss.backward()

        w.data = w.data - 1e-2 * w.grad.data  # 更新 w
        b.data = b.data - 1e-2 * b.grad.data  # 更新 b
        # print('epoch: {}, loss: {}, {}'.format(e, loss.item(), w.grad))
        print("epoch:{}, loss:{}, w:{}-{}, b:{}-{}".format(e, loss, w, w.grad, b, b.grad))
        # plt.figure()
        # plt.plot(x_train.data.numpy(), y_train.data.numpy(), 'bo', label='real')
        # plt.plot(x_train.data.numpy(), y_.data.numpy(), 'ro', label='estimated')
        # plt.show()
        # # plt.pause(0.5)
        # plt.close()
        # input("Press Enter to Continue")  # 之后改成识别输出

    y_ = linear_model(x_train)
    plt.figure()
    plt.plot(x_train.data.numpy(), y_train.data.numpy(), 'bo', label='real')
    plt.plot(x_train.data.numpy(), y_.data.numpy(), 'ro', label='estimated')
    plt.legend()
    print("w", w)
    print("b", b)
    print("-" * 10)


if __name__ == '__main__':
    fun1()

执行输出:

w: tensor([2.2691], requires_grad=True)
b: tensor([0.], requires_grad=True)
tensor(153.3520, grad_fn=<MeanBackward0>)
tensor([161.0043])
tensor([22.8730])
epoch:0, loss:3.135774850845337, w:tensor([0.4397], requires_grad=True)-tensor([21.9352]), b:tensor([-0.2576], requires_grad=True)-tensor([2.8870])
epoch:1, loss:0.3550890386104584, w:tensor([0.4095], requires_grad=True)-tensor([3.0163]), b:tensor([-0.2593], requires_grad=True)-tensor([0.1687])
epoch:2, loss:0.30295437574386597, w:tensor([0.4051], requires_grad=True)-tensor([0.4424]), b:tensor([-0.2573], requires_grad=True)-tensor([-0.2005])
epoch:3, loss:0.30131959915161133, w:tensor([0.4041], requires_grad=True)-tensor([0.0922]), b:tensor([-0.2548], requires_grad=True)-tensor([-0.2502])
...
epoch:93, loss:0.2522490322589874, w:tensor([0.3743], requires_grad=True)-tensor([0.0294]), b:tensor([-0.0477], requires_grad=True)-tensor([-0.2048])
epoch:94, loss:0.2518215477466583, w:tensor([0.3740], requires_grad=True)-tensor([0.0294]), b:tensor([-0.0456], requires_grad=True)-tensor([-0.2043])
epoch:95, loss:0.2513962388038635, w:tensor([0.3737], requires_grad=True)-tensor([0.0293]), b:tensor([-0.0436], requires_grad=True)-tensor([-0.2037])
epoch:96, loss:0.250973105430603, w:tensor([0.3734], requires_grad=True)-tensor([0.0292]), b:tensor([-0.0416], requires_grad=True)-tensor([-0.2032])
epoch:97, loss:0.250552237033844, w:tensor([0.3731], requires_grad=True)-tensor([0.0291]), b:tensor([-0.0395], requires_grad=True)-tensor([-0.2027])
epoch:98, loss:0.25013336539268494, w:tensor([0.3728], requires_grad=True)-tensor([0.0291]), b:tensor([-0.0375], requires_grad=True)-tensor([-0.2022])
epoch:99, loss:0.24971675872802734, w:tensor([0.3725], requires_grad=True)-tensor([0.0290]), b:tensor([-0.0355], requires_grad=True)-tensor([-0.2017])
w tensor([0.3725], requires_grad=True)
b tensor([-0.0355], requires_grad=True)
  1. 初始数据:
    1.png

  2. 第一次拟合:
    2.png

  3. 第二次拟合:
    3.png

  4. 最终效果:
    4.png

参考:

最小二乘法的本质是什么?
梯度下降算法原理讲解——机器学习
什么是梯度下降法?
梯度下降(Gradient Descent)小结