這里把按 [1] 推導的BP算法（Backpropagation）步驟整理一下。突然想整理這個的原因是知乎上看到了一個帥呆了的求矩陣微分的方法（也就是 [2]），不得不感嘆作者的功力。[1] 中直接使用矩陣微分的記號進行推導，整個過程十分簡潔。而且這種矩陣形式有一個非常大的優勢就是對照其進行編程實現時非常方便。
但其實用標量計算推導也有一定的好處，比如可以清楚地知道某個權重是被誰所影響的。
前向傳播過程：多層Logistic回歸
記號約定：
$L$：神經網絡的層數。輸入層不算。
$n^l$：第 $l$ 層神經元的個數。偏置神經元不算在內。
$W^{l}\in\mathbb R^{n^l\times n^{l-1}}$：第 $l-1$ 層到第 $l$ 層的權重矩陣。其中，$W_{ij}^{(l)}$ 表示第 $l-1$ 層第 $j$ 個神經元到第 $l$ 層第 $i$ 個神經元的權重。
$\textbf b^{(l)}\in\mathbb R^{n^l}$ ：第 $l$ 層的偏置。
$\textbf z^{(l)}\in\mathbb R^{n^l}$ ：第 $l$ 層各個神經元的輸入。
$f_l(\cdot)$ ：第 $l$ 層的激活函數。對於分類任務來說，最后一層為softmax。
$\textbf a^{(l)}\in\mathbb R^{n^l}$ ：第 $l$ 層各個神經元的輸出，也就是活性值。對於輸入層（第0層），$\textbf a^{(0)}=\textbf x$。
下圖給出了相應的圖示。
圖片來源：[1]
下面可以開始推導了。首先是前向傳播，計算網絡輸出：
$$\textbf z^{(l)}=W^{(l)}\textbf a^{(l-1)}+\textbf b^{(l)}$$
$$\textbf a^{(l)}=f_l(\textbf z^{(l)})=f_l(W^{(l)}\textbf a^{(l-1)}+\textbf b^{(l)})$$
$$\textbf a^{(L)}=f(\textbf x;W,\textbf b)$$
整個網絡實際上就可以看作是一個函數 $f$ ，1989年已經有人證明了單隱層網絡可以逼近任意函數。
既然是完成從輸入到輸出的映射，那么網絡除了可以看作分類器之外，還可以看作是特征提取器，再將網絡的輸出作為分類器的輸入，相當於把原特征 $\textbf x$ 映射成了新的特征 $\textbf a^{(L)}$ 。
如果把網絡看作分類器，那么最后一層的神經元個數 $n^L$ 應等於類別個數 $C$ ，且 $f_L(\cdot)=\text{softmax}(\cdot)$ ，歸一化成概率分布：
$$\hat{\textbf y}=\textbf a^{(L)}=\text{softmax}(\textbf z^{(L)})=\frac{\exp(\textbf z^{(L)})}{\textbf 1^{\top}\exp(\textbf z^{(L)})}\in\mathbb R^C$$
需要注意的是分子是列向量，分母是標量。$\hat{\textbf y}$ 的每一維 $\hat{\textbf y}_i$ 表示的是網絡 $f$ 給出的樣本 $\textbf x$ 屬於第 $i$ 類的概率。
給定樣本 $(\textbf x,\textbf y)$ ，其中 $\textbf y\in\mathbb R^C$ 是one-hot向量，那么使用交叉熵損失函數，得到網絡對一個樣本的損失為
$$\mathcal L(\textbf y,f(\textbf x;W,\textbf b))=-\textbf y^{\top}\ln\hat{\textbf y}$$
設訓練樣本數為 $N$ ，那么經驗風險為
$$R=\frac1N\sum_{i=1}^N\mathcal L(\textbf y_i,f(\textbf x_i;W,\textbf b))$$
如果使用權重矩陣的Frobenius范數作為正則項
$$||W||_F=\biggl(\sum_{l=1}^L\sum_{i=1}^{n^l}\sum_{j=1}^{n^{l-1}}(W_{ij}^{(l)})^2\biggr)^{\frac12}$$
那么結構風險為
$$R=\frac1N\sum_{i=1}^N\mathcal L(\textbf y_i,f(\textbf x_i;W,\textbf b))+\frac12\lambda||W||_F^2$$
一般來說，偏置是不加正則的。
現在考慮結構風險最小化，使用梯度下降法來更新權重矩陣和偏置。只需要求取模型對一個樣本的損失的梯度 $\dfrac{\partial \mathcal L(\textbf y,f(\textbf x;W,\textbf b))}{\partial W^{(l)}}$（因為對多個樣本的梯度就是對一個樣本的梯度的累加）：
$$W^{(l)}=W^{(l)}-\alpha \frac{\partial R}{\partial W^{(l)}}$$
$$\frac{\partial R}{\partial W^{(l)}}=\frac1N\sum_{i=1}^N\frac{\partial \mathcal L(\textbf y,f(\textbf x;W,\textbf b))}{\partial W^{(l)}}+\lambda W^{(l)}$$
對偏置的就不寫了。
反向傳播（Backpropagation，BP）算法
好了，墨跡了半天，重點終於來了。下面就是要求解 $\dfrac{\partial\mathcal L(\textbf y,f(\textbf x;W,\textbf b))}{\partial W^{(l)}}$ 、$\dfrac{\partial\mathcal L(\textbf y,f(\textbf x;W,\textbf b))}{\partial \textbf b^{(l)}}$ 。為了簡單起見，直接記作 $\dfrac{\partial \mathcal L}{\partial W^{(l)}}$ 、$\dfrac{\partial \mathcal L}{\partial \textbf b^{(l)}}$ 。
我所有的符號都是沿用的 [1] ，后面的推導過程也是。
首先擺上一些公式復習一下：
$$\frac{\partial A^{\top}\textbf x}{\partial\textbf x}=\frac{\partial \textbf x^{\top}A}{\partial\textbf x}=A$$
$$\frac{\partial \textbf y^{\top}\textbf z}{\partial\textbf x}=\frac{\partial \textbf y}{\partial\textbf x}\textbf z+\frac{\partial \textbf z}{\partial\textbf x}\textbf y$$
$$\frac{\partial \textbf y^{\top}A\textbf z}{\partial\textbf x}=\frac{\partial \textbf y}{\partial\textbf x}A\textbf z+\frac{\partial \textbf z}{\partial\textbf x}A^{\top}\textbf y$$
$$\frac{\partial y\textbf z}{\partial\textbf x}=\frac{\partial y}{\partial\textbf x}\textbf z^{\top}+y\frac{\partial \textbf z}{\partial\textbf x}$$
$$\frac{\partial \text{tr}AB}{\partial A}=B^{\top}\quad\quad\frac{\partial \text{tr}AB}{\partial A^{\top}}=B$$
$$\frac{\partial f(A)}{\partial A^{\top}}=(\frac{\partial f(A)}{\partial A})^{\top}$$
然后是鏈式法則（chain rule）：
$$\frac{\partial \textbf z}{\partial \textbf x}=\frac{\partial \textbf y}{\partial \textbf x}\frac{\partial \textbf z}{\partial \textbf y}$$
$$\frac{\partial z}{\partial X_{ij}}=(\frac{\partial z}{\partial\textbf y})^{\top}\frac{\partial\textbf y}{\partial X_{ij}}$$
$$\frac{\partial z}{\partial X_{ij}}=\text{tr}\biggl((\frac{\partial z}{\partial Y})^{\top}\frac{\partial Y}{\partial X_{ij}}\biggr)$$
以及逐元素計算的函數 $f$ （其導函數為 $f'$ ）及其梯度：
$$\frac{\partial f(\textbf x)}{\partial\textbf x}=\text{diag}(f'(\textbf x))$$
這里需要聲明一點，上面這套計算體系里，都是使用的分母布局，也就是說：維度為 $p$ 的列向量對維度為 $q$ 的列向量求導后得到的矩陣維數為 $q\times p$ 。關於矩陣求導，可以參考 [3] 。
又墨跡了半天，下面開始求 $\dfrac{\partial \mathcal L}{\partial W^{(l)}}$ 、$\dfrac{\partial \mathcal L}{\partial \textbf b^{(l)}}$ 。
BP的四個重要公式
下面這第一步在我認為是最重要的一步，因為它的思路很巧妙 —— 現在的目標是求標量對矩陣的導數，直接用chain rule拆成兩塊的話會發現：一塊是標量對向量的導數，好辦；還有一塊是向量對矩陣的導數，不好辦。所以 [1] 使用了一種很巧的方式，求標量對矩陣單個元素（標量）的導數，然后再歸納成矩陣形式。
根據鏈式法則，
$$\frac{\partial\mathcal L}{\partial W_{ij}^{(l)}}=(\frac{\partial\mathcal L}{\partial \textbf z^{(l)}})^{\top}\frac{\partial \textbf z^{(l)}}{\partial W_{ij}^{(l)}}$$
現在定義誤差項 $\delta^{(l)}=\dfrac{\partial\mathcal L}{\partial \textbf z^{(l)}}\in\mathbb R^{n^l}$ ，用來表征第 $l$ 層的神經元對於誤差的敏感程度。
先求 $\dfrac{\partial \textbf z^{(l)}}{\partial W_{ij}^{(l)}}$ ：
$$\frac{\partial \textbf z^{(l)}}{\partial W_{ij}^{(l)}}=\frac{\partial W^{(l)}\textbf a^{(l-1)}}{\partial W_{ij}^{(l)}}=(0,...,a_j^{(l-1)},...,0)^{\top}$$
其中只有第 $i$ 行的元素非零。既然這樣，那么
$$\frac{\partial\mathcal L}{\partial W_{ij}^{(l)}}=\delta_i^{(l)}a_j^{(l-1)}$$
從這個式子可以看出，如果神經元的輸出值 $a_j^{(l-1)}$ 很小，那么它連向下一層的權重 $W_{ij}^{(l)}$ 的更新就會很緩慢。
寫成矩陣的形式，就是下面的式子：
$$\frac{\partial\mathcal L}{\partial W^{(l)}}=\delta^{(l)}(\textbf a^{(l-1)})^{\top}$$
$$\frac{\partial\mathcal L}{\partial \textbf b^{(l)}}=\delta^{(l)}$$
下面的問題就是如何計算誤差項 $\delta^{(l)}$ 。對於誤差項的計算，需要分兩種情況考慮，一種情況是輸出層，另一種情況是隱層。
對於輸出層來說，
$$\begin{aligned}\delta^{(L)}&=\frac{\partial \mathcal L}{\partial \textbf z^{(L)}}\\&=\frac{\partial \textbf a^{(L)}}{\partial \textbf z^{(L)}}\frac{\partial \mathcal L}{\partial \textbf a^{(L)}}\\&=\text{diag}(f_L'(\textbf z^{(L)}))\frac{\partial \mathcal L}{\partial \textbf a^{(L)}}\\&=f_L'(\textbf z^{(L)})\odot \frac{\partial \mathcal L}{\partial \textbf a^{(L)}}\end{aligned}$$
對於隱層來說，
$$\begin{aligned}\delta^{(l)}&=\frac{\partial \mathcal L}{\partial \textbf z^{(l)}}\\&=\frac{\partial \textbf a^{(l)}}{\partial \textbf z^{(l)}}\frac{\partial\textbf z^{(l+1)}}{\partial \textbf a^{(l)}}\frac{\partial \mathcal L}{\partial \textbf z^{(l+1)}}\\&=\text{diag}(f_l'(\textbf z^{(l)}))(W^{(l+1)})^{\top}\delta^{(l+1)}\\&=f_l'(\textbf z^{(l)})\odot \bigl((W^{(l+1)})^{\top}\delta^{(l+1)}\bigr)\end{aligned}$$
從這個式子就可以看出，誤差項可以逐層往輸入的方向傳遞，這就是所謂的誤差反向傳播。
多提一句，這里出現了激活函數的導函數。有以下兩個常用關系：
$$\sigma'(\cdot)=\sigma(\cdot)\odot (\textbf 1-\sigma(\cdot))$$
$$\text{softmax}'(\cdot)=\text{softmax}(\cdot)\odot (\textbf 1-\text{softmax}(\cdot))$$
其中 $\sigma(\cdot)$ 代表logistic函數。
所以，BP算法歸結到最后其實只有四個式子比較重要：
$$\delta^{(L)}=f_L'(\textbf z^{(L)})\odot \frac{\partial \mathcal L}{\partial \textbf a^{(L)}}$$
$$\delta^{(l)}=f_l'(\textbf z^{(l)})\odot \bigl((W^{(l+1)})^{\top}\delta^{(l+1)}\bigr)$$
$$\frac{\partial \mathcal L}{\partial W^{(l)}}=\delta^{(l)}(\textbf a^{(l-1)})^{\top}$$
$$\frac{\partial \mathcal L}{\partial \textbf b^{(l)}}=\delta^{(l)}$$
BP算法與Delta規則的比較
1. 與二類感知器（只有一層，所以不需要看第二個式子）的Delta規則的比較：BP算法 $\delta^{(L)}$ 表達式中的 $f_L'(\textbf z^{(L)})$ 這項是額外多出來的。這也不難理解，因為二類感知器的最后是經過一個與0.5比較大小的閾值函數，這里則是Softmax函數（多分類）或Logistic函數（二分類）。
二類感知器中准則函數對參數的梯度就是用誤差（被錯誤分類的樣本為-1，被正確分類的樣本為0，對應到這里就是 $\dfrac{\partial \mathcal L}{\partial \textbf a^{(L)}}$ ）乘以樣本的特征（對應到這里就是上一層輸出 $\textbf a^{(L-1)}$ ）。
2. 另外可以看出，網絡的各個層的誤差項之間是相乘的關系，那么當層數很深的時候就可能會出現梯度消失／梯度爆炸的問題。
我也用標量的形式推導過，過程如下……
無論標量形式還是矩陣形式，關鍵就是一點 —— chain rule。
輸出層的誤差項 $\delta^{(L)}$求解
上面四個式子只是形式化地給出了反向傳播算法的梯度形式。下面求取：對於分類問題，也就是 $f_L(\cdot)=\text{softmax}(\cdot)$ ，並且使用交叉熵損失函數 $\mathcal L(\textbf y,f(\textbf x;W,\textbf b))=-\textbf y^{\top}\log\hat{\textbf y}$ 時，輸出層的誤差項 $\delta^{(L)}=\dfrac{\partial \mathcal L}{\partial \textbf z^{(L)}}$ 到底是個什么形式。
其實到這里問題已經很清楚了：這和Softmax回歸模型是一樣一樣的，這里的 $\delta^{(L)}$ 在Softmax回歸模型的求梯度過程中是一個中間結果。有三種方法：
（1）普通方法，按部就班地一步步推，可以參考我很早之前寫的一篇關於word2vec的博客；
（2）[1] 中的方法（第三章，在Softmax回歸那部分）；
（3）[2] 中的方法。
下面是基於 [2] 所給方法的完整求解過程。為了簡潔起見，所有的 $\textbf z^{(L)}$ 都記作 $\textbf z$：
[2]介紹的是這樣形式的求導：已知矩陣 $X$ ，函數 $f(X)$ 的函數值為標量，求 $\dfrac{\partial f}{\partial X}$ 。一種典型的例子就是求取損失對權重矩陣的導數。
對於一元微積分，$\text{d}f=f'(x)\text{d}x$ ；多元微積分，$\text{d}f=\sum_i\dfrac{\partial f}{\partial x_i}\text{d}x_i=(\dfrac{\partial f}{\partial \textbf x})^{\top}\text{d}\textbf x$；由此建立矩陣導數和微分的聯系：
$$\text{d}f=\sum_{i,j}\frac{\partial f}{\partial X_{ij}}\text{d}X_{ij}=\text{tr}((\frac{\partial f}{\partial X})^{\top}\text{d}X)$$
上式第二個等號成立是因為對於兩個同階方陣有 $\text{tr}(A^{\top}B)=\sum_{i,j}A_{ij}B_{ij}$ 。求解的流程就是，先求微分 $\text{d}f$ 表達式，然后再套上跡（因為標量的跡等於標量本身），然后再把表達式 $\text{tr}(\text{d}f)$ 和 $\text{tr}((\dfrac{\partial f}{\partial X})^{\top}\text{d}X)$ 進行比對，進而把 $\dfrac{\partial f}{\partial X}$ 給“挖”出來。
所以，問題就從求梯度轉化成了求微分。求微分當然少不了很多法則和技巧，下面隨着講隨着介紹。
首先求取 $\text{d} \mathcal L$ 。
$$\begin{aligned} \mathcal L&=-\textbf y^{\top}\ln\frac{\exp(\textbf z)}{\textbf 1^{\top}\exp(\textbf z)}\\&=-\textbf y^{\top}(\textbf z-\ln\begin{pmatrix}\textbf 1^{\top}\exp(\textbf z) \\ \textbf 1^{\top}\exp(\textbf z) \\ \vdots \\ \textbf 1^{\top}\exp(\textbf z)\end{pmatrix})\quad \textbf 1^{\top}\exp(\textbf z)\text{是標量}\\&=\ln(\textbf 1^{\top}\exp(\textbf z))-\textbf y^{\top}\textbf z\end{aligned}$$
根據法則 $\text{d}(g(X))=g'(X)\odot\text{d}X$ 、$\text{d}(XY)=(\text{d}X)Y+X(\text{d}Y)$，可得
$$\text{d}(\ln(\textbf 1^{\top}\exp(\textbf z)))=\frac{1}{\textbf 1^{\top}\exp(\textbf z)}\odot\text{d}(\textbf 1^{\top}\exp(\textbf z))$$
$$\text{d}(\textbf 1^{\top}\exp(\textbf z))=\textbf 1^{\top}\text{d}(\exp(\textbf z))=\textbf 1^{\top}(\exp(\textbf z)\odot\text{d}\textbf z)$$
所以
$$\text{d} \mathcal L=\frac{\textbf 1^{\top}(\exp(\textbf z)\odot\text{d}\textbf z)}{\textbf 1^{\top}\exp(\textbf z)}-\textbf y^{\top}\text{d}\textbf z$$
現在可以套上跡，根據恆等式 $\text{tr}(A^{\top}(B\odot C))=\text{tr}((A\odot B)^{\top}C)=\sum_{i,j}A_{ij}B_{ij}C_{ij}$ ，可得
$$\begin{aligned}\text{d} \mathcal L&=\text{tr}(\frac{(\textbf 1\odot \exp(\textbf z))^{\top}\text{d}\textbf z}{\textbf 1^{\top}\exp(\textbf z)})-\text{tr}(\textbf y^{\top}\text{d}\textbf z)\\&=\text{tr}(\biggl(\frac{(\exp(\textbf z))^{\top}}{\textbf 1^{\top}\exp(\textbf z)}-\textbf y^{\top}\biggr)\text{d}\textbf z)\\&=\text{tr}((\hat{\textbf y}-\textbf y)^{\top}\text{d}\textbf z)\\&=\text{tr}((\frac{\partial\mathcal L}{\partial\textbf z})^{\top}\text{d}\textbf z)\end{aligned}$$
這樣就得到了 $\delta^{(L)}$ 的形式：
$$\delta^{(L)}=\frac{\partial \mathcal L}{\partial \textbf z^{(L)}}=\hat{\textbf y}-\textbf y$$
這樣就不難看出為什么將 $\delta^{(l)}$ 稱為誤差項了。
學習率自適應
眾所周知，現在有很多算法用來改進SGD，以求更好的訓練效果（時間上 & 精度上）。下面這張圖總結了包括加動量、自適應學習率等眾多改進算法的更新公式，如果想詳細學習其中的理論，可以看 [1] ，也可以看Yoshua Bengio的那本最新出版的《Deep Learning》。對於這些算法，我認為沒有最好，只有最合適。在自己的任務上何者可以取得最滿意的效果，還是要多跑幾組實驗來確定。反正在TensorFlow里就是改個函數的事情，但是挑出最適合自己任務的那一個，還是挺不容易的。另外，諸如什么“先用 Adam 再用 SGD ”之類的我是沒有試過（我都是直接用Adam，逃。。。
圖片來源：http://weibo.com/1657470871/EACKVoLTI?type=comment#_rnd1490771184496
從這開始的內容與標題無關……
之前做過一個小實驗：裸寫一個單隱層神經網絡（實在是很簡單），隱層激活函數使用tanh、損失函數使用平方損失函數，樣本特征為三維、共三個類，每類10個點。只做訓練用，沒有測試集。
1. 改變隱層節點個數：使用批處理梯度下降（因為數據量小到可以忽略），迭代500輪后，訓練集上的loss變化（圖里iters改成epoch更合適）如下圖，可以看出來網絡表示能力大體上增強（由於隨機初始化的原因，每次實驗結果可能不同）。
2. 改變學習率：較大的學習率在早期會讓loss下降很快，但是后期可能會震盪。
參考資料：
[1] 《神經網絡與深度學習講義》
[2] 《矩陣求導術（上）》
[3] Matrix_calculus