TSP是什么？TSP全称Travelling salesman problem。中文名：旅行商问题。就是模拟退火中讲到的14个城市之间巡回旅行，求路径最短的问题。

为什么偏偏找他呢？因为这是一个衡量算法的“金标准”，而且他简单，用来介绍一个算法的基本思想再好不过。

下面我们用同样的方法看看蚁群算法。

参考资料：https://baike.baidu.com/item/TSP/2905216#viewPageContent

我尽量用白话去解释我所理解的蚁群算法：

这个是参数和概率计算公式的说明：（这个先大致看看，留个印象）

m = 200;                              % 蚂蚁数量

M=31                                %31个城市（大致模拟了中国省会城市的分布）

alpha =4;                            % 信息素重要程度因子

beta =1;                             % 启发函数重要程度因子

rho（就是ρ） = 0.1;                  % 信息素挥发因子

Q = 1;                               % 常系数

Eta（就是τ）= 1./D;                 % 启发函数

Tau（就是η） = ones(n,n);            % 信息素矩阵

Table = zeros(m,n);                    % 路径记录表

iter = 1;                             % 迭代次数初值

iter_max = 300;                       % 最大迭代次数

Route_best = zeros(iter_max,n);          % 各代最佳路径      

Length_best = zeros(iter_max,1);         % 各代最佳路径的长度 

Length_ave = zeros(iter_max,1);          % 各代路径的平均长度

以上是参数和公式，以下是网上介绍的过程，大致了解一下：

以上图片来源于：http://www.dataguru.cn/article-11200-1.html

我的理解：

蚂蚁在自然界里，知道怎么去选择路径，很大程度上是因为之前有蚂蚁在路径上留下了信息素去指导后面的蚂蚁运动，信息素在算法里就是τ。当然信息素作为一种化学物质，他是会挥发的，所以算法里面与之对应有一个ρ。

这个算法还加了一点什么呢？就是η，这是两个城市之间距离的倒数，意味着蚂蚁在城市A选择下一个城市的时候，会更有可能选择离城市A更近的城市。

开始，有200（m=200）只蚂蚁，但是这些蚂蚁，寿命很短，只够他们不重复爬一圈所有的城市，回到起点就死了，他们一生追求就是“世界那么大，我想去看看”，这一辈子就只想把31个城市看完然后回归故土。这些蚂蚁呢，我们一开始要为他们在这31(M=31)个城市中随机分布起点，也就是相对于出生了。之后他们一辈子还有三十个城市要走。然后，我们一次性要把200只蚂蚁接下来所有的行程安排好。怎么做呢？比如有一只蚂蚁1(1是该蚂蚁的序号,就是公式中的k)，他被分到了第6个城市出生，接下来怎么去安排他的第二个“旅游地点”呢？这就要请出我们的转移概率公式啦，就是那个在第一章图片中看上去很可怕的p（还有i,j作为下标，k作为上标，t作为输入参数），但是不用怕，i其实就是当前城市序号，j是下一步城市序号，k是蚂蚁序号，t是代数序号，就是算法迭代的代数。咱们只要把分子弄懂了，分母不过就是分子加起来的一个和。

我们要记得，城市不能重复，所以图片中才会有判断j从属关系的一个条件，如果j之前没有去过，就可以用公式算出概率，去参与概率竞争，否则概率就为0，就不去了。

重点是：分子是什么意思呢？一个是τ,一个是η，τ一开始全部是1，所有路径上的信息素都是1，在之后由于各个蚂蚁爬动产生的影响，以及信息素在每一次迭代之后都会挥发的效果，信息素在各个城市之间会差别越来越大，这是显然的。有一些路线由于很少有蚂蚁去爬动，导致挥发越来越多，渐渐就没有信息素了。其他那些有蚂蚁去爬动的路线上由于有信息素的补充，信息素能够维持在一个水平。信息素越大，选择的概率也就越大。

