CN111951277A - 一种基于cta图像的冠状动脉分割方法 - Google Patents

一种基于cta图像的冠状动脉分割方法 Download PDF

Info

Publication number
CN111951277A
CN111951277A CN202010740171.3A CN202010740171A CN111951277A CN 111951277 A CN111951277 A CN 111951277A CN 202010740171 A CN202010740171 A CN 202010740171A CN 111951277 A CN111951277 A CN 111951277A
Authority
CN
China
Prior art keywords
coronary artery
layer
region
formula
point
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010740171.3A
Other languages
English (en)
Other versions
CN111951277B (zh
Inventor
马骏
祝磊
薛凌云
张子恒
徐平
刘亦安
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hangzhou Dianzi University
Original Assignee
Hangzhou Dianzi University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hangzhou Dianzi University filed Critical Hangzhou Dianzi University
Priority to CN202010740171.3A priority Critical patent/CN111951277B/zh
Publication of CN111951277A publication Critical patent/CN111951277A/zh
Application granted granted Critical
Publication of CN111951277B publication Critical patent/CN111951277B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/11—Region-based segmentation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/136—Segmentation; Edge detection involving thresholding
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10072—Tomographic images
    • G06T2207/10081—Computed x-ray tomography [CT]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30004—Biomedical image processing
    • G06T2207/30101—Blood vessel; Artery; Vein; Vascular

Landscapes

  • Engineering & Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Apparatus For Radiation Diagnosis (AREA)

Abstract

本发明公开一种基于CTA图像的冠状动脉分割方法。首先，通过图像预处理有效抑制了非冠脉组织，提升了冠状动脉和背景的对比度；其次，结合光流法与心脏解剖结构先验知识检测带冠状动脉分叉的不规则升主动脉层，避免了手动初始化左右冠状动脉的起始点；最后，相比于传统区域生长法，提出的结合端点检测的自适应区域生长法对于灰度不均且拓扑结构复杂的细小分支有更好的分割能力和准确性。

Description

一种基于CTA图像的冠状动脉分割方法
技术领域
本发明属于医学图像处理领域，具体涉及一种非诊断目的的基于CTA图像的冠状动脉分割方法，是一种在冠状动脉CT血管造影(Computed Tomography Angiography，CTA)图像上左右冠状动脉支起点识别以及冠状动脉分割的方法。
背景技术
心脏是一个通过循环系统将血液泵送到全身的器官，它保证了全身各组织可通过血液获得足够的氧气和营养，并清除二氧化碳和其他废物，而冠状动脉则是包裹在心脏外周为心脏自身提供补给的血管，对心脏的正常工作起到关键性的作用。整个冠状动脉树拓扑结构复杂，分叉多且细小，其直径范围一般在2mm到7mm之间。从升主动脉根部出发，冠状动脉可以分为左主冠状动脉(Left Main Coronary Artery，LMCA)和右冠状动脉(RightCoronary Artery，RCA)两个主要分支，LMCA主要负责向左心房和左心室供应血液，其向下延伸还可再细分为左旋支(Left Circumflex Artery，LCX)和左前降支(Left AnteriorDescending Artery，LAD)，RCA主要负责向右心房、右心室、窦房结以及房室结提供血液，其继续向下延伸还可细分为右后降支(Posterior Descending Artery，PDA)和右边缘支(Right Marginal Artery，MA)。
CTA图像是一种通过静脉注射造影剂后，利用CT设备对患者心脏进行断层扫描形成的一个二维图像序列。随着成像分辨率的不断提升以及技术本身无创性和便捷性的优点，CTA逐渐成为临床上早期筛查冠状动脉血管疾病的有效方式。基于CTA的冠状动脉精准分割是后续冠脉狭窄定量分析和斑块分类等多种病灶识别的重要前置步骤。目前，心脏病专家的手动分割结果是最准确的，但手动分割非常耗时，并且不同专家的人为主观因素也会导致分割结果存在细微差异。近年来，各种计算机辅助诊断技术(Computer AidedDiagnosis，CAD)被用于冠状动脉分割，但是大部分仍属于半自动分割方法，需要手动设定左右冠状动脉各段分支起始种子点以及相应的阈值才能完成后续的冠状动脉分割，且容易遗漏冠状动脉的细小分支。少部分全自动方法很容易产生过分割，将非冠状动脉的血管也误分割为冠状动脉，从而降低了分割的准确度。
发明内容
本发明的目的在于提供一种基于CTA图像的全自动冠状动脉分割方法，能够完整、准确地分割并提取冠状动脉树，辅助医生提高心血管疾病诊断效率。
一种基于CTA图像的冠状动脉分割方法，其特征在于包括如下步骤：
步骤一、获取原始CTA心脏图像；
步骤二、分割冠状动脉感兴趣区域，抑制肺部血管等非心脏组织并去除大部分与冠状动脉灰度值相近的心房与心室等心脏组织，减少非冠状动脉组织对后续冠脉分割产生的影响：
2-1对CTA原始图像进行窗宽窗位调整，去除肺部血管；
获取Dicom头文件中的窗宽(Window Width，WW)和窗位(Window Level，WL)信息，截取CT值在
Figure BDA0002606470390000021
范围内的原始CTA图像并归一化到0～255的灰度级，CT值小于
Figure BDA0002606470390000022
的像素点设置为0，CT值大于
Figure BDA0002606470390000023
的像素点设置为255；映射过程如公式(1)所示：
Figure BDA0002606470390000024
式(1)中y表示归一化后每个像素点的灰度值,x表示原始CTA图像中每个像素点的CT值；
2-2对窗宽窗位调整后的CTA图像进行基于改进灰狼优化算法的多级阈值处理，用以去除与冠状动脉灰度值相近的心房与心室等组织，获得冠状动脉感兴趣区域；
2-2-1设定最大迭代次数K；随机初始化N组可行解向量
Figure BDA0002606470390000025
每组可行解向量
Figure BDA0002606470390000026
包含dim个取值范围在0到255之间的灰度阈值，可将0到255的图像灰度范围分为dim+1个阈值区间C0,C1…Cdim，其中上述每一组可行解向量
Figure BDA0002606470390000027
都为目标函数可行解；
2-2-2将上述N组可行解向量分别代入到目标函数，将目标函数值按照大小排序，取前三大目函数标值所对应的可行解向量分别赋予最优解
Figure BDA0002606470390000028
次优解
Figure BDA0002606470390000029
第三优解
Figure BDA00026064703900000210
多级阈值的目标函数：
Figure BDA00026064703900000211
式(2)中Hm是某可行解向量内灰度阈值所构建的第m个阈值区间内的熵函数；
Figure BDA0002606470390000031
式(3)中H0,H1,…Hm为每个阈值区间内的Kapur熵，t1,t2,…tm为每组可行解向量
Figure BDA0002606470390000032
中的灰度阈值，p0，p1…pi为阈值区间内每一个灰度值在图像中出现的概率，ω0,ω1…ωm为每一个阈值区间的概率，L＝256；
2-2-3根据公式(4-6)更新参数向量
Figure BDA0002606470390000033
和控制系数a：
Figure BDA0002606470390000034
Figure BDA0002606470390000035
a＝2(1-t/K)(6)
其中a根据迭代次数t从2线性递减至0，K为总的迭代次数，
Figure BDA0002606470390000036
是[0,1]内的随机向量,n＝1,2,3；
2-2-4根据当前迭代次数中三个等级依次递减的最优解向量
Figure BDA0002606470390000037
通过公式(7-10)更新其余解向量
Figure BDA0002606470390000038
Figure BDA0002606470390000039
式(7)中
Figure BDA00026064703900000310
代表每一组解向量
Figure BDA00026064703900000311
和三个最优解向量
Figure BDA00026064703900000312
之间的距离，
Figure BDA00026064703900000313
由公式(5)得到；
Figure BDA00026064703900000314
式(8)中
Figure BDA00026064703900000315
是每一组可行解向量
Figure BDA00026064703900000316
根据三个最优解解向量
Figure BDA00026064703900000317
计算得到的
Figure BDA00026064703900000318
候选可行解向量，
Figure BDA00026064703900000319
由公式(4)得到；
根据公式(2)分别计算
Figure BDA0002606470390000041
的目标函数值J1、J2、J3并根据公式(9)计算权重系数w1、w2、w3；
Figure BDA0002606470390000042
最后利用公式(10)根据权重系数w1、w2、w3更新每一组可行解向量
Figure BDA0002606470390000043
得到
Figure BDA0002606470390000044
Figure BDA0002606470390000045
2-2-5将更新后的N组可行解向量
Figure BDA0002606470390000046
分别代入到目标函数，将目标函数值按照大小排序，取前三大目函数标值所对应的可行解向量分别赋予最优解
Figure BDA0002606470390000047
次优解
Figure BDA0002606470390000048
第三优解
Figure BDA0002606470390000049
判断当前迭代次数是否满足最大迭代次数，若满足则结束，输出最优解
Figure BDA00026064703900000410
反之则返回2-2-3继续迭代；
2-2-6根据最优解
Figure BDA00026064703900000411
对窗宽窗位调整后的CTA图像保留最高一级灰度阈值区间Cdim内的像素，最终得到冠状动脉感兴趣区域；
步骤三、基于多尺度融合的血管特征提取，增强不同直径大小的冠状动脉；首先对CTA数据内所有体素点计算不同尺度σ下的Hessian矩阵H(A)，如式(11)所示；
Figure BDA00026064703900000412
式(10)中每一个子项Iij(A)表示图像I在点A处的二阶偏导，矩阵H(A)存在三个特征值λ1、λ2、λ3且满足|λ1|≤|λ2|≤|λ3|的关系，以及三个特征值各自对应的特征向量为
Figure BDA00026064703900000413
Frangi血管增强函数V(σ,A)可以利用多尺度Hessian矩阵计算所得的三个特征值λ1、λ2、λ3对血管状区域进行增强，对非血管区域进行抑制，血管相似性函数V(σ,A)如式(12)所示：
Figure BDA00026064703900000414
式(11)中，RA、RB用于量度管状、盘状和斑状等不同结构之间的特征差异，S是为对比度系数，用于描述目标和周围背景的亮度差异，S的值就越大，亮度差异越大；Frangi血管增强函数V(σ,A)的取值范围在0到1之间，V(σ,A)越接近于1表示该点处的组织与管状结构的相似性最大；α、β和γ为三个调参系数，通过影响RA、RB和S在血管相似性函数中所占的比重来调整式(12)的灵敏度；RA、RB和S定义分别如式(13)、式(14)和式(15)表示：
Figure BDA0002606470390000051
Figure BDA0002606470390000052
Figure BDA0002606470390000053
步骤四、升主动脉分割以及左右冠状动脉起始层识别；结合光流法检测相邻两层升主动脉间光流位移距离最大值的变化和心脏医学解剖结构知识判断左右冠状动脉的起始层；
4-1由于冠状动脉的根节点位于升主动脉根部，升主动脉的分割和提取是识别左右冠状动脉起始层的前提；采用霍夫变换检测类圆形的升主动脉区域，以被检测到的类圆形圆心作为种子点，利用区域生长法分割每一层升主动脉，最后得到完整的升主动脉；
所述的霍夫变换检测类圆形的升主动脉区域属于常规技术，故不详解；
4-2使用Lucas-Kanade光流法检测相邻两层升主动脉间的最大光流及所有像素中的最大位移距离；
在CTA图像序列中，大部份的升主动脉区域都呈现类圆形，但是当某一升主动脉层包含冠状动脉分支时，类圆形升主动脉区域会多一小段凸起的冠状动脉区域，该层相对于前一层在冠状动脉根结点处的光流位移距离会存在巨大的变化，通过对比相邻2层间光流位移距离最大值的变化即可初步筛选出冠状动脉的起始层；光流法假设相邻的2层CTA图像强度为I(x,y,t)和I(x+Δx,y+Δy,t+Δt)，并且相邻的2层图像强度不变，可以得到公式(16)
I(x+Δx,y+Δy,t+Δt)＝I(x,y,t)(16)
将公式(16)使用泰勒级数展开后可得：
Figure BDA0002606470390000054
其中
Figure BDA0002606470390000055
为高次项，大小可以忽略，为此结合公式(16)和(17)可得：
Figure BDA0002606470390000056
其中
Figure BDA0002606470390000057
和
Figure BDA0002606470390000058
分别为图像相对于X轴和Y轴的空间导数，
Figure BDA0002606470390000059
为时间导数，公式(18)可改写为
(Ix,Iy)(u,v)＝-It(19)
其中，u和v就是图像中每个像素分别沿着X轴和Y轴的光流位移距离，但是只有公式(19)一个约束方程无法求解2个未知变量，所以光流法假设目标像素的周围像素也存在相同的光流位移距离，即在目标像素周围m*m的窗口内，每个像素点都具有相似的光流位移距离u和v；简单的来说，就是通过增加约束方程来求解光流矢量u和v，如公式(20)所示；
Figure BDA0002606470390000061
式(19)中(Xn，Yn)表示m*m的窗口内第n个像素点；
通过公式(20)求解所得的u和v可以计算出每个像素点的位移距离，如公式(21)所示；
Figure BDA0002606470390000062
根据每一层Distances最大的点，筛选出最大的前10层做为候选冠脉起始层，然后结合三个心脏医学解剖结构先验知识进一步筛选出左右冠状动脉的起始层；
所述的三个心脏医学解剖结构先验知识具体如下：
第一，升主动脉根部只存在左右两个冠状动脉开口，在10层冠脉候选层中，位移距离最大的是左冠状动脉起始层Sleft；
第二，左右冠脉起始层间相隔10层以上，且左冠状动脉起始层位于右冠状动脉起始层的Z轴上方；
第三，左冠状动脉起始点OL位于该层质心的-90°～45°之间，右冠状动脉起始点OR位于该层质心的位于45°～135°之间；
最后在满足上述第二、三条件的剩余候选冠脉起始层中选择层号最小的层作为右侧冠状动脉起始层Sright；
步骤五、在步骤三增强后的CTA图像上，从步骤四获得的左右冠状动脉起始层开始，使用结合端点检测的区域生长法分割每一层CTA图像中的血管区域；具体是从步骤三获得的左右冠状动脉起始层Sleft和Sright开始，沿着Z轴向下或向上利用迭代区域生长法进行冠状动脉的初始分割；
所述的迭代区域生长法由步骤5-1种子点自动选择、步骤5-2迭代区域生长法生长准则、步骤5-3端点检测以及步骤5-4端点再生长四部分组成：
5-1种子点自动选择
利用血管上下层之间的空间连通性，即当前层Si与相邻层Si+1血管区域内的二维坐标存在部分重叠作为每层种子点的选择依据；由于三维的冠状动脉血管分叉众多，所以在二维的单层冠状动脉CTA图像中可能存在多个独立的血管连通区域，当Si层分割完成后，先检测血管连通区域个数c，然后分别在Si层的每个独立的血管区域中以一定的步长选取M个种子点作为相邻层Si+1血管区域生长的种子点；M为经验值；
5-2迭代区域生长法生长准则
利用迭代区域生长法分割每一层的冠状动脉血管区域；
5-2-1首先从5-1中选取的每一层c*M个种子点依次开始对该层血管区域进行分割，每一个种子点会将区域边缘像素点8领域内与种子点灰度值相差在1％以内的像素点不断归并到感兴趣区域中，直到不能生长为止；
5-2-2由于5-1中选取的种子点只有部分会落在血管区域内，所以每个种子点生长完成后，通过计算分割得到的血管感兴趣区域内的像素数量num判断该种子点是否在血管区域内：如果像素数num小于Nthresh，则表明该种子点位于血管区域内，保留该种子点分割得到的血管感兴趣区域；如果像素数num大于Nthresh，则说明该种子点落在血管区域之外，出现了严重的过分割现象，丢弃该种子点分割得到的区域，Nthresh为经验值；然后将所有保留下来的血管感兴趣区域之和作为该层最终的冠状动脉血管分割结果；
5-2-3判断每一层最终的冠状动脉血管分割结果像素数是否为0，如果等于0，则表明无法再继续分割血管，结束迭代，得到初始冠状动脉感兴趣区域；如果不等于0，则返回5-1继续执行下一层种子点的选择以及5-2迭代区域生长法；
5-3端点检测
对5-2中提取的初始冠状动脉感兴趣区域利用细化方法提取血管骨架线；具体是利用血管骨架线的26邻域单连通性，将骨架线上的体素分为端点P1、普通连接点P2、分叉点P3；从各个端点P1和分叉点P3开始，继续使用结合5-1的种子点选择机制和5-2的迭代区域生长法尝试Z轴向上和向下的双向生长，用来分割冠状动脉分支中拓扑结构复杂，上下起伏的血管部分；
所述端点P1在该体素的26邻域点中，只有1个邻域点位于骨架线上；所述普通连接点P2，在该体素的26邻域点中，只有2个邻域点属于骨架线上；所述分叉点P3，在该体素的26邻域点中，至少有3个邻域点位于骨架线上；
5-4端点再生长
判断结合端点检测的区域生长法继续分割后的冠状动脉区域是否有新的端点P1和分叉点P3，如果有，则返回5-4继续分割；如果没有，则结束分割。
步骤二中采用改进的灰狼优化算法，式(8)中采用离散的权重w1、w2、w3对每一次迭代中各组可行解向量
Figure BDA0002606470390000081
进行更新，相比于传统灰狼优化算法采用的平均权重
Figure BDA0002606470390000082
有助于提高最优解
Figure BDA0002606470390000083
的搜寻效率。
步骤三中采用基于多尺度融合的血管特征提取，增强不同直径大小的冠状动脉。由于整个冠状动脉树从升主动脉根部出发会变得越来越细，其直径范围一般在2mm到7mm之间。单一尺度σ下的Hessian矩阵H(A)特征值只能用于增强特定直径的血管，为此本专利计算多个不同尺度σ下的Hessian矩阵H(A)的特征值λ1、λ2、λ3用于提取不同直径血管的特征。最后，Frangi血管增强函数V(σ,A)可利用不同尺度σ下Hessian矩阵的特征值对不同直径的血管状区域进行增强，对非血管区域进行抑制。
本发明的有益效果如下：
本发明提出了一种全自动的冠状动脉分割方法，克服了半自动方法中需要手动初始化冠状动脉起始种子点的技术难点。本发明首先利用窗宽窗位调整和基于灰狼优化的多级阈值处理，有效去除原始CTA图像中非冠状动脉组织可能对后续冠状动脉分割产生的干扰；然后利用升主动脉和冠状动脉的空间位置关系，采用光流法识别冠状动脉的起始层；最后从左右冠状动脉的起始层开始，采用结合端点检测的区域生长法对冠状动脉进行分割，相比于传统区域生长法，对于拓扑结构复杂的细小分支有更好的分割能力和准确性。
附图说明
图1a为原始CTA图像；
图1b为原始CTA图像进行窗宽窗位调整和基于改进灰狼算法优化的多级阈值处理后的结果；
图2为基于多尺度融合的血管增强结果；
图3a为分割得到的升主动脉三维重构结果；
图3b为光流法检测得到的每层最大光流位移距离结果；
图3c为左右冠状动脉起始层中，具有最大光流位移距离的像素坐标和该层质心间的位置关系图；
图3d-e分别为光流法检测得到的左右冠状动脉起始层；
图4a为使用的冠状动脉骨架线示意图；
图4b为传统区域生长法分割的冠状动脉结果；
图4c为结合端点检测的区域生长法分割的冠状动脉结果；
图5为本发明方法流程图；
图6为基于改进灰狼优化算法的多级阈值处理流程图；
图7为基于结合端点检测的区域生长法冠状动脉分割流程图；
具体实施方式
以下结合附图对本发明作进一步说明。
如图5一种基于CTA图像的冠状动脉分割方法，包括如下步骤：
步骤一、获取原始CTA心脏图像；
步骤二、分割冠状动脉感兴趣区域：
2-1对CTA原始图像进行窗宽窗位调整，去除肺部血管；
获取Dicom头文件中的窗宽(Window Width，WW)和窗位(Window Level，WL)信息，截取CT值在
Figure BDA0002606470390000091
范围内的原始CTA图像并归一化到0～255的灰度级，CT值小于
Figure BDA0002606470390000092
的像素点设置为0，CT值大于
Figure BDA0002606470390000093
的像素点设置为255；映射过程如公式(1)所示：
Figure BDA0002606470390000094
式(1)中y表示归一化后每个像素点的灰度值,x表示原始CTA图像中每个像素点的CT值；
2-2对窗宽窗位调整后的CTA图像进行基于改进灰狼优化算法的多级阈值处理，用以去除与冠状动脉灰度值相近的心房与心室等组织，获得冠状动脉感兴趣区域如图6，具体是：
2-2-1设定最大迭代次数K；随机初始化N组可行解向量
Figure BDA0002606470390000095
每组可行解向量
Figure BDA0002606470390000096
包含dim个取值范围在0到255之间的灰度阈值，可将0到255的图像灰度范围分为dim+1个阈值区间C0,C1…Cdim，其中上述每一组可行解向量
Figure BDA0002606470390000097
都为目标函数可行解；
2-2-2将上述N组可行解向量分别代入到目标函数，将目标函数值按照大小排序，取前三大目函数标值所对应的可行解向量分别赋予最优解
Figure BDA0002606470390000098
次优解
Figure BDA0002606470390000099
第三优解
Figure BDA00026064703900000910
多级阈值的目标函数：
Figure BDA0002606470390000101
式(2)中Hm是某可行解向量内灰度阈值所构建的第m个阈值区间内的熵函数；
Figure BDA0002606470390000102
式(3)中H0,H1,…Hm为每个阈值区间内的Kapur熵，t1,t2,…tm为每组可行解向量
Figure BDA0002606470390000103
中的灰度阈值，p0，p1…pi为阈值区间内每一个灰度值在图像中出现的概率，ω0,ω1…ωm为每一个阈值区间的概率，L＝256；
2-2-3根据公式(4-6)更新参数向量
Figure BDA0002606470390000104
和控制系数a：
Figure BDA0002606470390000105
Figure BDA0002606470390000106
a＝2(1-t/K)(6)
其中a根据迭代次数t从2线性递减至0，K为总的迭代次数，
Figure BDA0002606470390000107
是[0,1]内的随机向量,n＝1,2,3；
2-2-4根据当前迭代次数中三个等级依次递减的最优解向量
Figure BDA0002606470390000108
通过公式(7-10)更新其余解向量
Figure BDA0002606470390000109
Figure BDA00026064703900001010
式(7)中
Figure BDA00026064703900001011
代表每一组解向量
Figure BDA00026064703900001012
和三个最优解向量
Figure BDA00026064703900001013
之间的距离，
Figure BDA0002606470390000111
由公式(5)得到；
Figure BDA0002606470390000112
式(8)中
Figure BDA0002606470390000113
是每一组可行解向量
Figure BDA0002606470390000114
根据三个最优解解向量
Figure BDA0002606470390000115
计算得到的
Figure BDA0002606470390000116
候选可行解向量，
Figure BDA0002606470390000117
由公式(4)得到；
根据公式(2)分别计算
Figure BDA0002606470390000118
的目标函数值J1、J2、J3并根据公式(9)计算权重系数w1、w2、w3；
Figure BDA0002606470390000119
最后利用公式(10)根据权重系数w1、w2、w3更新每一组可行解向量
Figure BDA00026064703900001110
得到
Figure BDA00026064703900001111
Figure BDA00026064703900001112
2-2-5将更新后的N组可行解向量
Figure BDA00026064703900001113
分别代入到目标函数，将目标函数值按照大小排序，取前三大目函数标值所对应的可行解向量分别赋予最优解
Figure BDA00026064703900001114
次优解
Figure BDA00026064703900001115
第三优解
Figure BDA00026064703900001116
判断当前迭代次数是否满足最大迭代次数，若满足则结束，输出最优解
Figure BDA00026064703900001117
反之则返回2-2-3继续迭代；
2-2-6根据最优解
Figure BDA00026064703900001118
对窗宽窗位调整后的CTA图像保留最高一级灰度阈值区间Cdim内的像素，最终得到冠状动脉感兴趣区域；
步骤三、基于多尺度融合的血管特征提取，增强不同直径大小的冠状动脉；首先对CTA数据内所有体素点计算不同尺度σ下的Hessian矩阵H(A)，如式(11)所示；
Figure BDA00026064703900001119
式(10)中每一个子项Iij(A)表示图像I在点A处的二阶偏导，矩阵H(A)存在三个特征值λ1、λ2、λ3且满足|λ1|≤|λ2|≤|λ3|的关系，以及三个特征值各自对应的特征向量为
Figure BDA00026064703900001120
Frangi血管增强函数V(σ,A)可以利用多尺度Hessian矩阵计算所得的三个特征值λ1、λ2、λ3对血管状区域进行增强，对非血管区域进行抑制，血管相似性函数V(σ,A)如式(12)所示：
Figure BDA0002606470390000121
式(11)中，RA、RB用于量度管状、盘状和斑状等不同结构之间的特征差异，S是为对比度系数，用于描述目标和周围背景的亮度差异，S的值就越大，亮度差异越大；Frangi血管增强函数V(σ,A)的取值范围在0到1之间，V(σ,A)越接近于1表示该点处的组织与管状结构的相似性最大；α、β和γ为三个调参系数，通过影响RA、RB和S在血管相似性函数中所占的比重来调整式(12)的灵敏度；RA、RB和S定义分别如式(13)、式(14)和式(15)表示：
Figure BDA0002606470390000122
Figure BDA0002606470390000123
Figure BDA0002606470390000124
步骤四、升主动脉分割以及左右冠状动脉起始层识别；结合光流法检测相邻两层升主动脉间光流位移距离最大值的变化和心脏医学解剖结构知识判断左右冠状动脉的起始层；
4-1由于冠状动脉的根节点位于升主动脉根部，升主动脉的分割和提取是识别左右冠状动脉起始层的前提；采用霍夫变换检测类圆形的升主动脉区域，以被检测到的类圆形圆心作为种子点，利用区域生长法分割每一层升主动脉，最后得到完整的升主动脉；
所述的霍夫变换检测类圆形的升主动脉区域属于常规技术，故不详解；
4-2使用Lucas-Kanade光流法检测相邻两层升主动脉间的最大光流及所有像素中的最大位移距离；
在CTA图像序列中，大部份的升主动脉区域都呈现类圆形，但是当某一升主动脉层包含冠状动脉分支时，类圆形升主动脉区域会多一小段凸起的冠状动脉区域，该层相对于前一层在冠状动脉根结点处的光流位移距离会存在巨大的变化，通过对比相邻2层间光流位移距离最大值的变化即可初步筛选出冠状动脉的起始层；光流法假设相邻的2层CTA图像强度为I(x,y,t)和I(x+Δx,y+Δy,t+Δt)，并且相邻的2层图像强度不变，可以得到公式(16)
I(x+Δx,y+Δy,t+Δt)＝I(x,y,t)(16)
将公式(16)使用泰勒级数展开后可得：
Figure BDA0002606470390000131
其中
Figure BDA0002606470390000132
为高次项，大小可以忽略，为此结合公式(16)和(17)可得：
Figure BDA0002606470390000133
其中
Figure BDA0002606470390000134
和
Figure BDA0002606470390000135
分别为图像相对于X轴和Y轴的空间导数，
Figure BDA0002606470390000136
为时间导数，公式(18)可改写为
(Ix,Iy)(u,v)＝-It(19)
其中，u和v就是图像中每个像素分别沿着X轴和Y轴的光流位移距离，但是只有公式(19)一个约束方程无法求解2个未知变量，所以光流法假设目标像素的周围像素也存在相同的光流位移距离，即在目标像素周围m*m的窗口内，每个像素点都具有相似的光流位移距离u和v；简单的来说，就是通过增加约束方程来求解光流矢量u和v，如公式(20)所示；
Figure BDA0002606470390000137
式(19)中(Xn，Yn)表示m*m的窗口内第n个像素点；
通过公式(20)求解所得的u和v可以计算出每个像素点的位移距离，如公式(21)所示；
Figure BDA0002606470390000138
根据每一层Distances最大的点，筛选出最大的前10层做为候选冠脉起始层，然后结合三个心脏医学解剖结构先验知识进一步筛选出左右冠状动脉的起始层；
所述的三个心脏医学解剖结构先验知识具体如下：
第一，升主动脉根部只存在左右两个冠状动脉开口，在10层冠脉候选层中，位移距离最大的是左冠状动脉起始层Sleft；
第二，左右冠脉起始层间相隔10层以上，且左冠状动脉起始层位于右冠状动脉起始层的Z轴上方；
第三，左冠状动脉起始点OL位于该层质心的-90°～45°之间，右冠状动脉起始点OR位于该层质心的位于45°～135°之间；
最后在满足上述第二、三条件的剩余候选冠脉起始层中选择层号最小的层作为右侧冠状动脉起始层Sright；
步骤五、在步骤三增强后的CTA图像上，从步骤四获得的左右冠状动脉起始层开始，使用结合端点检测的区域生长法分割每一层CTA图像中的血管区域；具体是从步骤三获得的左右冠状动脉起始层Sleft和Sright开始，沿着Z轴向下或向上利用迭代区域生长法进行冠状动脉的初始分割；
所述的迭代区域生长法由步骤5-1种子点自动选择、步骤5-2迭代区域生长法生长准则、步骤5-3端点检测以及步骤5-4端点再生长四部分组成，如图7具体是：
5-1种子点自动选择
利用血管上下层之间的空间连通性，即当前层Si与相邻层Si+1血管区域内的二维坐标存在部分重叠作为每层种子点的选择依据；由于三维的冠状动脉血管分叉众多，所以在二维的单层冠状动脉CTA图像中可能存在多个独立的血管连通区域，当Si层分割完成后，先检测血管连通区域个数c，然后分别在Si层的每个独立的血管区域中以一定的步长选取M(M＝15)个种子点作为相邻层Si+1血管区域生长的种子点；
举例：血管会有很多分叉，每个分叉在一层图像上就是一个区域，如果有三个分叉，那么这一层就有三个血管区域了，在每个血管区域中取15个种子点，一共会有15*3＝45个种子点；
5-2迭代区域生长法生长准则
利用迭代区域生长法分割每一层的冠状动脉血管区域；
5-2-1首先从5-1中选取的每一层c*15个种子点依次开始对该层血管区域进行分割，每一个种子点会将区域边缘像素点8领域内与种子点灰度值相差在1％以内的像素点不断归并到感兴趣区域中，直到不能生长为止；
5-2-2由于5-1中选取的种子点只有部分会落在血管区域内，所以每个种子点生长完成后，通过计算分割得到的血管感兴趣区域内的像素数量num判断该种子点是否在血管区域内：如果像素数num小于Nthresh，则表明该种子点位于血管区域内，保留该种子点分割得到的血管感兴趣区域；如果像素数num大于Nthresh，则说明该种子点落在血管区域之外，出现了严重的过分割现象，丢弃该种子点分割得到的区域，根据经验Nthresh＝2000；然后将所有保留下来的血管感兴趣区域之和作为该层最终的冠状动脉血管分割结果；
5-2-3判断每一层最终的冠状动脉血管分割结果像素数是否为0，如果等于0，则表明无法再继续分割血管，结束迭代，得到初始冠状动脉感兴趣区域；如果不等于0，则返回5-1继续执行下一层种子点的选择以及5-2迭代区域生长法；
5-3端点检测
对5-2中提取的初始冠状动脉感兴趣区域利用细化方法提取血管骨架线；具体是利用血管骨架线的26邻域单连通性，将骨架线上的体素分为端点P1、普通连接点P2、分叉点P3；从各个端点P1和分叉点P3开始，继续使用结合5-1的种子点选择机制和5-2的迭代区域生长法尝试Z轴向上和向下的双向生长，用来分割冠状动脉分支中拓扑结构复杂，上下起伏的血管部分；
所述端点P1在该体素的26邻域点中，只有1个邻域点位于骨架线上；所述普通连接点P2，在该体素的26邻域点中，只有2个邻域点属于骨架线上；所述分叉点P3，在该体素的26邻域点中，至少有3个邻域点位于骨架线上；
5-4端点再生长
判断结合端点检测的区域生长法继续分割后的冠状动脉区域是否有新的端点P1和分叉点P3，如果有，则返回5-4继续分割；如果没有，则结束分割。
为证明上述方法的可行性，下面采用具体的冠状动脉CTA图像进行举例说明。本发明验证所采用的CTA数据集来自山东大学齐鲁医学院，该数据集由年龄在51-85岁之间的14例男性和6例女性病例组成，均带有医生标注，每例数据都有大约200-300层包含冠脉血管的CTA轴向切片，每层轴向切片图像分辨率为512*512，像素间距dx，dy，dz为0.35mm，0.35mm，0.80mm。
实验过程如下：先对原始的CTA图像序列(如图1a所示)进行窗宽窗位调整和基于改进灰狼算法优化的多级阈值处理，分割得到的冠状动脉感兴趣区域如图1b所示，可有效抑制肺部血管等非心脏组织并去除大部分与冠状动脉灰度值相近的心房与心室等心脏组织，减少非冠状动脉组织对后续冠脉分割产生的影响；对预处理后的图1b进行基于多尺度融合的血管特征提取，结果如图2所示，不同直径的冠状动脉血管被增强；图3a为利用霍夫变换和区域生长分割得到的升主动脉三维重构图，图中升主动脉的左右凸起分叉就是左右冠状动脉的起始根结点，包含冠状动脉的升主动脉层不是普通的类圆形。本发明首先通过光流法计算每一层所有像素相对于前一层的光流位移距离，并保存每一层中的最大值，结果如图3b所示。在该例CTA数据的90层升主动脉中，普通类圆形升主动脉层，不会出现位移激增的现象，对应于图3b曲线中较平缓的区域，但是当遇到包含冠状动脉根节点的不规则升主动脉层时，类圆形的升主动脉层由于多了一段冠状动脉凸起，会出现最大位移距离激增的现象，所以可将图3b中5个最大峰值的层作为冠状动脉候选层，再结合冠状动脉的医学解剖结构先验知识进一步筛选得到左右冠状动脉的起始层。首先，5层冠状动脉候选层中位移距离最大峰即图3b中第53层就是左冠状动脉支的起始层，其次第53层中最大位移像素点的坐标(图3c中的黑圈)位于该层升主脉质心的-90°～45°之间，也符合心脏医学解剖结构的先验构知识，本例第53层左冠状动脉支的起始层结果如图3d所示。右冠状动脉起始层在本数据集中距离左冠状动脉支的起始层在10层以上即8mm左右，所以63层以后，峰值最大的74层最可能是右冠状动脉，其最大位移像素点的坐标(图3e中的黑圈)位于该层升主脉质心的45°～135°之间，也符合心脏医学解剖结构的先验构知识，本例第74层右冠状动脉支的起始层结果如图3e所示。本发明在20例CTA数据上对该算法进行了测试，结果如表1所示，左侧冠状动脉起始层检测成功19例，准确率达到95％，右侧冠状动脉起始层检测成功17例，准确率达到85％。
表1光流法识别冠状动脉起始层
Figure BDA0002606470390000161
最后从53层和74层分别使用结合端点检测的迭代区域生长法分割每一层CTA图像中的血管区域。图4b和4c所示分别为传统区域生长法与结合端点检测的区域生长法的冠状动脉分割结果对比，图4c所示结合端点检测的区域生长法对于冠状动脉末端拓扑结构复杂，上下起伏的细小分支有更好的分割能力。为了评估算法的分割结果精度，本发明使用Dice、Jaccard、MSD(Mean Squared Surface Distance)、MAXSD(Maximum SurfaceDistance)作为分割评价指标，结果如表2所示，采用本发明提出的结合端点检测的区域生长法得到的分割结果，Dice和Jaccard系数相比于传统区域生长法，分别平均提高了4％和8％，达到了0.70和0.57，MSD和MAXSD相比于传统区域生长法，分别平均降低了2％和4％，达到了0.40和2.66。
表2不同方法结果比较
Figure BDA0002606470390000162
本发明提出了一种面向CTA数据的冠状动脉自动分割方法。首先，通过图像预处理有效抑制了非冠脉组织，提升了冠状动脉和背景的对比度；其次，结合光流法与心脏解剖结构先验知识检测带冠状动脉分叉的不规则升主动脉层，避免了手动初始化左右冠状动脉的起始点；最后，相比于传统区域生长法，本文提出的结合端点检测的自适应区域生长法对于灰度不均且拓扑结构复杂的细小分支有更好的分割能力和准确性。

Claims (9)

1.一种基于CTA图像的冠状动脉分割方法，其特征在于包括如下步骤：
步骤一、获取原始CTA心脏图像；
步骤二、分割冠状动脉感兴趣区域；
步骤三、基于多尺度融合的血管特征提取，增强不同直径大小的冠状动脉；
步骤四、升主动脉分割以及左右冠状动脉起始层识别；
4-1采用霍夫变换检测类圆形的升主动脉区域，以被检测到的类圆形圆心作为种子点，利用区域生长法分割每一层升主动脉，最后得到完整的升主动脉；
4-2使用Lucas-Kanade光流法检测相邻两层升主动脉间的最大光流及所有像素中的最大位移距离，筛选出候选冠脉起始层，然后结合心脏医学解剖结构先验知识进一步筛选出左右冠状动脉的起始层；
步骤五、在步骤三增强后的CTA图像上，从步骤四获得的左右冠状动脉起始层开始，使用结合端点检测的区域生长法分割每一层CTA图像中的血管区域；从步骤三获得的左右冠状动脉起始层Sleft和Sright开始，沿着Z轴向下或向上利用迭代区域生长法进行冠状动脉的初始分割；具体是：
5-1种子点自动选择
利用血管上下层之间的空间连通性，即当前层Si与相邻层Si+1血管区域内的二维坐标存在部分重叠作为每层种子点的选择依据；由于三维的冠状动脉血管分叉众多，所以在二维的单层冠状动脉CTA图像中可能存在多个独立的血管连通区域，当Si层分割完成后，先检测血管连通区域个数c，然后分别在Si层的每个独立的血管区域中以一定的步长选取M个种子点作为相邻层Si+1血管区域生长的种子点；
5-2迭代区域生长法生长准则
利用迭代区域生长法分割每一层的冠状动脉血管区域；
5-3端点检测
对5-2中提取的初始冠状动脉感兴趣区域利用细化方法提取血管骨架线；具体是利用血管骨架线的26邻域单连通性，将骨架线上的体素分为端点P1、普通连接点P2、分叉点P3；从各个端点P1和分叉点P3开始，继续使用结合步骤5-1的种子点选择机制和步骤5-2的迭代区域生长法尝试Z轴向上和向下的双向生长，用来分割冠状动脉分支中拓扑结构复杂，上下起伏的血管部分；
5-4端点再生长
判断结合端点检测的区域生长法继续分割后的冠状动脉区域是否有新的端点P1和分叉点P3，如果有，则返回5-4继续分割；如果没有，则结束分割。
2.如权利要求1所述的一种基于CTA图像的冠状动脉分割方法，其特征在于步骤二具体是：
2-1对CTA原始图像进行窗宽窗位调整，去除肺部血管；
2-2对窗宽窗位调整后的CTA图像进行基于改进灰狼优化算法的多级阈值处理，用以去除与冠状动脉灰度值相近的心房与心室等组织，获得冠状动脉感兴趣区域。
3.如权利要求2所述的一种基于CTA图像的冠状动脉分割方法，其特征在于步骤2-1具体是：获取Dicom头文件中的窗宽(Window Width，WW)和窗位(Window Level，WL)信息，截取CT值在
Figure FDA0002606470380000021
范围内的原始CTA图像并归一化到0～255的灰度级，CT值小于
Figure FDA0002606470380000022
的像素点设置为0，CT值大于
Figure FDA0002606470380000023
的像素点设置为255；映射过程如公式(1)所示：
Figure FDA0002606470380000024
式(1)中y表示归一化后每个像素点的灰度值,x表示原始CTA图像中每个像素点的CT值。
4.如权利要求3所述的一种基于CTA图像的冠状动脉分割方法，其特征在于步骤2-2具体是：
2-2-1设定最大迭代次数K；随机初始化N组可行解向量
Figure FDA0002606470380000025
每组可行解向量
Figure FDA0002606470380000026
包含dim个取值范围在0到255之间的灰度阈值，可将0到255的图像灰度范围分为dim+1个阈值区间C0,C1…Cdim，其中上述每一组可行解向量
Figure FDA0002606470380000027
都为目标函数可行解；
2-2-2将上述N组可行解向量分别代入到目标函数，将目标函数值按照大小排序，取前三大目函数标值所对应的可行解向量分别赋予最优解
Figure FDA0002606470380000028
次优解
Figure FDA0002606470380000029
第三优解
Figure FDA00026064703800000210
多级阈值的目标函数：
Figure FDA0002606470380000031
式(2)中Hm是某可行解向量内灰度阈值所构建的第m个阈值区间内的熵函数；
Figure FDA0002606470380000032
式(3)中H0,H1,…Hm为每个阈值区间内的Kapur熵，t1,t2,…tm为每组可行解向量
Figure FDA0002606470380000033
中的灰度阈值，p0，p1…pi为阈值区间内每一个灰度值在图像中出现的概率，ω0,ω1…ωm为每一个阈值区间的概率，L＝256；
2-2-3根据公式(4-6)更新参数向量
Figure FDA0002606470380000034
和控制系数a：
Figure FDA0002606470380000035
Figure FDA0002606470380000036
a＝2(1-t/K) (6)
其中a根据迭代次数t从2线性递减至0，K为总的迭代次数，
Figure FDA0002606470380000037
是[0,1]内的随机向量,n＝1,2,3；
2-2-4根据当前迭代次数中三个等级依次递减的最优解向量
Figure FDA0002606470380000038
通过公式(7-10)更新其余解向量
Figure FDA0002606470380000039
Figure FDA00026064703800000310
式(7)中
Figure FDA00026064703800000311
代表每一组解向量
Figure FDA00026064703800000312
和三个最优解向量
Figure FDA00026064703800000313
之间的距离，
Figure FDA0002606470380000041
由公式(5)得到；
Figure FDA0002606470380000042
式(8)中
Figure FDA0002606470380000043
是每一组可行解向量
Figure FDA0002606470380000044
根据三个最优解解向量
Figure FDA0002606470380000045
计算得到的
Figure FDA0002606470380000046
候选可行解向量，
Figure FDA0002606470380000047
由公式(4)得到；
根据公式(2)分别计算
Figure FDA0002606470380000048
的目标函数值J1、J2、J3并根据公式(9)计算权重系数w1、w2、w3；
Figure FDA0002606470380000049
最后利用公式(10)根据权重系数w1、w2、w3更新每一组可行解向量
Figure FDA00026064703800000410
得到
Figure FDA00026064703800000411
Figure FDA00026064703800000412
2-2-5将更新后的N组可行解向量
Figure FDA00026064703800000413
分别代入到目标函数，将目标函数值按照大小排序，取前三大目函数标值所对应的可行解向量分别赋予最优解
Figure FDA00026064703800000414
次优解
Figure FDA00026064703800000415
第三优解
Figure FDA00026064703800000416
判断当前迭代次数是否满足最大迭代次数，若满足则结束，输出最优解
Figure FDA00026064703800000417
反之则返回2-2-3继续迭代；
2-2-6根据最优解
Figure FDA00026064703800000418
对窗宽窗位调整后的CTA图像保留最高一级灰度阈值区间Cdim内的像素，最终得到冠状动脉感兴趣区域。
5.如权利要求4所述的一种基于CTA图像的冠状动脉分割方法，其特征在于步骤三具体是：对CTA数据内所有体素点计算不同尺度σ下的Hessian矩阵H(A)，如式(11)所示；
Figure FDA00026064703800000419
式(10)中每一个子项Iij(A)表示图像I在点A处的二阶偏导，矩阵H(A)存在三个特征值λ1、λ2、λ3且满足|λ1|≤|λ2|≤|λ3|的关系，以及三个特征值各自对应的特征向量为
Figure FDA00026064703800000420
Figure FDA00026064703800000421
Frangi血管增强函数V(σ,A)可以利用多尺度Hessian矩阵计算所得的三个特征值λ1、λ2、λ3对血管状区域进行增强，对非血管区域进行抑制，血管相似性函数V(σ,A)如式(12)所示：
Figure 1
式(11)中，RA、RB用于量度管状、盘状和斑状等不同结构之间的特征差异，S是为对比度系数，用于描述目标和周围背景的亮度差异，S的值就越大，亮度差异越大；Frangi血管增强函数V(σ,A)的取值范围在0到1之间，V(σ,A)越接近于1表示该点处的组织与管状结构的相似性最大；α、β和γ为三个调参系数，通过影响RA、RB和S在血管相似性函数中所占的比重来调整式(12)的灵敏度；RA、RB和S定义分别如式(13)、式(14)和式(15)表示：
Figure 2
Figure 3
Figure 11
6.如权利要求5所述的一种基于CTA图像的冠状动脉分割方法，其特征在于步骤4-2具体是：在CTA图像序列中，大部份的升主动脉区域都呈现类圆形，但是当某一升主动脉层包含冠状动脉分支时，类圆形升主动脉区域会多一小段凸起的冠状动脉区域，该层相对于前一层在冠状动脉根结点处的光流位移距离会存在巨大的变化，通过对比相邻2层间光流位移距离最大值的变化即可初步筛选出冠状动脉的起始层；光流法假设相邻的2层CTA图像强度为I(x,y,t)和I(x+Δx,y+Δy,t+Δt)，并且相邻的2层图像强度不变，可以得到公式(16)
I(x+Δx,y+Δy,t+Δt)＝I(x,y,t) (16)
将公式(16)使用泰勒级数展开后可得：
Figure 6
其中θ2为高次项，大小可以忽略，为此结合公式(16)和(17)可得：
Figure 7
其中
Figure FDA0002606470380000057
和
Figure FDA0002606470380000058
分别为图像相对于X轴和Y轴的空间导数，
Figure FDA0002606470380000059
为时间导数，公式(18)可改写为
(Ix,Iy)(u,v)＝-It (19)
其中，u和v就是图像中每个像素分别沿着X轴和Y轴的光流位移距离，但是只有公式(19)一个约束方程无法求解2个未知变量，所以光流法假设目标像素的周围像素也存在相同的光流位移距离，即在目标像素周围m*m的窗口内，每个像素点都具有相似的光流位移距离u和v；通过增加约束方程来求解光流矢量u和v，如公式(20)所示；
Figure 8
式(19)中(Xn，Yn)表示m*m的窗口内第n个像素点；
通过公式(20)求解所得的u和v可以计算出每个像素点的位移距离，如公式(21)所示；
Figure 10
根据每一层Distances最大的点，筛选出最大的前10层作为候选冠脉起始层，然后结合心脏医学解剖结构先验知识进一步筛选出左右冠状动脉的起始层。
7.如权利要求6所述的一种基于CTA图像的冠状动脉分割方法，其特征在于所述的三个心脏医学解剖结构先验知识具体如下：
第一，升主动脉根部只存在左右两个冠状动脉开口，在10层冠脉候选层中，位移距离最大的是左冠状动脉起始层Sleft；
第二，左右冠脉起始层间相隔10层以上，且左冠状动脉起始层位于右冠状动脉起始层的Z轴上方；
第三，左冠状动脉起始点OL位于该层质心的-90°～45°之间，右冠状动脉起始点OR位于该层质心的位于45°～135°之间；
最后在满足上述第二、三条件的剩余候选冠脉起始层中选择层号最小的层作为右侧冠状动脉起始层Sright。
8.如权利要求6所述的一种基于CTA图像的冠状动脉分割方法，其特征在于所述的步骤5-2具体是：
5-2-1首先从5-1中选取的每一层c*15个种子点依次开始对该层血管区域进行分割，每一个种子点会将区域边缘像素点8领域内与种子点灰度值相差在1％以内的像素点不断归并到感兴趣区域中，直到不能生长为止；
5-2-2由于5-1中选取的种子点只有部分会落在血管区域内，所以每个种子点生长完成后，通过计算分割得到的血管感兴趣区域内的像素数量num判断该种子点是否在血管区域内：如果像素数num小于Nthresh，则表明该种子点位于血管区域内，保留该种子点分割得到的血管感兴趣区域；如果像素数num大于Nthresh，则说明该种子点落在血管区域之外，出现了严重的过分割现象，丢弃该种子点分割得到的区域；然后将所有保留下来的血管感兴趣区域之和作为该层最终的冠状动脉血管分割结果；
5-2-3判断每一层最终的冠状动脉血管分割结果像素数是否为0，如果等于0，则表明无法再继续分割血管，结束迭代，得到初始冠状动脉感兴趣区域；如果不等于0，则返回5-1继续执行下一层种子点的选择以及5-2迭代区域生长法。
9.如权利要求8所述的一种基于CTA图像的冠状动脉分割方法，其特征在于所述的步骤5-3中所述端点P1在该体素的26邻域点中，只有1个邻域点位于骨架线上；所述普通连接点P2，在该体素的26邻域点中，只有2个邻域点属于骨架线上；所述分叉点P3，在该体素的26邻域点中，至少有3个邻域点位于骨架线上。
CN202010740171.3A 2020-07-28 2020-07-28 一种基于cta图像的冠状动脉分割方法 Active CN111951277B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010740171.3A CN111951277B (zh) 2020-07-28 2020-07-28 一种基于cta图像的冠状动脉分割方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010740171.3A CN111951277B (zh) 2020-07-28 2020-07-28 一种基于cta图像的冠状动脉分割方法

Publications (2)

Publication Number Publication Date
CN111951277A true CN111951277A (zh) 2020-11-17
CN111951277B CN111951277B (zh) 2024-03-12

Family

ID=73338329

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010740171.3A Active CN111951277B (zh) 2020-07-28 2020-07-28 一种基于cta图像的冠状动脉分割方法

Country Status (1)

Country Link
CN (1) CN111951277B (zh)

Cited By (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113034451A (zh) * 2021-03-15 2021-06-25 北京医准智能科技有限公司 一种基于深度学习的胸部dr图像识别方法
CN113409333A (zh) * 2021-06-16 2021-09-17 青岛海信医疗设备股份有限公司 一种三维图像的切割方法及电子设备
CN113506277A (zh) * 2021-07-16 2021-10-15 推想医疗科技股份有限公司 图像处理方法及装置
CN113689388A (zh) * 2021-08-03 2021-11-23 慧影医疗科技（北京）有限公司 用于主动脉曲面重建的三维中心线起始点定位方法及装置
CN113888690A (zh) * 2021-10-19 2022-01-04 柏意慧心(杭州)网络科技有限公司 用于确定血管中的目标区段的方法、设备和介质
CN114299081A (zh) * 2021-12-16 2022-04-08 北京朗视仪器股份有限公司 一种上颌窦cbct图像分割方法、装置、存储介质及电子设备
CN114926700A (zh) * 2022-07-22 2022-08-19 浙江大学 冠状动脉类别确定方法、装置、电子设备及存储介质
CN116863013A (zh) * 2023-05-26 2023-10-10 新疆生产建设兵团医院 基于人工智能的扫描图像处理方法及其系统
CN117115186A (zh) * 2023-10-25 2023-11-24 高州市人民医院 基于区域生长的心血管分割方法
CN117197164A (zh) * 2023-11-08 2023-12-08 中国医学科学院北京协和医院 用于计算心肌血管供血区域的管道流域计算方法及系统
CN117274216A (zh) * 2023-10-09 2023-12-22 聆数医疗科技(苏州)有限公司 基于水平集分割的超声下颈动脉斑块检测方法及系统
CN117598782A (zh) * 2023-09-28 2024-02-27 杭州盛星医疗科技有限公司 用于经皮穿刺手术的手术导航方法、装置、设备及介质
CN117598782B (zh) * 2023-09-28 2024-06-04 苏州盛星医疗器械有限公司 用于经皮穿刺手术的手术导航方法、装置、设备及介质

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104978725A (zh) * 2014-04-03 2015-10-14 上海联影医疗科技有限公司 一种冠状动脉分割方法和装置
CN108335304A (zh) * 2018-02-07 2018-07-27 华侨大学 一种腹部ct扫描序列图像的主动脉瘤分割方法
CN110458847A (zh) * 2019-07-05 2019-11-15 心医国际数字医疗系统(大连)有限公司 基于cta影像的自动冠脉分割及中心线提取方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104978725A (zh) * 2014-04-03 2015-10-14 上海联影医疗科技有限公司 一种冠状动脉分割方法和装置
CN108335304A (zh) * 2018-02-07 2018-07-27 华侨大学 一种腹部ct扫描序列图像的主动脉瘤分割方法
CN110458847A (zh) * 2019-07-05 2019-11-15 心医国际数字医疗系统(大连)有限公司 基于cta影像的自动冠脉分割及中心线提取方法

Cited By (21)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113034451A (zh) * 2021-03-15 2021-06-25 北京医准智能科技有限公司 一种基于深度学习的胸部dr图像识别方法
CN113409333A (zh) * 2021-06-16 2021-09-17 青岛海信医疗设备股份有限公司 一种三维图像的切割方法及电子设备
CN113409333B (zh) * 2021-06-16 2022-07-22 青岛海信医疗设备股份有限公司 一种三维图像的切割方法及电子设备
CN113506277A (zh) * 2021-07-16 2021-10-15 推想医疗科技股份有限公司 图像处理方法及装置
CN113689388B (zh) * 2021-08-03 2024-02-06 慧影医疗科技(北京)股份有限公司 用于主动脉曲面重建的三维中心线起始点定位方法及装置
CN113689388A (zh) * 2021-08-03 2021-11-23 慧影医疗科技（北京）有限公司 用于主动脉曲面重建的三维中心线起始点定位方法及装置
CN113888690A (zh) * 2021-10-19 2022-01-04 柏意慧心(杭州)网络科技有限公司 用于确定血管中的目标区段的方法、设备和介质
CN113888690B (zh) * 2021-10-19 2022-08-12 柏意慧心(杭州)网络科技有限公司 用于确定血管中的目标区段的方法、设备和介质
CN114299081A (zh) * 2021-12-16 2022-04-08 北京朗视仪器股份有限公司 一种上颌窦cbct图像分割方法、装置、存储介质及电子设备
CN114299081B (zh) * 2021-12-16 2023-02-17 北京朗视仪器股份有限公司 一种上颌窦cbct图像分割方法、装置、存储介质及电子设备
CN114926700A (zh) * 2022-07-22 2022-08-19 浙江大学 冠状动脉类别确定方法、装置、电子设备及存储介质
CN116863013A (zh) * 2023-05-26 2023-10-10 新疆生产建设兵团医院 基于人工智能的扫描图像处理方法及其系统
CN116863013B (zh) * 2023-05-26 2024-04-12 新疆生产建设兵团医院 基于人工智能的扫描图像处理方法及其系统
CN117598782A (zh) * 2023-09-28 2024-02-27 杭州盛星医疗科技有限公司 用于经皮穿刺手术的手术导航方法、装置、设备及介质
CN117598782B (zh) * 2023-09-28 2024-06-04 苏州盛星医疗器械有限公司 用于经皮穿刺手术的手术导航方法、装置、设备及介质
CN117274216A (zh) * 2023-10-09 2023-12-22 聆数医疗科技(苏州)有限公司 基于水平集分割的超声下颈动脉斑块检测方法及系统
CN117274216B (zh) * 2023-10-09 2024-04-16 聆数医疗科技(苏州)有限公司 基于水平集分割的超声下颈动脉斑块检测方法及系统
CN117115186A (zh) * 2023-10-25 2023-11-24 高州市人民医院 基于区域生长的心血管分割方法
CN117115186B (zh) * 2023-10-25 2024-02-02 高州市人民医院 基于区域生长的心血管分割方法
CN117197164A (zh) * 2023-11-08 2023-12-08 中国医学科学院北京协和医院 用于计算心肌血管供血区域的管道流域计算方法及系统
CN117197164B (zh) * 2023-11-08 2024-03-08 中国医学科学院北京协和医院 用于计算心肌血管供血区域的管道流域计算方法及系统

Also Published As

Publication number Publication date
CN111951277B (zh) 2024-03-12

Similar Documents

Publication Publication Date Title
CN111951277B (zh) 一种基于cta图像的冠状动脉分割方法
CN107563983B (zh) 图像处理方法以及医学成像设备
JP5019825B2 (ja) ディジタル画像データにおけるコントラスト増強された血管の識別方法
CN110490040B (zh) 一种识别dsa冠状动脉图像中局部血管狭窄程度的方法
CN106228561B (zh) 血管提取方法
CN111047607B (zh) 一种自动分割冠状动脉的方法
CN109478327B (zh) 用于在任意视场计算断层摄影血管造影（cta）中自动检测全身动脉的方法
CN110648338B (zh) 图像分割方法、可读存储介质及图像处理装置
CN110223271B (zh) 血管图像的自动水平集分割方法及装置
JP2008521473A (ja) 複数要素血管セグメント化
CN102663824A (zh) 在心脏计算机断层成像容积中的心脏分离的方法和系统
CN109727242B (zh) 血管中心线提取方法、装置、计算机设备和存储介质
CN110910370B (zh) 一种cta图像冠脉狭窄检测的方法及装置
CN110570424B (zh) 一种基于cta动态影像的主动脉瓣半自动分割方法
CN111612743A (zh) 一种基于ct图像的冠状动脉中心线提取方法
CN105279759A (zh) 结合上下信息窄带约束的腹腔主动脉瘤外轮廓分割方法
US9042619B2 (en) Method and system for automatic native and bypass coronary ostia detection in cardiac computed tomography volumes
Vukadinovic et al. Segmentation of the outer vessel wall of the common carotid artery in CTA
Wang et al. Integrating automatic and interactive methods for coronary artery segmentation: let the PACS workstation think ahead
CN114998292A (zh) 一种基于残差双注意力机制的心血管钙化斑块检测系统
Karim et al. Left atrium segmentation for atrial fibrillation ablation
Khaled et al. Automatic fuzzy-based hybrid approach for segmentation and centerline extraction of main coronary arteries
Lavi et al. Single-seeded coronary artery tracking in CT angiography
Ma et al. A novel automatic coronary artery segmentation method based on region growing with annular and spherical sector partition
Garcia et al. Coronary vein extraction in MSCT volumes using minimum cost path and geometrical moments

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant