先验形状约束的超声前列腺图像分割方法

刘维平，杨新

0 引言

前列腺癌是我国男性中发病率呈增长趋势的恶性肿瘤，早期诊断对于及时治疗有重要的意义。医院确诊前列腺癌的主要手段是经直肠的超声活组织穿刺取样检验，而结合多种成像模式的穿刺活检能够提高前列腺癌的检出率。核磁共振成像对前列腺及其内部结构的成像效果远好于超声图像，便于医生观察以确定癌症位置。若将核磁共振成像中疑似癌症的位置映射到超声图像中作有针对性的穿刺，就需要对超声和核磁共振成像的前列腺分别分割，然后通过配准分割所得的轮廓使得超声和核磁共振成像前列腺的切片对应，方便与医生观察和比较。因此，超声前列腺分割对前列腺癌症的诊断有着重要的意义。除此之外，其在前列腺特异抗原计算、超声聚焦、超声引导的前列腺短程放疗等应用领域，同样具有极大的潜力。

近年来，基于先验知识的分割方法在超声前列腺分割中受到了广泛的关注。该方法的特点在于在传统的分割方法中引入了先验知识，对于待分割目标在复杂背景下产生的分割结果的不确定性，可依靠目标的先验知识得到解决，即利用先验知识对目标进行恢复和重构的过程。

Shen[1]提出处理局部形变的形状模型，采用旋转不变特征作为超声图像属性驱动前列腺形状模型的形变。为了避免一些优化过程中出现的局部能量最小，Shen采用多分辨率的策略，在低分辨率中粗略分割，在高分辨率中精确分割。

Xu[2]提出利用超声前列腺的先验直方图和先验形状模型实现分割。先验直方图和形状都是通过机器学习的方法从训练集中学习出来，通过一个手动调节的参数平衡这两项对曲线演化的影响。

Yan[3]采用稳定、正确的特征提取的结果来引导形状模型的形变。在特征提取的时候对特征进行初步的筛选，那些不好的或者不稳定的特征提取结果是不对形状模型的形变产生任何影响的。

本文采用基于先验形状的活动轮廓模型对超声前列腺图像分割。将前列腺的形状通过机器学习的方法求得形状模型，并将其作为先验知识引入活动轮廓模型，使得曲线满足先验形状的约束。

1 局部高斯拟合的活动轮廓

活动轮廓模型是一种图像分割方法，其主要思想是在图像域定义一个曲线或曲面，并且在由曲线或曲面自身相关的内力，以及由图像数据定义的外力的作用下移动。通常情况下，外力的定义是基于边缘特征或者区域特征的，它们分别采用了待分割图像的局部和全局信息。为了兼顾全局和局部特征对演化曲线的影响，近些年来，一些学者提出对背景和目标做局部统计的方法。最为代表性的思路是基于局部拟合方法的能量函数，用以解决灰度分布不均匀的图像分割问题。

Wang[4]提出基于局部高斯拟合的活动轮廓模型，对于图像上的某点，定义该点处的能量函数为

其中，Kσ(x)为截断的高斯核函数，σ为方差，pi,x(I(y ))，i= 1,2，为像素的灰度在区域Ω1，Ω2中的概率密度分布。点x处能量函数被限制在以该点为中心的某个窗口中，截断高斯核的方差决定了窗口ρ的大小，并且，图像中的点y对能量函数的贡献随着y与x距离的增大而减小，该截断高斯核函数如下所示：

其中，c是归一化的常量，通常情况下，令ρ=2σ。

在Wang提出的方法中，假设概率密度分布 pi,x(I(y ))满足高斯分布，即

其中， ui( x)和 σi(x)2分别为点x附近第i个局部区域灰度的均值和方差。

整个目标的边界通过最优化图像域Ω中所有局部高斯拟合能量函数得到，即

其中，φ为水平集函数，H(·)为 Heaviside函数，公式(5)是公式(4)的水平集表示。对该能量函数求极值并结合梯度下降方法可得到演化方程和某点处演化曲线内部和外部灰度值的局部拟合方程。该方法能够有效地克服图像亮度不一致的问题，并且在纹理图像分割中有着较好的效果。

2 先验形状模型约束的活动轮廓模型

