二进小波与扩散滤波结合的光学相干层析图像降噪

张 田，孙延奎，田小林

(1.清华大学计算机科学与技术系，北京100084;2.澳门科技大学资讯科技学院，澳门特别行政区)

其中q0为待处理图像中的斑点方差系数，ρ是常量，t=n×Δt，n是迭代次数，Δt是迭代步长。

与离散小波变换相比，二进小波变换具有平移不变性，因此，二进小波在图像降噪中具有重要应用［10］。二进小波变换具有尺度相关性，即图像重要边缘的小波系数会分布在多个尺度中，而图像噪声的小波系数主要集中在最低尺度中;特别的，图像重要边缘处的小波系数幅值随着尺度升高不会明显减小;而图像噪声的小波系数幅值随着尺度升高会急剧衰减。

光学相干层析成像技术［1］(OCT:Optical Coherence Tomography)是一种基于弱相干光的干涉原理，通过检测不同组织对入射的弱相干光的背向反射或散射信号得到生物组织的二维图像或三维结构的光学成像技术。因为其具有成像快、分辨率高、非侵入式等优点，自从上世纪90年代以来，OCT在医学领域获得了广泛应用，尤其在眼科疾病诊断方面，OCT系统获取的高清晰(μm量级)视网膜图像对于青光眼等眼科疾病的诊断提供了极大的帮助［2］。

但在OCT成像设备的实际应用中，光源、检测电路、振镜都会带来噪声，此外，光的多次散射也会产生很多散斑，这些都严重降低了图像质量，给OCT成像设备的临床应用带来了不利影响。OCT图像中存在的多种噪声中，由于使用的光源的相干作用产生的斑点噪声［3］是最主要的部分，这种噪声以乘法噪声的形式存在，其数学模型

乘法噪声的去除对降噪方法提出了新的要求，多年来国内外研究者对这个问题进行了大量研究，取得了众多成果。一类是对原图像进行一次对数运算，将乘法噪声转化为加法噪声后，再应用如BiShrink［4］，ProbShrink［5］，曲波(curvelet)系数收缩方法［6］等传统加法噪声去除方法。另一类是利用斑点噪声的特性开发的专用于去除乘法斑点噪声的方法，如利用归一化小波系数的小波扩散降斑方法［7］，斑点抑制各向异性扩散(SRAD: Speckle Reduction Anisotropic Diffusion)［8］等。SRAD算法是在传统各向异性扩散滤波方法［9］的基础上，通过设计一种专用于检测乘法斑点噪声的边缘检测算子，取得了较好的抑制斑点噪声的效果。但在对视网膜OCT图像的降噪实验中发现，SRAD算法虽然具有较强的去除斑点噪声的能力，但仍有很多较强的噪声会被误判为边缘从而被保留，使得SRAD降噪后的图像中仍存在较大量噪声。而小波变换具有很强的分离噪声和边缘的能力，而且图像的噪声和边缘的小波系数幅值在不同尺度间具有不同的衰减性。因此通过结合SRAD的降噪原理和小波系数尺度间的相关性，笔者提出一种改进的SRAD算法，在保持SRAD增强边缘的能力的同时，进一步提高其降噪能力。

1 算法描述

笔者利用OCT图像的对数图像的二进小波系数幅值的尺度相关性对边缘还是均质区域的指示作用，改进SRAD算法中的扩散系数计算公式，从而增强其保持边缘并去除斑点的能力。因为新方法结合了二进小波和SRAD算法的优点，故称为 DW－SRAD(Dyadic Wavelet Guided Speckle Reducing Anisotropic Diffusion)算法。

1.1 斑点抑制扩散滤波算法［8］

SRAD算法是在各向异性扩散滤波算法的基础上针对超声图像中普遍存在的乘性斑点噪声开发的。其具体过程如下:

对于一张能量有限的带有斑点噪声的图像I0(x，y)，令其支撑区间为Ω，则可根据下列方程计算其每一步迭代扩散后得到的新图像:

其中∂Ω代表支撑区间Ω的边界，n为∂Ω的外法向量，c(q)是扩散系数，

其中q(x，y;t)为方差系数。

q0(t)为斑点度量函数。

其中q0为待处理图像中的斑点方差系数，ρ是常量，t=n×Δt，n是迭代次数，Δt是迭代步长。

1.2 二进小波变换及其尺度相关性

与离散小波变换相比，二进小波变换具有平移不变性，因此，二进小波在图像降噪中具有重要应用［10］。二进小波变换具有尺度相关性，即图像重要边缘的小波系数会分布在多个尺度中，而图像噪声的小波系数主要集中在最低尺度中;特别的，图像重要边缘处的小波系数幅值随着尺度升高不会明显减小;而图像噪声的小波系数幅值随着尺度升高会急剧衰减。

