系统辨识大作业..docx

一、问题描述考虑仿真对象：Z(Z)+0.9z(A-1)+0.15z(Z-2)+0.02Z(A一3)=0.7“(A-1)-0.15(A-2)+贸衣)e()+1.0e(-1)+0.4le(-2)=Av(Jt),vN(OJ)式中，u(k)和z(k)是输入输出数据，v(k)是零均值、方差为1的不相关的随机噪声；u(k)采用与e(k)不相关的随机序列。1.设计实验，产生输入输出数据；2使用基本最小二乘法估计参数；3.考虑其他适用于有色噪声的辨识方法估计参数；4.模型验证.二、问题分析对于单输入单输出系统(SingleInputSingleOutput,SISO),如图1所示,将待辨识的系统看作"灰箱”，它只考虑系统的输入输出特性，而不强调系统的内部机理。图1中，输入u(k)和输出z(k)是可以测量的，G(ZT)是系统模型，用来描述系统的输入输出特性，y(k)是系统的实际输出。N(ZT)是噪声模型，v(k)是均值为零的不相关随机噪声，e(k)是有色噪声。图1SISO系统的"灰箱”结构对于SISo随机系统，被辨识模型G(Z)为：Xz)_bz'+b2z7+bnznm(z)+aiz+a2z'+>+anz其相应的差分方程为义幻=-aiy(k-i)+NbiU(k-z)若考虑被辨识系统或观测信息中含有噪声，被辨识模型可改写为2也)=-Xaj仕-£)+iu(-i)+v()E(=1式中，z(k)为系统输出量的第k次观测值;y(k)为系统输出量的第k次真值,y(kT)为系统输出量的第kT次真值，以此类推；u(k)为系统的第k个输入值，u(k-l)为系统的第k-1个输入值；V(k)为均值为O的不相关随机噪声。1.数据生成本部分需要生成系统的输入输出数据以及噪声数据。1)白噪声的生成辨识数据通常包含有噪声，如果该噪声相关性较弱或者强度很小，可近似看作白噪声。本次实验问题中vWN(0,1),即服从标准正态分布，可以将噪声看作为服从正态分布的白噪声过程，在Matlab中口1以由randn函数生成。2)输入数据的生成伪随机二进制序列(PseudoRandomBinarySequence,PRBS)是广泛应用的一种伪随机序列，所谓“二进制”是指序列中每个随机变量只有"0"或"1"两种逻辑状态。伪随机二进制序列可由多级线性反馈移位寄存器组成的随机信号发生器产生，其中具有最长循环周期的线性移位寄存器序列是伪随机二进制序列最常见的一种形式，简称M序列(MaximalLengthSequence)oM序列由于具有近似白噪声的性质，而且工程上易于实现，能筋保证较好的系统辨识精度，是普遍采用的一种辨识用输入信号。图2线性反馈移位寄存器产生伪随机二进制序列结构图以一个4级线性反馈移位寄存器产生伪随机二进制序列为例，如图2所示。假设4个移位寄存器/,q,%,%输出的初态非全零，移位寄存器的工作原理是:一个移位脉冲来到后，每级移位寄存器的输出移到下一级移位寄存器作为输入,最末一级移位寄存器的输出即为伪随机二进制序列。3)输出数据的生成根据给定的SISO系统，可以求出z(k)的表达式:Z（八）=-0.9z(-1)-0.15z(-2)-0.02z(-3)+0.7w(-1)0.15“(左一2)+NV（八）1.O(1)0.41e(k2)其理想系数值为l=0.9,2=0.15,3=0.02,b1=0.7,=-15,c1=1.0,c2=0.41.可以根据生成的白噪声序列和输入序列，以及必要的O初始值，带入表达式即可得到采样输出数据。2.差分模型阶检验在实际场景中，辨识模型的阶数和纯时延往往是未知的，在很多情况下仅仅依靠猜测。在模型的阶数和纯时延不确定时，设系统模型为y(0=-aiy(f£)+fh.y(t-i)+(t);=1J=其中n为模型的阶数,7(f)=C(ZT)e(f)模型的阶估计可以采用多种方法，本实验采用比较简单易行的损失函数检验法。定义预报误差(噪声方差的估计值)的平方和为损失函数，即当n从小增大时，JN应随之减小，当n增大到某一值时，人应近似白噪声过程。采用以下的检验原则：在n-1这一点，人最后一次出现陡峭的下降，此后就近似地保持不变或只有微小的下降，则取6="。3.参数辨识模型在系统辨识和参数估计领域中，最广泛的估计方法时最小二乘法和极大似然估计法。最小二乘法作为一种最基本的估计方法应用极为广泛，其他的大多数估计算法都与最小二乘法有关。它既可用于动态系统，也可用于静态系统；可用于线性系统，也可用于非线性系统；可用于离线估计，也可用于在线估计。在随机环境下利用最小二乘法时，无须知道观测数据的概率统计信息，而这种方法获得的结果，却有相当好的统计性质。最小二乘参数估计方法来源于数理统计的回归分析，它能提供一个在最小方差意义上与实验数据最好的你和的模型。该估计在一定条件下有最佳的统计特性，即它们是一致的，无偏的和有效的。最小二乘法时一个经典的方法，概念简明，适应范围广,在一些情况下，可得到与极大似然法一样好的统计效果，它能很方便地与其他辨识算法建立关系.(1)递推最小二乘算法当获得一批数据后，可一次求得相应的参数估计值，这样处理问题的方法称为一次完成算法或批处理算法。它在理论研究方面有许多方便之处，但当矩阵的维数增加时，矩阵求逆运算的计算量会急剧增加，将给计算机的计算速度和存储量带来负担，而且不适合在线辨识，无法跟踪参数随时间变化的情况。为了减少计算量，减少数据在计算机中所占的存储量，也为了实时地辨识出动态系统的特性，在用最小二乘法进行参数估计时，把它转化成参数递推的估计。参数递推估计是指被辨识的系统，每取得一次新的测量数据后，就在前一次估计结果的基础上，利用新引入的测量数据对前一次估计的结果进行修正，从而递推地得出新的参数估计值。这样，随着新测量数据的引入，一次接一次地进行参数估计，直到估计值达到满意的精确程度为止。最小二乘递推算法的基本思想可以概括为：当前的估计值0(外=上次估计值0伙-1)+修正项即新的估计值灰6是在旧的估计值灰%-1)的基础上，利用新的观测数据对旧的估计值进行修正而成的。考虑如下模型：A(ZT)Z伏)=B(zT)"(幻+v(Q,其中"伏)，Z（八）分别是系统的输入和输出；Wk)是均值为零，方差为一的不相关白噪声。且满足：AgT)=】+“/T+¾z-+anznB(z=biz-'+b2z-2+.+bmz-m人/力(&)=I-Z(&-1),.,-Z(A-),“(-1),."伏-m)=a,a2,.an,bl,b2,.m则使用加权最小二乘参数估计递推算法(RecursiveWeightedLeastSquares,RWLS)有：k)=灰Jt-1)+K(k)z(k)-hr(k)(k-1)jA伏)P(6=-K(k)h(k)P(k-1)伏)=z(1),.z(),u(k1),"伏m)M>Q 比 fit")J#小二*蒙数估计 递指整法图3最小二乘递推辨识参数估计过程中信息的变换可以看出，取A(Jl)=I的时候，加权最小二乘估计就退化成了最小二乘参数估计的递推算法(RecursiveLeastSquares,RLS)。加权参数L可以在(0,1范围内选择,A如果选一=1,所有的采样数据都是等同加权的，如果一<<1,则表示对新近获A(Z)A(L)得的数据给予充分大的加权因子，而削弱历史观测数据的影响。实际计算时，需要首先确定初始参数反0)和P(O)。般说来选取ZX为充分大实数£为充分小的向量对于这样的系统，使用最小二乘法参数估计的递推算法进行辨识可以得到无偏估计，但是如果噪声模型必须用C(ZT)V(Q表示时，此时就无法得到无偏估计了，因为该方法没有把噪声模型考虑进去。使用递推最小二乘法进行系统辨识的过程如图4所示。(2)增广最小二乘法当噪声均值为O时，最小二乘参数估计算法为无偏估计；当噪声的均值不为O时,最小二乘参数估计算法为有偏估计。为了解决最小二乘参数估计的有偏性，将噪声模型的辨识同时考虑进去，因此称为增广最小二乘法。该算法可以看成是对一般最小二乘参数辨识算法的简单推广或扩充，因此又称为扩充最小二乘算法。考虑如下模型：A(z')z()=B(z')u()+C(z')v(),其中(，Z(Z)分别是系统的输入和输出；v(k)是均值为零，方差为一的不相关白噪声。且满足:A(zi)1+izi+qz°+.+az'n,B(zl)blzi+b11.+bmzmC(Z)=1+CZ'+CjZ+.+Cl,Z'在模型阶次v%U已经确定的情况下，令力(k)=-z(k-1),一Z(A-“),“(-1),.“(-w),v(-l),.v(-v)r、6=,%,a”，仇,&，*”,，J，C2,，J由于v(6是白噪声，所以利用增广最小二乘法这一形式上的变换，即可获得参数e的无偏估计。不过由于数据向量6(左)中包含着不可测的噪声量v(t-l),v(2),.v(-v)f因此要用相应的估计值代替。当比W0,v（八）=0;当k>0,V（八）=Z（八）-必(&)0(A-I),即可得递推公式如下：(k)=(k-)+K(k)v(k)K（八）=P(k-)h(k)lr(k)P(k-1)力(&)+1,<P(k)=l-K(k)h(k)P(K-)v(k)=z(k)-h(k)(k-)h(k)=-z(k-),.z(k-n),-m(-1),.-w(-w),-v(-1),.-v(-v)7一般说来选取初值"°)=/伙O)=O增广最小二乘递推算法扩充了最小二乘法的参数向量8和数据向量/7(行的维数，把噪声模型的辨识同时考虑进去。最小二乘法的许多结论对它都是适用的，但最小二乘法只能获得模型的参数估计，对于有色噪声，也就是噪声模型必须用C(ZT)Ml)表示时，只能采用增广最小二乘递推算法，方可获得无偏估计。这是最小二乘参数估计的递推算法不可代替的。增广最小二乘算法的流程如图5所示。图4递推最小二乘参数辨识流程图图5增广最小二乘算法流程图三、MatIab仿真及运行结果在本次实验中，入取值为1,则系统辨识对象为：Z（八）=0.9Z(Z1)-0.15z(A2)0.02z(k3)+0.7<(左1)1.5(92)+%v(k)-1.0e(k-I)-0.4le(k-2)1.差分模型的阶检验在实际应用场景中，模型的阶数和纯时延往往是未知的，在此种情况下，-般采用如下方法来确定模型的阶数：固定一个阶数，当n从小增大时，定义损失函数JN应随之减小，当n增大到某一值时，人应近似白噪声过程。在nT这一点，JN最后7= H »使用一次出现陡峭的下降，此后就近似地保持不变或只有微小的下降，则取Matlab进行模型阶数确定代码如下：%魁%噪声方差的估计值和损失函数法定%¾clearcloseal1;1.=100:yl=l;y2=l;y3=l;y4=0;fori=l:L+5:xl=xor(y3ty4);%第一个移位寄存器的输入信号x2=yl;%第二个移位寄存器的输入信号x3=y2;%第三个移位寄存器的输入信号x4=y3;%第四个移位寄存器的输入信号yi)=y4;%第四个移位寄存器的输出，即M序列，幅值为"0”或"Iwify(i)>0+5:U=T；%M序列的值为1时，辨识的输入取TelseU=1;%M序列为0时，输入取1endyl=xl:y2=x2:y3r3;y4=x4;%移位寄存器的输出序列endv=randn(L+5,1);%ksi(k)=l三da*v(k)-ksi(k-l)-0,ll*ksi(k-2)%ksi(2)=Ojksi(I)=0:MakSi的前两个初始值为0Iawda=1;fork=3:L+5:ksi(k)=lamda*v(k)-ksi(k-l)-0.41*ksi(k-2);endz(3)=0;z(2)=0;z(l)=0;%取Z的前三个初始值为零fork二4:L+5;z(k)=-0.9*z(k-l)-0.15*z(k-2)-0.02*z(k-3)+0.7*u(k-l)-1.5*u(k-2)+ksi(k);%理想辨识输出采样信号end%模型阶次n=lfor1=1:1.ITKi.1)=2(i):Hl(i,2)=u(i):endestimatel=inv(Hl,*Hl)*Hl,*z(2=L+D,;Dl=(z(2±+0,-Hl*estimatel)'*(z(2:L+l)'-Hl*estimatel)L;%噪声方差的估计值AICl=Ltlog(Dl)+4*1;%模型阶次n=2for1=1:1.II(i,l)=s(itl);:H2(i,2)=z(i);IT(if3)=u(H-l):H2i,d=u；endestimate2=inv(H2'*112)*II2,*z(3:L+2)'D2=(z(3:L+2),-H2*estimate2),*(z(3:L+2),-H2*estimate2)/L;%噪声方差的估计值AIC2=L*log(D2)+4*2:%模型阶次n=3fori=l:LH3i,l)=zi+2);II3(i,2)=z(i+l);H3i,3)=zi;IB(i,4)=u(it2):I3(i,5)=u(H-l);H3.(i16)=u(i);endestimate3=inv(H3,*H3)*H3,*z(4:L+3),:D3=(z(4:L+3)'-H3*estimale3)'*(Z(4:L+3)'-H3*estimate3)/L;噪声方差的估计值AIC3=L*log(D3)+4*3;%模型阶次产4fori=l:LH4(itl)=z(i+3);H4(i,2)=z(i+2):II4(i,3)=z(i÷l);H4(i,4)=z(i);H4(i,5)=u(i÷3):H4(i,6)=u(i+2);H4(i,7)=u(i÷l);H4(i,8)=u:endestimate4=inv(H4,*IM)*II4,*z(5:L+4)*;04=(5:1+4)*-H4*estimate4)*(z(5:L+4)f-H4*estimate4)L;%嗓声方差的估计值AIC4=L*log(D4)+4*4;%模型阶次n=5fori=l:LH5(i,l)=z(i÷4);H5(i,2)=z(i*3);H5(i,3)=z(i+2);II5(i*4)=z(i÷l)：H5(i,5)=z(i);H5(i>6)=u(i+4):H5(i7)=d(i÷3);H5(i,8)=u(i÷2):H5(i,9)=u(i+1);H5(i,10)=u(i);endestimate5=inv(H5,*H5)*H5,*z(6:L+5)':D5=(z(6:L+5),-H5*estimate)*(z(6:L+5)*-H5*estimate5)L;%噪声方差的估计值AIC5=L*log(D5)+4*5;plot(l:5,LAIClAIC2AIC3AIC4AIC5)gridontitle(,损失函数法判断模型阶次)XlabelC阶次')ylabel('损失函数')程序运行结果如图6所示：Figure1。回区Eil Edit ?i“ Insert Xools D “let op xn<io* Klp已日J 4 '一、门卫/，<a 口目 损失函数法判断模型阶次180S523次 阶图6模型阶次的判定由图6可知，当系统阶次从1到5变化时，在n=3及以后，损失函数只有微小的下降，在n=4以后又出现上升的趋势，所以取n=3作为模型的阶次。2.参数辨识(1)递推最小二乘法根据前文所述递推最小二乘法算法流程，编写递推最小二乘法程序如下：%递推最小二乘辨识clcclearall%采用M序列方法产生输入数据1.=90;%90个数据namenda=l;yl=l;y2=l;y3=l;y4=0;u_f=2;fori=l:Lxl=xor(y3ty4);2=yl;*3书2;x4=y3;y(i)=y4;ify(i)>0.5u(i)=-u_f；elseu(i)=u_f;endyl=xl;y2=x2;y3=x3;y4=x4;endv=l*randn(L,1);eb(2)=0;eb(l)=0;fork=3:Leb(k)=-eb(k-1)-0.41*eb(k-2)+namenda*v(k);Cndz(3)=0;z(2)=0;z(1)=0;fork=4:L;Z(k)=-0.9*z(k-l)-0.15*z(k-2)-0.02*z(k-3)+0.7*(k-l)-l.5*u(k-2)+eb(k);endfigured)i=l:30;subplot(2,1,1)stem(u(i);title('输入M序列')gridonSUbPIot1,2)stem(z(i);titlc('输出序列')gridon;%cO=O.OOlO.OOlO.OOl0.0010.001'p0-103*eye(5,5);E=0.0000005;c=c,zeros(5,L-I)：e=zeros(5,L);%fork=4:Lhl=-z-l)z(k-2h-z-3)>u(k-l)1u(k-2)5;kl=pO*hl*inv(h*pO*hl+l);Cl=Co+kl*(z(k)-hl'*c);el=cl-c;e2=eLc;e(:>k)=e2;c=cl;c(:hk)=c1;pl(:,:,k)=p-kl*kl,*h*p*hl+l;=pl(:,:,k);ifabs(e2)<=EN=k;break;endN=k;endal_f=c(l,N)a2_f=c(2,N)a3_f=c(3,N)bl_f=c(4,N)b2_f=c(5,N)al=c(l,l:N);a2=c(2,l:N);a3=c(3,l:N);bl=c(4,1:N);b2=c(5,1:N);eal=e(l,:);ea2=e(2,:);ea3=e(3,:)ebl=e(4,:);eb2=e(5,:);i=l:N;figure(2)plot(i,al,'r',i,a2,'b',i,a3,'y',i,bl,'k',i,b2,'g');legend(*aa2ta3,11bT,'b2');gridon;title(,递推最小二乘参数估计)figure(3)i=l:30:plot(i,eal(k30),ri,ea2(l:30)/b,i,ea3(l:30)/y,i,ebl(1:30)/k,i,eb2(l:30)/g,);gridon;legend(参数al误差'，'参数a2误差'，'参数a3误差'，'参数bl误差','参数b2误差)；title(被估参数误差收敛情况)程序运行结果如图7至图9所示。图7输入及输出信号图9被估参数误差收敛情况使用递推最小二乘算法进行系统辨识的结果如表1所示。表1递推最小二乘算法的辨识结果参数.N真值0.90.150.020.77.5估计值0.93670.1591-0.01380.62757.3558从上表中可以看出，对于含有有色噪声的系统，采用递推最小二乘法进行辨识获得结果有较大误差，并不是无偏估计。所以，一般情况下，对于含有有色噪声的系统应采用增广递推最小二乘算法进行辨识。(2)增广最小二乘法根据前文所述增广最小二乘法算法流程，编写增广最小二乘法程序如下:Iw%l%+ft%嬲%增广最小二乘法的参数辨识算法%嬲WW%tww¾胱嬲closeall%M序列的产生%1.二60:%移位寄存器产生的M序列周期(取L:60)yl=l;y2=l;y3=l;y4=0;%四个移位寄存器输入的初值fori=l:L:xl=xor(y3,y4);%第一个移位寄存器输入的信号x2=yl:%第二个移位寄存器输入的信号x3=y2;%第三个移位寄存器输入的信号x4=y3:%第四个移位寄存器输入的信号y(i)=y4:%第四个移位寄存器的输出，即M序列,幅值为0或Iify(i)>0,5;u(i)=-l:/M序列的值为1时，辨识的输入为TelseU=1;为M序列的值为0时，辨识的输入为1endyl=xl;y2=x2;y3=x3;y4=M;/移位寄存器的输出序列,为下一次的输入信号做准备endfigured);subplot(3,1,D;stcm(u),gridon%画出V序列输入信号径线图并给图形加上网格title，输入信号(M序列)')%随机噪声的产生%v=randn(l,L);先产生一组L:60个正态分布的噪声subplot(3,1,2);plot(V),gridon%画出随机噪声信号title,随机噪声信号e(k)')%实际噪声的产生%ksi(k)=lamda*v(k)-ksi(k-l)-0.41*ksi(k2)%ksi(2) = 0si(l)二 0；%设ksi的前两个初始值为0Iamda=1:fork=3:L;ksi(k)=latnda*v(k)-ksi(k-l)-0.41*ksi(k2);endsubplot(3,1.3)plot(ksi),gridontitlc('实际噪声信号xi(k)')%z(3)=0;z(2)=0;z(l)=0; zm (3) =0; zm (2) =0;zm(l)=0; 时模型输出，模型输出的初值%输出采样，不考虑噪声时系统输出，不考虑噪声辨识过程c0=0.0010.0010.0010.0010.0010.0010.0010.001J,:%直接给出辨识被辨识参数的初始值,即一个充分Po=IO 6*eye(8, 8);小的实向量%直接给出初始向量P0,即一个充分大的实数单位矩阵c= c, zeros (8, Ll):M辨识参数参数矩阵的初始值及大小e=zeros(& L);*相对误差的初始值及大小%开始求Kfork=4:L;z(k)二-0,9*z(k-l)-0.15*z(k-2)-0.02*z(k-3)+0.7*u(k-l)-L5*u(k-2)+ksi(k);%系统在M序列下的输出采样信号hl=-z(k-l),-Z(k2),-Z(k3),u(k-l),u(k-2),v(k),ksi(k-l),ksi(k-2)*;1为求KCk)做准备x=h*p*hl+l;xl=inv(x)kl-p*hl*xl:%Kdl=z(k)-hl,*ccl=c÷kl*dl;%辨识参数cel=cl-c0;e2=el./c0;%求参数误差的相对变化e(：,k)=e2;c=cl:%给下一次用c(:,k)=cl;%把递推出的辨识参数c的列向量加入辨识参数矩阵pl=pO-kl*k*h*p*hl+l;%findp(k)p=pl:%给下次的计算使用end%整个循环结束c,c,先显示被辨识参数及参数收敛情况%分离变量al=c(l,:);a2=c(2,:);a3=c(3,:);bl=c(4,:);b2=c(5,:);%分离出al,a2ta3,bl,b2dl=c(6,:):d2=c(7,:);d3=c(8,:);%分离出的dl,d2,d3%分离变量的收敛情况Cal=C(L:);ea2=e(2,:):ca3=c(3,:);ebl=e(4,:):cb2=c(5,:);%分离出al,a2,bl,b2的收敛情况edl=e(6,:);ed2=e(7,:);ed3=c(8,:);%分离出dl,d2,d3的收敛情况figure(2);i=l：L;plot(i,al,'r',i,a2,'r:',i,a3,'r+',i,bl,*b,i,b2,'b:',i,dl,'g',i,d2,'g:',i,d3,'g+')gridonlegend(,aT,'a2','a3','bl','b2','cl','c2','c3')title(,辨识曲线')figure(3):i=l:L;plot(i,eal,r*,i,ea2,'r:',i,ca3,'r+',i,ebl,b,i,eb2,'b:',i,cdl,'g',i,cd2,'g:',i,ed3,'g+')%画出各个参数的收敛情况gridonlegend(*ea,tea2','ea3','ebl','eb2','ecl','ec2','ec3')title('辨识误差曲线')figure(4);i=l:L:PIot(i,z(i),'g'),gridon;title('辨识系统的输出响应')%画出被辨识系统的输出响应图10输入数据及噪声图12辨识误差曲线图13辨识的输出响应表2增广最小二乘算法的辨识结果参数.4q<真值0,90.150.020.7-1.5LO-EO-0.4估计值0.90.150.020.7-1.5Lo-1.0-0.4从上表中可以看出，对于含有有色噪声的系统，采用增广最小二乘法进行系统辨识可以取得无偏估计。四、模型与误差分析对系统进行辨识后，需要对辨识所得的系统进行进一步的模型验证，以观察所得的模型是否是对当前观测数据的最佳选择。对模型进行验证，可以采用事先生成的输入数据，分别输入到己知模型和经过系统辨识所得的模型中，检验所得模型是否存在问题。如果检验失败，可能存在的问题是：辨识所用的一组数据包含的信息部足或所选的模型类不合适。对本实验中的己知对象及经过系统辨识所得的模型，选取一组输入信号为：M=lt-tl,-tl,l1-l,Ll,-LLl,-l-l,-,l比较原参考对象的输出以及辨识模型的输出。具体代码如下：%w模型检验%w僦%clearallcloseallu=l,-l,1,-1,i,->1,1,-1,1,1,-1,-1,-1,1;%系统输入信号1.=20z=zeros(1,L);x=zeros(lrL);y=zeros(1,L);v=randn(1,L);%实际噪声的产生%ksi(k)=lamda*v(k)-ksi(k-l)-0.41*ksi(k2)%ksi(2)=0;ksi(l)=0;%设ksi的前两个初始值为0Iamda=1;fork=3:L;ksi(k)=lamda*v(k)-ksi(kl)-0.41*ksi(k-2);endforkX:16Z(k)=-0.9*z(k-l)-0.15*z(k-2)-0.02*z(k-3)+0.7*(k-l)-1.5*u(k-2)÷ksi(k)4真实模型的输出y(k)=-L2317*y(k-l)-0.4735*y(k-2)-0.0600*y(k-3)+0.8534*(k-l)-l.2373*(k-2)%递推最小二乘辨识模型输出X(k)=-0.9*x(k-l)-0.15*x(k-2)-0.02*x(k-3)+07*u(k-l)-1.5*u(k-2)÷ksi(k)%增广最小二乘辨识模型的输出endfigured);subplot(3,1,1)stem(u),gridontitleC系统输入信号)subplot(3,1,2)i二1:1:L;plot(i,ZXijy,'b'),gridonlegend(,真实模型输出'，'递推最小二乘辨识模型)titleC真实模型输出和递推最小二乘辨识模型输出值')subplot(3,1,3)i二1:1:L;plot(i,zri,x,'g,),gridonlegend(,真实模型输出'，'增广最小二乘辨识模型')titleC真实模型输出和增广最小二乘辨识模型输出值)figure(2);subplot(2,1,1)i=l;1:L;plot(i,z-y),gridon;titleC真实模型与递推最小二乘辨识模型的差值)subplot(2,1,2)i=1:1:1.PIot(i,z-),gridon;titleC真实模型与增广最小二乘辨识模型的差值)图14不同辨识模型输出和真实系统输出对比图15真实模型与辨识模型输出信号的差值由上图并结合表1和表2,分析可得表1递推最小二乘算法的辨识结果参数.%真值0.90+150,020+7-1,5估计值0.93670.1591-0.01380+6275-1.3558表2增广最小二乘算法的辨识结果参数.4q<真值0.90.150.020.7-1.51.O-1.0-0.1估计值0.90.150.020.7-1.51.O-1.0-0.1辨识给定的含有色噪声的系统时，最小二乘递推辨识算法的参数辨识结果是有偏差的，无法达到稳定值,但是这个误差是在-050.5范围内，还是可以接受的。而与最小二乘递推辨识算法相比，增广最小二乘递推辨识算法考虑了噪声模型，并且在递推到第九步时辨识结果就达到了稳定值，具有辨识速度快、辨识结果精确的特点。因此，对于含有有色噪声的系统，使用增广二乘辨识算法能够取得较好的效果。五、实验总结通过本次仿真实验，使我加深了对系统辨识算法的认识和了解，对于最小二乘参数辨识这种系统辨识领域中最基本最常用的方法有了较好的掌握，并通过MATLAB仿真熟悉了系统辨识的般内容和步骤，对系统辨识的基本过程有了一个较为清晰的认识口最后，王立琦老师这一学期给我们所带的系统辨识这门课程使我受益良多，对我的科研工作有很大的启发借鉴意义。在此十分感谢王老师的指导和帮助。