自适应滤波器对系统的辨识

1、自适应滤波器对系统的辨识1.自适应滤波器系统辨识原理对于一个真实的物理系统，人们主要关心其输入和输出特性，即对信号的传输特性，而不要求完全了解其内部结构。系统可以是一个或多个输入，也可以有一个或多个输出。通信系统的辨识问题是通信系统的一个非常重要的问题。所谓系统辨识，实质上是根据系统的输入和输出信号来估计或确定系统的特性以及系统的单位脉冲响应或传递函数。系统辨识和建模是一个非常广泛的概念，在控制、通信和信号处理等领域里都有重要意义。实际上，系统辨识和建模不仅局限于传统的工程领域，而且可以用来研究社会系统、经济系统和生物系统等。本节只讨论通信和信号处理中的系统辨识和建模问题。采用滤波器作为通信信

2、道的模型，并利用自适应系统辨识的方法对通信信道进行辨识，从而可以进一步地对通信信道进行均衡处理。如果把系统看成是一个“黑箱”，仅知道“黑箱”的输入和输出；以一个自适应滤波器作为这个“黑箱”的模型，并且使滤波器具有与“黑箱”同样的输入和输出。自适应滤波器通过调制自身的参数，使滤波器的输出与“黑箱”的输出相“匹配”。这里的“匹配”通常指最小二乘意义上的匹配(本节讨论的都以均方误差最小为准则的，即以LMS作为最优的准则)。这样，滤波器就模拟了系统对信号的传输行为。尽管自适应滤波器的结构和参数与真实的系统不一样，但是它们在输入、输出响应上保持高度一致。因此，在这个意义上，自适应滤波器就是这个未知“黑箱

3、”系统的模型。利用自适应滤波器辨识系统的原理下图所示：图1 自适应滤波器系统辨识原理假定未知信道为有限冲激响应(FIR)结构，构造一个FIR结构的自适应滤波器，如图1中所示。在图中，用一伪随机系列作为系统的输入信号x(n)，同时送入未知信道系统和自适应滤波器。调整自适应滤波器的系数，使误差信号e(n)的均方误差达到最小，则自适应滤波器的输出y(n)近似等于通信系统的输出d(n)。可以证明，加性噪声v(n)的存在并不影响自适应滤波器最终收敛到最优维纳解。可以认为，具有相同输入和相似输出的两个FIR系统，应该具有相似的特性。因此，可以采用自适应滤波器的特性或其单位脉冲响应来近似替代未知系统的特性或

4、单位脉冲响应。同理，将图1中的FIR系统换成IIR系统，辨识过程相同。2.辨识结果对FIR系统辨识的结果如下图所示，图1中蓝色曲线代表待辨识系统的输出，绿色代表自适应滤波器的输出，红色代表误差曲线，随着迭代次数的增加，误差越来越小，越来越平稳。对于辨识的FIR系统，为了确定其辨识的正确性，可以将辨识出自适应滤波器的系数和仿真时所用FIR系统的系数作对比，如图2中所示的真实值为实际仿真系统的滤波器系数，估计值为自适应滤波器辨识系统的系数。对于IIR系统，我们可以做出待辨识系统的幅频响应曲线和辨识出的自适应滤波器的幅频响应曲线来对比辨识系统的准确性。当然，对FIR系统也适用，图3是待对仿真的FIR

5、系统做出的幅频响应曲线，可以看出利用自适应滤波器辨识的FIR系统是准确的。 图2 FIR系统辨识结果图3 FIR辨识系统幅频响应下来我们对一个IIR系统利用自适应滤波器进行辨识，这个IIR系统本质上是一个巴特沃斯带通滤波器。图4中展示了对其辨识的结果。图4 IIR系统辨识结果在辨识的过程中，滤波器的长度、迭代次数和迭代因子对辨识的结果有影响，下面我们分别改变滤波器的长度，滤波器的迭代次数，滤波器的迭代因子对系统进行辨识。首先改变滤波器的长度：分别用长度L=10,20,30点的自适应滤波器对系统进行辨识，结果如图5，图6和图7所示。从图5,6,7中可以看出，滤波器选择的越长，系统辨识的误差收敛的

