小波变换应用(WT适用场合及其在R中的实现过程)

傅里叶变换如何工作的?

傅里叶变换(FFT)将一个信号从时域转换成频域。频谱中的峰值表示信号中出现次数最多的频率。频谱中锋越大,越尖,表示在信号中出现越普遍。频谱中峰值的位置(频率值)和高度(振幅)可以作为分类器的输入,如随机森林或梯度增强。

使用傅里叶变换能够很好的解决问题,这是因为信号不依赖时间变化而变化的;如果信号只频域信息,信号在任何时间段种都是相同的。这样越是非平稳信号,那傅里叶变换的结果越是糟糕。遗憾的是,生活中常见的大部分信号都是这种非平稳信号;例如ECG信号心电图eleCTrocardiogram,股票信息,设备或传感器的数据,更有效的方法就是使用小波变换来代替傅里叶变换。

小波变换应用(WT适用场合及其在R中的实现过程)(1)

傅里叶变换通过很多不同频率的正弦波叠加,就可以知道一个信号中存在哪些频率?如果信号和某个频率的正弦波之间的点积导致振幅很大,这意味着这两个信号之间有很多重叠,信号包含这个特定的频率。当然,这是因为点积是衡量两个向量/信号重叠程度的指标。

傅里叶变换的局限

傅里叶变换对频域变换有很好的解决方式,但是不能解决时域变换问题。这意味着,它能告诉我们信号包含哪些频段,但是不能告诉我们这个频段是在信号的哪个时间段出现的?

小波变换应用(WT适用场合及其在R中的实现过程)(2)

第一行的两个图表示的是一个信号包含4种频率,不因时间的变化而变化;第二行的两个图表示的是信号中的4种频率根据时间依次出现。图1中左上角的信号中包含4种不同的频率(4、30、60、90Hz) ,所有频率在这段时间内全部出现,右上角是这个信号的频谱。在下面的两个图中,同样能看到这4种频率,但是第一种出现在1/4时间上,之后的频率依次出现;并且能在右下部看到这个信号的频谱图。

值得注意的是两个频谱图都包含了4峰值,但是他们都不能明确表示出在信号的哪个时间段出现的?这就使得傅里叶变换不能区分这两个信号。

傅里叶变换处理非平稳信号有天生缺陷。它只能获取 一段信号总体上包含哪些频率的成分,但是 对各成分出现的时刻并无所知;因此时域相差很大的两个信号,可能频谱图一样。然而平稳信号大多是人为制造出来的,自然界的大量信号几乎都是非平稳的,所以在比如生物医学信号分析等领域的论文中,基本看不到单纯傅里叶变换这样naive的方法。

小波变换应用(WT适用场合及其在R中的实现过程)(3)

上图所示的是一个正常人的事件相关电位。对于这样的非平稳信号,只知道包含哪些频率成分是不够的,我们还想知道各个成分出现的时间。知道信号频率随时间变化的情况,各个时刻的瞬时频率及其幅值——这也就是时频分析。

傅里叶变换还存在吉布斯效应,即用无限长的三角函数怎么也拟合不好突变信号。

小波变换应用(WT适用场合及其在R中的实现过程)(4)


短时傅里叶变换的特点及适用条件

短时傅里叶变换(Short-time Fourier Transform, STFT),在频域和时域均具有中等分辨率,属于一种加窗的傅立叶变换。

小波变换应用(WT适用场合及其在R中的实现过程)(5)

时域上分成一段一段做FFT,不就知道频率成分随着时间的变化情况了吗!用这样的方法,可以得到一个信号的时频图了:

小波变换应用(WT适用场合及其在R中的实现过程)(6)

使用STFT存在一个问题,我们应该用多宽的窗函数?窗太宽太窄都有问题:窗太窄,窗内的信号太短,会导致频率分析不够精准,频率分辨率差。窗太宽,时域上又不够精细,时间分辨率低。

为什么不采用可变窗的STFT呢,是因为这样做冗余会太严重, STFT做不到正交化,这也是它的一大缺陷。

