一种基于凹形缝隙区域识别的终板自动分离方法与流程

本发明涉及医学图像分析及计算机辅助技术领域，特别是一种基于凹形缝隙区域识别的终板自动分离方法。

背景技术：

终板位于椎间盘与椎体中心之间，目前医学研究认为终板退化与椎间盘退变和骨质疏松存在紧密的相关性。为了深入分析椎间盘退变的具体机制，医学上常用小鼠动物模型，从高分辨率的小鼠脊椎图片中将终板组织分离出来，获取终板的精确形态用于病理分析。终板分离目前主要由医生根据经验进行手工分割。但小鼠脊椎图片量巨大，手工分离效率太低；同时小鼠终板是一块弯曲薄板，在尺寸上属于显微组织，且形态不规则，分割精度要求高；尤其在终板与椎体中心之间的狭窄空隙区域中，存在大量存在可见骨质连接。因此要从脊椎中分离出微型复杂的终板是一项具有挑战性的课题。

技术实现要素：

有鉴于此，本发明的目的是提出一种基于凹形缝隙区域识别的终板自动分离方法，能获取终板的精确边界信息和形态信息，鲁棒性好，自动化程度高，同时发明成果对研究椎间盘退变具体机制，以及终板和骨质疏松的相关性提供了科学依据。

本发明采用以下方案实现：一种基于凹形缝隙区域识别的终板自动分离方法，具体包括以下步骤：

步骤s1：对脊椎图像阈值化，获取图像中的骨骼对象；

步骤s2：采用自顶向下的方式确定终板感兴趣区域，去除终板无关数据；

步骤s3：自动识别终板与椎体中心相连的缝隙区域初始分界曲面；

步骤s4：对感兴趣区域骨架化，将骨组织之间的连接关系和拓扑进行压缩描述，且骨架化得到的图是后续图构建的基础和依据；

步骤s5：针对骨架化的图进行图构建，设计图割的代价函数，将分界曲面作为初始点，完成图割，识别出终板和椎体中心之间的相连骨小梁；

步骤s6：采用迭代切割线方法沿着相连骨小梁，搜索出最优位置和角度的线段，完全去除相连骨小梁，并进行三维重建获取终板模型。

进一步地，步骤s2具体包括以下步骤：

步骤s21：分离单个脊椎，由于图像的高分辨率，只要根据单个椎体内部骨骼的连通性，采用区域增长法能够很容易将脊椎与脊椎分离；

步骤s22：使用椎体两端投影法将脊椎背后的脊椎突去除，拟合出覆盖投影的圆柱形区域，其侧面的曲面形成分离曲面；

步骤s23：相比椎体中心内部的规律孔洞的松质骨，终板内部为致密组织，采用骨密度统计法统计椎体两端骨密度变化曲线，设置阈值获取密度最大且变化最大的区域作为终板出现的区域，即终板感兴趣区域。

进一步地，步骤s3具体包括以下步骤：

步骤s31：对终板感兴趣区域进行三维重建得到三角面片化的模型对象；

步骤s32：在模型表面构建调和场，调和场计算产生一系列候选等值曲线；

步骤s33：设计等值曲线分数计算函数，根据分数结果选出分数最高的候选等值线作为最终表面分界曲线；

步骤s34：结合终板远椎体中心端的外边缘曲线和表面分界曲线，采用三次b样条方法拟合出终板和椎体中间之间的初始分界曲面。

进一步地，步骤s34中，为了选出最贴合终板和椎体中心之间间隙区域的等值线，在候选等值线分数函数score(i)中引入平均梯度因子grd(i)和形状因子shape(i)设计，其中形状因子shape(i)以第i等值线为中心，沿着跨等值线方向，逐面片向左右两边扩展m步，分别得到左边和右边第j个扩展环，对应的长度分别为lij(left)和lij(right)，而对应的形状因子的函数为lmax是所有候选等值线长度最大值。

进一步地，其中形状因子计算使用了高斯多尺度卷积，可以保证对选择不敏感，将m设置为5，将σ设置为2。

进一步地，步骤s5具体包括以下步骤：

步骤s51：对骨架化的图像进行图构建，其中图的顶点代表骨架分叉结点，图中的两两结点之间的弧是反应骨小梁长度和方向的曲线路径；

步骤s52：根据相连骨小梁和椎体中心内部的骨小梁的形态差别特点，即相连骨小梁更接近椎体主方向且长度较短，因此将弧的权值函数设置为骨小梁长度与骨小梁和椎体主方向夹角绝对值之和，从而识别的相连骨小梁。

进一步地，步骤s6中，所述迭代切割线方法具体为：切割线从相连小梁中点开始，朝终板匀速靠近，当切割线到达终板表面时，切割线覆盖的骨小梁距离会急剧增加，当增加速度超过阈值20时，就代表切割直线接触到了终板表面，然后迭代调整切割线角度，将和调整差距最小的方向作为切割线最终方向。

进一步地，步骤s6具体包括以下步骤：

