基于统计形状模型的内耳三维水平集分割方法与流程
本发明涉及数字图像处理技术领域，特别涉及一种基于统计形状模型的内耳三维水平集分割方法。
背景技术：
内耳器官分割是人工耳蜗植入手术的术前影像学评估的重要环节，也是准确测量内耳结构的必要手段。内耳器官分割可以在术前提供可靠的解剖信息，为耳显微外科手术计划的制定提供依据，使耳显微外科疾病的诊断与治疗更加科学、更加直观，同时，基于标准空间坐标的内耳模型对于耳石病的理论研究和临床诊断具有重要意义。
如上所述，内耳的3d分割对临床研究和应用具有很高的价值。但是内耳图像的强度分布不均匀，结构精小，同时具有复杂的拓补管状结构，最小截面直径约为0.8mm，因此，内耳的分割仍然是非常具有挑战性的。
目前最常见的方法是半自动分割方法，如三维窄带水平集分割算法，再结合手动精细分割。但很多手工操作导致分割速度不理想。直到近几年，自动分割方法逐渐发展起来。在相关文献中，提出有一种基于k-means聚类的预分割结果作为可变形模型的初始轮廓来分割内耳磁共振图像(mri)的自动分割过程。但k-means聚类的方法忽略了内耳边界的模糊特征，容易引起边缘泄漏或感兴趣对象弱边缘分割不足的问题。在其他文献中还公开了，采用三维hessian矩阵对内耳的管状结构增强，然后通过模糊聚类方法对增强结果进行分割。最后，使用混合水平集分割方法来优化结果。然而，内耳不是完全管状的结构，hessian矩阵也可以抑制内耳的非管状结构。同时模糊聚类方法的收敛速度缓慢，影响分割速度，不是最好的预分割方法。
近年来，统计形状模型(ssm)已成为医学图像分割领域的一个热点。noble等人，使用基于ssm的方法在常规ct图像中分割内耳耳蜗。将尸体的高分辨率微米ct图像作为训练样本，建立统计形状模型，然后将模型投影到待分割的目标对象上获得初始轮廓，并通过非刚性配准的方法将初始轮廓驱动到更精准的目标边界处。但该方法用于建模的训练样本尸体耳蜗标本图像数据比较难获得，因此，ssm难以捕获足够的耳蜗变化模式去构建一个普适性好的形状模型。同时三维刚性配准的方法对于内耳的分割也比较耗时。除此之外，传统脑部ct图像只能够清晰地显示内耳的骨迷路结构，而脑部mri图像可以更好的显示膜迷路结构，膜迷路结构能够帮助诊断某些膜迷路疾病。目前，采用统计形状模型对内耳膜迷路进行分割的研究工作目前还比较少。
技术实现要素：
本发明所要解决的技术问题在于针对上述现有技术中的不足，提供一种基于统计形状模型的内耳三维水平集分割方法。
为解决上述技术问题，本发明采用的技术方案是：一种基于统计形状模型的内耳三维水平集分割方法，包括以下步骤：
步骤一：建立内耳的统计形状模型；
步骤二：通过体数据到体数据的刚性配准方法获取内耳的感兴趣区域，然后通过模型到体数据的刚性配准方法将统计形状模型的平均形状模型配准到内耳的感兴趣区域上，获得内耳的初始表面轮廓；
步骤三：利用基于阈值区域的三维水平集分割方法，以内耳的初始表面轮廓作为初始演化表面，进行水平集演化，得到目标轮廓表面，再结合原始脑部mri体数据最终计算得到需要的内耳轮廓。
优选的是，所述步骤一具体包括：
步骤1-1：采集足够数量的内耳体数据作为训练样本，再通过3d-slicer软件将其转换为均匀分布的三角形表面网格数据；
步骤1-2：任意选取训练样本中的一个作为参考表面网格数据，利用普氏分析法将其他训练样本中的训练形状与参考表面网格数据的形状对齐；
步骤1-3：利用3维刚性配准的方法建立所有经步骤1-2对齐处理后的训练样本之间的一对一的点对应关系，再任意选择其中一个训练样本作为参考体数据，其余的训练样本通过使用itk配准框架配准到参考体数据进行模型拟合，完成所有训练样本的最终对齐配准；
步骤1-4：定义由步骤1-3得到的训练样本si的样本总数为n个，其中i∈{0,1,2,......,n-1}，每个训练样本si都是一个点云分布的表面网格数据，表示为si＝{x0,y0,z0,x1,y1,z1,......,xk-1,yk-1,zk-1}，其中，x,y,z表示训练样本si的数据的坐标点，k表示训练样本si的数据的坐标点的总数；然后利用主成分分析法对训练样本si进行形状统计分析，得到内耳的统计形状模型v，其数学表达式如下：
其中，v(α1,......,αm)是个矢量，α表示每个空间位置点的坐标，
通过调节形状参数，以获得与训练样本集相似的新的形状个体。
优选的是，所述步骤二具体包括：
步骤2-1：采集原始脑部mri体数据作为待分割的测试体数据，从待分割的测试体数据中随机选择一个作为参考测试体数据，再通过手动交互式方法标定一个包含内耳的感兴趣区域矩形立方体掩膜，然后利用刚性配准的方法将该包含内耳的感兴趣区域矩形立方体掩膜自动的由参考测试体数据上传播到其他的待分割的测试体数据上，实现自动获得所有待分割体的测试体数据的内耳的感兴趣区域；
步骤2-2：通过由步骤一建立的内耳的统计形状模型产生多个新的形状个体，并计算得出产生的多个新的形状个体的平均形状模型，再利用模型到体数据的三维刚性配准方法将该平均形状模型配准到内耳的感兴趣区域，获得内耳的初始表面轮廓。其中，多个新的形状个体的平均形状模型的计算通过结合步骤1中的公式
优选的是，所述步骤三具体包括：
步骤3-1：构建高维水平集函数ψ(x,t)，驱使高维水平集函数以上述步骤二得到的初始表面轮廓作为初始演化表面进行表面演化，高维水平集函数ψ(x,t)的表面演化的数学表达式为：
其中▽是梯度算子，|▽ψ|为梯度算子的绝对值，
步骤3-2：在速度函数f的驱使下，使高维水平集函数ψ(x,t)逐步向目标轮廓表面的边界处传播，进行表面演化，到达边界处后停止，完成表面演化；再计算出高维水平集函数ψ(x,t)的零水平集，从而获得需要的目标轮廓表面γ(x,t)；
其中速度函数的数学表达式如下：
f＝αp(x)+βk(x)；
其中p(x)是传播项，k(x)是曲率项，α和β是与之对应的权重；
其中，目标轮廓表面γ(x,t)＝{ψ(x,t)＝0}；
步骤3-3：将由步骤3-2得到目标轮廓表面作为掩膜与由步骤二采集的原始脑部mri体数据做乘法计算，实现对采集的原始脑部mri体数据的分割，最终得到需要的内耳体数据。
本发明的有益效果是：本发明使用一些训练样本来构建统计形状模型，然后通过刚性配准自动将生成模型的平均形状模型匹配到待分割的目标体积上，最后，使用基于阈值区域的水平集分割方法来细化分割。由于内耳结构的复杂性，通常采用的球体或立方体等简单初始表面难以演化收敛到真实的内耳边界。然而，本发明的基于统计形状模型的内耳三维水平集分割方法可以弥补这一缺点。它包含了通过统计分析多个训练形状获得的先验形状信息，因此能更准确地接近感兴趣对象的实际边界以克服分割泄漏或分割不足的问题。本发明中使用的训练数据是脑部mri图像，容易获得大量的训练样本来建立更好的模型。本发明通过内耳训练样本的统计分析获得的ssm作为先验知识，已经十分接近内耳的真实边界，因此，基于统计形状模型的内耳三维水平集分割方法可以快速收敛到真正的内耳边界，最终的分割结果具有比基于区域生长的水平集分割方法更高的精度，本发明可以达到临床应用的精度。同时由于固有的先验信息，本发明的方法在对局部图像的伪影和扰动的抗干扰性方面明显要优于传统的低水平算法。并且水平集方法可以平滑曲面演化并得到一个光滑的内耳外表面。本发明方法的复杂度要低于通过各种预分割方法组合获得初始轮廓的方法，且本发明的方法在内耳分割中具有很高的鲁棒性和准确性。
附图说明
图1为本发明的基于统计形状模型的内耳三维水平集分割方法的一种实施例的流程简图；
图2为采用手工分割方法得到的最终内耳轮廓图的示意图；
图3为采用本发明的基于阈值区域的水平集分割方法得到的最终内耳轮廓示意图；
图4为采用手工分割方法得到的最终内耳轮廓示意图。
具体实施方式
下面结合实施例对本发明做进一步的详细说明，以令本领域技术人员参照说明书文字能够据以实施。
应当理解，本文所使用的诸如“具有”、“包含”以及“包括”术语并不排除一个或多个其它元件或其组合的存在或添加。
本实施例的一种基于统计形状模型的内耳三维水平集分割方法，其特征在于，包括以下步骤：
步骤一：建立内耳的统计形状模型。
步骤一具体包括：
步骤1-1：采集足够数量的内耳体数据作为训练样本，再通过3d-slicer软件将其转换为均匀分布的三角形表面网格数据；其中，可由经验丰富的耳科专家从脑部mri图像中手动分割出足够数量的内耳体数据作为训练样本；
步骤1-2：任意选取训练样本中的一个作为参考表面网格数据，利用普氏分析法将其他训练样本中的训练形状与参考表面网格数据的形状对齐；
步骤1-3：利用3维刚性配准的方法建立所有经步骤1-2对齐处理后的训练样本之间的一对一的点对应关系，再任意选择其中一个训练样本作为参考体数据，其余的训练样本通过使用itk配准框架配准到参考体数据进行模型拟合，完成所有训练样本的最终对齐配准；
步骤1-4：定义由步骤1-3得到的训练样本si的样本总数为n个，其中i∈{0,1,2,......,n-1}，每个训练样本si都是一个点云分布的表面网格数据，表示为si＝{x0,y0,z0,x1,y1,z1,......,xk-1,yk-1,zk-1}，其中，x,y,z表示训练样本si的数据的坐标点，k表示训练样本si的数据的坐标点的总数；然后利用主成分分析法对训练样本si进行形状统计分析，得到内耳的统计形状模型c，来简洁的描述和表达形状，可以通过这个统计形状模型来产生一些新的形状，或者来判断一个新的形状是否符合统计形状模型的描述，统计形状模型v数学表达式如下：
其中，v(α1,......,αm)为矢量，α表示每个空间位置点的坐标，
v为向量，通过调节形状参数，以获得与训练样本集相似的新的形状个体。
步骤二：通过体数据到体数据的刚性配准方法获取内耳的感兴趣区域，然后通过模型到体数据的刚性配准方法将统计形状模型的平均形状模型配准到内耳的感兴趣区域上，获得内耳的初始表面轮廓。本发明中的三维刚性配准需要两个体数据图像作为输入。一个是被选为固定图像的参考图像，另一个是待变换的浮动图像，利用互信息作为测度函数，使用梯度下降算法作为优化器。
所述步骤二具体包括：
步骤2-1：采集原始脑部mri体数据作为待分割的测试体数据，从待分割的测试体数据中随机选择一个作为参考测试体数据，再通过手动交互式方法标定一个包含内耳的感兴趣区域矩形立方体掩膜，然后利用刚性配准的方法将该包含内耳的感兴趣区域矩形立方体掩膜自动的由参考测试体数据上传播到其他的待分割的测试体数据上(体数据到体数据的刚性配准方法)，实现自动获得所有待分割体的测试体数据的内耳的感兴趣区域；之后的分割均是在roi(感兴趣区域)内进行，这样减少数据量，可以大大提高处理速度。
步骤2-2：通过由步骤一建立的内耳的统计形状模型产生多个新的形状个体，并计算得出产生的多个新的形状个体的平均形状模型，再利用模型到体数据的三维刚性配准方法将该平均形状模型配准到内耳的感兴趣区域，获得内耳的初始表面轮廓(模型到体数据)。其中，平均形状模型作为浮动图像，向目标体数据进行映射。
步骤三：利用基于阈值区域的三维水平集分割方法，以内耳的初始表面轮廓作为初始演化表面，进行水平集演化，使结果更靠近目标轮廓表面的边界，以得到目标轮廓表面，再结合原始脑部mri体数据最终计算得到需要的内耳体轮廓。
所述步骤三具体包括：
步骤3-1：构建高维水平集函数ψ(x,t)，水平集方法一种跟踪轮廓和表面运动的一种数字化方法，不直接对轮廓进行操作，而是将轮廓设置成一个高维函数ψ(x,t)的零水平集；再驱使高维水平集函数以上述步骤二得到的初始表面轮廓作为初始演化表面进行表面演化，高维水平集函数的运动(即表面演化)成为一个微分方程，高维水平集函数ψ(x,t)的表面演化的数学表达式为：
其中▽是梯度算子，|▽ψ|为梯度算子的绝对值，
步骤3-2：在速度函数f的驱使下，使高维水平集函数ψ(x,t)逐步向目标轮廓表面的边界处传播，进行表面演化，到达边界处后停止，完成表面演化；再计算高维水平集函数ψ(x,t)的零水平集，从而获得需要的目标轮廓表面γ(x,t)；
其中速度函数的数学表达式如下：
f＝αp(x)+βk(x)
其中p(x)是传播项，k(x)是曲率项，α和β是与之对应的权重；传播项与特征图像的强度成正比。在本方法中，特征图像是通过将原始图像中指定的阈值区域内的灰度图像映射到0至1范围，而阈值区间之外的图像灰度映射到-1至0范围内所获得。阈值范围内的图像作为前景，其余作为背景。所以本发明要根据经验指定出具体的待分割的目标对象的灰度分布阈值范围。如果图像像素值在阈值范围内，p(x)就取正值，反之则为负值。因为速度项p(x)能够实现在前景内扩张，在背景处收缩，直至收敛到内耳边界。本方法的曲率项k是采用常用的平均曲率计算方法，它可以控制演化表面的光滑度。此处的演化表面指的是嵌入到高维水平集函数中的零水平集，即目标轮廓表面。在这两项的驱动下，速度函数就会驱使表面在灰度均匀的区域以常量速度演化，而在靠近边界处减慢速度直至边界处停止演化。演化停止后，及代表高维水平集的轮廓的表面演化到达了其中嵌入的零水平集(演化表面)的目标边界处。然后，演化表面可以通过计算零水平集获得，其数学表达式如下：
γ(x,t)＝{ψ(x,t)＝0}；
该表面γ(x,t)即是我们需要的目标轮廓表面。
上述步骤3-1和3-2即为利用基于阈值区域的三维水平集分割方法对内耳进行分割的具体步骤。
步骤3-3：将由步骤3-2得到目标轮廓表面作为掩膜(实际上得到的是个二值化的图像，轮廓内的像素值设为1，轮廓外的额设为0)与由步骤二采集的原始脑部mri(核磁共振成像)体数据做乘法计算，实现对原始脑部mri体数据的分割，最终得到需要的内耳体数据，即最终的内耳轮廓。
如图1所示，还提供了本发明的基于统计形状模型的内耳三维水平集分割方法的一种实施例的流程简图。
图2、图3、图4给出了本发明的一种实施例中，采用手工分割方法、本发明的基于统计形状模型的内耳三维水平集分割方法和基于区域生长的水平集分割方法得到的最终分割结果的内耳轮廓的比较示意图。其中，图2为采用手工分割方法得到的最终内耳轮廓示意图，图3为采用本发明的基于阈值区域的水平集分割方法得到的最终内耳轮廓示意图，图4为采用手工分割方法得到的最终内耳轮廓示意图。
由图2、图3、图4可以看出，本方法能够清晰的分辨三个半规管，前庭和螺旋耳蜗结构，并且与手工分割结果比较有很好的相似性。从图中可以看出，基于区域生长的水平集分割方法没有清理干净内耳道部分，但是本发明的基于统计形状模型的内耳三维水平集分割方法成功处理了这个问题，并且获得完整精细的内耳结构。表1给出了本发明的基于统计形状模型的内耳三维水平集分割方法得到的最终内耳轮廓图和基于区域生长的水平集分割方法得到最终内耳轮廓图的定量分析的比较结果。
表1
表格中分割结果评价指标说明：
dsc(dicesimilaritycoefficient戴斯相似性系数)和mcc(matthew'scorrelationcoefficient马修斯相关系数)，通过计算像素重叠率来验证自动分割和手动分割之间的相似性。如果它们的值接近1，那么我们将得到更高的分割精度。fpr(假阳率)和fnr(假阴率)分别代表过分割和欠分割的程度，值越接近0，过度分割或欠分割的程度越低。
表格说明：这些评估指标是在每对采用算法分割和耳科医师手动分割(金标准)之间计算的。结果如表1所示，它们是9对(18个内耳)分割结果的平均值。
从表1可以看出，本发明的基于统计形状模型的内耳三维水平集分割方法的dsc平均值为0.9591，优于基于区域生长的水平集分割法的0.8745。对于大多数结构，通常认为0.8的dsc是很好的，因此，本发明的基于统计形状模型的内耳三维水平集分割方法的结果可以收敛到更加贴近解剖学正确地真实目标的边界。mcc也代表两个体数据的相似性。从结果可以看出，本发明的基于统计形状模型的内耳三维水平集分割方法与金标准有较高的相似性。本发明的基于统计形状模型的内耳三维水平集分割方法可以有效地提高体数据重叠率。本发明的方法的fpr和fnr的值较低，表明本发明相较于基于区域生长的分割方法的过分割和欠分割的程度更轻。
感兴趣目标的常见变体信息通过统计分析多个训练形状被包括在平均形状模型中，因此平均形状模型更接近感兴趣对象的实际边界，能为后面水平集精分割提供一个更加精准的初始表面轮廓，因此，本发明可以得到更令人满意的结果。
另外，分割结果仍然存在一些误差。主要有两个原因：第一，由于内耳结构微小，mri图像分辨率有限，手动分割可能有些粗糙，并不能完全代表“金标准”，只是近似值；第二，医学图像的固有噪声等因素引起的误差。
本发明通过配准的方法创新地将统计形状模型的平均形状模型作为初始表面轮廓添加到阈值区间的水平集分割演化过程中，通过内耳训练样本的统计分析获得的ssm作为先验知识，已经十分接近内耳的真实边界，因此，基于统计形状模型的内耳三维水平集分割方法可以快速收敛到真正的内耳边界，最终的分割结果具有比基于区域生长的水平集分割方法更高的精度，本发明可以达到临床应用的精度。同时由于固有的先验信息，本发明的方法在对局部图像的伪影和扰动的抗干扰性方面明显要优于传统的低水平算法。并且水平集方法可以平滑曲面演化并得到一个光滑的内耳外表面。本发明方法的复杂度要低于通过各种预分割方法组合获得初始轮廓的方法，且本发明的方法在内耳分割中具有很高的鲁棒性和准确性。
相较于之前已有的内耳分割方法，本文主要创造性的将统计形状模型与水平集分割结合到一起使用，同时目前内耳分割的相关研究还比较少，所以本发明极具探索价值。
本发明使用一些训练样本来构建统计形状模型，然后通过刚性配准自动将生成模型的平均形状模型匹配到待分割的目标体积上，最后，使用基于阈值区域的水平集分割方法来细化分割。由于内耳结构的复杂性，通常采用的球体或立方体等简单初始表面难以演化收敛到真实的内耳边界。然而，本发明的基于统计形状模型的内耳三维水平集分割方法可以弥补这一缺点。它包含了通过统计分析多个训练形状获得的先验形状信息，因此能更准确地接近感兴趣对象的实际边界以克服分割泄漏或分割不足的问题。本发明中使用的训练数据是脑部mri图像，容易获得大量的训练样本来建立更好的模型。本发明方法的复杂度要低于通过各种预分割方法组合获得初始轮廓的方法，且本发明的方法在内耳分割中具有很高的鲁棒性和准确性。
尽管本发明的实施方案已公开如上，但其并不仅仅限于说明书和实施方式中所列运用，它完全可以被适用于各种适合本发明的领域，对于熟悉本领域的人员而言，可容易地实现另外的修改，因此在不背离权利要求及等同范围所限定的一般概念下，本发明并不限于特定的细节。
- 还没有人留言评论。精彩留言会获得点赞！