deep learning 自編碼算法詳細理解與代碼實現（超詳細）

本文轉載自查看原文 2015-01-07 22:27 22594

在有監督學習中，訓練樣本是有類別標簽的。現在假設我們只有一個沒有帶類別標簽的訓練樣本集合 $\textstyle \{x^{(1)}, x^{(2)}, x^{(3)}, \ldots\}$ ，其中 $\textstyle x^{(i)} \in \Re^{n}$ 。自編碼神經網絡是一種無監督學習算法，它使用了反向傳播算法，並讓目標值等於輸入值，比如 $\textstyle y^{(i)} = x^{(i)}$ 。下圖是一個自編碼神經網絡的示例。通過訓練，我們使輸出 $\textstyle \hat{x}$ 接近於輸入 $\textstyle x$ 。當我們為自編碼神經網絡加入某些限制，比如限定隱藏神經元的數量，我們就可以從輸入數據中發現一些有趣的結構。舉例來說，假設某個自編碼神經網絡的輸入 $\textstyle x$ 是一張張8*8 圖像（共64個像素）的像素灰度值，於是 n=64，其隱藏層 $\textstyle L_2$ 中有25個隱藏神經元。注意，輸出也是64維的。由於只有25個隱藏神經元，我們迫使自編碼神經網絡去學習輸入數據的壓縮表示，也就是說，它必須從25維的隱藏神經元激活度向量中重構出64維的像素灰度值輸入 $\textstyle x$ 。如果網絡的輸入數據是完全隨機的，比如每一個輸入 $\textstyle x_i$ 都是一個跟其它特征完全無關的獨立同分布高斯隨機變量，那么這一壓縮表示將會非常難學習。但是如果輸入數據中隱含着一些特定的結構，比如某些輸入特征是彼此相關的，那么這一算法就可以發現輸入數據中的這些相關性。網絡訓練好以后，每一個輸入 $\textstyle x_i$ 對應的LayerL2 就相當於降維后的數據。(跟pca理解差不多，只是pca 是線性降維，這里因為有sigmoid非線性激活函數的原因，所以這里應該可看做是非線性的降維。）

現在我們實驗一下：我們有十張64*64的圖片,十張圖片如下：

我們從這十張圖片中隨機取出10000個小碎片，小碎片是8*8的小圖片（叫做patch，或者掩模），這10000個patch就是10000個樣本 $\textstyle x_i$ ，我們很容易理解這10000個樣本之間肯定有聯系。那么接下來我們就實現這個自編碼算法。

我們先得到這10000個樣本（注意樣本都約束在了[0.1,0.9]之間，不是說只要用神經網絡就要把數據歸一化把數據弄到0~1之間，而是因為這里是自編碼算法，神經網絡最后輸出的應該是原樣本，如果你樣本向量的元素有大於1的數，最后輸出怎么也不可能等於原樣本，因為sigmoid函數最大也才是趨近於1，所以這里需要先把數據約束到[0.1,0.9]之間。如果不是自編碼算法，那么是可以不用這樣歸一化到0~1之間的，你也不用擔心如果樣本的數值過大會使sigmoid函數輸出都是1，可能一開始你任意取的參數的時候會產生第二層很多神經元的激活值都是1，但是經過不斷的迭代，如果樣本值很大，那么參數w就會很小，最后肯定不會激活值都是1，並且由於加了稀疏項，第二層神經元激活值都會很小，這些會通過訓練調節w的大小還調節激活值的大小，進而調節最后輸出的大小。但注意最后神經網絡的輸出都在（0,1）之間）。

得到樣本的程序如下：

   1: function patches = sampleIMAGES()

   2: % sampleIMAGES

   3: % Returns 10000 patches for training

   4: load IMAGES;    % load images from disk IMAGES是一個512*512*10的數組。

   5: patchsize = 8;  % we'll use 8x8 patches

   6: numpatches = 10000;

   7: patches = zeros(patchsize*patchsize, numpatches);%這個patches是一個矩陣

   8: %每一列是一個64的列向量，是把每一個8*8的掩模變成這樣一個向量，然后這個矩陣有

   9: %10000個列。也就是說從10幅512*512的圖像中弄出來10000個掩模。現在初始化它。

  10: %% 從10幅512*512的圖像中隨機選出10000個8*8的掩模，然后制作 patches

  11: image_size=size(IMAGES);

  12: i=randi(image_size(1)-patchsize+1,1,numpatches);

  13: %在開區間（0，image_size(1)-patchsize+1）之間生成一個1行10000列的向量

  14: %也就是在一副512*512的圖像中的行號中（不包括最后的7行，因為選了也構不成8*8的掩模了）

  15: %隨意選出10000個數字組成向量A。下面這行是在列號中選10000個元素的向量B，

  16: %這樣每一次從A和B中選出一個數a，b，然后a+7，b+7，這樣a，b，a+7，b+7就組成了一個

  17: %8*8的掩模。（注意是選擇哪個地方的掩模之前，要先從十幅圖像中選一幅圖像)

  18: %然后依次這么做，就得到了10000個掩模，然后把每個掩模弄成8*8的向量，

  19: %依次排到 patches的每一列上就得到了 patches。

  20: j=randi(image_size(2)-patchsize+1,1,numpatches);

  21: k=randi(image_size(3),1,numpatches);

  22: for num=1:numpatches

  23:     patches(:,num)=reshape(IMAGES(i(num):i(num)+patchsize-1,j(num):j(num)+patchsize-1,k(num)),1,patchsize*patchsize);

  24: end

  25: %% ---------------------------------------------------------------

  26: % For the autoencoder to work well we need to normalize the data

  27: % Specifically, since the output of the network is bounded between [0,1]

  28: % (due to the sigmoid activation function), we have to make sure

  29: % the range of pixel values is also bounded between [0,1]

  30: patches = normalizeData(patches);

  31: display_network(patches(:,randi(size(patches,2),200,1)),8);%這個是隨機從那10000個patch里面選出來200個patch，看看是不是

  32: %弄對了。

  33: end %這個end是說明patches這個function結束了，因為下面還有一個function，在一個文件里，所以一個function完了必須寫end

