CN108961273A - 一种从ct影像中分割肺动脉和肺静脉的方法及系统 - Google Patents

Abstract

本发明涉及医学图像处理技术领域。本发明提出从CT影像中分割肺动脉和肺静脉的方法，包括：步骤S10：根据从CT影像中获取的肺血管数据和肺血管中心线图，确定肺血管中心线中的至少一个黏连点；并从所述肺血管数据中剔除所述至少一个黏连点，得到以肺血管中心线为根节点的多个血管子树；步骤S20：将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树；步骤S30：分别获取所述第一类血管子树和第二类血管子树对应的血管体积的数值，并确定血管体积数值较大的血管子树为肺静脉血管，血管体积数值较小的血管子树为肺动脉血管。本发明提供的分割方法的精度能够满足临床需求，且运算时间短。

Description

一种从CT影像中分割肺动脉和肺静脉的方法及系统

技术领域

本发明涉及医学图像处理技术领域，特别是一种从CT影像中分割肺动脉和肺静脉的方法及系统。

背景技术

肺动脉和肺静脉的分割对检测肺部疾病有重要意义，CT影像中的肺动脉和肺静脉分割是非常困难的，它的难度主要在于以下三点：一是，由于二者都是管状结构以及相似的灰度值，在CT影像中很难区分出肺动脉和肺静脉；二是，二者经常是缠绕的，并且接触的部分，边界不清晰；三是，血管结构是高度可变的，因人而异。

目前虽然有一些分割肺动脉和肺静脉的方法，但其分割准确度不够，难以满足临床应用需求。

上述缺陷是本邻域技术人员期望克服的。

发明内容

(一)要解决的技术问题

针对现有技术中从CT影像中分割肺动脉和肺静脉的方法准确度不够，难以满足临床应用需求的问题，本发明提供了一种基于肺血管中心线从CT影像中分割肺动脉和肺静脉的方法及系统，以提高从医学影像中从CT影像中分割肺动脉和肺静脉的准确度及处理效率。。

(二)技术方案

第一方面，为解决以上问题，本发明提供的技术方案为：

一种从CT影像中分割肺动脉和肺静脉的方法，包括：

步骤S10：根据从CT影像中获取的肺血管数据和肺血管中心线图，确定肺血管中心线中的至少一个黏连点；

并从所述肺血管数据中剔除所述至少一个黏连点，得到以肺血管中心线为根节点的多个血管子树；

步骤S20：将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树；

步骤S30：分别获取所述第一类血管子树和第二类血管子树对应的血管体积的数值，并确定血管体积数值较大的血管子树为肺静脉血管，血管体积数值较小的血管子树为肺动脉血管。

作为本发明的进一步改进，在所述步骤S10中，

所述根据从CT影像中获取的肺血管数据和肺血管中心线图，确定肺血管中心线中的黏连点，包括：

根据从CT影像中获取的肺血管数据和肺血管中心线图，确定几何图；

确定所述几何图中的至少一个黏连点，在所述几何图中，每个所述黏连点的周围有4个像素点，且所述4个像素点中属于同一血管分支级别的2个像素点与所述黏连点的连线张成的角度近似为180°。

作为本发明的进一步改进，在所述步骤S10中，

所述从所述肺血管数据中剔除所述至少一个黏连点，得到以肺血管中心线为根节点的多个血管子树，包括：

从所述肺血管数据中剔除所述至少一个黏连点，及每个所述黏连点26邻域的像素点。

作为本发明的进一步改进，所述步骤S20中，

将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树，包括：

根据获取的预先设定的距离阈值D_max，确定两棵子树ψ_i、ψ_j中匹配的叶子节点数n_match；

根据第一公式确定子树ψ_i、ψ_j间的匹配强度IMS_i,j，所述第一公式为：

其中，μ_i(i＝1,2)为权值，两者的和为1；N_i表示子树ψ_i中的叶子节点数，N_j表示子树ψ_j中的节点数；

根据获取的预先设定的匹配强度阈值IMS_max，根据第二公式判断子树ψ_i、ψ_j是否属于同一类血管子树，所述第二公式为：

其中，f(i,j)取值为“yes”，则ψ_i、ψ_j属于同一类血管；f(i,j)取值为“no”，则ψ_i、ψ_j不属于同一类血管。

作为本发明的进一步改进，所述步骤S20中，

将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树，还包括：

S201：在所有的血管子树中寻找匹配强度最小的两个血管子树ψ₁和ψ₂；

将ψ₁和ψ₂分别放在集合S₁和S₂中，则集合S₁和S₂内的血管子树分别为第一类血管子树和第二类血管子树；

