2022年 11月 9日

蚁群算法讲解python

简介

蚁群算法（Ant Clony Optimization， ACO）作为一个启发式群智能算法，它是由一群无智能或有轻微智能的个体通过相互协作而表现出智能行为，从而为求解复杂问题提供了一个新的可能性。

ACO是一种仿生学算法，是由自然界中蚂蚁觅食的行为而启发的。在自然界，蚂蚁觅食过程中，蚁群总能够按照寻找到一条从蚁巢和食物源的最优路径，这也就是蚁群算法的由来。

注：然而蚁群算法去做路径规划和优化智能算法或机器学习算法却不太一样。

路径规划的蚁群算法优化步骤（典型tsp问题）：

  1. 初始化一些参数。α——信息素重要程度因子；β——启发函数重要程度因子；ρ——信息素挥发因子；Q为常数；m为蚁群数量；城市规模num_city；迭代次数iter_max，信息素矩阵\tau _{ij}等。

  2. 计算概率。P_{ij}^{k}(t)=\begin{cases} \frac{\tau _{ij}(t)}{\sum_{s\in allowed_{k}}^{}\tau _{is}(t)}& \text{ if } j= allowed_{k} \\ 0& \text{ if } j\neq allowed_{k} \end{cases}，其中allowed_{k}为第k只蚂蚁后续可以走的城市集合；P_{ij}^{k}(t)为第k只蚂蚁在t时刻从城市i转移到城市j的概率。
  3. 更新信息素浓度。\begin{cases} \tau \left ( t+1 \right )=\left ( 1-\rho \right )\tau _{ij}\left ( t \right )+\sum_{k=1}^{m}\Delta\tau _{ij}^{k} & \text{ if }0< \rho < 1 \\ \Delta \tau _{ij}^{k}=\begin{cases} \frac{Q}{L_{k}}& \text{ if } L_{k}=distant(s_{ij}^{k}) \\ 0& \text{ if } else \end{cases}& \end{cases}，Lk为第k只蚂蚁曾经走过的路径。
  4. 判断是否达到最大迭代次数，达到从而结束算法。
  1. import random
  2. import math
  3. import numpy as np
  4. import matplotlib.pyplot as plt
  5.  
  6.  
  7. class ACO(object):
  8.     def __init__(self, num_city, data):
  9.         self.m = 50  # 蚂蚁数量
  10.         self.alpha = 0.5  # 信息素重要程度因子
  11.         self.beta = 5  # 启发函数重要因子
  12.         self.rho = 0.1  # 信息素挥发因子
  13.  
  14.         self.Q = 1  # 常量系数
  15.         self.num_city = num_city  # 城市规模
  16.         self.location = data  # 城市坐标
  17.         self.Tau = np.zeros([num_city, num_city])  # 信息素矩阵
  18.         self.Table = [[0 for _ in range(num_city)] for _ in range(self.m)]  # 生成的蚁群
  19.         self.iter = 1
  20.         self.iter_max = 500
  21.         self.dis_mat = self.compute_dis_mat(num_city, self.location)  # 计算城市之间的距离矩阵
  22.         self.Eta = 10. / self.dis_mat  # 启发式函数
  23.         self.paths = None  # 蚁群中每个个体的长度
  24.         # 存储存储每个温度下的最终路径，画出收敛图
  25.         self.iter_x = []
  26.         self.iter_y = []
  27.         # self.greedy_init(self.dis_mat,100,num_city)
  28.  
  29.     def greedy_init(self, dis_mat, num_total, num_city):
  30.         start_index = 0
  31.         result = []
  32.         for i in range(num_total):
  33.             rest = [x for x in range(0, num_city)]
  34.             # 所有起始点都已经生成了
  35.             if start_index >= num_city:
  36.                 start_index = np.random.randint(0, num_city)
  37.                 result.append(result[start_index].copy())
  38.                 continue
  39.             current = start_index
  40.             rest.remove(current)
  41.             # 找到一条最近邻路径
  42.             result_one = [current]
  43.             while len(rest) != 0:
  44.                 tmp_min = math.inf
  45.                 tmp_choose = -1
  46.                 for x in rest:
  47.                     if dis_mat[current][x] < tmp_min:
  48.                         tmp_min = dis_mat[current][x]
  49.                         tmp_choose = x
  50.  
  51.                 current = tmp_choose
  52.                 result_one.append(tmp_choose)
  53.                 rest.remove(tmp_choose)
  54.             result.append(result_one)
  55.             start_index += 1
  56.         pathlens = self.compute_paths(result)
  57.         sortindex = np.argsort(pathlens)       # argsort()是将X中的元素从小到大排序后，提取对应的索引index，然后输出到y
  58.         index = sortindex[0]
  59.         result = result[index]
  60.         for i in range(len(result)-1):
  61.             s = result[i]
  62.             s2 = result[i+1]
  63.             self.Tau[s][s2] = 1
  64.         self.Tau[result[-1]][result[0]] = 1
  65.         # for i in range(num_city):
  66.         #     for j in range(num_city):
  67.         # return result
  68.  
  69.     # 轮盘赌选择
  70.     def rand_choose(self, p):
  71.         x = np.random.rand()
  72.         for i, t in enumerate(p):
  73.             x -= t
  74.             if x <= 0:
  75.                 break
  76.         return i
  77.  
  78.     # 生成蚁群
  79.     def get_ants(self, num_city):
  80.         for i in range(self.m):
  81.             start = np.random.randint(num_city - 1)
  82.             self.Table[i][0] = start
  83.             unvisit = list([x for x in range(num_city) if x != start])
  84.             current = start
  85.             j = 1
  86.             while len(unvisit) != 0:
  87.                 P = []
  88.                 # 通过信息素计算城市之间的转移概率
  89.                 for v in unvisit:
  90.                     P.append(self.Tau[current][v] ** self.alpha * self.Eta[current][v] ** self.beta)
  91.                 P_sum = sum(P)
  92.                 P = [x / P_sum for x in P]
  93.                 # 轮盘赌选择一个一个城市
  94.                 index = self.rand_choose(P)
  95.                 current = unvisit[index]
  96.                 self.Table[i][j] = current
  97.                 unvisit.remove(current)
  98.                 j += 1
  99.  
  100.     # 计算不同城市之间的距离
  101.     def compute_dis_mat(self, num_city, location):
  102.         dis_mat = np.zeros((num_city, num_city))
  103.         for i in range(num_city):
  104.             for j in range(num_city):
  105.                 if i == j:
  106.                     dis_mat[i][j] = np.inf
  107.                     continue
  108.                 a = location[i]
  109.                 b = location[j]
  110.                 tmp = np.sqrt(sum([(x[0] - x[1]) ** 2 for x in zip(a, b)]))
  111.                 dis_mat[i][j] = tmp
  112.         return dis_mat
  113.  
  114.     # 计算一条路径的长度
  115.     def compute_pathlen(self, path, dis_mat):
  116.         a = path[0]
  117.         b = path[-1]
  118.         result = dis_mat[a][b]
  119.         for i in range(len(path) - 1):
  120.             a = path[i]
  121.             b = path[i + 1]
  122.             result += dis_mat[a][b]
  123.         return result
  124.  
  125.     # 计算一个群体的长度
  126.     def compute_paths(self, paths):
  127.         result = []
  128.         for one in paths:
  129.             length = self.compute_pathlen(one, self.dis_mat)
  130.             result.append(length)
  131.         return result
  132.  
  133.     # 更新信息素
  134.     def update_Tau(self):
  135.         delta_tau = np.zeros([self.num_city, self.num_city])
  136.         paths = self.compute_paths(self.Table)
  137.         for i in range(self.m):
  138.             for j in range(self.num_city - 1):
  139.                 a = self.Table[i][j]
  140.                 b = self.Table[i][j + 1]
  141.                 delta_tau[a][b] = delta_tau[a][b] + self.Q / paths[i]
  142.             a = self.Table[i][0]
  143.             b = self.Table[i][-1]
  144.             delta_tau[a][b] = delta_tau[a][b] + self.Q / paths[i]
  145.         self.Tau = (1 - self.rho) * self.Tau + delta_tau
  146.  
  147.     def aco(self):
  148.         best_lenth = math.inf           # math.inf返回浮点正无穷大
  149.         best_path = None
  150.         for cnt in range(self.iter_max):
  151.             # 生成新的蚁群
  152.             self.get_ants(self.num_city)  # out>>self.Table
  153.             self.paths = self.compute_paths(self.Table)
  154.             # 取该蚁群的最优解
  155.             tmp_lenth = min(self.paths)
  156.             tmp_path = self.Table[self.paths.index(tmp_lenth)]
  157.             # 可视化初始的路径
  158.             if cnt == 0:
  159.                 init_show = self.location[tmp_path]
  160.                 init_show = np.vstack([init_show, init_show[0]])
  161.             # 更新最优解
  162.             if tmp_lenth < best_lenth:
  163.                 best_lenth = tmp_lenth
  164.                 best_path = tmp_path
  165.             # 更新信息素
  166.             self.update_Tau()
  167.  
  168.             # 保存结果
  169.             self.iter_x.append(cnt)
  170.             self.iter_y.append(best_lenth)
  171.             print(cnt,best_lenth)
  172.         return best_lenth, best_path
  173.  
  174.     def run(self):
  175.         best_length, best_path = self.aco()
  176.         return self.location[best_path], best_length
  177.  
  178.  
  179. # 读取数据
  180. def read_tsp(path):
  181.     lines = open(path, 'r').readlines()
  182.     assert 'NODE_COORD_SECTION\n' in lines
  183.     index = lines.index('NODE_COORD_SECTION\n')
  184.     data = lines[index + 1:-1]
  185.     tmp = []
  186.     for line in data:
  187.         line = line.strip().split(' ')
  188.         if line[0] == 'EOF':
  189.             continue
  190.         tmpline = []
  191.         for x in line:
  192.             if x == '':
  193.                 continue
  194.             else:
  195.                 tmpline.append(float(x))
  196.         if tmpline == []:
  197.             continue
  198.         tmp.append(tmpline)
  199.     data = tmp
  200.     return data
  201.  
  202.  
  203. data = read_tsp('data/st70.tsp')
  204.  
  205. data = np.array(data)
  206. data = data[:, 1:]
  207. # 加上一行因为会回到起点
  208. show_data = np.vstack([data, data[0]])
  209.  
  210. aco = ACO(num_city=data.shape[0], data=data.copy())
  211. Best_path, Best = aco.run()
  212. print(Best)
  213. Best_path = np.vstack([Best_path, Best_path[0]])
  214. plt.plot(Best_path[:, 0], Best_path[:, 1])
  215. plt.title('st70:aco-tsp')
  216. plt.show()

