计算智能课程设计-粒子群优化算法求解旅行商问题-Matlab实现.doc

摘要：TSP是一个典型的NPC问题。本文首先介绍旅行商问题和粒子群优化算法的基本概念。然后构造一种基于交换子和交换序[1]概念的粒子群优化算法，通过控制学习因子和、最大速度，尝试求解旅行商问题。本文以中国31个省会城市为例，通过MATLAB编程实施对旅行商问题的求解,得到了一定优化程度的路径，是粒子群优化算法在TSP问题中运用的一次大胆尝试。
关键字：TSP问题；粒子群优化算法；MATLAB；中国31个城市TSP。
粒子群优化算法求解旅行商问题

1. 引言
1. 旅行商问题的概述
旅行商问题，即TSP问题（Traveling Salesman Problem）又译为旅行推销员问题货郎担问题，是数学领域中著名问题之一。假设有一个旅行商人要拜访n个城市，他必须选择所要走的路径，路径的限制是每个城市只能拜访一次，而且最后要回到原来出发的城市。路径的选择目标是要求得的路径路程为所有路径之中的最小值。
TSP问题是一个组合优化问题，其描述非常简单，但最优化求解非常困难，若用穷举法搜索，对N个城市需要考虑N！种情况并两两对比，找出最优，其算法复杂性呈指数增长，即所谓“指数爆炸”。所以，寻求和研究TSP问题的有效启发式算法，是问题的关键。
2. 粒子群优化算法的概述
粒子群优化算法(Particle Swarm optimization,PSO)又翻译为粒子群算法、微粒群算法、或微粒群优化算法。是通过模拟鸟群觅食行为而发展起来的一种基于群体协作的随机搜索算法。通常认为它是群集智能 (Swarm intelligence, SI) 的一种。它可以被纳入多主体优化系统 (Multiagent Optimization System, MAOS). 粒子群优化算法是由Eberhart博士和Kennedy博士发明。
PSO模拟鸟群的捕食行为。一群鸟在随机搜索食物，在这个区域里只有一块食物。所有的鸟都不知道食物在那里。但是他们知道当前的位置离食物还有多远。那么找到食物的最优策略是什么呢。最简单有效的就是搜寻目前离食物最近的鸟的周围区域。
PSO从这种模型中得到启示并用于解决优化问题。PSO中，每个优化问题的解都是搜索空间中的一只鸟。我们称之为“粒子”。所有的粒子都有一个由被优化的函数决定的适应值(fitnessvalue)，每个粒子还有一个速度决定他们飞翔的方向和距离。然后粒子们就追随当前的最优粒子在解空间中搜索。
PSO初始化为一群随机粒子(随机解)，然后通过迭代找到最优解，在每一次叠代中，粒子通过跟踪两个“极值”来更新自己。第一个就是粒子本身所找到的最优解，这个解叫做个体极值pBest，另一个极值是整个种群目前找到的最优解，这个极值是全局极值gBest。另外也可以不用整个种群而只是用其中一部分最优粒子的邻居，那么在所有邻居中的极值就是局部极值。
3. 粒子群优化算法求解旅行商问题
旅行商问题是一个典型的、易于描述却难于处理的组合优化问题,并且是一个NP完全难题,其可能的路径数目与城市数目n是成指数型增长的,对n个城市而言,可能的路径总(n-1)!。随着n的增加,路径数将按数率急剧增长,即所谓的“指数爆炸”,所以一般很难精确地求出其最优解,因而寻找出其有效的近似求解算法就具有重要的理论意义。而粒子群优化算法是解决复杂问题的有效方法,自然也能应用于解决旅行商问题。
PSO模拟鸟群的捕食行为。一群鸟在随机搜索食物，在这个区域里只有一块食物。所有的鸟都不知道食物在那里。但是他们知道当前的位置离食物还有多远。那么找到食物的最优策略是什么呢。最简单有效的就是搜寻目前离食物最近的鸟的周围区域。
2. 粒子群算法的基本思想[2]
1. 粒子群优化算法的基本原理
一个由个粒子(Particle)组成的群体(Swarm)在维搜索空间中以一定的速度飞行，每个粒子在搜索时，考虑到了自己搜索到的的历史最好点和群体内(或领域内)其它粒子的最好点，在此基础上进行位置(状态、也就是解)的变化。
第个粒子的位置表示为：
第个粒子的速度表示为：，
第个粒子所经历的历史最好点表示为：
群体内(或领域内)所有粒子所经历过的最好的点表示为：。
一般来说，粒子的位置和速度都是在连续的实数空间内进行取值，粒子的位置和速度根据如下方程进行变化

