【智能算法】协同进化算法

原创已于 2024-06-17 19:13:14 修改 · 3.5k 阅读

20 ·

CC 4.0 BY-SA版权

文章标签：

#人工智能 #机器学习 #网络 #算法

于 2024-06-16 23:27:12 首次发布

智能算法专栏收录该内容

243 篇文章

订阅专栏

该文章已生成可运行项目，

1、进化算法

自从达尔文的生物进化论被接受，基于自然界中生物优胜劣汰的生存规则发展起来的生物进化的理论研究得到了空前的发展。将生物遗传变异、优胜劣汰的生存机制应用到优化领域，就得到了进化计算(Evolutionary Computation, EC)。以种群形式存在的物种，想要生存下去，就必须通过遗传变异来适应环境，通过自身的不断完善来适应生存环境。遗传的目的在于将父代的优良性能传递给子代，让子代能更好更快地适应环境，而变异的目的则是使子代产生更好地更能适应当前环境的个体。适应能力好的个体被保留下来，而适应能力差的则被慢慢淘汰，这就是“物竞天择，适者生存”的原理。通过不断的自我完善和更新，物种的适应能力就得到了提升。

进化计算源于对生物进化论的研究，其基本原理如下：生物从简单到复杂、从低级到高级的进化过程本身就具有相当的自然性、自适应性、并行性和鲁棒性，其本身就是一个优化的过程；而进化计算正是模拟了这一优化过程，通过类似优胜劣汰遗传变异的操作来执行算法。

进化计算的提出是优化领域的一个里程碑。它比传统的优化算法具有更好的性能，主要表现为智能性和并行性，具体区别如下：

进化计算对变量进行编码，使所有被优化变量作为一个整体被操作。
适应度函数不受连续可微的限制，扩大了应用领域。
优化采用多点搜索，不同于传统算法的单点搜索，全局寻优性能得到提高。
采用启发式搜索，具有一定的隐含并行性，适用于非线性优化。
在适应度函数的引导下，采用搜索方向明确的随机化基本操作。

同时，进化计算也有一些自身的限制，比如：局部搜索能力较差、多样性有待提高、具体的遗传变异算子应结合不同情况调整等。进化计算的步骤可简要表述为图3.1。

依据不同的分类标准，不同的学者对进化算法做了不同的分类，进化算法的典型分类如下：遗传算法(GeneticAlgorithm，GA)、进化规划(Evolutionary Programming，EP)和进化策略(Evolutionary Strategies，ES)。这三种进化算法都是基于多点并行迭代优化算法提出的。

2、协同进化算法

达尔文通过对自然界生物进化过程的研究，提出了“适者生存，优胜劣汰”的生物进化论。生物进化论指出，自从生命诞生以来，以种群形式存在的物种经历了从简单到复杂，从低级到高级的进化过程。在这个过程中，越能适应环境并随环境改变的个体和种群，越能保持和生存下去；而不能随环境变化作出改变的物种则慢慢被淘汰，也就是说生物的进化是一个竞争的过程。但是生物的进化关系不仅仅只有竞争，同样存在着协同关系。个体是以种群存在的，而不同的种群构成了生态系统，种群是依赖生态系统存在的，也就是说生态系统中的各个种群不仅仅是竞争关系，同样存在着相互协同、共同发展的关系；任何种群想要独立存在都是不可能的。不同的种群相互竞争、同时又相互合作，使得各种群一代一代地发展下去。种群内存在着协同关系，比如，蜜蜂的分工合作；蚂蚁的相互配合。同时种群间也存在着协同关系，比如，蜜蜂和花粉的合作，鳄鱼和鸟的合作等。

达尔文进化论经过长期的发展，逐渐变得完善。Jazen认为协同进化是指两个种群中个体的特性为适应对方而相互进化；罗森则认为协同进化有联系的休戚与共的种群之间的联合进化。卡罗尔认为，广义上的协同进化是生物之间的交互变化。不论怎么定义，可以看出，协同进化的主体是种群，种群之间通过发生一定的关系来达到共同进化的目的。

协同进化算法的生物学基础就是生物之间的协同进化发展。各个不同种群通过或合作或竞争来达到一起进化的目的。基于此提出的协同进化算法已经被广泛应用到了各个领域，比如机器学习等，同时可以跟神经网络、心理学等结合使用。

协同进化的研究范围包括：竞争物种之间的协同进化；捕食者和猎物系统的协同进化；寄生物和寄主系统的协同进化；拟态的协同进化；互利作用的协同进化等等。

3、协同进化算法基本原理

协同进化算法和普通进化算法的最大的区别在于协同进化算法利用了种群之间相互协同、共同进化的关系，而不仅仅是只考虑竞争的关系。协同进化还有一些自己本身的特点：

