chapter_linear-networks：线性回归的从零开始实现

线性回归的从零开始实现

🏷sec_linear_scratch

在了解线性回归的关键思想之后，我们可以开始通过代码来动手实现线性回归了。
在这一节中，(我们将从零开始实现整个方法，
包括数据流水线、模型、损失函数和小批量随机梯度下降优化器)。
虽然现代的深度学习框架几乎可以自动化地进行所有这些工作，但从零开始实现可以确保你真正知道自己在做什么。
同时，了解更细致的工作原理将方便我们自定义模型、自定义层或自定义损失函数。
在这一节中，我们将只使用张量和自动求导。
在之后的章节中，我们会充分利用深度学习框架的优势，介绍更简洁的实现方式。

%matplotlib inline
import random
import torch
from d2l import torch as d2l

X = torch.normal(0, 1, (1000, 2))
X

tensor([[ 0.3820, -0.3786],
        [-1.0111, -1.2613],
        [-0.8860, -0.1149],
        ...,
        [ 0.2313,  0.7420],
        [ 0.2294,  0.7855],
        [-0.6714,  2.0202]])

生成数据集

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

(**我们使用线性模型参数$\mathbf{w}=[2,-3.4{]}^{\top }$、$b=4.2$
和噪声项$ϵ$生成数据集及其标签：

$$\mathbf{y}=\mathbf{X}\mathbf{w}+b+ϵ.$$
**)

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

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)

注意，[features中的每一行都包含一个二维数据样本，
labels中的每一行都包含一维标签值（一个标量）]。

print('features:', features[0],'\nlabel:', labels[0])

features: tensor([0.0297, 0.0525]) 
label: tensor([4.0805])

通过生成第二个特征features[:, 1]和labels的散点图，
可以直观观察到两者之间的线性关系。

d2l.set_figsize()
d2l.plt.scatter(features[:, (1)].detach().numpy(), labels.detach().numpy(), 1); # detach()分离计算

读取数据集

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

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

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)]) # min是为了防止超出index
        yield features[batch_indices], labels[batch_indices] # yield是返回迭代器(单步与return类似)，顺序迭代返回，直到数据没有迭代结束

通常，我们利用GPU并行运算的优势，处理合理大小的“小批量”。
每个样本都可以并行地进行模型计算，且每个样本损失函数的梯度也可以被并行计算。
GPU可以在处理几百个样本时，所花费的时间不比处理一个样本时多太多。

我们直观感受一下小批量运算：读取第一个小批量数据样本并打印。
每个批量的特征维度显示批量大小和输入特征数。
同样的，批量的标签形状与batch_size相等。

batch_size = 10

for X, y in data_iter(batch_size, features, labels):
    print(X, '\n', y)
    break

tensor([[ 0.6194, -0.1824],
        [ 1.0520,  0.0790],
        [ 0.3021, -0.6188],
        [ 0.3905,  0.6211],
        [-0.6219,  1.4667],
        [ 1.4290,  1.3414],
        [-0.6279,  0.9669],
        [-1.6481,  0.9980],
        [ 1.1159,  2.4694],
        [ 0.4971, -1.2591]]) 
 tensor([[ 6.0515],
        [ 6.0413],
        [ 6.9232],
        [ 2.8701],
        [-2.0101],
        [ 2.5135],
        [-0.3240],
        [-2.4966],
        [-1.9816],
        [ 9.4869]])

当我们运行迭代时，我们会连续地获得不同的小批量，直至遍历完整个数据集。
上面实现的迭代对于教学来说很好，但它的执行效率很低，可能会在实际问题上陷入麻烦。
例如，它要求我们将所有数据加载到内存中，并执行大量的随机内存访问。
在深度学习框架中实现的内置迭代器效率要高得多，
它可以处理存储在文件中的数据和数据流提供的数据。

初始化模型参数

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

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

