本文着重讲述经典BP神经网络的数学推导过程，并辅助一个小例子。本文不会介绍机器学习库(比如sklearn, TensorFlow等)的使用。 欲了解卷积神经网络的内容，请参见我的另一篇博客一文搞定卷积神经网络——从原理到应用。

  本文难免会有叙述不合理的地方，希望读者可以在评论区反馈。我会及时吸纳大家的意见，并在之后的chat里进行说明。

本文参考了一些资料，在此一并列出。

• http://neuralnetworksanddeeplearning.com/chap2.html
• https://www.deeplearning.ai/ coursera对应课程视频讲义
• coursera华盛顿大学机器学习专项
• 周志华《机器学习》
• 李航《统计学习方法》
• 张明淳《工程矩阵理论》

0. 什么是人工神经网络？

  首先给出一个经典的定义：“神经网络是由具有适应性的简单单元组成的广泛并行互连网络，它的组织能够模拟生物神经系统对真实世界物体所作出的交互反应”[Kohonen, 1988]。

  这种说法虽然很经典，但是对于初学者并不是很友好。比如我在刚开始学习的时候就把人工神经网络想象地很高端，以至于很长一段时间都不能理解为什么神经网络能够起作用。类比最小二乘法线性回归问题，在求解数据拟合直线的时候，我们是采用某种方法让预测值和实际值的“偏差”尽可能小。同理，BP神经网络也做了类似的事情——即通过让“偏差”尽可能小，使得神经网络模型尽可能好地拟合数据集。

1. 神经网络初探

1.1 神经元模型

  神经元模型是模拟生物神经元结构而被设计出来的。典型的神经元结构如下图1所示：

【图1 典型神经元结构 （图片来自维基百科）】

  神经元大致可以分为树突、突触、细胞体和轴突。树突为神经元的输入通道，其功能是将其它神经元的动作电位传递至细胞体。其它神经元的动作电位借由位于树突分支上的多个突触传递至树突上。神经细胞可以视为有两种状态的机器，激活时为“是”，不激活时为“否”。神经细胞的状态取决于从其他神经细胞接收到的信号量，以及突触的性质（抑制或加强）。当信号量超过某个阈值时，细胞体就会被激活，产生电脉冲。电脉冲沿着轴突并通过突触传递到其它神经元。（内容来自维基百科“感知机”）

  同理，我们的神经元模型就是为了模拟上述过程，典型的神经元模型如下：

【图2 典型神经元模型结构 （摘自周志华老师《机器学习》第97页）】

  这个模型中，每个神经元都接受来自其它神经元的输入信号，每个信号都通过一个带有权重的连接传递，神经元把这些信号加起来得到一个总输入值，然后将总输入值与神经元的阈值进行对比（模拟阈值电位），然后通过一个“激活函数”处理得到最终的输出（模拟细胞的激活），这个输出又会作为之后神经元的输入一层一层传递下去。

1.2 神经元激活函数

  本文主要介绍2种激活函数，分别是$sigmoid$和$relu$函数，函数公式如下：

  做函数图如下：

$sigmoid(z)$

$relu(z)$
【图3 激活函数】

补充说明
【补充说明的内容建议在看完后文的反向传播部分之后再回来阅读，我只是为了文章结构的统一把这部分内容添加在了这里】

  引入激活函数的目的是在模型中引入非线性。如果没有激活函数，那么无论你的神经网络有多少层，最终都是一个线性映射，单纯的线性映射无法解决线性不可分问题。引入非线性可以让模型解决线性不可分问题。

  一般来说，在神经网络的中间层更加建议使用$relu$函数，两个原因：

• $relu$函数计算简单，可以加快模型速度；
• 由于反向传播过程中需要计算偏导数，通过求导可以得到$sigmoid$函数导数的最大值为0.25，如果使用$sigmoid$函数的话，每一层的反向传播都会使梯度最少变为原来的四分之一，当层数比较多的时候可能会造成梯度消失，从而模型无法收敛。

