附錄4 遞推最小二乘算法在白噪聲和有色噪聲下工作的源程序

1、nsr=0.1;v=randn(1,1024);v1=nsr*v;b2=1 0.5;a2=1 -0.9 0.95;v2=filter(b2,a2,v1);y1=x+v1;y2=x+v2;附錄4遞推最小二乘算法在白噪聲和有色噪聲下工作的源程序%定義噪信比的取值范圍%產(chǎn)生指定方差的白噪聲信號%在此基礎(chǔ)上生成有色噪聲%分別將白噪聲和有色噪聲疊加在理想輸出上c10=0.001 0.001 0.001 0.001;c20=0.001 0.001 0.001 0.001;p10=10A6*eye(4);p20=10A6*eye(4);E=0.000000005;c1=c10,zeros(4,1023);c

2、2=c20,zeros(4,1023);e1=zeros(4,1024);e2=zeros(4,1024);for k=3:1024h1=-y1(k-1) -y1(k-2) u(k-1) u(k-2);h2=-y2(k-1) -y2(k-2) u(k-1) u(k-2);z1=h1*p10*h1+1; z11=inv(z1);z2=h2*p20*h2+1; z21=inv(z2);k1=p10*h1*z11;k2=p20*h2*z21;d1=y1(k)-h1*c10; c11=c10+k1*d1;d2=y2(k)-h2*c20; c21=c20+k2*d2;p11=p10-k1*k1*h1*p

3、10*h1+1;p21=p20-k2*k2*h2*p20*h2+1;e11=(c11-c10)./c10;e21=(c21-c20)./c20;e1(:,k)=e11;e2(:,k)=e21;c1(:,k)=c11;c2(:,k)=c21;c10=c11;c20=c21;p10=p11;p20=p21;if (e11E)&(e21E) j=k,break;end%準(zhǔn)備求k的值%求出k的值%求出theta的估計值%求出p的值%求出兩次估計值之間的相對誤差%保存相對誤差%保存估計值%將本次估計值作為下次估計的初值enda11=c1(1,:);a12=c1(2,:);b11=c1(3,:);b12=

4、c1(4,:);a21=c2(1,:);a22=c2(2,:);b21=c2(3,:);b22=c2(4,:);ea11=e1(1,:);ea12=e1(2,:);eb11=e1(3,:);eb12=e1(4,:);ea21=e2(1,:);ea22=e2(2,:);eb21=e2(3,:);eb22=e2(4,:);figure% 繪圖i=1:j; plot(i,a11(i),r,i,a12(i),b,i,b11(i),m,i,b12(i),k);title(NSR=0.1時白噪聲影響下各參數(shù)辨識結(jié)果隨遞推次數(shù)的變化)； xlabel(遞推次數(shù));ylabel(各參數(shù)辨識結(jié)果)；axis t

5、ight;grid on;figureplot(i,a21(i),r,i,a22(i),b,i,b21(i),m,i,b22(i),k);title(NSR=0.1時有色噪聲影響下各參數(shù)辨識結(jié)果隨遞推次數(shù)的變化);xlabel(遞推次數(shù));ylabel(各參數(shù)辨識結(jié)果)；axis tight;grid on;figureplot(i,ea11(i),r,i,ea12(i),b,i,eb11(i),m,i,eb12(i),k);title(NSR=0.1時白噪聲影響下各參數(shù)辨識誤差隨遞推次數(shù)的變化)； xlabel(遞推次數(shù));ylabel (各參數(shù)辨識誤差)；axis tight;grid o

6、n;figureplot(i,ea21(i),r,i,ea22(i),b,i,eb21(i),m,i,eb22(i),k);title(NSR=0.1時有色噪聲影響下各參數(shù)辨識誤差隨遞推次數(shù)的變化);xlabel(遞推次數(shù));ylabel (各參數(shù)辨識誤差)；axis tight;grid on;nsr1=0.4;v1=nsr1*v;v2=nsr1/nsr*v2;y1=x+v1;y2=x+v2;%重新定義噪信比的取值范圍%產(chǎn)生指定方差的白噪聲信號%將原有色噪聲擴(kuò)大4倍%分別將白噪聲和有色噪聲疊加在理想輸出上c10=0.001 0.001 0.001 0.001;c20=0.001 0.001

7、0.001 0.001;p10=10A6*eye(4);p20=10A6*eye(4);E=0.000000005;c1=c10,zeros(4,1023);c2=c20,zeros(4,1023);e1=zeros(4,1024);e2=zeros(4,1024);for k=3:1024h1=-y1(k-1) -y1(k-2) u(k-1) u(k-2);h2=-y2(k-1) -y2(k-2) u(k-1) u(k-2);z1=h1*p10*h1+1; z11=inv(z1);z2=h2*p20*h2+1; z21=inv(z2);% 準(zhǔn)備求 k 的值k1=p10*h1*z11;k2=p

8、20*h2*z21;% 求出 k 的值d1=y1(k)-h1*c10; c11=c10+k1*d1;%求出theta的估計值%求出p的值%求出兩次估計值之間的相對誤差%保存相對誤差%保存估計值%將本次估計值作為下次估計的初值d2=y2(k)-h2*c20; c21=c20+k2*d2; p11=p10-k1*k1*h1*p10*h1+1;p21=p20-k2*k2*h2*p20*h2+1;e11=(c11-c10)./c10;e21=(c21-c20)./c20;e1(:,k)=e11;e2(:,k)=e21;c1(:,k)=c11;c2(:,k)=c21;c10=c11;c20=c21;p1

9、0=p11;p20=p21;if (e11E)&(e21E) j1=k,break;endenda11=c1(1,:);a12=c1(2,:);b11=c1(3,:);b12=c1(4,:);a21=c2(1,:);a22=c2(2,:);b21=c2(3,:);b22=c2(4,:);ea11=e1(1,:);ea12=e1(2,:);eb11=e1(3,:);eb12=e1(4,:);ea21=e2(1,:);ea22=e2(2,:);eb21=e2(3,:);eb22=e2(4,:);figure% 繪圖i=1:j1;plot(i,a11(i),r,i,a12(i),b,i,b11(i)

10、,m,i,b12(i),k);title(NSR=0.4時白噪聲影響下各參數(shù)辨識結(jié)果隨遞推次數(shù)的變化)；xlabel(遞推次數(shù));ylabel(各參數(shù)辨識結(jié)果)；axis tight; grid on;figureplot(i,a21(i),r,i,a22(i),b,i,b21(i),m,i,b22(i),k);title(NSR=0.4時有色噪聲影響下各參數(shù)辨識結(jié)果隨遞推次數(shù)的變化);xlabel(遞推次數(shù));ylabel(各參數(shù)辨識結(jié)果)；axis tight; grid on;figureplot(i,ea11(i),r,i,ea12(i),b,i,eb11(i),m,i,eb12(i),k);title(NSR=0.4時白噪聲影響下各參數(shù)辨識誤差隨遞推次數(shù)的變化)；xlabel(遞推次數(shù));ylabel (各參數(shù)辨識誤差)；axis ti