[摘要]改进了一种寻找匹配于特定信号的自适应正交小波基的方法,并将其应用于一些实用信号与原方法以及Daubechies小波进行对比。实验结果表明:此方法具有明显的优点,小波的设计更加灵活,逼近误差更小,而且在尺度滤波器长度不变的情况下,可以灵活设计消失矩,是一种寻找最优正交小波基的有效方法。

[关键词]Mallat算法;Daubechies小波;尺度滤波器;消失矩

[中图分类号]O174.2[文献标识码]A

对正交小波而言,我们希望它是有限支撑的,以使Mallat算法[1]更加快捷;希望它是光滑的,以便高精度地模拟和分析信号;希望它的时域和频域的局部化是强劲的以便在信号处理中发挥突出作用。Daubechies小波[2]为此做出了杰出的贡献。

我们更希望在用小波进行信号处理时,对一般信号或一类信号,能够自适应地找到一组较好的小波基,以便更好地逼近信号。一种可选的方案就是寻找匹配于特定信号的自适应小波基。近些年,人们对寻找自适应小波的有效算法已经取得一定成果。在这一领域,早期Gopinath[3],Tewfik[4]等人提出了几种匹配于特定信号的正交小波基算法。Jian-kangZhang[5]在文献[3]的基础上利用内点算法[6],凸最优[7]等方法寻找最优的正交小波基,使得逼近误差有了明显改善。本文针对文献[5]中的算法,通过对其算法进行改进,进一步优化了算法,减小了误差。

在本文中,连续函数或信号用如f(t)的形式表示,其傅里叶变换(CTFT)我们用f^(ω)的形式表示;离散的函数或信号用如h[l]的形式表示,其傅里叶变换(DTFT)我们用h^(ω)(其中h^(ω)=∑lh[l]e-iωl)的形式表示;矩阵X的元素(i,j)我们用[X]ij表示,半正定矩阵A我们用A≥0表示。

1预备知识

设h[l]为0≤l≤L-1上的有限脉冲响应(FIR)滤波器的脉冲响应。那么,构成正交小波基的必要条件包括:

(1)双尺度方程:

其中,φ(t)为尺度函数,ψ(t)为小波函数;

(2)正交条件:

然而,仅仅满足以上两个条件是不够的,下面给出了构成正交小波基的充分条件:

长为L的尺度滤波器h[l]满足正交条件(3),并且存在某个整数N>0使得

h^(ω)=(1+e-jω2)Nq^(ω),(4)

其中,q^(ω)=∑

Lq-1

l=0

q(l)e-jωl且q^(0)=2;如果

|q^(ω)|<2N-1/2,ω∈(-π,π),(5)

那么,由双尺度方程(1)式定义的尺度函数是正交的,而且由(2)式定义的集合{ψj,n(t)=2j/2ψ(2jt-

n)}j,n∈Z构成L2(R)的一组正交基。

由引理1可知,对任意信号f(t)∈L2(R),我们有f(t)=∑j,n∈Zd(j,n)ψj,n(t),其中d[j,n]=

∫+∞-∞f(t)ψj,n(t)dt,由多分辨分析的知识[1]可知,如果空间VJ=span{2J/2φ(2Jt-n)}n∈Z,那么空间VJ=

span{2j/2ψ(2jt-n)}j<J,n∈Z。因此信号f(t)在VJ上的投影可以表示为PJ(f,φ)=∑na[J,n]φJ,n(t)=

∑j<J∑nd[j,n]ψj,n(t),其中a[J,n]=∫+∞-∞f(t)φJ,n(t)dt。即信号f(t)=PJ(f,φ)+∑j>J∑nd[j,

n]ψj,n(t),而信号f(t)在VJ上的投影的最小平方误差定义为:

εJ(f,φ)=‖f(t)-PJ(f,φ)‖2L2(R).(6)

对于带限于[-π,π]上的信号f(t),因为Nyquist的抽样间隔为1s,所以我们取J=0,由文献[3]可

以把(6)式转化为

ε0(f,φ)=12π∫π-π|f^(ω)|2(1-|φ^(ω)|2)dω,(7)

