超声弹性图像去噪方法研究

邵党国，刘 帆，相 艳，马 磊，易三莉，贺建峰

（昆明理工大学 信息工程与自动化学院和云南省计算机应用重点实验室，云南 昆明 650500）

超声弹性图像去噪方法研究

邵党国，刘 帆，相 艳，马 磊，易三莉，贺建峰

（昆明理工大学 信息工程与自动化学院和云南省计算机应用重点实验室，云南 昆明 650500）

超声弹性成像能直观展现各组织的弹性信息，已经成为肿瘤检测的重要工具。由于成像系统中的固有特性，使得图像遍布着大量伪影噪声，降低了图像的可辨读性。为了获取伪影噪声较少且组织边缘结构清晰的弹性图像，本文提出一种分数阶偏微分图像去噪方法，并通过各向异性算法有效平衡了梯度阈值k与微分阶数v的关系。人体体模弹性图像和仿真弹性图像的测试结果表明，该算法能有效改善弹性图像质量，增强边缘结构信息，且去噪后图像的信噪比及噪声比例差得到显著提升。

超声弹性图像；分数阶微分；各向异性扩散；去噪

0 引言

Ophir等人提出的超声弹性成像技术在近些年得到迅速发展，已经成为现代医学成像的新模式[1]。该技术主要是利用生物组织与相关病变组织间的弹性差异进行成像，原理如下：首先给待检测组织施加一个内部或外部的动态/静态压力，受到同等压力的软组织块相比硬组织块产生的形变程度更大；分别采集组织块形变前后的回波信号，对其进行分析，计算出各组织的位移分布，最终得到待检测组织的应变分布及其弹性系数分布并进行成像[2]。

由于超声成像系统的固有特性，使得弹性成像和其他超声成像模式产生的图像一样，必然会受到噪声影响，进而降低图像的可辨读性。弹性图像的主要噪声来源有两个：解相关误差和幅度调制误差。解相关误差是由压缩、变形、平滑移动或其他运动而导致运动跟踪不符合常规的时间延迟估计模型；幅度调制误差是由信号的随机波动，使得位移估计的位置不在分析窗口内。各种图像去噪技术，如直方图匹配算法[3]、双边滤波算法[4]、小波去噪算法[5]以及复合方法[6]等已经被运用到超声弹性图像的去噪研究中，并取得了一定的研究成果。噪声抑制对图像分析、边缘检测、图像分割、图像重建等都具有重要意义[7-8]。

分数阶微积分算法是近几年图像去噪领域的重要工具，其原理是利用分数阶微积分在提高信号高频成分的同时能够非线性保留信号中低频成分的特殊性质，将其应用于数字图像领域，达到抑制噪声并保留边缘的目的[9-10]。尽管分数阶微分算法在数字图像中得到了较好的应用，但在超声图像中的研究仍较少，特别是弹性图像。为了获取边缘结构清晰且噪声区域平滑的超声弹性图像，本文提出了一种基于分数阶微分的各向异性扩散算法（FAD）。

1 分数阶微分

分数阶微分至今为止并没有统一的定义，从不同角度推导可得到不同的定义表达式，目前较常用的是 Grünwld-Letnikov定义，Riemann-Liouville定义和 Caputo 定义[11-13]。由于 Grünwld-Letnikov 定义是将连续可导的微分阶数由整数阶推导至分数阶，通过对整数阶微分的差分近似递推式求极限推衍而来，因而该定义更适合应用于图像处理领域[14]，其推导过程如下：

可微函数f(x)一阶导数定义：

其二阶导数定义：

推导至函数f(x)的n阶导数定义：

根据上述推导，引入 Gamma函数概念，可以得到Grünwld-Letnikov定义表达式如下：

式中，一元信号f(x)的持续期t为[a, t]，h=(t-a)/n为步长，v为微分阶数。若将一元信号f(x)的持续期按单位h=1进行等分，则n=[(t-a)/h]h=1=(t-a)，可以得到一元信号f(x)的差分表达式如下：

1.1 滤波模板的构建