S202：从剩余的血管子树中寻找匹配强度最小的两个血管子树ψ₃和ψ₄，并将ψ₃和ψ₄，分别归入集合S₁或S₂：

S2021：分别确定ψ₃与S₁中所有子树的距离之和Sum₁，ψ₃与S₂中所有子树的距离之和Sum₂；

S2022：若ψ₃与S₁中所有子树的距离之和Sum₁小于与S₂中所有子树的距离之和Sum₂，则ψ₃属于S₁；否则，ψ₃属于S₂；

S2023：若ψ₃属于S₁，则ψ₄属于S₂；若ψ₃属于S₂，则ψ₄属于S₁；

S203：重复步骤S202，直到所有的血管子树被归入集合S₁或S₂。

作为本发明的进一步改进，所述根据获取的预先设定的距离阈值D_max，确定两棵子树ψ_i、ψ_j中匹配的叶子节点数n_match，包括：

根据第三公式确定子树ψ_i、ψ_j中匹配的叶子节点数n_match：

其中，g_current表示本次运算时，确定的子树ψ_i、ψ_j匹配的节点数；g_previous表示上一次运算时，确定的子树ψ_i、ψ_j匹配的节点数；n_i,s表示子树ψ_i中第s个节点，n_j,t表示子树ψ_j中第t个节点，d(n_i,s,n_j,t)表示节点n_i,s和n_j,t的欧氏距离，D_max为预先设定的距离阈值。

作为本发明的进一步改进，所述步骤S30中，

所述分别获取所述第一类血管子树和第二类血管子树对应的血管体积的数值，包括：

以剔除黏连点后的所述肺血管中心线上的像素点为种子点，在剔除所述黏连点后的所述肺血管数据中，分别确定所述第一类血管子树和第二类血管子树对应的血管数据；

分别获取所述第一类血管子树和第二类血管子树对应的血管数据的总体积V及对应的血管中心线的总长度L，根据第四公式确定所述第一类血管子树和第二类血管子树对应的血管体积的数值V_{normalization}，所述第四公式为：

第二方面，为解决以上问题，本发明提供的技术方案为：

一种从CT影像中分割肺动脉和肺静脉的系统，包括：

存储器、处理器、总线以及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现如第一方面中说明的方法的步骤。

第三方面，为解决以上问题，本发明提供的技术方案为：

一种计算机存储介质，其上存储有计算机程序，所述程序被处理器执行时实现第一方面中说明的方法的步骤。

(三)有益效果

与传统的采用灰度阈值分割法分割肺动脉和肺静脉不同，本发明提供的分割方法根据基于在肺血管中心线上，肺动脉和肺静脉的黏连点及与其周围相邻结点所组成的平面角度近似为平角这一特性，找到肺血管中心线中发生黏连的部分并剔除这些黏连的部分，得到多个血管子树；然后根据肺动静脉的血管子树中，叶子节点在CT影像上彼此靠近这一原理，设计了基于叶子节点欧氏距离的子树匹配算法将肺血管中心线分为两类；随后以肺血管中心线上的像素点为种子点，利用2D区域生长法，以肺血管中心线上的像素点为种子点，将两类血管分别生长出来；并根据肺动脉的体积比肺静脉的体积小这一医学常识确定两类血管分别所属的血管类别。

本发明提供的分割方法的精度能够满足临床需求，且运算时间短。

附图说明

图1为本发明实施例的肺动静脉分割流程示意图；

图2为肺血管图及标识有顶点、边界点及结点的肺中心血管图的示意图；

图3为肺血管中心线图转化为折线图的示意图；

图4为肺血管图及中心线及黏连点示意图；

图5为本发明实施例的黏连点检测与剔除流程示意图；

图6为两个血管子树及确定出的2个匹配节点的示意图；

图7为本发明实施例中基于叶子节点匹配将血管子树分为两类的流程示意图；

图8为本发明实施例中的3D肺血管图及标识顶点、边界点及结点后的3D肺血管图的示意图；

图9为本发明实施例中的黏连点检测结果示意图；

图10为本发明实施例中的划分得到的两类血管子树的示意图；

图11为本发明实施例中的肺动脉和肺静脉分割后的示意图；

图12为本发明实施例的分割方法与金标准中的分割方法的分割结果的对比示意图。

具体实施方式

为了更好的解释本发明，以便于理解，下面结合附图，通过具体实施方式，对本发明作详细描述。

