现代控制论文-系统辨识 贝叶斯

基于贝叶斯方法的系统辨识、贝叶斯、最大后验概率估计叶斯的起源和基本原理对于统计决策函数、统计推断、统计的估算等做出了贡献.1763年发表了这方面。贝叶斯(ReverendThomasBayes1702-1761)学派奠基性的作品是贝叶斯的的数学家拉普拉(Laplace,P.S)用贝叶斯的方法导出了重要的“相继律”，贝叶斯意大利的菲纳特(B.deFinetti)及其英国的杰弗(Jeffreys,H.)都对贝叶斯学派的第二次世界大战后，瓦尔德(Wald,A.)提出了统计的决策理论，在这一理1958年英国最悠久的统计杂志Biometrika全文重新刊登了贝叶斯的论文，20世RobbinsH为代表，提出了经验贝叶斯方法和经典方法P(BA)P(A)P(BA,C)P(AC)P(BP(BA)P(A)P(BA,C)P(AC)P(B)P(BC)由已知的先验概率P(BA)求出P(AB)Zk=z(k),u(k),z(k-1),u(k-1),,z(1),u(1),z(0),u(0)}Zkz(k-1),u(k-1),z(k-2),u(k-2),,z(1),u(1),z(0),u(0)}Zk=z(k),u(k),Zk-1}p(9|Zk)=p(9|z(k),Zk-1)=p(9|Z1)=p(9|z(1),z0)=原则上说，根据(1)式可求得θ的后验概率密度函数，但实际上这是有困难有可能得到(1)式的解析解。求得参数θ的后验概率密度后，就可以利用它进一参数估计方法2.1极大后验参数估计方法极大后验参数估计方法就是把后验概率密度函数p(9|Zk)，达到极大值作为Maxp(9|Zk)orMaxp(9|Zk)or?9 kalogp|Zk-1)9ˆ=0(4)MAPMAP若让式(4)第一项取0，则对应的估计就是极大似然估计。可见，极大后2.2条件期望参数估计方法3.1最小二乘模型的Bayes参数辨识kzk,-zT(k-n),u(k-1),u(k-m)]TP(9Zk)=(6)Zi，均值为，协方差阵为 P0lJ P0lJ(9)p(z(k)|9,Z)=1exp|-p(9|Z)=Norm.exp|-1「z(k)-QNorm：Norm=(2)-N/2pzkZk将(11)式两边取对数后分别对θ求导，得：?9Lk」k-1Lk-1」?logp(9|Zk)=1Q「z(k)-Q?9Lk」k-1Lk-1」由(13)式计算得：TT-1ˆ1T-1TII-2Qk令1Kk=PkQk，则得出：k22k2k-2logp(9|Z)=Const+2(k)-(k)0T9-1k22k2k2kkk-1k-1k-1k-1k-1k-1k-1k-1+19T02kkk-1k-1k-1k-1k-1k-1k-1k-1PPTP-9T0kz(k)+9ˆ1k-1-0kz(k)+9ˆ1k-1]|T9=Const，+「|9-1P0z(k)-PP-19ˆ]|TP-1|「9-1P0z(k)-PP-19ˆ]|L|kkkk-1k-1」|kL|2kkkk-1k-1」|3.2贝叶斯估计收敛性k2kk利用K=1k2kkz(k)-0T9ˆ=z(k)-0T{9ˆ+K[z(k)-0T9ˆ]}=(1-0T)1(k)k-0T]kk-1kkkk-1=(1-0kTPk0k)[z(k)-0T9ˆ]2kk-1E(2)=E[D(z(k))]=E[0T9-0T9ˆ]2=E[0T(9-9ˆ)]2kkkkkkkkkkkkkk=E[0T(9-9ˆ)(9k-)T]k=0Tkkkkkkkkkkkkkk4.1.仿真模型V(k)e(k)++u(u(k)z(k)yz(k)(|A(z-1)=1-1.5z-1+0.7z-2=C(z-1)zzu(k)G(z-1)+v(k)N(z-1)=z(k)u(k)+v(k)=z(k)u(k)B(z-1)+v(k)D(z-1)=z(k)A(z-1)ukzzv(k)(1-z-1+0.2z-2)=z(k)(1-1.5z-1+0.7z-2)1.0u(k)z-1+0.5u(k)z-2+v(k)-v(k)z-1+0.2v(k)z-2=z(k)-1.5z(k)z-1+0.7z(k)z-24.2程序流程图应及模型敛情况NY应及模型敛情况NY给被识别的e和P赋初值参数bb1bb 2dd1dd2dd3值析态，即a1=-1.5000,a2=0.7000,b1=1.0000,b2=0.5000,d1=1.0000,例的Bayes最小二乘递推算法的辨识结果与增广最小二乘递推算法的辨识结果sy1;y2=1;y3=1;y4=0;%四个移位积存器的输出初始值xxory,y4);%第一个移位积存器的输入信号xy%第二个移位积存器的输入信号xy%第三个移位积存器的输入信号xy%第四个移位积存器的输入信号y(i)=y4;%第四个移位积存器的输出信号，幅值"0"和"1"ify(i)>0.5,u(i)=-1;%M序列的值为"1"时,辨识的输入信号取“-1”elseu(i)=1;%M序列的值为"0"时,辨识的输入信号取“1”endyx1;y2=x2;y3=x3;y4=x4;%为下一次的输入信号作准备figure形subplot(2,1,1);%画第一个图形的第一个子图stem(u),gridon%画出M序列输入信号v=randn(1,60);%产生一组60个正态分布的随机噪声subplot(2,1,2);%画第一个图形的第二个子图plot(v),gridon;%画出随机噪声信号R=corrcoef(u,v);%计算输入信号与随机噪声信号的相关系数r=R(1,2)%取出互相关系数rv=std(v)*std(v)%计算随机噪声的方差u%显示输入信号v%显示噪声信号zzeros0);zm=zeros(1,60);%定义输出采样矩阵与模型输出矩阵的大小z;z(1)=0;zs(2)=0;zs(1)=0;%输出采样、系统实际输出、模型输出赋初值zm;zm(1)=0;%模型输出赋初值sc0=[0.0010.0010.0010.0010.0010.0010.001]';%直接给出被辨识参数的初始值,p0=10^6*eye(7,7);%直接给出初始状态P0，即一个充分大的实数单位矩阵E=5.e-9;%相对误差E=0.000000005c=[c0,zeros(7,59)];%被辨识参数矩阵的初始值及大小e=zeros(7,60);%相对误差的初始值及大小fork=3:60;%开始求Kxh*p0*h1+rv;xinvx);k1=p0*h1*x1;%Kccce2=e1./c0;%求参数的相对变化e(:,k)=e2;c0=c1;%给下一次用c(:,k)=c1;%把辨识参数c列向量加入辨识参数矩阵p1=p0-k1*k1'*x;%findp(k)p0=p1;%给下次用ife2<=Ebreak;%若收敛情况满足要求，终止计算end%判断结束end%循环结束c,e,z,zs,zm%显示被辨识参数、误差情况、输出采样值、实际输出值、模型输出值%分离赋值a1=c(1,:);a2=c(2,:);b1=c(3,:);b2=c(4,:);%分离出a1、a2、b1、b2figure形plot(i,a1,'r',i,a2,'r:',i,b1,'b',i,b2,'b:',i,d1,'g',i,d2,'g:',i,d3,'g+')%画出各个被辨识参数lePa