草庐IT

如何对时间序列进行小波分析,得出其周期?

哥廷根数学学派 2023-09-20 原文

从信号处理角度进行分析 简单的时间序列直接做各种谱分析(频谱,包络谱,平方包络谱,功率谱,倒谱等等) 比如一些简单的旋转机械振动时间序列信号

​如果频谱不好分析,那可以分析如下图所示的时间序列的时频谱

​给个简单的模拟信号的例子 t = 0:1/2000:1-1/2000; dt = 1/2000; x1 = sin(50*pi*t).*exp(-50*pi*(t-0.2).^2); x2 = sin(50*pi*t).*exp(-100*pi*(t-0.5).^2); x3 = 2*cos(140*pi*t).*exp(-50*pi*(t-0.2).^2); x4 = 2*sin(140*pi*t).*exp(-80*pi*(t-0.8).^2); x = x1+x2+x3+x4; figure; plot(t,x) title('Superimposed Signal')

​其连续小波变换时频谱如下

​一个模拟的轴承内圈故障振动信号,带有明显的周期性

​相应的频谱如下,红色虚线代表故障特征频率及相应的倍频

包络谱如下

​看一下相应的CWT时频谱,很明显能看出冲击性

​还可以试试小波相干与交叉小波分析 小波相干、交叉小波,可以很好地反映两个不同时间序列变化之间的“相关性”。小波相干分析,一般反映序列间周期性“变化趋势”的一致性,但不直接反映变化周期的强度关系;交叉小波分析,一般反映序列间“共有周期”的强度。

此外,如果时频谱线能量发散,时频脊线模糊,还可以试试同步压缩之类的算法 时间序列信号处理系列-基于Python的同步压缩变换 - 哥廷根数学学派的文章 - 知乎 时间序列信号处理系列-基于Python的同步压缩变换 - 知乎 当时间序列信号中噪声较大时,为了有利于周期性分析,不可避免地要进行降噪前处理 比如K-SVD降噪

样条框架降噪

Morlet小波降噪

当待分析的时间序列过于复杂时,那可能要引入模态分解(多分辨分析),比如小波分解,经验模态分解及其变体,变分模态分解,经验小波变换,局部均值分解,辛几何模态分解,各种各样的自适应分解算法 基于小波脊线的时间序列分解

好多同学都对各种模态分解方法的时间序列处理感兴趣,那就随便说一下 实际上,时间序列通常由多个具有物理意义的分量组成,在很多时候,为了更容易的研究信号,我们希望在与原始数据相同的时间尺度上单独研究这些分量中的一个或多个,理想情况下,我们希望这些经MRA分解到的多个分量在物理上是有意义的,可容易解释的。多分辨率分析MRA通常与小波或小波包相关联,但诸如经验模态分解EMD,变分模态分解VMD等模态分解方法也可以构成MRA。 先给一个简单的合成信号,信号以1000Hz的频率采样1秒钟。 Fs = 1e3; t = 0:1/Fs:1-1/Fs; comp1 = cos(2*pi*200*t).*(t>0.7); comp2 = cos(2*pi*60*t).*(t>=0.1 & t<0.3); trend = sin(2*pi*1/2*t); rng default wgnNoise = 0.4*randn(size(t)); x = comp1+comp2+trend+wgnNoise; plot(t,x) xlabel('Seconds') ylabel('Amplitude') title('Synthetic Signal')

该信号由3个主要分量组成:频率为 60Hz的时间局部振荡分量、频率为 200 Hz的时间局部振荡分量和趋势项分量。趋势项分量为正弦曲线,频率为0.5Hz。60Hz的振荡分量发生在 0.1到 0.3 秒之间,而 200Hz的振荡分量发生在 0.7 到 1 秒之间。 但这些分量从时域波形中无法分辨,因此进行频域变换。 xdft = fft(x); N = numel(x); xdft = xdft(1:numel(xdft)/2+1); freq = 0:Fs/N:Fs/2; plot(freq,20*log10(abs(xdft))) xlabel('Cycles/second') ylabel('dB') grid on

