神經網絡學習筆記_1(BP網絡分類雙螺旋線)

這個也是本學期模式識別的一個小作業題。

實驗目的：采用神經網絡的方法對二維雙螺旋樣本點進行分類。

實驗數據：由於本次的實驗數據是雙螺旋，需要用數學公式產生，其產生方法見參考文獻[1].

　　         即由下面的公式產生：

　　

 

　　實驗數據分為2部分：訓練數據和測試數據。由於雙螺旋曲線有自己的方程表達式，我們產生的測試數據不能與訓練數據一樣，否則訓練出來的網絡沒有說服力，因此我們在程序中2個樣本集在雙螺旋曲線上的采樣間隔不同，且起始位置不同，這樣保證了測試樣本和訓練樣本的數據集沒有重疊，訓練出來的模型的分類能力更具有說服力。

　　下圖是試驗過程中產生的訓練樣本和測試樣本分布圖：

　　

 

　　其中訓練樣本200個，每條螺旋線100個；測試樣本數據也為200個，每天螺旋線100個。

試驗方法：本次試驗我采用了2種方法，第一種是利用matlab自帶的ANN工具箱函數，來訓練出一個BP神經網絡；另一種方法是不采用matlab自帶的工具箱，而自己用代碼實現 　　　　  　　       一個BP神經網絡，最后給出了2者的對比。

 

方法一：采用matlab自帶神經網絡工具箱

 

實驗說明：采用的是matlab自帶的BP神經網絡，隱含層只有1層，該層有8個神經元。程序代碼中的studytra_data, studytra_label為訓練數據，分別是2*200的矩陣。　　   　　   　　      predict_test_data, predict_test_label 為測試數據，也是2*200的矩陣。網絡形成用的是feedforwardnet函數，BP網絡訓練采用的是train函數。

實驗結果：

　　實驗預測的結果如下所示：

　　

　　即第一條螺旋線的預測精度為92%,第二條螺旋線的預測精度為68%.

程序代碼(matlab)：

%% 產生雙螺旋數據

studytra_num=100;  %訓練樣本數目
predict_test_number=100;    %測試樣本數目

i=(1:1:studytra_num)';
%雙螺旋數據點的產生方程
alpha1=pi*(i-1)/25;
beta=0.4*((105-i)/104);
x0=0.5+beta.*sin(alpha1);
y0=0.5+beta.*cos(alpha1);
z0=zeros(studytra_num,1);
x1=0.5-beta.*sin(alpha1);
y1=0.5-beta.*cos(alpha1);
z1=ones(studytra_num,1);

%% 事實證明BP神經網絡在訓練數據時與輸入數據正負樣本的順序是有關系的
% 如果一開始的一半數據都是正(負)樣本，后面的全是負(正)樣本，則訓練出來的
%效果不好，所以這里需要隨機打亂
k=rand(1,2*studytra_num);
[m,n]=sort(k);

studytra=[x0 y0 z0;x1,y1,z1]; %1條螺旋線數據點,200*3的矩陣
trian_label1=studytra(n(1:2*studytra_num),end)';   %訓練數據類別，1*200的行向量
studytra_data1=studytra(n(1:2*studytra_num),1:end-1)'; %訓練數據屬性，2*200的矩陣

%把1維的輸出變成2維的輸出,studytra_labe2為200*2的矩陣
for i=1:2*studytra_num
    switch trian_label1(i)
        case 0
            studytra_label2(i,:)=[1 0];
        case 1
            studytra_label2(i,:)=[0 1];
    end
end

[studytra_data,studytra_datas]=mapminmax(studytra_data1);%studytra_data為2*200的矩陣
studytra_label=studytra_label2'; %studytra_label為2*200的矩陣
         
plot(x0,y0,'r+');
hold on;
plot(x1,y1,'go');


%% 產生雙螺旋測試數據

i=(1.5:1:predict_test_number+0.5)';    %每類51個樣本
%雙螺旋數據點的產生方程
alpha2=pi*(i-1)/25;
beta2=0.4*((105-i)/104);
m0=0.5+beta2.*sin(alpha2);
n0=0.5+beta2.*cos(alpha2);
s0=zeros(predict_test_number,1);
m1=0.5-beta2.*sin(alpha2);
n1=0.5-beta2.*cos(alpha2);
s1=ones(predict_test_number,1);