协同进化算法的本质是多种群同时进化，跟普通进化算法的单种群进化不同，协同进化算法的优势就在于种群间的相互作用。
因为协同进化算法具有很大的搜索空间，所以避免了早熟收敛。
协同进化算法具有并行性，但它主要是依赖种群间的协同模式来实现并行性的；而普通的并行进化算法是利用多个处理器来实现的。

根据协同方式的不同，协同进化算法大致可以分为以下几种：

（1）合作型协同进化；

（2）竞争型协同进化；

（3）辅助式共生协同进化；

（4）基于混合策略的协同进化；

（5）基于种群密度的协同进化；

目前，合作型协同进化和竞争型协同进化应用比较广泛，基于混合策略的协同进化思想跟别的协同进化思想略有不同，它不强调多种群，而是将解的个体和不同的进化策略结合，通过调整不同的策略来达到协同的目的。下面对这几种协同进化算法做一简单介绍：

（1）合作型协同进化

遗传算法进化采用单种群进化，每个个体表示优化问题的一个完整的解，而合作型协同进化算法采用多种群进化，同时每个种群中的个体只维持整体解的一部分；同时个体适应度的计算不再只依赖于本种群，需要联合别的种群中的个体进行适应度评价。合作型协同进化算法将复杂问题分解为多个子问题，每个子问题维持一个种群，各子种群进行独立进化，进化过程中个体适应度的评价需要和其他子种群进行信息交互，而不仅仅依赖于当前子种群；进化完成各子种群的解组合成为问题的最终解。

合作型协同进化算法步骤如下：

①编码：根据实际问题，采用不同的编码方式。

②种群初始化：依据实际问题，分为多个种群，对各种群随机初始化，每个种群维持解

的一部分。

③适应度评价：对于每个种群中的待评价个体，需要和别的种群中的个体联合来进行适

应度的评价：其他联合个体可以选取其余种群中的当前最优个体，也可以选取随机个体。

④进化操作：每个种群执行或相同或不同的进化流程，比如选择交叉变异等遗传流程，生成新的子种群。

⑤判断是否结束：如果满足各种群组成的完整解满足终止条件，或者达到进化代数等，

则结束进化；否则返回③继续。

（2）竞争型协同进化

竞争型协同进化算法是基于自然界中捕食和被捕食的关系提出的。它们之间是相互竞争的关系，被捕食者必须通过各种方式来提升自己的防御能力，才能保证自己生存下去；而捕食者也必须加强自己的攻击力，来维持生存。

在竞争过程中，一个种群中个体的适应度通过和另一个种群中的个体竞争来计算得到的。两个种群根据待评价个体所属的范围，依次作为主种群和从种群。计算出适应度以后，再依赖这个适应度来进行选择交叉编译等遗传操作。随着进化的不断进行，各个种群将进化出不同的特性。

基于上述的描述，竞争型协同进化的步骤一般如下：

初始化两个种群：种群1和种群2。
种群1中的个体通过与种群2的子集竞争，得到适应度。
利用此适应度，对种群1执行遗传操作；直到满足其终止条件。
种群2中的个体和种群1的子集竞争，得到适应度。
基于上述适应度，对种群2执行遗传操作；直到满足其终止条件。
判断是否满足算法的总终止条件，满足就结束算法；否则转回②继续。

（3）基于混合策略的协同进化

基于混合策略的协同进化不同于上述两种协同进化采用多种群的形式，它是由董红斌等人为了解决数值优化问题而提出的一种不同思想的协同进化[39]。它实际上是各种不同的进化策略和种群中个体结合的一种算法。它将个体和对应的进化策略绑定，个体随进化过程改变，相应的进化策略也随之改变，通过个体和进化策略的协同来达到进化的目的。

协同进化算法更精细地描述了生物进化过程，所以比传统的进化算法具有更多的优势，能解决更多的问题，尤其是复杂问题。目前，协同进化算法已经被广泛应用到了各个领域并取得了很好的效果，大致可以分为如下几个领域：

（1）大规模的复杂优化问题。这类型的优化问题一般都由多个子域组成，合作型协同进化算法正是适用于解决此类问题，它可以将各个子域或者各个子问题分而治之，对每个子问题维持一个种群，通过多个种群联合进化来得到问题的最终解。各个子种群通过个体适应度的评价相互联系，相互带动，合作型协同进化通过分解问题的策略，极大地提高了算法的效率。