影响超声前列腺图像分割的主要原因有如下几个方面：1）斑点噪声是超声图像中常见的干扰，导致分割困难。2）在超声前列腺区域中，灰度或者纹理分布是不均匀的，有些部位灰度级较低，而另外一些部位灰度值较高；在超声前列腺外部，有些地方的灰度分布和前列腺区域内部的分布一样，故难以采用全局的纹理或者灰度特征区别前列腺的内外部分。3）超声伪影也会对前列腺分割产生严重的影响。因此，超声图像的分割迄今仍然是一个难题。

在超声前列腺图像分割中，单纯采用基于边缘或者区域的传统活动轮廓模型的方法很难得到较好的分割结果，因此引入前列腺的先验性形状作为传统活动轮廓模型的约束。对于先验形状约束的图像分割方法，首先要将待分割目标形状的训练集配准，然后将配准后的形状采用机器学习的方法提取形状特征并且对形状建模，最后将形状模型嵌入到传统活动轮廓模型当中。

主成分分析方法(Principle Component Analysis, PCA)常常用于提取形状模型的形状特征，而形状形变能力主要由形状模型形变参数的概率密度估计表示。当目标的形状变化不是很大的情况下，可以假设形变参数的概率密度函数满足高斯分布，而当形变很大的情况下，就需要采用非参数的方法估计形变参数[5]。但是超声前列腺形状的形变不是严重，因此针对这个具体问题，采用高斯分布的形状模型。

2.1 形状配准

形状配准是机器视觉中的一个非常复杂的问题，其基本思想可以描述如下：给定目标形状、当前要配准的形状和相似性度量，求解某变换使得变换后两个形状的相似性度量最小。Paragios[6]采用符号距离函数描述形状，通过最优化两形状的平方误差和从而实现形状的刚性配准。本章采用这种方法将形状训练集中的图像配准。

假设目标形状为ψT(X , Y )，当前形状为ψC(U , V )，目的是找到两者间的刚性变换

其中，F(·)为前向映射；s，θ，Tx，Ty分别为从坐标系(U, V)前向变换到坐标系(X, Y)的尺度参数，旋转参数，x轴方向的平移参数，y轴方向的平移参数，使得两个形状的平方误差和最小，即最优化如下能量函数

其中，Φ为坐标系(X, Y)空间；B(·)为后向映射，即为 F(·)的逆映射

采用梯度下降的方法对公式(12)求最优，得到配准参数的迭代方程如下

通过迭代公式(13)～(16)直至稳定，可将两个以符号距离函数表示的形状配准。在训练形状模型的时候，给定形状的训练集，选取其中任意一个形状的作为目标形状，将其余的形状与该形状采用上述方法配准即可得到配准后的训练集。

2.2 训练超声前列腺的形状模型

其中，U的列向量描述了形状模型变化的正交系，Λ是相应奇异值的对角矩阵。新形状可以用k个主分量和k维列向量系数γ计算如下：

其中，Uk是d×k维矩阵，包含了k个最大特征值对应的U的列向量。

一个确定k的方法是选择满足公式

最小的k。假设参数向量γj满足高斯分布，并且满足-3≤γj≤3，γj是 γ的第j个分量，λj是(1 n) MMT的第 j个特征值。通过调整γ使形状模型产生形变。

2.3 基于先验形状模型的前列腺分割

定义基于先验形状模型的活动轮廓的能量函数如下：

其中，ETotal为总的能量函数，ELGDF为局部高斯拟合的活动轮廓模型的能量函数，EShape为形状模型的能量函数。分割的目的是通过最优化以上的总的能量函数来驱动水平集曲线演化。λ为权值系数，用以平衡两个能量项对曲线演化的影响。

在基于先验形状的活动轮廓模型中，需要最优化形状模型与演化曲线的平方误差和，使得形状模型与演化曲线最匹配。因此，整个系统中有两套参数需要被优化：配准参数和形变参数。构造的能量函数如下所示：

其中，B(·)为后向映射，包括尺度s，旋转θ，x和y方向上的平移参数分别为Tx和Ty。配准参数的迭代方程如下

其中F(·)为前向映射。通过梯度下降方法得到形变参数的迭代方程：

其中，jγ为形状模型的第j个形变参数，j=1,…,k。

对能量函数求最优解并且利用梯度先将方法可以得到如下公式

其中，δ(φ)为Dirac函数。公式 (23) 为曲线演化方程，公式 (24)～ (27) 为计算局部高斯拟合方法的参数计算公式。

其中，α=0.1为设定的参数。