在初始化参数之后，我们的任务是更新这些参数，直到这些参数足够拟合我们的数据。
每次更新都需要计算损失函数关于模型参数的梯度。
有了这个梯度，我们就可以向减小损失的方向更新每个参数。
因为手动计算梯度很枯燥而且容易出错，所以没有人会手动计算梯度。
我们使用 :numref:sec_autograd中引入的自动微分来计算梯度。

定义模型

接下来，我们必须[定义模型，将模型的输入和参数同模型的输出关联起来。]
回想一下，要计算线性模型的输出，
我们只需计算输入特征$\mathbf{X}$和模型权重$\mathbf{w}$的矩阵-向量乘法后加上偏置$b$。
注意，上面的$\mathbf{X}\mathbf{w}$是一个向量，而$b$是一个标量。
回想一下 :numref:subsec_broadcasting中描述的广播机制：
当我们用一个向量加一个标量时，标量会被加到向量的每个分量上。

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

[定义损失函数]

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

def squared_loss(y_hat, y):  #@save
    """均方损失"""
    return (y_hat - y.reshape(y_hat.shape)) ** 2 / 2 # 返回向量，没有相加也没有均值

(定义优化算法)

正如我们在 :numref:sec_linear_regression中讨论的，线性回归有解析解。
尽管线性回归有解析解，但本书中的其他模型却没有。
这里我们介绍小批量随机梯度下降。

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

def sgd(params, lr, batch_size):  #@save
    """小批量随机梯度下降"""
    with torch.no_grad(): # 不需要计算梯度，因为更新的时候不参与梯度计算
        for param in params:
            param -= lr * param.grad / batch_size # batch_size作平均loss，小批量中数据样本的平均误差的梯度
            param.grad.zero_()

训练

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

概括一下，我们将执行以下循环：

• 初始化参数
• 重复以下训练，直到完成
  • 计算梯度$\mathbf{g}\leftarrow {\partial }_{(\mathbf{w},b)}\frac{1}{|\mathcal{B}|}{\sum }_{i\in \mathcal{B}}l({\mathbf{x}}^{(i)},y^{(i)},\mathbf{w},b)$
  • 更新参数$(\mathbf{w},b)\leftarrow (\mathbf{w},b)-\eta \mathbf{g}$

在每个迭代周期（epoch）中，我们使用data_iter函数遍历整个数据集，
并将训练数据集中所有样本都使用一次（假设样本数能够被批量大小整除）。
这里的迭代周期个数num_epochs和学习率lr都是超参数，分别设为3和0.03。
设置超参数很棘手，需要通过反复试验进行调整。
我们现在忽略这些细节，以后会在 :numref:chap_optimization中详细介绍。

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的小批量损失，在sgd中做平均误差
        # 因为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}') # 输出平均误差

epoch 1, loss 0.050640
epoch 2, loss 0.000224
epoch 3, loss 0.000051

因为我们使用的是自己合成的数据集，所以我们知道真正的参数是什么。
因此，我们可以通过[比较真实参数和通过训练学到的参数来评估训练的成功程度]。
事实上，真实参数和通过训练学到的参数确实非常接近。

print(f'w的估计误差: {true_w - w.reshape(true_w.shape)}') # SubBackward0：相减
print(f'b的估计误差: {true_b - b}')

w的估计误差: tensor([ 0.0004, -0.0008], grad_fn=<SubBackward0>)
b的估计误差: tensor([0.0006], grad_fn=<RsubBackward1>)

# 练习2
import torch
import random
from d2l import torch as d2l
#生成数据集
def synthetic_data(r, b, num_examples):
    I = torch.normal(0, 1, (num_examples, len(r)))
    u = torch.matmul(I, r) + b
    u += torch.normal(0, 0.01, u.shape) # 噪声
    return I, u.reshape((-1, 1)) # 标量转换为向量

true_r = torch.tensor([20.0])
true_b = 0.01
features, labels = synthetic_data(true_r, true_b, 1000)

#读取数据集
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]
        
batch_size = 10
# 初始化权重
r = torch.normal(0,0.01,size = ((1,1)), requires_grad = True)
b = torch.zeros(1, requires_grad = True)

