小波去噪及其MATLAB中的函数

小波去噪及其MATLAB中的函数一、小波去噪中信号阈值的估算信号去噪是信号处理领域的经典问题之一。传统的去噪方法主要包括线性滤波方法和非线性滤波方法，如中值滤波和wiener滤波等。传统去噪方法的不是在于使信号变换后的熵增高、无法刻画信号的非平稳特性并且无法得到信号的相关性。为了克服上述缺点，人们开始使用小波变换解决信号去噪问题。小波变换具有下列良好特性：(1)低熵性：小波系数的稀疏分布，使信号变换后的熵降低；(2)多分辨率特性：可以非常妤地刻画信号的非平稳特性，如边缘、尖峰、断点等；(3)去相关性：可取出信号的相关性，且噪声在小波变换后有白化趋势，所以比时域更利于去噪；(4)选基灵活性：由于小波变换可以灵活选择基函数，因此可根据信号特点和去噪要求选择适合小波。小波在信号去噪领域已得到越来越广泛的应用。阈值去噪方法是一种实现简单、效果较好的小波去噪方法。阈值去噪方法的思想就是对小波分解后的各层系数中模大于和小于某阈值的系数分别处理，然后对处理完的小波系数再进行反变换，重构出经过去噪后的信号。下面从阈值函数和阈值估计两方面对阈值去噪方法进行介绍。1．阈值函数常用的阈值函数主要是硬阈值函数和软阈值函数。(1)硬阈值函数。表达式为)|(|)(TwwIw，如图4.18所示，其中横坐标表示信号的原始小波系数，纵坐标表示阈值化后的小波系数。(2)软阈值函数。表达式为)|(|))sgn(()(TwIT，如图4.19所示，其中横坐标表示信号的原始小波系数，纵坐标表示阈值化后的小波系数。一般来说，硬阈值方法可以很好地保留信号边缘等局部特征，软阈值处理相对要平滑，但会造成边缘模糊等失真现象。为了克服上述缺陷，最近提出了一种半软阈值函数，如图4.20所示。它可以兼顾软阈值和硬阈值方法的优点，其表达式为)|(|)||()|(|)sgn()(2211212TwwITwTITTTwTww其中。210TT在软阈值的基础上，可以对其改进使其具有更高阶，如图4.21所示。可以看出它在噪声(小波系数)与有用信号(小波系数)之间存在一个平滑过渡区，更符合自然信号/图像的连续特性。其表达式为TwkTTw12TwkTTwTwwTkwkk12||)12(1)(1222．阈值估计Donoho在1994年提出了VisuShrink方法(或称统一阈值去噪方法)。它是态变量联合分布，在维数趋向无穷时得出的结论，在最小最大估针对多维独立正计的限制下得出的最优阈值。阈值的选择满足：NTnln2其中n，是噪声标准方差，N是信号的长度。Donoho给出了证明这种估计在信号属于Besov集时，在大量风险函数下获得近似理想的去噪风险。Donoho的统一阈值方法在实际应用中效果欠理想，产生过扼杀现象，1997年Janse提出了基于无偏估计的阈值计算法。风险函数定义为：NfftR/ˆ)(2。由于小波变换的正交性，风险函数可NXYtRt/)()(2以同样在小波域中写成形式。设NYYtRt/)()(2，则YXXYYXXYENYVENtERYYNtETtttnt,)(2)(1)(,2)()(1)(2222最后可得到风险函数的表达式：tYINtYNtYINtETtERiNinniNiiNinn||2||1||2)()(12221122其中是示性函数，I为两数取小。于是，最佳的阈值选择可以通过最小化风险函数得到，即)(minarg0*tERttMATLAB中实现了信号的阈值去噪，主要包括阈值获取和阈值去噪两方而。下面对它们进行说明。二、小波去噪在MATLAB中的函数1)阈值获取MATLAB中实现信号阈值获取的函数有ddencmp、thselect、wbmpen和wdcbm，下面对它们的用法进行简单的说明。Ddencmp的调用格式有以下三种(1)[THR，SORH，KEEPAPP，CRIT]=ddencmp(IN1，IN2，X)(2)[THR，SORH，KEEPAPP，CRIT]=ddencmp(IN1，'wp'，X)(3)[THR，SORH，KEEPAPP]=ddencmp(IN1，'wv'，X)函数ddencmp用于获取在消噪或压缩过程中的默认阈值。输入参数X为一维或二维信号；IN1取值为'den'或'crop'，den表示进行去噪，crop表示进行压缩；IN2取值为'wv'或'wp'，wv表示选择小波，wp表示选择小波包。返回值THR是返回的阈值；SORH是软阈值或硬阈值选择参数；KEEPAPP表示保存低频信号；CRIT是熵名(只在选择小波包时用)。函数thselect的调用格式如下：THR=thselect(X，TPTR)THR=thselect(X，TPTR)根据字符串TPTR定义的阈值选择规则来选择信号X的自适应阈值。自适应阈值选择规则包括下面四种：·TPTR='rigrsure'，自适应阈值选择使用Stein的无偏风险估计原理。·TPTR='heursure'，使用启发式阈值选择。·TPTR='sqtwolog'，阈值等于sqrt(2*log(1ength(X)))。·TPTR='minimaxi'，用极大极小原理选择阈值。阈值选择规则基于模型etfy)(，是高斯A噪声N(O，1)。e函数wbmpen的调用格式如下：THR=wbmpen(C，L，SIGMA，ALPHA)THR=wbmpen(C，L，SIGMA，ALPHA)返回去噪的全局阈值THR。THR通过给定的一种小波系数选择规则计算得到，小波系数选择规则使用Birge-Massart的处罚算法。[C，L]是进行去噪的信号或图像的小波分解结构；SIGMA是零均值的高斯白噪声的标准偏差；ALPHA用于处罚的调整参数，它必须是一个大于1的实数，一股取ALPHA=2。设t*是crit(t)=-sum(c(k)^2，k=t)+2*SIGMA^2*t*(ALPHA+log(n/t))的最小值，其中c(k)是按绝对值从大到小排列的小波包系数，n是系数的个数，则THR=c(t*)。wbmpen(C，L，SIGMA，ALPHA，ARG)计算阈值并画出三条曲线。·2*SIGMA^2*t*(ALPHA+10g(n/t))·Sum(c(k)^2，k=t)·crit(t)函数wdcbm的调用格式有以下两种：(1)[THR，NKEEP]=wdcbm(C，L，ALPHA)(2)[THR，NKEEP]=wdcbm(C，L，ALPHA，M)函数wdcbm用于使用Birge-Massart算法获取一维小波变换的阈值。返回值THR是与尺度无关的阈值，NKEEP是系数的个数。[C，L]是要进行消噪或压缩的信号在j=length(L)-2层的分解结构；ALPHA和M必须是大于1的实数；THR是关于j的向量，THR(i)是第i层的阈值；NKEEP也是关于j的向量，NKEEP(i)是i层的系数个数。一般压缩时ALPHA取1.5，去噪时ALPHA取3。2)信号的阈值去噪MATLAB中实现信号的阈值去噪的函数有wden、wdencmp、wthresh、wthcoef、wpthcoef以及wpdencmp。下面对它们的用法进行简单的介绍。函数wden的调用格式有以下两种：(1)[XD，CXD，LXD]=wden(X，TPTR，SORH，SCAL，N，'wname')(2)[XD，CXD，LXD]=wden(C，L，TPTR，SORH，SCAL，N，'wname')函数wden用于一维信号的自动消噪。X为原始信号，[C，L]为信号的小波分解，N为小波分解的层数。TPTR为阈值选择规则，TPTR的取值有以下四种：·TPTR='rigrsure'，采用Stein无偏似然估计。·TPTR='heursure'，采用启发式阈值选择。·TPTR='sqtwolog'，取通用阈值)log(2。·TPTR='minimaxi'，采用极大极小值进行阈值选择。SORH是软阈值或硬阈值的选择(分别对应's'和'h')。SCAL指所使用的阈值是否需要重新调整，包含下而三种：·SCAL='one'，不调整。·SCAL='sln'，根据第一层的系数进行噪声层的估计来调整阈值。·SCAL='mln'，根据不同层的噪声估计来调整阈值。XD为消噪后的信号，[CXD，LXD]为消噪后信号的小波分解结构。格式(1)返回对信号X经过N层分解后的小波系数进行阈值处理后的消噪信号XD和信号XD的小波分解结构[CXD，LXD]。格式(2)返回参数与格式(1)相同，但其结构是由直接对信号的小波分解结构[C,L]进行阈值处理得到的。函数wdencmp的调用格式有下面三种：(1)[XC，CXC，LXC，PERF0，PERFL2]=wdenemp('gbl’，X，'wname'，N，THR，SORH,KEEPAPP)(2)[XC，CXC，LXC，PERF0，PERFL2]=wdencmp('1vd'，X，'wname'，N，THR，SORH)(3)[XC，CXC，LXC，PERF0，PERFL2]=wdencmp('1vd’，C，L，'wname'，N，THR，SORH)函数wdencmp用于一维或二维信号的消噪或压缩。wname是所用的小波函数，gbl(global的缩写)表示每层都采用同一个阈值进行处理，lvd表示每层用不同的阈值进行处理，N表示小波分解的层数，THR为阈值向量，对于格式(2)和(3)每层部要求有一个阈值，因此阈值向量THR的长度为N，SORH表示选择软阈值或硬阈值(分别取值为's'和'h)，参数KEEPAPP取值为1时，则低频系数不进行阈值量化，反之，则低频系数要进行阈值量化。XC是消噪或压缩后的信号，[CXC,LXC]是XC的小波分解结构，PERF0和PERFL2是恢复和压缩的范数百分比。如果[C，L]是X的小波分解结构，则PERFL2=100*(CXC向量的范数/C向量的范数)2；如果X是一维信号，小波wname是一个正交小波，则PERFL2=2L22/100XXC。函数wthresh的调用格式如下：Y=wthresh(X，SORH，T)Y=wthresh(X，SORH，T)返回输入向量或矩阵X经软阈值(如果SORH='s')或硬阈值(如果SORH='h')处理后的信号。T是阈值。Y=wthresh(X，'s'，T)返回的是Y=)|(|)(TXXSING，即把信号的绝对值与阈值进行比较，小于或等于阈值的点变为0，大于阈值的点变为该点值与阈值的差值。Y=_wthresh(X，'h'，T)返回的是Y=)|(|1TXX，即把信号的绝对值与阈值比较，小于或等于阈值的点变为0，大于阈值的点保持不变。一股来说，用硬阈值处理后的信号比软阈值处理后的信号更粗糙。函数wpthcoef的调用格式如下：T=wpthcoef(T，KEEPAPP，SORH，THR)NT=wpthcoef(T，KEEPAPP，SORH，THR)通过对小波包树T的系数进行阈值处理后返回一个新的小波包树NT。如果KEEPAPP=1，则细节信号的系数不进行阈值处理；否则，就要进行阈值处理。如果SORH='s'，使用软阈值，如果SORH='h'，则使用硬阈值。THR是阈值。函数wthcoef的调用格式有下面四种：(1)NC=wthcoef('d'，C，L，N，P)(2)NC=wthcoef('d'，C，L，N)(3)NC=wthcoef('a'，C，L)(4)NC=wthcoef('t'，C，L，N，T，SORH)函数wthcoef用于一维信号小波系数的阈值处理。格式(1)返回小波分解结构[C，L]经向量N和P定义的压缩率处理后的新的小波分解向量NC，[NC，L]即构成一个新的小波分解结构。N包含要被压缩的细节向量，P是把较小系数置0的百分比信息的向量。N和P的长度必须相同，向量N必须满足1≤N(i)≤length(L)-