智能優化算法——蝙蝠算法（PythonMatlab實現）
?錄
1 前?
蝙蝠算法是2010年楊教授基于群體智能提出的啟發式搜索算法，是?種搜索全局最優解的有效?法。該算法基于迭代優化，初始化為?組隨機解，然后迭代搜尋最優解，且在最優解周圍通過隨機飛?產?局部新解，加強局部搜索速度。該算法具有實現簡單、參數少等特點。
該算法主要?于?標函數尋優，基于蝙蝠種群利?產?的聲波搜索獵物和控制飛??向的特征來實現函數的尋優。以?只蝙蝠作為基本單元，且每只蝙蝠都有?個適應值來對函數解空間進?優化。每只蝙蝠可以調
整??發射聲波的響度、頻率等對空間進?搜索，使整個種群的活動逐步由?序變為有序。但蝙蝠算法在尋優末段容易陷?局部的極值，本?引?變速度權重因?修正系數(和粒?群參數調整的?適應過程相似)，盡可能避免局部極值的困境，從?達到全局尋優的效果。
2 蝙蝠算法原理細講
?先在?變量范圍內產?蝙蝠的隨機位置；
然后每個蝙蝠向周圍某處移動，飛?靠聲波反饋來更變?向，?向周圍移動是隨機的，蝙蝠下?時刻移動到某?位置有?定的出現頻率，所以本?在運動公式中加了聲波和頻率兩個因素。每?時刻的位移看作是?次飛?，每次飛?的距離很短（距離長短反映出搜素精度）。
每只蝙蝠初始位置是隨機的，在初始位置時刻其中?只蝙蝠對應的函數最值，算法開始會記錄該位置，然后所有蝙蝠逐漸向該位置靠近，飛??向?致向當前最值?向（即該?向的位置該蝙蝠的出現頻率更?），但是在飛??向也是隨機飛?的，相當于逐步搜索過去。蝙蝠飛?過程中不會因為當前飛?出現的最?值位置?改變種群的飛?趨勢，這樣可以盡量避免算法陷?局部極值。
算法每飛??次就記錄?次種群所處位置的最值，最后找出記錄中的最值。
算法有兩個參數可以影響最終的結果：種群數量和飛?次數。其中種群數量的影響是最?的。
3 詳細步驟
3.1 初始化相關參數
蝙蝠的位置為Xi，飛?速度Vi，聲?響度為Ai，頻率yi范圍，設有?標函數為 ：
w(t)其中為時刻變速慣性權重因?，作?是使蝙蝠的前期搜索對后期搜索提供參照，wmax為w(t)的最?值、wmin為w(t)的最?值；
，?般取2，Tmax為最?迭代次數；為當前位置最優解；y(i)為頻率滿?正態均勻分布的?個隨機數，β是?個隨機變量，且。開始運?時，蝙蝠在隨機進?頻率分配。
3.4 通過蝙蝠多次飛?產?多個新解，進?全局搜索，若得到的新解
3.5 排列所有蝙蝠的位置，并找出當前最優值及對應的位置
4 Python實現
4.1 代碼
#=========導?相關庫===============
import numpy as np
from numpy.random import random as rand
#========參數設置==============
# objfun:?標函數
# N_pop: 種群規模，通常為10到40
# N_gen: 迭代數
# A: 響度（恒定或降低）
# r: 脈沖率（恒定或減?）
# 此頻率范圍決定范圍
# 如有必要，應更改這些值
# Qmin: 頻率最?值
# Qmax: 頻率最?值
# d: 維度
# lower: 下界
# upper: 上界
def bat_algorithm(objfun, N_pop=20, N_gen=1000, A=0.5, r=0.5,
Qmin=0, Qmax=2, d=10, lower=-2, upper=2):
N_iter = 0 # Total number of function evaluations
#=====速度上下限================
Lower_bound = lower * np.ones((1,d))
Upper_bound = upper * np.ones((1,d))
Q = np.zeros((N_pop, 1)) # 頻率
v = np.zeros((N_pop, d)) # 速度
S = np.zeros((N_pop, d))
#=====初始化種群、初始解=======
# Sol = np.random.uniform(Lower_bound, Upper_bound, (N_pop, d)) # Fitness = objfun(Sol)
Sol = np.zeros((N_pop, d))
Fitness = np.zeros((N_pop, 1))
for i in range(N_pop):
Sol[i] = np.random.uniform(Lower_bound, Upper_bound, (1, d)) Fitness[i] = objfun(Sol[i])
#====找出初始最優解===========
fmin = min(Fitness)
Index = list(Fitness).index(fmin)
best = Sol[Index]
#======開始迭代=======
for t in range(N_gen):
#====對所有蝙蝠/解決?案進?循環 ======
for i in range(N_pop):
# Q[i] = Qmin + (Qmin - Qmax) * np.random.rand
Q[i] = np.random.uniform(Qmin, Qmax)
v[i] = v[i] + (Sol[i] - best) * Q[i]
S[i] = Sol[i] + v[i]
#===應?簡單的界限/限制====
Sol[i] = simplebounds(Sol[i], Lower_bound, Upper_bound)
# Pul rate
if rand() > r:
# The factor 0.001 limits the step sizes of random walks
S[i] = best + 0.001*np.random.randn(1, d)
#====評估新的解決?案 ===========
# print(i)
Fnew = objfun(S[i])
#====如果解決?案有所改進，或者聲?不太?，請更新====
if (Fnew <= Fitness[i]) and (rand() < A):
Sol[i] = S[i]
Fitness[i] = Fnew
Fitness[i] = Fnew
#====更新當前的最佳解決?案======
if Fnew <= fmin:
best = S[i]
fmin = Fnew
N_iter = N_iter + N_pop
print('Number of evaluations: ', N_iter)
print("Best = ", best, '\n fmin = ', fmin)
return best
def simplebounds(s, Lower_bound, Upper_bound):
Index = s > Lower_bound
s = Index * s + ~Index * Lower_bound
Index = s < Upper_bound
s = Index * s + ~Index * Upper_bound
return s
#====?標函數=============
def test_function(u):
a = u ** 2
return a.sum(axis=0)
if __name__ == '__main__':
# print(bat_algorithm(test_function))
bat_algorithm(test_function)
4.2 結果
5 Matlab實現
5.1 代碼
clear
wmax=0.9;%慣性權重最?值
wmin=0.4;%慣性權重最?值
n=10000; % 群體??
A=rand(1,n); % 聲?響度 (不變或者減?)
%% 頻率范圍
Qmin=0; % 最低頻率
Qmax=1; % 最?頻率
d=2;% 搜索變量的維數（即頻率和速度）
%% 初始矩陣
Q=zeros(n,1); % 頻率矩陣初始化
v=zeros(n,d); % 速度矩陣初始化，初始化意義就是產??個初始矩陣
%% x?變量范圍
u=-3;
o=12.1;
% y?變量范圍
p=4.1;
l=5.8;
%% 初始化群體/解
for i=1:n
Sol(i,1)=-3+(12.1+3)*rand(1,1);%x?變量范圍【-3，12.1】
Sol(i,2)=4.1+(5.8-4.1)*rand(1,1);%y?變量【4.1，5.8范圍】
%將隨機?成的兩個?變量帶?函數式
Fitness(i)=Fun(Sol(i,:));%函數值
end
%% 尋找當前最優解
[fmax,I]=max(Fitness);
best=Sol(I,:);
T=100;%飛?次數
%% 開始飛?
for t=1:T
for i=1:n,
Q(i)=Qmin+(Qmin-Qmax)*rand;%rand均勻分布的隨機數
%v(i,:)=v(i,:)+(Sol(i,:)-best)*Q(i);（原速度）
w=(wmax-wmin)*exp(-2*(t/T)^2)+wmin;%慣性權重因?
v(i,:)=w*v(i,:)+(Sol(i,:)-best)*A(i)*Q(i);%更改后的速度
S(i,:)=Sol(i,:)+v(i,:);%位置移動
%% 邊界問題，如果下次飛?超出?變量范圍外了，那么下次飛?的位置為投影在的邊界上的位置 %x軸
if S(i,1)>o
S(i,1)=o;
end
if S(i,1)<u
S(i,1)=u;
end
%y軸
if S(i,2)>l
S(i,2)=l;
end
if S(i,2)<p
S(i,2)=p;
end
%% 評估該次飛?后產?的新解
Fnew(i)=Fun(S(i,:));
end
[Fmax,Z]=max(Fnew);%找出該次飛?后產?的最?值
C(t,:)=S(Z,:);
FFnew(t)=Fmax;
end
[Ffmax,N]=max(FFnew);%找出整個飛?過程中的最?值
M=C(N,:)
Ffmax
