螞蟻幾乎沒有視力,但他們卻能夠在黑暗的世界中找到食物,而且能夠找到一條從洞穴到食物的最短路徑。它們是如何做到的呢?

簡介
由來
蟻群算法是一種用來尋找優化路徑的機率型算法。它由Marco Dorigo于1992年在他的博士論文中提出,其靈感來源于螞蟻在尋找食物過程中發現路徑的行為。
這種算法具有分布計算、資訊正回報和啟發式搜尋的特征,本質上是進化算法中的一種啟發式全局優化算法。
思想
将蟻群算法應用于解決優化問題的基本思路為:用螞蟻的行走路徑表示待優化問題的可行解,整個螞蟻群體的所有路徑構成待優化問題的解空間。路徑較短的螞蟻釋放的資訊素量較多,随着時間的推進,較短的路徑上累積的資訊素濃度逐漸增高,選擇該路徑的螞蟻個數也愈來愈多。最終,整個螞蟻會在正回報的作用下集中到最佳的路徑上,此時對應的便是待優化問題的最優解。
應用
蟻群算法應用廣泛,如旅行商問題(traveling salesman problem,簡稱TSP)、指派問題、Job-shop排程問題、車輛路徑問題(vehicle routing problem)、圖着色問題(graph coloring problem)和網絡路由問題(network routing problem)等等,下面以TSP的求解為例。
scikit-opt庫實作
import numpy as np
from scipy import spatial
import pandas as pd
import matplotlib.pyplot as plt
num_points = 25
points_coordinate = np.random.rand(num_points, 2) # generate coordinate of points
distance_matrix = spatial.distance.cdist(points_coordinate, points_coordinate, metric=\'euclidean\')
def cal_total_distance(routine):
num_points, = routine.shape
return sum([distance_matrix[routine[i % num_points], routine[(i + 1) % num_points]] for i in range(num_points)])
from sko.ACA import ACA_TSP
aca = ACA_TSP(func=cal_total_distance, n_dim=num_points,
size_pop=50, max_iter=200,
distance_matrix=distance_matrix)
best_x, best_y = aca.run()
# %% Plot
fig, ax = plt.subplots(1, 2)
best_points_ = np.concatenate([best_x, [best_x[0]]])
best_points_coordinate = points_coordinate[best_points_, :]
ax[0].plot(best_points_coordinate[:, 0], best_points_coordinate[:, 1], \'o-r\')
pd.DataFrame(aca.y_best_history).cummin().plot(ax=ax[1])
plt.show()
View Code
Python實作
# -*- coding: utf-8 -*-
import random
import copy
import time
import sys
import math
import tkinter #//GUI子產品
import threading
from functools import reduce
# 參數
\'\'\'
ALPHA:資訊啟發因子,值越大,則螞蟻選擇之前走過的路徑可能性就越大
,值越小,則蟻群搜尋範圍就會減少,容易陷入局部最優
BETA:Beta值越大,蟻群越就容易選擇局部較短路徑,這時算法收斂速度會
加快,但是随機性不高,容易得到局部的相對最優
\'\'\'
(ALPHA, BETA, RHO, Q) = (1.0,2.0,0.5,100.0)
# 城市數,蟻群
(city_num, ant_num) = (50,50)
distance_x = [
178,272,176,171,650,499,267,703,408,437,491,74,532,
416,626,42,271,359,163,508,229,576,147,560,35,714,
757,517,64,314,675,690,391,628,87,240,705,699,258,
428,614,36,360,482,666,597,209,201,492,294]
distance_y = [
170,395,198,151,242,556,57,401,305,421,267,105,525,
381,244,330,395,169,141,380,153,442,528,329,232,48,
498,265,343,120,165,50,433,63,491,275,348,222,288,
490,213,524,244,114,104,552,70,425,227,331]
#城市距離和資訊素
distance_graph = [ [0.0 for col in range(city_num)] for raw in range(city_num)]
pheromone_graph = [ [1.0 for col in range(city_num)] for raw in range(city_num)]
#----------- 螞蟻 -----------
class Ant(object):
# 初始化
def __init__(self,ID):
self.ID = ID # ID
self.__clean_data() # 随機初始化出生點
# 初始資料
def __clean_data(self):
self.path = [] # 目前螞蟻的路徑
self.total_distance = 0.0 # 目前路徑的總距離
self.move_count = 0 # 移動次數
self.current_city = -1 # 目前停留的城市
self.open_table_city = [True for i in range(city_num)] # 探索城市的狀态
city_index = random.randint(0,city_num-1) # 随機初始出生點
self.current_city = city_index
self.path.append(city_index)
self.open_table_city[city_index] = False
self.move_count = 1
# 選擇下一個城市
def __choice_next_city(self):
next_city = -1
select_citys_prob = [0.0 for i in range(city_num)] #存儲去下個城市的機率
total_prob = 0.0
# 擷取去下一個城市的機率
for i in range(city_num):
if self.open_table_city[i]:
try :
# 計算機率:與資訊素濃度成正比,與距離成反比
select_citys_prob[i] = pow(pheromone_graph[self.current_city][i], ALPHA) * pow((1.0/distance_graph[self.current_city][i]), BETA)
total_prob += select_citys_prob[i]
except ZeroDivisionError as e:
print (\'Ant ID: {ID}, current city: {current}, target city: {target}\'.format(ID = self.ID, current = self.current_city, target = i))
sys.exit(1)
# 輪盤選擇城市
if total_prob > 0.0:
# 産生一個随機機率,0.0-total_prob
temp_prob = random.uniform(0.0, total_prob)
for i in range(city_num):
if self.open_table_city[i]:
# 輪次相減
temp_prob -= select_citys_prob[i]
if temp_prob < 0.0:
next_city = i
break
# 未從機率産生,順序選擇一個未通路城市
# if next_city == -1:
# for i in range(city_num):
# if self.open_table_city[i]:
# next_city = i
# break
if (next_city == -1):
next_city = random.randint(0, city_num - 1)
while ((self.open_table_city[next_city]) == False): # if==False,說明已經周遊過了
next_city = random.randint(0, city_num - 1)
# 傳回下一個城市序号
return next_city
# 計算路徑總距離
def __cal_total_distance(self):
temp_distance = 0.0
for i in range(1, city_num):
start, end = self.path[i], self.path[i-1]
temp_distance += distance_graph[start][end]
# 回路
end = self.path[0]
temp_distance += distance_graph[start][end]
self.total_distance = temp_distance
# 移動操作
def __move(self, next_city):
self.path.append(next_city)
self.open_table_city[next_city] = False
self.total_distance += distance_graph[self.current_city][next_city]
self.current_city = next_city
self.move_count += 1
# 搜尋路徑
def search_path(self):
# 初始化資料
self.__clean_data()
# 搜素路徑,周遊完所有城市為止
while self.move_count < city_num:
# 移動到下一個城市
next_city = self.__choice_next_city()
self.__move(next_city)
# 計算路徑總長度
self.__cal_total_distance()
#----------- TSP問題 -----------
class TSP(object):
def __init__(self, root, width = 800, height = 600, n = city_num):
# 建立畫布
self.root = root
self.width = width
self.height = height
# 城市數目初始化為city_num
self.n = n
# tkinter.Canvas
self.canvas = tkinter.Canvas(
root,
width = self.width,
height = self.height,
bg = "#EBEBEB", # 背景白色
xscrollincrement = 1,
yscrollincrement = 1
)
self.canvas.pack(expand = tkinter.YES, fill = tkinter.BOTH)
self.title("TSP蟻群算法(n:初始化 e:開始搜尋 s:停止搜尋 q:退出程式)")
self.__r = 5
self.__lock = threading.RLock() # 線程鎖
self.__bindEvents()
self.new()
# 計算城市之間的距離
for i in range(city_num):
for j in range(city_num):
temp_distance = pow((distance_x[i] - distance_x[j]), 2) + pow((distance_y[i] - distance_y[j]), 2)
temp_distance = pow(temp_distance, 0.5)
distance_graph[i][j] =float(int(temp_distance + 0.5))
# 按鍵響應程式
def __bindEvents(self):
self.root.bind("q", self.quite) # 退出程式
self.root.bind("n", self.new) # 初始化
self.root.bind("e", self.search_path) # 開始搜尋
self.root.bind("s", self.stop) # 停止搜尋
# 更改标題
def title(self, s):
self.root.title(s)
# 初始化
def new(self, evt = None):
# 停止線程
self.__lock.acquire()
self.__running = False
self.__lock.release()
self.clear() # 清除資訊
self.nodes = [] # 節點坐标
self.nodes2 = [] # 節點對象
# 初始化城市節點
for i in range(len(distance_x)):
# 在畫布上随機初始坐标
x = distance_x[i]
y = distance_y[i]
self.nodes.append((x, y))
# 生成節點橢圓,半徑為self.__r
node = self.canvas.create_oval(x - self.__r,
y - self.__r, x + self.__r, y + self.__r,
fill = "#ff0000", # 填充紅色
outline = "#000000", # 輪廓白色
tags = "node",
)
self.nodes2.append(node)
# 顯示坐标
self.canvas.create_text(x,y-10, # 使用create_text方法在坐标(302,77)處繪制文字
text = \'(\'+str(x)+\',\'+str(y)+\')\', # 所繪制文字的内容
fill = \'black\' # 所繪制文字的顔色為灰色
)
# 順序連接配接城市
#self.line(range(city_num))
# 初始城市之間的距離和資訊素
for i in range(city_num):
for j in range(city_num):
pheromone_graph[i][j] = 1.0
self.ants = [Ant(ID) for ID in range(ant_num)] # 初始蟻群
self.best_ant = Ant(-1) # 初始最優解
self.best_ant.total_distance = 1 << 31 # 初始最大距離
self.iter = 1 # 初始化疊代次數
# 将節點按order順序連線
def line(self, order):
# 删除原線
self.canvas.delete("line")
def line2(i1, i2):
p1, p2 = self.nodes[i1], self.nodes[i2]
self.canvas.create_line(p1, p2, fill = "#000000", tags = "line")
return i2
# order[-1]為初始值
reduce(line2, order, order[-1])
# 清除畫布
def clear(self):
for item in self.canvas.find_all():
self.canvas.delete(item)
# 退出程式
def quite(self, evt):
self.__lock.acquire()
self.__running = False
self.__lock.release()
self.root.destroy()
print (u"\n程式已退出...")
sys.exit()
# 停止搜尋
def stop(self, evt):
self.__lock.acquire()
self.__running = False
self.__lock.release()
# 開始搜尋
def search_path(self, evt = None):
# 開啟線程
self.__lock.acquire()
self.__running = True
self.__lock.release()
while self.__running:
# 周遊每一隻螞蟻
for ant in self.ants:
# 搜尋一條路徑
ant.search_path()
# 與目前最優螞蟻比較
if ant.total_distance < self.best_ant.total_distance:
# 更新最優解
self.best_ant = copy.deepcopy(ant)
# 更新資訊素
self.__update_pheromone_gragh()
print (u"疊代次數:",self.iter,u"最佳路徑總距離:",int(self.best_ant.total_distance))
# 連線
self.line(self.best_ant.path)
# 設定标題
self.title("TSP蟻群算法(n:随機初始 e:開始搜尋 s:停止搜尋 q:退出程式) 疊代次數: %d" % self.iter)
# 更新畫布
self.canvas.update()
self.iter += 1
# 更新資訊素
def __update_pheromone_gragh(self):
# 擷取每隻螞蟻在其路徑上留下的資訊素
temp_pheromone = [[0.0 for col in range(city_num)] for raw in range(city_num)]
for ant in self.ants:
for i in range(1,city_num):
start, end = ant.path[i-1], ant.path[i]
# 在路徑上的每兩個相鄰城市間留下資訊素,與路徑總距離反比
temp_pheromone[start][end] += Q / ant.total_distance
temp_pheromone[end][start] = temp_pheromone[start][end]
# 更新所有城市之間的資訊素,舊資訊素衰減加上新疊代資訊素
for i in range(city_num):
for j in range(city_num):
pheromone_graph[i][j] = pheromone_graph[i][j] * RHO + temp_pheromone[i][j]
# 主循環
def mainloop(self):
self.root.mainloop()
#----------- 程式的入口處 -----------
if __name__ == \'__main__\':
print (u"""
--------------------------------------------------------
程式:蟻群算法解決TPS問題程式
作者:許彬
日期:2015-12-10
語言:Python 2.7
--------------------------------------------------------
""")
TSP(tkinter.Tk()).mainloop()
View Code
C實作
#include<bits/stdc++.h>
using namespace std;
const int INF = 0x3f3f3f3f;
#define sqr(x) ((x)*(x))
#define eps 1e-8
string file_name;
int type;// type == 1 全矩陣, type == 2 二維歐拉距離
int N;//城市數量
double **dis;//城市間距離
double **pheromone;//資訊素
double **herustic;//啟發式值
double **info;// info = pheromone ^ delta * herustic ^ beta
double pheromone_0;//pheromone初始值,這裡是1 / (avg * N)其中avg為圖網中所有邊邊權的平均數。
int m;//種群數量
int delta, beta;//參數
double alpha;
int *r1, *s, *r;//agent k的出發城市,下一個點,目前點。
int MAX, iteration;//最大疊代次數,疊代計數變量
set<int> empty, *J;
struct vertex{
double x, y;// 城市坐标
int id;// 城市編号
int input(FILE *fp){
return fscanf(fp, "%d %lf %lf", &id, &x, &y);
}
}*node;
typedef pair<int, int> pair_int;
struct Tour{//路徑
vector<pair_int> path;//path[i],存儲一條邊(r,s)
double L;
void clean(){
L = INF;
path.clear();
path.shrink_to_fit();
}//清空
void calc(){
L = 0;
int sz = path.size();
for (int i = 0; i < sz; i ++){
L += dis[path[i].first][path[i].second];
}
}//計算長度
void push_back(int x, int y){
path.push_back(make_pair(x, y));
}
int size(){
return (int)path.size();
}
int r(int i){
return path[i].first;
}
int s(int i){
return path[i].second;
}
void print(FILE *fp){
int sz = path.size();
for (int i = 0; i < sz; i ++){
fprintf(fp, "%d->", path[i].first + 1);
}
fprintf(fp, "%d\n", path[sz - 1].second + 1);
}
bool operator <(const Tour &a)const{
return L < a.L;
}//重載
} *tour, best_so_far;
double EUC_2D(const vertex &a, const vertex &b){
return sqrt(sqr(a.x - b.x) + sqr(a.y - b.y));
}
void io(){//輸入
printf("input file_name and data type\n");
cin >> file_name >> type;
FILE *fp = fopen(file_name.c_str(), "r");
fscanf(fp, "%d", &N);
node = new vertex[N + 5];
dis = new double*[N + 5];
double tmp = 0;
int cnt = 0;
if (type == 1){
for (int i = 0; i < N; i ++){
dis[i] = new double[N];
for (int j = 0; j < N; j ++){
fscanf(fp, "%lf", &dis[i][j]);
tmp += i != j ? dis[i][j] : 0;// i == j的時候dis不存在,是以不考慮。
cnt += i != j ? 1 : 0;// i == j的時候dis不存在,是以不考慮。
}
}
}else{
for (int i = 0; i < N; i ++)
node[i].input(fp);
for (int i = 0; i < N; i ++){
dis[i] = new double[N];
for (int j = 0; j < N; j ++){
dis[i][j] = EUC_2D(node[i], node[j]);// 計算距離
tmp += i != j ? dis[i][j] : 0;// i == j的時候 dis不存在,是以不考慮。
cnt += i != j ? 1 : 0;// i == j的時候dis不存在,是以不考慮。
}
}
}
pheromone_0 = (double)cnt / (tmp * N);//pheromone初始值,這裡是1 / (avg * N)其中avg為圖網中所有邊邊權的平均數。
fclose(fp);
return;
}
void init(){//初始化
alpha = 0.1;//evaporation parameter,揮發參數,每次資訊素要揮發的量
delta = 1;
beta = 6;// delta 和 beta分别表示pheromones和herustic的比重
m = N;
pheromone = new double*[N + 5];
herustic = new double*[N + 5];
info = new double*[N + 5];
r1 = new int[N + 5];
r = new int[N + 5];
s = new int[N + 5];
J = new set<int>[N + 5];
empty.clear();
for (int i = 0; i < N; i ++){
pheromone[i] = new double[N + 5];
herustic[i] = new double[N + 5];
info[i] = new double[N + 5];
empty.insert(i);
for (int j = 0; j < N; j ++){
pheromone[i][j] = pheromone_0;
herustic[i][j] = 1 / (dis[i][j] + eps);//加一個小數eps,防止被零除
}
}
best_so_far.clean();
iteration = 0;
MAX = N * N;
}
double power(double x, int y){//快速幂,計算x ^ y,時間複雜度O(logn),感興趣可以百度
double ans = 1;
while (y){
if (y & 1) ans *= x;
x *= x;
y >>= 1;
}
return ans;
}
void reset(){
tour = new Tour[m + 5];
for (int i = 0; i < N; i ++){
tour[i].clean();
r1[i] = i;//初始化出發城市,
J[i] = empty;
J[i].erase(r1[i]);//初始化agent i需要通路的城 市
r[i] = r1[i];//目前在出發點
}
for (int i = 0; i < N; i ++)
for (int j = 0; j < N; j ++){
info[i][j] = power(pheromone[i][j], delta) * power(herustic[i][j], beta);
}//選擇公式
}
int select_next(int k){
if (J[k].empty()) return r1[k]; //如果J是空的,那麼傳回出發點城市
double rnd = (double)(rand()) / (double)RAND_MAX;//産生0..1的随機數
set<int>::iterator it = J[k].begin();
double sum_prob = 0, sum = 0;
for (; it != J[k].end(); it ++){
sum += info[r[k]][*it];
}//計算機率分布
rnd *= sum;
it = J[k].begin();
for (; it != J[k].end(); it ++){
sum_prob += info[r[k]][*it];
if (sum_prob >= rnd){
return *it;
}
}//依照機率選取下一步城市
}
void construct_solution(){
for (int i = 0; i < N; i ++){
for (int k = 0; k < m; k ++){
int next = select_next(k);//選擇下一步的最優決策
J[k].erase(next);
s[k] = next;
tour[k].push_back(r[k], s[k]);
r[k] = s[k];
}
}
}
void update_pheromone(){
Tour now_best;
now_best.clean();//初始化
for (int i = 0; i < m; i ++){
tour[i].calc();
if (tour[i] < now_best)
now_best = tour[i];//尋找目前疊代最優解
}
if (now_best < best_so_far){
best_so_far = now_best;//更新最優解
}
for (int i = 0; i < N; i ++)
for (int j = 0; j < N; j ++)
pheromone[i][j] *= (1 - alpha);//資訊素揮發
int sz = now_best.size();
for (int i = 0; i < sz; i ++){
pheromone[now_best.r(i)][now_best.s(i)] += 1. / (double)now_best.L;
pheromone[now_best.s(i)][now_best.r(i)] = pheromone[now_best.r(i)][now_best.s(i)];// 對稱
}//更新資訊素含量
}
int main(){
srand((unsigned) time(0));//初始化随機種子
io();
init();
double last = INF;
int bad_times = 0;
for (; iteration < MAX; iteration ++){
if (bad_times > N) break;//進入局部最優
reset();//初始化agent資訊
construct_solution();//對于所有的agent構造一個完整的tour
update_pheromone();//更新資訊素
printf("iteration %d:best_so_far = %.2f\n", iteration, best_so_far.L);
if (last > best_so_far.L)
last = best_so_far.L, bad_times = 0;
else bad_times ++;//記錄目前未更新代數,若疊代多次未更新,認為進入局部最優
}
printf("best_so_far = %.2f\n", best_so_far.L);// 輸出目标值
best_so_far.print(stdout);//輸出路徑
}
算例示範
例一 滿秩矩陣式(type = 1)
輸入檔案格式為:
File_name File_type
salesman.in 1
5
0 1 2 2 3
2 0 3 4 2
3 2 0 4 1
3 4 5 0 5
2 4 1 4 0
輸出結果為:
opt_solution:
11
例二 二維坐标式(type = 2)
輸入檔案格式為:
File_name File_type
KroA100.tsp 2
100
1 1380 939
2 2848 96
3 3510 1671
4 457 334
5 3888 666
6 984 965
7 2721 1482
8 1286 525
9 2716 1432
10 738 1325
11 1251 1832
12 2728 1698
13 3815 169
14 3683 1533
15 1247 1945
16 123 862
17 1234 1946
18 252 1240
19 611 673
20 2576 1676
21 928 1700
22 53 857
23 1807 1711
24 274 1420
25 2574 946
26 178 24
27 2678 1825
28 1795 962
29 3384 1498
30 3520 1079
31 1256 61
32 1424 1728
33 3913 192
34 3085 1528
35 2573 1969
36 463 1670
37 3875 598
38 298 1513
39 3479 821
40 2542 236
41 3955 1743
42 1323 280
43 3447 1830
44 2936 337
45 1621 1830
46 3373 1646
47 1393 1368
48 3874 1318
49 938 955
50 3022 474
51 2482 1183
52 3854 923
53 376 825
54 2519 135
55 2945 1622
56 953 268
57 2628 1479
58 2097 981
59 890 1846
60 2139 1806
61 2421 1007
62 2290 1810
63 1115 1052
64 2588 302
65 327 265
66 241 341
67 1917 687
68 2991 792
69 2573 599
70 19 674
71 3911 1673
72 872 1559
73 2863 558
74 929 1766
75 839 620
76 3893 102
77 2178 1619
78 3822 899
79 378 1048
80 1178 100
81 2599 901
82 3416 143
83 2961 1605
84 611 1384
85 3113 885
86 2597 1830
87 2586 1286
88 161 906
89 1429 134
90 742 1025
91 1625 1651
92 1187 706
93 1787 1009
94 22 987
95 3640 43
96 3756 882
97 776 392
98 1724 1642
99 198 1810
100 3950 1558
輸出結果為:
best_known_solution: 21282
View Code
參考連結:
1. 百度百科-蟻群算法
2. 掘金-10分鐘搞懂蟻群算法
3. scikit-opt官方文檔-ACA部分
4. CSDN-fanxin_i-蟻群算法原理及實作(Python)
5. 知乎-資料魔法師-十分鐘快速get蟻群算法(附C代碼)