參考與評述

參考書目《Deep Learning》Lan Goodfellow.
經典的深度學習框架是以計算圖&梯度下降方法實現對前饋網絡的有監督學習。
這裏復現了前饋計算圖的梯度計算實現。

一、前饋計算圖實現

1. 前向與梯度計算

• 結果數組 （保存輸入節點與計算節點的輸出值，能夠反映節點在計算方向的拓撲排序）
• 梯度數組 （保存輸入節點與計算節點的梯度，能夠反映節點在計算方向的拓撲排序）
• 連接圖 （反映每個節點的父節點）
• 輸出函數集合 （反映每個計算節點如何根據其輸入得到輸出）
• 梯度函數集合 （反映每個計算節點如何根據輸入和它的梯度計算對其任意父節點的梯度）

即可組成一個完整前饋計算圖。
我們以下圖所示全連接神經網絡為例構建計算圖（我們將每個神經元看作一個節點）。

(1)\---/(3)\---/(5)     X       X      (7)--(8)-->
(2)/---\(4)/---\(6)/

其中(1)(2)為輸入點，(3)(4)(5)(6)為隱層，(7)為輸出層，(8)為計算Loss的輸出。這樣，除了(1)(2)外其它都是計算節點。我們可以為該網絡定義一個包含8個元素的數組，以及計算關系。

數組及連接圖如下：

//結果數組與梯度數組,我們把輸入節點放到計算節點之前
list=[1] [2] [3] [4] [5] [6] [7] [8]
grad=[1] [2] [3] [4] [5] [6] [7] [8]
//父節點集合
Par[1]=null  Par[2]=null
Par[3]=1,2   Par[4]=1,2   //註意由於是前饋網絡，Pa(i)=j則i>j
Par[5]=3,4   Par[6]=3,4   
Par[7]=5,6   Par[8]=7
 

計算節點輸出函數集合如下：（以下輸出均為標量，輸入\(x\)為向量）

\[ Fun(i,x)=sigmoid(b_{i}+\sum_j w_{ij}*x_{j}),i\in [3,7] \Fun(8,x)=0.5 * (x-E)^2 \]

計算節點梯度函數集合如下：
註： dcFun(i,x,p) ，即求節點i在輸入x下的輸出對其父節點p的輸出的偏導數。

\[ dcFun(i,x,p)=sigmoid‘(b_{i}+\sum_j {w_{ij}*x_{j}})*w_{ip}, i \in [3,7] \dcFun(8,x,p)=(x-E) \]

先解釋下下面所用到的一些方法。

Par[i]->get_output_array(); //返回節點i的所有父節點的輸出列表，即i節點輸入向量
Par[i]->get_index_array();  //返回節點i的所有父節點的索引列表
L.has_element();            //列表中還有值存在
L.get_next_element();       //返回列表下一個值
put_input_in(list,i,j);     //向list的i到j索引處輸入值
 

前饋實現：

//前饋計算
put_input_in(list,1,2);
for(i=3;i<=8;i++)
    list[i]=Fun(i,Par[i]->get_output_array());

前饋傳播之後，為了最小化最終輸出list[8] ，即由(8)定義的損失函數輸出，我們需要計算每一節點的梯度。
反饋實現：

//反饋梯度計算
for(i=1;i<8;i++) grad[i]=0; //清空梯度數組
grad[8]=1;
for(i=8;i>=3;i--) //叠代每個計算節點，累加其各個父節點梯度
{
    input=Par[i]->get_output_array();
    par_array=Par[i]->get_index_array();
    while(par_array.has_element()) //叠代本節點的所有父節點
    {
        par_index=par_array.get_next_element();
        grad[par_index]+=grad[i]*dcFun(i,input,par_index); 
        //這裏dcFun(c,x,p)即求節點c在輸入x下的輸出對其父節點p的輸出的偏導數
    }
}

以計算（4）節點的輸出梯度為例。我們先得到其子節點（5）和（6）的梯度。則
\(grad[4] = grad[5]*(5對4輸出的偏導) + grad[6]*(6對4輸出的偏導)\)
如何求5對4輸出的偏導？我們由公式推導。
\(out5=sigmoid(out3*w_{51}+out4*w_{52}+b_{5})\)
\(我們令a5=out3*w_{51}+out4*w_{52}+b_{5}\)
\(\frac{\partial{out5}}{\partial{out4}}=sigmoid‘(a5)*w_{52}\)
這樣我們可以得到grad[4]的表達式。
\(grad[4] = grad[5]*(sigmoid‘(a5)*w_{52}) + grad[6]*(sigmoid‘(a6)*w_{62})\)
在反饋計算叠代到（5）（6）節點時，都會累加grad[4]這個值。
我們根據任意節點輸出的梯度，以及其輸入，就能調整這個節點的一些參數 。

2. 參數更新

在上一步中，我們得到了每個節點的輸出，以及Loss對每個節點輸出的梯度。
參數可以放在如上所示節點的內部，也可以單獨作為一個節點的輸出。
如果參數在節點內部：

\[ \frac{\partial L}{\partial A_W}= \frac{\partial L}{\partial A_{out}} \frac{\partial A_{out}}{\partial A_W} \]

由於A的輸出對其參數W的導數只由節點A決定，因此這些操作可以在計算所有梯度後並行執行。
如果參數由單個節點定義：

X  --A--->
    /     //參數放在W節點輸出中，這樣A就只是一個計算形式，需要計算材料X與W
   W      //這種形式在如Tensorflow框架中出現

這種形式中，我們將W當作一個輸入節點，這樣，在梯度計算時我們的A將會直接算出Loss對於W輸出的梯度。

二、全連接神經網絡

1. MLP前向計算

一個全連接神經網絡(MLP)可以當作一個整體作為計算圖中的一個計算節點，它有它的依賴，輸出方法，以及求父節點梯度的計算方法，權值更新方法。為了方便易用，也可以每一層當作一個計算節點。(PyTorch)
我們還可以將權值放到某個輸入節點中，為了區分它和輸入，把它定義成變量節點。(Tensorflow)

Require： 網絡深度l

