零基础学习机器学习(七)线性回归的基本原理和简单实现

本节一起跟随李沐老师的课程动手学习线形回归的基本思想,并从零实现。

1. 线性回归的基本原理

回归(regression)是能为一个或多个自变量与因变量之间关系建模的一类方法。 在自然科学和社会科学领域,回归经常用来表示输入和输出之间的关系。

在机器学习领域中的大多数任务通常都与预测(prediction)有关。 当我们想预测一个数值时,就会涉及到回归问题。 常见的例子包括:预测价格(房屋、股票等)、预测住院时间(针对住院病人等)、 预测需求(零售销量等)。 但不是所有的预测都是回归问题。 在后面的章节中,我们将介绍分类问题。分类问题的目标是预测数据属于一组类别中的哪一个。

1.1 线性回归的基本元素

线性回归(linear regression)可以追溯到19世纪初, 它在回归的各种标准工具中最简单而且最流行。 线性回归基于几个简单的假设: 首先,假设自变量x和因变量y之间的关系是线性的, 即y可以表示为x中元素的加权和,这里通常允许包含观测值的一些噪声; 其次,我们假设任何噪声都比较正常,如噪声遵循正态分布。

为了解释线性回归,我们举一个实际的例子: 我们希望根据房屋的面积(平方英尺)和房龄(年)来估算房屋价格(美元)。 为了开发一个能预测房价的模型,我们需要收集一个真实的数据集。 这个数据集包括了房屋的销售价格、面积和房龄。 在机器学习的术语中,该数据集称为训练数据集(training data set) 或训练集(training set)。 每行数据(比如一次房屋交易相对应的数据)称为样本(sample), 也可以称为数据点(data point)或数据样本(data instance)。 我们把试图预测的目标(比如预测房屋价格)称为标签(label)或目标(target)。 预测所依据的自变量(面积和房龄)称为特征(feature)或协变量(covariate)。

1.1.1 线性模型

线性假设是指目标(房屋价格)可以表示为特征(面积和房龄)的加权和,如下面的式子:

price=warea×area+wage×age+bprice = w_{area} \times area + w_{age} \times age + b

上式中的wareaw_{area}wagew_{age} 称为权重(weight),权重决定了每个特征对我们预测值的影响。b称为偏置(bias)、偏移量(offset)或截距(intercept)。 偏置是指当所有特征都取值为0时,预测值应该为多少。 即使现实中不会有任何房子的面积是0或房龄正好是0年,我们仍然需要偏置项。 如果没有偏置项,我们模型的表达能力将受到限制。 严格来说, 上式是输入特征的一个 仿射变换(affine transformation)。 仿射变换的特点是通过加权和对特征进行线性变换(linear transformation), 并通过偏置项来进行平移(translation)。

给定一个数据集,我们的目标是寻找模型的权重w和偏置b, 使得根据模型做出的预测大体符合数据里的真实价格。 输出的预测值由输入特征通过线性模型的仿射变换决定,仿射变换由所选权重和偏置确定。

而在机器学习领域,我们通常使用的是高维数据集,建模时采用线性代数表示法会比较方便。 当我们的输入包含d个特征时,我们将预测结果y^\hat{y}(通常使用“尖角”符号表示的估计值)表示为:

y^=w1x1++wdxd+b\hat{y} = w_1x_1 + \ldots + w_dx_d + b

将所有特征放到向量xRdx \sub R^d中, 并将所有权重放到向量wRdw \sub R^d中, 我们可以用点积形式来简洁地表达模型

y^=wTx+b\hat{y} = w^Tx+b

用符号表示的矩阵X可以很方便地引用我们整个数据集的n个样本。 其中,X的每一行是一个样本,每一列是一种特征。

对于特征集合X,预测值y^\hat{y} 可以通过矩阵-向量乘法表示为:

y^=Xw+b\hat{y} = Xw + b

