CN100466974C - 使用n-MIP光线投射技术的自动冠状动脉分离 - Google Patents

Abstract

提出一种检测心脏的冠状动脉和其他外周血管的新方法。在通过分割方法、例如图形理论分割方法找到心肌位置后，该方法使用双轴椭球对心脏建模。对于椭球的每一个点，计算与表面垂直的强度集。该集然后被滤波以检测心血管结构。最后，使用血管跟踪技术，血管中心线点被检测并被连接在一起以形成完整的冠状动脉树。

Description

使用n-MIP光线投射技术的自动冠状动脉分离

相关申请的对照

本申请要求2003年11月17日申请的美国临时申请60/520789的优先权，该临时申请全部包括在这里作为参考。

技术领域

本发明通常涉及医学成像领域，尤其涉及使用n-MIP光线投射技术的自动冠状动脉分离(coronary isolation)。

背景技术

在美国心脏病是主要的死亡原因。一种心脏病是冠心病。冠心病是向心脏供血的冠状动脉狭窄或堵塞。这种狭窄或堵塞潜在减缓了血流量，这就减少了送往心脏的氧和营养物质。这可能导致胸痛(即，心绞痛)，或者，如果一条或更多的冠状动脉完全堵塞，则可能导致心脏病发作。其他由冠状动脉病变导致的问题包括充血性心力衰竭(下文中的“CHF”)。

由于动脉粥样斑在动脉血管壁上的沉积，冠状动脉变得狭窄或堵塞。动脉粥样斑可以包括血流中的过量胆固醇和其他物质，例如引起炎症的细胞，蛋白质和钙。动脉粥样斑在动脉中的沉积被称为动脉硬化症，也被称为“动脉硬化”。

冠状动脉的动脉粥样斑、动脉瘤和异常结构的早期检测是诊断和治疗冠状动脉病变的重要因素。超声心动描记术是检测冠状动脉的动脉粥样斑、动脉瘤和异常结构的主要工具。超声心动描记术使用高频声波显示心脏的组织、结构和功能。超声心动描记术的缺点包括低分辨率和噪声。通常，不能用超声心动描记术显示冠状动脉。

最近，计算机断层摄影术(“CT”)已在各种医学诊断领域、例如癌症检测中，得到了广泛地应用。CT提供亚毫米空间分辨率和低噪声级。但是，使用CT进行心脏外周血管检测的技术数量极其有限。

磁共振(“MR”)提供了和由CT提供的相同的分辨率和低噪声级，此外MR完全是非侵入性的。使用MR的人不需要象CT所要求的那样注射任何造影剂。另外，CT使用X射线，长期照射X线会有潜在的健康风险，而MR使用人体磁场。但是，起搏器、动脉瘤夹和其他植入物的使用者不能使用MR。

发明内容

在本发明的一个方面，提供一种自动获得图像中的心脏外周血管中心线的方法。该方法包括对图像中的心脏进行分割以产生三维(“3D”)分割的心肌；由三维分割的心肌计算3D距离映像；向3D距离映像中心脏表面上的每一个体素投射法向最大强度投影光线；对投影光线滤波以获得血管候选点；和精选(refining)候选点以获得中心线。

在本发明的另一个方面，提供一种自动获得图像中的心脏外周血管中心线的方法。该方法包括对图像中的心脏进行分割以产生三维(“3D”)分割的心肌；由三维分割的心肌计算3D距离映像；使用3D距离映像，用回转椭球对心壁建模；对回转椭球上的每一个体素投射法向最大强度投影光线；对投影光线滤波以获得血管候选点；和精选血管候选点以获得中心线。

在本发明的另外一个方面，提供一种其中已经存储指令的机器可读介质，通过处理器执行介质中已存储的指令来执行自动获得图像中心脏外周血管中心线的方法。该方法包括对图像中的心脏进行分割以产生三维(“3D”)分割的心肌；由三维分割的心肌计算3D距离映像；向3D距离映像中心脏表面上的每一个体素投射法向最大强度投影光线；对投影光线滤波以获得血管候选点；和精选血管候选点以获得中心线。

在本发明的另外一个方面，提供一种其中已经存储指令的机器可读介质，通过处理器执行介质中已存储指令来执行自动获得图像中的心脏外周血管中心线的方法。该方法包括对图像中的心脏进行分割以产生三维(“3D”)分割的心肌；由三维分割的心肌计算3D距离映像；使用3D距离映像，用回转椭球对心壁建模；对回转椭球上的每一个体素投射法向最大强度投影光线；对投影光线滤波以获得血管候选点；和精选血管候选点以获得中心线。

