時間：2019.8.12
老師：祁永強(qiáng)
內(nèi)容：蟻群算法基本原理
個性：“無中生有，自圓其說”，戴準(zhǔn)帽子（套框架），優(yōu)化算法必須能用matlab仿真，狠勁，煉獄式成長，書由厚到薄（按頁）由薄到厚（），撕書式成長，佯裝努力結(jié)果是沒有差別的，為人師在于喚人醒，“日鋸其半，萬日不竭”，20頁最佳
感想：從原理上了解了蟻群算法，所找的一份ACATSP代碼成功復(fù)現(xiàn)，并理解了全過程。唯獨(dú)在下一個訪問地點(diǎn)的概率計算與選擇上稍有存疑，寫在了最后。歡迎交流探討觀點(diǎn)，e-mail：03170908@cumt.edu.cn
1.概述概念
起源：意大利M.Dorigo等，通過模擬自然界螞蟻搜索路徑的行為，提出的新型模擬進(jìn)化算法。
原理：螞蟻在運(yùn)動過程中，能感知信息素，指導(dǎo)自己的運(yùn)動方向（正反饋）。信息素越多越好，為什么，因?yàn)槌跏紩r，信息素越多說明往返越多（同樣時間時），即路徑最短。
核心：信息素＋正反饋（所有群體算法的核心）
適用：求解TSP問題、分配問題、job-shop調(diào)度問題
優(yōu)勢：求解復(fù)雜優(yōu)化問題（特別是離散優(yōu)化問題）
比較優(yōu)勢：與大多數(shù)基于梯度的優(yōu)化算法想比
1.無集中控制約束，不會因個別的故障影響整個問題的求解，保證魯棒性
2.非直接信息交換，確保系統(tǒng)的擴(kuò)展性
3.并行分布式
4.對問題定義的連續(xù)性無特殊要求
5.算法實(shí)現(xiàn)簡單
2.蟻群算法基本思想
信息素要按一定比率揮發(fā)（避免局部最優(yōu)）
每只螞蟻的一部轉(zhuǎn)移概率由圖中的每條邊上的兩類參數(shù)決定：
1.信息素值也稱信息素痕跡
2.可見度，即先驗(yàn)值
信息素的更新方式：
1 揮發(fā) ：所有路徑上的信息素以一定比率進(jìn)行減少
2 增強(qiáng) ：給螞蟻評價，好螞蟻加強(qiáng)
基于圖的蟻群算法（GBAS）
3.基本蟻群算法
表示
時刻螞蟻
從城市i到城市
的轉(zhuǎn)移概率，計算公式如下：
：信息素
: 信息素增強(qiáng)（啟發(fā)式的體現(xiàn)，和距離負(fù)相關(guān)）
: 推薦1，5
信息素?fù)]發(fā)系數(shù)：推薦0.5
信息素增強(qiáng)=
，Q推薦1，L表示距離長度
規(guī)模和終止條件：
1 給定終止輪數(shù)
2 最優(yōu)解了連續(xù)次相同
3 給定目標(biāo)值
4.【進(jìn)階】收斂
馬于爾科夫收斂定義，依概率收斂于1
證明收斂，才能說算法是可以實(shí)現(xiàn)的
蟻群算法分支：MAX-MIN螞蟻，精英螞蟻，蟻周算法，蟻密算法，蟻量算法
5.代碼實(shí)例
蟻群算法針對TSP應(yīng)用（超詳細(xì)分析）
ant-colony-algorithm-for-TSP(MATLAB可直接運(yùn)行代碼)
原始出處：解放軍信息工程大學(xué)老師
下面代碼注釋為本人針對上面的資料進(jìn)行學(xué)習(xí)之后添加，歡迎交流。轉(zhuǎn)移矩陣部分務(wù)必參考上面的公式
function [R_best,L_best,L_ave,Shortest_Route,Shortest_Length]=ACATSP(C,NC_max,m,Alpha,Beta,Rho,Q)
%%-------------------------------------------------------------------------
%% 主要符號說明
%% C n個城市的坐標(biāo)，n×2的矩陣
%% NC_max 最大迭代次數(shù)
%% m 螞蟻個數(shù)
%% Alpha 表征信息素重要程度的參數(shù) 1
%% Beta 表征啟發(fā)式因子重要程度的參數(shù) 5
%% Rho 信息素蒸發(fā)系數(shù) 0.5
%% Q 信息素增加強(qiáng)度系數(shù) 1
%% R_best 各代最佳路線
%% L_best 各代最佳路線的長度
%%=========================================================================
%%第一步：變量初始化
n=size(C,1);%n表示問題的規(guī)模（城市個數(shù)）
D=zeros(n,n);%D表示完全圖的賦權(quán)鄰接矩陣
for i=1:n
for j=1:n
if i~=j
D(i,j)=((C(i,1)-C(j,1))^2+(C(i,2)-C(j,2))^2)^0.5;
else
D(i,j)=eps; %i=j時不計算，應(yīng)該為0，但后面的啟發(fā)因子要取倒數(shù)，用eps（浮點(diǎn)相對精度）表示
end
D(j,i)=D(i,j); %對稱矩陣
end
end
Eta=1./D; %Eta為啟發(fā)因子，這里設(shè)為距離的倒數(shù)
Tau=ones(n,n); %Tau為信息素矩陣
Tabu=zeros(m,n); %存儲并記錄路徑的生成
NC=1; %迭代計數(shù)器，記錄迭代次數(shù)
R_best=zeros(NC_max,n); %各代最佳路線
L_best=inf.*ones(NC_max,1); %各代最佳路線的長度，先賦值inf無窮大，后面更新
L_ave=zeros(NC_max,1); %各代路線的平均長度
while NC<=NC_max %停止條件之一：達(dá)到最大迭代次數(shù)，停止
%%第二步：將m只螞蟻放到n個城市上，全放完需要[m/n]輪
Randpos=[]; %#隨機(jī)存取，此處初始起點(diǎn)，可以固定設(shè)置
for i=1:(ceil(m/n))
Randpos=[Randpos,randperm(n)];%randperm產(chǎn)生1到n的隨機(jī)序列
end %使用這種隨機(jī)分配方法可以保證初始去往各個城市的螞蟻數(shù)量均勻，比全部隨機(jī)要更好
Tabu(:,1)=(Randpos(1,1:m))'; %取出前m個，賦給Tabu路徑矩陣，即給出了m只螞蟻的初始去往地點(diǎn)
%%第三步：m只螞蟻按概率函數(shù)選擇下一座城市，對應(yīng)公式中的轉(zhuǎn)移矩陣
for j=2:n %所在城市不計算，j表示即將訪問的下一座城市
for i=1:m
visited=Tabu(i,1:(j-1)); %記錄已訪問的城市，避免重復(fù)訪問
J=zeros(1,(n-j+1)); %待訪問的城市，總數(shù)n-j+1
P=J; %待訪問城市的選擇概率分布
Jc=1; %待訪問的城市個數(shù)，初始為1（用于修改矩陣J時做下標(biāo)）
for k=1:n
if isempty(find(visited==k, 1))%這里改善了尋找的性能，只要能找到一個等于K就非空
J(Jc)=k;
Jc=Jc+1; %待訪問的城市個數(shù)自加1
end
end %循環(huán)完成后將更新完待訪問矩陣J
%下面計算從當(dāng)前城市visted(end)轉(zhuǎn)移到待選城市k的概率分布，參照公式的分子
%信息素的作用也僅體現(xiàn)在這里，計算轉(zhuǎn)移矩陣概率
for k=1:length(J)
P(k)=(Tau(visited(end),J(k))^Alpha)*(Eta(visited(end),J(k))^Beta);
end
P=P/(sum(P));
%按概率原則選取下一個城市
Pcum=cumsum(P); %cumsum，元素向上（或向前）累加求和
Select=find(Pcum>=rand); %若計算的概率大于rand的就可以選擇這條路線
to_visit=J(Select(1)); %選擇可選城市列表Select中的第一個作為將訪問的下一個
Tabu(i,j)=to_visit;
end
end
if NC>=2
Tabu(1,:)=R_best(NC-1,:);
%NC=2開始，第一只螞蟻路徑使用上一輪的最佳路線，參與下面的信息素更新
%#所以其實(shí)上面i=1那只螞蟻的路徑推演是沒必要的，可以用if去優(yōu)化掉
end
%%第四步：記錄本次迭代最佳路線
L=zeros(m,1); %開始距離為0，m*1的列向量
for i=1:m
R=Tabu(i,:);
for j=1:(n-1)
L(i)=L(i)+D(R(j),R(j+1)); %原距離加上第j個城市到第j+1個城市的距離
end
L(i)=L(i)+D(R(1),R(n)); %#一輪下來后走過的距離，此處假定了最終回到起點(diǎn)，可以修改
end
L_best(NC)=min(L); %最佳距離取最小
pos=find(L==L_best(NC));
R_best(NC,:)=Tabu(pos(1),:); %此輪迭代后的最佳路線
L_ave(NC)=mean(L); %此輪迭代后的平均距離
NC=NC+1; %迭代繼續(xù)
%%第五步：更新信息素
Delta_Tau=zeros(n,n); %開始時信息素增量為n*n的0矩陣
for i=1:m
for j=1:(n-1)
Delta_Tau(Tabu(i,j),Tabu(i,j+1))=Delta_Tau(Tabu(i,j),Tabu(i,j+1))+Q/L(i);
%此次循環(huán)在路徑（i，j）上的信息素增量
end
Delta_Tau(Tabu(i,n),Tabu(i,1))=Delta_Tau(Tabu(i,n),Tabu(i,1))+Q/L(i);
%#此次循環(huán)在從n會到起點(diǎn)上的信息素增量，對應(yīng)要求，單程可刪去
end
Tau=(1-Rho).*Tau+Delta_Tau; %考慮信息素?fù)]發(fā)掉Rho，更新后的信息素,加上Delta
%%第六步：路徑表清零
Tabu=zeros(m,n); %%直到最大迭代次數(shù)
end
%%第七步：輸出結(jié)果
Pos=find(L_best==min(L_best)); %找到最佳路徑（非0為真）
Shortest_Route=R_best(Pos(1),:) %最大迭代次數(shù)后最佳路徑
Shortest_Length=L_best(Pos(1)) %最大迭代次數(shù)后最短距離
subplot(1,2,1) %繪制第一個子圖形
DrawRoute(C,Shortest_Route) %畫路線圖的子函數(shù)
subplot(1,2,2) %繪制第二個子圖形
plot(L_best)
hold on %保持圖形
plot(L_ave,'r')
title('平均距離和最短距離') %標(biāo)題
function DrawRoute(C,R)
%%=========================================================================
%% DrawRoute.m
%% 畫路線圖的子函數(shù)
%% C Coordinate 節(jié)點(diǎn)坐標(biāo)，由一個N×2的矩陣存儲
%% R Route 路線
%%=========================================================================
N=length(R);
scatter(C(:,1),C(:,2));
hold on
plot([C(R(1),1),C(R(N),1)],[C(R(1),2),C(R(N),2)],'g')%#這里繪制了起點(diǎn)到終點(diǎn)的直連線，非循環(huán)路徑可刪去
hold on
for i=2:N
plot([C(R(i-1),1),C(R(i),1)],[C(R(i-1),2),C(R(i),2)],'g')
hold on
end
title('旅行商問題優(yōu)化結(jié)果 ')
關(guān)于代碼的一點(diǎn)說明，結(jié)合蟻群算法的原理和主要的轉(zhuǎn)移概率公式，理解代碼不難。上面的代碼是針對完全圖的回路所計算的，其中有多處可以加以改動，在注釋號之后加上了#號，體現(xiàn)我的一點(diǎn)個人見解。
此外，一個困惑：上面在計算待訪問城市之后，關(guān)于下一個城市的選擇，使用的是向前累加，取判斷大于某個隨機(jī)值rand，這樣體現(xiàn)的是輪盤賭算法嗎？
