系統(tǒng)辨識(shí)最小二乘參數(shù)估計(jì)matlab

1、最小二乘參數(shù)估計(jì)摘要：最小二乘的一次性完成辨識(shí)算法（也稱批處理算法），他的特點(diǎn)是直接利用已經(jīng)獲得的所有（一批）觀測(cè)數(shù)據(jù)進(jìn)行運(yùn)算處理。這種算法在使用時(shí)，占用內(nèi)存大，離線辨識(shí)，觀測(cè)被辨識(shí)對(duì)象獲得的新數(shù)據(jù)往往是逐次補(bǔ)充到觀測(cè)數(shù)據(jù)集合中去的。在應(yīng)用一次完成算法時(shí)，如果要求在每次新增觀測(cè)數(shù)據(jù)后，接著就估計(jì)出系統(tǒng)模型的參數(shù)，則需要每次新增數(shù)據(jù)后要重新求解矩陣方程。最小二乘辯識(shí)方法在系統(tǒng)辯識(shí)領(lǐng)域中先應(yīng)用上已相當(dāng)普及，方法上相當(dāng)完善，可以有效的用于系統(tǒng)的狀態(tài)估計(jì)，參數(shù)估計(jì)以及自適應(yīng)控制及其他方面。關(guān)鍵詞：最小二乘（Least-squares），系統(tǒng)辨識(shí)（System Identification）目錄： 1.

2、目的22.設(shè)備23引言23.1 課題背景24數(shù)學(xué)模型的結(jié)構(gòu)辨識(shí)35 程序45.1 M序列子函數(shù)45.2主程序56實(shí)驗(yàn)結(jié)果：77參考文獻(xiàn)：71.目的1.1掌握系統(tǒng)辨識(shí)的理論、方法及應(yīng)用1.2熟練Matlab下最小二乘法編程1.3掌握M序列產(chǎn)生方法2.設(shè)備 PC機(jī)1臺(tái)（含Matlab軟件） 3引言3.1 課題背景最小二乘理論是有高斯（，它使各次實(shí)際觀測(cè)值和計(jì)算值之間的差值的平方乘以度量其精度的數(shù)值以后的和最小?！边@就是最小二乘法的最早思想。最小二乘辨識(shí)方法提供一個(gè)估算方法，使之能得到一個(gè)在最小方差意義上與實(shí)驗(yàn)數(shù)據(jù)最好擬合的數(shù)學(xué)模型。遞推最小二乘法是在最小二乘法得到的觀測(cè)數(shù)據(jù)的基礎(chǔ)上，用新引入的數(shù)據(jù)

3、對(duì)上一次估計(jì)的結(jié)果進(jìn)行修正遞推出下一個(gè)參數(shù)估計(jì)值，直到估計(jì)值達(dá)到滿意的精確度為止。4數(shù)學(xué)模型的結(jié)構(gòu)辨識(shí)4.1根據(jù)漢格爾矩陣估計(jì)模型的階次設(shè)一個(gè)可觀可控的SISO過(guò)程的脈沖響應(yīng)序列為個(gè)g(1),g(2),g(L),可以通過(guò)漢格爾（Hankel）矩陣的秩來(lái)確定系統(tǒng)的階次。令Hankel陣為：，其中決定陣地維數(shù)，k可在1至間任意選擇。則有。如果（過(guò)程的真實(shí)階次），那么Hankel陣的秩等于。因此可以利用Hankel陣的奇異性來(lái)確定系統(tǒng)的階次。4.2根據(jù)殘差平方和估計(jì)模型的階次SISO過(guò)程的差分方程模型的輸出殘差為，數(shù)據(jù)長(zhǎng)度L，為階時(shí)的數(shù)據(jù)矩陣，為階時(shí)的參數(shù)的估計(jì)量，為模型階次估計(jì)值，為真實(shí)階次，則殘