1.3 神经网络结构

  我们使用如下神经网络结构来进行介绍，第0层是输入层（3个神经元）， 第1层是隐含层（2个神经元），第2层是输出层：

【图4 神经网络结构（手绘）】

  

  现在，我们使用矩阵形式重写这个公式：

  

  同理，

  更一般地，我们可以把前向传播过程表示：

  到这里，我们已经讲完了前向传播的过程，值得注意的是，这里我们只有一个输入样本，对于多个样本同时输入的情况是一样的，只不过我们的输入向量不再是一列，而是m列，每一个都表示一个输入样本。

  多样本输入情况下的表示为：

其中，此时

每一列都表示一个样本，从样本1到m

之后我们的叙述都是先讨论单个样本的情况，再扩展到多个样本同时计算。

2. 损失函数和代价函数

  说实话，**损失函数（Loss Function）和代价函数（Cost Function）**并没有一个公认的区分标准，很多论文和教材似乎把二者当成了差不多的东西。

  为了后面描述的方便，我们把二者稍微做一下区分（这里的区分仅仅对本文适用，对于其它的文章或教程需要根据上下文自行判断含义）：

  损失函数主要指的是对于单个样本的损失或误差；代价函数表示多样本同时输入模型的时候总体的误差——每个样本误差的和然后取平均值。

  举个例子，如果我们把单个样本的损失函数定义为：

$$L(a,y)=-[y\cdot log(a)+(1-y)\cdot log(1-a)]$$

  那么对于m个样本，代价函数则是：

3. 反向传播

  反向传播的基本思想就是通过计算输出层与期望值之间的误差来调整网络参数，从而使得误差变小。

  反向传播的思想很简单，然而人们认识到它的重要作用却经过了很长的时间。后向传播算法产生于1970年，但它的重要性一直到David Rumelhart，Geoffrey Hinton和Ronald Williams于1986年合著的论文发表才被重视。

  事实上，人工神经网络的强大力量几乎就是建立在反向传播算法基础之上的。反向传播基于四个基础等式，数学是优美的，仅仅四个等式就可以概括神经网络的反向传播过程，然而理解这种优美可能需要付出一些脑力。事实上，反向传播如此之难，以至于相当一部分初学者很难进行独立推导。所以如果读者是初学者，希望读者可以耐心地研读本节。对于初学者，我觉得拿出1-3个小时来学习本小节是比较合适的，当然，对于熟练掌握反向传播原理的读者，你可以在十几分钟甚至几分钟之内快速浏览本节的内容。

3.1 矩阵补充知识

  对于大部分理工科的研究生，以及学习过矩阵论或者工程矩阵理论相关课程的读者来说，可以跳过本节。

  本节主要面向只学习过本科线性代数课程或者已经忘记矩阵论有关知识的读者。

  总之，具备了本科线性代数知识的读者阅读这一小节应该不会有太大问题。本节主要在线性代数的基础上做一些扩展。（不排除少数本科线性代数课程也涉及到这些内容，如果感觉讲的简单的话，勿喷）

3.1.1 求梯度矩阵

  

  即

  同理，一个输入是向量（向量一般指列向量，本文在没有特殊声明的情况下默认指的是列向量）的函数

，则有：

  根据定义，很容易得到以下性质：

  有了上述知识，我们来举个例子：

3.1.2 海塞矩阵

  定义一个输入为$n$维向量，输出为实数的函数

，那么海塞矩阵（Hessian Matrix）定义为多元函数$f$的二阶偏导数构成的方阵：

  由上式可以看出，海塞矩阵总是对称阵。

3.1.3 总结

  根据3.1.1和3.1.2小节的内容很容易得到以下等式：

 

  比如对于矩阵

矩 阵

  矩阵相乘

  对应元素相乘使用符号$\odot$表示：

3.3 梯度下降法原理

  通过之前的介绍，相信大家都可以自己求解梯度矩阵（向量）了。

  那么梯度矩阵（向量）求出来的意义是什么？从几何意义讲，梯度矩阵代表了函数增加最快的方向，因此，沿着与之相反的方向就可以更快找到最小值。如图5所示：

