遗传算法（Genetic Algorithm）实现TSP问题

一、遗传算法的基本思想

遗传算法（Genetic Algorithm, GA）是模拟达尔文生物进化论的自然选择和遗传学机理的生物进化过程的计算模型，是一种通过模拟自然进化过程搜索最优解的方法。

遗传算法（Genetic Algorithm, GA）起源于对生物系统所进行的计算机模拟研究。它是模仿自然界生物进化机制发展起来的随机全局搜索和优化方法，借鉴了达尔文的进化论和孟德尔的遗传学说。其本质是一种高效、并行、全局搜索的方法，能在搜索过程中自动获取和积累有关搜索空间的知识，并自适应地控制搜索过程以求得最佳解。

达尔文进化论的主要观点是：物竞天择，适者生存。遗传算法（Genetic Algorithm）的基本思想就是模仿自然进化过程，通过对群体中具有某种结构形式的个体进行遗传操作，从而生成新的群体，逐渐逼近最优解。在求解过程中设定一个固定规模的种群，种群中的每个个体都表示问题的一个可能解，个体适应环境的程度用适应度函数判断，适应度差的个体被淘汰，适应度好的个体得以继续繁衍，繁衍的过程中可能要经过选择、交叉、变异，形成新的族群，如此往复，最后得到更多更好的解。

二、遗传算法的主要步骤

（1）编码：将问题的候选解用染色体表示，实现解空间向编码空间的映射过程。遗传算法不直接处理解空间的决策变量，而是将其转换成由基因按一定结构组成的染色体。编码方式有很多，如二进制编码、实数向量编码、整数排列编码、通用数据结构编码等等。本文将采用二进制编码的方式，将十进制的变量转换成二进制，用0和1组成的数字串模拟染色体，可以很方便地实现基因交叉、变异等操作。

（2）种群初始化：产生代表问题可能潜在解集的一个初始群体（编码集合）。种群规模设定主要有以下方面的考虑：从群体多样性方面考虑，群体越大越好，避免陷入局部最优；从计算效率方面考虑，群体规模越大将导致计算量的增加。应该根据实际问题确定种群的规模。产生初始化种群的方法通常有两种：一是完全随机的方法产生；二是根据先验知识设定一组必须满足的条件，然后根据这些条件生成初始样本。

（3）计算个体适应度：利用适应度函数计算各个个体的适应度大小。适应度函数（Fitness Function）的选取直接影响到遗传算法的收敛速度以及能否找到最优解，因为在进化搜索中基本不利用外部信息，仅以适应度函数为依据，利用种群每个个体的适应程度来指导搜索。

（4）进化计算：通过选择、交叉、变异，产生出代表新的解集的群体。选择（selection）：根据个体适应度大小，按照优胜劣汰的原则，淘汰不合理的个体；交叉（crossover）：编码的交叉重组，类似于染色体的交叉重组；变异（mutation）：编码按小概率扰动产生的变化，类似于基因突变。

（5）解码：末代种群中的最优个体经过解码实现从编码空间向解空间的映射，可以作为问题的近似最优解。这是整个遗传算法的最后一步，经过若干次的进化过程，种群中适应度最高的个体代表问题的最优解，但这个最优解还是一个由0和1组成的数字串，要将它转换成十进制才能供我们理解和使用。
以上文字原文链接：https://blog.csdn.net/u011125673/article/details/93393041

三、几点说明

1.问题的表示

路径表示是表示旅程对应的基因编码的最自然，最简单的表示方法。它在编码，解码，存储过程中相对容易理解和实现。例如：旅程（5-1-7-8-9-4-6-2-3）可以直接表示为（5 1 7 8 9 4 6 2 3）

2.轮盘赌选择算法（Roulette Wheel Selection）

