批處理最小二乘算法

1、理論問題：1)白噪聲信號是什么？ Matlab 中如何產(chǎn)生？2)逆M序列是什么？ Matlab中如何產(chǎn)生？答：( 1)白噪聲( white noise) 系統(tǒng)辨識(shí)中所用到的數(shù)據(jù)通常都含有噪聲， 從工程實(shí)際出發(fā)， 這 種噪聲往往可以視為具有理譜密度的平穩(wěn)隨機(jī)過程。 白噪聲是一種最 簡單的隨機(jī)過程， 是由一系列不相關(guān)的隨機(jī)變量組成的理想化隨機(jī)過 程。白噪聲的數(shù)學(xué)描述如下：如果隨機(jī)過程 (t) 均值為 0、自相關(guān)系 數(shù)為 2 ( ) ，即R( )2 (t)式中， (t) 為單位脈沖函數(shù)，即,0(t) 0, 00,且- (t)d 1則稱該隨機(jī)過程是白噪聲。Matlab 中產(chǎn)生代碼為：xi=randn(

2、L,1);%隨機(jī)產(chǎn)生均值為 0，方差為 1 的高斯隨機(jī)白噪聲序列( 2)逆 M 序列譜分析表明， M 序列含有直流成分， 將造成對辨識(shí)系統(tǒng)的 “凈擾 動(dòng)”，這通常不是所希望的。而逆 M 序列將克服這一缺點(diǎn)，是一種比 M 序列更為理想的偽隨機(jī)碼序列。設(shè)M(k)是周期為Npbit、元素的取值為0或1的M序列，S(k)是 周期為 2bit 、元素依次取值為 0 或 1 的方波序列，將這兩序列按位進(jìn) 行異或運(yùn)算，得到的符合序列就是周期為 2Npbit 、元素取值為 0 或 1 的逆 M 序列，記作 IM(k) ，即有IM(k)M(k)S(k)上述逆 M 序列的邏輯值“ 0”或“1”分別變換為 -1 或

3、 1，此時(shí)逆 M 序列的均值為 0。Matlab中產(chǎn)生代碼為：IM=xor(S,x4);if IM=0u(k)=-1;elseu(k)=1;End代碼中S為元素依次取0或1的方波序列。題目：考慮如下系統(tǒng)：y(k) 1.5y(k 1) 0.7y(k 2) u(k 3) 0.5u(k 4)(k)式中， (k) 為方差為 1 的白噪聲。選用幅值為 1 的逆 M 序列作為輸入信號 u(k) ，利用 LS 算法進(jìn)行 參數(shù)估計(jì)(數(shù)據(jù)長度 L=100)Matlab 代碼如下：>> clear allL=100;d=3;%數(shù)據(jù)長度為 100，純遲延為 3T,uk=zeros(4,1);yk=zer

4、os(2,1);%初始輸入 u(0)=u(-1)=u(-2)=u(-3)=0%初始輸出 y(0)=y(-1)=0x1=1; x2=1; x3=1; x4=0; S=1;xi=randn(L,1);%數(shù)據(jù)長度為 100，方差為 1的白噪聲theta=-1.5;0.7;1.0;0.5;%a1=-1.5,a2=0.7,b0=1,b1=0.5for k=1:Lphi(k,:)=-yk;uk(3:4)'%100X4矩陣phi第k行對應(yīng)的 y(k-1),y(k-2),u(k-3),u(k-4)y(k)=phi(k,:)*theta+xi(k);%10CK1 矩陣 y第k行對應(yīng)的 y(k)IM=xo

5、r(S,x4);%把M序列和方波進(jìn)行異或運(yùn)算得到逆 M序列if IM=0u(k)=-1;elseu(k)=1;endS=not(S); M=xor(x3,x4);%繼續(xù)產(chǎn)生M序列，為制造逆M序列做準(zhǔn)備x4=x3; x3=x2; x2=x1; x1=M;%更新數(shù)據(jù)for i=4:-1:2uk(i)=uk(i-1);enduk(1)=u(k);%更新輸入值 ukfor i=2:-1:2yk(i)=yk(i-1);endyk(1)=y(k);%更新輸出值 ykendthetae=inv(phi'*phi)*phi'*y'%最小二乘估計(jì)公式，計(jì)算參數(shù)估計(jì)值結(jié)果thetae =-

6、1.54890.74791.12200.5626分析如下系統(tǒng)方程為y(k)1.5y(k1) 0.7y(k2) u(k3)0.5u(k4)(k)由CAR模型：y(k)1.5y(k1) 0.7y(k2) u(k3)0.5u(k4)(k)可以得到：A(z 1)11.5z 10.7z21B(z 1)z 3(1 0.5z 1)我們可以知道純遲延為 3T，na=2，nb=1， d=3，phi(k,:)=-yk;uk(3:4)' 100 4矩陣 ph第k行對應(yīng)的 y(k-1),y(k-2),u(k-3), u(k-4) y(k)=phi(k,:)*theta+xi(k);輸出等于矩陣phi與真值矩陣

7、相乘加上白噪聲。這是循環(huán)的主體，接下來一開始我們賦給初值給輸入以及輸出初始輸入 u(0)=u(-1)=u(-2)=u(-3)=0 ;初始輸出 y(0)=y(-1)=0 ；所以y(1) 1.5y(0) 0.7 y( -1) u(-2) 0.5u(-3) (1) 經(jīng)過白噪聲的更新代碼后， 以及輸入輸出值得更新后我們有， 也就是把計(jì)算出來 的y(1)給到下一組y(2)的計(jì)算中區(qū)去：y(2) 1.5y(1) 0.7 y( 0) u(-1) 0.5u(-2) (2)y(99) 1.5 y (98) 0.7y(97) u(96) 0.5u(95) (99) 如此不停的循環(huán)直到滿足循環(huán)次數(shù) L：L次數(shù)越大，也就是我們的估計(jì)結(jié)果會(huì)越接近真值。最后我們獲得矩陣 也就是我們代碼中所說的 phi：0000T(1)0.8404000T(2)