自适应滤波算法原理及其应用.docx

1、自适应滤波算法原理及其应用自适应滤波算法原理与应用经典的滤波算法包括，维纳滤波，卡尔曼滤波，自适应滤波。维纳滤波与卡尔 曼滤波能够满足一些工程问题的需求，得到较好的滤波效果。但是他们也存在局限 性，对于维纳滤波来说，需要得到足够多的数据样本时，才能获得较为准确的自相 关函数估计值，一旦系统设计完毕，滤波器的长度就不能再改变，这难以满足信号 处理的实时性要求；对于卡尔曼滤波，需要提前对信号的噪声功率进行估计，参数 估计的准确性直接影响到滤波的效果。在实际的信号处理中，如果系统参数能够随 着输入信号的变化进行自动调整，不需要提前估计信号与噪声的参数，实现对信号 的自适应滤波，这样的系统就是自适应滤

2、波系统。1基本自适应滤波算法自适应滤波算法的基本思想是根据输入信号的特性自适应调整滤波器的系数, 实现最优滤波。图1自适应滤波结构框图若自适应滤波的阶数为M,滤波器系数为W,输入信号序列为X,则输出为:A/-13)=艺“，(加兀(一7) ( 1)加= d (口) 一 y(n) (2)其中d()为期望信号，巩“)为误差信号。Af-1 MyW =工 w(in)x(n 一 m) T yj =工、內 (3)m-0 i-l令 W = %，吗，叫_$,色=勺,切,，兀J ( 4)则滤波器的输出可以写成矩阵形式：刀=X：W = WTX, ( 5)勺=d厂兀=dj -X；W = dj -心j ( 6)定义代价

3、函数:(7)J(j) = Ee2. = Edry.r= Edj-W,Xj)1当使上式中的代价函数取到最小值时，认为实现最优滤波，这样的自适应滤波 成为最小均方自适应滤波（LMS）o对于最小均方自适应滤波，需要确定使得均方误差最小的滤波器系数，一般使 用梯度下降法求解这类问题。滤波器系数向量的迭代公式为：(8)咋严叱+”（V厶）式中，为步长因子，人为代价函数的梯度。(9)因为瞬时梯度-2X丿勺为真实梯度值的无偏佔讣，实际应用中可使用瞬时梯度代替真实梯度，即有：X旳(11)通过逐步迭代，即可得到最优的滤波器系数，实现对输入信号的自适应滤波。2自适应滤波的工程应用为了比较不同滤波算法的滤波效果，这里

4、仍然采用前面用到的二维圆周运动轨 迹追踪的问题作为工程背景。自适应滤波算法的程序设计思路如图2所示。生成现测佶兮丄初鮒化滤波系数与步长1计算误差1更新淀波系数工得到淀波结果ng 心 TV输出淀波佶专图2自适应滤波算法流程图迭代步长/ = 0.2时，得到的滤波结果为：3000 r2000 -1000A I l 0 50 100 150 200 250 300 350 400 450 500图3 X方向自适应滤波结果基本自适应滤波/ = 0.2图4 Y方向自适应滤波结果基本自适应滤波“ =02从X与Y方向上的位移变化曲线与方差变化曲线上可以看出，滤波结果出现了 发现，最终得到的结果并没有达到最优解