6、也越快，误差也越小，对系统的辨识就越准确，着从幅度频率响应曲线中就可以看出来。但滤波器选择的长，相应的计算量也将会增加。图5 L=10点时对IIR系统的辨识图6 L=20点时对IIR系统的辨识图7 L=30点时对IIR系统的辨识接下来改变滤波器的迭代次数：分别为N=200,500,1000点的自适应滤波器对系统进行辨识，结果如图8，图9和图10所示。从图中可以看到，迭代次数越多，对系统的辨识越准确，辨识的误差也越小。图8 N=200点时候对IIR系统辨识结果图9 N=500点时候对IIR系统辨识结果图10 N=1000点时候对IIR系统辨识结果最后改变滤波器的迭代因子：分别为u=0.08,0.

7、008,0.0008，对系统进行辨识，结果如图11，图12和图13所示。迭代因子的选择对辨识的结果有重要的影响。从图11,12和13中可以看出，迭代因子选的合适，系统辨识收敛的速度快，辨识结果准确，如图12。迭代因子如果选择的太大的话，此时自适应滤波器可能不收敛，就会发散，如图11所示。相反，如果迭代因子选择的过小，系统的收敛就会很慢，需要较多的迭代才能收敛到稳定的值，计算耗时。如图12所示，由于迭代因子选择过小，因此在迭代500次后还没有收敛，因此用自适应滤波器辨识的系统和待辨识系统的幅度频率响应误差比较大。如果增加迭代次数，最终也会收敛到相应的最优值，但这需要耗费较长时间。图11 u=0.

8、08时对IIR系统的辨识图12 u=0.008时对IIR系统的辨识图13 u=0.0008时对IIR系统的辨识3.辨识程序代码Part1:%本程序实现利用自适应滤波器进行系统辨识;close all;clear;clc;%辨识系统%辨识FIR系统;x = randn(1,500); %输入信号;b = fir1(31,0.5); %待辨识的FIR系统;n = 0.1*randn(1,500); %噪声信号;d = filter(b,1,x); %参考信号; mu = 0.008; %下降速率;ha = adaptfilt.lms(32,mu);y,e = filter(ha,x,d);figu

9、re(1);subplot(2,1,1);plot(1:500,d;y;e);title('FIR系统的辨识');legend('参考值','输出值','误差');xlabel('时间值');ylabel('信号值');subplot(2,1,2); stem(b.',ha.coefficients.');legend('实际值','估计值');xlabel('系数');ylabel('系数值'); figure(2);

10、H1,W1=freqz(b,1);Hf1=abs(H1);plot(W1,Hf1);hold on;H2,W2=freqz(ha.coefficients,1);Hf2=abs(H2);plot(W2,Hf2,'r');title('待辨识系统和经辨识后自适应滤波器的幅频响应曲线');legend('待辨识系统幅频响应曲线','辨识系统的自适应滤波器的幅频响应曲线');%辨识IIR系统;x1=randn(1,1000); %输入信号N1,Wc=buttord(0.4,0.7,0.3,0.8,2,40);b1,a1=butter(N

11、1,Wc); % 需要辨识的FIR系统，此系统为一个巴特沃斯带通滤波器;d1=filter(b1,a1,x1); % 参考信号mu=0.008; %下降速率;ha1 = adaptfilt.lms(40,mu);y1,e1 = filter(ha1,x1,d1);figure(3);subplot(2,1,1);plot(1:1000,d1;y1;e1);title('IIR系统的辨识');legend('参考值','输出值','误差');xlabel('时间值');ylabel('信号值');su

12、bplot(2,1,2); H1,W1=freqz(b1,a1);Hf1=abs(H1);plot(W1,Hf1);hold on;H2,W2=freqz(ha1.coefficients,1);Hf2=abs(H2);plot(W2,Hf2,'r');title('待辨识系统和经辨识后自适应滤波器的幅频响应曲线');legend('待辨识系统幅频响应曲线','辨识系统的自适应滤波器的幅频响应曲线');%step1:产生用于辨识系统的信号，白噪声声或者M序列都可以用来辨识系统，这里选用白噪声作为辨识信号;N=1000; %用于辨识