【图5 梯度下降法 图片来自百度】

  反向传播的过程就是利用梯度下降法原理，慢慢的找到代价函数的最小值，从而得到最终的模型参数。梯度下降法在反向传播中的具体应用见下一小节。

3.4 反向传播原理（四个基础等式）

  反向传播能够知道如何更改网络中的权重$w$ 和偏差$b$ 来改变代价函数值。最终这意味着它能够计算偏导数

​ 和

  

  

  

  于是，每一层的误差向量可以表示为：

  下面开始正式介绍四个基础等式【确切的说是四组等式】

  **注意：**这里我们的输入为单个样本(所以我们在下面的公式中使用的是损失函数而不是代价函数)。多个样本输入的公式会在介绍完单个样本后再介绍。

• 等式1 ：输出层误差

  其中，$L$表示输出层层数。以下用$\partial L$ 表示

  写成矩阵形式是：

  【注意是对应元素相乘，想想为什么？】

  

• 等式2： 隐含层误差
  

  写成矩阵形式：

 

• 等式3：参数变化率

  写成矩阵形式：

  

• 等式4：参数更新规则

  这应该是这四组公式里最简单的一组了，根据梯度下降法原理，朝着梯度的反方向更新参数：

  写成矩阵形式：

  这里的$\alpha$指的是学习率。学习率指定了反向传播过程中梯度下降的步长。

3.5 反向传播总结

  我们可以得到如下最终公式：

3.5.1 单样本输入公式表

3.5.2 多样本输入公式表

  多样本：需要使用代价函数，如果有m个样本，那么由于代价函数有一个 1\m的常数项，因此所有的参数更新规则都需要有一个1\m的前缀。

  多样本同时输入的时候需要格外注意维度匹配，一开始可能觉得有点混乱，但是不断加深理解就会豁然开朗。

  通过前面的介绍，相信读者可以发现BP神经网络模型有一些参数是需要设计者给出的，也有一些参数是模型自己求解的。

  那么，哪些参数是需要模型设计者确定的呢？

  比如，学习率 α ，隐含层的层数，每个隐含层的神经元个数，激活函数的选取，损失函数（代价函数）的选取等等，这些参数被称之为超参数。

  其它的参数，比如权重矩阵$w$和偏置系数$b$在确定了超参数之后是可以通过模型的计算来得到的，这些参数称之为普通参数，简称参数。

  超参数的确定其实是很困难的。因为你很难知道什么样的超参数会让模型表现得更好。比如，学习率太小可能造成模型收敛速度过慢，学习率太大又可能造成模型不收敛；再比如，损失函数的设计，如果损失函数设计不好的话，可能会造成模型无法收敛；再比如，层数过多的时候，如何设计网络结构以避免梯度消失和梯度爆炸……

  神经网络的程序比一般程序的调试难度大得多，因为它并不会显式报错，它只是无法得到你期望的结果，作为新手也很难确定到底哪里出了问题（对于自己设计的网络，这种现象尤甚，我目前也基本是新手，所以这些问题也在困扰着我）。当然，使用别人训练好的模型来微调看起来是一个捷径……

  总之，神经网络至少在目前来看感觉还是黑箱的成分居多，希望通过大家的努力慢慢探索吧。

4. 是不是猫？

  本小节主要使用上述公式来完成一个小例子，这个小小的神经网络可以告诉我们一张图片是不是猫。本例程参考了coursera的作业，有改动。

  在实现代码之前，先把用到的公式列一个表格吧，这样对照着看大家更清晰一点(如果你没有2个显示器建议先把这些公式抄写到纸上，以便和代码对照)：

  那先看结果吧，我们的代价函数是：

  其中，$m$是样本数量；

  下面简单介绍一下这个代价函数是怎么来的（作者非数学专业，不严谨的地方望海涵）。