34:

  35: %% ---------------------------------------------------------------

  36: function patches = normalizeData(patches)

  37: % Remove DC (mean of images).

  38: patches = bsxfun(@minus, patches, mean(patches));

  39: % mean(patches)對patches每一列算出一個平均值。也就是每一組樣本算出一個平均值。

  40: %mean(patches)是一個1*10000的行向量，這里 patches是一個64*10000的矩陣，按理說

  41: %不能相減，但這里 用到了bsxfun（@minus,），那么mean(patches)會自動復制64行，

  42: %然后也變成一個64*10000的矩陣，然后相減。

43:

  44: % Truncate to +/-3 standard deviations and scale to -1 to 1

  45: pstd = 3 * std(patches(:));%patches(:)是把這個矩陣變成一個向量。然后std，也就是

  46: %求所有10000個掩模所有像素值的標准差，這里如果你不是很明白為什么要把所以樣本一塊求標准差，

  47: %而不是一個樣本組一個樣本組的求標准差。並且還最后把所有像素的標准差*3。解釋一下，以前求

  48: %數據都是一個樣本組一個樣本組進行歸一化，那是把數據進行歸一化，也就是均值為0，方差為1，但是

  49: %這樣並不能保證所有的數據都在[-1,1]之間，只是方差為1。但是我們知道我們把所有數據減去均值

  50: %的話，那么所有數據的99.7%都會落在[-3*標准差，3*標准差]之間，所以我們只需要把剩下的0.03%的數據都

  51: %置成-3*標准差或3*標准差即可。這樣所有數據都在[-3*標准差，3*標准差]之間，然后我們除以

  52: %3*標准差，那么素有數據都會在[-1,1]之間了。注意這里是將所有數據一塊處理，而不需要一組一組樣本處理，

  53: %因為這里並不需要處理后的每個樣本滿足標准正態分布:均值為0，方差為1那種關系，也不可能滿足，所以

  54: %數據一塊處理即可。

  55: patches = max(min(patches, pstd), -pstd) / pstd;

  56: %min(patches, pstd)意思是patches這個矩陣中元素如果大於pstd，那么這個元素就是prsd。

  57: %也就是上限是pstd，然后這里再max，也就是矩陣元素的下限是-petd。

58:

  59: % Rescale from [-1,1] to [0.1,0.9]

  60: patches = (patches + 1) * 0.4 + 0.1;

  61: end

隨機選出200個樣本/pitch（也就是從那10幅圖像中隨機取出的8*8的小碎片）可以看一下：

接下來要做的就是求出cost function函數，那么這里的cost function是什么那？我們已經從上一節中知道cost function（加了L2正則項）：

那么我們通過梯度下降法就可以得到收斂的參數。那么關鍵就是求cost function的偏導數。

並且我們知道代價函數 $\textstyle J(W,b)$ 對每一個參數的偏導數可通過如下公式求：

……………………………………………….………(1)

而每一個樣本的 $\textstyle \frac{\partial}{\partial W_{ij}^{(l)}} J(W,b; x, y)$ 和 $\textstyle \frac{\partial}{\partial b_{i}^{(l)}} J(W,b; x, y)$ 可以通過如下公式求出：

………………………………………………………………………………………………….(2)

而即為每個單元的殘差。為每個單元的激活值，所以說，利用前向算法和后向算法就可以求出每個單元的激活值和殘差。

那么任務就完成了。

但是這里我們在剛才cost function中還想添加一項。正如我們添加的L2正則項為了防止過擬合產生一樣，我們添加的這一項是為了讓隱藏層的激活值大多數都在0附近。（注意L2正則項會使大部分參數在0附近，而不是激活值）這樣做有什么好處呢？我們這里定義如果當神經元的輸出接近於1的時候我們認為它被激活，而輸出接近於0的時候認為它被抑制，這里我們限制讓隱藏層的激活值大部分都在0附近，也就是說訓練好以后，當輸入一個 $\textstyle x_i$ 后，進行前向算法大部分的隱藏層神經元都是處於被抑制狀態，只是很少的神經元被激活（你可能會說現在的輸入不都限制在[0.1,0.9]之間了嗎，而我們知道sigmoid函數的圖像如下圖，[0.1,0.9]對應的函數值不可能接近1啊？你忽略了一點，這里sigmoid輸入是前一層所有神經元的線性疊加，而不是只有一個，[0.1,0.9]的值輸入，所有激活函數肯定是可以接近1的。）

