螢火蟲算法

1.1 螢火蟲的生物現象

在熱帶的夏夜中，螢火蟲會聚集在一起產生短暫而有節奏的光，不同種類的螢火蟲的閃光模式往往是不同的。這種閃光的基本功能有三個：①是吸引交配的異性伙伴 ②則是吸引潛在的獵物 ③保護性的警告機制，告訴捕獵者其毒性或苦味

在光源遠處r距離的光強服從平方反比定律: 光強I隨著距離r的增加而減小，I∝1/r^2
并且，空氣也會吸收光，隨著距離增大，光會變得越來越弱。

這兩個綜合的因素使大多數螢火蟲視覺到一個有限的距離，通常在晚上幾百米，但這也足以讓螢火蟲交流。

閃爍光的表述可以與要優化的目標函數相關聯，這使得制定新的優化算法成為可能。

2.1 設計螢火蟲算法

我們可以由螢火蟲的一些閃爍特性，從而開發螢火蟲啟發的算法。為了簡單地描述螢火蟲算法(FA)，現在使用以下三個理想的規則：
①所有的螢火蟲雌雄同體的，所以一只螢火蟲有可能會被其他任意的螢火蟲所吸引（只要這只螢火蟲光強更大），無論它們的性別。
②吸引力與它們的亮度成正比，因此對于任何兩只閃爍的螢火蟲，越不亮的螢火蟲就會向越亮的螢火蟲移動。吸引力與亮度成正比，它們都隨著距離的增加而減小。如果沒有比某只螢火蟲更亮的螢火蟲，它會隨機移動；
③螢火蟲的亮度受到目標函數的影響。也就是說，螢火蟲本身的亮度是由目標函數決定的，而螢火蟲對別的螢火蟲的吸引是由別的螢火蟲看到他的光強決定的，這個光強會隨著距離增大而減少。

2.2 偽代碼 Pseudo code of the firefly algorithm (FA)

基于以上三條規則
螢火蟲算法(FA)的基本步驟可以總結為如圖所示的偽代碼。

3.1 光的強度和吸引力

在螢火蟲算法中，有兩個重要的問題：光強和吸引力。為簡單起見，我們可以假設螢火蟲的吸引力是由它的亮度決定的，而其亮度是由目標函數決定的。

螢火蟲在特定位置x的亮度I可以選擇為I(x)∝f(x)。然而，吸引力是相對的，它應該由其他螢火蟲判斷，并且隨著螢火蟲i和螢火蟲j之間的距離的變化而變化。距離越大，光強被空氣吸收的越多，吸引力也因此而減少。

①平方反比定律和吸收的聯合效應可以近似為以下高斯形式：

其中 i₀是初始的光強，γ為固定光吸收系數，r為兩個螢火蟲之間的距離。

②由于螢火蟲的吸引力與相鄰螢火蟲看到的光強度成正比，我們現在可以定義螢火蟲的吸引力β為

其中，β₀是r=0的吸引力，β₀一般取為1。任意兩個螢火蟲i和j在xi和xj之間的距離分別為笛卡爾距離。

③并且為了增大跳出局部最優的能力，還加入了隨機性因素

其中α為隨機化參數，α ∈ [0, 1]， 是由高斯分布或均勻分布中抽取的隨機數向量，在最簡單的二維空間下，使用rand ? 1/2進行表示，即此向量的每個分量是[-0.5,0.5]之間的隨機數。

螢火蟲i被另一個更有吸引力（更亮）的螢火蟲j所吸引，i向j移動的過程綜合上述三條，表示為以下公式：
其中，坐標xj-坐標xi意味著i的移動方向為j。

值得指出的是，這是一種偏向于向著更明亮的螢火蟲的隨機行走。如果β₀=0，它將變成一個簡單的隨機游走。此外，隨機化項可以很容易地擴展到其他分布，如L‘evy flight。

參數γ現在表征了吸引力的變化，它的值在決定收斂速度和FA算法的行為方面至關重要。理論上，γ∈[0，∞)，但在實踐中，對于大多數應用程序，它通常從0.1到10不等。