其中，为惯性权重。和称为学习因子(Learning Factor)或加速系数(Acceleration Coefficient)，一般为正常数。学习因子使粒子具有自我总结和向群体中优秀个体学习的能力，从而向自己的历史最优点以及群体内或领域内的历史最优点靠近。和通常等于2。，，是在区间内均匀分布的伪随机数。粒子的速度被限制在一个最大的范围内。
当把群体内所有粒子都作为领域成员时，得到粒子群优化算法的全局版本；当群体内部分成员组成领域时得到粒子群优化算法的局部版本。局部版本中，一般有两种方式组成领域，一种是索引号相邻的粒子组成领域，另一种是位置相邻的粒子组成领域。粒子群优化算法的领域定义策略又可以称为粒子群的领域拓扑结构。[1]
2. 粒子群优化算法的流程

3. 粒子群优化算法的主要构成要素
1. 群体大小
是个整型参数。当很小的时候，陷入局优的可能性很大。然而，群体过大将导致计算时间大幅增加。并且当群体树木增长至一定水平时，再增长将不再有显著的作用。当=1时，PSO算法变成基于个体搜索的技术，一旦陷入局优，将不可能跳出。当m很大时，PSO的优化能力很好，可是收敛速度将非常慢。
2. 学习因子和
学习因子使粒子具有自我总结和向群体中优秀个体学习的能力，从而向群体内或领域内最优点靠近。和通常都等于2，不过也可能有其他取值。但是一般等于，并且范围在0和4之间。
3. 最大速度：
最大速度决定粒子在一次迭代中最大的移动距离。较大，探索能力增强，但是粒子容易飞过最好的解。较小时，开发能力增强，但是容易陷入局优。有分析和实验表明，设定的作用可以通过惯性权重的调整来实现。所以现在的实验基本上使用进行初始化，将设定为每维变量的变化范围，而不必进行细致的选择与调节。
4. 惯性权重
智能优化方法的运行是否成功，探索能力和开发能力的平衡是非常关键的。对于粒子群优化算法来说，这两种能力的平衡就是靠惯性权重来实现的。较大的惯性权重使粒子在自己原来的方向上具有更大的速度，从而在原方向上飞行更远，具有更好的搜索能力；较小的惯性权重使粒子继承了较少的原方向上的速度，从而飞行较近，具有更好的开发能力。通过调节惯性权重能够调节粒子群的搜索能力。
5. 领域拓扑结构
全局版本粒子群优化算法将整个群体作为粒子的领域，速度快，不过有时会陷入局优；局部版本粒子群优化算法将索引号相近或者位置相近的个体作为粒子的领域，收敛速度慢一点，不过很难陷入局部最优。显然，全局版本的粒子群优化算法可以看作局部版本粒子群优化算法的一个特例，即将整个群体都作为领域。
6. 停止准则
一般使用最大迭代次数或可以接受的满意解作为停止准则。
7. 粒子空间的初始化
较好地选择粒子的初始化空间，将大大缩短收敛时间。这是问题依赖的。
3. 粒子群算法求解旅行商问题的流程
粒子群优化算法虽然成功地应用于连续优化的问题中，而在离散上的研究和应用还很少，尤其是用PSO解决TSP问题是一个新的研究方向[2]。
最初的PSO是用来解决连续空间的问题的，为了适合求解TSP问题，人们对算法进行了各种改进。本文采用王岚[3]重新定义PSO的运算符号和规则，引入交换子和交换序的概念，构造一种特殊的PSO用于求解TSP的方法。先对这种改进的PSO进行简略介绍：
定义1 设n个城市的TSP问题的解序列为S=（ai），i=1,2，…，n.定义交换子SO（i1，i2）为交换解S中的点ai1和ai2，则S’=S+ SO（i1，i2）为解S经算子SO（i1，i2）操作后的新解，这里为符号“+”赋予了新的含义.
例1 有一个有5个城市的TSP问题，其解为S=（1,3,5,2,4），交换算子为SO（1，2），则S’=S+SO(1,2)=(1,3,5,2,4)+SO(1,2)=(3,1,5,2,4).
定义2 一个或多个交换子的有序队列就是交换序，记作SS。其中，SS=（SO1，SO2，…，SOn），SO1，SO2，…，SOn是交换子，它们之间的顺序是有意义的。
交换序作用于一个TSP解上意味着这个交换序中的所有交换子依次作用于该解上，即S’=S+SS=S+（SO1，SO2，…，SOn）=[(S+SO1)+SO2]+…+SOn.
定义3 不同的交换序作用于同一解上可能产生相同的新解，所有有相同效果的交换序的集合称为交换序的等价集。
定义4 若干个交换序可以合并成一个新的交换序，定义为两个交换序的合并算子。
例2 设两个交换序SS1和SS2按先后顺序作用于解S上，得到新解S’。假设另外有一个交换序SS’作用于同一解上，能够得到相同的解S’，可定义SS’=SS1+SS2。SS’和SS1+SS2属于同一等价集。一般来说，SS’不唯一。
定义5 在交换序等价集中，拥有最少交换子的序称为该等价集的基本交换序。可按如下方法构造一个基本交换序。
设给定两个解路径A和B，需要构造一个基本交换序SS，使得B+SS=A。
A：（1,2,3,4,5）； B：（2,3,1,5,4）
可以看出，A(1)=B(3)=1，所以第一个交换子是SO（1,3），B1=B+SO（1,3），得到B1（1,3,2,5,4）；A(2)=B1(3)=2，所以第二个交换子是SO（2,3），B2=B1+SO（2,3），得到B2=（1,2,3,5,4）；同理，第三个交换子是SO（4,5），得到B3=B2+SO（4,5）=A。这样就得到一个基本交换序：
SS=A-B=（SO（1,3），SO（2,3），SO（4,5））。
基本PSO算法中的速度算式在基于以上交换子和交换序的概念下各符号有了新的定义。其中、仍为学习因子，、仍为 [0,1]内 的随机数；表示基本交换序以的概率保留，表示基本交换序以的概率保留。由此可以看出的值越大，保留的交换子就越多，的影响就越大；同理，的值越大， 的影响也就越大。这也正是学习因子的含义所在。
求解TSP的PSO算法步骤描述如下：
（1）. 初始化粒子群，即给群体中每一个粒子赋一个随机的初始解和交换序；
（2）. 如果满足条件，转步骤（5）；
（3）. 根据粒子的当前位置，计算下一个位置，即新解：
i. 计算与之间的差A，A=，其中A是一个基本交换序表示A作用于得到；
ii. 计算B=，其中B也是一个基本交换序；
iii. 由计算出速度，并将交换序转换为一个基本交换序；
iv. 如果找到一个更好的解，则更新
（4）. 如果群体找到一个更好的解，则更新，转步骤（2）；
（5）. 输出结果。
4. 实例验证
表一:中国31个省会城市的坐标
城市序号
1
2
3
4
5
6
7
8
横坐标
1304
3639
4177
3712
3488
3326
3238
4196
纵坐标
2312
1315
2244
1399
1535
1556
1229
1044
城市序号
9
10
11
12
13
14
15
16
横坐标
4312
4386
3007
2562
2788
2381
1332
3715
纵坐标
790
570
1970
1756
1491
1676
695
1678
城市序号
17
18
19
20
21
22
23
24
横坐标
3918
4061
3780
3676
4029
4263
3429
3507
纵坐标
2179
2370
2212
2578
2838
2931
1908
2376
城市序号
25
26
27
28
29
30
31

