主要内容

连续小波分析的实用介绍

此示例演示如何执行和解释连续小波分析。此示例帮助您回答常见问题,例如:连续小波分析和离散小波分析之间的区别是什么?为什么连续小波分析中的频率或尺度是对数间隔的?在哪些类型的信号分析问题中是连续的美国的小波技术特别有用吗?

连续小波分析与离散小波分析

人们经常问的一个问题是:什么是不断的关于连续小波分析?毕竟,你大概对在电脑上做小波分析感兴趣而不是用纸和笔,而在电脑上,没有什么是真正连续的。当这个词连续小波分析用于科学计算环境,意味着一个小波分析技术,具有多个小波的每个速度,或频率加倍,以及在时间之间的偏移是一个样本。这提供了由此产生的连续小波变换(CWT)具有两个在应用中非常有用的属性:

  • 信号的频率内容比以离散小波技术更精细地捕获。

  • CWT具有与每个频带中的原始数据相同的时间分辨率。

此外,在CWT的大多数应用中,复值小波比实值小波更有用。其主要原因是复值小波包含相位信息。看比较两种信号的时变频率含量例如,相位信息非常有用。本教程中的示例专门使用复数小波。

参见[1为小波信号处理的详细处理,包括复值小波的连续小波分析。

每八度的滤镜或声音

通常用来表示每个八度的小波滤波器个数的术语是声音每倍频程.为了说明这一点,构造一个默认的连续小波滤波器组。

fb = cwtfilterbank
fb = cwtfilterbank with properties: VoicesPerOctave: 10 Wavelet: 'Morse' SamplingFrequency: 1 SamplingPeriod: [] PeriodLimits: [] SignalLength: 1024 FrequencyLimits: [] TimeBandwidth: 60 Wavelet parameters: [] Boundary: 'reflection'

如果检查滤波器组的属性,默认情况下,每倍频程有10个小波滤波器(或VoicesPerOctave).绘制小波滤波器的频率响应。

freqz(神奇动物)

图中包含一个轴对象。标题为CWT Filter Bank的axis对象包含71个类型为line的对象。

小波滤波器在连续分析中占有重要的地位常q所有小波滤波器的性质,即它们在频率或带宽上的扩展与它们的中心频率成正比。换句话说,小波滤波器在高频时比在低频时更宽。由于时间和频率之间的互反关系,这意味着高频小波比低频小波更好地在时间上本地化(有更小的时间支持)。万博1manbetx为了看到这一点,从滤波器组中提取中心频率和时域小波。绘制两个小波的实部,一个高频和一个低频,以作比较。

PSI =小波(FB);F = centerFrequencies (fb);图绘图(Psi(9,:))ylabel(“振幅”) yyaxisplot(真实(psi(psi(ex-16,:)))ylabel(“振幅”)轴S1 ='具有中心频率的小波';S2 = [num2str (F (9),' % 1.2 f '”和“num2str (F (end-16),' % 1.2 f '“周期/样本”];标题({S1;S2})包含('样本')传说(“0.25 c / s波”“0.01 c / s波”

图中包含一个轴对象。以中心频率为0.25和0.01 cycles/sample的小波为标题的轴对象包含2个线型对象。这些对象代表0.25 c/s小波,0.01 c/s小波。

在频率响应图中,中心频率不超过0.45个循环/样本。找出滤波器组在一个倍频程或倍频中有多少个小波滤波器。确认号码等于VoicesPerOctave

numFilters10 = nnz(F >= 0.2 & F<= 0.4)
numFilters10 = 10

创建过滤器组时,可以指定不同数量的过滤器。创建一个每倍频程8个声音的滤波器组,并绘制频率响应图。

fb = cwtfilterbank (“VoicesPerOctave”,8); 图Z(fb)

图中包含一个轴对象。标题为CWT Filter Bank的axis对象包含57个类型为line的对象。

确认滤波器组中每个八度有八个滤波器。

F = centerFrequencies (fb);numFilters8 = nnz(F >= 0.2 & F<= 0.4)
numfilters8 = 8.

你可以阅读连续小波分析和离散小波分析之间的更详细的解释连续和离散小波变换

对数间隔中心频率

人们可能会发现,小波分析的一个令人困惑的方面是滤波器的对数间距。

绘制小波滤波器组的中心频率。

次要情节(2,1,1)情节(F) ylabel (“周期/样本”)标题(“小波中心频率”网格)子图(2,1,2)半机(F)网格ylabel (“周期/样本”)标题(“半对数刻度”

图中包含2个轴对象。轴对象1带标题小波中心频率包含类型线的对象。轴对象2具有标题半志量表包含类型线的对象。

图表明,小波中心频率不是线性间隔,这是一般情况下,与其他滤波器组。具体来说,中心频率呈指数递减,高中心频率之间的步长大于低中心频率之间的步长。为什么这对小波有意义?回想一下小波是常数q滤波器,这意味着它们的带宽与中心频率成比例。如果你有一个滤波器组每个滤波器都有相同的带宽,比如短时间傅里叶变换,或者光谱图,滤波器组,那么要保持滤波器之间的频率重叠你需要一个恒定的步长。然而,对于小波,步长应该与频率成正比,就像带宽一样。对于连续小波分析,最常见的间距是基数2^(1/)NV), 在哪里NV是每个octave的过滤器数量,提升到整数的力量。为了比较,专门用于离散小波分析的间距是升高到整数功率的基础2。

参见[2]下表总结了离散小波技术和连续小波技术之间的主要异同点。对于离散小波技术,括号中提供了MATLAB中代表性算法的名称。

时频分析

由于小波同时在时间和频率上进行局部化,因此它们在许多应用中都很有用。对于连续小波分析,最常见的应用领域是时频分析。此外,连续小波变换的以下特性使它对某些类型的信号特别有用。

  • 小波滤波器的恒Q特性,即高频小波的持续时间较短,低频小波的持续时间较长。

  • 连续分析中小波滤波器间的单样本时间偏移

为了理解哪类信号是连续小波分析的好候选,考虑双曲线性调频信号。

负载hyperbolicChirp图绘制(t, hyperbolchirp)轴Xlabel(‘时间’) ylabel (“振幅”)标题('双曲线唧唧'

图中包含一个轴对象。标题为“双曲啁啾”的轴对象包含一个类型为line的对象。

双曲啁啾是 1 5. π 0. 8. - T. 5. π 0. 8. - T. .第一组分从0.1至0.68秒有效,第二个组分从0.1到0.75秒。这些组件的每个瞬间或瞬时频率的频率是 7. 5. 0. 8. - T. 2 2 5. 0. 8. - T. 2 分别循环/秒。这意味着瞬时频率近乎很低 T. = 0. 并随着时间接近0.8秒迅速增加。

获得啁啾信号的连续小波变换,并将连续小波变换与真实的瞬时频率绘制成虚线。

Fs = 1 /意味着(diff (t));[cfs f] = cwt (hyperbolchirp Fs);helperHyperbolicChirpPlot (cfs f t)

图中包含一个轴对象。标题为“幅值比例图”的轴对象包含3个类型为“表面”、“线”的对象。

我们看到CWT显示了一种时频表示,它准确地捕获了双曲啁啾的瞬时频率。现在用短时间傅里叶变换对相同的信号进行分析,短时间傅里叶变换采用恒定带宽滤波器。使用默认的窗口大小和重叠。

pspectrum(hyperbolchirp,2048,的谱图

图中包含一个轴对象。标题为Fres = 41.0679 Hz, Tres = 62.5 ms的轴对象包含一个类型为image的对象。

请注意,当频率较低时,谱图无法区分两个不同的瞬时频率,或者当两个瞬时频率非常接近时,谱图无法区分两个不同的瞬时频率。根据输入长度2048,谱图默认使用128个样本的窗长进行计算。在2048hz的采样率下,这相当于62.5 msec的时间分辨率。基于默认的Kaiser窗口,这导致频率分辨率为41 Hz。这太大了,无法区分t=0附近的瞬时频率。如果我们为了更好地分辨低频而将频率分辨率降低大约1/2会发生什么?

pspectrum(hyperbolchirp,2048,的谱图“FrequencyResolution”, 20岁,...'重叠的',90)

图中包含一个轴对象。标题为Fres = 20.0637 Hz, Tres = 127.9297 ms的axis对象包含一个类型为image的对象。

结果并不好,因为我们在较低频率处获得的任何频率分离都会导致较高的瞬时频率在时间上重叠在一起。谱图是一种强大的时频分析技术,实际上在许多应用中都是最优的,但就像所有的技术一样,它也有局限性。这种特殊类型的信号对固定带宽滤波器提出了主要挑战。理想情况下,您希望在低频率下有长时间响应,频率支持较窄;在高频率下有短时间响应,频率支持较宽。万博1manbetx小波变换恰好提供了这一点,因此在这种情况下非常成功。一个等价的说法是小波变换在较高的频率有更好的时间分辨率,在较低的频率有更好的频率分辨率。

在真实世界的信号中,高频事件通常持续时间较短,而低频事件持续时间较长。例如,加载并绘制人体心电图(ECG)信号。数据以180hz采样。

负载wecgtm = 0:1/180:元素个数(wecg) * 1/180-1/180;wecg情节(tm)网格标题(“人体心电图”)包含(‘时间’) ylabel (“振幅”

图中包含一个Axis对象。标题为Human ECG的Axis对象包含一个line类型的对象。

这些数据由缓慢变化的部分组成,这些部分被高频瞬变信号打断,这些高频瞬变信号代表着通过心脏传播的电脉冲。如果我们对这些脉冲事件在时间上的精确局部化感兴趣,我们希望分析窗口(过滤器)是短的。为了更精确的时间定位,我们愿意牺牲频率分辨率。另一方面,为了识别低频振荡,最好使用较长的分析窗口。

[cfs f] = cwt (wecg, 180);显示亮度图像(tm、f、abs (cfs))包含(‘时间’) ylabel (‘频率(Hz)’)轴xycaxis([0.025 - 0.25])标题(心电数据的CWT

图中包含一个轴对象。心电数据以CWT为标题的轴对象包含一个类型为图像的对象。

CWT显示近30hz和37hz的稳态振荡以及表示心跳的瞬态事件(QRS复合体)。通过这种分析,您可以在相同的时频表示中同时分析这两种现象。

作为慢慢变化组件的高频瞬态的另一个例子,在1993年12月4日至12:00 UT上,将每小时在南极南极上记录的太阳磁场大小的时间序列。1994年5月24日。见[2]第218-220页,获取该数据的完整描述。

负载SolarmfMagnitudes.绘图(SM_dates,SM)轴标题(《太阳磁场强度》)包含('日期') ylabel ('震级'

图中包含一个axes对象。标题为“太阳磁场幅值”的axes对象包含一个line类型的对象。

原始时间数据显示了具有几种冲动事件或冲击波的整体大致线性趋势。当快速的太阳风或快速冠状质量喷射超过太阳风时产生冲击波。在时间序列中,我们看到大约在以下日期发生的显着冲击波结构:2月11日,2月26日,3月10日,4月3日和4月23日。

对数据进行连续小波分析,并将CWT幅度与时间序列数据绘制在同一图上。

sm_dates helperSolarMFDataPlot (sm)

图中包含一个轴对象。轴对象包含2个类型图像的对象,线。

CWT在脉冲事件在时间序列中发生的同时捕获脉冲事件。然而,CWT也揭示了隐藏在时间序列中的数据的低频特征。例如,从1994年2月11日的激波结构之前开始并延伸到激波结构之外,出现了一个低频稳态事件(接近0.04周期/天)。这是由于日冕物质抛射(CME)事件可能发生在不同的太阳区域,与更近的激波同时到达尤利西斯宇宙飞船。

在心电图和太阳磁场量级的例子中,我们有两个来自非常不同的机制的时间序列产生相似的数据:既有短时瞬变信号,也有长时低频振荡信号。用相同的时频表示来描述两者是有利的。

本地化瞬变

时间序列数据中的瞬态事件通常是最具信息丰富的事件之一。瞬态事件可以指示数据生成机制中的突然变化,您希望及时检测和本地化。例子包括机械的故障,传感器问题,经济时序中的财务“冲击”以及许多其他人。例如,考虑以下信号。

rng默认的;dt = 0.001;t = 0: dt: 1.5 dt;addNoise = 0.025 * randn(大小(t));x = cos(2 *π* 150 * t) * (t > = 0.1 & t < 0.5) +罪(2 *π* 200 * t) * 0.7 (t > & t < = 1.2);x = x + addNoise;X ([222 800]) = X ([222 800]) +[-2 2];图;情节(t。* 1000 x);Xlabel(的毫秒);ylabel (“振幅”);

图中包含一个axes对象。axes对象包含一个line类型的对象。

该信号由两个正弦分量组成,在222毫秒和800毫秒时突然开启和关闭,并伴有额外的“缺陷”。获得信号的小波变换,并仅绘制最小尺度,或最高中心频率,小波系数。创建第二个轴并绘制原始时间数据。

cfs = cwt (x, 1000,“爱”);情节(t、abs (cfs (1:)),“线宽”(2) ylabel'震级')set(gca,“克斯蒂克”,[0.1 0.222 0.5 0.7 0.8 1.2]) yyaxis情节(t x,'--') ylabel (“振幅”)包含('秒')举行

图中包含一个axes对象。axes对象包含两个line类型的对象。

最好的CWT系数本地化数据中的所有突发更改。CWT系数显示正弦分量打开和关闭的峰,以及将正弦曲线的缺陷定位在222毫秒和800毫秒。因为CWT系数具有与数据的相同时间分辨率,所以使用CWT精细尺度大小通常是有用的,以便检测数据中的瞬态变化。为了最准确地本地化这些瞬态,使用最精细的尺度(最高频率)系数。

提高时频分析

任何使用滤波器的时频变换,如CWT中的小波,或短时间傅里叶变换中的调制窗口,都必然会在时间和频率上抹去信号的图像。信号能量在时间和频率上的定位的不确定性来自滤波器在时间和频率上的扩展。同步压缩是一种试图通过沿频率轴“压缩”变换来补偿这种模糊的技术。为了用小波来说明这一点,获得由调幅调频(AM-FM)信号和18赫兹正弦信号组成的信号的CWT。调幅-调频信号由公式定义

2 + 余弦 4. π T. 2 π 2 3. 1 T. + 9. 0. 3. π T.

负载多重压缩sig = sig1 + sig2;类(团体、sampfreq“爱”

图中包含一个轴对象。具有标题幅度标题标题尺度标题的轴对象包含3个类型图像,行,区域的对象。

虽然CWT同时显示18 Hz正弦波和AM-FM分量,但我们可以看到18 Hz正弦波在频率上明显分散。现在使用同步压缩来锐化频率估计。

海温,[f] =墓场(sig, sampfreq);图helperSSTPlot (sst、t、f)

图中包含一个轴对象。轴对象包含类型图像的对象。

SyschroSqueezed变换在很大程度上补偿了小波滤波器引入的时间频率的扩展。SynchroSqueezing在时间频率上挤压相对宽的峰值缩小.你可以提取这些脊然后从时频脊中重建单个的分量。

(冰箱,iridge) = wsstridge (sst 5 f,“NumRidges”2);图;轮廓(t、f、abs (sst));网格;标题(“带脊线的同步压缩变换”);Xlabel(的时间(秒));ylabel (“赫兹”); 持有;情节(t,冰箱,'k-'“线宽”2);持有

图形包含一个轴对象。标题为“带脊线的轴对象”的轴对象包含3个等高线、直线类型的对象。

最后根据时频脊重构出分量的近似,并与原始分量进行比较。

iridge xrec = iwsst (sst);helperPlotComponents (xrec sig1 sig2 t)

图中包含4个轴对象。具有标题重建模式的轴对象1包含线条类型的对象。具有标题原始组件的轴对象2包含线条类型的对象。具有标题重建模式的轴对象3包含线条类型的对象。具有标题原始组件的轴对象4包含线条类型的对象。

比较两种信号的时变频率含量

通常你会有两种信号它们在某种程度上是相关的。一个信号可能决定另一个信号的行为,或者由于对两个信号的外部影响,这些信号可能只是相互关联。如果你得到两个信号的CWT,你可以使用这些变换来获得小波相干性,衡量时变相关性的衡量标准。

小波的每一时刻和中心频率的相干性产生一个0到1之间的相干值,量化了两个信号之间的相关性强度。由于使用了复数小波,所以也有相位信息,可以用来推断信号之间的超前滞后关系。以下面两个信号为例。

t = 0:0.001:2;X = cos(2*pi*10*t).*(t>=0.5 & t<1.1)+...cos(2*pi*50*t)。*(t>=0.2&t<1.4)+0.25*randn(大小(t));Y=sin(2*pi*10*t)。*(t>=0.6&t<1.2)+...罪(2 *π* 50 * t)。*(t>= 0.4 & t<1.6)+ 0.35*randn(size(t)));图subplot(2,1,1) plot(t,X) ylabel(“振幅”)子地块(2,1,2)图(t,Y)ylabel(“振幅”)包含(‘时间’

图中包含2个轴对象。axis对象1包含一个类型为line的对象。axis对象2包含一个类型为line的对象。

两个信号都是由两个正弦波(10hz和50hz)组成的白噪声。正弦波的时间支持略有不同。万博1manbetx注意,10赫兹和50赫兹的分量y相对于中相应的组分滞后1/4个周期X.绘制两个信号的小波相干性。

1000年图wcoherence (X, Y,,'Phasedisplaythreshold', 0.7)

图中包含一个轴对象。以小波相干为标题的轴对象包含类型为图像、线、片的157个对象。

小波相干估计捕捉到信号在50和10赫兹附近高度相关,并且只在正确的时间间隔内。图上的暗箭头表示信号之间的相位关系。这里的箭头指向上方,这表明了两个信号之间的90度相位关系。这正是由信号中的cos -sin项所决定的关系。

作为一个真实的例子,考虑厄尔尼诺地区3的数据和从1871年到2003年底的季节性全印度降雨量指数。这些数据每月抽样一次。尼诺3时间序列是在赤道太平洋西经90度至150度,北纬5度至南纬5度记录的每月摄氏度海面温度异常记录。全印度降雨量指数代表了印度的平均降雨量,以毫米为单位,剔除了季节性成分。

负载ninoairdata图subplot(2,1,1) plot(datyear,nino) title(“厄尔尼诺3区——海温异常”) ylabel (“摄氏度”)包含('年')轴次要情节(2,1,2)情节(datayear、空气)轴标题(“全印度降水指数季节化”) ylabel (“毫米”)包含('年'

图中包含2个轴对象。带有标题EL Nino Region 3  -  SST异常的轴对象1包含类型线的对象。轴对象2带有标题叛过的所有印度雨降雨指数包含了类型线的对象。

计算两个时间序列之间的小波相干性。将相位显示阈值设置为0.7。要显示以年为单位的周期,指定采样间隔为一年的1/12。

图wcoherence(尼诺、空气、年(1/12),'Phasedisplaythreshold', 0.7);

图中包含一个轴对象。以小波相干为标题的轴对象包含图像、线、片等类型的61个对象。

这幅图显示了强相干的时间局部区域,其发生周期与典型的2至7年厄尔尼诺周期相对应。该图还显示,在这些时期,两个时间序列之间大约有3/8到1/2个周期延迟。这表明,与南美洲海岸记录的厄尔尼诺相一致的海洋变暖时期,与17000公里以外的印度的降雨量相关,但这种影响大约延迟了1/2个周期(1至3.5年)。

结论

在这个例子中,你学到:

  • 连续小波分析与离散小波分析的区别。

  • CWT的恒Q和单采样时移特性如何允许您同时分析非常不同的信号结构。

  • 一种使用连续小波变换的时频重分配技术。

  • 如何使用CWT比较两个信号中的频率内容。

时频分析有关连续小波分析的更多突出主题和特色示例。

参考文献

马拉特,s.g.。信号处理的小波之旅:稀疏方法第三版阿姆斯特丹 ; 波士顿:爱思唯尔/学术出版社,2009年。

Percival, Donald B.和Andrew T. Walden。时间序列分析的小波方法剑桥:剑桥大学出版社,2000年。https://doi.org/10.1017/CBO9780511841040.

也可以看看

职能

应用程序

相关的话题