H5(i,8)=u(i+2); H5(i,9)=u(i+1); H5(i,10)=u(i); end
estimate5=inv(H5'*H5)*H5'*z(6:L+5)';
D5=(z(6:L+5)'-H5*estimate5)'*(z(6:L+5)'-H5*estimate5)/L; %噪声方差的估计值 AIC5=L*log(D5)+4*5;
plot(1:5,[AIC1 AIC2 AIC3 AIC4 AIC5]) grid on
title('损失函数法判断模型阶次') xlabel('阶次') ylabel('损失函数')
程序运行结果如图6所示:
图6 模型阶次的判定
由图6可知,当系统阶次从1到5变化时,在n=3及以后,损失函数只有微小的下降,在n=4以后又出现上升的趋势,所以取n=3作为模型的阶次。
2. 参数辨识
(1)递推最小二乘法
根据前文所述递推最小二乘法算法流程,编写递推最小二乘法程序如下:
% 递推最小二乘辨识 clc clear all
%采用M序列方法产生输入数据 L=90;?个数据 namenda=1;
y1=1;y2=1;y3=1;y4=0; u_f=2; for i=1:L
x1=xor(y3,y4); x2=y1; x3=y2; x4=y3; y(i)=y4; if y(i)>0.5 u(i)=-u_f; else
u(i)=u_f; end
y1=x1;y2=x2;y3=x3;y4=x4; end
v=1*randn(L,1); eb(2)=0;eb(1)=0; for k=3:L
eb(k)=-eb(k-1)-0.41*eb(k-2)+namenda*v(k); end
z(3)=0;z(2)=0;z(1)=0; for k=4:L;
z(k)=-0.9*z(k-1)-0.15*z(k-2)-0.02*z(k-3)+0.7*u(k-1)-1.5*u(k-2)+eb(k); end figure(1) i=1:30; subplot(2,1,1) stem(u(i));
title('输入M序列') grid on subplot(2,1,2) stem(z(i)); title('输出序列') grid on; %
c0=[0.001 0.001 0.001 0.001 0.001]'; p0=10^3*eye(5,5); E= 0.0000005; c=[c0,zeros(5,L-1)]; e=zeros(5,L); %
for k=4:L
h1=[-z(k-1),-z(k-2),-z(k-3),u(k-1),u(k-2)]'; k1=p0*h1*inv(h1'*p0*h1+1); c1=c0+k1*(z(k)-h1'*c0); e1=c1-c0; e2=e1./c0; e(:,k)=e2; c0=c1; c(:,k)=c1;
p1(:,:,k)=p0-k1*k1'*[h1'*p0*h1+1]; p0=p1(:,:,k); if abs(e2)<=E N=k; break; end N=k; end
% 递推最小二乘仿真图 a1_f=c(1,N) a2_f=c(2,N) a3_f=c(3,N) b1_f=c(4,N) b2_f=c(5,N)
a1=c(1,1:N); a2=c(2,1:N);a3=c(3,1:N);b1=c(4,1:N); b2=c(5,1:N); ea1=e(1,:); ea2=e(2,:); ea3=e(3,:);eb1=e(4,:); eb2=e(5,:); i=1:N; figure(2)
plot(i,a1,'r',i,a2,'b',i,a3,'y',i,b1,'k',i,b2,'g'); legend('a1','a2','a3','b1','b2'); grid on;
title('递推最小二乘参数估计') figure(3) i=1:30;
plot(i,ea1(1:30),'r',i,ea2(1:30),'b',i,ea3(1:30),'y',i,eb1(1:30),'k',i,eb2(1:30),'g'); grid on;
legend('参数a1误差','参数a2误差','参数a3误差','参数b1误差','参数b2误差'); title('被估参数误差收敛情况')
程序运行结果如图7至图9所示。
图7 输入及输出信号
百度搜索“77cn”或“免费范文网”即可找到本站免费阅读全部范文。收藏本站方便下次阅读,免费范文网,提供经典小说综合文库系统辨识大作业(3)在线全文阅读。
相关推荐: