小波阈值滤波算法（软阈值、硬阈值、改进的阈值）

最后更新于：2024-04-07 21:29:30

一、代码运行环境

MATLAB2018b及更新版本。

二、程序介绍

注：图标代表该m文件为脚本文件，可以直接运行；图标代表函数文件，在没有输入变量的情况下无法直接运行。更详细的解释可以看这里。

1.demoWaletThFilter.m文件

小波阈值滤波函数filterWaveletTh.m的测试脚本文件可以直接运行。程序运行完成后，在MATLAB的工作区，双击fil变量，可以查看滤波得到的具体数值。

同时程序将会画出滤波前后对比图：

同时在命令行窗口将会打印出滤波结果评价指标：

2.filterWaveletTh.m文件

使用小波阈值法实现滤波的主函数。

阈值滤波的不同参数设置全部封装在此函数，只需要输入待滤波数据、小波名称、阈值类型、小波分解水平和阈值选择规则，实现一行代码完成小波阈值滤波。

function s = filterWaveletTh(data,wname,SORH,lev,tptr)

% 使用小波阈值法实现滤波，本文件是为了形式统一方便调用的一层封装，具体实现在kwden.m、kwdencmp.m、kwthresh.m、kthselect.m文件中
% 这4个文件均为经本人调整过后的函数文件，原函数分别为wden,wdencmp,wthresh,thselect，所以关于这几个函数更多介绍可以查看对应的官方帮助文档
% 输入：
% data：待滤波数据
% wname：小波名称，可选范围参考这里：https://ww2.mathworks.cn/help/wavelet/ref/wfilters.html?searchHighlight=wname&s_tid=srchtitle_wname_2#d123e130597
% SORH：  阈值类型 可以选择：
%       -'s' 软阈值
%       -'h' 硬阈值
% 
%       SORH ='a1'时，采用改进方法1，改进后的软阈值：
%                     参考论文：孙万麟, 王超. 基于改进的软阈值小波包网络的电力信号消噪[J]. 海军工程大学学报, 2019(4).
%       SORH = 'a2'时，采用改进方法2：
%                     参考论文：刘冲, 马立修, 潘金凤,等. 联合VMD与改进小波阈值的局放信号去噪[J]. 现代电子技术, 2021, 44(21):6.
%       SORH = 'a3'时，采用改进方法4：
%                     参考论文：基于改进小波阈值-CEEMDAN算法的ECG信号去噪研究
%       SORH = 'a4'时，采用改进方法3：
%                     参考论文：基于改进小波阈值函数的语音增强算法研究
%                     采用该方法需要输入两个调节因子，分别是a4_alpha和a4_gamma，取a4_alpha>0，0<a4_gamma<1
%       SORH = 'a5'时，采用改进方法3：
%                     参考论文：《基于VMD与改进小波阈值的地震信号去噪方法研究》
%                     采用该方法需要输入两个调节因子，分别是alpha和beta
% lev：  小波分解水平
% tptr： 阈值选择规则，有可选类型：
%       '' — Adaptive threshold selection using the principle of Stein's Unbiased Risk Estimate (SURE).
%       'sqtwolog' — Fixed-form threshold is sqrt(2*log(length(X))).
%       'heursure' — Heuristic variant of 'rigrsure' and 'sqtwolog'.
%       'minimaxi' — Minimax thresholding.
%       'visushrink' - 通用阈值去噪方法
% options：部分阈值函数需要补充的参数设置，即某些调节因子，需要用结构体方式调用幅值，如options.a4_alpha=2。如果不需要设置，可以赋值为options=[]。具体如下：
%          -a4_alpha：改进方法4的alpha值设置 a4_alpha>0
%          -a4_gamma：改进方法4的gamma值设置 0<a4_gamma<1
%          -a5_alpha：改进方法5的alpha值设置
%          -a5_beta：改进方法5的beta值设置
% 输出：
% s：经过滤波后得到的数据

3.FilterEffectEvaluation.m文件

对滤波效果进行评估。目前包含的评估指标有：信噪比SNR，均方差MSE，波形相似系数NCC。

function [SNR,MSE,NCC] = FilterEffectEvaluation(ori,fil)

% 对滤波效果进行评估
% 目前包含的评估指标有：信噪比SNR，均方差MSE，波形相似系数NCC
% 输入：
% ori：无噪声的原始数据，一维序列
% fil：滤波后的数据，一维序列
% 输出：
% SNR：信噪比
% MSE：均方误差
% NCC：波形相似系数

4.kthselect.m文件

经店长修改后的小波阈值滤波阈值选择函数（原函数为thselect），加入了 VisuShrink 方法，即统一阈值去噪方法。

function thr = kthselect(x,tptr,lev,wname)


% 经 Mr.看海 修改后的小波阈值滤波阈值选择函数（原函数为thselect）
% 加入了 VisuShrink 方法，即通用阈值去噪方法
% 该方法的实现主要参考：https://core.ac.uk/download/pdf/217399263.pdf （Denoising techniques - a comparison）
% 如果要使用VisuShrink方法，需要令TPTR =
% 'visushrink'，注意大小写，并指定lev和wname值。此时输出thr为lev*1维度的数组
% 其中lev：  小波分解水平
%     wname：小波名称

5.kwden.m文件

经店长修改后的小波阈值滤波函数（原函数为kwden），加入了 visushrink 方法，即统一阈值去噪方法。

function [xd,cxd,lxd,thrs] = kwden(in1,in2,in3,in4,in5,in6,in7)

% 经 Mr.看海 修改后的小波阈值滤波函数（原函数为wden）
% 加入了 visushrink 方法，即统一阈值去噪方法
% 加入了改进阈值函数，具体见kwthresh.m文件
% 加入了为满足改进阈值函数额外的调节因子的入口参数