小波变换是一种权衡,在与时间相关的特征数据上,它在时域具有高分辨率;在与频率相关的特征数据上,它在频域具有高分辨率。具体的特点:

  • 对于小频率值,频域分辨率高,时域分辨率低。
  • 对于大频率值,频域分辨率低,时域分辨率高。
  • 傅立叶变换就相当于: 你只能在远距离观察油画;
  • 加窗傅立叶变换相当于: 你只能在固定的距离观察油画;
  • 而小波变换相当于,你可以在任意的距离观察油画。

小波变换应用(WT适用场合及其在R中的实现过程)(7)

小波变换是如何工作的?

傅里叶变换是通过一系列频率不同的正弦波来拟合信号;也就是说,信号是通过正弦波的线性组合来表示的。

小波变换的概念是由法国从事石油信号处理的工程师J.Morlet在1974年首先提出的,通过物理的直观和信号处理的实际经验的需要建立了反演公式;小波函数源于多分辨分析,其基本思想是将函数f(t)表示为一系列逐次逼近表达式, 其中每一个都是f(t)动经过平滑后的形式,它们分别对应不同的分辨率。多分辨分析又称多尺度分析,是建立在函数空间概念基础上的理论,其思想的形成来源于工程。

小波变换的出发点和STFT还是不同的。 STFT是给信号加窗,分段做FFT;而小波直接把傅里叶变换的基给换了——将 无限长的三角函数基换成了 有限长的会衰减的小波基。这样 不仅能够获取频率,还可以 定位到时间了。

小波变换使用一系列称为小波的函数,每个函数具有不同的尺度;小波这个词的意思就是很小的波。小波中的正弦波是不限延长的,小波则是在一个时间点上的局域波。

小波变换应用(WT适用场合及其在R中的实现过程)(8)

既然小波是在时间上局域的,就可以叠加不同时间段的小波;从信号开始的时间点慢慢向信号结尾移动小波,这个方法也就是所说的卷积。在对原始(母)小波做完这些,就可以改变小波尺寸,使得它变大然后继续上面的操作。

小波变换应用(WT适用场合及其在R中的实现过程)(9)

由上图可知,经过小波变换信号从1维变成2维;这种小波变换的2维输出是以分布图的形式通过时间尺度变换表示的信号。

那么这个scale是什么维度?由于傅里叶变换保留了频率项,所以小波变换通常用尺度来表示;这就是为什么二维分布图有时间和尺度。对于那些认为频率比尺度更直观的人来说,用这个方程可以把尺度转换成伪频率

越是高的比例(长的小波)对应一个较小的频率,因此在时间维度通过修改小波的比例因子可以在频域分析更小的频率(获得更高的分辨率,但是包含的频率信息就少了)。反之一样,通过使用更小的比例因子可以在时间维度上获得更多的细节,所以尺度基本上是频率的倒数。

然而衰减的小波对于突变信号的处理就不一样了:

小波变换应用(WT适用场合及其在R中的实现过程)(10)


小波变换的特点及适用场合

傅里叶变换和小波变换的另一个不同的地方是小波中有很多不同的族。小波族之间是不同的,因为每个小波族在小波的紧凑性和平滑性方面都做出了不同的权衡。这意味着我们可以选择一个特定的小波族,它最适合在信号中寻找的特征。

不同种类的小波有着不同的形状,平滑的和紧密的,用于不同的目的;由于小波只需要满足两个数学条件,因此很容易产生一种新的小波;这两个数学条件就是所谓的归一化约束和正交化约束

小波变换应用(WT适用场合及其在R中的实现过程)(11)

第一行包含四个离散小波,第二行包含四个连续小波

小波变换有两种不同的风格——连续小波变换和离散小波变换

小波变换应用(WT适用场合及其在R中的实现过程)(12)

el-Nino数据(最上),经过傅里叶变换后(中间),连续小波变换(底部)

在标度图中可以看到,在1920年之前有很多波动,而在1960 - 1990年之间没有那么多波动;而且随着时间的推移,从较短的周期向较长的周期转变。这是一种信号的动态行为可以用小波变换来表示但不能用傅里叶变换来表示。


利用连续小波变换和卷积神经网络对信号进行分类

小波变换对于机器学习的作用很大。小波变换在el-Nino数据集上的应用,它不仅能告诉最大振荡的周期是多少?还能告诉这些振荡是什么时候出现的,什么时候没有?