步骤s61：在单个相连小梁上构建一根切割线，切割线初始位置位于相连小梁骨架的中点，角度和相连小梁垂直，并记录下切割线与相连小梁区域相交形成的初始端点距离d1；

步骤s62：切割线沿着骨架超终板方向匀速移动，每次平移i，获得相应的端点距离di，且该距离会随着平移发生变化；当切割线到达终板表面时端点距离会急剧增加，当超过阈值20时，便停止前进，获得最终切割线位置；

步骤s63：切割线停止移动，重新调整其方向，每次旋转一定角度，获得新的端点距离，当新的端点距离和旋转前端点距离差最小时，当前角度即为终板表面相切角度；

步骤s64：按照步骤s63确定的最终的切割线位置和角度，终板上的相连小梁得以完全切除；

步骤s65：针对终板数据，使用marchingcube算法进行三维重建，得到终板的三维模型。

与现有技术相比，本发明有以下有益效果：

1、本发明无需用户交互，包括分界曲面的提取和相连小梁的识别在内所有操作都能实现自动化。

2、本发明的自顶向下的终板兴趣区域提取去除了终板无关数据，有效地提升系统效率。

3、本发明对兴趣区域骨架化提取，在保证骨骼内部的连通性和拓扑完整表达的前提下，压缩了目标数据，极大地减少了图割方法的计算时间。

4、本发明能获取终板的精确边界信息和形态信息，鲁棒性好，实现方便。

5、本发明对研究椎间盘退变具体机制，以及终板和骨质疏松的相关性提供了科学依据。

附图说明

图1为本发明实施例的方法流程示意图。

具体实施方式

下面结合附图及实施例对本发明做进一步说明。

应该指出，以下详细说明都是示例性的，旨在对本申请提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本申请所属技术领域的普通技术人员通常理解的相同含义。

需要注意的是，这里所使用的术语仅是为了描述具体实施方式，而非意图限制根据本申请的示例性实施方式。如在这里所使用的，除非上下文另外明确指出，否则单数形式也意图包括复数形式，此外，还应当理解的是，当在本说明书中使用术语“包含”和/或“包括”时，其指明存在特征、步骤、操作、器件、组件和/或它们的组合。

如图1所示，本实施例提供了一种基于凹形缝隙区域识别的终板自动分离方法，具体包括以下步骤：

步骤s1：对脊椎图像阈值化，获取图像中的骨骼对象；

步骤s2：由于原始图像中常包含好几个脊椎，因此采用自顶向下的方式确定终板感兴趣区域，去除终板无关数据，能够提高处理效率；

步骤s3：自动识别终板与椎体中心相连的缝隙区域初始分界曲面；

步骤s4：对感兴趣区域骨架化(得到宽度为一个像素的骨架)，将骨组织之间的连接关系和拓扑进行压缩描述，且骨架化得到的图是后续图构建的基础和依据；

步骤s5：针对骨架化的图进行图构建，设计图割的代价函数，将分界曲面作为初始点，完成图割，识别出终板和椎体中心之间的相连骨小梁；

步骤s6：采用迭代切割线方法沿着相连骨小梁，搜索出最优位置和角度的线段，完全去除相连骨小梁，并进行三维重建获取终板模型。

较佳的，在本实施例中，步骤s1中，小鼠脊椎原始图像除了骨骼，还包括周围包裹的肌肉等软组织。根据软组织和骨骼像素对应的灰度值差异较大，设置图像灰度阈值为220(图像像素最大灰度值为255)，利用阈值法可以顺利得到图像中的骨骼像素。

在本实施例中，所述步骤s2中利用自动向下的方式完成图像感兴趣区域，脊椎由椎体和背部的脊椎突组成，其中椎体又包含两端上下终板以及中间的椎体中心。一组小鼠脊椎图像通常包含2到3个脊椎，两两脊椎通过背后的脊椎突上下相连。步骤s2具体包括以下步骤：

步骤s21：分离单个脊椎，由于图像的高分辨率，只要根据单个椎体内部骨骼的连通性，利用区域增长法可以顺利将脊椎与脊椎分离；

步骤s22：使用椎体两端投影法将脊椎背后的脊椎突去除，拟合出覆盖投影的圆柱形区域，其侧面的曲面形成分离曲面；

步骤s23：针对脊椎体，利用骨密度分析方法，获得终板所处的图像帧数范围。即相比椎体中心内部的规律孔洞的松质骨，终板内部为致密组织，两者存在较大的骨密度差别，因此统计椎体两端骨密度变化曲线，设置阈值获取密度最大且变化最大的区域作为终板出现的区域，即终板感兴趣区域。

在本实施例中，步骤s3具体包括以下步骤：

步骤s31：对终板感兴趣区域使用marchingcubes算法进行三维重建得到三角面片化的模型对象；

步骤s32：在模型表面构建调和场，调和场计算产生一系列候选等值曲线；其中，在调和场的构建时，分别设置调和场约束和权值函数，具体实施时约束选三个终板种子点设置为1，而三个椎体中心种子点设置为0，而在权值函数时添加凹度因子，具体来说，如果模型上的点是凹点，其凹度因子为0.001，否则，其凹度因子为0，这使每个具有高曲率凹度较大的点权值将大大降低，使调和场产生的一系列候选等值线在缝隙区域更稠密；

