遗传算法MATLAB程序设计.docx
《遗传算法MATLAB程序设计.docx》由会员分享，可在线阅读，更多相关《遗传算法MATLAB程序设计.docx（17页珍藏版）》请在冰豆网上搜索。
遗传算法MATLAB程序设计
摘自Matlab在数学建模中的应用，北航出版社，2011.4
4.2遗传算法MATLAB程序设计
4.2.1程序设计流程及参数选取
4.2.1.1遗传算法程序设计伪代码
BEGIN
t=0;%GenerationsNO.
初始化P（t）;%InitialPopulationorChromosomes
计算P（t）的适应值;
while（不满足停止准则）do
begin
t=t+1;
从P（t-1）中选择P（t）;%Selection
重组P（t）;%CrossoverandMutation
计算P（t）的适应值;
end
END
4.2.1.2遗传算法的参数设计原则
在单纯的遗传算法当中，也并不总是收敛，即使在单峰或单调也是如此。
这是因为种群的进化能力已经基本丧失，种群早熟。
为了避免种群的早熟，参数的设计一般遵从以下原则[5]：
（1）种群的规模：
当群体规模太小时，很明显会出现近亲交配，产生病态基因。
而且造成有效等位基因先天缺乏，即使采用较大概率的变异算子，生成具有竞争力高阶模式的可能性仍很小，况且大概率变异算子对已有模式的破坏作用极大。
同时遗传算子存在随机误差（模式采样误差），妨碍小群体中有效模式的正确传播，使得种群进化不能按照模式定理产生所预测的期望数量；种群规模太大，结果难以收敛且浪费资源，稳健性下降。
种群规模的一个建议值为0~100。
（2）变异概率：
当变异概率太小时，种群的多样性下降太快，容易导致有效基因的迅速丢失且不容易修补；当变异概率太大时，尽管种群的多样性可以得到保证，但是高阶模式被破坏的概率也随之增大。
变异概率一般取0.0001~0.2。
（3）交配概率：
交配是生成新种群最重要的手段。
与变异概率类似，交配概率太大容易破坏已有的有利模式，随机性增大，容易错失最优个体；交配概率太小不能有效更新种群。
交配概率一般取0.4~0.99。
（4）进化代数：
进化代数太小，算法不容易收敛，种群还没有成熟；代数太大，算法已经熟练或者种群过于早熟不可能再收敛，继续进化没有意义，只会增加时间开支和资源浪费。
进化代数一般取100~500。
（5）种群初始化：
初始种群的生成是随机的；在初始种群的赋予之前，尽量进行一个大概的区间估计，以免初始种群分布在远离全局最优解的编码空间，导致遗传算法的搜索范围受到限制，同时也为算法减轻负担。
4.2.1.3适应度函数的调整
（1）在遗传算法运行的初期阶段
群体中可能会有少数几个个体的适应度相对其他个体来说非常高。
若按照常用的比例选择算子来确定个体的遗传数量时，则这几个相对较好的个体将在下一代群体中占有很高的比例，在极端情况下或当群体现模较小时，新的群体甚至完全由这样的少数几个个体所组成。
这时交配运算就起不了什么作用，因为相同的两个个体不论在何处发生交叉行为都永远不会产生新的个体。
这样就会使群体的多样性降低，容易导致遗传算法发生早熟现象（或称早期收敛），使遗传算法所求到的解停留在某一局部最优点上。
因此，我们希望在遗传算法运行的初期阶段，算法能够对一些适应度较高的个体进行控制，降低其适应度与其他个体适应度之间的差异程度，从而限制其复制数量，以维护群体的多样性。
（2）在遗传算法运行的后期阶段
群体中所有个体的平均适应度可能会接近于群体中最佳个体的适应度。
也就是说，大部分个体的适应度和最佳个体的适应度差异不大，它们之间无竞争力，都会有以相接近的概率被遗传到下一代的可能性，从而使得进化过程无竞争性可言，只是一种随机的选择过程。
这将导致无法对某些重点区域进行重点搜索，从而影响遗传算法的运行效率。
因此，我们希望在遗传算法运行的后期阶段，算法能够对个体的适应度进行适当的放大，扩大最佳个体适应度与其他个体适应度之间的差异程度，以提高个体之间的竞争性。
4.2.2MATLAB遗传算法工具箱
4.2.2.1GA工具箱版本
为了省略艰深难懂的遗传算法数学理论和降低程序开发的难度，MATLAB软件已经将遗传算法命令进行了封装，做成专门的遗传算法工具箱（GAToolbox），方便用户调用。
关于遗传算法工具箱，需要说明三点：
①目前基于MATLAB环境下遗传算法工具箱的版本较多，各版本的功能和用法也不完全相同，需要加以区分。
倘若想要使用某个工具箱，但是MATLAB没有自带，则用户需要自行安装。
目前GAToolbox主要有下列几种：
ØMATLAB7.xVersion
GeneticAlgorithmandDirectSearchToolbox2.0.1；
ØMATLAB6.x（or7.x）Version+GAOT
GAOT：
GeneticAlgorithmOptimizationToolbox；
美国NorthCarolinaStateUniversity开发；
ØMATLAB6.x（or7.x）Version+GEATbx
GEATbx:
GeneticandEvolutionaryAlgorithmToolbox。
英国TheUniversityofSheffield开发，目前此工具箱功能最完善性能最优越。
②少数非常规、非主流的盗版MATLAB软件，其自身所带的遗传算法工具箱可能会没有。
原因在于盗版软件在破解安装注册序列号过程中，主程序可以正常使用，但是部分工具箱功能在破解过程中被丢失。
不同的破解版本，工具箱丢失的种类和程度会有所差异。
③遗传算法工具箱已经将常用的遗传运算命令进行了集成，用户使用很方便。
但是封装的工具箱的内部命令不能根据特殊需要进行相应调整和修改。
从这个角度上说，具有人工智能性质的GAToolbox是一种“傻瓜式的智能”或“智能式的傻瓜”。
作者建议尽量不要使用工具箱，尝试编写遗传算法源程序去解决问题。
本章基于MATLAB7.0Version自带的遗传算法与直接搜索工具箱（GeneticAlgorithmandDirectSearchToolbox）来优化目标函数。
使用的方法是：
在命令窗口输入>>gatool指令然后在弹出的GUI界面上填写相关参数和函令句柄或使用M文件进行程序设计。
遗传算法与直接搜索工具箱有ga和gaoptimset两个核心函数。
4.2.2.2工具箱核心函数用法
（1）ga函数
【语法格式】[x,fval,reason]=ga（@fitnessfun,nvars,options）
其中，x——经过遗传进化以后自变量最佳染色体返回值；
fval——最佳染色体的适应度；
reason——算法停止的原因；
@fitnessfun——适应度句柄函数；
nvars——目标函数自变量的个数；
options——算法的属性设置，该属性是通过函数gaoptimset赋予的。
【实现功能】对目标函数进行遗传运算。
（2）gaoptimset函数
【语法格式】options=gaoptimset（‘PropertyName1’,’PropertyValue1’,
‘PropertyName2’,’PropertyValue2’,‘PropertyName3’,’PropertyValue3’……）
【实现功能】设置遗传算法的参数和句柄函数，表4-4介绍常用的11种属性。
由于遗传算法本质上是一种启发式的随机运算，算法程序经常重复运行多次才能得到理想结果。
鉴于此，可以将前一次运行得到的最后种群作为下一次运行的初始种群，如此操作会得到更好的结果：
[x,fval,reason,output,final_pop]=ga（@fitnessfcn,nvars）;
表4-4gaoptimset函数设置属性
序号
属性名
默认值
实现功能
1
PopInitRange
[0;1]
初始种群生成区间
2
PopulationSize
20
种群规模
3
CrossoverFraction
0.8
交配概率
4
MigrationFraction
0.2
变异概率
5
Generations
100
超过进化代数时算法停止
6
TimeLimit
Inf
超过运算时间限制时算法停止
7
FitnessLimit
-Inf
最佳个体等于或小于适应度阈值时算法停止
8
StallGenLimit
50
超过连续代数不进化则算法停止
9
StallTimeLimit
20
超过连续时间不进化则算法停止
10
InitialPopulation
[]
初始化种群
11
PlotFcns
[]
绘图函数，可供选择有@gaplotbestf@gaplotbestindiv等
最后一个输出变量final_pop返回的就是本次运行得到的最后种群。
再将final_pop作为ga函数的初始种群，语法格式为：
options=gaoptimset（'InitialPopulation',finnal_pop）;
[x,fval,reason,output,finnal_pop2]=ga（@fitnessfcn,nvars,options）;
4.2.2.3GeneticAlgorithmandDirectSearchToolbox适应度函数设计
Case1遗传算法和直接搜索工具箱中的ga函数是求解目标函数的最小值，所以求目标函数最小值的问题，可直接令目标函数为适应度函数（注意这里是求最小值，与上节求最大的例子恰好相反）。
编写适应度函数，语法格式如下：
functionf=fitnessfcn（x）%x为自变量向量
f=f（x）;
Case2如果有约束条件（包括自变量的取值范围），对于求解函数的最小值问题，可以使用如下语法格式：
functionf=fitnessfcn（x）
if（x<=-1）|x>3）%表示有约束x>-1和x<=3，其它约束条件类推
f=inf;
else
f=f（x）;
Case3如果有约束条件（包括自变量的取值范围），对于求解函数的最大值问题，可以使用如下语法格式：
functionf=fitnessfcn（x）
if（x<=-1）|x>3）
f=inf;
else
f=-f（x）;%注意这里f=-f（x）而不是f=f（x）
若目标函数作为适应度函数，则最终得到的目标函数值为-fval而不是fval。
4.2.2.4程序设计范例
求解目标函数
图4-4目标函数图像
此函数在[-30,30]区间上都只有一个最小值，这个理论最小值为：
图4-4为N=1和N=2时，目标函数的图像。
现在利用GA工具箱求解上述函数N=10时的最小值，也就是有10个变量，
10个约束条件
，具体求解过程如下。
先编写M文件：
%if和else后面的语句不必分行书写
functionf=lbw（x）%存储M文件名为lbw.m
if（x
（1）>30|x
（1）<-30|x
（2）>30|x
（2）<-30|x（3）>30|x（3）<-30|x（4）>30|x（4）<-30|x（5）>30|x（5）<-30|x（6）>30|x（6）<-30|x（7）>30|x（7）<-30|x（8）>30|x（8）<-30|x（9）>30|x（9）<-30|x（10）>30|x（10）<-30）
f=300;
else
f=-2*pi*exp（-0.2*sqrt（1/10*（（x
（1））.^2+（x
（2））.^2+（x（3））.^2+（x（4））.^2+（x（5））.^2+（x（6））.^2+（x（7））.^2+（x（8））.^2+（x（9））.^2+（x（10））.^2）））-exp（1/10*（cos（2*pi*x
（1））+cos（2*pi*x
（2））+cos（2*pi*x（3））+cos（2*pi*x（4））+cos（2*pi*x（5））+cos（2*pi*x（6））+cos（2*pi*x（7））+cos（2*pi*x（8））+cos（2*pi*x（9））+cos（2*pi*x（10））））+2*pi;
end
然后在命令窗口输入：
>>clearall
>>options=gaoptimset（'Generations',800,'StallGenLimit',100,'PlotFcns',@gaplotbestf;
>>[x,f]=ga（@lbw,10,options）
程序输出的结果为：
Optimizationterminated:
stalltimelimitexceeded.
x=
0.00610.00460.00330.0116-0.0051-0.00350.0012-0.01180.00560.0053
f=
-2.7075
计算的结果表明，与理论值相比，误差只有0.3966%，不到1%，误差很小。
图4-5为遗传算法在进化过程中最佳个体适应度值和种群平均适应度值变化趋势图，从图4-5可以看出种群没有进化到800代就停止了，原因是种群连续100代最佳个体的适应度值都没有变化，说明种群足够成熟了，算法自行停止。
图4-5适应度变化趋势
4.3遗传算法应用案例
4.3.1案例一无约束目标函数最大值遗传算法求解策略
遗传算法虽然有工具箱，但是工具箱并不是万能的，很多情况下需要具体问题具体对待。
另外，使用工具箱无助于从本质上理解遗传算法的精义所在。
以下所举案例使用的MATLAB代码是用源程序编写，读者尝试逐条读懂语句，方能针对复杂的多约束非线性规划问题设计出对口的程序。
本案例待求解
问题：
用MATLAB画出该函数在[-2,2]区间的图形，如图4-6所示。
利用MATLAB自带工具“DataCursor”菜单探测出
在区间
上最大值大概位于（1.534，185.1）坐标点附近，具体求解程序如P4-1。
图4-6目标函数最大值示意图
本程序取运算精度precision=0.0001，初始种群popsize=50，进化代数Generationnmax=12，交配概率pcrossover=0.90，变异概率pmutation=0.09，经过12次迭代，程序输出的结果为：
Bestpopulation=
1.5205
Besttargetfunvalue=
185.1242
程序编号
P4-1
文件名称
GA401.m
说明
遗传算法源程序求解最大值
%主程序：
用遗传算法求解y=200*exp（-0.05*x）.*sin（x）在[-22]区间上的最大值
clc;
clearall;
closeall;
globalBitLength
globalboundsbegin
globalboundsend
bounds=[-22];%一维自变量的取值范围
precision=0.0001;%运算精度
boundsbegin=bounds（:
1）;
boundsend=bounds（:
2）;
%计算如果满足求解精度至少需要多长的染色体
BitLength=ceil（log2（（boundsend-boundsbegin）'./precision））;
popsize=50;%初始种群大小
Generationnmax=12;%最大代数
pcrossover=0.90;%交配概率
pmutation=0.09;%变异概率
%产生初始种群
population=round（rand（popsize,BitLength））;
%计算适应度,返回适应度Fitvalue和累积概率cumsump
[Fitvalue,cumsump]=fitnessfun（population）;
Generation=1;
whileGenerationforj=1:
2:
popsize
%选择操作
seln=selection（population,cumsump）;
%交叉操作
scro=crossover（population,seln,pcrossover）;
scnew（j,:
）=scro（1,:
）;
scnew（j+1,:
）=scro（2,:
）;
%变异操作
smnew（j,:
）=mutation（scnew（j,:
）,pmutation）;
smnew（j+1,:
）=mutation（scnew（j+1,:
）,pmutation）;
end
population=smnew;%产生了新的种群
%计算新种群的适应度
[Fitvalue,cumsump]=fitnessfun（population）;
%记录当前代最好的适应度和平均适应度
[fmax,nmax]=max（Fitvalue）;
fmean=mean（Fitvalue）;
ymax（Generation）=fmax;
ymean（Generation）=fmean;
%记录当前代的最佳染色体个体
x=transform2to10（population（nmax,:
））;
%自变量取值范围是[-22],需要把经过遗传运算的最佳染色体整合到[-22]区间
xx=boundsbegin+x*（boundsend-boundsbegin）/（power（（boundsend）,BitLength）-1）;
xmax（Generation）=xx;
Generation=Generation+1
end
Generation=Generation-1;
Bestpopulation=xx
Besttargetfunvalue=targetfun（xx）
%绘制经过遗传运算后的适应度曲线。
一般地，如果进化过程中种群的平均适应度与最大适
%应度在曲线上有相互趋同的形态，表示算法收敛进行得很顺利，没有出现震荡；在这种前
%提下，最大适应度个体连续若干代都没有发生进化表明种群已经成熟。
figure
（1）;
hand1=plot（1:
Generation,ymax）;
set（hand1,'linestyle','-','linewidth',1.8,'marker','*','markersize',6）
holdon;
hand2=plot（1:
Generation,ymean）;
set（hand2,'color','r','linestyle','-','linewidth',1.8,...
'marker','h','markersize',6）
xlabel（'进化代数'）;ylabel（'最大/平均适应度'）;xlim（[1Generationnmax]）;
legend（'最大适应度','平均适应度'）;
boxoff;holdoff;
%子程序：
新种群交叉操作,函数名称存储为crossover.m
functionscro=crossover（population,seln,pc）;
BitLength=size（population,2）;
pcc=IfCroIfMut（pc）;%根据交叉概率决定是否进行交叉操作，1则是，0则否
ifpcc==1
chb=round（rand*（BitLength-2））+1;%在[1,BitLength-1]范围内随机产生一个交叉位
scro（1,:
）=[population（seln
（1）,1:
chb）population（seln
（2）,chb+1:
BitLength）];
scro（2,:
）=[population（seln
（2）,1:
chb）population（seln
（1）,chb+1:
BitLength）];
else
scro（1,:
）=population（seln
（1）,:
）;
scro（2,:
）=population（seln
（2）,:
）;
end
%子程序：
计算适应度函数,函数名称存储为fitnessfun
function[Fitvalue,cumsump]=fitnessfun（population）;
globalBitLength
globalboundsbegin
globalboundsend
popsize=size（population,1）;%有popsize个个体
fori=1:
popsize
x=transform2to10（population（i,:
））;%将二进制转换为十进制
%转化为[-2,2]区间的实数
xx=boundsbegin+x*（boundsend-boundsbegin）/（power（（boundsend）,BitLength）-1）;
Fitvalue（i）=targetfun（xx）;%计算函数值，即适应度
end
%给适应度函数加上一个大小合理的数以便保证种群适应值为正数
Fitvalue=Fitvalue'+230;
%计算选择概率
fsum=sum（Fitvalue）;
Pperpopulation=Fitvalue/fsum;
%计算累积概率
cumsump
（1）=Pperpopulation
（1）;
fori=2:
popsize
cumsump（i）=cumsump（i-1）+Pperpopulation（i）;
end
cumsump=cumsump';
%子程序：
新种群变异操作，函数名称存储为mutation.m
functionsnnew=mutation（snew,pmutation）;
BitLength=size（snew,2）;
snnew=snew;
pmm=IfCroIfMut（pmutation）;%根据变异概率决定是否进行变异操作，1则是，0则否
ifpmm==1
chb=round（rand*（BitLength-1））+1;%在[1,BitLength]范围内随机产生一个变异位
snnew（chb）=abs（snew（chb）-1）;
end
%子程序：
判断遗传运算是否需要进行交叉或变异,函数名称存储为IfCroIfMut.m
functionpcc=IfCroIfMut（mutORcro）;
test（1:
100）=0;
l=round（100*mutORcro）;
test（1:
l）=1;
n=round（rand*99）+1;
pcc=test（n）;
%子程序：
新种群选择操作,函数名称存储为selection.m
functionseln=selection（population,cumsump）;
%从种群中选择两个个体
fori=1:
2
r=rand;%产生一个随机数
prand=cumsump-r;
j=1;
whileprand（j）<0
j=j+1;
end
seln（i）=j;%选中个体的序号
end
%子程序：
将2进制数转换为10进制数,函数名称存储为transform2to10.m
functionx=transform2to10（Population）;
BitLength=size（Population,2）;
x=Population（BitLength）;
fori=1:
BitLe
升级会员 