predict_test=[m0 n0 s0;m1,n1,s1]; %1條螺旋線數據點,3*102的矩陣
predict_test_label1=predict_test(:,end)';   %測試數據類別，1*102的行向量
predict_test_data1=predict_test(:,1:end-1)'; %測試數據屬性，2*102的矩陣

%把1維的輸出變成2維的輸出,studytra_labe2為200*2的矩陣
for i=1:2*predict_test_number
    switch predict_test_label1(i)
        case 0
            predict_test_label2(i,:)=[1 0];
        case 1
            predict_test_label2(i,:)=[0 1];
    end
end

predict_test_label=predict_test_label2'; %predict_test_label為2*102的矩陣
         
%%  畫出測試數據雙螺旋曲線
plot(m0,n0,'c+');
hold on;
plot(m1,n1,'yo');
legend('訓練數據螺旋線1','訓練數據螺旋線2','測試數據螺旋線1','測試數據螺旋線2');

predict_test_data=mapminmax('apply',predict_test_data1,studytra_datas);

%% 至此，訓練和測試數據如下所示：
studytra_data  ;    %2*200的矩陣
studytra_label  ;   %2*200的矩陣
predict_test_data   ;    %2*200的矩陣
predict_test_label   ;   %2*200的矩陣

% %%  調用matlab的神經網絡工具箱函數用來訓練,方法1
% net=newff(minmax(studytra_data),[10,2],{'tansig','purelin','studytragdm'})
% inputWeights=net.IW{1,1};
% inputbias=net.b{1};
% layerWeights=net.LW{2,1};
% layerbias=net.b{2};
% 
% net.studytraParam.show=50;
% net.studytraParam.lr=0.05;
% net.studytraParam.mc=0.9;
% net.studytraParam.epochs=1000;
% net.studytraParam.goal=1e-2;
% 
% [net,tr]=studytra(net,studytra_data,studytra_label);


%% 調用matlab的神經網絡工具箱函數用來訓練,方法2
net=feedforwardnet(8);
%注意此時進入train函數的樣本每一列為一個樣本的屬性，列數為樣本數，且命名時最好不要含有train,test字樣，否則matlab會報錯
net=train(net,studytra_data,studytra_label);
%view(net);
predict_label=sim(net,predict_test_data);


% %% 用訓練到的模型預測數據
% for i=1:2*predict_test_number
%     for j=1:midnum
%         I(j)=predict_test_data(:,i)'*w1(j,:)'+b1(j);
%         Iout(j)=1/(1+exp(-I(j)));%Iout為1*3的行向量
%     end
%     predict_test(:,i)=w2'*Iout'+b2;%predict_test為2*102的矩陣
% end


%% 預測結果分析
for i=1:2*predict_test_number
    output_pred(i)=find(predict_label(:,i)==max(predict_label(:,i)));    %out_pred為1*102的矩陣
end

error=output_pred-predict_test_label1-1;    %


%% 計算出每一類預測錯誤的個數總和
k=zeros(1,2); %k=[0 0]
for i=1:2*predict_test_number
    if error(i)~=0    %matlab中不能用if error(i)！=0 
        [b c]=max(predict_test_label(:,i));
        switch c
            case 1
                k(1)=k(1)+1;
            case 2
                k(2)=k(2)+1;
        end
    end
end


%% 求出每一類總體的個數和
kk=zeros(1,2); %k=[0 0]
for i=1:2*predict_test_number
    [b c]=max(predict_test_label(:,i));
    switch c
        case 1
            kk(1)=kk(1)+1;
        case 2
            kk(2)=kk(2)+1;
    end
end


%% 計算每一類的預測正確率
accuracy=(kk-k)./kk

 

方法二：自己實現BP神經網絡的訓練和測試函數

 

實驗說明：實現BP神經網絡的方法參考的是文獻[2]，主要分為下面6個步驟：

1. 網絡的初始化
2. 隱含層輸出計算
3. 輸出層輸出計算
4. 誤差計算
5. 權值更新
6. 閾值更新

　　同樣，訓練和測試樣本都為2*200的矩陣，且訓練過程中設置的循環次數為200次，只有1個隱含層，且神經元個數為5個。

　　實驗結果：

　　實驗預測的精度如下所示：

　　

　　即第1條螺旋線的預測精度為94%,，第2條螺旋線的預測精度為37%。

實驗代碼(matlab)：

%% 產生雙螺旋數據,每類100個樣本點，共200個樣本
train_num=100;
train_circle_number=5000;
test_number=100;
i=(1:1:train_num)';