（2）问题候选解需要通过交互测试来进行性能评价的问题。竞争型协同进化可以解决此类问题。它维持一个主种群和一个测试种群，分别放置候选解和测试样本。候选解的适应度是通过其在测试种群中的表现来给定的，主种群正是基于这种适应度的评价方式来进化的；而测试种群中个体的适应度则是基于它给予候选解的测试信息量来评定的。通过这样的方式使得当前候选解评价总是能选取到合适的测试样本，避免了大测试集效率低的问题和固定测试集不准确的问题，使得进化更快收敛，更有针对性地收敛。

（3）一些结构非常复杂的优化问题。对于这类问题，普通进化算法很难在它极不规则的搜索空间中、在缺乏相关信息的情况下找到最优解。而协同进化算法的多种群机制正可以解决此类问题，比如可以通过构建辅助种群，采用竞争型协同进化来解决。

虽然协同进化算法具有如此多的优点，但它也有不足的地方，比如虽然协同进化算法适应度评价采用联合其他种群来评价的方式，是适应度评价的一个突破；不过这种评价方式也存在一些缺点：不稳定、容易出现进化动态性等。

数值优化作为测试算法性能的一种最基本方法，是算法用于实际生活应用领域之前的必要步骤，近年来许多学者对协同进化数值优化作出了大量的研究，其中利用精英作为协同进化的主体，提出的精英协同进化算法在数值优化领域取得了非常好的效果。

4、协同进化算法matlab代码解析

function m_main()
%author Xing Peng
%@copyright reserved
G=[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 1 1 1 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0;
   0 0 1 1 1 0 0 1 0 0 0 0 0 0 1 1 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0;
   0 0 0 1 1 0 0 1 0 0 0 0 0 0 0 0 0 0 1 1;
   0 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0;
   1 1 1 1 1 0 0 1 1 1 0 0 0 0 1 0 0 0 0 0;
   1 1 1 1 1 0 0 1 1 1 0 0 0 1 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 1 0 0 0 0 0 1 1 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0;
   0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 1 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0;];
MM=size(G,1);      % G 地形图为01矩阵，如果为1表示障碍物
Tau=ones(MM*MM,MM*MM);% Tau 初始信息素矩阵（认为前面的觅食活动中有残留的信息素）
Tau=8.*Tau;
K=100;          % K 迭代次数（指蚂蚁出动多少波）
M=50;         % M 蚂蚁个数（每一波蚂蚁有多少个）
S=1 ;         % S 起始点（最短路径的起始点）
E=MM*MM;         % E 终止点（最短路径的目的点）
Alpha=1;    % Alpha 表征信息素重要程度的参数
Beta=7;     % Beta 表征启发式因子重要程度的参数
Rho=0.3 ;     % Rho 信息素蒸发系数
Q=1;          % Q 信息素增加强度系数
minkl=inf;
mink=0;
minl=0;

D=G2D(G);
N=size(D,1);%N表示问题的规模（象素个数）
a=1;%小方格象素的边长
Ex=a*(mod(E,MM)-0.5);%终止点横坐标
if Ex==-0.5
Ex=MM-0.5;
end
Ey=a*(MM+0.5-ceil(E/MM));%终止点纵坐标
Eta=zeros(N);%启发式信息，取为至目标点的直线距离的倒数
%下面构造启发式信息矩阵
for i=1:N
 ix=a*(mod(i,MM)-0.5);
   if ix==-0.5
   ix=MM-0.5;
   end
iy=a*(MM+0.5-ceil(i/MM)); 
   if i~=E
   Eta(i)=1/((ix-Ex)^2+(iy-Ey)^2)^0.5;
   else
   Eta(i)=100;
   end
end
ROUTES=cell(K,M);%用细胞结构存储每一代的每一只蚂蚁的爬行路线
PL=zeros(K,M);%用矩阵存储每一代的每一只蚂蚁的爬行路线长度
%%-----------启动K轮蚂蚁觅食活动，每轮派出M只蚂蚁--------------------
for k=1:K
   for m=1:M
    %% 第一步：状态初始化
    W=S;%当前节点初始化为起始点
    Path=S;%爬行路线初始化
    PLkm=0;%爬行路线长度初始化
    TABUkm=ones(N);%禁忌表初始化
    TABUkm(S)=0;%已经在初始点了，因此要排除
    DD=D;%邻接矩阵初始化
    %% 第二步：下一步可以前往的节点
    DW=DD(W,:);
    DW1=find(DW);
    for j=1:length(DW1)
        if TABUkm(DW1(j))==0
            DW(DW1(j))=0;
        end
    end
    LJD=find(DW);
    Len_LJD=length(LJD);%可选节点的个数

   %% 觅食停止条件：蚂蚁未遇到食物或者陷入死胡同
    while W~=E&&Len_LJD>=1
    %% 第三步：转轮赌法选择下一步怎么走
         PP=zeros(Len_LJD);
         for i=1:Len_LJD
              PP(i)=(Tau(W,LJD(i))^Alpha)*((Eta(LJD(i)))^Beta);
         end
        sumpp=sum(PP);
        PP=PP/sumpp;%建立概率分布
        Pcum(1)=PP(1);
        for i=2:Len_LJD
              Pcum(i)=Pcum(i-1)+PP(i);
        end
        Select=find(Pcum>=rand);
        to_visit=LJD(Select(1));
      %% 第四步：状态更新和记录
        Path=[Path,to_visit];%路径增加
        PLkm=PLkm+DD(W,to_visit);%路径长度增加
        W=to_visit;%蚂蚁移到下一个节点
        for kk=1:N
            if TABUkm(kk)==0
                 DD(W,kk)=0;
                 DD(kk,W)=0;
            end
        end
       TABUkm(W)=0;%已访问过的节点从禁忌表中删除
       DW=DD(W,:);
       DW1=find(DW);
       for j=1:length(DW1)
           if TABUkm(DW1(j))==0
               DW(j)=0;
           end
       end
       LJD=find(DW);
       Len_LJD=length(LJD);%可选节点的个数
     end
