基于Hermite微分算子的水平集MR图像分割方法

1、基于Hermite微分算子的水平集MR图像分割方法 作者：杨丰；金大年；陈武凡；罗敏(南方医科大学生物医学工程学院，广东  广州  510515)摘要：本文介绍一种基于Hermite微分算子用于MR图像分割的水平集方法。该方法采用Hermite微分算子来代替传统一阶差分算子，即在水平集曲线演化时函数微分用像素点的二阶邻域差值求得，而不是传统方法由一阶邻域决定。实验结果表明，对于相同的分割过程，运用了Hermite微分算子的水平集方法，其分割结果更加精确。尤其是对于由噪声等因素所引起的退化图像，其分割效果明显优于传统方法，而运算速度与传统方法相差无几。关键

2、词：MRI；图像分割；水平集；有限差分；Hermite中图分类号：R318.5  文献标识码：A  文章编号：1673-4254(2006)01-0036-05A level set method based on Hermite derivative filter for segmentation of magnetic resonance imagesYANG Feng；JIN Da-nian；CHEN Wu-fan；LUO MinCollege of Biomedical Engineering, Southern Medical Universi

3、ty, Guangzhou 510515, ChinaAbstract: A level set segmentation algorithm based on Hermite derivative filter is proposed for segmentation of human magnetic resonance images (MRI). Instead of utilizing the traditional first difference, Hermite derivative filter was used to calculate the differential co

4、efficients in the course of level set interface evolution, so that the differential coefficients were no longer decided by the first neighbor, but by the second neighbor of the examined pixel. Results of the segmentation tests proved that to the same segmentation process, the level set method utiliz

5、ing Hermite derivative filter produced a more accurate result. The proposed method showed especial superiority over the conventional method for images with interferences by noise. At the same time, the new algorithm did not increase the time for the segmentation.Key words: magnetic resonance imaging

6、; image segmentation; level set; finite difference; Hermite收稿日期：2005-06-12基金项目：“973”国家重点基础研究发展规划项目(2003CB716100)；广东省自然科学基金项目(32891)Supported by Natural Science Foundation of Guangdong Province(32891)and National "973" Development Project in High-tech Research(2003CB716100)作者简介：杨丰(1965-)，男，

7、副教授，研究方向：医学图像处理和模式识别，E-mail:yangf     MR图像分割是计算机对人体组织进行定性或定量分析所必不可少关键步骤。基于几何变形模型的水平集图像分割方法123，将二维(三维)的闭合曲线(曲面)的演化问题转化为高维空间中水平集函数曲面演化的隐含方式来求解，适应于对拓扑结构变化的处理，其计算精度高，算法稳定，因而广泛应用于二维、三维医学图像分割。    传统水平集方法在进行数值求解过程中常用一阶差分和中点差分来近似代替函数的导数值，易受由噪声等因素引起的图像退化影响，造成MR图像分割结果不理

8、想，主要表现为边界泄漏现象，虽然通过在停止力项加入基于区域信息的策略可以在一定程度上减小噪声对分割过程的影响，但仍不能完全克服边界泄漏现象。针对此问题，本文介绍一种基于Hermite微分算子的水平集图像分割方法，该方法采用Hermite微分算子来代替传统一阶差分算子，即在水平集曲线演化时函数微分用像素点的二阶邻域差值求得，而不是传统方法由一阶邻域决定，使得分割过程更加精确，能较好地抵抗噪声对分割过程所产生的干扰。    本文第一部分简要介绍水平集方法基本思想，其中包括函数的数值求解过程；第二部分介绍改进的基于Hermite微分算子的水平集方法；第三部分是

9、对实验结果的分析和讨论。    1  水平集方法简介    1.1  基本思想    以二维情况为例，水平集方法(level set method, LSM)的主要思想23是把二维平面的闭合曲线(t)视为三维空间连续函数曲面(水平集函数)的一个=0的零水平层。这样，遵循一定的关系，二维曲线的演化就转换为三维水平集函数曲面的演化。通过截取不同时刻三维空间连续函数曲面的水平层，我们就可以得到当前时刻二维平面的闭合曲线的位置和形状。 

10、;   界面(t)在力场(又称速度函数)F作用下，沿切线方向运动的欧氏方程为：    (x,y,t=0)=±d  (1)    其中：d为从象素(x,y)到界面(t)的有符号距离，且规定象素点(x,y)在初始界面外部为正，内部为负。因此，关于初始界面函数为：    (t=0)=(x, y)|(x, y, t=0)=0  (2)    在给定(x, y, t=

11、0)值，关于函数(x, y, t)的水平集运动方程为：      (3)    上式(3)是一个标准的Hamilton-Jacobi方程3。方程中的Vext是与轮廓线法线方向的速度相关的某种函数，由于通常速度只在轮廓线的位置有意义，也即定义在零水平集，因此常需要将零水平集的速度扩展到整个函数去曲面定义域，称为扩展速度场Vext。    1.2  近似数值解法    由于水平集函数在演化过程中始终保