肺血管分割对于肺部疾病检测有着重要的意义，利用计算机辅助诊断系统检测肺结节以及肺栓塞等疾病都涉及到肺血管的精确分割。肺血管中，肺动脉和肺静脉经常缠绕并接触。而在胸部CT图像中，肺动脉和肺静脉的灰度值是相近的，采用阈值分割法很难将二者分割开来。

由于肺血管自身结构复杂和CT影像上的噪声影响，用传统的影像分割方法很难达到理想的分割结果，本发明提供的分割方法根据基于在肺血管中心线上，肺动脉和肺静脉的黏连点及与其周围相邻结点所组成的平面角度近似为平角这一特性，找到肺血管中心线中发生黏连的部分并剔除(或去除)这些黏连的部分，得到多个血管子树；然后根据肺动静脉的血管子树中，叶子节点在CT影像上彼此靠近这一原理，设计了基于叶子节点欧氏距离的子树匹配算法将肺血管中心线分为两类；随后以肺血管中心线上的像素点为种子点，利用2D区域生长法，以肺血管中心线上的像素点为种子点，将两类血管分别生长出来；并根据肺动脉的体积比肺静脉的体积小这一医学常识确定两类血管分别所属的血管类别。

本发明提供的分割方法的精度能够满足临床需求，且运算时间短。

本发明提供的方法从CT影像中分割肺动脉和肺静脉经过以下三个步骤：首先，针对肺动脉和肺静脉发生黏连的部分，需要检测到黏连点并剔除；然后，基于叶子节点的欧氏距离及子树间的匹配度，将剔除黏连点后的多个血管子树分为两类；最后，基于归一化体积，确定每一类血管子树对应的肺部血管是属于肺动脉还是肺静脉。

具体地，如图1所示，本发明一个实施例的从CT影像中分割肺动脉和肺静脉的方法，包括：

步骤S10：根据从CT影像中获取的肺血管数据和肺血管中心线图，确定肺血管中心线中的至少一个黏连点；

并从所述肺血管数据中剔除所述至少一个黏连点，得到以肺血管中心线为根节点的多个血管子树；

步骤S20：将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树；

步骤S30：分别获取所述第一类血管子树和第二类血管子树对应的血管体积的数值，并确定血管体积数值较大的血管子树为肺静脉血管，血管体积数值较小的血管子树为肺动脉血管。

优选地，在所述步骤S10中，

所述根据从CT影像中获取的肺血管数据和肺血管中心线图，确定肺血管中心线中的黏连点，包括：

根据从CT影像中获取的肺血管数据和肺血管中心线图，确定几何图；

确定所述几何图中的至少一个黏连点，在所述几何图中，每个所述黏连点的周围有4个像素点，且所述4个像素点中属于同一血管分支级别的2个像素点与所述黏连点的连线张成的角度近似为180°。

需要说明的是，这里的几何图可以为图2所示的折线图。

优选地，在所述步骤S10中，

所述从所述肺血管数据中剔除所述至少一个黏连点，得到以肺血管中心线为根节点的多个血管子树，包括：

从所述肺血管数据中剔除所述至少一个黏连点，及每个所述黏连点26邻域的像素点。

优选地，所述步骤S20中，

将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树，包括：

根据获取的预先设定的距离阈值D_max，确定两棵子树ψ_i、ψ_j中匹配的叶子节点数n_match；

根据第一公式确定子树ψ_i、ψ_j间的匹配强度IMS_i,j，所述第一公式为：

其中，μ_i(i＝1,2)为权值，两者的和为1；N_i表示子树ψ_i中的叶子节点数，N_j表示子树ψ_j中的结点数；

根据获取的预先设定的匹配强度阈值IMS_max，根据第二公式判断子树ψ_i、ψ_j是否属于同一类血管子树，所述第二公式为：

其中，f(i,j)取值为“yes”，则ψ_i、ψ_j属于同一类血管；f(i,j)取值为“no”，则ψ_i、ψ_j不属于同一类血管。

至此，已经可以判断出任意的两个血管子树是否为同一类血管，但还需要进一步确定这两个血管子树分别属于第一类血管还是第二类血管。

具体编程实现时，所述步骤S20中，

将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树，可以包括：

S201：在所有的血管子树中寻找匹配强度最小的两个血管子树ψ₁和ψ₂；

将ψ₁和ψ₂分别放在集合S₁和S₂中，则集合S₁和S₂内的血管子树分别为第一类血管子树和第二类血管子树；

S202：从剩余的血管子树中寻找匹配强度最小的两个血管子树ψ₃和ψ₄，并将ψ₃和ψ₄，分别归入集合S₁或S₂：

S2021：分别确定ψ₃与S₁中所有子树的距离之和Sum₁，ψ₃与S₂中所有子树的距离之和Sum₂；

S2022：若ψ₃与S₁中所有子树的距离之和Sum₁小于与S₂中所有子树的距离之和Sum₂，则ψ₃属于S₁；否则，ψ₃属于S₂；

S2023：若ψ₃属于S₁，则ψ₄属于S₂；若ψ₃属于S₂，则ψ₄属于S₁；

S203：重复步骤S202，直到所有的血管子树被归入集合S₁或S₂。

需要说明的是，以上S2022和S2023中以ψ₃作为示例；应该理解为，以ψ₄作为示例时的步骤和逻辑是相似的，这里不再赘述。

优选地，所述根据获取的预先设定的距离阈值D_max，确定两棵子树ψ_i、ψ_j中匹配的叶子节点数n_match，包括：

根据第三公式确定子树ψ_i、ψ_j中匹配的叶子节点数n_match：

其中，g_current表示本次运算时，确定的子树ψ_i、ψ_j匹配的节点数；g_previous表示上一次运算时，确定的子树ψ_i、ψ_j匹配的节点数；n_i,s表示子树ψ_i中第s个节点，n_j,t表示子树ψ_j中第t个节点，d(n_i,s,n_j,t)表示节点n_i,s和n_j,t的欧氏距离，D_max为预先设定的距离阈值。

优选地，所述步骤S30中，

所述分别获取所述第一类血管子树和第二类血管子树对应的血管体积的数值，包括：

以剔除黏连点后的所述肺血管中心线上的像素点为种子点，在剔除所述黏连点后的所述肺血管数据中，分别确定所述第一类血管子树和第二类血管子树对应的血管数据；

分别获取所述第一类血管子树和第二类血管子树对应的血管数据的总体积V及对应的血管中心线的总长度L，根据第四公式确定所述第一类血管子树和第二类血管子树对应的血管体积的数值V_{normalization}，所述第四公式为：

如图1所示，本发明实施例的分割肺动脉和肺静脉的方法主要包括：黏连点判断与剔除、基于叶子节点欧氏距离的子树匹配和基于体积的肺动脉和肺静脉分割。

(1)、黏连点判断与剔除

在构成肺血管中心线的全部像素点中，可以划分为以下3类：

a、顶点(即端点)，顶点在其26领域内只有1个相邻点；

b、结点，结点在其26领域内有至少3个相邻点；

c、边界点，边界点在其26邻域中有且只有2个相邻点。

在编程实现时，可以将全部的肺血管中心线数据经过26邻域中相邻点的数目分别划分到顶点、结点和边界点3个数据集合内；也可以将全部的肺血管数据经过26邻域中相邻点的数目分别划分到顶点、结点和边界点3个数据集合内。

将图2(a)所示的肺血管的中心线上的顶点、边界点及结点分别用图形符号标识出；如图2(b)所示，其中，空心圆圈为顶点或结点，顶点的灰度值较大(具体地，在采用彩色配色方案时，为蓝色)，在其26领域内只有1个相邻点；结点的灰度值较小(具体地，在采用彩色配色方案时，为黄色)，在其26领域内具有3个以上的相邻点(也即该处的血管分支数大于等于3)；灰度值较大的线条(具体地，在采用彩色配色方案时，为黑色)由边界点组成，每个边界点在其26邻域中有且只有2个相邻点。

如图2所示，将相邻的边界点依次连接起来，可以得到一条直线；将多个不同延伸方向的直线利用节点连接起来，可以得到一组折线。

进一步地，如图3(a)所示，将肺血管中心线中每一个像素点用单元正方形来表示，每个单元正方形内的数字代表该像素点的相邻点个数，从而得到肺血管中心线图。

应该理解为，肺血管中心线图中的单元正方形的个数与肺血管中心线中的像素点个数是相同的；而单元正方形则可以由空心圆圈等其他便于显示的方式来替代。

更进一步地，为了更加直观地展示这三类像素点内在的几何关系，将肺血管中心线上的顶点和结点用实心圆点来表示，将多个相邻的边界点表示为两个顶点或结点之间的一条连线，即“边”edge，可以得到与肺血管中心线对应的折线图。

需要说明的是，以上肺血管中心线图和肺血管图均为三维3D图；在附图中，为了方便展示，可以选择在一个便于观察的剖面内的投影，即2D平面图形内进行直观、简洁地展示。

如图3(b)所示的肺血管中心线折线图中，具有7个实心圆点，其中，包括4个顶点和3个结点，n_i(i＝1,2,...，7)；具有6条边，也就是6条折线。