%%第五步：记下每一代每一只蚂蚁的觅食路线和路线长度
    ROUTES{k,m}=Path;
    if Path(end)==E
          PL(k,m)=PLkm;
         if PLkm<minkl
             mink=k;minl=m;minkl=PLkm;
         end
    else
          PL(k,m)=0;
    end
   end
  %% 第六步：更新信息素
  Delta_Tau=zeros(N,N);%更新量初始化
  for m=1:M
     if PL(k,m) 
        ROUT=ROUTES{k,m};
        TS=length(ROUT)-1;%跳数
        PL_km=PL(k,m);
        for s=1:TS
          x=ROUT(s);
          y=ROUT(s+1);
          Delta_Tau(x,y)=Delta_Tau(x,y)+Q/PL_km;
          Delta_Tau(y,x)=Delta_Tau(y,x)+Q/PL_km;
        end
     end
  end
       Tau=(1-Rho).*Tau+Delta_Tau;%信息素挥发一部分，新增加一部分
  end
%%---------------------------绘图--------------------------------
  plotif=1;%是否绘图的控制参数
  if plotif==1
%绘收敛曲线
      minPL=zeros(K);
      for i=1:K
          PLK=PL(i,:);
          Nonzero=find(PLK);
          PLKPLK=PLK(Nonzero);
          minPL(i)=min(PLKPLK);
      end
   figure(1)
   plot(minPL);
   text(K/2,minPL(K/2)+0.5, 'rebot1');
   hold on
   grid on
      title('收敛曲线');   %title('收敛曲线（最小路径长度）');
   xlabel('迭代次数');
      ylabel('适应度值'); %ylabel('路径长度');
   %绘爬行图
   figure(2)
   axis([0,MM,0,MM])
   for i=1:MM
       for j=1:MM
            if G(i,j)==1
                 x1=j-1;y1=MM-i;
                 x2=j;y2=MM-i;
                 x3=j;y3=MM-i+1;
                 x4=j-1;y4=MM-i+1;
                fill([x1,x2,x3,x4],[y1,y2,y3,y4],[0.2,0.2,0.2]);
                hold on
            else
                 x1=j-1;y1=MM-i;
                 x2=j;y2=MM-i;
                 x3=j;y3=MM-i+1;
                 x4=j-1;y4=MM-i+1;
                 fill([x1,x2,x3,x4],[y1,y2,y3,y4],[1,1,1]);
                 hold on
            end
       end
   end
   hold on
   ROUT=ROUTES{mink,minl};
   LENROUT=length(ROUT);
   Rx=ROUT;
   Ry=ROUT;
   for ii=1:LENROUT
        Rx(ii)=a*(mod(ROUT(ii),MM)-0.5);
        if Rx(ii)==-0.5
             Rx(ii)=MM-0.5;
        end
        Ry(ii)=a*(MM+0.5-ceil(ROUT(ii)/MM));
   end
   plot(Rx,Ry)
   plot(Rx(1),Ry(1),'r.','MarkerSize',50);
   plot(Rx(LENROUT),Ry(LENROUT),'g.','MarkerSize',50);
   text(Rx(1)-0.2,Ry(1)+0.3,'rebot1 start');
   text(Rx(LENROUT)-1.5,Ry(LENROUT)-0.3,'rebot1 end');
   