所以通过看比例图可以区分出一个坏的马达和一个工作的马达,一个健康的人和一个生病的人,一个上楼的人和一个下楼的人等等。但是如果你像我一样懒,您可能不想手动筛选数千个标度图。实现这一过程自动化的一种方法是建立卷积神经网络,该网络可以自动检测每个尺度图所属的类并对其进行分类。

在图像处理方面,小波变换可以对数据分析前对数据进行分层处理;将数据分为高频和低频,由母小波经过平移和尺度伸缩表示原本的数据的高频和低频。

网络信息的数据挖掘方面,小波变换可以进行数据特征的提取。小波分析的应用领域还包括:信号分析、图像处理;量子力学、理论物理;军事电子对抗与武器的智能化;计算机分类与识别;音乐与语言的人工合成;医学成像与诊断;地震勘探数据处理;大型机械的故障诊断等方面;在信号分析方面的滤波、去噪声、压缩、传递等;在图像处理方面的图像压缩、分类、识别与诊断;在医学成像方面的减少B超、CT、核磁共振成像的时间,提高分辨率等。

每个小波变换都会有一个mother wavelet,称之为母小波,同时还有一个father wavelet,就是scaling function。而该小波的basis函数其实就是对这个母小波和父小波缩放和平移形成的。缩放倍数都是2的级数,平移的大小和当前其缩放的程度有关。小波变换的三个特点,就是小波级数是二维的,能定位时域和频域,计算很快。

小波函数是由尺度函数构造的,尺度函数的性质决定了小波函数的性质;尺度函数从滤波器的角度看是低通滤波器,而小波函数是高通滤波器。尺度函数又称为小波父函数,根据双尺度方程,可以由尺度函数生成小波。父小波呢,主要是为了方便做多解析度分析(multiresolution analysis, MRA);而小波变换的强大,就体现在这个多解析度上。

进行信号处理时,先用尺度函数对信号进行分解;尺度函数的频带与待分析信号的频带相同,然后将逼近函数分别在尺度空间和小波空间中进行分解.就得到了信号的低频粗略部分和高频细节部分.此时新的尺度函数频带是原信号频带的一半.小波函数的频带是另一半(高频部分).由此实现了对原信号的按频带分解!

时间序列分类是根据已标注的时间序列建立一个分类模型,然后使用分类模型预测未标记时间序列的类别;从时间序列中抽取出新特征有助于提高分类模型的性能。特征提取技术有奇异值分解(SVD)、离散傅立叶变换(DFT)、离散小波变换(DWT)、分段积累近似法(PAA)、连续重要点(PIP)、分段线性表示,以及符号表示。

小波变换必须掌握的若干概念

真正理解透小波变换,至少要知道有 “尺度函数”的存在,它是构造“小波函数”的关键,并且是它和小波函数一起才构成了小波多分辨率分析,理解了它才有可能利用小波做一些数字信号处理;另外,还要理解离散小波变换DWT、正交小波变换、二维小波变换、小波包

为什么说小波能实现正交化是优势?简单说,如果采用正交基,变换域系数会没有冗余信息,变换前后的信号能量相等,等于是用最少的数据表达最大的信息量,利于数值压缩等领域。JPEG2000压缩就是用正交小波变换。比如典型的正交基:二维笛卡尔坐标系的(1,0)、(0,1),用它们表达一个信号显然非常高效,计算简单。而如果用三个互成120°的向量表达,则会有信息冗余,有重复表达。但是并不意味着正交一定优于不正交。比如如果是做图像增强,有时候反而希望能有一些冗余信息,更利于对噪声的抑制和对某些特征的增强。

小波变换的不足

A.作为图像处理方法,和多尺度几何分析方法(超小波)比:对于图像这种二维信号的话,二维小波变换只能沿2个方向进行,对图像中点的信息表达还可以,但是对线就比较差。而图像中最重要的信息恰是那些边缘线,这时候ridgelet(脊波), curvelet(曲波)等多尺度几何分析方法就更有优势了。

B. 作为时频分析方法,和HHT比:

1998年,Norden E. Huang(黄锷:中国台湾海洋学家)等人提出了经验模态分解方法,并引入了Hilbert谱的概念和Hilbert谱分析的方法,美国国家航空和宇航局(NASA)将这一方法命名为Hilbert-Huang Transform,简称HHT,即希尔伯特-黄变换。

