- 1 从布谷鸟的育雏到布谷鸟算法
- 2 布谷鸟算法
- 3 萊维飞行与公式(1)的深层含义
- 4 附:CS算法求解函数最小值代码
- 5 源码下载
- 6 参考文献
1 从布谷鸟的育雏到布谷鸟算法
布谷鸟不会做窝,也不会育雏,在春末夏初,向北飞,趁别的鸟(宿主鸟)外出觅食时,将卵蛋产在宿主鸟窝里,让宿主鸟抚养自己孩子 。当然,布谷鸟在产卵前,为了不被宿主鸟发现鸟窝的异常,会把宿主的卵移走。而一旦靠养母孵化的雏鸟,也有将宿主鸟本身的雏鸟推出巢穴的本性,并且会模仿其他鸟的行为来增大不被宿主鸟发现的概率1
2009年,Xin-She Yang2 与Suash Deb在《Cuckoo Search via Levy Flights》一文中提出了布谷鸟算法(简称CS)。假设每只布谷鸟一次只产一枚卵 ,并且宿主鸟发现外来鸟蛋后,就舍弃该鸟窝,另寻他地建造新的鸟窝 ,那么可以认为 :鸟窝=卵蛋=解,卵蛋是否能够成功被宿主鸟孵化并茁长成长是衡量解好坏的唯一标准 。布谷鸟寻找鸟窝下蛋的过程就是在D维空间中寻找解的过程 ,而鸟窝的好坏象征着解的好坏。
2 布谷鸟算法
布谷鸟算法是布谷鸟育雏行为和萊维飞行结合的一种算法 。
在CS算法中,有两个路径(或者说成是两个位置的更新)备受关注:
一个是布谷鸟寻找鸟窝下蛋的寻找路径是采用早已就有的萊维飞行3,如上图所示,无敌的走位是一种长步长与短步长相间的走位,这其实就是萊维飞行的主要特点,学者们也证实了自然界中很多鸟类的飞行也遵从萊维飞行,这也是最有效寻找目标的方法之一 。所以采用萊维飞行更新鸟窝位置的公式被定义如下:
Xt+1=Xt+α⨂Levy(β)
X
t
+
1
=
X
t
+
α
⨂
L
e
v
y
(
β
)
, 公式(1)
其中 ,
α
α
是步长缩放因子,
Levy(β)
L
e
v
y
(
β
)
是萊维随机路径,
⨂
⨂
就是
.∗
.
∗
运算
另一个是宿主鸟以一定概率Pa发现外来鸟后重新建窝的位置路径,这个路径可以用萊维飞行或者随机方式4,(本文采用随机) , 除此之外,这个位置普遍采用偏好随机游动的方式,即利用了其他鸟窝的相似性5。所以新建的鸟窝的位置的公式被定义如下:
Xt+1=Xt+r⨂Heaviside(Pa−ϵ)⨂(Xi−Xj)
X
t
+
1
=
X
t
+
r
⨂
H
e
a
v
i
s
i
d
e
(
P
a
−
ϵ
)
⨂
(
X
i
−
X
j
)
, 公式(2)
其中,
r,ϵ
r
,
ϵ
是服从均匀分布的随机数,
Heaviside(x)
H
e
a
v
i
s
i
d
e
(
x
)
是跳跃函数(x>0,=1;x<0,=0) ,
Xi,Xj
X
i
,
X
j
是其他任意的连个鸟窝。
CS算法的执行过程如下:
![这里写图片描述](https://img-blog.csdn.net/20180703215829945?watermark/2/text/aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2c0MjU2ODA5OTI=/font/5a6L5L2T/fontsize/400/fill/I0JBQkFCMA==/dissolve/70)
3 萊维飞行与公式(1)的深层含义
从数学的发展史上说,早在1937年, P. Levy6确定了对称Levy稳定分布的积分形式为
Levy(s)=1π∫+∞0exp(−β|k|λ)cos(ks)dk
L
e
v
y
(
s
)
=
1
π
∫
0
+
∞
e
x
p
(
−
β
|
k
|
λ
)
c
o
s
(
k
s
)
d
k
,但是该积分并没有明确的解析,要生成一个服从该分布的随机数是难上加难的问题,不过当
s≫s0>0,即s→∞
s
≫
s
0
>
0
,
即
s
→
∞
时,
Levy(s)≈λβΓ(λ)sin(πλ2)π.1s1+λ
L
e
v
y
(
s
)
≈
λ
β
Γ
(
λ
)
s
i
n
(
π
λ
2
)
π
.
1
s
1
+
λ
,通常
β=1
β
=
1
。这个近似的分布呈现幂律行为(重尾或长尾巴),这个行为类似于二八原则[^6],或者说少部分人集中了世界大部分的财富,正如下图所示的,这个分布总是有一个长尾巴或者称之为重尾巴,有时也叫做一个翼。
![这里写图片描述](https://img-blog.csdn.net/20180703215842697?watermark/2/text/aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2c0MjU2ODA5OTI=/font/5a6L5L2T/fontsize/400/fill/I0JBQkFCMA==/dissolve/70)
萊维飞行的方差随时间呈现指数的关系,即
σ2(t)~t3−β,1≤β≤3
σ
2
(
t
)
~
t
3
−
β
,
1
≤
β
≤
3
,所以萊维飞行比布朗运动更加的出色。
此后,不少学者根据这个近似部分提出很多用于生成服从萊维分布的随机数的实现方法,其中就包含了Mantegna7在1994年提出的一种用正太分布求解随机数的方法,有时也叫Mantegna方法,生成服从萊维分布的随机步长的方法如下:
s=u|v|1β
s
=
u
|
v
|
1
β
其中,
u~N(0,σ2),v~N(0,1)
u
~
N
(
0
,
σ
2
)
,
v
~
N
(
0
,
1
)
,
σ={Γ(1+β)sin(πβ2)βΓ(1+β2)2β−12}1β
σ
=
{
Γ
(
1
+
β
)
s
i
n
(
π
β
2
)
β
Γ
(
1
+
β
2
)
2
β
−
1
2
}
1
β
在matlab中用Mantegna方法模拟二维平面萊维飞行:
x = [0,0];
y = [0,0];
beta = 1.5;
sigma_u = (gamma(1+beta)*sin(pi*beta/2)/(gamma((1+beta)/2)*beta*2^((beta-1)/2)))^(1/beta);
sigma_v = 1;
for i=1:1000
u = normrnd(0,sigma_u);
v = normrnd(0,sigma_v);
s = u/(abs(v))^(1/beta);
x(:,1) = x(:,2);
x(:,2) = x(:,1)+1*s;
u = normrnd(0,sigma_u);
v = normrnd(0,sigma_v);
s = u/(abs(v))^(1/beta);
y(:,1) = y(:,2);
y(:,2) = y(:,1)+1*s;
plot(x,y);
hold on;
end
axis square;
![这里写图片描述](https://img-blog.csdn.net/20180703215858214?watermark/2/text/aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2c0MjU2ODA5OTI=/font/5a6L5L2T/fontsize/400/fill/I0JBQkFCMA==/dissolve/70)
从模拟上来看,图形的路径确实符合萊维飞行的长短相间的特征,Mantegna用正太分布实现了生成服从萊维分布随机步长的方法是可靠的 。
时间到了2009年,Xin-She Yang 与Suash Deb提出了布谷鸟算法,同时,Yang把Levy分布函数经过简化和傅立叶变换后得到其幂次形式的概率密度函数8 :
Levy~u=t−β,1≤β≤3
L
e
v
y
~
u
=
t
−
β
,
1
≤
β
≤
3
。并把萊维飞行用在了鸟窝位置的更新上,于是产生了公式(1)
Xt+1=Xt+α⨂Levy(β)
X
t
+
1
=
X
t
+
α
⨂
L
e
v
y
(
β
)
。这个计算式其实就是
Xt+1=Xt+αS
X
t
+
1
=
X
t
+
α
S
,
S
S
就是服从Levy分布Levy~u=t−β,1≤β≤3的随机步长,考虑到具体怎么计算时,Yang采用的正是1994年的Mantegna方法 。
所以在布谷鸟算法中,我们可以用下面的具体计算公式来计算鸟窝的更新位置:
Xt+1=Xt+αS=Xt+α.∗服从N(0,σ2)的随机数|服从N(0,1)的随机数|1β
X
t
+
1
=
X
t
+
α
S
=
X
t
+
α
.
∗
服
从
N
(
0
,
σ
2
)
的
随
机
数
|
服
从
N
(
0
,
1
)
的
随
机
数
|
1
β
其中,
σ={Γ(1+β)sin(πβ2)βΓ(1+β2)2β−12}1β
σ
=
{
Γ
(
1
+
β
)
s
i
n
(
π
β
2
)
β
Γ
(
1
+
β
2
)
2
β
−
1
2
}
1
β
,通常,
β=1.5
β
=
1.5
这在matlab等一些编程工具中都是可以计算的。
值得一提的是,
α
α
是步长缩放因子,通常
α=1
α
=
1
,在之后的布谷鸟算法发展中,针对
α
α
有各种各样的变种,如Yang[^8]为了让算法适应不同的解,让
α=α0(Xi−Xj)
α
=
α
0
(
X
i
−
X
j
)
,
Xi,Xj
X
i
,
X
j
为任意不同的解 。
4 附:CS算法求解函数最小值代码
求函数
f(x)=∑ni=1x2i,(−20≤x≤20,n=10)
f
(
x
)
=
∑
i
=
1
n
x
i
2
,
(
−
20
≤
x
≤
20
,
n
=
10
)
最小值
clear all ;
close all ;
clc ;
N = 25;
D = 10 ;
T = 200 ;
Xmax = 20 ;
Xmin = -20 ;
Pa = 0.25 ;
nestPop = rand(N,D)*(Xmax-Xmin)+Xmin ;
for t=1:T
levy_nestPop = func_levy(nestPop,Xmax,Xmin) ;
nestPop = func_bestNestPop(nestPop,levy_nestPop);
rand_nestPop = func_newBuildNest(nestPop,Pa,Xmax,Xmin);
nestPop = func_bestNestPop(nestPop,rand_nestPop) ;
[~,index] = max(func_fitness(nestPop)) ;
trace(t) = func_objValue(nestPop(index,:)) ;
end
figure
plot(trace);
xlabel('迭代次数') ;
ylabel('适应度值') ;
title('适应度进化曲线') ;
function [ result ] = func_levy( nestPop,Xmax,Xmin)
[N,D] = size(nestPop) ;
beta = 1.5 ;
alpha = 1 ;
sigma_u = (gamma(1+beta)*sin(pi*beta/2)/(beta*gamma((1+beta)/2)*2^((beta-1)/2)))^(1/beta) ;
sigma_v = 1 ;
u = normrnd(0,sigma_u,N,D) ;
v = normrnd(0,sigma_v,N,D) ;
step = u./(abs(v).^(1/beta)) ;
nestPop = nestPop+alpha.*step ;
nestPop(find(nestPop>Xmax)) = Xmax ;
nestPop(find(nestPop<Xmin)) = Xmin ;
result = nestPop ;
end
function [ nestPop ] = func_bestNestPop( nestPop,new_nestPop )
index = find(func_fitness(nestPop)<func_fitness(new_nestPop)) ;
nestPop(index,:) = new_nestPop(index,:) ;
end
function [ nestPop ] = func_newBuildNest( nestPop ,Pa ,Xmax,Xmin)
[N,D] = size(nestPop) ;
nestPop = nestPop+rand.*heaviside(rand(N,D)-Pa).*(nestPop(randperm(N),:)-nestPop(randperm(N),:));
nestPop(find(nestPop>Xmax)) = Xmax ;
nestPop(find(nestPop<Xmin)) = Xmin ;
end
![这里写图片描述](https://img-blog.csdn.net/20180703215912847?watermark/2/text/aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2c0MjU2ODA5OTI=/font/5a6L5L2T/fontsize/400/fill/I0JBQkFCMA==/dissolve/70)
5 源码下载
https://download.csdn.net/download/g425680992/10517545
6 参考文献
- 布谷鸟搜索算法研究综述 , 兰少峰 ↩
- Cuckoo search via Lévy flights. Yang XS ↩
- Cuckoo search via Lévy flights. Yang XS ↩
- 逐维改进的布谷鸟搜索算法 , 王李进 ↩
- Cuckoo search for inverse problems and simulated-driven shape optimization ,Yang XS ↩
- P. Levy, Theoric de l’Addition des Variables Aleatoires ↩
- Nature-Inspired Metaheuristic Algorithms Second Edition,Yang XS ↩
- Nature-Inspired Metaheuristic Algorithms Second Edition,Yang XS ↩
本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)