针对各向同性扩散模型去噪时出现的边缘模糊问题[15]，Perona和Malik于1990年提出了一种非线性偏微分方程模型（简称 P-M 模型），有效的保护了图像边缘，该模型定义如下[16]：

式中，div为散度算子，∇为梯度算子，∇I表示图像I的梯度，为扩散系数方程，I0（x,y）为原始图像，t为表示降噪过程与扩散时间有关的时间算子。Perona等人给出了扩散系数方程，如下所示：

式中，k表示扩散阈值，可阻止方程进一步扩散，图像梯度∇I可视为边缘检测器。若∇I远远大于扩散门限k值，趋于0，扩散得到抑制；反之，若∇I远远小于扩散门限k值，）趋于1，扩散得到加强，此时可视为高斯滤波。

由一元信号的差分表达式（5）可知，表达式左侧前n项的系数总和为一个非零数值，这是分数阶微分与整数阶微分显著特性区别之一。图像相邻像素间存在着高度的自相关性，为了使滤波模板具有抗旋转性，本文构建了具有 x轴、y轴正负方向及两对角线这8个方向的滤波模板；同时对滤波模板进行了部分修改以得到更好的去噪效果，修改后的模板如图1所示：

图1 微分掩模算子Fig.1 The fractional order differentiation mask

其中v为微分阶数，用上述滤波模板对图像进行卷积运算，结果发现：当k值不变，v逐渐增加时，图像的保边能力逐渐增强；当v值不变，k逐渐增加时，图像的平滑能力逐渐增强。可见扩散阈值 k和微分阶数v是图像平滑与图像保边的重要参数。因此，如何平衡k、v之间的关系是本文重点解决的问题。

1.2 平衡k、v关系

用 FAD算法对体模弹性图像和数字仿真弹性图像进行处理，分别求取 k为[1,35],v为[0.1,5.5]时

各图像所对应的SNRe及CNRe[17-18]，结果如图2、图3所示。其中，横坐标为v值（范围为[0.1,5.5]，间隔0.01），纵坐标为k值（范围为[1,35]，间隔1）。其中每个（a）为SNRe的结果，每个（b）为CNRe的结果，图中颜色越亮意味着SNRe值或CNRe值越大；反之则越小。

图2 人体体模弹性图像去噪后的SNRe、CNRe值，横坐标为v（范围为0.1-5.5）纵坐标为k（范围为1-35）（a）、（b）分别为SNRe和CNRe值Fig.2 The de-noising values of SNRe and CNRe for the phantom ultrasonic elastography.The horizontal axis represents v, and k values on the vertical. (a) and (b) are the SNRe and CNRe values respectively

图3 数字仿真弹性图像去噪后的SNRe、CNRe值，横坐标为v（范围为0.1-5.5）纵坐标为k（范围为1-35）(a)、(b)分别为SNRe和CNRe值Fig.3 The de-noising values of SNRe and CNRe for the digital simulated ultrasonic elastography.The horizontal axis represents v, and k values on the vertical. (a) and (b) are the SNRe and CNRe values respectively.

由图2、图3可以看出k、v之间存在一定的关系，根据实验所得的SNRe、CNRe值，选取处理效果较好的k、v值进行曲线拟合，结果如图4所示，得到二者间的关系表达式如下：

式中，v为微分阶数，k为扩散门限值。表达式均方根误差为0.09721，误差平方和为0.567，拟合系数为0.9943。

2 实验结果

图4 弹性图像k、v拟合关系图Fig.4 The curve fitting relation between k and v for elastography

为了体现FAD算法的有效性，将FAD算法与传统PM算法进行结果对比。在PM算法中，除了滤波模板和微分阶数 v有所不同，其他参数均与FAD算法保持一致（即积分常数设置为 0.2，传导系数函数为1，迭代次数设置为100次）。图5-6中，每个（a）为原始超声弹性图像；每个（b）为 PM算法k=6的结果；每个（c）为FAD算法k=6,v=0.56的结果；每个（d）为PM算法k=9的结果； 每个（e）为FAD算法k=9,v=0.85的结果。

