CN115661467A - 脑血管图像分割方法、装置、电子设备及存储介质 - Google Patents

Abstract

本发明涉及图像分割技术领域，提供一种脑血管图像分割方法、装置、电子设备及存储介质，其中方法包括：获取待分割的三维图像；对三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；基于海森矩阵血管滤波方法，得到血管增强图像；基于初始血管分割结果，确定分割阈值；基于分割阈值和血管增强图像所对应的血管响应函数，确定第一能量函数；基于初始血管分割结果，确定第二能量函数；基于第一能量函数和第二能量函数确定各体素的类先验概率；基于隶属度矩阵，以及类先验概率，确定血管类别的后验概率和背景类别的后验概率，并基于后验概率，确定血管分割结果。本发明提供的方法、装置、电子设备及存储介质，提高了分割质量。

Description

脑血管图像分割方法、装置、电子设备及存储介质

技术领域

本发明涉及图像分割技术领域，尤其涉及一种脑血管图像分割方法、装置、电子设备及存储介质。

背景技术

TOF-MRA（Time of flight-MR angiography）是一种非侵入性的脑血管***成像技术，拥有良好的空间分辨率，相比于计算机断层扫描血管造影（Computed TomographicAngiography，CTA）图像，TOF-MRA受骨骼影响小，已广泛应用于临床。从TOF-MRA图像可以准确检测脑血管，以支持血管相关疾病的早期诊断和治疗，对于许多临床应用至关重要。

传统脑血管结构提取与识别主要依赖于经验丰富的专家手动操作，具有主观性强、复现性差、效率低及复杂性高等缺点。而基于深度学习的脑血管图像分割算法需要基于大量高质量的标注数据学习才能获得精准血管分割模型。然而，对于专家来说，脑血管的手动标注也是一项具有挑战性的任务，获取良好的脑血管标注需要耗费专家大量的时间与精力，成本高昂，且存在明显的主观差异，不符合实际临床应用中的迫切需求。

发明内容

本发明提供一种脑血管图像分割方法、装置、电子设备及存储介质，用以解决现有技术中脑血管结构提取与识别主要依赖于经验丰富的专家手动操作，具有主观性强、复现性差、效率低及复杂性高的缺陷。

本发明提供一种脑血管图像分割方法，包括：

获取待分割的三维图像；

对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；

基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像；

基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值；

基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数；

基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数；

基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率；

基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

根据本发明提供的一种脑血管图像分割方法，所述对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果，包括：

基于K-Means算法，对基准图像的灰度直方图进行聚类，得到初始聚类中心；

基于所述初始聚类中心，以及所述三维图像的灰度直方图，对所述三维图像进行聚类处理，得到所述初始血管分割结果。

根据本发明提供的一种脑血管图像分割方法，所述基于所述初始聚类中心，以及所述三维图像的灰度直方图，对所述三维图像进行聚类处理，得到所述初始血管分割结果，包括：

基于所述初始聚类中心，确定当前聚类中心；

基于所述当前聚类中心以及所述三维图像的灰度直方图，确定所述三维图像中全体体素的隶属度矩阵；

基于所述灰度直方图和所述全体体素的隶属度矩阵，对所述当前聚类中心进行更新，得到更新聚类中心；

基于所述当前聚类中心和所述更新聚类中心之间的差异程度，确定目标函数；

若所述目标函数大于等于预设阈值，则将所述更新聚类中心作为当前聚类中心，返回对所述当前聚类中心进行更新，否则将所述更新聚类中心下全体体素的隶属度矩阵，确定所述初始血管分割结果。

根据本发明提供的一种脑血管图像分割方法，所述基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，包括：

基于如下公式，构建马尔科夫随机场：

其中，

表示马尔科夫随机场，

表示第一能量函数，

表示第二能量函数，

表示当前处理体素，

表示当前处理体素的邻域体 素，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体素的初始分割结 果，

和

表示权重系数，其中α+β= 1。

根据本发明提供的一种脑血管图像分割方法，所述基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数，包括：

基于如下公式，确定所述三维图像中各体素的第一能量函数：

其中，

表示第一能量函数，

表示当前处理体素，

表示当前处理体 素的邻域体素，

表示当前处理体素的初始分割结果，

表示血管响应函数，L _v 表示血管 点，L _b 表示非血管点，T表示血管增强图像的分割阈值。

根据本发明提供的一种脑血管图像分割方法，所述基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数，包括：

基于如下公式，确定所述三维图像中各体素的第二能量函数：

其中，

表示第二能量函数，

表示当前处理体素，

表示当前处理 体素的邻域体素，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体素的 初始分割结果。

根据本发明提供的一种脑血管图像分割方法，所述三维图像的灰度直方图的确定步骤包括：

获取所述三维图像的原始灰度直方图；

基于基准灰度直方图，对所述原始灰度直方图进行组映射，得到所述三维图像的灰度直方图；