5、。分析其原因，可能是迭代步长太大，将 迭代步长减小之后，取 =01,得到较为理想的滤波结果，示于图5和6.图5X方向自适应滤波结果基本自适应滤波“ =01图6 丫方向自适应滤波结果-基本自适应滤波“ =0.1可以看出，减小步长因子之后，两个方向上的滤波轨迹与期望的轨迹之间的误 差明显减小，证明了自适应滤波的有效性。3自适应滤波的收敛性分析在上一节的讨论中，迭代步长选择对于算法的收敛性具有决定性作用，步长值 的微小改变即可对算法的收敛效果产生明显影响，因此如何确定合适的步长值是自 适应滤波算法中重要的内容。Ee； = E0-兀戸=Ed 訂一 2E(d ,XjW + WyE X, X ； W (

6、12)=Edj - 2R +V.=2/?xrW-2/?rfv ( 13)系统的最小均方误差最小时，有：V,=0则下式成立： ”=心 (14)对于滤波器系数的迭代过程，有：EWjJ = EWj + EejXj= EWj + E(dj-W,Xj)XJ= EIVJ + /E(/?dv-/?JV,) ( 15)=(1-“心)EW+“心对自相关矩阵进行分解，即： RZZQ ( 16)则相邻两次迭代过程的滤波器系数之间满足关系式：(17)QHEWM-Wopi= (I-)Q,EWj -%= (I-yQnEW.-W(pl咼怙=除 +Q(I-yQll(W0-Vi)t) ( is)当迭代次数为无穷大时，理论上可以

7、实现最优滤波，即迭代步长应该满足：lim/-/M?=O II一“& Ivl 21,2,W ( 19)从而有：九 (20)式20即为确保算法收敛迭代步长应满足的条件。得到步长的收敛性条件，即可在满足要求的范圉内调整步长因子，选择最佳的 步长，在确保算法收敛的前提下，提高收敛速度。对于二维轨迹追踪问题，取步长 因子为“ = 0.6“，得到的滤波结果如图7至9所示。图7 X方向自适应滤波结果-基本自适应滤波“ =图8 丫方向自适应滤波结果-基本自适应滤波“ =0.6“期冬丸迹LMS自适应遥浪回周运动轨述启踪图9二维圆周运动轨迹滤波结果-基本自适应滤波“ =0.6“亦从X方向，Y方向上的滤波结果可以看

8、出，滤波轨迹在起初的一段时间内与期 望轨迹存在较大的误差，但随着迭代次数增加，两者的误差逐渐减小，最终得到误 差的最小值。二维轨迹图上也能得到类似的结论。4变步长自适应滤波在满足收敛性条件的要求下选择迭代步长，可以确保最终得到收敛的结果，但 是这一步长在整个过程中是固定的。然而，更为理想的情况是在滤波的初始阶段， 误差值很大时，迭代步长可以取较大的值，以取得较快的收敛速度，随着误差减小， 逐渐接近最优L1标时，迭代步长也相应减小，从而得到较好的收敛精度，这就是变 步长自适应滤波算法。变步长的自适应滤波算法已经有了较长时间的发展，前人发展了很多有效的变 步长算法，这里仅选择两种常用的方法。（1）

9、 归一化变步长自适应滤波算法，= 0+X兀咲严咒+ “旳X,其中0为常数，且满足0va0,030“,皿。Sigmod函数使用滤波器的输出误差对迭代步长进行控制，从表达式中可以看出,误差较大时，步长因子的值较大，误差减小时，步长因子的值也会相应减小。图10变步长自适应滤波算法程序设计流程图釆用变步长的自适应滤波算法对二维圆周运动的轨迹进行追踪，滤波结果示于 图至13。其中参数a = l2,0 = 2。图13二维圆周运动轨迹滤波结果-变步长自适应滤波从X方向与Y方向上的滤波曲线可以看岀，变步长的自适应滤波输出结果与期 望信号之间的误差更小，固定步长时起始阶段的大幅度波动也消失了，对运动轨迹 的追踪

10、效果也更好。5解相关自适应滤波当输入信号之间具有较强的相关性时，自适应滤波的效果并不理想，因此改进 自适应滤波算法的一个方法就是消除相邻输入信号序列的相关性，称为解相关自适 应滤波。解相关自适应滤波算法的实现过程为：Zj = XjiXji该算法通过求解相邻两个输入信号序列的相关系数，在当前输入信号中减去与 上一输入信号的相关部分，作为当前的输入信号，实现解相关的自适应滤波。图14 给出了解相关自适应滤波算法的程序设计流程。生成观测信号1初始化滤波系数更新步长与波波器系数变换至时域得到滤波结果_输出滤波信号图18变换域自适应滤波算法流程图参考文献李方伟,张浩.一种新的变步长LMS自适应滤波算法及

11、其仿真J.重庆邮电大学学 报(自然科学版),2009,(05):591-594.2齐林，周丽晓.变换域自适应滤波算法的研究J.郑州大学学报(理学 版),2007,(01):61-66.冯存前，张永顺.变步长频域快速自适应收发隔离算法研究电子对抗技 术,2004,(05):22-25+45.4Deherty J, Porayath R A robust echo canceler for acoustic environments. I EEE Trans. Circuits and Systems, llz 1997, 44:389-398.张贤达现代信号处理(第三版).北京:清华大学出版社,

12、2015.高西全，丁玉美.数字信号处理-时域离散随机信号处理.西安:西安电子科技大学 出版社，2002.代码：自适应滤波算法：%该程序实现对二维圆周运动轨迹的自适应滤波%该程序为主函数，调用不同的子函数实现不同的滤波方法%子函数1： funplms_filter2-固定步长自适应滤波%子函数2： fun_cplms_filter2-变步长自适应滤波%子函数3： fun_lms_filter_der2-解相关自适应滤波clearclose allN=2000;theta=linspace(0,2*pi,N); e_x=cos(theta);e_y二sin (theta);%极坐标参数%x, y方

13、向上的期望信号%高斯白噪声%观测信号no_x=normrnd(0,sqrt(0.08)/l,N); no_y=normrnd(0,sqrt(0.12),1,N); m_x=e_x+no_x;m_y 二 e_y+no_y;%fixed step% Err_x,f_x=funplms_filter2(e_x,m_x,N,10);% Err_y,f_y=fun_fplms_filter2(e_y/m_y,N/10); %changed step% Err_xzf_x=fun_cplms_filter2(e_x7m_xzN,10);% Err_y,f_y =fun_cplmsilter2(e_y,m_

14、y,N JO); %decorrelati on Err_x,f_x=funJms_filter_der2(e_x,m_x,N,10);Err_y,f_y=fun _lms_f ilter_der2(e_y,m_y,N,10)； figureplot(e_x/e_y,k,linewidth,/2)hold onplot(m_x,m_y,b,)hold onplot(f_x,f_y,W)title(LMS自适应滤波圆周运动轨迹追踪) legendC期望轨迹丁观测轨迹T滤波轨迹*) figuresubplot(211)plot(e_x/lk,)hold onplotfmx/b1)hold onpl

15、ot(f_x,r)title(x方向上信号滤波效果对比1legend/期望信号T观测信号T滤波信号,4)subplot(212)plot(Err_x,k)titlefx方向上滤波方差变化曲线)figuresubplot(211)Plot(e_y；k-)hold onplot(m_y；b)hold onplot(f_y；r)title(y方向上信号滤波效果对比*)legend(期望信号丁观测信号T滤波信号,4)subplot(212)plotfErry/k)title(y方向上滤波方差变化曲线*)function SE,x_f=fun_fplms_filter2(x0/xm/n,m)%this