从频率中可以更容易地辨别振荡分量的频率,但时间局部性信号却丢失。为了同时定位时间和频率信息,使用连续小波变换进行分析。

从CWT时频谱图中可以看出60Hz和200Hz分量的时间范围,但没有发现趋势项分量。为了分离出信号的分量并单独进行分析,接下来使用多分辨分析,直接在时域中进行相关操作。 多分辨分析通过将信号分成不同分辨率的分量进而缩小分析范围,而提取不同分辨率的信号分量相当于分解数据在不同时间尺度上的变化,或等效地在不同频带上进行分析。首先,采用离散小波变换的变体最大重叠离散小波变换对信号进行多分辨分析,分解层数为8。关于最大重叠离散小波变换的相关内容,请查看如下文献。

最大重叠离散小波变换的8层多分辨分析分解如下:

如果从上向下看,会看到所分解的分量变得越来越平滑,即分量频率越来越低。回想一下,原始信号包含3个主要分量,一个 200 Hz 的高频振荡成分、一个 60 Hz 的低频振荡成分和一个趋势成分,它们都被加性噪声破坏了。 从D2 图中可以看出时间局部化的高频分量被分解出来,而下面的两个图包含较低频率的振荡分量,这是多分辨率分析的一个重要方面,最后S8子图中包含了趋势项分量。 除了小波多分辨分析,经验模态分解 (EMD) 是一种所谓的数据自适应多分辨技术。 EMD 在不使用固定基函数的情况下递归地从数据中提取不同的分辨率成分,关于EMD相关文献浩如烟海,不做赘述了。EMD的多分辨分析分解如下所示:

​虽然MRA分解分量的数目不同,但 EMD MRA和小波 MRA会产生相似的信号波形,在 EMD MRA分解中,高频振荡成分位于第1个本征模态函数中 (IMF1),低频振荡成分主要位于IMF2和IMF3中,IMF6 中的趋势项分量与小波技术提取的趋势分量非常相似。 自适应多分辨分析的另一种技术是变分模态分解 (VMD),VMD 从信号中提取固有模式函数或振荡模式,并不使用固定基函数进行分析。EMD在时域上递归,以逐步提取低频IMF分量,而VMD 首先识别频域中的信号峰值并同时提取所有模式,相关文献如下: Dragomiretskiy, Konstantin, and Dominique Zosso. “Variational Mode Decomposition.” IEEE Transactions on Signal Processing 62, no. 3 (February 2014): 531–44. https://doi.org/10.1109/TSP.2013.2288675. VMD的多分辨分析分解如下所示:

​由上图可知,与小波和EMD类似,VMD将3个分量基本分离了出来。 还有一种数据自适应多分辨分析技术:经验小波变换 (EWT) ,EWT根据分析信号的频率构造 Meye小波进而进行自适应小波,之前写过EWT相关的内容: 经验小波变换在信号处理及轴承故障诊断中的应用 - 哥廷根数学学派的文章 - 知乎https://zhuanlan.zhihu.com/p/53 EWT的多分辨分析分解如下所示:

​与之前的EMD和小波MRA类似,EWT分解出了相关的振荡分量,用于执行分析的滤波器及其通带信息如下:

下面考虑一段神户地震信号,源于1995 年 1 月 16 日在澳大利亚霍巴特的塔斯马尼亚大学记录,从 20:56:51 (GMT) 开始,以 1 秒的间隔持续 51分钟。 figure plot(T,kobe) title('Kobe Earthquake Seismograph') ylabel('Vertical Acceleration (nm/s^2)') xlabel('Time') axis tight grid on

以最大重叠离散小波变换为例,其8层MRA分解如下:

​从D4和D5子图中可以看出初级与延迟次级波分量,地震波中的分量以不同的速度传播,初级波比次级(剪切)波传播的更快。 将信号分解为若干分量的目的通常是去除某些分量以减轻对信号分析的影响,MRA技术的关键是重建原始信号的能力,如下:

​每种方法的最大重建误差约为10^(-12) 或更小,表明它们可以完美地对信号进行重建。在很多研究中我们对趋势项不感兴趣,由于趋势项一般位于最后一个 MRA 分解分量中,因此只需将该分量去除,然后进行重建。

​此外,再删除第1个MRA分解分量(看起来主要是噪声)

​在前面我们将趋势项删除,然而在许多应用中,趋势项可能是我们的主要研究部分,因此可视化几种MRA方法所提取的趋势项分量。

根据以上的分析,小波MRA技术可以更平滑且最准确地提取趋势项,EMD提取了一个平滑的趋势项,但它相对于真实趋势幅度发生了偏移,而 VMD似乎比小波和EMD更偏向于提取振荡分量。 在前面的示例中,强调了多分辨分析在检测数据中的振荡分量和总体趋势中的作用,然而MRA还可以定位和检测信号中的瞬态成分。为了说明这一点,以1947年第一季度至 2011 年第四季度美国实际国内生产总值 (GDP) 数据,垂直的黑线标志着“大缓和”的开始,标志着从 1980 年代中期开始,美国宏观经济波动性减弱的时期,很难从原始数据中辨别出来。

​部分代码链接见评论区,面包多下载,内容包括现代信号处理,机器学习,深度学习,故障诊断等

