小波阈值去噪MATLAB函数集：软硬阈值与5种改进方法封装

简介：面向信号处理与图像去噪方向的学习者与工程实践者，这份资源提供了一套可直接调用的小波阈值去噪 MATLAB 实现，覆盖软阈值、硬阈值以及五种改进阈值方法，适合已具备小波变换基础、希望快速验证不同阈值策略效果的中高级用户。作者在改造 wden、thselect、wthresh 三个函数的基础上进一步封装为 filterWaveletTh，延续专栏一贯风格，配合少量参数设置即可完成一次完整去噪流程。压缩包共 7 个文件，约 18KB，以 p 文件为主承载核心函数，另含 1 个 m 示例脚本与 1 份 docx 代码说明，便于对照理解调用方式与参数含义。已有 1324 人学习，说明该封装思路在同类需求中具有一定参考价值。读者可获得现成的阈值函数集合、可复用的封装接口以及效果评估与示例脚本，用于对比软硬阈值与改进方法的差异，并快速嵌入自己的滤波实验。

1. 小波阈值去噪：一份能直接调用的 MATLAB 函数集，把软硬阈值和 5 种改进方法都封装好了

做信号处理的人大概都遇到过这种场景：采集到的一段振动信号、心电信号或者语音波形，噪声混在里面，直接做频谱分析全是毛刺。这时候小波阈值去噪几乎是绕不开的一步——它不像傅里叶变换那样把噪声和信号摊在同一频域里，而是利用小波多分辨率特性，把噪声集中到高频细节系数上，再通过阈值把它们压掉。道理不复杂，但真到写代码的时候，软阈值、硬阈值、半软阈值、Garrote 阈值、改进的指数阈值……每种方法的阈值函数都不一样，每次换方法就得重写一遍循环，时间全耗在重复劳动上。

这份 MATLAB 代码解决的就是这个痛点：把软阈值、硬阈值以及 5 种改进阈值方法全部封装成独立函数，输入信号和参数就能直接出结果。适合正在做去噪实验、需要横向对比不同阈值方法效果的研究生和工程师，也适合刚接触小波分析、想先跑通再理解原理的新手。代码结构清晰，函数之间解耦，改一个阈值规则不影响其他部分，二次开发的门槛很低。

2. 项目结构与技术选型：为什么用 MATLAB 函数封装而不是脚本堆叠

2.1 整体目录与文件组织

这份代码没有复杂的工程结构，核心就是一组.m函数文件加一个演示脚本。常见的组织方式是这样：

wavelet_threshold/ ├── main_demo.m % 演示脚本，加载信号、调用去噪函数、画对比图 ├── wdenoise_custom.m % 主入口函数，负责小波分解、阈值处理、重构 ├── threshold_soft.m % 软阈值函数 ├── threshold_hard.m % 硬阈值函数 ├── threshold_garrote.m % Garrote 阈值函数 ├── threshold_semisoft.m % 半软阈值函数 ├── threshold_improved1.m % 改进方法一（指数型阈值） ├── threshold_improved2.m % 改进方法二（自适应阈值） ├── threshold_improved3.m % 改进方法三（分段阈值） ├── threshold_improved4.m % 改进方法四（对数型阈值） ├── threshold_improved5.m % 改进方法五（混合阈值） └── utils/ ├── calc_threshold.m % 通用阈值计算（通用阈值、SUREShrink、Minimax 等） └── plot_compare.m % 结果对比绘图辅助函数

这种「一个方法一个文件」的做法看起来有点笨，但实际用起来很顺手。因为每种阈值函数的数学表达式不同，参数个数也可能不一样，拆成独立文件后，函数签名可以各自定制，调用时通过函数句柄或者 switch-case 分发，互不干扰。我一般会把阈值函数统一成y = threshold_xxx(x, T, varargin)的形式，x是小波系数，T是阈值，后面的可变参数留给改进方法用，比如调节因子、分段点之类的。

提示：如果你拿到的代码里没有utils文件夹，阈值计算逻辑可能直接写在主函数里，不影响运行，但改起来会稍微麻烦一点。

2.2 技术选型：MATLAB 小波工具箱的依赖边界

这份代码依赖 MATLAB 的 Wavelet Toolbox，核心用到的是wavedec、waverec、detcoef、appcoef这几个函数。wavedec做多尺度分解，返回各层细节系数和最后一层近似系数；waverec负责重构。阈值处理只作用于细节系数，近似系数一般保留不动——因为噪声主要落在高频细节上，近似系数代表信号的低频轮廓，动了容易丢信息。

选 MATLAB 而不是 Python，主要是因为这个场景下 MATLAB 的小波工具箱接口更直接，wavedec一行就能拿到系数元胞数组，重构也只要一行。Python 的 PyWavelets 当然也能做，但如果你手头的实验数据、绘图流程都在 MATLAB 里，来回导数据反而增加出错概率。代码本身不依赖任何第三方工具箱之外的包，只要 MATLAB 版本带 Wavelet Toolbox 就能跑。

有一个细节值得注意：wavedec返回的系数结构是[cA_n, cD_n, cD_{n-1}, ..., cD_1]，顺序是从最粗尺度到最细尺度。很多新手第一次写阈值处理时容易搞反，把cD_1当成最粗的系数去处理，结果去噪效果完全不对。封装函数里一般会显式注释这个顺序，调用时按索引逐层处理即可。

2.3 主入口函数的调用约定

主函数wdenoise_custom的典型签名是这样的：

function [x_denoised, coeffs_out] = wdenoise_custom(x, wname, level, method, rule, varargin) % x 输入信号（列向量） % wname 小波基名称，如 'db4'、'sym8' % level 分解层数 % method 阈值方法：'soft'、'hard'、'garrote'、'semisoft'、'imp1'~'imp5' % rule 阈值规则：'sqtwolog'、'minimax'、'heursure'、'rigrsure' % varargin 改进方法的额外参数

调用示例：

load noisdopp; % MATLAB 自带示例信号 x = noisdopp; xd = wdenoise_custom(x, 'sym8', 5, 'soft', 'sqtwolog'); plot(x); hold on; plot(xd, 'r', 'LineWidth', 1.2); legend('原始信号', '去噪后');

这里level的选择有讲究。层数太少，噪声压不干净；层数太多，信号边缘容易被过度平滑。常见做法是level = fix(log2(length(x)))再减 1 到 2 层，或者直接用wmaxlev算最大允许层数再取个折中值。我一般会先跑 3 到 5 层看效果，再微调。

3. 核心模块拆解：从分解到重构的完整链路

3.1 小波分解与系数提取

整个去噪流程的第一步是把信号拆成多尺度系数。核心代码不长，但每一步都有对应的参数含义：

% 多尺度小波分解 [C, L] = wavedec(x, level, wname); % 提取第 1 层到第 level 层的细节系数 cD = cell(1, level); for k = 1:level cD{k} = detcoef(C, L, k); end % 提取最后一层近似系数 cA = appcoef(C, L, wname, level);

wavedec返回的C是所有系数拼接成的一个长向量，L是各段长度记录。detcoef(C, L, k)按层号取出细节系数，k=1对应最细尺度，k=level对应最粗尺度。这个层号顺序和很多人直觉相反，第一次用容易翻车——如果你把k=1当成最粗层去加大阈值，结果就是把信号的高频细节全削掉了，重构出来波形发钝。

appcoef取的是第level层的近似系数，也就是低频轮廓。这一步通常不做阈值处理，原样保留。有些改进方法会对近似系数也做轻微收缩，但那属于进阶操作，默认流程里不动它。

3.2 阈值计算：通用阈值、SUREShrink 与 Minimax

阈值T怎么定，直接决定去噪效果。代码里封装的calc_threshold一般支持四种经典规则：

function T = calc_threshold(detail_coef, rule, sigma) % detail_coef 某一层细节系数 % rule 阈值规则 % sigma 噪声标准差估计 n = length(detail_coef); switch lower(rule) case 'sqtwolog' % 通用阈值（Universal Threshold） T = sigma * sqrt(2 * log(n)); case 'minimax' % Minimax 阈值 T = sigma * (0.3936 + 0.1829 * log2(n)); case 'rigrsure' % SUREShrink 阈值，基于 Stein 无偏风险估计 T = sure_threshold(detail_coef, sigma); case 'heursure' % 启发式阈值，在 sqtwolog 和 rigrsure 之间选 T = heursure_threshold(detail_coef, sigma); otherwise error('未知阈值规则: %s', rule); end end

sigma的估计通常用median(abs(detail_coef)) / 0.6745，这是基于噪声系数服从高斯分布的中位数估计法，对异常值不敏感，比直接算标准差稳。通用阈值sqrt(2*log(n))偏保守，容易把信号也压掉一部分；SUREShrink 更自适应，但计算量稍大。实际对比时，我一般先用sqtwolog跑基线，再用rigrsure看有没有提升。

注意：sigma如果估错，后面所有阈值方法的效果都会偏。建议在每层细节系数上单独估sigma，而不是全局估一个，尤其是非平稳信号。

3.3 软阈值与硬阈值的实现差异

软阈值和硬阈值是最基础的两个函数，代码都很短，但行为差异很大：

function y = threshold_soft(x, T) % 软阈值：系数绝对值小于 T 的置零，大于 T 的向零收缩 y = sign(x) .* max(abs(x) - T, 0); end function y = threshold_hard(x, T) % 硬阈值：系数绝对值小于 T 的置零，大于 T 的原样保留 y = x .* (abs(x) >= T); end

软阈值的输出在|x| = T处连续，但导数不连续，重构信号更平滑，代价是幅度有偏——大系数也被统一减了T，信号边缘会变模糊。硬阈值保留了大系数的原始幅度，边缘更锐利，但在T附近不连续，容易产生伪吉布斯振荡。这两种方法的取舍是去噪实验里最常被讨论的点，没有绝对优劣，看你的信号更怕模糊还是更怕振铃。

代码里sign(x) .* max(abs(x) - T, 0)这个写法比用if判断向量化程度高，MATLAB 里对长信号处理速度快很多。硬阈值的x .* (abs(x) >= T)同理，逻辑索引直接乘，避免循环。

3.4 五种改进阈值函数的统一接口

改进方法的核心思路都是在软硬阈值之间找折中，或者让阈值函数随系数幅度自适应变化。代码里封装的 5 种方法，函数签名统一成y = threshold_impN(x, T, params)，方便主函数用switch分发：

function y = apply_threshold(x, T, method, varargin) % 统一阈值分发入口 switch lower(method) case 'soft' y = threshold_soft(x, T); case 'hard' y = threshold_hard(x, T); case 'garrote' y = threshold_garrote(x, T); case 'semisoft' y = threshold_semisoft(x, T, varargin{:}); case 'imp1' y = threshold_improved1(x, T, varargin{:}); case 'imp2' y = threshold_improved2(x, T, varargin{:}); case 'imp3' y = threshold_improved3(x, T, varargin{:}); case 'imp4' y = threshold_improved4(x, T, varargin{:}); case 'imp5' y = threshold_improved5(x, T, varargin{:}); otherwise error('未知阈值方法: %s', method); end end

以 Garrote 阈值和半软阈值举例，它们的表达式分别是：

function y = threshold_garrote(x, T) % Garrote 阈值：y = x - T^2/x，|x| > T；否则为 0 y = zeros(size(x)); idx = abs(x) > T; y(idx) = x(idx) - T^2 ./ x(idx); end function y = threshold_semisoft(x, T, alpha) % 半软阈值：在 [T, alpha*T] 区间内做线性过渡 if nargin < 3, alpha = 2; end y = zeros(size(x)); absx = abs(x); idx1 = absx < T; idx2 = absx >= T & absx < alpha * T; idx3 = absx >= alpha * T; y(idx2) = sign(x(idx2)) .* (absx(idx2) - T) * alpha / (alpha - 1) / alpha; y(idx3) = x(idx3); end

Garrote 阈值在|x|较大时接近硬阈值，在|x|接近T时平滑过渡到零，兼顾了边缘保留和连续性。半软阈值的alpha参数控制过渡区间宽度，alpha越大过渡越缓，去噪后信号越平滑但细节损失也越多。这 5 种改进方法的具体公式各有不同，但接口一致，换方法只需要改method字符串，不用动主流程。

3.5 系数回写与信号重构

所有细节系数处理完之后，要把它们按原顺序拼回C向量，再调waverec重构：

% 将处理后的细节系数写回 C C_new = C; offset = L(1); % 跳过近似系数段 for k = level:-1:1 len = L(level - k + 2); C_new(offset + 1 : offset + len) = cD_processed{k}; offset = offset + len; end % 重构信号 x_denoised = waverec(C_new, L, wname);