13、系统的信号的长度;x2=randn(1,N);z(1)=-1;z(2)=0;for i=3:N z(i)=1.5*z(i-1)-0.7*z(i-2)+x2(i-1)+0.5*x2(i-2)+x2(i);end%step2: 利用自适应滤波器开始辨识系统;for M=10:10:30 %自适应滤波器的长度;分别给出的M=10,20,30时候的情况; w=zeros(1,M); y=zeros(1,N); u=0.005; for i=1:N if i<M for j=1:i, y(i)=y(i)+w(j)*x2(i+1-j); % 用卷积求y(i) end e(i)=z(i)-y(i);

14、for j=1:i, w(j)=w(j)+2*u*e(i)*x2(i+1-j); % 修正w(n) end else for j=1:M, y(i)=y(i)+w(j)*x2(i+1-j); % 用卷积求y(i) end e(i)=z(i)-y(i); for n=1:M; w(n)=w(n)+2*u*e(i)*x2(i+1-n); % 修正w(n) end end end figure(4); subplot(2,1,1); if M=20 plot(1:N,z;y;e); title('自适应滤波器长度M=20时候系统辨识'); legend('参考值',&

15、#39;输出值','误差'); hold off; end %绘制辨识系统的幅频响应曲线和自适应滤波器的幅频响应曲线; B=1,1,0.5; A=1,-1.5,0.7; H1,W1=freqz(B,A); Hf1=abs(H1); b=w; a=1; H2,W2=freqz(b,a); Hf2=abs(H2); if M=10 W10=W2; HF10=Hf2; elseif M=20 W20=W2; HF20=Hf2; elseif M=30 W30=W2; HF30=Hf2; endendsubplot(2,1,2);hold on;plot(W1,Hf1,W10,

16、HF10,W20,HF20,W30,HF30);title('待辨识系统和经辨识后自适应滤波器的幅频响应曲线');legend('待辨识系统幅频响应曲线','M=10自适应滤波器的幅频响应曲线','M=20自适应滤波器的幅频响应曲线','M=30自适应滤波器的幅频响应曲线');Part2:%辨识系统%辨识FIR系统;x = randn(1,500); %输入信号;b = fir1(31,0.5); %待辨识的FIR系统;n = 0.1*randn(1,500); %噪声信号;d = filter(b,1,x); %