它模拟博彩游戏的轮盘赌。一个轮盘被划分为N个扇形表示种群的一个染色体，而每个扇形的面积与它所表示的染色体的适应值成正比，为了选择种群中的个体，设想有一个指针指向轮盘，转动轮盘，当轮盘停止后，指针所之乡的染色体被选择。因为一个染色体的适应值越大表示该染色体的扇形面积就越大，因此它被选择的可能性也就越大。

3.顺序交叉概念（后面代码会用到）

 上图原文链接：遗传算法解决TSP问题_里昂科科的博客-CSDN博客

4.变异操作

在旅行路线上随机选择2个城市进行交换

5.算法的局限性

 当问题规模n比较小时，得到的一般都是最优解；当规模比较大时，一般只能得到近似解。这时可以通过增大种群大小和增加最大遗传代数或者增大变异几率使得优化值更接近最优解。

6.关于遗传算法所涉及的变量，在代码中均有说明，关于遗传算法更多的细节和背景说明，可以看其他博主，这里主要用代码展示其实现过程

 四、完整代码

%使用遗传算法实现TSP问题
clc,clear all;close all;
%随机化城镇
N=20;    %城镇数量
C=floor(rand(N,2).*100.+1);%随机化城市位置（小概率的两城市重合不考虑）
distance=dis_arr(C,N);  %城市距离矩阵
distance=distance+distance';
srt=[];%初始旅行路线
for i=1:N
    srt=[srt i];
end
subplot(221);
grid on;
title('原始TSP图像');
xlabel('x');
ylabel('y');
 netplot(C,srt,N);
%遗传算法数据初始化;
max_Nc=1000;%进化代数
sizepop=60;%种群规模
crs_p=0.8;  %交叉概率
var_p=N/150; %变异概率
individual=struct('fitness',[],'chrom',[]);%存储适应值（路线长度）和城市路线，这里适应值越小越好
best_fitness=[max_Nc,1];
%best_chrom=[max_Nc,1];
%----------------种群初始化---------------------
for i=1:sizepop
   individual.chrom(i,:)=randperm(N); 
   tra_path=individual.chrom(i,:);
   individual.fitness(i)=cpt_path_len(distance,tra_path,N);
end
[min_fitness,min_index]=min(individual.fitness); %此处将路线总路程当作适应值，按最小看
min_path=individual.chrom(min_index,:);
%---------------种群开始进化---------------------
for Nc=1:max_Nc
    %种群选择
    individual=select(individual,sizepop);
    [min_fitness,min_index]=min(individual.fitness); %此处将路线总路程当作适应值，按最小看
    min_path=individual.chrom(min_index,:);
    % 种群交叉
    individual=cross(individual,sizepop,crs_p,N,distance);  
    %种群变异
    individual=variation(individual,sizepop,var_p,N,distance);
    %记录每一次迭代的适应值和旅行路程
    
    [min_fitness1,min_index1]=min(individual.fitness); %此处将路线总路程当作适应值，按最小看
    %将遗传前最好的个体留下来
    if(min_fitness<min_fitness1)
        individual.fitness(min_index1)=min_fitness;
        individual.chrom(min_index1,:)=min_path;
    end
    min_path=individual.chrom(min_index1,:);  
    best_fitness(Nc,1)=min_fitness;
end
subplot(222);
netplot(C,min_path,N);
grid on
title('TSP问题优化结果');
xlabel('x');
ylabel('y');
subplot(223);
grid on
title('每次迭代总路程');
plot(best_fitness);
xlabel('迭代次数');
ylabel('路线长度');
%---------------函数定义-----------------------
function var_individual=variation(individual,sizepop,var_p,N,distance)
    var_individual=individual;
    [~,index]=min(individual.fitness);
    for i=1:sizepop
        if i==index
            continue;
        end
        if(var_p>rand)
            p1=floor(rand*N+1);
            p2=floor(rand*N+1);
            while p1==p2
                p1=floor(rand*N+1);
                p2=floor(rand*N+1);
            end      %随机选2个位置进行交换      
            tmp=var_individual.chrom(i,p1);
            var_individual.chrom(i,p1)=var_individual.chrom(i,p2);            
            var_individual.chrom(i,p2)=tmp;
            var_individual.fitness(i)=cpt_path_len(distance,var_individual.chrom(i,:),N); %计算适应度           
        end
    end
