第八讲基于群体智能的参数辨识

基于群体智能算法的参数辨识遗传算法(GeneticAlgorithm)优胜劣汰遗传过程常用术语

个体(染色体)

群体

复制

交叉

变异

选择

适应度

适应值优化目标目标函数适应度函数优化目标函数适应度函数适应值通过适应度函数来决定染色体的优劣程度遗传编码二进制编码（编码，译码）对每个参数确定它的变化范围，并用一个二进制数来表示，若参数a的变化范围为(amin,amax)，用m位二进制数b来表示，则两者之间满足实数编码选择方法轮盘赌交叉操作(二进制)单点交叉：交叉点k的范围为[1，N-1]，N为个体变量数目，在该点分界相互交换变量多点交叉，m个交叉位置可无重复随机地选择，在交叉点之间的变量间续地相互交换，产生两个新的后代均匀交叉将每个点都作为潜在的交叉点，随机地产生与个体等长的0、1掩码，掩码的片断表明了哪个父体向子体提供变量值交叉操作（实值）子个体1=a*父个体1+(1-a)*父个体2子个体2=(1-a)*父个体1+a*父个体2a为1到0之间的随机数变异（二进制）对于二进制编码个体而言，变异意味着变量的对某个二进制位执行求反运算变异前：101100111001101变异后：101101111001101变异（实数）新个体=旧个体+(rand-0.5)*(amax-amin)

clc;clearall;closeall;Ge=100;%代数m=80;%一次选择个体数目Pc=0.9;%交叉概率Pm=0.05;%变异概率ts=0.5;%采样周期%用于记忆最优个体对应适应值%对应参数值bestQQ=10000;k_s=0;T1_s=0;T2_s=0;%被拟合对象参数st=70;lp=st/ts;yout0=0;yout=0;k1=4;t11=5;t22=9;timef=0;x00=0;x01=0;r=1.0;%阶跃输入a00=exp(-ts/t11);b00=1-a00;a01=exp(-ts/t22);b01=1-a01;forii=1:lpx00=a00*x00+k1*b00*r;x01=a01*x01+1*b01*x00;yout0(ii)=x01;timef(ii)=ii*ts;endforkk=1:Gex1l=3.8;x1r=4.2;%kx2l=4;x2r=6;%T1x3l=7;x3r=9;%T2fori=1:m%产生初始种群

pop1(i)=rand(1)*(x1r-x1l)+x1l;pop2(i)=rand(1)*(x2r-x2l)+x2l;pop3(i)=rand(1)*(x3r-x3l)+x3l;end

%按照误差平方的最小值计算目标函数

fori=1:mk=pop1(i);T1=pop2(i);T2=pop3(i);x=[0,0];yout=0;error=0;Q=0;a1=exp(-ts/T1);b1=1-a1;a2=exp(-ts/T2);b2=1-a2;forj=1:lpx(1)=a1*x(1)+k*b1*r;x(2)=a2*x(2)+1*b2*x(1);yout(j)=x(2);error(j)=yout(j)-yout0(j);Q=Q+(error(j)^2);endQQ(i)=Q;end%计算适应度，由于是直接求最小值问题，所以直接将Q取为个体的适应度

sumQQ=sum(QQ);cfitness=QQ/sumQQ;fori=1:mcfitness(i)=(1-cfitness(i))/(m-1);end%根据适应度大小选择下一代，（比例选择）

fori=2:mcfitness(i)=cfitness(i-1)+cfitness(i);end%选择运算fori=1:mp=rand(1);index=1;while(p>cfitness(index))index=index+1;endnewpop1(i)=pop1(index);newpop2(i)=pop2(index);newpop3(i)=pop3(index);end%放回原来数组里fori=1:m

pop1(i)=newpop1(i);pop2(i)=newpop2(i);pop3(i)=newpop3(i);endfori=1:m%交叉运算

p=rand(1);index1(i)=floor(rand(1)*(m-1))+1;index2(i)=floor(rand(1)*(m-1))+1;index3(i)=floor(rand(1)*(m-1))+1;if(p<Pc)A=rand(1);pop1(i)=pop1(i)*A+(1-A)*pop1(index1(i));B=rand(1);pop2(i)=pop2(i)*B+(1-B)*pop2(index2(i));C=rand(1);pop3(i)=pop3(i)*C+(1-C)*pop3(index3(i));endend

%变异运算fori=1:mp=rand(1);ifp<Pmpop1(i)=pop1(i)+Leng1*(rand(1)-0.5)*(x1r-x1l);pop2(i)=pop2(i)+Leng2*(rand(1)-0.5)*(x2r-x2l);pop3(i)=pop3(i)+Leng3*(rand(1)-0.5)*(x3r-x3l);if(pop1(i)>x1r)pop1(i)=x1r;endif(pop1(i)<x1l)pop1(i)=x1l;endif(pop2(i)>x2r)pop2(i)=x2r;endif(pop2(i)<x2l)pop2(i)=x2l;endif(pop3(i)>x3r)pop3(i)=x3r;endif(pop3(i)<x3l)pop3(i)=x3l;endendend

