《武汉工程大学学报》 2013年04期
61-65
出版日期：2013-04-30
ISSN:1674-2869
CN:42-1779/TQ
改进的二维经验模式分解方法
0引言经验模式分解（Empirical mode decomposition, EMD\[1\]）是由Norden Huang等人于1998年提出的，它是一种全新的用于处理非平稳和非线性信号的方法.EMD可以将信号分解为一些内在模式函数（IMF）和表示信号变化趋势的残余量，其中IMF需要满足两个条件：（1）在完整信号序列中，信号的极值点（包括极大值和极小值）的数量必须和过零点的数量相等或是相差不多于一个；（2）在信号的每个点，它形成的局部极大值和局部极小值的包络均值为零.这一序列的IMF包含了信号从高频到低频的不同频率段的成分.J.C.Nunes等人在2003年将一维EMD扩展到二维，提出了二维经验模式分解（Bidimensional Empirical mode decomposition, BEMD\[2\]）方法，并应用在图像纹理分析中.二维经验模式分解被国内外学者广泛地应用于图像处理中，比如应用在图像压缩、图像去噪\[3\]、图像缩放、图像轮廓提取等方面.由于经验模式分解在算法方面还不够完善，学者们对它进行了一系列的改进，比如对边界效应的处理、插值方法的选取、极值点的查找以及筛分停止条件的选取等等.黄大吉等人提出了镜像闭合延拓法\[4\]，通过镜像对称对信号进行延拓，使信号不存在端点，从而避免了边界效应，但是此方法用于二维图像信号时，使图像的数据量扩大8倍，使算法的时间复杂度增大.本文的实现过程是针对边界效应做出的改进方法.1BEMD的实现过程二维经验模式分解被广泛应用于图像处理中，它可以将将一副图像分解为若干表示图像不同频率的内在模式函数和表示图像变化趋势的残余量.BEMD的实现步骤如下：（1） 初始化：输入二维图像f(x,y)，令r1,1(x,y)=f(x,y)，r1,1(x,y)作为待处理的图像.（2） 分别计算rl,m(x,y)的局部极大值和局部极小值，其中l表示由图像中分解出的第l个IMF，m表示第m次分解（m=1…M）.（3） 分别对局部极大值点和局部极小值点进行曲面插值包络，得到上包络面emaxl,m(x,y)和下包络面eminl,m(x,y).（4） 确定这些包络的均值曲面meanl,m(x,y)：meanl,m(x,y)=emaxl,m(x,y)+eminl,m(x,y)2（5） 计算输入图像rl,m(x,y)与包络均值meanl,m的差值得到hl,m(x,y)：hl,m(x,y)=rl,m(x,y)-meanl,m(x,y)（6） 判断hl,m(x,y)是否符合IMF的筛分终止条件SD，若符合则判定hl,m(x,y)是第l个IMF；若不满足，则令rl,m+1(x,y)=hl,m(x,y)，返回步骤（2），重复执行以上过程，直到经过M次分解后图像信号满足条件为止，令m=M，故筛分出第l个IMF为Cl(x,y)=hl,M(x,y)，其中SD为SD=∑x=Ix=0∑y=Jy=0|rl,m(x,y)-rl,m+1(x,y)|2r2l,m(x,y)（7） 计算残余量rl+1,1(x,y)：rl+1,1(x,y)=rl,1(x,y)-Cl(x,y)（8） 若残余量中含有两个以上的极值点，继续返回步骤（1）分解图像并把残余量rl+1,1(x,y)作为新的输入图像，直到图像分解出的残余量rL(x,y)没有极值点为止，结束整个BEMD过程.可以将原图像信号表示为f(x,y)=∑l=Ll=1Cl(x,y)+rL(x,y)2算法改进2.1BEMD的边界效应在BEMD分解中首先求局部极值点，然后基于Delaunay三角剖分利用三次样条函数，通过插值得到极大值包络线和极小值包络线，由于端点部分不一定是极值点，在端点处不满足三次样条插值函数的条件，使得分解出的内在模式函数在端点处出现虚假的成分，并且随着筛分过程的进行，这些虚假的成分会逐渐向内影响整个数据，最终使得分解结果失真，即产生了“边界效应”，信号的频率越小，其边界效应越明显，对结果的危害性越大.针对边界效应这一问题，很多学者提出了解决的方法，其中Huang 本人提出了特征波法\[5\],还有支持矢量回归机法\[6\]，自回归模型法\[7\]，多项式拟合法\[8\],波形镜像法\[9\]等改进算法的提出.其中文献\[10\]比较了几种边界处理的方法，表明镜像法是相对较好的边界处理方法.第4期张彦铎，等：改进的二维经验模式分解方法武汉工程大学学报第35卷2.2对边界效应的处理处理方法是使用图像数据信号的边缘延拓，即在信号的边界处增加一部分数据，这些数据是相对于边界数据向外的平移，并不像镜像法那样将信号扩大8倍，而是只增加了边界处的极少数据，这样相当于在图像的周围形成保护圈把图像边界保护起来，让边界效应不至于污染到原始图像信号，这样就避免了边界效应的产生.此方法操作简单，而且增加的数据量不大，所以可以提高计算速度.信号数据延拓的方法如图1、图2所示.图1原数据Fig.1Source data图2延拓后的数据Fig.2Data of extension2.3筛分停止准则Huang提出了确定筛分过程是否结束的准则，即使通过限制标准差的大小来实现，标准差(Standard Deviation,SD)通过两次连续的处理结果计算得出：SD=∑x=Ix=0∑y=Jy=0|rl,m(x,y)-rl,m+1(x,y)|2r2l,m(x,y)SD称为筛分门限，它确定筛分过程是否截止，故它的取值直接影响BEMD分解的质量.SD的阈值选择还没有严格的理论依据，都是建立在大量实验的基础上取得经验值，一般取值在0.2和0.3之间，而且对同一个图像信号都是取得相同的SD阈值.若阈值设置过大，这就表示结束条件宽松而使得提取的IMF质量差；若SD设置过小，这就表示结束条件过于严格，筛分次数增多，就会出现图像过度筛分的情况.文献\[11\]提出分解每个IMF的合适阈值并不相同，它的解决方法是从筛分过程中提取极值点曲面的空间位置信息，根据这些信息判断是否应该停止筛分过程.此方法比固定阈值有一定优势，但是实现较复杂.针对SD阈值相对小从而导致图像分解过度的情况提出了解决方法，即是对筛分处理过度的IMF加一定的补偿量，每个IMF加上不同的补偿量，根据实验经验得出这个补偿量的值是每个IMF最小值的相反数，即imf=imf-min(imf)；具体实现是先根据算法步骤完成每个IMF的筛分过程，在筛分截止后对每个IMF进行处理，计算出它的最小值，然后将这个最小值的绝对值作为补偿量加到原来的IMF中，即是对IMF的图像信息作了一部分的还原.3实现过程和图像处理结果将Lena（256×256）灰度图像作为原始数据源，首先对图像按以上方法进行边界延拓，其次利用八邻域法求图像的局部极值点，然后使用基于Delaunay三角剖分的三次样条函数进行插值形成上下包络曲面，接下来做差值，判断如果满足条件，则第一个IMF分解完成，最后对IMF进行值的调整，下一个图像分解的输入图像要去掉已分解完成的上一个IMF.按以上的步骤循环下去，直到图像残余量的极值点少于两个就全部结束.实验结果如图3~5所示，图3是IMF1、图4是IMF2、图5是IMF3，图6是残余量、图7是重构图像，图8是原图像，重构图像是3个IMF与残余量相加得到的，从图7和图8可以看出重构图和原图吻合得很好，用图像标准差\[12\]来计算重构图与原图的相似程度,标准差S为S=1I×J∑Ii=1∑Jj=1\[fres(i,j)-fn(i,j)\]212其中I表示图像像素的行数，J表示图像像素的列数，fres(i,j)表示重构图像矩阵的第i行j列像素值，fn(i,j)表示原图像矩阵的第i行j列像素值.通过计算可以得到标准差为1541×10-6.所以重构图像与原始图像之间的灰度波动很小，由此可见，改进的方法能解决边界问题和分解过度问题.分别使用镜像延拓法和基于Delaunay三角剖分对称法来解决边界效应，Lena（256×256）灰度图像作为原始数据源，得到的结果分析如表1所示.表1标准差对比Table 1The comparison of standard deviation方法镜像延拓法基于Delaunay三角剖分对称法本文方法标准差5473×10-56667×10-61541×10-6 图3第一层Fig.3IMF1图4第二层Fig.4IMF2图5第三层Fig.5IMF3图6残余量Fig.6Residual图7重构图Fig.7Reconstruction maps图8Lena灰度图Fig.8Lena gray image4结语以上对图像处理中二维经验模式分解算法提出了改进，在分解过程中需要对像素极值点的选取以及对极值点进行插值形成包络曲面，这些步骤会造成边界效应，本文应用了一种简单的边界延拓方法，避免了边界效应的产生，与其它处理边界效应问题的方法比较最显著的特点是简单易行，附加的图像边缘数值少且计算量小，所以计算使用的时间相对来说较短，时间在5 s到6 s之间，虽然还达不到实时处理系统的要求，但是速度也有了提高，今后的努力方向是进一步地提高速度；采用给每个IMF加补偿量的方法来解决筛分停止阈值小的情况下产生的分解过度问题，计算出来的图像标准差很小，所以证实了重构图与原始图像吻合得很好.因此此方法能在基于BEMD的图像处理［13］中起到积极的作用. 致谢本文在撰写的过程中，得到了闵峰博士在搜集资料和实验数据采集方面的大力支持，特此向他表示感谢