给出先验形状约束的活动轮廓模型的具体计算步骤：

1.给定演化方程的终止条件

2.先验形状约束的前列腺分割算法步骤：

3.初始化活动轮廓和形状模型；

4.将形状模型与演化曲线配准，公式(17)～(20)；

5.优化形变参数，公式(22)；

6.迭代曲线演化方程，公式(23)～(27)；

7.判断曲线是否满足终止条件(28)，如果不满足则转到第4步；如果满足则结束算法。

3 实验结果

超声数据库有12个病例，每个病例有6个2维超声前列腺图像的切片，每个切片的大小为 700×400。这些超声图像数据来自于新华医院超声科，并且由医生手动分割完毕。采用留一法，即每个病历留下一个作为分割结果的测试集，其余的图像的手动分割结果作为形状模型的训练集。因此，在这个模型当中，采用60幅前列腺的形状作为形状模型的训练集，12幅图像作为算法的测试。

将配准后的前列腺形状作为训练集来训练前列腺的形状模型。通过分别调整公式(13) 中的jγ，j=1,…,5，得到前列腺形状模型的形变结果，如图1所示：

图1 前列腺形状模型的形变

图1可以看出，这种形状模型能够较好地模拟前列腺形变。当然，有学者提出更为复杂的形状模型，用以提取形变较大的形状，例如基于核密度估计的形状模型。但是在二维前列腺分割中，基于PCA和高斯分布的形状模型足以应付前列腺的形变，因此，针对特定问题采用了该种形状模型。

在超声前列腺的分割中，通过人为在前列腺图像中选取上下左右四点拟合成一个闭合区域作为活动轮廓模型演化曲线的初始化，如图2(a)所示。先将λ赋予较小的值，再将该值在每次迭代曲线演化方程的时候乘以一个大于 1的参数。这样使得形状约束逐渐作用于曲线演化，避免了能量函数的局部最小。图2(b)～(d)为曲线演化和形状模型配准、形变的过程，可以看到曲线最终收敛到形状模型附近，得到了较好的分割效果，如图2(a) ～(d)所示：

图2 前列腺分割过程

从另外一些超声前列腺的分割结果中可以看出，尽管超声图像中的噪声、伪影、灰度分布不均匀等问题很严重，提出的方法仍然可以较好地提取前列腺的边缘，如图3所示：

图3 其它分割结果

为了对分割算法的准确程度有一个客观的评测，医学图像中，普遍采用的是和临床医生的手工分割的标准轮廓进行比较。采用 Dice相似性系数 (Dice Similarity Coefficient,DSC) 作为衡量标准。Dice相似性系数定义如下

其中，As和Ag分别是算法分割结果和医生手动分割结果的面积。算法的分割结果与医生手动分割结果的Dice相似性系数，如表1所示：

表1 算法分割结果和医生手动分割结果的Dice相似性系数

4 结论

提出一种半监督的、基于先验形状约束和高斯拟合活动轮廓的分割模型，主要应用于超声前列腺图像分割。通过将分割结果与医生手动分割结果的对比，验证了所提出方法的有效性。

[1]Dinggang Shen, Yiqiang Zhan, Christos Davatzikos.Segmentation of prostate boundaries from ultrsound images using statistical shape model [J].IEEE Transactions on Medical Imaging, 2003, 22 (4): 539-551.

[2]Robert S.Xu, Oleg V.Michailovich, Igor Solovey, etc.A probability tracking approach to segmentation of ultrasound prostate images using weak shape priors [C].In Proceeding of SPIE Medical Imaging, 2010.

[3]Pingkun Yan, Sheng Xu, Baris Turbey, etc.Discrete deformable model guided by partial active shape model for TRUS image segmentation [J].IEEE Transactions on Biomedical Engineering, 2010, 57 (5): 1158-1166.

[4]Li Wang, Lei He, Arabinda Mishra, etc.Active contours driven by local Gaussian distribution fitting energy [J].Signal Processing, 2009, 89 (12): 2435-2447.

[5]Weiping Liu, Yanfeng Shang, Xin Yang, etc.A shape prior constraint for implicit active contours, Pattern Recognition Letters [J].2011, 32(15): 1937-1947

[6]Nikos Paragios, Mikael Rousson, Visvanathan Ramesh.Non-rigid registration using distance functions [J].Computer Vision and Image Understanding, 2003, 89(2-3): 142-165.