蚁群算法优化智能算法（以优化函数为例）： 

  1. 初始化种群pop_size=100 ，最大迭代次数NGEN=50，信息素挥发因子rou=0.8，Q常量为1。
  2. 随机产生蚂蚁初始位置，计算适应度函数值，设为初始信息素，计算信息素转移概率 
    pi[i] = (t_max - t[i]) / t_max

  3. 依据概率计算下一时刻位置 
    1.  if pi[i] < np.random.uniform(0, 1):
    2.                     self.pop_x[i][j] = self.pop_x[i][j] + np.random.uniform(-1, 1) * lamda
    3.                 else:
    4.                     self.pop_x[i][j] = self.pop_x[i][j] + np.random.uniform(-1, 1) * (
    5.                             self.bound[1][j] - self.bound[0][j]) / 2

  4. 信息素更新。 
    t[i] = (1 - rou) * t[i] + Q * self.fitness(self.pop_x[i])

  5. 反复迭代，直到达到最大迭代次数

完整python代码：

  1. import numpy as np
  2. import matplotlib.pyplot as plt
  3.  
  4.  
  5. class ACO:
  6.     def __init__(self, parameters):
  7.         # 初始化
  8.         self.NGEN = parameters[0]  # 迭代的代数
  9.         self.pop_size = parameters[1]  # 种群大小
  10.         self.var_num = len(parameters[2])  # 变量个数
  11.         self.bound = []  # 变量的约束范围
  12.         self.bound.append(parameters[2])
  13.         self.bound.append(parameters[3])
  14.  
  15.         self.pop_x = np.zeros((self.pop_size, self.var_num))  # 所有蚂蚁的位置
  16.         self.g_best = np.zeros((1, self.var_num))  # 全局蚂蚁最优的位置
  17.  
  18.         # 初始化第0代初始全局最优解
  19.         temp = -1
  20.         for i in range(self.pop_size):
  21.             for j in range(self.var_num):
  22.                 self.pop_x[i][j] = np.random.uniform(self.bound[0][j], self.bound[1][j])
  23.             fit = self.fitness(self.pop_x[i])
  24.             if fit > temp:
  25.                 self.g_best = self.pop_x[i]
  26.                 temp = fit
  27.  
  28.     def fitness(self, ind_var):
  29.         """
  30.         个体适应值计算
  31.         """
  32.         x1 = ind_var[0]
  33.         x2 = ind_var[1]
  34.         x3 = ind_var[2]
  35.  
  36.         y = 4*x1 ** 2 + 2*x2 + x3 ** 3
  37.         return y
  38.  
  39.     def update_operator(self, gen, t, t_max):
  40.         """
  41.         更新算子：根据概率更新下一时刻的位置
  42.         """
  43.         rou = 0.8  # 信息素挥发系数
  44.         Q = 1  # 信息释放总量
  45.         lamda = 1 / gen
  46.         pi = np.zeros(self.pop_size)
  47.         for i in range(self.pop_size):
  48.             for j in range(self.var_num):
  49.                 pi[i] = (t_max - t[i]) / t_max
  50.                 # 更新位置
  51.                 if pi[i] < np.random.uniform(0, 1):
  52.                     self.pop_x[i][j] = self.pop_x[i][j] + np.random.uniform(-1, 1) * lamda
  53.                 else:
  54.                     self.pop_x[i][j] = self.pop_x[i][j] + np.random.uniform(-1, 1) * (
  55.                             self.bound[1][j] - self.bound[0][j]) / 2
  56.                 # 越界保护
  57.                 if self.pop_x[i][j] < self.bound[0][j]:
  58.                     self.pop_x[i][j] = self.bound[0][j]
  59.                 if self.pop_x[i][j] > self.bound[1][j]:
  60.                     self.pop_x[i][j] = self.bound[1][j]
  61.             # 更新t值
  62.             t[i] = (1 - rou) * t[i] + Q * self.fitness(self.pop_x[i])
  63.             # 更新全局最优值
  64.             if self.fitness(self.pop_x[i]) > self.fitness(self.g_best):
  65.                 self.g_best = self.pop_x[i]
  66.         t_max = np.max(t)
  67.         return t_max, t
  68.  
  69.     def main(self):
  70.         popobj = []
  71.         best = np.zeros((1, self.var_num))[0]
  72.         for gen in range(1, self.NGEN + 1):
  73.             if gen == 1:
  74.                 tmax, t = self.update_operator(gen, np.array(list(map(self.fitness, self.pop_x))),
  75.                                                np.max(np.array(list(map(self.fitness, self.pop_x)))))
  76.             else:
  77.                 tmax, t = self.update_operator(gen, t, tmax)
  78.  
  79.             print('############ Generation {} ############'.format(str(gen)))
  80.             print(self.g_best)
  81.             print(self.fitness(self.g_best))
  82.             if self.fitness(self.g_best) > self.fitness(best):
  83.                 best = self.g_best.copy()
  84.             popobj.append(self.fitness(best))
  85.             print('最好的位置：{}'.format(best))
  86.             print('最大的函数值：{}'.format(self.fitness(best)))
  87.         print("---- End of (successful) Searching ----")
  88.  
  89.         plt.figure()
  90.         plt.title("Figure1")
  91.         plt.xlabel("iterators", size=14)
  92.         plt.ylabel("fitness", size=14)
  93.         t = [t for t in range(1, self.NGEN + 1)]
  94.         plt.plot(t, popobj, color='b', linewidth=2)
  95.         plt.show()
  96.  
  97.  
  98. if __name__ == '__main__':
  99.     NGEN = 100
  100.     popsize = 50
  101.     low = [1, 1, 1]
  102.     up = [30, 30, 30]
  103.     parameters = [NGEN, popsize, low, up]
  104.     aco = ACO(parameters)
  105.     aco.main()

 效果不好时，大家记得调节参数。

总结：

蚁群算法缺点：

  1. 收敛速度慢
  2. 易于陷入局部最优

参数作用：

  1. 蚂蚁数量：种群数量越多，得到的最优解就越精确,但是会有大量蚂蚁重复路过同一路径,运行时间增多。
  2. alpha  信息素重要程度因子：alpha值过大,蚂蚁选择之前走过的路径可能性就越大,蚁群搜索路径的随机性减弱。alpha值过小,蚁群搜索范围就会减少,易陷入局部最优解。
  3. beta 启发函数因子：beta值增大,蚁群更容易选择局部较短路径,能加快算法的运行速度,但是可能陷入局部最优解。
  4. ruo 信息素挥发因子：当ruo很小的时候，每条路径的残留信息很多，就会反复进行搜索，运算时间增加；ruo很大的时候,会放弃搜索很多有效路径,可能会丢失最优值。