图5 人体体模超声弹性图像去噪效果对比Fig.5 The contrast of denoising in phantom ultrasonic elastography

图6 数字仿真超声弹性图像去噪效果对比Fig,6 The contrast of denoising in digital simulated ultrasonic elastography

由图5-6，可以明显看到两种算法均能够平滑同质区域的伪影噪声，但FAD算法相比PM算法在边缘信息保留效果上占有明显优势。SNRe及 CNRe值可以定量分析算法的平滑能力。为了减少分数阶微分估计误差，k、v值应满足关系式：v = 0.001509*k2+0.07144*k + 0.05085，k [6, 16]任意整数。体模弹性图像中取一个 25*25的区域测试 SNRe及 CNRe值；数字仿真弹性图像中取一个60*60的区域测试SNRe及CNRe值。

表1和表2为图5、图6各图像的SNRe、CNRe及运算时间。通过实验分析可知：当k值小于10时，如图5、图6的（b）（d），PM算法出现过度保边现象，进而降低了SNRe、CNRe；当k值大于10时，P-M算法由于平滑过度而丢失了部分边缘信息，显然这并不符合图像去噪的目的。FAD算法平衡了抑制噪声与边缘保留的关系，可得到较理想的去噪结果。且当k[6, 16]时，FAD算法在抑制噪声及边缘信息保留两个方面具有明显效果。若 k值小于 6，FAD算法的平滑噪声区域能力下降，出现过度保边的现象；若k值大于16，FAD算法的边缘保边能力下降，出现过度平滑，图像部分边缘细节信息被丢失。综上所述，可知梯度阈值k和微分阶数v之间的平衡是图像平滑与保边的两个重要参数。大量实验数据表明，基于FAD算法所提出的k、v关系式，去噪后弹性图像噪声抑制效果明显，且边缘信息也得到了较好的保留。

表1 SNRe、CNRe以及运算时间Table 1 SNRe、CNRe and the operation times

表2 SNRe、CNRe以及运算时间Table 2 SNRe、CNRe and the operation times

3 结论

超声弹性图像弥补了其他医学成像技术无可比拟的优点，即能提供组织的硬度信息，成为临床检测的重要工具。然而超声成像系统固有的噪声极大地降低了图像的质量。为了平滑伪影噪声并保留特征区域边缘信息，进而改善图像质量，本文提出一种基于分数阶微分的图像去噪方法，并通过各向异性扩散算法平衡扩散阈值k与微分阶数v之间的关系。人体体模弹性图像和数字仿真弹性图像的实验结果显示，本文采用的FAD算法能够有效平滑弹性图像伪影噪声，提高SNRe及CNRe值；与此同时，特征区域边缘结构信息也得到了较好的保留。相比于传统PM算法，FAD算法占有明显优势。

[1] 陈蒙. 医学X光成像中图像去噪算法研究[J]. 软件,2013(11): 62-63.

[2] 王春光, 刘文, 刘东权. 一种基于编码激励的超声弹性成像噪声抑制方法[J]. 计算机应用研究, 2013, 30(5): 1596-1600.

[3] 何颖妮, 邵党国, 刘东权. 直方图匹配算法在超声弹性成像上的应用研究[J]. 计算机应用研究, 2013, 30(4): 1266-1270.

[4] Zhou X, Liu W, Zhang C, et al. Speckle Reduction of Ultrasound Elastography with Bilateral Filter[J]. Advanced Science Letters, 2013, 11(4): 3025-3027.

[5] 顾伟, 杜景琦. 基于小波阈值的图像去噪[J]. 软件, 2017,38(7): 107-110.

[6] Cui S and Liu DC. Noise reduction for ultrasonic elastography using transmit-side frequency compounding: a preliminary study[J]. IEEE Transactions on Ultrasonics, Ferroelectrics,and Frequency Control, 2011, 58(3): 509-516.

[7] 钟彩. 边缘检测算法在图像预处理中的应用[J]. 软件,2013, 34(1): 158-159.