回答一下為什么要加這種洗屬性約束呢？因為這個網絡是模擬的人腦，而人腦也是一個神經元分層工作的系統，進來一個圖像后，首先第一層的神經元開始工作，因為神經元被激活是需要能量的，所以只有很少的神經元被激活，工作，大部分都處於被抑制狀態。

然后，這個稀疏約束項是什么呢？首先對於隱藏層的一個神經元，我們使用 $\textstyle a^{(2)}_j(x)$ 來表示在給定輸入為 $\textstyle x$ 情況下，自編碼神經網絡隱藏神經元 $\textstyle j$ 的激活度，隱藏神經元 $\textstyle j$ 的平均活躍度就是：

，m是樣本的個數。

現在我們想把這個隱藏神經元 $\textstyle j$ 的平均激活度限制在比如=0.2，那么我們怎么讓等於或接近呢？我們使用這個公式：

作為這個隱藏神經元 $\textstyle j$ 的懲罰因子。為什么選擇這個作為懲罰項呢？我們看這個公式的圖像為：

是自變量，=0.2，所以從圖像可以當到為0.2附近的時候，這個懲罰項最小，而與0.2離的越遠，懲罰項值越大。所以如果把這個當這個隱藏神經元的懲罰項加到cost function中，那么最小化這個cost function，就可以使等於或接近。當然我們不能只是限制隱藏層的這一個神經元，而是所有的神經元都必須約束，即：

為了書寫方面，我們令，

加入稀疏懲罰項以后，我們的cost function 就變成了：

$\textstyle \beta$ 控制稀疏性懲罰因子的權重。

現在我們該如何對新的cost function 求各參數的偏導呢？

我們之前在沒有加稀疏項時，我們是一個樣本一個樣本進行前向算法求各單元激活值，進行后項算法求各單元殘差，然后根據公式（2），求每個樣本的小cost function 對各參數的偏導，然后利用公式(1)求出大的 cost function對各參數的偏導（是在這一步對L2正則項求導的）。

加入稀疏項后，我們不改變之前的步驟和程序，我們只需要將前面在后向傳播算法中計算的隱藏層的殘差

$\begin{align} \delta^{(2)}_i = \left( \sum_{j=1}^{s_{2}} W^{(2)}_{ji} \delta^{(3)}_j \right) f'(z^{(2)}_i), \end{align}$

換成

$\begin{align} \delta^{(2)}_i = \left( \left( \sum_{j=1}^{s_{2}} W^{(2)}_{ji} \delta^{(3)}_j \right) + \beta \left( - \frac{\rho}{\hat\rho_i} + \frac{1-\rho}{1-\hat\rho_i} \right) \right) f'(z^{(2)}_i) . \end{align}$

就可以了。

有一個需要注意的地方就是我們需要知道 $\textstyle \hat\rho_i$ 來計算這一項更新。所以在計算每個樣本的后向傳播之前，你需要對所有的訓練樣本計算一遍前向傳播（注意是所有訓練樣本都計算一遍前向傳播，因為要獲取平均激活度，而不像之前沒加稀疏項的時候，只需要一個樣本一個樣本的計算前向后向即可。），從而獲取平均激活度。如果你的訓練樣本可以小到可以整個存到內存之中（對於編程作業來說，通常如此），你可以方便地在你所有的樣本上計算前向傳播並將得到的激活度存入內存並且計算平均激活度。然后你就可以使用事先計算好的激活度來對所有的訓練樣本進行后向傳播的計算。如果你的數據量太大，無法全部存入內存，你就可以掃過你的訓練樣本並計算一次前向傳播，然后將獲得的結果累積起來並計算平均激活度 $\textstyle \hat\rho_i$ （當某一個前向傳播的結果中的激活度 $\textstyle a^{(2)}_i$ 被用於計算平均激活度 $\textstyle \hat\rho_i$ 之后就可以將此結果刪除）。然后當你完成平均激活度 $\textstyle \hat\rho_i$ 的計算之后，你需要重新對每一個訓練樣本做一次前向傳播從而可以對其進行后向傳播的計算。對於后一種情況，你對每一個訓練樣本需要計算兩次前向傳播，所以在計算上的效率會稍低一些。