步骤s33：设计等值曲线分数计算函数(在选取最优等值线时，结合当前等值线所在局部区域的强度和形状变化，两者乘积结果作为每条等值线的分数)，根据分数结果选出分数最高的候选等值线作为最终表面分界曲线。

步骤s34：结合终板远椎体中心端的外边缘曲线和表面分界曲线(在得到的分界曲线的基础上，均匀选取终板图像中外边缘曲线中的10个点作为曲线拟合控制点)，采用三次b样条方法拟合出终板和椎体中间之间的初始分界曲面。

在本实施例中，步骤s34中，为了选出最贴合终板和椎体中心之间间隙区域的等值线，在候选等值线分数函数score(i)中引入平均梯度因子grd(i)和形状因子shape(i)设计，其中形状因子shape(i)以第i等值线为中心，沿着跨等值线方向，逐面片向左右两边扩展m步，分别得到左边和右边第j个扩展环，对应的长度分别为lij(left)和lij(right)，而对应的形状因子的函数为lmax是所有候选等值线长度最大值。

在本实施例的步骤s34中，形状因子计算使用了高斯多尺度卷积，可以保证对选择不敏感，将m设置为5，将σ设置为2。

较佳的，在本实施中，步骤s4采用cardenes提出的骨架化方法完成图像的骨架化，获得一个像素宽度的骨架线状网络结构图。在实施骨架化过程中，采用前景8连通，背景4连通的方式，考虑到后续终板分割属于序列图像分割，因此骨架化实现采用二维方式，骨架化后的图像中通常存在些许细小不连通的骨架，使用保留最大连通域的方法去除，从而获得最终一个像素宽度的线性网状骨架。

在本实施例中，步骤s5具体包括以下步骤：

步骤s51：对骨架化的图像进行图构建，其中图的顶点代表骨架分叉结点，图中的两两结点之间的弧是反应骨小梁长度和方向的曲线路径(不是传统图中的直线)；

步骤s52：根据相连骨小梁和椎体中心内部的骨小梁的形态差别特点，即相连骨小梁更接近椎体主方向且长度较短，因此将弧的权值函数设置为骨小梁长度与骨小梁和椎体主方向夹角绝对值之和，从而识别的相连骨小梁。

较佳的，图构建实施时，图中顶点对应着骨架分叉结点，图中两结点点间的弧对应着相连骨组织的曲线路径，两顶点间弧长长度是沿着弧长路径统计弧上像素总数。图割方法中，选择终板上的点作为源点，选择椎体中心上的点作为汇点。在图割目标函数的数据项中，设置规律为：若图中顶点离初始分界曲线越远，则顶点标签可信度越大，若顶点离分界曲线越远，则标签可信度越低。而在图割的光滑项中，设置规律为：具有较短长度且椎体主方向接近的弧具有低权值，即弧权值为弧长度与弧和椎体主方向之间夹角绝对值的和。通过实施图割方法，便可以顺利识别出所有骨小梁中属于终板和椎体中心缝隙区域内的相连骨小梁。

在本实施例中，步骤s6中，所述迭代切割线方法具体为：切割线从相连小梁中点开始，朝终板匀速靠近，当切割线到达终板表面时，切割线覆盖的骨小梁距离会急剧增加，当增加速度超过阈值20时，就代表切割直线接触到了终板表面，然后迭代调整切割线角度，将和调整差距最小的方向作为切割线最终方向。

在本实施例中，步骤s6具体包括以下步骤：

步骤s61：在单个相连小梁上构建一根切割线，切割线初始位置位于相连小梁骨架的中点，角度和相连小梁垂直，并记录下切割线与相连小梁区域相交形成的初始端点距离d1；

步骤s62：切割线沿着骨架超终板方向匀速移动，每次平移i，获得相应的端点距离di，且该距离会随着平移发生变化；当切割线到达终板表面时端点距离会急剧增加，当超过阈值20时，便停止前进，获得最终切割线位置；

步骤s63：切割线停止移动，重新调整其方向，每次旋转一定角度，获得新的端点距离，当新的端点距离和旋转前端点距离差最小时，当前角度即为终板表面相切角度；

步骤s64：按照步骤s63确定的最终的切割线位置和角度，终板上的相连小梁得以完全切除；

步骤s65：针对终板数据，使用marchingcube算法进行三维重建，得到终板的三维模型。

本领域内的技术人员应明白，本申请的实施例可提供为方法、系统、或计算机程序产品。因此，本申请可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、cd-rom、光学存储器等)上实施的计算机程序产品的形式。

本申请是参照根据本申请实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。

这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。

这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。

以上所述，仅是本发明的较佳实施例而已，并非是对本发明作其它形式的限制，任何熟悉本专业的技术人员可能利用上述揭示的技术内容加以变更或改型为等同变化的等效实施例。但是凡是未脱离本发明技术方案内容,依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与改型，仍属于本发明技术方案的保护范围。