ARTICLE DETAIL

资讯详情

深耕网站建设、视觉设计与SEO优化的一线实战洞察。

MATLAB小波分析:从尺度函数到信号去噪的工程实践

MATLAB小波分析:从尺度函数到信号去噪的工程实践 简介面向MATLAB使用者的小波分析学习资源适合信号处理、图像分析初学者及需要快速上手小波工具箱的工程人员。内容围绕小波函数和尺度函数的核心性质展开重点演示如何在MATLAB中调用wavemngr、wavfun等函数生成小波并绘制对应曲线帮助读者建立从构造、计算到可视化的整体认知。压缩包共2个文件包含1个.m脚本和1个.fig图形文件整体大小仅54KB脚本用于套用常见小波基完成函数取值计算fig文件则可直接查看小波函数与尺度函数波形便于对照运行结果深化理解。已有1199人学习过该资源。借助包内案例可快速认识Haar、Daubechies、Morlet等常用小波基的差异与特点通过逐行运行脚本也能直观理解尺度函数的构造过程为后续在信号去噪、图像压缩、特征提取及多尺度分析等场景中的应用提供清晰的实践入口。1. 小波函数与尺度函数先看这对搭档再打开MATLAB小波变换和傅里叶变换最大的区别在于它同时保留时域和频域信息。要做到这一点靠的并不是一个函数而是一对搭档。尺度函数负责抓信号里缓慢变化的趋势小波函数负责抓突变和细节两者通过多分辨率分析互相补充共同构成小波变换的骨架。很多人能熟练使用dwt和wavedec却不太清楚wfilters返回的四组滤波器到底是干什么的也不知道wavefun画出来的波形代表什么分辨率。这部分内容会把“小波函数、尺度函数”这两个概念落到MATLAB的具体函数和可执行代码上。你不必先把多分辨率分析的推导全部啃完按顺序跑一遍wfilters、wavefun、dwt再回头看二尺度方程那些公式符号会发现就是手边的Lo_D和Hi_D。适合正在做信号去噪、故障诊断或时间序列预测的工程师也适合准备把小波系数作为特征输入给神经网络的读者。2. 多分辨率分析尺度函数定骨架小波函数补细节2.1 尺度空间与小波空间两个子空间如何拼出完整信号多分辨率分析的出发点是构造一列嵌套子空间记作V_j。V_j由尺度函数的整数平移族φ_{j,k}(t)2^{j/2}φ(2^j t-k)张成V_j随着j增大而扩大且满足嵌套关系V_j⊂V_{j1}。单纯用尺度函数逼近信号能拿到越来越精细的近似但两个相邻尺度之间还有“差”的部分没有表达这部分就是由小波函数生成的子空间W_j。形式化地看V_{j1}V_j⊕W_j也就是下一层尺度空间被拆成上一层尺度空间与细节空间的直和。直和的含义是信号中的信息在每一层被划分成两个正交的分量属于V_j的低频趋势以及属于W_j的高频细节。把这种拆分由深到浅逐层做下去就得到常见的小波分解树。工程上关心的是上述抽象的“子空间逼近”最终都会被转化为两组离散滤波器系数一组给尺度函数一组给小波函数这决定了你在MATLAB里操作的对象其实是滤波器向量而非连续函数。所以理解尺度函数和小波函数重点并不是无穷迭代的数学定义而是这三件事两个函数对应两个子空间、两个子空间互不重叠、两层之间存在二尺度递推关系。下一节从这个递推关系直接跳进MATLAB的系数世界。2.2 二尺度方程与滤波器系数先用wfilters抓出h和g两个子空间之间的递推关系由二尺度方程描述φ(t)√2 Σ_n h(n) φ(2t-n) ψ(t)√2 Σ_n g(n) φ(2t-n)其中h是尺度函数对应的低通系数g是小波函数对应的高通系数。MATLAB的wfilters返回的正是这两种系数以及它们的重构版本。下面用db2小波做一次完整读取% 读取 db2 小波的四组滤波器系数 [Lo_D, Hi_D, Lo_R, Hi_R] wfilters(db2); disp(分解低通 Lo_D:); disp(Lo_D); disp(分解高通 Hi_D:); disp(Hi_D);Lo_D和Hi_D分别对应二尺度方程里的h和g用于分解信号Lo_R和Hi_R用于重构是把信号从系数域拼回时间域时使用的对偶滤波器。Lo表示低通D是decompositionR是reconstruction命名规则直接告诉你在哪一个环节使用Lo和Hi不能对调D和R也不能混用。正交小波中Lo_R恰好是Lo_D的时间反转Hi_D也可以由Lo_R推导出来因此只需要存一半系数。公式里的√2来自归一化约定MATLAB的wfilters返回系数已经吸收了这部分因子不需要自己再乘√2。如果从文献里手写系数要先把约定核对清楚再交给dwt这是新手很容易踩的坑。二尺度方程里的系数与MATLAB变量一一对应如下方程符号含义MATLAB变量用在哪个环节h(n)尺度函数的低通系数Lo_D分解出近似系数cAg(n)小波函数的高通系数Hi_D分解出细节系数cDh(n)重构端低通系数Lo_R由cA重建低频分量g(n)重构端高通系数Hi_R由cD重建高频分量2.3 用wavefun算出来看一眼两个函数长什么样% 查看 db2 小波的尺度函数和小波函数时域波形 [phi, psi, xval] wavefun(db2, 10); figure; subplot(2,1,1); plot(xval, phi, LineWidth, 1.5); % 尺度函数波形 title(db2 尺度函数); grid on; subplot(2,1,2); plot(xval, psi, LineWidth, 1.5); % 小波函数波形 title(db2 小波函数); grid on;wavefun的第二个参数是迭代次数。每次迭代相当于把内部生成函数的采样点加密一倍数字越大输出波形越接近真正的连续函数。工程画图用8~10次迭代即可论文插图想要平滑曲线可以取12但计算量也随之增大。输出xval是等间距采样横轴phi和psi的纵轴幅度没有统一标准haar这类特殊小波才是规整的±1db族大部分值落在-2到2之间直接用plot画不需要手工缩放。到这里可以直观看到尺度函数是带支撑区的低通形状小波函数是带振荡的高通形状两者在时域上互补。这个认识会直接影响第4章里小波基的选择细节信号明显时优先关心小波函数的形状和消失矩。3. MATLAB的小波函数与尺度函数入口wfilters和wavefun的参数细节3.1 wfilters返回的四组系数分解与重构别混用wfilters除了接受db2这样的名称外也接受sym4、coif3、bior3.5等。下面的代码把sym4四组系数画成离散序列能更清楚看到四个滤波器的长度和支持范围差异wname sym4; [Lo_D, Hi_D, Lo_R, Hi_R] wfilters(wname); figure; subplot(2,2,1); stem(Lo_D); title(Lo\_D); subplot(2,2,2); stem(Hi_D); title(Hi\_D); subplot(2,2,3); stem(Lo_R); title(Lo\_R); subplot(2,2,4); stem(Hi_R); title(Hi\_R);对正交小波Lo_R是Lo_D的时间反转Hi_D与Lo_R之间满足调制关系所以四组系数并不独立。双正交小波则不同四组滤波器都要单独保存这也是wfilters在bior族上返回更特殊结构的原因。如果写错小波名先执行waveinfo(db)、waveinfo(bior)查看支持列表再核对括号里的族名和序号。3.2 wavefun的调用方式正交小波和双正交小波返回值不同正交小波调用wavefun时用三个输出[phi, psi, xval]。双正交小波在分解端和重构端各有一套尺度函数和小波函数因此要接收五个输出% 双正交小波返回两套尺度函数 [phi1, psi1, phi2, psi2, xval] wavefun(bior3.5, 8); figure; plot(xval, phi1, b, xval, phi2, r, LineWidth, 1.2); legend(分解端 \phi, 重构端 \phi); title(bior3.5 的两套尺度函数);phi1和psi1对应分解端的尺度函数与小波函数phi2和psi2对应重构端。因为双正交小波不再要求两个子空间严格正交而是要求分解和重构这对对偶系统满足完全重构条件所以会出现两套函数。图像处理里bior族更常见正是因为它有线性相位能避免重构时边缘相位畸变。3.3 小波族选型速查db、sym、bior各管一摊小波族wname示例正交性线性相位典型场景haarhaar正交有入门演示、快速验证daubechiesdb4正交无通用信号去噪、故障诊断symletssym4正交近似线性需要减少相位失真时coifletscoif3正交无平衡支撑长度与光滑性biorthogonalbior3.5双正交有图像处理、信号重构选型顺序我一般这样走先看是否需要线性相位需要就去bior族不需要再比较消失矩和支撑长度db4到db8是常用的折中区间最后用同样的阈值去噪流程跑一遍比较输出信噪比而不是只看曲线形状。MATLAB里通过waveinfo可以查看每个族消失矩、支撑长度的官方说明比在网上找二手对比表更可靠。4. 小波分解重构与去噪从dwt到waverec的参数设定4.1 dwt/idwt一层分解先看系数怎么拆、怎么拼用一个叠加了50 Hz正弦、300 Hz正弦和随机噪声的信号做单层分解fs 1000; t (0:999) / fs; x sin(2*pi*50*t) 0.2*sin(2*pi*300*t) 0.05*randn(1, 1000); % 单层小波分解cA 是近似cD 是细节 [cA, cD] dwt(x, db4); % 直接重构验证信息是否完整 xr idwt(cA, cD, db4); fprintf(最大重构误差: %.3e\n, max(abs(x - xr)));dwt的输出cA和cD长度约为输入信号的一半精确长度由信号长度、滤波器长度和边界延拓模式共同决定不要假设一定是floor(N/2)。idwt能把系数还原成和x等长的序列。如果重构误差明显偏离机器精度优先检查边界模式和信号长度而不是怀疑小波函数选错了。4.2 wavedec/waverec多层分解层数和每层频带怎么对应N 4; [C, L] wavedec(x, N, db4); % 取出第 4 层近似和第三层细节的时间序列 A4 wrcoef(a, C, L, db4, 4); D3 wrcoef(d, C, L, db4, 3);wavedec把所有层的系数拼在一维数组C里L记录每段分界点通常不需要手工解析。wrcoef按类型和层号直接把系数恢复成与x等长的时间序列A4是第4层近似D3是第3层细节。N的上限由信号长度决定理论上不能超过floor(log2(length(x)))工程上取到目标频段所在层即可。以fs1000为例各层分量对应的频带近似如下分量频带fs1000 HzD1250 ~ 500 HzD2125 ~ 250 HzD362.5 ~ 125 HzD431.25 ~ 62.5 HzA40 ~ 31.25 Hz这是理想滤波器组下的划分真实小波滤波器存在过渡带跨层会有轻微重叠。机械振动诊断里常用这个表反推知道故障特征频率落在哪一段就取对应层数分解再做包络谱。4.3 小波阈值去噪thselect和wthresh怎么配合thr thselect(C, sqtwolog); % 选择阈值 Cd wthresh(C, s, thr); % 软阈值收缩 xd waverec(Cd, L, db4); % 重构去噪thselect的四种阈值规则适用于不同噪声强度最小化风险规则适合弱噪声固定阈值规则最常见但噪声强时容易过度平滑。sorh参数决定收缩方式s软阈值会把所有系数向零压缩h硬阈值保留超过阈值的原值硬阈值重构的信号更锐利但会在断点处留下毛刺。工程上可以先跑一组不同规则的对比用去噪后信号的均方根误差和峰度两个指标一起判断。规则含义rigrsureStein无偏风险估计适合弱噪声heursure启发式混合规则兼顾强、弱噪声sqtwolog固定阈值噪声较强时常用minimaxi极小极大准则保留更多细节小波基和阈值的搭配没有固定答案。同一个信号db4加sqtwolog可能比sym8加rigrsure更合适。建议把阈值规则、小波族、分解层数三个变量做成循环在验证集上选最优组合而不是依赖经验拍板。5. 进阶自定义滤波器组、边界模式与Elman网络特征提取5.1 自定义滤波器组检查正交条件后再交给dwtMATLAB的小波函数库不是封闭的dwt允许直接传入滤波器系数。系数必须满足正交条件最常用的一组验证是直流增益sum(h)/sqrt(2)接近1且偶数位自相关接近0。下面用db2的系数做示例实际使用时替换成你自己设计的hh [0.482962913144534, 0.836516303737808, ... 0.224143868042013, -0.129409522551260]; Lo_D h; Hi_D fliplr(h) .* (-1).^(0:length(h)-1); Lo_R fliplr(Lo_D); Hi_R fliplr(Hi_D); fprintf(直流增益 %.6f\n, sum(h) / sqrt(2)); [cA, cD] dwt(x, Lo_D, Hi_D); xr idwt(cA, cD, Lo_R, Hi_R);高通由低通翻转调制得到这是正交小波的通用构造方式。若直流增益明显偏离1重构误差会立刻放大。先跑重构误差误差接近机器精度再把这些系数用于后续特征提取。5.2 边界模式用dwtmode和最大重构误差验证边界延拓会直接影响重构误差。dwtmode(per)做周期延拓对本来就近似周期的信号效果好dwtmode(sym)做对称延拓对一般实测信号更常用。对比代码dwtmode(per); [ca_p, cd_p] dwt(x, db4); xr_p idwt(ca_p, cd_p, db4); dwtmode(sym); [ca_s, cd_s] dwt(x, db4); xr_s idwt(ca_s, cd_s, db4); fprintf(per 误差: %.3e\n, max(abs(x - xr_p))); fprintf(sym 误差: %.3e\n, max(abs(x - xr_s)));注意dwtmode修改的是全局设置调试完要恢复默认模式否则后续脚本都可能被影响。5.3 把分解系数做成Elman网络特征四个验证点小波分解配合Elman网络做时间序列预测是一类常见组合特征通常取各层细节系数的能量统计量[C, L] wavedec(x, 4, db4); feat zeros(1, 5); for k 1:4 d wrcoef(d, C, L, db4, k); feat(k) sum(d.^2) / length(d); end feat(5) sum(wrcoef(a, C, L, db4, 4).^2);接模型之前做四个验证一是重构误差基准用rng(0)、randn(1,512)固定输入把最大重构误差压到1e-12以下二是特征归一化到同一量级三是按时间段切分训练集和测试集避免随机抽样破坏时间相关性四是对比做与不做阈值去噪两种特征的效果小波去噪带来的提升有时在特征层面已经体现。前两步过关后再讨论网络结构才有意义。本文还有配套的精品资源点击获取
返回列表