所述基准灰度直方图是对基准图像的灰度直方图进行融合平均得到的。

本发明还提供一种脑血管图像分割装置，包括：

获取单元，用于获取待分割的三维图像；

聚类单元，用于对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；

血管增强单元，用于基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像；

确定分割阈值单元，用于基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值；

确定第一能量函数单元，用于基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数；

确定第二能量函数单元，用于基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数；

确定类先验概率单元，用于基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率；

血管分割单元，用于基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

本发明还提供一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现如上述任一种所述脑血管图像分割方法。

本发明还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现如上述任一种所述脑血管图像分割方法。

本发明还提供一种计算机程序产品，包括计算机程序，所述计算机程序被处理器执行时实现如上述任一种所述脑血管图像分割方法。

本发明提供的脑血管图像分割方法、装置、电子设备及存储介质，对三维图像的灰度直方图进行聚类处理，得到初始血管分割结果，相较于传统的聚类算法基于原始图像像素或体素值进行计算，提高了聚类处理的效率，节省了大量的时间成本；并且，基于初始血管分割结果中全体体素的隶属度矩阵，以及类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于血管类别的后验概率、背景类别的后验概率以及最大后验分类准则，确定血管分割结果，提高了分割质量，也可以避免耗费专家大量的时间与精力，降低了时间成本，且不存在明显的主观差异，符合实际临床应用中的迫切需求。

附图说明

为了更清楚地说明本发明或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作一简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明提供的脑血管图像分割方法的流程示意图之一；

图2是本发明提供的确定三维图像的灰度直方图的流程示意图；

图3是本发明提供的脑血管图像分割方法的流程示意图之二；

图4是本发明提供的脑血管图像分割装置的结构示意图；

图5是本发明提供的电子设备的结构示意图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合本发明中的附图，对本发明中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明的说明书和权利要求书中的术语“第一”、“第二”等是用于区别类似的对象，而不用于描述特定的顺序或先后次序。应该理解这样使用的数据在适当情况下可以互换，以便本申请的实施例能够以除了在这里图示或描述的那些以外的顺序实施，且“第一”、“第二”等所区分的对象通常为一类。

相关技术中，TOF-MRA图像中的脑血管的准确分割提取主要受以下几个方面的制约：脑血管形状的复杂性和可变性、血管灰度的动态变化性、部分血管结构的细小性、血管空间分布复杂性、图像噪声和局部容积效应明显等。这些特性直接导致对TOF-MRA图像中的脑血管的分割困难。

目前用于脑血管分割的方法包括：1）基于区域相似性模型的分割方法；2）基于活动轮廓模型的分割方法；3）基于统计模型的分割方法；4）基于深度学习的分割方法。

1）基于区域相似性模型的分割方法：该类方法将图像分割成具有一定相似性的体素区域，进一步可以通过两类方法去分割，即阈值分割法和区域增长法。Wang 等人利用Otsu阈值法从 TOF-MRA 图像中提取包含灰白质、眼睛以及脑血管的前景区域和包含脑脊液、骨骼和空气的背景区域，再对前景区域进行阈值分割提取脑血管。区域增长方法从血管内的手工种子点出发，逐渐地把领域范围内满足生长条件的体素点加入集合中以实现分割，与单纯的阈值分割法相比，能够使得血管分割具有较好的连通性并保留一定的几何拓扑结构。

2）基于活动轮廓模型的分割方法：该类方法使用连续曲线来表达目标边缘，并定义能量泛函使得其自变量包括边缘曲线从而将分割过程转变为求解能量泛函最小值的过程。根据轮廓的表示方法，该方法进一步可划分为参数化活动轮廓模型以及几何活动轮廓模型。其中，参数化活动轮廓模型中典型的模型可以是Snake模型，几何活动轮廓模型中典型的方法可以是水平集方法。Snake模型将轮廓曲线S以参数形式表示，内部力Fin与S形状有关，被用来控制S的弹性形变，保证曲线的平滑连续性。外部力Fout用于促进S趋近血管边缘。Lee等人提出了结合卡尔曼滤波的参数化活动轮廓模型分割方法。他们把体素灰度和梯度作为外部力Fout，对第一层CTA图像切片进行分割，然后使用卡尔曼滤波跟踪其他切片中的血管位置从而完成血管分割。水平集方法将轮廓曲线设置为一个高维水平集函数的零水平集，然后通过对水平集函数进行微分从而提取零水平集来得到物体的轮廓。Wang等人提出了一种基于自适应阈值的水平集方法用于脑血管分割。他们设计了包含两个区域项、一个边界项和一个惩罚项的新型能量方程。其中第一个区域项利用最大密度投影获得全局阈值，第二个区域项使用动态灰度阈值保证细小血管的提取，边界项用来控制轮廓曲线趋于血管边界，惩罚项用来避免重新初始化水平集函数。

3）基于统计模型的分割方法：该类方法主要通过聚类分析实现图像分割。Wilson等人假设TOF-MRA图像直方图由两个高斯分布和一个均匀分布加权组合而成的有限混合模型（FMM）表示，利用期望最大算法（Expectation-Maximization algorithm，EM）对有限混合模型中的参数进行估计从而实现脑血管分割。Zhou等人则假设TOF-MRA图像直方图由一个瑞利分布和三个高斯分布加权组合而成，通过增大模型的复杂度提高直方图曲线拟合精度，并进一步通过马尔科夫随机场（Markov Random Field，MRF）考虑血管体素之间的邻域关系，从而优化分割结果。Zhang等人提出了一种设备无关的FMM-MRF脑血管分割模型。他们通过基于血管知识的EM算法对三个高斯分布组成的FMM进行参数估计，并进一步结合多尺度血管滤波增强图像以及血管体素邻域信息，提出了双能量约束的马尔科夫随机场用于分割结果的优化。

4）基于深度学***衡问题。Chen等人结合UNet和Transformer，提出了一种基于生成一致性的半监督模型用于脑血管分割。他们通过图像重建的一致性来提高模型对纹理的表征能力，同时使用语义融合的Transformer网络来增强模型的语义特征表达能力。

由于容易受图像噪声影响以及不同个体的TOF-MRA图像血管灰度差异较大等因素，基于手工设置阈值的全局或局部阈值法或生长种子点的区域生长法难以获得鲁棒准确的分割结果。此外，基于深度学习的脑血管图像分割算法需要基于大量高质量的标注数据学习才能获得精准血管分割模型。然而，即使对于专家来说，脑血管的手动标注也是一项具有挑战性的任务，获取良好的脑血管标注需要耗费专家大量的时间与精力，成本高昂，且存在明显的主观差异，不符合实际临床应用中的迫切需求。

相比之下，基于统计模型的分割方法鲁棒性较强，自动化程度较高，提取的血管网络较为完整，相较于区域相似性模型的分割方法或基于深度学习的分割方法具有明显优势。然而，现有的基于统计模型的分割方法仍然未能充分利用脑血管特征，抗噪能力低，在图像对比度低，血管灰度动态变化大等情况下分割出的脑血管会出现较多的断裂情况或假阳性情况，迫切需要更为精确、鲁棒、高度自动化的脑血管分割方法。

基于上述问题，本发明提供一种脑血管图像分割方法，图1是本发明提供的脑血管图像分割方法的流程示意图之一，如图1所示，该方法包括：

步骤110，获取待分割的三维图像。

具体地，可以获取待分割的三维图像，此处的待分割的三维图像是指后续用于图像分割的图像，此处的三维图像可以是磁共振血管造影（TOF-MRA）图像，也可以是CTA图像（Computed Tomography Angiography，计算机断层摄影血管造影），本发明实施例对此不作具体限定。

步骤120，对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果。

具体地，在获取到待分割的三维图像后，可以对三维图像进行去颅骨、去噪声和偏置场校正的预处理操作。其中，可以使用磁共振影像分割工具BET（Brain ExtractionTool）或DeepBET对三维图像进行去颅骨操作，得到脑组织区域，提高血管信噪比，减轻计算负担。然后，可以使用开源软件ANTs（Advanced Net Tool Suite）对去颅骨后的血管图像进行去噪声和偏置场校正。

其中，可以使用空间自适应非局部均值滤波器（Non-Local Means，NL-means）对图像进行去噪声，去除的噪声可以包括高斯噪声和莱斯噪声。由于磁场的变化，图像中相同组织在视觉上可能会变得更亮或更暗，这些变化的映射称为偏置场。由于这些信号灰度的变化并不是解剖学上的差异，偏置场将给分类器或算法的性能带来一定的影响。因此，可以使用ANTs的N4BiasFieldCorrection算法对去噪声后的三维图像进行偏置场校正，有利于更加精准地提取脑血管。

在对三维图像进行预处理后，可以得到三维图像的灰度直方图，然后，可以对三维图像的灰度直方图进行直方图规定化处理后再进行聚类处理，从而得到初始血管分割结果。此处的聚类处理可以使用模糊C均值聚类算法（Fuzzy C-Means，FCM）。

此处的直方图规定化是指在获取到原始灰度直方图后，基于基准灰度直方图，对原始灰度直方图进行组映射（Group Mapping Law，GML），得到三维图像的灰度直方图的过程。此处的基准灰度直方图是对基准图像的灰度直方图进行融合平均得到的，此处的基准图像是具有相同的灰度直方图峰值最多的那些图像。

可以理解的是，在对三维图像的灰度直方图进行直方图规定化处理后，可以得到变换后的三维图像及其三维图像的灰度直方图。

此处，可以使用K-Means算法作用于基准灰度直方图，得到初始聚类中心。

