MATLAB神经网络(三):遗传算法优化BP

遗传算法是1962年由Holland提出的模拟自然界遗传机制和生物进化论而成的一泓并行随机搜索最优化方法。它把自然界优胜劣汰、适者生存的生物进化原理引入优化参数形成的编码串联群体中,按照所选择的适应度函数并通过遗传中的选择、交叉和变异对个体进行筛选,使适应度值好的个体被保留,适应度较差的个体被淘汰,新的群体既继承了上一代的信息,又优于上一代。这样反复循环,直至满足条件。遗传算法基本的操作为:
1.选择操作: 从旧群体中以一定概率选择个体到新群体中,个体被选中的概率跟适应度值有关,个体适应度值越好,被选中的概率越大。

2。 交叉操作:从个体中选择两个个体,通过两个染色体的交换组合,来产生新的优秀个体。交叉过程为从群体中任选两个染色体,随机选择一点或多点染色体位置进行交换。交叉操作如图。

3.变异操作:从群体中任选一个个体,选择染色体中的一点进行变异以产生更优秀的个体。

遗传算法具有高效启发式搜索、并行计算等特点,目前已经应用在函数优化、组合优化以及生产调度等方面。

遗传算法的基本要素:染色体编码方法、适应度函数、遗传操作和运行参数。染色体编码方法是指个体的编码方法,目前包括二进制法、实数法等。二进制法把个体编码成为一个二进制串,实数法是指把个体编码成为一个实数串。适应度函数是指根据进化目标编写的计算个体适应度值的函数。遗传操作是指选择操作、交叉操作和变异操作。

运行参数是算法在初始化时确定的参数,主要包括群体大小M,遗传代数G,交叉概率P和变异概率Pm。

我们通过遗传算法优化最基础的BP神经网络,以拟合非线性函数:

[公式]

本例中,由于拟合非线性函数有2个输入参数,1个输出参数,所以设置的BP神经网络结构为2-5-1,共有2*5+5*1个权值,5+1=6个阈值,所以遗传算法个体编码长度为16+5=21。

适应度函数用训练数据训练BP神经网络,并把训练数据预测误差作为个体适应度。

function error=fun(x,inputnum,hiddennum,outputnum,net,inputn,outputn)
%该函数用来计算适应度值
%x input个体
%inputnum input 输入层节点数
%hiddennum input 隐含层节点数
%outputnum input 输出层节点数
%net input 网络
%inputn input 训练输入数据
%outputn input 训练输出数据

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=newff(inputn,outputn,hiddennum);
net.trainParam.epochs=20;
net.trainParam.lr=0.1;
net.trainParam.goal=0.00001;
net.trainParam.show=100;
net.trainParam.showWindow=0;

%神经网络训练
net=train(net,inputn,outputn);
%BP神经网络预测
an=sim(net,inputn);
%预测误差和作为个体适应度值
error=sum(abs(an-outputn));

选择操作采用轮盘赌法从种群中选择适应度好的个体组成新种群。

function ret=select(individuals,sizepop)
%该函数用于进行选择操作
%individuals input 种群信息
%sizepop input 种群规模
%ret output 选择后的新种群
%求适应度值倒数
fitness1=10./individuals.fitness; %individuals.fitness为个体适应度值
%个体选择概率
sumfitness=sum(fitness1);
sumf=fitness1./sumfitness;

%采用轮盘赌法选择新个体
index=[];
for i=1:sizepop %sizepop为种群数
 pick=rand;
 while pick==0
       pick=rand;
  end
  for i=1:sizepop
       pick=pick-sumf(i);
       if pick<0
            index=[index i];
          break;
     end
   end
end
%新种群
individuals.chrom=individuals.chrom(index,:); %individuals.chrom为种群中个体
individuals.fitness=individuals.fitness(index);
ret=individuals;

交叉操作:从种群中选择两个个体,按一定概率交叉得到新个体。

function ret=Cross(pcross,lenchrom,chrom,sizepop,bound)
%该函数用于进行交叉操作
%pcross 交叉概率
%lenchrom 个体长度
%chrom 种群个体
%sizepop 种群规模
%ret 交叉后的新种群

for i=1:sizepop
      %随机选择两个个体
    pick=rand(1,2);
     while prod(pick)==0
       pick=rand(1,2);
    end

    index=ceil(pick,*sizepop)

判断是否交叉
  pick=rand;
   while pick==0
    pick=rand;
  end
    if pick>pcross 
      continue;
  end
 flag=0;
   while flag==0
    %选择交叉位置
     pick=rand;
        while pick==0
          pick=rand;
     end
     pos=ceil(pick,*sum(lenchrom)); %lenchrom为个体长度