.
  代价函数的确定用到了统计学中的**“极大似然法”**，既然这样，那就不可避免地要介绍一下“极大似然法”了。极大似然法简单来说就是“在模型已定，参数未知的情况下，根据结果估计模型中参数的一种方法"，换句话说，极大似然法提供了一种给定观察数据来评估模型参数的方法。

  举个例子（本例参考了知乎相关回答），一个不透明的罐子里有黑白两种球（球仅仅颜色不同，大小重量等参数都一样）。有放回地随机拿出一个小球，记录颜色。重复10次之后发现7次是黑球，3次是白球。问你罐子里白球的比例？

  相信很多人可以一口回答“30%”，那么，为什么呢？背后的原理是什么呢？

  

  极大似然法有一个重要的假设：

假设所有样本独立同分布！！！

  好了，现在来看看我们的神经网络模型。

  最后一层我们用sigmoid函数求出一个激活输出a，如果a大于0.5，就表示这个图片是猫（$y=1$），否则就不是猫（$y=0$）。因此:
$$P(y=1|x;\theta )=a$$
$$P(y=0|x;\theta )=1-a$$

公式解释：
上述第一个公式表示，给定模型参数$\theta$和输入$x$，是猫的概率是$P(y=1|x;\theta )=a$

  把两个公式合并成一个公式，即

这里的$\theta$指的就是我们神经网络的权值参数和偏置参数。

  那么似然函数

  变成对数形式：

  所以我们的目标就是最大化这个对数似然函数，也就是最小化我们的代价函数：

  其中，$m$是样本数量；

  好了，终于可以开始写代码了，码字手都有点酸了，不得不说公式真的好难打。

由于代码比较简单就没有上传github。本文代码和数据文件可以在这里下载: https://pan.baidu.com/s/1q_PzaCSXOhRLOJVF5-vy2Q，密码： d7vx

其他下载源：
https://drive.google.com/file/d/0B6exrzrSxlh3TmhSV0ZNeHhYUmM/view?usp=sharing

4.1 辅助函数

  辅助函数主要包括激活函数以及激活函数的反向传播过程函数：
其中，激活函数反向传播代码对应公式4和9.

def sigmoid(z):
    """
    使用numpy实现sigmoid函数
    
    参数：
    Z numpy array
    输出：
    A 激活值（维数和Z完全相同）
    """
    return 1/(1 + np.exp(-z))

def relu(z):
    """
    线性修正函数relu
    
    参数：
    z numpy array
    输出：
    A 激活值（维数和Z完全相同）
    
    """
    return np.array(z>0)*z

def sigmoidBackward(dA, cacheA):
    """
    sigmoid的反向传播
    
    参数：
    dA 同层激活值
    cacheA 同层线性输出
    输出：
    dZ 梯度
    
    """
    s = sigmoid(cacheA)
    diff = s*(1 - s)
    dZ = dA * diff
    return dZ

def reluBackward(dA, cacheA):
    """
    relu的反向传播
    
    参数：
    dA 同层激活值
    cacheA 同层线性输出
    输出：
    dZ 梯度
    
    """
    Z = cacheA
    dZ = np.array(dA, copy=True) 
    dZ[Z <= 0] = 0
    return dZ

  另外一个重要的辅助函数是数据读取函数和参数初始化函数：

def loadData(dataDir):
    """
    导入数据
    
    参数：
    dataDir 数据集路径
    输出：
    训练集，测试集以及标签
    """
    train_dataset = h5py.File(dataDir+'/train.h5', "r")
    train_set_x_orig = np.array(train_dataset["train_set_x"][:]) # your train set features
    train_set_y_orig = np.array(train_dataset["train_set_y"][:]) # your train set labels

    test_dataset = h5py.File(dataDir+'/test.h5', "r")
    test_set_x_orig = np.array(test_dataset["test_set_x"][:]) # your test set features
    test_set_y_orig = np.array(test_dataset["test_set_y"][:]) # your test set labels

    classes = np.array(test_dataset["list_classes"][:]) # the list of classes
    
    train_set_y_orig = train_set_y_orig.reshape((1, train_set_y_orig.shape[0]))
    test_set_y_orig = test_set_y_orig.reshape((1, test_set_y_orig.shape[0]))
    
    return train_set_x_orig, train_set_y_orig, test_set_x_orig, test_set_y_orig, classes