有关如何对时间序列进行小波分析,得出其周期?的更多相关文章

  1. ruby - 如何使用 Nokogiri 的 xpath 和 at_xpath 方法 - 2

    我正在学习如何使用Nokogiri,根据这段代码我遇到了一些问题:require'rubygems'require'mechanize'post_agent=WWW::Mechanize.newpost_page=post_agent.get('http://www.vbulletin.org/forum/showthread.php?t=230708')puts"\nabsolutepathwithtbodygivesnil"putspost_page.parser.xpath('/html/body/div/div/div/div/div/table/tbody/tr/td/div

  2. ruby - 如何从 ruby​​ 中的字符串运行任意对象方法? - 2

    总的来说,我对ruby​​还比较陌生,我正在为我正在创建的对象编写一些rspec测试用例。许多测试用例都非常基础,我只是想确保正确填充和返回值。我想知道是否有办法使用循环结构来执行此操作。不必为我要测试的每个方法都设置一个assertEquals。例如:describeitem,"TestingtheItem"doit"willhaveanullvaluetostart"doitem=Item.new#HereIcoulddotheitem.name.shouldbe_nil#thenIcoulddoitem.category.shouldbe_nilendend但我想要一些方法来使用

  3. ruby-on-rails - 使用 Ruby on Rails 进行自动化测试 - 最佳实践 - 2

    很好奇,就使用ruby​​onrails自动化单元测试而言,你们正在做什么?您是否创建了一个脚本来在cron中运行rake作业并将结果邮寄给您?git中的预提交Hook?只是手动调用?我完全理解测试,但想知道在错误发生之前捕获错误的最佳实践是什么。让我们理所当然地认为测试本身是完美无缺的,并且可以正常工作。下一步是什么以确保他们在正确的时间将可能有害的结果传达给您? 最佳答案 不确定您到底想听什么,但是有几个级别的自动代码库控制:在处理某项功能时,您可以使用类似autotest的内容获得关于哪些有效,哪些无效的即时反馈。要确保您的提

  4. python - 如何使用 Ruby 或 Python 创建一系列高音调和低音调的蜂鸣声? - 2

    关闭。这个问题是opinion-based.它目前不接受答案。想要改进这个问题?更新问题,以便editingthispost可以用事实和引用来回答它.关闭4年前。Improvethisquestion我想在固定时间创建一系列低音和高音调的哔哔声。例如:在150毫秒时发出高音调的蜂鸣声在151毫秒时发出低音调的蜂鸣声200毫秒时发出低音调的蜂鸣声250毫秒的高音调蜂鸣声有没有办法在Ruby或Python中做到这一点?我真的不在乎输出编码是什么(.wav、.mp3、.ogg等等),但我确实想创建一个输出文件。

  5. ruby-on-rails - 如何验证 update_all 是否实际在 Rails 中更新 - 2

    给定这段代码defcreate@upgrades=User.update_all(["role=?","upgraded"],:id=>params[:upgrade])redirect_toadmin_upgrades_path,:notice=>"Successfullyupgradeduser."end我如何在该操作中实际验证它们是否已保存或未重定向到适当的页面和消息? 最佳答案 在Rails3中,update_all不返回任何有意义的信息,除了已更新的记录数(这可能取决于您的DBMS是否返回该信息)。http://ar.ru

  6. ruby-on-rails - 按天对 Mongoid 对象进行分组 - 2

    在控制台中反复尝试之后,我想到了这种方法,可以按发生日期对类似activerecord的(Mongoid)对象进行分组。我不确定这是完成此任务的最佳方法,但它确实有效。有没有人有更好的建议,或者这是一个很好的方法?#eventsisanarrayofactiverecord-likeobjectsthatincludeatimeattributeevents.map{|event|#converteventsarrayintoanarrayofhasheswiththedayofthemonthandtheevent{:number=>event.time.day,:event=>ev

  7. ruby-on-rails - 'compass watch' 是如何工作的/它是如何与 rails 一起使用的 - 2

    我在我的项目目录中完成了compasscreate.和compassinitrails。几个问题:我已将我的.sass文件放在public/stylesheets中。这是放置它们的正确位置吗?当我运行compasswatch时,它不会自动编译这些.sass文件。我必须手动指定文件:compasswatchpublic/stylesheets/myfile.sass等。如何让它自动运行?文件ie.css、print.css和screen.css已放在stylesheets/compiled。如何在编译后不让它们重新出现的情况下删除它们?我自己编译的.sass文件编译成compiled/t

  8. ruby - 如何将脚本文件的末尾读取为数据文件(Perl 或任何其他语言) - 2

    我正在寻找执行以下操作的正确语法(在Perl、Shell或Ruby中):#variabletoaccessthedatalinesappendedasafileEND_OF_SCRIPT_MARKERrawdatastartshereanditcontinues. 最佳答案 Perl用__DATA__做这个:#!/usr/bin/perlusestrict;usewarnings;while(){print;}__DATA__Texttoprintgoeshere 关于ruby-如何将脚

  9. ruby - 如何指定 Rack 处理程序 - 2

    Rackup通过Rack的默认处理程序成功运行任何Rack应用程序。例如:classRackAppdefcall(environment)['200',{'Content-Type'=>'text/html'},["Helloworld"]]endendrunRackApp.new但是当最后一行更改为使用Rack的内置CGI处理程序时,rackup给出“NoMethodErrorat/undefinedmethod`call'fornil:NilClass”:Rack::Handler::CGI.runRackApp.newRack的其他内置处理程序也提出了同样的反对意见。例如Rack

  10. ruby - 如何每月在 Heroku 运行一次 Scheduler 插件? - 2

    在选择我想要运行操作的频率时,唯一的选项是“每天”、“每小时”和“每10分钟”。谢谢!我想为我的Rails3.1应用程序运行调度程序。 最佳答案 这不是一个优雅的解决方案,但您可以安排它每天运行,并在实际开始工作之前检查日期是否为当月的第一天。 关于ruby-如何每月在Heroku运行一次Scheduler插件?,我们在StackOverflow上找到一个类似的问题: https://stackoverflow.com/questions/8692687/

随机推荐