图4中示出了部分肺血管及这部分肺血管中心线对应的折线图。如图4所示，结点n₁₁周围有4个其他结点：n₄，n₈，n₁₅，n₁₆，其中，结点n₁₆和n₁₁之间的边与结点n₁₁和n₈之间的边组成的角∠n₁₆n₁₁n₈约等于180°，结点n₁₅和n₁₁之间的边与结点n₁₁和n₄之间的边所组成的角∠n₁₅n₁₁n₄也约等于180°。进一步结合图4中血管的形态，可以判断n₁₁是黏连点。

需要说明的是，“当肺血管发生黏连时，在黏连点周围将具有4个像素点(这4个像素点可以是顶点、结点或边界点，但不能全部是是顶点)，并且这4个像素点中具有相同数目相邻点(属于同一分支)的2个像素点与黏连点组成的角度在3D空间近似为180°”是满足一个统计学意义上的自然规律，并非是人为指定的规则。

从图4中可知，当肺血管发生黏连时，在黏连点周围将具有4个像素点(这4个像素点可以是顶点、结点或边界点，但不能全部是是顶点)，并且这4个像素点中具有相同数目相邻点(属于同一分支)的2个像素点与黏连点组成的角度在3D空间近似为180°。

如图5所示，在编程实现时，利用肺血管中心线的对应的折线图来检测黏连点的步骤如下：

首先确定相邻点数目为4的全部像素点；

然后逐一判断该像素点是否为黏连点：

A1、周围具有4个像素点；

A2、且这4个像素点中，有2个像素点的相邻点数目相同；

A3、且这2个像素点分别与该像素点的连线在空间组成的角度近似等于180°。

若以上A1-A3均成立，则确定该像素点为黏连点；否则，该像素点不是黏连点。

针对每个确定为黏连点的像素点，将其从肺血管中心线数据中或肺血管数据中剔除，从而以将肺血管中心线分割为多个互相不连通的血管子树。

剔除了黏连点之后的肺血管中心线图或肺血管图，从拓扑结构中看，为多个互相不连通的血管子树。

优选地，为了避免把黏连点剔除后，相邻的血管子树还可能是连通的，将黏连点以及黏连点26邻域的像素点全部从肺血管中心线上剔除。

记黏连点剔除后形成的血管子树集合为S＝{ψ_i}(i＝1,2,...)，其中ψ_i为第i个血管子树，每个血管子树中包括根节点和叶子节点。

确定血管子树的节点总数及叶子节点数以及位置与常规的树及子树的结构参数确定方式是相同的，且为本领域技术人员所公知，这里不再赘述。

(2)、利用叶子节点的欧氏距离将血管子树分类

鉴于肺血管中只有肺动脉和肺静脉这两类，因此依据血管子树中的叶子节点的距离可以将血管子树集合分为两类，其中一类为肺动脉血管子树，另一类为肺静脉血管子树。

这里的事实前提是：静脉血管末梢与静脉血管末梢的欧氏距离是小的；静脉血管末梢与动脉血管末梢的欧氏距离是小的。

1)、基于欧氏距离的叶子节点匹配

利用叶子节点欧氏距离将血管子树分类的步骤包括：通过两棵子树ψ_i、ψ_j中叶子节点的欧氏距离寻找子树ψ_i、ψ_j中匹配的叶子节点的数目，进而确定子树间匹配强度；进而利用该匹配强度来判断ψ_i、ψ_j是否属于同一类血管。

首先，根据公式(1.1)确定子树ψ_i、ψ_j中匹配的叶子节点数n_match：

其中，g_current表示本次运算时，确定的子树ψ_i、ψ_j匹配的节点数，g_previous表示上一次运算时，确定的子树ψ_i、ψ_j匹配的节点数；n_i,s表示子树ψ_i中第s个节点，n_j,t表示子树ψ_j中第t个节点，d(n_i,s,n_j,t)表示节点n_i,s和n_j,t的欧氏距离，D_max为预先设定的距离阈值。

也即，如果子树ψ_i中第s个节点与子树ψ_j中第t个节点的欧氏距离d(n_i,s,n_j,t)小于预先设定的距离阈值，则判断子树ψ_i中第s个节点与子树ψ_j中第t个节点匹配；如果子树ψ_i中第s个节点与子树ψ_j中第t个节点的欧氏距离d(n_i,s,n_j,t)等于或大于预先设定的距离阈值，则判断子树ψ_i中第s个节点与子树ψ_j中第t个节点不匹配。

在最后一次运算完成时确定的子树ψ_i、ψ_j匹配的节点数g_current即子树ψ_i、ψ_j中匹配的叶子节点数n_match。

