系统辨识大作业(4)

图8 递推最小二乘参数估计

图9 被估参数误差收敛情况

使用递推最小二乘算法进行系统辨识的结果如表 1所示。

表 1 递推最小二乘算法的辨识结果

参数

0.9 0.9367

0.15 0.1591

0.02 -0.0138

0.7 0.6275

-1.5 -1.3558

真值 估计值

从上表中可以看出，对于含有有色噪声的系统，采用递推最小二乘法进行辨

识获得结果有较大误差，并不是无偏估计。所以，一般情况下，对于含有有色噪声的系统应采用增广递推最小二乘算法进行辨识。

（2）增广最小二乘法

根据前文所述增广最小二乘法算法流程，编写增广最小二乘法程序如下：

%%%%%%%%%%%%%%%%%%增广最小二乘法的参数辨识算法%%%%%%%%%%%%%%%%%%% clear all

close all

%--------------------------M序列的产生-----------------------% L=60 ; %移位寄存器产生的M序列周期（取L=60） y1=1;y2=1;y3=1;y4=0; %四个移位寄存器输入的初值 for i=1:L;

x1=xor(y3,y4); %第一个移位寄存器输入的信号 x2=y1; %第二个移位寄存器输入的信号 x3=y2; %第三个移位寄存器输入的信号 x4=y3; %第四个移位寄存器输入的信号

y(i)=y4; %第四个移位寄存器的输出，即M序列,幅值为0或1 if y(i)>0.5;

u(i)=-1; %M序列的值为1时，辨识的输入为-1 else

u(i)=1; %M序列的值为0时，辨识的输入为1 end

y1=x1;y2=x2;y3=x3;y4=x4; %移位寄存器的输出序列,为下一次的输入信号做准备 end figure(1);

subplot(3,1,1);

stem(u),grid on %画出M序列输入信号径线图并给图形加上网格 title('输入信号（M序列）')

%-----------------------------随机噪声的产生--------------------% v=randn(1,L); %产生一组L=60个正态分布的噪声 subplot(3,1,2);

plot(v),grid on %画出随机噪声信号 title('随机噪声信号e(k)')

%---------------------------实际噪声的产生-----------------------% %--------- ksi(k)=lamda*v(k)-ksi(k-1)-0.41*ksi(k-2)--------------% ksi(2) = 0;ksi(1) = 0; %设ksi的前两个初始值为0 lamda = 1; for k=3:L;

ksi(k) = lamda*v(k) - ksi(k-1) - 0.41*ksi(k-2); end

subplot(3,1,3) plot(ksi), grid on

title('实际噪声信号\\xi(k)')

%--------------------------辨识过程----------------------------%

z(3)=0;z(2)=0;z(1)=0;zm(3)=0;zm(2)=0;zm(1)=0; %输出采样，不考虑噪声时系统输出，不考虑噪声时模型输出，模型输出的初值

c0=[0.001 0.001 0.001 0.001 0.001 0.001 0.001 0.001]'; %直接给出辨识被辨识参数的初始值，即一个充分小的实向量

p0=10^6*eye(8,8); %直接给出初始向量P0，即一个充分大的实数单位矩阵 c=[c0,zeros(8,L-1)]; %被辨识参数参数矩阵的初始值及大小 e=zeros(8,L); %相对误差的初始值及大小 for k=4:L; %开始求K

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) + ksi(k); %系统在M序列下的输出采样信号

h1=[-z(k-1),-z(k-2),-z(k-3),u(k-1),u(k-2),v(k),ksi(k-1),ksi(k-2)]'; %为求K（k）做准备

x=h1'*p0*h1+1;x1=inv(x);k1=p0*h1*x1; %K d1=z(k)-h1'*c0;c1=c0+k1*d1; e1=c1-c0;e2=e1./c0; e(:,k)=e2;

c0=c1;

%给下一次用

%辨识参数c

%求参数误差的相对变化

c(:,k)=c1; %把递推出的辨识参数c的列向量加入辨识参数矩阵

p1=p0-k1*k1'*[h1'*p0*h1+1]; % find p(k)

p0=p1; %给下次的计算使用 end %整个循环结束 c,e, % 分离变量

a1=c(1,:);a2=c(2,:);a3=c(3,:);b1=c(4,:);b2=c(5,:); %分离出a1,a2,a3,b1,b2 d1=c(6,:);d2=c(7,:);d3=c(8,:); %分离出的d1,d2,d3 % 分离变量的收敛情况

ea1=e(1,:);ea2=e(2,:);ea3=e(3,:);eb1=e(4,:);eb2=e(5,:); %分离出a1,a2,b1,b2的收敛情况 ed1=e(6,:);ed2=e(7,:);ed3=e(8,:); %分离出d1,d2,d3的收敛情况 figure(2); i=1:L;

plot(i,a1,'r',i,a2,'r:',i,a3,'r+',i,b1,'b',i,b2,'b:',i,d1,'g',i,d2,'g:',i,d3,'g+') grid on

legend('a1','a2','a3','b1','b2','c1','c2','c3') title('辨识曲线')

figure(3); i=1:L;

plot(i,ea1,'r',i,ea2,'r:',i,ea3,'r+',i,eb1,'b',i,eb2,'b:',i,ed1,'g',i,ed2,'g:',i,ed3,'g+') %画出各个参数的收敛情况 grid on

legend('ea1','ea2','ea3','eb1','eb2','ec1','ec2','ec3') title('辨识误差曲线') figure(4);

i=1:L;plot(i,z(i),'g'), grid on;

title('辨识系统的输出响应') %画出被辨识系统的输出响应

%显示被辨识参数及参数收敛情况

程序运行结果如图10至图13所示。

图10输入数据及噪声

百度搜索“77cn”或“免费范文网”即可找到本站免费阅读全部范文。收藏本站方便下次阅读，免费范文网，提供经典小说综合文库系统辨识大作业(4)在线全文阅读。