需要注意的是：在算法的實際代碼中，螢火蟲之間的吸引并不是只有最亮那個才能對別的螢火蟲產生吸引，而是對于螢火蟲i而言，其他所有的螢火蟲都會比較自己與螢火蟲i的光強，如果i更亮，那么別的螢火蟲再去計算他們之間的吸引力，然后根據吸引力向著螢火蟲i進行移動，也就是說，只有一個螢火蟲的光強大于別的任何螢火蟲，就會對別的螢火蟲產生吸引，微觀上來看是最亮那個吸引力最強而已。

4.1 參數的理解

當γ→0和γ→∞時分別有兩種重要的極限情況。

①對于γ→0來說：吸引力是恒定的β=β0這相當于說在一個理想化的天空中，光的強度不會降低。因此，在該領域的任何地方都可以看到一只閃爍的螢火蟲。因此，可以很容易地達到一個單一的（通常是全局的）最佳狀態。如果我們去掉FA偽代碼中j的內環（inner loop），并將xj替換為當前的全局最佳g?，那么螢火蟲算法就成為粒子群優化算法(PSO)的特殊情況。因此，這種特殊情況的效率與PSO相同。

② 極限情況γ→∞ 這意味著在其他螢火蟲的視線中，吸引力幾乎為零。這相當于螢火蟲在一個非常厚的有霧的區域隨機漫游的情況。看不到其他的螢火蟲，每只螢火蟲都以一種完全隨機的方式漫游。因此，這與完全隨機搜索方法相對應。

由于螢火蟲算法通常處于這兩個極端之間，因此可以調整參數γ和α，使其性能優于隨機搜索和PSO。實際上，FA可以同時有效地找到全局最優和局部最優。

5.1 改進措施

基礎的螢火蟲算法是非常有效的，但我們可以預料到，隨著迭代次數增加，生成的解逐漸接近最優解，但是解決方案仍在變化并且可能由于隨機項的突然增大反而偏離最優解。這種情況下，我們可以通過隨著迭代次數增加，隨機性減少這種方法來提高解決方案的質量。

有兩種方案可供選擇：
①

其中t_max為最大迭代次數。α₀為初始隨機化參數，α_∞為最終值。
很明顯，算法迭代的開始階段的α是接近于α₀的，之后便逐漸減少，直到趨近于α_∞

②

這是類似于模擬退火算法的方法，θ∈(0,1]是隨機性還原常數，也即是一個0~1的隨機數。
很明顯，算法迭代的開始階段的α是接近于α₀的，之后便逐漸減少，直到趨近于0。

6.1 實際應用及代碼

為了證明全局最優和局部最優可以同時找到，我們現在使用下面的四峰函數，在實現過程中，參數的值分別為α=0.2、γ=1和β₀=1，delta=0.97并且x、y的取值范圍都是[-5,5]

由Yang xin she教授給出的Matlab代碼：