Require： \(W_{i}\ ,\ i\in{1,...,l} （W_i的每一列表示i層某個神經元的全部權值）\)

Require： \(b_{i}\ ,\ i\in{1,...,l} （b_i表示i層各個神經元的偏置）\)

Require： X，程序輸入 （X每一行為一個輸入樣本，行數為多少批樣本，應用SGD）

Require： Y，目標輸出 （Y每一行為一個樣本對應標簽，行數為多少批樣本標簽）

H_{0}=X

for k=1:l do

\[ A_{k}=b_{k}+H_{k-1}W_{k} （A_k行數為批次數目，列數為k層神經元數目，這裏加法為行行加）\H_{k}=f(A_{k}) （對A_k逐元素計算非線性函數） \]

end for

$E=H_{l} $
\(L=Loss(E,Y)+\lambda\Omega(\theta) （\lambda\Omega(\theta)為正則項，用於優化學習）\)

2. 對向量偏導的定義：

\[ \frac{\partial y_m}{\partial x_n}= \begin{pmatrix} \frac{\partial y_1}{\partial x_1} & \frac{\partial y_1}{\partial x_2} & ... & \frac{\partial y_1}{\partial x_n}\ \frac{\partial y_2}{\partial x_1} & \frac{\partial y_2}{\partial x_2} & ... & \frac{\partial y_2}{\partial x_n}\ ...&...&...&... \ \frac{\partial y_m}{\partial x_1} & \frac{\partial y_m}{\partial x_2} & ... & \frac{\partial y_m}{\partial x_n}\\end{pmatrix}_{m \times n} \]

\[ \frac{\partial z}{\partial y_m}=(\frac{\partial z}{\partial y_1},\frac{\partial z}{\partial y_2},...,\frac{\partial z}{\partial y_m})^T \]

\[ \frac{\partial z}{\partial x_n}= (\frac{\partial y_m}{\partial x_n})^T(\frac{\partial z}{\partial y_m})= \begin{pmatrix} \frac{\partial y_1}{\partial x_1}\frac{\partial z}{\partial y_1} +\frac{\partial y_2}{\partial x_1}\frac{\partial z}{\partial y_2} +...\ \frac{\partial y_1}{\partial x_2}\frac{\partial z}{\partial y_1} +\frac{\partial y_2}{\partial x_2}\frac{\partial z}{\partial y_2} +...\ ...\ \frac{\partial y_1}{\partial x_n}\frac{\partial z}{\partial y_1} +\frac{\partial y_2}{\partial x_n}\frac{\partial z}{\partial y_2} +... \end{pmatrix}_n= \begin{pmatrix} \sum_{k=1}^m{\frac{\partial z}{\partial y_k}\frac{\partial y_k}{\partial x_1}} \\sum_{k=1}^m{\frac{\partial z}{\partial y_k}\frac{\partial y_k}{\partial x_2}} \... \\sum_{k=1}^m{\frac{\partial z}{\partial y_k}\frac{\partial y_k}{\partial x_n}} \\end{pmatrix}_n \]

3. 線性單元梯度計算：

已知 \(AB=C\) ， \(\frac{\partial L}{\partial C}=G\)，求 \(\frac{\partial L}{\partial A}與\frac{\partial L}{\partial B}\) ：（ \(L\)對某個矩陣\(X\)的偏導\(G\)的形式與\(X\)一模一樣）

\[ \begin{pmatrix} a_{11} & a_{12} & a_{13} \ a_{21} & a_{22} & a_{23} \\end{pmatrix}_{2 \times 3} \times \begin{pmatrix} b_{11} & b_{12} \ b_{21} & b_{22} \ b_{31} & b_{32} \\end{pmatrix}_{3 \times 2}= \begin{pmatrix} a_{11} b_{11}+a_{12} b_{21}+a_{13} b_{31} & a_{11} b_{12}+a_{12} b_{22}+a_{13} b_{32} \ a_{21} b_{11}+a_{22} b_{21}+a_{23} b_{31} & a_{21} b_{12}+a_{22} b_{22}+a_{23} b_{32} \\end{pmatrix}_{2 \times 2} \]

\[ \frac{\partial L}{\partial A}= \begin{pmatrix} g_{11} b_{11}+g_{12} b_{12} & g_{11} b_{21}+g_{12} b_{22} & g_{11} b_{31}+g_{12} b_{32} \ g_{21} b_{11}+g_{22} b_{12} & g_{21} b_{21}+g_{22} b_{22} & g_{21} b_{31}+g_{22} b_{32} \\end{pmatrix}_{2 \times 3}=G \times B^T \]

\[ \frac{\partial L}{\partial B}= \begin{pmatrix} g_{11} a_{11}+g_{21} a_{21} & g_{12} a_{11}+g_{22} a_{21} \ g_{11} a_{12}+g_{21} a_{22} & g_{12} a_{12}+g_{22} a_{22} \ g_{11} a_{13}+g_{21} a_{23} & g_{12} a_{13}+g_{22} a_{23} \\end{pmatrix}_{3 \times 2}=A^T \times G \]

對於偏置：已知 \(A+b=C,\frac{\partial L}{\partial C}=G\) 求 \(\frac{\partial L}{\partial A},\frac{\partial L}{\partial b}\)

\[ \begin{pmatrix} a_{11} & a_{12} & a_{13} \ a_{21} & a_{22} & a_{23} \\end{pmatrix}_{2 \times 3}+ \begin{pmatrix} b_{1} & b_{2} & b_{3} \\end{pmatrix} = \begin{pmatrix} a_{11}+b_1 & a_{12}+b_2 & a_{13}+b_3 \ a_{21}+b_1 & a_{22}+b_2 & a_{23}+b_3 \\end{pmatrix}_{2 \times 3} \]

\[ \frac{\partial L}{\partial A}= \begin{pmatrix} g_{11} & g_{12} & g_{13} \ g_{21} & g_{22} & g_{23} \\end{pmatrix}_{2 \times 3}=G \]

\[ \frac{\partial L}{\partial b}= \begin{pmatrix} g_{11}+g_{21} & g_{12}+ g_{22} & g_{13}+ g_{23} \\end{pmatrix}=G的行和 \]

