摘要

本文利用人工神經網路實現對垃圾郵件的識別和分類。神經網路型別選擇常用的BP神經網路，程式設計環境為Matlab。為了加深對神經網路的理解，不使用已有的神經網路工具包，而是根據底層原理編寫實現完整的正向輸出過程和引數調整過程。此外，根據實驗中獲取的經驗，引入了一種學習速度動態調整的方法，可以加快前期學習速度，在相同學習次數下獲得更小的累計誤差。作為效能評價的一部分，將BP神經網路與另一種可行分類方法——貝葉斯分類進行橫向對比，顯示了神經網路的優越效能

1.BP神經網路

1.1原理

如圖1所示，BP神經網路由輸入層，隱含層，輸出層組成。輸入層節點數量為n，序號用i表示，即每個節點表示為\({x_i}{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} i = 1,2 \cdots n\)；隱含層的節點數量為p，序號用j表示，隱含層各節點的輸入值表示為\(h{I_j}{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} j = 1,2 \cdots p\)，輸出值為；輸出層節點數為q，序號用k表示，所以輸出層各節點的輸入值為\(y{I_k}{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} k = 1,2 \cdots q\)，輸出值為\(y{O_k}{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} k = 1,2 \cdots q\). \({w_{ij}}\)是從輸入節點xi到隱含節點hj的權值，\({w_{jk}}\)是從隱含層節點到輸出層節點的權值。\(b{H_j}\)是隱含層各節點的啟用閾值，\(b{Y_k}\)是輸出層各節點的啟用閾值。

Fig. 1 BP神經網路結構

這樣，各節點之間的關係就可以如下表示：

隱含層輸入：\(h{I_j} = \sum\limits_{i = 1}^n {{x_i}{w_{ij}} - } b{H_j}\)

隱含層輸出：\(h{O_j} = f(h{I_j})\)

輸出層輸入：\(y{I_k} = \sum\limits_{j = 1}^p {{w_{jk}}h{I_j} - b{Y_k}} \)

輸出層輸出：\(y{O_k} = f\left( {y{I_k}} \right)\)

其中，\(f\left(  \cdot  \right)\)是啟用函式，選用單極性s型函式可以簡化計算，因為其倒數可以用其自身表示。單極性的s型啟用函式表示為：\(f\left( x \right) = \frac{1}{{1 + {e^{ - bx}}}}\)，b為大於零的係數。

對於一個訓練樣本\({{\bf{x}}_r} = ({x_1},{x_2} \ldots {x_n})\)，及其對應的期望輸出\({{\bf{d}}_r} = ({d_1},{d_2} \ldots {d_q})\)，神經網路輸出的均方誤差定義為\({E_r} = \frac{1}{2}\sum\limits_{k = 1}^q {{{\left( {{d_k} - y{O_k}} \right)}^2}} \)：，這裡乘以常係數1/2是為了推導計算方便。

神經網路訓練的過程就是為了使單次訓練誤差\({E_r}\)最小化，而最小化的方法是根據誤差調整係數w和閾值b. 一種方法是以誤差的梯度下降方向調整係數。\({E_r}\)分別對wij和wjk求導，得：$$\frac{{\partial {E_r}}}{{\partial {w_{jk}}}} = \frac{{\partial {E_r}}}{{\partial yI}}\frac{{\partial yI}}{{\partial {w_{jk}}}} = {\delta _r}{\kern 1pt}  \cdot h{O_j}$$ ，其中 $${\delta _k}{\kern 1pt}  = ({d_k} - y{O_k})f\left( {y{I_k}} \right)\left( {1 - f\left( {y{I_k}} \right)} \right)$$ 所以\(Wjk\)的更新為\({w_{jk}}^{N + 1} = {w_{jk}}^N + \eta  \cdot {\delta _k} \cdot h{O_j}\)， 閾值的更新：$$b{Y_k}^{N + 1} = b{Y_k}^N - \eta  \cdot {\delta _k}$$

類似地，可以得到隱含層權值的更新為：

$${w_{ij}}^{N + 1} = {w_{ij}}^N + \eta  \cdot {g_j} \cdot {x_i}$$