# 定义模型
def linreg(I, r, b):
    return torch.matmul(I, r) + b
# 损失函数
def square_loss(u_hat, u):
    return (u_hat - u.reshape(u_hat.shape)) ** 2/2
# 优化算法
def sgd(params, lr, batch_size):
    with torch.no_grad():
        for param in params:
            param -= lr * param.grad/batch_size
            param.grad.zero_()

lr = 0.03
num_epochs = 10
net = linreg
loss = square_loss

for epoch in range(num_epochs):
    for I, u in data_iter(batch_size, features, labels):
        l = loss(net(I, r, b), u)
        l.sum().backward()
        sgd([r, b], lr, batch_size)
    with torch.no_grad():
        train_l = loss(net(features, r, b), labels)
        print(f'epoch {epoch + 1}, loss {float(train_l.mean()):f}')
print(r)
print(b)
print(f'r的估计误差: {true_r - r.reshape(true_r.shape)}')
print(f'b的估计误差: {true_b - b}')

epoch 1, loss 0.396170
epoch 2, loss 0.000816
epoch 3, loss 0.000053
epoch 4, loss 0.000051
epoch 5, loss 0.000051
epoch 6, loss 0.000051
epoch 7, loss 0.000051
epoch 8, loss 0.000051
epoch 9, loss 0.000051
epoch 10, loss 0.000051
tensor([[20.0004]], requires_grad=True)
tensor([0.0103], requires_grad=True)
r的估计误差: tensor([-0.0004], grad_fn=<SubBackward0>)
b的估计误差: tensor([-0.0003], grad_fn=<RsubBackward1>)

注意，我们不应该想当然地认为我们能够完美地求解参数。
在机器学习中，我们通常不太关心恢复真正的参数，而更关心如何高度准确预测参数。
幸运的是，即使是在复杂的优化问题上，随机梯度下降通常也能找到非常好的解。
其中一个原因是，在深度网络中存在许多参数组合能够实现高度精确的预测。

小结

• 我们学习了深度网络是如何实现和优化的。在这一过程中只使用张量和自动微分，不需要定义层或复杂的优化器。
• 这一节只触及到了表面知识。在下面的部分中，我们将基于刚刚介绍的概念描述其他模型，并学习如何更简洁地实现其他模型。

练习

1. 如果我们将权重初始化为零，会发生什么。算法仍然有效吗？
   答：将权重初始化为零，算法依然有效。但网络层数加深后，在全连接的情况下，由于权重的对称性会导致出现隐藏神经元的对称性，由于多个值相同的隐藏神经元就好比一个神经元，在反向传播时，虽然有多个参数但实际效果是只更新一个参数，严重影响算法效果。
2. 假设你是乔治·西蒙·欧姆，试图为电压和电流的关系建立一个模型。你能使用自动微分来学习模型的参数吗?
   答：数学公式：U = R*I。步骤与本章从零开始实现一样。
3. 您能基于普朗克定律使用光谱能量密度来确定物体的温度吗？
   答：略。
4. 如果你想计算二阶导数可能会遇到什么问题？你会如何解决这些问题？
   答：二阶导数可能不存在，或无法得到显式的一阶导数，可以在求一阶导数是通过retain_graph=True参数保存计算图，进而求二阶导。
5. 为什么在squared_loss函数中需要使用reshape函数？
   答：y与y_hat的维度可能不一样，可能分别用行向量和列向量表示。我们需要将真实值y的形状转换为和预测值y_hat的形状相同。
6. 尝试使用不同的学习率，观察损失函数值下降的快慢。
   答：学习率越小收敛越慢，学习率越大收敛越快，但学习率过大可能会导致左右横跳直到NAN无法收敛。
7. 如果样本个数不能被批量大小整除，data_iter函数的行为会有什么变化？
   答：①最后一批未整除样本单独作为一个批次
   ②最后一批未整除样本扔掉处理
   ③从下一个epoch借一部分合成一个完整的批次。

Discussions