读入数据
通过如下代码读入数据,了解下波士顿房价的数据集结构,数据存放在本地目录下housing.data文件中。
# 导入需要用到的package
import numpy as np
import json
# 读入训练数据
datafile = './work/housing.data'
data = np.fromfile(datafile, sep=' ')
data
array([6.320e-03, 1.800e+01, 2.310e+00, ..., 3.969e+02, 7.880e+00, 1.190e+01])
数据形状变换
由于读入的原始数据是1维的,所有数据都连在一起。因此需要我们将数据的形状进行变换,形成一个2维的矩阵,每行为一个数据样本(14个值),每个数据样本包含13个XXX(影响房价的特征)和一个YYY(该类型房屋的均价)。
# 读入之后的数据被转化成1维array,其中array的第0-13项是第一条数据,第14-27项是第二条数据,以此类推....
# 这里对原始数据做reshape,变成N x 14的形式
feature_names = [ 'CRIM', 'ZN', 'INDUS', 'CHAS', 'NOX', 'RM', 'AGE','DIS',
'RAD', 'TAX', 'PTRATIO', 'B', 'LSTAT', 'MEDV' ]
feature_num = len(feature_names)
data = data.reshape([data.shape[0] // feature_num, feature_num])
# 查看数据
x = data[0]
print(x.shape)
print(x)
(14,) [6.320e-03 1.800e+01 2.310e+00 0.000e+00 5.380e-01 6.575e+00 6.520e+01 4.090e+00 1.000e+00 2.960e+02 1.530e+01 3.969e+02 4.980e+00 2.400e+01]
在本案例中,我们将80%的数据用作训练集,20%用作测试集,实现代码如下。通过打印训练集的形状,可以发现共有404个样本,每个样本含有13个特征和1个预测值。
ratio = 0.8
offset = int(data.shape[0] * ratio)
training_data = data[:offset]
training_data.shape
(404, 14)
数据归一化处理
对每个特征进行归一化处理,使得每个特征的取值缩放到0~1之间。这样做有两个好处:一是模型训练更高效;二是特征前的权重大小可以代表该变量对预测结果的贡献度(因为每个特征值本身的范围相同)。
# 计算train数据集的最大值,最小值,平均值
maximums, minimums, avgs = \
training_data.max(axis=0), \
training_data.min(axis=0), \
training_data.sum(axis=0) / training_data.shape[0]
# 对数据进行归一化处理
for i in range(feature_num):
#print(maximums[i], minimums[i], avgs[i])
data[:, i] = (data[:, i] - minimums[i]) / (maximums[i] - minimums[i])
封装成load data函数
将上述几个数据处理操作封装成load data函数,以便下一步模型的调用,实现方法如下。
def load_data():
# 从文件导入数据
datafile = './work/housing.data'
data = np.fromfile(datafile, sep=' ')
# 每条数据包括14项,其中前面13项是影响因素,第14项是相应的房屋价格中位数
feature_names = [ 'CRIM', 'ZN', 'INDUS', 'CHAS', 'NOX', 'RM', 'AGE', \
'DIS', 'RAD', 'TAX', 'PTRATIO', 'B', 'LSTAT', 'MEDV' ]
feature_num = len(feature_names)
# 将原始数据进行Reshape,变成[N, 14]这样的形状
data = data.reshape([data.shape[0] // feature_num, feature_num])
# 将原数据集拆分成训练集和测试集
# 这里使用80%的数据做训练,20%的数据做测试
# 测试集和训练集必须是没有交集的
ratio = 0.8
offset = int(data.shape[0] * ratio)
training_data = data[:offset]
# 计算训练集的最大值,最小值,平均值
maximums, minimums, avgs = training_data.max(axis=0), training_data.min(axis=0), \
training_data.sum(axis=0) / training_data.shape[0]
# 对数据进行归一化处理
for i in range(feature_num):
#print(maximums[i], minimums[i], avgs[i])
data[:, i] = (data[:, i] - minimums[i]) / (maximums[i] - minimums[i])
# 训练集和测试集的划分比例
training_data = data[:offset]
test_data = data[offset:]
return training_data, test_data
# 获取数据
training_data, test_data = load_data()
x = training_data[:, :-1]
y = training_data[:, -1:]
# 查看数据
print(x[0])
print(y[0])
[0. 0.18 0.07344184 0. 0.31481481 0.57750527 0.64160659 0.26920314 0. 0.22755741 0.28723404 1. 0.08967991] [0.42222222]
模型设计
模型设计是深度学习模型关键要素之一,也称为网络结构设计,相当于模型的假设空间,即实现模型“前向计算”(从输入到输出)的过程。
如果将输入特征和输出预测值均以向量表示,输入特征xxx有13个分量,yyy有1个分量,那么参数权重的形状(shape)是13×113\times113×1。假设我们以如下任意数字赋值参数做初始化:
w=[0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,−0.1,−0.2,−0.3,−0.4,0.0]w=[0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, -0.1, -0.2, -0.3, -0.4, 0.0]w=[0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,−0.1,−0.2,−0.3,−0.4,0.0]
w = [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, -0.1, -0.2, -0.3, -0.4, 0.0]
w = np.array(w).reshape([13, 1])
取出第1条样本数据,观察样本的特征向量与参数向量相乘的结果。
x1=x[0]
t = np.dot(x1, w)
print(t)
[0.69474855]
完整的线性回归公式,还需要初始化偏移量bbb,同样随意赋初值-0.2。那么,线性回归模型的完整输出是z=t+bz=t+bz=t+b,这个从特征和参数计算输出值的过程称为“前向计算”。
b = -0.2
z = t + b
print(z)
[0.49474855]
将上述计算预测输出的过程以“类和对象”的方式来描述,类成员变量有参数www和bbb。通过写一个forward函数(代表“前向计算”)完成上述从特征和参数到输出预测值的计算过程,代码如下所示。
class Network(object):
def __init__(self, num_of_weights):
# 随机产生w的初始值
# 为了保持程序每次运行结果的一致性,
# 此处设置固定的随机数种子
np.random.seed(0)
self.w = np.random.randn(num_of_weights, 1)
self.b = 0.
def forward(self, x):
z = np.dot(x, self.w) + self.b
return z
基于Network类的定义,模型的计算过程如下所示。
net = Network(13)
x1 = x[0]
y1 = y[0]
z = net.forward(x1)
print(z)
[2.39362982]
训练配置
模型设计完成后,需要通过训练配置寻找模型的最优值,即通过损失函数来衡量模型的好坏。训练配置也是深度学习模型关键要素之一。
通过模型计算x1x_1x1表示的影响因素所对应的房价应该是zzz, 但实际数据告诉我们房价是yyy。这时我们需要有某种指标来衡量预测值zzz跟真实值yyy之间的差距。对于回归问题,最常采用的衡量方法是使用均方误差作为评价模型好坏的指标,具体定义如下:
Loss=(y−z)2Loss = (y - z)^2Loss=(y−z)2
上式中的LossLossLoss(简记为: LLL)通常也被称作损失函数,它是衡量模型好坏的指标。在回归问题中,均方误差是一种比较常见的形式,分类问题中通常会采用交叉熵作为损失函数,在后续的章节中会更详细的介绍。对一个样本计算损失函数值的实现如下:
Loss = (y1 - z)*(y1 - z)
print(Loss)
[3.88644793]
因为计算损失函数时需要把每个样本的损失函数值都考虑到,所以我们需要对单个样本的损失函数进行求和,并除以样本总数NNN。
Loss=1N∑i=1N(yi−zi)2Loss= \frac{1}{N}\sum_{i=1}^N{(y_i - z_i)^2}Loss=N1i=1∑N(yi−zi)2
在Network类下面添加损失函数的计算过程如下:
class Network(object):
def __init__(self, num_of_weights):
# 随机产生w的初始值
# 为了保持程序每次运行结果的一致性,此处设置固定的随机数种子
np.random.seed(0)
self.w = np.random.randn(num_of_weights, 1)
self.b = 0.
def forward(self, x):
z = np.dot(x, self.w) + self.b
return z
def loss(self, z, y):
error = z - y
cost = error * error
cost = np.mean(cost)
return cost
使用定义的Network类,可以方便的计算预测值和损失函数。需要注意的是,类中的变量xxx, www,bbb, zzz, errorerrorerror等均是向量。以变量xxx为例,共有两个维度,一个代表特征数量(值为13),一个代表样本数量,代码如下所示。
net = Network(13)
# 此处可以一次性计算多个样本的预测值和损失函数
x1 = x[0:3]
y1 = y[0:3]
z = net.forward(x1)
print('predict: ', z)
loss = net.loss(z, y1)
print('loss:', loss)
predict: [[2.39362982] [2.46752393] [2.02483479]] loss: 3.384496992612791
训练过程
上述计算过程描述了如何构建神经网络,通过神经网络完成预测值和损失函数的计算。接下来介绍如何求解参数www和bbb的数值,这个过程也称为模型训练过程。训练过程是深度学习模型的关键要素之一,其目标是让定义的损失函数LossLossLoss尽可能的小,也就是说找到一个参数解www和bbb,使得损失函数取得极小值。
我们先做一个小测试:如 图5 所示,基于微积分知识,求一条曲线在某个点的斜率等于函数在该点的导数值。那么大家思考下,当处于曲线的极值点时,该点的斜率是多少?
图5:曲线斜率等于导数值
这个问题并不难回答,处于曲线极值点时的斜率为0,即函数在极值点的导数为0。那么,让损失函数取极小值的www和bbb应该是下述方程组的解:
∂L∂w=0\frac{\partial{L}}{\partial{w}}=0∂w∂L=0
∂L∂b=0\frac{\partial{L}}{\partial{b}}=0∂b∂L=0
将样本数据(x,y)(x, y)(x,y)带入上面的方程组中即可求解出www和bbb的值,但是这种方法只对线性回归这样简单的任务有效。如果模型中含有非线性变换,或者损失函数不是均方差这种简单的形式,则很难通过上式求解。为了解决这个问题,下面我们将引入更加普适的数值求解方法:梯度下降法。
梯度下降法
在现实中存在大量的函数正向求解容易,但反向求解较难,被称为单向函数,这种函数在密码学中有大量的应用。密码锁的特点是可以迅速判断一个密钥是否是正确的(已知xxx,求yyy很容易),但是即使获取到密码锁系统,无法破解出正确的密钥是什么(已知yyy,求xxx很难)。
这种情况特别类似于一位想从山峰走到坡谷的盲人,他看不见坡谷在哪(无法逆向求解出LossLossLoss导数为0时的参数值),但可以伸脚探索身边的坡度(当前点的导数值,也称为梯度)。那么,求解Loss函数最小值可以这样实现:从当前的参数取值,一步步的按照下坡的方向下降,直到走到最低点。这种方法笔者称它为“盲人下坡法”。哦不,有个更正式的说法“梯度下降法”。
训练的关键是找到一组(w,b)(w, b)(w,b),使得损失函数LLL取极小值。我们先看一下损失函数LLL只随两个参数w5w_5w5、w9w_9w9变化时的简单情形,启发下寻解的思路。
L=L(w5,w9)L=L(w_5, w_9)L=L(w5,w9)
这里我们将w0,w1,...,w12w_0, w_1, ..., w_{12}w0,w1,...,w12中除w5,w9w_5, w_9w5,w9之外的参数和bbb都固定下来,可以用图画出L(w5,w9)L(w_5, w_9)L(w5,w9)的形式。
net = Network(13)
losses = []
#只画出参数w5和w9在区间[-160, 160]的曲线部分,以及包含损失函数的极值
w5 = np.arange(-160.0, 160.0, 1.0)
w9 = np.arange(-160.0, 160.0, 1.0)
losses = np.zeros([len(w5), len(w9)])
#计算设定区域内每个参数取值所对应的Loss
for i in range(len(w5)):
for j in range(len(w9)):
net.w[5] = w5[i]
net.w[9] = w9[j]
z = net.forward(x)
loss = net.loss(z, y)
losses[i, j] = loss
#使用matplotlib将两个变量和对应的Loss作3D图
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
fig = plt.figure()
ax = Axes3D(fig)
w5, w9 = np.meshgrid(w5, w9)
ax.plot_surface(w5, w9, losses, rstride=1, cstride=1, cmap='rainbow')
plt.show()
<Figure size 640x480 with 1 Axes>
对于这种简单情形,我们利用上面的程序,可以在三维空间中画出损失函数随参数变化的曲面图。从图中可以看出有些区域的函数值明显比周围的点小。
需要说明的是:为什么这里我们选择w5w_5w5和w9w_9w9来画图?这是因为选择这两个参数的时候,可比较直观的从损失函数的曲面图上发现极值点的存在。其他参数组合,从图形上观测损失函数的极值点不够直观。
观察上述曲线呈现出“圆滑”的坡度,这正是我们选择以均方误差作为损失函数的原因之一。图6 呈现了只有一个参数维度时,均方误差和绝对值误差(只将每个样本的误差累加,不做平方处理)的损失函数曲线图。
图6:均方误差和绝对值误差损失函数曲线图
由此可见,均方误差表现的“圆滑”的坡度有两个好处:
- 曲线的最低点是可导的。
- 越接近最低点,曲线的坡度逐渐放缓,有助于通过当前的梯度来判断接近最低点的程度(是否逐渐减少步长,以免错过最低点)。
而绝对值误差是不具备这两个特性的,这也是损失函数的设计不仅仅要考虑“合理性”,还要追求“易解性”的原因。
现在我们要找出一组[w5,w9][w_5, w_9][w5,w9]的值,使得损失函数最小,实现梯度下降法的方案如下:
- 步骤1:随机的选一组初始值,例如:[w5,w9]=[−100.0,−100.0][w_5, w_9] = [-100.0, -100.0][w5,w9]=[−100.0,−100.0]
- 步骤2:选取下一个点[w5′,w9′][w_5^{'} , w_9^{'}][w5′,w9′],使得L(w5′,w9′)<L(w5,w9)L(w_5^{'} , w_9^{'}) < L(w_5, w_9)L(w5′,w9′)<L(w5,w9)
- 步骤3:重复步骤2,直到损失函数几乎不再下降。
如何选择[w5′,w9′][w_5^{'} , w_9^{'}][w5′,w9′]是至关重要的,第一要保证LLL是下降的,第二要使得下降的趋势尽可能的快。微积分的基础知识告诉我们,沿着梯度的反方向,是函数值下降最快的方向,如 图7 所示。简单理解,函数在某一个点的梯度方向是曲线斜率最大的方向,但梯度方向是向上的,所以下降最快的是梯度的反方向。
图7:梯度下降方向示意图
计算梯度
上面我们讲过了损失函数的计算方法,这里稍微改写,为了使梯度计算更加简洁,引入因子12\frac{1}{2}21,定义损失函数如下:
L=12N∑i=1N(yi−zi)2L= \frac{1}{2N}\sum_{i=1}^N{(y_i - z_i)^2}L=2N1i=1∑N(yi−zi)2
其中ziz_izi是网络对第iii个样本的预测值:
zi=∑j=012xij⋅wj+bz_i = \sum_{j=0}^{12}{x_i^{j}\cdot w_j} + bzi=j=0∑12xij⋅wj+b
梯度的定义:
𝑔𝑟𝑎𝑑𝑖𝑒𝑛𝑡=(∂L∂w0,∂L∂w1,...,∂L∂w12,∂L∂b)𝑔𝑟𝑎𝑑𝑖𝑒𝑛𝑡 = (\frac{\partial{L}}{\partial{w_0}},\frac{\partial{L}}{\partial{w_1}}, ... ,\frac{\partial{L}}{\partial{w_{12}}} ,\frac{\partial{L}}{\partial{b}})gradient=(∂w0∂L,∂w1∂L,...,∂w12∂L,∂b∂L)
可以计算出LLL对www和bbb的偏导数:
∂L∂wj=1N∑i=1N(zi−yi)∂zi∂wj=1N∑i=1N(zi−yi)xij\frac{\partial{L}}{\partial{w_j}} = \frac{1}{N}\sum_{i=1}^N{(z_i - y_i)\frac{\partial{z_i}}{\partial{w_j}}} = \frac{1}{N}\sum_{i=1}^N{(z_i - y_i)x_i^{j}}∂wj∂L=N1i=1∑N(zi−yi)∂wj∂zi=N1i=1∑N(zi−yi)xij
∂L∂b=1N∑i=1N(zi−yi)∂zi∂b=1N∑i=1N(zi−yi)\frac{\partial{L}}{\partial{b}} = \frac{1}{N}\sum_{i=1}^N{(z_i - y_i)\frac{\partial{z_i}}{\partial{b}}} = \frac{1}{N}\sum_{i=1}^N{(z_i - y_i)}∂b∂L=N1i=1∑N(zi−yi)∂b∂zi=N1i=1∑N(zi−yi)
从导数的计算过程可以看出,因子12\frac{1}{2}21被消掉了,这是因为二次函数求导的时候会产生因子222,这也是我们将损失函数改写的原因。
下面我们考虑只有一个样本的情况下,计算梯度:
L=12(yi−zi)2L= \frac{1}{2}{(y_i - z_i)^2}L=21(yi−zi)2
z1=x10⋅w0+x11⋅w1+...+x112⋅w12+bz_1 = {x_1^{0}\cdot w_0} + {x_1^{1}\cdot w_1} + ... + {x_1^{12}\cdot w_{12}} + bz1=x10⋅w0+x11⋅w1+...+x112⋅w12+b
可以计算出:
L=12(x10⋅w0+x11⋅w1+...+x112⋅w12+b−y1)2L= \frac{1}{2}{({x_1^{0}\cdot w_0} + {x_1^{1}\cdot w_1} + ... + {x_1^{12}\cdot w_{12}} + b - y_1)^2}L=21(x10⋅w0+x11⋅w1+...+x112⋅w12+b−y1)2
可以计算出LLL对www和bbb的偏导数:
∂L∂w0=(x10⋅w0+x11⋅w1+...+x112⋅w12+b−y1)⋅x10=(z1−y1)⋅x10\frac{\partial{L}}{\partial{w_0}} = ({x_1^{0}\cdot w_0} + {x_1^{1}\cdot w_1} + ... + {x_1^{12}\cdot w_12} + b - y_1)\cdot x_1^{0}=({z_1} - {y_1})\cdot x_1^{0}∂w0∂L=(x10⋅w0+x11⋅w1+...+x112⋅w12+b−y1)⋅x10=(z1−y1