$$b{H_j}^{N + 1} = b{H_j}^N - \eta  \cdot {g_j}$$ ，其中$${g_j} = f\left( {h{I_j}} \right)\left( {1 - f\left( {h{I_j}} \right)} \right)\sum\limits_{k = 1}^q {{\delta _k}{w_{jk}}} $$

1.2演算法描述

上一小節中詳細說明了BP神經網路的結果，並從降低均方誤差出發，推匯出了梯度下降法的權值調整方法。利用BP神經網路解決問題首先利用已知樣本訓練神經網路，權係數收斂到一定程度之後，就可以輸入未知樣本，正向計算輸出了。由此可見，關鍵在於如何訓練網路的引數。

輸入訓練集的的累計誤差可以作為是否結束訓練的一個依據，定義累計誤差：$$E = \frac{1}{m}\sum\limits_{r = 1}^m {{E_k} = } \frac{1}{m}\sum\limits_{r = 1}^m {\left( {\frac{1}{2}\sum\limits_{k = 1}^q {{{\left( {{d_k} - y{O_k}} \right)}^2}} } \right)} $$ 當累計誤差下降到小於預設的閾值時就可以認為訓練已經達到了要求。但是當閾值設定不當時，也許永遠也降不到預設的累計誤差，這時可以設定最大訓練次數M，訓練次數達到這個上界時即結束訓練。

訓練BP神經網路的過程可以用以下流程圖表示：

Fig. 2 訓練流程圖

2.實現細節與優化

2.1 資料樣本

2.1.1 資料說明

本實驗利用BP神經網路實現對垃圾郵件的識別。實驗所用的資料樣本spambase.data是一個垃圾郵件的資料庫，來自於惠普公司的Hewlett Packard Labs實驗室。該資料庫中包含了4601個樣本，其中1813條為垃圾郵件(spam)，佔比約40%；每封郵件根據文字內容提取出57個屬性，並且具有是否為垃圾郵件的標記。在資料檔案中，每一行為一個樣本，58個屬性按順序排列，使用”,”分隔。各屬性的簡單說明表所述。

表 1樣本屬性說明

2.1.2 資料初步處理

歸一化：

資料中第1-57維是訓練輸入，第58維是期望的輸出。由於輸入樣本的範圍各不相同，首先應當歸一化到相同的範圍，這樣權值W才能最終收斂在相對集中的範圍內。又由於這些資料大多數為較小值，例如1-48維的取值範圍是[0,100]但大部分都是接近於0的數值，整體上表現為泊松分佈。為了讓歸一化後的屬性值大部分分佈在0值附近，應當取上下界不等的歸一化方法，上界的絕對值大於下界絕對值。用數學式表示為：x ＝ (x - min)/(max - min). 根據資料特點，具體取max=5，min=-1.

期望輸出轉換：

樣本類別為2類，第58維是類別標記，原始資料中用0表示正常郵件，1表示垃圾郵件。在神經網路中，c類可以用c個輸出節點表示，因此對期望輸出轉化為：正常郵件：{1,0}, 垃圾郵件：{0,1}.

亂序：

原始樣本中按照正常郵件，垃圾郵件排序，無論從訓練的效能角度，還是測評的科學性出發，都應當先將資料打亂。具體的打亂方法可以這樣實現：1.利用rand函式產生4601個隨機資料。2.對這些資料進行排序，輸出每個資料對應的順序，就得到了亂序的索引。3.利用這個索引從原始資料中取出樣本，重新排列。

2.2 動態調整的學習速率優化

在引數調整方法(梯度下降)確定之後，學習速率\(\eta \)將會是影響訓練速度和效能的重要因素。速率過大，在前期雖然具有較快的累計誤差下降速度，但後期會產生震盪，難以收斂；速率過小則訓練速度慢，所需的訓練次數加大。顯然，學習速率\(\eta \)應當隨著訓練次數的增長而逐漸下降，進一步說，應當是隨著累計誤差的下降而調低。一種明顯的解決方法是將速率\(\eta \)設定為累計誤差的函式，在誤差較大時速率較大，而誤差降低之後\(\eta \)也同步縮小。設定最小控制速率\({\eta _{\min }}\)和最大控制速率\({\eta _{\max }}\)，累計誤差\({E _{accu }}\)，累計誤差閾值\({E_{thr}}]\)，則學習速率可以這樣設定：

