基于监督学习的前列腺MR/TRUS图像分割和配准方法

丛 明，吴 童，刘 冬✉，杨德勇，杜 宇

1) 大连理工大学机械工程学院，大连 116024 2) 大连医科大学附属第一医院泌尿外科，大连 116024 3) 大连大华中天科技有限公司，大连 116024

目前，前列腺癌的诊断主要依靠前列腺特异性抗原（Prostate specific antigen, PSA），并继以穿刺活检，由于PSA特异性较差，故临床上将穿刺活检作为确诊前列腺癌的金标准[1]. 经直肠超声引导下的前列腺穿刺活检因实时性、低成本、易操作等优点，成为了临床上最为普遍的诊断前列腺癌的方法. 然而，由于直肠超声图像低的成像质量，因此很难从图像上准确定位恶性肿瘤区域. 经直肠超声引导下前列腺6针法穿刺活检的假阴性率高达30%[2]. 另一方面，多参数核磁共振图像（Multiparametric magnetic resonance imaging, mpMRI）是目前公认的诊断前列腺癌的最佳影像技术，可精确定位可疑病灶区，从而达到靶向穿刺的目的[3]. 为了提高前列腺癌的检出率，核磁/超声（MR/TRUS）融合引导靶向穿刺技术应运而生[4-5]，它旨在对术前MR图像和实时TRUS图像进行配准融合从而提高穿刺精度. 然而，由于MR扫描过程中线圈或直肠内气体的膨胀、病人的呼吸和不自觉的运动、超声探头在直肠中对前列腺腺体的挤压，使得术前MR与实时TRUS下前列腺的形态并不一致.为了补偿形状变化的影响，通常需要对MR/TRUS图像进行变形配准.

由于超声图像低的信噪比及MR/TRUS图像上复杂的灰度变化关系，仅依据图像灰度特征很难对MR/TRUS图像上对应的结构特征进行定位配准，因此前列腺MR/TRUS的配准融合通常是基于已分割的前列腺表面进行变形配准. 现有的配准方法主要有基于灰度的配准方法[6-7]和基于分割表面[8-9]的配准方法. 前者主要针对MR/TRUS图像复杂的灰度变化下构造区域相似性的度量标准，后者则在于对MR/TRUS图像进行分割，并通过特定的标记点或直接对整个前列腺切片，进行配准融合. 然而，由于前列腺的自动分割是一项极具挑战的任务，因此现有的融合方法主要是基于手动或半自动的方式对前列腺MR/TRUS图像进行分割，并基于手动放置或特殊的生理特征点进行刚体或非刚体配准[10-11].

为了充分发挥前列腺MR/TRUS图像在肿瘤诊断、引导穿刺过程中的优势，对前列腺病人术前的MR图像进行手动分割并标记穿刺区域，在前列腺TRUS图像实时引导穿刺的过程中，将术前标注的MR图像信息配准到TRUS图像上[12]，即采用术前MR图像手动分割的前列腺区域到术中TRUS图像自动分割的前列腺区域的配准方式进行图像融合引导穿刺. 考虑到前列腺TRUS图像在分割过程中形状尺寸的变化及纹理灰度的影响，采用基于形状和灰度变化的混合方法建立前列腺TRUS图像的活动表观模型[13]，用于图像的自动分割，同时采用薄板样条（Thin plate spline, TPS）[14]对前列腺MR/TRUS图像上对应的轮廓区域进行配准融合.

本文针对活动表观模型对轮廓进行自动分割时姿态参数的初始化问题，提出一种基于监督学习的超声图像自动分割方法，并基于随机森林[15]的分类结果建立了姿态参数边界驱动的数学模型. 在配准过程中，针对Kuhn-Munkres（KM）算法进行改进以降低局部匹配代价，同时引入轮廓特征点的定位误差[16]作为薄板样条的正则因子进行图像配准. 最后，本文将分割结果和标准轮廓进行对比，同时与其他方法比较了配准精度.

1 基于监督学习的图像分割方法

应用活动表观模型进行前列腺轮廓分割的过程中，需要解决姿态参数的初始化问题，考虑基于随机森林实现姿态参数的自动估计，以完成前列腺超声图像的自动分割.

1.1 前列腺活动表观模型

活动表观模型是一种由形状模型和灰度模型建立的组合模型，用于对可变形状进行数学描述.

本文采用多幅前列腺超声图像轮廓上的40个轮廓标签点训练前列腺活动表观模型，建立形状模型和灰度模型：

对于一幅新的超声图像来说，采用活动表观模型进行分割的过程即是求解图像过程，因此需要估计式（3）中对应的组合参数. 在估计组合参数的过程中，Cootes建立了多元线性回归模型[13]，如式（4），认为在前后两次迭代更新的过程中，组合参数的更新量与 灰度变化量呈线性关系：

1.2 前列腺超声图像随机森林模型

图1 初始姿态 对收敛结果的影响. （a～d）初始姿态参数过大无法收敛；（e～f）初始姿态满足收敛条件Fig.1 Effect of the initial position parameter on the results of convergence: (a-d) large initial position parameter resulted error convergence; (e-f) initial position parameter met convergence conditions

为了准确的估计一张新的超声图像上前列腺区的姿态参数，建立前列腺超声图像的随机森林模型，对图像上的灰度特征进行二分类，从而实现前列腺图像的预分割. 在图像数据的训练过程中，为了减轻超声图像姿态变化和不同病人图像的灰度差异的影响，采用增强系数（Enhanced correlation coefficient, ECC）算法[17]对所有的超声图像进行通道对齐，同时进行像素灰度的归一化，区间值为[0,1]. 对每一张训练图像提取特征标签P=(X,V)，其中X=(xima,yima)表示像素点的图像坐标，V=(d1,d2,a1,···,a9)表示以为中心的3×3邻域的均值、标准差和邻域灰度值(a1,···,a9).

对训练数据集有放回的随机抽样，并随机选择特征子集，以信息增益作为特征的选择指标生产每一棵子树，当信息增益达到最小值或者达到最大指定深度时停止生长. 在随机森林中，每一棵子树的结构都是不同的，因此在分类过程中，每一个测试样本都会在不同的子树中得到不同分类结果，统计出现次数最多的决策结果作为随机森林的最终分类结果. 前列腺随机森林模型的训练和预测过程如图2所示，对于一张新的超声图像，随机森林分类器将输出前列腺区域预测结果的二值图像.

1.3 估计形状模型的姿态参数

图2 前列腺随机森林模型的训练和预测过程Fig.2 Training and prediction of the prostate random forest model

图3 均值形状的初始姿态Fig.3 Initial position parameters of the mean shape model

图4 姿态 下的坐标变化关系Fig.4 Coordinate transformation relationship at position

图5 生成最小外接矩形Fig.5 Generation of the minimum enclosing rectangle

图6 前列腺超声图像的自动分割过程. （a）前列腺TRUS图像；（b）参数寻优结果；（c）图像分割结果；（d）分割对比结果Fig.6 Automatic segmentation process of prostate TRUS images:(a) prostate TRUS image; (b) parameters optimization result; (c) segmentation result; (d) image segmentation comparison results

2 基于形状特征的图像配准方法

将由第一部分自动分割的TRUS图像与术前分割的MR图像在引入定位误差的TPS框架下进行融合配准，为了自动提取特征点的配对信息，构造了形状特征算法，并基于改进的KM算法实现求解.

2.1 提取轮廓形状特征点

图7为一组待配准的前列腺MR/TRUS图像，其中前列腺MR图像轮廓为术前手动分割，前列腺TRUS图像轮廓由第1部分方法自动分割，为了提取轮廓的ShapeContext特征[20]，需要对轮廓离散采样，并建立每一点的ShapeContext描述子作为轮廓点的特征矢量.

图7 待配准的前列腺MR/TRUS图像. （a）MR图像轮廓点；（b）TRUS图像轮廓点Fig.7 Prostate MR/TRUS images to be registered: (a) contour points on MR image; (b) contour points on TRUS image

通常采用KM算法[21]求解上述指派问题，由于KM算法是针对完备匹配问题进行求解，因此部分奇异的轮廓点会“误匹配”，从而产生较大的局部误差. 本文对KM算法进行改进，即仅考虑待匹配点的最优和次优解的增广路径，如下所示，并对改进前后的算法的匹配结果进行对比，如图9所示，改进后的KM算法消除了部分轮廓奇异点，提高了匹配点对的局部精度.

图8 形状描述符的建立Fig.8 Construction of the shape descriptor

图9 改进KM算法的对比结果Fig.9 Results compared with the improved KM algorithm

Step 2: Search corresponding matching point forand the objects considered are onlyandwherehas higher priority than. If the augmented paths produced fromandhave been considered, give up the matching of current point.

Step 3: Output the remaining answers of matching of

2.2 引入各向异性误差薄板样条配准方法

薄板样条作为一种常用于生物图像变形配准的插值工具，在医学图像配准方面有着广泛的应用. 使用薄板样条对变形场进行估计时，需要指定多组配对特征点进行约束. 考虑到轮廓分割中产生的人为误差，引入配对特征点基于图像梯度的各向异性误差作为正则因子，进行图像配准.

式中：m表示特征点为中心的邻域像素，为高斯核的方差. 对邻域m上的像素计算高斯第一阶导数，记作，且满足：

图10 图像配准融合结果. （a）MR图像的变换结果；（b）MR/TRUS图像融合结果Fig.10 Registration results: (a) transformation result of MR image;(b) registration result of MR/TRUS images

3 实验结果

