ex4

文章框架

大体的思路
这部分主要是课程上的内容为主
算法整体思路

下面不以作业的标题而以作业中我们要实现的程序模块为线索来为主要的框架，中间用到哪些课程上讲到的matlab技巧把他们结合在一起

其实实现的难点有两点一个是matlab本身的技巧，还有一个就是公式的向量化。

1. 神经网络训练算法（coursera版）整体描述

1.1 不同点1——代价函数的形式不同

参考下面的2.3.1

1.2 不同点2——最后一层的误差表示方式不同

这里为了简化输出层的误差就用样本的值和输出值相减得到，其实比较严谨的计算方法是样本值和输出值相减后再乘以输出层的偏导数。
参考我之前写的两篇文章。

2. 在完成ex4中过程中的值得提及的实现细节

2.1 数据集

其实这个数据集中暗含了一种解决问题的思路，这就是怎么样把一个现实的问题转化成一个计算机可以计算的问题，我们的目标是要利用计算机识别手写的数字，而这个问题是如何建模的？其实一张图片在计算机里本身也是以数字的方式存储的，在现实中具备某一种特征的物体（比如汽车或者猫）其图片应该也蕴含着它的特征信息（因为人看到这张图片之后我们能够识别出这是什么东西，比如说我们看到汽车的图片不仅能知道它是汽车还能知道它的牌子），这张图片在计算机中是以数字的方式存储的，那么这个物体的特征信息也一定是由数字的方式体现，或者说这个存储这个图片的数字信息中蕴含这某种数字特征，这个数字特征就是这个物体的特征在计算机空间中的数字表现。
（
图片在计算机里的存储原理很像我们在国庆阅兵上看到的组字方队，图片放大了来看都是马赛克，其实一张图片就是一个一个马赛克构成的，对于一张黑白的图片，一个点是从白色到黑色之间的一个颜色（相当于把从白色到黑色之间的渐变分成了255的梯度，深浅，）上是一个0~255的。如果是彩色的，如果是RGB的方式存储那就是三个值（R,G,B），还有类似与CMKY
一个图片的大小是如何决定的呢？组字方队的人数（像素点的多少），以及每一个点容量的上线
）

2.2 本次作业用到网络结构

网络结构分三层，第一层输入层有400个节点，正好和表示每个样本的向量的分量数（或者叫做样本的特征数）一样；第二层隐藏层有25个节点，第三层输出层有10个节点分别对应1到10这10种分类情况。

其实这个结构严格上来说是不太好的，因为如果按照教程里的内容，隐藏层的单元数量最好要大于输入层的单元数量

2.3 代价函数的实现

2.3.0 整体描述

代价函数（nnCostFunction.m）的实现几乎是这次作业的全部，从得分值就可以看出来，这次作业提交了四次代价函数的实现程序，第一次要实现前馈算法以及代价函数（未正则化），第二次要实现正则化的代价函数，第三次要实现反向传播算法（未正则化），第四次要实现正则化的反向传播算法，因此与其说它是求神经网络的代价值不如说是一个将中间执行过程中cost返回的利用批量梯度下降来训练神经网络的函数。
其实也可以这么理解你要求得神经网络的代价，要求得神经网络的代价，你要用前馈传播算法求得输出层的计算值，而要实现出前馈传播算法首先你的神经网络的参数要齐备吧，对于一个参数为知的神经网络，那么我们首先要用批量梯度下降法训练出神经网络的参数，而使用这个方法的过程中需要利用反向传播算法来求得梯度，当我们把所有的前置条件实现了，代价也可以求出来了。
注意求得神经网络的代价和实现代价函数是两回事，代价函数就是一行公式，实现起来很容易，但是代价是把具体值带入到代价函数得到的结果

2.3.1 代价函数的形式问题

在UFDL描述的NN代价函数与coursera课程里代价函数不同，到底哪个是对的？两种都是对的，可以认为这个可以是同一个东西从不同侧面描述的结果，代价函数有两种常见形式，一种是“均方差形式”，还有一种是“交叉熵”形式。UFDL中对的代价函数属于“均方差的形式”，coursera课程的代价函数属于“交叉熵”形式。其实logistic回归所用的代价函数就是交叉熵形式，而logistic回顾可以看成神经网络的一种特殊形式（K=1）的情况，交叉熵的结果和用极大似然估计推导出来的结果是一样的。
本次作业中采用交叉熵的形式。

【机器学习】代价函数（cost function）
【机器学习详解】线性回归、梯度下降、最小二乘的几何和概率解释

2.3.2 关于正则化的一些问题

在这一章的学习过程中，在讲到代价函数的时候，有包含正则项的代价函数，

$$\begin{array}{c}J(\mathrm{\Theta })=-\frac{1}{m}\sum _{i=1}^m\sum _{k=1}^K[y_k^{(i)}\log ((h_{\mathrm{\Theta }}(x^{(i)}){)}_k)+(1-y_k^{(i)})\log (1-(h_{\mathrm{\Theta }}(x^{(i)}){)}_k)]+\frac{\lambda }{2m}\sum _{l=1}^{L-1}\sum _{i=1}^{s_l}\sum _{j=1}^{s_{l+1}}({\mathrm{\Theta }}_{j,i}^{(l)}{)}^2\end{array}$$

在讲到反向传播算法的时候，梯度下降的公式里也有偏导项，

$$D_{i,j}^{(l)}:={\displaystyle \frac{1}{m}}({\mathrm{\Delta }}_{i,j}^{(l)}+\lambda {\mathrm{\Theta }}_{i,j}^{(l)})$$

那么我们把正则化引入到神经网络参数训练的过程中的时候是在代价函数里加入还是在梯度下降的过程是加入，这两个是只实现一个就行，还是两个都要实现呢？
之前我们讲正则化的时候，只要在代价函数里定义了正则项就足够了，因为在实现梯度下降的过程中会对代价函数求偏导，这样梯度下降的公式也会带有正则项了。而神经网络比较特殊，它的梯度不是由代价函数直接推导而来的，它的梯度是由反向传播算法计算而来的，反向传播算法本身是不带正则项的，所以需要我们再梯度下降的公式中加上一个正则项。

作业中首先会让我们实现一个不包含正则项的代价函数，之后利用已经提供好的参数和样本，来计算出这个网络状态下的代价值
之后会让我们实现一个包含正则项的代价函数，之后利用已经提供的参数和样本计算出这个网络下的代价值
之后会让我们实现一个不包含正则项的反向传播算法（比较严谨的说法是“对一个非正则化的神经网络实现反向传播算法以计算梯度”），
当我们验证在非正则化的情况是计算结果是正确的，我们会实现包含正则项的反向传播算法（对一个正则化的神经网络实现反向传播算法以计算梯度）
这其实是step by step的引导你实现最终的整体，分步完成的好处就是也是减少出错方便debug的一种方法

那么这种情况又怎么办呢？有没有可能代价函数实现了正则项，但是求梯度的时候我不带正则项呢？因为我们第三步实现不带正则项的梯度的时候，中间计算过程会不会用到代价函数的计算，那么这样不是代价函数包含正则项，求梯度不包含正则项了吗？这种情况是否有意义，该怎么理解呢？如何解决这个问题，我们还是从引入正则项的目的开始谈起，引入正则项的目的主要是为了弱化某些参数在整个模型中起到的作用，从数学操作上来说就是把这个参数的值缩小，我们就拿简单的线性回归的正则化代价函数为例说明这个问题：

$$J(\theta )=\frac{1}{2m}[\sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)}{)}^2+\lambda \sum _{j=1}^n{\theta }_j^2]$$

为了方便思考，我们再假设一些理想的条件，假设这个函数是一个凸函数，存在唯一的极值，那么这个极值所对应的自变量的值也是固定的，比如说$[{\theta }_1,{\theta }_2,{\theta }_3,...{\theta }_n]$在极值点的取值是$[a_1,a_2,a_3,...,a_n]$,假设我对${\theta }_1$进行正则化，前面加上了一个系数，变成$\lambda {\theta }_1$，这个时候整个代价函数的极值点还是不变的，这点很重要，因为虽然代价函数的样子发生了变化，但是其实本质是不变的，为什么这么说呢，在未正则化的时候，我们假设原来的代价函数可以配方成其中一项为$({\theta }_1-a_1{)}^2$的形式，等到${\theta }_1$被正则化了，我们把$\lambda {\theta }_1$看成一个新的参数，比如说就叫${{\theta }_1}^{{}^{\prime }}$，代价函数这一项就变成$({{\theta }_1}^{{}^{\prime }}-a_1{)}^2$的形式，本质不变，${{\theta }_1}^{{}^{\prime }}=a_1$ 即$\lambda {\theta }_1=a_1$，那么${\theta }_1=\frac{1}{\lambda }{a}_1$（其实严格上上来说，并不是所有代价函数都能通过配方来求极值，我这个是用一个最简单的情况来说明问题）
如果说在利用梯度下降的过程中，代价函数带正则项，而梯度下降的公式不带正则项，那又会出现什么样的问题呢？
我们仍然以线性回归作为例子，

$${\theta }_j:={\theta }_j(1-\alpha \frac{\lambda }{m})-\alpha \frac{1}{m}\sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)})x_j^{(i)}$$

可以看到增加了正则化项的梯度更新公式相对于原来的梯度更新公式来说，每次迭代都在原有更新规则的基础上减少了一个额外的值，也就是达到同样的目标每次步子迈的大了一点。
如果没有这个正则项，每次迈步子还是原来的大小，如果代价函数带上了正则项同时梯度下降也带上正则项，那么不仅$\theta$不用迭代到像原来那么大的值（比如$\lambda =1000$的时候，$\theta$需要迭代到的值只有原来的千分之一），而且迭代的步子也大，总的来说就是更快的收敛。如果代价函数带上了正则项但是梯度下降没有带上正则项，那么就是收敛的不像都带上正则项那么快，但是比都不带上正则项要快，就是这样的区别。
其实我当时思考这个问题的时候，核心的灵感主要是两点，一个是极值点不变，另一个是把带上正则项的参数看成一个整体

2.3.3 偏置处理的问题

怎么样处理偏置项，偏置是放在神经元里面还是放在神经元外面，其实作业上已经给出了说明，就是“Model representation”这一节，在输出层之外的每一层都多加了一个输出始终为“1”的神经元。
具体的做法是在每个$\mathrm{\Theta }$矩阵前面加一列1

2.3.4 nnCostFunction.m第一次实现——实现前馈传播算法、不带正则化的代价函数实现，并且计算出代价

2.3.4.1 想到细节1

这四次用的同一个函数头部，参数列表是一样的，那么必然有些参数是在这个分步实现的过程是用不到的，

2.3.4.2 想到的细节2：能否向量化

很多时候代码的实现难点在于向量化，这次能不能进行向量化，其实我认为是可以的。我们先从不带正则化的公式开始看起

$$\begin{array}{c}J(\mathrm{\Theta })=-\frac{1}{m}\sum _{i=1}^m\sum _{k=1}^K[y_k^{(i)}\log ((h_{\mathrm{\Theta }}(x^{(i)}){)}_k)+(1-y_k^{(i)})\log (1-(h_{\mathrm{\Theta }}(x^{(i)}){)}_k)]\end{array}$$

这个式子比较复杂，为了方便理解我们做一些简化，我们假设m=2，k=3当m=1时，也就是第一个样本$(x^{(1)},y^{(1)})$来说，在$J(\mathrm{\Theta })$展开式对应的部分为：

$$\begin{array}{cc}[y_1^{(1)}\log ((h_{\mathrm{\Theta }}(x^{(1)}){)}_1)+(1-y_1^{(1)})\log (1-(h_{\mathrm{\Theta }}(x^{(1)}){)}_1)]+& \text{对应m=1,k=1的情况}\\ [y_2^{(1)}\log ((h_{\mathrm{\Theta }}(x^{(1)}){)}_2)+(1-y_2^{(1)})\log (1-(h_{\mathrm{\Theta }}(x^{(1)}){)}_2)]+& \text{对应m=1,k=2的情况}\\ [y_3^{(1)}\log ((h_{\mathrm{\Theta }}(x^{(1)}){)}_3)+(1-y_3^{(1)})\log (1-(h_{\mathrm{\Theta }}(x^{(1)}){)}_3)]& \text{对应m=1,k=3的情况}\end{array}$$

我们看这一部分