具体编程实现时，将全部的血管子树两两组合分别运算，即可以得到任意两个子树ψ_i与子树ψ_j中匹配的节点数。

2)、基于匹配强度进行子树分类

进一步地，根据公式(1.2)确定子树间匹配强度IMS_i,j：

其中，N_i表示子树ψ_i中的叶子节点数，N_j表示子树ψ_j中的叶子节点数，n_match表示子树ψ_i、ψ_j中匹配的叶子节点数，μ_i(i＝1,2)分别表示子树ψ_i和子树ψ_j的权值，两个权值之和恒等于1。

具体编程实现时，将全部的血管子树两两组合，可以得到任意两个子树ψ_i与子树ψ_j的子树间匹配强度IMS_i,j。

接下来，根据公式(1.3)判断子树ψ_i、ψ_j是否属于同一类血管：

其中，IMS_max为预先设置的匹配强度阈值，IMS_i,j为子树ψ_i和ψ_j的间匹配强度。

即，若IMS_i,j大于匹配强度阈值IMS_max，则f(i,j)取值为“yes”，则ψ_i、ψ_j属于同一类血管；若IMS_i,j不大于IMS_max，则f(i,j)取值为“no”，则ψ_i、ψ_j不属于同一类血管。

图6中示出了两个血管子树的示意图及确定出的2个匹配节点。如图6所示，两棵子树间相匹配的叶子节点有2个，假设IMS_max＝0.70，μ₁＝0.5，μ₂＝0.5，则由公式(1.2)可得IMS_i,j为0.5；由公式(1.3)可得IMS_i,j小于IMS_max，所以图6中所示两棵子树属于同一类血管。

基于叶子节点欧氏距离及子树匹配度将全部的血管子树分为两类的流程示意图如图7所示。

具体实施时，可以采用以下算法：

S100：首先在所有的血管子树中寻找匹配强度最小两个血管子树ψ₁和ψ₂；然后将这两个血管子树分别放在集合S₁和S₂中，则集合S₁和S₂内的血管可以认为是互斥的两类，也即第一类和第二类；

S200：接着从剩余的血管子树集合中继续寻找匹配强度最小的两个血管子树ψ_i和ψ_j，并将ψ_i和ψ_j分别归入集合S₁或S₂：

S201：分别确定ψ_i与S₁中所有子树的距离之和Sum₁，ψ_i与S₂中所有子树的距离之和Sum₂；

S202：利用属于同一类的子树的叶子节点距离相近这一准则，若ψ_i与S₁中所有子树的距离之和Sum₁小于与S₂中所有子树的距离之和Sum₂，则ψ_i属于S₁；否则，ψ_i属于S₂；

因为ψ_i和ψ_j属于两个不同的类别；若ψ_i所属的类别确定了，则ψ_j所属的类别也就确定了。

相应地，若ψ_i属于S₁，则ψ_j属于S₂；若ψ_i属于S₂，则ψ_j属于S₁。

以上在不知道每个血管子树所属的血管类别的情况下，利用子树叶子节点的欧氏距离将肺血管树分为两类。

3)、基于归一化体积进行肺动脉和肺静脉分类

至此，已经可以将肺血管数据或肺血管中心线数据分为两个子树集合S₁和S₂，但是还无法确定这两个子集合中哪一个是肺动脉，哪一个是肺静脉。

结合以下医学事实：肺动脉的血压高于肺静脉的血压；而肺血管中的血流量是一定的，即肺动脉中的血流量和肺静脉中的血流量是相等的，则可以确定肺动脉的总体积比肺静脉的总体积要小。

因此，可以根据这两类血管的总体积来对他们进行分类。

优选地，在确定血管总体积时，为了不受血管子树长度变化影响，根据公式(1.4)归一化每一类血管的体积：

V_{normalization}为归一化后的体积，V为某一类的血管子树的总体积，L为该类血管子树的血管中心线的总长度。

具体实施时，采用以下算法确定肺血管中每一类血管子树的总体积：

首先，以剔除了黏连点后的中心线上的像素点为种子点，以这些种子点在已经分割出肺血管的CT图像上进行2D区域生长，从而避免因动静脉黏连而发生的生长错误；

接着，统计每一类的血管子树的总体积，以及中心线的总长度，从而计算出归一化后的总体积；

最后，依据肺动脉的归一化的总体积小于肺静脉的归一化的总体积这一原理，确定这两类血管是肺动脉还是肺静脉。

以下使用多组数据对本发明实施例所提出的肺动脉和肺静脉分割方法进行验证。这些数据包括已经分割出的肺血管数据和已提取出的肺血管中心线数据。