4. 非線性單元梯度計算：

已知 \(f(A)=H,\frac{\partial L}{\partial H}=G\) 求 \(\frac{\partial L}{\partial A}\)

\[ H= \begin{pmatrix} f(a_{11}) & f(a_{12}) \ f(a_{21}) & f(a_{22}) \\end{pmatrix}_{2 \times 2}, \frac{\partial L}{\partial A}= \begin{pmatrix} g_{11}f‘(a_{11}) & g_{12}f‘(a_{12}) \ g_{21}f‘(a_{21}) & g_{22}f‘(a_{22}) \\end{pmatrix}_{2 \times 2}= G⊙f‘(A) \]

5. MLP梯度計算：

使用SGD算法，我們的輸入是一次一批的，標簽也是一次一批。

\(G\gets\nabla_{E}L\) （得到矩陣G，E為最後一層的輸出）

for \(k=l:1\) do ：

\[ G\gets\nabla_{A_k}L=G⊙f‘(a_k) ，（得到對線性單元輸出的梯度）\\nabla_{b_k}L=G行和，（得到對偏置的梯度） \\nabla_{W_k}=G \times H_{k-1}^T ，（得到對權值的梯度）\G\gets\nabla_{H_{k-1}}L=W_k^T \times G，（梯度傳播到父節點）\\]

end for

6. 優化線性單元：

為了增加代碼局部性來提高CPU運算速度，我們優化了計算方式。
與之前的矩陣乘法不同，權值矩陣行數為每個樣本輸出特征個數，列數與輸入特征個數相同。

\[ H=Linear( \begin{pmatrix} a_{00} & a_{01} \ a_{10} & a_{11} \ a_{20} & a_{21} \\end{pmatrix}_{3 \times 2} , \begin{pmatrix} w_{00} & w_{01} \ w_{10} & w_{11} \ w_{20} & w_{21} \ w_{30} & w_{31} \\end{pmatrix}_{4 \times 2}, \begin{pmatrix} b_{0} & b_{1} & b_{2} & b_{3} \\end{pmatrix} ) \= \begin{pmatrix} a_{00}w_{00}+a_{01}w_{01}+b_0 & a_{00}w_{10}+a_{01}w_{11}+b_1 & a_{00}w_{20}+a_{01}w_{21}+b_2 & a_{00}w_{30}+a_{01}w_{31}+b_3\ a_{10}w_{00}+a_{11}w_{01}+b_0 & a_{10}w_{10}+a_{11}w_{11}+b_1 & a_{10}w_{20}+a_{11}w_{21}+b_2 & a_{10}w_{30}+a_{11}w_{31}+b_3\ a_{20}w_{00}+a_{21}w_{01}+b_0 & a_{20}w_{10}+a_{21}w_{11}+b_1 & a_{20}w_{20}+a_{21}w_{21}+b_2 & a_{20}w_{30}+a_{21}w_{31}+b_3\\end{pmatrix}_{3 \times 4} \]

\[ \frac {\partial L}{\partial A}= \begin{pmatrix} g_{00}w_{00}+g_{01}w_{10}+g_{02}w_{20}+g_{03}w_{30} &g_{00}w_{01}+g_{01}w_{11}+g_{02}w_{21}+g_{03}w_{31} \ g_{10}w_{00}+g_{11}w_{10}+g_{12}w_{20}+g_{13}w_{30} &g_{10}w_{01}+g_{11}w_{11}+g_{12}w_{21}+g_{13}w_{31} \ g_{20}w_{00}+g_{21}w_{10}+g_{22}w_{20}+g_{23}w_{30} &g_{20}w_{01}+g_{21}w_{11}+g_{22}w_{21}+g_{23}w_{31} \\end{pmatrix}\=G \times W \]

\[ \frac {\partial L}{\partial W}= \begin{pmatrix} g_{00}a_{00}+g_{10}a_{10}+g_{20}a_{20} &g_{00}a_{01}+g_{10}a_{11}+g_{20}a_{21} \ g_{01}a_{00}+g_{11}a_{10}+g_{21}a_{20} &g_{01}a_{01}+g_{11}a_{11}+g_{21}a_{21} \ g_{02}a_{00}+g_{12}a_{10}+g_{22}a_{20} &g_{02}a_{01}+g_{12}a_{11}+g_{22}a_{21} \ g_{03}a_{00}+g_{13}a_{10}+g_{23}a_{20} &g_{03}a_{01}+g_{13}a_{11}+g_{23}a_{21} \\end{pmatrix} \=G^T \times A \]

\[ \frac {\partial L}{\partial b}= \begin{pmatrix} g_{00}+g_{10}+g_{20} & g_{01}+g_{11}+g_{21} & g_{02}+g_{12}+g_{22} & g_{03}+g_{13}+g_{23} \\end{pmatrix}\=G列和 \]

7. Batch Normalization：

引入BN層的作用是為了加速學習，試想如下網絡：

x   A-->    //二維坐標輸入(x,y)經過一個線性變換得到輸出（二分類）
y/

如果(x,y)的整個數據集是偏離原點很遠的一些點。由於在二分類問題中我們需要找到一條分割不同類點的直線，而初始化的bias(約等於0)表示的直線在靠近原點的地方，因此需要很多步叠代才能達到準確位置。同時，如果這些點非常密集，那麽得到的梯度就會非常微小，造成梯度彌散。

• 我們希望將這些點保留相對位置地移動到原點附近（也相當於將原點移動到這些點中心位置）。
  因此我們就把每個x減去x的平均值，每個y減去y的平均值。

• 我們還希望讓這些點有著歸一化的縮放大小（保持相對位置縮放到一個標準窗口內）
  因此我們就在讓上述處理後的每個x‘,y‘除以他們各自的均方差根（強行拉成標準正態分布）

• 還希望在這個歸一化標準基礎上，增加一點靈活性，以適應非線性函數
  因此我們加上一個乘性可學習參數、一個加性可學習參數。

因此，BN層一般都是加在非線性層之前，線性層之後，且BN層之前的線性層可以不要偏置（由於BN含偏置）。
每個神經元輸出都套一個BN塊，這個BN塊搜集一個batch內該神經元所有輸出，得到均值和均方差根，並且計算出標準化後的值到下一層。