虽然我们相信给定x预测y的最佳模型会是线性的, 但我们很难找到一个有
个样本的真实数据集,其中对于所有的1in1 \le i \le ny(i)y^{(i)}完全等于wTx(i)+bw^Tx^{(i)} + b。 无论我们使用什么手段来观察特征X和标签y, 都可能会出现少量的观测误差。 因此,即使确信特征与标签的潜在关系是线性的, 我们也会加入一个噪声项来考虑观测误差带来的影响。

在开始寻找最好的模型参数(model parameters)w和b之前, 我们还需要两个东西: (1)一种模型质量的度量方式; (2)一种能够更新模型以提高模型预测质量的方法。

1.1.2损失函数

在我们开始考虑如何用模型拟合(fit)数据之前,我们需要确定一个拟合程度的度量。 损失函数(loss function)能够量化目标的实际值与预测值之间的差距。 通常我们会选择非负数作为损失,且数值越小表示损失越小,完美预测时的损失为0。 回归问题中最常用的损失函数是平方误差函数。 当样本i的预测值为y^(i)\hat{y}^(i),其相应的真实标签为y(i)y^(i)时, 平方误差可以定义为以下公式:

l(i)(w,b)=12(y(i)y^(i))2l^{(i)}(w,b) = \frac{1}{2} (y^{(i)} - \hat{y}^{(i)})^2

常数12\frac{1}{2} 不会带来本质的差别,但这样在形式上稍微简单一些 (因为当我们对损失函数求导后常数系数为1)。 由于训练数据集并不受我们控制,所以经验误差只是关于模型参数的函数。 为了进一步说明,来看下面的例子。 我们为一维情况下的回归问题绘制图像

由于平方误差函数中的二次方项, 估计值y^(i)\hat{y}^{(i)}和观测值y(i)y^{(i)}之间较大的差异将导致更大的损失。 为了度量模型在整个数据集上的质量,我们需计算在训练集n个样本上的损失均值(也等价于求和)。

L(w,b)=1ni=1nl(i)(w,b)=1ni=1n12(wTx(i)+by(i))2L(w,b) = \frac{1}{n} \sum_{i=1}^n l^{(i)}(w,b) = \frac{1}{n} \sum_{i=1}^n\frac{1}{2}(w^Tx^{(i)} + b - y^{(i)})^2

在训练模型时,我们希望寻找一组参数(w,b)(w^*, b^*), 这组参数能最小化在所有训练样本上的总损失。如下式

w,b=argminL(w,b)w^*,b^* = argmin L(w,b)

1.1.3 解析解

线性回归刚好是一个很简单的优化问题。 与我们将在本书中所讲到的其他大部分模型不同,线性回归的解可以用一个公式简单地表达出来, 这类解叫作解析解(analytical solution)。 首先,我们将偏置b合并到参数w中,合并方法是在包含所有参数的矩阵中附加一列。 我们的预测问题是最小化yXw2||y - Xw|| ^2。 这在损失平面上只有一个临界点,这个临界点对应于整个区域的损失极小点。 将损失关于w的导数设为0,得到解析解:

w=(XTX)1XTyw* = (X^TX)^{-1}X^Ty

像线性回归这样的简单问题存在解析解,但并不是所有的问题都存在解析解。 解析解可以进行很好的数学分析,但解析解对问题的限制很严格,导致它无法广泛应用在深度学习里。

1.1.4 随机梯度下降

即使在我们无法得到解析解的情况下,我们仍然可以有效地训练模型。 在许多任务上,那些难以优化的模型效果要更好。 因此,弄清楚如何训练这些难以优化的模型是非常重要的。

本书中我们用到一种名为梯度下降(gradient descent)的方法, 这种方法几乎可以优化所有深度学习模型。 它通过不断地在损失函数递减的方向上更新参数来降低误差。

梯度下降最简单的用法是计算损失函数(数据集中所有样本的损失均值) 关于模型参数的导数(在这里也可以称为梯度)。 但实际中的执行可能会非常慢:因为在每一次更新参数之前,我们必须遍历整个数据集。 因此,我们通常会在每次需要计算更新的时候随机抽取一小批样本, 这种变体叫做小批量随机梯度下降(minibatch stochastic gradient descent)。