$$y_1^{(1)}\log ((h_{\mathrm{\Theta }}(x^{(1)}){)}_1)+y_2^{(1)}\log ((h_{\mathrm{\Theta }}(x^{(1)}){)}_2)+y_3^{(1)}\log ((h_{\mathrm{\Theta }}(x^{(1)}){)}_3)$$

，
这一部分可以看成向量$y^{(1)}$和向量$log(h_{\theta }(x^{(1)}))$的内积，而$h_{\theta }(x^{(1)})$是由第一个样本$x^{(1)}$带入到这个神经网络所代表的公式中计算出来的，而这里的计算用的“前馈传播算法”。同理，

$$(1-y_1^{(1)})\log (1-(h_{\mathrm{\Theta }}(x^{(1)}){)}_1)+(1-y_2^{(1)})\log (1-(h_{\mathrm{\Theta }}(x^{(1)}){)}_2)+(1-y_3^{(1)})\log (1-(h_{\mathrm{\Theta }}(x^{(1)}){)}_3)$$

这部分可以可以看成向量$1-y^{(1)}$（在matlab里面一个向量加减标量相当于这个向量里的每一个分量加减这个向量）以及向量$log(1-h_{\theta }(x^{(1)}))$的内积，就相当于上面那个公式里$y$替换成$1-y$，$h_{\theta }(x^{(i)})$替换成$1-h_{\theta }(x^{(i)})$。所以针对第一个样本来说，它对整体代价的贡献用向量形式来表达为：

$$y^{(1)}\cdot log(h_{\theta }(x^{(1)}))+(1-y^{(1)})\cdot log(1-h_{\theta }(x^{(1)}))$$

那么对于第$i$个样本，它对整体代价的贡献部分的向量表达形式为：

$$y^{(i)}\cdot log(h_{\theta }(x^{(i)}))+(1-y^{(i)})\cdot log(1-h_{\theta }(x^{(i)}))$$

那么对于全部样本的代价函数来说，我只要循环m次把结果累加就可以了。

2.3.4.3 想到的细节3：能否在同一个循环里实现cost值的计算和前馈传播算法

在上面讨论向量化的过程中，其中$h_{\theta }(x^{(i)})$的实现是需要利用前馈传播算法的，那么我是先对所有的样本跑一边前馈神经算法把所有的预测值（也就是$h_{\theta }(x^{(i)})$）计算出来，准备好，之后再计算cost值的时候，需要用到预测值的时候直接拿过来用？还是说我直接计算cost值（因为我的目标就是计算cost值）在需要用到的时候再去调用前馈神经算法计算需要的预测值。从累加的思想来看，其实每一个样本对应的预测值和其它样本没有关系嘛，我需要用的时候调用前馈神经算法计算一下行了吗，那么其实我走一遍循环（m次）就可以了，每次循环的过程中，我先调用前馈神经算法计算当前这个样本的预测值，之后在带入$y^{(i)}\cdot log(h_{\theta }(x^{(i)}))+(1-y^{(i)})\cdot log(1-h_{\theta }(x^{(i)}))$这个公式中计算出这个样本对整体代价的贡献程度，之后再累加起来。
这样只用一个循环，cost只用一个变量就可以了。到这里我们可以把计算代价函数的整体算法思想用伪代码的形式表现出来为：

J = 0
For i =1 to m
  计算出$h_{\theta }(x^{(i)})$
  $J=J+y^{(i)}\cdot log(h_{\theta }(x^{(i)}))+(1-y^{(i)})\cdot log(1-h_{\theta }(x^{(i)}))$
$J=J\times -\frac{1}{m}$

2.3.4.4 想到的细节4：

一开始我们是用提供好的参数值计算输出层的值，但是在后面几步，我们的目标是训练神经网络的参数，这个时候没有参数值了，那么前馈传播算法是不是不能用了，或者要修改一下代码呢？其实不用，并不是没有参数了，而只是没有“最优”参数了，每一步其实都有参数，最开始的参数是随即初始化的。

2.3.4.5 想到的细节5：自创的一些表达方法

我用 ${\phantom{a}}^{(1)}{a}^{(3)}$ 表示用第一个样本计算出神经网络第三层输出值，左上标表示样本的序号，右上角表示层数的序号。

2.3.4.6 nnCostFunction.m的具体落地

其实这个部分的程序只要上一次作业中 predict.m 里面修改就可以了。
这里不需要PredictMatrix这个数据结构，这个是为了存储这个5000个样本分别对应与0到9的概率，其实就是保存每一次循环输出层的结果了，但是这里我们是拿输出层的结果来计算cost值，所以不用保存。

接下来我们还发现一个问题就是：原有每个样本中的目标值（也就是$(x^{(i)},y^{(i)})$中的$y^{(i)}$），需要将其向量化。比如说 $y=1$ 对应 $y=(1,0,0,0,0,0,0,0,0,0{)}^T$这个向量。
那么如何进行向量化呢，需不需要提前把所有样本对应的向量准备出来呢？也没有必要，因为我们每次计算也是利用当前值，和之前与之后的都没有关系，那么这个问题就就是找到一个标量和一个向量之间的映射关系，这个关系就是：如果y=i，那么向量中的第i个分量就为1，其它都为0。那么我们可以先生成一个全零的向量，之后把当前的y值作为向量的索引，把向量中的第y个分量置为1。
用代码表示就是：

Y = zeros(10,1);
Y(y(i)) = 1;

综合以上的结论给出全部的具体代码是：

X = [ones(m, 1) X];