12、持为一个函数，因此，我们可以利用离散网格来表达水平集函数(x, y, t)。设离散网格的间隔分别为xy，且在n时刻，结点(i, j)处的水平集函数为，则(ix,jy,nt)。这里t是时间步长，则方程(3)可以离散化为：      (4)    上式中 (x,y,t)的有限数值差分 必须采取适当的形式，以避免常量演化中导致的水平集函数奇异性问题。正确处理此问题必须基于波前传播的“熵守恒”理论1，应用在偏微分方程的求解上，就是偏微分方程的迟滞解方法。于是方程(4)的近似数值表达式为： &#

13、160;      其中，和分别表示像素点(i, j)在当前与下一步的水平集函数；而          (6)    这种差分方法称为逆向有限差分法(Upwind Finite Difference)3。    分别表示不同方向的差分值，为中点差分，它们的表达式为：        

14、;  (7)    而kij是水平集函数在(i,y)的曲率，由下式定义：      (8)    可以由中点差分计算出来。    2  基于Hermite微分算子的水平集方法    2.1  Hermite微分算子    如上所述，在求解水平集偏微分方程时，之前的几乎所有的水平集

15、方法基本都是采用式(7)所示的一阶差分或中点差分微分算子来近似函数的导数。采用一阶差分微分算子来求水平集偏微分方程近似解的好处是计算简便，耗时较少，但缺点是误差较大，尤其当图像质量退化比较严重的时候，往往得不到理想的分割结果(见实验结果图1C)，有时候甚至会导致分割失败。    我们引入了一种基于Hermite多项式逼近的微分算子4。    由于幂函数系xk, k=0, 1,在任何区间上线性无关，所以，可采用Gram-Schmidt(克莱姆施密特)正交化方法由幂函数系产生在指定区间a, b上带指定权函数(x)4的正

16、交多项式系，正交化方法如下：        (9)    其中，    若直接采用1, x, x2, ,xn作基底，当权函数(x)=1时，虽然()=(xk,xj)容易计算，但由此形成的法方程系数矩阵在n稍大时是病态矩阵，用单字长在计算机上求解方程，其结果往往不太可靠。为避免上述弊端，采用正交基底。    Hermite多项式插值是带导数条件的代数插值，即，它满足如下条件：  &

17、#160;     (10)    用Gram-Schmidt正交化方法计算Hermite正交多项式Hn(x)，它满足递推关系：        (11)    以Hn(x)系为基底，将原离散图像函数表示为p*(x)=。    用最常用也最容易实现的最佳平方逼近5为标准来度量逼近误差。    求最佳平方逼近元素

18、p*(x)=，就是求它所含的系数ck(k=0, 1, )。由6可知，      (12)    离散形式为：      (13)    其中xi是区间(-, )上带权的n次Hermite多项式：    的零点即Hermite多项式的插值结点，见表1。     我们取n=5，因为其插值结点为±2.0201828

19、705，±0.9585724646，0，可以通过计算将k表示为相邻五个象素的函数。    计算得出Hermite微分算子：    f'=0.100232f(x-2)-f(x+2)-0.432275f(x-1)-f(x+1)  (14)    上式就是二邻域的Hermite微分算子的表达式。    我们再来考察Hermite微分算子的截断误差，用泰勒级数公式展开后可以求得：  &#

20、160;   (15)    当我们取n=5时，可以看到Hermite微分算子的截断误差对于函数是2×n=10阶的，而一阶差分的截断误差是一阶的，中点差分的截断误差是二阶的。因此，运用Hermite微分算子的来逼近函数的导数比一阶差分微分算子更为准确。    2.2  基于Hermite微分算子的水平集差值求解    在基于Hermite微分算子的的算法中，我们把式(7)中(i, j)点沿x轴方向的向后一阶差分微分改

21、为(i-1, j)点的Hermite微分，向前一阶差分微分改为(i+1, j)点的Hermite微分；同理，沿y轴方向的向后、向前一阶差分微分也分别由后一像素点和前一像素点的Hermite微分代替。这样，我们把式(7)改为：        而在计算式(8)的曲率kij时，也把原来的中点差分改为Hermite微分。这样，对于水平集数值求解过程中的导数，我们把原来的一阶差分微分算子和中点差分微分算子都用Hermite微分算子来代替。原方法每个像素点的导数仅由该点一阶邻域决定，误差比较大，而且当临近像素点出现噪声时，往往