def iniPara(laydims):
    """
    随机初始化网络参数
    
    参数：
    laydims 一个python list
    输出：
    parameters 随机初始化的参数字典（”W1“，”b1“，”W2“，”b2“, ...）
    """
    np.random.seed(1)
    parameters = {}
    for i in range(1, len(laydims)):
        parameters['W'+str(i)] = np.random.randn(laydims[i], laydims[i-1])/ np.sqrt(laydims[i-1])
        parameters['b'+str(i)] = np.zeros((laydims[i], 1))
    return parameters

4.2 前向传播过程

对应公式1和2.

def forwardLinear(W, b, A_prev):
    """
    前向传播
    """
    Z = np.dot(W, A_prev) + b
    cache = (W, A_prev, b)
    return Z, cache

def forwardLinearActivation(W, b, A_prev, activation):
    """
    带激活函数的前向传播
    """
    Z, cacheL = forwardLinear(W, b, A_prev)
    cacheA = Z
    if activation == 'sigmoid':
        A = sigmoid(Z)
    if activation == 'relu':
        A = relu(Z)
    cache = (cacheL, cacheA)
    return A, cache

def forwardModel(X, parameters):
    """
    完整的前向传播过程
    """
    layerdim = len(parameters)//2
    caches = []
    A_prev = X
    for i in range(1, layerdim):
        A_prev, cache = forwardLinearActivation(parameters['W'+str(i)], parameters['b'+str(i)], A_prev, 'relu')
        caches.append(cache)
        
    AL, cache = forwardLinearActivation(parameters['W'+str(layerdim)], parameters['b'+str(layerdim)], A_prev, 'sigmoid')
    caches.append(cache)
    
    return AL, caches

4.3 反向传播过程

线性部分反向传播对应公式5和6。

def linearBackward(dZ, cache):
    """
    线性部分的反向传播
    
    参数：
    dZ 当前层误差
    cache （W, A_prev, b）元组
    输出：
    dA_prev 上一层激活的梯度
    dW 当前层W的梯度
    db 当前层b的梯度
    """
    W, A_prev, b = cache
    m = A_prev.shape[1]
    
    dW = 1/m*np.dot(dZ, A_prev.T)
    db = 1/m*np.sum(dZ, axis = 1, keepdims=True)
    dA_prev = np.dot(W.T, dZ)
    
    return dA_prev, dW, db

非线性部分对应公式3、4、5和6 。

def linearActivationBackward(dA, cache, activation):
    """
    非线性部分的反向传播
    
    参数：
    dA 当前层激活输出的梯度
    cache （W, A_prev, b）元组
    activation 激活函数类型
    输出：
    dA_prev 上一层激活的梯度
    dW 当前层W的梯度
    db 当前层b的梯度
    """
    cacheL, cacheA = cache
    
    if activation == 'relu':
        dZ = reluBackward(dA, cacheA)
        dA_prev, dW, db = linearBackward(dZ, cacheL)
    elif activation == 'sigmoid':
        dZ = sigmoidBackward(dA, cacheA)
        dA_prev, dW, db = linearBackward(dZ, cacheL)
    
    return dA_prev, dW, db

完整反向传播模型：

def backwardModel(AL, Y, caches):
    """
    完整的反向传播过程
    
    参数：
    AL 输出层结果
    Y 标签值
    caches 【cacheL, cacheA】
    输出：
    diffs 梯度字典
    """
    layerdim = len(caches)
    Y = Y.reshape(AL.shape)
    L = layerdim
    
    diffs = {}
    
    dAL = - (np.divide(Y, AL) - np.divide(1 - Y, 1 - AL))
    
    currentCache = caches[L-1]
    dA_prev, dW, db =  linearActivationBackward(dAL, currentCache, 'sigmoid')
    diffs['dA' + str(L)], diffs['dW'+str(L)], diffs['db'+str(L)] = dA_prev, dW, db
    
    for l in reversed(range(L-1)):
        currentCache = caches[l]
        dA_prev, dW, db =  linearActivationBackward(dA_prev, currentCache, 'relu')
        diffs['dA' + str(l+1)], diffs['dW'+str(l+1)], diffs['db'+str(l+1)] = dA_prev, dW, db
        
    return diffs