對於一個batch中神經元A輸出的所有x

\[ m=\frac 1n \sum x_i \std=\frac 1n \sum (x_i-m)^2 \x‘_i \gets \frac {x_i-m}{\sqrt {std+eps}} \y_i \gets \gamma x‘_i +\beta \]

其中eps為一個很小值（1e-5）來防止分母為零
以下舉例對一個MLP輸出矩陣求BN
（每一行一個batch每一列一個特征,每一列就是一個神經元輸出）

\[ A= \begin{pmatrix} a_{00} & a_{01} \ a_{10} & a_{11} \ a_{20} & a_{21} \\end{pmatrix} \m_0=(a_{00}+a_{10}+a_{20})/3 \m_1=(a_{01}+a_{11}+a_{21})/3 \std_0=((a_{00}-m_0)^2+(a_{10}-m_0)^2+(a_{20}-m_0)^2)/3 \std_1=((a_{01}-m_1)^2+(a_{11}-m_1)^2+(a_{21}-m_1)^2)/3 \y_{ij}= \gamma_j a‘_{ij}+\beta_j =\gamma_j \frac{ (a_{ij}-m_j)}{\sqrt{std_j+eps} } +\beta_j \]

下面求梯度，使用鏈式法則:

\[ \frac {\partial L}{\partial std_j}=\frac{-\gamma_j }{2}(std_j+eps)^{-3/2} \sum_i g_{ij}(a_{ij}-m_j) \\frac {\partial L}{\partial m_j}=\sum_i g_{ij} \gamma_j \frac {-\sqrt{std_j+eps}-(a_{ij}-m_j)(1/2)(std_j+eps)^{-1/2}\frac{\partial std_j}{\partial m_j}}{std_j+eps} \=\frac{-\gamma_j}{\sqrt{std_j+eps}}\sum_i g_{ij}-\frac{\partial L}{\partial std_j}\frac1n \sum_i 2(a_{ij}-m_j) \\frac {\partial L}{\partial a_{ij}}=\sum_i g_{ij} \gamma_j \frac {(1-\frac 1n) \sqrt {std_j+eps}-(a_{ij}-m_j)\frac 12(std_j+eps)^{-1/2}\frac 2n(a_{ij}-m_j)(1-\frac 1n)}{std_j+eps} \=\frac 1n \frac{\partial L}{\partial m_j}+\frac {2(a_{ij}-m_j)}{n} \frac{\partial L}{\partial std_j} +\frac {g_{ij} \gamma_{j}}{\sqrt {std_j+eps}} \\frac {\partial L}{\partial \gamma_j}=\sum_i g_{ij}a‘_{ij} \\frac {\partial L}{\partial \beta_j}=\sum_i g_{ij} \]

三、損失函數

1. Softmax單元：

在分類問題中，神經網絡的輸出層屬於One-Hot類型。
比如手寫數字識別中只需識別0~9數字，那麽神經網絡的輸出層一共有10個神經元（對應10分類問題）。我們需要將這些輸出表示為對應概率，因此和必須=1且每個輸出大於0。這就需要通過Softmax層處理。

\[ S_i=\frac{e^i}{\sum_{k=0}^n e^k} \]

\[ \frac{\partial S_i}{\partial i}=\frac{e^i \sum e^k - e^i e^i}{(\sum e^k)^2} =\frac{e^i ( \sum e^k - e^i)}{(\sum e^k)^2}=S_i(1-S_i) \]

\[ \frac{\partial S_i}{\partial j}_{i \neq j}= -\frac{e^i e^j}{(\sum e^k)^2} = -S_iS_j \]

已知 \(S=softmax(A),\frac{\partial L}{\partial S}=G\) 求 \(\frac{\partial L}{\partial A}\)：
註意：A的每一行為一個樣本，因此Softmax是對每一行進行操作

\[ softmax \begin{pmatrix} a_{11} & a_{12} & a_{13} \ a_{21} & a_{22} & a_{23} \\end{pmatrix}_{2 \times 3} = \begin{pmatrix} \frac{e^{a11}}{e^{a11}+e^{a12}+e^{a13}} &\frac{e^{a12}}{e^{a11}+e^{a12}+e^{a13}} & \frac{e^{a13}}{e^{a11}+e^{a12}+e^{a13}} \ \frac{e^{a21}}{e^{a21}+e^{a22}+e^{a23}} &\frac{e^{a22}}{e^{a21}+e^{a22}+e^{a23}} & \frac{e^{a23}}{e^{a21}+e^{a22}+e^{a23}} \\end{pmatrix}_{2 \times 3} = \begin{pmatrix} s_{11} & s_{12} & s_{13} \ s_{21} & s_{22} & s_{23} \\end{pmatrix} \]

\[ \frac{\partial L}{\partial A}= \begin{pmatrix} s_{11}(g_{11}(1-s_{11})-g_{12}s_{12}-g_{13}s_{13}) & s_{12}(g_{11}s_{11}-g_{12}(1-s_{12})-g_{13}s_{13}) & ... \ s_{21}(g_{21}(1-s_{21})-g_{22}s_{22}-g_{23}s_{23}) & s_{21}(g_{21}s_{21}-g_{22}(1-s_{22})-g_{23}s_{23}) & ... \\end{pmatrix} \]

2. Cross-Entropy單元：

Cross-Entropy是分類問題常用的損失函數，配合Softmax使用，其輸入為一個表示不同類別識別概率的向量，以及One-Hot類型的標簽向量，輸出代價。

\[ L_{Cross-Entropy}=- \frac 1n \sum y \log a + (1-y) \log (1-a) \]

比如我們得到最終Softmax輸出的分類概率為 [0.8 , 0.1 , 0.1] One-Hot標簽為[1 , 0 , 0]。
\(Loss=-\frac 13 (\log 0.8+\log 0.9+\log 0.9)\)
對於一個batch中我們可以這樣計算：

\[ L=CrossEntropy( \begin{pmatrix} a_{11} & a_{12} & a_{13} \ a_{21} & a_{22} & a_{23} \\end{pmatrix}, \begin{pmatrix} 1 & 0 & 0 \ 0 & 1 & 0 \\end{pmatrix} ) \=-\frac 16(\log a_{11} + \log (1-a_{12}) + \log (1-a_{13}) \+ \log (1-a_{21}) + \log a_{22} + \log (1-a_{23})) \]

對於Cross-Entropy我們僅需對其輸入求導，對標簽矩陣無需求導。
求導十分方便，對矩陣A每個位置的求導，僅與該位置的原始數據和Y有關。
\(\frac{\partial L}{\partial a_{ij}}=y_{ij}是否為1? 是-\frac 1{na_{ij}}:否\frac 1{n(1-a_{ij})}\)

3. 分類問題的SCE單元：

\[ L=SCE( \begin{pmatrix} a_{00} & a_{01} & a_{02} \ a_{10} & a_{11} & a_{12} \\end{pmatrix}, \begin{pmatrix} 0 \ 1 \\end{pmatrix} ) \=-\frac 12 (\ln \frac {e^{a_{00}}}{e^{a_{00}}+e^{a_{01}}+e^{a_{02}}} + \ln \frac {e^{a_{11}}}{e^{a_{10}}+e^{a_{11}}+e^{a_{12}}}) \= \frac 12(\ln s_{00}+\ln s_{11}) \= -\frac 12( a_{00}-\ln (e^{a_{00}}+e^{a_{01}}+e^{a_{02}})+a_{11}-\ln (e^{a_{10}}+e^{a_{11}}+e^{a_{12}}) ) \\]

\(\frac{\partial L}{\partial a_{ij}}=y_{i}是否為j? 是\frac {s_{ij}-1}{n}:否\frac {s_{ij}}{n}\)

四、卷積神經網絡

1. 單通道圖像卷積：

我們用一個例子對圖A進行二維卷積運算，圖像為單通道，使用一個卷積核。

\[ Conv2d( \begin{pmatrix} a_{00} & a_{01} & a_{02} \ a_{10} & a_{11} & a_{12} \ a_{20} & a_{21} & a_{22} \\end{pmatrix}, \begin{pmatrix} w_{00} & w_{01} \ w_{10} & w_{11} \\end{pmatrix},b)\= \begin{pmatrix} a_{00}w_{00}+a_{01}w_{01}+a_{10}w_{10}+a_{11}w_{11}+b & a_{01}w_{00}+a_{02}w_{01}+a_{11}w_{10}+a_{12}w_{11}+b\ a_{10}w_{00}+a_{11}w_{01}+a_{20}w_{10}+a_{21}w_{11}+b & a_{11}w_{00}+a_{12}w_{01}+a_{21}w_{10}+a_{22}w_{11}+b\\end{pmatrix} \]

如果圖A高度Ah寬度Aw，卷積核高度Wh寬度Ww。
輸出圖高：Ah-Wh+1
輸出圖寬：Aw-Ww+1

for(i=0;i<Ah-Wh+1;i++)
    for(j=0;j<Aw-Ww+1;j++)
    {
        sum=0;
        for(di=0;di<Wh;di++)
            for(dj=0;dj<Ww;dj++)
                sum+=A[i+di][j+dj]*W[di][dj];
        Out[i][j]=sum;
    }

2. 單通道卷積梯度計算：

上面例子中，對W的偏導如下：

\[ \frac {\partial L}{\partial W}= \begin{pmatrix} g_{00}a_{00}+g_{01}a_{01}+g_{10}a_{10}+g_{11}a_{11} & g_{00}a_{01}+g_{01}a_{02}+g_{10}a_{11}+g_{11}a_{12} \ g_{00}a_{10}+g_{01}a_{11}+g_{10}a_{20}+g_{11}a_{21} & g_{00}a_{11}+g_{01}a_{12}+g_{10}a_{21}+g_{11}a_{22} \\end{pmatrix} \=Conv2d(G,A,0)，（G形狀與輸出一樣） \]

上面例子中，對A的偏導如下：

\[ \frac {\partial L}{\partial A}= \begin{pmatrix} g_{00}w_{00} & g_{00}w_{01}+g_{01}w_{00} & g_{01}w_{01} \ g_{00}w_{10}+g_{10}w_{00} & g_{00}w_{11}+g_{01}w_{10}+g_{10}w_{01}+g_{11}w_{00} & g_{01}w_{11}+g_{11}w_{01}\ g_{10}w_{10} & g_{10}w_{11}+g_{11}w_{10} & g_{11}w_{11} \\end{pmatrix} \]

看似十分淩亂，這裏運用了一個小技巧：
先對G做padding變為G‘：
高度方向：上下各增加Wh-1行0。
寬度方向：左右各增加Ww-1列0。
再對W做元素逆置：
W矩陣的一維排列：W00,W01,W10,W11 變換成： W11,W10,W01,W00 再整合成原形狀矩陣。
最終結果為：
\[ \frac {\partial L}{\partial A}=Conv2d(G‘, \begin{pmatrix} w_{11} & w_{10} \ w_{01} & w_{00} \\end{pmatrix},0) \]

來驗證一下尺寸是否正確：
Gw=Aw-Ww+1
G‘w=Aw-Ww+1+2(Ww-1)=Aw+Ww-1
Aw=G‘w-Ww+1

\[ \frac {\partial L}{\partial b}=G各元素和 \]

上面表明，我們同樣可以用卷積操作來求卷積核的梯度，十分方便。

3. Conv2d-Padding：

如果需要輸入圖像與輸出圖像形狀一樣，可以在進行卷積運算之前將輸入圖周圍補0。
為了盡可能將圖像放到中間。
如果Ww為奇數，Whfw=Ww/2（向下取整），圖像左右各加Whfw列0，行操作類似。
如果Ww為偶數，Whfw=Ww/2，圖像左加Whfw-1列0，右加Whfw列0，行操作類似。
在反傳梯度的時候，原來Padding過的位置不需要傳梯度。

4. 多通道多圖卷積：

