差分进化算法解决多目标优化问题--内附matlab代码
一 多目标优化所谓的多目标优化是尝试同时最小化K个独立的目标函数。其目标是:
求
由此可知只有当有一个单一的向量同时取得K个目标函数的最小值,及可以说明其是多目标的解。故此解向量x满足条件:
往往这k个目标是相互冲突的,想要找到一个最优解满足k目标这是比较困难的,所以需要解决这一问题,排除彼此冲突的情况。最优解是一个折中的方案,取决于哪个目标更重要。如果每个目标函数能被指定一个权重来衡量其相对重要性,那么最优解就更加明显。
二 差分进化算法解决多目标问题
在了解多目标之后,现在就利用差分进化算法(DE)实战应用在其上。这里采用基准测试集(ZDT测试函数,来源于文献:Comparison of Multiobjective Evolutionary Algorithms: Empirical Results;SCH,FON,POL,KUR,来源于文献:A Fast and Elitist Multiobjective Genetic Algorithm:
NSGA-II;DTLZ,来源于文献:Scalable Multi-Objective Optimization Test Problems;MOP,来源于文献:Decomposition of a Multiobjective Optimization Problem into a Number of Simple Multiobjective Subproblems)
测试函数文件:
%%%%%%%%%%%%%%%%%%多目标测试函数%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function y=muti_fun(x,type)
global D f_num
if type==1%SCH----维度1---边界范围[-10^3,10^3]
y(1)=x.^2;
y(2)=(x-2).^2;
elseif type==2%FON-------维度3----边界范围[-4,4]
L=3;
sumFON1=0;
sumFON2=0;
for i=1:L
sumFON1=sumFON1+(x(i)-1/sqrt(3)).^2;
sumFON2=sumFON2+(x(i)+1/sqrt(3)).^2;
end
y(1)=1-exp(-sumFON1);
y(2)=1-exp(-sumFON2);
elseif type==3 %POL-----维度2-----边界范围[-pi,pi]
A1=0.5*sin(1)-2*cos(1)+sin(2)-1.5*cos(2);
A2=1.5*sin(1)-cos(1)+sin(2)-0.5*cos(2);
B1=0.5*sin(x(1))-2*cos(x(1))+sin(x(2))-1.5*cos(x(2));
B2=1.5*sin(x(1))-cos(x(1))+2*sin(x(2))-0.5*cos(x(2));
y(1)=1+(A1-B1).^2+(A2-B2).^2;
y(2)=(x(1)+3).^2+(x(2)+1).^2;
elseif type==4%KUR-----------维度3-----边界范围[-5,5]
sumKUR1=0;
sumKUR2=0;
for i=1:D-1
sumKUR1=sumKUR1+(-10*exp(-0.2*sqrt(x(i).^2+x(i+1).^2)));
sumKUR2=sumKUR2+abs(x(i)).^(0.8)+5*sin(x(i).^3);
end
y(1)=sumKUR1;
y(2)=sumKUR2;
elseif type==5%ZDT1----维度30----边界范围
y(1)=x(1);
y(2)=1-sqrt(y(1)/g(x));
elseif type==6%ZDT2----维度30----边界范围
y(1)=x(1);
y(2)=1-(y(1)/g1(x)).^2;
elseif type==7%ZDT3----维度30----边界范围
y(1)=x(1);
y(2)=1-sqrt(y(1)/g1(x))-(y(1)/g1(x))*sin(10*pi*y(1));
elseif type==8%ZDT4----维度10----边界范围x1=,x(2,...,10)=[-5,5]------PF distribution of objection space have problems
y(1)=x(1);
% y(2)=g2(x)*(1-(y(1)/g2(x)).^2);
y(2)=1-sqrt(y(1)/g2(x));
elseif type==9%ZDT6----维度10----边界范围
y(1)=1-exp(-4*x(1))*sin(6*pi*x(1)).^6;
y(2)=1-(y(1)/g3(x)).^2;
elseif type==10%MOP1----维度10----边界范围
y(1)=(1+g4(x))*x(1);
y(2)=(1+g4(x))*(1-sqrt(x(1)));
elseif type==11 %MOP2----维度10----边界范围
y(1)=(1+g5(x))*x(1);
y(2)=(1+g5(x))*(1-x(1).^2);
elseif type==12 %MOP3----维度10----边界范围
y(1)=(1+g6(x))*cos((pi*x(1))/2);
y(2)=(1+g6(x))*sin((pi*x(1))/2);
elseif type==13%MOP4----维度10----边界范围
y(1)=(1+g7(x))*x(1);
y(2)=(1+g7(x))*(1-x(1).^2*cos(2*pi*x(1)).^2);
elseif type==14 %MOP5----维度10----边界范围
y(1)=(1+g8(x))*x(1);
y(2)=(1+g8(x))*(1-sqrt(x(1)));
elseif type==15 %MOP6----三目标----维度10----边界范围
y(1)=(1+g9(x))*x(1)*x(2);
y(2)=(1+g9(x))*x(1)*(1-x(2));
y(3)=(1+g9(x))*(1-x(1));
elseif type==16%MOP7----三目标----维度10----边界范围
y(1)=(1+g10(x))*cos((x(1)*pi)/2)*cos((x(2)*pi)/2);
y(2)=(1+g10(x))*cos((x(1)*pi)/2)*sin((x(2)*pi)/2);
y(3)=(1+g10(x))*sin((x(1)*pi)/2);
elseif type==17 %DTLZ1----三目标----维度2----边界范围
for i=1:f_num
y(1)=(1/2)*x(1)*x(2)*(1+g11(x));
y(2)=(1/2)*x(1)*(1-x(2))*(1+g11(x));
y(3)=(1/2)*(1-x(1))*(1+g11(x));
end
% y(1)=1/2*x(1)*(1+g11(x));
% y(2)=1/2*(1-x(1))*(1+g11(x));
elseif type==18 %DTLZ2----三目标----维度2----边界范围
y(1)=(1+g12(x))*cos((pi/2)*x(1))*cos((pi/2)*x(2));
y(2)=(1+g12(x))*sin((pi/2)*x(1))*sin((pi/2)*x(2));
y(3)=(1+g12(x))*sin((pi/2)*x(1));
elseif type==19 %DTLZ3----
y(1)=(1+g11(x))*cos((pi/2)*x(1));
y(2)=(1+g11(x))*sin((pi/2)*x(1));
elseif type==20 %DTLZ4----
a=100;
y(1)=(1+g12(x))*cos((pi/2)*x(1).^a);
y(2)=(1+g12(x))*sin((pi/2)*x(1).^a);
elseif type==21 %DTLZ7----
y(1)=x(1);
y(2)=(1+g13(x))*h(x);
end
end
% function f=fun1(x,a,b)
% if x>a && x<=b
% f=exp(-2*((x-0.1)/0.8).^2*log(2))*abs(sin(5*pi*x)).^0.5;
% else
% f=exp(-2*((x-0.1)/0.8).^2*log(2))*abs(sin(5*pi*x)).^6;
% end
% end
%
% function f=fun2(x,e)
% f=sin(5*pi*(x.^e-0.05)).^6;
% end
function f=g(x)
global D
f=1+9*(sum(x(2:D))/(D-1));
end
function f=g1(x)
global D
f=1+9*(sum(x(2:D))/(D-1));
end
function f=g2(x)
global D
D=10;
f=0;
for i=2:D
f=f+x(i).^2-10*cos(4*pi*x(i));
end
f=f+1+10*(D-1);
end
function f=g3(x)
global D
f=1+9*(sum(x(:,2:D))/(D-1)).^0.25;
end
function f=g4(x) %MOP1
global D
f=0;
for i=2:D
t(i)=x(i)-sin(0.5*pi*x(1));
f=-0.9*t(i).^2+abs(t(i)).^0.6+f;
end
f=f*2*sin(pi*x(1));
end
function f=g5(x) %MOP2
global D
f=0;
for i=2:D
t(i)=x(i)-sin(0.5*pi*x(1));
f=(abs(t(i))/(1+exp(5*abs(t(i)))))+f;
end
f=f*10*sin(pi*x(1));
end
function f=g6(x) %MOP3
global D
f=0;
for i=2:D
t(i)=x(i)-sin(0.5*pi*x(1));
f=(abs(t(i))/(1+exp(5*abs(t(i)))))+f;
end
f=f*10*sin(pi*x(1)*(1/2));
end
function f=g7(x) %MOP4
global D
f=0;
for i=2:D
t(i)=x(i)-sin(0.5*pi*x(1));
f=(abs(t(i))/(1+exp(5*abs(t(i)))))+f;
end
f=f*10*sin(pi*x(1));
end
function f=g8(x) %MOP5
global D
f=0;
for i=2:D
t(i)=x(i)-sin(0.5*pi*x(1));
f=-0.9*t(i).^2+abs(t(i)).^0.6+f;
end
f=f*2*abs(cos(pi*x(1)));
end
function f=g9(x) %MOP6
global D
f=0;
for i=2:D
t(i)=x(i)-x(1)*x(2);
f=-0.9*t(i).^2+abs(t(i)).^0.6+f;
end
f=f*2*sin(pi*x(1));
end
function f=g10(x) %MOP7
global D
f=0;
for i=3:D
t(i)=x(i)-x(1)*x(2);
f=-0.9*t(i).^2+abs(t(i)).^0.6+f;
end
f=f*2*sin(pi*x(1));
end
function f=g11(x)%DTLZ1 DTLZ3
global D
f=0;
for i=1:D
f=(x(i)-0.5).^2-cos(20*pi*(x(i)-0.5))+f;
end
f=(f+sqrt(sum(x.^2)))*100;
end
function f=g12(x) %DTLZ2 DTLZ4
global D
f=0;
for i=1:D
f=f+(x(i)-0.5).^2;
end
end
function f=g13(x) %DTLZ7
f=1+(9/sqrt(sum(x.^2)))*sum(x);
end
function f=h(x)
global f_num
f=f_num-(x(1)/(1+g13(x)))*(1+sin(3*pi*x(1)));
end测试主m文件
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%DE解决多目标优化问题-----传统的差分进化算法
%时间:2021.10.27
%开发者:gouhuipeng
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clc
clear
close all
Np=100;%种群数
D=30; %个体数(问题的维度)
F=0.5; %突变因子(变异因子)
Cr=0.3;%交叉因子
%xmax(1)=1; %问题上界
%xmin(1)=0; %问题下届
for i=1:D
xmin(:,i)=0;
xmax(:,i)=1;
end
type=5; %测试问题的类型
K=500; %迭代次数
for j=1:D
x(:,j)=xmin(:,j)+(xmax(:,j)-xmin(:,j))*rand(Np,1);%产生初始值------第一代----对于这里的初始值可不可以进行给定一个最优初始值
end
best=x(1,:);%全局最优个体 ---之后不断更新--令其等于我们的初始值的第一个值
%end
%plot(x(:,1),x(:,2),&#39;.&#39;);
for i=2:Np
if(muti_fun(x(i,:),type)<muti_fun(best(1,:),type)) %编写函数有问题,下去思考
best=x(i,:); %选择最优个体
end
end
muti_y=muti_fun(best,type); %保证个体全局最优
%%进入循环直到满足精度要求或者迭代次数达到
for n=1:K %做迭代
time(n)=n;
%F=0.9-(n*(0.9-0.5)./K); %采用线性递减的方式做
%第二步 变异
for i=1:Np
for j=1:D %检查是否越界
r1=1;r2=1;r3=1;%使得个体满足变异条件.此处与Java有点不一样,他是从1开始
while(r1==r2||r1==r3||r2==r3||r1==i||r2==i||r3==i)
r1=ceil(Np*rand(1)); %保持其中的r1,r2,r3互异,这样做的目的是为了防止种群的单一性
r2=ceil(Np*rand(1));
r3=ceil(Np*rand(1));
end
v(i,j)=x(r1,j)+F*(x(r2,j)-x(r3,j));
%做一个防止越界
if v(i,j)<xmin(:,j)
v(i,j)=xmin(:,j);
elseif x(i,j) >xmax(:,j)
x(i,j)=xmax(:,j);
end
end
%交叉
for j=1:D
temper=rand;%随机产生一个数,用于进行多点交叉
if(temper<Cr)
u(i,j)=v(i,j);
else
u(i,j)=x(i,j);
end
end
mutifun_Uresult(i,:)=muti_fun(u(i,:),type);
mutifun_Xresult(i,:)=muti_fun(x(i,:),type);
%选择----------修改选择方法
% if(muti_fun(u(i,:),type)<=muti_fun(x(i,:),type))
% x(i,:)=u(i,:);
% end
% plot(i,x,&#39;bo&#39;);
ifmutifun_Uresult(i,:)<= mutifun_Xresult(i,:)
x(i,:)=u(i,:);
end
fun_muti_x(i,:)=muti_fun(x(i,:),type);
% plot( fun_muti_x(i,1), fun_muti_x(i,2),&#39;b--&#39;);
if fun_muti_x(i,:)< muti_y
%fi=fun_DE1(x(i,:));
best=x(i,:);
end
Best_f(n,:)=muti_fun(best,type);
end
end
%ZDT1理论前沿面
x1=0:0.01:1;
f2=1-sqrt(x1);
plot(x1,f2,&#39;k-&#39;);title(&#39;ZDT1&#39;)
%ZDT2理论前沿面
% x1=0:0.01:1;
% f2=1-x1.^2;
% plot(x1,f2,&#39;k-&#39;);
%ZDT3理论前沿面
%ZDT4理论前沿面
%ZDT6理论前沿面
hold on
plot( fun_muti_x(:,1), fun_muti_x(:,2),&#39;k+&#39;);
xlabel(&#39;f-1&#39;);ylabel(&#39;f-2&#39;);zlabel(&#39;f-3&#39;);
% legend(&#39;optimal PF&#39;,&#39;DE&#39;)
%optvalue=Best_f;
%fprintf(&#39;最优解结果为%f,%f,%f,%f,%f,%f,%f,%f,%f,%f&#39;,best);
disp([&#39;最小函数值为:&#39;, num2str(Best_f(n,:))]);三 输出结果
对于想学习matlab和差分进化算法的,推荐如下几本书,这里建议大家做一下机器学习和进化算法的结合。
欢迎关注知乎,留言咨询!
页:
[1]