D2L自学日志02：线性回归
线性回归
线性假设是指目标可以表示为特征的加权和。
举个例子：我想做一个预测房屋售价的模型。
房屋售价相关的数量特征有：房屋面积、房龄。那么我们就可以得到这样的式子（虽然不一定有道理）：
$$ price = w_{area} \cdot area + w_{age} \cdot age + b $$其中$w_{area}$和$w_{age}$称为权重，代表这个特征对于结果的影响大小，比如在这个例子中，房屋面积和房龄对于售价的影响力可能不同；$b$代表偏置，用来拓展模型的表达能力（不加偏置模型只能进行线性变换，加了偏置模型可以实现仿射变换）。
在机器学习领域中，我们处理的对象往往是高维数据集，即特证数很可能非常多。我们可以将所有的相关特征全部塞进一个向量中。比如，我们的目标可以用$x_1,x_2,x_3,...$来描述，那么我们可以这样来表示我们的预测值$\hat{y}$：
$$ \hat{y} = \sum_{i = 1}^{n}{w_ix_i} + b $$将所有特征塞进向量$\vec{x}$后，我们可以得到：
$$ \hat{y} = w^Tx + b $$其中$w^T$表示权重向量，每个元素表示对应的特征的权重大小。这样，我们就用点积的形式简洁的表示了模型。
通过这样的模型，我们已经可以实现将一个样本的特征与其预测值实现映射。但是通常情况下，这样做并不足以解决问题，因为我们多数时候不会只有一个样本。我们理想的效果是这样：通过同一套映射规则，将多个样本的特征分别映射到各自的预测值。非常巧的是，对于线性回归模型，我们恰好有一个强大的工具可以实现这一功能：矩阵。
对于一个样本的特征$x_{1},x_{2},x_{3},...$，我们可以把它塞进一个行向量$\vec{x}^T$中，组成一个特征行向量；随后，让这些行向量从上至下排列，我们就得到了一个由特征数值组成的矩阵。
比如，对于样本a，特征行向量为$[x_{a_{1}}, x_{a_{2}}, x_{a_{3}}, ... ,x_{as}]$;对于样本b，我们有特征行向量$[x_{b_{1}},x_{b_{2}},...]$，以此类推，我们可以最后得到每一个样本的特征行向量，随后我们将它拼在一起，得到如下矩阵：
$$ X = \begin{pmatrix} x_{a_{1}}&x_{a_{2}}&x_{a_{3}}&{...}&x_{as} \\ x_{b_{1}}&x_{b_{2}}&x_{b_{3}}&{...}&x_{bs} \\ \vdots \\ \end{pmatrix} $$这样，我们就将整个数据集的所有样本的所以特征汇总在了一个矩阵中。我们同样将预测值从上至下排列为一个列向量$\hat{y}$，这样我们就得到了：
$$ \hat{y}= Xw + b $$对于给定训练数据特征$\mathbf{X}$和对应的已知标签$\mathbf{y}$，线性回归的目标是找到一组权重向量$\mathbf{w}$和偏置$b$：当给定从$\mathbf{X}$的同分布中取样的新样本特征时，这组权重向量和偏置能够使得新样本预测标签的误差尽可能小。说白了，我们认为$\hat{y}$和$X$之间可能有线性关系，而我们所做的其实是希望让机器通过学习已有的样本数据来自行归纳这种潜在的线性关系。
但是，一个很现实的问题是：以上的所有假设都建立在目标与特征真的存在线性关系这一前提，然而现实世界中我们几乎不可能找到严格线性的数据集（现实数据的浮点数误差，你可以这么理解）。
所以，既然这种误差不可避免，一个很天才的想法是：我们引入一个噪声项来观察误差带来的影响。这就像听歌的时候音质不好，所以你干脆放点噪音让音质的缺陷显得非常合理。
总而言之，我们的根本目标是寻找描述$X$和$\hat{y}$之间映射关系最好的一组$w,b$。但是在正式开始寻找之前，我们需要回答两个颇具哲理的问题：
- 什么是最好？
- 如何寻找？
他们分别对应于：
- 模型质量的度量方式。
- 一种能更新模型进而提高模型质量的方式。
度量质量：损失函数
对于第一个问题，我们希望通过我们的模型，最好输入$x$就能直接精确预测到$y$，换言之，什么样的映射能让$\hat{y}$和$y$之间差距最小，这样的映射就是我们所在寻找的。因此，我们可以通过量化模型的误差进而体现模型的质量。而这种量化的方式就是所谓的损失函数（loss funciton）。
对于一个通常的误差函数，我们希望它的输出值为非负数，这样他的大小即直接展现损失大小。在回归问题中，一个常用的损失函数是平方误差函数。我们记样本$i$的预测值为$\hat{y}^{(i)}$，其对应的真实标签为$y^{(i)}$，那么平方误差就可以定义为如下的式子：
$$ l^{(i)}(w, b)=\frac{1}{2}(\hat{y}^{(i)}-y^{(i)})^2 $$这里，之所以有一个$\frac{1}{2}$，只是因为这样做后续求导系数为1，稍微方便一点。
你会注意到这个式子的作用对象是一个标量，即它只反应了一个样本中的损失，但是我们希望能作用于一个数据集的所有样本。于是有一个很简单粗暴的方法：让我们算出每一个样本的损失，然后求平均。
$$ L(w,b) = \frac{1}{n}\sum_{i = 1}^{n}{l^{(i)}(w,b)} $$这样，我们就将我们训练模型时的目的明确为：
寻找一组参数$w^*,b^*$使得这组参数可以最小化所有训练样本中的总损失。
解析解
不同于我们后续会涉及到的大多数模型，线性回归的解可以用一个公式简单表达出来，即：在给定了一个数据集的时候，我们可以直接算出来最优的参数$w^*,b^*$，而无需复杂的训练过程。
首先，对于偏置量$b$来说，我们没必要去专门考察它在什么时候最优。我们可以对公式$\hat{y} = Xw + b$做如下调整：我们在X最左侧强行加上一列全为1的向量，然后在$w$中强行插入一个$b$，变为$\beta = [b, w_1,w_2,...]$。这样，我们可以让b的地位和其他权重保持一致。我们将问题简化为：求使得$\hat{y} = X\beta$误差最小的权重向量$\beta$。
对于解答这个问题，我们可以使用线性代数优美的几何本质。$X$是一个已知的矩阵，它可以把一个空间中的所有向量映射到另一个空间（即$X$的列空间）中，假如$y$就在$X$的列空间中，这就意味着我们可以用一个向量通过X变换严丝合缝的映射到$y$，这种情况就对应与数据集中的特征与目标呈现严格的线性关系。但是现实世界中我们找不到如此理想的数据集，这意味着$y$多数时候并不在$X$的列空间中。
此时我们回顾前面的误差函数定义，我们会发现其实我们所计算的误差是$\hat{y}$和$y$对应坐标之间的欧几里得距离（平方和，只是没开根号）。我们可以将问题简化如下：
已知$\hat{y}$被约束在$X$的列空间中，$y$不一定在$X$的列空间中，什么情况下$\hat{y}$和$y$之间的欧几里得距离会最小？
如果这还是很抽象，我们来看一个三维空间下的表述：
已知$\hat{y}$在$X$列向量张成的平面或直线上，而$y$并不一定在这个平面或直线上，什么情况下$\hat{y}$和$y$之间的欧几里得距离最小？
答案很显然：$\hat{y}$与$y$之间的差向量垂直于$X$的列空间。
我们令$e = y - X\beta$表示误差向量，我们希望$e$可以垂直于$X$的每一个列向量，这就意味着$e$与每个$X$列向量的点积都为0。也就是说：
$$ X^Te = 0 $$将$e$的表达式代入，我们就可以直接求得最佳的$\beta$：
$$ \beta = (X^TX)^{-1}X^Ty $$这样我们就一次性求得了线性回归中的最优解$w^*,b^*$。
但是，虽然像线性回归这样的简单问题存在解析解，但是对于更加复杂的问题，我们大概率找不到一个解析解。事实上，很多问题找不到解析解并不是因为难算，而是数学意义上的无解，即我们很可能无法将特征与目标的最佳映射写成一个简单的数学公式。因此，我们没办法通过求解析解来处理深度学习的所有问题。
随机梯度下降
即使我们不知道解析解长成什么样，但是有一点我们是确定的：特征与目标之间是存在映射关系的，即特征和目标的分布之间存在着某种规律。
虽然我还没有更深入地学习机器学习，但是我隐约感觉这有点像人类物理学家对物理规律的总结。比如当我们遇到一些从来没有人解释过的现象时，我们会尝试提出各种各样的猜想，并尝试依据猜想，从相关因素去预测结果。最后，我们只会留下预测最精准的猜想作为我们的科学理论。
但是，人类会改变思考的角度，机器不会。比如在尝试分析天体问题时，过去人们站在地球的视角用本轮和均轮描述天体运动，但是后来发现似乎从太阳的视角分析问题得到的结论会更简单。但是机器学习的原理就类似于无限地调整本轮均轮尺寸来让预测值和实际值靠近。我们只能输入数据来训练机器，让它构建一套精密的本轮均轮模型来进行预测，但是我们没办法直接从机器中提取简洁的规律，因为我们能从机器中得到的也只是那一套本轮均轮模型罢了。
对于几乎所有的深度学习模型，有这样一种普适的优化方法：梯度下降。它的本质是让我们在损失函数递减的方向上更新参数（即各种权重）来降低误差。
我们已经构建起了由特征值到预测值的映射，对于给定的数据集，我们已知了特征和对应的真实值，而预测值可以通过我们构建的映射得到，这就意味着在给定数据集中我们实际上构建了由特征值到误差的映射，而我们要做的事情实际就是调整这套映射来让误差尽可能小。
我们可以用一个经典的比喻来理解梯度下降的工作原理：我们将参数视作经纬度，而误差视作一座山的海拔。对于每一个样本，我们都有这样的一座山，但是每一个样本心中的山的形状并不一样。比如，假如我们的模型有两个参数$a,b$，在样本A中可能${(a = 3, b = 2)}$的时候预测的误差最小，但是对于样本B，可能$(a = 4, b =5)$的时候误差才最小。不过，尽管山的形状大概率不同，我们仍然认为不同样本的山的形状存在一定的分布规律。
对于一个数据集，我们很可能有成千上万个样本，每个样本的山都堆叠在相同的经纬度，而我们一次只能存在于一个经纬度。现在我们希望找到这样的一个经纬度能让我们在每一座山上的海拔的平均值最低。
梯度下降法认为我们可以这样实现目的：首先，我们不再去想象山的全貌，而是只去想象我们所在的经纬度的山的切面。山是一个曲面，在我们所在经纬度上与这个曲面相切的平面就是我们要找的切面。
然后，我们想象分别在每一个参数对应的方向上移动一个微小的位移（比如1），这样每个样本的误差山的切面上高度都会改变，有的高度会增大，有的高度会减小，我们将所有高度变化求平均值，就得到了这个参数对于平均误差变化的贡献值。
这个贡献值为正数时，意味着随着参数变大，平均误差值也将会变大；当这个贡献值为负数时，参数变大，平均误差值则会随之变小。贡献值的绝对值越大，平均误差值随参数变化的速度也会越快。这个贡献值在数学上实际就是误差关于参数的偏导数 。
我们求得所有参数对误差的贡献值，然后让每个参数的贡献值乘以这个参数对应方向的单位向量，最后做矢量和，我们就得到了让平均参数上升最快的参数改变方向，而它的反向则是使平均参数下降最快的参数改变方向。
在物理上，当一个标量随着空间位置而发生改变，我们称之为标量场。比如，一个中心天体附近的引力势符合一定的空间分布。在经典力学的范畴下，我们可以这样描述一个中心天体附近的引力势分布：
$$ \phi = -\frac{GM}{\sqrt{ x^2 + y^2 + z^2} } $$其中$x,y,z$是在以天体为中心的坐标轴下的三维坐标。对于在这个坐标系下坐标为$(x,y,z)$处质量为$m$的小天体，它的引力势能等于：
$$ E_p = -\frac{GMm}{\sqrt{ x^2+y^2+z^2 }} $$这时候，你会思考：沿着什么方向这个引力势能减小的速度最大呢？如果你有一定物理知识，你会知道等势面与保守力方向垂直，保守力的方向就是势下降最快的方向。回忆起高中物理知识，我们知道引力的表达式（不考虑相对论效应）为：
$$ \vec{F} = -\frac{GMm}{(x^2+y^2+z^2)^\frac{5}{2}}(x\hat{i} + y \hat{j} + z \hat{k}) $$如果你注意力惊人，你会发现这个表达式实际上还有另一种写法：
$$ \vec{F} = -(\hat{i}\frac{\partial E_p}{\partial x}+\hat{j}\frac{\partial E_p}{\partial y}+\hat{k}\frac{\partial E_p}{\partial z}) $$进一步缩写：
$$ \vec{F} = -\nabla E_p $$而我们知道，误差$L=\sum w^Tx + b=F(w_{1},w_{2},w_{3},...)$，其中$L,w_{1},w_{2},w_{3},...$均为标量；引力势$\phi=F(x_{1},x_{2},x_{3})$，其中$\phi,x_{1},x_{2},x_{3}$均为标量。一个合理的猜想是：如果我们想像引力那样沿最快改变方向修改参数来降低误差，我们只需要沿着：
$$ -\nabla L $$的方向修改参数即可。
事实上，这样做既在数学上是正确的，同时也是工程上已经成熟使用的方法。比如在机器学习中我们可以算出平均误差之后通过反向传播算出平均误差关于各个参数的梯度，继而通过让参数减去梯度向量的特定倍数来更新参数进而降低误差。这就是梯度下降的含义。
但是，在机器学习中我们往往需要调整大量的参数并计算数以万计的样本，假如每一次更新参数都必须遍历所有的样本，这样做就必然会导致训练效率极其底下。
假如我们仔细思考，不难发现，其实如果我们只是确定一个大致的方向，我们用不着遍历所有的样本。理论上我们只需要从数据集中随机抽取几个样本，这几个样本同样可以一定程度反应整个数据集的特征。这就是目前工程上使用更广泛的随机梯度下降法。如果我们更新参数时只从全数据集中抽取一个固定的小数量的样本，这就是所谓的小批量随机梯度下降法。
在每次迭代时，我们随机抽取一个固定大小的小批量样本$\mathcal{B}$，然后我们计算这个小批量$\mathcal{B}$的平均损失，并求它对于参数的梯度，进而更新参数。
我们可以用下面这一公式来表示这一过程：
$$(\mathbf{w},b) \leftarrow (\mathbf{w},b) - \frac{\eta}{|\mathcal{B}|} \sum_{i \in \mathcal{B}} \partial_{(\mathbf{w},b)} l^{(i)}(\mathbf{w},b).$$从零实现线性回归
包括：
- 数据流水线
- 模型
- 损失函数
- 小批量随机梯度下降优化器
step01：生成数据集
为了训练模型，我们需要准备训练用的数据集。因为我们将要实现的是线性回归模型，理所当然的，我们希望我们的数据在一定误差范围内符合线性关系。所以，我们生成数据集的思路是在一个线性模型的基础上添加随机噪声。
在本例中，我们构建一个这样的线性模型：一个标签，两个特征。其中参数$w = [2,-3.4]^T,b = 4.2$，此外还有一个随机噪声$\epsilon$，服从均值为0，标准差为0.01的正态分布：
$$ y = xw + b + \epsilon $$实现如下：
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))
torch.normal()：三个参数依次是平均值、标准差、形状。
torch.matmul()：顾名思义，矩阵乘法。
true_w = torch.Tensor([2,-3.4])
true_b = 4.2
features, labels = synthetic_data(true_w, true_b, 1000)
按照计划设置数据集参数。如果一切正常，features将会是一个形状为(1000, 2)的tensor，labels将会是一个形状为(1000, 1)的tensor。
我们尝试打印某次生成的数据的前10个：
print("features: ", features[0:10], "\nlabels: ", labels[0:10])
结果：
features: tensor([[ 0.6174, 0.7957],
[ 1.5489, -0.3956],
[ 0.0109, 1.5542],
[ 0.5645, 1.4945],
[ 0.7462, 0.4104],
[ 0.8710, -0.4334],
[-0.2848, 0.1450],
[-1.0261, -1.5976],
[ 1.0452, 0.4699],
[ 0.0985, -0.5582]])
labels: tensor([[ 2.7369],
[ 8.6462],
[-1.0495],
[ 0.2449],
[ 4.2774],
[ 7.4146],
[ 3.1246],
[ 7.5757],
[ 4.6824],
[ 6.2991]])
简单计算可以验证这组数据大致符合我们的需求。
step02：读取数据集
成功建立数据集之后我们需要读取数据集。在前文中我们介绍了小批量随机梯度下降，此处我们尝试来实现这一模型优化策略。我们定义一个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)])
yield features[batch_indices], labels[batch_indices]
此处使用了python独有的生成器语法。简单来说，调用这个函数我们并不直接得到小批量样本，而是得到一个迭代器。比如iter = data_iter(...)后，iter并不是一个小批量样本张量。想要从中提取一组样本，我们需要调用next(iter)，或者使用for循环遍历。这样做的好处是节省内存，做到流式内容读取，适用于处理大数据集。（ps：如果学过CS61A，你会很熟悉这个语法。）
如上的采样实现原理是：先生成一个打乱的索引列表，然后从中每batch_size个索引组成一个batch，并返回对应的特征以及标签。
BATCH_SIZE = 10
counter = 0
for X, y in data_iter(BATCH_SIZE, features, labels):
print(X, "\n", y)
counter += 1
if counter >= 4:
break
data_iter本质上一次性生成了所有的batch。在如上的检验代码中，我们设定一个batch大小为10，并从生成的batch中读取前4个。在某次读取后结果如下：
tensor([[-1.1099, 1.6687],
[ 0.3618, 0.1041],
[-0.0529, -0.2535],
[-1.2494, -0.7070],
[ 2.1251, -2.0001],
[ 1.4389, -0.8865],
[ 0.2365, -0.4293],
[-1.3268, -0.2949],
[-1.1298, -0.2571],
[ 2.1921, 0.2929]])
tensor([[-3.7069],
[ 4.5750],
[ 4.9620],
[ 4.0994],
[15.2621],
[10.0946],
[ 6.1369],
[ 2.5313],
[ 2.8041],
[ 7.5802]])
tensor([[-1.6883, 1.4908],
[ 1.0874, -0.0492],
[ 0.4867, -0.2907],
[-1.6112, -0.7551],
[-0.6087, 0.1846],
[ 0.5941, -1.2927],
[ 1.2502, -2.2800],
[-1.9777, 0.4598],
[-0.2728, 1.0429],
[ 0.2092, 1.1009]])
tensor([[-4.2361],
[ 6.5448],
[ 6.1642],
[ 3.5481],
[ 2.3488],
[ 9.7791],
[14.4587],
[-1.3148],
[ 0.1061],
[ 0.8728]])
tensor([[-0.1582, -0.5111],
[ 0.0206, -1.8523],
[-0.6847, 0.0036],
[ 1.1984, 0.6702],
[-0.5487, 1.7227],
[ 3.1342, -1.0611],
[-0.8778, 0.0670],
[-2.2774, 1.6131],
[-0.4362, 0.2129],
[ 0.4503, -0.7968]])
tensor([[ 5.6233],
[10.5495],
[ 2.8155],
[ 4.3185],
[-2.7354],
[14.0680],
[ 2.2155],
[-5.8575],
[ 2.6139],
[ 7.8082]])
tensor([[-0.1396, -0.8441],
[ 0.4309, -0.8643],
[-0.3181, -0.7736],
[ 1.1347, 0.1844],
[ 0.3201, -0.1031],
[ 0.1408, 0.5805],
[ 0.2353, -0.6993],
[ 0.6425, 1.3151],
[ 0.4438, -1.0837],
[-2.6588, -0.6255]])
tensor([[6.7791],
[7.9966],
[6.1920],
[5.8500],
[5.1947],
[2.5164],
[7.0426],
[1.0002],
[8.7780],
[0.9943]])
每两个tensor构成一组feature和label。
step03：初始化模型参数
随机创建一组初始参数作为初始值。在下面的代码中，我们通过均值为0、标准差为0.01的正态分布采样来初始化权重，并将偏置初始化为0。
w = torch.normal(0, 0.01, size=(2,1), requires_grad=True)
b = torch.zeros(1, requires_grad=True)
在初始化参数之后，我们的任务是更新这些参数，直到这些参数足够拟合我们的数据。每次更新都需要计算损失函数关于模型参数的梯度。有了这个梯度，我们就可以向减小损失的方向更新每个参数。
step04：定义模型
我们已经获得了数据集，得到了参数张量，现在我们想要计算得到预测值，所以我们需要按照我们之前规划的参数、特征与目标之间的映射关系建立模型。
def linreg(X, w, b):
"""线性回归模型"""
return torch.matmul(X, w) + b
此处使用了pytorch以及很多其他数据科学库中常见的广播机制，b是一个数字，但是前面的torch.matmul()会得到一个向量。向量+数字，结果是向量中的每个元素都加上数字。
模型的地位是我们希望产出的产品。我们希望通过训练来得到参数最优的模型。因此，训练的过程本质就是不断的让模型进行计算，每次计算后根据结果更新参数，参数达到最优后，这套嵌入了最优参数的模型就可以用作实际用途，你可以输入特征让它预测结果。
step05：定义损失函数
按照计划，我们将预测值与真实值之差的平方值的$\frac{1}{2}$作为损失值。
def squared_loss(y_hat, y):
"""均方损失"""
return (y_hat - y.reshape(y_hat.shape)) ** 2 / 2
得到的是一个损失向量，每个元素代表对应的样本在当前模型计算下的损失值。
step06：定义优化算法
尽管线性回归存在解析解，此处我们仍然使用小批量梯度下降法来尝试找到最优解。
def sgd(params, lr, batch_size):
"""小批量随机梯度下降"""
with torch.no_grad():
for param in params:
param -= lr * param.grad / batch_size
param.grad.zero_()
with torch.no_grad()明确在此代码块中进行的计算不会被记入torch计算图中。我们优化参数不免要利用一些和参数相关的计算，但是这些计算并不在模型的计算过程内，更新参数的计算过程和误差的产生毫无关联，因此我们需要这样做来与计算图隔离。
你可能会注意到，我们前面定义的损失函数并不是一个平均值，而是单纯的平方和，这就意味着损失大小会随着batch大小改变而改变。因此在更新参数的时候我们需要做一步额外的标准化操作：让param的改变值除以batch_size，这样对于相同的样本，不管batch_size大小如何，我们都可以得到一样的改变量。
最后，每次更新完一个参数，我们将参数的梯度归零，方便进行下一次训练。
step07：训练
训练的过程中我们不断地会经历正向传播、反向传播、更新参数的循环，周而复始直到完成一个数据集全部数据的训练，我们称这样一个过程为一个迭代周期。我们大致要做这样的逻辑：
- 初始化参数
- 重复以下训练直到完成所有数据集的训练：
- 计算梯度
- 更新参数
在每个迭代周期中，我们将使用data_iter遍历数据集中的每个样本。在此过程中迭代周期数量以及学习率都是所谓的超参数，它们并不是我们通过训练可以优化的，只能通过实验经验来调整优化。
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}')
因为数据集是我们自己生成的，我们知道最优参数的标准答案为$w = [2,-3.4]^T,b = 4.2$，所以我们可以据此来衡量模型训练的成效。
print(f'w的估计误差: {true_w - w.reshape(true_w.shape)}')
print(f'b的估计误差: {true_b - b}')
完整代码：
import random
import torch
from d2l import torch as d2l
def synthetic_data(w: torch.Tensor, b: float, num_examples: int):
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)
# print("features: ", features[0:10], "\nlabels: ", labels[0:10])
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]
w = torch.normal(0, 0.01, size=(2, 1), requires_grad=True)
b = torch.zeros(1, requires_grad=True)
def linreg(X, w, b):
"""线性回归模型"""
return torch.matmul(X, w) + b
def squared_loss(y_hat, y):
"""均方损失"""
return (y_hat - y.reshape(y_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 = 3
net = linreg
loss = squared_loss
batch_size = 10
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}")
print(f"w的估计误差: {true_w - w.reshape(true_w.shape)}")
print(f"b的估计误差: {true_b - b}")
某次输出：
epoch 1, loss 0.036131
epoch 2, loss 0.000128
epoch 3, loss 0.000046
w的估计误差: tensor([-0.0004, -0.0010], grad_fn=<SubBackward0>)
b的估计误差: tensor([0.0002], grad_fn=<RsubBackward1>)
训练了3个迭代周期，最后得到的参数与最优参数之间误差在e-04数量级。
总结来说，使用梯度下降法无法得到最优参数的真实值，我们也并不关心如何恢复真实参数。我们真正希望实现的是高度准确的预测参数。