end
function crs_individual=cross(individual,sizepop,crs_p,N,distance) %交叉运算
    crs_individual=individual;
    rot=floor(sizepop/2)*2;  
    for i=1:2:rot        
       pick=rand;
       if(crs_p>pick)
           start=randi([1,N-ceil(N/2)]);
           num=randi([0,ceil(N/2)]);
           cut1=crs_individual.chrom(i,start:num+start);
           cut2=crs_individual.chrom(i+1,start:num+start);
           crschrom1=crs_individual.chrom(i,:);
           crschrom2=crs_individual.chrom(i+1,:);
           schrom1=[];
           schrom2=[];
           for j=1:N
              if (length(find(cut1==crschrom2(j)))==0) %采用顺序交叉，看文章图示说明
                  schrom2=[schrom2 crschrom2(j)];              
              end
              if (length(find(cut2==crschrom1(j)))==0) 
                  schrom1=[schrom1 crschrom1(j)];                 
              end
           end           
           n=1;m=1;
           newchrom1=[];
           newchrom2=[];
           while(n<=N)
               if(n==start)
                  newchrom1=[newchrom1 cut2];
                  newchrom2=[newchrom2 cut1];
                  n=n+num+1;
               else
                   newchrom1=[newchrom1 schrom1(m)];
                   newchrom2=[newchrom2 schrom2(m)];
                   m=m+1;
                   n=n+1;
               end              
           end
           crs_individual.chrom(i,:)=newchrom1;
           crs_individual.fitness(i)=cpt_path_len(distance,newchrom1,N);
           crs_individual.chrom(i+1,:)=newchrom2;
           crs_individual.fitness(i+1)=cpt_path_len(distance,newchrom2,N);
       end     
    end
end
function slt_individual=select(individual,sizepop) %选择运算
    p1=sum(1./individual.fitness);
    p=1./individual.fitness./p1;
    p=cumsum(p);%以下开始采用轮盘赌选择法
    slt_individual=individual;    
    for i=1:sizepop
        select1=find(p>rand);        
        slt_individual.fitness(i)=individual.fitness(select1(1));
        slt_individual.chrom(i,:)=individual.chrom(select1(1),:);        
    end
    %将上次适应度前5%的精英留下
    num=ceil(sizepop/20);
    for i=1:num
       [last_minfit,index]=min(individual.fitness );
       last_minpath=individual.chrom(index,:);
       individual.fitness(index)=inf;
       [~,maxindex]=max(slt_individual.fitness);
       slt_individual.fitness(maxindex)=last_minfit;
        slt_individual.chrom(maxindex,:)=last_minpath;
    end    
end
function path_len=cpt_path_len(distance,path,N) %计算路径长度
    path_len=0;
    for i=1:N-1
        path_len=path_len+distance(path(i),path(i+1));
    end
    path_len=path_len+distance(path(N),path(1));
end
function len_arr=dis_arr(C,N) %求距离矩阵
len_arr=zeros(N);
    for i=1:N-1
       for j=i+1:N
           len_arr(i,j)=sqrt((C(i,1)-C(j,1))^2+(C(i,2)-C(j,2))^2);
       end
    end
end
function netplot(C,path,N) %画图
    hold on;
    for i=1:N-1
        k=path(i);
        m=path(i+1);
        plot(C(k,1),C(k,2),'ro');
        line([C(k,1) C(m,1)],[C(k,2) C(m,2)],'Color','green');
    end
    k=path(N);
    m=path(1);
    plot(C(k,1),C(k,2),'ro');
    line([C(k,1) C(m,1)],[C(k,2) C(m,2)],'Color','green');
end

五、结果展示

 代码为自己独立编写完成，可能还有不少值得优化的地方哦~