主要内容

小波去噪与非参数函数估计

小波工具箱™提供了许多函数,用于在噪声中估计未知函数(信号或图像)。您可以使用这些函数去噪信号,并作为非参数函数估计的方法。

最通用的一维模型是

(n) =f(n) + σe(n)

哪里n=0,1,2,…N-1。这个e(n)高斯随机变量分布为N(0, 1)。σ的方差e(n)是σ吗2.

在实践中,(n)通常是一个离散时间信号,其等时间步长被加性噪声破坏,您正试图恢复该信号。

更一般地说,您可以查看(n)作为N维随机向量

( f ( 0 ) + σ e ( 0 ) f ( 1 ) + σ e ( 1 ) f ( 2 ) + σ e ( 2 ) f ( N 1 ) + σ e ( N 1 ) ) = ( f ( 0 ) f ( 1 ) f ( 2 ) f ( N 1 ) ) + ( σ e ( 0 ) σ e ( 1 ) σ e ( 2 ) σ e ( N 1 ) )

在一般情况下,去噪和回归之间的关系是清楚的。

你可以用n × m的随机矩阵替换n × 1的随机向量来获得被加性噪声损坏的图像的恢复问题。

您可以使用以下代码获得该模型的一维示例。

负载cuspamax;y = cuspamax + 0.5 * randn(大小(cuspamax));情节(y);持有在;绘图(cuspamax,“r”,“线宽”轴紧;传奇(f (n) + \σe (n) ',‘f(n)’,“位置”,“西北”);

对于具有一定平滑特性的一大类函数(信号、图像),小波技术对于函数恢复来说是最优或接近最优的。

具体来说,该方法对函数族是有效的f只有少数几个非零小波系数。这些函数具有稀疏小波表示。例如,一个几乎处处光滑的函数,只有很少的突变,就具有这样的性质。

一般的基于小波的去噪和非参数函数估计方法是将数据变换到小波域,对小波系数进行阈值化,然后进行逆变换。

您可以将这些步骤总结为:

  1. 分解

    选择一个小波和一个级别N.计算信号的小波分解到水平N

  2. 阈值的细节系数

    对于从1到1的每个级别N,对细节系数设置阈值。

  3. 重建

    利用原始水平近似系数计算小波重构N以及从1到3个级别的修改细节系数N

去噪方法

小波工具箱支持多种去噪方法。文中实现了四万博1manbetx种去噪方法thselect.每个方法对应一个tptr命令中的选项

用力推= thselect (y, tptr)

返回阈值。

选项

去噪方法

“rigrsure”

使用Stein无偏风险估计原则进行选择(SURE)

“sqtwolog”

固定形式(通用)阈值等于

2 ln ( N )

具有N信号的长度。

“heursure”

使用前两个选项的混合选择

“minimaxi”

选择使用极大极小原理

  • 选项“rigrsure”软阈值估计器使用基于Stein无偏风险估计(二次损失函数)的阈值选择规则。您可以获得特定阈值的风险估计t.将风险降到最低t给出阈值的选择。

  • 选项“sqtwolog”使用一个固定形式的阈值产生的极小极大的性能乘以一个小的因素成正比日志(长度()).

  • 选项“heursure”是前两个选项的混合。因此,如果信噪比非常小,则当然,估计是非常嘈杂的。因此,如果检测到这种情况,将使用固定形式的阈值。

  • 选项“minimaxi”使用选定的固定阈值来生成对均方误差的极小极大性能的理想程序。在统计学中,极小极大原理被用来设计估计量。由于去噪信号可以被同化为未知回归函数的估计器,极小极大估计器是实现在给定函数集上的最大均方误差最小的选项。

下面的例子展示了一个1000 × 1 N(0,1)向量的去噪方法。这里的信号是

f ( n ) + e ( n ) e ( n ) ~ N ( 0 , 1 )

具有f (n)= 0.

rng违约; sig=randn(1e3,1);thr_rigrsure=thselect(信号,“rigrsure”)thr_univthresh=thselect(信号,“sqtwolog”) thr_heursure = thselect(sig,“heursure”) thr_minimaxi = thselect(sig,“minimaxi”)直方图(sig);h=findobj(gca,“类型”,“补丁”);集(h,“FaceColor”,[0.7 0.7 0.7],“EdgeColor”,“w”);持有在;绘图([thr_rigrsure thr_rigrsure],[0 300],“线宽”2);Plot ([thr_univthresh thr_univthresh], [0 300],“r”,“线宽”2);Plot ([thr_minimaxi thr_minimaxi], [0 300],“k”,“线宽”2);Plot ([-thr_rigrsure -thr_rigrsure], [0 300],“线宽”,2); 绘图([-thr_univthresh-thr_univthresh],[0 300],“r”,“线宽”,2); 绘图([-thr_minimaxi-thr_minimaxi],[0 300],“k”,“线宽”2);

对于Stein的无偏风险估计(SURE)和minimax阈值,保留大约3%的系数。对于通用阈值,拒绝所有值。

我们知道细节系数向量是系数的叠加f以及e,以及e导致细节系数,这是标准的高斯白噪声。

你使用后thselect要确定阈值,您可以设置每个级别的阈值wthcoef,直接处理原始信号的小波分解结构

软阈值或硬阈值

硬阈值和软阈值是收缩规则。确定阈值后,必须决定如何将该阈值应用于数据。

最简单的方案是坚硬的阈值。让T表示阈值x您的数据。硬性阈值是

η ( x ) = { x | x | T 0 | x | < T

软阈值为

η ( x ) = { x T x > T 0 | x | T x + T x < T

可以使用硬规则或软规则应用阈值wthresh

y = linspace (1100);用力推= 0.4;ythard = wthresh (y,“h”,用力推);ytsoft = wthresh (y,'s',用力推);次要情节(131);情节(y);标题(“原始数据”);次要情节(132);情节(ythard‘*’);标题(“硬阈值”);次要情节(133);情节(ytsoft‘*’);标题(“软阈值”);

处理无比例噪声和非白噪声

通常在实践中,基本模型不能直接使用。我们在这里检查可用于处理主要去噪功能中的模型偏差的选项wdenoise

最简单的用法wdenoise

sd=wdenoise(s)
哪个返回去噪版本sd原始信号的使用默认的参数设置,包括小波,去噪方法,软或硬阈值。任何默认设置都可以更改:
sd=wdenoise(s,n,'DenoisingMethod',tptr,'Wavelet',wav,'ThresholdRule',sorh,'NoiseEstimate',scal)
哪个返回去噪版本sd原始信号的使用tptr去噪方法。需要的其他参数有sorh,鳞片,wname.参数sorh指定在级别上分解的细节系数的阈值n通过小波wav.其余的参数鳞片是要具体说明的。它对应于估计数据中噪声方差的方法。

选项

噪声估计方法

“水平独立”

“水平独立”基于最细尺度(最高分辨率)小波系数估计噪声的方差。

“与级别相关”

“与级别相关”基于每个分辨率的小波系数估计噪声的方差。

对于更一般的程序,wdencmp函数执行小波系数阈值化,用于去噪和压缩目的,同时直接处理一维和二维数据。它允许您在中定义自己的阈值化策略

xd = wdencmp(选择x, wav, n,用力推,sorh, keepapp);

哪里

  • 选择= ' gbl ('苏氨酸是统一阈值的正实数。

  • 选择= ' lvd '苏氨酸是与级别相关的阈值的向量。

  • keepapp = 1保持近似系数,如前所述和

  • keepapp = 0允许近似系数阈值化。

  • x是否要对信号进行去噪,以及wav,n,sorh和上面一样。

小波去噪的作用

我们从Donoho和Johnstone的第一个示例开始一维去噪方法的示例。

块信号阈值

首先设置信噪比(SNR)并设置随机种子。

sqrt_snr=4;初始值=2055615866;

产生原始信号xref还有一个嘈杂的版本x加入标准高斯白噪声。同时绘制信号。

[xref x] = wnoise (sqrt_snr 1, 11日,init);次要情节(2,1,1)情节(xref)轴紧标题(原始信号的)子图(2,1,2)图(x)轴紧标题(噪声信号的)

图中包含2个轴对象。标题为“原始信号”的轴对象1包含线条类型的对象。标题为“噪音信号”的轴对象2包含线条类型的对象。

对小波分解得到的细节系数进行软启发式阈值去噪x使用sym8小波。使用的默认设置wdenoise查看其余参数。与原始信号进行比较。

xd=wdenoise(x,“小波”,“sym8”,“DenoisingMethod”,“确定”,“阈值规则”,“软的”);图subplot(2,1,1) plot(xref)轴紧标题(原始信号的)子地块(2,1,2)绘图(xd)轴紧标题(的去噪信号)

图中包含2个轴对象。标题为原始信号的轴对象1包含line类型的对象。标题为去噪信号的轴对象2包含line类型的对象。

由于只有少量的大系数表征原始信号,因此该方法的性能非常好。

电信号去噪

当您怀疑非白噪声时,必须根据电平噪声的电平相关估计重新调整阈值。作为第二个例子,让我们在电信号的高扰动部分尝试这种方法。

首先加载电信号并从中选择一个片段。绘制部分。

负载leleccumindx=2000:3450;x=leleccum(indx);图形绘制(indx,x)轴紧标题(原始信号的)

图中包含一个轴对象。标题为“原始信号”的轴对象包含类型为line的对象。

去噪信号使用db4小波和三层小波分解和软固定形式阈值。为了处理非白噪声,使用电平相关噪声大小估计。与原始信号进行比较。

xd = wdenoise (x 3“小波”,“db4”,...“DenoisingMethod”,“UniversalThreshold”,...“阈值规则”,“软的”,...“噪音估计”,“与级别相关”); 图子地块(2,1,1)图(indx,x)轴紧标题(原始信号的)子地块(2,1,2)绘图(indx,xd)轴紧标题(的去噪信号)

图中包含2个轴对象。标题为原始信号的轴对象1包含line类型的对象。标题为去噪信号的轴对象2包含line类型的对象。

尽管在时间2410前后传感器故障开始前后噪声性质的时间异质性,结果还是相当好的。

扩展到图像去噪

针对一维情况描述的去噪方法也适用于图像,并且很好地适用于几何图像

( , j ) = f ( , j ) + σ e ( , j )

哪里 e 是单位方差的高斯白噪声。

二维去噪过程采用相同的三个步骤,使用二维小波工具代替一维小波工具。对于阈值的选择,产品(尺寸)被用来代替长度(s)如果使用固定形式阈值。

请注意,除“自动”1-D去噪情况外,2-D去噪和压缩使用wdencmp.为了说明二维去噪,加载一个图像并创建一个有噪声的版本。为了重现性,设置随机种子。

初始值=205565866;rng(初始);负载女人img=X;imgNoisy=img+15*randn(大小(img));

使用ddencmp在这种情况下,使用固定形式的阈值对水平噪声进行估计,阈值是软的,并且保留了近似系数。

[thr,sorh,KEEPAP]=ddencmp(“书房”,西弗吉尼亚州的, imgNoisy);苏氨酸
thr=107.9838

苏氨酸等于estimated_sigma *√日志(prod(大小(img))))

使用全局阈值选项去噪噪声图像。显示结果。

imgDenoised=wdencmp(gbl(的,不知道,“sym4”,2,thr,sorh,keepapp);figure colormap(粉红色(255))sm=尺寸(图1);子地块(2,2,1)图像(wcodemat(img,sm))标题(原始图像的)子地块(2,2,2)图像(wcodemat(imgNoisy,sm))标题(“噪音图像”)子地块(2,2,3)图像(wcodemat(imgDenoised,sm))标题(”“去噪图像)

图中包含3个轴对象。标题为原始图像的轴对象1包含图像类型的对象。标题为噪波图像的轴对象2包含图像类型的对象。标题为去噪图像的轴对象3包含图像类型的对象。

去噪后的图像与原始图像对比良好。

一维小波方差自适应阈值

其思想是逐级定义时间相关阈值,然后提高去噪策略处理非平稳方差噪声模型的能力。

更准确地说,模型假设(如前所述)观测值等于叠加在噪声上的感兴趣信号

( n ) = f ( n ) + σ e ( n )

但噪声方差会随时间而变化。在几个时间间隔上有几个不同的方差值。值和间隔都是未知的。

让我们关注估计变化点或等效间隔的问题。所使用的算法基于Marc Lavielle关于使用动态规划检测变化点的原著(参见中的[Lav99])参考文献).

让我们从固定设计回归模型中生成一个信号,其中两个噪声方差变化点位于位置200和600处。为了再现性,设置随机种子。

init = 2055615866;rng (init);x = wnoise (10);bb = randn(1、长度(x));cp1 = 200;cp2 = 600;x = x + (bb (1: cp1), bb (cp1 + 1: cp2) / 4, bb (cp2 + 1:结束)];情节(x)标题(噪声信号的)

图中包含一个轴对象。标题为Noise Signal的轴对象包含一个line类型的对象。

这个例子的目的是从信号中恢复两个变化点x

步骤1. 通过抑制近似值来恢复噪声信号。首先使用db4小波。然后在第1级重建细节。

wname =“db4”;列弗= 1;[c、l] = wavedec (x,列弗,wname);侦破= wrcoef (“d”,c,l,wname,1);地物图(det)标题(“1级细节”)

图中包含一个轴对象。标题为Level 1 Detail的axes对象包含一个类型为line的对象。

在这一阶段恢复的1级重建细节几乎没有信号。如果信号的感兴趣部分具有稀疏小波表示,则它从变化点检测的角度捕获噪声的主要特征。为了去除几乎所有信号,我们用平均值替换最大值。

步骤2。要去除几乎所有的信号,用平均值替换最大值的2%。

x = (abs(检波器));v2p100 = x(修复(长度(x) * 0.98));印第安纳州=找到(abs(检波器)> v2p100);依据(印第安纳州)=意味着(依据);

第三步。使用wvarchg函数使用以下参数估计变化点:

  • 两个改变点之间的最小延迟为d = 10。

  • 最大更改点数为5。

[cp_est,kopt,t_est]=wvarchg(第5段)
cp_est =1×2259 611
kopt = 2
东东=6×61024 0 0 0 0 0 612 1024 0 0 0 0 259 611 1024 0 0 0 198 259 611 1024 0 0 198 235 259 611 1024 0 198 235 260 346 611 1024

提出了两个变化点和三个间隔。由于噪声的三个区间方差是非常不同的,优化程序很容易检测到正确的结构。估计的变化点与真实的变化点很接近:200和600。

第四步。(可选)替换估计的更改点。

两个 6.t_est(1张):包含方差变化点的瞬间,以及自科普特是建议的变更点数量,然后

cp_est=t_est(kopt+1,1:kopt);

您可以通过计算来替换估计的更改点:

对于cp_New = t_est(k+1,1:k)结束
cp_New=612
新科大=1×2259 611
新科大=1×3198 259 611
新科大=1×4198 235 259 611
新科大=1×5198 235 260 346 611

小波去噪分析

下面的测量和设置对于分析小波信号和图像是有用的:

  • M S E-均方误差(MSE)是数据与信号或图像近似之间的差除以元素数的平方范数。

  • 最大的错误-信号或图像近似中的最大绝对平方偏差。

  • L2范数比-信号或图像近似的平方l2范数与输入信号或图像的比率。对于图像,图像在取l2范数之前被重塑为列向量

  • P S N R-峰值信噪比(PSNR),单位为分贝。PSNR仅对以每样本位或每像素位编码的数据有意义。

  • B P P-比特/像素比(BPP),即压缩比(压缩比)乘以8,假设每个像素一个字节(8比特)。

  • 补偿比-压缩比,即压缩图像中的元素数除以原始图像中的元素数,以百分比表示。