网站编程 外包类型,爱奇艺会员做任务送十天网站,佛山网站建设公司排行,四川省商投建设公司官网目录1.算法1.1.原理1.2.性能比较1.3.步骤2.代码2.1.源码及注释2.2.执行与效果1.算法
1.1.原理 \qquad建议没接触过粒子群算法的朋友先看较为基础的全局粒子群算法原理及介绍#xff0c;以下博文链接有详细的讲解、代码及其应用举例#xff1a; 【Simulink】粒子群算法#…
目录1.算法1.1.原理1.2.性能比较1.3.步骤2.代码2.1.源码及注释2.2.执行与效果1.算法
1.1.原理
\qquad建议没接触过粒子群算法的朋友先看较为基础的全局粒子群算法原理及介绍以下博文链接有详细的讲解、代码及其应用举例 【Simulink】粒子群算法PSO整定PID参数附代码和讲解 \qquad这里就介绍一下全局粒子群算法和混合粒子群算法的区别。全局粒子群算法General PSO将粒子速度矢量影响因子分为了三类惯性因子、自身最优因子、社会因子。每个粒子的初速度定为0即v00v_00v00第jjj个粒子1≤j≤m1≤j≤m1≤j≤m的下一次迭代的速度v(j)v^{(j)}v(j)由以下三部分组成 v(j)w⋅v0c1⋅rand⋅(P(j)−X(j))c2⋅rand⋅(PG−X(j))v^{(j)}w\cdot v_0c_1\cdot rand\cdot (P^{(j)}-X^{(j)})c_2\cdot rand\cdot(P_G-X^{(j)}) v(j)w⋅v0c1⋅rand⋅(P(j)−X(j))c2⋅rand⋅(PG−X(j)) 注rand是(0,1)的随机数v0v_0v0代表上一次粒子的速度。 第一部分为自身惯性因子因为下一次的迭代次数保留了上一次的速度信息 第二个部分为自身最优因子P(j)\color{blue}{P^{(j)}}P(j)为第jjj个因子在之前所有迭代次数中自适应度最高的位置可以理解为历史上自身离最优解最近的位置。 第三部分为社会因子PG\color{blue}{P_G}PG为种群在之前所有迭代次数中自适应度最高的位置可以理解为历史上所有粒子离最优解最近的位置中的最近的一个。 \qquad而混合粒子群算法Hybrid PSO与全局粒子群算法唯一的区别就是将社会因子分解为了两部分——局部社会因子和全局社会因子。即v(j)w⋅v0c1⋅rand⋅(P(j)−X(j))c2[q⋅rand⋅(PG−X(j))(1−q)⋅rand⋅(PL(j)−X(j))]v^{(j)}w\cdot v_0c_1\cdot rand\cdot (P^{(j)}-X^{(j)})c_2[q\cdot rand\cdot(P_G-X^{(j)})(1-q)\cdot rand\cdot(P_L^{(j)}-X^{(j)})]v(j)w⋅v0c1⋅rand⋅(P(j)−X(j))c2[q⋅rand⋅(PG−X(j))(1−q)⋅rand⋅(PL(j)−X(j))] 一般0q1,即两部分社会因子都为正。与上面不同的是PL(j)P_L^{(j)}PL(j)是第jjj个粒子的局部最优解的坐标。 \qquad全局社会因子使得粒子的速度矢量有向目前的全局最优解趋向的趋势而局部社会因子则使速度矢量包含向粒子的目前的局部最优解趋向的趋势。全局最优解好理解局部最优解就是该粒子一定范围内粒子中的最优解一般这个范围是固定的而包含的粒子数量如果为1则是粒子本身不确定。 \qquad全局粒子群算法收敛速度快但是对于多极值问题容易收敛到一个局部最优解混合粒子群算法由于要计算局部最优解的原因收敛速度慢但可以很好地解决全局粒子群算法收敛到局部最优解的问题。
1.2.性能比较
\qquad对粒子群算法感兴趣的朋友可以读一下面这两部分解释。
1.为什么全局粒子群算法容易收敛到局部最优解呢 \qquad我们不妨假设仅有3个粒子目前有2个极小值点A和B假设第一次迭代位置随机初速度为0一个粒子a恰好就位于其中一个极小值点A附近另2个粒子b和c彼此较近但离B较远函数值f(A)f(C)f(B)f(A)f(C)f(B)f(A)f(C)f(B)。初速度为0惯性因子为0。又因为a就是当前全局最优点自身最优因子也是本身因为迭代次数为1所以后两项都是0a不会动。粒子b和c的历史最优因子就是它现在的位置而社会因子使它向粒子a移动并且每移动一次它的自身最优因子都比上次的好即不会对它的目前的速度造成任何影响。
\qquad最终两个粒子都会收敛到极值点A。即使极值点B的目标函数值更小一些但由于极值点A先被粒子群“找到”因此就像一个极小值点黑洞一样把所有的粒子都“吸”了过去。所以简单地来说全局最优解就像一个吸引其他粒子的“黑洞”如果在找到第二个“黑洞”之前所有粒子都被吸引到了第一个“黑洞”的周围那就只能找到第一个局部最优解了。
123
2.那为什么混合粒子群算法就可以避免这个问题呢 \qquad我们不妨假设一个和上文一样的场景。目标函数只有2个极小值点A和B第一次迭代粒子a就位于A的附近粒子b和c不位于B的附近但彼此离得较近离A较远。而且目标函数值都比a要大。初速度为0因此三者的惯性因子为0。a的自身最优因子就是它本身同时也是全局最优解因此a不会动。而b、c的自身最优因子也是本身然而全局社会最优因子是a即二者都有向a偏移的趋势。
\qquad现在考虑局部社会最优因子b和c距离较小会受到局部社会因子的影响而二者都离a较远不会受到a的影响同样a也不会受到b和c的影响。假设c的函数值比b小因此b是c的局部最优解同时也是它自身的。在局部社会因子的作用下b有向c移动的趋势c和a均不受影响c和a的局部最优解都是本身。假设局部最优因子的权重合适b应该向极小值B移动c向a点移动。因为b有惯性因子很可能在第二次迭代的过程中离极小值B的距离会比c要短满足f(B)f(C)f(A)f(B)f(C)f(A)f(B)f(C)f(A)。那么在接下来的几步迭代中三者就会顺利收敛到极小值B附近而不是更大的一个极小值A进而找到全局最优解。
1234
1.3.步骤
\qquad步骤和全局粒子群算法基本一致下面将不一样的步骤标红看过全局粒子群算法那篇博文的读者可以重点关注一下。 【step0】确定参数维度N\color{blue}{N}N、参数范围即每只鸟的初始位置确定惯性系数c1,c2,wc_1,c_2,wc1,c2,w确定种群规模m迭代次数n以及局部因子作用半径RRR。
每个粒子的初始位置是随机的设输入向量x(x1,x2,...,xN)Tx(x_1,x_2,...,x_N)^Tx(x1,x2,...,xN)T必须满足参数范围限制为 {xmin(1)x1xmax(1)xmin(2)x2xmax(2)...xmin(N)x1xmax(N)\begin{cases}x_{min}^{(1)}x_1x_{max}^{(1)} \\x_{min}^{(2)}x_2x_{max}^{(2)} \\... \\x_{min}^{(N)}x_1x_{max}^{(N)} \end{cases}⎩⎪⎪⎪⎨⎪⎪⎪⎧xmin(1)x1xmax(1)xmin(2)x2xmax(2)...xmin(N)x1xmax(N) 记Xmin(xmin(1),xmin(2),...xmin(N)),Xmax(xmax(1),xmax(2),...xmax(N))\color{blue}X_{min}(x_{min}^{(1)},x_{min}^{(2)},...x_{min}^{(N)}),X_{max}(x_{max}^{(1)},x_{max}^{(2)},...x_{max}^{(N)})Xmin(xmin(1),xmin(2),...xmin(N)),Xmax(xmax(1),xmax(2),...xmax(N)) 则输入向量xxx应满足XminxXmaxX_{min}xX_{max}XminxXmax 【step1】计算每个粒子的速度 每个粒子的初速度定为0即v00v_00v00第jjj个粒子1≤j≤m1≤j≤m1≤j≤m的下一次迭代的速度v(j)v^{(j)}v(j)由三部分组成 v(j)w⋅v0c1⋅rand⋅(P(j)−X(j))c2[q⋅rand⋅(PG−X(j))(1−q)⋅rand⋅(PL(j)−X(j))]v^{(j)}w\cdot v_0c_1\cdot rand\cdot (P^{(j)}-X^{(j)})c_2[q\cdot rand\cdot(P_G-X^{(j)})(1-q)\cdot rand\cdot(P_L^{(j)}-X^{(j)})]v(j)w⋅v0c1⋅rand⋅(P(j)−X(j))c2[q⋅rand⋅(PG−X(j))(1−q)⋅rand⋅(PL(j)−X(j))]
注rand是(0,1)的随机数v0v_0v0代表上一次粒子的速度。 第一部分为自身惯性因子因为下一次的迭代次数保留了上一次的速度信息 第二个部分为自身最优因子P(j)\color{blue}{P^{(j)}}P(j)为第jjj个因子在之前所有迭代次数中自适应度最高的位置可以理解为历史上自身离食物最近的位置。 第三部分为社会因子由全局社会因子和局部社会因子组成。
PG\color{blue}{P_G}PG为全局最优解是种群在之前所有迭代次数中自适应度最高的位置可以理解为历史上所有粒子离食物最近的位置中的最近的一个。 一般情况下取w1,c1c22w1,c_1c_22w1,c1c22当种群规模较大时可以让www的值随迭代次数减小以增加收敛速度。 PL(j)\color{blue}{P_L^{(j)}}PL(j)为局部最优解每个粒子都有各自的局部最优解可以理解为以该粒子为中心R\color{blue}{R}R为半径的超球体内包含的所有粒子中自适应度最高的粒子的位置。q是全局社会因子的占比比例越高局部社会因子所占的权重越小。R\color{blue}{R}R称为局部因子作用半径和粒子分布平均密度有关。
【step2】按照step1的速度公式计算下一次的速度并计算下一次的粒子位置。对于第jjj个粒子第k1k1k1次迭代第k1k1k1代的位置Xk1(j)\color{blue}{X_{k1}^{(j)}}Xk1(j)与第kkk次迭代的位置XK(j)\color{blue}{X_K^{(j)}}XK(j)与速度vk(k1)\color{blue}{v_k^{(k1)}}vk(k1)关系为 Xk1(j)Xk(j)vk(j1)⋅dtX^{(j)}_{k1}X^{(j)}_{k}v_k^{(j1)}\cdot dtXk1(j)Xk(j)vk(j1)⋅dt 其中dt\color{blue}{dt}dt是仿真间隔一般情况下可以取1如果希望仿真得慢一点搜索平滑一点可以适当减小dt\color{blue}{dt}dt。 【step3】计算每个粒子的自适应度Fk1(j)\color{blue}{F^{(j)}_{k1}}Fk1(j)为了计算出step1公式中的P(j)、PG\color{blue}{P^{(j)}、P_G}P(j)、PG和PL(j)\color{red}{P_L^{(j)}}PL(j)。为了方便起见我们记前k次计算得到了的PGP_GPG为PG(k)P_G^{(k)}PG(k)则第k1次迭代中我们将适应度最高的粒子位置记为PG(k1)P_G^{(k1)}PG(k1)则最终的PGP_GPG为 PG{PG(k)F(PG(k))F(PG(k1))PG(k1)F(PG(k))F(PG(k1))P_G\begin{cases}P_G^{(k)}\qquad F(P_G^{(k)})F(P_G^{(k1)}) \\[2ex]P_G^{(k1)}\quad F(P_G^{(k)})F(P_G^{(k1)}) \end{cases}PG⎩⎨⎧PG(k)F(PG(k))F(PG(k1))PG(k1)F(PG(k))F(PG(k1)) 同样我们记前k次计算得到的第jjj个粒子的位置为Pk(j)P^{(j)}_{k}Pk(j)第k1次计算得到的第jjj个粒子的位置为Pk1(j)P^{(j)}_{k1}Pk1(j)则最终的P(j)P^{(j)}P(j)为 P(j){Pk(j)F(Pk(j))F(Pk1(j))Pk1(j)F(Pk(j))F(Pk1(j))P^{(j)}\begin{cases}P_{k}^{(j)}\quad F(P_{k}^{(j)})F(P_{k1}^{(j)}) \\[2ex]P_{k1}^{(j)}\quad F(P_{k}^{(j)})F(P_{k1}^{(j)}) \end{cases}P(j)⎩⎨⎧Pk(j)F(Pk(j))F(Pk1(j))Pk1(j)F(Pk(j))F(Pk1(j)) 我们计第k次迭代的PL(j)P_L^{(j)}PL(j)为PLk(j)P_{L_k}^{(j)}PLk(j)第kkk次迭代的第iii粒子的位置为Pk(i)P_k^{(i)}Pk(i)则 PLk1(j)avg[min(F(Pk(i)))],∣∣Pk(i)−Pk(j)∣∣≤RP_{L_{k1}}^{(j)}avg[min(F(P_k^{(i)}))],||P_k^{(i)}-P_k^{(j)}||≤RPLk1(j)avg[min(F(Pk(i)))],∣∣Pk(i)−Pk(j)∣∣≤R 即所有距离第jjj个粒子不超过R的粒子中函数值最小的那个。 \qquad注意PL(j)P_L^{(j)}PL(j)不需要保存为数组也不需要和上次的值比较它的值依赖于每次迭代粒子的位置也会跟着来粒子位置分布的变化为变化。
【step4】更新每个粒子的信息。 由于我们在step2的位置迭代中已经更新过粒子的位置信息在step1中的速度迭代中已经更新过粒子的速度信息而在step3中又更新了每个粒子的历史最优位置P(j)P^{(j)}P(j)及种群最优位置PGP_GPG因此实际上如果仅需要知道最优解的情况下我们不需要做这一步。 但如果需要作出参数随迭代次数的改变的图的话每次迭代产生最优解PG(k)P_G^{(k)}PG(k)及最优适应度F(PG(k))F(P_G^{(k)})F(PG(k))需要用数组保存下来。
2.代码
2.1.源码及注释
\qquad混合粒子群算法与全局粒子群算法的唯一区别就在一多出一个局部社会因子计算局部社会因子就需要计算每个粒子为中心R为半径的“局部”内自适应度最高的粒子。因此判断粒子与粒子间的距离成为了一大难题。将这部分代码写成并行形式可以大大降低算法时间复杂度和空间复杂度。 \qquad我们定义一个邻接矩阵 D[d11d12...d1nd21d22...d2n⋮⋮⋱⋮dn1dn2...dnn]D\begin{bmatrix} d_{11} d_{12} ... d_{1n}\\ d_{21} d_{22} ... d_{2n}\\ \vdots \vdots \ddots \vdots\\ d_{n1} d_{n2} ... d_{nn} \end{bmatrix}D⎣⎢⎢⎢⎡d11d21⋮dn1d12d22⋮dn2......⋱...d1nd2n⋮dnn⎦⎥⎥⎥⎤ 其中dijd_{ij}dij代表第iii个粒子距离第jjj个粒子的距离dii0d_{ii}0dii0。 每次迭代过程中都计算邻接矩阵DDD计算第jjj粒子的局部社会因子时即参考DDD矩阵的第jjj行或列的数据即可。 计算邻接矩阵DDD的函数书写如下(Distance_Matrix.m)
function DMDistance_Matrix(X)
msize(X,2);%样本数
nsize(X,1);%维数
DMzeros(n,n);%距离矩阵
for t1:mDifones(n,1)*X(t,:)-X;%坐标差矩阵Dissum(Dif.^2,2);%第t个样本距离其他所有样本的距离矩阵DM(t,:)Dis;
end注意这里X是样本矩阵样本是以行向量的形式组成矩阵的。
主程序代码如下所示(Hybrid_POS.m)
function [Pg,fmin]Hybrid_PSO(f,dimension,n,m,xmax,xmin,vmax,vmin,w,c1,c2,R,p,eta,ifdraw)
%混合粒子群算法,f为目标函数,dimension为维度,n为代数,m为种群规模
%w,c1,c2分别为惯性权重因子历史因素因子社会因素因子
%位置限幅为[xmin,xmax]速度限幅为[vmin,vmax]
%R为局部粒子群算法作用范围
%p为局部粒子群占的比重Max1
%eta为w的衰减因子范围[0,1]0代表不衰减1代表epoch全跑完w0
%ifdraw: True时为绘制粒子位置变化图False时不绘制
% 以下四行为位置限幅和速度限幅对应的四个参数均为行向量
Xmaxones(m,1)*xmax;
Xminones(m,1)*xmin;
Vmaxones(m,1)*vmax;
Vminones(m,1)*vmin;
Savefzeros(n1,1);%记录自适应度F的数组
if ifdrawSaveDatazeros(m,dimension,10);%记录11代种群的位置只需结果时不需要此数组
end
dt0.3;%位移仿真间隔
vzeros(m,dimension);%初始速度为0
X(Xmax-Xmin).*rand(m,dimension)Xmin;%初始位置满足(-xmax,xmax)内均匀分布
PX;%P为每个粒子每代的最优位置
Local_Pgzeros(m,dimension);%局部最优解矩阵
last_ff(X);
[fmin,min_i]min(last_f);%Pg为所有代中的最优位置
PgX(min_i,:);%获取这个最优位置的坐标
Savef(1)fmin;N0;
for iter1:nDMDistance_Matrix(X);%距离矩阵for i1:mXR_ifind(DM(i,:)R);[~,min_i]min(last_f(XR_i));%按行筛选小于R的点Local_Pg(i,:)X(XR_i(min_i),:);%局部最优解的坐标endrandM rand(m,dimension);vw*vc1*rand*(P-X)c2*randM.*((1-p)*(ones(m,1)*Pg-X)p*(Local_Pg-X));v(vVmax).*Vmax(vVmin vVmax).*v(vVmin).*Vmin;XXv*dt;X(XXmax).*Xmax(XXmin XXmax).*X(XXmin).*Xmin;new_ff(X);%新的目标函数值update_jfind(new_flast_f);P(update_j,:)X(update_j,:);%修正每个粒子的历史最优值[new_fmin,min_i]min(new_f);new_PgX(min_i,:);Pg(new_fminfmin)*new_Pg(new_fminfmin)*Pg;last_fnew_f;%保存当前的函数值fminmin(new_fmin,fmin);%更新函数最小值Savef(iter)fmin;if mod(iter,floor(n/10))0%每10代记录一次种群参数NN1;if ifdrawSaveData(:,:,N)X;endendww-iter/n*eta*w;
end
if ifdrawfor j1:10figure(j)plot(SaveData(:,1,j),SaveData(:,2,j),o);xlim([xmin(1),xmax(1)])ylim([xmin(2),xmax(2)])end
end
figure % 另起一个图
plot(Savef,b-)
title(混合粒子群优化)
disp(Pg)
disp(fmin)
end注代码注释部分是绘制数次迭代粒子分布位置要解除需要一起解除CtlrT,否则会报错
2.2.执行与效果
\qquad我们首先来尝试使用RosenBrock函数香蕉函数对其进行试验。香蕉函数的极小值点为(1,1)极小值为0。书写为匿名函数如下 f(x)((1-x(:,1)).^2100*(x(:,2)-x(:,1).^2).^2);完整的测试例程M文件
close all
clear
f(x)((1-x(:,1)).^2100*(x(:,2)-x(:,1).^2).^2);
Pg Hybrid_PSO(f,2,200,50,[10,10],[-10,-10],[3,3],[-3,-3],2,1,1,3,0.3,0.8,true);
disp([Pg:,num2str(Pg)])
disp([fmin:,num2str(f(Pg))])运行窗口输出 start_HPSO1.0000 1.00005.9764e-12Pg:1 1
fmin:5.9764e-12注这里的冒号不可以省略因为是样本为多维行向量组成的矩阵所以函数也应返回同维向量而不是一个值。
我们设定一个稍大的范围速度限幅设置为位置限幅约1/3的大小局部社会因子取0.3局部社会因子影响半径R取位置限幅的1/6。仿真命令及结果如下 目标函数的下降趋势如下
粒子群的变化趋势如下图
①②③
可以发现粒子群算法顺利收敛到了(1,1)附近。
希望本文对您有帮助谢谢阅读