如果一個圖像是三通道的，那麽每一個卷積核也應該為三通道（一個卷積核偏置只有一個）。
一個卷積核輸出：該卷積核各通道分別卷積圖像各通道，各通道輸出圖疊加，再每個位置加一個偏置。
對於100張2通道4*4圖形成的batch，用6個3*3卷積核操作的形狀如下：

\[ Conv2d( A_{100 \times 2 \times 4 \times 4},W_{6 \times 2 \times 3 \times 3},b_6)= O_{100 \times 6 \times 2 \times 2} \]

A1 A2 A3           K11 K12 K13      A^1 A^2 //A的3個通道分別和K1_,K2_的3個通道卷積
          Conv2d                =  
B1 B2 B3           K21 K22 K23      B^1 B^2 //B的3個通道分別和K1_,K2_的3個通道卷積
//A^1這樣形成：A1卷積K11+A2卷積K12+A3卷積K13+偏置1（標量）
//A^2這樣形成：A1卷積K21+A2卷積K22+A3卷積K23+偏置2（標量）  
//以上加法為矩陣對應位置相加

為了使輸入與輸出圖形狀一樣，對輸入圖padding，高寬分別補上2*Whfh與2*Whfw個0。

5. 多通道多圖卷積梯度計算：

用一個例子來說明：

//以下每一個元素都是一張圖
    I00 I01 I02
I=  I10 I11 I12  //四張三通道圖
    I20 I21 I22
    I30 I31 I32

W=  W00 W01 W02  //兩個三通道卷積核
    W10 W11 W12
  
    H00 H01
H=  H10 H11  //輸出四張兩通道圖（每一個通道對應一個卷積核輸出）
    H20 H21
    H30 H31

    G00 G01  //該矩陣每一個元素都是一張圖
G=  G10 G11  //矩陣高度為一個batch的樣本數
    G20 G21  //矩陣寬度為一個樣本的輸出特征圖多少，也即該層卷積核數目
    G30 G31

下面用類似 矩陣乘法 的操作來求多圖多通道卷積的梯度：

\[ Conv2dGrad\_I( \begin{pmatrix} G_{00} & G_{01} \ G_{10} & G_{11} \ G_{20} & G_{21} \\end{pmatrix}, \begin{pmatrix} W_{00} & W_{01} & W_{02} \ W_{10} & W_{11} & W_{12} \\end{pmatrix}) \\ = \begin{pmatrix} G_{00}W_{00}+G_{01}W_{10} & G_{00}W_{01}+G_{01}W_{11} & G_{00}W_{02}+G_{01}W_{12} \ G_{10}W_{00}+G_{11}W_{10} & G_{10}W_{01}+G_{11}W_{11} & G_{10}W_{02}+G_{11}W_{12} \ G_{20}W_{00}+G_{21}W_{10} & G_{20}W_{01}+G_{21}W_{11} & G_{20}W_{02}+G_{21}W_{12} \\end{pmatrix}) \\ =G \times W \\（乘法為求單圖梯度，加法為圖對應位置疊加） \]

\[ Conv2dGrad\_W( \begin{pmatrix} G_{00} & G_{01} \ G_{10} & G_{11} \ G_{20} & G_{21} \\end{pmatrix}, \begin{pmatrix} I_{00} & I_{01} & I_{02} \ I_{10} & I_{11} & I_{12} \ I_{20} & I_{21} & I_{22} \\end{pmatrix}) \\ = \begin{pmatrix} G_{00}I_{00}+G_{10}I_{10}++G_{20}I_{20} & G_{00}I_{01}+G_{10}I_{11}++G_{20}I_{21} & G_{00}I_{02}+G_{10}I_{12}++G_{20}I_{22} \ G_{01}I_{00}+G_{11}I_{10}++G_{21}I_{20} & G_{01}I_{01}+G_{11}I_{11}++G_{21}I_{21} & G_{01}I_{02}+G_{11}I_{12}++G_{21}I_{22} \\end{pmatrix} \\ = G^T \times I \\（乘法為求單卷積核梯度，加法為對應位置疊加） \]

\[ Conv2dGrad\_b( \begin{pmatrix} G_{00} & G_{01} \ G_{10} & G_{11} \ G_{20} & G_{21} \\end{pmatrix})= \begin{pmatrix} \sum G_{00} + \sum G_{10} + \sum G_{20} \ \sum G_{01} + \sum G_{11} + \sum G_{21} \\end{pmatrix} \（求和符號為所有元素累加） \]

6. 多步長卷積 ：

為了在卷積的同時對圖像進行降維，可以指定步長進行卷積。

stride=2; //要指定的步長
//得到輸出圖的形狀
Out_h=(Ah-Wh+1)%stride==0?(Ah-Wh+1)/stride:(Ah-Wh+1)/stride+1; //(Ah-Wh+1)/stride向上取整
Out_w=(Aw-Ww+1)%stride==0?(Aw-Ww+1)/stride:(Aw-Ww+1)/stride+1;
Out_1d=(float*)Out; //變成一維數組操作
index=0;
for(i=0;i<Ah-Wh+1;i+=stride)
    for(j=0;j<Aw-Ww+1;j+=stride) {
        sum=0;
        for(di=0;di<Wh;di++)
            for(dj=0;dj<Ww;dj++)
                sum+=A[i+di][j+dj]*W[di][dj];
        Out_1d[index++]=sum;
    }

7. 多步長卷積梯度計算 ：

對卷積核W梯度：

//CPU優化代碼如下
G_1d=(float*)G; //變成一維數組操作
for(i=0;i<Wh;i++)
    for(j=0;j<Ww;j++) {
        sum=0;
        index=0;
        for(di=0;di<Grad_h*stride;di+=stride)
            for(dj=0;dj<Grad_w*stride;dj+=stride)
                sum+=G_1d[index++]*A[i+di][j+dj]; //相當於卷積核G內間隔stride卷積
        Wgrad[i][j]+=sum;
    }
/*如需要並行操作，可以沿用Conv2d操作，具體如下：
1.創造一個Ah-Wh+1行Aw-Ww+1列零矩陣G‘
2.將G‘用G間隔stride-1行stride-1列填充（G平均分散成G‘且G‘(0,0)=G(0,0)）
3.對W梯度=Conv2d(A,G‘,0)*/