4、差平方和函數(shù)：殘差平方和有這樣的性質(zhì)：當(dāng)L足夠大時(shí)，隨著增加先是顯著地下降，當(dāng)時(shí)，值顯著下降的現(xiàn)象就終止。這就是損失函數(shù)法來(lái)定階的原理。圖1-1 損失函數(shù)法原理。 4.3根據(jù)AIC準(zhǔn)則估計(jì)模型的階次.圖1-2 具體的定階用法是：對(duì)不同階次首先用極大似然法估計(jì)參數(shù)，然后計(jì)算似然函數(shù)值及值，找到使的作為。5 程序%待辨識(shí)系統(tǒng) z(k)=0.1*z(k-1)-0.5*z(k-2)+0.6*z(k-3)+u(k-1)+0.5*u(k-2)-0.2*u(k-3)+v(k)/800%clcclear %清理工作間變量L=300; % M序列的周期x1=1;x2=1;x3=1;x4=0;x5=1;x6=0;

5、 %四個(gè)移位積存器的輸出初始值for k=1:L; %開(kāi)始循環(huán)，長(zhǎng)度為L(zhǎng) u(k)=xor(x3,x4); %第一個(gè)移位積存器的輸入是第3個(gè)與第4個(gè)移位積存器的輸出的“或” x6=x5;x5=x4;x4=x3;x3=x2;x2=x1;x1=u(k);end %大循環(huán)結(jié)束，產(chǎn)生輸入信號(hào)u plot(u) %繪圖M序列v=randn(300,1); %隨機(jī)誤差干擾z=zeros(1,300);for k=4:300 z(k)=0.10*z(k-1)+0.55*z(k-2)+0.40*z(k-3)+0.70*u(k-1)+0.90*u(k-2)-3.50*u(k-3)+v(k)/400; %用理想輸

6、出值作為觀測(cè)值endH=zeros(300,6); %定義一個(gè)H“0”矩陣for i=4:300 H(i,:)=-z(i-1) -z(i-2) -z(i-3) u(i-1) u(i-2) u(i-3);%用循環(huán)產(chǎn)生H矩陣 z1(i,:)=z(i); %用循環(huán)產(chǎn)生z矩陣end%計(jì)算參數(shù)%c=inv(H*H)*H*z1%帶入公式書(shū)上a1=c(1),a2=c(2),a3=c(3),b1=c(4),b2=c(5),b3=c(6)%辨識(shí)出參數(shù)%系統(tǒng)階次辨識(shí)AIC算法%bb=zeros(5,1);n=1; %假設(shè)為1階for i=2:300 H1(i,:)=-z(i-1) u(i-1); zz1(i,:)

7、=z(i);endaa1=inv(H1*H1)*H1*zz1bb(1)=(zz1-H1*aa1)*(zz1-H1*aa1)/L;AIC(1)=L*log(bb(1)+4*n;n=2; %假設(shè)為2階for i=3:300 H2(i,:)=-z(i-1) -z(i-2) u(i-1) u(i-2); zz2(i,:)=z(i);endaa2=inv(H2*H2)*H2*zz2bb(2)=(zz2-H2*aa2)*(zz2-H2*aa2)/L;AIC(2)=L*log(bb(2)+4*n;n=3; %假設(shè)為3階for i=4:300 H3(i,:)=-z(i-1) -z(i-2) -z(i-3) u

8、(i-1) u(i-2) u(i-3); zz3(i,:)=z(i);endaa3=inv(H3*H3)*H3*zz3bb(3)=(zz3-H3*aa3)*(zz3-H3*aa3)/L;AIC(3)=L*log(bb(3)+4*n;n=4; %假設(shè)為4階for i=5:300 H4(i,:)=-z(i-1) -z(i-2) -z(i-3) -z(i-4) u(i-1) u(i-2) u(i-3) u(i-4); zz4(i,:)=z(i);endaa4=inv(H4*H4)*H4*zz4bb(4)=(zz4-H4*aa4)*(zz4-H4*aa4)/L;AIC(4)=L*log(bb(4)+4*n;n=5; %假設(shè)為5階for i=6:300 H5(i,:)=-z(i-1) -z(i-2) -z(i-3) -z(i-4) -z(i-5) u(i-1) u(i-2) u(i-3) u(i-4) u(i-5); zz5(i,:)=z(i);endaa5=inv(H5*H5)*H5*zz5bb(5)=(zz5-H5*aa5)*(zz5-H5*aa5)/L;AIC(5)=L*log(bb(5)+4*n;x=min(AIC)for i=1:5 if(AIC(i)=x) N=i %所辨識(shí)出的階次N endendplot(1:5,AIC(1) AIC(2) AIC(3)