在得到三维图像的灰度直方图和初始聚类中心后，可以使用模糊C均值聚类算法 进行聚类处理，得到初始血管分割结果。例如，可以对三维图像的灰度直方图进行聚类处 理，得到全体体素的隶属度矩阵

，对

进行去模糊化，即根据最大隶属度原则赋予 每个灰度值最大隶属度所在的类为其聚类结果。将所有灰度值所属类映射回原始图像空间 当中，获得三维图像中每个体素所属类别，取四个类别中的最后两类所包含的体素作为血 管体素，从而得到初始血管分割结果。

传统的聚类算法是基于原始图像像素或体素值进行计算，增加了大量的时间成本以及增加了大量的内存。尽管对于二维图像这些缺陷可以忽略，但对于大型三维图像这些缺陷无法忽略，且会明显影响算法的执行效率，因此，本发明实施例对三维图像的灰度直方图进行聚类处理，可以提高聚类处理的效率，节约了大量的时间成本。

步骤130，基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像。

具体地，当三维图像是TOF-MRA图像时，由于TOF-MRA是非对比剂增强血管成像技术，基于血液流入增强效应，使得血管呈现高亮信号。然而，层流引起的血流速度变化会导致血液中自旋信号灰度不均，导致细小血管灰度低，粗大血管灰度高的情况。为了降低血管造影图像中这些难以避免的灰度变化，抑制非血管结构或组织和图像噪声，而血管增强滤波可以作为一个重要的预处理技术，其中，可以基于海森矩阵血管滤波方法，对三维图像进行血管增强，得到血管增强图像，具体步骤如下：

第一步，可以计算图像在x点处的海森矩阵

：

其中，

表示x处的灰度值，σ表示血管尺度，G(σ, x)表示高斯函数。

第二步，对海森矩阵进行特征值分解：

，并对特 征值的绝对值进行排序。

第三步，对于管状结构的血管，海森矩阵特征值存在关系

。 针对正常管状结构的血管，定义血管响应函数F如下：

其中，

，

，

为海森矩阵特征值。

第四步，为了改善低对比度下的血管显影，在尺度σ下对

重新定义如下：

其中τ是0到1之间的阈值，控制对噪声目标的响应。

此外，为增强椭圆血管结构，修改血管响应函数F如下：

由此，可以基于血管响应函数，对三维图像进行血管增强，得到血管增强图像。

步骤140，基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值。

具体地，在得到初始血管分割结果后，可以基于初始血管分割结果，确定血管比 例，此处的血管比例可以是血管所占三维图像中颅内空间比例，可以用

表示。

在确定血管比例后，可以基于血管比例，从血管增强图像中确定分割阈值，分割阈值T满足如下公式：

（1）

其中，N表示体素总数，T表示分割阈值，F表示血管增强图像所对应的血管响应函 数，公式（1）表示血管增强图像中满足F≥T的体素总数占比为

。

步骤150，基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数。

具体地，在确定分割阈值后，可以基于分割阈值和血管增强图像所对应的血管响应函数，确定三维图像中各体素的第一能量函数，此处的第一能量函数可以是基于血管增强图像的单基团能量函数。

步骤160，基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数。

具体地，在得到初始血管分割结果后，可以基于初始血管分割结果，确定三维图像中各体素的第二能量函数，此处的第二能量函数可以是基于初始血管分割结果的双基团六邻域能量函数。

步骤170，基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率。

具体地，初始血管分割结果中包含了许多噪声点、孤立血管点或者血管片段，需要结合血管几何知识以及空间邻域信息对初始血管分割结果进一步优化。

因此，在得到血管增强图像、初始血管分割结果、第一能量函数和第二能量函数后，可以基于第一能量函数和第二能量函数构建马尔科夫随机场（Markov Random Field，MRF），并基于马尔科夫随机场，得到各体素的类先验概率。此处的马尔可夫随机场可用于解决离散的优化问题，其每个标签是独立的，非连续的。此处的类先验概率是指根据以往经验和分析得到的概率，如全概率公式，它往往作为“由因求果”问题中的“因”出现的概率。

各体素的类先验概率的公式如下：

其中，

表示类先验概率，

表示尔科夫随机场，

是血管点L _v 和非血管点L _b 的总和，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体 素的初始分割结果。

步骤180基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

具体地，在得到类先验概率后，可以基于初始血管分割结果中全体体素的隶属度矩阵，以及类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率。即，三维图像在进行去颅骨、去噪声和偏置场校正的预处理操作后，再对三维图像的灰度直方图进行聚类处理，可以将三维图像中的目标分为血管类别和背景类别两类。此处的迭代条件模式是一种确定性算法，用于获得马尔可夫随机场联合概率的局部最大值的配置，迭代条件模式通过迭代最大化以其余变量为条件的每个变量的概率来实现。

例如，可以将初始血管分割结果中全体体素的隶属度矩阵看作类条件概率，则血管类别隶属度和背景类别隶属度的公式如下：

其中，

表示血管类别隶属度，

表示背景类别隶属度，

表示当前处理体 素

的第

类隶属度，

表示初始分割结果中第k类体素占所有体素的比例。

因此，血管类别的后验概率可以表示为：

背景类别的后验概率可以表示为：

其中，

表示整体能量函数，

表示血管类别隶属度，

表示正比符 号。

在确定血管类别的后验概率和背景类别的后验概率后，可以基于血管类别的后验概率、背景类别的后验概率以及最大后验分类准则，确定血管分割结果。

公式如下：

即，可以以血管类别的后验概率大于背景类别的后验概率作为最大后验分类准则，确定血管分割结果。

由此，可以从三维图像中，将血管清晰地分割出来，方便后续的研究以及临床应用。

本发明实施例提供的方法，对三维图像的灰度直方图进行聚类处理，得到初始血管分割结果，相较于传统的聚类算法基于原始图像像素或体素值进行计算，提高了聚类处理的效率，节省了大量的时间成本；并且，基于初始血管分割结果中全体体素的隶属度矩阵，以及类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于血管类别的后验概率、背景类别的后验概率以及最大后验分类准则，确定血管分割结果，提高了分割质量，也可以避免耗费专家大量的时间与精力，降低了时间成本，且不存在明显的主观差异，符合实际临床应用中的迫切需求。

基于上述实施例，步骤120包括：

步骤121，基于K-Means算法，对基准图像的灰度直方图进行聚类，得到初始聚类中心；

步骤122，基于所述初始聚类中心，以及所述三维图像的灰度直方图，对所述三维图像进行聚类处理，得到所述初始血管分割结果。

具体地，可以基于K-Means算法，对基准图像的灰度直方图进行聚类，得到初始聚类中心。此处的基准图像是具有相同的灰度直方图峰值最多的那些图像。

在得到初始聚类中心后，可以基于初始聚类中心，以及三维图像的灰度直方图，对三维图像进行聚类处理，得到初始血管分割结果。

基于上述实施例，步骤122包括：

基于所述初始聚类中心，确定当前聚类中心；

基于所述当前聚类中心以及所述三维图像的灰度直方图，确定所述三维图像中全体体素的隶属度矩阵；

基于所述灰度直方图和所述全体体素的隶属度矩阵，对所述当前聚类中心进行更新，得到更新聚类中心；

基于所述当前聚类中心和所述更新聚类中心之间的差异程度，确定目标函数；

若所述目标函数大于等于预设阈值，则将所述更新聚类中心作为当前聚类中心，返回对所述当前聚类中心进行更新，否则将所述更新聚类中心下全体体素的隶属度矩阵，确定所述初始血管分割结果。

具体地，可以基于初始聚类中心，确定当前聚类中心，例如，可以以初始聚类中心 作为当前聚类中心，当前聚类中心可以用

表示。

然后，可以基于当前聚类中心以及三维图像的灰度直方图，确定三维图像中全体 体素的隶属度矩阵

：

其中，

表示体素灰度值，

，M表示图像的最大灰度值，l表示类别，c表 示类别总数，m表示模糊系数，

表示

的范数，

表示

的范数，Z 表示聚类中心，l表示遍历所有类别，

表示第s类到第k类中心的距离。

再基于灰度直方图和全体体素的隶属度矩阵，对当前聚类中心进行更新，得到更新聚类中心。

其中，

为更新聚类中心，

为隶属度矩阵，s表示灰度，

表示灰度为s的体 素频率，

表示体素灰度值，M表示图像的最大灰度值。

进一步，可以基于当前聚类中心和更新聚类中心之间的差异程度，确定目标函数

：

其中，

为更新聚类中心，

为当前聚类中心。

最后，若目标函数

大于等于预设阈值

，则将更新聚类中心作为当前聚类中心， 返回对当前聚类中心进行更新，否则，即

，将更新聚类中心下各体素的隶属度矩阵，确 定初始血管分割结果。

此外，还可以通过测试发现聚类中心为4的时候图像分割的效果最佳，观察图像分割结果发现前两类分别集中于脑脊液和灰白质，后两类集中于细长的、灰度值较低的动、静脉血管以及粗大的、灰度值较高的动脉血管，因此，还可以预设所有三维图像的聚类中心为4。

本发明实施例提供的方法，基于初始聚类中心，以及三维图像的灰度直方图，对三维图像进行聚类处理，得到初始血管分割结果，由此得到的初始血管分割结果的分割效果更好，有利于后续确定血管分割结果。

基于上述实施例，步骤170包括：

基于如下公式，构建马尔科夫随机场：

其中，

表示马尔科夫随机场，

表示第一能量函数，

表示第二能量函数，

表示当前处理体素，

表示当前处理体素的邻域体素，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体素的初始分割结果，

和

表示权重系数，其中α+β= 1。

基于上述实施例，步骤150包括：

基于如下公式，确定所述三维图像中各体素的第一能量函数：

其中，

表示第一能量函数，

表示当前处理体素，

表示当前处理体素 的邻域体素，

表示当前处理体素的初始分割结果，

表示血管响应函数，L _v 表示血管 点，L _b 表示非血管点，T表示血管增强图像的分割阈值。

基于上述实施例，步骤160包括：

基于如下公式，确定所述三维图像中各体素的第二能量函数：

其中，

表示第二能量函数，

表示当前处理体素，

表示当前处理体 素的邻域体素，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体素的初 始分割结果。

基于上述实施例，图2是本发明提供的确定三维图像的灰度直方图的流程示意图，图2中，左侧为三维图像的原始灰度直方图，中间为基准灰度直方图，右侧为三维图像的灰度直方图。如图2所示，所述三维图像的灰度直方图的确定步骤包括：

步骤210，获取所述三维图像的原始灰度直方图；

步骤220，基于基准灰度直方图，对所述原始灰度直方图进行组映射，得到所述三维图像的灰度直方图；

所述基准灰度直方图是对基准图像的灰度直方图进行融合平均得到的。

具体地，可以获取三维图像的原始灰度直方图，此处的三维图像的原始灰度直方图是对三维图像进行去颅骨、去噪声和偏置场校正的预处理操作得到的。

为了减小从不同站点，不同扫描仪获取的三维图像之间的差异，以提高算法的泛化能力以及自动化程度，本发明实施例基于基准灰度直方图，对原始灰度直方图进行组映射（Group Mapping Law，GML），作为额外的预处理步骤。

在获取到原始灰度直方图后，可以基于基准灰度直方图，对原始灰度直方图进行组映射，得到三维图像的灰度直方图，此处的基准灰度直方图是对基准图像的灰度直方图进行融合平均得到的，此处的基准图像是具有相同的灰度直方图峰值最多的那些图像。

首先，从所有待分割的三维图像中选取具有相同的灰度直方图峰值最多的那些图像，即基准图像，将基准图像的直方图进行融合平均，并进行腐蚀、闭运算等形态学操作，得到接近正态分布的直方图作为基准灰度直方图。此处的形态学操作可以包括腐蚀操作和闭运算操作。

然后，获取三维图像的原始灰度直方图，基于基准灰度直方图，对原始灰度直方图进行组映射，得到三维图像的灰度直方图。为增加组映射后的血管体素灰度值连续性，在对原始灰度直方图进行组映射之前采用3×3×3的高斯核对图像进行滤波处理。具体步骤如下：

第一步，可以获取三维图像的原始灰度直方图

，以及累计直方图

，其中，

是三维图像的概率密度函数。

第二步，获取基准灰度直方图

，其中，

是基准灰度直 方图对应的概率密度函数。

第三步，令

，使用GML（Group Mapping Law）组映射规则，从基准灰度 直方图映射到原始灰度直方图，为每一个基准灰度直方图的概率密度找到一组映射，依据 该映射求得

，即由原始灰度直方图的灰度值r得到的基准灰度直方图对应的 灰度值z。

基于上述任一实施例，图3是本发明提供的一种脑血管图像分割方法的流程示意图之二，如图3所示，该方法包括：

第一步，可以获取待分割的三维图像。

第二步，可以对三维图像进行去颅骨、去噪声、偏置场校正等预处理操作，得到原始灰度直方图。

第三步，可以选取基准图像，对基准图像的灰度直方图进行融合平均得到基准灰度直方图，然后，可以基于基准灰度直方图，对原始灰度直方图进行组映射，得到三维图像的灰度直方图。

第四步，可以对三维图像的灰度直方图进行聚类处理，得到初始血管分割结果。

第五步，可以基于海森矩阵血管滤波方法，对三维图像进行血管增强，得到血管增强图像。

第六步，可以基于血管增强图像和初始血管分割结果，构建马尔科夫随机场，并基于马尔科夫随机场，得到类先验概率。

第七步，基于初始血管分割结果中全体体素的隶属度矩阵，以及类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于血管类别的后验概率、背景类别的后验概率以及最大后验分类准则，确定血管分割结果。

下面对本发明提供的脑血管图像分割装置进行描述，下文描述的脑血管图像分割装置与上文描述的脑血管图像分割方法可相互对应参照。

基于上述任一实施例，本发明提供一种脑血管图像分割装置，图4是本发明提供的脑血管图像分割装置的结构示意图，如图4所示，该装置包括：

获取单元410，用于获取待分割的三维图像；

聚类单元420，用于对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；

血管增强单元430，用于基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像；

确定分割阈值单元440，用于基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值；

确定第一能量函数单元450，用于基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数；

确定第二能量函数单元460，用于基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数；

确定类先验概率单元470，用于基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率；