附图说明

通过参考下列结合附图的描述可以理解本发明，其中相同的参考数字表示相同的元件，其中：

根据本发明的一个典型实施方案，图1描述了一种从3D空间中的点到2D空间中的点的典型映射；

根据本发明的一个典型实施方案，图2描述了不同患者心脏外周血管的典型展开的n-MIP视图；

根据本发明的一个实施方案，图3描述了在心脏表面两个不同位置上的典型轮廓曲线；

根据本发明的一个典型实施方案，图4描述了典型血管候选点；

根据本发明的一个典型实施方案，图5描述了在连接在一起的血管中心线上的典型点；以及

根据本发明的一个典型实施方案，图6描述了心脏外周血管分割方法的流程图。

具体实施方式

下面描述本发明的说明性实施方案。为了清楚起见，在本说明书中没有描述实际实施中的所有特征。当然应当理解，在任何这种实际实施方案的开发中，要达到开发人员的特定目的，例如遵守系统相关和商业相关的限制，必须做出许多特定实施决定，这些特定目的随着实施不同而发生变化。此外，应当理解这种开发计划可能是复杂和耗时的，但对于那些得益于本发明的本领域普通技术人员来说不过是常规任务。

尽管对本发明容易进行各种修改和替换，它的特定实施方案已经以实例的方式显示在附图中，并在这里进行详细描述。但是，应当理解这里关于特定实施方案的描述并不意味着将本发明限制在已公开的特殊形式上，相反，本发明涵盖了落入由后附权利要求限定的本发明精神和范围内的所有修改、等价和替代形式。

应当理解这里描述的系统和方法可以以各种形式的硬件、软件、固件、专用处理机或它们的组合实现。尤其是，至少本发明的一部分作为包括程序指令的应用程序优选实现，这些程序指令被确实地包含在一个或多个程序存储设备里(例如，硬盘，磁性软盘，RAM，ROM，CD ROM等)并且可由任何包括适当结构的设备或机器、例如一种具有处理器、内存和输入/输出接口的通用数字计算机执行。还应当理解，由于附图中描述的一些系统组成部件和处理步骤优选地以软件实现，所以系统模块(或方法步骤的逻辑流程)之间的连接可能不同，这取决于本发明的编程方式。根据这里给出的指导，相关领域的普通技术人员将能够设计出本发明的这些及类似的实施方式。

本发明提出一种心脏外周血管重建方法，该方法假定存在代表心肌的分割体。这种方法基于以下事实：血管通常和心脏表面平行。因此，对心壁分割可以提供关于血管的重要信息。这里描述的方法能够用于MR和CT模态图像。

本发明考虑一种自动冠状血管分割的新颖的、多阶段方法。

(1)分割心脏并获取心壁外形。

(2)通过简单的几何体、例如回转椭球对心脏表面建模。

(3)光线从回转椭球表面上的每一个点投射，形成三维(“3D”)光线集，并记录分割强度。3D光线集用作预分割工具。和光线相交的每一条血管在光线轮廓曲线(profile curve)上产生一个强度峰。这是一种检测属于血管的体素的简单技术。然后，在这个体素空间里的高强度管状结构能够用于检测血管。

(4)使用血管跟踪技术和最小生成树，建立完整的血管树。

心脏分割和距离映像计算

本发明考虑一种从带有形状约束的图形优化技术衍生的分割方法。这种图形分割的原理是，根据最小权重的分割使图形中定义的能量函数最小化。上述能量被记作两项之和：对分割映像施加平滑约束的E_smooth(f)，以及衡量标号f如何适应数据的E_data(f)：

E(f)＝E_smooth(f)+E_data(f)，       (1)

$E_{\mathrm{smooth}}\left(f\right)=\underset{p,q\∈\mathrm{neighbors}}{\mathrm{\Σ}}{V}_{p,q}(f\left(p\right),f\left(q\right)),---\left(2\right)$

$E_{\mathrm{data}}\left(f\right)=\underset{p\∈P}{\mathrm{\Σ}}{D}_p\left(f\left(p\right)\right)---\left(3\right)$

方程(2)中的V_p，q是一对相邻像素{p，q}之间的相互作用函数。方程(3)中的D_p衡量标号f和像素p强度的接近程度。通常，这种方法提供一种针对二进制值f(p)情况的整体最优解决方案。

上面提供的心脏分割方法仅是示例。应当理解，如本领域技术人员所设想的，可以使用本领域进行分离的其它方法，例如基于模型的分割和基于水平集方法的分割方法。