下面我們來實現一個函數，函數功能：給定神經網絡的結構，網絡參數值，L2正則項、稀疏約束項的權重，稀疏性參數。這個函數就會給出這個神經網絡的cost function的值，和cost function對各個參數的偏導。代碼如下：

   1: function [cost,grad] = sparseAutoencoderCost(theta, visibleSize, hiddenSize, ...

   2:                                              lambda, sparsityParam, beta, data)

   3: %這個sparseAutoencoderCost函數的作用是：你給定神經網絡的結構，網絡參數值，L2正則項、稀疏約束項的權重，稀疏性參數。

   4: %那么這個函數就會給出這個神經網絡的cost function的值，和cost function對各個參數的偏導。注意這個函數不優化各參數的

   5: %值，那個是優化minfun的活，這里是你給定這些參數，然后給出cost function 和偏導。如果你用優化算法比如梯度下降法優化參數

   6: %時，每一次迭代，你都要通過給出當下的參數，然后通過這個函數得到偏導，然后利用優化算法得到新的參數值，再帶入這個函數

   7: %得到偏導，不停迭代，直到收斂。

   8: %注意visibleSize 應該是data的行數，既然data矩陣有了，就應該不需要再輸入visibleSize了，確實

   9: %是這樣，但是這里因為需要使用 visibleSize這個量，所以就又給出了一遍。你輸入的visibleSize要保證等於

  10: %size（data，1）。

  11: % visibleSize: the number of input units (probably 64)

  12: % 輸入層神經元個數，不包括始終為1的偏置節點

  13: % hiddenSize: the number of hidden units (probably 25) 這里神經元網絡隱藏層是1層，

  14: %hiddenSize是這層神經元的個數,同樣不包括始終為1的偏置節點

  15: % lambda: weight decay parameter L2正則項權重系數

  16: % sparsityParam: 就是稀疏性參數，一個標量，一般取為0.2等等

  17: % beta: weight of sparsity penalty term 稀疏項權重系數

  18: % data: Our 64x10000 matrix containing the training data.  So, data(:,i) is

  19: % the i-th training example. 樣本數據

  20: % The input theta is a vector (because minFunc expects the parameters to be a vector).

  21: % We first convert theta to the (W1, W2, b1, b2) matrix/vector format, so that this

  22: % follows the notation convention of the lecture notes. W1指第一層的參數W2第二層參數

  23: %輸入的theta是個向量，向量元素的個數是神經元網絡所有參數的個數。在這里是按這個順序給參數的，theta的前

  25: %hiddenSize*visibleSize是神經網絡第一層和W有關的參數，如果W下標的編號和

  26: %http://www.cnblogs.com/happylion/p/4193527.html第一個圖中的w下標的編號一樣，再結合下面reshape“按列”

  27: %進行重組成矩陣的該函數特性，這里把theta向量的前hiddenSize個元素給了W1矩陣的第一列，而W1矩陣的第一列各元素的含義是

  28: %輸入層第一個神經元分別對應第二層隱藏層的各單元（不包括第二層偏置節點）的參數。W1矩陣的列數就是輸入層神經元的個數

  31: W1 = reshape(theta(1:hiddenSize*visibleSize), hiddenSize, visibleSize);%W1是一個hiddenSize*visibleSize的矩陣

  32: W2 = reshape(theta(hiddenSize*visibleSize+1:2*hiddenSize*visibleSize), visibleSize, hiddenSize);
%因為這個神經網絡是稀疏編碼，所以輸出層的神經元個數和輸入相同，所以W1和為W2參數個數相同

  33: b1 = theta(2*hiddenSize*visibleSize+1:2*hiddenSize*visibleSize+hiddenSize);

  34: b2 = theta(2*hiddenSize*visibleSize+hiddenSize+1:end);

  35:  %注意這里不是對神經網絡的參數初始化，這里的參數都是在每一次迭代中特定給出的，每一次迭代不能每次都進行初始化，
%只有第一次進行初始化。這里只是把給定的特定參數分配給網絡

  38: % Cost and gradient variables (your code needs to compute these values).

  39: % Here, we initialize them to zeros.

  40: cost = 0;

  41: W1grad = zeros(size(W1)); %W1grad 也是一個矩陣，跟W1的行列數是一樣的。里面的元素與神經網絡里面的參數對應的情況和W1

  42: %里元素與神經網絡參數對應的情況是一樣的。第一列各元素的含義是cost function 對輸入層第一個神經元分別對應第二層隱藏層的

  43: %各單元（不包括第二層偏置節點）的參數的偏導。

  44: W2grad = zeros(size(W2));

  45: b1grad = zeros(size(b1));

  46: b2grad = zeros(size(b2));

