BP神經網絡 MATLAB源程序

和之前的習慣同樣，只舉典例，而後給程序和運行結果進行說明。

問題背景是：

給定某地區20年的數據，6列，21行，第一列值爲年份，第二列爲人數，第三列爲機動車數量，第四列爲公路面積，第五列爲公路客運量，第六列爲公路貨運量，這20年是1990年到2009年，如今給咱們2010和2011年，第2、三和四列的數據，讓咱們用BP網絡預測該地區2010年和2011年公路的客運量和公路貨運量，也就是第五和六列的值。

問題背景詳見《MATLAB在數學建模中的應用》卓金武第二版 134頁 或者聯繫我。

其中，CUMCM 2006 B題 艾滋病治療最佳停藥時間的肯定也可用此法，程序改改數據，還有神經網絡輸入輸出層和隱含層的部分便可，詳情多看書。

代碼：

function main()
clc             % 清屏
clear all;      % 清除內存以便加快運算速度
close all;      % 關閉當前全部figure圖像
warning off;    % 屏蔽沒有必要的警告
SamNum=20;      % 輸入樣本數量爲20
TestSamNum=20;  % 測試樣本數量也是20
ForcastSamNum=2;% 預測樣本數量爲2
HiddenUnitNum=8;% 中間層隱節點數量取8
InDim=3;        % 網絡輸入維度爲3
OutDim=2;       % 網絡輸出維度爲2

% 原始數據
% 人數
sqrs=[20.55 22.44 25.37 27.13 29.45 30.10 30.96 34.06 36.42 38.09 39.13 39.99 41.93 44.59 47.30 52.89 55.73 56.76 59.17 60.63];
% 機動車數量
sqjdcs=[0.6 0.75 0.85 0.9 1.05 1.35 1.45 1.6 1.7 1.85 2.15 2.2 2.25 2.35 2.5 2.6 2.7 2.85 2.95 3.1];
sqglmj=[0.09 0.11 0.11 0.14 0.20 0.23 0.23 0.32 0.32 0.34 0.36 0.36 0.38 0.49 0.56 0.59 0.59 0.67 0.69 0.79];
% 公路客運量
glkyl=[5126 6217 7730 9145 10460 11387 12353 15750 18304 19836 21024 19490 20433 22598 25107 33442 36836 40548 42927 43467];
% 公路貨運量
glhyl=[1237 1379 1385 1399 1663 1714 1834 4322 8132 8936 11099 11203 10524 11115 13320 16762 18673 20724 20803 21804];

p=[sqrs;sqjdcs;sqglmj];  % 輸入數據矩陣
t=[glkyl;glhyl];         % 目標數據矩陣
[SamIn,minp,maxp,tn,mint,maxt]=premnmx(p,t); % 原始樣本對（輸入和輸出）初始化

rand('state',sum(100*clock));   % 依據系統時鐘種子產生隨機數
NoiseVar=0.01;                  % 噪聲強度爲0.01（添加噪聲的目的是爲了防止網絡過分擬合）
Noise=NoiseVar*randn(2,SamNum); % 生成噪聲
SamOut=tn+Noise;                % 將噪聲添加到輸出樣本上

TestSamIn=SamIn;                % 這裏取輸入樣本與測試樣本相同，由於樣本容量偏少
TestSanOut=SamOut;              % 也取輸出樣本與測試樣本相同

MaxEpochs=50000;                % 最多訓練次數爲50000
lr=0.035;                       % 學習速率爲0.035     
E0=0.65*10^(-3);                % 目標偏差爲0.65*10^(-3)
W1=0.5*rand(HiddenUnitNum,InDim)-0.1;% 初始化輸入層與隱含層之間的權值
B1=0.5*rand(HiddenUnitNum,1)-0.1;% 初始化輸入層與隱含層之間的權值
W2=0.5*rand(OutDim,HiddenUnitNum)-0.1;% 初始化輸出層與隱含層之間的權值
B2=0.5*rand(OutDim,1)-0.1;% 初始化輸出層與隱含層之間的權值

ErrHistory=[];  % 給中間變量預先佔據內存
for i=1:MaxEpochs          
    HiddenOut=logsig(W1*SamIn+repmat(B1,1,SamNum)); % 隱含層網絡輸出
    NetworkOut=W2*HiddenOut+repmat(B2,1,SamNum);   %輸出層網絡輸出
    Error=SamOut-NetworkOut;  % 實際輸出與網絡輸出之差
    SSE=sumsqr(Error);   % 能量函數(偏差平方和)
    ErrHistory=[ErrHistory SSE];
    if SSE<E0,break,end  % 若是達到偏差要求則跳出學習循環
    
    % 如下6行是BP網絡最核心的程序
    % 它們是權值（闕值）依據能量函數負梯度降低原理所作的每一步動態調整
    Delta2=Error;
    Delta1=W2'*Delta2.*HiddenOut.*(1-HiddenOut);
    % 對輸出層與隱含層之間的權值和闕值進行修正
    dW2=Delta2*HiddenOut';
    dB2=Delta2*ones(SamNum,1);
    % 對輸入層與隱含層之間的權值和闕值進行修正
    dW1=Delta1*SamIn';
    dB1=Delta1*ones(SamNum,1);
    
    W2=W2+lr*dW2;
    B2=B2+lr*dB2;
    
    W1=W1+lr*dW1;
    B1=B1+lr*dB1;
    
end

HiddenOut=logsig(W1*SamIn+repmat(B1,1,TestSamNum));  % 隱含層輸出最終結果
NetworkOut=W2*HiddenOut+repmat(B2,1,TestSamNum);     % 輸出層輸出最終結果
a=postmnmx(NetworkOut,mint,maxt);                    % 還原網絡輸出層的結果
x=1990:2009;                                         % 時間軸刻度
newk=a(1,:);                                         % 網絡輸出客運量
newh=a(2,:);                                         % 網絡輸出貨運量
figure;
subplot(2,1,1);plot(x,newk,'r-o',x,glkyl,'b--+');  % 繪製公路客運量對比圖
legend('網絡輸出客運量','實際客運量');
xlabel('年份'); ylabel('客運量/萬人');
title('源程序神經網絡客運量學習和測試對比圖');

subplot(2,1,2);plot(x,newh,'r-o',x,glhyl,'b--+');  % 繪製公路貨運量對比圖
legend('網絡輸出貨運量','實際貨運量');
xlabel('年份'); ylabel('貨運量/萬人');
title('源程序神經網絡貨運量學習和測試對比圖');

% 利用訓練好的數據進行預測
% 當用訓練好的網絡對新數據pnew進行預測時，也應作相應的處理
pnew=[73.39 75.55 
    3.9635 4.0975
    0.9880 1.0268];   % 2010年和2011年的相關數據
pnewn=tramnmx(pnew,minp,maxp);  %利用原始輸入數據的歸一化參數對新數據進行歸一化
HiddenOut=logsig(W1*pnewn+repmat(B1,1,ForcastSamNum));  % 隱含層輸出預測結果
anewn=W2*HiddenOut+repmat(B2,1,ForcastSamNum);          % 輸出層輸出預測結果
% 把網絡預測獲得的數據還原爲原始的數量級
format short
anew=postmnmx(anewn,mint,maxt)

運行結果：

anew =
  1.0e+004 *
    4.1324    4.0941
    2.1521    2.1440

實際樣本與網絡輸出的對比圖以下所示：

版權聲明：本文爲博主原創文章，未經博主容許不得轉載。