$$\eta  = {\eta _{\min }} + {\eta _{\max }}{\left( {\frac{{{E_{accu}} - {E_{thr}}}}{{{E_{accu}} + {E_{thr}}}}} \right)^2}$$ .

2.3 引數設定

前述所需設定的引數設定如下，累計誤差閾值取0.02；最大訓練輪數100；由於輸入樣本維數是57，輸出結果為2維，因此輸入節點數57，輸出節點2. 隱含層的確定方法根據經驗公式hiddenNode=sqrt(inputNode+outputNode) + a, 0<a<10. 因此，隱含層節點數量取15.

學習速率動態調整方面，最低值控制在etaMin=0.01, 最大值控制在etaMax=0.4, 初始時速率為0.05, 在第二輪時eta就會根據誤差大小進行調整。所有的引數設定如下：

errorThr=0. 02; %最大容許誤差

M=100; %最大訓練輪數

inputNode=57;

outputNode=2;

hideNode=15; %sqrt(inoutNode+outputNode)+a , 0<a<10

eta=0.05;%本次訓練速度

etaMin=0.01;

etaMax=0.4;

3. 效能測試

3.1 分類正確性

樣本資料庫中共有4601條資料，選擇其中約2/3用作訓練，1/3用作測試。由於在第3節預處理過程中已經將資料亂序排列了，因此可以直接取前3000條作為訓練資料，剩下的作為測試資料。

由於資料的選擇、權值的初始化都是隨機的，因此訓練結束後，利用剩下的測試資料對所得模型進行正確性檢測，每次結果都是稍有差異的。並且顯然，識別效能是與累計誤差相關的，多次執行次程式，設定不同的訓練次數，考察分類效能，訓練與檢驗結果記錄在表2中。

表 2 訓練與檢驗統計結果

3.2 累計誤差下降速度

訓練過程實質上是累計誤差下降的過程，因此每一輪訓練結束後的累計誤差E_accu可以展現訓練的動態過程。本文中使用了一種學習速度自動調整的動態設定方法，可以在初始時使用較大的學習速度，隨後慢慢降低。為了證實這種方法的優點，訓練100輪，分佈採用固定學習速率和採用這種動態調整方法，比較在累計誤差下降曲線上的差異。

圖3是學習過程中累計誤差的下降曲線與訓練輪數M的關係圖。圖中共同展示了本文中提出的學習速率動態調整方法與固定速率方法所得的結果。使用較大的學習速率，在初始時能更快降低累計誤差，但隨後會出現震盪；較低的學習速率雖然能穩定降低誤差，但初始學習速率過慢。動態速率方法(x符號表示的曲線)起初快速下降，隨著誤差的降低，學習速率也減小，因此能穩定地逼近最優點。從圖3中可見，係數動態調整方法在初始時的下降速率與eta=0.4相當，最終穩定下降，在所有固定速率取值方案中得到了最低的累計誤差，且沒有產生震盪。

Fig. 3 本文中使用的學習速率動態調整方法與固定學習速率對比，動態調整方法初始時具有較快的下降速率，而後能平穩降低累計誤差

Fig. 4 學習速率動態調整方法中，eta值逐漸減小

3.3 與貝葉斯分類器橫向對比

垃圾郵件的樣本是從實際郵件中統計得到的，對垃圾郵件的識別是一個二分類問題。貝葉斯分類器是一種適合用於該場景下的分類器，並且在現實中曾經採用過這種方案進行垃圾郵件的識別。

貝葉斯分類是一種有監督的分類方法。其原理如下：首先統計各類的先驗概率\(P({\omega _i})\)，以及類條件概率分佈\(P({\bf{x}}|{\omega _i})\)，然後通過貝葉斯公式