16、fun ction con ducts the adaptive filteri ng with fixed step lengthx_e=xO;x_mO=xm;N=n;M=m;x_f二x_mO; %order of filter and initial weight valuesw=zeros(l,M);SE=zeros(l,N);% fundmental LMS adptive filter Modern SP Zxd P183rxx=xcorr(x_mO)/N;Rxx=toeplitz(rxx(N:e nd);mui_max=l/max(eig(Rxx);% trace(Rxx)% mu

17、i_max=l/trace(Rxx); %convergenee conditionmui=0.6*mui_max; %initial step length% normallized LMS Modern SP Zxd P183% alpha=0.8;beta=2;% the iterative filterx_m=zeros(l,M) x_mO;for i=l:Nun=x_m(M+i:-l:i+l);u_out=sum(u.*w); err=x_e(i)-u_out;% mui=alpha/(beta+sum(un42);% mui=0.06;w=w+mui*u_i nerr;x_f(i)

18、=u_out;se=x_e-x_f;SE(i)=sum(se.A2)/N;end function SE,x_f=fun_cplms_filter2(xO/xm/n,m)%this fun ction con ducts the adaptive filteri ng with vari ng length x_e=xO; %parameter in function modex_mO=xm;N=n;M=m;x_f=x_mO; %order of filter and initial weight valuesw=zeros(l,M);SE=zeros(l,N);% fundmental LM

19、S adptive filter Modern SP Zxd P183 % rxx=xcorr(x_mO)/N;% Rxx=toeplitz(rxx(N:end);% mui_max=l/max(eig(Rxx);% % trace(Rxx)%convergenee condition%initial step length% % mui_max=l/trace(Rxx);% mui=0.6*muLmax;% normallized LMS Modern SP Zxd P183% alpha=0.8;beta=2;alpha=-6.6;beta=0.18;a=2; %Sigmoid fucnt

20、ion% the iterative filterx_m=zeros(l,M) x_mO;for i=l:Nun=x_m(M+i:-l:i+l);u_out=sum(unw);err二 x_e(i)u_out;mui=beta*(l-exp(alpha*errAa);% mui=alpha/(beta+sum(un42);w=w+mui*u_i nerr;x_f(i)=u_out;se=x_ex_f;SE(i)=sum(se.A2)/N;end% x_f(M:M+3)=x_m(M:M+3);% to improve the initial steps of filter% Err=x_f-x_

21、e;% Re=Err./x_e;% find(abs(Re)=max(abs(Re)% plot the filter outputsfigureplot(x_e,k,)hold onplot(x_mO/b)hold onplotfxj/r)% legend(,expected,z,measured1,filtered)legend/期望信号T观测信号T滤波信号)% figure% plotfSE/k)% title)自适应滤波方差变化曲线1function SE,x_f=funms_filter_der2(x0,xm,n,m)%this fun ction con ducts the ada

22、ptive filteri ng using de-corellatio n method N=n;M=m;x_e=xO;x_mO=xm;x_f=x_mO;SE=zeros(l,N);%order of filter and initial weight valuesw=zeros(l,M);% fundmental LMS adptive filter Modern SP Zxd P183% rxx=xcorr(x_m)/N;% Rxx=toeplitz(rxx(N:end);% mui_max=2/trace(Rxx); %convergenee condition% mui=0.4*mu

23、i_max; %initial step length% normallized LMS Modern SP Zxd P183% alpha=0.5;beta=l;% decorrelation LMS Modern SP Zxd P183%to cancle the correlation between u(n) and u(n-l)rho=0.8;% the iterative filter x_m=zeros(l,M) x_mO;u_inO=x_m(M+l:-l:2);SE(l)=sum(x_e-x_f).A2)/N;for i=2:Nunt=x_m(M+i:-l:i+l);err=x_e(i)-sum(u _in t.*w);corr=sum(ut.*u_inO)/sum(un0卢un0); %correlation cofficientv=ut-corr*u_inO;mui=rho*err/sum(u_int.*v);w=w+mui*v;u_out=sum(u_int.*w);x_f(i)=u_out;u_inO=u_i nt;se=x_ex_f;SE(i)=sum(se.A2)/N;end