η是两个城市之间距离的倒数，也就是有更大的概率去到那些离自己更近的城市。这个η称为启发式因子。η=1/d称为启发函数。

α，β是信息素和启发式因子的权重因子，值越大在概率公式中也越重要。

回归正题，蚂蚁到了第6个城市后，会得到的一个禁忌列表（Taboo），就是他以前走过的城市的列表，总的城市列表减去这个禁忌列表就会得到允许列表即（Allow列表），然后就根据允许列表中的所有可能（比如6-7,6-9,6-25,6-10等等有三十个“数对”）产生他们所对应的η，τ，然后算出他们mij=（τ^α）*（η^β），上面举了四个“数对”,就有四个不同的m67，m69,m625,m610然后每一个的概率pij=mij/,分母是所有“可能”的m之和，因为Taboo列表里面m为0。这样第1只蚂蚁走完了第一个城市，又产生第二个城市的概率，那么电脑根据概率随机生成下一个城市，又用同样的方法随机生成第三，第四第五直到第三十一个城市。紧接着第2只蚂蚁重复上述过程，直到全部200只蚂蚁全部安排完31个城市旅行顺序。

这里需要注意的是，在同一批次，即同一代的蚂蚁，互相之间是不会受到之前蚂蚁留下信息素浓度的影响，他们计算概率时所用到的τ是上一代蚂蚁留下的信息素。（第一代蚂蚁用的信息素用的是初始值，初始值是所有城市之间信息素浓度一样。）

接下来是更新信息素了。更新信息素学术界有三种方法

第一种是：ant cycle system

△τ=Q/L

第二种是：ant quantity system

△τ=Q/dij

第三种是：ant density system

△τ=Q

首先我们脑海里要先有一个31*31的矩阵（矩阵是关于对角线对称的，对角线上我们取0.0001，不取0因为到时候这个矩阵要全部取倒数，0.0001还是0.1都没有关系，反正这个数据不会用到，只是为了分母不为零），用来描述这个三十一城市互相之间的信息素浓度。一开始都是1。上面所说的第一代蚂蚁爬过之后信息素会怎么变化呢？

首先，所有τ=τ*（1-rho（就是ρ）），先挥发掉一部分。挥发操作的频率是每一代完成之后进行一次。这是环境因素使信息素变少。

之后再讲生物因素，这里我介绍第一种方法ant cycle system。“在每一次迭代之后”（注意频率），对于蚂蚁1来说，算出他旅行的总长度L，这是分母。然后分子Q在程序中永远是一个常数，用来描述蚂蚁释放信息素的能力，然后相除算出△τ，之后在每一个通过的路径都加上这个△τ，于是我们脑海中的那个矩阵就会有62项加上了这个△τ（因为这个矩阵是关于对角线对称的）。比如6-7-1-2-3-。。。。。（省略）-25-30-6，那么6-7,7-1,1-2等等等一直到30-6都会加上这个相同的△τ，那么紧接着再更新蚂蚁2走过的所有路径的τ，加上蚂蚁2算出来的△τ（与蚂蚁1的△τ区分开来，不同蚂蚁的L不同，△τ自然也不同），紧接着蚂蚁3,4,5。。。200,的信息素全部更新完成。

至此，第一代完成，第一代蚂蚁全部死亡。

程序随机产生第二代200只蚂蚁，随机初始出生地，再重复上述规则。完成

直到达到最大迭代次数。iter_max。（iteration 迭代 ）

下面再说说第二种生物因素更新规则，ant quantity system,我们注意到在第一种方法中，对于同一只蚂蚁，所有的路径上△τ是一样的，△τ的计算频率是每一代计算一次，第二种规则就不同了，△τ的分母是两个城市之间的距离，所以是每经过一个城市就计算一次，对于同一个蚂蚁，不同“城市数对”之间的信息素增量△τ是不同的。每一个城市数对单独计算信息素增量，再在矩阵中单独加上去，规则一计算了200个△τ，规则一计算了6200个△τ。