为了验证所提方法在超声图像自动分割、MR/TRUS图像配准方面的精度，进行了两组对比实验. 第一组实验采用监督学习的初始化方法并基于前列腺活动表观模型完成图像分割，并与影像专家分割的前列腺轮廓进行对比，同时计算轮廓之间豪斯多夫距离[24]（Hausdorff distance, HD）. 第二组实验对比分析基于形状轮廓配准方法与传统方法的实验结果，比较的指标包括配准结果的区域重叠戴斯系数DSC[25]及尿道特征点的配准定位精度（Target error, TE）. 实验图像取自美国影像引导治疗中心提供的数据集：https://zenodo.org/record/16396#.XPJ6GaIzZdh，收集了前列腺患者的MR/TRUS影像图片.

3.1 图像分割实验

模型的数据标签借助T.F.Cootes的工具箱完成，图像的训练和收敛过程借助OpenCV视觉库，采用C++语言实现，遗传算法的优化工作借助了Python下的Geatpy包完成.

训练采用的数据集为3例前列腺病人中层切片的MR/TRUS影像图片，其中18张作为训练集，5张作为测试集，图像的分割结果如图11所示，其中图11（a1～a5）为待分割的超声图像，图11（b1～b5）为随机森林的预分割结果，图11（c1～c5）为轮廓收敛的迭代过程. 采用豪斯多夫距离对轮廓相似度进行评估.

在算法的实现过程中，首先基于随机森林模型对测试图像进行预分割，如图11（b1～b5），其中预分割的区域用蓝色轮廓标出，标准区域用黑色轮廓标出. 从图11（b1～b5）可以看出预分割的轮廓边界比较粗糙，但能初步定位前列腺的区域.图11（c1～c5）表示预分割轮廓基于遗传算法的优化结果逐步收敛的过程. 其中蓝色轮廓表示前列腺形状均值模型的基于预分割图像的初始姿态，红色轮廓是采用遗传算法进行姿态优化后的结果，粉色轮廓是基于活动表观模型的收敛结果，黑色轮廓为手动绘制的标准轮廓. 从图11（c1～c5）的分割结果可以看出，轮廓的最终收敛结果与标准轮廓具有高的轮廓相似度，其HD均值分别为1.20，1.04，2.48，1.40, 1.57，可以得出分割轮廓与标准轮廓之间的豪斯多夫距离的均值较低为1.54，证明了分割方法的有效性.

3.2 图像配准实验

在评估配准方法精度的过程中，使用超声图像标准轮廓区域与术前核磁图像轮廓区域进行配准. 实验选取前列腺病人9组MR/TRUS影像进行图像配准，同时，与文献方法[26]进行对比，并评估配准区域重叠度DSC系数，其结果如图12和图13所示.

图12第1行和第2行分别为前列腺核磁和超声图像，第3行和第4行是本文方法计算出来的配准结果和融合结果，第5行和第6行是采用文献方法[26]得出的配准和融合结果. 图13对比了两种配准结果的区域重叠度DSC，结果表明本文方法具有良好的配准精度.

图11 前列腺超声图像分割结果. （a1～a5）待分割的TRUS图像；（b1～b5）随机森林预分割结果；（c1～c5）轮廓分割收敛过程Fig.11 Segmentation results of prostate US images: (a1-a5) initial TRUS images; (b1-b5) pre-segmentation results of random forest; (c1-c5) convergence processes of contour segmentation

图12 核磁超声图像配准结果对比Fig.12 Results of MR/TRUS images registration

图13 本文方法与传统方法对比结果Fig.13 Results compared with the traditional method

为了进一步验证配准方法的目标定位精度，即核磁超声图像上对应特征点之间的配准精度，选择5组前列腺MR/TRUS切片，以尿道口作为参考点进行评估，其结果如图14所示，图像红色点表示尿道口处对应特征点，蓝色点和黄色点分别为本文方法和传统方法定位结果. 使用软件3Dslicer测量尿道口配准定位误差，并进行对比平均定位精度AP和标准差，结果如表1.

从表1中可以看出，本文方法针对尿道特征点的定位精度为1.59±0.23 mm，临床对可疑肿瘤区域的穿刺误差在3 mm以内均可视为有效穿刺，认为配准结果可以满足临床需求.

图14 尿道特征点定位结果Fig.14 Location results of urethral points

表1 尿道特征点定位结果对比Table 1 Comparison of location results of urethral points

4 结论

（1）针对前列腺核磁超声融合引导穿刺手术中前列腺的变形配准问题，提出一种基于监督学习的超声图像自动分割方法，并与术前核磁图像进行图像配准.

（2）新方法通过引入随机森林分类器建立了边界驱动的姿态估计模型，实现前列腺超声图像的自动分割，与专家轮廓对比具有高的精度.

（3）在图像配准方面，使用了形状矢量来构建薄板样条的配对特征点同时引入各向异性误差作为正则因子，配准结果表明，与传统方法相比，新方法在图像的配准方面具有较高的精度，在前列腺核磁超声图像融合引导方面具有临床应用价值，是一种准确、稳定的图像配准方法.