血管分割单元480，用于基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

本发明实施例提供的装置，对三维图像的灰度直方图进行聚类处理，得到初始血管分割结果，相较于传统的聚类算法基于原始图像像素或体素值进行计算，提高了聚类处理的效率，节省了大量的时间成本；并且，基于初始血管分割结果中全体体素的隶属度矩阵，以及类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于血管类别的后验概率、背景类别的后验概率以及最大后验分类准则，确定血管分割结果，提高了分割质量，也可以避免耗费专家大量的时间与精力，降低了时间成本，且不存在明显的主观差异，符合实际临床应用中的迫切需求。

基于上述任一实施例，聚类单元具体包括：

确定初始聚类中心单元，用于基于K-Means算法，对基准图像的灰度直方图进行聚类，得到初始聚类中心；

确定初始血管分割结果单元，用于基于所述初始聚类中心，以及所述三维图像的灰度直方图，对所述三维图像进行聚类处理，得到所述初始血管分割结果。

基于上述任一实施例，确定初始血管分割结果单元具体用于：

基于所述初始聚类中心，确定当前聚类中心；

基于所述当前聚类中心以及所述三维图像的灰度直方图，确定所述三维图像中各体素的隶属度矩阵；

基于所述灰度直方图和所述各体素的隶属度矩阵，对所述当前聚类中心进行更新，得到更新聚类中心；

基于所述当前聚类中心和所述更新聚类中心之间的差异程度，确定目标函数；

若所述目标函数大于等于预设阈值，则将所述更新聚类中心作为当前聚类中心，返回对所述当前聚类中心进行更新，否则将所述更新聚类中心下各体素的隶属度矩阵，确定所述初始血管分割结果。

基于上述任一实施例，确定类先验概率单元具体包括：

基于如下公式，构建马尔科夫随机场：

其中，

表示马尔科夫随机场，

表示第一能量函数，

表示第二能量函数，

表示当前处理体素，

表示当前处理体素的邻域体素，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体素的初始分割结果，

和

表示权重系数，其中α+β= 1。

基于上述任一实施例，确定第一能量函数单元具体用于：

基于如下公式，确定所述三维图像中各体素的第一能量函数：

其中，

表示第一能量函数，

表示当前处理体素，

表示当前处理体 素的邻域体素，

表示当前处理体素的初始分割结果，

表示血管响应函数，L _v 表示血管 点，L _b 表示非血管点，T表示血管增强图像的分割阈值。

基于上述任一实施例，确定第二能量函数单元具体用于：

基于如下公式，确定所述三维图像中各体素的第二能量函数：

其中，

表示第二能量函数，

表示当前处理体素，

表示当前处理体 素的邻域体素，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体素的初 始分割结果。

基于上述任一实施例，所述三维图像的灰度直方图的确定步骤包括：

获取所述三维图像的原始灰度直方图；

基于基准灰度直方图，对所述原始灰度直方图进行组映射，得到所述三维图像的灰度直方图；

所述基准灰度直方图是对基准图像的灰度直方图进行融合平均得到的。

图5示例了一种电子设备的实体结构示意图，如图5所示，该电子设备可以包括：处理器（processor）510、通信接口（Communications Interface）520、存储器（memory）530和通信总线540，其中，处理器510，通信接口520，存储器530通过通信总线540完成相互间的通信。处理器510可以调用存储器530中的逻辑指令，以执行脑血管图像分割方法，该方法包括：获取待分割的三维图像；对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像；基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值；基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数；基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数；基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率；基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

此外，上述的存储器530中的逻辑指令可以通过软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备（可以是个人计算机，服务器，或者网络设备等）执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：U盘、移动硬盘、只读存储器（ROM，Read-Only Memory）、随机存取存储器（RAM，Random Access Memory）、磁碟或者光盘等各种可以存储程序代码的介质。

另一方面，本发明还提供一种计算机程序产品，所述计算机程序产品包括计算机程序，计算机程序可存储在非暂态计算机可读存储介质上，所述计算机程序被处理器执行时，计算机能够执行上述各方法所提供的脑血管图像分割方法，该方法包括：获取待分割的三维图像；对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像；基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值；基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数；基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数；基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率；基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

又一方面，本发明还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现以执行上述各方法提供的脑血管图像分割方法，该方法包括：获取待分割的三维图像；对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像；基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值；基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数；基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数；基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率；基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

以上所描述的装置实施例仅仅是示意性的，其中所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。本领域普通技术人员在不付出创造性的劳动的情况下，即可以理解并实施。

通过以上的实施方式的描述，本领域的技术人员可以清楚地了解到各实施方式可借助软件加必需的通用硬件平台的方式来实现，当然也可以通过硬件。基于这样的理解，上述技术方案本质上或者说对现有技术做出贡献的部分可以以软件产品的形式体现出来，该计算机软件产品可以存储在计算机可读存储介质中，如ROM/RAM、磁碟、光盘等，包括若干指令用以使得一台计算机设备（可以是个人计算机，服务器，或者网络设备等）执行各个实施例或者实施例的某些部分所述的方法。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (10)