end
plotif2=0;%绘各代蚂蚁爬行图
if plotif2==1
    figure(3)
    axis([0,MM,0,MM])
    for i=1:MM
         for j=1:MM
              if G(i,j)==1
                   x1=j-1;y1=MM-i;
                   x2=j;y2=MM-i;
                   x3=j;y3=MM-i+1;
                   x4=j-1;y4=MM-i+1;
                   fill([x1,x2,x3,x4],[y1,y2,y3,y4],[0.2,0.2,0.2]);
                   hold on
              else
                   x1=j-1;y1=MM-i;
                   x2=j;y2=MM-i;
                   x3=j;y3=MM-i+1;
                   x4=j-1;y4=MM-i+1;
                   fill([x1,x2,x3,x4],[y1,y2,y3,y4],[1,1,1]);
                   hold on
              end
         end
    end
    for k=1:K
         PLK=PL(k,:);
         minPLK=min(PLK);
         pos=find(PLK==minPLK);
         m=pos(1);
         ROUT=ROUTES{k,m};
         LENROUT=length(ROUT);
         Rx=ROUT;
         Ry=ROUT;
         for ii=1:LENROUT
             Rx(ii)=a*(mod(ROUT(ii),MM)-0.5);
             if Rx(ii)==-0.5
                  Rx(ii)=MM-0.5;
             end
             Ry(ii)=a*(MM+0.5-ceil(ROUT(ii)/MM));
         end
        plot(Rx,Ry)
        hold on
      end
end

%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%  2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
G=[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 1 1 1 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0;
   0 0 1 1 1 0 0 1 0 0 0 0 0 0 1 1 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0;
   0 0 0 1 1 0 0 1 0 0 0 0 0 0 0 0 0 0 1 1;
   0 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0;
   1 1 1 1 1 0 0 1 1 1 0 0 0 0 1 0 0 0 0 0;
   1 1 1 1 1 0 0 1 1 1 0 0 0 1 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 1 0 0 0 0 0 1 1 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0;
   0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 1 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0;];
MM=size(G,1);      % G 地形图为01矩阵，如果为1表示障碍物
Tau=ones(MM*MM,MM*MM);% Tau 初始信息素矩阵（认为前面的觅食活动中有残留的信息素）
Tau=8.*Tau;
K=100;          % K 迭代次数（指蚂蚁出动多少波）
M=50;         % M 蚂蚁个数（每一波蚂蚁有多少个）
S=MM ;         % S 起始点（最短路径的起始点）
E=MM*10+1;         % E 终止点（最短路径的目的点）
Alpha=1;    % Alpha 表征信息素重要程度的参数
Beta=7;     % Beta 表征启发式因子重要程度的参数
Rho=0.3 ;     % Rho 信息素蒸发系数
Q=1;          % Q 信息素增加强度系数
minkl=inf;
mink=0;
minl=0;

D=G2D(G);
N=size(D,1);%N表示问题的规模（象素个数）
a=1;%小方格象素的边长
Ex=a*(0.5);%终止点横坐标
if Ex==-0.5
Ex=MM-0.5;
end
Ey=a*(5.5);%终止点纵坐标
Eta=zeros(N);%启发式信息，取为至目标点的直线距离的倒数
%下面构造启发式信息矩阵
for i=1:N
 ix=a*(mod(i,MM)-0.5);
   if ix==-0.5
   ix=MM-0.5;
   end
iy=a*(MM+0.5-ceil(i/MM)); 
   if i~=E
   Eta(i)=1/((ix-Ex)^2+(iy-Ey)^2)^0.5;
   else
   Eta(i)=100;
   end
end
ROUTES=cell(K,M);%用细胞结构存储每一代的每一只蚂蚁的爬行路线
PL=zeros(K,M);%用矩阵存储每一代的每一只蚂蚁的爬行路线长度
%% -----------启动K轮蚂蚁觅食活动，每轮派出M只蚂蚁--------------------
for k=1:K
   for m=1:M
    %% 第一步：状态初始化
    W=S;%当前节点初始化为起始点
    Path=S;%爬行路线初始化
    PLkm=0;%爬行路线长度初始化
    TABUkm=ones(N);%禁忌表初始化
    TABUkm(S)=0;%已经在初始点了，因此要排除
    DD=D;%邻接矩阵初始化
    %% 第二步：下一步可以前往的节点
    DW=DD(W,:);
    DW1=find(DW);
    for j=1:length(DW1)
        if TABUkm(DW1(j))==0
            DW(DW1(j))=0;
        end
    end
    LJD=find(DW);
    Len_LJD=length(LJD)%可选节点的个数

   %% 觅食停止条件：蚂蚁未遇到食物或者陷入死胡同
    while W~=E&&Len_LJD>=1
    %% 第三步：转轮赌法选择下一步怎么走
         PP=zeros(Len_LJD);
         for i=1:Len_LJD
              PP(i)=(Tau(W,LJD(i))^Alpha)*((Eta(LJD(i)))^Beta);
         end
        sumpp=sum(PP);
        PP=PP/sumpp;%建立概率分布
        PP
        Pcum(1)=PP(1);
        for i=2:Len_LJD
              Pcum(i)=Pcum(i-1)+PP(i);
        end
      %  Select=zeros(Len_LJD);
      ra=rand
      Pcum
        Select=find(Pcum>=ra);
        Select
        to_visit=LJD(Select(1));
      %% 第四步：状态更新和记录
        Path=[Path,to_visit];%路径增加
        PLkm=PLkm+DD(W,to_visit);%路径长度增加
        W=to_visit;%蚂蚁移到下一个节点
        for kk=1:N
            if TABUkm(kk)==0
                 DD(W,kk)=0;
                 DD(kk,W)=0;
            end
        end
       TABUkm(W)=0;%已访问过的节点从禁忌表中删除
       DW=DD(W,:);
       DW1=find(DW);
       for j=1:length(DW1)
           if TABUkm(DW1(j))==0
               DW(j)=0;
           end
       end
       LJD=find(DW);
       Len_LJD=length(LJD);%可选节点的个数
     end