$$P({\omega _i}|{\bf{x}}) = \frac{{P({\bf{x}}|{\omega _i}) \cdot P({\omega _i})}}{{P({\bf{x}})}} = \frac{{p({\bf{x}}|{\omega _i})d{\bf{x}} \cdot P({\omega _i})}}{{p({\bf{x}})d{\bf{x}}}} = \frac{{p({\bf{x}}|{\omega _i}) \cdot P({\omega _i})}}{{p({\bf{x}})}}$$ 將觀測量轉換得到後驗概率\(P({\omega _i}|{\bf{x}})\) . 對於兩類分類問題，根據後驗概率的大小判決所屬分類，即：

$$P({\omega _1}|{\bf{x}}){\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \matrix     >   \\      <   \\ 

 \endmatrix {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} P({\omega _2}|{\bf{x}}){\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt}  \Rightarrow {\bf{x}} \in \matrix    {{\omega _1}}  \\     {{\omega _2}}  \\ 

 \endmatrix $$

上式中由於分母都一樣，可以得到變換形式的決策式：

$$p({\bf{x}}|{\omega _1}){\kern 1pt} {\kern 1pt} P({\omega _1}){\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \matrix     >   \\      <   \\ 

 \endmatrix {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} P({\bf{x}}|{\omega _2}){\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} P({\omega _2}){\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt}  \Rightarrow {\bf{x}} \in \matrix    {{\omega _1}}  \\     {{\omega _2}}  \\ 

 \endmatrix $$

筆者曾經利用這些資料進行過貝葉斯分類器的設計和測試。同樣取樣2/3資料樣本作為訓練集，1/3資料作為測試集，分類效果如下。這裡涉及到將離散屬性近似到連續的概率密度函式的問題，因此需要引用“頻數量化區間的概念”。總之，利用貝葉斯分類器，在不同統計量化區間下，識別率如下表。

表 3 橫向對比，利用貝葉斯分類器分類的效能

使用貝葉斯分類器，雖然正常郵件的識別率也很高，但垃圾郵件的檢出率最大在85.5%左右。對比BP神經網路的效能，正常郵件的識別率與垃圾郵件的檢出率都比貝葉斯分類器高，垃圾郵件的識別率在90%以上，長時間訓練後(300輪)，垃圾郵件檢出率達到93.5%，正常郵件識別率達到92.9%. 可見神經網路方法在充分訓練之後，得到的兩類分類正確率都較高且平衡。

在計算量方面，訓練神經網路需要耗費大量的算力，但一旦訓練完成，對未知樣本正向計算時並不複雜，只是兩次矩陣乘法運算。

1. 總結

本實驗根據BP神經網路原理從底層編碼實現了正向計算和引數反向調整的過程。通過動手操作，我進一步加深了對神經網路的理解。

神經網路訓練的目標是要降低累計誤差，因而一種顯而易見的方法是沿梯度下降方向調整權值。

神經網路具通用性，對不同問題的解決具有極大的普適性。特別是包含隱含層的BP神經網路，可以描述非線性問題。對於模型未知或難以描述的情形，可以直接進行訓練，在忽略其模型的情況下得到解決方案。

1. Reference

1. 附錄-程式程式碼

%2018年4月29日

%使用自適應調整的學習速率
%方法：eta=etaMin+( (errNow-errTh)/(errNow+errTh))^3*etaMax

clear all;

tic;

%%%----引數設定----
b=1; %激發函式f(x)中的係數
f = @(x) 1./(1+exp(-b*x)); %s型激發函式,其導數為b*f(x)*(1-f(x))


errorThr=0.02; %最大容許誤差
M=100; %最大訓練輪數



inputNode=57;
outputNode=2;
hideNode=15; %sqrt(inoutNode+outputNode)+a , 0<a<10



eta=0.05;%本次訓練速度
etaMin=0.01;
etaMax=0.4;
%連線權值和閾值
wij=rand(inputNode,hideNode)/5;
wjk=rand(hideNode,outputNode )/5;
bH=rand(1,hideNode)/4;
bY=rand(1,outputNode)/4;

%%%----引數設定結束



data_all=dlmread('spambase.data', ',');
%讀取原始資料，58維，1-48：特定單詞的頻率[0,100]的實數；49-54：特定字元的頻率[0,100]

%55；大寫字母遊程平均長度； 56:大寫字母最大遊程；57：最大字母總數； 58：類別標記



%打亂資料，劃分訓練集和測試集

data_all_samples=mapminmax(data_all(:,1:57)',-1,10);%樣本歸一化處理
% data_all_samples=data_all(:,1:57)';
data_all_label=data_all(:,58)'; %提取標籤



%打亂原始資料，並分為訓練集和測試集

k=rand(1,size(data_all,1));
[m,n]=sort(k);%隨機打亂原始資料，n作為序號
 data_train_input=data_all_samples(:,n(1:3000));
 data_train_output=data_all_label(:,n(1:3000));
data_test_input=data_all_samples(:, n(3001:4601));
data_test_target=data_all_label(:, n(3001:4601),:);



%對訓練集的標記擴充套件為2維，作為輸出

data_train_target=zeros(2,size(data_train_output,2));

for i=1:size(data_train_output,2)
 switch data_train_output(i)
        case 0
            data_train_target(:,i)=[1;0]; %非垃圾郵件
        case 1
            data_train_target(:,i)=[0;1]; %垃圾郵件
    end
end



errorAccuRound=[];%每輪訓練結束後的累計誤差
etaRecord=[];
for round =1:M
    errorAccu=0;%單輪訓練內的誤差
for index=1:3000;

   

    %樣本index;的輸出
    hI=wij'* data_train_input(:,index ) -bH';
    hO=f(hI);
    yI=(wjk'* hO -bY');
    yO=f(yI);



  %更新系數：

    delta_O= (data_train_target(:,index)-yO).*f(yI).*(1-f(yI));
    %更新Wjk和bY

    wjkold=wjk;
    for j=1:hideNode
        for k=1:outputNode
            wjk(j,k)=wjk(j,k)+eta*delta_O(k)*hO(j);
        end
    end

    bY=(bY' - eta*delta_O)';


    % hI=wij'* data_train_input(:,1) -bH';
    % hO=f(hI);
    % yI=(wjk'* hO -bY');
    % yO=f(yI)



    %更新Wij和bH

    delta_H=sum(repmat(delta_O',hideNode,1).*wjkold,2).*f(hI).*(1-f(hI));
    for i=1:inputNode

        for j=1:hideNode
            wij(i,j)=wij(i,j)+eta*delta_H(j)*data_train_input(i,index);

        end

    end

    bH=(bH' - eta*delta_H)';

   
     %計算誤差：
    errorAccu=errorAccu+0.5*(data_train_target(:,index)-yO)'*(data_train_target(:,index)-yO);



end

   

    %根據累計誤差判斷是否結束訓練
    errorAccuRound=[errorAccuRound, errorAccu/3000];%儲存本輪的累計誤差
    if errorAccu/3000 < errorThr %累計誤差小於要求時結束訓練

        break

    end

   

    %學習速度的動態調整 
    eta=etaMin + etaMax*((errorAccu/3000 - errorThr) / (errorThr+errorAccu/3000) )^2
    etaRecord=[etaRecord eta];

end

%%

count0=0;
count1=0;
for index_test=1:1601;

   

  

    hI=wij'* data_test_input(:,index_test) -bH';
    hO=f(hI);
    yI=(wjk'* hO -bY');
    yO=f(yI);
    if (yO(1)>yO(2))
        class=0;
    else

        class=1;

    end

       if class==0 && data_test_target(index_test)==0
           count0=count0+1;
        elseif class==1 && data_test_target(index_test)==1
               count1=count1+1;
        end

end

fprintf('%s%.3d\n','非垃圾郵件正確識別率',count0/(1601-sum(data_test_target)))
fprintf('%s%.3d\n','垃圾郵件正確識別率',count1/sum(data_test_target))
% data_train_target(:,index_test);

toc
plot(errorAccuRound)