【神經網絡學習筆記】遺傳算法優化BP神經網絡-非線性函數拟合

我們知道，在建立神經網絡的時候，需要給各個節點和連接配接賦予門檻值和權值，但一般我們都是随機賦予，然後讓系統在一次次的訓練中需找到最小值。這種方法很有一種碰運氣的成分在裡面，而用遺傳算法來優化BP神經網絡則可以讓神經網絡從一個較接近的水準開始訓練。

所謂用遺傳算法優化，優化的是神經網絡訓練前各節點的初始值，這些初始值不再取随機值，而是用我們通過遺傳算法得到的值來代替。

我們把每個物種賦予一個DNA序列，這個序列包含着各個門檻值和權值，假設我們的網絡是一個2-5-1的網絡，那麼DNA的長度就可以是2*5+5+5*1+1 = 21，每位表示一個門檻值或權值。

再設定好這個種群的大小，即有多少個個體，然後就可以開始進化這個物種啦。

進化過程不短重複這幾步

1.計算個體的适應度，即代入它的基因（門檻值權值）計算所得到的與期望值的內插補點。

2.記錄适應度最好即內插補點最小的個體。

3.輪盤法選擇個體，适應度好的更容易被選到。

4.把選擇出來的個體進行交叉和變異，相當于生物的有性繁殖和基因突變。

到了設定的遺傳代數之後，得到的最優個體基因就比較接近最終的門檻值權值啦，我們用這個個體的基因來初始化BP神經網絡，再進行訓練和預測，可以得到更好的效果。

下面是主函數的代碼，其中的子函數請下載下傳資源。

%% 該代碼為基于遺傳算法神經網絡的預測代碼
% 清空環境變量
clc
clear
% 
%% 網絡結建構立
%讀取資料
load data input output

%節點個數
inputnum=2;
hiddennum=5;
outputnum=1;

%訓練資料和預測資料
input_train=input(1:1900,:)';
input_test=input(1901:2000,:)';
output_train=output(1:1900)';
output_test=output(1901:2000)';

%選連樣本輸入輸出資料歸一化
[inputn,inputps]=mapminmax(input_train);
[outputn,outputps]=mapminmax(output_train);

%建構網絡
net=newff(inputn,outputn,hiddennum);

%% 遺傳算法參數初始化
maxgen=20;                         %進化代數，即疊代次數
sizepop=10;                        %種群規模
pcross=[0.2];                       %交叉機率選擇，0和1之間
pmutation=[0.1];                    %變異機率選擇，0和1之間

%節點總數
numsum=inputnum*hiddennum+hiddennum+hiddennum*outputnum+outputnum;

lenchrom=ones(1,numsum);        
bound=[-3*ones(numsum,1) 3*ones(numsum,1)];    %資料範圍

%------------------------------------------------------種群初始化--------------------------------------------------------
individuals=struct('fitness',zeros(1,sizepop), 'chrom',[]);  %将種群資訊定義為一個結構體
avgfitness=[];                      %每一代種群的平均适應度
bestfitness=[];                     %每一代種群的最佳适應度
bestchrom=[];                       %适應度最好的染色體
%初始化種群
for i=1:sizepop
    %随機産生一個種群
    individuals.chrom(i,:)=Code(lenchrom,bound);    %編碼（binary和grey的編碼結果為一個實數，float的編碼結果為一個實數向量）
    x=individuals.chrom(i,:);
    %計算适應度
    individuals.fitness(i)=fun(x,inputnum,hiddennum,outputnum,net,inputn,outputn);   %染色體的适應度
end
FitRecord=[];
%找最好的染色體
[bestfitness bestindex]=min(individuals.fitness);
bestchrom=individuals.chrom(bestindex,:);  %最好的染色體
avgfitness=sum(individuals.fitness)/sizepop; %染色體的平均适應度
% 記錄每一代進化中最好的适應度和平均适應度
trace=[avgfitness bestfitness]; 
 
%% 疊代求解最佳初始閥值和權值
% 進化開始
for i=1:maxgen
    i
    % 選擇
    individuals=Select(individuals,sizepop);
    avgfitness=sum(individuals.fitness)/sizepop;
    %交叉
    individuals.chrom=Cross(pcross,lenchrom,individuals.chrom,sizepop,bound);
    % 變異
    individuals.chrom=Mutation(pmutation,lenchrom,individuals.chrom,sizepop,i,maxgen,bound);
    
    % 計算适應度 
    for j=1:sizepop
        x=individuals.chrom(j,:); %解碼
        individuals.fitness(j)=fun(x,inputnum,hiddennum,outputnum,net,inputn,outputn);   
    end
    
  %找到最小和最大适應度的染色體及它們在種群中的位置
    [newbestfitness,newbestindex]=min(individuals.fitness);
    [worestfitness,worestindex]=max(individuals.fitness);
    % 代替上一次進化中最好的染色體
    if bestfitness>newbestfitness
        bestfitness=newbestfitness;
        bestchrom=individuals.chrom(newbestindex,:);
    end
    individuals.chrom(worestindex,:)=bestchrom;
    individuals.fitness(worestindex)=bestfitness;
    
    avgfitness=sum(individuals.fitness)/sizepop;
    
    trace=[trace;avgfitness bestfitness]; %記錄每一代進化中最好的适應度和平均适應度
    FitRecord=[FitRecord;individuals.fitness];
end

%% 遺傳算法結果分析 
figure(1)
[r c]=size(trace);
plot([1:r]',trace(:,2),'b--');
title(['适應度曲線  ' '終止代數＝' num2str(maxgen)]);
xlabel('進化代數');ylabel('适應度');
legend('平均适應度','最佳适應度');
disp('适應度                   變量');

%% 把最優初始閥值權值賦予網絡預測
% %用遺傳算法優化的BP網絡進行值預測
w1=x(1:inputnum*hiddennum);
B1=x(inputnum*hiddennum+1:inputnum*hiddennum+hiddennum);
w2=x(inputnum*hiddennum+hiddennum+1:inputnum*hiddennum+hiddennum+hiddennum*outputnum);
B2=x(inputnum*hiddennum+hiddennum+hiddennum*outputnum+1:inputnum*hiddennum+hiddennum+hiddennum*outputnum+outputnum);

net.iw{1,1}=reshape(w1,hiddennum,inputnum);
net.lw{2,1}=reshape(w2,outputnum,hiddennum);
net.b{1}=reshape(B1,hiddennum,1);
net.b{2}=B2;

%% BP網絡訓練
%網絡進化參數
net.trainParam.epochs=100;
net.trainParam.lr=0.1;
%net.trainParam.goal=0.00001;

%網絡訓練
[net,per2]=train(net,inputn,outputn);

%% BP網絡預測
%資料歸一化
inputn_test=mapminmax('apply',input_test,inputps);
an=sim(net,inputn_test);
test_simu=mapminmax('reverse',an,outputps);
error=test_simu-output_test;
           

END