第三种规则，△τ只有一种，对于所有蚂蚁的所有路径都是一模一样的Q常数。

还说一下，Q越大，那么没有蚂蚁爬过的路径和有蚂蚁爬过的路径之间的信息素差距就越大，也就是说这个时候蚂蚁爬过起到的作用就越大，相当于在增大α。换句话说α起到的作用也是使路径有没有蚂蚁爬过所产生的影响越来越大。

Q,α太大会使得早期一些碰巧被蚂蚁爬过的路径迅速占据大部分的概率，使得模型快速局部收敛，削弱了模型的全局搜索能力。

代码如下：参考自《matlab在数学建模中的应用》

tic%计时

%% I. 清空环境变量

%clear all

%clc

% data2=0;

for circle=11:20%将每一次的一些数据专门保存在data2矩阵里面

%% II. 导入数据

load citys_data.mat

%##

1304         2312

3639         1315

4177         2244

3712         1399

3488         1535

3326         1556

3238         1229

4196         1004

4312         790

4386         570

3007         1970

2562         1756

2788         1491

2381         1676

1332         695

3715         1678

3918         2179

4061         2370

3780         2212

3676         2578

4029         2838

4263         2931

3429         1908

3507         2367

3394         2643

3439         3201

2935         3240

3140         3550

2545         2357

2778         2826

2370         2975

%##

%% III. 计算城市间相互距离

n = size(citys,1);%n是城市的个数

D = zeros(n,n);

 for i = 1:n

    for j = 1:n

        if i ~= j

            D(i,j) = sqrt(sum((citys(i,:) – citys(j,:)).^2));

        else

            D(i,j) = 1e-4;     

        end

    end   

 end

%% IV. 初始化参数

m = 100;                              % 蚂蚁数量

alpha =4;                           % 信息素重要程度因子

beta =2;                            % 启发函数重要程度因子

rho = 0.1;                           % 信息素挥发因子

Q = 1;                               % 常系数

Eta = 1./D;                          % 启发函数

Tau = ones(n,n);                     % 信息素矩阵

Table = zeros(m,n);                  % 路径记录表

iter = 1;                            % 迭代次数初值

iter_max = 120;                      % 最大迭代次数

Route_best = zeros(iter_max,n);      % 各代最佳路径      

Length_best = zeros(iter_max,1);     % 各代最佳路径的长度 

Length_ave = zeros(iter_max,1);      % 各代路径的平均长度 

best_length =zeros(iter_max,1);

iter_limit=1;

