系统辨识各类最小二乘法汇总

1、最小二乘方法的总结与比较Y30150638控制二班王思明（作业二）最小二乘法大约是1795年高斯在他那著名的星体运动轨迹研究工作中提出来的。后来，最小二乘法就成了理论的奠基石。由于最小二乘原理简单，编程程序也不困难，所以它颇受人们重视，应用相当广泛。最小二乘的基本结果有两种形式，一种是经典的一次完成算法，另一种是现代的递推算法。本次作业只讨论递推算法和其他扩展的最小二乘方法比较。一、最小二乘法最小二乘法具有两个缺陷：1、当噪声模型是有色噪声时，最小二乘参数估计不是你无偏、一致估计。2、随着数据的增长，最小二乘法出现数据饱和现象，以致算法慢慢失去修正能力。下面是其分别在白噪声和有色噪声时的参数估

2、计曲线。最小二乘的程序见附录。图1-1RLS白噪声时参数估计变化曲线QFigure4fileEditViewInsertTeoDesktopWindowHelpClUUJ氐1%踏工/氢E国图1-2RLS有色噪声时参数估计变化曲线从图中可以很明显的看出，当噪声变成有色噪声时，RLS算法是有偏的，始终不能得出参数的真值，而且可以看出随着数据的增多，算法对参数的修正能力明显下降，到最后参数值几乎不变，这说明RLS算法进入了数据饱和状态。二、遗忘因子法和限定记忆法这两种算法（RFERFM）是基于RLS的数据饱和现象而提出来的一种辨识算法。RFF很简单，在老数据加上遗忘因子以降低老数据提供的信息量，增加

3、新数据的信息量。仅仅在RLS兑法内部加入遗忘因子系数，这里因为篇幅有限，不给出其程序以及最后的参数变化曲线。RFM就比较特殊，其限定有效数据的长度L,即始终只有长度L的数据在影响参数估计，可以从根本上解决数据饱和现象。下面是其实验的参数估计变化曲线。图2-1RFM白噪声时的参数估计曲线从图中可以看出，相比于RLS的白噪声时的参数估计，变化曲线围绕真值震荡变化，这说明没有数据饱和现象，辨识结果良好。三、偏差补偿法当噪声是有色噪声时，我们提出偏差补偿最小二乘法，因为从图12我们可以看出，最后的参数估计是有偏的，所以我们在算法中每次都减去这个偏差量,使得最后的结果能够逼近真值。我们称这种算法为RCL

4、So但是，既然我们加入了偏差量，如果我们的噪声是白噪声时，那这个偏差量就变成了误差量，使得参数估计轨迹偏真值。下图是其参数估计曲线，程序见附录。图3-1RCLS有色噪声时的参数变化曲线图,2RCLS白噪声时的参数估计曲线从图中可以很清楚的看出，RCLS可以对有色噪声情况下，对参数进行无偏一致估计，但是在白噪声的情况下，却不能得到无偏估计，其原因上面也提到了,这里就不在赘述了。四、增广最小二乘和广义最小二乘针对不同的有色噪声模型，我们用不同的辨识方法，并能顺带求出噪声的模型模型，通过扩充参数向量和数据向量:的维数，把噪声模型的辨识同时考虑进去。当噪声模型：e(k)=D(Z-1)*v(k)(v(k

5、)为白噪声)时，我们采用增广最小二乘方法。其实验结果如下，其程序见附录。图4-1RELS参数估计曲线当噪声模型：e(k)=i*v(k)(v(k)为白噪声)时，我们采用广义最C'Z)小二乘方法。其实验结果如下，其程序见附录。图4-2RGLS参数估计曲线从上图可以看出，针对不同的参数模型，RELS和RGLS分别能发挥不同的作用，并且给出参数的无偏一致估计，并顺带求出了噪声参数，其辨识结果良好，符合要求.五、辅助变量法在有些问题中，我们并不能事先知道参数的模型，而且我们也不需要求出噪声参数，这时普通的RLS算法已经不能满足我们的要求，我们提出了辅助变量法，RIV.下图给出其实验结果，具体程序

6、见附录。DFigure4-deFileEditViewInsertJoolsDesktopVtfndowHelp引JGIiblQI0AO电娓/|。|。目|口参数变化曲线图5-1RIV有色噪声时参数估计变化曲线从图中可以看出，RIV算法对有色噪声的情况，也能给出参数的无偏估计,其结果与RLS（图12）比较可以看出区别。六、总结RLS：当噪声模型是有色噪声时，最小二乘参数估计不是你无偏、一致估计。而且随着数据的增长，最小二乘法出现数据饱和现象，以致算法慢慢失去修正能力。RFF.RFM：RFF很简单，在老数据加上遗忘因子以降低老数据提供的信息量，增加新数据的信息量，RFM就比较特殊，其限定有效数据的

