NSGA-II matlab 源码

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % MATLAB Code for % % % % Non-dominated Sorting Genetic Algorithm II (NSGA-II) % % Version 1.0 - April 2010 % % % % Programmed By: S. Mostapha Kalami Heris % % % % e-Mail: sm.kalami@gmail.com % % kalami@ee.kntu.ac.ir % % % % Homepage: http://www.kalami.ir % % % % nsga2.m : main file of the program % % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clc; clear; close all; % for m=1:20 CostFunction=@MyCost5; TestProblem=5; npop=100; % main population size nvar=10; % number of unknown variables nobj=numel(CostFunction(zeros(1,nvar))); % number of objective functions maxit=20; % maximum number of iterations 最大迭代数 VarMin=0; % lower bound of unknown variables VarMax=1; % upper bound of unknown variables VarRange=[VarMin VarMax]; % range of unknown variables pc=0.8; % crossover ratio 交叉率 nc=round(pc*npop/2)*2; % number of parents (also offsprings) pm=0.3; % mutation ratio 变异率 nm=round(pm*npop); % number of mutants 变异个体数 % initialization pop=CreateEmptyIndividuals(npop);% 初始化种群 %disp(['number of population size is ' num2str(numel(pop))]); for i=1:numel(pop) pop(i).Position=unifrnd(VarMin,VarMax,[1 nvar]); pop(i).Cost=CostFunction(pop(i).Position); end % main loop for it=1:maxit %非支配排序 选出所有非支配解 [pop F]=NonDominatedSorting(pop); pop=CalcCrowdingDistance(pop,F); disp(['Iteration ' num2str(it) ': Number of 1st Front Individuals = ' num2str(numel(F{1}))]); if it==maxit break; end pop2=CreateEmptyIndividuals(nc); for k=1:nc/2 p1=BinaryTournamentSelection(pop); p2=BinaryTournamentSelection(pop); ch=Crossover(p1,p2); ch(1).Cost=CostFunction(ch(1).Position); ch(2).Cost=CostFunction(ch(2).Position); pop2(2*k-1)=ch(1); pop2(2*k)=ch(2); end pop3=CreateEmptyIndividuals(nm); for k=1:nm p=BinaryTournamentSelection(pop); q=Mutate(p,VarRange); q.Cost=CostFunction(q.Position); pop3(k)=q; end pop=[pop pop2 pop3]; [pop F]=NonDominatedSorting(pop); pop=CalcCrowdingDistance(pop,F); pop=SortPopulation(pop); pop=pop(1:npop); end rep_costs=GetCosts(pop); Np=500; P1=0:1/Np:(Np-1)/Np; switch TestProblem case 1 P2=1-sqrt(P1); %ZDT1 and 4 case 2 P2=1-P1.^2; %for ZDT2 case 3 P2=1-sqrt(P1)-P1.*sin(10*pi*P1);%for ZDT3 Pd=P2(2:Np)-P2(1:Np-1); ind=find(Pd>0); P1(ind)=[];P2(ind)=[]; Ps=size(P2,2); Pd=P2(2:Ps)-P2(1:Ps-1); ind=find(Pd>0); pss=size(ind,2); P2(ind(pss)+1:Ps)=[]; P1(ind(pss)+1:Ps)=[]; for j=pss-1:-1:1 indJ=find(P2<P2(ind(j))); P2(ind(j)+1:indJ(1)-1)=5; end ind=find(P2==5); P1(ind)=[]; P2(ind)=[]; case 4 P2=1-sqrt(P1); %ZDT1 and 4 case 5 P1=1-exp(-4*P1).*sin(6*pi*P1).^6;%for ZDT6 P2=1-P1.^2;%for ZDT6 end P=[P1;P2]; % Calculate Closeness to P* (Generational Distance) of the Pareto front in % the final population.计算外部文档中的非支配解与函数真正的帕雷托解之间的平均距离 M=size(rep_costs,2); N=size(P,2); for i=1:M for j=1:N Dq1=0; for k=1:nobj Dq1=Dq1+(rep_costs(k,i)-P(k,j))^2; end Dq2(j)=sqrt(Dq1); end Dq(i)=min(Dq2); end GD=sqrt(sum(Dq.^2))/M; %计算当前文档中接的分布平均度 for i=1:size(rep_costs,2) d1=[]; d2=[]; for k=1:size(rep_costs,2) for j=1:nobj d1(j)=abs(rep_costs(j,i)-rep_costs(j,k)); end d2(k)=sum(d1); end s=find(d2==0); d2(s)=[]; d(i)=min(d2); end d_average=mean(d); eta=sqrt(sum((d_average-d).^2)/(M-1)); % Calculate diversity metrics of the Pareto front in the final population. [rep_costs(1,:),IF1]=sort(rep_costs(1,:)); rep_costs(2,:)=rep_costs(2,IF1); for k=1:nobj dm1(k)=(rep_costs(k,1)-P(k,1))^2; dm2(k)=(rep_costs(k,M)-P(k,N))^2; end DM1=sqrt(sum(dm1)); DM2=sqrt(sum(dm2)); DM=DM1+DM2; for i=1:M-1 for j=1:nobj di1(j)=(rep_costs(j,i)-rep_costs(j,i+1))^2; end di(i)=sqrt(sum(di1)); end Davg=mean(di); Ds=(DM+sum(abs(di-Davg)))/(DM+Davg*(M-1)); %画出非支配解的分布图（在适应度空间上）即帕雷托前沿的分布图 figure; hold on plot(P(1,:),P(2,:),'b-'); plot(rep_costs(1,:),rep_costs(2,:),'rx'); % end