17、参考信号; mu = 0.008; %下降速率;ha = adaptfilt.lms(32,mu);y,e = filter(ha,x,d);figure(1);subplot(2,1,1);plot(1:500,d;y;e);title('FIR系统的辨识');legend('参考值','输出值','误差');xlabel('时间值');ylabel('信号值');subplot(2,1,2); stem(b.',ha.coefficients.');legend('实际值

18、','估计值');xlabel('系数');ylabel('系数值'); figure(2);H1,W1=freqz(b,1);Hf1=abs(H1);plot(W1,Hf1);hold on;H2,W2=freqz(ha.coefficients,1);Hf2=abs(H2);plot(W2,Hf2,'r');title('待辨识系统和经辨识后自适应滤波器的幅频响应曲线');legend('待辨识系统幅频响应曲线','辨识系统的自适应滤波器的幅频响应曲线');%辨识IIR系

19、统;%改变迭代次数辨识效果N=200;x1=randn(1,N); %输入信号N1,Wc=buttord(0.4,0.7,0.3,0.8,2,40);b1,a1=butter(N1,Wc); % 需要辨识的IIR系统，此系统为一个巴特沃斯带通滤波器;d1=filter(b1,a1,x1); % 参考信号M=30; %滤波器的长度;mu=0.008; %下降速率(收敛因子);ha1 = adaptfilt.lms(M,mu);y1,e1 = filter(ha1,x1,d1);figure(1);subplot(2,1,1);plot(1:N,d1;y1;e1);title('N=200

20、点IIR系统的辨识');legend('参考值','输出值','误差');xlabel('时间值');ylabel('信号值');subplot(2,1,2); H1,W1=freqz(b1,a1);Hf1=abs(H1);plot(W1,Hf1);hold on;H2,W2=freqz(ha1.coefficients,1);Hf2=abs(H2);plot(W2,Hf2,'r');title('待辨识系统和经辨识后自适应滤波器的幅频响应曲线');legend('待

21、辨识系统幅频响应曲线','辨识系统的自适应滤波器的幅频响应曲线');N=500;x1=randn(1,N); %输入信号N1,Wc=buttord(0.4,0.7,0.3,0.8,2,40);b1,a1=butter(N1,Wc); % 需要辨识的IIR系统，此系统为一个巴特沃斯带通滤波器;d1=filter(b1,a1,x1); % 参考信号M=30; %滤波器的长度;mu=0.008; %下降速率(收敛因子);ha1 = adaptfilt.lms(M,mu);y1,e1 = filter(ha1,x1,d1);figure(2);subplot(2,1,1);pl

22、ot(1:N,d1;y1;e1);title('N=500点IIR系统的辨识');legend('参考值','输出值','误差');xlabel('时间值');ylabel('信号值');subplot(2,1,2); H1,W1=freqz(b1,a1);Hf1=abs(H1);plot(W1,Hf1);hold on;H2,W2=freqz(ha1.coefficients,1);Hf2=abs(H2);plot(W2,Hf2,'r');title('待辨识系统和经辨识后

23、自适应滤波器的幅频响应曲线');legend('待辨识系统幅频响应曲线','辨识系统的自适应滤波器的幅频响应曲线');N=1000;x1=randn(1,N); %输入信号N1,Wc=buttord(0.4,0.7,0.3,0.8,2,40);b1,a1=butter(N1,Wc); % 需要辨识的IIR系统，此系统为一个巴特沃斯带通滤波器;d1=filter(b1,a1,x1); % 参考信号M=30; %滤波器的长度;mu=0.008; %下降速率(收敛因子);ha1 = adaptfilt.lms(M,mu);y1,e1 = filter(ha1,

24、x1,d1);figure(3);subplot(2,1,1);plot(1:N,d1;y1;e1);title('N=1000点IIR系统的辨识');legend('参考值','输出值','误差');xlabel('时间值');ylabel('信号值');subplot(2,1,2); H1,W1=freqz(b1,a1);Hf1=abs(H1);plot(W1,Hf1);hold on;H2,W2=freqz(ha1.coefficients,1);Hf2=abs(H2);plot(W2,Hf2

25、,'r');title('待辨识系统和经辨识后自适应滤波器的幅频响应曲线');legend('待辨识系统幅频响应曲线','辨识系统的自适应滤波器的幅频响应曲线');%改变迭代因子系统辨识N=500;x1=randn(1,N); %输入信号N1,Wc=buttord(0.4,0.7,0.3,0.8,2,40);b1,a1=butter(N1,Wc); % 需要辨识的IIR系统，此系统为一个巴特沃斯带通滤波器;d1=filter(b1,a1,x1); % 参考信号M=30; %滤波器的长度;mu=0.08; %下降速率(收敛因子);h

26、a1 = adaptfilt.lms(M,mu);y1,e1 = filter(ha1,x1,d1);figure(4);subplot(2,1,1);plot(1:N,d1;y1;e1);title('N=500点、u=0.08时IIR系统的辨识');legend('参考值','输出值','误差');xlabel('时间值');ylabel('信号值');subplot(2,1,2); H1,W1=freqz(b1,a1);Hf1=abs(H1);plot(W1,Hf1);hold on;H2,W

27、2=freqz(ha1.coefficients,1);Hf2=abs(H2);plot(W2,Hf2,'r');title('待辨识系统和经辨识后自适应滤波器的幅频响应曲线');legend('待辨识系统幅频响应曲线','辨识系统的自适应滤波器的幅频响应曲线');%改变迭代因子系统辨识N=500;x1=randn(1,N); %输入信号N1,Wc=buttord(0.4,0.7,0.3,0.8,2,40);b1,a1=butter(N1,Wc); % 需要辨识的IIR系统，此系统为一个巴特沃斯带通滤波器;d1=filter(b1,a1,x1); % 参考信号M=30; %滤波器的长度;mu=0.008; %下降速率(收敛因子);ha1 = adaptfilt.lms(M,mu);y1,e1 = filter(ha1,x1,d1);figure(5);subplot(2,1,1);plot(1:N,d1;y1;e1);title('N=500、u=0.008时点IIR系统的辨识');legend('参考值','输出值','误差');xlabel(&#