在每次迭代中,我们首先随机抽样一个小批量B\mathcal{B}, 它是由固定数量的训练样本组成的。 然后,我们计算小批量的平均损失关于模型参数的导数(也可以称为梯度)。 最后,我们将梯度乘以一个预先确定的正数η\eta,并从当前参数的值中减掉。

我们用下面的数学公式来表示这一更新过程

(w,b)=(w,b)ηBiB(w,b)l(i)(w,b)(w,b) = (w,b) - \frac{\eta}{|\mathcal{B}|} \sum_{i \sub \mathcal{B}} \partial_{(w,b)}l^{(i)}(w,b)

公式 (3.1.10)中的w和x都是向量。 在这里,更优雅的向量表示法比系数表示法(如w1,w2,,wdw_1,w_2,\ldots, w_d)更具可读性。B|\mathcal{B}|表示每个小批量中的样本数,这也称为批量大小(batch size)。η\eta表示学习率(learning rate)。 批量大小和学习率的值通常是手动预先指定,而不是通过模型训练得到的。 这些可以调整但不在训练过程中更新的参数称为超参数(hyperparameter)。 调参(hyperparameter tuning)是选择超参数的过程。 超参数通常是我们根据训练迭代结果来调整的, 而训练迭代结果是在独立的验证数据集(validation dataset)上评估得到的。

在训练了预先确定的若干迭代次数后(或者直到满足某些其他停止条件后), 我们记录下模型参数的估计值,表示为w^,b^\hat{w}, \hat{b}。 但是,即使我们的函数确实是线性的且无噪声,这些估计值也不会使损失函数真正地达到最小值。 因为算法会使得损失向最小值缓慢收敛,但却不能在有限的步数内非常精确地达到最小值。

线性回归恰好是一个在整个域中只有一个最小值的学习问题。 但是对像深度神经网络这样复杂的模型来说,损失平面上通常包含多个最小值。 深度学习实践者很少会去花费大力气寻找这样一组参数,使得在训练集上的损失达到最小。 事实上,更难做到的是找到一组参数,这组参数能够在我们从未见过的数据上实现较低的损失, 这一挑战被称为泛化(generalization)。

2. 从零实现线性回归

我们使用pytorch实现

2.1 生成数据集

为了简单起见,我们将根据带有噪声的线性模型构造一个人造数据集。 我们的任务是使用这个有限样本的数据集来恢复这个模型的参数。 我们将使用低维数据,这样可以很容易地将其可视化。 在下面的代码中,我们生成一个包含1000个样本的数据集, 每个样本包含从标准正态分布中采样的2个特征。 我们的合成数据集是一个矩阵XR1000×2X \sub R^{1000 \times 2}

我们使用线性模型参数w=[2,3.4]w = [2,-3.4]、b = 4.2和噪声项ϵ\epsilon 生成数据集及其标签:

y=Xw+b+ϵy = Xw + b + \epsilon

可以视为模型预测和标签时的潜在观测误差。 在这里我们认为标准假设成立,即
服从均值为0的正态分布。 为了简化问题,我们将标准差设为0.01。 下面的代码生成合成数据集。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
import random
import torch
from d2l import torch as d2l

def synthetic_data(w, b, num_examples): #@save
"""生成y=Xw+b+噪声"""
X = torch.normal(0, 1, (num_examples, len(w)))
y = torch.matmul(X, w) + b
y += torch.normal(0, 0.01, y.shape)
return X, y.reshape((-1, 1))

true_w = torch.tensor([2, -3.4])
true_b = 4.2
features, labels = synthetic_data(true_w, true_b, 1000)

我们此时可以可视化观察下生成数据的样子

1
2
d2l.set_figsize()
d2l.plt.scatter(features[:, (1)].detach().numpy(), labels.detach().numpy(), 1)

2.2 批量读取数据集

回想一下,训练模型时要对数据集进行遍历,每次抽取一小批量样本,并使用它们来更新我们的模型。 由于这个过程是训练机器学习算法的基础,所以有必要定义一个函数, 该函数能打乱数据集中的样本并以小批量方式获取数据。

在下面的代码中,我们定义一个data_iter函数, 该函数接收批量大小、特征矩阵和标签向量作为输入,生成大小为batch_size的小批量。 每个小批量包含一组特征和标签。

1
2
3
4
5
6
7
8
9
def data_iter(batch_size, features, labels):
num_examples = len(features)
indices = list(range(num_examples))
# 这些样本是随机读取的,没有特定的顺序
random.shuffle(indices)
for i in range(0, num_examples, batch_size):
batch_indices = torch.tensor(
indices[i: min(i + batch_size, num_examples)])
yield features[batch_indices], labels[batch_indices]

2.3 初始化参数模型

在我们开始用小批量随机梯度下降优化我们的模型参数之前, 我们需要先有一些参数。 在下面的代码中,我们通过从均值为0、标准差为0.01的正态分布中采样随机数来初始化权重, 并将偏置初始化为0。

1
2
w = torch.normal(0, 0.01, size=(2,1), requires_grad=True)
b = torch.zeros(1, requires_grad=True)

2.4 定义模型

接下来,我们必须定义模型,将模型的输入和参数同模型的输出关联起来。 回想一下,要计算线性模型的输出, 我们只需计算输入特征X和模型权重w的矩阵-向量乘法后加上偏置b。 注意,上面的Xw是一个向量,b是一个标量。

1
2
3
def linreg(X, w, b):  #@save
"""线性回归模型"""
return torch.matmul(X, w) + b

2.5 定义损失函数

因为需要计算损失函数的梯度,所以我们应该先定义损失函数。 这里我们使用 3.1节中描述的平方损失函数。 在实现中,我们需要将真实值y的形状转换为和预测值y_hat的形状相同。

1
2
3
def squared_loss(y_hat, y):  #@save
"""均方损失"""
return (y_hat - y.reshape(y_hat.shape)) ** 2 / 2

2.6 定义优化算法

尽管线性回归有解析解,但本书中的其他模型却没有。 这里我们介绍小批量随机梯度下降。

在每一步中,使用从数据集中随机抽取的一个小批量,然后根据参数计算损失的梯度。 接下来,朝着减少损失的方向更新我们的参数。 下面的函数实现小批量随机梯度下降更新。 该函数接受模型参数集合、学习速率和批量大小作为输入。每 一步更新的大小由学习速率lr决定。 因为我们计算的损失是一个批量样本的总和,所以我们用批量大小(batch_size) 来规范化步长,这样步长大小就不会取决于我们对批量大小的选择。

1
2
3
4
5
6
def sgd(params, lr, batch_size):  #@save
"""小批量随机梯度下降"""
with torch.no_grad():
for param in params:
param -= lr * param.grad / batch_size
param.grad.zero_()

2.7 训练

现在我们已经准备好了模型训练所有需要的要素,可以实现主要的训练过程部分了。 理解这段代码至关重要,因为从事深度学习后, 相同的训练过程几乎一遍又一遍地出现。 在每次迭代中,我们读取一小批量训练样本,并通过我们的模型来获得一组预测。 计算完损失后,我们开始反向传播,存储每个参数的梯度。 最后,我们调用优化算法sgd来更新模型参数。

在每个迭代周期(epoch)中,我们使用data_iter函数遍历整个数据集, 并将训练数据集中所有样本都使用一次(假设样本数能够被批量大小整除)。 这里的迭代周期个数num_epochs和学习率lr都是超参数,分别设为3和0.03。 设置超参数很棘手,需要通过反复试验进行调整。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
lr = 0.03
num_epochs = 3
net = linreg
loss = squared_loss

for epoch in range(num_epochs):
for X, y in data_iter(batch_size, features, labels):
l = loss(net(X, w, b), y) # X和y的小批量损失
# 因为l形状是(batch_size,1),而不是一个标量。l中的所有元素被加到一起,
# 并以此计算关于[w,b]的梯度
l.sum().backward()
sgd([w, b], lr, batch_size) # 使用参数的梯度更新参数
with torch.no_grad():
train_l = loss(net(features, w, b), labels)
print(f'epoch {epoch + 1}, loss {float(train_l.mean()):f}')