四、用简单神经网络识别手写数字（内含代码详解及订正）

本博客主要内容为图书《神经网络与深度学习》和National Taiwan University (NTU)林轩田老师的《Machine Learning》的学习笔记，因此在全文中对它们多次引用。初出茅庐，学艺不精，有不足之处还望大家不吝赐教。

1. 前期准备

1.1 数据集

  MNIST数据集是基于NIST（美国国家标准与技术研究院）收集的两个数据集合。为了构建MNIST，NIST数据集合被Yann LeCun，Corinna Cortes和Christopher J. C. Burges拆分放入一个更方便的格式，本文所使用的数据集是在一种更容易在Python中加载和操纵MNIST数据的形式。从蒙特利尔大学的LISA机器学习实验室获得了这个特殊格式的数据。这个数据可以通过克隆这本书的代码仓库获得数据；如果你不使用git，那么你能够在这里下载数据和代码。需要注意的是在这里的代码的版本是基于python2的，如果使用python3的用户直接运行这个代码会出现错误，在本文的后续给出经过经过修改后的基于python3 的代码版本。
  在这个数据集中包含60,000个训练图像和10,000个测试图像，在训练的过程中将测试集（testing data）保持不变，将训练集(training data)分为由50000个图像组成的训练集以及由剩下的10000个图像组成的验证集（validation set）。验证集对于解决如何去设置神经网络中的超参数（hyper-parameter）是十分有用的。因此以后提到「MNIST训练数据」指的不是原始的60,000图像数据集，而是我们的50,000图像数据集。

1.2 软件的配置

  除了MNIST数据之外，我们还需要一个叫做Numpy的用于处理快速线性代数的Python库。考虑到今后会更加深入的研究神经网络，因此本文给出配置Keras的方法，其中几乎已经涵盖了今后会用到的绝大多数的库。在这里的方式是默认这是一台全新的电脑从头配置Keras，已经配置过的用户可以跳过此步。

2. 搭建分类用神经网络

  为了识别数字，我们将会使用一个如图1的三层神经网络

图1. 用于手写数字分类的神经网络

  这个网络的输入层是用于训练数据的 28×28 的手写数字位图，因此我们的输入层包含了 28×28=784 个神经元。为了方便起见，在图1中没有完全画出 784 个输入神经元。输入的像素点是其灰度值，0.0 代表白色，1.0 代表黑色，中间值表示不同程度的灰度值。
  在网络的第二层，即隐层中设置了 n 个神经元， n 是一个需要通过实验确定的值，在这里取 n 的值为15。
  网络的输出层包含了10个神经元，把输出层神经元依次标记为 0 到 9，找到拥有最高的激活值的神经元，将它的标记作为神将网络的结果进行输出。例如 6 号神经元有最高值，那么我们的神经网络预测输入数字是 6，对其它的神经元也如此。
  为什么用 10 个输出神经元，如果采用二进制编码的形式只需要4个神经元即可。最终的判断是基于经验主义的：我们可以实验两种不同的网络设计，结果证明对于这个特定的问题而言，10 个输出神经元的神经网络比 4 个的识别效果更好。理解这一现象的原因可以参考第一篇博客中第二小节“迈向深度学习”的内容，神经网络是通过一部分一部分学习的，因此把数字的最高有效位和数字的形状联系起来并不是一个简单的问题。很难想象出有什么恰当的历史原因，一个数字的形状要素会和一个数字的最高有效位有什么紧密联系。

3. 手写数字分类的主要代码

  这里是用python搭建神经网络识别手写数字的全部代码，实际上程序只包含74行非空行、非注释代码。在这里面主体是一个叫做Network的类，在这个类里面定义了正向传播、反向传播、随机梯度下降等方法，其次是一些辅助函数，即sigmoid函数与sigmoid导数计算函数，具体方法如下所示

def init(self, sizes)：初始化
def feedforward(self, a ): 前馈
def SGD(self, training_data, epochs, mini_batch_size, eta,test_data=None) : 随机梯度下降函数
def update_mini_batch(self, mini_batch, eta) : 小块训练集更新参数
def backprop(self, x, y) : 反向传播
def cost_derivative(self, output_activations, y) : 损失函数在输出层的导数（微商）
def evaluate(self, test_data) : 返回神经网络分类正确的总个数
def sigmoid(z) : sigmoid函数值
def sigmoid_prime(z) : sigmoid导数值

在上述这些方法中存在调用关系，将他们的调用关系绘制出来如图2所示

图2. 神经网络中函数相互调用关系

  因此我们从左至右一次对这就个方法进行详尽的分析

3.1 辅助函数

#### Miscellaneous functionsdef sigmoid(z):    """The sigmoid function."""    return 1.0/(1.0+np.exp(-z))def sigmoid_prime(z):    """Derivative of the sigmoid function."""    return sigmoid(z)*(1-sigmoid(z))

  这两个函数是计算sigmoid函数及导数值的方法，十分容易理解，因此不做分析