% ======================================================== % % Files of the Matlab programs included in the book: % % Xin-She Yang, Nature-Inspired Metaheuristic Algorithms, % % Second Edition, Luniver Press, (2010). www.luniver.com % % ======================================================== % % =========================================================% % Firefly Algorithm by X S Yang (Cambridge University) % % Usage: firefly_simple([number_of_fireflies,MaxGeneration]) % eg: firefly_simple([12,50]); % % ======================================================== % % This is a demo for 2D functions; for higher dimenions, % % you should use fa_ndim.m or fa_mincon.m % % Parameters choice: % Gamma should be linked with scales. Otherwise, the FA % % the efficiency will be significantly reduced because % % the beta term may be too small. % % Similarly, alpha should also be linked with scales, % % the steps should not too large or too small, often % % steps are about 1/10 to 1/100 of the domain size. % % In addition, alpha should be reduced gradually % % using alpha=alpha_0 delta^t during eteration t. % % Typically, delta=0.9 to 0.99 will be a good choice. % % ======================================================== %function [best]=firefly_simple(instr) % n=number of fireflies % MaxGeneration=number of pseudo time steps if nargin<1 %%不知道這個nargin是干嘛的 也不知道是這個if是干嘛的，應該有無都無所謂 確定的是instr=[12 50]中12是螢火蟲個數，50是迭代次數instr=[12 50]; end n=instr(1); MaxGeneration=instr(2); % Show info help firefly_simple.m %%沒啥用 rand('state',0); % Reset the random generator 產生一個隨機數，"state"的意思有下標地創建隨機數，如果換一個下標比如從0換成1則重新隨機，如果不換下標則隨機數無論怎么rand都還是這個不會變化 % ------ Four peak functions --------------------- str1='exp(-(x-4)^2-(y-4)^2)+exp(-(x+4)^2-(y-4)^2)'; str2='+2*exp(-x^2-(y+4)^2)+2*exp(-x^2-y^2)'; funstr=strcat(str1,str2); % Converting to an inline function f=vectorize(inline(funstr));%%inline把串轉化為函數 vectorize將函數的標量運算轉換為函數運算：如矩陣相乘*轉換為.*點乘 ^2轉換為.^2 前者是矩陣乘矩陣 后者是對矩陣中每個元素進行平方 % range=[xmin xmax ymin ymax]; range=[-5 5 -5 5];% ------------------------------------------------ alpha=0.2; % Randomness 0--1 (highly random) gamma=1.0; % Absorption coefficient delta=0.97; % Randomness reduction (similar to % an annealing schedule) % ------------------------------------------------ % Grid values are used for display only Ngrid=100; dx=(range(2)-range(1))/Ngrid dy=(range(4)-range(3))/Ngrid [x,y]=meshgrid(range(1):dx:range(2),...range(3):dy:range(4));%%生成一個三維圖像 并對x、y軸-5~5進行分格 間隔為0.1 再加上零點所以一共101個點 z=f(x,y);%z為計算此四峰函數在不同分隔點上對應的函數值 % Display the shape of the objective function 作出此四峰函數的圖像 figure(1); surfc(x,y,z);% ------------------------------------------------ % generating the initial locations of n fireflies 對12個螢火蟲的初始位置進行初始化 [xn,yn,Lightn]=init_ffa(n,range);%%使用init_ffa函數進行初始化并返回值，xn為螢火蟲的x軸初始坐標、yn為螢火蟲的y軸初始坐標，lightn為 % Display the paths of fireflies in a figure with % contours of the function to be optimizedfigure(2); % Iterations or pseudo time marching for i=1:MaxGeneration, %%%%% start iterations 開始迭代 % Show the contours of the functioncontour(x,y,z,15); hold on;%創建等高線圖 15條等高線，函數x軸為x矩陣，y軸為y矩陣，函數值大小為z矩陣，顏色由藍變淺，越淺說明函數值越大 % Evaluate new solutions 計算xn，yn坐標下的函數值zn zn=f(xn,yn);% Ranking the fireflies by their light intensity 比較不同螢火蟲的亮度 [Lightn,Index]=sort(zn);%lightn中存儲的是zn進行升序排序之后的矩陣，Index存儲的是進行升序排序之后lightn中每個元素對應于排序前zn中每個元素的位置如A={19,23,12} 排序后lightn={12,19,23},Index={3,1,2} xn=xn(Index); yn=yn(Index);%排序之后xn=xn(Index); yn=yn(Index);意味著xn中存儲的是對應函數值升序排序的x坐標，比如第一個xn(1)存的就是函數值zn最小的那個螢火蟲的x坐標，y同理 xo=xn; yo=yn; Lighto=Lightn;%使用新的xo、yo、lightno存儲，作用是作為參數傳入ffa_move來進行螢火蟲每次迭代的吸引與隨機移動工作，移動之后的新坐標再傳回xn、yn % Trace the paths of all roaming fireflies plot(xn,yn,'.','markersize',10,'markerfacecolor','g');%畫螢火蟲 表示為一個小綠點 % Move all fireflies to the better locations [xn,yn]=ffa_move(xn,yn,Lightn,xo,yo,Lighto,alpha,gamma,range);%通過ffa_move函數進行吸引與隨機移動 drawnow; % Use "hold on" to show the paths of fireflieshold off;% Reduce randomness as iterations proceed 隨著迭代次數增加減少隨機性alpha alpha=newalpha(alpha,delta);end %%%%% end of iterations %迭代完成之后，獲取最終的答案 %best中存的就是螢火蟲（x,y）處對應的函數值第一列存的是x，第二列存的y，第三列存的函數值 best(:,1)=xo'; best(:,2)=yo'; best(:,3)=Lighto';end %%end for firfly_simple% ----- All subfunctions are listed here --------- % The initial locations of n fireflies 初始化螢火蟲坐標 function [xn,yn,Lightn]=init_ffa(n,range) xrange=range(2)-range(1); yrange=range(4)-range(3); xn=rand(1,n)*xrange+range(1); yn=rand(1,n)*yrange+range(3); Lightn=zeros(size(yn)); end% Move all fireflies toward brighter ones function [xn,yn]=ffa_move(xn,yn,Lightn,xo,yo,...Lighto,alpha,gamma,range) ni=size(yn,2); nj=size(yo,2); for i=1:ni, % The attractiveness parameter beta=exp(-gamma*r)for j=1:nj, r=sqrt((xn(i)-xo(j))^2+(yn(i)-yo(j))^2); if Lightn(i)<Lighto(j), % Brighter and more attractive beta0=1; beta=beta0*exp(-gamma*r.^2); xn(i)=xn(i).*(1-beta)+xo(j).*beta+alpha.*(rand-0.5); yn(i)=yn(i).*(1-beta)+yo(j).*beta+alpha.*(rand-0.5); endend % end for j end % end for i [xn,yn]=findrange(xn,yn,range); end% Reduce the randomness during iterations function alpha=newalpha(alpha,delta) alpha=alpha*delta; end % Make sure the fireflies are within the range function [xn,yn]=findrange(xn,yn,range) for i=1:length(yn),if xn(i)<=range(1), xn(i)=range(1); endif xn(i)>=range(2), xn(i)=range(2); endif yn(i)<=range(3), yn(i)=range(3); endif yn(i)>=range(4), yn(i)=range(4); end end end % ============== end =====================================