上述分割产生3D掩模(即，像素标记的对象和背景)。在恢复外周血管期间，从3D掩模表面的距离映像能够提供重要的约束。例如，当外周血管和心脏表面保持平行时，在距离映像中它们的距离(即，血管到心脏表面的距离)变化平稳。

通过两次解析3D掩模计算距离映像，所述两次解析中一次在一个方向，一次在另一方向。距离映像基本上是与每个点和参照物(这里是心脏表面)之间的距离相关的体积。和本领域技术人员公知的一样，距离映像中的每一个体素在边缘处(对象—背景)被3D斜面掩模M滤波。如下面将要详细描述的，通过一个更简单的几何对象，例如椭球，已滤波的体素被用于对心脏建模，以便使用制图学方法容易地平滑心脏表面。

用回转椭球对外形建模

距离映像用于依据椭球或双轴回转椭球对心壁建模。尽管可以考虑更加精确的模型，例如三轴椭球，但主要的缺陷是不存在3D—2D投影的简单数学方法。

双轴椭球投影是制图学中的公知技术，但在心脏的医疗成像中却不是公知技术。双轴椭球减少了使用球体时可能发生的变形。对于半轴长度为a和b的双轴椭球，表面方程为

$\frac{x^2}{a^2}+\frac{y^2}{a^2}+\frac{x^2}{b^2}=1---\left(4\right)$

或者，用参数形式：x＝acos(λ)sin(φ)，y＝asin(λ)sin(φ)，z＝bcos(φ)，其中λ∈[0，2π]且φ∈[0，π]。和通常在制图学中使用的类似，λ和φ被分别称为经度和纬度。椭球中心被计算作为位于距离映像等值面0(心脏外形)上的所有点的重心G。长轴Z是向量GM，其中M是等值面上使得长度‖GM‖最大的点。同样，短轴X是向量GN，其中N是等值面上使得长度‖GN‖最小的点。推导出轴Y具有直(direct)正交基，B＝(G，X，Y，Z)。

应当理解建模质量不依赖于分割质量和抗噪声性能，这使这里描述的分割方法和上述分割方法的选择无关。此外，和其他器官不同，心脏天然具有椭球外形。光线从椭球表面投射以计算心脏表面的二维(“2D”)视图。

光线投射和N-MIP投影

然后，从椭球上投射光线，并采集心壁内外预定方向上和范围内的体素强度。通过下面的简单变换，计算用于投射光线的椭球上的结点分布。现在参考图1，对于3D空间内(椭球上)的每一个点P(λ，φ)105，根据下式计算2D点X(μ，v)110(平面上)

如果光线直接从距离映像自身被投射，光线将不够均匀而不能被利用。实际上，由于心脏表面不规则，这种光线会遗漏许多血管。另一种解决方法可以根据双轴椭球的横轴墨卡托投影(transverse Mercator projection)。这种方法的缺点是该体积必须是回转椭球，这就意味着短半轴位于极面上，长半轴位于赤道面上。

现在参考图2，用反色显示了不同患者心脏外周血管展开的n-MIP(即，法向最大强度投影)视图的实例。因为为了使视图看上去与血管造影片相似，强度已经被反转，所以冠状动脉同心脏其他部分比起来显得更暗。2D视图中的每一个像素I(x，y)是每一个光线强度I_n的标准化之和

$I(x,y)=\frac{1}{N}\underset{0<n<N+1}{\mathrm{\&Sigma;}}{I}_n---\left(5\right)$

轮廓匹配技术用于检测光线是否和血管结构相交。现在参考图3，光线轮廓曲线起始于心壁，心壁被表示为大的、平坦、高强度区域。图3显示了心脏表面上两个不同位置处的轮廓曲线。心肌区域用阴影表示。血管结构显示为距离心包一定距离的高强度、有限直径的峰。

下面，对这些峰进行处理，并根据下面的标准确定它们是否是血管点的合适候选点：

1.它关系到特定的强度属性(强度值和峰形)

2.它到心壁的距离在一定范围内

被加到局部最大值(脊检测)的这些强度峰允许光线轮廓曲线上的血管的公平检测。为了检测这些峰，本发明使用具有下述核：[-12-1]的拉普拉斯算子的零交点。

现在参考图4，一旦这个过程完成，在光线上检测到的峰就和3D点联系在一起。如图4所示，这些3D点可能不是准确地位于血管中心，但是紧邻或位于血管腔内。这种情况源于光线投射效果，因为光线在整个空间是均匀分布的。下一步将这些候选点集中，滤除噪声，并跟踪检测到的血管。