對輸入圖A梯度：

/*無明顯優化算法，統一用Conv2d操作，步驟如下：
1.創造一個Ah-Wh+1行Aw-Ww+1列零矩陣G‘
2.將G‘用G間隔stride-1行stride-1列填充（G平均分散成G‘且G‘(0,0)=G(0,0)）
3.對G‘進行如下Padding：高度方向上下各增加Wh-1行零，寬度方向左右各增加Ww-1列零，記G‘‘。
4.將卷積核W元素逆置變為W‘
5.對A梯度=Conv2d(G‘‘,W‘,0)*/
//對偏置b的梯度依然是G內所有元素相加

8. 池化層 ：

將圖像切成一塊一塊，求每塊最大值Max或平均值Avg做為輸出圖的相應點。
MaxPool作用：降維，減少對識別物體的平移敏感度。
AvgPool作用：降維，保持原有平移敏感度。

MaxPool 前向計算：

kernel_size=2; //池化方框邊長
stride=2;//池化步長
//得到輸出圖的形狀
Out_h=(Ah-kernel_size+1)%stride==0?(Ah-kernel_size+1)/stride:(Ah-kernel_size+1)/stride+1; 
Out_w=(Aw-kernel_size+1)%stride==0?(Aw-kernel_size+1)/stride:(Aw-kernel_size+1)/stride+1;
Out_1d=(float*)Out; //變成一維數組操作
Out_Max_Y=(int*)malloc(Out_h*Out_w*sizeof(int)); //保存最大值坐標方便梯度計算
Out_Max_X=(int*)malloc(Out_h*Out_w*sizeof(int));
index=0;
for(i=0;i<Ah-kernel_size+1;i+=stride)
    for(j=0;j<Aw-kernel_size+1;j+=stride) {
        maxi=i;maxj=j;maxm=A[i][j];
        for(pi=i+1;pi<i+kernel_size;pi++)
            for(pj=j+1;pj<j+kernel_size;pj++)
                if(A[pi][pj]>maxm){
                    maxi=pi;
                    maxj=pj;
                    maxm=A[pi][pj];
                }
        Out_Max_Y[index]=maxi;
        Out_Max_X[index]=maxj;
        Out_1d[index++]=maxm;
    }

MaxPool 梯度計算：

//梯度矩陣G與輸出矩陣大小相同，因此是被降維過的。
//要計算輸入矩陣A的梯度
//先得到每塊中哪個是最大值
//G中某個位置的值直接傳給對應塊最大值的位置，該塊中其他的梯度不傳遞值
G_1d=(float*)G; //變成一維度數組操作
index=0;
for(i=0;i<Ah-kernel_size+1;i+=stride)
    for(j=0;j<Aw-kernel_size+1;j+=stride) {
        Agrad[Out_Max_Y[index]][Out_Max_X[index]]+=G_1d[index]; //註意反向梯度計算是累加的
        index++;
    }
free(Out_Max_Y); //如果是動態內存需要回收
free(Out_Max_X);

AvgPool 前向計算：

kernel_size=2; //池化方框邊長
stride=2;//池化步長
//得到輸出圖的形狀
Out_h=(Ah-kernel_size+1)%stride==0?(Ah-kernel_size+1)/stride:(Ah-kernel_size+1)/stride+1; 
Out_w=(Aw-kernel_size+1)%stride==0?(Aw-kernel_size+1)/stride:(Aw-kernel_size+1)/stride+1;
Out_1d=(float*)Out; //變成一維數組操作
Out_Max_Y=(int*)malloc(Out_h*Out_w*sizeof(int)); //保存最大值坐標方便梯度計算
Out_Max_X=(int*)malloc(Out_h*Out_w*sizeof(int));
index=0;
for(i=0;i<Ah-kernel_size+1;i+=stride)
    for(j=0;j<Aw-kernel_size+1;j+=stride) {
        sum=0;
        for(pi=i;pi<i+kernel_size;pi++)
            for(pj=j;pj<j+kernel_size;pj++)
                sum+=A[pi][pj]
        Out_1d[index++]=sum/(kernel_size*kernel_size);
    }

AvgPool 梯度計算：

//G中某個位置的值，平均分配給對應塊各個位置
G_1d=(float*)G; //變成一維數組操作
index=0;
for(i=0;i<Ah-kernel_size+1;i+=stride)
    for(j=0;j<Aw-kernel_size+1;j+=stride) {
        for(pi=i;pi<i+kernel_size;pi++)
            for(pj=0;pj<j+kernel_size;pj++)
                Agrad[pi][pj]+=G_1d[index]/(kernel_size*kernel_size); //平均分配
        index++;
    }

9. Conv2d-BN-Layer ：

由於Batch Normalization是套在每個神經元輸出後面的，一個卷積核即一個神經元，因此操作如下：

A1 A2 //假設輸出為三張圖，每張圖兩個通道
B1 B2 //由於A1,B1,C1都是卷積核W1的輸出，A2,B2,C2都為卷積核W2的輸出
C1 C2 //因此有兩個BN節點，分別對應W1與W2
//1.計算(A1,B1,C1)上所有像素的平均值和均方差根，求出對應各像素BN輸出
//2.計算(A2,B2,C2)上所有像素的平均值和均方差根，求出對應各像素BN輸出
//3.計算梯度參考前面，操作集合確定後方法是一樣的。

五、深度殘差網絡

如果要訓練很深的網絡，會面臨一個微妙的問題：由於每次梯度傳播都是按順序從最高層到最底層，如果某個很低的層需要很高層的某個信息，那麽在此之間所有層都應該學習去保留這些信息，這增加了學習的復雜度。ResNet為了解決低層無法充分利用高層的有用信息(梯度)這一問題。

1. 殘差網絡基本單元 ：

X--A1--relu--A2-- add--relu-->
 \               /
  ---------------
//A1、A2為單層線性單元，relu非線性激活函數，add為對應元素相加
//註：可以用任意深度的非線性塊替代上圖中A1-relu-A2的位置（最後為線性層）