1.一种脑血管图像分割方法，其特征在于，包括：

获取待分割的三维图像；

对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；

基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像；

基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值；

基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数；

基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数；

基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率；

基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

2.根据权利要求1所述的脑血管图像分割方法，其特征在于，所述对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果，包括：

基于K-Means算法，对基准图像的灰度直方图进行聚类，得到初始聚类中心；

基于所述初始聚类中心，以及所述三维图像的灰度直方图，对所述三维图像进行聚类处理，得到所述初始血管分割结果。

3.根据权利要求2所述的脑血管图像分割方法，其特征在于，所述基于所述初始聚类中心，以及所述三维图像的灰度直方图，对所述三维图像进行聚类处理，得到所述初始血管分割结果，包括：

基于所述初始聚类中心，确定当前聚类中心；

基于所述当前聚类中心以及所述三维图像的灰度直方图，确定所述三维图像中全体体素的隶属度矩阵；

基于所述灰度直方图和所述全体体素的隶属度矩阵，对所述当前聚类中心进行更新，得到更新聚类中心；

基于所述当前聚类中心和所述更新聚类中心之间的差异程度，确定目标函数；

若所述目标函数大于等于预设阈值，则将所述更新聚类中心作为当前聚类中心，返回对所述当前聚类中心进行更新，否则将所述更新聚类中心下全体体素的隶属度矩阵，确定所述初始血管分割结果。

4.根据权利要求1所述的脑血管图像分割方法，其特征在于，所述基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，包括：

基于如下公式，构建马尔科夫随机场：

其中，

表示马尔科夫随机场，

表示第一能量函数，

表示第二能量函数，

表示当前处理体素，

表示当前处理体素的邻域体 素，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体素的初始分割结 果，

和

表示权重系数，其中α+β= 1。

5.根据权利要求1所述的脑血管图像分割方法，其特征在于，所述基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数，包括：

基于如下公式，确定所述三维图像中各体素的第一能量函数：

其中，

表示第一能量函数，

表示当前处理体素，

表示当前处理体素的 邻域体素，

表示当前处理体素的初始分割结果，

表示血管响应函数，L _v 表示血管点，L _b 表示非血管点，T表示血管增强图像的分割阈值。

6.根据权利要求1所述的脑血管图像分割方法，其特征在于，所述基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数，包括：

基于如下公式，确定所述三维图像中各体素的第二能量函数：

其中，

表示第二能量函数，

表示当前处理体素，

表示当前处理体素 的邻域体素，

表示当前处理体素的初始分割结果，

表示当前体素的邻域体素的初 始分割结果。

7.根据权利要求1至6中任一项所述的脑血管图像分割方法，其特征在于，所述三维图像的灰度直方图的确定步骤包括：

获取所述三维图像的原始灰度直方图；

基于基准灰度直方图，对所述原始灰度直方图进行组映射，得到所述三维图像的灰度直方图；

所述基准灰度直方图是对基准图像的灰度直方图进行融合平均得到的。

8.一种脑血管图像分割装置，其特征在于，包括：

获取单元，用于获取待分割的三维图像；

聚类单元，用于对所述三维图像的灰度直方图进行聚类处理，得到初始血管分割结果；

血管增强单元，用于基于海森矩阵血管滤波方法，对所述三维图像进行血管增强，得到血管增强图像；

确定分割阈值单元，用于基于所述初始血管分割结果，确定血管比例，并基于所述血管比例，确定所述血管增强图像的分割阈值；

确定第一能量函数单元，用于基于所述分割阈值和所述血管增强图像所对应的血管响应函数，确定所述三维图像中各体素的第一能量函数；

确定第二能量函数单元，用于基于所述初始血管分割结果，确定所述三维图像中各体素的第二能量函数；

确定类先验概率单元，用于基于所述第一能量函数和所述第二能量函数构建马尔科夫随机场，并基于所述马尔科夫随机场，确定各体素的类先验概率；

血管分割单元，用于基于所述初始血管分割结果中全体体素的隶属度矩阵，以及所述类先验概率，在迭代条件模式下，确定血管类别的后验概率和背景类别的后验概率，并基于所述血管类别的后验概率和所述背景类别的后验概率，确定血管分割结果。

9.一种电子设备，包括存储器、处理器及存储在所述存储器上并可在所述处理器上运行的计算机程序，其特征在于，所述处理器执行所述程序时实现如权利要求1至7任一项所述脑血管图像分割方法。

10.一种非暂态计算机可读存储介质，其上存储有计算机程序，其特征在于，所述计算机程序被处理器执行时实现如权利要求1至7任一项所述脑血管图像分割方法。

Images