
简介本资源是一套基于Matlab实现的小波相干性Wavelet Coherence分析完整代码包面向本科及硕士阶段的信号处理、地球物理、气候时序分析等方向的学习者与科研人员用于量化两组非平稳时间序列在时频域内的协同变化特征。压缩包共含109个文件主体为35个Matlab函数.m与15篇Markdown格式说明文档.md辅以12个文本参数配置.txt、11个示例数据.sample及11张结果可视化图.png整体大小3.08MB结构清晰、模块分明便于理解算法原理与复现实验流程。已有140人学习下载内含Matlab 2014a/2019a双版本兼容代码、详细注释、运行截图及典型数据集如sst_nino3.dat支持开箱即用对运行异常提供常见排错指引适合作为课程设计、毕业论文或科研预研的可靠技术支撑。1. 小波相干性到底是什么为什么值得花时间研究收到这份wavelet-coherence matlab代码.zip的时候我第一反应是这哥们儿要么在做信号处理相关的课题要么就是被导师逼着分析两串时间序列之间的关系。不管你是哪种情况我都能确定一件事你正在研究的东西研究对了方向但前面的路有不少坑在等着你。小波相干性Wavelet Coherence简称 WTC是时频分析领域里非常实用的一项技术。一句话说透它它能告诉你两个信号在哪个时间段、哪个频率带上存在稳定的相关性以及它们的相位关系是怎样的。普通的相关性分析只能给你一个笼统的数字比如“相关系数0.8”但小波相干性可以告诉你这个0.8是发生在第5秒到第10秒之间的高频段还是发生在第30秒之后整个频谱都高度同步。这种精细程度对很多实际问题来说是决定性的。小波相干性最典型的应用场景包括地球物理与气象研究分析海表温度与气压指数、降水量与季风活动之间的关系在哪个年代际尺度上耦合最强。神经科学分析脑电信号不同通道之间的功能连接比如睡眠纺锤波频率12-16Hz上的同步性。金融时间序列研究股票指数与宏观经济变量在不同周期短期、中期、长期上的联动效应。机械故障诊断比较振动信号与噪声信号定位故障特征频率判断故障传播路径。这套MATLAB代码能直接帮你省掉从头造轮子的时间。你拿到手之后只要把数据整理成两个等长的时间序列跑一遍示例脚本就能得到一张典型的“小波相干谱图”——横轴是时间纵轴是频率或周期颜色代表相干性强弱箭头代表相位差。这张图背后承载的信息量远超传统的互相关分析。适合谁来参考我已经默认拿到这份代码的人至少会基础的MATLAB操作会读矩阵会跑脚本。如果你只是刚接触MATLAB连plot都写不利索那建议你先把数据整理好再找一段现成的示例数据跑通之后再去逐步改参数。下面我详细拆解这套代码的核心逻辑和实际操作中会遇到的问题。2. 小波相干性背后的数学原理以及MATLAB实现的核心思路2.1 从普通相干性到小波相干为什么傅里叶不够用我见过很多做信号处理的人一上来就想着用互功率谱密度算相干函数相干系数Cxy(f)公式是用两个信号的互功率谱密度除以各自自功率谱密度的平方根。这个方法的局限非常明显它假设信号是平稳的也就是统计特性不随时间变化。现实世界里的信号几乎没有平稳的——脑电信号会有事件相关电位振动信号会因为转速变化而改变频率成分气象信号更是有各种尺度的周期叠加。小波变换的出现解决的就是这个问题。它用一组伸缩平移的小波基函数去匹配信号的局部特征相当于在时间-频率平面上铺了一张网每个网格点都代表时间位置和频率尺度上的局部能量。小波相干性的本质就是在小波变换的基础上对两个信号的互小波谱进行归一化得到0到1之间的相干值。如果你问为什么不用短时傅里叶变换STFT做同样的事答案在于窗函数。STFT的窗口宽度固定一旦选定时间分辨率和频率分辨率就锁死了也就是所谓的不确定性原理限制。而小波变换的窗口是可变的——分析高频成分时窗口自动变窄获得更好的时间精度分析低频成分时窗口自动变宽获得更好的频率精度。这个自适应特性在处理宽频带非平稳信号时优势极其明显。2.2 小波相干性的计算公式拆解直接看公式可能有点劝退但你花十分钟理解它后面调参数心里就有底了。小波相干性R²(s, τ)的定义是R²(s, τ) |S(s⁻¹ W_xy(s, τ))|² / ( S(s⁻¹ |W_x(s, τ)|²) × S(s⁻¹ |W_y(s, τ)|²) )其中s 是尺度参数对应于频率尺度越大频率越低。τ 是时间平移参数。W_x(s, τ) 和 W_y(s, τ) 分别是两个信号x(t)和y(t)的连续小波变换结果。W_xy(s, τ) 是互小波谱Cross-Wavelet Spectrum等于 W_x 乘以 W_y 的共轭。S 表示平滑算子。这个公式要表达的意思其实和普通相干系数类似——互谱的模平方除以各自功率谱的乘积。但关键在于那个平滑算子 S它必须在时间维和尺度维上分别做平滑否则分子分母会出现零的比值导致任何两个信号算出来的相干性都恒等于1那这指标就没意义了。套用我们行业里常用的说法平滑窗口就是小波相干的“灵魂参数”一不对结论全崩。2.3 MATLAB实现时的关键参数选择不管是你自己写的代码还是基于包里的现成脚本以下几个参数必须搞清楚小波基的选择。绝大多数小波相干分析用的都是复数Morlet小波核心参数ω₀中心频率默认取6。为什么是6因为当ω₀6时Morlet小波的时间分辨率与频率分辨率乘积接近最优而且小波函数的形状接近高斯包络下的正弦波物理意义清晰。如果你把ω₀调小时间分辨率变好但频率分辨率变差反过来ω₀调大频率分辨率变好但时间分辨率变差。我建议新手不要乱动这个参数等把整套流程跑通了再根据具体信号特征去尝试。尺度频率范围的设置。代码里通常会用一个对数分布的尺度向量最常见的是 s s₀ × 2^(j × dj)其中 s₀ 是最小尺度dj 是尺度间隔j 从0到J。dj一般取0.1或者0.125对应频率轴大概每倍频程有几条线间隔太密计算量暴增间隔太疏图看起来像马赛克。J的选择则决定了最低频率一般用最大频率除以最小频率的对数关系推导出来代码里通常已经有默认计算不用手动写死。平滑窗口参数。这一点我要重点强调。平滑操作分为时间维和尺度维。时间维的平滑窗口宽度一般取与尺度成正比的高斯窗这样低频段的平滑范围大、高频段平滑范围小符合信号处理的直觉。尺度维的平滑通常用Boxcar窗口宽度在0.6左右相对合适。你的代码包里如果没有暴露这个参数那大概率写死了如果你想微调要找到对应的函数把平滑窗口宽度改成一个输入变量方便后期做参数敏感性分析。显著性检验的设定。很多代码会在图上叠加细黑线表示该区域的相干性通过了95%置信水平检验。显著性检验一般基于蒙特卡洛模拟生成大量红噪声或白噪声序列对计算它们的小波相干性分布取95%分位数作为阈值。模拟次数太少比如小于100次置信线位置不稳模拟次数太多计算时间成倍增长。我的实测经验是300次是一个性价比不错的值既能得到稳定的噪声背景分布又不会让你在笔记本上等到怀疑人生。3. 实操过程如何跑通这套小波相干MATLAB代码3.1 先看代码包的文件结构拿到zip包之后第一步不是直接双击运行而是先看目录结构。通常一份完整的小波相干代码至少要包含以下几类文件主脚本demo.m或example.m跑通全流程的入口包含数据生成或数据加载、参数设置、函数调用、出图命令。核心函数wave_coherence.m或类似名称实现小波相干计算的核心逻辑。辅助函数小波变换计算、平滑函数、显著性检验函数、相位角计算函数。测试数据文件.mat或.csv方便你直接看到预期效果。用MATLAB打开主脚本按F5之前仔细看一下文件路径设置。我遇到过太多人报错“未定义函数或变量”最后发现是当前文件夹不对或者数据文件不在搜索路径里。用which 函数名命令可以快速确认函数是否被MATLAB识别。3.2 数据准备与导入这一步做不好后面全是白费代码的核心输入是两组等长的一维时间序列。我刚拿到这类代码时会先做一个“造假数据”的验证实验生成一段采样率1000Hz、时长5秒的模拟信号第一段包含一个10Hz的振荡分量第二段前半部分包含10Hz的振荡但与第一段相位差90°后半部分换成20Hz的振荡。然后用这份代码跑一遍看相干谱图能不能在预期的时间-频率位置出现高相干区以及箭头方向是否和预设相位差一致。如果你用已知答案的数据来验证代码后面分析真实数据时的信心就完全不一样了。如果你的数据时间宽度不同采样率不同需要先处理成相同的长度和采样频率。MATLAB中常用的手段是resample函数或者interp1插值。注意插值会引入虚假的低频成分千万别在插值之后再去做小波分析否则低频端出现的高相干完全可能是插值产生的假象。最好就是裁剪到相同时间段然后统一采样率。数据里如果有NaN值必须先处理。MATLAB的快速傅里叶变换遇到NaN直接罢工小波变换虽然不一定会报错但计算出来的相干谱在NaN附近可能出现条带状异常。常用的处理方法是线性插值填充或者把NaN对应的数据段直接标注出来在解读结果时剔除。3.3 主函数调用与关键参数调整以我常用的代码路径为例核心调用格式大致是% 输入数据x, y 两个行向量或列向量dt为采样间隔秒 % 其中dt 1/fsfs为采样率 dt 1/1000; % 核心调用输出wtc为小波相干矩阵period为周期向量coi为锥形影响区域边界 [wtc, period, coi, phase] wavelet_coherence(x, y, dt); % 绘制相干谱图 figure; imagesc(t, log2(period), abs(wtc).^2);这里我加了个.^2代表取相干系数R的平方——很多论文里画的是R²你自己出图的时候要明确画的是哪个量别混着用。phase输出是每一点的相位差单位是弧度表示y相对于x在对应时频点上的相位领先/滞后。这一步我会顺手做一次数据长度上的检查。小波变换要求数据长度至少能够容纳最低频率的几个完整周期否则低频端的边界效应会非常夸张整个图的下半部分基本不能看。经验公式是数据长度至少是最低周期的4到8倍。如果你的数据只有100个点却想分析0.01Hz的周期成分那结果没有任何可信度。3.4 绘图细节与图形解读技巧出图之后第一眼看到的是下面这样的信息分层横轴是时间纵轴是周期注意小波分析里习惯用对数纵轴标注周期值而不是频率因为这样低频部分能看得更清楚。颜色从深蓝到深红表示相干性从0到1。白色或浅色区域表示相干性低。图中的细黑线表示95%显著性水平黑线内部的区域才是置信的“有效强相干区”。颜色的设置很有讲究。MATLAB自带的jet色图虽然好看但在低频端容易出现色带过大、层次不分明的问题。我一般用parula或viridis这类感知均匀的色图避免颜色深浅误导强弱判断。由于实际信号的标准谱图通常需要更直观地区分低相干和高相干我会手动调整caxis的范围比如默认的相干性由caxis([0 1])控制但如果你的数据普遍偏低可以缩到[0.4 1]高亮有效区。箭头方向是相位信息——右向箭头表示y与x同相左向箭头表示反相180°上箭头表示y领先x90°也就是x滞后y90°等等。在解读箭头时千万注意相位差的解释周期依赖在同一个时间点如果你看的是高频成分和低频成分它们的箭头方向可能截然不同代表着不同频段上的领先-滞后关系完全不同。4. 常见问题与排查技巧实录4.1 图形上出现整条或整块的“红色珊瑚礁”是怎么回事如果你跑出来的相干谱图整体通红几乎所有区域都大于0.9这大概率不是信号本身真的那么强相关而是平滑参数设置出了问题。我遇到过的场景往往是把时间维平滑窗口设成了常数而不是与尺度成比例结果高频段的平滑窗口相对过大相当于把相邻时频点的值糊成一团相干性被“糊”高了。另外还有一种情况数据本身包含趋势项或均值不为零。小波变换对常数偏移非常敏感低频端会出现一条高相干宽条带。处理办法是对数据做去均值必要时做趋势去除detrend函数甚至可以加一个带通滤波预处理把无关的极低频漂移滤掉。4.2 低频段整片被锥形阴影覆盖是不是没法用了锥形影响区域Cone of InfluenceCOI在图上以白色半透明区域标出表示该区域内的值受到边界效应污染严重不能当作真实值来解读。很多新手拿到图之后发现低频段一大半都在COI里顿时觉得分析没法做了。我的处理习惯是低频端的结论宁可保守也不要强行解读。如果你关心的频段正好落在COI内说明数据长度不足以在该频率上给出可靠结论这种情况下你应该考虑采集更长的数据而不是尝试用代码去修补边界效应。如果只是为了参考趋势COI内部还是可以看的但任何结论都要在论文或报告里明确标注置信度存疑。4.3 计算速度慢得离谱怎么优化小波相干计算涉及两个信号各自的小波变换、互谱计算、两次平滑、显著性检验的蒙特卡洛模拟数据长度一上去超过几万个点耗时骤增。几个实操优化手段如果只是为了快速预览效果先把显著性检验关掉或把模拟次数从300降到50出图确认形态之后再跑正式版。适当增加尺度间隔dj从0.1改成0.2频率轴上的分辨率降低一半但计算量几乎减半。用tic/toc记录各段耗时定位瓶颈函数。如果瓶颈在蒙特卡洛循环考虑把它改成parfor并行循环 —— 前提是你装了Parallel Computing Toolbox。不少代码包为了通用性在函数内部绘制了大量辅助图。正式分析时关闭多余图窗只保留核心相干谱能省下不少时间。4.4 箭头方向乱糟糟的没有规律性怎么解读箭头方向看似随机实际上它反映的是相位差在时频面上的局部变化。如果你预期某个频段上两个信号有稳定的相位差比如y始终领先x约90°但图上箭头忽左忽右先别急着下“无相关性”的结论。检查一下这个频段是否在COI内以及该区域的相干性是否达到显著水平如果相关性本身就低于0.5相位箭头参考意义有限。如果想提取特定区域的平均相位差不要直接对箭头角度做算术平均因为角度有周期性179°和-179°的平均值应该是180°而不是0°。正确的做法是先把所有相位差转换为单位复数exp(1i*angle)求复数平均之后再转换回角度这样才能得到正确的循环平均相位。4.5 工具箱版本不同导致的兼容性问题MATLAB官方自R2016b开始内置了wcoherence函数需要Wavelet Toolbox它可以直接替代第三方的自写代码而且绘图样式更规范。如果你用的是官方工具箱版的代码请注意wcoherence的输入变量和输出格式跟Grinsted等第三方工具包不完全一致——尤其是phase输出的单位是弧度还是角度、period列向量还是行向量这两个最容易踩坑。如果你拿到的zip包是第三方代码比如经典的Grinsted工具包在较新版本的MATLABR2022b及以后上可能会出现不兼容警告常见原因是老的绘图函数被官方废弃了。解决办法是把报错行里的plot、imagesc等老旧语法改成新语法或者干脆转用官方wcoherence。5. 进阶扩展从双变量走向偏小波相干性5.1 为什么需要偏小波相干双变量小波相干分析只能告诉你x和y之间在何时何频上相关但它无法排除第三个变量z的影响。举个具体例子你想研究气温与降水量之间的关系但它们都受到季节的影响。双变量相干分析会显示出强烈的年周期相干这并不令人意外但你需要知道的是扣除季节因素之后气温与降水是否仍有额外的同步性。这就需要用偏小波相干Partial Wavelet CoherencePWC分析。偏小波相干的思路和非偏版本类似但需要在三个变量的互谱之间做偏相关运算。数学上它是在复数域上做类似“偏相关系数”的运算R²_xy|z |R_xy - R_xz × conj(R_yz)|² / ( (1 - |R_xz|²) × (1 - |R_yz|²) )其中 R_xy、R_xz、R_yz 分别是两两之间的复相关系数由平滑互谱计算得到conj表示共轭。这个公式和统计学里的偏相关系数公式长得一模一样只是全部换成了复数运算。如果你的zip包里不包含偏小波相干的代码自己在原基础上扩展也很快。关键点在于将之前算好的互小波谱输出复用而不是重新做变换。这样能省掉大量的重复计算。扩展时要注意平滑窗口参数保持一致否则偏相干的结果和双变量结果之间不具备可比性。5.2 时间-频率-空间多维度的网络分析思路现在的信号分析早就不局限于两个通道了。比如脑电有几十个通道你需要分析通道之间的同步网络看哪些脑区在特定频段耦合紧密。此时单靠小波相干矩阵会生成海量的图你需要把相干值压缩成网络指标节点强度、特征路径长度、聚类系数等。这些指标可以随时间和频段变化动态呈现大脑功能连接的演变。在实现上最直接的做法是分层处理先用小波相干计算每个通道对在目标频段上的平均相干值。设定一个阈值比如显著性水平的相干值只保留超过阈值的连接。将连接矩阵导入图分析工具MATLAB的Graph and Network Algorithms工具箱甚至可以直接用graph函数计算各种网络指标。这种做法的好处是结论的可视化非常直观——你不再需要一张张翻相干图而是看到网络拓扑结构随时频演变。缺点是要想结果稳定数据量要求很高通道数和时间点数都会影响网络指标的方差需要足够长的数据才能得到可靠结论。5.3 小波相位同步方法的替代小波相干算的是幅度归一化后的互谱关系它隐含了一种假设两个信号在局部的幅度变化模式一致。但有些场景下你更关心的是相位同步性——比如两个振荡器虽然有相同的振荡频率但彼此的相位始终保持着稳定的锁定关系此时无论幅度变化如何你都想捕捉这种耦合。这种情况下可以用相位锁定值Phase Locking ValuePLV。它直接由瞬时相位差计算PLV |mean( exp(1i × (φ_x(t) - φ_y(t))) )|在MATLAB中实现的时候你可以复用已经算好的小波变换结果提取特定频带的相位序列。平滑处理时PLV的计算窗口长度决定了最终的结果稳定性窗口太短PLV虚高因为点数少相位差的随机波动也能得到较大的统计值窗口太长时间分辨率损失。如果数据点很少建议用复杂的全局PLV而不是滑动窗口PLV。这里有一个反直觉的经验小波相干性高并不代表相位锁定强相位锁定强也不代表相干性一定高——因为相干性考虑了幅度的归一化但PLV不考虑幅度的时变特性。两种方法结合着看能更全面地描述两个非线性系统之间的耦合关系。最后补一句实操心得这份代码拿到手之后别急着拿自己最宝贵的数据直接跑先用造出来的简单信号验证一遍代码行为是否符合理论预期再分析真实数据。我自己测试代码至少花掉一整个晚上但这一步省下的排查时间远超过投入。遇到看不懂的中间变量就动手打印出来看尺寸、看数值范围比翻手册效率高得多。小波相干分析不是那种“输进去就能出结论”的工具它的每一步参数选择都要对数据特点有清醒的认识。祝你好运跑图顺利。本文还有配套的精品资源点击获取