for i = 1 : 1 : m
    a1 = X(i,:);
    a2 = sigmoid(Theta1 * a1');
    a3 = sigmoid(Theta2 * [1;a2]);
    Y = zeros(10,1);
    Y(y(i)) = 1;
    J = J + dot(Y,log(a3))+dot((1-Y),log(1-a3));
end

但是在运行的时候出现这样的报错：

错误使用 dot (line 33)
A 和 B 的大小必须相同。

难道是dot这个函数的参数一个是行向量另外一个是列向量导致的？后来排除了这种可能，都是列向量。
我也不知道怎么灵机一动，我尝试用ex4.m这个大脚本来调用这个函数试试，结果能运行出来，但是结果是错的，结果是“-1438.145”，后来我重新看了一下公式忘了除以$-\frac{1}{m}$，最后加上一句变成

X = [ones(m, 1) X];

for i = 1 : 1 : m
    a1 = X(i,:);
    a2 = sigmoid(Theta1 * a1');
    a3 = sigmoid(Theta2 * [1;a2]);
    Y = zeros(10,1);
    Y(y(i)) = 1;
    J = J + dot(Y,log(a3))+dot((1-Y),log(1-a3));
end
J = J * (-1) / m ;

再次利用ex4.m调用就对了。这就奇怪了，为什么利用大脚本调用是正确的，而利用submit脚本调用是错的，于是看sumbit.m这个函数，发现这个是submit是为了验证程序正确而专门设计的一个测试用例，但是我还是觉得好像没有问题啊，于是单步调试，发现向量化过程中犯了了一个错误，就是把向量的纬度给写死了，我们大作业里面由于是要区分10个数字，所以向量的纬度是10，但是这个submit函数用例并没有那么复杂，他神经网络的架构是(2:4:4)这么一个结构，也就是测试用例结果向量只要4个纬度就行了，这样就会出现问题，所以在向量化的过程中要用变量来初始化一个全零向量。更改后的程序如下：

X = [ones(m, 1) X];

for i = 1 : 1 : m
    a1 = X(i,:);
    a2 = sigmoid(Theta1 * a1');
    a3 = sigmoid(Theta2 * [1;a2]);
    Y = zeros(num_labels,1);
    Y(y(i)) = 1;
    J = J + dot(Y,log(a3))+dot((1-Y),log(1-a3));
end
J = J * (-1) / m ;

这次submit就通过了。

2.3.5 nnCostFunction.m第二次实现：代价函数增加正则项

2.3.5.1 大体的思路

这一步的主要目的就是把正则项部分求出来之后和上面cost值相加。如何求正则项呢？

$$\begin{array}{c}J(\mathrm{\Theta })=-\frac{1}{m}\sum _{i=1}^m\sum _{k=1}^K[y_k^{(i)}\log ((h_{\mathrm{\Theta }}(x^{(i)}){)}_k)+(1-y_k^{(i)})\log (1-(h_{\mathrm{\Theta }}(x^{(i)}){)}_k)]+\frac{\lambda }{2m}\sum _{l=1}^{L-1}\sum _{i=1}^{s_l}\sum _{j=1}^{s_{l+1}}({\mathrm{\Theta }}_{j,i}^{(l)}{)}^2\end{array}$$

这个公式里指示了正则项的一般表达方式，如果仅就当前大作业的网络结构而言，正则项可以具体化为

$$\frac{\lambda }{2m}[\sum _{j=1}^{25}\sum _{k=1}^{400}({\mathrm{\Theta }}_{j,k}^{(1)}{)}^2+\sum _{j=1}^{10}\sum _{k=1}^{25}({\mathrm{\Theta }}_{j,k}^{(2)}{)}^2]$$

这样一看就非常容易找到规律了（这也是一种利用“具体化”的方式辅助思考从而找到规律的思路），正则项就是神经网络所有参数的平方和。了解到这点了下面就很容易了，就把参数遍历一遍，在这个过程中每个参数求一个平方，之后再进行累加。而${\mathrm{\Theta }}^{(1)}$和${\mathrm{\Theta }}^{(2)}$题目一开始就已经给了。于是我们立刻就可以写出循环遍历之后累加的方案。

    T1 = 0;
    T2 = 0;
    for i =1:1:hidden_layer_size
        for j = 1:1:input_layer_size
            T1 = T1 + Theta1(j,i) * Theta1(j,i);
        end
    end
    for i = 1:1:num_labels
        for j = 1:1:hidden_layer_size
            T2 = T2 + Theta2(j,i) * Theta2(j,i);
        end
    end
    J = J + (T1+T2)*lambda/(2*m);

2.3.5.2 向量化

我们看代码的一个片段

    for i =1:1:hidden_layer_size
        for j = 1:1:input_layer_size
            T1 = T1 + Theta1(j,i) * Theta1(j,i);
        end
    end

这个不就是第一个参数矩阵(${\mathrm{\Theta }}^{(1)}$)按行进行数递增的顺序把矩阵的每一个元素拿出来进行处理之后再相加么，用了matlab一段时间之后你再碰到循环语句就会自然而然的考虑有没有批量处理的替代方案（向量化），因为一个一个处理数据实在是太低效了嘛，我直接以矩阵为最小单位进行运算行不行，可以啊，这里面涉及到两种操作，一个是元素平方，可以理解为两个相同的矩阵进行hardmard乘积；另一个是把矩阵里的元素全部相加，可以利用matlab里面的sum操作（sum这个函数也是之前参考作业ex2网友的代码得到的启发，sum函数作用于向量的时候是吧向量里面的元素都加起来，sum函数作用于矩阵的时候，它会按列把每列元素进行求和得到一个行向量）
于是结果向量化的结果为

(sum(sum(Theta1 .* Theta1))+sum(sum(Theta2 .* Theta2))) * lambda / (2 * m );

2.3.5.3 代码的通用性

作业文档中提到代码要适应任意尺寸的${\mathrm{\Theta }}^{(1)}$和${\mathrm{\Theta }}^{(2)}$，我估计是为了用于测试的代码。但是其实不用它说我们都应该写一个通用的代码，而且不仅仅是三层的，甚至是多层的。
其实这个是非常好实现的，只要神经网络的结构是确定的那么实现一个通用代码是很容易的。我们不需要为每一层参数矩阵的计算结果准备一个变量，只需要一个变量通过累加的方式就可以了。层数用作最外层的循环即可。

2.3.5.4 其它要注意的细节

还是偏置的问题，作业文档中提到，在正则的时候需要把偏置项去掉，具体来说就是参数矩阵Theta1和Theta2的第一列去掉。因此最后的代码为：

T1 = Theta1(:,2:end);
T2 = Theta2(:,2:end);
J = J + (sum(sum(T1 .* T1))+sum(sum(T2 .* T2))) * lambda / (2 * m);

一些测试

>> a = [1,2,3,4;5,6,7,8;9,10,11,12]

a =

     1     2     3     4
     5     6     7     8
     9    10    11    12

>> b = a(2:end,:)

b =

     5     6     7     8
     9    10    11    12

>> c = a(:,2:end)

c =

     2     3     4
     6     7     8
    10    11    12

2.3.5.5 代码现状

X = [ones(m, 1) X];

for i = 1 : 1 : m
    a1 = X(i,:);
    a2 = sigmoid(Theta1 * a1');
    a3 = sigmoid(Theta2 * [1;a2]);
    Y = zeros(num_labels,1);
    Y(y(i)) = 1;
    J = J + dot(Y,log(a3))+dot((1-Y),log(1-a3));
end
J = J * (-1) / m ;

T1 = Theta1(:,2:end);
T2 = Theta2(:,2:end);
J = J + (sum(sum(T1 .* T1))+sum(sum(T2 .* T2))) * lambda / (2 * m);

2.4 sigmoid函数梯度

这一部分比较简单，公式已经给你了 $g^{\prime }(z)=g(z)(1-g(z))$，我们主要实现起向量形式。难点主要是想象g'(z)的向量形式是什么样的。我的绝招还是用一个小的数值看一下规律。
设$a$ = $\left(\begin{array}{c}1\\ 2\\ 3\end{array}\right)$，那么g(a) = $\left(\begin{array}{c}\frac{1}{1+e^{-1}}\\ \frac{1}{1+e^{-2}}\\ \frac{1}{1+e^{-3}}\end{array}\right)$，那么1-g(a) = $\left(\begin{array}{c}1-\frac{1}{1+e^{-1}}\\ 1-\frac{1}{1+e^{-2}}\\ 1-\frac{1}{1+e^{-3}}\end{array}\right)$，
设b =1 ，那么g(b)=$\frac{1}{1+e^{-1}}$ ,g'(b) = $\frac{1}{1+e^{-1}}(1-\frac{1}{1+e^{-1}})$ 那么 g'(a) = $\left(\begin{array}{c}\frac{1}{1+e^{-1}}(1-\frac{1}{1+e^{-1}})\\ \frac{1}{1+e^{-2}}(1-\frac{1}{1+e^{-2}})\\ \frac{1}{1+e^{-3}}(1-\frac{1}{1+e^{-3}})\end{array}\right)$ = g(a) .* (1-g(a))，也就是 g(a)和(1 - g(a))是“.*”的关系。
所以代码就是

g = sigmoid(z) .* (1 - sigmoid(z));

2.5 神经网络训练

ex4上说这部分叫做反向传播，其实是包括梯度下降的整个训练过程

其实整体的思路是非常简单的，首先先得到最后一层的误差，之后这个误差再通过神经网络反向往前传递（直观但是不严谨的理解方式是用$l$层的误差乘以$l$层的参数矩阵的转置得到上一层（$l-1$）的误差，但是还要乘以点东西），得到每一层的误差以后怎么办呢？我们干吗要得到每一层的误差，这是由于代价函数在当前这个样本对于参数的梯度可以用误差来表达，得到了单个样本的的梯度，把每个样本得到的梯度累加起来一平均就是整个代价函数对参数的梯度。

代码也很容易就写出来，在写代码的时候比较难的一点主要是在于如何给执行过程中所需要的变量进行命名，我认为应该把公式推导中所利用到的一些量都应该命名，这样代码比较好理解。
需要命名的变量有这些
首先${\delta }^{(l+1)}(a^{(l)}{)}^T$ 这个东西叫啥，这个就是代价函数在单个样本的情况下对${\mathrm{\Theta }}^{(l)}$的偏导，或者梯度（因为我都是向量操作），用数学公式表示就是${\mathrm{\nabla }}_{{\mathrm{\Theta }}^{(l)}}J(\mathrm{\Theta };x^{(t)},y^{(t)})$，这个东西在吴恩达这版教程里直接作为一个累加公式的一部分一笔带过了，我觉得这个东西还是很重要的，有必要单拿出来说，所以我给这部分起的名字叫做——SingleGradTheta，${\mathrm{\nabla }}_{{\mathrm{\Theta }}^{(1)}}J(\mathrm{\Theta };x^{(t)},y^{(t)})$为SingleGradTheta1，${\mathrm{\nabla }}_{{\mathrm{\Theta }}^{(1)}}J(\mathrm{\Theta };x^{(t)},y^{(t)})$为SingleGradTheta2，
另外${\mathrm{\Delta }}^{(l)}$如何命名，这个就是SingleGradTheta的累加啊，如果用首字母大写的Delta也是可以的，但是总觉得和误差delta会混淆，所以我这里用SumGradTheta来命名，其中对${\mathrm{\Theta }}_{(1)}$梯度的累和为SumGradTheta1，其它依次类推，那么这样写出来的代码就是

m = size(X, 1)
J = 0;
Theta1_grad = zeros(size(Theta1));
Theta2_grad = zeros(size(Theta2));
X = [ones(m, 1) X];

SumGradTheta1 = zeros(hidden_layer_size,input_layer_size);
SumGradTheta2 = zeros(num_labels,hidden_layer_size);


for i = 1 : 1 : m
    a1 = X(i,:);
    z2 = Theta1 * a1'; a2 = sigmoid(z2);
    z3 = Theta2 * [1;a2] ; a3 = sigmoid(z3); 
    Y = zeros(num_labels,1);
    Y(y(i)) = 1;

    delta_3 = a3 - Y;
    delta_2 =  Theta2' * delta_3 .* sigmoidGradient(z2);
    delta_2 = delta_2(2:end);

    SingleGradTheta2 = delta_3 * a2';
    SumGradTheta2 = SumGradTheta2 + SingleGradTheta2;
    SingleGradTheta1 = delta_2 * a1(2:end);  
    SumGradTheta1 = SumGradTheta1 + SingleGradTheta1;

    J = J + dot(Y,log(a3))+dot((1-Y),log(1-a3));
end
    Theta1_grad = SumGradTheta1 / m ;    
    Theta2_grad = SumGradTheta2 / m ;

J = J * (-1) / m ;
T1 = Theta1(:,2:end);
T2 = Theta2(:,2:end);
J = J + (sum(sum(T1 .* T1))+sum(sum(T2 .* T2))) * lambda / (2 * m);

整体代码很快就写出来了，但是从写完整体代码到能正确的输出结果确经过了比较长的调试过程，刚一开始主要的问题就是报错向量预算维度对不上的问题，这个就只能一点一点的trace看看到底哪里对不上了嘛，
我们看维度的变化
起始状态

变量	维度	备注
${\mathrm{\Theta }}^{(1)}$	25 x 401	原始数据 ‘ex4weights.mat’自带
${\mathrm{\Theta }}^{(2)}$	10 x 26	原始数据 ‘ex4weights.mat’自带

FP过程

变量	维度	备注
$a^{(1)}$	401 x 1	X = [ones(m, 1) X] X矩阵多了一列
$a^{(2)}$	25 x 1 ? 26 x 1 ?	如果单从$a^{(2)}=g(z^{(2)})$来看应该是25 x 1 ， 但是之后计算$z^{(3)}={\mathrm{\Theta }}^{(2)}{a}^{(2)}$确需要$a^{(2)}$维度为26 x 1
$a^{(3)}$	10 x 1	这个没有什么异议，因为最后输出就是10个节点

BP过程

变量	维度	备注
${\delta }^{(3)}=a^{(3)}-Y$	10 x 1	这个没有什么异议 输出层是最好确定的
${\delta }^{(2)}=({\mathrm{\Theta }}^{(2)}{)}^T{\delta }^{(3)}.\ast g^{\prime }(z^{(2)})$	$({\mathrm{\Theta }}^{(2)}{)}^T{\delta }^{(3)}$这部分 $(10\times 26{)}^T(10\times 1)=26\times 1$ $z^{(2)}={\mathrm{\Theta }}^{(1)}{a}^{(1)}$ 这部分 (25 x 401)(401 x 1)= 25 x 1	两部分维度不一样，如何统一维度？如果不统一就这个就没法进行下去
${\delta }^{(1)}$	没有

从这里可以看到在BP的过程中，再求${\delta }^{(2)}$的时候出现了计算公式里的两个构成因子的维度出现了不一致的情况，如果说在FP过程中，再求$a^{(3)}$的时候还可以利用“z3 = Theta2 * [1;a2]”的方法绕过到底要不要更改$a^{(2)}$的问题，但是在BP的过程里这个问题实在是绕不过去了。其实仔细观察这两个因子，$g^{\prime }(z^{(2)})$这部分能改吗？好像也没有改的空间，因为你无论怎么改都说不过去。
那么这个问题如何解决，解决的方式就是“扩枝和减枝”，也就是在FP的过程中，把$a^{(2)}$增加一个分量（a2 = [1;a2]），这个分量的值为1，这个叫做“扩枝”，因为从拓扑结构来看，这一层增加了一个叶子节点嘛；而在BP的过程中，计算${\delta }^{(2)}=({\mathrm{\Theta }}^{(2)}{)}^T{\delta }^{(3)}.\ast g^{\prime }(z^{(2)})$这个公式的时候，在计算$({\mathrm{\Theta }}^{(2)}{)}^T{\delta }^{(3)}$的时候要把这个计算结果向量的第一个分量去掉（delta_2_intermediate = delta_2_intermediate(2:end)），这个叫做“剪枝”，因为从拓扑结构来看，如同把一个叶子节点从树上摘掉了。
之所以这么做就是因为“把偏置单拎出来”这种处理方法对整个神经网络拓扑结构产生了微小的调整，进而“向前传播”以及“向后传播造”传递机制也需要根据这种拓扑结构产生相应的调整。
你看原来网络中（除了输入层和输出层）的任何一个节点既跟上一层的所有节点联系，又跟下一层的所有节点联系。而我增加了这个偏置节点之后，这个偏置节点只跟下一层的所有节点联系，而不跟上一层的节点发生联系。假设一个信息从输出层的某一个节点出发向输入层传播，那么当信息传递到偏置节点的时候，就进入了一个死胡同，无法继续向前了。同时信息从输入层往输出层方向传递的时候，如果仅用原有网络上的节点作为输入的节点是不够的，还要再加上偏置节点，才是完整的输入节点。
这样一个特征反映在算法实现的过程中就是——在实现ForwardPropagation的时候，在每一层准备发往下一层的输出数据的时候，要扩充一个输入为1的节点，如果用向量化运算的方式的话，就是这一层的输出向量要增加一个分量1，才是这一层完整意义的输出向量。而在实现backpropagation的时候，上一层($l+1$)的误差通过参数矩阵反向传递给下一层（$l$）的时候，在传递完成的时候（也就是$({\mathrm{\Theta }}^{(l)}{)}^T{\delta }^{(l+1)}$ 这个公式计算完成的时候）,需要把无法继续往下传导到下一层的那个节点去掉，也就是把向量中代表偏置节点的那个分量去掉，才能进行下一步的操作。
其实在作业描述文档的图片里已经暗示了这点操作，比如FP过程里hidden layer里有一小的注脚就是（add $a_0^{(2)}$）其实它暗示的就是 a2 = [1;a2]

而在BP那部分提示的配图中，在hidden layer中也有一个小脚注（reomve ${\delta }^{(2{)}_0}$）,其实我注意到这点，但是我以为的是求完了整个公式后再做这个操作，但是经过上面的分析发现在其中的一个部分，也就是$({\mathrm{\Theta }}^{(2)}{)}^T{\delta }^{(3)}$这个部分就要更新。

经过上述的分析，下面的就是修改代码，修改的地方主要有两处，一个是FP部分，在得到a2之后，增加代码a2 = [1;a2];之后要在BP过程中在计算delta_2中对$({\mathrm{\Theta }}^{(2)}{)}^T{\delta }^{(3)}$这部分进行更新，那么最好对$({\mathrm{\Theta }}^{(2)}{)}^T{\delta }^{(3)}$这部分起个名字吧，想了半天，我把它称之为“误差中间变量”，对因为它只是整个计算过程中的一环，所以是中间步骤，而又是一个提现了某种含义的变量，就是误差反向扩散的那一步，拿我之前的图来说明，就是绿圈的那个。这块我起名为delta_2_intermediate

多说一句，认识深入的过程是伴随这一系列“起名字”的过程的，其实这里确大量的知识细节需要补充，我称这些待补充的知识细节为“知识台阶”，你要越过这些台阶才能获得认识，而我“起名字”的过程就是为这些台阶建立脚手架的过程，起名字首先就把不同的事物区别开来，事物首先区别开来才能进一步深入认识，这就是分析的第一步，否则东西浑然一体，如何研究，这是一个学习心得。

改进后的代码如下：

m = size(X, 1)
J = 0;
Theta1_grad = zeros(size(Theta1));
Theta2_grad = zeros(size(Theta2));
X = [ones(m, 1) X];

SumGradTheta1 = zeros(hidden_layer_size,input_layer_size);
SumGradTheta2 = zeros(num_labels,hidden_layer_size);


for i = 1 : 1 : m
    a1 = X(i,:);
    z2 = Theta1 * a1'; a2 = sigmoid(z2); a2 = [1;a2];
    z3 = Theta2 * a2 ; a3 = sigmoid(z3); 
    Y = zeros(num_labels,1);
    Y(y(i)) = 1;

    delta_3 = a3 - Y;
    delta_2_intermediate = Theta2' * delta_3;
    delta_2_intermediate = delta_2_intermediate(2:end);
    delta_2 = delta_2_intermediate .* sigmoidGradient(z2);

    SingleGradTheta2 = delta_3 * a2';
    SumGradTheta2 = SumGradTheta2 + SingleGradTheta2;
    SingleGradTheta1 = delta_2 * a1(2:end);  %这个公式在第一层的时候要注意，不用再转置了，主要也取决于输入向量是列向量还是行向量，因为输入向量本来就是行向量
    SumGradTheta1 = SumGradTheta1 + SingleGradTheta1;

    J = J + dot(Y,log(a3))+dot((1-Y),log(1-a3));
end
    Theta1_grad = SumGradTheta1 / m ;    %  这块两个矩阵的形状不一样，居然能赋值他通过也是醉了
    Theta2_grad = SumGradTheta2 / m ;

J = J * (-1) / m ;
T1 = Theta1(:,2:end);
T2 = Theta2(:,2:end);
J = J + (sum(sum(T1 .* T1))+sum(sum(T2 .* T2))) * lambda / (2 * m);

这次倒是不报错了，但就是正确性检测通不过去，这是怎么回事。我这里是参考了其它人的代码，发现我用来保存梯度累和的变量的维度弄错了。也就是SumGradTheta1和SumGradTheta2的维度弄错了，SumGradTheta1和SumGradTheta2的维度分别等于权重矩阵${\mathrm{\Theta }}^{(1)}$和${\mathrm{\Theta }}^{(2)}$是对的，但是由于我们当前的模型是“把偏置单拎出来”的，所以权重矩阵${\mathrm{\Theta }}^{(1)}$和${\mathrm{\Theta }}^{(2)}$要在原来的基础上增加一列全1的一列。而其实程序一开始就把权重矩阵${\mathrm{\Theta }}^{(1)}$和${\mathrm{\Theta }}^{(2)}$为我们弄好了，就是在原始数据'ex4weights.mat'中就是这样。于是把SumGradTheta1和SumGradTheta2的size直接等于${\mathrm{\Theta }}^{(1)}$和${\mathrm{\Theta }}^{(2)}$就行了。
另外还有一个问题，原来我设的SumGradTheta1维度是25 x 400 ,而Theta1_grad的维度是 25 x 401，而

Theta1_grad = SumGradTheta1 / m ;

这句居然能过去。说明就是在matlab中，不同维度的矩阵可以相互赋值。下面这个实验可以验证这点

A = zeros(10,26);
B = ones(10,25);
A = B；

这次经过修改后的代码为：

X = [ones(m, 1) X];

SumGradTheta1 = zeros(size(Theta1_grad));
SumGradTheta2 = zeros(size(Theta2_grad));


for i = 1 : 1 : m
    a1 = X(i,:);
    z2 = Theta1 * a1'; a2 = sigmoid(z2); a2 = [1;a2];
    z3 = Theta2 * a2 ; a3 = sigmoid(z3); 
    Y = zeros(num_labels,1);
    Y(y(i)) = 1;

    delta_3 = a3 - Y;
    delta_2_intermediate = Theta2' * delta_3;
    delta_2_intermediate = delta_2_intermediate(2:end);
    delta_2 = delta_2_intermediate .* sigmoidGradient(z2);

    SingleGradTheta2 = delta_3 * a2';
    SumGradTheta2 = SumGradTheta2 + SingleGradTheta2;
    SingleGradTheta1 = delta_2 * a1;  
    SumGradTheta1 = SumGradTheta1 + SingleGradTheta1;

    J = J + dot(Y,log(a3))+dot((1-Y),log(1-a3));
end
    Theta1_grad = SumGradTheta1 / m ;
    Theta2_grad = SumGradTheta2 / m ;

J = J * (-1) / m ;
T1 = Theta1(:,2:end);
T2 = Theta2(:,2:end);
J = J + (sum(sum(T1 .* T1))+sum(sum(T2 .* T2))) * lambda / (2 * m);

这次终于通过了

最后我再总结一下，在本例中实现BP算法的难点在于“把偏置单拎出来”这个处理方法对于整个实现过程的影响，具体的来说它产生了两个后果
1. 它将对于权重以及偏置这两类参数的训练统一为对于权重的训练，简化了一个步骤
2. 它对整个网络的拓扑结构进行了微调，因此FP以及BP的传播方式也受到了影响，我们在用程序实现传播的过程中需要进行微调。
在之前代码报错的分析中我们已经详细的描述了第二点，现在我们补充描述一下“把偏置单拎出来”是如何把对偏置的训练转化为对网络参数的训练的。我们还是用一个简单的小模型来描述一下。
假设原始的模型是这样的

这个图形所代表的权重矩阵为$\mathrm{\Theta }$ = $\left(\begin{array}{ccc}{\theta }_{11}& {\theta }_{12}& {\theta }_{13}\\ {\theta }_{21}& {\theta }_{22}& {\theta }_{23}\end{array}\right)$。
当我把第二层偏置单拎出放到第一层上，就变成了下面这种结构，

第0号节点就是偏置节点，这个节点上的输出永远是1，由于增加的这个节点导致拓扑形态发生了变化，
（1）对于下层的输入来说，上一层的输出增加了一个分量1，a = [1;a]，
（2）权重矩阵需要增加一列，权重矩阵变为：$\mathrm{\Theta }$ = $\left(\begin{array}{cccc}{\theta }_{10}& {\theta }_{11}& {\theta }_{12}& {\theta }_{13}\\ {\theta }_{20}& {\theta }_{21}& {\theta }_{22}& {\theta }_{23}\end{array}\right)$。从新的权重矩阵来看，比原来的权重矩阵多了两个参数，其实这两个参数就是偏置以权重存在的表现，但是这个时候我们可以通过仅训练权重从而得到偏置。

2.6 梯度检验

代码并不需要自己实现，我简述一下我对于原理的理解，其实原理就是用一种近似计算来验证偏导的计算是否正确，我们知道导数的定义是通过极限来定义的，函数f(x)在$x_0$处的导数定义为

$$f^{\prime }(x_0)=\underset{\mathrm{\Delta }x\to 0}{lim}\frac{f(x_0+\mathrm{\Delta }x)-f(x_0)}{\mathrm{\Delta }x}$$

当$\mathrm{\Delta }x$特别小的时候(${10}^{(-4)}$) 极限可以近似的用算术表达式来替代：

$$f^{\prime }(x_0)=\frac{f(x_0+\mathrm{\Delta }x)-f(x_0-\mathrm{\Delta }x)}{2\mathrm{\Delta }x}$$

你要主要上面极限式是从两个方向逼近$x_0$所以下面的算术式是从$x_0$的两边两个增量端点进行计算的。

2.7 正则化神经网络

在神经网络训练中增加正则化和代价函数增加正则不一样，不要弄混，代价函数的正则化是在代价函数的计算公式里加上一个正则项，而神经网络的正则化是在计算梯度的公式中增加一个正则项。
我非常不喜欢coursera里对于一些表述方法比如

$$\frac{\mathrm{\partial }J(\mathrm{\Theta })}{\mathrm{\partial }{\mathrm{\Theta }}_{ij}^{(l)}}=D_{ij}^{(l)}=\frac{1}{m}{\mathrm{\Delta }}_{ij}^{(l)}\text{for j =0}$$

$$\frac{\mathrm{\partial }J(\mathrm{\Theta })}{\mathrm{\partial }{\mathrm{\Theta }}_{ij}^{(l)}}=D_{ij}^{(l)}=\frac{1}{m}{\mathrm{\Delta }}_{ij}^{(l)}+\frac{\lambda }{m}{\mathrm{\Theta }}_{ij}^{(l)}\text{for j }\ge \text{ 1}$$

中的$D_{ij}^{(l)}$和${\mathrm{\Delta }}_{ij}^{(l)}$这都是一些意味不明的表达方式，如果是我的话，我会用

$$\frac{\mathrm{\partial }J(\mathrm{\Theta })}{\mathrm{\partial }{\mathrm{\Theta }}_{ij}^{(l)}}=\frac{1}{m}\sum _{t=1}^m\frac{\mathrm{\partial }J(\mathrm{\Theta };x^{(t)},y^{(t)})}{\mathrm{\partial }{\mathrm{\Theta }}^{(l)}}+\frac{\lambda }{m}{\mathrm{\Theta }}^{(l)}$$

或者

$$\frac{\mathrm{\partial }J(\mathrm{\Theta })}{\mathrm{\partial }{\mathrm{\Theta }}_{ij}^{(l)}}=\frac{1}{m}\sum _{t=1}^m{\mathrm{\nabla }}_{{\mathrm{\Theta }}^{(l)}}J(\mathrm{\Theta };x^{(t)},y^{(t)})+\frac{\lambda }{m}{\mathrm{\Theta }}^{(l)}$$

其中${\mathrm{\Theta }}^{(l)}$是是吧第一列全部置0后的权重矩阵，因为bias不参加正则化。如果是分量形式的话，还要分$j=0$以及$j\ge 1$两种情况，如果利用向量形式，就把相应的权重矩阵的第一列全部置为0即可。
以上就是原理。那么写代码非常简单，就是在之前代码中的

    Theta1_grad = SumGradTheta1 / m ;
    Theta2_grad = SumGradTheta2 / m ;

这两句后面加上正则项就行了，改完的代码为：

    Theta1_grad = SumGradTheta1 / m + lambda / m * [zeros(size(Theta1,1),1) Theta1(:,2:end)];
    Theta2_grad = SumGradTheta2 / m + lambda / m * [zeros(size(Theta2,1),1) Theta2(:,2:end)];

难点主要是如何把原来的权重矩阵的第一列置为0，思路就是把原来权重矩阵第一列去掉（实际的实现是拿出从第二列到最后一列的数据Theta(:,2:end)），之后生成一列行数和原矩阵一样的全零的向量和去掉第一列的矩阵合在一起。
以Theta1为例
拿出从第二列到最后一列的数据利用Theta1(:,2:end)，
生成一列行数和矩阵一样的全零的向量zeros(size(Theta1,1),1),
两个矩阵合在一起 [zeros(size(Theta1,1),1) Theta1(:,2:end)]
最后的完整代码为：

X = [ones(m, 1) X];

SumGradTheta1 = zeros(size(Theta1_grad));
SumGradTheta2 = zeros(size(Theta2_grad));


for i = 1 : 1 : m
    a1 = X(i,:);
    z2 = Theta1 * a1'; a2 = sigmoid(z2); a2 = [1;a2];
    z3 = Theta2 * a2 ; a3 = sigmoid(z3); 
    Y = zeros(num_labels,1);
    Y(y(i)) = 1;

    delta_3 = a3 - Y;
    delta_2_intermediate = Theta2' * delta_3;
    delta_2_intermediate = delta_2_intermediate(2:end);
    delta_2 = delta_2_intermediate .* sigmoidGradient(z2);

    SingleGradTheta2 = delta_3 * a2';
    SumGradTheta2 = SumGradTheta2 + SingleGradTheta2;
    SingleGradTheta1 = delta_2 * a1;   
    SumGradTheta1 = SumGradTheta1 + SingleGradTheta1;

    J = J + dot(Y,log(a3))+dot((1-Y),log(1-a3));
end
    Theta1_grad = SumGradTheta1 / m + lambda / m * [zeros(size(Theta1,1),1) Theta1(:,2:end)];
    Theta2_grad = SumGradTheta2 / m + lambda / m * [zeros(size(Theta2,1),1) Theta2(:,2:end)];

J = J * (-1) / m ;
T1 = Theta1(:,2:end);
T2 = Theta2(:,2:end);
J = J + (sum(sum(T1 .* T1))+sum(sum(T2 .* T2))) * lambda / (2 * m);

来自为知笔记(Wiz)