%% 第五步：记下每一代每一只蚂蚁的觅食路线和路线长度
    ROUTES{k,m}=Path;
    if Path(end)==E
          PL(k,m)=PLkm;
         if PLkm<minkl
             mink=k;minl=m;minkl=PLkm;
         end
    else
          PL(k,m)=0;
    end
   end
  %% 第六步：更新信息素
  Delta_Tau=zeros(N,N);%更新量初始化
  for m=1:M
     if PL(k,m) 
        ROUT=ROUTES{k,m};
        TS=length(ROUT)-1;%跳数
        PL_km=PL(k,m);
        for s=1:TS
          x=ROUT(s);
          y=ROUT(s+1);
          Delta_Tau(x,y)=Delta_Tau(x,y)+Q/PL_km;
          Delta_Tau(y,x)=Delta_Tau(y,x)+Q/PL_km;
        end
     end
  end
       Tau=(1-Rho).*Tau+Delta_Tau;%信息素挥发一部分，新增加一部分
  end
%% ---------------------------绘图--------------------------------
  plotif=1;%是否绘图的控制参数
  if plotif==1
%绘收敛曲线
      minPL=zeros(K);
      for i=1:K
          PLK=PL(i,:);
          Nonzero=find(PLK);
          PLKPLK=PLK(Nonzero);
          minPL(i)=min(PLKPLK);
         
      end
   figure(1)
   plot(minPL);
    text(K/2,minPL(K/2)-0.5, 'rebot2');
   hold on
   grid on
   title('收敛曲线（最小路径长度）');
   xlabel('迭代次数');
   ylabel('路径长度');
   %绘爬行图
   figure(2)
   hold on
   ROUT=ROUTES{mink,minl};
   LENROUT=length(ROUT);
   Rx=ROUT;
   Ry=ROUT;
   for ii=1:LENROUT
        Rx(ii)=a*(mod(ROUT(ii),MM)-0.5);
        if Rx(ii)==-0.5
             Rx(ii)=MM-0.5;
        end
        Ry(ii)=a*(MM+0.5-ceil(ROUT(ii)/MM));
   end
   plot(Rx,Ry)
    plot(Rx(1),Ry(1),'r.','MarkerSize',50);
    plot(Rx(LENROUT),Ry(LENROUT),'g.','MarkerSize',50);
   text(Rx(1)-0.2,Ry(1)+0.3,'rebot2 start');
   text(Rx(LENROUT)-0.2,Ry(LENROUT)-0.3,'rebot2 end');
   
end
plotif2=0;%绘各代蚂蚁爬行图
if plotif2==1
    figure(3)
    axis([0,MM,0,MM])
    for i=1:MM
         for j=1:MM
              if G(i,j)==1
                   x1=j-1;y1=MM-i;
                   x2=j;y2=MM-i;
                   x3=j;y3=MM-i+1;
                   x4=j-1;y4=MM-i+1;
                   fill([x1,x2,x3,x4],[y1,y2,y3,y4],[0.2,0.2,0.2]);
                   hold on
              else
                   x1=j-1;y1=MM-i;
                   x2=j;y2=MM-i;
                   x3=j;y3=MM-i+1;
                   x4=j-1;y4=MM-i+1;
                   fill([x1,x2,x3,x4],[y1,y2,y3,y4],[1,1,1]);
                   hold on
              end
         end
    end
    for k=1:K
         PLK=PL(k,:);
         minPLK=min(PLK);
         pos=find(PLK==minPLK);
         m=pos(1);
         ROUT=ROUTES{k,m};
         LENROUT=length(ROUT);
         Rx=ROUT;
         Ry=ROUT;
         for ii=1:LENROUT
             Rx(ii)=a*(mod(ROUT(ii),MM)-0.5);
             if Rx(ii)==-0.5
                  Rx(ii)=MM-0.5;
             end
             Ry(ii)=a*(MM+0.5-ceil(ROUT(ii)/MM));
         end
        plot(Rx,Ry)
        hold on
      end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%   3
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
G=[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 1 1 1 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0;
   0 0 1 1 1 0 0 1 0 0 0 0 0 0 1 1 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0;
   0 0 0 1 1 0 0 1 0 0 0 0 0 0 0 0 0 0 1 1;
   0 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0;
   1 1 1 1 1 0 0 1 1 1 0 0 0 0 1 0 0 0 0 0;
   1 1 1 1 1 0 0 1 1 1 0 0 0 1 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 1 0 0 0 0 0 1 1 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0;
   0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 0 1 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0;
   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0;];