精选血管候选点以找出中心线

血管候选体素能够假定为紧邻或位于血管腔内。然而，如图4所示，体素仍不在血管中心线上。此外，同样如图4所示，峰检测可能对噪声敏感。因此，需要额外的处理将这些点集中和消除噪声。精选操作包括在每个3D点对应的局部强度空间近邻中对它进行特征值分析。有

$A_{\mathrm{ij}}=\frac{{\&PartialD;}^{2}I}{\&PartialD;i\&PartialD;j},---\left(6\right)$

其中i和j等于x，y或z，和图像强度函数I。在空间中的点P，海赛矩阵(Hessian matrix)M定义为

M(P)＝[A_ij]_{[i＝x..z，j＝x..z]}。

因为矩阵M定义为对称正矩阵，因此特征值的计算就一目了然了。在血管中心，3个特征值λ₁，λ₂和λ₃验证了下面的方程。

由特征向量给出两个主方向。和λ₁相关的向量V1对应局部血管走向，同时V2和V3定义血管切面。

对于每一个候选点P，由

和

定义的切面上的强度函数通过方差σ递增的高斯函数G_σ滤波。G_σ最小化(I-G_σ)的L²范数，G_σ被当作模型。然后，势函数V_σ被建立并被最小化。V_σ表示在点P处高斯血管模型和方差σ的匹配程度

$A=\frac{\left|{\mathrm{\&lambda;}}_2\right|}{\left|{\mathrm{\&lambda;}}_3\right|}$

$B=\frac{\left|{\mathrm{\&lambda;}}_2\right|}{\sqrt{\left|{\mathrm{\&lambda;}}_3{\mathrm{\&lambda;}}_2\right|}}$

$S=\sqrt{{\mathrm{\&lambda;}}_1+{\mathrm{\&lambda;}}_2+{\mathrm{\&lambda;}}_3}$

$V_{\mathrm{\&sigma;}}\left(p\right)=(1-e^{-\frac{A^2}{2{\mathrm{\&alpha;}}^2}})e^{-\frac{B^2}{2{\mathrm{\&beta;}}^2}}(1-e^{-\frac{A^S}{2{\mathrm{\&gamma;}}^2}})+\&Integral;{|G_{\mathrm{\&sigma;}}-I|}^2$

其中α，β和γ用于控制V_σ对于不同比例A，B和S的灵敏度。这种方法基于强度分析。例如，其他方法可能包括求邻近区域候选点的平均值。但是，噪声和依赖于光线密度的点之间的大间距使得这些方法不适合本申请。为了解决血管点之间的大间距，本发明考虑一种增加点的数量并将这些点连接起来产生完整的冠状动脉树的血管跟踪技术。

在特征分析之后，假定已将这些点正确地放置，即使这些点的数量还不足以跟踪血管中心线。因此，为了放置更多的点，在线性强度变化、到被分割的心壁距离不变且方向不变的情况下，对每一个3D候选点执行简单的血管跟踪操作。换句话说，从由特征分析产生的点，本发明建立另一个点并和前一个点连接起来，从而跟踪血管。速度向量

是两个向量的加权和：速度向量和源于局部血管走向的局部方向向量

$\stackrel{\&RightArrow;}{v_t}=\mathrm{\&lambda;}\stackrel{\&RightArrow;}{v_{t-1}}+(1-\mathrm{\&lambda;})\stackrel{\&RightArrow;}{v_l}$

其中λ∈[0，1]。该速度向量，强度以及当前候选点到心脏的距离必须均匀。本发明通过拉普拉斯算子分析检测边缘：

$\&Integral;\left|\stackrel{\&RightArrow;}{\&dtri;}I\stackrel{\&RightArrow;}{\mathrm{dl}}\right|<I_{\mathrm{coronary}}-I_{\mathrm{exterior}}$

如果这种不等式不适用于距离函数D和图像强度I，则边缘已经被检测到，且跟踪操作停止。否则，该跟踪操作被执行直到到达血管末端，或直到能够将所有的链连接到另一个点上。在这个阶段，候选点被集中并在血管中心线上均匀分布。如图5所示，通过最小生成树算法，这些候选点被连到一起，该最小生成树算法使P₁和P₂两点间的下述价值函数(cost function)最小化：

