机器学习Hello World：波士顿房价预测（线性回归模型）

虽然之前接触过一些深度学习的知识，了解了一些理论知识，但是其实做出来的东西都是一些浅层的应用（只会把现成的模型进行迁移学习），而且很多时候都是知其然不知其所以然。现在决定好好学一学。

波士顿房价预测模型可以说是机器学习的hello world教程，使用了线性回归模型。

下面的表格列出了影响波士顿地区的房价的因素以及房价的中位数。

数据是一个仅仅以空格进行分隔的文件，打开来可以看到长这个样子：

那么我们需要构建一个单层的神经网络来解决这个问题：

在解决这个问题的时候，我们假设房价中位数与各因素之间的关系可以用线性关系表示：

模型的求解就是通过给定的数据，拟合出每个w_j和b。其中，w_j和b分别表示该线性模型的权重和偏置。一维情况下，w_j和b 是直线的斜率和截距。

我们使用均方误差（MSE）作为损失函数来表示预测值和真实值之间的差距。

下面讲解几个概念：

数据集拆分：对于给定的数据，我们需要把一些用于训练，叫做训练集。另一些用于验证，称为测试集。

数据归一化：我们需要把所有的特征值都缩放到0~1之间，因为这样才能在后期进行梯度下降的时候，统一步长。

这样做有两个好处：一是模型训练更高效；二是特征前的权重大小可以代表该变量对预测结果的贡献度（因为每个特征值本身的范围相同）

损失函数：作为评价预测值和真实值之间差距的标准，我们需要让损失函数尽可能地小。

梯度下降法：求损失函数的梯度，然后反向传播，沿着梯度的反方向，更新权重w和偏置b，使得各神经元的权重能够使得损失函数更小。

随机梯度下降法：和梯度下降法类似，只是在每一轮训练(epoch)前，把训练集进行乱序处理，并且使用mini_batch进行训练，从而避免靠近训练集末尾的数据对训练结果影响过大的问题，并且能够提高训练速度。

为什么要使用均方误差作为损失函数，而不是绝对值误差？

由此可见，均方误差表现的“圆滑”的坡度有两个好处：

曲线的最低点是可导的。
越接近最低点，曲线的坡度逐渐放缓，有助于通过当前的梯度来判断接近最低点的程度（是否逐渐减少步长，以免错过最低点）。

而绝对值误差是不具备这两个特性的，这也是损失函数的设计不仅仅要考虑“合理性”，还要追求“易解性”的原因。

求解梯度

为了使得计算过程更简洁易懂，我们改写损失函数，改成：

其中，z_i是网络对第i个样本的预测值

求出损失函数对w0的偏导数：

对于其他权重，同理。从而能得出L对于w和b的梯度。

每一步训练，我们都需要将w和b向梯度的反方向移动一小段距离。经过大量训练，就能使得loss降低。

代码：

import numpy as np
import json


# 读入数据的函数
def load_data():
    data_file = 'housing.data'
    # 从文件读入数据，并指定分隔符为空格
    data = np.fromfile(data_file, sep=' ')
    # 此时data.shape为(7084,)

    # 每条数据包含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])
    # 于是，data.shape变成了(506, 14)

    # 将数据集拆分为训练集和测试集
    # 这里使用80%为训练集， 20%为测试集
    # 训练集和测试集必须没有交集

    ratio = 0.8
    offset = int(data.shape[0] * ratio)
    training_data = data[: offset]


    # 计算训练集的最大值、最小值、平均值     形状为（14,）
    maximums = training_data.max(axis=0)
    minimums = training_data.min(axis=0)
    avgs = training_data.sum(axis=0) / training_data.shape[0]

    # 对数据进行归一化
    for i in range(feature_num):
        data[:, i] = (data[:, i] - minimums[i]) / (maximums[i] - minimums[i])

    training_data = data[: offset]
    test_data = data[offset:]
    return training_data, test_data


class NetWork(object):
    def __init__(self, num_of_weights):
        # 随机产生w的初始值
        # 线性回归模型
        # w的形状是(13, 1)
        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
        num_samples = error.shape[0]
        cost = error * error
        # 把所有样本的cost相加，求平均
        cost = np.sum(cost)/num_samples
        return cost

    # 梯度
    def gradient(self, x, y):
        z = self.forward(x)
        gradient_w = (z-y) * x
        gradient_w = np.mean(gradient_w, axis=0)
        gradient_w = gradient_w[:, np.newaxis]

        gradient_b = z-y
        gradient_b = np.mean(gradient_b)
        return gradient_w, gradient_b

    # 更新梯度
    def update(self, gradient_w, gradient_b, eta=0.01):
        self.w -= eta * gradient_w
        self.b -= eta * gradient_b

    # 训练函数
    def train(self, training_data, num_epochs, batch_size=10, eta=0.01):
        losses = []
        n = len(training_data)
        for epoch_id in range(num_epochs):
            # 在每轮迭代开始之前，将训练数据的顺序随机打乱
            # 然后再按每次取batch_size条数据的方式取出
            np.random.shuffle(training_data)

            # 将训练数据拆分
            mini_batches = [training_data[k: k+ batch_size] for k in range(0, n, batch_size)]

            for iter_id, mini_batch in enumerate(mini_batches):
                x = mini_batch[:, :-1]
                y = mini_batch[:, -1:]

                a = self.forward(x)
                L = self.loss(a, y)
                gradient_w, gradient_b = self.gradient(x, y)
                self.update(gradient_w, gradient_b, eta)
                losses.append(L)
                print('Epoch {:3d} / iter {:3d}, loss = {:.4f}'.
                      format(epoch_id, iter_id, L))

        return losses


if __name__ == '__main__':
    traing_data, test_data = load_data()
    print(traing_data)
    # x的形状是(404, 13)
    # y的形状是(404, 1)
    x = traing_data[:, : -1]
    y = traing_data[:, -1:]

    net = NetWork(13)
    num_epoches = 10000
    losses = net.train(traing_data, num_epochs=num_epoches, batch_size=100, eta=0.01)
    # 训练结果可视化
    import matplotlib.pyplot as plt
    plot_x = np.arange(len(losses))
    plot_y = np.array(losses)
    plt.plot(plot_x, plot_y)
    plt.show()

经过10000个epoch的训练，最终输出的loss为：0.0032。（在不使用随机梯度下降的情况下，经过实验，训练1000000次，loss趋近于0.11）。并且在同样训练10000个epoch的情况下，随机梯度下降更快。

在训练过程中的损失函数值的图表：