最小二乘法在系统辨识中的应用(包含相关的三种算法)
合集下载
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
)]1()()()[() 1()(kkkzkkkfffhK) 1()]()([)(kkkkffffPhKIP1] 1)() 1()()[() 1()(kkkkkkffffffhPhhPK)]1()()( e)[() 1()(kkkkkkeeeeehK) 1()]()([)(kkkkeeeePhKIP1] 1)() 1()()[() 1()(kkkkkkeeeeeehPhhPK利用上述公式即可求得参数的估计值并能求出噪声模型的估计。
但是，数据向量)(kfh中的变量均需要按照（3.2.2）式计算，然而噪声模型)(1zC并不知道。
为此需要用迭代的方法来估计)(1zC。
令)()(1z)(1kvCke（3.2.5）置,2[k],)()]3(),2(),1([)(31ccckekekekeபைடு நூலகம்h（3.2.6）就把噪声模型（3.2.5）也化成最小二乘格式)()()(kvkkeeeh由于上式的噪声已是白噪声，所以再次利用最小二乘法可获得噪声模型参数e的无偏估计。
通过极小化（1.1.4）式来计算的方法称作最小二乘法，未知模型参数最可能的值是在实际观测值与计算值之累次误差的平方和达到最小处所得到的，这种模型输出能最好地接近实际过程的输出。
2、辨识原理考虑模型（1.1.2）式的辨识问题，其中)(kz和)(kh都是可观测的数据，是待估计参数，准则函数取（1.1.4）根据（1.1.3）的定义，准则函数)(J可写成二次型的形式)()()(HzHzllllJ（1.2.1）显然上式中的Hl代表模型的输出，或者说是过程的输出预报值。
)()() 1()()()()()(1111kkkkkiikkihhPhhhhP（2.2.3）令：
] ) 1(),2 (z),1 (z[1kzkz则：
] )i()()[1()() 1(1111111ikkkkkzikkhPzHHH于是有i111)()() 1() 1(kizikkhP（2.2.4）令] )k(),2 (z),1 (z[zkz利用（2.2.3）和（2.2.4）式，可得)]1()()()[()() 1()}()() 1()]()()(){[()]()() 1() 1()[(] )i()()[()()(1111ikkkzkkkkzkkkkkkkzkkkkzikkkkkkkhhphhhPPhPPhPzHHH（2.2.5）引进增益矩阵)(kK，定义为)()()(kkkhPK（2.2.6）则（2.2.5）式写成)]1()()()[() 1()(kkkzkkkhK（2.2.7）进一步把（2.2.3）式写成11)]()() 1([)(kkkkhhPP（2.2.8）为了避免矩阵求逆运算，利用矩阵反演公式可将（2.2.8）式演变成) 1()]()([)]()() 1([)(11kkkkkkkPhKIhhPP（2.2.9）将（2.2.9）式代入（2.2.6）式，整理后有1] 1)() 1()()[() 1()(kkkkkkhPhhPK（2.2.10）综合（2.2.7）、（2.2.9）、（2.2.10）式便得到最小二乘参数估计递推算法。
但是，数据向量)(ke h依然包含着不可测得噪声量) 3(),2(),1(kekeke，它可用相应的估计值代替，置)]3( e),2( e),1( e[)(kkkke h（3.2.7）其中0, 0)( ekk，当0k时，按照式子)()()()( ekkhkzk计算，式子中)]3(),2(),1(),3(),2(),1([)(kukukukzkzkzkh综上分析，广义最小二乘递推算法可归纳成：
因此)(J可以看作来衡量模型输出与实际过程输出的接近情况。
极小化)(J，求得参数的估计值llllzHHH1)(（1.2.2）将使模型的输出最好的预报过程的输出。
3、辨识结果]0.00000.00000.00000.00000.00001.0000- [][321321bbbaaa二、递推最小二乘法（RLS）1、数学模型在第一部分中建立了最小二乘法，并用一次完成算法进行了计算。
cbaanncnbnnnzczczczCzbzbzbzBzazazazA2211122111221111)()(1)(在本题中，在本题中，a n =bn =c n =3.即广义最小二乘发的基本思想是基于对数据先进行一次滤波预处理，然后利用普通最小二乘法对滤波后的数据进行辨识。
21bbb)33211332211133221111)((1)(zczczczCzzzzBzazazazA（3.1.2）2、辨识原理由（3.1.1）式可得)()()()()()()(1111kvkuzCzBkzzCzA（3.2.1）令)()()()()()(11kuzCkukzzCkzff，（3.2.2）及,,],,,[)]3(),2(),1(),3(),2(),1([)(321321bbbaaakukukukzkzkzkfffffffh（3.2.3）可将模型（3.1.1）式化成最小二乘格式)()()(kvkkzffh（3.2.4）由于上式)(kv是白噪声，所以利用最小二乘法即可获得参数的无偏估计。
)]1()()()[() 1()(kkkzkkkhK) 1()]()([)]()() 1([)(11kkkkkkkPhKIhhPP1] 1)() 1()()[() 1()(kkkkkkhPhhPK利用上述公式即可求得参数的估计值3、辨识结果]0.0001-0.0005-0.00000.0002-0.0007-1.0005- [][321321bbbaaa三、广义最小二乘法（GLS）1、数学模型设时不变SISO动态过程的数学模型为)()(1)()()()(111kvzCkuzBkzzA（3.1.1）其中，)(ku和)(kz分别表示过程的输入和输出，)(kv是均值为零的不相关随机噪声，多项式)(1zA、)(1zB和)(1zC为：
bbb1nbnnnzzzzBzaz2az1azAbaa2122111)(1)(在本题中，a n =bn =3.即bbb1)3322113322111(1)(zzzzBzazazazA将此模型写成最小二乘格式)()()(knkkzh（1.1.2）其中，)(kz是过程的输出量；)(kh是可观测的数据向量；)(kn是均值为零的随机噪声。
但是由于具体使用时占用内存量大，而且不能用于在线辨识。
所以引入了最小二乘参数估计的递推算法。
递推算法的基本思想可以概括成：
新的估计值)(k=老的估计值) 1(k+修正项（2.1.1）在此算法中，时不变SISO动态过程的数学模型仍与最小二乘法的一样。
模型的最小二乘格式也相同，只是计算方法不同，具体计算方法，在递推最小二乘法的辨识原理中描述。
式中k,],,,,[)]3(),2(),1(),3(),2(),1([)(321321bbbaaakukukukzkzkzkh对于L, 2 , 1，方程式（1.1.2）构成一个线性方程组，可以把它写成) 3() 2() 1() 3() 2() 1(0) 0 (u) 1 (u0) 0 (z) 1 (z00) 0 (u00) 0 (z)() 2 () 1 ()](, ),2 (n),1 (n[)](,),2 (z),1 (z[lulululzlzlzllnnlzzlllhhhh（1.1.3）利用数据序列)}({kz和)}({kh，极小化准则函数kLkkzJ12])()([)(h（1.1.4）使min)(J的估计值记作，称为参数的最小二乘估计值。
最小二乘法在系统辨识中的应用(包含相关的三种算法)
2008级硕士研究生《系统建模理论》试卷已知一个三阶线性离散系统的输入、输出数据，共有40个采样值，试分别用：
最小二乘法（LS）、递推最小二乘法（RLS）、广义最小二乘法（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.2705 k 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.0786 k 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.1568 k 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动态过程的数学模型为)()()()()(11knkuzBkzzA（1.1.1）其中，)(ku为过程的输入量，)(kz为过程的输出量，)(kn是噪声，多项式)(1zA和)(1zB为：
2、辨识原理首先将1.2.2式一次完成算法写成] )i()([] )i()([)()(1111iillllllllziizlzhhhHPHHH（2.2.1）定义ik1111111)()() 1()()()(kkkikkiikPiikPHHhhHHhh（2.2.2）其中：
)() 2 () 1 (kkhhhH，) 1() 2 () 1 (1kkhhhH由（2.2.2）式可得：
通过极小化114式来计算的方法称作最小二乘法知模型参数最可能的值是在实际观测值与计算值之累次误差的平方和达到最小处所得到的这种模型输出能最好地接近实际过程的输辨识原理考虑模型112式的辨识问题其中kz和kh都是可观测的数据是待估计参数准则函数取114根据113的定义准则函数j可写成二次型的形121显然上式中的hl代表模型的输出或者说是过程的输出预报值
3、估计结果123123[][-1.0005-0.00080.0002-0.0000-0.0005-0.0001]aaabbb噪声模型的估计值为1.0e-003 *[ -0.0000 0.4115 0.0001 ]即1(C z)1附录：
MATLAB程序1、一次完成的最小二乘法%%一次完成的最小二乘法clear s=(1: 40) ; 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] ; 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. 0027 2. 1145 1. 1157]; h=zeros(40, 6) ; h(1, : ) =[-z(1) 0 0 u(1) 0 0] ; h(2, : ) =[-z(2) -z(1) 0 u(2) u(1) 0]; for i=3: 1:40 h(i, : ) =[-z(i) -z(i-1) -z(i-2) u(i) u(i-1) u(i-2) ] ; end o=inv(h’ *h) *h’ *z’ %误差分析J=0; for i=1: 1:40 e(i) =[z(i) -h(i, :) *o] ; J=J+e(i) _; end J plot(s, e) ; 2、递推的最小二乘法%%递推的最小二乘法clear s=(1: 40) ; 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] ; 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. 0027 2. 1145 1. 1157]; h=zeros(40, 6) ; h(1, : ) =[-z(1) 0 0 u(1) 0 0] ; h(2, : ) =[-z(2) -z(1) 0 u(2) u(1) 0]; for i=3: 1:40 h(i, : ) =[-z(i) -z(i-1) -z(i-2) u(i) u(i-1) u(i-2) ] ; end p0=10_; o0=0. 001; k(1, : ) =[p0*h(1, :) ‘ *inv(h(1, :) *p0*h(1, : ) ‘ +1) ] ‘ ; K=zeros(6, 40) ; K(: , 1) =k(1, : ) ‘ ; p=zeros(6, 6, 40) ; p(: , : , 1) =[eye(6) -K(: , 1) *h(1, : ) ]*p0; o=zeros(6, 40) ‘ ; o(1, : ) =o0+K(:, 1) *[z(1) -h(1) *o0]; for i=2: 1:40 k(i, : ) =[p(: , : , i-1) *h(i, :) ‘ *inv(h(i, :) *p(: , : , i-1) *h(i, :) ‘ +1) ]’ ; p(:, : , i) =[eye(6) -k(i, : ) ‘ *h(i, :) ] *p(:, : , i-1) ; o(i, : ) =[o(i-1, :) ‘ +k(i, : ) ‘ *[z(i) -h(i, : ) *o(i-1, : ) ‘ ] ] ‘ ; end %辨识结果o(40, : ) %误差分析J=0; for i=2: 1:40 e(i) =[z(i) -h(i, :) *o(i-1, : ) ‘ ]; J=J+e(i) _; end plot(s, e) ; J 3、广义最小二乘法%%广义预测的最小二乘法clear s=(1: 40) ; 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] ; 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.
但是，数据向量)(kfh中的变量均需要按照（3.2.2）式计算，然而噪声模型)(1zC并不知道。
为此需要用迭代的方法来估计)(1zC。
令)()(1z)(1kvCke（3.2.5）置,2[k],)()]3(),2(),1([)(31ccckekekekeபைடு நூலகம்h（3.2.6）就把噪声模型（3.2.5）也化成最小二乘格式)()()(kvkkeeeh由于上式的噪声已是白噪声，所以再次利用最小二乘法可获得噪声模型参数e的无偏估计。
通过极小化（1.1.4）式来计算的方法称作最小二乘法，未知模型参数最可能的值是在实际观测值与计算值之累次误差的平方和达到最小处所得到的，这种模型输出能最好地接近实际过程的输出。
2、辨识原理考虑模型（1.1.2）式的辨识问题，其中)(kz和)(kh都是可观测的数据，是待估计参数，准则函数取（1.1.4）根据（1.1.3）的定义，准则函数)(J可写成二次型的形式)()()(HzHzllllJ（1.2.1）显然上式中的Hl代表模型的输出，或者说是过程的输出预报值。
)()() 1()()()()()(1111kkkkkiikkihhPhhhhP（2.2.3）令：
] ) 1(),2 (z),1 (z[1kzkz则：
] )i()()[1()() 1(1111111ikkkkkzikkhPzHHH于是有i111)()() 1() 1(kizikkhP（2.2.4）令] )k(),2 (z),1 (z[zkz利用（2.2.3）和（2.2.4）式，可得)]1()()()[()() 1()}()() 1()]()()(){[()]()() 1() 1()[(] )i()()[()()(1111ikkkzkkkkzkkkkkkkzkkkkzikkkkkkkhhphhhPPhPPhPzHHH（2.2.5）引进增益矩阵)(kK，定义为)()()(kkkhPK（2.2.6）则（2.2.5）式写成)]1()()()[() 1()(kkkzkkkhK（2.2.7）进一步把（2.2.3）式写成11)]()() 1([)(kkkkhhPP（2.2.8）为了避免矩阵求逆运算，利用矩阵反演公式可将（2.2.8）式演变成) 1()]()([)]()() 1([)(11kkkkkkkPhKIhhPP（2.2.9）将（2.2.9）式代入（2.2.6）式，整理后有1] 1)() 1()()[() 1()(kkkkkkhPhhPK（2.2.10）综合（2.2.7）、（2.2.9）、（2.2.10）式便得到最小二乘参数估计递推算法。
但是，数据向量)(ke h依然包含着不可测得噪声量) 3(),2(),1(kekeke，它可用相应的估计值代替，置)]3( e),2( e),1( e[)(kkkke h（3.2.7）其中0, 0)( ekk，当0k时，按照式子)()()()( ekkhkzk计算，式子中)]3(),2(),1(),3(),2(),1([)(kukukukzkzkzkh综上分析，广义最小二乘递推算法可归纳成：
因此)(J可以看作来衡量模型输出与实际过程输出的接近情况。
极小化)(J，求得参数的估计值llllzHHH1)(（1.2.2）将使模型的输出最好的预报过程的输出。
3、辨识结果]0.00000.00000.00000.00000.00001.0000- [][321321bbbaaa二、递推最小二乘法（RLS）1、数学模型在第一部分中建立了最小二乘法，并用一次完成算法进行了计算。
cbaanncnbnnnzczczczCzbzbzbzBzazazazA2211122111221111)()(1)(在本题中，在本题中，a n =bn =c n =3.即广义最小二乘发的基本思想是基于对数据先进行一次滤波预处理，然后利用普通最小二乘法对滤波后的数据进行辨识。
21bbb)33211332211133221111)((1)(zczczczCzzzzBzazazazA（3.1.2）2、辨识原理由（3.1.1）式可得)()()()()()()(1111kvkuzCzBkzzCzA（3.2.1）令)()()()()()(11kuzCkukzzCkzff，（3.2.2）及,,],,,[)]3(),2(),1(),3(),2(),1([)(321321bbbaaakukukukzkzkzkfffffffh（3.2.3）可将模型（3.1.1）式化成最小二乘格式)()()(kvkkzffh（3.2.4）由于上式)(kv是白噪声，所以利用最小二乘法即可获得参数的无偏估计。
)]1()()()[() 1()(kkkzkkkhK) 1()]()([)]()() 1([)(11kkkkkkkPhKIhhPP1] 1)() 1()()[() 1()(kkkkkkhPhhPK利用上述公式即可求得参数的估计值3、辨识结果]0.0001-0.0005-0.00000.0002-0.0007-1.0005- [][321321bbbaaa三、广义最小二乘法（GLS）1、数学模型设时不变SISO动态过程的数学模型为)()(1)()()()(111kvzCkuzBkzzA（3.1.1）其中，)(ku和)(kz分别表示过程的输入和输出，)(kv是均值为零的不相关随机噪声，多项式)(1zA、)(1zB和)(1zC为：
bbb1nbnnnzzzzBzaz2az1azAbaa2122111)(1)(在本题中，a n =bn =3.即bbb1)3322113322111(1)(zzzzBzazazazA将此模型写成最小二乘格式)()()(knkkzh（1.1.2）其中，)(kz是过程的输出量；)(kh是可观测的数据向量；)(kn是均值为零的随机噪声。
但是由于具体使用时占用内存量大，而且不能用于在线辨识。
所以引入了最小二乘参数估计的递推算法。
递推算法的基本思想可以概括成：
新的估计值)(k=老的估计值) 1(k+修正项（2.1.1）在此算法中，时不变SISO动态过程的数学模型仍与最小二乘法的一样。
模型的最小二乘格式也相同，只是计算方法不同，具体计算方法，在递推最小二乘法的辨识原理中描述。
式中k,],,,,[)]3(),2(),1(),3(),2(),1([)(321321bbbaaakukukukzkzkzkh对于L, 2 , 1，方程式（1.1.2）构成一个线性方程组，可以把它写成) 3() 2() 1() 3() 2() 1(0) 0 (u) 1 (u0) 0 (z) 1 (z00) 0 (u00) 0 (z)() 2 () 1 ()](, ),2 (n),1 (n[)](,),2 (z),1 (z[lulululzlzlzllnnlzzlllhhhh（1.1.3）利用数据序列)}({kz和)}({kh，极小化准则函数kLkkzJ12])()([)(h（1.1.4）使min)(J的估计值记作，称为参数的最小二乘估计值。
最小二乘法在系统辨识中的应用(包含相关的三种算法)
2008级硕士研究生《系统建模理论》试卷已知一个三阶线性离散系统的输入、输出数据，共有40个采样值，试分别用：
最小二乘法（LS）、递推最小二乘法（RLS）、广义最小二乘法（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.2705 k 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.0786 k 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.1568 k 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动态过程的数学模型为)()()()()(11knkuzBkzzA（1.1.1）其中，)(ku为过程的输入量，)(kz为过程的输出量，)(kn是噪声，多项式)(1zA和)(1zB为：
2、辨识原理首先将1.2.2式一次完成算法写成] )i()([] )i()([)()(1111iillllllllziizlzhhhHPHHH（2.2.1）定义ik1111111)()() 1()()()(kkkikkiikPiikPHHhhHHhh（2.2.2）其中：
)() 2 () 1 (kkhhhH，) 1() 2 () 1 (1kkhhhH由（2.2.2）式可得：
通过极小化114式来计算的方法称作最小二乘法知模型参数最可能的值是在实际观测值与计算值之累次误差的平方和达到最小处所得到的这种模型输出能最好地接近实际过程的输辨识原理考虑模型112式的辨识问题其中kz和kh都是可观测的数据是待估计参数准则函数取114根据113的定义准则函数j可写成二次型的形121显然上式中的hl代表模型的输出或者说是过程的输出预报值
3、估计结果123123[][-1.0005-0.00080.0002-0.0000-0.0005-0.0001]aaabbb噪声模型的估计值为1.0e-003 *[ -0.0000 0.4115 0.0001 ]即1(C z)1附录：
MATLAB程序1、一次完成的最小二乘法%%一次完成的最小二乘法clear s=(1: 40) ; 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] ; 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. 0027 2. 1145 1. 1157]; h=zeros(40, 6) ; h(1, : ) =[-z(1) 0 0 u(1) 0 0] ; h(2, : ) =[-z(2) -z(1) 0 u(2) u(1) 0]; for i=3: 1:40 h(i, : ) =[-z(i) -z(i-1) -z(i-2) u(i) u(i-1) u(i-2) ] ; end o=inv(h’ *h) *h’ *z’ %误差分析J=0; for i=1: 1:40 e(i) =[z(i) -h(i, :) *o] ; J=J+e(i) _; end J plot(s, e) ; 2、递推的最小二乘法%%递推的最小二乘法clear s=(1: 40) ; 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] ; 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. 0027 2. 1145 1. 1157]; h=zeros(40, 6) ; h(1, : ) =[-z(1) 0 0 u(1) 0 0] ; h(2, : ) =[-z(2) -z(1) 0 u(2) u(1) 0]; for i=3: 1:40 h(i, : ) =[-z(i) -z(i-1) -z(i-2) u(i) u(i-1) u(i-2) ] ; end p0=10_; o0=0. 001; k(1, : ) =[p0*h(1, :) ‘ *inv(h(1, :) *p0*h(1, : ) ‘ +1) ] ‘ ; K=zeros(6, 40) ; K(: , 1) =k(1, : ) ‘ ; p=zeros(6, 6, 40) ; p(: , : , 1) =[eye(6) -K(: , 1) *h(1, : ) ]*p0; o=zeros(6, 40) ‘ ; o(1, : ) =o0+K(:, 1) *[z(1) -h(1) *o0]; for i=2: 1:40 k(i, : ) =[p(: , : , i-1) *h(i, :) ‘ *inv(h(i, :) *p(: , : , i-1) *h(i, :) ‘ +1) ]’ ; p(:, : , i) =[eye(6) -k(i, : ) ‘ *h(i, :) ] *p(:, : , i-1) ; o(i, : ) =[o(i-1, :) ‘ +k(i, : ) ‘ *[z(i) -h(i, : ) *o(i-1, : ) ‘ ] ] ‘ ; end %辨识结果o(40, : ) %误差分析J=0; for i=2: 1:40 e(i) =[z(i) -h(i, :) *o(i-1, : ) ‘ ]; J=J+e(i) _; end plot(s, e) ; J 3、广义最小二乘法%%广义预测的最小二乘法clear s=(1: 40) ; 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] ; 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.