首先进行黏连点判断与剔除。因为左右两个肺在进行动静脉分割时是相互独立的，为了不相互影响，将中心线转化为几何图时，将一个肺的中心线分别单独地进行转化。

将一组512*512*386的肺血管中心线数据转化为几何图的结果如图8所示。图8(a)为左肺中心线转几何图；图8(b)为右肺中心线转几何图。其中，转化后的几何图中分别示出了在26邻域中只有两个相邻像素点的边界点；只有1个相邻像素点的顶点；在26邻域中相邻像素点多于两个的结点。

接下来利用几何图中节点的属性及几何关系进行黏连点的检测与剔除。在检测黏连点时也是分别针对左右侧的每一个肺中的中心线单独地进行，检测得到黏连点如图9所示，其中，图9(a)为左肺，图9(b)为右肺。

图9中用到的这组CT数据为512*512*403，其中，浅色线条为肺血管中心线，深色的局部点为黏连点。

检测出黏连点后，将黏连点及其周围的26邻域的像素点剔除掉，这时，根据肺血管中心线图及数据就可以得到互不联通的多个血管子树。

采用前述的方法将这些子树分为两类，以某组数据为例，两类血管子树的分类结果如图10所示。

为最终确定血管所属类别，需要沿着肺血管中心线通过2D区域生长法先将肺血管生长出来，然后比较两类的体积大小以最终确认血管子树所属的类别。

以某三组数据为例，统计了两类肺血管的体积大小以及中心线的大小，如表1所示。

表1某三组数据的左右肺血管归一化体积对照表

根据归一化的肺血管体积来判断两类血管的所属类别，其中，数值较大的归一化血管体积属于肺静脉，较小的归一化血管体积属于肺动脉。其中，某组数据的最终的分类结果如图11所示，其中，图11中较深的颜色(红色)示出的属于肺动脉，较深的颜色(蓝色)示出的的属于肺静脉。

图11只是定性地显示了肺动脉和肺静脉的分割结果。为了能够定量地分析肺动脉和肺静脉分割的准确性，将本发明实施例分割的结果与专业人员手工标记的标准子树结果进行比较，以得到本发明实施例的算法的分割准确率。以某组数据为例，两种方法分割后的对比结果如图12所示，其中，多个椭圆框示出的区域显示了本发明实施例算法与金标准的之间的差异。

采用本发明实施例算法的分割结果与金标准结果取交集，然后剔除金标准结果得到分割准确率，以此来评价本发明实施例分割算法的准确性，并对某5组数据分别进行了统计，结果如表2所示。

表2某五组数据的肺动脉和肺静脉分割准确率

由表2可知，无论是左肺的肺动脉和肺静脉的平均分割准确率还是右肺的肺动脉和肺静脉的平均分割准确率都在80％左右，达到了预期的分割要求。

相应地，统计上述5组数据在进行肺动静脉分割时的时间如表3所示，由表3可知在1分钟内就可以自动地将肺动脉和肺静脉分割完毕，这种量级的计算时间可以满足临床需求。

表3对5组CT影像数据分割肺动脉和肺静脉所需的运算时间

综上，本发明实施例的分割方法根据在肺血管中心线上，肺动脉和肺静脉的黏连点及与其周围相邻结点所组成的平面角度近似为平角这一特性，找到肺血管中心线中发生黏连的部分并剔除这些黏连的部分，得到多个血管子树；然后根据肺动静脉的血管子树中，叶子节点在CT影像上彼此靠近这一原理，设计了基于叶子节点欧氏距离的子树匹配算法将肺血管中心线分为两类；随后以肺血管中心线上的像素点为种子点，利用2D区域生长法，以肺血管中心线上的像素点为种子点，将两类血管分别生长出来；并根据肺动脉的体积比肺静脉的体积小这一医学常识确定两类血管分别所属的血管类别。

经过与金标准对比，本发明实施例的分割方法的分割精度能够满足临床需求，且运算时间短。

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

Claims (9)

1.一种从CT影像中分割肺动脉和肺静脉的方法，其特征在于，包括：

步骤S10：根据从CT影像中获取的肺血管数据和肺血管中心线图，确定肺血管中心线中的至少一个黏连点；

并从所述肺血管数据中剔除所述至少一个黏连点，得到以肺血管中心线为根节点的多个血管子树；

步骤S20：将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树；

步骤S30：分别获取所述第一类血管子树和第二类血管子树对应的血管体积的数值，并确定血管体积数值较大的血管子树为肺静脉血管，血管体积数值较小的血管子树为肺动脉血管。

