实验二 自适应滤波信号一、实验目的：1．利用自适应LMS 算法实现FIR 最佳维纳滤波器。

2．观察影响自适应LMS算法收敛性，收敛速度以及失调量的各种因素，领会自适应信号处理方法的优缺点。

3．通过实现AR 模型参数的自适应估计，了解自适应信号处理方法的应用。

二、实验原理及方法自适应滤波是一种自适应最小均方误差算法（LMS ），这种算法不像维纳滤波器需要事先知道输入和输出信号的自相关和互相关矩阵，它所得到的观察值,滤波器等价于自动“学习”所需要的相关函数，从而调整FIR 滤波器的权系数，并最终使之收敛于最佳值，即维纳解。

)(n y 下面是自适应FIR 维纳滤波器的LMS 算法公式：（2-1）)()()(0^^m n y n h n x Mm m -=∑= （2-2） ^)()()(n x n x n e -=M m m n y n e n h n h m m ⋯=-∙+=+,1)()(2)()1(^^μ （2-3）其中FIR 滤波器共有M+1个权系数，表示FIR 滤波器第m 个权系数在第n 步的估计值。

),0)((^M m n h m ⋯=因此，给定初始值)M ,0(),0(⋯=m h m ，每得到一个样本，可以递归得到一组新的滤波器权系数，只要步长)(n y μ满足max10λμ<< （2-4）其中max λ为矩阵R 的最大特征值，当∞→n 时，)M ,0(),0(⋯=m h m 收敛于维纳解。

现在我们首先考察只有一个权系数h 的滤波器，如图2.1所示。

假如信号由下式确定：)(n y )()()(y n w n s n += （2-5） )()(n hx n s = （2-6） 其中h 为标量常数，与互不相关，我们希望利用和得到)(n x )(n w )(n y )(n x )(n s图1利用公式（2-1），（2-2），（2-3），我们可以得到下面的自适应估计算法：（2-7） )()()(^^n x n h n s = （2-8） )())()()((2)()1(^^^n x n x n h n y n h n h -+=+μ其框图如图所示。

图2选择的初始值为，对式2-8取数学期望可得^)(n h ^)0(h （2-9） ))0(()21(])([^^h h R h n h E n--+=μ其中（2-10） ])()([R T n x n x E =因此，只要满足R10<<μ （2-11） 的条件，总归可以收敛于最佳值，从而也逐渐收敛于。

^)(n h h ^)(n s )(n s 自适应信号处理的一个重要应用是用来进行参数估计。

下面是利用LMS 算法实现AR 模型参数的估计。

如果信号为一个M 阶的AR 模型，即)(n y )()()2()1()(21n w M n y a n y a n y a n y M +--⋯-----= （2-12） 通过解Yule-Walker 方程可以得到AR 模型的参数估计，同样，利用LMS 算法，我们也可以对AR 模型的参数估计进行自适应估计，其算法如下：（2-13）∑=--=Mm m m n y n a n y 1^^)()()( (2-14) ^)()()(n y n y n e -= (2-15)M m m n y n e a n a m m ≤≤-∙-=+1)()(2)1(^^μ这种算法的实现框图如图2.3所示。

图3同样可以证明，只要步长μ值选择合适，当∞→n 时，上述自适应算法得到的)(^n a M =μ也收敛于AR 模型的参数。

m a ([h E 三、实验结果及分析1. 时8.0,1,03.02-==h w σ)]n 和)(n h 的比较：图4 自适应滤波器)]([n h E 和)(n h的比较图1是时100=L )]([n h E 和)(n h的比较，由图可以看出的均值刚开始一直在变化，然后逐渐趋于稳定，最后收敛于最佳值。

^)(n h2. 自适应滤波器的效果和)(n s )(n s的比较：图5 和)(n s )(n s的比较图2是时和100 L )(n s )(n s的比较比较，由于滤波器逐渐收敛于最佳值，所以随着L 的增大滤波效果越好。

3.利用实验一中的维纳滤波器估计的h 和估计的到的)(n s )(n s与2中的估计进行比较：图6 自适应滤波器和维纳滤波器的估计比较图3是利用卷积得到的100=L 时的)(n s 和有2中得到的)(n s进行比较，可得出自适应滤波算法大大改善了滤波性能。

4、的情况8.0,1,01.02-===h w σμ图701.0=μ时自适应滤波器的性能8.0,1,1.02-===h w σμ的情况图81.0=μ时自适应滤波器的性能8.0,1,12-===h w σμ的情况图91=μ时自适应滤波器的性能由图4-图6可知，当μ增大时，自适应滤波器的性能随之降低。

当1=μ时滤波器失配。

5、时：8.0,01.0,03.02-===h w σμ图10时自适应滤波器的性能8.0,01.0,03.02-===h w σμ由图7可知，其他条件不变的情况下，随着 μ的减小，收敛的更快，失调量越大。

)(n h 6. 自适应AR 模型：当时和的估计情况：1,01.0,8.0,3.1,100,2,2221===-====w a a L p M σμ1a 2a图11和的估计情况1a 2a7. 通过实验一解Yule-Walker 方程的方法解的和的估计值和通过自适应AR 模型得到的估计值的比较1a 2a图12 Yule-Walker 方程得到的估计值和自适应AR 模型得到的估计值的比较8、时自适应AR 的估计情况：01.0,01.0,8.0,3.1,100,2,2221===-====w a a L p M σμ图13时a1和a2的估值01.02=w σ图10是时的a1和a2的估计值，由图可以看出，虽然a1和a2的估计值都逐渐收敛，但是已经产生了很大的误差，收敛速度比时的情况慢，失调量也变大。

01.02=w σ12=w σ四、思考题1、s(n)和w(n)的方差越大，收敛速度越快，失调量也越大；在噪声方差越大时，随即起伏项越大，表现为权值振动也就越大，为了追求更快的收敛性，往往增大μ值，但μ值越大，就使得失调系数也就越大，在有限次起动迭代过程中，也就很难收敛到较稳态值，所以必须选择合适的μ值。

2、定义权矢量A=[]T ^^2^1)(,,)(),(n a n a n a m ⋯且X(n)=[X(n-1),X(n-2),⋯,X(n-m)] T 误差信号为：（3-1） An X n y n X A n y n y n y n e T T )()()()()()()(^-=-=-=误差平方和为（3-2） A n X n X A A n x n y n y n e T T T )()()()(2)()(22+-=上式两边去数学期望得均方误差为（3-3） A n X n X E A A n X n y E n y E n e E T T T )]()([)]()([2)]([)]([22+-=定义互相关函数行向量：XY T R （3-4） )]()([n X n y E R T XY T =和自相关函数矩阵（3-5） )]()([n X n X E R T XX =则均方误差（3-3）式可表示为（3-6） A R A A R n y E n e E XX T XY T +-=2)]([)]([22这说明均方误差是权系数向量A 的二次函数，它是一个中间向上凹的抛物形曲面，是具有唯一最小值的函数。

调节权系数使均方误差为最小，相当于沿抛物形曲面下降找最小值。

可以用梯度来求该最小值。

将（3-6）对权系数A 求导，得到均方误差函数的梯度AR R n a n e E n a n e E n e E n XX XY m 22)()]([,,)()]([)]([)(^2^122+-=⎥⎥⎦⎤⎢⎢⎣⎡∂∂⋯∂∂=∇=∇T （3-7） 令，即可求出最佳权系数向量0)(=∇n （3-8） XY XX opt R R A 1-=将A 带入（3-3）式得最小均方误差opt XY T A R n Y E n e E -=)]([)]([2min 2 （3-9） 根据最陡下降法，“下一时刻”权系数向量应该等于“现在时刻”权系数向量加上一个负均方误差梯度^)1(+n a m )n ^)(n a m (∇-的比例项，即 （3-10） )()()1(^^n n a n a m m ∇-=+μ式中，μ是一个控制收敛速度与稳定性的常数，称之为收敛因子。

一种的近似计算方法是：直接取作为均方误差的估计值，即 )(n ∇)(2n e )]([2n e E （3-11） )]([)(2)]([)(2^n e n e n e n ∇=∇=∇式中 )]([n e ∇为 （3-12） )()]()([)]([n X n X A n y n e T -=-∇=∇将式（3-12）带入式（3-11）得梯度估值)()(2)(^n X n e n -=∇于是Widrow-Hoff LMS 算法最终为 )()(2)()1(^^n X n e n a n a m m μ-=+其中X （n ）相当于y(n)的时刻以前的输入，即y （n-m ）故 )()(2)()1(^^m n y n e n a n a m m --=+μ五、程序源代码a ）自适应滤波器：%% 数字信号处理： 自适应信号滤波器 （1）clear;clc;%% input parameterL=input('样本个数L= ');u=input('\n步长u= ');sigam=input('\n方差sigam=');h=input('\n一阶自适应滤波器的参数h=');h_=zeros(1,L);h_(1)=input('\n估计之前估计一阶自适应滤波器的参数h_='); %%xn=randn(1,L);wn=sigam*randn(1,L);sn=h*xn;yn=sn+wn;figure(1)stem(sn,'r')hold onstem(yn,'b')title('原信号（红色）和观测信号（蓝色）')R=mean(xn.*xn);% R=1/L*sum((xn).*(xn))sn_=zeros(1,L);% E_h_=zeros(1,L);for ii=1:Lsn_(ii)=h_(ii)*xn(ii);% E_h_(ii)=h+(h_(1)-h)*(1-2*u*R)^ii;h_(ii+1)=h_(ii)+2*u*(yn(ii)-h_(ii)*xn(ii))*xn(ii);endE_h_=h+(h_(1)-h)*(1-2*u*R).^(1:L);stem(E_h_,'r')hold onstem(h_(2:L+1),'b')title('自适应滤波器估计h期望（红色)和其估计值h（蓝色）')figure(3)stem(sn,'r')hold onstem(sn_,'b')title('原信号（红色）滤波后的信号（蓝色）')sn_n=conv(xn,h)figure(4)stem(sn_,'b')hold onstem(sn_n,'r')title('自适应滤波器得到的原信号（绿色）维纳滤波器得到的原信号（红色）')b）自适应AR模型：%% 数字信号处理：自适应滤波器（2）%%%%%close all;clear;clc;%% input parameterL=100;sigam=0.01;M=2;p=2;u=0.01;a1=-1.3;a2=0.8;a=zeros(M,L);% a1_=zeros(1,L);% a2_=zeros(1,L);wn=sqrt(sigam)*randn(1,L);A=[1,a1,a2];yn=filter(1,A,wn);% yn=zeros(1,L);% yn(1)=wn(1);% yn(2)=-a1*yn(1)+wn(2);% for kk=3:L% yn(kk)=-a1*yn(kk-1)-a2*yn(kk-2)+wn(kk);% endyn_=zeros(1,L);en=zeros(1,L);for n=1:Lfor m=1:Mif n-m>0yn_(n)=yn_(n)-a(m,n)*yn(n-m);%滤波器输出endenden(n)=yn(n)-yn_(n);%计算误差for m=1:Mif n-m>0a(m,n+1)=a(m,n)-2*u*en(n)*yn(n-m);%更新权值endendend%%%for in=1:L-1% if in==1% yn_(in)=-a1_(in)*0;% en(in)=yn(in)-yn_(in);% a1_(in+1)=a1_(in)-2*u*en(in)*0;% a2_(in+1)=0;% elseif in==2% yn_(in)=-a1_(in)*yn(in-1)-a2_(in)*0;% en(in)=yn(in)-yn_(in);% a1_(in+1)=a1_(in)-2*u*en(in)*yn(in-1);% a2_(in+1)=a2_(in)-2*u*en(in)*0;% else% yn_(in)=-a1_(in)*yn(in-1)-a2_(in)*yn(in-1); % en(in)=yn(in)-yn_(in);% a1_(in+1)=a1_(in)-2*u*en(in)*yn(in-1);% a2_(in+1)=a2_(in)-2*u*en(in)*yn(in-2); % end%endfigure(1)%stem(a1,'r')%hold onstem(a(1,1:L),'b')title('a1的值（red）及其迭代估计值（blue）') figure(2)%stem(a2,'r')%hold onstem(a(2,1:L),'b')title('a2的值（red）及其迭代估计值（blue）')%% 维纳滤波器Wn=randn(L,1);AA=[1,a1,a2];sn=filter(1,AA,Wn);r_ss = xcorr(sn,'unbiased');p=p+1;R_ss=zeros(p);for i=1:pR_ss(i,:)=r_ss(L-i+1:L+p-i);endI=zeros(1,p);for i=1:pif i==1I(i)=sigam;elseI(i)=0;endendA=inv(R_ss)*I';A1=[a(1,100),a(2,100)];figure(3),clfstem(A(2:p),'*r')hold onstem(A1,'*b')axis([-1 5 -2 2 ])title('自适应AR估计的a（红色）和Yule-Walker方程法估计的a（蓝色）')。