以 TSP 问题为例谈对 遗传算法 GA 的理解【MATLAB代码】
遗传算法，顾名思义，是一种仿生学算法，原理就是中学学过的达尔文定律。
用8个字概括就是：
物竞天择，适者生存。
根据我的亲身体验，刚接触智能算法的时候总是想刨根问底，实际上这是没有很大意义的，需要的是先对整个算法的大框架了解并理解就行了。
当你了解了遗传算法后，我想很多人和我最开始学习遗传算法一样，觉得框架很简单，但就是不能自己来完成，这就需要大量的敲代码了，从借鉴别人的代码开始，慢慢学习。然后你会发现遗传算法也不过如此，当然这是对你见过的遗传算法的问题错觉，要想完全掌握遗传算法，核心是能自己设计遗传算法的 基因编码 、 基因变异 和 基因重组 过程以及 个体适应度函数 的选择。
下面我先简单介绍一下我对遗传算法的理解：
实现遗传算法，举一个通俗易懂的例子，就相当于是袁隆平爷爷做杂交水稻实验一样，需要先建立一个大规模水稻筛选群（在遗传算法中称为种群），当然这些水稻个体间的基因一般是不同的，这样才有筛选的必要，（每个水稻的基因这个在遗传算法中也称为个体的基因），然后就是观察不同水稻的产量情况（这个产量在遗传算法中称为适应度），产量高的就留下来（在遗传算法中就是选出来适应度较高的个体），当然产量高的水稻也不一定能全部活下来，如果一不小心把水稻踩了一脚给踩死了，虽然是小概率事件，但这毕竟是事实，完成了这一步这样就算是完成了第一次的筛选。然后后面的事情就是一直重复这样的过程，但是想要产生更加高产的水稻，无非两种办法，一个是基因重组，一个是基因突变，在遗传算法中确实也是这样实现的。通过这两种方式更新种群，其中重组和变异的发生概率可以设置，更新种群后再继续筛选，一直这样，慢慢的就能筛选出最好的或者次好的个体了。但是需要这样更新多少代就要取决于种群的规模、变异和重组的方式。原因很简单：数量多，变异和重组的方式才更多，不同的重组方式和变异方式，重组概率，变异概率都影响新型个体的数目。
下面我将引入一个典型的 TSP 问题作为讲解：
假设有50个城市，已知各个城市之间的距离，现在要寻找一条能走完所有城市的最短路径，那么这个问题肯定事不能用暴力枚举的，因为使用暴力枚举确实太恐怖了，50！估计要把电脑枚举炸裂，因此就用解决该问题比较简单的一种方式------> 遗传算法 试一试吧。后面我将再介绍一种 蚁群算法 的解决方式。当然我认为这两种算法是差不多的，但 蚁群算法 的适应度能做到定点更小，或者说我认为 蚁群算法 在某种程度上包含了遗传算法，而遗传算法在某一些方面比 蚁群算法 运用面更广。
我对实现 TSP 的整个遗传算法的流程理解如下：
1、 选择一个合适的基因编码方式，这里的话就随机生成一个一维数组吧，数组长度和城市数目一样，这个数组中的数字就是城市的序号，数组元素的顺序就表示访问城市的顺序。
2、 选择一个大小合适的种群，也就是说有多少个不同的存储城市访问顺序的数组。
3、 有了编码，有了种群后，就设计一个适应度函数，明确我们筛选的原则，这里的话就按照总路程的大小来衡量吧。
4、 当然我们还要考虑怎么来筛选，但咱也不能让适应度不是很好的个体全死掉吧，比较还是要救一下吧，但生死有命富贵在天，能不能救得活就看天意吧，而对于适应度高的个体咱也不能保证他就能长命百岁吧，万一他突发心肌梗塞，只能说他活下来的概率大一点。
那么，怎么来模拟这个过程呢？让适应度高的大概率存活，让适应度低的小概率存活。
这里就用最简单的轮盘赌吧，不要嫌弃这个方式简单，其实这是用的最多的方式，当然如果你还不知道轮盘赌是怎么回事的话，就想象一下幸运大转盘吧，转到奖品的概率和转到谢谢惠顾的概率和其所占轮盘的面积是成正比的，轮盘赌就是这样的。好了，你现在想象完了后基本明白了，一定在想怎么用程序实现。这个轮盘赌的实现方式是这样的：用一个随机数做比较就可以实现了啊？
可能会有人这样想：将所有适应度归一化为0——1之间的数，然后每去筛选一个个体就生成一个0——1的随机数，如果其对应的适应度值小于这个随机数就认为他被Game Over 了。
但是这样其实是不正确的，因为这样并不是轮盘赌，幸运大转盘上的总面积是由所有的谢谢惠顾、一等奖、二等奖、三等奖构成的，因此在做轮盘赌的时候这样才是正确的：将所有的适应度值累加，适应度的总和就相当于是幸运大转盘的总面积，而由于使用了累加，那么累加数组内的任意两个个体位置之间的差值一定是和这个个体的适应度成正比的，这样就构成了一个正确的轮盘，然后就是进行筛选了。比如一共有100个个体，我要筛选出80%的个体，那么就需要在0——适应度总和之间生成80个随机数，当这个随机数落入到哪里说明这个个体就被选中了，这样的方式是一定可以选出刚好80个个体的。
5、 完成选择操作后，一部分个体不幸死了，然后就需要补充一部分进来，这一部分补充的个体一般是通过上一代种群中的优秀个体重组和变异得到的。
6、 然后就重复上面的步骤，继续选择，继续变异…
7、 当进化到你设定的代数后程序就截止了。也就可以得到一个最优解或者次优解。
下面直接上MATLAB代码吧【摘录自 智能算法30个例子》----北航出版社 】，看了代码相信你一定能明白：
clear all;
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('-------------------------------------------------------------')
%% 计算两两城市之间的距离
%输入 a 各城市的位置坐标
%输出 D 两两城市之间的距离
function D=Distanse(a)
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
%% 画路径函数
%输入
% Chrom 待画路径
% X 各城市坐标位置
function DrawPath(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
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)
%% 适配值函数
%输入：
%个体的长度（TSP的距离）
%输出：
%个体的适应度值
function FitnV=Fitness(len)
FitnV=1./len;
%% 初始化种群
%输入：
% NIND：种群大小
% N： 个体染色体长度（这里为城市的个数）
%输出：
%初始种群
function Chrom=InitPop(NIND,N)
Chrom=zeros(NIND,N);%用于存储种群
for i=1:NIND
Chrom(i,:)=randperm(N);%随机生成初始种群
end
%% 变异操作
%输入：
%SelCh 被选择的个体
%Pm 变异概率
%输出：
% SelCh 变异后的个体
function SelCh=Mutate(SelCh,Pm)
[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
%% 输出路径函数
%输入：R 路径
function p=OutputPath(R)
R=[R,R(1)];
N=length(R);
p=num2str(R(1));
for i=2:N
p=[p,'—>',num2str(R(i))];
end
disp(p)
%% 计算各个体的路径长度
% 输入：
% D 两两城市之间的距离
% Chrom 个体的轨迹
function len=PathLength(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
%% 交叉操作
% 输入
%SelCh 被选择的个体
%Pc 交叉概率
%输出：
% SelCh 交叉后的个体
function SelCh=Recombin(SelCh,Pc)
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
%% 重插入子代的新种群
%输入：
%Chrom 父代的种群
%SelCh 子代种群
%ObjV 父代适应度
%输出
% Chrom 组合父代与子代后得到的新种群
function Chrom=Reins(Chrom,SelCh,ObjV)
NIND=size(Chrom,1);
NSel=size(SelCh,1);
[~,index]=sort(ObjV);
Chrom=[Chrom(index(1:NIND-NSel),:);SelCh];
%% 进化逆转函数
%输入
%SelCh 被选择的个体
%D 个城市的距离矩阵
%输出
%SelCh 进化逆转后的个体
function SelCh=Reverse(SelCh,D)
[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,:);
%% 选择操作
%输入
%Chrom 种群
%FitnV 适应度值
%GGAP：选择概率
%输出
%SelCh 被选择的个体
function SelCh=Select(Chrom,FitnV,GGAP)
NIND=size(Chrom,1);
NSel=max(floor(NIND*GGAP+.5),2);
ChrIx=Sus(FitnV,NSel);
SelCh=Chrom(ChrIx,:);
% 输入:
%FitnV 个体的适应度值
%Nsel 被选择个体的数目
% 输出:
%NewChrIx 被选择个体的索引号
function NewChrIx = Sus(FitnV,Nsel)
% 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
7万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