因为(1)式在频域等价于φ^(ω)=1

2h^(ω2)φ^(ω2),由文献[1]可知

φ^(ω)=∏∞k=1(1/2)h^(ω/2k),

通过有限去逼近无限,有

^φ(ω)≈φ~(K)(ω)=∏Kk=112h^(ω2k),(8)

将(8)式代入(7)式,得

ε0(f,φ)≈ε0~K(f,h)=12π∫π-π|f^(ω)|2(1-|φ^(K)(ω)|2)dω,(9)

显然(9)式实际上已经转化为关于尺度滤波器h[l]的关系式。

由文献[2b]我们还能得到,满足(4)式的N是0<N<L/2的,且N是小波ψ(t)的消失矩。

2算法

由以上这些知识,现在我们来引出关于如何设计尺度滤波器h[l]问题:

已知带限于[-π,π]上的信号f(t),其尺度滤波器的长度为L,0<N<L/2,怎样寻找使得ε0(f,φ)最小

且满足(3),(4),(5)式的正交的尺度滤波器h[l]?

下面,我们对此问题进行分析:

如果给定了逼近误差ε0(f,φ),从(9)式是可以解得尺度滤波器h(l),文献[3]就是这样做的,但计

算很复杂,而且要求复杂的局部最小值管理,逼近效果也不理想。因为它是非线性的(非凸的),难于寻

找局部最小值,其次随着K的增大,计算会变得更加复杂。为此,文献[5]提出了一种解决方案,将目标

函数(9)式设法转化为线性函数,再利用线性规划来解决此问题。本文在满足引理1的前提下,通过消

去文献[5]中算法的一个参数M优化其算法,简化了过程,而且在尺度滤波器长度、消失矩一定的情况

下,算得的尺度滤波器系数是[5]中方法里最优的一组。步骤如下:

)式中|q^(ω)|<2N-1/2,所以必定存在一个大于零的参数λ,使得

|q^(ω)|2≤λ<22(N-1/2),(10)

那么由(10)式和文献[5]我们可以将(9)式转化为时域中的等价形式,得到

ε0(f,φ)≈minh[l],λb1[0]2+∑Lq-1k=0(-1)kb1[k]rh[k]+βN(f)λ,(11)

其中

b1[k]=12∑m∈Zf(m2)f(m+k2),(12)

rh[k]=∑m∈Zh[m]h[m+k],(13)

βN(f)=12π∫π-πω2N|f^(ω)|2dω24N+1(22N-1),(14)

因为由(12)式可知b1[0]=12∑m∈Zf2(m2)为常数,所以可以将问题简化为下列数学模型:

模型1已知带限于[-π,π]上的信号f(t),其尺度滤波器的长度为L,0<N<L/2,寻找满足(3),

(4),(10)式的正交的尺度滤波器h[l]使得

minh[l],λ∑Lq-1k=0(-1)kb1[k]rh[k]+βN(f)λ.(15)

第二步:由(4)式在时域的等价形式为

h[l]=2-N∑Nn=0cnNq[l-n],(16)

将(16)式代入(13)式,令rq[k]=∑m∈Zq[m]q[m+k](q[l]由(5)式定义),得

rn[k]=2-2N∑Nn=-Ncn+N2Nrq[k-n],(17)

那么由(17)式可以分别把(3),(4),(10),(15)式最终转化为关于rq[k]的关系式。下面我们分别进行分析:

(1)约束条件(3),(4)式根据文献[5]可以分别转化为:

正交条件:

rq[0]+2∑Lq-1k=1rq[k]=2,(18)

正则性条件:

m[2k]rq[0]+∑Lq-1n=1(m[2k+n]+m[2k-n])rq[n]=22Nδ[k],k=0,1,2,…,[L-12];(19)

(2)由文献[5],[8]可将(10)式转化为关于rq[k]与参数λ的关系式:

∑n[X1]n,n+k=rq[k],0≤k≤Lq-1,(20)

∑n[X2]n,n+k=λδ[k]-rq[k],0≤k≤Lq-1,(21)

其中[X1],[X2]为Lq×Lq的对称半正定矩阵,而由(10)式参数λ满足0<λ<22(N-1/2);(22)