图像I(x，y)连续做2次二维二进小波变换后，可得到一个低频图像和4个高频图像(x，y)(x，y)(x，y)，(x，y)。在二维小波域，每一个点(x0，y0)在第j尺度都有一个二维影响锥(COI:Cone of Influence(x0，y0)［10－12］。令N2j/(x0，y0)表示影响锥(x0，y0)中各小波系数绝对值之和，即

定义:

称τ(x0，y0)为点(x0，y0)处二进小波系数尺度间的相关系数。二进小波变换的尺度相关性表明，较大的τ(x0，y0)意味着点(x0，y0)更有可能位于重要边缘上;较小的τ(x0，y0)意味着点(x0，y0)很可能位于均质区域。这表明，τ(x0，y0)能够起到指示点(x0，y0)是位于图像的重要边缘上还是位于噪声的均质(homogenous)区域的作用。利用τ值的特征，可改进SRAD中的扩散系数计算公式，从而进一步提高SRAD降低斑点噪声并保持OCT图像边缘的能力。

1.3 算法流程

笔者提出了一个结合各向异性扩散滤波算法和二进小波系数尺度间相关性的OCT图像降斑降噪方法。该方法利用OCT图像的对数图像在二进小波域的小波系数尺度间的相关系数，设计一种更适合OCT图像的扩散系数计算公式，从而起到更好的保留OCT图像重要边缘，同时有效去除斑点噪声的效果。具体操作过程如下:

1)计算二进小波系数尺度间的相关性系数。

对原图像I做一次对数运算，得到对数图像Ilog;对Ilog进行2次二维二进小波变换，得到一个低频图像S22(x，y)和4个高频图像(x，y)，(x，y)，(x，y)，(x，y)。再根据式(7)、式(8)计算每点处二进小波系数尺度间的相关性系数τ(x0，y0)。

2)设计新的扩散系数计算公式。

利用τ具有区分噪声和边缘的能力，引入一个扩散系数修正因子s(τ):定义s(τ)为τ的单调递减函数，使得在图像的边缘上，τ相对较大时s(τ)较小;而在图像的均质区域，τ较小时s(τ)较大。例如，可定义s(τ)的计算公式为:

其中K为一个常量，由实验给出。

由s(τ)修正SRAD的扩散系数计算公式c(q)得到一个新的扩散系数计算公式。

3)利用cnew(q，τ)代替SRAD中的扩散系数计算公式c(q)进行扩散滤波

与c(q)相比，cnew(q，τ)具有更好的定位边缘与噪声的能力，因而将cnew(q，τ)应用于图像的迭代扩散降噪过程，会具有更好的降斑降噪能力，即在避免对边缘区域过平滑的同时更有效地抑制斑点噪声。

2 实验结果

为验证所提出算法的有效性，笔者利用视网膜OCT图像进行实验。数据采集系统的有关参数为:光源中心波长 830 nm，有效探测深度3.4 mm，信噪比为51 dB，纵向分辨率8.5μm［13］。实验中，迭代步长Δt为0.25，ρ取1，迭代次数为50次，常数K为20。对15张视网膜OCT图像进行了降噪实验，其中一幅图像及其降噪后图像的对比结果如图1所示。

图1 降噪前后的效果对比图Fig.1 The denoising results of a retina OCT image by DW_SRAD method

图1(a)为原OCT图像，绿色框中为背景区域，红色框中为前景区域，用于计算4种降噪指标;图1(b)为图1(a)中蓝线标注的A－scan中的信号;图1(c)为DW_SRAD降噪后的图像;图1(d)为图1(c)中蓝线标注的A－scan中的信号。从图1中可看出，原图中大量噪声被有效的降低了，包括黑色背景区域中的离散噪声和高亮前景区域中的噪声。图1中(b)和(d)的对比能清楚地看到重要边缘信号得到了较好的保留。

为定量地评价算法降噪效果，该次实验引入了几个国际通行的指标:

1)信噪比SNR(Signal－to－Noise Ratio)［14］

用于衡量图像中信号和噪声强度:

其中μ为信号强度的期望值;σ为噪声的标准差。

2)等效视数 ENL(Effective Number of Looks)［14］

用于度量图像中由于噪声引入的扰动，即一块本应该是平滑区域但被噪声污染后的区域平滑性度量

其中t代表图像中的前景区域，即图1中红色框中的区域，μt为该区域中的像素均值，σt为该区域中像素标准差。

3)对比度噪声比率 CNR(Contrast－to－Noise Ratio)［14］

提供图像特征区域(感兴趣区域)与背景区域(例如空气区域)的有效对比度的一个客观度量，计算公式:

其中符号t，μt和σt的含义与式(12)中的相同。r表示图像中的背景参考区域，即图1中的绿色框中的区域，μr和σr分别为该区域中的信号均值与标准差。