3.2 cost_derivative（）

def cost_derivative(self, output_activations, y):     """Return the vector of partial derivatives \partial C_x /     \partial a for the output activations."""    return (output_activations-y)

  这个方法反回了损失函数对输出层激活函数的值的导数。这是十分易于理解的，因为在这里采用了二次函数作为损失函数，对该损失进行求导即为上述的用激活值减去实际值的计算方式。

3.3 feedforward（）

def feedforward(self, a):    """Return the output of the network if ``a`` is input."""    for b, w in zip(self.biases, self.weights):         a = sigmoid(np.dot(w, a)+b)    return a

  这个方法返回输入值为a的时候的输出值。在这里需要注意的是zip()函数，其作用是将对应的元素组成一个元组，将元组构成列表，在此段程序中是将biases与weights的值组合起来构成一个列表。

3.4 evaluate（）

def evaluate(self, test_data):    """Return the number of test inputs for which the neural    network outputs the correct result. Note that the neural    network's output is assumed to be the index of whichever    neuron in the final layer has the highest activation."""    test_results = [(np.argmax(self.feedforward(x)), y)                    for (x, y) in test_data]    return sum(int(x == y) for (x, y) in test_results)

  返回神经网络对测试集样本分类正确的数量，应该注意的是这里假定神经网络的输出是输出层中具有最高激活值的神经元的索引值。其中test_results中存储的值是由元组（x, y）构成的列表，其中x代表输出层激活值最大的神经元的索引值，因此test_results中实际存储的是有神经网络预测出的类别与实际的类别。

3.5 backprop（）

def backprop(self, x, y):    """Return a tuple ``(nabla_b, nabla_w)`` representing the    gradient for the cost function C_x.  ``nabla_b`` and    ``nabla_w`` are layer-by-layer lists of numpy arrays, similar    to ``self.biases`` and ``self.weights``."""    nabla_b = [np.zeros(b.shape) for b in self.biases]    nabla_w = [np.zeros(w.shape) for w in self.weights]    # feedforward    activation = x    activations = [x] # list to store all the activations, layer by layer    zs = [] # list to store all the z vectors, layer by layer    for b, w in zip(self.biases, self.weights):        z = np.dot(w, activation)+b        zs.append(z)        activation = sigmoid(z)        activations.append(activation)    # backward pass    delta = self.cost_derivative(activations[-1], y) * \        sigmoid_prime(zs[-1])    nabla_b[-1] = delta    nabla_w[-1] = np.dot(delta, activations[-2].transpose())    # Note that the variable l in the loop below is used a little    # differently to the notation in Chapter 2 of the book.  Here,    # l = 1 means the last layer of neurons, l = 2 is the    # second-last layer, and so on.  It's a renumbering of the    # scheme in the book, used here to take advantage of the fact    # that Python can use negative indices in lists.    for l in range(2, self.num_layers):        z = zs[-l]        sp = sigmoid_prime(z)        delta = np.dot(self.weights[-l+1].transpose(), delta) * sp        nabla_b[-l] = delta        nabla_w[-l] = np.dot(delta, activations[-l-1].transpose())    return (nabla_b, nabla_w)

  返回一个元组(nabla_b, nabla_w)用来代表损失函数C_x的梯度，其中英文单词nabla[‘næblə] 的含义是向量微分算子。nabla_b和nabla_w是逐层列出NumPy数组，与self.biases和self.weights十分相似。np.zeros(b.shape)中的shape是读取矩阵的形状或者长度，因此该行及下一行是对梯度矩阵的初始化。
  根据第三节神经网络的训练中“反向传播算法整体描述”中所讲，要计算每一层网络的权值输入 z 和激活值 activation，因此在# feedforward部分首创建了用于储存这两个值的列表，计算他们并存储起来。
  在# backward pass部分，self.cost_derivative(activations[-1], y)是计算损失函数对于输出层激活值的导数值，即 ∇aC，其中activations[-1]指的是activations中的最后一层，即输出层的值；而sigmoid_prime(zs[-1])计算的是输出层激活值对于权值输入的导数，即 σ′(z)，因此delta实际上是上一节中提到的反向传播四个基本等式中的第一个等式δL=∇aC⨀σ′(zL)；之后按照四个等式中的后两个 ∂C∂blj=δj l，∂C∂wjkl=al−1kδlj 计算对于最后一层的梯度。其中需要转置transpose()的原因在于这里的np.dot()是真正意义上的矩阵运算，因此activations[-2]需要添加转置，才能通过计算得到 w 矩阵
  需要注意的是下面循环中的变量L与本书之前介绍的神经网络位于第几层中的符号有点不同。在这里 l=1 是神经网络的最后一层， l=2 是神经网络的倒数第二层。这只是一种编号方案，为了利用Python可以在列表中使用负指数的特性。
  需要注意的是，在python3中不再提供xrange()函数，因此本文采用range()，它整合了xrange()的所有功能。假设包括输入层、输出层在内一共有四层神经网络，则range(2, self.num_layers)=（2，3），而z = zs[-l]相当于去的是z = zs[-2]，即权值输入的倒数第二层，因此对应了注释中这里的变量 l 是与之前的定义不同的。剩下的与之前讲的相同，是根据基本等式进行计算的过程。