47:

  48: %% 求cost function 對各參數的偏導數

  49: %  Instructions: Compute the cost/optimization objective J_sparse(W,b) for the Sparse Autoencoder,

  50: %                and the corresponding gradients W1grad, W2grad, b1grad, b2grad.

  51: %

  52: % W1grad, W2grad, b1grad and b2grad should be computed using backpropagation.

  53: % Note that W1grad has the same dimensions as W1, b1grad has the same dimensions

  54: % as b1, etc.  Your code should set W1grad to be the partial derivative of J_sparse(W,b) with

  55: % respect to W1.  I.e., W1grad(i,j) should be the partial derivative of J_sparse(W,b)

  56: % with respect to the input parameter W1(i,j).  Thus, W1grad should be equal to the term

  57: % [(1/m) \Delta W^{(1)} + \lambda W^{(1)}] in the last block of pseudo-code in Section 2.2

  58: % of the lecture notes (and similarly for W2grad, b1grad, b2grad).

  59: %

  60: % Stated differently, if we were using batch gradient descent to optimize the parameters,

  61: % the gradient descent update to W1 would be W1 := W1 - alpha * W1grad, and similarly for W2, b1, b2.

  62: %

  63: %1.forward propagation

  64: data_size=size(data);

  65: biasPara_1=repmat(b1,1,data_size(2));%b1是第1層（輸入層）總是1的那個偏置節點

  66: %與第二層各單元（不包括第二層偏置節點）對應的各個參數，b1這個向量元素的個數是第二層單元數-1個第二層偏置節點數

  67: %所以這里的biasPara_1這個矩陣是hiddenSize*樣本個數的矩陣，並且矩陣每一列都是一樣的。為什么要這么做，復制那么多列，

  68: %因為這里是所有樣本一塊處理的，不是一個樣本一個樣本的處理。因為data也是所有樣本的數據，所以，你就能看出來下面的

  69: %W1*data+biasPara_1，所有樣本是一塊處理的：W1*data得到的矩陣第一列第一個元素是第一個樣本各子元素在第二層

  70: %第一個元素上的線性疊加，所以還要加上b1的第1個元素；W1*data得到的矩陣第一列第2個元素是第一個樣本各子元素在第2層

  71: %第1個元素上的線性疊加，所以還要加上b1的第2個元素；以此類推得到完整第一列，就是處理完第一個樣本，然后處理第二列，

  72: %就是處理第二個樣本，以此類推..所以W1*data+biasPara_1這個矩陣第m列第n行的元素就是第m個樣本，對第二層第n單元

  73: %（不包括偏置節點，因為偏置節點不需要什么sigmoid函數，也不需要有輸入，它始終是1，可以認為激活值始終是1）

  74: %的sigmoid函數的輸入Z。

  75: biasPara_2=repmat(b2,1,data_size(2));

  76: active_value2=sigmoid(W1*data+biasPara_1);

  77: active_value3=sigmoid(W2*active_value2+biasPara_2);

  78: %W2*active_value2+biasPara_2這個矩陣第m列第n行的元素就是第m個樣本，對第3層第n單元的sigmoid函數的輸入Z。

  79: %active_value2這個矩陣每一列是每一個樣本前向算法中第二層各單元的激活值（不包括偏置節點，因為偏置節點

  80: %相當於激活值始終為1）。

  81: %active_value3是第三層各單元的激活值。