7、长度L,即始终只有长度L的数据在影响参数估计，可以从根本上解决数据饱和现象RCLS：当噪声是有色噪声时，我们提出偏差补偿最小二乘法，在算法中每次都减去这个偏差量，使得最后的结果能够逼近真值。但却不适合噪声为白噪声的时候。RELS：当噪声模型：e(k)=D(Z-1)*v(k)(v(k)为白噪声)时，我们采用增广最小二乘方法。能辨识出参数(包括噪声参数)的无偏估计。RGLS：当噪声模型：e(k)=-4r-*v(k)(v(k)为白噪声)时，我们采用广义最C(Zx)小二乘方法。能辨识出参数(包括噪声参数)的无偏估计。RIV：当不确定参数模型时，我们可以直接采用RIV模型得到参数的无偏一致估计。七、附录

8、1、RLS实验程序：clc%产牛高斯白噪声clear%yl=randii(l2500)：%yl=yl/std(yl);%yl=yl-mean(yl);%a=0;%b=sqrt(0.1);%yl=a+b*yl;%figure(l);%plot(yl);%title(高斯白噪声);fori=l:2500%产生有色噪声yi(i)=i;endN=400;uk=rand(l,N);fori=l:N%产生随机输入uk(i)=uk(i)*(-l)A(i-l);endfigure。)；plot(uk);titlef随机输入曲线);yk=zeros(l,N);fork=3:N%定义输出yk(k)=l.5*yk(

9、k-l)-0.7*yk(k-2)+uk(k-l)+0.5*ukQc-2)+yl(k);endfigure。)；plot(yk);title,对应输出曲线)；theta=O;O;O;O;p=10A6*eye(4);fort=3:N%迭代计算（Is算法）h=(-yk(t-l);-yk(t-2);uk(t-l);uk(t-2);x=l+h,*p*h;p=(p-p*h*l/x*h，+p);theta=theta+p*h*(yk(t)-hr*theta);alt(t)=theta(l);a2t(t)=theta(2);blt(t)=theta(3);b2t(t)=tlieta(4);endthetat=

10、l:N;%绘制参数变化曲线figure(4);titleC参数变化曲线)；text(300,1.5,'al);text(300,0.7/a2r);text(300,l/bl);text(300,0.5,1529;2、RFM实验程序（部分）只给出算法更改部分：fors=3:100;h=(yk(s-l);yk(s-2);uk(s-l);uk(s-2);p=pp*h*(inv(l+lT*p*h)*h*p;theta=theta+p*h*(yk(s)-ir*tlieta);endsfort=3:N-100m=s-H;h=(yk(m-l);yk(m-2);uk(m-l);uk(m-2);p=p-

11、p*h*(inv(l+li,*p*h)*h*p;tlieta=theta-Fp*h*(yk(m)-h,*tlieta);h=(yk(t-l);yk(t-2);uk(t-l);uk(t-2);p=p+p*h*(inv(l-h,*p*h)*h,*p;theta=theta+p*h*(yk(t)-ir*theta);alt(t)=theta(l);a2t(t)=theta(2);blt(t)=tlieta(3);b2t(t)=tlieta(4);end3、RCLS实验程序(部分)theta=O;O;O;O;p=10A6*eye(4);j=0；thetc=O;O;O;O;fort=3:Nh=(-yk(

12、t-l);-yk(t-2);uk(t-l);uk(t-2);x=l+l产p*h;j=j+(yk(t)-h*theta)A2)/(l-Hi,*p*h);p=(p-p*h*l/x+h,*p);theta=theta+p*h*(yk-htheta)d=l000;0100;0000;0000;w=j/(t*(14-thetcr*d*theta);thetc=theta-t*w*p*d*thetc;alt(t)=thetc(l);a2t(t)=thetc(2);blt(t)=tlietc(3);b2t(t)=thetc(4);end4、RELS实验程序(部分)fort=3:Nx=l+iV*p*h;p=(

13、p-p*h*l/x*hr*p);theta=theta+p*h*(yk(t)-h,*theta);alt(t)=theta(l);a2t(t)=theta(2);blt(t)=tlieta(3);b2t(t)=tlieta(4);dlt(t)=theta(5);d2t(t)=theta(6);end5、RGLS试验程序(部分)fort=3:Nhe=(-e(t-l);-e(t-2);xe=H-he，+pe+he;pe=(pe-pe*he*l/xe*he,*pe);thete=thete+pe*he*(e(t)-he*thete);clKt)=thete(l)；c2t(t)=thete(2);endfort=3:Nykf(t)=yk(t)+clt(t)*yk(t-l)4-c2t(t)*yk(t-2);ukf(t)=iik(t)+clt(t)*uk(t-l)+c2t(t)*uk(t-2);h=(-ykf(t-l);-ykf(t-2);ukf(t-l);ukf(t.2);x=l+hr*p*h;p=(p-p*h*l/x*h，+p);theta=tlieta+p*h*(ykf(t)-h'*theta);alt(t)=theta(l);a2t(t)=theta(2);blt(t)=theta(3);b2t(t)=tlieta(4);