%个体交叉
pick=rand;
v1=chrom(index(1),pos);
v2=chrom(index(2),pos);
chrom(index(1),pos)=pick*v2+(1-pick)*v1;
chrom(index(2),pos)=pick*v1+(1-pick)*v2;

%测试新个体是否满足约束要求
flag1=test(lenchrom,bound,chrom(index(1),:));
flag2=test(lenchrom,bound,chrom(index(2),:));
    if flag1*flag2==0
     flag=0;
      else flag=1;
     end
    end
   end
ret=chrom

变异操作:从种群中随机选择一个个体,按一定概率变异得到新个体。

function ret=Mutation(pmutation,lenchrom,chrom,sizepop,num,maxgen,bound)
%pmutation 变异概率
%lenchrom 个体长度
%chrom 种群个体
%sizepop 种群规模
%bound 个体上界和下界
%maxgen 最大迭代次数
%num 当前迭代次数
%ret 交叉后的新种群

for i=1:sizepop
     pick=rand;
     while pick==0
          pick=rand;
     end
     index=ceil(pick*sizepop);
     %判断是否变异
     pick=rand;
    if pick>pmutation %pmutation 为变异概率
        continue;
      end
     flag=0;
     while flag==0
          %随机选择变异位置
          pick=rand;
         while pick==0
               pick=rand;
          end
         pos=ceil(pick*sum(lenchrom)); %lenchrom为个体长度
         %变异操作
         v=chrom(i,pos);
         v1=v-bound(pos,1);
         v2=bound(pos,2)-v;
        pick=rand;
        fg=(rand*(1-num/maxgen)^2;
        if pick>0.5
             chrom(i,pos)=chrom(i,pos)+(bound(pos,2)-chrom(i,pos))*fg;
        else
             chrom(i,pos)=chrom(i,pos)-(chrom(i,pos)-bound(pos,1))*fg;
        end
        flag=test(lenchrom,bound,chrom(i,:));   %新个体是否满足约束要求
     end
end
ret=chrom

遗传算法主函数

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)';

%数据归一化
intput_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(intputn,outputn,hiddennum);

%遗传算法参数初始化
maxgen=50; %迭代次数
sizepop=10;  %种群规模
pcross=[0.4]; %交叉概率
pmutation=[0.2];   %变异概率
numsum=inputnum*hiddennum+hiddennum+hiddenum*outputnum+outputnum;
lenchrom=cones(1,numsum);    %个体长度
bound=[-3*ones(numsum,1)3*ones(numsum,1)]; %个体范围

%种群信息定义为结构体
indivisuals=struct('fitness',zeros(1,sizepop),'chrom',[]);
avgfitness=[]; %每代平均适应度值
bestfitness=[]; %每代最佳适应度值
bestchrom=[];  %最优个体

%计算个体适应度值

    for i=1:sizepop
             %个体初始化
             individuals.chrom(i,:)=Code(lenchrom,bound);
                
            %计算个体适应度值
            x=individuals.chrom(i,:);
            individuals.fitness(i)=fun(x,inputnum,hiddennum,outputnum,net,inputn,outputn);
    end
    %迭代寻优
    for i=1:maxgen
            
           %选择操作
           individuals=Select(individuals,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);
           [newworstfitness,newworstindex]=max(individuals,fitness);
           %最优个体更新
           if bestfitness>newbestfitness
                    bestfitness=newbestfitness;
                    bestchrom=individuals.chrom(newbestindex,:);
           end
           individuals.chrom(newworstindex,:)=bestchrom
           individuals.fitness(newworstindex)=bestfitness;
           
           %记录最优个体适应度值和平均适应度值
avgfitness=sum(individuals.fitness)/sizepop;
trace=[trace;avgfitness bestfitness];
end
%把最优个体赋给BP神经网络,用该网络拟合非线性函数
x=bestchrom
w1=x(1:inputnum*hiddennum);
B1=x(inputnum*hiddennum+1:inputnum*hiddennum+hiddennum);
w2=x(inputnum*hiddennum+hiddennum+1:inputnnum*hiddennum+hiddenum+hiddennum*outputnum);
B2=x(inputnnum*hiddennum+hiddennum+hiddennum*outputnum+1:inputnnum*hiddennum+hiddennum+hiddennum*outputnum+outputnum);

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

%神经网络参数
net.trainParam.epochs=100;
net,trainParam.lr=0.1;
%net.trainParam.goal=0.00001;

%BP神经网络训练
[net,per2]=train(net,inputn,outputn);

%BP神经网络预测
inputn_test=mapminmax('apply',input)test,inputps);
an=sim(net,input_test);

%预测结果反归一化
test_simu=mapminmax('reverse',an,outputps);