82:

  83: %2.computing cost function

  84: ave_square=sum(sum((active_value3-data).^2)./2)/data_size(2);

  85: %這是cost

  86: %function的均方誤差項，因為這里data的label還相當於是data（自編碼），首先說一下一個矩陣“.^2”是矩陣里每個元素都平方。

  87: %sum((active_value3-data).^2）得到的是一個行向量，每一個元素是每一個樣本通過神經網絡得到的輸出層單元激活值（是一個向量）

  88: %與每一個樣本的label（也就是它自身）（也是一個向量）的誤差平方（計算是通過兩個向量對應的各個元素分別對應做差然后平方和）

  89: %然后sum((active_value3-data).^2）得到的這個向量各元素都乘以1/2（也就是每個樣本的誤差平方乘以一個系數1/2），然后再sum

  90: %一下，就是把所有樣本的誤差平方累加起來，最后再除以樣本的個數data_size(2)，就是cost function的均方誤差項。

  91: weight_decay=lambda/2*(sum(sum(W1.^2))+sum(sum(W2.^2)));

  92: %這是const function 中的L2正則項

  93: active_average=sum(active_value2,2)./data_size(2);%sum(active_value2,2)這里是active_value2按行相加，也就是對隱藏層每一個

  94: %神經元，將所有樣本關於這個神經元的激活值累加起來。sum（，2）是按行相加，如果只是sum（），就是默認的按列相加。matlab中

  95: %很多函數默認的都是按列進行某種操作。因此這里得到的是一個列向量。然后點除data_size(2)就是向量每個元素都除以樣本個數

  96: %這樣得到的列向量active_average的每個元素就是隱藏層每一個神經元對於所有樣本的平均激活度。

  97: p_para=repmat(sparsityParam,hiddenSize,1);

  98: %因為下面那行代碼中有p_para./active_average所以要弄成和active_average一樣size的列向量，才能與active_average里的元素點除。

  99: sparsity=beta.*sum(p_para.*log(p_para./active_average)+(1-p_para).*log((1-p_para)./(1-active_average)));

 100: %這是稀疏項

 101: cost=ave_square+weight_decay+sparsity;%cost function就完成了。

 102:

 103: %求殘差

 104: delta3=(active_value3-data).*(active_value3).*(1-active_value3);

 105: %首先說博客中求殘差的公式是按一個樣本的一個神經元來定位的，而這里用matlab求的時候是所有的樣本（一列是一個樣本），

 106: %一層所有的神經元一塊運算的，這是matlab 矩陣運算的特點，避免了for循環。所以這里跟公式中略有不同。雖然是一起運算的，

 107: %但每一個神經元殘差的計算還是要根據公式來，所以這里用到了點乘：矩陣對應元素之間運算。

 108: %公式在http://www.cnblogs.com/happylion/p/4193527.html可以看到。

 109: %這里得到的delta3是一個矩陣，每一列是一個樣本輸入神經網絡后，第三層各個神經元的殘差，列數代表樣本數。

 110: active_average_repmat=repmat(sum(active_value2,2)./data_size(2),1,data_size(2));

 111: %sum(active_value2,2)./data_size(2)就是active_average，然后repmat，因為active_average是一個列向量，再平鋪

 112: %data_size(2)列，那么active_average_repmat就是一個矩陣，弄成矩陣的的形式是為了后面所有樣本第二層所有神經元

 113: %一起運算。

 114: default_sparsity=repmat(sparsityParam,hiddenSize,data_size(2));

 115: %default_sparsity也是一個hiddenSize*data_size(2)矩陣，弄成這個形式也是為了下面所有樣本第二層所有神經元一起的矩陣

 116: %運算。

 117: sparsity_penalty=beta.*(-(default_sparsity./active_average_repmat)+((1-default_sparsity)./(1-active_average_repmat)));

 118: delta2=(W2'*delta3+sparsity_penalty).*((active_value2).*(1-active_value2));

 119: %求第二層的殘差的時候跟輸出層略有不同，因為第二層是隱藏層，又因為稀疏項中是關於隱藏層的內容，所以這里把稀疏項的內容

 120: %融入到殘差里面。這樣由殘差計算偏導的公式就不用改變了。同樣這里是所有樣本第二層所有神經元的一起計算的矩陣計算，

 121: %所以也用到了點乘，矩陣之間對應元素的運算就是博客中的公式。((active_value2).*(1-active_value2))這一項是f'(z(i))

 122: %W2'*delta3得到的是一個矩陣，矩陣的每一列是一個樣本輸入神經網絡后，第二層各個神經元的殘差公式的第一項部分（可以看

 123: %http://www.cnblogs.com/happylion/p/4193527.html中第二個圖就會很明白了），列數代表樣本數。

 124:

 125: %求各參數偏導

 126: W2grad=delta3*active_value2'./data_size(2)+lambda.*W2;

 127: %delta3*active_value2'得到的是一個矩陣，矩陣的第一個元素是所有樣本對應的 小cost function

 128: %對第二層的W11的偏導 的累加形成的總cost function對第二層的W11的偏導。矩陣第一行第二個元素是

 129: %總cost function對第二層的W21的偏導。所以你看這里的W下標的順序和前面的W矩陣是一樣的。然后這個矩陣的元素都除以樣本

 130: %個數，得到的就是cost function中只有均方誤差項時，對第二層各個參數的偏導。然后這個矩陣每個元素加上lambda倍的W2這個

 131: %矩陣對應的元素，得到的W2grad矩陣的元素是總cost function對第二層各參數的偏導。

 132: W1grad=delta2*data'./data_size(2)+lambda.*W1;

 133: b2grad=sum(delta3,2)./data_size(2);%b2grad向量是總cost function對第二層偏置節點各參數的偏導

 134: b1grad=sum(delta2,2)./data_size(2);

 135:

 136: %-------------------------------------------------------------------

 137: % After computing the cost and gradient, we will convert the gradients back

 138: % to a vector format (suitable for minFunc).  Specifically, we will unroll

 139: % your gradient matrices into a vector.

 140:

 141: grad = [W1grad(:) ; W2grad(:) ; b1grad(:) ; b2grad(:)];%正好矩陣（:）這種形式變成一個向量也是按列的，正好和

 142: %reshape（）按列的方式一樣。所以元素與網絡參數對應的順序沒有變化。

 143:

 144: end

 145:

 146: %-------------------------------------------------------------------

 147: % Here's an implementation of the sigmoid function, which you may find useful

 148: % in your computation of the costs and the gradients.  This inputs a (row or

 149: % column) vector (say (z1, z2, z3)) and returns (f(z1), f(z2), f(z3)).

 150: %若果x是個矩陣，也是可以的，igmoid(x)也是一個矩陣，里面的每個元素就是對應x矩陣對應的元素進行sigmoid計算。

 152: function sigm = sigmoid(x)

 154:     sigm = 1 ./ (1 + exp(-x));

 155: end

當我們得到這個神經網絡的cost function 對各參數的偏導數時候，我們想知道這個求得的偏導數是不是正確，我們可以用求導的定義公式來求得，具體內容可以參考這里。具體的代碼如下：

   1: function numgrad = computeNumericalGradient(J, theta)

   2: %% computeNumericalGradient函數的功能是：給定一個方程J(theta)，給定的方式是@你想使用的方程，

   3: %然后你再給一個theta的值，那么這個函數就可以利用求導定義的那個公式來求出這個J方程在theta這點的導數，如果

   4: %theta是一個向量，那么就是求出J方程在theta這點所有未知量的偏導數，因此得出的結果是一個向量，以 numgrad

   5: %給出。

   6: numgrad = zeros(size(theta));%初始化

   7: EPSILON=0.0001;

   8: for i=1:size(theta)  %每一個i是J方程對theta的第i個變量求偏導數。

   9:     theta_plus=theta;

  10:     theta_minu=theta;

  11:     theta_plus(i)=theta_plus(i)+EPSILON;%用定義公式求偏導，那么theta這點的很小的增量或減量中只有第i個

  12:     %變量增加或減少了EPSILON

  13:     theta_minu(i)=theta_minu(i)-EPSILON;

  14:     numgrad(i)=(J(theta_plus)-J(theta_minu))/(2*EPSILON);

  15: end

  16: end

現在我們就檢驗一下：

   1: %% 當我們檢查sparseAutoencoderCost函數寫的對不對或這個函數得到的偏導數是不是跟導數的定義一致，

   2: %我們隨機給定一個theta 即可，如果這個隨機給出的theta使上述兩個函數的偏導數求得的一致，那么我們就認為任何theta

   3: %都會使之相應的一致。並且我們不必把所有樣本都拿來檢驗，只需要把其中的10個做為樣本，把隱藏層單元弄成2，3個即可。

   4: %但輸入層單元個數不能改變，因為它等於樣本矩陣的行數。最后我們用norm(numgrad-grad)/norm(numgrad+grad)

   5: %進行判定numgrad、grad兩個向量的相似程度，如果結果 <1e-9,那么就說他們相似到我們希望的程度了，即檢查結果滿意！

   6: visibleSize = 8*8;

   7: hiddenSize = 3;

   8: sparsityParam = 0.01;

   9: lambda = 0.0001;

  10: beta = 3;

  11: patches = sampleIMAGES;

  12: theta = initializeParameters(hiddenSize, visibleSize);

  13: [cost, grad] = sparseAutoencoderCost(theta, visibleSize, hiddenSize, lambda, ...

  14:                                      sparsityParam, beta, patches(:,1:10));

  15: numgrad = computeNumericalGradient( @(x) sparseAutoencoderCost(x, visibleSize, ...

  16:                                                   hiddenSize, lambda, ...

  17:                                                   sparsityParam, beta, ...

  18:                                                   patches(:,1:10)), theta);

  19: disp([numgrad grad]);

  20: diff = norm(numgrad-grad)/norm(numgrad+grad);

  21: disp(diff);

檢驗結果為：3.8100e-11

所以檢驗滿意！

下面我們要做的就是尋找最優的神經網絡（隱藏層是25，不是測驗時候的3個了）參數。我們要先初始化一個任意的參數向量。我們可以選定好一個優化算法，比如梯度下降法：

$K3T2{(1}X62{%HL9__%2G45$

那么我們就將這個任意的參數向量代入sparseAutoencoderCost函數，得到這個參數向量下的cost function的值和cost function 對這個參數向量各參數的偏導數。然后利用上面的這個梯度下降公式得到新的參數向量。然后再帶入sparseAutoencoderCost函數等得到新的cost function的值和cost function 對這個參數向量各參數的偏導數。不停的迭代cost function 值不斷下降，直到幾乎不下降為止，那么這時候得到的參數向量就是我們要求的神經網絡的參數。

我們這里不用梯度下降法，而是用一種比較流行的擬牛頓法 L-BFGS進行優化。而上面說的梯度下降法工作模式類似，這里也是依靠sparseAutoencoderCost得到的cost function的值和cost function 對這個參數向量各參數的偏導數，並且這個L-BFGS算法代碼（minfun文件代碼里面的一種優化算法，下載地址： minfun下載，並把它放到上面那些代的.m文件的同一文件夾）還會將cost function和偏導數加工得到更高階的偏導數，這樣收斂速度更快。

代碼如下：

   1: visibleSize = 8*8;

   2: hiddenSize = 25;

   3: sparsityParam = 0.01;

   4: lambda = 0.0001;

   5: beta = 3;

   6: patches = sampleIMAGES;

   7: theta = initializeParameters(hiddenSize, visibleSize);

   8: %  Use minFunc to minimize the function

   9: addpath minFunc/

  10: options.Method = 'lbfgs'; % Here, we use L-BFGS to optimize our cost

  11:                           % function. Generally, for minFunc to work, you

  12:                           % need a function pointer with two outputs: the

  13:                           % function value and the gradient. In our problem,

  14:                           % sparseAutoencoderCost.m satisfies this.

  15: options.maxIter = 400;      % Maximum number of iterations of L-BFGS to run

  16: options.display = 'on';

  17: [opttheta, cost] = minFunc( @(p) sparseAutoencoderCost(p, ...

  18:                                    visibleSize, hiddenSize, ...

  19:                                    lambda, sparsityParam, ...

  20:                                    beta, patches), ...

  21:                               theta, options);

  22: %因為每次迭代在sparseAutoencoderCost函數只需要更改參數網絡，網絡結構不變，所以用到了

  23: %@(p) sparseAutoencoderCost(p, visibleSize, hiddenSize, lambda, sparsityParam, beta, patches), theta, options)

運行結果cost function 在0.1248 附近得到收斂，部分如圖：

最后程序中的opttheta就是我們訓練好的神經網絡參數。

我們既然訓練好了，我們就想好好我們的訓練學習好了什么東西、怎么看呢？

我們之前說過，第二層神經元與第一層的關系相當於Pca一樣，第二層神經元值的向量以原始的樣本向量相當於降維以后的結果。只是pca是線性的，而這里因為有logistic函數，所以是非線性的降維。那么隱藏層那一層的每一個神經元的值代表的含義是什么呢？我們仍然可以類比pca，如圖所示：

我們看到pca是將原來的數據在新坐標系下表示，並且把某些坐標軸去掉（把區分數據不明顯的坐標軸去掉），而使數據降維的過程。所以降維后樣本每一維的數據可以看成：在新的基下或新的特征下，樣本在這個特征表達多少或程度多少。

所以類別一下：隱藏層那一層的每一個神經元的值代表的就是在新選擇的某一個特征下（因為隱藏層有25個神經元，所以新選了25個特征，然后將樣本用這25個特征表示出來（應該不是線性疊加），新的數據就有25維，每一維的數值就是樣本在這個特征的表達式多少。），原樣本數據在這個樣本表達了多少。那么我們應該如何找到這個新特征是什么呢？肯定不是的值，因為它是樣本在某個特征下的表達程度，因為樣本可以用這25個特非線性性疊加起來。所以這些特征應該和原樣本pitch大小相同，類似於特征量一樣。那么應該怎么求呢？做法是對下面這個方程求最大值。

因為參數是已知的，如果我們求隱藏層某一個神經元的激活值求最大值，那么就可以找到使這個神經元激活的輸入是多少，那么這個輸入就是其中的一個特征。一般來說使其中的一個神經元激活，另外其他的就不一定激活了，因為我們已經加入了激活值的稀疏項。還有一個問題需要注意：如果你的樣本xi的值都很大的話，得到平凡解，所以必須給 $\textstyle x$ 加約束，若假設輸入有范數約束 $\textstyle ||x||^2 = \sum_{i=1}^{100} x_i^2 \leq 1$ ，則令隱藏單元 $\textstyle i$ 得到最大激勵的輸入應由下面公式計算的像素 $\textstyle x_j$ 給出（共需計算100個像素，j=1,…,100）：

$\begin{align} x_j = \frac{W^{(1)}_{ij}}{\sqrt{\sum_{j=1}^{100} (W^{(1)}_{ij})^2}}. \end{align}$

具體的推導過程見下一個博客

當我們用上式算出各像素的值、把它們組成一幅圖像、並將圖像呈現在我們面前之時，那么這個特征我們就可以看到了。並且我們看到的這個圖像是帶有有界范數的圖像。注意我們這是為了向可視化特征才給輸入x加了范數約束，因為不加約束會得到平凡解，（平凡解就是說根本不用解誰都知道某個量是這個方程的解，比如AX=0，誰都知道X=0是這個方程的解，所以X=0就是平凡解）但是一般的情況你的輸入樣本是不需要約束 $\textstyle ||x||^2 = \sum_{i=1}^{100} x_i^2 \leq 1$ 這種的約束。不過樣本應該也要進行歸一化，把數據歸一到[0,1]之間。

可視化特征的程序為：

W1 = reshape(opttheta(1:hiddenSize*visibleSize), hiddenSize, visibleSize);
display_network(W1', 12); %display_network函數中已經包括對輸入加了范數約束。

那么25個特征如圖所示：

我們可以看出來這25個特征都是一些“邊緣”特征，為什么會是邊緣特征呢？

因為人的神經系統也是分層的，第一層的神經層也是對一幅圖像的邊緣特征處於激活。詳細並且很經典的講解特征的是這一篇博客，可以參考這里。

參考：http://deeplearning.stanford.edu/wiki/index.php

http://www.cnblogs.com/hrlnw/archive/2013/06/08/3127162.html

http://deeplearning.stanford.edu/wiki/index.php/Exercise:Sparse_Autoencoder

免責聲明！

本站轉載的文章為個人學習借鑒使用，本站對版權不負任何法律責任。如果侵犯了您的隱私權益，請聯系本站郵箱yoyou2525@163.com刪除。

猜您在找 Deep Learning學習隨記（一）稀疏自編碼器 Deep Learning 學習筆記（8）：自編碼器( Autoencoders ) 超詳細的EM算法理解【deep learning精華部分】稀疏自編碼提取高階特征、多層微調完全解釋及代碼逐行詳解超詳細的RNN代碼實現(tensorflow) Deep Learning 閱讀筆記：Convolutional Auto-Encoders 卷積神經網絡的自編碼表達【優化算法】Greedy Randomized Adaptive Search算法超詳細解析，附代碼實現TSP問題求解 TextRank算法詳細講解與代碼實現（完整）基本算法思想Java實現的詳細代碼 Keras代碼超詳細講解LSTM實現細節