Hilbert变换的表达式实际上就是将原始信号和一个信号做卷积的结果;Hilbert变换可以看成是将原始信号通过一个滤波器,或者一个系统,这个系统的冲击响应为h(t)。HHT在每一时刻得到的只有一个值,而不是像小波之类的得到一系列的值(多尺度分析)。所以我们从其时频分布图上看到的是一条线,而不是一幅图。

小波变换应用(WT适用场合及其在R中的实现过程)(13)

复数信号是完整的,而实信号只是在复平面的实轴上的一个投影;解析信号可以计算包络(瞬时振幅)和瞬时相位。在上图中可以看到,实际上计算的包络就是黑色的线围成的立体图形的边界在实部的投影,而计算这个边的投影也很简单,就是在复平面上的螺旋线中的每一个点的模值,也就是A(t) = sqrt(x^2(t) Hilbert(x(t))^2),而瞬时相位就是虚部(Hilbert变换后的)和实部(原始信号)在某一时间点的比值的arctan,瞬时频率就是它的导数。

相比于HHT等时频分析方法,小波依然没脱离海森堡测不准原理的束缚,某种尺度下,不能在时间和频率上同时具有很高的精度;以及小波是非适应性的,基函数选定了就不改了。


小波变换的R语言实现

小波变换应用(WT适用场合及其在R中的实现过程)(14)

左图为原始数据,右图为小波分析加显著性检测结果。

小波变换应用(WT适用场合及其在R中的实现过程)(15)

install.packages('WaveletComp');library(WaveletComp);series.length=6*128*24

periodic.series {WaveletComp} 计算线性变化周期的(确定性)周期性时间序列;它计算并返回一个指定长度的正弦曲线,它具有给定的初始相位,以及从给定周期长度开始到最后给定长度的线性变化周期;

x1<-periodic.series(start.period=1*24,length=series.length)

x2<-periodic.series(start.period=8*24,length=series.length)

x3<-periodic.series(start.period=32*24,length=series.length)

x4<-periodic.series(start.period=128*24,length=series.length)

x<-x1 x2 x3 x4

plot(ts(x,start=0,frequency=24),type='l',xlab='time(days)',ylab='hourly data',main='周期为1、8、32、128天的系列数据')

小波变换应用(WT适用场合及其在R中的实现过程)(16)

my.data<-data.frame(x=x)

my.wt=analyze.wavelet(my.data,'x',loess.span=0,dt=1/24,dj=1/20,lowerPeriod=1/4,make.pval=T,n.sim=10)

wt.image(my.wt,color.key='interval',legend.params=list(lab='wavelet power levels')) #单个时间序列的小波功率谱图像图.

wt.avg(my.wt,siglvl=0.05,sigcol='red')


包2:waveslim;包3:wavelets

使用哈尔过滤器提取DWT系数,wavelets包用于实现离散小波变换,其中函数dwt(x,filter,n.levels,...)计算DWT系数,x为一个单变量或多变量的时间序列,filter指明所使用的小波过滤器,n.levels指明需要分解的层数,该函数返回一个dwt对象,其中槽W包包含小波系数,V包包含相似系数。原始时间序列还可以使用wavelets包中的idwt()函数通过逆离散小波变换进行重构。

不管是信号压缩,滤波,还是别的方式处理,只要是用小波变换,都按照下面的流程:

1. 选取合适的wavelet function和scaling function,从已有的信号中,反算出系数c和d。

2. 对系数做对应处理

3. 从处理后的系数中重新构建信号。

如何构建自己的scaling function,选取合适的系数集h[k],并由此构建自己的wavelet function。为什么要强调尺度函数和小波函数组成一个标准正交基orthonormal basis呢?计算方便是一方面,还有一个原因是,如果他们满足这个性质,就满足瑞利能量定理,也就是说,信号的能量,可以完全用每个频域里面的展开部分的能量,也就是他们的展开系数表示。

,

免责声明:本文仅代表文章作者的个人观点,与本站无关。其原创性、真实性以及文中陈述文字和内容未经本站证实,对本文以及其中全部或者部分内容文字的真实性、完整性和原创性本站不作任何保证或承诺,请读者仅作参考,并自行核实相关内容。文章投诉邮箱:anhduc.ph@yahoo.com

    分享
    投诉
    首页