本文详细介绍了遗传算法的概念、特点和基本操作，包括种群初始化、适应度函数、选择、交叉和变异等步骤。重点讲解了MATLAB中遗传算法工具箱的使用，如谢菲尔德遗传算法工具箱的安装配置和功能。此外，还探讨了遗传算法在非线性规划、TSP问题、多目标优化及车间调度等领域的应用实例。
1.遗传算法的介绍
遗传算法（Genetic Algorithm,GA）关键是把问题参数编码成染色体，再利用迭代方式进行选择、交叉、变异等运算交换种群中的染色体信息，根据适应度的好坏选择一定数量优秀的染色体，经过若干代的进化后，最终收敛于最好的符合优化目标染色体。
特点：
- 全局搜索能力强，局部搜索能力弱
- 一般得到的是次优解，而非最优解（问题规模小可得最优解，问题规模大一般只能得近似解）
术语：
- 个体(individual)
- 适应度函数(fitness function)、适应度函数值(fitness value)
- 种群(population)、种群大小(population size)
- 代(generation)
- 父代(parents)、子代(children)
- 选择(select)、交叉(crossover)、变异(mutation)
- 精英数目(eliteCount)、交叉后代比例(crossover fraction)
步骤：
- 初始化种群：位串编码、Grey编码、实数编码、多级参数编码、有序串编码、结构式编码
- 适应度函数：由目标函数变换得到，函数值越大/小，个体越优
- 选择操作：个体适应度函数值越高，被选中概率越大，轮盘赌法、锦标赛法
- 交叉操作：随机抽取种群两个个体，通过染色体的交换组合，把父串的优秀特征传给子串，从而产生新的优秀个体
- 变异操作：维持种群多样性，随机抽取种群单个体，选择染色体的一点变异
- 收敛寻优：通过迭代产生的新种群为下次计算的初始值，经过众多次迭代，算法收敛于“最优解”
- 结果统计：多次测试，根据需求，取最优或平均
参数：
- 种群规模：100
- 进化次数：30
- 交叉概率：0.6
- 变异概率：0.1
改进：
- 精英策略：子代种群中的最优个体永远不会比父代最优的个体差，使得父代好的个体不至于交叉、变异而丢失
- 进化逆转：TSP解的关键在于码串的顺序，交叉一方面可以维持种群的多样性，避免陷入局部最优，另一方面它不利于子代继承亲代的较多信息，特别是进化到了后期，群体充斥大量高适应度个体，交叉对亲代的较优基因破坏很大，使子代难以继承到亲代的优良基因，从而使交叉算子的搜索能力大大降低；单向进化逆转，即只接受朝着好的方向的逆转，搜索最优解的能力强于交叉算子
2.遗传算法工具箱
2.1 GA工具箱
- 英国谢菲尔德大学的遗传算法工具箱
- 美国北卡罗来纳州立大学的遗传算法最优化工具箱
- MATLAB自带的遗传算法工具箱与直接搜索工具箱
2.2 谢菲尔德遗传算法工具箱
2.2.1 安装配置到 Matlab
- 将 gatbx 文件夹复制粘贴到 Matlab 软件安装目录下的 toolbox 目录中
- 将所在文件夹添加到 Matlab 的搜索路径：str=[matlabroot,'\toolbox\gatbx']; addpath(str)
- 测试是否安装成功：v = ver('gatbx')
2.2.2 工具箱相关函数功能介绍
2.2
2.3 MATLAB自带的遗传算法工具箱与直接搜索工具箱
Genetic Algorithm and Direct Search Toolbox，GADST
2.3.1 GADST GUI
命令行输入：
optimtool('ga')
注意：GADST是对目标函数去最小值进行优化，若要最大值优化，将适应度函数*-1
3 非线性规划的函数寻优
1.初始化种群：实数编码：不必进行数值转换，可以直接在解的表现型上进行遗传算法操作
function ret=Code(lenchrom, bound)
% lenchrom input : 染色体长度
% bound input : 变量的取值范围
% ret output: 染色体的编码值
flag=0;
while flag==0
pick=rand(1,length(lenchrom));
ret=bound(:,1)'+(bound(:,2)-bound(:,1))'.*pick; % 线性插值
flag=test(lenchrom,bound,ret); % 检验染色体的可行性
end
2.适应度函数：求函数最小值，将函数值的倒数作为个体的适应度值，即函数值越小的个体越优
function y = fun(x)
y=-5*sin(x(1))*sin(x(2))*sin(x(3))*sin(x(4))*sin(x(5))-sin(5*x(1))*sin(5*x(2))*sin(5*x(3))*sin(5*x(4))*sin(5*x(5))+8;
3.选择操作：轮盘赌：基于适应度比例的选择策略
function ret=Select(individuals,sizepop)
% individuals input : 种群信息
% sizepop input : 种群规模
% opts input : 选择方法的选择
% ret output : 经过选择后的种群
individuals.fitness= 1./(individuals.fitness);
sumfitness=sum(individuals.fitness);
sumf=individuals.fitness./sumfitness;
index=[];
for i=1:sizepop % 转sizepop次轮盘
pick=rand;
while pick==0
pick=rand;
end
for j=1:sizepop
pick=pick-sumf(j);
if pick<0
index=[index j];
break; % 寻找落入的区间，此次转轮盘选中了染色体i，注意：在转sizepop次轮盘的过程中，有可能会重复选择某些染色体
end
end
end
individuals.chrom=individuals.chrom(index,:);
individuals.fitness=individuals.fitness(index);
ret=individuals;
4.交叉操作：实数交叉
function ret=Cross(pcross, lenchrom, chrom, sizepop, bound)
% pcorss input : 交叉概率
% lenchrom input : 染色体的长度
% chrom input : 染色体群
% sizepop input : 种群规模
% ret output : 交叉后的染色体
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));
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),:)); % 检验染色体1的可行性
flag2=test(lenchrom,bound,chrom(index(2),:)); % 检验染色体2的可行性
if flag1*flag2==0
flag=0;
else flag=1;
end % 如果两个染色体不是都可行，则重新交叉
end
end
ret=chrom;
5.变异操作
function ret=Mutation(pmutation, lenchrom, chrom, sizepop, pop, bound)
% pcorss input : 变异概率
% lenchrom input : 染色体长度
% chrom input : 染色体群
% sizepop input : 种群规模
% pop input : 当前种群的进化代数和最大的进化代数信息
% ret output : 变异后的染色体
for i=1:sizepop
% 随机选择一个染色体进行变异
pick=rand;
while pick==0
pick=rand;
end
index=ceil(pick*sizepop);
% 变异概率决定该轮循环是否进行变异
pick=rand;
if pick>pmutation
continue;
end
flag=0;
while flag==0
% 变异位置
pick=rand;
while pick==0
pick=rand;
end
% 随机选择了染色体变异的位置，即选择了第pos个变量进行变异
pos=ceil(pick*sum(lenchrom));
v=chrom(i,pos);
v1=v-bound(pos,1);
v2=bound(pos,2)-v;
pick=rand; % 变异开始
if pick>0.5
delta=v2*(1-pick^((1-pop(1)/pop(2))^2));
chrom(i,pos)=v+delta;
else
delta=v1*(1-pick^((1-pop(1)/pop(2))^2));
chrom(i,pos)=v-delta;
end % 变异结束
flag=test(lenchrom,bound,chrom(i,:)); % 检验染色体的可行性
end
end
ret=chrom;
6.收敛寻优
基本GA寻优
%% 清空环境
clc
clear
%% 遗传算法参数
maxgen=30; % 进化代数
sizepop=100; % 种群规模
pcross=[0.6]; % 交叉概率
pmutation=[0.01]; % 变异概率
lenchrom=[1 1 1 1 1]; % 变量字串长度
bound=[0 0.9*pi;0 0.9*pi;0 0.9*pi;0 0.9*pi;0 0.9*pi]; % 变量范围
%% 个体初始化
individuals=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); % 个体适应度
end
% 找最好的染色体
[bestfitness bestindex]=min(individuals.fitness);
bestchrom=individuals.chrom(bestindex,:); % 最好的染色体
avgfitness=sum(individuals.fitness)/sizepop; % 染色体的平均适应度
% 记录每一代进化中最好的适应度和平均适应度
trace=[];
%% 进化开始
for i=1:maxgen
% 选择操作
individuals=Select(individuals,sizepop);
avgfitness=sum(individuals.fitness)/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);
end
% 找到最小和最大适应度的染色体及它们在种群中的位置
[newbestfitness,newbestindex]=min(individuals.fitness);
[worestfitness,worestindex]=max(individuals.fitness);
% 代替上一次进化中最好的染色体
if bestfitness>newbestfitness
bestfitness=newbestfitness;
bestchrom=individuals.chrom(newbestindex,:);
end
individuals.chrom(worestindex,:)=bestchrom;
individuals.fitness(worestindex)=bestfitness;
avgfitness=sum(individuals.fitness)/sizepop;
% 记录每一代进化中最好的适应度和平均适应度
trace=[trace;avgfitness bestfitness];
end %进化结束
%% 结果显示
[r c]=size(trace);
figure
plot([1:r]',trace(:,1),'r-',[1:r]',trace(:,2),'b--');
title(['函数值曲线 ' '终止代数＝' num2str(maxgen)],'fontsize',12);
xlabel('进化代数','fontsize',12);ylabel('函数值','fontsize',12);
legend('各代平均值','各代最佳值');
disp('函数值 变量');
ylim([1.5 8])
%xlim([1,size(trace,1)])
grid on
% 窗口显示
disp([bestfitness x]);
工具箱fmincon函数优化
function ret = nonlinear(chrom,sizepop)
for i=1:sizepop
x=fmincon(inline('-5*sin(x(1))*sin(x(2))*sin(x(3))*sin(x(4))*sin(x(5))-sin(5*x(1))*sin(5*x(2))*sin(5*x(3))*sin(5*x(4))*sin(5*x(5))'),chrom(i,:)',[],[],[],[],[0 0 0 0 0],[2.8274 2.8274 2.8274 2.8274 2.8274]);
ret(i,:)=x';
end
%% 清空环境
clc
clear
warning off
%% 遗传算法参数
maxgen=30; % 进化代数
sizepop=100; % 种群规模
pcross=[0.6]; % 交叉概率
pmutation=[0.01]; % 变异概率
lenchrom=[1 1 1 1 1]; % 变量字串长度
bound=[0 0.9*pi;0 0.9*pi;0 0.9*pi;0 0.9*pi;0 0.9*pi]; %变量范围
%% 个体初始化
individuals=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); % 个体适应度
end
%找最好的染色体
[bestfitness bestindex]=min(individuals.fitness);
bestchrom=individuals.chrom(bestindex,:); % 最好的染色体
avgfitness=sum(individuals.fitness)/sizepop; % 染色体的平均适应度
% 记录每一代进化中最好的适应度和平均适应度
trace=[];
%% 进化开始
for i=1:maxgen
% 选择操作
individuals=Select(individuals,sizepop);
avgfitness=sum(individuals.fitness)/sizepop;
% 交叉操作
individuals.chrom=Cross(pcross,lenchrom,individuals.chrom,sizepop,bound);
% 变异操作
individuals.chrom=Mutation(pmutation,lenchrom,individuals.chrom,sizepop,[i maxgen],bound);
if mod(i,10)==0
individuals.chrom=nonlinear(individuals.chrom,sizepop);
end
% 计算适应度
for j=1:sizepop
x=individuals.chrom(j,:);
individuals.fitness(j)=fun(x);
end
% 找到最小和最大适应度的染色体及它们在种群中的位置
[newbestfitness,newbestindex]=min(individuals.fitness);
[worestfitness,worestindex]=max(individuals.fitness);
% 代替上一次进化中最好的染色体
if bestfitness>newbestfitness
bestfitness=newbestfitness;
bestchrom=individuals.chrom(newbestindex,:);
end
individuals.chrom(worestindex,:)=bestchrom;
individuals.fitness(worestindex)=bestfitness;
avgfitness=sum(individuals.fitness)/sizepop;
% 记录每一代进化中最好的适应度和平均适应度
trace=[trace;avgfitness bestfitness];
end % 进化结束
%% 结果显示
figure
[r c]=size(trace);
plot([1:r]',trace(:,1),'r-',[1:r]',trace(:,2),'b--');
title(['函数值曲线 ' '终止代数＝' num2str(maxgen)],'fontsize',12);
xlabel('进化代数','fontsize',12);ylabel('函数值','fontsize',12);
legend('各代平均值','各代最佳值','fontsize',12);
ylim([1.5 8])
disp('函数值 变量');
% 窗口显示
disp([bestfitness x]);
grid on
7.结果
基本GA优化迭代到第30次函数值收敛于2，工具箱优化迭代到第10次函数值收敛于2
4 基于遗传算法的TSP算法
TSP：从一个城市出发，访问每个城市且只去一次，最后回到出发城市，使路径最短（一条遍历n个城市的最短曼哈顿回路）
| 城市编号 | X坐标 | Y坐标 | 城市编号 | X坐标 | Y坐标 |
| 1 | 16.47 | 96.1 | 8 | 17.2 | 96.29 |
| 2 | 16.47 | 94.44 | 9 | 16.3 | 97.38 |
| 3 | 20.09 | 92.54 | 10 | 14.05 | 98.12 |
| 4 | 22.39 | 93.37 | 11 | 16.53 | 97.38 |
| 5 | 25.23 | 97.24 | 12 | 21.52 | 95.59 |
| 6 | 22 | 96.05 | 13 | 19.41 | 97.13 |
| 7 | 20.47 | 97.02 | 14 | 20.09 | 92.55 |
1.编码
对14个城市的TSP问题：|1|14|2|4|8|5|6|12|9|7|11|13|3|10| 可作为一条合法染色体
2.初始化种群：依据城市规模而定 [50, 200]
function Chrom=InitPop(NIND, N)
% 输入：
% NIND：种群大小
% N： 个体染色体长度（这里为城市的个数）
Chrom=zeros(NIND,N); % 用于存储种群
for i=1:NIND
Chrom(i,:)=randperm(N); % 随机生成初始种群
end
3.适应度函数：适应度函数值越大，个体值越优
function FitnV=Fitness(len)
% 输入：
% 个体的长度（TSP的距离）
% 输出：
% 个体的适应度值
FitnV=1./len;
4.选择
function SelCh=Select(Chrom,FitnV,GGAP)
% 输入
% Chrom 种群
% FitnV 适应度值
% GGAP：选择概率
% 输出
% SelCh 被选择的个体
NIND=size(Chrom,1);
NSel=max(floor(NIND*GGAP+.5),2);
ChrIx=Sus(FitnV,NSel);
SelCh=Chrom(ChrIx,:);
function NewChrIx = Sus(FitnV,Nsel)
% 输入:
% FitnV 个体的适应度值
% Nsel 被选择个体的数目
% 输出:
% NewChrIx 被选择个体的索引号
% Identify the population size (Nind)
[Nind,ans] = size(FitnV);
% Perform stochastic universal sampling
cumfit = cumsum(FitnV);
trials = cumfit(Nind) / Nsel * (rand + (0:Nsel-1)');
Mf = cumfit(:, ones(1, Nsel));
Mt = trials(:, ones(1, Nind))';
[NewChrIx, ans] = find(Mt < Mf & [ zeros(1, Nsel); Mf(1:Nind-1, :) ] <= Mt);
% Shuffle new population
[ans, shuf] = sort(rand(Nsel, 1));
NewChrIx = NewChrIx(shuf);
% End of function
5.交叉：部分映射杂交，确定交叉操作的父代并将其两两分组，每组重复
以10个城市为例演示交叉操作
function SelCh=Recombin(SelCh,Pc)
% 输入
% SelCh 被选择的个体
% Pc 交叉概率
% 输出：
% SelCh 交叉后的个体
NSel=size(SelCh,1);
for i=1:2:NSel-mod(NSel,2)
if Pc>=rand % 交叉概率Pc
[SelCh(i,:),SelCh(i+1,:)]=intercross(SelCh(i,:),SelCh(i+1,:));
end
end
% 输入：
% a和b为两个待交叉的个体
% 输出：
% a和b为交叉后得到的两个个体
function [a,b]=intercross(a,b)
L=length(a);
r1=randsrc(1,1,[1:L]);
r2=randsrc(1,1,[1:L]);
if r1~=r2
a0=a;b0=b;
s=min([r1,r2]);
e=max([r1,r2]);
for i=s:e
a1=a;b1=b;
a(i)=b0(i);
b(i)=a0(i);
x=find(a==a(i));
y=find(b==b(i));
i1=x(x~=i);
i2=y(y~=i);
if ~isempty(i1)
a(i1)=a1(i);
end
if ~isempty(i2)
b(i2)=b1(i);
end
end
end
% %交叉算法采用部分匹配交叉%交叉算法采用部分匹配交叉
% function [a,b]=intercross(a,b)
% L=length(a);
% r1=ceil(rand*L);
% r2=ceil(rand*L);
% r1=4;r2=7;
% if r1~=r2
% s=min([r1,r2]);
% e=max([r1,r2]);
% a1=a;b1=b;
% a(s:e)=b1(s:e);
% b(s:e)=a1(s:e);
% for i=[setdiff(1:L,s:e)]
% [tf, loc] = ismember(a(i),a(s:e));
% if tf
% a(i)=a1(loc+s-1);
% end
% [tf, loc]=ismember(b(i),b(s:e));
% if tf
% b(i)=b1(loc+s-1);
% end
% end
% end
6.变异：随机选两个点交换位置
以10个城市为例演示交叉操作
function SelCh=Mutate(SelCh,Pm)
% 输入：
% SelCh 被选择的个体
% Pm 变异概率
% 输出：
% SelCh 变异后的个体
[NSel,L]=size(SelCh);
for i=1:NSel
if Pm>=rand
R=randperm(L);
SelCh(i,R(1:2))=SelCh(i,R(2:-1:1));
end
end
7.进化逆转：逆转算子的单方向性，即只有经过逆转后，适应度值有提高的才保留，否则逆转无效
以10个城市为例演示交叉操作
function SelCh=Reverse(SelCh,D)
% 输入
% SelCh 被选择的个体
% D 个城市的距离矩阵
% 输出
% SelCh 进化逆转后的个体
[row,col]=size(SelCh);
ObjV=PathLength(D,SelCh); % 计算路径长度
SelCh1=SelCh;
for i=1:row
r1=randsrc(1,1,[1:col]);
r2=randsrc(1,1,[1:col]);
mininverse=min([r1 r2]);
maxinverse=max([r1 r2]);
SelCh1(i,mininverse:maxinverse)=SelCh1(i,maxinverse:-1:mininverse);
end
ObjV1=PathLength(D,SelCh1); % 计算路径长度
index=ObjV1<ObjV;
SelCh(index,:)=SelCh1(index,:);
9.寻优
clear
clc
close all
X =[16.47,96.10
16.47,94.44
20.09,92.54
22.39,93.37
25.23,97.24
22.00,96.05
20.47,97.02
17.20,96.29
16.30,97.38
14.05,98.12
16.53,97.38
21.52,95.59
19.41,97.13
20.09,92.55]; % 个城市坐标位置,可以换成load CityPosition1.mat
NIND=100; % 种群大小
MAXGEN=200;
Pc=0.9; % 交叉概率
Pm=0.05; % 变异概率
GGAP=0.9; % 代沟(Generation gap)
D=Distanse(X); % 生成距离矩阵
N=size(D,1); % (34*34)
%% 初始化种群
Chrom=InitPop(NIND,N);
%% 在二维图上画出所有坐标点
% figure
% plot(X(:,1),X(:,2),'o');
%% 画出随机解的路线图
DrawPath(Chrom(1,:),X)
pause(0.0001)
%% 输出随机解的路线和总距离
disp('初始种群中的一个随机值:')
OutputPath(Chrom(1,:));
Rlength=PathLength(D,Chrom(1,:));
disp(['总距离：',num2str(Rlength)]);
disp('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~')
%% 优化
gen=0;
figure;
hold on;box on
xlim([0,MAXGEN])
title('优化过程')
xlabel('代数')
ylabel('最优值')
ObjV=PathLength(D,Chrom); %计算路线长度
preObjV=min(ObjV);
while gen<MAXGEN
%% 计算适应度
ObjV=PathLength(D,Chrom); %计算路线长度
% fprintf('%d %1.10f\n',gen,min(ObjV))
line([gen-1,gen],[preObjV,min(ObjV)]);pause(0.0001)
preObjV=min(ObjV);
FitnV=Fitness(ObjV);
%% 选择
SelCh=Select(Chrom,FitnV,GGAP);
%% 交叉操作
SelCh=Recombin(SelCh,Pc);
%% 变异
SelCh=Mutate(SelCh,Pm);
%% 逆转操作
SelCh=Reverse(SelCh,D);
%% 重插入子代的新种群
Chrom=Reins(Chrom,SelCh,ObjV);
%% 更新迭代次数
gen=gen+1 ;
end
%% 画出最优解的路线图
ObjV=PathLength(D,Chrom); %计算路线长度
[minObjV,minInd]=min(ObjV);
DrawPath(Chrom(minInd(1),:),X)
%% 输出最优解的路线和总距离
disp('最优解:')
p=OutputPath(Chrom(minInd(1),:));
disp(['总距离：',num2str(ObjV(minInd(1)))]);
disp('-------------------------------------------------------------')
% 计算两两城市之间的距离
function D=Distanse(a)
% 输入 a 各城市的位置坐标
% 输出 D 两两城市之间的距离
row=size(a,1);
D=zeros(row,row);
for i=1:row
for j=i+1:row
D(i,j)=((a(i,1)-a(j,1))^2+(a(i,2)-a(j,2))^2)^0.5;
D(j,i)=D(i,j);
end
end
% 输出路径函数
function p=OutputPath(R)
% 输入：R 路径
R=[R,R(1)];
N=length(R);
p=num2str(R(1));
for i=2:N
p=[p,'—>',num2str(R(i))];
end
disp(p)
% 计算各个体的路径长度
function len=PathLength(D,Chrom)
% 输入：
% D 两两城市之间的距离
% Chrom 个体的轨迹
[row,col]=size(D);
NIND=size(Chrom,1);
len=zeros(NIND,1);
for i=1:NIND
p=[Chrom(i,:) Chrom(i,1)];
i1=p(1:end-1);
i2=p(2:end);
len(i,1)=sum(D((i1-1)*col+i2));
end
%重插入子代的新种群
function Chrom=Reins(Chrom,SelCh,ObjV)
% 输入：
% Chrom 父代的种群
% SelCh 子代种群
% ObjV 父代适应度
% 输出
% Chrom 组合父代与子代后得到的新种群
NIND=size(Chrom,1);
NSel=size(SelCh,1);
[TobjV,index]=sort(ObjV);
Chrom=[Chrom(index(1:NIND-NSel),:);SelCh];
function varargout = dsxy2figxy(varargin)
if length(varargin{1}) == 1 && ishandle(varargin{1}) ...
&& strcmp(get(varargin{1},'type'),'axes')
hAx = varargin{1};
varargin = varargin(2:end);
else
hAx = gca;
end;
if length(varargin) == 1
pos = varargin{1};
else
[x,y] = deal(varargin{:});
end
axun = get(hAx,'Units');
set(hAx,'Units','normalized');
axpos = get(hAx,'Position');
axlim = axis(hAx);
axwidth = diff(axlim(1:2));
axheight = diff(axlim(3:4));
if exist('x','var')
varargout{1} = (x - axlim(1)) * axpos(3) / axwidth + axpos(1);
varargout{2} = (y - axlim(3)) * axpos(4) / axheight + axpos(2);
else
pos(1) = (pos(1) - axlim(1)) / axwidth * axpos(3) + axpos(1);
pos(2) = (pos(2) - axlim(3)) / axheight * axpos(4) + axpos(2);
pos(3) = pos(3) * axpos(3) / axwidth;
pos(4) = pos(4) * axpos(4 )/ axheight;
varargout{1} = pos;
end
set(hAx,'Units',axun)
function DrawPath(Chrom,X)
% 输入
% Chrom 待画路径
% X 各城市坐标位置
R=[Chrom(1,:) Chrom(1,1)]; % 一个随机解(个体)
figure;
hold on
plot(X(:,1),X(:,2),'o','color',[0.5,0.5,0.5])
plot(X(Chrom(1,1),1),X(Chrom(1,1),2),'rv','MarkerSize',20)
for i=1:size(X,1)
text(X(i,1)+0.05,X(i,2)+0.05,num2str(i),'color',[1,0,0]);
end
A=X(R,:);
row=size(A,1);
for i=2:row
[arrowx,arrowy] = dsxy2figxy(gca,A(i-1:i,1),A(i-1:i,2)); % 坐标转换
annotation('textarrow',arrowx,arrowy,'HeadWidth',8,'color',[0,0,1]);
end
hold off
xlabel('横坐标')
ylabel('纵坐标')
title('轨迹图')
box on
10.结果
初始种群中的一个随机值:
2—>8—>14—>1—>10—>12—>5—>6—>11—>3—>7—>13—>9—>4—>2
总距离：63.9542
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
最优解:
14—>3—>4—>5—>6—>12—>7—>13—>8—>11—>9—>10—>1—>2—>14
总距离：29.3405
----------------------------------------------------------------
5 多种群遗传算法的函数优化
未成熟收敛：
表现：群体中所有个体都趋于同一状态而停止进化，算法得不到令人满意的结果
问题：传统GA通常采用轮盘赌选择、单点交叉、位点变异
- 当种群中存在个别超常个体（适应度函数值远远高于其他个体）时，该个体在选择算子作用下，将会被多次选中，下一代群体很快被该个体控制，群体失去竞争性，导致群体停滞不前，陷入局部最优
- 交叉和变异的概率变化：全局搜索能力和局部搜索能力的均衡
- 群体规模：群体规模小，多样性差，个体间竞争弱，随着进化的进行，群体很快出现单一化，群体更新只能靠变异操作；群体规模大，计算量增加，计算效率不高
- 终止依据：迭代次数。次数少，进化不充分
多种群遗传算(Multiple Population Genetic Algorithm，MPGA)
- 引入多个种群同时优化搜索，不同的种群以不同的控制参数，实现不同的搜索目的
- 多种群间通过移民算子联系，实现多种群的协同进化
- 通过人工选择算子保存各种群每代中的最优个体，并作为算法收敛的依据
二元函数求最值
移民算子
| 变量名 | 类型 | 说明 | |
| 输入参数 | Chrom | Cell | 每个元胞单元为一个种群的编码 |
| ObjV | Cell | 每个元胞为一个种群所有个体的目标值 | |
| 输出参数 | Chrom | Cell | 每个元胞单元为一个种群的编码 |
| ObjV | Cell | 每个元胞为一个种群所有个体的目标值 |
%移民算子
function [Chrom,ObjV]=immigrant(Chrom,ObjV)
MP=length(Chrom);
for i=1:MP
[MaxO,maxI]=max(ObjV{i}); % 找出第i种群中最优的个体
next_i=i+1; % 目标种群（移民操作中）
if next_i>MP;next_i=mod(next_i,MP);end
[MinO,minI]=min(ObjV{next_i}); % 找出目标种群中最劣的个体
% 目标种群最劣个体替换为源种群最优个体
Chrom{next_i}(minI,:)=Chrom{i}(maxI,:);
ObjV{next_i}(minI)=ObjV{i}(maxI);
end
人工选择算子
| 变量名 | 类型 | 说明 | |
| 输入参数 | Chrom | Cell | 每个元胞单元为一个种群的编码 |
| ObjV | Cell | 每个元胞为一个种群所有个体的目标值 | |
| MaxObjV | Double | 各个种群当前最优个体的目标值 | |
| MaxChrom | Double | 各个种群当前最优个体的编码 | |
| 输出参数 | MaxObjV | Double | 各个种群当前最优个体的目标值 |
| MaxChrom | Double | 各个种群当前最优个体的编码 |
% 人工选择算子
function [MaxObjV,MaxChrom]=EliteInduvidual(Chrom,ObjV,MaxObjV,MaxChrom)
MP=length(Chrom); % 种群数
for i=1:MP
[MaxO,maxI]=max(ObjV{i}); % 找出第i种群中最优个体
if MaxO>MaxObjV(i)
MaxObjV(i)=MaxO; % 记录各种群的精华个体
MaxChrom(i,:)=Chrom{i}(maxI,:); % 记录各种群精华个体的编码
end
end
目标函数
function obj=ObjectFunction(X)
% X的每行为一个个体
col=size(X,1);
for i=1:col
obj(i,1)=21.5+X(i,1)*sin(4*pi*X(i,1))+X(i,2)*sin(20*pi*X(i,2));
end
主函数
多种群遗传：
clear;
clc
close all
NIND=40; % 个体数目
NVAR=2; % 变量的维数
PRECI=20; % 变量的二进制位数
GGAP=0.9; % 代沟
MP=10; % 种群数目
FieldD=[rep(PRECI,[1,NVAR]);[-3,4.1;12.1,5.8];rep([1;0;1;1],[1,NVAR])]; % 译码矩阵
for i=1:MP
Chrom{i}=crtbp(NIND, NVAR*PRECI); % 创建初始种群
end
pc=0.7+(0.9-0.7)*rand(MP,1); % 在【0.7,0.9】范围内随机产生交叉概率
pm=0.001+(0.05-0.001)*rand(MP,1); % 在【0.001,0.05】范围内随机产生变异概率
gen=0; % 初始遗传代数
gen0=0; % 初始保持代数
MAXGEN=10; % 最优个体最少保持代数
maxY=0; % 最优值
for i=1:MP
ObjV{i}=ObjectFunction(bs2rv(Chrom{i}, FieldD)); % 计算各初始种群个体的目标函数值
end
MaxObjV=zeros(MP,1); % 记录精华种群
MaxChrom=zeros(MP,PRECI*NVAR); % 记录精华种群的编码
while gen0<=MAXGEN
gen=gen+1; % 遗传代数加1
for i=1:MP
FitnV{i}=ranking(-ObjV{i}); % 各种群的适应度
SelCh{i}=select('sus', Chrom{i}, FitnV{i},GGAP); % 选择操作
SelCh{i}=recombin('xovsp',SelCh{i}, pc(i)); % 交叉操作
SelCh{i}=mut(SelCh{i},pm(i)); % 变异操作
ObjVSel=ObjectFunction(bs2rv(SelCh{i}, FieldD)); % 计算子代目标函数值
[Chrom{i},ObjV{i}]=reins(Chrom{i},SelCh{i},1,1,ObjV{i},ObjVSel); % 重插入操作
end
[Chrom,ObjV]=immigrant(Chrom,ObjV); % 移民操作
[MaxObjV,MaxChrom]=EliteInduvidual(Chrom,ObjV,MaxObjV,MaxChrom); % 人工选择精华种群
YY(gen)=max(MaxObjV); % 找出精华种群中最优的个体
if YY(gen)>maxY % 判断当前优化值是否与前一次优化值相同
maxY=YY(gen); % 更新最优值
gen0=0;
else
gen0=gen0+1; % 最优值保持次数加1
end
end
%% 进化过程图
plot(1:gen,YY)
xlabel('进化代数')
ylabel('最优解变化')
title('进化过程')
xlim([1,gen])
%% 输出最优解
[Y,I]=max(MaxObjV); %找出精华种群中最优的个体
X=(bs2rv(MaxChrom(I,:), FieldD)); %最优个体的解码解
disp(['最优值为：',num2str(Y)])
disp(['对应的自变量取值：',num2str(X)])
标准遗传：
clear
clc
pc=0.7;
pm=0.05;
% 定义遗传算法参数
NIND=40; % 个体数目
MAXGEN=500; % 最大遗传代数
NVAR=2; % 变量的维数
PRECI=20; % 变量的二进制位数
GGAP=0.9; % 代沟
trace=zeros(MAXGEN,1);
FieldD=[rep(PRECI,[1,NVAR]);[-3,4.1;12.1,5.8];rep([1;0;1;1],[1,NVAR])]; % 建立区域描述器
Chrom=crtbp(NIND, NVAR*PRECI); % 创建初始种群
gen=0; % 代计数器
maxY=0; % 最优值
ObjV=ObjectFunction(bs2rv(Chrom, FieldD)); % 计算初始种群个体的目标函数值
while gen<MAXGEN % 迭代
FitnV=ranking(-ObjV); % 分配适应度值(Assign fitness values)
SelCh=select('sus', Chrom, FitnV, GGAP); % 选择
SelCh=recombin('xovsp', SelCh, pc); % 重组
SelCh=mut(SelCh,pm); % 变异
ObjVSel=ObjectFunction(bs2rv(SelCh, FieldD)); % 计算子代目标函数值
[Chrom ObjV]=reins(Chrom, SelCh, 1, 1, ObjV, ObjVSel); % 重插入
gen=gen+1; % 代计数器增加
if maxY<max(ObjV)
maxY=max(ObjV);
end
trace(gen,1)=maxY;
end
%% 进化过程图
plot(1:gen,trace(:,1));
%% 输出最优解
[Y,I]=max(ObjV);
X=bs2rv(Chrom, FieldD);
disp(['最优值为：',num2str(Y)])
disp(['对应的自变量取值：',num2str(X(I,:))])
结果
6.基于遗传算法的多目标优化
Pareto最优解(Pareto optima)：某一目标函数的提高需要以另一个目标函数的降低为代价
带精英策略的快速非支配排序遗传算法(Nondominated Sorting Genetic Algorithm II, NSGA-II)
概念：
- 支配(dominate)、非劣(non-inferior)：若个体p至少有一个目标必个体q好，且个体p的所有目标都不必个体q的差，则个体p支配个体q，个体p非劣与个体q
- 序值(rank)、前端(front)：如果个体p支配个体q，则个体p的序值比个体q的序值低；如果个体p和q互不支配，则p和q拥有相同序值。序值为1的个体属于第一前端，序值为2的个体属于第二前端，以此类推。显然，在当前种群中，第一前端是完全不受支配的，第二前端受第一前端中个体的支配。通过排序，可以将种群中的个体分到不同的前端
- 拥挤距离(Crowding Distance)：计算某前端中的某个体与该前端中其他个体之间的距离，用以表征个体间的拥挤程度。显然拥挤距离的值越大，个体间就越不拥挤，种群的多样性就越好。需要注意，只有处于同一前端的个体间才需要计算拥挤距离，不同前端之间的个体计算拥挤距离是没有意义的
- 最优前端个体系数(ParetoFraction)：最优前端中的个体在种群中所占的比例
通过GUI方式调用gamultiobj
optimtool('gamultiobj')
通过命令行方式调用gamultiobj
[x,fval] = gamultiobj(fitnessfcn, nvars, A, b, Aeq, beq, lb, ub, options)
选择：stepgamultiobj
function parents = selectiontournament(expectation,nParents,options,tournamentSize)
%SELECTIONTOURNAMENT Each parent is the best of a random set.
% PARENTS = SELECTIONTOURNAMENT(EXPECTATION,NPARENTS,OPTIONS,TOURNAMENTSIZE)
% chooses the PARENTS by selecting the best TOURNAMENTSIZE players out of
% NPARENTS with EXPECTATION and then choosing the best individual
% out of that set.
%
% Example:
% Create an options structure using SELECTIONTOURNAMENT as the selection
% function and use the default TOURNAMENTSIZE of 4
% options = optimoptions('ga','SelectionFcn',@selectiontournament);
%
% Create an options structure using SELECTIONTOURNAMENT as the
% selection function and specify TOURNAMENTSIZE to be 3.
%
% tournamentSize = 3;
% options = optimoptions('ga','SelectionFcn', ...
% {@selectiontournament, tournamentSize});
%
% (Note: If calling gamultiobj, replace 'ga' with 'gamultiobj')
% Copyright 2003-2017 The MathWorks, Inc.
% How many players in each tournament?
if nargin < 4 || isempty(tournamentSize)
tournamentSize = 4;
end
% Choose the players
playerlist = ceil(size(expectation,1) * rand(nParents,tournamentSize));
% Play tournament
parents = tournament(playerlist,expectation);
function champions = tournament(playerlist,expectation)
%tournament between players based on their expectation
playerSize = size(playerlist,1);
champions = zeros(1,playerSize);
% For each set of players
for i = 1:playerSize
players = playerlist(i,:);
% For each tournament
winner = players(1); % Assume that the first player is the winner
for j = 2:length(players) % Winner plays against each other consecutively
score1 = expectation(winner,:);
score2 = expectation(players(j),:);
if score2(1) > score1(1)
winner = players(j);
elseif score2(1) == score1(1) && ... % score(1) is an integer rank
(length(score1) > 1 && score2(2) > score1(2))
% socre(2) may not be present for single objective problems
winner = players(j);
end
end
champions(i) = winner;
end
% 锦标赛选择函数
function parents = selectiontournament(expectation, nParents, options, tournamentSize)
% 输入
% expectation 在函数 stepgamultiobj 中赋值，赋值为[-rank,Distance]
% 函数 gamultiobj 将 tournamentSize 默认赋值为2
if nargin<4 || isempty(tournamentSize)
tournamentSize = 4;
end
playerlist = ceil(size(expectation,1)) * rand(nParents, tournamentSize);
% 产生 nParents 行 tournamentSize 矩阵，其元素为 [0,size(expectation,1)] 间的整数
% 选择操作是锦标赛，该矩阵可理解为运动员名单，该名单有 nParents 组，每组有 tournamentSize 个运动员
parents = tournament(playerlist, expectation);
function champions = tournament(playerlist,expectation)
playerSize = size(playerlist,1);
champions = zeros(1,playerSize);
for i=1:playerSize
players = playerlist(i,:);
winner = players(1);
for j=2:length(players)
score1 = expectation(winner,:); % 第一个运动员成绩
score2 = expectation(players(j),:); % 第 j 个运动员成绩
if score2(1) > score1(1)
winner = players(j); % 序值小的获胜
elseif score2(1) == score1(1)
try
if score2(2) > score1(2) % 序值相同，拥挤距离大的获胜
winner = players(j);
end
catch
end
end
end
champions(i) = winner;
end
交叉、变异、产生子种群、父子种群合并
nextPopulation = [xoverKids,mutateKids] % 生成子种群
nextPopulation = [eliteKids,xoverKids,mutateKids] % 生成子种群
population = [population,nextPopulation] %合并父子种群
计算序值和拥挤距离
% 非适配排序函数
function nondominatedRank = nonDominatedRank(score, nParent)
popSize = size(score,1); % 计算需要排序的个体数目
if nargin < 2
nParent = popSize;
end
rankedIndiv = false(popSize,1); % 表征个体是否已经被赋予序值的矩阵，默认每个个体都没排序
nondominatedRank = inf(popSize,1); % 函数返回的序值矩阵，初始序值为inf
rankCounter = 1; % 序值从1开始，依次加1
while ~all(rankedIndiv) && nnz(isfinite(nondominatedRank))<=nParent
front = false(popSize,1);
for p=1:popSize % 对每个个体排序
if rankedIndiv(p) % 已经排序的跳出循环
continue;
end
dominates = true; % 表征某个体是否被支配
for q=1:popSize % 个体p依次与个体q比较
if rankedIndiv(q) || p==q % 个体q已排序或p,q为同一体，跳出循环
continue;
end
if any(score(q,:)<score(p,:)) && all(score(q,:)<=score(p,:)) % 检查个体q是否支配个体p
dominates = false; % 个体q支配个体p，标志dominates位false
break;
end
end
if dominates
nondominatedRank(p) = rankCounter; % 个体p序值为true，则将当前序值 rankCounter赋予个体p，个体p被排序
front(p) = true; % 个体p序值为默认的inf，要等到下一个rnakCounter再被赋予序值
end
end
rankedIndiv(front) = true; % 更新 rankCounter
rankCounter = rankCounter+1;
end
function crowdingDistance = distancecrowding(pop,score,~,space)
%DISTANCECROWDING Assign local crowding distance to each individual
% CROWDINGDISTANCE = DISTANCECROWDING(POP,SCORE,OPTIONS,SPACE) Calculates
% crowding distance for each individuals on a non-dominated front. The
% fourth argument SPACE can be 'phenotype' or 'genotype' for distance to
% be in function space or decision variable space respectively.
%
% Example:
% Create an options structure using DISTANCECROWDING as the distance
% function in decision variable space
% options = optimoptions('gamultiobj','DistanceMeasureFcn',{@distancecrowding,'genotype'});
% Reference: Kalyanmoy Deb, "Multi-Objective Optimization using
% Evolutionary Algorithms", John Wiley & Sons ISBN 047187339, pg: 245 -
% 253
% Copyright 2007-2015 The MathWorks, Inc.
if nargin < 4
space = 'phenotype';
end
% Score should be finite to work with 'phenotype' distance
if strcmpi(space,'phenotype') && nnz(~isfinite(score)) == 0
y = score;
else % if strcmpi(space,'genotype')
y = pop;
end
popSize = size(y,1);
numData = size(y,2);
crowdingDistance = zeros(popSize,1);
for m = 1:numData
data = y(:,m);
% Normalize obective before computing distance
data = data./(1 + max(abs(data(isfinite(data)))));
[sorteddata,index] = sort(data);
% The best and worst individuals are at the end of Pareto front and
% they are assigned Inf distance measure
crowdingDistance([index(1),index(end)]) = Inf;
% Distance measure of remaining individuals
i = 2;
while i < popSize
crowdingDistance(index(i)) = crowdingDistance(index(i)) + ...
min(Inf, (data(index(i+1)) - data(index(i-1))));
i = i+1;
end
end
function [pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas,maxLinInfeas] = ...
trimPopulation(pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas, ...
maxLinInfeas,popSize,nScore,nParents,ParetoFraction)
%trimPopulation Trim population to get nParents with specified distribution
%among all fronts
% Copyright 2007-2014 The MathWorks, Inc.
if nScore > 1 % front diversity make sense for multi-objective
[pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas,maxLinInfeas,popSize] = ...
frontDiversity(pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas, ...
maxLinInfeas,popSize,nParents,ParetoFraction);
end
% trim individuals to nParents size based on distance
[pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas,maxLinInfeas] = distanceDiversity( ...
pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas,maxLinInfeas,nParents,popSize);
%-------------------------------------------------------------
function [pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas,maxLinInfeas] = distanceDiversity( ...
pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas,maxLinInfeas,nParents,popSize)
% Remove indivduals from last front if length(nonDomRank) is greater than nParents
to_remove = length(nonDomRank) - nParents;
if to_remove > 0
[~,I] = sortrows([nonDomRank, crowdingDistance],[1,-2]);
indiv_to_remove = I(popSize:-1:(popSize - to_remove + 1));
pop(indiv_to_remove,:) = [];
score(indiv_to_remove,:) = [];
nonDomRank(indiv_to_remove) = [];
crowdingDistance(indiv_to_remove) = [];
c(indiv_to_remove,:) = [];
ceq(indiv_to_remove,:) = [];
isFeas(indiv_to_remove) = [];
maxLinInfeas(indiv_to_remove) = [];
end
%-------------------------------------------------------------------------------
function [pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas,maxLinInfeas,popSize] = ...
frontDiversity(pop,score,nonDomRank,crowdingDistance,c,ceq,isFeas,maxLinInfeas, ...
popSize,nParents,ParetoFraction)
numRank = unique(nonDomRank); % numRank will be sorted
% Trim front size based on a pre-determined count
totalNumOfRank = length(numRank);
retainIndiv = zeros(totalNumOfRank,1); % Individuals to be retained in each front
availableIndiv = zeros(totalNumOfRank,1); % Available individuals in each front
% First rank is treated differently for number of desired individuals
availableIndiv(1) = nnz(nonDomRank == 1);
allowedIndiv = ceil(nParents*ParetoFraction);
if allowedIndiv <= availableIndiv(1)
retainIndiv(1) = allowedIndiv;
else
retainIndiv(1) = availableIndiv(1);
end
% If there are more than one fronts then we calculate the maximum number of
% individuals we can retain
if totalNumOfRank > 1
gpRatio = 0.8; % Geometric progression ratio
carryover = 0;
for i = numRank(2:end)'
availableIndiv(i) = nnz(nonDomRank == i);
% allowedIndiv redcuces as geometric progression
allowedIndiv = ceil(nParents*gpRatio^(i-1)*((1-gpRatio)/(1-gpRatio^totalNumOfRank))) + carryover;
if allowedIndiv <= availableIndiv(i)
retainIndiv(i) = allowedIndiv;
else
retainIndiv(i) = availableIndiv(i);
carryover = allowedIndiv - availableIndiv(i);
end
end
end
% Adjust retainIndiv if required
front = totalNumOfRank; increase = 0.10;
while sum(retainIndiv) <= ceil(nParents)
if retainIndiv(front) < availableIndiv(front)
retainIndiv(front) = min(availableIndiv(front),ceil(retainIndiv(front)*(1+increase)));
increase = increase*0.95;
end
if front == 1
front = totalNumOfRank;
increase = 0.10;
else
front = front - 1;
end
end
for i = numRank'
popIndices = 1:popSize;
index = (nonDomRank == i);
if retainIndiv(i) < nnz(index) % We need to trim the front
if retainIndiv(i) == 0
continue;
end
% Select individuals after a tournament selection
% Choose the players for each round
playerlist = popIndices(index);
losers = tournament(playerlist,crowdingDistance,retainIndiv(i));
pop(losers,:) = [];
score(losers,:) = [];
nonDomRank(losers) = [];
crowdingDistance(losers) = [];
c(losers,:) = [];
ceq(losers,:) = [];
isFeas(losers) = [];
maxLinInfeas(losers) = [];
popSize = size(score,1);
else
% No need to modify the front; do nothing
end
end
%-------------------------------------------------------------
function losers = tournament(playerlist,expectation,allowed)
%tournament between players based on their distance values
[~,sortedIndex] = sort(expectation(playerlist),'descend');
losers = playerlist(sortedIndex(allowed+1:end));
7. 基于多层编码遗传算法的车间调度
对于较为复杂的优化问题，染色体难以充分表达问题的解。多层编码把个体编码分为多层，每层含义不同，多层编码共同完整表达了问题的解，从而用一条染色体准确表达复杂问题的解。
车间调度是指根据产品制造的合理需求分配加工车间顺序，从而达到合理产生产品制造资源，提高企业经济效益的目的
1.个体编码
染色体对应工件的加工顺序，染色体前半部分表示所有工件在机器上的加工顺序，后半部分表示工件每道工序的加工机器序号
2.适应度值：时间最短，适应度值越高
function P=calp(S,PNumber)
% 功能说明： 根据基因S,计算调度工序P
% 输入参数：
% S 为基因
% PNumber 为工件个数
% 输出参数:
% P 为输出的调度工序
% 比如数字1002表示工件10的工序03
WNumber=length(S);%工序总个数
WNumber=WNumber/2;
%取工序基因，取基因的一半
S=S(1,1:WNumber);
%初始化
temp=zeros(1,PNumber);
P=zeros(1,WNumber);
%解码生成调度工序
for i=1: WNumber
%工序加+1
temp(S(i))=temp(S(i))+1;
P(i)=S(i)*100+temp(S(i));
end
function PVal=caltime(S,P,JmNumber,T,Jm)
% 功能说明： 根据调度工序,计算出调度工序时间
% 输入参数：
% P 为调度工序
% JmNumber 为机器个数
% T 为各工件各工序的加工时间
% Jm 为各工件各工序使用的机器
% 输出参数:
% PVal 为调度工序开始加工时间及完成时间
% 工件个数 工序个数
[PNumber MNumber]=size(Jm);
%取机器基因，取基因的一半
M=S(1,PNumber*MNumber+1:PNumber*MNumber*2);
%工序总个数
WNumber=length(P);
%初始化
TM=zeros(1,JmNumber);
TP=zeros(1,PNumber);
PVal=zeros(2,WNumber);
%计算调度工序时间
for i=1: WNumber
% 取机器号
val= P(1,i);
a=(mod(val,100)); %工序
b=((val-a)/100); %工件
Temp=Jm{b,a};
m=Temp(M(1,i));
%取加工时间
Temp=T{b,a};
t=Temp(M(1,i));
%取机器加工本工序的开始时间和前面一道工序的完成时间
TMval=TM(1,m);
TPval=TP(1,b);
%机器加工本工序的开始时间 大于前面一道工序的完成时 ，取机器加工本工序的开始时间
if TMval>TPval
val=TMval;
%取前面一道工序的完成时间
else
val=TPval;
end
%计算时间
PVal(1,i)=val;
PVal(2,i)=val+t;
%记录本次工序的机器时间和工序时间
TM(1,m)=PVal(2,i);
TP(1,b)=PVal(2,i);
end
3.选择操作：轮盘赌
% SELECT.M (universal SELECTion)
%
% This function performs universal selection. The function handles
% multiple populations and calls the low level selection function
% for the actual selection process.
%
% Syntax: SelCh = select(SEL_F, Chrom, FitnV, GGAP, SUBPOP)
%
% Input parameters:
% SEL_F - Name of the selection function
% Chrom - Matrix containing the individuals (parents) of the current
% population. Each row corresponds to one individual.
% FitnV - Column vector containing the fitness values of the
% individuals in the population.
% GGAP - (optional) Rate of individuals to be selected
% if omitted 1.0 is assumed
% SUBPOP - (optional) Number of subpopulations
% if omitted 1 subpopulation is assumed
%
% Output parameters:
% SelCh - Matrix containing the selected individuals.
% Author: Hartmut Pohlheim
% History: 10.03.94 file created
function SelCh = select(SEL_F, Chrom, FitnV, GGAP, SUBPOP);
% Check parameter consistency
if nargin < 3, error('Not enough input parameter'); end
% Identify the population size (Nind)
[NindCh,Nvar] = size(Chrom);
[NindF,VarF] = size(FitnV);
if NindCh ~= NindF, error('Chrom and FitnV disagree'); end
if VarF ~= 1, error('FitnV must be a column vector'); end
if nargin < 5, SUBPOP = 1; end
if nargin > 4,
if isempty(SUBPOP), SUBPOP = 1;
elseif isnan(SUBPOP), SUBPOP = 1;
elseif length(SUBPOP) ~= 1, error('SUBPOP must be a scalar'); end
end
if (NindCh/SUBPOP) ~= fix(NindCh/SUBPOP), error('Chrom and SUBPOP disagree'); end
Nind = NindCh/SUBPOP; % Compute number of individuals per subpopulation
if nargin < 4, GGAP = 1; end
if nargin > 3,
if isempty(GGAP), GGAP = 1;
elseif isnan(GGAP), GGAP = 1;
elseif length(GGAP) ~= 1, error('GGAP must be a scalar');
elseif (GGAP < 0), error('GGAP must be a scalar bigger than 0'); end
end
% Compute number of new individuals (to select)
NSel=max(floor(Nind*GGAP+.5),2);
% Select individuals from population
SelCh = [];
for irun = 1:SUBPOP,
FitnVSub = FitnV((irun-1)*Nind+1:irun*Nind);
ChrIx=feval(SEL_F, FitnVSub, NSel)+(irun-1)*Nind;
SelCh=[SelCh; Chrom(ChrIx,:)];
end
% End of function
selectJm
function SelS=selectJm(S,S_T)
Num=length(S_T);
MaxVal=0;
for i=1:Num
if S_T(i)>MaxVal;
MaxVal=S_T(i);
end
end
MaxVal=2*MaxVal;
for i=1:Num
S_T(i)=MaxVal-S_T(i);
end
eVal=0;
for i=1:Num
eVal=eVal+S_T(i);
end
for i=1:Num
P(i)=S_T(i)/eVal;
end
for i=2:Num
P(i)=P(i)+P(i-1);
end
num=rand;
SelS=1;
while(num>P(SelS))
SelS=SelS+1;
end
4.交叉操作
function NewChrom=across(Chrom,XOVR,Jm,T)
% Chrom=[1 3 2 3 1 2 1 3 2;
% 1 1 2 3 3 1 2 3 2;
% 1 3 2 3 2 2 1 3 1;
% 1 3 3 3 1 2 1 2 2;
% ];
% XOVR=0.7;
[NIND,WNumber]=size(Chrom);
WNumber=WNumber/2;
NewChrom=Chrom;%初始化新种群
[PNumber MNumber]=size(Jm);
Number=zeros(1,PNumber);
for i=1:PNumber
Number(i)=1;
end
%随机选择交叉个体(洗牌交叉)
SelNum=randperm(NIND);
Num=floor(NIND/2);%交叉个体配对数
for i=1:2:Num
if XOVR>rand;
Pos=unidrnd(WNumber);%交叉位置
while Pos==1
Pos=unidrnd(WNumber);
end
%取两交叉的个体
S1=Chrom(SelNum(i),1:WNumber);
S2=Chrom(SelNum(i+1),1:WNumber);
S11=S2;S22=S1; %初始化新的个体
%新个体中间片断的COPY
S11(1:Pos)=S1(1:Pos);
S22(1:Pos)=S2(1:Pos);
%比较S11相对S1,S22相对S2多余和缺失的基因
S3=S11;S4=S1;
S5=S22;S6=S2;
for j=1:WNumber
Pos1=find(S4==S3(j),1);
Pos2=find(S6==S5(j),1);
if Pos1>0
S3(j)=0;
S4(Pos1)=0;
end
if Pos2>0
S5(j)=0;
S6(Pos2)=0;
end
end
for j=1:WNumber
if S3(j)~=0 %多余的基因
Pos1=find(S11==S3(j),1);
Pos2=find(S4,1);%查找缺失的基因
S11(Pos1)=S4(Pos2);%用缺失的基因修补多余的基因
S4(Pos2)=0;
end
if S5(j)~=0
Pos1=find(S22==S5(j),1);
Pos2=find(S6,1);
S22(Pos1)=S6(Pos2);
S6(Pos2)=0;
end
end
% 保存交叉前的机器 基因
S1=Chrom(SelNum(i),:);
S2=Chrom(SelNum(i+1),:);
for k=1:WNumber
Pos1=Find(S11(k),S1);
S11(WNumber+k)=S1(WNumber+Pos1);
S1(Pos1)=0;
Pos1=Find(S22(k),S2);
S22(WNumber+k)=S2(WNumber+Pos1);
S2(Pos1)=0;
end
%生成新的种群
NewChrom(SelNum(i),:)=S11;
NewChrom(SelNum(i+1),:)=S22;
end
end
5.变异操作
main
%% 清空环境
clc;clear
%% 下载数据
load scheduleData Jm T JmNumber
%工序 时间
%% 基本参数
NIND=40; %个体数目
MAXGEN=50; %最大遗传代数
GGAP=0.9; %代沟
XOVR=0.8; %交叉率
MUTR=0.6; %变异率
gen=0; %代计数器
%PNumber 工件个数 MNumber 工序个数
[PNumber MNumber]=size(Jm);
trace=zeros(2, MAXGEN); %寻优结果的初始值
WNumber=PNumber*MNumber; %工序总个数
%% 初始化
Number=zeros(1,PNumber); % PNumber 工件个数
for i=1:PNumber
Number(i)=MNumber; %MNumber工序个数
end
% 代码2层，第一层工序，第二层机器
Chrom=zeros(NIND,2*WNumber);
for j=1:NIND
WPNumberTemp=Number;
for i=1:WNumber
%随机产成工序
val=unidrnd(PNumber);
while WPNumberTemp(val)==0
val=unidrnd(PNumber);
end
%第一层代码表示工序
Chrom(j,i)= val;
WPNumberTemp(val)=WPNumberTemp(val)-1;
%第2层代码表示机器
Temp=Jm{val,MNumber-WPNumberTemp(val)};
SizeTemp=length(Temp);
%随机产成工序机器
Chrom(j,i+WNumber)= unidrnd(SizeTemp);
end
end
%计算目标函数值
[PVal ObjV P S]=cal(Chrom,JmNumber,T,Jm);
%% 循环寻找
while gen<MAXGEN
%分配适应度值
FitnV=ranking(ObjV);
%选择操作
SelCh=select('rws', Chrom, FitnV, GGAP);
%交叉操作
SelCh=across(SelCh,XOVR,Jm,T);
%变异操作
SelCh=aberranceJm(SelCh,MUTR,Jm,T);
%计算目标适应度值
[PVal ObjVSel P S]=cal(SelCh,JmNumber,T,Jm);
%重新插入新种群
[Chrom ObjV] =reins(Chrom, SelCh,1, 1, ObjV, ObjVSel);
%代计数器增加
gen=gen+1;
%保存最优值
trace(1, gen)=min(ObjV);
trace(2, gen)=mean(ObjV);
% 记录最佳值
if gen==1
Val1=PVal;
Val2=P;
MinVal=min(ObjV);%最小时间
STemp=S;
end
%记录 最小的工序
if MinVal> trace(1,gen)
Val1=PVal;
Val2=P;
MinVal=trace(1,gen);
STemp=S;
end
end
% 当前最佳值
PVal=Val1; %工序时间
P=Val2; %工序
S=STemp; %调度基因含机器基因
%% 描绘解的变化
figure(1)
plot(trace(1,:));
hold on;
plot(trace(2,:),'-.');grid;
legend('解的变化','种群均值的变化');
%% 显示最优解
figure(2);
MP=S(1,PNumber*MNumber+1:PNumber*MNumber*2);
for i=1:WNumber
val= P(1,i);
a=(mod(val,100)); %工序
b=((val-a)/100); %工件
Temp=Jm{b,a};
mText=Temp(MP(1,i));
x1=PVal(1,i);
x2=PVal(2,i);
y1=mText-1;
y2=mText;
plotRec(x1,x2,mText);
plotRec(PVal(1,i),PVal(2,i),mText);
hold on;
fill([x1,x2,x2,x1],[y1,y1,y2,y2],[1-1/b,1/b,b/PNumber]);
text((x1+x2)/2,mText-0.25,num2str(P(i)));
end
aberranceJm
function ChromNew=aberranceJm(Chrom,MUTR,Jm,T)
%初始化
[NIND,WNumber]=size(Chrom);
WNumber=WNumber/2;
ChromNew=Chrom;
[PNumber MNumber]=size(Jm);
Number=zeros(1,PNumber);
for i=1:PNumber
Number(i)=1;
end
for i=1:NIND
%取一个个体
S=Chrom(i,:);
WPNumberTemp=Number;
for j=1:WNumber
JMTemp=Jm{S(j), WPNumberTemp(S(j))};
SizeTemp=length(JMTemp);
%是否变异
if MUTR>rand;
% 选择机器（随机选择）
% S(j+WNumber)=unidrnd(SizeTemp);
%选择机器（ 加工时间少的选择几率大）
if SizeTemp==1
S(j+WNumber)=1;
else
S(j+WNumber)=selectJm(S(j++WNumber),T{S(j),WPNumberTemp(S(j))});
end
end
WPNumberTemp(S(j))=WPNumberTemp(S(j))+1;
end
%数据放入新群
ChromNew(i,:)=S;
end
cal
function [PVal ObjV P S]=cal(Chrom,JmNumber,T,Jm)
% 功能说明： 根据基因群,计算出个群中每个个体的调度工序时间，
% 保存最小时间的调度工序和调度工序时间
% 输入参数：
% Chrom 为基因种群
% T 为各工件各工序使用的时间
% Jm 为各工件各工序使用的机器
% 输出参数:
% PVal 为最佳调度工序时间
% P 为最佳输出的调度工序
% ObjV 为群中每个个体的调度工序时间
% S 为最佳输出的调度基因
%初始化
NIND=size(Chrom,1);
ObjV=zeros(NIND,1);
% 工件个数 工序个数
[PNumber MNumber]=size(Jm);
for i=1:NIND
%取一个个体
S=Chrom(i,:);
%根据基因，计算调度工序
P= calp(S,PNumber);
%根据调度工序，计算出调度工序时间
PVal=caltime(S,P,JmNumber,T,Jm);
%取完成时间
MT=max(PVal);
TVal=max(MT);
%保存时间
ObjV(i,1)=TVal;
%初始化
if i==1
Val1=PVal;
Val2=P;
MinVal=ObjV(i,1);
STemp=S;
end
%记录 最小的调度工序时间、最佳调度工序时间 最佳输出的调度工序
if MinVal> ObjV(i,1)
Val1=PVal;
Val2=P;
MinVal=ObjV(i,1);
STemp=S;
end
end
%最佳调度工序时间 最佳输出的调度工序
PVal=Val1;
P=Val2;
S=STemp;
Find
function Pos=Find(FindVal,S)
% S=[1 3 2 3 1 2 1 3 2];
% FindVal=3;
[m n]=size(S);
Pos=-1;
for i=1:n
if FindVal==S(i)
Pos=i;
break;
end
end
PlotRec
function PlotRec(mPoint1,mPoint2,mText)
%mPoint1=10;
%mPoint2=15;
%mText=1;
vPoint =zeros(4,2);
vPoint(1,:)=[mPoint1,mText-0.5];
vPoint(2,:)=[mPoint2,mText-0.5];
vPoint(3,:)=[mPoint1,mText];
vPoint(4,:)=[mPoint2,mText];
plot([vPoint(1,1),vPoint(2,1)],[vPoint(1,2),vPoint(2,2)])
hold on;
plot([vPoint(1,1),vPoint(3,1)],[vPoint(1,2),vPoint(3,2)])
plot([vPoint(2,1),vPoint(4,1)],[vPoint(2,2),vPoint(4,2)])
plot([vPoint(3,1),vPoint(4,1)],[vPoint(3,2),vPoint(4,2)])
end
ranking
% RANKING.M (RANK-based fitness assignment)
%
% This function performs ranking of individuals.
%
% Syntax: FitnV = ranking(ObjV, RFun, SUBPOP)
%
% This function ranks individuals represented by their associated
% cost, to be *minimized*, and returns a column vector FitnV
% containing the corresponding individual fitnesses. For multiple
% subpopulations the ranking is performed separately for each
% subpopulation.
%
% Input parameters:
% ObjV - Column vector containing the objective values of the
% individuals in the current population (cost values).
% RFun - (optional) If RFun is a scalar in [1, 2] linear ranking is
% assumed and the scalar indicates the selective pressure.
% If RFun is a 2 element vector:
% RFun(1): SP - scalar indicating the selective pressure
% RFun(2): RM - ranking method
% RM = 0: linear ranking
% RM = 1: non-linear ranking
% If RFun is a vector with length(Rfun) > 2 it contains
% the fitness to be assigned to each rank. It should have
% the same length as ObjV. Usually RFun is monotonously
% increasing.
% If RFun is omitted or NaN, linear ranking
% and a selective pressure of 2 are assumed.
% SUBPOP - (optional) Number of subpopulations
% if omitted or NaN, 1 subpopulation is assumed
%
% Output parameters:
% FitnV - Column vector containing the fitness values of the
% individuals in the current population.
%
% Author: Hartmut Pohlheim (Carlos Fonseca)
% History: 01.03.94 non-linear ranking
% 10.03.94 multiple populations
function FitnV = ranking(ObjV, RFun, SUBPOP)
% Identify the vector size (Nind)
[Nind,ans] = size(ObjV);
if nargin < 2, RFun = []; end
if nargin > 1, if isnan(RFun), RFun = []; end, end
if prod(size(RFun)) == 2,
if RFun(2) == 1, NonLin = 1;
elseif RFun(2) == 0, NonLin = 0;
else error('Parameter for ranking method must be 0 or 1'); end
RFun = RFun(1);
if isnan(RFun), RFun = 2; end
elseif prod(size(RFun)) > 2,
if prod(size(RFun)) ~= Nind, error('ObjV and RFun disagree'); end
end
if nargin < 3, SUBPOP = 1; end
if nargin > 2,
if isempty(SUBPOP), SUBPOP = 1;
elseif isnan(SUBPOP), SUBPOP = 1;
elseif length(SUBPOP) ~= 1, error('SUBPOP must be a scalar'); end
end
if (Nind/SUBPOP) ~= fix(Nind/SUBPOP), error('ObjV and SUBPOP disagree'); end
Nind = Nind/SUBPOP; % Compute number of individuals per subpopulation
% Check ranking function and use default values if necessary
if isempty(RFun),
% linear ranking with selective pressure 2
RFun = 2*[0:Nind-1]'/(Nind-1);
elseif prod(size(RFun)) == 1
if NonLin == 1,
% non-linear ranking
if RFun(1) < 1, error('Selective pressure must be greater than 1');
elseif RFun(1) > Nind-2, error('Selective pressure too big'); end
Root1 = roots([RFun(1)-Nind [RFun(1)*ones(1,Nind-1)]]);
RFun = (abs(Root1(1)) * ones(Nind,1)) .^ [(0:Nind-1)'];
RFun = RFun / sum(RFun) * Nind;
else
% linear ranking with SP between 1 and 2
if (RFun(1) < 1 | RFun(1) > 2),
error('Selective pressure for linear ranking must be between 1 and 2');
end
RFun = 2-RFun + 2*(RFun-1)*[0:Nind-1]'/(Nind-1);
end
end;
FitnV = [];
% loop over all subpopulations
for irun = 1:SUBPOP,
% Copy objective values of actual subpopulation
ObjVSub = ObjV((irun-1)*Nind+1:irun*Nind);
% Sort does not handle NaN values as required. So, find those...
NaNix = isnan(ObjVSub);
Validix = find(~NaNix);
% ... and sort only numeric values (smaller is better).
[ans,ix] = sort(-ObjVSub(Validix));
% Now build indexing vector assuming NaN are worse than numbers,
% (including Inf!)...
ix = [find(NaNix) ; Validix(ix)];
% ... and obtain a sorted version of ObjV
Sorted = ObjVSub(ix);
% Assign fitness according to RFun.
i = 1;
FitnVSub = zeros(Nind,1);
for j = [find(Sorted(1:Nind-1) ~= Sorted(2:Nind)); Nind]',
FitnVSub(i:j) = sum(RFun(i:j)) * ones(j-i+1,1) / (j-i+1);
i =j+1;
end
% Finally, return unsorted vector.
[ans,uix] = sort(ix);
FitnVSub = FitnVSub(uix);
% Add FitnVSub to FitnV
FitnV = [FitnV; FitnVSub];
end
% End of function
REINS
% REINS.M (RE-INSertion of offspring in population replacing parents)
%
% This function reinserts offspring in the population.
%
% Syntax: [Chrom, ObjVCh] = reins(Chrom, SelCh, SUBPOP, InsOpt, ObjVCh, ObjVSel)
%
% Input parameters:
% Chrom - Matrix containing the individuals (parents) of the current
% population. Each row corresponds to one individual.
% SelCh - Matrix containing the offspring of the current
% population. Each row corresponds to one individual.
% SUBPOP - (optional) Number of subpopulations
% if omitted or NaN, 1 subpopulation is assumed
% InsOpt - (optional) Vector containing the insertion method parameters
% ExOpt(1): Select - number indicating kind of insertion
% 0 - uniform insertion
% 1 - fitness-based insertion
% if omitted or NaN, 0 is assumed
% ExOpt(2): INSR - Rate of offspring to be inserted per
% subpopulation (% of subpopulation)
% if omitted or NaN, 1.0 (100%) is assumed
% ObjVCh - (optional) Column vector containing the objective values
% of the individuals (parents - Chrom) in the current
% population, needed for fitness-based insertion
% saves recalculation of objective values for population
% ObjVSel - (optional) Column vector containing the objective values
% of the offspring (SelCh) in the current population, needed for
% partial insertion of offspring,
% saves recalculation of objective values for population
%
% Output parameters:
% Chrom - Matrix containing the individuals of the current
% population after reinsertion.
% ObjVCh - if ObjVCh and ObjVSel are input parameter, than column
% vector containing the objective values of the individuals
% of the current generation after reinsertion.
% Author: Hartmut Pohlheim
% History: 10.03.94 file created
% 19.03.94 parameter checking improved
function [Chrom, ObjVCh] = reins(Chrom, SelCh, SUBPOP, InsOpt, ObjVCh, ObjVSel);
% Check parameter consistency
if nargin < 2, error('Not enough input parameter'); end
if (nargout == 2 & nargin < 6), error('Input parameter missing: ObjVCh and/or ObjVSel'); end
[NindP, NvarP] = size(Chrom);
[NindO, NvarO] = size(SelCh);
if nargin == 2, SUBPOP = 1; end
if nargin > 2,
if isempty(SUBPOP), SUBPOP = 1;
elseif isnan(SUBPOP), SUBPOP = 1;
elseif length(SUBPOP) ~= 1, error('SUBPOP must be a scalar'); end
end
if (NindP/SUBPOP) ~= fix(NindP/SUBPOP), error('Chrom and SUBPOP disagree'); end
if (NindO/SUBPOP) ~= fix(NindO/SUBPOP), error('SelCh and SUBPOP disagree'); end
NIND = NindP/SUBPOP; % Compute number of individuals per subpopulation
NSEL = NindO/SUBPOP; % Compute number of offspring per subpopulation
IsObjVCh = 0; IsObjVSel = 0;
if nargin > 4,
[mO, nO] = size(ObjVCh);
if nO ~= 1, error('ObjVCh must be a column vector'); end
if NindP ~= mO, error('Chrom and ObjVCh disagree'); end
IsObjVCh = 1;
end
if nargin > 5,
[mO, nO] = size(ObjVSel);
if nO ~= 1, error('ObjVSel must be a column vector'); end
if NindO ~= mO, error('SelCh and ObjVSel disagree'); end
IsObjVSel = 1;
end
if nargin < 4, INSR = 1.0; Select = 0; end
if nargin >= 4,
if isempty(InsOpt), INSR = 1.0; Select = 0;
elseif isnan(InsOpt), INSR = 1.0; Select = 0;
else
INSR = NaN; Select = NaN;
if (length(InsOpt) > 2), error('Parameter InsOpt too long'); end
if (length(InsOpt) >= 1), Select = InsOpt(1); end
if (length(InsOpt) >= 2), INSR = InsOpt(2); end
if isnan(Select), Select = 0; end
if isnan(INSR), INSR =1.0; end
end
end
if (INSR < 0 | INSR > 1), error('Parameter for insertion rate must be a scalar in [0, 1]'); end
if (INSR < 1 & IsObjVSel ~= 1), error('For selection of offspring ObjVSel is needed'); end
if (Select ~= 0 & Select ~= 1), error('Parameter for selection method must be 0 or 1'); end
if (Select == 1 & IsObjVCh == 0), error('ObjVCh for fitness-based exchange needed'); end
if INSR == 0, return; end
NIns = min(max(floor(INSR*NSEL+.5),1),NIND); % Number of offspring to insert
% perform insertion for each subpopulation
for irun = 1:SUBPOP,
% Calculate positions in old subpopulation, where offspring are inserted
if Select == 1, % fitness-based reinsertion
[Dummy, ChIx] = sort(-ObjVCh((irun-1)*NIND+1:irun*NIND));
else % uniform reinsertion
[Dummy, ChIx] = sort(rand(NIND,1));
end
PopIx = ChIx((1:NIns)')+ (irun-1)*NIND;
% Calculate position of Nins-% best offspring
if (NIns < NSEL), % select best offspring
[Dummy,OffIx] = sort(ObjVSel((irun-1)*NSEL+1:irun*NSEL));
else
OffIx = (1:NIns)';
end
SelIx = OffIx((1:NIns)')+(irun-1)*NSEL;
% Insert offspring in subpopulation -> new subpopulation
Chrom(PopIx,:) = SelCh(SelIx,:);
if (IsObjVCh == 1 & IsObjVSel == 1), ObjVCh(PopIx) = ObjVSel(SelIx); end
end
% End of function
RWS
% RWS.m - Roulette Wheel Selection
%
% Syntax:
% NewChrIx = rws(FitnV, Nsel)
%
% This function selects a given number of individuals Nsel from a
% population. FitnV is a column vector containing the fitness
% values of the individuals in the population.
%
% The function returns another column vector containing the
% indexes of the new generation of chromosomes relative to the
% original population matrix, shuffled. The new population, ready
% for mating, can be obtained by calculating
% OldChrom(NewChrIx, :).
% Author: Carlos Fonseca, Updated: Andrew Chipperfield
% Date: 04/10/93, Date: 27-Jan-94
function NewChrIx = rws(FitnV,Nsel);
% Identify the population size (Nind)
[Nind,ans] = size(FitnV);
% Perform Stochastic Sampling with Replacement
cumfit = cumsum(FitnV);
trials = cumfit(Nind) .* rand(Nsel, 1);
Mf = cumfit(:, ones(1, Nsel));
Mt = trials(:, ones(1, Nind))';
[NewChrIx, ans] = find(Mt < Mf & ...
[ zeros(1, Nsel); Mf(1:Nind-1, :) ] <= Mt);
结果：
扫一扫
4873
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