MM=size(G,1);      % G 地形图为01矩阵，如果为1表示障碍物
Tau=ones(MM*MM,MM*MM);% Tau 初始信息素矩阵（认为前面的觅食活动中有残留的信息素）
Tau=8.*Tau;
K=100;          % K 迭代次数（指蚂蚁出动多少波）
M=50;         % M 蚂蚁个数（每一波蚂蚁有多少个）
S=MM*(MM-1)+1 ;         % S 起始点（最短路径的起始点）
E=MM;         % E 终止点（最短路径的目的点）
Alpha=1;    % Alpha 表征信息素重要程度的参数
Beta=7;     % Beta 表征启发式因子重要程度的参数
Rho=0.3 ;     % Rho 信息素蒸发系数
Q=1;          % Q 信息素增加强度系数
minkl=inf;
mink=0;
minl=0;

D=G2D(G);
N=size(D,1);%N表示问题的规模（象素个数）
a=1;%小方格象素的边长
Ex=a*(mod(E,MM)-0.5);%终止点横坐标
if Ex==-0.5
Ex=MM-0.5;
end
Ey=a*(MM-0.5);%终止点纵坐标
Eta=zeros(N);%启发式信息，取为至目标点的直线距离的倒数
%下面构造启发式信息矩阵
for i=1:N
 ix=a*(mod(i,MM)-0.5);
   if ix==-0.5
   ix=MM-0.5;
   end
iy=a*(MM+0.5-ceil(i/MM)); 
   if i~=E
   Eta(i)=1/((ix-Ex)^2+(iy-Ey)^2)^0.5;
   else
   Eta(i)=100;
   end
end
ROUTES=cell(K,M);%用细胞结构存储每一代的每一只蚂蚁的爬行路线
PL=zeros(K,M);%用矩阵存储每一代的每一只蚂蚁的爬行路线长度
%% -----------启动K轮蚂蚁觅食活动，每轮派出M只蚂蚁--------------------
for k=1:K
   for m=1:M
    %% 第一步：状态初始化
    W=S;%当前节点初始化为起始点
    Path=S;%爬行路线初始化
    PLkm=0;%爬行路线长度初始化
    TABUkm=ones(N);%禁忌表初始化
    TABUkm(S)=0;%已经在初始点了，因此要排除
    DD=D;%邻接矩阵初始化
    %% 第二步：下一步可以前往的节点
    DW=DD(W,:);
    DW1=find(DW);
    for j=1:length(DW1)
        if TABUkm(DW1(j))==0
            DW(DW1(j))=0;
        end
    end
    LJD=find(DW);
    Len_LJD=length(LJD);%可选节点的个数

   %% 觅食停止条件：蚂蚁未遇到食物或者陷入死胡同
    while W~=E&&Len_LJD>=1
    %% 第三步：转轮赌法选择下一步怎么走
         PP=zeros(Len_LJD);
         for i=1:Len_LJD
              PP(i)=(Tau(W,LJD(i))^Alpha)*((Eta(LJD(i)))^Beta);
         end
        sumpp=sum(PP);
        PP=PP/sumpp;%建立概率分布
        Pcum(1)=PP(1);
        for i=2:Len_LJD
              Pcum(i)=Pcum(i-1)+PP(i);
        end
        Select=find(Pcum>=rand);
        to_visit=LJD(Select(1));
      %% 第四步：状态更新和记录
        Path=[Path,to_visit];%路径增加
        PLkm=PLkm+DD(W,to_visit);%路径长度增加
        W=to_visit;%蚂蚁移到下一个节点
        for kk=1:N
            if TABUkm(kk)==0
                 DD(W,kk)=0;
                 DD(kk,W)=0;
            end
        end
       TABUkm(W)=0;%已访问过的节点从禁忌表中删除
       DW=DD(W,:);
       DW1=find(DW);
       for j=1:length(DW1)
           if TABUkm(DW1(j))==0
               DW(j)=0;
           end
       end
       LJD=find(DW);
       Len_LJD=length(LJD);%可选节点的个数
     end