%% V. 迭代寻找最佳路径

 while iter <= iter_max

     % 随机产生各个蚂蚁的起点城市

      start = zeros(m,1);

      for i = 1:m

          temp = randperm(n);

          start(i) = temp(1);

      end

      Table(:,1) = start;

      citys_index = 1:n;

      % 逐个蚂蚁路径选择

      for i = 1:m

          % 逐个城市路径选择

         for j = 2:n

             tabu = Table(i,1:(j – 1));           % 已访问的城市集合(禁忌表)

             allow_index = ~ismember(citys_index,tabu);

             allow = citys_index(allow_index);  % 待访问的城市集合

             P = allow;

             % 计算城市间转移概率

             for k = 1:length(allow)

                 P(k) = Tau(tabu(end),allow(k))^alpha * Eta(tabu(end),allow(k))^beta;

             end

             P = P/sum(P);

             % 轮盘赌法选择下一个访问城市

             Pc = cumsum(P);    

            target_index = find(Pc >= rand);

            target = allow(target_index(1));

            Table(i,j) = target;

         end

      end

      % 计算各个蚂蚁的路径距离

      Length = zeros(m,1);

      for i = 1:m

          Route = Table(i,:);

          for j = 1:(n – 1)

              Length(i) = Length(i) + D(Route(j),Route(j + 1));

          end

          Length(i) = Length(i) + D(Route(n),Route(1));

      end

      % 计算最短路径距离及平均距离

      if iter == 1

          [min_Length,min_index] = min(Length);

          Length_best(iter) = min_Length; 

          Length_ave(iter) = mean(Length);

          Route_best(iter,:) = Table(min_index,:);

          best_length = min_Length;

      else

          [min_Length,min_index] = min(Length);

          Length_best(iter) = min(Length_best(iter – 1),min_Length);

          Length_ave(iter) = mean(Length);

          best_length(iter)=min_Length;

          if Length_best(iter) == min_Length

              Route_best(iter,:) = Table(min_index,:);

              iter_limit=iter;

          else

              Route_best(iter,:) = Route_best((iter-1),:);

          end

      end

      % 更新信息素

      Delta_Tau = zeros(n,n);

      % 逐个蚂蚁计算

      for i = 1:m

          % 逐个城市计算

          for j = 1:(n – 1)

              Delta_Tau(Table(i,j),Table(i,j+1)) = Delta_Tau(Table(i,j),Table(i,j+1)) + Q/Length(i);

          end

          Delta_Tau(Table(i,n),Table(i,1)) = Delta_Tau(Table(i,n),Table(i,1)) + Q/Length(i);

      end

      Tau = (1-rho) * Tau + Delta_Tau;

    % 迭代次数加1，清空路径记录表

    iter = iter + 1;

    Table = zeros(m,n);

 end

toc

disp([‘运行时间’,num2str(toc)])

%% VI. 结果显示

[Shortest_Length,index] = min(Length_best);

Shortest_Route = Route_best(index,:);

disp([‘最短距离:’ num2str(Shortest_Length)]);

disp([‘最短路径:’ num2str([Shortest_Route Shortest_Route(1)])]);

%% VII. 绘图

figure(1)

plot([citys(Shortest_Route,1);citys(Shortest_Route(1),1)],…

     [citys(Shortest_Route,2);citys(Shortest_Route(1),2)],’o-‘);

grid on

 for i = 1:size(citys,1)

    text(citys(i,1),citys(i,2),[‘   ‘ num2str(i)]);

 end

text(citys(Shortest_Route(1),1),citys(Shortest_Route(1),2),’       起点’);

text(citys(Shortest_Route(end),1),citys(Shortest_Route(end),2),’       终点’);

xlabel(‘城市位置横坐标’)

ylabel(‘城市位置纵坐标’)

title([‘蚁群算法优化路径(最短距离:’ num2str(Shortest_Length) ‘)’])

figure(2)

plot(1:iter_max,Length_best,’b’,1:iter_max,Length_ave,’r:’,1:iter_max,best_length,’*-‘)

legend(‘最短距离’,’平均距离’,’每一代最短距离’)

xlabel(‘迭代次数’)

ylabel(‘距离’)

title(‘各代最短距离与平均距离对比’)

text(iter_limit,best_length(iter_limit)-200,[‘最后一次收敛点:’,num2str(iter_limit),’   ‘,num2str(best_length(iter_limit))]);

text(80,best_length(end)+400,[‘收敛距离’,num2str(best_length(end))]);

differ=best_length(end)-best_length(iter_limit);

text(50,best_length(end)+2000,[[‘m:  ‘,num2str(m),’ ‘],[‘a: ‘,num2str(alpha),’ ‘],[‘b:  ‘,num2str(beta)]…

,[‘  differ: ‘,num2str(differ)]]);

data2(circle,[1,2,3,4,5,6,7,8])=[circle,m,alpha,beta,best_length(end),best_length(iter_limit),iter_limit,differ]

end

图像展示：

第一组（alpha=2,beta=1,2,3）

第二组：（alpha=4,beta=1,2,3）

第三组（alpha=6,beta=1,2,3）

第四组（alpha=8,beta=1,2,3）:

一些研究alpha,beta关系的实验数据：

☆平均数：

小结：当alpha增大的时候，全局搜索能力下降，局部搜索上升，收敛速度加快。刚开始几代随机性很强，那么先撒下信息素的路径就会由于alpha而产生巨大优势，从而以后很多代都会围绕着他们收敛，向外发散能力减弱，能很快收敛。收敛代数会越来越小，如数据。

在同样alpha的情况下，增大beta可以使模型加快收敛，因为alpha和beta的扩大都会增加数据差异所带来的收益，差一点点可能就在概率上差很多，那么有一些路径信息素和η总体来讲都不错的话，就会在巨大的alpha和beta下很快拉开差距，从而加快收敛。

第三，图像中还反映了一个现象，最优路径往往在50代左右就出现， 但是后面由于波动，最终收敛于一个较大距离的路径。为什么？这正体现了概率是由信息素和η两方面去决定的，我们反过来想想，当一个路径最短时，毫无以为信息素浓度浓度是最大的，但是你得有蚂蚁去爬啊，蚂蚁选择你这条路径一是由于信息素，二是他与下一个城市的距离，我们反过来想想，总路径最小，等不等价于每一次选择城市时都是最近的呢？显然不一定吧，那么η产生的概率就会夺走那些最短路径的优势，而有些“较短路径”，虽然不是最短的，信息素分布也还比较多，但是总距离上的牺牲换来的却是η从总体来讲的提高。所以笑到最后的，一定是那些“文武双全”的路径，两方面制衡下达到了一个比较好的值，从而称霸。

第四，我们再去看看第一组的图像。很明显，当beta=1时，最终收敛的距离与最优距离相差无几，但是随着beta增大，差距越来越大，因为beta=1时，η占的权重比较小，这个时候信息素是决定性因素，那么随后收敛的时候偏来偏去，总不会离最优距离太远，最优距离的优势是很明显的，但是beta增大，最优距离的信息素优势被“冲淡”了，就偏得越来越远。这佐证了第三点。

第五，再从数据的波动性看第一组的图像。达到最优距离点之后，beta=1时，波动比较缓慢，几乎很平滑。beta=2，波动明显剧烈了一点，beta=3,波动更剧烈了，这就说明beta破坏了最短距离所建造的“平衡体系”，向外发散的可能性大大增加，越大，效果越明显。

第六，增大m，蚂蚁的数量，会一定程度上增加全局搜索能力。样本多了嘛。

第七，很重要的一点，以上的大小都是“相对”的，并且很多都是基于实验结果的经验之谈，算法分析不好定量分析，大多数情况下只能定性分析。

第八，书上给出的参考值：m/M(城市数)=1.5，较稳妥，alpha属于[1,4],beta属于[3,4.5],rho属于[0.2,0.5],Q属于[10,1000]。

迭代次数选择，一般来讲，先取200，执行程序看算法收敛轨迹，再判断。

第九，蚁群算法对于TSP及其类似问题有很好的适应性，对于其他问题还有待考虑，有一定局限性，这个“城市”“距离”怎么抽象出来嫁接到另一个情况下的实体是关键问题。

语法总结：

##

1、[A,B]=max(P),当P为一列或一行的时候，对应的行数（A）,列数（B）会变成值。因为此时行数和列数肯定是1，没有必要再说行数和列数了。

##

2、mean()函数总结

参考链接：https://zhidao.baidu.com/question/51932148.html

matlab中的mean函数函数功能是求数组的平均数或者均值。
使用方法如下：
M = mean(A)
返回沿数组中不同维的元素的平均值。
如果A是一个向量，mean(A)返回A中元素的平均值。
如果A是一个矩阵，mean(A)将其中的各列视为向量，把矩阵中的每列看成一个向量，返
M = mean(A,dim)
返回A中沿着标量dim指定的维数上的元素的平均值。对于矩阵，mean(A,2)就是包含每一行的平均值的列向量。
比如:
A = [1 2 3; 3 3 6; 4 6 8; 4 7 7];
mean(A)：相当于mean(A,1),求每一列所有行的平均值。这个1表示是以行为对象去求平均值。
ans =
       3.0000 4.5000 6.0000