4)边缘保持参数β［15］

用于度量降噪后的图像边缘锐度降低的程度。

其中I和Id分别为原OCT图像和降噪后的OCT图像，Γ(X，Y)=∑XijYij，ΔI和ΔId分别为对原图像和降噪后图像进行拉普拉斯滤波的结果和为图像像素均值，该实验中使用3×3拉氏算子。

对15张视网膜OCT图像利用SARD与DW_ SRAD两种算法进行实验后，得到的以上4种指标的平均结果如表1所示。

表1 DW_SRAD与SRAD降噪性能对比Table 1 Performance indexes comparison between the DW_SRAD and SRAD method

从表1可看出，新方法在没有明显加重对原图像的模糊(用于衡量降噪后图像锐度或者与原图像相比的失真程度的指标β值保持基本不变)的前提下，其它3种衡量降噪效果的指标都有明显的提高。

3 结论

利用相邻尺度间的二进小波系数的相关性，提出了一种新的扩散系数计算公式，对斑点抑制各向异性扩散滤波方法SRAD进行了改进。改进后的算法在保持经典SRAD算法的增强边缘能力的同时，显著提高了其降低噪声的能力，信噪比、等效视觉等客观量化指标都有所提高。利用小波系数的尺度相关性可以改善扩散滤波效果。该方法对一般斑点图像噪声降噪有普适性。

［1］Zysk A M，Nguyen F T，Oldenburg A L，et al.Optical coherence tomography:a review of clinical development from bench to bedside［J］.Journal of Biomedical Optics，2007，12(5):12－16.

［2］Jindahra P，Hedges T R，Mendoza－Santiesteban C E，et al.Optical coherence tomography of the retina:applications in neurology［J］.Current Opinion In Neurology，2010，23(1):16－23.

［3］Karamata B，Hassler K，Laubscher M，etal.Speckle statistics in optical coherence tomography［J］.Journal of the Optical Society of America A－Optics Image Science and Vision，2005，22(4):593－596.

［4］Sendur L，Selesnick IW.Bivariate shrinkage functions for wavelet－based denoising exploiting interscale dependency［J］.IEEE Transactions on Signal Processing，2002，50(11):2744－2756.

［5］Pizurica A，Philips W.Estimating the probability of the presence of a signal of interest inmultiresolution singleand multiband image denoising［J］.IEEE Transactions on Image Processing，2006，15(3):654－665.

［6］Jian Zhong－ping，Yu Zhao－xia，Yu Ling－feng，et al. Speckle attenuation in optical coherence tomography by curvelet shrinkage［J］.Optics Letters，2009，34(10): 1516－1518.

［7］Yong Yue，Croitoru M M，Bidani A，et al.Nonlinear multiscale wavelet diffusion for speckle suppression and edge enhancement in ultrasound images［J］.IEEE Transactions on Medical Imaging，2006，25(3):297－311.

［8］Yu Y J，Acton S T.Speckle reducing anisotropic diffu－sion［J］.IEEE Transactions on Image Processing，2002，11(11):1260－1270.

［9］Perona P，Malik J.Scale－Space and edge－detection using anisotropic diffusion［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1990，12(7):629－639.

［10］Mallat S，Zhong S.Singularity detection and processing with wavelets［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1992，14(7):710－732.

［11］Junmei Zhong，Ruola Ning，David Conover.Image denoising based on multiscale singularity detection for cone beam CT breast imaging［J］.IEEE Transactions on Medical Imaging，2004，23(6):696－703.

［12］Zhong Jun－mei，Sun Hui－fang.Wavelet－based multiscale anisotropic diffusion with Adaptive Statistical analysis for image restoration［J］.IEEE Transactions on Circuits and Systems—I:Regular Papers，2008，55(9):2716－2725.

［13］段炼，何永红，朱锐，等.三维谱域光学相干层析成像系统的研制［J］.中国激光，2009，36(10):2528－2533.

Duan Lian，He Yong－hong，Zhu Rui，et al.Development of a spectrum domain 3D optical coherence tomography system［J］.Chinese Journal of Lasers，2009，36(10): 2528－2533.

［14］Salinas H M，Fernandez D C.Comparison of PDE－based nonlinear diffusion approaches for image enhancement and denoising in optical coherence tomography［J］. IEEE Transactions on Medical Imaging，2007，26(6): 761－771.

［15］邓菊香，梁艳梅.光学相干层析图像的小波降噪方法研究［J］.光学学报，2009，29(8):2138－2141.

Deng Ju－xiang，Liang Yan－mei.Noise reduction with wavelet transform in optical coherence tomographic images［J］.Acta Optica Sinica，2009，29(8):2138－2141.