執行結果：

7.1 帶約束的優化問題——Spring Design Problem

一個帶約束實際優化問題 彈簧設計問題：拉開彈簧和壓縮彈簧的設計是一個著名的基準優化問題。其主要目的是根據偏轉、應力、喘振頻率和幾何形狀的限制來最小化重量。它涉及三個設計變量：線圈直徑x₁，線圈直徑x₂和線圈的數量/長度x₃。這個問題可以概括為
受以下約束條件的約束

MATLAB代碼：from Dr.Yang

% ======================================================== % % Files of the Matlab programs included in the book: % % Xin-She Yang, Nature-Inspired Metaheuristic Algorithms, % % Second Edition, Luniver Press, (2010). www.luniver.com % % ======================================================== % % -------------------------------------------------------- % % Firefly Algorithm for constrained optimization using % % for the design of a spring (benchmark) % % by Xin-She Yang (Cambridge University) Copyright @2009 % % -------------------------------------------------------- %function fa_mincon % parameters [n N_iteration alpha betamin gamma] para=[40 500 0.5 0.2 1]; format long;help fa_mincon.m % This demo uses the Firefly Algorithm to solve the % [Spring Design Problem as described by Cagnina et al., % Informatica, vol. 32, 319-326 (2008). ]% Simple bounds/limits disp('Solve the simple spring design problem ...'); Lb=[0.05 0.25 2.0]; Ub=[2.0 1.3 15.0];% Initial random guess 隨機猜一個解（x1，x2，x3） u0=Lb+(Ub-Lb).*rand(size(Lb));[u,fval,NumEval]=ffa_mincon(@cost,@constraint,u0,Lb,Ub,para);% Display results bestsolution=u bestojb=fval total_number_of_function_evaluations=NumEval%%% Put your own cost/objective function here --------%%% %% Cost or Objective functionfunction z=cost(x) z=(2+x(3))*x(1)^2*x(2);% Constrained optimization using penalty methods % by changing f to F=f+ \sum lam_j*g^2_j*H_j(g_j) % where H(g)=0 if g<=0 (true), =1 if g is false%%% Put your own constraints here --------------------%%% function [g,geq]=constraint(x) % All nonlinear inequality constraints should be here % If no inequality constraint at all, simple use g=[]; g(1)=1-x(2)^3*x(3)/(71785*x(1)^4); % There was a typo in Cagnina et al.'s paper, % the factor should 71785 insteady of 7178 ! tmpf=(4*x(2)^2-x(1)*x(2))/(12566*(x(2)*x(1)^3-x(1)^4)); g(2)=tmpf+1/(5108*x(1)^2)-1; g(3)=1-140.45*x(1)/(x(2)^2*x(3)); g(4)=x(1)+x(2)-1.5;% all nonlinear equality constraints should be here % If no equality constraint at all, put geq=[] as follows 如果沒有等式約束就用geq=[]; geq=[];%%% End of the part to be modified -------------------%%%%%% --------------------------------------------------%%% %%% Do not modify the following codes unless you want %%% %%% to improve its performance etc %%% % ------------------------------------------------------- % ===Start of the Firefly Algorithm Implementation ====== % Inputs: fhandle => @cost (your own cost function, % can be an external file ) % nonhandle => @constraint, all nonlinear constraints % can be an external file or a function % Lb = lower bounds/limits % Ub = upper bounds/limits % para == optional (to control the Firefly algorithm) % Outputs: nbest = the best solution found so far % fbest = the best objective value % NumEval = number of evaluations: n*MaxGeneration % Optional: % The alpha can be reduced (as to reduce the randomness) % ---------------------------------------------------------% Start FA function [nbest,fbest,NumEval]...=ffa_mincon(fhandle,nonhandle,u0, Lb, Ub, para) % Check input parameters (otherwise set as default values) %nargin為關鍵字：代表本函數輸入的參數個數 if nargin<6, para=[20 50 0.25 0.20 1]; end if nargin<5, Ub=[]; end if nargin<4, Lb=[]; end if nargin<3, disp('Usuage: FA_mincon(@cost, @constraint,u0,Lb,Ub,para)'); end% n=number of fireflies % MaxGeneration=number of pseudo time steps % ------------------------------------------------ % alpha=0.25; % Randomness 0--1 (highly random) % betamn=0.20; % minimum value of beta % gamma=1; % Absorption coefficient % ------------------------------------------------ n=para(1); MaxGeneration=para(2); alpha=para(3); betamin=para(4); gamma=para(5);% Total number of function evaluations NumEval=n*MaxGeneration;% Check if the upper bound & lower bound are the same size if length(Lb) ~=length(Ub),disp('Simple bounds/limits are improper!');return end% Calcualte dimension d=length(u0);% Initial values of an array zn=ones(n,1)*10^100; % ------------------------------------------------ % generating the initial locations of n fireflies [ns,Lightn]=init_ffa(n,d,Lb,Ub,u0);%ns為初始化后的各個螢火蟲坐標% Iterations or pseudo time marching for k=1:MaxGeneration, %%%%% start iterations% This line of reducing alpha is optionalalpha=alpha_new(alpha,MaxGeneration);% Evaluate new solutions (for all n fireflies) for i=1:n,zn(i)=Fun(fhandle,nonhandle,ns(i,:));Lightn(i)=zn(i); end% Ranking fireflies by their light intensity/objectives [Lightn,Index]=sort(zn);% lightn是排序后的函數值 ns_tmp=ns; for i=1:n,ns(i,:)=ns_tmp(Index(i),:);%對ns也按照lightn排序后的順序重新對應索引 end%% Find the current best nso=ns; Lighto=Lightn; nbest=ns(1,:); Lightbest=Lightn(1);%nbest為最佳的螢火蟲的坐標 Lightbest為最佳的螢火蟲的適應度函數值% For output only fbest=Lightbest;% Move all fireflies to the better locations [ns]=ffa_move(n,d,ns,Lightn,nso,Lighto,nbest,...Lightbest,alpha,betamin,gamma,Lb,Ub);7end %%%%% end of iterations% ------------------------------------------------------- % ----- All the subfunctions are listed here ------------ % The initial locations of n fireflies function [ns,Lightn]=init_ffa(n,d,Lb,Ub,u0)%n螢火蟲個數，d維度，u0隨機猜的一個解，ns初始解空間，lightn是排序后的函數值% if there are bounds/limits,有上下限就在上下限之間進行生成隨機解 if length(Lb)>0,for i=1:n,ns(i,:)=Lb+(Ub-Lb).*rand(1,d);end else% generate solutions around the random guess 無上下限就在隨機解附近去生成所以隨機解for i=1:n,ns(i,:)=u0+randn(1,d);end end% initial value before function evaluations Lightn=ones(n,1)*10^100;% Move all fireflies toward brighter ones function [ns]=ffa_move(n,d,ns,Lightn,nso,Lighto,...nbest,Lightbest,alpha,betamin,gamma,Lb,Ub) % Scaling of the system scale=abs(Ub-Lb);% Updating fireflies for i=1:n, % The attractiveness parameter beta=exp(-gamma*r)for j=1:n,r=sqrt(sum((ns(i,:)-ns(j,:)).^2));% Update moves if Lightn(i)>Lighto(j), % Brighter and more attractivebeta0=1; beta=(beta0-betamin)*exp(-gamma*r.^2)+betamin;%這里beta為什么是需要計算的呢，而且與firefly_simple不同？tmpf=alpha.*(rand(1,d)-0.5).*scale;ns(i,:)=ns(i,:).*(1-beta)+nso(j,:).*beta+tmpf;%移動endend % end for jend % end for i% Check if the updated solutions/locations are within limits [ns]=findlimits(n,ns,Lb,Ub);% This function is optional, as it is not in the original FA % The idea to reduce randomness is to increase the convergence, % however, if you reduce randomness too quickly, then premature % convergence can occur. So use with care. function alpha=alpha_new(alpha,NGen) % alpha_n=alpha_0(1-delta)^NGen=10^(-4); % alpha_0=0.9 delta=1-(10^(-4)/0.9)^(1/NGen); alpha=(1-delta)*alpha;% Make sure the fireflies are within the bounds/limits function [ns]=findlimits(n,ns,Lb,Ub) for i=1:n,% Apply the lower boundns_tmp=ns(i,:);I=ns_tmp<Lb;ns_tmp(I)=Lb(I);% Apply the upper boundsJ=ns_tmp>Ub;ns_tmp(J)=Ub(J);% Update this new movens(i,:)=ns_tmp; end% ----------------------------------------- % d-dimensional objective function function z=Fun(fhandle,nonhandle,u) % Objective z=fhandle(u);%直接求出適應度函數值% Apply nonlinear constraints by the penalty method % Z=f+sum_k=1^N lam_k g_k^2 *H(g_k) where lam_k >> 1 z=z+getnonlinear(nonhandle,u);%應用約束function Z=getnonlinear(nonhandle,u) %%這里有約束的情況其實就是把約束帶上，直接求整體的函數值，因為約束前面的系數極大10^15所以如果要整體最小的話，那么必然要讓約束趨近于0才可以，那么當整體最小的時候，約束也必然被滿足，這就是懲罰法 Z=0; % Penalty constant >> 1 lam=10^15; lameq=10^15; % Get nonlinear constraints [g,geq]=nonhandle(u)%g是非等式約束求出的約束值 geq是什么等式約束求出的約束值% Apply inequality constraints as a penalty function %%加上不等式約束條件 for k=1:length(g),Z=Z+ lam*g(k)^2*getH(g(k))%lam應該是對應的懲罰法中的ui（u），lameq對應vi（v） end % Apply equality constraints (when geq=[], length->0) %加上等式約束條件 for k=1:length(geq),Z=Z+lameq*geq(k)^2*geteqH(geq(k)); end% Test if inequalities hold % H(g) which is something like an index function function H=getH(g) if g<=0,H=0; elseH=1; end% Test if equalities hold function H=geteqH(g) if g==0,H=0; elseH=1; end %% ==== End of Firefly Algorithm implementation ======

總結

以上是生活随笔為你收集整理的《Nature-Inspired Metaheuristic Algorithms》——萤火虫算法 FIREFLY ALGORITHM的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。