mean(A,2)：
ans =
       2.0000
       4.0000
       6.0000
       6.0000
 

##
3、在MATLAB中如何计算程序运行时间

参考链接：https://blog.csdn.net/colourful_sky/article/details/70437001

1、tic和toc组合

tic

%代码块

toc

%disp([‘运行时间: ‘,num2str(toc)]);

程序遇到tic时Matlab自动开始计时，运行到toc时自动计算此时与最近一次tic之间的时间

2、etime()与clock组合

t1=clock;

%代码块

t2=clock;

etime(t2,t1)

计算t1，t2之间的时间差，它是通过调用windows系统的时钟进行时间差计算得到运行时间的。

3、cputime函数

t0=cputime

%代码块

t1=cputime-t0

上以上三种方法，都是可以进行程序运行时间的计算，但是Matlab官方推荐使用tic，toc组合。但是使用tic/toc的时候一定要注意，toc计算的是与最后一次（即离它最近）运行的tic之间的时间。（和if else对应关系一样的）

##

4、matlab代码有“节”这个概念，调试的时候在编辑器中有运行当前节，运行下一节的选项。建立节的方法是用%%加换行的方式：

%%

正文

%%

正文

##

1. 如何更改图像坐标轴的属性

参考链接：https://jingyan.baidu.com/article/948f5924231210d80ff5f9db.html

先用get(gca)得到坐标图的各个属性再去设置，熟悉之后就可以不看了。

set(gca,’属性名’，’属性值’)

set(gca,’box’,’on’,’xlim’,[02*pi],’YDir’,’reverse’)后，图形变成下图所示，出现坐标边界（box），x轴显示坐标范围缩小（xlim）,y轴方向反转（ydir）

m/M约为1.5时，比较合理

时间已过 25.728734 秒。

运行时间25.7289

最短距离:15601.9195

最短路径:29   1  15  14  12  13  11  23  16   5   6   7   2   4   8   9  10   3  18  17  19  24  25  20  21  22  26  28  27  30  31  29

##

1. Legend是图线的图标

x = 0:pi/50:2*pi;

y = sin(x);

y1=cos(x);

figure  

plot(x,y,’r’)

hold on 

plot(x,y1,’b’)

legend(‘a’,’b’);

这里的对应关系是用位置关系来传达的。即参数的先后关系。

1. 加了分号可以抑制变量在窗口栏的输出

x = 0:pi/50:2*pi;

x = 0:pi/50:2*pi （会输出x）

8、轮盘赌法是个很常见很常用的方法。

先算出各个情况的概率，然后依次相加。

0.00144420638923246  0.0122002605472786    0.0681045118239203         0.0814283585181726    0.0977425110049918    0.111676411900531       0.118714416724746         0.124147333734136       0.127994554381832       0.130889461704665       0.144406233096476         0.149649194261840       0.155546015598805       0.159497127626474       0.160566366422886         0.191211982730686       0.631695720512809       0.717028980967632       0.778285709428712         0.797824613344900       0.809645018666816       0.850774702749089       0.940758458789103         0.967250090309183       0.975353382420096       0.980361409312193       0.984392722029262         0.990128055314803       0.996549676670898       1.00000000000000

这是其中某只蚂蚁的三十个城市的“轮盘赌”概率。第一个数据为p1，第二个数据为p1+p2，第三个数据为p1+p2+p3，。。。所以最后一项为1。第一个区间是0到第一个数据，第二个区间是第一个数据到第二个数据，以此类推。

然后用rand产生随机数（0,1），随机数在第k个区间就去第k个城市。