蚁群算法

蚁群算法（Ant Colony Algorithm）最初于1992年由意大利学者M.Dorigo等人提出，它是一种模拟自然界中真实蚁群觅食行为的仿生优化算法。研究发现：每只蚂蚁觅食时在走过的路线上会留下一种称为信息素的物质，蚂蚁之间靠感知这种物质的浓度进行信息传递。蚂蚁在选择路径时总是倾向于朝信息索浓度高的方向移动，而距离短的路径上走过的蚂蚁多，留下的信息素也多，后续蚂蚁选择它的概率也会越大；其他路径上的信息素会随着时间的推移不断挥发，这样就形成了一种正反馈机制，最后整个蚁群聚集到最短路径上。

人工蚁群算法模拟了这一过程。每只蚂蚁在解空间独立地搜索可行解，解越好留下的信息素越多，随着算法推进，较优解路径上的信息素增多，选择它的蚂蚁也随之增多，最终收敛到最优或近似最优的解上。

一、算法原理

蚂蚁系统是最早的蚁群系统，它采用 $τij(t)\tau_{ij}(t)$ 来模仿t时刻路径 $i$ 到 $j$ 上面的信息残留量，即信息素浓度。类似于蚂蚁觅食过程，每条路径上面的信息素会挥发，如果有蚂蚁经过的时候，信息素的浓度会相应增加。因此，蚂蚁系统中的信息素浓度的更新公式为：

τij(t+n)=ρ⋅τij(t)+Δτij\tau_{ij}(t+n)=\rho·\tau_{ij}(t)+\Delta\tau_{ij}

式中， $ρ\rho$ 是一个 $0$ 到 $1$ 的数字， $(1−ρ)(1-\rho)$ 为挥发因子。另外， $Δτij\Delta\tau_{ij}$ 表示一次旅行（遍历完所有城市）后，所有路径i到j的蚂蚁留下的信息素总量，即：

Δτij=∑k=1mΔτijkr\Delta\tau_{ij}=\sum_{k=1}^m\Delta\tau_{ij}^kr

式中， $Δτijk\Delta\tau_{ij}^k$ 表示第k只蚂蚁在路径 $i$ 到 $j$ 上面留下的信息素量。如果第 $k$ 只蚂蚁经过路径 $i$ 到 $j$ ，则：

Δτijk=QLkr\Delta\tau_{ij}^k=\frac{Q}{L_kr}

式中， $Q$ 为一个常数， $L_k$ 为蚂蚁已经走过路径的总长度。否则，第 $k$ 只蚂蚁在 $i$ 到 $j$ 上面留下的信息素量为 $0$ 。

一般来说有了信息素浓度的更新公式，就可以直接给出蚂蚁对每条路径的选择概率了。然而，为了更好的利用TSP问题自身的性质，M.Dorigo等引入了一个启发项： $η=1dij\eta=\frac{1}{d_{ij}}$ 。通过结合信息素浓度和启发因子，可以得到蚂蚁选择路径 $i$ 到 $j$ 的概率为：