2.根据权利要求1所述的方法，其特征在于，在所述步骤S10中，

所述根据从CT影像中获取的肺血管数据和肺血管中心线图，确定肺血管中心线中的黏连点，包括：

根据从CT影像中获取的肺血管数据和肺血管中心线图，确定几何图；

确定所述几何图中的至少一个黏连点，在所述几何图中，每个所述黏连点的周围有4个像素点，且所述4个像素点中属于同一血管分支级别的2个像素点与所述黏连点的连线张成的角度近似为180°。

3.根据权利要求1所述的方法，其特征在于，在所述步骤S10中，

所述从所述肺血管数据中剔除所述至少一个黏连点，得到以肺血管中心线为根节点的多个血管子树，包括：

从所述肺血管数据中剔除所述至少一个黏连点，及每个所述黏连点26邻域的像素点。

4.根据权利要求1所述的方法，其特征在于，所述步骤S20中，

将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树，包括：

根据获取的预先设定的距离阈值D_max，确定两棵子树ψ_i、ψ_j中匹配的叶子节点数n_match；

根据第一公式确定子树ψ_i、ψ_j间的匹配强度IMS_i,j，所述第一公式为：

其中，μ_i(i＝1,2)为权值，两者的和为1；N_i表示子树ψ_i中的叶子节点数，N_j表示子树ψ_j中的节点数；

根据获取的预先设定的匹配强度阈值IMS_max，根据第二公式判断子树ψ_i、ψ_j是否属于同一类血管子树，所述第二公式为：

其中，f(i,j)取值为“yes”，则ψ_i、ψ_j属于同一类血管；f(i,j)取值为“no”，则ψ_i、ψ_j不属于同一类血管。

5.根据权利要求4所述的方法，其特征在于，所述步骤S20中，

将所述多个血管子树按照匹配度划分为两类，分别为第一类血管子树和第二类血管子树，还包括：

S201：在所有的血管子树中寻找匹配强度最小的两个血管子树ψ₁和ψ₂；

将ψ₁和ψ₂分别放在集合S₁和S₂中，则集合S₁和S₂内的血管子树分别为第一类血管子树和第二类血管子树；

S202：从剩余的血管子树中寻找匹配强度最小的两个血管子树ψ₃和ψ₄，并将ψ₃和ψ₄，分别归入集合S₁或S₂：

S2021：分别确定ψ₃与S₁中所有子树的距离之和Sum₁，ψ₃与S₂中所有子树的距离之和Sum₂；

S2022：若ψ₃与S₁中所有子树的距离之和Sum₁小于与S₂中所有子树的距离之和Sum₂，则ψ₃属于S₁；否则，ψ₃属于S₂；

S2023：若ψ₃属于S₁，则ψ₄属于S₂；若_ψ3属于S₂，则ψ₄属于S₁；

S203：重复步骤S202，直到所有的血管子树被归入集合S₁或S₂。

6.根据权利要求4所述的方法，其特征在于，

所述根据获取的预先设定的距离阈值D_max，确定两棵子树ψ_i、ψ_j中匹配的叶子节点数n_match，包括：

根据第三公式确定子树ψ_i、ψ_j中匹配的叶子节点数n_match：

其中，g_current表示本次运算时，确定的子树ψ_i、ψ_j匹配的节点数；g_previous表示上一次运算时，确定的子树ψ_i、ψ_j匹配的节点数；n_i,s表示子树ψ_i中第s个节点，n_j,t表示子树ψ_j中第t个节点，d(n_i,s,n_j,t)表示节点n_i,s和n_j,t的欧氏距离，D_max为预先设定的距离阈值。

7.根据权利要求1所述的方法，其特征在于，所述步骤S30中，

所述分别获取所述第一类血管子树和第二类血管子树对应的血管体积的数值，包括：

以剔除黏连点后的所述肺血管中心线上的像素点为种子点，在剔除所述黏连点后的所述肺血管数据中，分别确定所述第一类血管子树和第二类血管子树对应的血管数据；

分别获取所述第一类血管子树和第二类血管子树对应的血管数据的总体积V及对应的血管中心线的总长度L，根据第四公式确定所述第一类血管子树和第二类血管子树对应的血管体积的数值V_{normalization}，所述第四公式为：

8.一种从CT影像中分割肺动脉和肺静脉的系统，其特征在于，包括：

存储器、处理器、总线以及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现如权利要求1-7任意一项的方法的步骤。

9.一种计算机存储介质，其上存储有计算机程序，其特征在于，所述程序被处理器执行时实现如权利要求1-7任意一项的方法的步骤。