22、得到错误的结果，从而可能导致整个分割过程的失败。而新的方法该像素点的导数通过该点二阶邻域求得，在计算上更加精确，即使周围出现噪声点，其造成的影响也相对原方法要减小了很多。    3  实验结果与分析    图1是传统LSM方法、BLSM方法和HBLSM方法分别对一加了高斯噪声的脑部MR图像进行白质分割的结果。其中BLSM方法7是基于贝叶斯分类的水平集方法(Level Set Method Based on Bayesian classification)，该方法在计算水平集函数的停止力时加入了基

23、于贝叶斯分类的图像区域信息；HBLSM是将BLSM的数值求解过程用上述基于Hermite微分算子的算法代替。    其中图1 A为原图像；图1 B为传统LSM方法分割结果，可以看到泄漏现象非常严重，整个分割过程完全失败，主要原因是传统水平集方法仅采用图像局部的梯度信息作为轮廓线停止运动的停止力，当边界比较模糊而且噪声比较严重的情况下，轮廓线无法正确识别边界而继续向边界外运动，从而造成广泛的泄漏；图1 C为BLSM方法分割结果，可以看到此方法由于在识别边界时结合了贝叶斯分类的区域信息，因此较好地解决了边界泄漏问题，结果明显好于传统LSM方法，但仍在个别边界

24、发生了泄漏，这主要是因为某些噪声点正好位于本来就比较模糊的边界附近，使得轮廓线不能正确地停止在边界；图1 D为HBLSM方法分割结果，因为该方法使得各个像素点的水平集函数的导数通过该点二阶邻域求得，而不是传统方法由一阶邻域决定，从而极大减弱了噪声对水平集函数数值求解过程的影响，因此可以看到泄漏现象已经不存在了，分割结果比较理想。    我们采用Harvard大学附属医院Massachusetts General脑部形态分析中心网络图像标准分割库IBSR提供的20组正常人的脑部MR图像8作为实验数据，利用HBLSM方法和BLSM方法对这20组 “金标准”数

25、据进行分割比较。表2是分别采用BLSM方法和HBLSM方法对这20组数据关于白质的重叠率和平均重叠率统计。从表中可以看到，对于相同的分割过程，运用了Hermite微分算子的水平集方法，具有更优的分割结果。反映在BLSM方法上，平均重叠率为76.44%，较BLSM法的74.47%提高了1.97%。虽然对于大多数的数据组，HBLSM方法重叠率与BLSM方法重叠率比较接近，提高幅度不大，有个别组甚至是一样的结果，这是因为这20组标准脑部MR图像数据中的大部分图像质量比较理想，因此在水平集函数数值求解过程中，对导数采用一阶差分微分算子和Hermite微分算子差别不大。但对于数据编号为15_3、16_3

26、和5_8的这几组图像，提高幅度比较大，原因是这几组图像不同程度的具有图像质量退化，尤其以5_8这组添加了噪声的图像数据提高幅度最大，达到8.7个百分点。上述结果显示了基于Hermite微分算子的水平集数值求解方法对于由噪声等原因引起的图像退化，分割结果具有理想的改进。    在分割结果有所改进的同时，我们再来考察基于Hermite微分算子的水平集方法的运算速度。BLSM方法和HBLSM方法对上述20组图像数据进行分割的平均用时分别为15.4和15.5 s，算法都以VC+6.0为平台，在AMD1.8 GHz处理器、512 MB内存的微机上实现。可以看到两种

27、方法用时相差无几。因为Hermite微分算子计算每个导数需要2次乘法，3次减法，比一阶差分的1次减法相比，多了2次乘法和2次减法；而与中点差分的1次减法和1次除法相比，多了2次乘法和2次减法，但是也少了1次除法。我们都知道，微机对加法、减法和乘法的计算速度较快，只是对除法比较慢。而Hermite微分算子并没有除法，所以运算速度与传统方法相差不大。     图1  3种方法对脑部噪声MR图像的分割结果    Fig.1  Segmentation results of brai

28、n MRI with noise for the three different algorithms    4  结论    本文介绍了一种基于Hermite微分算子的水平集方法。该方法在水平集函数数值求解过程中，利用了Hermite微分算子来近似代替函数的导数，较好地克服了传统水平集在计算过程中仅利用一阶差分和中点差分微分算子所带来的误差，使水平集数值求解更加精确和稳定，尤其对包含噪声的图像分割更加具有其独特的优势。    参考文献：    1Sethian JA. An analysis of flame propagationD. Berkeley: Univ California, 1982.    2Sethian JA. Level set methods and fast marching methods: evolving inte