[8] 徐红娟, 韩玉合, 刘文洲, 等. 一种基于Ray Casting算法的医学影像重建系统[J]. 软件, 2012, 33(9): 118-118.

[9] 黄果, 许黎, 蒲亦非. 分数阶微积分在图像处理中的研究综述[J]. 计算机应用研究, 2012, 29(2): 414-420.

[10] 王福兴, 蒲亦非, 周激流. 分数阶微积分的应用研究[J].无线互联科技, 2011, (8): 12-14.

[11] Scherer R, Kalla S L, Tang Y, et al. The Grünwald-Letnikov method for fractional differential equations[J]. Computers &Mathematics with Applications, 2011, 62(3): 902-917.

[12] 陈庆利, 蒲亦非, 黄果, 等. 数字图像的0～1阶Riemann-Liouville分数阶微分增强模板[J]. 电子科技大学学报,2011, 40(5): 772-776.

[13] 陈庆利, 黄果, 张秀琼, 等. 数字图像的Caputo分数阶微分增强[J]. 计算机辅助设计与图形学学报, 2013, 25(4):519-525.

[14] 杨柱中, 周激流, 晏祥玉, 等. 基于分数阶微分的图像增强[J]. 计算机辅助设计与图形学学报, 2008, 20(3): 343-348.

[15] 苏永利, 张博, 张书玲. 基于各向同性逆扩散方程的指纹图像增强方法[J]. 西北大学学报自然科学版, 2009, 39(1):9-12.

[16] 陈婷婷, 沈民奋, 杨金耀. 一种复合各向异性扩散的图像去噪算法[J]. 数据采集与处理, 2014, 29(5): 757-763.

[17] ]Paul Liu and Dong Liu. Filter-Based Compounded Delay Estimation With Application to Strain Imaging[J]. IEEE Transactions on Ultrasonics Ferroelectrics & Frequence Control,2011, 58(10): 2078-2095.

[18] Shao Dangguo, Liu Paul, Liu DongC. Characteristic matchingbased adaptive fast bilateral filter for ultrasound speckle reduction[J]. Pattern Recognition Letters, 2013, 34(5): 463-9.

Reseach of Elastography Denoising Methods

SHAO Dang-guo, LIU Fan, XIANG Yan, MA Lei, YI San-li, HE Jian-feng
(Information Engineering and Automation, Kunming University of Science and Technology and Key Lab of Computer Technologies Application of Yunnan Province, 650500, Kunming, Yunnan)

Ultrasonic elastography which can directly show the elasticity of the organization's information has become an important tool for tumor detection. The inherent characteristics of the imaging system made the elastography distributed with a large number of artifacts. And then reduced the readability of the elastography. In order to obtain a high quality elastography with few artifacts and a clear edge structure, a fractional partial differential image denoising model was constructed, which can effectively balance k (the gradient threshold) and v (the differential order) by anisotropic algorithm. The results show that the algorithm can effectively improve the elastography quality and enhance the edge structure information. The signal to noise ratio and the noise ratio difference of the elastography after denoising are significantly improved.

Ultrasonic elastography; Fractional order differentiation; Anisotropic diffusion; Artifact denoising

TP391

A

10.3969/j.issn.1003-6970.2017.12.015

本文著录格式：邵党国，刘帆，相艳. 超声弹性图像去噪方法研究[J]. 软件，2017，38（12）：81-86

中国博士后科学基金：(2016M592894XB)；云南省科技厅面上项目(KKS0201703015)

邵党国(1979-)，男，博士后，硕士生导师，主要研究方向为医学图像处理、机器学习等；刘帆(1995-)，女，硕士研究生，主要研究方向为医学信号图像处理；相艳(1979-)，女，讲师，博士研究生，主要研究方向为图像处理、自然语言处理、数据挖掘、机器学习；马磊，男，博士研究生，主要研究方向为数据挖掘与分析，机器学习；易三莉(1977-)，女，讲师，硕士生导师，主要研究方向为医学图像处理；贺建峰(1965-)，男，教授，硕士生导师，主要研究方向为医学图像处理、数据分析等。