$E_{P_1{P}_2}={|D\left(P_1\right)-D\left(P_2\right)|}^2+{|I\left(P_1\right)-I\left(P_2\right)|}^2+{\left|\right|\stackrel{\&RightArrow;}{P_1{P}_2}\left|\right|}^2+\left|\stackrel{\&RightArrow;}{v_1}\left(P_1\right)\stackrel{\&RightArrow;}{v_1}\left(P_2\right)\right|,$

并形成外周血管中心线。更具体而言，图5显示了连接在一起的血管中心线上的点。

现在参考图6，和上面较详细的描述的一样，图中显示了根据本发明一个典型实施方案的基于光线滤波方法的心脏外周血管分割方法。在逐步的方法中，首先例如使用根据图形分割和测量表面的方法对心脏进行分割(在605)。从分割输出中计算3D距离映像(在610)。然后，从该距离映像，计算对心脏建模的双轴椭球(在615)。从椭球表面向心脏的方向投射光线(在620)。然后，沿着n-MIP光线轮廓，检测可能属于血管的候选点(在625)。最后，从候选点跟踪血管中心线(在630)。

例如，这种方法能够用于显示整个外围血管树和检测动脉粥样斑和动脉瘤。一旦被分割，血管结构能够被展开，而且管道狭窄和动脉瘤的数量将一目了然。通过允许自动生成整体视图显示、用于虚拟内窥镜的引导、产生沿血管的弯曲MPR或横截面图形，分割结果支持高效率的报告。此外，和应用到CT数据集一样，这种方法检测管腔结构对于工业应用来说足够地快(对于900MHz、1GBRAM双处理器而言需要60秒)。另外，这种方法与依赖前向传播(frontpropagation)的技术相比具有很强的抗噪声性能。

以上公开的特定实施例仅是说明性的，因为对受益于本发明的教导的本领域技术人员而言显而易见的是，可以用不同的但等同的方式修改和实现本发明。而且，本发明不受这里所描述的具体结构或设计的限制，而是由下面的权利要求所限制，因此显而易见的是上面所公开的特定实施例可以被修改或改变，并且所有这些改变都视为在本发明的精神和范围内。因此，这里所请求保护的范围由下面的权利要求限定。

Claims (14)

1.一种自动获得图像中的心脏外周血管中心线的方法，包括：

对图像中的心脏分割以产生三维分割的心肌；

从三维分割的心肌计算三维距离映像；

对三维距离映像中的心脏表面上的每一个体素投射法向最大强度投影光线；

对光线滤波以获得血管候选点；和

精选血管候选点以获得中心线。

2.如权利要求1所述的方法，其中心脏分割步骤包括使用图形分割方法对心脏分割。

3.如权利要求1所述的方法，其中心脏分割步骤包括使用基于模型的分割方法对心脏分割。

4.如权利要求1所述的方法，其中心脏分割步骤包括使用水平集方法对心脏分割。

5.如权利要求1所述的方法，其中每一条光线是与一直线相交的体素强度的集合。

6.如权利要求1所述的方法，其中对光线滤波以获得血管候选点的步骤包括：检测高强度的峰作为血管候选点。

7.如权利要求6所述的方法，其中检测高强度的峰作为血管候选点的步骤包括：如果峰在预定强度范围内，则检测高强度的峰作为血管候选点。

8.如权利要求6所述的方法，其中检测高强度的峰作为血管候选点的步骤包括：如果峰到心壁为预定距离，则检测高强度的峰作为血管候选点。

9.如权利要求6所述的方法，其中检测高强度的峰作为血管候选点的步骤包括：如果峰在预定范围内，则检测高强度的峰作为血管候选点。

10.如权利要求9所述的方法，其中如果峰在预定范围内则检测高强度的峰作为血管候选点的步骤包括：如果峰在预定的1mm到10mm的范围内，则检测高强度的峰作为血管候选点。

11.如权利要求6所述的方法，其中检测高强度的峰作为血管候选点的步骤包括：使用局部最大值检测方法。

12.如权利要求1所述的方法，其中精选血管候选点以获得中心线的步骤包括：计算针对每一个血管候选点的海赛矩阵和特征值。

13.如权利要求1所述的方法，其中图像是CT图像和MR图像中的一种。

14.一种自动获得图像中的心脏外周血管中心线的方法，包括：

对图像中的心脏分割以产生三维分割的心肌；

从三维分割的心肌计算三维距离映像；

使用三维距离映像，按照回转椭球对心壁建模；

对回转椭球上的每一个体素投射法向最大强度投影光线；

对光线滤波以获得血管候选点；和

精选血管候选点以获得中心线。

Images