《系统建模理论与自适应控制》试卷一、已知三阶线性离散系统的输入、输出数据，共有40个采样值，分别用最小二乘法（LS ）的一次完成算法和递推算法、广义最小二乘法（GLS ）进行参数估计，并阐述各种算法的辨识原理，给出各种算法的程序流程图及程序注释。

k 1 2 3 4 5 6 7 8 9 10 u (k ) 0.8251 0.0988 0.4628 -0.9168 2.2325 0.0777 2.3654 0.3476 1.1473 -1.9035 y (k ) 1.5333 -1.0680 1.0666 -0.5284 -0.5835 3.1471 -3.7185 6.2149 -6.3026 7.2705k 11 12 13 14 15 16 17 18 19 20 u (k ) -0.9229 1.6400 -0.8410 0.7599 -0.4739 -0.1784 -1.7760 -1.6722 1.2959 -0.0591 y (k ) -9.0552 8.1735 -5.9004 3.9870 -2.2486 0.9525 -0.5325 -1.5227 0.4200 1.0786k 21 22 23 24 25 26 27 28 29 30 u (k ) -1.0576 -1.0071 1.1342 -0.0740 0.6759 0.5221 0.9954 0.5271 -1.7656 0.4936 y (k ) -1.5579 0.6640 -1.4222 2.6444 -2.9572 3.6340 -3.1281 3.8334 -3.2542 1.1568k 31 32 33 34 35 36 37 38 39 40 u (k ) 1.4810 0.9591 -3.1293 -0.3604 -0.4251 0.4185 -0.6728 -0.0027 2.1145 1.1157 y (k ) 0.0615 0.9120 -0.0692 -3.2731 3.7486 -4.3194 4.7230 -5.2781 5.1507 -2.7235一、最小二乘法（LS ） 1、数学模型设时不变SISO 动态过程的数学模型为)()()()()(11k n k u z B k z z A +=-- （1.1.1）其中，)(k u 为过程的输入量，)(k z 为过程的输出量，)(k n 是噪声，多项式)(1-z A 和)(1-z B 为：⎪⎩⎪⎨⎧+++=++++=--------nbn n n z b z b z b z B z a z a z a z A b aa 2211122111)(1)( 在本题中，a n =b n =3.即⎩⎨⎧++=+++=--------33221113322111)(1)(z b z b z b z B z a z a z a z A 将此模型写成最小二乘格式)()()(k n k k z +=θh τ （1.1.2） 其中，是过程的输出量；)(k τh 是可观测的数据向量；是均值为零的随机噪声。

式中⎩⎨⎧=---------=ττ],,,,,[)]3(),2(),1(),3(),2(),1([)(321321b b b a a a k u k u k u k z k z k z k θh 对于L k ,2,1=，方程式（1.1.2）构成一个线性方程组，可以把它写成⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡------------=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡===)3()2()1()3()2()1(0)0()1(0)0()1(00)0(00)0()()2()1()](,),2(),1([)](,),2(),1([l u l u l u l z l z l z u u z z u z l l n n n n l z z z z l l l τττττh h h h （1.1.3） 利用数据序列({k z 和)}({k h ，极小化准则函数∑=-=Lk k k z J 12])()([)(θh θτ （1.1.4）使min )(=θJ 的θ估计值记作θˆ，称为参数θ的最小二乘估计值。

通过极小化（1.1.4）式来计算θˆ的方法称作最小二乘法， 未知模型参数θ最可能的值是在实际观测值与计算值之累次误差的平方和达到最小处所得到的，这种模型输出能最好地接近实际过程的输出。

2、辨识原理考虑模型（1.1.2）式的辨识问题，其中(k z 和)(k h 都是可观测的数据，θ是待估计参数，准则函数取（1.1.4）根据（1.1.3）的定义，准则函数)(θJ 可写成二次型的形式)()()(θH z H z θl l l l J --=τθ （1.2.1） 显然上式中的θH l 代表模型的输出，或者说是过程的输出预报值。

因此)(θJ 可以看作来衡量模型输出与实际过程输出的接近情况。

极小化)(θJ ，求得参数θ的估计值l l l l z ττH H H θ1)(ˆ-= （1.2.2） 将使模型的输出最好的预报过程的输出。

3、程序流程图图1 一次完成算法流程图4、程序及辨识结果程序： clc clear% 给定系统输入输出数据u=[ 0.8251 0.0988 0.4628 -0.9168 2.2325 0.0777 2.3654 0.3476 1.1473 -1.9035 ...-0.9229 1.6400 -0.8410 0.7599 -0.4739 -0.1784 -1.7760 -1.6722 1.2959 -0.0591 ...-1.0576 -1.0071 1.1342 -0.0740 0.6759 0.5221 0.9954 0.5271 -1.7656 0.4936 ...1.4810 0.9591 -3.1293 -0.3604 -0.4251 0.4185 -0.6728 -0.00272.1145 1.1157 ];开始为输入信号赋初值定义输出观测值的长度和系统的输出值给样本矩阵L H 和L z 赋值计算参数∧θ从∧θ中分离并显示出被辨识参数123123,,,,,a a a b b b结束z=[ 1.5333 -1.0680 1.0666 -0.5284 -0.5835 3.1471 -3.7185 6.2149 -6.3026 7.2705...-9.0552 8.1735 -5.9004 3.9870 -2.2486 0.9525 -0.5325 -1.5227 0.4200 1.0786 ...-1.5579 0.6640 -1.4222 2.6444 -2.9572 3.6340 -3.1281 3.8334 -3.2542 1.1568 ...0.0615 0.9120 -0.0692 -3.2731 3.7486 -4.3194 4.7230 -5.2781 5.1507 -2.7235 ];%对可观测矩阵H 赋值 H=zeros(40,6); H(1,:)=[0 0 0 0 0 0];H(2,:)=[-z(1) 0 0 u(1) 0 0];H(3,:)=[-z(2) -z(1) 0 u(2) u(1) 0]; for i=4:1:40H(i,:)=[-z(i-1) -z(i-2) -z(i-3) u(i-1) u(i-2) u(i-3)]; end%计算参数CS=inv(H'*H)*H'*z' %误差分析计算 J=0;for k=1:1:40e(k)=z(k)-H(k,:)*CS; J=J+e(k)^2; end J辨识结果：]0.00000.00000.00000.00000.00001.0000-[][ˆ321321==b b b a a a θ二、递推最小二乘法（RLS ） 1、数学模型在第一部分中建立了最小二乘法，并用一次完成算法进行了计算。

但是由于具体使用时占用内存量大，而且不能用于在线辨识。

所以引入了最小二乘参数估计的递推算法。

递推算法的基本思想可以概括成：新的估计值)(ˆk θ=老的估计值)1(ˆ-k θ+修正项 （2.1.1） 在此算法中，时不变SISO 动态过程的数学模型仍与最小二乘法的一样。

模型的最小二乘格式也相同，只是计算方法不同，具体计算方法，在递推最小二乘法的辨识原理中描述。

2、辨识原理首先将1.2.2式一次完成算法写成])()([])()([)(ˆ)(ˆ1111∑∑=-=-===li li ll l l l l i z i i i z l z h h h H P H H H θττττ （2.2.1）定义⎪⎪⎩⎪⎪⎨⎧==-==∑∑-=---=-1111111ˆ)()()1(ˆ)()()(k i k k ki k k i i k P i i k P H H h h H H h h ττττ （2.2.2） 其中：⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡=)()2()1(k k τττh h h H ，⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡-=-)1()2()1(1k k τττh h h H 由（2.2.2）式可得：)()()1()()()()()(1111k k k k k i i k k i τττh h P h h h h P +-=+=--=-∑ （2.2.3）令：τ])1(),2(),1([1-=-k z z z k z则：])()()[1()()1(ˆ1111111∑-=------==-k i k k k k i z i k k h P z H H H θττ 于是有∑-=-=--111)()()1(ˆ)1(k i i z i k k h θP （2.2.4） 令τ])(),2(),1([k z z z k =z 利用（2.2.3）和（2.2.4）式，可得)]1(ˆ)()()[()()1(ˆ)}()()1(ˆ)]()()(){[()]()()1(ˆ)1()[(])()()[()()(ˆ1111--+-=+--=+--===--=-∑k k k z k k k k z k k k k k k k z k k k k i z i k k ki kk k k θh h p θh θh h P P h θP P h P z H H H θττττ （2.2.5）引进增益矩阵)(k K ，定义为)()()(k k k h P K = （2.2.6）则（2.2.5）式写成)]1(ˆ)()()[()1(ˆ)(ˆ--+-=k k k z k k k θh K θθτ （2.2.7） 进一步把（2.2.3）式写成11)]()()1([)(--+-=k k k k τh h P P （2.2.8） 为了避免矩阵求逆运算，利用矩阵反演公式可将（2.2.8）式演变成)1()]()([)]()()1([)(11--=+-=--k k k k k k k τP h K I h h P P τ （2.2.9） 将（2.2.9）式代入（2.2.6）式，整理后有1]1)()1()()[()1()(-+--=k k k k k k h P h h P K τ （2.2.10） 综合（2.2.7）、（2.2.9）、（2.2.10）式便得到最小二乘参数估计递推算法。