%雙螺旋數據點的產生方程
alpha1=pi*(i-1)/25;
beta=0.4*((105-i)/104);
x0=0.5+beta.*sin(alpha1);
y0=0.5+beta.*cos(alpha1);
z0=zeros(train_num,1);
x1=0.5-beta.*sin(alpha1);
y1=0.5-beta.*cos(alpha1);
z1=ones(train_num,1);

%% 事實證明BP神經網絡在訓練數據時與輸入數據正負樣本的順序是有關系的
% 如果一開始的一半數據都是正(負)樣本，后面的全是負(正)樣本，則訓練出來的
%效果不好，所以這里需要隨機打亂
k=rand(1,2*train_num);
[m,n]=sort(k);

train=[x0 y0 z0;x1,y1,z1]; %1條螺旋線數據點,200*3的矩陣
trian_label1=train(n(1:2*train_num),end)';   %訓練數據類別，1*200的行向量
train_data1=train(n(1:2*train_num),1:end-1)'; %訓練數據屬性，2*200的矩陣

%把1維的輸出變成2維的輸出,train_labe2為200*2的矩陣
for i=1:2*train_num
    switch trian_label1(i)
        case 0
            train_label2(i,:)=[1 0];
        case 1
            train_label2(i,:)=[0 1];
    end
end

train_label=train_label2'; %train_label為2*200的矩陣
         
plot(x0,y0,'r+');
hold on;
plot(x1,y1,'go');
%legend();

%% BP神經網絡結構的初始化
%網絡結構，2個輸入，3個神經元，2個輸出
innum=2;
midnum=5;
outnum=2;

[train_data,train_datas]=mapminmax(train_data1);

%輸入輸出取值閾值隨機初始化
%w1矩陣表示每一行為一個隱含層神經元的輸入權值
w1=rands(midnum,innum);%rands函數用來初始化神經元的權值和閾值是很合適的,w1為3*2的矩陣
b1=rands(midnum,1);%b1為3*1的矩陣
%w2矩陣表示每一列為一個輸出層神經元的輸入權值
w2=rands(midnum,outnum);%w2為3*2的矩陣
b2=rands(outnum,1);%b2為2*1的矩陣

%用來保存上一次的權值和閾值，因為后面的更新方差是遞歸的，要用到
w1_1=w1;w1_2=w1_1;
b1_1=b1;b1_2=b1_1;
w2_1=w2;w2_2=w2_1;
b2_1=b2;b2_2=b2_1;

%學習率的設定
alpha=0.05;

%訓練10次就ok了，而不管訓練后的結果如何
for train_circle=1:train_circle_number  ;
    for i=1:2*train_num; %200個訓練樣本
       %% 輸入層的輸出
        x=train_data(:,i);%取出第i個樣本，x(i)為2*1的列向量
        %% 隱含層的輸出
        for j=1:midnum;
            I(j)=train_data(:,i)'*w1(j,:)'+b1(j);  %I(j)為1*1的實數
            Iout(j)=1/(1+exp(-I(j)));   %Iout(j)也為1*1的實數
        end     %Iout為1*3的行向量   
        %% 輸出層的輸出
         yn=(Iout*w2)'+b2;   %yn為2*1的列向量，因此此時的傳函為線性的，所以可以一步到位，不必上面
        
        %% 計算誤差
        e=train_label(:,i)-yn; %e為2*1的列向量，保存的是誤差值
        
        %計算權值變換率
        dw2=e*Iout; %dw2為2*3的矩陣，每一行表示輸出接點的輸入權值變化率
        db2=e'; %e為1*2的行向量
        
        for j=1:midnum
            S=1/(1+exp(-I(j)));
            FI(j)=S*(1-S);  %FI(j)為一實數，FI為1*3的行向量
        end
        
        for k=1:1:innum
            for j=1:midnum
                dw1(k,j)=FI(j)*x(k)*(e(1)*w2(j,1)+e(2)*w2(j,2));    %dw1為2*3的矩陣
                db1(j)=FI(j)*(e(1)*w2(j,1)+e(2)*w2(j,2));   %db1為1*3的矩陣
            end
        end
        
        %% 權值更新方程
        w1=w1_1+alpha*dw1'; %w1仍為3*2的矩陣
        b1=b1_1+alpha*db1'; %b1仍為3*1的矩陣
        w2=w2_1+alpha*dw2'; %w2仍為3*2的矩陣
        b2=b2_1+alpha*db2'; %b2仍為2*1的矩陣
        
        %% 保存上一次的權值和閾值
        w1_2=w1_1;w1_1=w1;
        b1_2=b1_1;b1_1=b1;
        w2_2=w2_1;w2_1=w2;
        b2_2=b2_1;b2_1=b2;
    end