这里最容易出错的是索引偏移。L的结构是[length(cA_level), length(cD_level), ..., length(cD_1)]，所以细节系数在C里的排列顺序是从粗到细，而cD{k}里k=1是最细层。回写时循环要从level递减到 1，偏移量按L对应段累加。写错顺序的话，重构出来的信号长度对但内容全乱，波形看着像噪声——这个坑我踩过不止一次，后来每次都在回写后加一句长度断言assert(length(C_new) == length(C))才放心。

4. 跑起来：从加载信号到出对比图的完整步骤

4.1 环境准备与依赖检查

先把 MATLAB 打开，确认 Wavelet Toolbox 可用：

% 检查小波工具箱是否安装 ver('wavelet')

如果输出里有 Wavelet Toolbox 的版本信息，说明环境没问题。没有的话需要在附加功能管理器里装一下。代码本身不依赖其他工具箱，plot、median、log2这些都是基础函数。

把代码文件夹加到路径里：

addpath(genpath('wavelet_threshold'));

genpath会把子文件夹一起加进来，utils里的辅助函数也能直接调用。

4.2 用自带信号跑通第一个例子

main_demo.m里一般会加载 MATLAB 自带的noisdopp或者leleccum信号做演示。直接运行：

main_demo

脚本内部大致做这几件事：加载信号、加噪（或者用自带含噪信号）、调用wdenoise_custom分别用软阈值和硬阈值去噪、画三张子图对比。第一次跑通之后，你会看到原始信号、软阈值结果、硬阈值结果并排显示，直观感受两种方法的差异。

如果想用自己的数据，把加载那行换掉就行：

x = load('your_signal.txt'); % 假设是列向量文本 x = x(:); % 强制列向量 xd = wdenoise_custom(x, 'db4', 4, 'garrote', 'rigrsure');

x(:)这个操作别省，行向量传进去某些函数会按行处理，结果维度对不上，报错信息还不一定明显。

4.3 批量对比 7 种阈值方法的效果

单次调用只能看一种方法，做实验时通常要横向对比。写个循环批量跑：

methods = {'soft','hard','garrote','semisoft','imp1','imp2','imp3'}; results = cell(1, length(methods)); for i = 1:length(methods) results{i} = wdenoise_custom(x, 'sym8', 5, methods{i}, 'sqtwolog'); end % 计算信噪比和均方根误差 snr_vals = zeros(1, length(methods)); rmse_vals = zeros(1, length(methods)); for i = 1:length(methods) snr_vals(i) = 20 * log10(norm(x) / norm(x - results{i})); rmse_vals(i) = sqrt(mean((x - results{i}).^2)); end

这里x是干净信号，results{i}是去噪结果。SNR 越高、RMSE 越低说明去噪效果越好。实际含噪场景下没有干净信号做参考，可以用相邻层系数的一致性或者视觉判断来辅助。批量跑完把snr_vals和rmse_vals打成表格，哪种方法在你的数据上表现好一目了然。

提示：semisoft和部分改进方法需要额外参数，批量调用时如果没传，函数内部要有默认值兜底，否则会报「输入参数不足」。检查一下各函数的nargin处理。

4.4 参数调优：小波基、层数与阈值规则的组合

三个关键参数需要调：小波基wname、分解层数level、阈值规则rule。常见组合经验是：

调参时不要一次全改，固定两个调一个。我一般先定wname和level，用sqtwolog跑基线，再换rigrsure和heursure看阈值规则的影响，最后在表现最好的规则下试改进阈值方法。这样每一步的变化都能归因，不会出现「改了一堆参数不知道哪个起了作用」的情况。

5. 避坑与排查：跑这份代码时最容易翻车的几个地方

5.1 重构信号长度和原始信号对不上

现象：waverec出来的x_denoised比输入x长一两个点，画图时尾部多出一截。

原因：小波分解对信号长度有要求，wavedec内部会做边界延拓，重构时如果C_new的长度和L记录不一致，waverec会按延拓后的长度输出。

解决：重构后截断到原始长度x_denoised = x_denoised(1:length(x));，或者在分解前把信号补零到2^level的整数倍。前者简单，后者更规范。我一般用截断，因为补零会引入边界效应。

5.2 阈值处理后系数全为零

