PSO(Particle Swarm Optimization,粒子族群最佳化)
developed in 1995 by Kennedy and Eberhart
inspired by social behavior of bird flocking
原理:一群particle分散在一個soluction space移動,每個particle會記錄自己的position和velocity
概念
1在pso中,每個particle代表一個solution
2判斷solution的好壞,是根據每個particle的fitness,該值依不同position會重新計算
3fitness的產生方式,依據解不同的問題,方法也會不同
ex:NOP(Numerical Optimization),MNO(Multiobjective Numerical Optimization),BIN(Binary Problem),clustering
4每個particle會在solution space中根據velocity不斷移動position,一直到停止條件成立
5velocity由3個值組成:1上次的velocity,2該particle本身的思考,3所有particle相互作用下的思考
pso資源
http://www.adaptivebox.net/research/bookmark/psocodes_link.html
http://www.particleswarm.info/
http://clerc.maurice.free.fr/pso/
資料結構
在D維空間下,
第i個粒子代表Xi,也等於(Xi1,Xi2,…XiD)
pbest:第i個粒子經歷過最好適應值的位置 Pi,也等於(Pi1,Pi2,…PiD)
gbest,所有第g群的粒子經歷過最好適應值的位置,Pg = (Pg1, Pg2, …, PgD),
第i個粒子的速度(移動的距離)代表Vi
在第d個維度的移動公式為
Vid_new=w*Vid+c1*rand()*(Pid-Xid)+c2*rand()*(Pgd-Xid)
ps:為了限制velocity在solution space,ex:若Vid_new > MAX_Vid,則Vid=MAX_Vid
Xid_new = Xid + Vid_new
ps:為了限制position在solution space,可檢查xid是否在該維度最大值與最小值之間
公式說明
w*Vid 粒子之前速度的慣性,值小容易local exploration有助於收斂,值大則相反
c1*rand()*(Pid-Xid)為粒子本身的思考
c2*rand()*(Pgd-Xid)為粒子相互影響下的思考
參數說明
Vid 粒子上次的速度,也就是粒子上次移動的距離
Vmax 粒子最大速度,太高會超過好解,太低則探索不足,常見的設定值是該維度最大及最小的距離
w:Inertia Weight (慣性權重)
c1,c2:acceleration constants
rand():0到1之間的隨機數
步驟大致如下
1random initial,隨機初始化Vid,Xid
2evaluate:計算每個粒子的fitneess(適應值)
3record pbest and gbest
對每個粒子,若該粒子現在位置的適應值優於pbest位置的適應值,則pbest位置=現在位置
對每個粒子,若該粒子現在位置的適應值優於gbest位置的適應值,則gbest位置=現在位置
4update Vid and Xid,根據移動公式調整粒子位置
5若達到結束條件(),否則返回步驟2
ps:結束條件設定通常為:設定固定世代數,或解一直沒有改善
演算法大致如下
//
initialize //step1
//
for t=1 to t_max do //step5
For each particle i do
assign zp to cluster
calculate fitness //step2
update pid(local best) and pgd(global best) //step3
//
For each particle i do
For each data vector zp
vid = w*vid+c1*rand()*(pid-xid)+c2*Rand()*(pgd-xid) //step4
if(vid>vmax)vid=vmax
xid = xid+vid
…………………………………………………….
運用在分群的做法
1
維度部份
在目標資料有D維空間下,若要分N群
粒子要維護每群的群中心,所以粒子Xi等於(Xi1,Xi2,…,XiN)
將每一群拆開後如下
(Xi11,Xi12,…,Xi1D
,Xi21,Xi22,…,Xi2D
,…
,XiN1,XiN2,…,XiND)
因此會對N*D個維度做計算
2
fitness部份
公式有以下兩種:
第一種:sum(j,sum_distance(Zp,Mj)/count(Cij) ) / N
第二種:sum(j,sum_distance(Zp,Mj))
j,表示第j群
count(Cij),第i個粒子的第j個群的資料筆數
Zp,第p個資料向量
Mj,第j個群中心
sum_distance(Zp,Mj),將組成Mj的所有Zp,拿來計算與Mj的距離並加總
…………………………………………………………………………………….
pso clustering for php example
/*********************************************
pso example
$data=psodata_fromfile(data_iris.txt);
$para['gene']=200;
$para['in']=10;
$para['jn']=3;
$para['w']=0.72;
$para['c1']=1.49;
$para['c2']=1.49;
$pso=pso($data['z'],$para);
print_r($pso['x']);
while($c=each($pso['c'])){
while($j=each($c[1])){
echo 'cluster'.$c[0].'=';
while($p=each($data['full'][$j[1]])){
echo $p[1].' : ';
}
echo "n";
}
}
echo '['.$pso['fitness'].'][time:'.$pso['time_total'].']';
***********************************************/
///////////////////////////////////////////////////////////////////////////////////////////////////////
///////////////////////////////////////////////////////////////////////////////////////////////////////
function pso($z,$para){
$time_start=caclutime();
/*the structure as follow
$para['gene']=200;
$para['in']=10;
$para['jn']=3;
$para['w']=0.72;
$para['c1']=1.49;
$para['c2']=1.49;
$p[i][1-D] 粒子i記錄的pbest
$pg[1-D] 記錄gbest
d=all atribute
p=all data
$z[1-D][1-p] 資料p的各維度數據
$z[attribute_1][data_1]
$z[attribute_1][data_2]
$z[attribute_1][data_p]
$z[attribute_2][data_1]
$z[attribute_2][data_2]
$z[attribute_2][data_p]
$z[attribute_3][data_1]
$z[attribute_3][data_2]
$z[attribute_3][data_p]
$x[i][1-N][1-D] 粒子i要維護的N個群中心
$x:record center of cluster
$x[patical_1][cluster_1][attribute_1]
$x[patical_1][cluster_1][attribute_d]
$x[patical_1][cluster_j][attribute_1]
$x[patical_1][cluster_j][attribute_d]
$c[i][j][1-count(Cij)] 在粒子i下,屬於第j群的所有資料
$c:record data of cluster
$c[patical_1][cluster_1][0]=sample1
$c[patical_1][cluster_1][1]=sample3
$c[patical_1][cluster_1][n]=sample6
$c[patical_1][cluster_2][0]=sample2
$c[patical_1][cluster_2][1]=sample4
$c[patical_1][cluster_2][n]=sample5
$c[patical_1][cluster_3][0]=sample7
*/
$pn=count($z[1]);
$dn=count($z);
for($d=1;$d<=$dn;$d++){
arsort($z[$d]);
$zdmax=each($z[$d]);
asort($z[$d]);
$zdmin=each($z[$d]);
$xida['max'][$d]=$zdmax[1];
$xida['min'][$d]=$zdmin[1];
$vmax[$d]=abs($zdmax[1]-$zdmin[1]);
}
for($i=1;$i<=$para['in'];$i++){
$fitness_pbest[$i]=9999999;
}
$fitness_gbest=9999999;
//init_clustercenter
reset($z);
$x=init_clustercenter($z,$para,$dn,$pn);
/////////
while($t<$para['gene']){
$t++;
//evaluate fitness
reset($x);
while($i=each($x)){
$state=asigncluster($i[1],$z,$dn,$pn);
$c[$i[0]]=$state['c_i'];
$f[$i[0]]=$state['d_total'];
$fitness=$f[$i[0]];
$fitness_gbest_history2_tmp[$i[0]]=$fitness;
$fitness_pbest_history2[$i[0]][$t]=$fitness;
if($fitness_pbest[$i[0]] > $fitness){
$pbest[$i[0]]=$x[$i[0]];
$fitness_pbest[$i[0]]=$fitness;
if($fitness_gbest > $fitness_pbest[$i[0]]){
$fitness_gbest=$fitness;
$particle_gbest=$i[0];
$gbest=$x[$i[0]];
$bestc=$c[$i[0]];
}
}
$fitness_pbest_history[$i[0]][$t]=$fitness_pbest[$i[0]];
}
asort($fitness_gbest_history2_tmp);
$tmp_localgbest=each($fitness_gbest_history2_tmp);
$fitness_gbest_history2[$t]=$tmp_localgbest[1];
$fitness_gbest_history[$t]=$fitness_gbest;
for($j=1;$j<=$para['jn'];$j++){
$bestc_count[$j][$t]=count($bestc[$j]);
}
//update
reset($x);
while($i=each($x)){
while($j=each($i[1])){
while($d=each($j[1])){
$pid=$pbest[$i[0]][$j[0]][$d[0]];
$pgd=$gbest[$j[0]][$d[0]];
$xid=$d[1];
$c1rand=rand(1,1000)/1000;
$c2rand=rand(1,1000)/1000;
$vid_new[$i[0]][$d[0]]=$para['w']*$vid[$i[0]][$d[0]]+$para['c1']*$c1rand*($pid-$xid)+$para['c2']*$c2rand*($pgd-$xid);
if($vid_new[$i[0]][$d[0]]>$vmax[$d[0]])$vid_new[$i[0]][$d[0]]=$vmax[$d[0]];
if($vid_new[$i[0]][$d[0]]<(0-$vmax[$d[0]]))$vid_new[$i[0]][$d[0]]=(0-$vmax[$d[0]]);
$vid_abs[$i[0]][$d[0]]=abs($vid_new[$i[0]][$d[0]]);
$new_xid=$vid_new[$i[0]][$d[0]]+$xid;
if($new_xid>$xida['max'][$d[0]])$new_xid=$xida['max'][$d[0]];
if($new_xid<$xida['min'][$d[0]])$new_xid=$xida['min'][$d[0]];
$x[$i[0]][$j[0]][$d[0]]=$new_xid;
}
}
$vidsum_history[$i[0]][$t]=array_sum($vid_abs[$i[0]]);
}
$vid=$vid_new;
}
////////////////////end
$pso['c']=$c[$particle_gbest];
$pso['fitness']=$fitness_gbest;
$pso['fitness_history']=$fitness_gbest_history;
$pso['fitness_history2']=$fitness_gbest_history2;
$pso['fitness_pbest_history']=$fitness_pbest_history;
$pso['fitness_pbest_history2']=$fitness_pbest_history2;
$pso['x']=$gbest;
$pso['time_total']=caclutime()-$time_start;
$pso['vidsum_history']=$vidsum_history;
$pso['bestc_count']=$bestc_count;
return $pso;
}
///////////////////////////////////////////////////////////////////////////////////////////////////////
function asigncluster($x_i,$z,$dn,$pn){
$d_total=0;
for($p=1;$p<=$pn;$p++){
reset($x_i);
while($j=each($x_i)){
$distance_j=0;
$distance_pj=0;
for($d=1;$d<=$dn;$d++){
$data_p_dime=$z[$d][$p];
$clustercenter_j_dime=$j[1][$d];
$distance_pj+=pow($data_p_dime-$clustercenter_j_dime,2);
}
$distancepj[$j[0]]=sqrt($distance_pj);
}
asort($distancepj);
$min=each($distancepj);
$c_i[$min[0]][]=$p;
$d_total+=$min[1];
$check[$min[0]]=1;
}
for($j=1;$j<=$para['jn'];$j++){
if($check[$j]!=1)$c_i[$j]='';
}
$state['c_i']=$c_i;
$state['d_total']=$d_total;
return $state;
}
//////////////////////////////////////////////////////////////////////////////////////////
function init_clustercenter($z,$para,$dn,$pn){
/*
$x[1][1][1]=5.006;
$x[1][1][2]=3.418;
$x[1][1][3]=1.464;
$x[1][1][4]=0.244;
$x[1][2][1]=5.883606557377;
$x[1][2][2]=2.7409836065574;
$x[1][2][3]=4.3885245901639;
$x[1][2][4]=1.4344262295082;
$x[1][3][1]=6.8538461538462;
$x[1][3][2]=3.0769230769231;
$x[1][3][3]=5.7153846153846;
$x[1][3][4]=2.0538461538462;
*/
$len=ceil($pn/$para['jn']);
for($i=2;$i<=$para['in'];$i++){
for($j=1;$j<=$para['jn'];$j++){
$len_e=$len*$j;
$len_s=$len*($j-1)+1;
if($len_e>=$pn)$len_e=$pn;
$center=rand($len_s,$len_e);
for($d=1;$d<=$dn;$d++){
$x[$i][$j][$d]=$z[$d][$center];
}
}
}
return $x;
}
//////////////////////////////////////////////////////////////////////////////
function caclutime(){
$time = explode( " ", microtime());
$usec = (double)$time[0];
$sec = (double)$time[1];
return $sec + $usec;
}
/////////////////////////////////////////////////////////////////////////////
function psodata_fromfile($sourcedata){
/*
data fromat description
name1,attribute1_value,attribute2_value,...attributeN_value
name2,attribute1_value,attribute2_value,...attributeN_value
...
nameN,attribute1_value,attribute2_value,...attributeN_value
ex:
Iris-setosa,5.1,3.5,1.4,0.2
Iris-setosa,4.9,3,1.4,0.2
Iris-setosa,4.7,3.2,1.3,0.2
Iris-setosa,4.6,3.1,1.5,0.2
Iris-setosa,5,3.6,1.4,0.2
ps:attribute must is continue,not discrete
*/
$data=file($sourcedata);
$dime=each($data);
$dime_a=explode(',',$dime[1]);
$dn=count($dime_a)-1;
reset($data);
while(list($key,$dot_id)=each($data)){
$i++;
$dot_id = trim($dot_id, " n.");
$dot_id_a=explode(',',$dot_id);
$fullz[$i][0]=$dot_id_a[0];
for($d=1;$d<=$dn;$d++){
$z[$d][$i]=$dot_id_a[$d];
$fullz[$i][$d]=$dot_id_a[$d];
}
}
$data['z']=$z;
$data['full']=$fullz;
return $data;
}