輸入X經過一個旁路直接與A2輸出疊加。
理論上，一個A1--relu--A2塊能夠逼近任意函數，姑且稱它為萬能塊。萬能塊並不好用，比如讓它逼近函數H(X)=X十分困難。並且如果要通過萬能塊傳梯度到X，勢必會經過一定衰減。
如果任務要求讓不加殘差網絡的萬能塊逼近H(X)，那麽殘差網絡的萬能塊轉化為逼近H(X)-X。
有如下兩個優勢：

• 因為X一般都是正的元素（如圖像像素，前面經過relu情況），並且H(X)中真正有用的信息必然擁有正的響應值（relu的特性使上層只把梯度傳給A2輸出為正的神經元，A2輸出為負的成為死神經元），因此H(X)-X在有用信息的方向上更接近原點，因此更容易學習。

• 傳到X的梯度為萬能塊與後面層梯度之和，如果萬能塊貢獻很小，那麽其傳下去的梯度幾乎不影響前面各層，該萬能層也會被漸漸忽略，即逼近H(X)，也即網絡有明顯的忽略一些層的能力。如果萬能塊與後面層貢獻相當，那麽前面層能學習到兩者共同的有用信息。

上述殘差網絡塊需要輸入與輸出形狀一致，事實上很難保證，我們可以通過線性變化讓其一致。

2. 殘差卷積單元 ：

為了將殘差結構應用到卷積神經網絡上，我們可以做成一個特殊的單元，在應用時疊加這些單元，下述兩個單元結構引用了論文 Deep Residual Learning for Image Recognition 中的結構。

對於較淺的層數<=34：

Node* ResidualCNNBlock(Node* X,int inchannel,int outchannel,int stride) {
    //輸入圖X(batch,inchannel,h,w)
    out=Conv2d(X,inchannel,outchannel,3,stride,padding=1,bias=False);//stride步3*3卷積
    out=BatchNorm2d(out);
    out=ReLU(out);
    //(batch,outchannel,h/stride,w/stride)
    out=Conv2d(out,outchannel,outchannel,3,1,padding=1,bias=False);
    out=BatchNorm2d(out);
    //萬能塊輸出(batch,outchannel,h/stride,w/stride)
    if(stride!=1||inchannel!=outchannel) { //輸入與輸出形狀不一致，需要線性變換
        sc_out=Conv2d(X,inchannel,outchannel,1,stride,padding=0,bias=False);
        sc_out=out=BatchNorm2d(sc_out);
        out=Add(out,sc_out);
    } else {
        out=Add(out,X);
    }
    out=ReLU(out); //最後一個非線性映射
    return out; //輸出(batch,outchannel,h/stride,w/stride)
}

對於較深的層數：

//即便輸入通道比較多，也可以用比較少的卷積核卷積，來節省計算
Node* ResidualCNNBlock(Node* X,int inchannel,int outchannel,int stride,int knum) {
    //(batch,inchannel,h,w)
    out=Conv2d(X,inchannel,knum,1,1,padding=0,bias=False);//1步1*1卷積降維
    out=BatchNorm2d(out);
    out=ReLU(out);
    //(batch,knum,h,w)
    out=Conv2d(out,knum,knum,3,stride,padding=1,bias=False);//stride步3*3卷積
    out=BatchNorm2d(out);
    out=ReLU(out);
    //(batch,knum,h/stride,w/stride)
    out=Conv2d(out,knum,outchannel,1,1,padding=0,bias=False);//1步1*1卷積升維
    out=BatchNorm2d(out);
    //(batch,outchannel,h/stride,w/stride)
    if(stride!=1||inchannel!=outchannel) { //輸入與輸出形狀不一致，需要線性變換
        sc_out=Conv2d(X,inchannel,outchannel,1,stride,padding=0,bias=False);
        sc_out=out=BatchNorm2d(sc_out);
        out=Add(out,sc_out);
    } else {
        out=Add(out,X);
    }
    out=ReLU(out); //最後一個非線性映射
    return out; //輸出(batch,outchannel,h/stride,w/stride)
}

3. 深度殘差卷積網絡案例 ：

由於當圖像尺寸減半時卷積核數目加倍，我們將相同卷積核數目的殘差網絡層綁定在一起，設計成一個塊：

Node* ResNetBlocks(Node* X,int inchannel,int outchannel,int stride,int depth,int knum) {
    //(batch,inchannel,h,w)
    out=ResidualCNNBlock(X,inchannel,outchannel,stride,knum);
    //(batch,outchannel,h/stride,w/stride)
    for(i=1;i<depth;i++)
        out=ResidualCNNBlock(out,outchannel,outchannel,1,knum);
    return out; //(batch,outchannel,h/stride,w/stride)
}

以下為官方的幾種設計案例：

//34-Layer:
out=Conv2d_Padding_BN_ReLU(X,Xchannel,kernel_size=7,outchannel=cnl);//預處理
out=ResNetBlocks(out,inchannel=cnl,outchannel=64,stride=2,depth=3);
out=ResNetBlocks(out,inchannel=64,outchannel=128,stride=2,depth=4);
out=ResNetBlocks(out,inchannel=128,outchannel=256,stride=2,depth=6);
out=ResNetBlocks(out,inchannel=256,outchannel=512,stride=2,depth=3);
out=GlobalAvgPool(out); //每張圖像尺寸降成1
out=MLP(infeature=512,outfeature=分類數目);
//50-layer,101-layer,152-layer:   (差別僅關註depth參數)
out=Conv2d_Padding_BN_ReLU(X,Xchannel,kernel_size=7,outchannel=cnl);//預處理
out=ResNetBlocks(out,inchannel=cnl,outchannel=256,stride=2,depth=3,knum=64);
out=ResNetBlocks(out,inchannel=256,outchannel=512,stride=2,depth=4,knum=128);//4,4,8
out=ResNetBlocks(out,inchannel=512,outchannel=1024,stride=2,depth=6,knum=256);//6,23,56
out=ResNetBlocks(out,inchannel=1024,outchannel=2048,stride=2,depth=3,knum=512);
out=GlobalAvgPool(out); //每張圖像尺寸降成1
out=MLP(infeature=2048,outfeature=分類數目);

深度學習計算框架實現