(3)将(15)式转化为关于rq[k]的关系式为

minrq[k],λ,X1,X2≥0∑Lq-1k=0s[k]rq[k]+βN(f)λ,(23)

s[0]=2-2N∑L-1n=1(-1)nm[n]b1[n],s[k]=2-2N∑L-1n=1(-1)n(m[n+k]+m[n-k])b1[n],k=1,2,…,Lq-1,m[k]=cn+N2N,k=0,±1,…,±N,0,其他.

可以看出,本文中的约束条件(18),(19),(20),(21)式跟[5]中的formulation3的前4个约束条件是一样的,唯独的一点区别就是λ的取值范围变成了0<λ<22(N-1/2)(相当于[5]中的M=0),也就是说λ的取值范围变大了。由以上这些知识,我们可以将模型1转化下面的数学模型:

模型2已知带限于[-π,π]上的信号f(t),其尺度滤波器的长度为L,0<N<L/2,令Lq=L-N,

寻找满足(18),(19),(20),(21),(22)式的点列rq[k],以及参数λ和Lq×Lq的对称半正定矩阵X1,

X2,使得

第三步:将数学模型2转化为关于rq[k]与参数λ的线性规划问题,其形式如下:

mincx,约束条件:Ax=b,

其中A为矩阵,b为列向量,c为横向量,x为所求列向量。

从上面的形式看本文的改进是显而易见的。因为对线性规划来说,约束条件越宽松,目标函数cx才可能取得越小,也就是说对信号的逼近误差才可能越小。又因为本文的算法是由引理1出发的,所以我们能够保证算得的尺度滤波器,是共轭镜像滤波器,得到的小波基是正交小波基。

因为此线性规划的约束条件中含有对称半正定矩阵,所以它并不是一般的线性规划问题,要用到内点法[6]和SDP(半正定规划)[9]的知识,不过现在建立在MATLAB基础上的SeDuMi程序包[10]已经可以解决此类问题。因此,通过此种方法就得到了使得逼近误差最小的最优解rq[k],而rq[k]=∑m∈Zq[m]q[m+k],对此我们可以用谱因式分解方法[11]进行分解,解出q[l]。再由(16)式就得到我们所要找的尺度滤波器h[l]。

3实验结果及分析

在这一部分本文将举一些例子与[5]中的方法以及Daubechies尺度滤波器进行对比。

例1为了进行对比,本文仍以信号f1(t)=sin(πt)/πt为例,此信号在频域的[-π,π]上有着单位矩形谱,因此是带限于[-π,π]上的。表1是消失矩N=6与尺度滤波器长度L分别取16,18时,用我们的方法得到对信号f1(t)的逼近误差与[5]中M分别取1,2时对信号f1(t)的逼近误差对比。ε0(f,φ)在这里我们用(9)式的ε0~(K)(f,φ)来逼近,其中K=10.从表1可以看出用我们的方法得到的逼近误差要比[5]中M=1,2时得到的逼近误差要小。实验证明,我们的方法是正确的,得到的滤波器系数是最优的。

例2我们以最常用的高斯信号f2(t)=e-t2/4/(2π)为例,其傅里叶变换为f^2(ω)=e-ω2,仍为高斯函数,是几乎带限于[-π,π]上的。表2是我们取尺度滤波器长度L=20时,与长为20的标准波器的逼近误差对比(表2中,Daubechies尺度滤波器的逼近误差为3.9584×10-7),ε0(f,φ)在这里我们用(9)式的ε0~(K)(f,φ)来逼近,其中K=10.可以发现此种方法的逼近误差相对于长为20的Daubechies尺度滤波器的逼近误差的改善程度是非常大的,尤其是N=4时,改善了将近一倍。

我们选取表2中消失矩N=4的尺度滤波器来进行详细说明。图1是我们通过双尺度方程得到的尺度函数与长为20的Daubechies滤波器的尺度函数对比,可以发现我们的尺度函数与Daubechies的有着相似的形状。图2是相应的尺度滤波器h[l]的频率响应对比,可以看到我们的尺度滤波器在带通与带阻,低通与高通之间转化更快,更能显示低通特性。

