小波变换与尺度函数
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
小波分析里,很容易混淆的一个概念就是小波函数(wavelet function)和尺度函数(scaling function)的关系。
本文将不涉及小波分析的由来及发展历史,也不谈小波分析应用,本文主要目标仅是试着解释清楚小波函数和尺度函数两者的关系,同时也解释一些小波分析中的其他必要相关概念。
当然,要更好理解小波分析,一些傅里叶变换的知识是必要的。
我们知道,傅里叶变换分三种不同但又紧密相连的形式:
1,积分傅里叶变换,时域频域都连续;
2,傅里叶级数展开,时域连续,频域离散;
3,离散傅里叶变换,时域频域都离散。
同样,在小波分析中,也有三种类似的形式。
积分(连续)小波变换(CWT),小波级数展开,以及离散小波变换(DWT)。
先看看连续小波变换,连续小波正变换为[1]:
(1)
逆变换为:
(2)
其中*号表示复共轭,为小波基函数(basis function)。
不同小波基函数,都是由同一个基本小波(basic wavelet)ψ(t),经缩放和平移生成,即:
(3)
傅里叶变换把一个信号f(t)分解为一系列不同频率正弦型信号的叠加,而傅里叶变换系数就代表不同正弦型信号的幅值。
其中,所有正弦型基函数都由傅里叶基函数生成。
类似于傅里叶基函数,所有小波基函数也由同一个基本小波生成[2]。
不同的是,傅里叶基函数是固定的正弦型信号,而基本小波并未指定,需要根据实际的信号形式,在满足基本小波约束条件下进行设计。
可以看到,连续小波变换采用积分形式,而实际应用中,我们计算的都是采样后的信号,也需要通过离散形式来处理和表达,所以更加有用的是时域频域都离散的DWT,离散小波变换。
但是离散小波变换的计算将引入三个问题:
1,数据冗余。
观察式(1),可以看到,小波变换将一个一维信号变换为二维小波系数。
同样,若信号是二维,变换后将得到三维小波系数。
这反映了小波变换的优点,变换不仅具有傅里叶变换的频域分辨率,同时具有了时域或空域分辨率。
但是一维信号用二维系数来表达,这就意味着必然有很大的冗余性。
2,与数据冗余紧密相连的就是CWT中无限数量小波的问题。
从式(3)中看出,即使我们把平移量τ 和缩放量s都离散化,仍将是一个无限的序列,无法实际应用。
矩阵形式表达这个问题就是y=Wx,W为小波基函数矩阵,x为小波系数,y为离散信号向量。
这种冗余性表现之一就是W中列数远多于行数。
相比较,傅里叶变换中的W为一个正交归一矩阵。
3,对大多数信号来说,小波变换得不到解析解,所以只能通过数值算法得到。
这样,就需要一个快速计算方法,来进行小波信号分解(decomposition)。
第一个问题,可以通过引入二进小波(Dyadic wavelet)来解决。
二进小波,把由基本小波生成小波基函数的方式表达为:
(4)
式中j决定缩放,k决定平移幅度。
这样得到的二进形式小波基函数。
但是问题并没有解决,这样的小波基函数,仍然在j和k两个维度上无限延伸。
所以我们进一步引入紧支二进小波(compact dyadic wavelets)的概念。
如果把函数f(t)和基本小波限制为仅在[0,1]区间内有非零值,在[0,1]区间外为零的函数(解释:如果f(t)不满足此条件,可以通过对f(t)的缩放平移,使其满足,然后再设计相应的基本小波。
关于为什么强加这样一个区间限制,将在后面讲述)。
紧支二进小波表达为:
(5)
其中,。
j是满足的最大整数。
比如n=3时对应着j=1且
k=1。
相应的函数的表达改为:
(6)
这样在小波变换后就不会出现增维现象。
但是仍然没有解决无限系数的问题。
于是在解决后两个问题的过程中,引入了尺度函数。
不过在讨论尺度函数之前,还需要解释其他几个相关概念,包括滤波器族(Filter Bank),子带编码(Subband Coding),多分辨率分析(Multi-Resolution Analysis,MRA)。
滤波器族:
下图是一系列带通滤波器的频域图。
图(1)
一个信号离散信号x(n)经过这一系列带通滤波器滤波后,将得到一组系数
Vi(n)。
如下图所示:
图(2)
这样,我们就把一个信号分解成了不同频率的分量。
只要这些带通滤波器的频率能够覆盖整个原信号x(n)的频谱范围,反变换时,把这些不同频率信号,按其分量大小组合起来,就可得到原信号x(n)。
这样一组带通滤波器就称为滤波器族。
滤波器族能实现将信号分为不同频率分量,从而实现分解信号并分析信号的目的。
但是在滤波器族的计算中,我们需要指定频域分割方式。
研究者们给出了一种分割方式,即均分法,从而引出了子带编码的概念。
子带编码:
子带编码通过使用均分频域的滤波器,将信号分解为若干个子带[2]。
这样是可以实现无冗余且无误差地对数据分解及重建目的。
但是Mallat在1989年的研究表明,如果只分为2个子带,可以实现更高效的分解效率。
从而引入了多分辨率分析(MRA)。
多分辨率分析:
如果子带编码时将信号带宽先对分为高通(实际为带通)和低通两个部分,对应于两个滤波器。
然后对低通部分继续等分。
下图为子带编码示意图。
图(3)
从图中看出,每次分割保留高通部分的滤波结果,因为这里已经是信号的细节了,而且通常我们分析的信号,其绝大部分能量都在低频部分。
所以高频部分的分割可以到此为止,但是低通部分仍然有更多的细节可以划分划分出来,所以将低通部分继续等分。
分割迭代进行。
这样做的优点是,我们只需要设计两个滤波器,然后每次迭代将其对分。
缺点是,频域的分割方式确定。
对于某些信号来说,这样的划分并不是最优的。
同时,可以看到,我们上面说的二进小波,使用的正是这样一种多分辨率分析方式。
根据傅里叶变换的相似性定理:
(7)
则随着式(5)中n的增大,缩放尺度j在增大,此时时域的函数将被压窄,同时对应的频域傅氏变换将带宽加倍同时中心频谱位置也加倍。
可由下图表明。
图中为j在逐渐增大过程中,对应频谱发生的改变[1]。
图(4)
这里仍然有个问题。
每次都将频谱分为剩下的一半,那实际上,我们永远也取不到整个频段。
就好比一杯水,每次都只许喝一半,那将永远无法把它完全喝完(除非喝到了原子级别,无法再分为一半)。
所以,这样分割后的函数仍然是无限多的。
为解决这个问题,终于引出了我们最初想讨论的尺度函数的概念。
尺度函数:
在上图中,我们对频域进行分割,当分割到某个频率j时,不再继续分割了,剩下的所有低频部分由一个低通滤波器来表示,这就可以实现对信号频谱的完整分割。
这个剩余低通滤波器就是尺度函数。
事实上,很容易看出,尺度函数无非就是某级多分辨率分析中的低通滤波器。
也就是图(3)中最下面一级的LP。
那么现在可以看看能不能确定缩放系数j的范围。
如果信号被采样为N个样本(设我们采样时,把N定为2的幂指次方),那么应满足j<log2(N)。
原因是由于上面在紧支二进小波那里所讲的[0,1]区间限制。
由于基本小波和函数f(t)在[0,1]区间外取零,所以假设j=log2(N),那么由于紧支二进小波中规定
,则k=0,小波函数非零值区间应满足。
而信号均匀采样后,两个样本间的距离也是1/N,这就是说,低通滤波器在时域的宽度将等于采样数据间距。
所以,用这样一个低通滤波器同信号卷积(也就是滤波)时,计算将永远只涉及一个采样点,从而失去了低通滤波器对多个采样点取加权平均的作用。
当j>log2(N)时,低通滤波器宽度将小于采样间距,
同样无法实现低通滤波的作用。
所以我们就得到了尺度j的取值上限。
这也解释了为什么紧支二进小波要对信号和基本小波有取值区间限制。
而尺度j的下限,应如上述引入尺度函数处所讲,取为同尺度函数宽度相同的值。
从而实现小波函数和尺度函数共同分割信号频域的目的。
这样,从CWT 到DWT的前两个问题都得到了解决。
现在,任何一个信号的频谱都可由一系列小波函数ψ(t)和一个尺度函数φ(t)实现完全分割。
那相应的信号重建,即式(6)可表达为[3]:
(8)
也就是说,信号分解为尺度函数系数和小波函数系数两部分。
还有最后一个问题,就是快速的小波分解算法。
也就是说,怎么才能根据尺度函数和小波函数,得到相应分解系数。
没有一个快速的算法,小波分析仍然无法应用在实际中。
关于这个问题,Mallat在1989年的MRA建立过程中,给出了一个快速算法,鱼骨算法(Herringbone algorithm)。
同时,我们也要引入尺度向量(Scaling Vector)的概念。
首先,我们明确一下尺度函数和小波函数索引j的关系。
如图(4)中j=n-1处,我们把包括j=n-1及更大频率部分都用小波函数表达,以下部分用一个尺度函数表达。
这时,尺度函数和小波函数都称为在尺度j下的函数。
也就是说,尺度j的尺度函数和尺度j及其以上的所有小波函数,可以覆盖整个频域。
这样,如果一个函数f(t)的频谱范围在频域尺度j-1内,那该函数可表达为j-1尺度下尺度函数(低频)和小波函数(高频)滤波(卷积)的输出,即:
(9)
其中λ和γ分别为尺度函数系数和小波函数系数。
类似于傅里叶级数展开,它们由函数分别同尺度函数及小波函数内积得到,
(10)
从图(4)中可以看到,如j=n-3下的尺度函数频谱等同于j=n-2下尺度函数频谱对分后的低频部分。
推而广之,我们可以说,某个尺度j下的尺度函数,可以表达为j+1下尺度函数与一个低通滤波器作用(卷积)后的结果。
即[1]:
(11)
这个关系表达了两个相邻尺度下,尺度函数间的关系,因此称为双尺度关系(two-scale relation),或多分辨率公式(multiresolution formulation)。
同理,还可从图(4)中看出,尺度j下的小波函数可由尺度j+1下的尺度函数经一个高通滤波器作用(卷积)后得到,
(12)
如果我们把式(10)中的φ和ψ分别用式(11)和式(12)中的表达进行替换(推导过程省略,感兴趣的请参阅[1]),我们将得到一个重要表达式[1]:
(13)
(14)
这个表达式说明,尺度j-1下的小波函数系数和尺度函数系数都可以通过上一个尺度j下的尺度函数系数分别经过低通滤波器h(k)和高通滤波器g(k)得到。
注意式(13,14)中的离散性质,说明h(k)和g(k)仅仅是一个固定的离散的数列,这个数列h(k)就是传说中神奇的尺度向量(scaling vector)。
式(13,14)的过程如下图所示:
图(5)
也就是说,式(13,14)中,每当k增加一个单位时,滤波器h和g将平移两个单位。
由于λ是上一级低通滤波器的输出系数,所以也只在一定范围内非零。
这样,h或g将以两倍的速度移出λ的非零取值区。
换句话说,式(13,14)为对上一级j下λ的间隔采样。
计算后j-1级的λ和γ将只有j级时一半数量的采样点。
如果我们把采样后的信号f(i Δt)看作最高一级j下的尺度函数系数λj,当我们已经得到函数h(k)和g(k)时,就可以利用式(13,14)对信号进行分解,迭代进行,于是得到不同频率下的小波系数。
这样,知道了h(k)和g(k),就不需要明确求出小波函数和尺度函数,离散小波分解过程也变成了简单的求内积(式(13,14))。
这就是Mallat的鱼骨算法(因为这个过程,层层分解,像图(3)右鱼骨头一样),一个快速小波系数分解计算过程。
于是,到这里,我们的第三个问题也得到了解答。
总结一下,从CWT到DWT的三个问题,
第1个问题:维数冗余,通过紧支二进小波解决。
第2个问题:无限延伸,通过引入尺度函数解决。
第3个问题:快速算法,通过多分辨率分析解决。
这样,在计算DWT时,我们就有两种方法:
1,构造出尺度函数和所有小波函数,通过式(8)直接计算。
2,设计低通滤波器h(k),利用鱼骨型算法迭代计算。
其中方法2中,只需要设计出一个低通滤波器h(k),而方法1需要计算出所有小波函数。
而且方法2计算过程明确,所以,实际中常常应用方法2来计算离散小波变换。
到这里,我们必须讲讲以上两种计算方法的联系,尤其是尺度向量h(k)和尺度函数的关系。
事实上,如我们可以预料到的,知道尺度向量h(k),我们就可以
根据h(k)构造出g(k),同时也也能构造出尺度函数φ(x),并利用尺度函数构造出基本小波ψ(x)。
同样地,如果已经知道尺度函数φ(x),也可以构造出其他所有信息。
这里顺便说一下,小波函数通常叫做小波族,英文能更好的反映这个称呼,叫做wavelet family。
其中,尺度函数叫做father wavelet,基本小波叫做mother wavelet。
正说明了这种,通过尺度函数和基本小波,可以得到所有小波函数序列的关系。
接下来,我们先说说由尺度向量构造尺度函数的过程。
这里,需要引入哈尔尺度函数(Haar Scaling Function)。
Haar尺度函数是数学家Haar在1910年提出来的,非常简单,是一个[0,1]区间内为1,区间外取0的函数,表达如下
(15)
而尺度函数应该满足三个条件:
如果尺度向量也满足一定约束条件时[4](page 216),假设φ0(x)是Haar尺度函数,则令
(16)
当n→∞时,φn(x)收敛至φ(x),且φ(x)满足尺度函数条件,也就是说,该极限值就是我们要求的尺度函数。
一个很好的例子就是鼎鼎大名的Daubechies小波。
她计算出的一组尺度向量值为:
利用这组尺度向量值,并根据上述所讲的过程求尺度函数,得到Daubechies 尺度函数,如图(貌似图中的是Dau 12,不是Dau 4,不过Dau 4图也类似):
图(6)
同时,也可以根据尺度函数h(k)计算出g(k)。
(17)
具体推导过程,参看[2]。
而有了尺度函数和高通滤波器g(k)后,就可以计算出基本小波。
从图(4)中看出,基本小波是由尺度函数频谱向右平移整个尺度函数带宽得到。
计算结果:
(18)
这样就得到了所有信息。
从这里以及以上Daubechies 的尺度向量可以看出,仅仅区区几个值(Daubechies 4中为4个值),我们就可以算出尺度函数,小波函数等所有信息。
这就是为什么之前说尺度向量是传说中神奇的尺度向量。
当已知尺度函数,求尺度向量时,可根据下式进行[2]:
(19)
到这里,我觉得我想讲的已经都解释清楚了。
关于一些公式的推导,本文中没有涉及,若感兴趣需要自己去资料里找。