%% 第五步：记下每一代每一只蚂蚁的觅食路线和路线长度
    ROUTES{k,m}=Path;
    if Path(end)==E
          PL(k,m)=PLkm;
         if PLkm<minkl
             mink=k;minl=m;minkl=PLkm;
         end
    else
          PL(k,m)=0;
    end
   end
  %% 第六步：更新信息素
  Delta_Tau=zeros(N,N);%更新量初始化
  for m=1:M
     if PL(k,m) 
        ROUT=ROUTES{k,m};
        TS=length(ROUT)-1;%跳数
        PL_km=PL(k,m);
        for s=1:TS
          x=ROUT(s);
          y=ROUT(s+1);
          Delta_Tau(x,y)=Delta_Tau(x,y)+Q/PL_km;
          Delta_Tau(y,x)=Delta_Tau(y,x)+Q/PL_km;
        end
     end
  end
       Tau=(1-Rho).*Tau+Delta_Tau;%信息素挥发一部分，新增加一部分
  end
%% ---------------------------绘图--------------------------------
  plotif=1;%是否绘图的控制参数
  if plotif==1
%绘收敛曲线
      minPL=zeros(K);
      for i=1:K
          PLK=PL(i,:);
          Nonzero=find(PLK);
          PLKPLK=PLK(Nonzero);
          minPL(i)=min(PLKPLK);
      end
   figure(1)
   plot(minPL);
    text(K/2,minPL(K/2)-0.5, 'rebot3');
   hold on
   grid on
  title('收敛曲线');% title('收敛曲线（最小路径长度）');
   xlabel('迭代次数');
   ylabel('适应度值'); %ylabel('路径长度');
   %绘爬行图
   figure(2)
 
   hold on
   ROUT=ROUTES{mink,minl};
   LENROUT=length(ROUT);
   Rx=ROUT;
   Ry=ROUT;
   for ii=1:LENROUT
        Rx(ii)=a*(mod(ROUT(ii),MM)-0.5);
        if Rx(ii)==-0.5
             Rx(ii)=MM-0.5;
        end
        Ry(ii)=a*(MM+0.5-ceil(ROUT(ii)/MM));
   end
   plot(Rx,Ry)
   plot(Rx(1),Ry(1),'r.','MarkerSize',50);
   plot(Rx(LENROUT),Ry(LENROUT),'g.','MarkerSize',50);
   text(Rx(1)-0.2,Ry(1)+0.3,'rebot3 start');
   text(Rx(LENROUT)-1.5,Ry(LENROUT)+0.3,'rebot3 end');
   
end
plotif2=0;%绘各代蚂蚁爬行图
if plotif2==1
    figure(3)
    axis([0,MM,0,MM])
    for i=1:MM
         for j=1:MM
              if G(i,j)==1
                   x1=j-1;y1=MM-i;
                   x2=j;y2=MM-i;
                   x3=j;y3=MM-i+1;
                   x4=j-1;y4=MM-i+1;
                   fill([x1,x2,x3,x4],[y1,y2,y3,y4],[0.1,0.1,0.1]);
                   hold on
              else
                   x1=j-1;y1=MM-i;
                   x2=j;y2=MM-i;
                   x3=j;y3=MM-i+1;
                   x4=j-1;y4=MM-i+1;
                   fill([x1,x2,x3,x4],[y1,y2,y3,y4],[1,1,1]);
                   hold on
              end
         end
    end
    for k=1:K
         PLK=PL(k,:);
         minPLK=min(PLK);
         pos=find(PLK==minPLK);
         m=pos(1);
         ROUT=ROUTES{k,m};
         LENROUT=length(ROUT);
         Rx=ROUT;
         Ry=ROUT;
         for ii=1:LENROUT
             Rx(ii)=a*(mod(ROUT(ii),MM)-0.5);
             if Rx(ii)==-0.5
                  Rx(ii)=MM-0.5;
             end
             Ry(ii)=a*(MM+0.5-ceil(ROUT(ii)/MM));
         end
        plot(Rx,Ry)
        hold on
      end
end

function D=G2D(G)
l=size(G,1);
D=zeros(l*l,l*l); 
for i=1:l
    for j=1:l
        if G(i,j)==0
            for m=1:l
                for n=1:l
                    if G(m,n)==0
                        im=abs(i-m);jn=abs(j-n);
                        if im+jn==1||(im==1&&jn==1)
                        D((i-1)*l+j,(m-1)*l+n)=(im+jn)^0.5;
                        end
                    end
                end
            end
        end
    end
end

本文章已经生成可运行项目