图1关于信号f2(t)的尺度函数对比

图2关于信号f2(t)的尺度滤波器的频响对比

例3此例本文用我们的方法与[5]中的方法以及Daubechies的方法得到的重构信号进行对比。图3a是我们从MATLAB中读取的一个波形函数leleccum;图3b是用我们的方法取尺度滤波器L=8,消失矩N=4时的尺度滤波器得到的重构信号;图3c是我们用[5]中的方法L=8,N=4,M=1时得到的重构信号;图3d是长为8的标准Daubechies尺度滤波器得到的重构信号。通过对比可以发现我们的方法得到信号逼近效果比其他两种逼近效果要好。

图3信号重构对比

4讨论

本文以构成交小波基的充分条件(引理1)为直接出发点进行推导,通过将问题转化为线性规划

的方式寻求最优的紧支撑正交小波基,并将其应用于一些实用的信号。实验证明,通过此种方法设计的正交小波基是灵活有效的,不仅能够灵活设计小波消失矩,而且使得信号的逼近误差得到明显减小;并且由于我们是建立在线性规划基础上的,所以我们得到的最优正交小波基一定是全局最优的。实例证明此种方法是一种寻找最优的正交小波基的有效方法。

[参考文献]

[1]StephaneMallat.信号处理的小波导引[M].杨力华译.第2版.北京:机械工业出版,2002.173.

[2]DaubechiesI.TenLecturesonWavelets[M].PA:SIAM,Philadelphia,1992.a:182,b:171.

[3]GopinathRA,OdegardJE,BurrusCS.OptimalWaveletRepresentationofSignalsandtheWaveletSamplingTheorem[J].IEEETrans.CircuitsSystems,1994,41(4):262—277.

[4]TewfikAH,SinhaD,JorgensenP.Ontheoptimalchoiceofawaveletforsignalrepresentation[J].IEEEInform.The-ory,1992,38(2):747—765.

[5]Jian-kangZhang,TimothyNDavidson,KonMaxWong.Efficientdesignoforthonormalwaveletbasesforsignalrepre-sentation[J].IEEEtrans.signalprocessing,2004,52(6):1983—1996.

[6]InteriorYYe.PointAlgorithms:TheoryandAnalysis[M].NewYork:Wiley,1997.

[7]AkkarakaranS,VaidyanathanPP.Filterbankoptimizationwithconvexobjectivesandtheoptimalityofprincipalcompo-nentforms[J].IEEETrans.SignalProcessing,2001,49(1):100—114.

[8]DavidsonTN,LuoZQ,SturmJF.LinearmatrixinequalityformulationofspectralmaskconstraintswithapplicationstoFIRfilterdesign[J].IEEETrans.SignalProcessing,2002,50(11):2701—2715.

[9]VandenbergheL,BoydS.Semidefiniteprogramming[J].SIAMRev.,1996,31(5):49—95.

[10]SturmJF.UsingSeDuMi1.02,aMATLABtoolboxforoptimizationoversymmetriccones[J].Optim.Meth.Software,1999,11—12:625—653.

[11]GoodmanTNT,MiccchelliCA,RodriguezF,eta.lSpectralfactorizationofLaurentpolynomials[J].-pute.Math.,1997,7(4):429—454.

SignalprocessingmethodbasedonadaptivewaveletsLIWan-she,RENJun-wei(SchoolofMathematicsandInformationScience,ShaanxiNormalUniversity,Xi′an710062,China)

Abstract:Thispaperimprovedasortofanefficientmethodforselectinganorthonormalwaveletbasewhichmatchestoagivensigna.lThismethodwasappliedtosomepracticalsignalsandcomparedwithoriginalmethodandDaubechieswavelets.Theexperimentsproved:Severalobviousadvantagescanbediscovered,suchasdesigningismoreflexible,theapproximateerrorissmallerandthedesignofzeromomentsunderthesamelengthofscalingfilterismoreflexible.Soitisanefficientmethodtofindoptimalwaveletbases.

Keywords:Mallatalgorithm;Daubechieswavele;tscalingfilter;zeromoment