4.4 测试结果

  打开你的jupyter notebook，运行我们的BP.ipynb文件，首先导入依赖库和数据集，然后使用一个循环来确定最佳的迭代次数大约为2000：

【图6】

  最后用一个例子来看一下模型的效果——判断一张图片是不是猫：

【图7】

好了，测试到此结束。你也可以自己尝试其它的神经网络结构和测试其它图片。

5. 本文小结

  本文主要叙述了经典的全连接神经网络结构以及前向传播和反向传播的过程。通过本文的学习，读者应该可以独立推导全连接神经网络的传播过程，对算法的细节烂熟于心。另外，由于本文里的公式大部分是我自己推导的，瑕疵之处，希望读者不吝赐教。

  虽然这篇文章实现的例子并没有什么实际应用场景，但是自己推导一下这些数学公式并用代码实现对理解神经网络内部的原理很有帮助，继这篇博客之后，我还计划写一个如何自己推导并实现卷积神经网络的教程，如果有人感兴趣，请继续关注我！

  本次内容就到这里，谢谢大家。

订正与答疑：

前向传播过程比较简单，我就不再赘述了。

这里主要针对反向传播过程中可能会出现的问题做一个总结：

1. 具体解释一下公式1里面的“堆砌”是什么意思？

有读者对这里不太理解，这其实是因为，我们的输出层不一定是只有一个神经元，可能有好多个神经元，因此损失函数是每个输出神经元“误差”之和，因此才会出现这种$\sum$的形式，然后每个输出神经元的误差函数与其它神经元没有关系，所以只有$k=j$的时候值不是0.

另外，这里说的“堆砌”指的就是：

2. 公式2写成矩阵形式为什么系数矩阵会有转置？自己没搞懂。

这里可能有一点绕，有的读者感觉我的推导不是很明白，所以有必要详细说明一下。

很多读者不明白，写成矩阵形式的时候

里面的“系数矩阵转置”是怎么来的。这里就主要说明一下：

相信大家都已经理解了下面这个前向传播公式：

求偏导这里在原文中有一点错误，应该是：

为了大家有一个直观的感受，来一个具体的例子：

第 1 层的系数矩阵比方是：

那么，

那么，根据之前介绍的求解梯度向量的定义：

这就解释了，为什么会出现转置了。

然后排布成矩阵形式：

也可以解释为什那么会变成转置。

写成矩阵形式，注意检查一下维度匹配的问题。

3. 公式3能具体讲一下矩阵形式是怎么来的吗？

4. 为什么会损失函数不用最小二乘法？

有的读者问，为什么使用极大似然估计，而不用最小二乘法？

其实在线性回归模型中，损失函数也是使用极大似然法来估计的，http://www.jianshu.com/p/0d25be8901c9，只不过线性回归模型中，我们假设残差是高斯分布，因此最终使用极大似然法和最小二乘法的结果是一样的。

如果你直接使用误差平方和最小的话，也不是不可以，但是效果可能会比较差。因为他是非凸的，可能会收敛到局部最优解。

而使用对数似然函数作为损失函数，是凸函数。

5. 能不能顺便介绍一下Python环境的管理？

一般我用anaconda管理Python环境，另外IDE推荐spyder，因为你可以像使用MATLAB的workspace一样直接观察中间变量，对初学者检查自己的程序的正确性很有帮助。

6. 为什么w的初始化使用随机初始化，而b参数的初始化全部初始化为0？

首先，b参数也可以用随机初始化。

为什么不把w全部初始化为0呢？因为这样的话，每次学习所有的隐含神经元学到的东西都一样，最终会导致，和层一个神经元没有区别，所以我们需要随机初始化开打破这种局面，让每个神经元都更好地“学习“。