现象：去噪结果是一条直线，所有细节系数都被压没了。

原因：阈值T算得过大，或者sigma估计偏大。常见于sqtwolog规则下信号长度很短的情况，sqrt(2*log(n))虽然随n增长慢，但sigma如果用了全局标准差而不是中位数估计，会被大幅值系数拉高。

解决：检查sigma是否用median(abs(cD))/0.6745估计；换rigrsure规则试试；或者对阈值乘一个小于 1 的缩放因子（常见 0.5～0.8）做保守处理。

5.3 改进阈值函数报「未定义参数」

现象：调用imp1到imp5时报错，提示某个变量未定义。

原因：改进方法通常有额外参数（调节因子、分段系数），主函数分发时varargin没传进去，或者函数内部默认值没写全。

解决：检查apply_threshold里varargin{:}是否透传；在各threshold_improvedN.m开头加if nargin < 3, param = default_value; end兜底。改完用nargin打印一下实际传入参数个数确认。

5.4 不同 MATLAB 版本函数行为差异

现象：同一份代码在旧版本能跑，新版本报错或者结果不同。

原因：wavedec的边界处理模式在不同版本间有过调整，detcoef的返回维度也可能从行向量变成列向量。

解决：在分解时显式指定边界模式wavedec(x, level, wname, 'mode', 'symmetric')，避免默认值变化；取完系数后统一cD{k} = cD{k}(:)强制列向量。这样跨版本行为一致。

5.5 绘图对比时信号重叠看不清

现象：原始信号和去噪信号画在同一张图上，颜色接近，分不清哪条是哪条。

原因：plot默认颜色循环，多条线容易撞色。

解决：显式指定颜色和线宽，plot(x, 'b', 'LineWidth', 1); hold on; plot(xd, 'r--', 'LineWidth', 1.2);，图例用legend('原始','去噪')。对比多种方法时用subplot分行画，比叠在一起清楚得多。

6. 二次开发：换阈值规则、加新方法、做批量实验的改造思路

拿到这份代码之后，最常见的改造需求是加一种自己的阈值函数，或者把阈值规则换成更适合特定信号的版本。加新方法的流程很固定：新建threshold_my.m，签名对齐y = threshold_my(x, T, varargin)，然后在apply_threshold的switch里加一个case 'my'，最后在wdenoise_custom的参数校验里把'my'加进合法方法列表。三步做完就能调用，不用动分解和重构的代码。

换阈值规则也是类似的路子。calc_threshold里加一个case，实现你的阈值公式，返回标量T。如果新规则需要额外输入（比如先验噪声水平），通过varargin传进来，别改函数签名，否则所有调用点都得跟着改。我一般会在函数开头写一段注释，标明这个规则的参考文献和适用条件，过几个月回来看还能想起当时为什么这么写。

做批量实验时，建议把参数组合和结果指标抽成一个独立脚本，用表格驱动：

wnames = {'db4','sym8','coif3'}; levels = [3 4 5]; rules = {'sqtwolog','rigrsure','heursure'}; methods = {'soft','hard','garrote','semisoft','imp1','imp2','imp3','imp4','imp5'}; results = []; for w = 1:length(wnames) for l = 1:length(levels) for r = 1:length(rules) for m = 1:length(methods) xd = wdenoise_custom(x, wnames{w}, levels(l), methods{m}, rules{r}); snr = 20*log10(norm(x)/norm(x-xd)); results = [results; {wnames{w}, levels(l), rules{r}, methods{m}, snr}]; end end end end

跑完按 SNR 降序排一下，最优组合直接出来。这种全组合跑起来耗时不小，但一次跑完比手动试几十次靠谱。注意results用元胞数组存混合类型，后面取数值列时记得cell2mat转换。

还有一个容易被忽略的点：改进阈值方法里的调节参数，默认值往往是对某类信号调好的，换到你的数据上不一定最优。二次开发时可以把这些参数也纳入网格搜索，但组合数会爆炸，建议先用默认值筛出前几名方法，再对这几个方法单独调参。从那以后我每次加新阈值函数，都强制先在一个已知含噪的合成信号上验证——干净信号加高斯白噪声，看 SNR 提升是否符合预期，确认函数本身没写错，再上真实数据。这个习惯帮我省了不少排查时间。希望帮到你。

本文还有配套的精品资源，点击获取