[τik(t)]α⋅[ηik]β,j∈allowedk0,elsep_{ij}^{k}(t)=\left\{\begin{array}{rcl}\frac{[\tau_{ij}(t)]^\alpha·[\eta_{ij}]^{\beta}}{\sum_{k= \in allowed \: [\tau_{ik}(t)]^\alpha·[\eta_{ik}]^{\beta} } }, && {j \in allowed_k}\\ 0, & & {else}\end{array} \right.

式中， $α\alpha$ 和 $β\beta$ 是调节因子，用于调节 $τij(t)\tau_{ij}(t)$ 和 $ηij\eta_{ij}$ 之间的作用。此外 $allowed_k$ 表示蚂蚁 $k$ 还没有走过的路径（用禁忌表存储已经走过的路径），通过这种存储可以保证所有解的逻辑可行。如果路径 $i$ 到 $j$ 上的信息浓度越大 $τij(t)\tau_{ij}(t)$ 的值就越大，该路径被选择的概率就越大；同样，如果该路径长度越短，则 $η=1dij\eta=\frac{1}{d_{ij}}$ 越大，该路径被选择的概率也越大。

求解TSP问题的蚁群算法中的人工蚂蚁具有以下特点：

1）他们概率性地选择下一条路径，该概率与路径长度和路径上的信息素浓度有关；

2）为了保证解的逻辑可行，蚂蚁不允许选择已经走过的路径（通过禁忌表实现）；

3）蚂蚁走过一条路径时会在该路径上面分泌一种叫做信息素的物质。

二、代码实现

import random
import numpy as np
import math



#==========================================
#对称矩阵，计算任意两个城市之间的距离
def distance_p2p_mat():
    dis_mat=[]
    for i in range(num_city):
        dis_mat_each=[]
        for j in range(num_city):
            dis=math.sqrt(pow(location[i][0]-location[j][0],2)+pow(location[i][1]-location[j][1],2))
            dis_mat_each.append(dis)
        dis_mat.append(dis_mat_each)
   # print(dis_mat)
    return dis_mat

#计算所有寻找到的路径对应的距离
def cal_newpath(dis_mat,path_new):
    dis_list=[]
    for each in path_new:
        dis=0
        for j in range(num_city-1):
            dis=dis_mat[each[j]][each[j+1]]+dis
        dis=dis_mat[each[num_city-1]][each[0]]+dis#回家
        dis_list.append(dis)
    return dis_list
#==========================================






location=np.loadtxt('./city_location.txt')
num_ant=200 #蚂蚁个数
num_city=30 #城市个数

alpha=1 #信息素影响因子
beta=5  #期望影响因子
info=0.1 #信息素的挥发率
Q=1 #常数

count_iter = 0  #迭代计数器
iter_max = 30   #迭代次数





dis_list=distance_p2p_mat()     #计算任意两个城市间的距离
dis_mat=np.array(dis_list)      #将list转化为矩阵

e_mat_init=1.0/(dis_mat+np.diag([10000]*num_city))      #期望矩阵。加对角阵是原矩阵对角线为0，而除数不能是0，所以先用一个比较大的数垫一下
diag=np.diag([1.0/10000]*num_city)          #上一步生成的num_city*num_city维的对角线为1000的对角矩阵
e_mat=e_mat_init-diag           #已经做过除法了，所以让对角线元素复原，再减去那个对角矩阵

pheromone_mat=np.ones((num_city,num_city))    #信息浓度矩阵。初始化每条边的信息素浓度，全1矩阵

path_mat=np.zeros((num_ant,num_city)).astype(int)     #蚂蚁的路径矩阵。初始化每只蚂蚁路径，都从0城市出发.(如果不加数据转化，则默认生成的是float类型的0，即0.0)



while count_iter < iter_max:    #最外层迭代
    for ant in range(num_ant):  #对每一只蚂蚁进行分析
        visit=0     #都从0城市出发
        unvisit_list=list(range(1,30))  #未访问的城市。这个语句生成一个[1..29]的数组。再加上统一的出发点0，共30个城市
        for j in range(1,num_city):  #j代表第ant个蚂蚁的第j步
            
            trans_list=[]
          
            trans=0
            for k in range(len(unvisit_list)):  #第ant个蚂蚁的第j步取哪个城市
                trans +=np.power(pheromone_mat[visit][unvisit_list[k]],alpha)*np.power(e_mat[visit][unvisit_list[k]],beta)  #计算第ant个蚂蚁由visit位置向k位置走的概率。这里要注意：直接累加
                trans_list.append(trans)      #将 每一步 累加的结果保存到一个数组中
             
            #轮盘法选择下一个城市
            rand=random.uniform(0,trans)#产生随机数
            for t in range(len(trans_list)):
                if(rand <= trans_list[t]):  #因为之前就已经累加了，trans_list[t]一定是一个递增数组，所以可以直接与trans_list[t]相比较
                    visit_next=unvisit_list[t]  #选择下标为t的城市作为这个蚂蚁下一步的方向
                    break
      
            path_mat[ant,j]=visit_next  #装填这只蚂蚁的路径矩阵

            unvisit_list.remove(visit_next) #在未走的城市列表中删去这个结点。这个操作，就会使unvisit_list这个数组变成断断续续的
            visit=visit_next    #更新这只蚂蚁的当前位置
    

    #所有蚂蚁的路径表填满之后，算每只蚂蚁的总距离
    dis_allant_list=cal_newpath(dis_mat,path_mat)



    #选取拥有最短路径的蚂蚁的路径
    # 注意：这里其实可以不用写成if-else结构，但那样对于每一次迭代都需要一次min和max计算，而有时候后一次迭代结果并不一定比前一次迭代结果更优，就会造成冗余计算。
    #       使用if-else结构后，每一次迭代只是一个判断，并不需要每次都进行min和max计算
    if count_iter == 0:
        dis_new=min(dis_allant_list)
        path_new=path_mat[dis_allant_list.index(dis_new)].copy()      
    else:
        if min(dis_allant_list) < dis_new:
            dis_new=min(dis_allant_list)
            path_new=path_mat[dis_allant_list.index(dis_new)].copy() 



    # 为了避免残留信息素过多而淹没启发式信息，所以要及时的更新信息素矩阵
    pheromone_change=np.zeros((num_city,num_city))
    for i in range(num_ant):
        for j in range(num_city-1):
            pheromone_change[path_mat[i,j]][path_mat[i,j+1]] += Q/dis_mat[path_mat[i,j]][path_mat[i,j+1]]
        pheromone_change[path_mat[i,num_city-1]][path_mat[i,0]] += Q/dis_mat[path_mat[i,num_city-1]][path_mat[i,0]] #最后一个结点到起点
    
    pheromone_mat=(1-info)*pheromone_mat + pheromone_change

    count_iter += 1 #迭代计数+1，进入下一次迭代


        
print('最短距离：',dis_new)
print('最短路径：',path_new)

city_location.txt内容表示城市的横纵坐标，如下所示：