3.6 update_mini_batch（）

def update_mini_batch(self, mini_batch, eta):    """Update the network's weights and biases by applying    gradient descent using backpropagation to a single mini batch.    The ``mini_batch`` is a list of tuples ``(x, y)``, and ``eta``    is the learning rate."""    nabla_b = [np.zeros(b.shape) for b in self.biases]    nabla_w = [np.zeros(w.shape) for w in self.weights]    for x, y in mini_batch:        delta_nabla_b, delta_nabla_w = self.backprop(x, y)        nabla_b = [nb+dnb for nb, dnb in zip(nabla_b, delta_nabla_b)]        nabla_w = [nw+dnw for nw, dnw in zip(nabla_w, delta_nabla_w)]    self.weights = [w-(eta/len(mini_batch))*nw                    for w, nw in zip(self.weights, nabla_w)]    self.biases = [b-(eta/len(mini_batch))*nb                   for b, nb in zip(self.biases, nabla_b)]

  利用反向传播梯度下降法对单个批处理网络的权值和偏差进行更新。其中的”mini_batch”是元组(x, y)的一个列表，而eta是学习率。按照梯度下降法公式更新参数。

3.7 def SGD（）

def SGD(self, training_data, epochs, mini_batch_size, eta,        test_data=None):    """Train the neural network using mini-batch stochastic    gradient descent.  The ``training_data`` is a list of tuples    ``(x, y)`` representing the training inputs and the desired    outputs.  The other non-optional parameters are    self-explanatory.  If ``test_data`` is provided then the    network will be evaluated against the test data after each    epoch, and partial progress printed out.  This is useful for    tracking progress, but slows things down substantially."""    if test_data: n_test = len(test_data)    n = len(training_data)    for j in range(epochs):        random.shuffle(training_data)        mini_batches = [            training_data[k:k+mini_batch_size]            for k in xrange(0, n, mini_batch_size)]        for mini_batch in mini_batches:            self.update_mini_batch(mini_batch, eta)        if test_data:            print "Epoch {0}: {1} / {2}".format(                j, self.evaluate(test_data), n_test)        else:            print "Epoch {0} complete".format(j)

  使用随机梯度下降法训练神经网络，training_data是元组（x, y）的列表，其中x代表输入，而y代表期望输出。其他的非优化参数是自我解释的。如果测试数据被提供给，网络将在每一轮训练完成后使用训练数据进行评估，并把评估值打印输出出来。这对于训练过程是有用的，但是会明显降低训练速度。
  len(test_data)的返回值是测试数据test_data的长度；random.shuffle()将序列的所有元素随机排序，并按块取训练数据进行参数更新。有必要的时候输出测试集结果。

4. 手写数字分类实验及结果分析

  先加载MNIST数据。我将用下面所描述的一小 段辅助程序 mnist_loader.py 来完成。我们在一个Python shell中执行下面的命令，

>>> import  mnist_loader >>> training_data,  validation_data,    test_data   =   \ ...mnist_loader.load_data_wrapper()

  当然，这也可以被做成一个单独的Python程序，但在Python shell执行最方便。在加载完MNIST数据之后，我们将设置一个有30个隐层神经元的 Network 。我们在导入如上 所列的名字为 network 的Python程序后做，

>>> import  network >>> net = network.Network([784,30,10])

  最后，我们将使用随机梯度下降来从MNIST training_data 学习超过30次迭代，迷你块大小 为10，学习率 η=3.0，

>>> net.SGD(training_data,  30, 10, 3.0,    test_data=test_data)

  最后神经网络输出的正确率大概在95%左右，但是每个人的结果可能不尽相同，因为随机初始化的权重和偏置的值可能是不一样的。但是当增加隐层的神经元之后效果可能提高或者降低，我们将在下一章采用其他技术使得增加神经元的个数对实验结果的影响没有这样剧烈。

4.1 关于学习率对于实验结果的影响分析

  在训练神经网络的过程中主要有迭代次数、mini-batch的大小和学习率 η 三个超参数。如果我们设定学习率为 η=0.001

>>> net =   network.Network([784,   100,    10]) >>> net.SGD(training_data,  30, 10, 0.001,  test_data=test_data)

下面是实验结果，可以看到效果不是很好

Epoch 0: 1139 / 10000
Epoch 1: 1136 / 10000
Epoch 2: 1135 / 10000
…
Epoch 27: 2101 / 10000
Epoch 28: 2123 / 10000
Epoch 29: 2142 / 10000

  然而，你能够看到网络的性能随着时间的推移在缓慢的变好。这意味着着我们应该增大学习率，逐渐增大学习率直至得到一个较为理想的结果。其他超参数的调节将在后面继续进行介绍。