platEMO：一款强大的多目标优化工具（MATLAB）的扩展操作-操作实例

platEMO：⼀款强⼤的多⽬标优化⼯具（MATLAB）的扩展操作——操作实例platEMO：⼀款强⼤的多⽬标优化⼯具（MATLAB）的扩展操作——操作实例前⾔增加问题每⼀个多⽬标优化问题在platEMO中由⼀个.m⽂件表⽰，该⽂件应该放在⽂件夹Problems中。即我们⾸先需要在Problems这个⽂件夹新建⼀个我们问题的⽂件夹。我们定义⼀个多⽬标优化问题⼀般需要考虑他的问题个数Global.M，变量的个数Global.D，变量的上界Global.upper，变量的下界Global.lower，以及变量的编码⽅式Global.encoding，最后就是⽬标函数的形式(包括约束条件)了。以下我们以带约束的DTLZ1()函数为例说明问题是如何增加的。classdefC1_DTLZ1<PROBLEM%<problem><DTLZvariant>%如果你想在GUI界⾯去调⽤你的问题，这⾥<>必须起⼀个你问题⾃⼰的名字%ConstrainedDTLZ1%如果你想在GUI界⾯去调⽤你的问题，这⾥<>必须起⼀个你问题⾃⼰的名字methods%%Initializationfunctionobj=C1_DTLZ1()ifisempty(obj.Global.M)obj.Global.M=3;%定义你优化问题的⽬标个数endifisempty(obj.Global.D)obj.Global.D=obj.Global.M+4;%定义你优化问题的变量个数，在DTLZ1中为⽬标个数加4endobj.Global.lower=zeros(1,obj.Global.D);定义你的优化问题的下界obj.Global.upper=ones(1,obj.Global.D);定义你优化问题的上界obj.Global.encoding='real';%定义你优化问题的编码⽅式，例如实数编码码'real',platEMO提供三种编码⽅式：'binary'⼆进制编码；'permutation'排列编码以及其他实数编码，根据你的优化问题的实际搜索空间可以⾃⾏选择end%%Calculateobjectivevalues计算⽬标函数值functionPopObj=CalObj(obj,PopDec)[N,D]=size(PopDec);%⾸先我们需要得到种群⼤⼩N和⾃变量个数DM=obj.Global.M;%得到⽬标个数g=100*(D-M+1+sum((PopDec(:,M:end)-0.5).^2-cos(20.*pi.*(PopDec(:,M:end)-0.5)),2));PopObj=0.5*repmat(1+g,1,M).*fliplr(cumprod([ones(N,1),PopDec(:,1:M-1)],2)).*[ones(N,1),1-PopDec(:,M-1:-1:1)];%给出⾃⼰问题的表达式，如果是⿊盒优化问题，这⾥需要将⿊盒问题抽象成包含输⼊和输出的接⼝形式，注意输⼊的矩阵PopDec每⼀⾏代表种群中的⼀条染⾊体，⿊盒问题需要可以接受这样的输⼊形式，如果⿊盒优化的函数只能⼀个个染⾊体计算，那么需要⽤for循环算出每⼀个染⾊体在不同⽬标上对应的⽬标函数值，再合成矩阵PopObj。end%%Calculateconstraintviolations计算⽬标函数的约束值，没有约束的话，可以不写这个接⼝函数functionPopCon=CalCon(obj,PopDec)PopObj=obj.CalObj(PopDec);PopCon=PopObj(:,end)/0.6+sum(PopObj(:,1:end-1)/0.5,2)-1;%值得注意的是，约束必须要写成g(x)>=0的形式哦，这⾥计算的g(x)的值end%%SamplereferencepointsonParetofront%这⾥是优化问题的最优PF，在实际问题中我们是得不到最优的PF的，因此实际问题时直接省略就可以的，这⾥C1_DTLZ1函数是学者们构造出来⽤来测试算法性能的，因此我们可以得到PF的理论解。functionP=PF(obj,N)P=UniformPoint(N,obj.Global.M)/2;endendend⾸先由于C1_DTLZ1是PROBLEM的⼀个⼦类，每个操作都是⼀个重载⽅法。以上每条代码都给出来详细的注释，参考这些就可以加⼊⾃⼰的⽅法啦。⼀下有两点值得注意的事项：1.决策变量在三种情况下可能是不可⾏的。⾸先，对于连续MOPs，它可能⼤于全局上界（GLOBAL.upper）或⼩于下界（GLOBAL.lower）。在这种情况下，它将被INDIVIDUAL的类设置为边界值，⽽MOP类不需要处理这种情况。其次，它可能没有满⾜约束(⼀个正的约束违背表明这个约束没有被满⾜)，在这种情况下，约束违背应该通过重载⽅法PROBLEM.CalCon()来计算。最后，对于组合MOPs，它可能是⼀个⾮法字符，在这种情况下，它应该通过重载⽅法PROBLEM.CalDec()来修复，例如在MOKP.m中：2.functionPopDec=CalDec(obj,PopDec)3.C=sum(obj.W,2)/2;4.[~,rank]=sort(max(obj.P./obj.W));5.fori=1:size(PopDec,1)6.whileany(obj.W*PopDec(i,:)’>C)7.k=find(PopDec(i,rank),1);8.PopDec(i,rank(k))=0;9.end10.end11.end2.问题类还可以接受⽤户的参数设置，⽐如我们⽬标函数中有⼀个参数可能需要经常改变，我们可以把他定义为该问题类的参数，例如：WFG1中的参数KclassdefWFG1<PROBLEM%<problem><WFG>%BenchmarkMOPproposedbyWalkingFishGroup%K------Thepositionparameter,whichshouldbeamultipleofM-1%如果你想在GUI界⾯显⽰，这⾥必须写上properties(Access=private)K;%Positionparameterendmethods%%Initializationfunctionobj=WFG1()ifisempty(obj.Global.M)obj.Global.M=3;endifisempty(obj.Global.D)obj.Global.D=obj.Global.M+9;endobj.K=obj.Global.ParameterSet(obj.Global.M-1);obj.Global.lower=zeros(1,obj.Global.D);obj.Global.upper=2:2:2*obj.Global.D;obj.Global.encoding='real';end%%CalculateobjectivevaluesfunctionPopObj=CalObj(obj,PopDec)[N,D]=size(PopDec);M=obj.Global.M;K=obj.K;L=D-K;D=1;S=2:2:2*M;A=ones(1,M-1);z01=PopDec./repmat(2:2:size(PopDec,2)*2,N,1);t1=zeros(N,K+L);t1(:,1:K)=z01(:,1:K);t1(:,K+1:end)=s_linear(z01(:,K+1:end),0.35);t2=zeros(N,K+L);t2(:,1:K)=t1(:,1:K);t2(:,K+1:end)=b_flat(t1(:,K+1:end),0.8,0.75,0.85);t3=zeros(N,K+L);t3=zeros(N,K+L);t3=b_poly(t2,0.02);t4=zeros(N,M);fori=1:M-1t4(:,i)=r_sum(t3(:,(i-1)*K/(M-1)+1:i*K/(M-1)),2*((i-1)*K/(M-1)+1):2:2*i*K/(M-1));endt4(:,M)=r_sum(t3(:,K+1:K+L),2*(K+1):2:2*(K+L));x=zeros(N,M);fori=1:M-1x(:,i)=max(t4(:,M),A(i)).*(t4(:,i)-0.5)+0.5;endx(:,M)=t4(:,M);h=convex(x);h(:,M)=mixed(x);PopObj=repmat(D*x(:,M),1,M)+repmat(S,N,1).*h;end%%SamplereferencepointsonParetofrontfunctionP=PF(obj,N)M=obj.Global.M;P=UniformPoint(N,M);c=ones(size(P,1),M);fori=1:size(P,1)forj=2:Mtemp=P(i,j)/P(i,1)*prod(1-c(i,M-j+2:M-1));c(i,M-j+1)=(temp^2-temp+sqrt(2*temp))/(temp^2+1);endendx=acos(c)*2/pi;temp=(1-sin(pi/2*x(:,2))).*P(:,M)./P(:,M-1);a=0:0.0001:1;E=abs(temp*(1-cos(pi/2*a))-1+repmat(a+cos(10*pi*a+pi/2)/10/pi,size(x,1),1));[~,rank]=sort(E,2);fori=1:size(x,1)x(i,1)=a(min(rank(i,1:10)));endP=convex(x);P(:,M)=mixed(x);P=repmat(2:2:2*M,size(P,1),1).*P;endendendfunctionOutput=s_linear(y,A)Output=abs(y-A)./abs(floor(A-y)+A);endfunctionOutput=b_flat(y,A,B,C)Output=A+min(0,floor(y-B))*A.*(B-y)/B-min(0,floor(C-y))*(1-A).*(y-C)/(1-C);Output=roundn(Output,-6);endfunctionOutput=b_poly(y,a)Output=y.^a;endfunctionOutput=r_sum(y,w)Output=sum(y.*repmat(w,size(y,1),1),2)./sum(w);endfunctionOutput=convex(x)Output=fliplr(cumprod([ones(size(x,1),1),1-cos(x(:,1:end-1)*pi/2)],2)).*[ones(size(x,1),1),1-sin(x(:,end-1:-1:1)*pi/2)];endendfunctionOutput=mixed(x)Output=1-x(:,1)-cos(10*pi*x(:,1)+pi/2)/10/pi;end这样定义之后，我们就可以直接在外部通过下⾯命令接受⽤户设定的参数了：1.obj.K=obj.Global.ParameterSet(obj.Global.M-1);增加算法接下来，我们主要讲如何增加⾃⼰设计的多⽬标进化算法(MOEA)了。MOEA函数由PlatEMO中的.m⽂件表⽰，该⽂件应该放在Algorithm⽂件夹中。即我们⾸先需要在Algorithm这个⽂件夹新建⼀个我们算法的⽂件夹。我们知道⼀个进化算法最重要的就是五个部分，初始化，函数评价，交配选择，交配产解，环境选择，接下来我们以NSGAII为例,说明如何在PlatEMO中增加⼀个多⽬标进化算法。⾸先我们给出NSGAII的主流程代码NSGAII.m：functionNSGAII(Global)%<algorithm><N>%NondominatedsortinggeneticalgorithmII%如果你想在GUI界⾯去调⽤你的算法，这⾥<>必须起⼀个你算法的名字%%GeneraterandompopulationPopulation=Global.Initialization();%这⾥是对种群进⾏初始化，Global类中已经完成了，不需要我们⾃⼰编写，如果你想设计特殊的初始化⽅式，可以在Population的基础上进⾏⾃⼰的操作即可。值得注意的是，在platEMO中，个体被初始化或者利⽤遗传算⼦产⽣之后（INDIVIDUAL⽬标被实例化）就会被⾃动评价，这是不需要我们⼈为操作的。[~,FrontNo,CrowdDis]=EnvironmentalSelection(Population,Global.N);%这⾥是对种群中的个体进⾏选择，详情见后⾯的分析。%%OptimizationwhileGlobal.NotTermination(Population)%这⾥是对算法停⽌条件的判断，⾥⾯包含了函数评估次数的检查，以及⼀些数据保存和画图等操作。如果需要修改算法的停⽌条件，需要⾃⼰修改Global，个⼈建议不要修改，因为如果对架构不了解有可能改错。函数评估次数是⽬前多⽬标优化算法中最常⽤的停⽌条件了，毕竟MOEA的收敛性并不好最为停⽌条件去度量。MatingPool=TournamentSelection(2,Global.N,FrontNo,-CrowdDis);%这⾥是进⾏交配选择，即给出交配的⽗代，platEMO中给出了常⽤的选择⽅式,⽐如⼆进制锦标赛和轮盘赌选择，他们在Public这个⽂件夹中，我们可以直接调⽤，需要⾃⼰设计选择⽅式的话，按照这些⽂件修改即可。Offspring=GA(Population(MatingPool));%这⾥是进⾏⼦代产⽣，paltEMO提供了很多⼦代产⽣⽅式，包含：DE,FEP,GA,PSO等，我们可以直接调⽤，这些函数位于Operater中。对于函数GA()，如果输⼊是INDIVIDUAL对象的数组，那么输出也是INDIVIDUAL对象的数组。如果输⼊是决策变量矩阵，那么输出也是决策变量矩阵。[Population,FrontNo,CrowdDis]=EnvironmentalSelection([Population,Offspring],Global.N);%这⾥是进⾏⾃然选择，在MOEA中，如何选择解以兼顾收敛性和多样性是学者们研究的重点之⼀，这⾥给出的是NSGAII中最经典的⽅法.endendEnvironmentalSelection.mfunction[Population,FrontNo,CrowdDis]=EnvironmentalSelection(Population,N)%TheenvironmentalselectionofNSGA-II%%Non-dominatedsorting[FrontNo,MaxFNo]=NDSort(Population.objs,Population.cons,N);%这⾥进⾏⾮⽀配排序，这是MOEA中最常⽤的⼀种⽅法Next=FrontNo<MaxFNo;%%Calculatethecrowdingdis