end


%% 產生雙螺旋測試數據
%% 產生雙螺旋數據,每類100個樣本點，共200個樣本
i=(1.5:1:test_number+0.5)';    %每類51個樣本

%雙螺旋數據點的產生方程
alpha2=pi*(i-1)/25;
beta2=0.4*((105-i)/104);
m0=0.5+beta2.*sin(alpha2);
n0=0.5+beta2.*cos(alpha2);
s0=zeros(test_number,1);
m1=0.5-beta2.*sin(alpha2);
n1=0.5-beta2.*cos(alpha2);
s1=ones(test_number,1);

test=[m0 n0 s0;m1,n1,s1]; %1條螺旋線數據點,3*102的矩陣
test_label1=test(:,end)';   %測試數據類別，1*102的行向量
test_data1=test(:,1:end-1)'; %測試數據屬性，2*102的矩陣

%把1維的輸出變成2維的輸出,train_labe2為200*2的矩陣
for i=1:2*test_number
    switch test_label1(i)
        case 0
            test_label2(i,:)=[1 0];
        case 1
            test_label2(i,:)=[0 1];
    end
end

test_label=test_label2'; %test_label為2*102的矩陣
         
%%  畫出測試數據雙螺旋曲線
plot(m0,n0,'c+');
hold on;
plot(m1,n1,'yo');
legend('訓練數據螺旋線1','訓練數據螺旋線2','測試數據螺旋線1','測試數據螺旋線2');

test_data=mapminmax('apply',test_data1,train_datas);

% %% 用訓練到的模型對訓練數據本身進行預測
% for i=1:102
%     for j=1:midnum
%         I(j)=train_data(:,i)'*w1(j,:)'+b1(j);
%         Iout(j)=1/(1+exp(-I(j)));%Iout為1*3的行向量
%     end
%     predict(:,i)=w2'*Iout'+b2;%predict為2*102的矩陣
% end
% 
% test_data=mapminmax('apply',train_data1,train_datas);
% test_label=train_label;
% test_label1=trian_label1;

%% 用訓練到的模型預測數據
for i=1:2*test_number
    for j=1:midnum
        I(j)=test_data(:,i)'*w1(j,:)'+b1(j);
        Iout(j)=1/(1+exp(-I(j)));%Iout為1*3的行向量
    end
    predict(:,i)=w2'*Iout'+b2;%predict為2*102的矩陣
end

%% 預測結果分析
for i=1:2*test_number
    output_pred(i)=find(predict(:,i)==max(predict(:,i)));    %out_pred為1*102的矩陣
end

error=output_pred-test_label1-1;    %


%% 計算出每一類預測錯誤的個數總和
k=zeros(1,2); %k=[0 0]
for i=1:2*test_number
    if error(i)~=0    %matlab中不能用if error(i)！=0 
        [b c]=max(test_label(:,i));
        switch c
            case 1
                k(1)=k(1)+1;
            case 2
                k(2)=k(2)+1;
        end
    end
end


%% 求出每一類總體的個數和
kk=zeros(1,2); %k=[0 0]
for i=1:2*test_number
    [b c]=max(test_label(:,i));
    switch c
        case 1
            kk(1)=kk(1)+1;
        case 2
            kk(2)=kk(2)+1;
    end
end


%% 計算每一類的正確率
accuracy=(kk-k)./kk
               

 

實驗總結：

　　可以看出，自己實現的BP網絡其預測效果教matlab自帶的差些，這主要是其中一些算法沒有進行優化，比如說初始權值的選取，自己實現時是隨機產生的，這樣的話不同的權值在相同迭代次數的情況下效果顯然不同。另外本次試驗的訓練樣本數才200個，也比較少，且中間層神經網絡只有5個，進一步影響了試驗的效果。

　　所以下一步的工作是，不管是matlab自帶的工具還是自己寫的函數，對參數的尋優過程要做一定的優化。

 

參考文獻：

1. W Zhao, DS Huang，The structure optimization of radial basis probabilistic neural networks based on genetic algorithms.
2. Matlab中文論壇，MATLAB神經網絡30個案例分析。

大家有更好的解法歡迎提出並交流。