横坐标
3394
3439
2935
3140
2545
2778
2370

纵坐标
2643
3201
3240
3550
2357
2826
2975

(表中数据来自网站
本文以中国31个省会城市的旅行商问题作为验证标准，验证以上改进的粒子群算法的实用效果。
编写基于以上思想的求解中国31个城市旅行商问题的matlab程序（代码见附件），当主要参数为
编号
snn
vmax
c1
c2
gnmax
1
100
30
2
2
100
（其中snn为最初粒子数，vmax为最大速度， c1为粒子基于自己当前位置与自己历史最优位置的自我学习因子， c2 为粒子基于自己当前位置与群体最优位置的群体学习因子，gnmax为最大迭代次数。）
时，进行十次实验，得到
次数
1
2
3
4
5
路径长度
31997
31653
31804
31946
31056
次数
6
7
8
9
10
路径长度
32298
33198
31335
32416
32964
可知第二次得到的结果最好，其对应的遍历路径为L: [ 6 4 11 13 29 21 5 16 14 30 15 28 26 27 25 17 24 20 19 3 22 18 10 2 7 12 1 31 23 8 9 ]、相应最佳粒子变化趋势及最好路径图为：
此结果与其他算法得到的最优值相差较大，但是可以看出，粒子确实在往好的方面变化。因此可以尝试改变参数，多次试验，看能否得到更优值。
5. 算法的改进
多次实验结果如下：
编号
snn
vmax
c1
c2
gnmax
最优值
1
100
30
1
3
100
31556
2
100
25
1
3
100
32211
3
100
20
1
3
100
31648
4
100
15
1
3
100
32649
5
100
10
1
3
100
31787
6
100
5
1
3
100
30590
7
100
5
1
2
100
32491
8
100
5
1
1
100
33169
9
100
5
2
3
100
32863
10
100
5
3
3
100
32109
11
150
5
1
3
100
31178
12
200
5
1
3
100
32106
13
250
5
1
3
100
31940
14
300
5
1
3
100
30050
15
350
5
1
3
100
30148
16
400
5
1
3
100
31069
17
300
5
1
3
50
32222
18
300
5
1
3
150
29746
19
300
5
1
3
200
29803
20
300
5
1
3
250
30598
21
300
5
1
3
300
29042
由以上实验结果可知，当初始化粒子总数snn为300，最大速度vmax为5，学习因子c1为1，c2为3，最大迭代次数为300时（即编号为21的实验），得到最短遍历长度为29042，其对应的遍历路径为L: [ 21 22 16 29 11 5 13 4 1 15 24 17 18 12 14 7 6 23 8 9 10 2 3 27 25 28 20 ]、相应最佳粒子变化趋势及最好路径图为：

6. 结论
改变各个参数值后求得的最优值差别不大，并且所得的最优值与其他算法求得的最优值相比更劣（遗传算法求解同样的问题得到的最优值为15404[3]），但是路径确实优化了，这种基于交换子和交换序的改进的粒子群优化算法求解TSP问题是有一定成效的。
由最佳粒子变化趋势图可知最优值的变化依赖于突变，并且这种突变发生效率很低，然而通过修改学习因子c1、c2、以及最大速度vmax并不能明显改变这种突变率。我想是因为这种基于交换子和交换序的改进的粒子群优化算法的性质所决定的。当然也有可能由于编码不恰当造成的，但是目前我尚未找到更好的编码方法。
这种源于大自然群体智能的粒子群优化算法一定具有一定的通用性的，我相信通过恰当的改进一定能够用于现实生活中的绝大多数优化问题。本文所用的交换子和交换序的概念以及对应的编码方法也一定能够改进的。
7. 参考文献
[1] 黄岚，王康平，周春光，庞魏，董龙江，彭利.粒子群优化算法求解旅行商问题.吉林大学学报（理学版），2003,4:477-480.
[2] 汪定伟，王俊伟，王洪峰等.智能优化方法.高等教育出版社.2007
[3] 丛爽，贾亚军.进化策略与蚁群算法融合的求解旅行商问题.控制工程.2011,v.8,no.1.
10