6.kwthresh.m文件

经店长修改过后的阈值函数（原函数为wthresh）。除软阈值、硬阈值外，添加了改进阈值函数类型。

function y = kwthresh(x,sorh,t)

% 经 Mr.看海 修改过后的阈值函数（原函数为wthresh）
% 添加了改进阈值函数类型
% 输入：
% x：阈值函数的输入值，不需修改
% sorh：阈值函数类型
%       sorh ='a1'时，采用改进方法1，改进后的软阈值：
%                     参考论文：孙万麟, 王超. 基于改进的软阈值小波包网络的电力信号消噪[J]. 海军工程大学学报, 2019(4).
%       sorh = 'a2'时，采用改进方法2：
%                     参考论文：刘冲, 马立修, 潘金凤,等. 联合VMD与改进小波阈值的局放信号去噪[J]. 现代电子技术, 2021, 44(21):6.
%       sorh = 'a3'时，采用改进方法4：
%                     参考论文：基于改进小波阈值-CEEMDAN算法的ECG信号去噪研究
%       sorh = 'a4'时，采用改进方法3：
%                     参考论文：基于改进小波阈值函数的语音增强算法研究
%                     采用该方法需要输入两个调节因子，分别是alpha和gamma，取alpha>0，0<gamma<1
%       sorh = 'a5'时，采用改进方法3：
%                     参考论文：《基于VMD与改进小波阈值的地震信号去噪方法研究》
%                     采用该方法需要输入两个调节因子，分别是alpha和beta
% t：阈值
% 输出：
% y：阈值函数输出值
% 测试阈值函数图线可以用下述代码：
% figure('color','w');plot(-5:0.1:5,kwthresh(-5:0.1:5,'a1',1));
% 测试阈值函数图线与软硬阈值图线对比图可用下述代码：
% figure('color','w');plot(-5:0.1:5,kwthresh(-5:0.1:5,'a1',1));hold on
% plot(-5:0.1:5,kwthresh(-5:0.1:5,'s',1))
% plot(-5:0.1:5,kwthresh(-5:0.1:5,'h',1))
% legend('改进阈值方法','软阈值','硬阈值')

硬阈值、软阈值、改进阈值方法1、改进阈值方法2的阈值函数分别如下图：

软阈值

硬阈值

改进阈值方法1

改进阈值方法2

改进阈值方法3

改进阈值方法4

改进阈值方法5

三、快速开始

1.运行测试脚本

先在MATLAB里打开下载好的文件夹，然后运行demoWaletThFilter.m程序，程序运行完毕后如果没报错，且正常画出滤波对比图像，则说明运行环境正常，程序正确。

2.修改仿真数据/导入数据

复制一个demoWaletThFilter.m的文件副本，在这个副本里做如下修改：

（1）第一种情况，你可能想要对你自己要研究的仿真数据进行滤波测试，此时你需要对 demoWaletThFilter.m 脚本文件中的第1小节内容进行修改替换即可。需要注意的是，请最好保持变量名的一致，即将无噪声信号命名为x，添加噪声后的信号命名为sig，x和sig都应该是一维数据。此时数据替换完成。

（2）第二种情况，你可能是想对一段真实采集的数据进行滤波，此时需要根据你的文件类型的不同（excel,txt,csv等），将数据导入MATLAB的方法有所不同。同学们可以看博主针对常用文件的导入方法的这个教程，教程上未包含的数据类型，大家可以再参考这个文档。导入完成后请讲这个待滤波信号命名为sig。此时你可能没有纯净信号（在大多数实际应用中是没有纯净信号的），这不影响你是用小波阈值滤波进行信号去噪，但是你将无法计算滤波评价指标。

3.实现小波阈值滤波

参照 demoWaletThFilter .m文件中第二步，根据需要调整小波名称、小波分解水平、阈值类型、阈值选择规则。

（需要注意，对于使用真实采集的数据进行滤波，且无纯净信号的情况下，需要删除代码中的第3节滤波评价指标计算，并将第4节滤波对比图像按照注释进行修改）

运行程序即可。

四、关于完整版与公开版代码

五、获取公开版程序（需使用电脑浏览器打开）

注：公开版代码需使用MATLAB2022a及以上版本

六、获取完整版程序（使用电脑浏览器或者手机浏览器打开）

获取通道一（淘宝）：点击此处获取完整版程序

获取通道二（本页面）：点击下面“立即支付”按钮，付款后获取完整版代码下载链接和售后联系方式~本通道处于测试阶段，使用该通道可以额外优惠（仅需56元）。付款完成后刷新一下本页面即可看到下载链接。

（注意支付跳转失败的话，请使用浏览器打开本页面）

七、完整版代码重要更新

20230524 修复已知bug。

20221207 补充了改进方法3~改进方法5，增加了如何修改添加新的阈值函数代码的视频教程（完整版独有）。

20220818 完成初版代码

八、常见问题

1.使用改进的阈值方法的滤波效果一定比软、硬阈值要好吗？

不是的。对于同一段信号，软、硬阈值的方法算出来信噪比等评估参数效果比改进的阈值方法要好也是有可能的。具体的阈值函数的选择很多时候是要和研究对象相结合的。

2.为什么在使用改进方法时，改变调节因子，滤波的结果不发生变化？

这可能是阈值过大或者过小（通常是过大）导致的，如果阈值过大，那么大部分甚至所有的小波系数都可能小于这个阈值，此时调整调节因子滤波结果也不会发生变化。根据经验，将tptr设置为’visushrink’可以较少发生这种情况，如果tptr已经设置为’visushrink’，则建议使用别的阈值选择规则试一试哦！