蟻群演算法簡介
自1991年由義大利學者 M. Dorigo,V. Maniezzo 和 A. Colorni 通過模擬蟻群覓食行為提出了一種基於種群的模擬進化演算法——蟻群優化。該演算法的出現引起了學者們的極大關注,蟻群演算法的特點:
① 其原理是一種正反饋機制或稱增強型學習系統; 它通過【最優路徑上螞蟻數量的增加→資訊素強度增加→後來螞蟻選擇概率增大→最優路徑上螞蟻數量更大增加】達到最終收斂於最優路徑上L
② 它是一種通用型隨機優化方法, 它吸收了螞蟻的行為特(內在搜尋機制) , 它是使用人工螞蟻模擬(也稱螞蟻系統) 來求解問題L但人工螞蟻決不是對實際螞蟻的一種簡單模擬, 它融進了人類的智慧L人工螞蟻有一定的記憶; 人工螞蟻不完全是瞎的; 人工螞蟻生活的時空是離散的L
③ 它是一種分散式的優化方法, 不僅適合目前的序列計算機, 而且適合未來的平行計算機L
④ 它是一種全域性優化的方法, 不僅可用於求解單目標優化問題, 而且可用於求解多目標優化問題L
⑤ 它是一種啟發式演算法, 計算複雜性為o (Nc*n2*m) , 其中Nc 是迭代次數, m 是螞蟻數目, n 是目的節點數目L
蟻群發現最短路徑的原理和機制[1]
下面用圖 1解釋蟻群發現最短路徑的原理和機制。
如圖 1(a)所示,在蟻巢和食物源之間有兩條道路 Nest-A-B-D-Food 和Nest-A-C-D-Food,其長度分別為 4 和 6。單位時間內螞蟻可移動一個單位長度的距離。開始時所有路徑上都沒有外激素。
如圖 1(b),在 t=0 時刻,20 只螞蟻從蟻巢出發移動到 A。由於路徑上沒有外激素,它們以相同概率選擇左側或右側道路,因此平均有 10 只螞蟻走左側,另外 10 只走右側。
如圖 1(c),在 t=4 時刻,第一組先到達食物源的螞蟻將折回。
如圖 1(d),在 t=5 時刻,兩組螞蟻將在 D 點相遇。此時 BD 上的外激素數量與 CD 上的相同,因此返回的 10 只螞蟻中有 5 只選擇 BD 而另 5 只選擇 CD。
如圖 1(e),在 t=8 時刻,前 5 個螞蟻將返回巢穴,而在 AC、CD 和 AB 上各有 5 個螞蟻。
如圖 1(f),在 t=9 時刻,前 5 個螞蟻又回到 A 並且再次面對往左還是往右的選擇。這時,AB 上的軌跡數是 20 而 AC 上是 15,因此將有較為多數的螞蟻選擇往右,從而增強了 AB 上外激素的量。隨著該過程的繼續,兩條道路上外激素數量的差距將越來越大,直至絕大多數螞蟻都選擇了最短的路徑。正是由於一條道路要比另一條道路短,因此,在相同的時間間隔內,短的路線會有更多的機會被選擇。
根據仿生學家的研究結果,螞蟻憑藉路徑尋優的能力能夠找到蟻巢與食物之間的最短路徑,其原理在於:螞蟻在所經過的路徑上留下一種揮發性分泌物(pheromone,以下稱為資訊素),資訊素隨著時間的推移會逐漸揮發消失.螞蟻在覓食過程中能夠感知這種物質的存在及其強度,並以此來指導自己的運動方向,傾向於朝著這種物質強度高的方向移動,即選擇該路徑的概率與當時這條路徑上該物質的強度成正比.資訊素強度越高的路徑,選擇它的螞蟻就越多,則在該路徑上留下的資訊素的強度就更大,而強度大的資訊素又吸引更多的螞蟻,從而形成一種正反饋.通過這種正反饋,螞蟻最終可以發現最佳路徑,導致大部分的螞蟻都會走此路徑.
以求解n個城市的TSP旅行商問題為例說明ACA模型.
設蟻群中螞蟻的數量為m,dij (i,j=1,2,…,n)表示城市i和城市j之間的距離,bi(t)表示t時刻位於城市i的螞蟻的個數,則有 表示t時刻在城市i,j連線上殘留的資訊量.初始時刻,各條路徑上資訊量相等,設τij(0)=C(C為常數).螞蟻k(k=1,2,…,m)在運動過程中,根據各條路徑上的資訊量決定轉移方向. 表示在t時刻螞蟻k由城市i轉移到城市j的概率.
(1)
殘留資訊的重要程度;β——啟發資訊的重要程度;tabuk——記錄螞蟻k當前所走過的城市,稱為記憶列表,k=1,2,…,m,集合tabuk隨著進化過程作動態調整.經過n個時刻,所有螞蟻都完成了一次遍歷.此時,計算每一隻螞蟻所走過的路徑Lk,並儲存最短路徑Lmin=min{Lk︱k=1,2,…,m}.在螞蟻完成一次迴圈以後,各路徑上的資訊量進行如下調整
τij(t+1)=(1-ρ)τij(t)+Δτij (2)
式中ρ∈(0,1),表示資訊素τij(t)隨時間的推移而衰減的程度.所以1-ρ為資訊素殘留因子,開始時Δτij(0)=0,
資訊素增量Δτij可表示 (3)
式中Δτkij為螞蟻k在本次迴圈中在城市i和j之間留下的資訊量,它的計算公式根據具體問題而定.Dorigo曾給出Δτkij3種不同的模型,分別稱為Ant-Cycle模型、Ant-Quantity模型、Ant-Density模型,它們的區別就在於資訊素的更新機制,即其差別在於Δτkij
在Ant-Cycle模型中:
(4) 式中,Q表示資訊素強度,它在一定程度上影響演算法的收斂速度;Lk表示第K只螞蟻在本次迴圈中所奏路徑的總長度。
在Ant-Quantity模型中:
(5) 式中,Q表示資訊素強度,它在一定程度上影響演算法的收斂速度;dij表示第K只螞蟻在t和t+1之間經過的( i, j )
在Ant-Density模型中:
(6) 區別:式(5)式(6)中利用的是區域性資訊,即螞蟻完成一步後更新路徑上的資訊素;而式(4)中利用的是整體資訊,即螞蟻完成一個迴圈後所有路徑上的資訊素。經過大量試驗總結研究,採用式(4)效能較好,所以 Ant-Cycle模型是最優的。
以上說明了資訊素殘留因子1-ρ、資訊啟發式因子α、期望啟發式因子β、資訊素強度Q、螞蟻數目M等都是非常重要的引數,其選區方式和選區原則直接影響到蟻群演算法的全域性收斂性和求解效率。我們學習到這種“三步走”[2]選擇蟻群演算法最優組合引數的有效方法:
(1) 確定螞蟻數目M,根據 城市規模 / 螞蟻數目 ≈1.5的選擇策略來確定螞蟻的總數目。
(2) 引數粗調,即調整數值範圍較大的資訊啟發式因子α、期望啟發式因子β、資訊素強度Q等引數,已得到較理想的解。
(3) 引數微調,即調整數值範圍較小的資訊素殘留因子1-ρ。
2 目前蟻群演算法的應用
雖然對蟻群演算法的研究時間不長, 但是初步研究已顯示出它在求解複雜優化問題方面具有很大的優勢, 特別是1998 年在比利時布魯塞爾專門召開了第一屆螞蟻優化國際研討會後, 現在每兩年召開一次這樣的螞蟻優化國際研討會。這標誌著蟻群演算法的研究已經得到了國際上的廣泛支援,使得這種新興的智慧進化仿生演算法展現出了勃勃生機[3]。
以蟻群演算法為代表的群體智慧已成為當今分散式人工智慧研究的一個熱點,許多源於蜂群和蟻群模型設計的演算法已越來越多地被用於企業的運轉模式的研究。美國五角大樓正在資助關於群體智慧系統的研究工作--群體戰略(SWARM STRATEGY),它的一個實戰用途是通過運用成群的空中無人駕駛飛行器和地面車輛來轉移敵人的注意力,讓自己的軍隊在敵人後方不被察覺地安全行進。英國電信公司和美國世界通訊公司以電子螞蟻為基礎,對新的電信網路管理方法進行了試驗。群體智慧還被應用於工廠生產計劃的制定和運輸部門的後勤管理。美國太平洋西南航空公司採用了一種直接源於螞蟻行為研究成果的運輸管理軟體,結果每年至少節約了1000萬美元費用開支。英國聯合利華公司已率先利用群體智慧技術改善其一家牙膏廠的運轉狀況。美國通用汽車公司,法國液氣公司,荷蘭公路交通部和美國一些移民事務機構也都採用這種技術來改善其運轉的機能。又如美國MCIWorld.com公司一直研究人工螞蟻,並用於管理公司的電話網,對使用者記賬收費等工作。另外,還設計“人工螞蟻”打算用於因特網的路由管理。鑑於群體智慧廣闊的應用前景,美國和歐洲聯盟均於近幾年開始出資資助基於群體智慧模擬的相關研究專案, 關在一些院校開設群體智慧的相關課程.牛津大學出版社1999年版的E.Bonabeau和M.Dorigo等人編寫的專著《群體智慧:從自然到人工系統》(Swarm Intelligence:From Natural to Artificial System),以及2001年出版的J.Kennedy和R.Eberhart編著的《群體智慧》(Swarm Intelligence)進一步擴大了群體智慧的影響.IEEE進化計算會刊也於2002年8月出版了蟻群優化算特刊。國內也有研究者用螞蟻演算法求解全國144個城市的最短迴路問題,求得的解同其它方法求到得解一樣精確,這說明螞蟻演算法不但是求解組合優化問題的可行方法,而且是一種很有競爭力的演算法。國家自然科學基金"十五"期間學科交叉類優先資助領域中的認知科學及其資訊處理的研究內容中也明確列出了群體智慧領域的進化,自適應與現場認知主題[4]。而且從1999年開始,幾乎每年都會有幾項相關專案獲得資助。蟻群演算法是一種新型的模擬進化演算法,其在資料探勘中的應用正逐步引起人們的關注。目前,人工蟻群在知識發現的過程中主要用於發掘聚類模型和分類模型。
2.1蟻群演算法在資料探勘中的應用
聚類是將一組物件分成若干個群體,每個群體構成一個簇,使得簇內的物件儘可能具有最大的相似性,不同簇之間的物件儘可能有最大的相異性。目前,聚類方法主要有K均值法,模糊聚類、神經網路聚類、基於遺傳演算法的聚類、小波變換聚類以及將這些演算法有效結合而形成的改進方法。隨著蟻群演算法研究的興起,人們發現在某些方面採用蟻群模型進行聚類更加接近實際的聚類問題。將蟻群演算法用於聚類分析,靈感源於螞蟻堆積他們的屍體和分類他們的幼體。基於蟻群演算法的聚類方法從原理上可分為兩種:一種是基於蟻堆形成原理來實現資料聚類,另一種是運用螞蟻覓食的原理,利用資訊來實現聚類分析。
而資料是資料探勘的另一個重要主題,它是在資料庫物件集合中尋找屬性,並根據分類模式將其劃分為不同類別的過程。分類過程利用歷史資料記錄自動推匯出對給定資料的分類樹。分類器構造方法有統計學方法、機器學習法、神經網路、決策樹等。從知識發現的觀點來看,分類規則的表達方式形如if<條件>then<類>規則前件(if 部分)包含一組條件集合,一般由邏輯連線符連線;規則結論(then部分)定義了樣本的預測類,這些樣本的預測屬性滿足規則前件所定義的所有條件[5]。將蟻群演算法引入分類規則的發現,是利用蟻群覓食原理在資料庫中進行搜尋,對隨機產生的一組規則進行選擇優化,直到資料庫能被該組規則覆蓋,從而挖掘出隱含在資料庫中的規則,建立最優的分類模型。蟻群演算法搜尋的初始條件為發現規則的集合為空,且訓練集包含所有的訓練樣本。螞蟻搜尋一次要完成規則生成、規則剪枝、資訊素更新三個任務。一次搜尋生成一條規則,並且將這條規則加入發現規則集合,同時將該條規則所覆蓋的訓練樣本從訓練集中刪除。如果未覆蓋訓練樣本的數目大於使用者定義的閾值,即最大未覆蓋樣本數,就反覆執行上述過程,最終演算法將得到一組最優分類規則集合[5]。
最早在這一領域開展工作的是Deneubourg 等[6],他們根據資料物件與其周圍物件的相似性,讓螞蟻隨機地移動、拾起或放下資料物件,以達到聚類資料的目的,這個基本模型已成功地應用於機器人領域。Lumer 等首先改進此演算法,提出了LF演算法。Wu 等、Ramos等、Yang等[7]從不同角度對LF演算法進行了改進,在用蟻群演算法進行聚類分析方面取得了一定成效。近幾年,學者在這方面的研究從來沒有間斷過,也取得了一定的研究成果。
2.2 結論
不過,將蟻群演算法運用於資料發掘還存在一些問題,需要進一步研究:
(1)如何將現實的挖掘任務轉換成蟻群求解的問題空間,並用適當的方式表達。如何定義“人工螞蟻”以及螞蟻間的非直接通訊方式(如路徑上的資訊素、物件的分佈狀態等)的選擇。
(2)如何建立正反饋機制,定義啟發函式,遞增地進行問題求解,並且使得到的解與問題定義中現實世界的情況相對應。
(3)基於蟻群的演算法要初始化大量的引數,這些引數的選擇會對演算法的效能產生較大的影響,但其選取的方法和原則目前尚無理論上的依據,只能通過多次實驗調優,因此引數的最佳設定原則還有待進一步研究。
(4)蟻群演算法的搜尋時間較長,如何將蟻群演算法與遺傳演算法、免疫演算法等優化演算法相結合,改善和提高演算法效能,以適應海量資料庫的知識發現。
所以如何在資料探勘中運用蟻群演算法快速、高效地獲得高質量的知識越來越受到人們的關注,逐漸成為近期的研究熱點[5]。
以下是解放軍資訊工程大學一個老師編的matlab程式,請尊重原作者勞動,引用時請註明出處。
我經過修改增加了註釋,已經執行過,無誤,
function [R_best,L_best,L_ave,Shortest_Route,Shortest_Length]=ACATSP(C,NC_max,m,Alpha,Beta,Rho,Q)
%%-------------------------------------------------------------------------
%% 主要符號說明
%% C n個城市的座標,n×2的矩陣
%% NC_max 最大迭代次數
%% m 螞蟻個數
%% Alpha 表徵資訊素重要程度的引數
%% Beta 表徵啟發式因子重要程度的引數
%% Rho 資訊素蒸發係數
%% Q 資訊素增加強度係數
%% R_best 各代最佳路線
%% L_best 各代最佳路線的長度
%%=========================================================================
%%第一步:變數初始化
n=size(C,1);%n表示問題的規模(城市個數)
D=zeros(n,n);%D表示完全圖的賦權鄰接矩陣
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時不計算,應該為0,但後面的啟發因子要取倒數,用eps(浮點相對精度)表示
end
D(j,i)=D(i,j); %對稱矩陣
end
end
Eta=1./D; %Eta為啟發因子,這裡設為距離的倒數
Tau=ones(n,n); %Tau為資訊素矩陣
Tabu=zeros(m,n); %儲存並記錄路徑的生成
NC=1; %迭代計數器,記錄迭代次數
R_best=zeros(NC_max,n); %各代最佳路線
L_best=inf.*ones(NC_max,1); %各代最佳路線的長度
L_ave=zeros(NC_max,1); %各代路線的平均長度
while NC<=NC_max %停止條件之一:達到最大迭代次數,停止
%%第二步:將m只螞蟻放到n個城市上
Randpos=[]; %隨即存取
for i=1:(ceil(m/n))
Randpos=[Randpos,randperm(n)];
end
Tabu(:,1)=(Randpos(1,1:m))'; %此句不太理解?
%%第三步:m只螞蟻按概率函式選擇下一座城市,完成各自的周遊
for j=2:n %所在城市不計算
for i=1:m
visited=Tabu(i,1:(j-1)); %記錄已訪問的城市,避免重複訪問
J=zeros(1,(n-j+1)); %待訪問的城市
P=J; %待訪問城市的選擇概率分佈
Jc=1;
for k=1:n
if length(find(visited==k))==0 %開始時置0
J(Jc)=k;
Jc=Jc+1; %訪問的城市個數自加1
end
end
%下面計算待選城市的概率分佈
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); %若計算的概率大於原來的就選擇這條路線
to_visit=J(Select(1));
Tabu(i,j)=to_visit;
end
end
if NC>=2
Tabu(1,:)=R_best(NC-1,:);
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)); %一輪下來後走過的距離
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 %迭代繼續
%%第五步:更新資訊素
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);
%此次迴圈在路徑(i,j)上的資訊素增量
end
Delta_Tau(Tabu(i,n),Tabu(i,1))=Delta_Tau(Tabu(i,n),Tabu(i,1))+Q/L(i);
%此次迴圈在整個路徑上的資訊素增量
end
Tau=(1-Rho).*Tau+Delta_Tau; %考慮資訊素揮發,更新後的資訊素
%%第六步:禁忌表清零
Tabu=zeros(m,n); %%直到最大迭代次數
end
%%第七步:輸出結果
Pos=find(L_best==min(L_best)); %找到最佳路徑(非0為真)
Shortest_Route=R_best(Pos(1),:) %最大迭代次數後最佳路徑
Shortest_Length=L_best(Pos(1)) %最大迭代次數後最短距離
subplot(1,2,1) %繪製第一個子圖形
DrawRoute(C,Shortest_Route) %畫路線圖的子函式
subplot(1,2,2) %繪製第二個子圖形
plot(L_best)
hold on %保持圖形
plot(L_ave,'r')
title('平均距離和最短距離') %標題
function DrawRoute(C,R)
%%=========================================================================
%% DrawRoute.m
%% 畫路線圖的子函式
%%-------------------------------------------------------------------------
%% C Coordinate 節點座標,由一個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')
hold on
for ii=2:N
plot([C(R(ii-1),1),C(R(ii),1)],[C(R(ii-1),2),C(R(ii),2)],'g')
hold on
end
title('旅行商問題優化結果 ')