%循环结束，求出最优个体输出a1=exp(-ts/T1);b1=1-a1;a2=exp(-ts/T2);b2=1-a2;x=[0,0];yout=0;forj=1:lpx(1)=a1*x(1)+k*b1*r;x(2)=a2*x(2)+1*b2*x(1);yout(j)=x(2);endfigure;plot(timef,yout,‘b--’,timef,yout0,’k’);figure;Plot(errorhistor);%一次迭代结束，%找出最优个体存储下来[f,I]=min(QQ);iff<bestQQbestQQ=f;k_s=pop1(I);T1_s=pop2(I);T2_s=pop3(I);enderrorhistory=[errorhistorybestQQ];end%TPSO(ParticleSwarmOptimization)ParticleSwarmOptimizationTom:WhereshouldIflyto?Tom:Thisisthebestplaceinmysearchingfoodhistory.John:Thisisbestplaceforourswarmnow.PSO算法采用速度一位置搜索模型。每个粒子代表解空间的一个候选解，解的优劣程度由适应度函数决定，而适应度函数根据优化目标定义。PSO算法随机初始化一群粒子，然后通过迭代找到最优解。在每一次迭代中，粒子通过跟踪两个“极值”来更新自己。一个是粒子本身所找到的最优解，即个体极值pbest；，另一个是整个种群目前找到的最优解，称之为全局极值gbest。PSO的初始化（1）设定群体规模m，确定搜索空间[xmin,xmax]和单步前进的最大值vmax（一般应小于最大搜索空间的20%）；（2）对任意粒子，在[xmin,xmax]内服从均匀分布产生个体xi；（3）对任意粒子，在[-vmax,vmax]内服从均匀分布产生个体vin；（4）对任意粒子，设Pi=Xi。%对一个系统进行仿真，并且计算绝对误差的矩积分，评价函数functionQ=jianmoobj(k,T1,T2);globaltimefyoutyout0rints=0.5;st=100;lp=st/ts;rin=1.0;x0=[0,0];x=[0,0];yout0=0;yout=0;error=0;Q=0;k0=1.5;t01=6;t02=16;a01=exp(-ts/t01);b01=1-a01;a02=exp(-ts/t02);b02=1-a02;a11=exp(-ts/T1);b11=1-a11;a12=exp(-ts/T2);b12=1-a12;fori=1:lptimef(i)=i*ts;r(i)=rin;x(1)=x(1)*a11+k*b11*r(i);x(2)=x(2)*a12+1*b12*x(1);yout(i)=x(2);x0(1)=x0(1)*a01+k0*b01*r(i);x0(2)=x0(2)*a02+1*b02*x0(1);yout0(i)=x0(2);error(i)=yout(i)-yout0(i);Q=Q+(error(i)^2);end%绘制最终的输出曲线functionQ=jianmoout(k,T1,T2);globaltimefyoutyout0rints=0.5;rin=1.0;x0=[0,0];x=[0,0];yout0=0;yout=0;error=0;st=100;lp=st/ts;k0=1.5;t01=6;t02=16;a01=exp(-ts/t01);b01=1-a01;a02=exp(-ts/t02);b02=1-a02;a11=exp(-ts/T1);b11=1-a11;a12=exp(-ts/T2);b12=1-a12;fori=1:lptimef(i)=i*ts;r(i)=rin;x(1)=x(1)*a11+k*b11*r(i);x(2)=x(2)*a12+1*b12*x(1);yout(i)=x(2);x0(1)=x0(1)*a01+k0*b01*r(i);x0(2)=x0(2)*a02+1*b02*x0(1);yout0(i)=x0(2);endfigure(2);plot(timef,rin,'r',timef,yout,'b',timef,yout0,'g');xlabel('Time(s)');ylabel('rin,yout,yout0');%PSO算法closeall;clearall;globaltimefyoutyout0rin%随机初始化Gp=50;%循环代数errorhistory=[];%初始化各个粒子x=rand(30,3);[n,m]=size(x);%初始化各个粒子的第一步前进的速度v=rand(size(x));dd=0.2;xmin=[1514];xmax=[2718];vmax=dd.*(xmax-xmin);fori=1:nforj=1:mx(i,j)=rand(1)*(xmax(1,j)-xmin(1,j))+xmin(1,j);v(i,j)=rand(1)*vmax(1,j);endendfori=1:ny(i)=jianmoobj(x(i,1),x(i,2),x(i,3));%计算原始粒子的目标值

py=zeros(1,n);endp=x;%原始粒子的历史最优位置既是本身ymin=inf;%原始粒子中的最优值k=1;%找出未迭代前所有粒子中的最优值fori=1:nif((ymin)>y(i))ymin=y(i);k=i;endendg=x(k,:);%未迭代前原始粒子中的最优值c1=1.5;%常数c2=1.5;y0=y;%保留各个粒子目标值的前一时刻值gy=inf;%全局目标最优值fork=1:gp%惯性因子

%w=1-0.8*k/lp;w=0.73;%个体历史最优值

fori=1:nif(y0(i)>=y(i))forj=1:mp(i,j)=x(i,j);%p为个体最优值

endpy(i)=y(i);%py为个体最优目标值

endendy0=y;%全局最优值

fori=1:nif((gy)>(py(i)))g=p(i,:);gy=py(i);endend%计算个体前进速度，并限制最大值

fori=1:nforj=1:mv(i,j)=w*v(i,j)+c1*rand(1)*(p(i,j)-x(i,j))+c2*rand(1)*(g(j)-x(i,j));end%速度最大值限制

forj=1:m