CN103339653A - 勾画组织损伤的方法 - Google Patents

勾画组织损伤的方法 Download PDF

Info

Publication number
CN103339653A
CN103339653A CN2011800654880A CN201180065488A CN103339653A CN 103339653 A CN103339653 A CN 103339653A CN 2011800654880 A CN2011800654880 A CN 2011800654880A CN 201180065488 A CN201180065488 A CN 201180065488A CN 103339653 A CN103339653 A CN 103339653A
Authority
CN
China
Prior art keywords
image
damage
frisket
algorithm
ttp
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
CN2011800654880A
Other languages
English (en)
Other versions
CN103339653B (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.)
Aarhus Universitet
Original Assignee
Aarhus Universitet
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 Aarhus Universitet filed Critical Aarhus Universitet
Publication of CN103339653A publication Critical patent/CN103339653A/zh
Application granted granted Critical
Publication of CN103339653B publication Critical patent/CN103339653B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/002D [Two Dimensional] image generation
    • G06T11/003Reconstruction from projections, e.g. tomography
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/0002Inspection of images, e.g. flaw detection
    • G06T7/0012Biomedical image inspection
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T1/00General purpose image data processing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/149Segmentation; Edge detection involving deformable models, e.g. active contour models
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/155Segmentation; Edge detection involving morphological operators
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10072Tomographic images
    • G06T2207/10088Magnetic resonance imaging [MRI]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20036Morphological image processing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30004Biomedical image processing
    • G06T2207/30016Brain

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Medical Informatics (AREA)
  • Quality & Reliability (AREA)
  • Radiology & Medical Imaging (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Software Systems (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)
  • Image Processing (AREA)
  • Apparatus For Radiation Diagnosis (AREA)

Abstract

本发明涉及一种用于估算生物组织中的半暗带的半暗带尺寸的度量,其中分别通过灌注加权成像(PWI)和弥散加权成像(DWI)分析获得的第一和第二图像;且其中第一图像的分析包括水平集方法的应用,第二图像的分析包括灰度形态运算的应用。在本发明进一步的实施例中,连通组件标记算法可以应用在第一和第二图像中的任何一个上。本发明进一步涉及***、计算机程序产品和相应方法的使用。

Description

勾画组织损伤的方法
技术领域
本发明涉及一种组织损伤的计算机辅助勾画的方法,尤其涉及一种卒中患者脑部组织损伤的计算机辅助勾画。本发明进一步涉及***、计算机程序产品和相应方法的使用。
背景技术
在病人发生急性脑卒中的情况下,供血的减少引起一部分的脑组织局部遭受缺血,例如遭受氧气和营养物供应不足的病患,这部分被称为灌注损伤。灌注损伤的子部分是不可逆地损伤,这部分被称为核心损伤或弥散损伤。灌注损伤的其余部分,所谓缺血半暗带,如果供血在短时间内恢复,例如几个小时内,可被抢救。在决定每个卒中病人的治疗策略中,识别半暗带的体积和位置很关键。由于治疗方法如使用血栓溶解剂,必须在短时间内开始,快速评估很关键。
核心和灌注损伤的大小和位置典型地是通过人工检查弥散加权和灌注加权磁共振成像(DWI/PWI)形态来确定的。这需要技术精湛的专家介入,例如卒中神经系放射学家。这样的人工确定DWI/PWI不匹配,需要耗费大量的时间和成本。
WO2007/058632A1描述了在计算机辅助下,处理和分析的脑图像以用于卒中诊断或决定治疗决策,其中多个卒中相关的图像是叠加的。
由于它们不能区别损伤体积和损伤之外的体积,这样的自动图像处理方法是不可靠的。
为了自动确定低灌注组织,已尝试使用简单阈值算法。然而,这种方法,问题在于缺乏已建的适用每个人的阈值,如灌注阈值。此外,研究表明阈值取决于许多因素,包括症状出现后的时间和后期处理的方法。
因此,一种改进的计算机辅助勾画组织损伤的方法是有益的,尤其更快和/或更可靠的方法是有益的。
发明内容
本发明的目的是提供一种组织损伤的计算机辅助勾画的方法,该方法解决了上述的现有技术中的方法的不快速和不可靠的问题。
此外,本发明的一个优势是本发明能给出可重现的结果。另一个优势在于实施该方法以便以全自动方式操作。
本发明的进一步的目的是对现有技术提供一个可替代的选择。
因此,上述目标和几个其他目标旨在获得本发明的第一个方面,通过提供估算生物组织中半暗带的半暗带尺寸的方法,该方法包括步骤:
—接收第一图像,所述第一图像包括表示空间分辨灌注参数的数据;
—接收第二图像,所述第二图像包括表示空间分辨弥散参数的数据;
—识别所述第一图像上的第一区域,所述第一区域对应于灌注损伤;
—识别所述第二图像上的第二区域,所述第二区域对应于弥散损伤;
—基于所述第一和第二区域估算半暗带尺寸的度量;
其中所述方法进一步包括:
—在所述灌注加权图像上应用第一图像处理方法,所述第一图像处理方法包括水平集方法;
—在所述弥散加权图像上应用第二图像处理方法,第二图像处理方法包括灰度形态运算。
本发明特别地,但不局限地,在改善速度和/或可再现性和/或极小化用户干预上是有益的。该方法产生具有平滑边界的损伤估算。本发明提供了一种勾画第一区域和第二区域的方法,以此勾画潜在可挽救的组织的位置和区域。本发明的一个优势在于它使得快速估算生物组织中半暗带的半暗带尺寸的度量成为可能。另一个潜在的优势可在于它使得可重现地估算生物组织中半暗带的半暗带尺寸的度量成为可能。另一个优势在于它使得估算生物组织中半暗带的半暗带尺寸的度量成为可能,其中度量与人工提供的估算关联得很好,例如由经过培训的专家人工提供的估算。
半暗带尺寸可以基于第一和第二区域确定,例如半暗带是以第一区域和第三区域之间的差异而确定的,第三区域是包含在第一区域和第二区域中的一个区域。在另一个实施例中,半暗带是以第一区域与第二区域之间的差异而得到的。
图像被理解为一组表示空间分辨参数的数据点,例如一组空间分辨值,其中每个数据点对应于处于某一个位置的参数的值。可以理解,多个位置可以被包含在一个平面内,对应于二维图像,或者这些位置可以跨更多个维度分布,例如三维。可以进一步理解,尽管位置被描述为空间数学上理想的点,每个数据点可以对应一个有限的面积或体积,例如有限的面积或体积被指定给每一个数据点。在本申请中,“图”与“图像”可交换使用。
可以理解,即使在二维图像中识别一个区域的情况下,相应的体积可以被推导出,因为可以理解的是:为二维图像指定一个延拓,例如深度,在与两个图像维度正交的第三空间维度中。换句话说,图像中每一个像素表示一个体素。因此,可以理解,一个面积可以对应于一个体积。
图像可以表示生物组织的切片,例如人的大脑可以用大约20幅图像来表示。
第一图像被理解为包括表示空间分辨灌注参数的数据点的图像,例如一组空间分辨值,其中每个值对应于给定位置的灌注参数的值。第一图像可以是通过灌注加权成像(PWI)获得的图像。PWI在现有技术中是众所周知的。
灌注参数理解为一种度量,该度量量化、表示或相关于通过一个元素,例如生物组织如生物器官的路径,尤其是动脉血到毛细血管的输送过程。
第二图像理解为包括表示空间分辨弥散参数的数据点的图像,例如一组空间分辨值,其中每个值对应于给定位置的弥散参数的值。第二图像可以是通过弥散加权成像(DWI)获得的图像。DWI在现有技术中是众所周知的。
第一图像和/或第二图像可以借助于核磁共振NMR成像提供。可选地,应用其他成像技术例如正电子发射断层扫描(PET)、超音波或X射线计算机断层扫描(CT)以获得第一和/或第二图像。
弥散参数理解为度量,该度量量化、表示或相关于分子在它们的微环境中的机动性,例如水分子的动态移动。
第一区域理解为第一图像中对应于灌注损伤的区域。
灌注损伤理解为一个严重缺血、功能受损具有梗塞风险的组织区域(或复数个区域)。
第二区域理解为第二图像中对应于弥散损伤的区域。
弥散损伤理解为一个(或复数个)组织区域,表示受到不可逆损伤的组织。
损伤理解为对生物组织的急性损害。梗塞理解为生物组织的永久损害。
半暗带理解为包含在第一区域内,但不包含在第二区域内的区域。因此半暗带区域表示严重缺血、功能受损具有梗塞风险的组织;如果在不可逆损伤之前再灌注,半暗带区域会被挽救;否则半暗带区域逐步增添到弥散损伤中,直到达到最大梗塞扩展。
半暗带尺寸理解为半暗带的大小,例如,如果半暗带以对应于图像中若干数据点的区域被给出,半暗带的尺寸可以对应于指定给这些数据点中的每一个的大小的和的大小被给出。
图像处理方法理解为一种方法,其接收输入图像,例如第一图像或第二图像,作为输入,并输出一个新的图像,例如处理过的图像。处理过的图像可不同于输入图像,不同之处在于一个或更多个值被改变和/或一个或更多个位置被改变。图像处理方法可以把单个图像作为输入,输出单个图像;然而图像处理方法也可以把复数个图像作为输入,输出一个或更多个图像。例如,图像处理方法可以涉及3D核心程序,因此包括来自其他图像的信息,例如表示生物组织中相邻的切片的图像。
水平集方法(LSM)理解为一种方法,其隐式地通过一个点集表示一个图像中区域间的边界,该点集为函数具有固定值的点集一个定义在图像上的被称为水平集的函数具有固定值,该函数定义在图像上,被称为水平集函数。这可以是零水平集。此外,该方法可包括描述随水平集函数时间变化的微分方程。LSM可用作极小化能量泛函的迭代技术,这样图像被优化地分割为若干个区域。例如,我们可以设法将一个图像分成两个区域,一种用于此目的的方法包括:
识别将图像分成两个区域的曲线,其中曲线或表面的位置和形状可对应于极小化能量泛函的位置和形状,能量泛函表示对应于差值的和的误差,该差值为在两个区域中的每一个区域中的每一个数据点的值分别与每个区域中的值的平均值之间的差值。
水平集方法通常是根据微分方程(最初Jacobi-Hamilton类型)用于演化曲线的过程。演化的“驱动力”或者是本身等值线(例如局部曲率)的内在属性,或者是图像分割的值,可方便地依赖局部图像属性。力的一个选择可以是(相反)图像梯度,例如,在物体边界演化放缓。这要求图像梯度的计算,或者更普遍地,边缘检测(该领域中已知的例子,除图像梯度之外,包括Canny边缘检测,Sobel边缘检测)。在一些图像中,例如在MTT图像中,常常在缺血和正常组织之间,只有一个模糊的梯度,这使得边缘检测效率很低。在一个实施例中,为了规避(Mumford-Shah类型:异质性极小化)直接引用边缘,使能“模糊”边界的识别(有时称为认知边界),使用变分公式。
形态灰度开运算和闭运算在现有技术中是众所周知的,指的是灰度腐蚀
Figure BDA00003530197200043
和膨胀⊕的组合。
形态灰度重构理解为开重构运算和闭重构运算的任何组合。重构是在原始图像中标记图像的连续(灰度)膨胀或腐蚀。设h1=标记图像:
hk+1=(hk⊕B)∩g,
直到稳定,
hk+1=hk
其中B是结构元素,g是原始图像。
在本发明的一个实施例中,标记图像可以通过灰度腐蚀或灰度膨胀获得,或者从/向输入图像减去或加上一个常数。
形态灰度运算可以是灰度腐蚀、膨胀、开运算、闭运算、开重构运算和闭重构运算的任何组合。
形态灰度重构是有益的,因为它独特地截去图像中的高峰,如具有特别高的数据值的像素,而同时易于保留相关区域原始边界,正如由内核大小决定的。形态灰度重构描述在“图像分析中的形态灰度重构:应用软件及高效算法”,文森特,L,IEEE图像处理汇刊,第2卷,第2部,第176-201页,1993年4月,通过引用包含于本文中。
达峰时间(TTP)图像在现有技术中众所周知,TTP图像理解为一个图像,该图像中空间分辨参数对应于时间间隔的长度,该时间间隔为从开始时间到一个对应于在给定位置测量的组织浓度曲线的最大值的时间。
平均通过时间(MTT)图像在现有技术中众所周知,MTT图像理解为一个图像,该图像中空间分辨参数对应于流体,如血液,通过生物组织的毛细血管的平均前置时间。这可以由用动脉输入函数去卷积的组织浓度曲线的最大高度与脑血容量的比确定,脑血容量用组织浓度曲线下的面积确定。
连通区域标记理解为一种方法,其中为了识别数据点与至少一个其他数据点在值和位置方面连通而处理图像。CCL在参考文献C算法中描述,第3版,塞奇威克,罗伯特,艾迪生-韦斯利,1998,第11-20页,通过引用包含于本文中。
蒙片理解为一组勾画图像中一个区域的数据点。在当前的上下文中,蒙片理解为包括关于勾画轮廓的信息,包括蒙片是一组表示一条线的数据,这条线本身表示轮廓,或蒙片是一组表示勾画的区域的数据,或者一组表示包围区域的数据。进一步包括蒙片是描述勾画边界或勾画的区域或体积的数学表达式。
开运算、闭运算或填充理解为众所周知的图像处理方法,例如形态二值和/或灰度形态运算的各种组合,例如腐蚀和膨胀。
正则化理解为对图像强度强加的罚项,以及随后基于在图像强度上强加的罚项的修正。这可产生在图像强度上的软约束。
种子点理解为包含在损伤内的一个数据点。例如,第一图像中的种子点理解为包含在灌注损伤内的一个数据点。例如,第二图像中的种子点理解为包含在弥散损伤内的一个数据点。进一步理解为,种子点可以是定义为被包含在损伤中的一个数据点。
表观弥散系数(ADC)图像在现有技术中众所周知,ADC图像理解为一个图像,在图像中空间分辨参数对应于一种度量,这种度量量化、表示或涉及分子在它们的微环境中的机动性,尤其水分子的动态移动。ADC图像强度是绝对标度。
CSF理解为脑脊液,脑脊液是占据了大脑中脑室***和脊髓的液体。
动态磁敏感对比(DSC)MRI(核磁共振成像)在现有技术中众所周知,指的是对非弥散、顺磁性造影剂的小圆团通过大脑的路径的度量,典型地使用T2/T2*加权MRI。
脑血流量(CBF)在现有技术中众所周知,指的是血液输送到组织的速度。在DSC MRI中,典型地为每个体积元素(体素)通过去卷积的组织曲线的最大函数值计算CBF。
脑血量(CBV)在现有技术中众所周知,指的是组织区域中的血液的容积率。在DSC MRI中,CBV通常以造影剂浓度曲线下的面积计算,且由动脉输入函数下的面积规格化。
弥散加权成像在现有技术中众所周知,指的是对分子在它们的微环境中机动性的度量。
弥散加权图像(DWI)在现有技术中众所周知,指的是用不同实验方法的、所谓b因子获得的图像,其依赖于应用的梯度脉冲振幅和时序,具有向弥散特性(弥散加权)改变图像灵敏度的效果。
表观弥散系数(ADC)在现有技术中众所周知,指的是弥散系数的度量,对照简单获得弥散加权图像。典型地使用若干具有不同b因子的弥散加权图像获得ADC。
MCA理解为源于颈内动脉的大脑中动脉,向额骨、颅顶骨和颞叶提供血液供给,通常被称为MCA区域。
STM理解为根据预先确定的阈值,将灰度或RGB图像转换为二值图像的图像处理技术。
结构图理解为区别(通过图像强度)大脑中不同组织成分的图像,例如白质、灰质骨骼。
体素在现有技术中众所周知,可理解为表示给定体积的实体,例如生物组织内的体积。体素可以被指定一个参数的值,例如灌注参数或弥散参数。一个或更多个体素可以构成图像。
“算法”被理解为对应于“方法”,并因此与“方法”可以互换使用。可以进一步理解,算法可以在计算机程序产品中实现,例如适合实现一种特别的方法的程序。
在本发明的一个实施例中,提供了一种方法,该方法进一步包括在从图像组中选择的任何一个图像上应用连通区域标记,该图像组包括第一图像和第二图像。使用CCL的优势可在于相比于例如阈值转换对噪声更不敏感。STM不限制在某个体素群集中根据种子点的搜索,但认为所有超过某个阈值的体素是损伤体素,不管它们的解剖学位置或关于它们的实际损伤位置的邻近位置。
在本发明的另一个实施例中,提供了一种方法,其中该方法进一步包括在从图像组中选择的任何一个图像上应用形态二值重构,该图像组包括第一图像和第二图像。形态二值重构可适用于确定哪些体素被连接到种子点。使用形态重构的优势可在于相比于例如阈值转换对噪声更不敏感。STM不限制在某个像素群集中根据种子点的搜索,但认为所有超过某个阈值的体素是损伤体素,不管它们的解剖学位置或关于它们的实际损伤位置的邻近位置。形态重构的另一个可的优势是它是快速的。
“形态二值重构”理解为涉及两个图像和一个结构元素的形态转换。一个图像(标记)是转换的起点,另一个图像(蒙片),限制转换。使用的结构元素定义连通性。为提取连接部件目的的形态二值重构,在参考文献“图像分析中的形态灰度重构:应用程序和高效算法”中进一步描述,文森特,L.,IEEE图像处理汇刊,1993.2(2):第176-201页,完整地通过引用包含于本文中。
在本发明的一个实施例中,提供了一种方法,其中第二图像处理方法进一步包括接收由第一图像处理方法确定的参数,且其中识别第二区域的步骤基于该参数。这可在用于某些图像形态的自动分割算法上是有益的,例如MTT;DWI可以方便地把一个参数当作输入,该参数从另一个图像形态的分析中获得,例如TTP图像。参数可以是感兴趣的区域或蒙片。接收这样的参数对于减少运行算法需要的时间和避免假阳性损伤是有益的。
在本发明的另一个实施例中,提供了一种方法,其中第一图像处理方法包括:
-在初级第一图像上识别初级蒙片;
-基于次级第一图像和初级蒙片,在所述第一图像上识别第一区域,作为次级第一图像上的一个区域;
其中所述初级第一图像没有对于动脉输入函数被修正,而所述次级第一图像对于动脉输入函数被修正。
然而,在本发明的另一个实施例中,提供了一种方法,其中步骤:
-在所述初级第一图像上识别初级蒙片;
包括在所述初级第一图像上应用水平集方法。
在本发明的另一个实施例中,提供了一种方法,其中第一图像处理方法包括在初级第一图像上识别初级蒙片,且其中在第一图像上识别第一区域的步骤,把次级第一图像和初级蒙片作为输入,其中初级和次级第一图像是不相似的,次级第一图像已基于动脉输入函数被处理过。可以理解,初级的第一图像可以是未对动脉输入函数修正的第一图像,而次级第一图像可以是对动脉输入函数修正的第一图像。这样做的一个优势在于它有助于分割过程中的步骤的快速收敛,例如水平集迭代。另一个优势可在于极小化噪声的影响,这种影响典型地在次级第一图像中,比在初级第一图像中更明显。在特定的实施例中,初级第一图像是TTP图像。在特定的实施例中,次级第一图像是MTT图像。TTP损伤可以被认为灌注损伤的合理的最初估算,在这个“第一次猜测”时,寻找非正常MTT代表的损伤应因此被启动。非知情的MTT算法初始化可造成轮廓被吸引到远离感兴趣区域的人工高信号区域的风险。
在本发明的另一个实施例中,提供了一种方法,其中第一图像处理方法进一步包括应用灰度形态运算。应用灰度形态运算的优势可在于图像中的噪声可被降低,而损伤边缘可被保留。
在本发明的另一个实施例中,提供了一种方法,其中第二图像处理方法进一步包括应用水平集方法。ADC损伤和DWI损伤不必在空间上匹配,因此使用水平集方法调整DWI损伤边界是有益的。
在本发明的另一个实施例中,提供了一种方法,其中第一图像处理方法的水平集方法是Mumford-Shah分割的水平集实现。这个实施例可以自适应地将一个图像分割成两个具有不同平均图像强度,例如参数的值的区域,这样在每个区域中的这些值,从紧紧围绕各自的平均值变化的意义上来说是(最佳)同质的。因此Mumford-Shah分割包括异质性减到最小程度。在另一种实施例中,曲率罚项确保边界是平滑曲线。这个实施例的优势在于它符合损伤体积的专家轮廓。上面提出的问题,例如,将图像分割成两个具有不同平均图像强度,例如参数的值的区域的问题,这样在每个区域中的这些值,从紧紧围绕各自的平均值变化的意义上来说是(最佳)同质的,可与成本函数(或能量)极小化有关,该成本函数由两个区域中覆盖所有可轮廓的强度均匀度的总和具体表达。这实际上是一个变分法的问题,正如已经指出可以通过根据水平集函数重新计算该问题和计算相关的Euler-Lagrange方程,有效的迭代解法可以用公式表示。在特别的实施例中,提供了一种方法,其中第一图像处理方法的水平集方法是Chan-Vese模型的实现。Chan-Vese模型在现有技术中众所周知,且可被看作Mumford-Shah模型水平集实现的一种简化,其中分段常数函数是估计的。
人工延长平均通过时间的区域,例如脑室,是全自动勾画MTT损伤的一个挑战。在另一个实施例中,在没有额外的以用户介入或速度的形式成本的情况下,并发DWI图像可以被用来作为方便的排除CSF体素的方法。作为一种选择,脑室可以在最初的、原始的T2加权DSC图像上被识别,或者通过使用CBV标准被识别。
在本发明的另一个实施例中,提供了一种方法,该方法进一步包括在从第一图像和第二图像选择的任何一个图像上应用二值形态运算。二值形态运算在现有技术中是众所周知的。应用二值形态运算的一个优势可在于,在第一图像上识别第一区域的步骤,和/或在第二图像上识别第二区域的步骤可以改善,例如执行更快,例如导致关于第一和/或第二区域的专家识别量化更好的结果。然而,在另一个实施例中,该方法进一步包括应用二值形态重构。二值形态重构在现有技术中是众所周知的。
二值形态运算的实例可以包括:填充损伤内由信号脱离引起的“小孔”,删除显示与损伤连通的小伪影,以及平滑损伤边界,例如形态开运算和闭运算。为了调解边界不规则和/或减小假阳性损伤的数量,应用二值形态运算可以是有益的。
在本发明的另一个实施例中,提供了一种方法,其中二值形态运算包括从开运算、闭运算、填充中选择的任何一项图像处理技术。
在本发明的另一个实施例中,提供了一种方法,该方法进一步包括在从第一图像和第二图像组成的组中选择的任何一个图像上识别种子点。优势是通过结合CCL与确定种子点,例如通过自动或人工确定,使得图像分割能够不包括假阳性对象。通过在CCL中实施来自区域增长的种子点特征,获得一种快速和更低噪声的TTP损伤分割。
在本发明的另一个实施例中,提供了一种方法,其中识别种子点的步骤包括在图像中识别主要点,主要点对应图像中的一个点,该点关于图像中的其他点具有更高或更低的数据值。这个实施例的优势在于它使种子点自动确定成为可能。
在本发明的另一个实施例中,提供了一种方法,其中识别种子点的步骤包括损伤侧性的检测,例如损伤侧性的检测包括:在单成像形态中为每个半球确定平均体素值;例如损伤侧性的检测包括:从多个不同成像形态中为每个半球提取一个或更多个特征;例如不同形态是TTP图像和ADC图像。种子点可被限制在损伤侧性中。在另一个特别的实施例中,提供一种方法,其中识别种子点的步骤包括损伤侧性的检测,其中损伤侧性的确定包括:基于第一成像形态为每个半球确定第一特征的步骤,和基于第二成像形态为每个半球确定第二特征的步骤;第二成像形态不同于第一成像形态。在另一个特别的实施例中,损伤侧性可以基于一个或更多个特征被确定,例如在过滤后的TTP图上任何一个半球上提取的特征,例如特征:1)左和右半球的平均强度值;2)在TTP图上最大强度体素的侧性;以及3)具有强度<600*10-6mm2/sec的ADC体素的数量。对于每一个特征,“获胜的”半球被指定值1,最终具有至少两个正数结果的半球可以被认为是同侧半球。从多个不同成像形态中为每个半球确定一个或更多个特征的优势,可在于损伤侧性的确定因而以更高的确定性完成,例如选择错误半球的风险更低。因此,当限定种子点在确定的半球内,种子点在正确的半球内的确定性更大,例如具有更小的选择错误半球的风险。
根据本发明的第二个方面,本发明进一步涉及一个***,该***估算受测对象的半暗带的半暗带尺寸的度量,该***包括:
—用于接收第一图像的输入,第一图像包括表示空间分辨灌注参数的数据;
—用于接收第二图像的输入,第二图像包括表示空间分辨弥散参数的数据;
—图像处理单元设置为在第一图像上识别第一区域,第一区域对应于灌注损伤;
—图像处理单元设置为在第二图像上识别第二区域,第二区域对应于弥散损伤;
—处理器设置为估算半暗带尺寸的度量,作为第一区域和第三区域之间的区别,第三区域是包含在第一区域和第二区域两个区域中的一个区域;
其中该方法进一步包括:
—图像处理单元设置为在灌注加权图像上应用第一图像处理方法,第一图像处理方法包括水平集方法;
—图像处理单元设置为在弥散加权图像上应用第二图像处理方法,第二图像处理方法包括灰度形态运算。
在另一个实施例中,提供一个***,该***估算受测对象的半暗带的半暗带尺寸的度量,该***包括成像单元,例如核磁共振扫描器。
根据本发明的第三个方面,本发明进一步涉及计算机程序产品,该计算机程序产品能够实现根据第一方面的方法。通过使用计算机程序产品,所述方法可以在最少的或无用户介入情况下,又快又可靠地实施。
根据本发明的第四个方面,本发明进一步涉及根据第一方面的方法的使用,所述方法用于估算生物组织中的半暗带的半暗带尺寸的度量,其中生物组织是脑组织。
本发明的第一、第二、第三和第四方面的每一个可以与其他方面中的任何一个组合。本发明的这些和其他方面可从通过以下相关具体实施例的阐明而显而易见。
附图说明
现在结合附图,更详细地描述根据本发明的计算机辅助勾画组织损伤的方法。附图示出实现本发明的一种方式,不应将其理解为对权利要求的范围内的其他可能的实施方式限制。
图1显示了根据本发明的一个实施例的流程图;
图2显示了最佳阈值时根据例1的算法、STM和人工勾画轮廓的体积对比;
图3显示了TTP>4秒时根据例1的算法与专家相比(斯皮尔曼相关R2=0.87)低估损伤体积;
图4显示了根据例1的算法与STM的Dice系数的箱形图;
图5显示了通过根据例1的算法和STM产生的分割的实例;
图6:以若干迭代获得的估算的损伤边界的实例;
图7比较了水平集损伤体积与专家损伤体积;
图8:由根据例2的算法估算的损伤边界的实例,以及对于六个具有典型损伤的患者的专家间一致性;
图9显示了四个专家共识区域与水平集算法之间的一致性;
图10显示了四个专家共识区域与阈值转换之间的一致性;
图11显示DWI图的预处理;
图12显示DWI阈值的确定;
图13:通过不同方法估算的DWI损伤的体积比较;
图14:通过不同方法估算的半暗带的体积比较;
图15:分类性能;
图16:通过不同方法勾画的DWI损伤轮廓的实例;
图17:通过不同方法勾画的半暗带轮廓的实例;
图18:TTP损伤的体积比较;
图19:Dice系数的箱形图;
图20:分割实例;
图21:初始蒙片到包围DWI损伤的最终蒙片的演化;
图22:DWI损伤的体积比较;
图23:半暗带的体积比较;
图24:分类性能;
图25:半暗带实例;
图26:TTP损伤分割中的步骤的示意图;
图27:TTP损伤的体积比较;
图28:Dice系数的箱形图;
图29:通过APS自动生成的TTP损伤蒙片的实例;
图30:DWI损伤分割算法中步骤的示意图;
图31:通过自动方法和人工勾画轮廓估算的DWI损伤的体积比较;
图32:通过自动方法和人工勾画轮廓估算的半暗带的体积比较;
图33:总结半暗带分类的性能度量。
具体实施方式
图1显示了根据本发明的一个实施例的方法的流程图。图1显示了接收输入100的步骤,例如从图形用户界面(GUI)接收输入,向控制该过程的处理器提供输入。图形用户界面可以是主图形用户,它控制分割过程,例如第一图像和第二图像和/或它们的子类型的分割,它包括与损伤体积相关的、由分割过程返回的数据。由于输入100,一个或更多个分割过程(例如TTP、MTT、DWI和FU类型的图像的分割)可单独启动,如虚线箭头102所示,或者顺序启动(例如以时间顺序TTP、MTT和DWI),如箭头101、103、104所示。在图1中,每一列中的步骤表示为某个图像类型执行的步骤。箭头103、104表明用于一种图像类型的方法步骤中获得的数据,可以用作用于另一个图像类型的方法步骤的输入。
对于TTP图像,方法步骤包括确定112种子点,应用114连通区域标记和水平集方法用于产生TTP图像上的蒙片。所述方法也包括由用户的调整116,例如通过GUI,和在先步骤112、114、116中获得的体积发送118给处理器。
对于MTT图像,方法步骤包括在TTP图像中确定122种子点,在TTP图像上的应用124种子点确定和连通区域标记以及水平集方法,用以产生TTP图像上的蒙片。TTP蒙片于是用作在MTT图像上应用125实现水平集方法的算法的输入。该方法也包括由用户调整126,例如通过GUI,和在先步骤122、124、125、126中获得的统计数据发送128给处理器。
对于DWI图像,方法步骤包括应用133规则化和灰度形态运算;在TTP图像上应用134种子点确定和连通区域标记以及水平集方法,以生成TTP蒙片;确定135DWI阈值;应用137连通区域标记以获得DWI蒙片。所述方法也包括由用户调整136,例如通过GUI;和在先步骤133、134、135、137、136中获得的统计数据发送138给处理器。
如以上描述的,每个图像类型可以单独地被分析,然而,为了估算生物组织中半暗带的半暗带尺寸的度量,有益的是:循序地分析图像,并使用从一个图像类型获得的统计数据,例如蒙片作为另一个图像类型分析的输入,例如如箭头103所示,与TTP蒙片相关的统计数据,以及用户的调整被提供给在MTT图像上应用水平集方法的步骤作为输入。
根据这个实施例的方法还具有特征:获得后续蒙片的步骤,例如卒中后相对长的时间段获得的图像上的蒙片,例如3个月,其表示永久的梗塞。这包括步骤:在结构图中应用143规则化和灰度形态运算;在TTP图像上应用144CCL和水平集方法;使用TTP蒙片以确定145结构图中的阈值;在结构图上应用147连通区域标记。所述方法也包括由用户调整146,例如通过GUI,和在先步骤143、144、145、147、146中获得的统计数据发送148给处理器。
实例
例1-包括TTP数据的图像的分割
材料与方法
形态灰度重构
在灌注损伤内寻找一个体素,所谓种子点,用眼睛容易识别;然而搜索图像强度的最大值经常导致假阳性种子点,例如对应于无关区域的图像区域中的高强度点,例如眼睛、脑脊液(CSF)、大脑外强度等。为自动检测种子点,我们使用形态灰度重构,以截去小刺突,同时保留TTP损伤的边缘,并均质化TTP图像的背景。与常规平滑相比,例如使用高斯核,形态灰度重构对刺突敏感度更低,由于这个原因更适合于种子点检测。在当前的上下文中,刺突理解为与损伤相比相对较小的区域,且刺突包括具有极端值的数据点,例如极端值不同于图像中相同方向的平均值,例如极性,作为损伤的数据点。换句话说,如果损伤包括亮像素,刺突可被看作亮像素的小区域。
连通区域标记
连通区域标记(CCL)算法检测图像,以寻找像素群集,其中每个群集中的像素被连接到给定的连通路径。在一个典型的实施例中,CCL算法包括四步:为简单起见,考虑具有4-连通搜索路径的二值图像。步骤1:该算法从左到右遍历所有像素;步骤2:当它到达一个强度=1的像素(p)时,它检查已经遇到的邻近的像素(对于4-连通性[-1,0;0,1]):i)如果两个邻近的像素=0,给p一个新的标记;ii)如果邻近的像素中的一个=1,指定那个标记给p;iii)如果两个邻近的像素=1,指定其中的一种标记给p,且记录几个像素等值;步骤3:等值对在等值群集中排序;步骤4:CCL算法遍历所有指定新标号的像素。在这项研究中,我们在三维中执行CCL,且使用6-连通路径。
根据例1的算法
自动生成的CSF蒙片删除DWI图中的脑室,CSF蒙片通过使用群集分析获得。k-均值群集的输入图像,是由结构图(T2-加权)和由弥散加权成像(表观弥散系数,ADC)获得的图像组合的图像。向PWI图共同注册MNI(蒙特利尔神经学研究所)模板,自动生成左和右半球蒙片。非脑组织(例如眼睛、骨骼、空气等)和静脉(横窦和乙状窦)在PWI体积中的劣质切片中,对噪声成分有贡献,并增大假阳性率。因此,我们通过分析每一切片中的强度变化,删除劣质切片。
TTP图通过用于损伤侧性和种子点自动检测的形态灰度重构过滤。通过在过滤的TTP图上估算任何一个半球上的平均强度,确定TTP损伤的侧性。种子点定义为在过滤的TTP图中具有最高强度的体素。
未过滤的TTP图的值,在对侧的半球中关于平均TTP延时规格化。对侧半球,是未受到影响的一侧,具有“正常”达峰时间。同侧TTP图设置阈值>4秒,且连通区域通过CCL算法找到。连通区域,包含种子点,被认为是TTP损伤。最后,小的对象(<5mL)从损伤蒙片中删除,小孔被填充且损伤边界被平滑。
根据例1的TTP分割算法,在
Figure BDA00003530197200131
R2010a中,由
Figure BDA00003530197200132
实现,但本领域技术人员可以理解根据例1的分割算法也可以用其他编程语言实现。
可以理解,根据例1的算法可以用相应的方法代替。
患者和图像获取
这项研究包括119名患者(♀=49),且他们准许进入四个不同国家(英国、法国、德国、西班牙和丹麦)的医院。这些医院是I-KNOW联盟的成员。患者的中值年龄是71岁,从症状出现到最初的MRI扫描的中值时间是148分钟(附录1中的表I)。
表I描述了淤血患者特征,范围表示第一和第三四分位数。
急性TTP损伤中值是107mL,国家卫生研究院的卒中尺度(NIHSS)中值是11。I-KNOW协议定义了七个卒中子类,其中心源性脑梗塞是主要的卒中子类,随后是带有明显颈动脉狭窄的大血管疾病和未确定的子类。标准动态磁敏感对比MRI在不同型号的扫描器上完成(GE Signa Excite1.5、GE Signa Excite3、GE Signa HDx1.5、Siemens TrioTim 3、Siemens Avanto 1.5、Siemens Sonata 1.5、PhilipsGyroscan NT 1.5以及Philips Intera 1.5)。
根据例1的算法与人工勾画TTP损伤轮廓的比较
具有丰富的临床工作经验的神经系放射学家,使用内部开发的软件半自动勾画TTP损伤轮廓。软件使用标准阈值转换方法,其大致确定TTP损伤,随后神经系放射学家可以人工调整损伤的边界,神经系放射学家看不到任何强度数值范围和其他脑图。我们进行优化试验,通过极小化整个数据集的平方误差的总和,为当前的神经系放射学家估算最佳TTP阈值。最后,TTP损伤>4秒也使用标准阈值方法STM估算。为量化根据例1的算法/STM和人工勾画的损伤之间的几何和体积相似度,我们使用Dice系数,Dice系数是0和1之间的比值,1值表示完全重叠,0值表示无重叠。敏感性和特异性在同侧半球中确定。
结果
体积比较
根据例1的算法的最佳性能,在TTP>2.8秒时获得。人工勾画轮廓与根据例1方法的算法获得的TTP损伤体积之间的斯皮尔曼相关系数是R2=0.89,参见图2A。
图2显示了在最佳阈值时,依据例1的算法、STM和人工勾画轮廓的体积对比;
图2A:当前神经系放射学家的最佳阈值估计为>2.8秒(R2=0.89)。图2A显示了根据例1的算法的数据,第一个轴表示人工确定的蒙片体积[mL],而第二个轴表示通过根据例1的算法确定的损伤体积[mL]。
图2B:与人工勾画轮廓相比,STM高估TTP损伤(R2=0.7)。图2B显示了STM方法的数据,第一个轴表示人工确定的蒙片体积[mL],而第二个轴表示通过STM方法确定的损伤体积[mL]。
图2C:人工方法与根据例1的算法相比的Bland-Altman图以及图2D:人工方法与STM相比的Bland-Altman图。实线表示两个比较的方法之间的平均值,而虚线表示平均差±2标准差(SD)。第一个轴表示平均损伤体积[mL],第二个轴表示人工确定的体积和根据例1的算法和STM方法确定的体积之间的差值[mL],分别在图2C和图2D中。
在最佳TTP阈值转换,人工勾画轮廓与根据例1的算法之间的平均差为1.6mL(±2SD:74.5mL),参见图2C。
在TTP>4秒时,与人工勾画轮廓相比,根据例1的算法略微低估损伤体积。在TTP>4秒时,使用人工和根据例1的算法获得的TTP损伤体积之间的相关是R2=0.87,且在TTP>4秒时人工勾画轮廓和根据例1的算法之间的平均差是35.3mL(±2SD:76.1mL),参见图3。
图3显示了TTP>4秒时根据例1的算法与专家相比(斯皮尔曼相关R2=0.87)低估损伤体积。
STM在TTP>2.8秒时,与人工勾画轮廓相比,高估了TTP损伤体积,参见图2B。人工勾画轮廓与STM之间的平均差为-47.5mL(±2SD:115mL),参见图2D,且相关是R2=0.70。
图5显示了通过根据例1的算法和STM产生的分割的实例。在上面的行中TTP图像和随后行中表明:a)人工蒙片;b)根据例1的算法获得的蒙片;以及c)在TTP图上覆盖STM获得的蒙片(白色)。两个非人工图在TTP>2.8秒时获得。通常,与阈值转换方法相比,根据例1的算法估算边界为平滑曲线。此外,根据例1的算法能够很大程度避免假阳性损伤(图5)。
分类性能
在最佳TTP阈值,根据例1的算法的中值灵敏度(中值:80%和25th-75th四分位数:64-88%)没有显著(p=0.22)高于STM的中值灵敏度(中值:74%和25th-75th四分位数:57-83),然而,根据例1的算法和STM(p<0.001)的特异性之间存在明显差异。对于根据例1的算法的中值特异性是95%(25th-75th四分位数:92-98%)和STM的中值特异性是91%(25th-75th四分位数:88-95%),参见附录1中的表II。
表II描述了在最佳TTP阈值时,根据例1的算法和STM的灵敏度和特异性,用百分比表示。括号中给出了第一和第三四分位数,另外也给出了秩和检验的统计数据。
此外,Dice系数分析表明由根据例1的算法产生的蒙片(中值Dice系数:0.73),与通过STM生成的蒙片(中值Dice系数:0.57)相比,具有明显(p<0.001)更高的与人工勾画的TTP损伤蒙片相似度,参见图4。
图4显示了根据例1的算法与STM的Dice系数的箱形图。
例2-类型MTT数据的图像的分割
材料和方法
水平集
这个方法从假设中出现,该假设为平均MTT值,Mhypo,在低灌注的区域中,比比正常MTT值的平均值Mnorm更高。在低灌注的区域中,MTT值在Mhypo附近变化,而在正常的组织中,它们接近Mnorm。这对应于图像MTT值之间低平方误差,MTT(x,y),和各自的平均值。识别损伤,例如缺血性损伤,于是可以被表述为寻找一条平滑、封闭曲线C,其极小化总变差(等式(1))
E=∫Inside C(MTT(x,y)-Mhypo)2dxdy+∫Outside C(MTT(x,y)-Mnorm)2dxdy
由于缺血区域不一定是连在一起的,例如在闭塞引起分水岭区的损伤的情况下,C可以表示多条曲线是很重要的。这可以通过隐式定义C为一个函数φ的零水平集来保证,其给平面中的每一个点指定一个值,C={φ(x,y)=0}。该函数φ被称为水平集函数。
如果我们假设H(z)表示Heaviside函数,除了z≥0时该函数为零,且它的导函数δ(z)=dH(z)/dz,那么(1)可以写成(等式(2))
E=λ1∫(MTT(x,y)-Mhypo)2H(φ)dxdy+λ2∫(MTT(x,y)-Mnorm)2(1-H(φ))dxdy
+α∫δ(φ)|▽φ|dxdy
最后一项测量轮廓的长度,实际上控制平滑度。我们也可以引入参数λ1和λ2,控制两个误差项的权重。可以表明E可以通过迭代求解微分方程极小化(等式(3))
d&phi; dt = | &dtri; &phi; | [ &alpha; &dtri; &CenterDot; ( &dtri; &phi; | &dtri; &phi; | ) + &lambda; 1 ( MTT ( x , y ) - M hypo ) 2 - &lambda; 2 ( MTT ( x , y ) - M norm ) 2 ]
因此,通过演化根据等式(3)的最初的轮廓,我们得到将MTT图像分成两个同质区域的平滑边界。Mhypo和Mnorm的值,在每次迭代中作为轮廓内和轮廓外的平均值,被自适应地估算。
根据例2的算法
为极小化个体之间在总体图像强度上的变化,在应用水平集过程之前,图像必须被标准化。在这个方法中,我们勾画了一个小的对侧的正常显示白质的区域,减去相应的所有体素的平均MTT值。为了避免假的脑室的轮廓(一些后处理算法生成具有对应于CSF体素的长MTT的区域),基于获得的DWI作为卒中协议的部分,假的脑室的轮廓被自动分割。遵循PWI和DWI图像的线性配准,脑室体素通过它们的表观弥散系数(ADC)被识别,这样具有超过系数为ADCV=2.9的平均正常白质ADC的ADC值的体素,在规格化的MTT图像中设置为0。用于根据例2的水平集算法的最初损伤估计,定义为具有小于ADCV的相对ADC和超过1.2的相对MTT的区域。
最初评估之后,等式(3-4)中的加权参数固定为λ1=1.0,λ2=1.25和α=0.5。这对应于与损伤中相比,对于正常组织中MTT异质性稍高罚项。灵敏度在加权参数中的变化描述如下(也参见图6)。
图6显示了估算的损伤边界的实例,所述估算的损伤边界在按照等式(3)的演化过程中,对于平滑参数α和异质性罚项λ2不同值,在若干迭代次数下获得的。首行显示了随用于这项研究的参数的演化。在初始化时,规格化的MTT图像是经过阈值转换的,在主损伤内产生小孔,并识别出许多小的“假的”损伤。使用该研究参数,这些在10次迭代和100次迭代后解决,轮廓与损伤接近得很好。600次迭代后仅观察到较小的轮廓变化,这表明了收敛性。应用较低的平滑度(α=0.2),估算的轮廓对于图像中沿损伤边界的不规则更敏感,然而边界长度(α=0.8)的进一步限制导致更多的圆形轮廓。减小λ2在非损伤体素中允许更大的异质性,一些高强度体素可因此被包含在这个组中,导致更小的损伤估算(λ2=1.00)。相反地,高的值(λ2=1.50)导致更大和更多异质的损伤估算。
在每一切片中完成总数为600次的迭代。所有算法由MathWorks
Figure BDA00003530197200171
在Matlab
Figure BDA00003530197200172
R2007a中实现。等式(3)的实现改编自www.shawnlankton.com/2007/05/active-contours/。详细资料可以在参考文献中找到。参考文献“无边缘的活动轮廓”,Chan TF,Vese LA著,IEEE图像处理汇刊,第10卷,第2期,第266-277页,2001,在此被包含引用。
可以理解,根据例2的算法可以用相应的方法代替。
患者和图像获取
我们包括14个表现急性缺血卒中症状的患者(NIHSS 12.5±5.64)。标准动态敏感性对比MRI,在发病3小时内完成,使用一台3.0T MR扫描器(SignaExcite HDx,通用电气医疗***,密尔沃基,威斯康辛州,美国)。平均通过时间由循环分块奇异值分解(oSVD)计算,具有自动动脉输入函数搜索算法,使用灌注分析软件PENGUIN(www.cfin.au.dk/software/penguin)。循环分块SVD进行:使用动脉输入函数的观察结果,形成循环分块
Figure BDA00003530197200173
矩阵;基于振荡指数对特征值进行阈值转换(如参考文献Wu O,
Figure BDA00003530197200174
L,Weisskoff RM,BennerT,Rosen BR,Sorensen AG。示踪剂到达时间不敏感技术,使用带循环分块去卷积矩阵的奇异值分解,用于估算在mr灌注加权成像中的流。磁共振医学:磁共振医学学会官方杂志/磁共振医学学会。2003:50:164-174,通过引用包含于本文中)。
人工勾画损伤轮廓
四位具有丰富的临床评估MTT图经验的神经系放射学家,使用免费提供的图像分析软件(“MRIcro”作者:Chris Rorden,2005,www.sph.sc.edu/comd/rorden/mricro.html),人工勾画急性MTT损伤。读者被要求基于他们最佳临床判断,勾画低灌注组织的范围。为了使勾画损伤轮廓反映“真实”独立临床专家之间评估人之间的一致性,在这项研究之前没有给定进一步的标准和说明。评估人根据自己的喜好,自由调整窗口/水平设置。未使用自动预处理,例如阈值转换。读者看不到所有临床和其他成像数据。
使用根据例2的算法与人工损伤勾画的轮廓的比较
以前的研究表明:在勾画损伤轮廓中评估人之间的可变性值得考虑。在没有损伤勾画轮廓的黄金法则的情况下,因此我们比较建议的、根据例2的算法与损伤轮廓置信度的程度,作为由专家间共识的可变程度的评估。因此,生成四个损伤评估,最大程度与由至少一个专家分类为低灌注的体素相一致,且最多地保留由所有专家分类为低灌注缺血的体素。我们于是可以比较由根据例2的算法估算的损伤轮廓,与这四个分等级的共识体积。
最后,我们比较水平集技术与标准阈值转换技术的性能,如果MTT超过在正常出现对侧的白质的平均值多于2个标准偏差,标记组织体素缺血。至于根据例2的水平集算法,这项技术适用于经过预处理删除脑室体素的相关MTT图。
结果
收敛
图6表明最初轮廓向患者的损伤边界的收敛,在右半球的两个分水区域具有提高的MTT。最初估计(Init.),如图6所示,显示了若干不连通的区域,在扩展的损伤内具有“正常”组织的小区域,正如基本阈值转换方法预期的。10次迭代后(Iter.)(对应于0.2秒),两个大的粘连的区域出现,而大多数小的区域消失。由于等式(3)中的曲线长度罚项,进一步迭代主要导致估算的轮廓平滑(100次迭代,1.8秒)。在100次和600次迭代(10.0秒)仅观察到较小的变化,表明快速收敛。在这项研究中,使用600次迭代,在所有患者中单切片计算时间的中值是10.0秒[9.9,10.1]。
与专家相比
图7是比较水平集损伤体积与每位专家的损伤体积,显示了人工估算的损伤体积与根据例2的算法的勾画MTT损伤轮廓之间的一致性。尽管我们注意到专家体积众有相当大的变化,未观察到与标识线有***偏差,这表明手动估算损伤体积和根据实施例2的算法估算的损伤体积之间有良好的整体一致性。专家之间,可以观察到趋势,评估人2个和3个,比评估人1个和4个典型地勾画更小的低灌注体积。
图7还说明了由四个专家估算的损伤体积的大的变化。在所有患者中,专家估算的最大和最小的体积之间的差异,大于平均损伤大小的10%,且14例中有6例这种差异超过30%。这种变化是由损伤边界轮廓的位置控制的。
图8显示了由根据例2(红色轮廓线)的算法估算的损伤边界的实例,以及六个具有典型损伤的患者的专家间一致性。一个在后部MCA区域中具有MTT损伤的患者显示在图8(a)中,图8(b)表明专家之间以及根据例2的算法和专家之间具有极好的一致性。图8(c)显示了一个实例,其中正常和低灌注组织之间后部边界的准确位置很难定义。根据例2的算法,然而,勾画区域的轮廓与专家相差很小,图8(d)。图8也显示了根据例2的算法识别离散损伤的能力,在这种情况下,在中脑动脉区域内具有低灌注区域的患者(图8(e-f)),在前颈动脉区域内具有损伤的患者(图8(g))。图8(i,k)再次表示低灌注的和未受影响的组织之间的边界不清楚的情况,导致在专家损伤勾画中的变化。在这两种情况中,根据例2的算法匹配具有多于2个专家间的共识的区域。我们也注意到根据例2的算法正确地避免了在脑室中的MTT高强度信号,图8(d,f,h,j)。在所有的情况中,根据例2的算法估算边界为一条平滑曲线,避免了使用阈值转换观察到的分散的小的“假的”损伤(参见图6)。我们注意到根据例2的算法识别边界,所述边界匹配多于两个专家之间共识的损伤,即使在低灌注的和正常组织之间没有清楚边界的情况下,例如(c-d),(i-j),(k-l)。然而,在一些情况下,根据例2的算法沿着大脑的边缘错误地检测小的损伤(j,l),这些小的损伤可用进一步的后处理删除。
图9显示了由根据例2的算法与专家根据他们的共识估算的损伤体积之间的相关性。估算的体积,与至少两个专家之间的共识确定的损伤,相关良好。与所有专家在大的损伤上的共识相比,观察到一些高估。相比之下,由根据例2的算法检测到的损伤体积,比所有专家勾画的轮廓包围的所有体素的并集定义的损伤更小。只有一个专家识别为低灌注的组织产生的损伤,比根据例2的算法确定的MTT损伤更大。然而,根据例2的算法,与2个或多个专家之间共识定义的损伤之间观察到良好的一致性。
我们使用人工和根据例2的算法之间的体积差,勾画的MTT共识损伤作为一致性的度量。表1显示了水平集技术和阈值转换的体积差的平均和标准偏差。在水平集技术和至少三个专家之间共识的损伤之间,观察到最小差-9.3±45.2ml(人工-根据例2的算法)。阈值转换方法,只有在损伤定义为所有专家选择的体素的并集的情况下,比水平集技术具有更低的偏离;但在所有情况下显示出更大的标准偏差。单个数据点落在95%一致性限制之外,正如从附录2中表III计算得到的,对于每一个共识程度和技术(然而对应于不同的患者)。
表III显示了专家共识损伤与STM方法以及对应于根据例3的算法的水平集方法之间的平均体积差(用ml表示)。相对于两个或三个专家之间共识的区域,根据例2的水平集技术算法显示出低偏离。STM(阈值转换)方法,只有在损伤定义为所有专家选择的体素的并集的情况下,比水平集技术具有更低的偏离;但在所有情况下显示出更大的标准偏差。
图10显示了四个专家共识损伤和阈值转换技术之间的相关性,其中当MTT值超过对侧白质中的平均值两个标准偏差时,体素识别为低灌注。请注意:相对于共识区域,阈值转换***地高估低灌注的组织体积。与水平集技术(参见图9)形成对比,阈值转换方法显示出与由专家选择的所有体素的并集获得的损伤的某些一致性;而与两个或更多个专家之间共识的区域相比,阈值转换方法始终高估体积。
在检测全共识区域时,根据例2的算法的中值灵敏度是0.74[0.60;0.89]。在小的损伤中观察到更低的灵敏度,这种测量数值上不稳定。考虑最大的损伤的75%(22ml以上),灵敏度是0.82[0.73;0.90]。如图8(l)所示,假的高强度可通过根据例2的算法,错误地检测为缺血组织,总的特异性是0.95[0.93;0.98],当忽略小于22ml的损伤时,总的特异性是0.96[0.91;0.98]。
讨论
人工延长MTT的区域,例如脑室,代表全自动勾画MTT损伤的关键挑战。我们建议使用并发DWI成像,作为方便的排除CSF体素的方法,而没有额外的以用户介入或速度的形式的成本。作为一种选择,脑室可以在最初的(无限TR)原始的T2加权DSC图像上被识别,或者通过使用CBV标准被识别。
尽管在CSF中大多数具有高MTT值的体素,在预处理过程中被删除,假的“高MTT”体素保留,典型地靠近脑桥和沿着皮质表面(参见图8(k)和8(l))。这样的区域可被建议的根据例2的算法错误地分割,从而降低特异性。进一步的后处理可被应用于排除例如,与主损伤对侧的体素。增大等式(3-4)中λ1的项导致更小和更同质的体素组,这些体素被勾画为低灌注。然而,在极小化假阳性率的同时,具有对MTT提高更敏感体素可丢失,增大假阴性分类率。作为一种选择,隔离的假阳性体素容易识别,且有限的用户介入要求为精确体积量化,删除这样的区域。
例3-通过DWI获得的图像的分割
材料与方法
形态灰度重构
DWI损伤容易被眼睛识别作为高信号区域,尽管弥散损伤包含伪影高信号;且由于身体结构和噪声伪影,高信号遍及图像。因此,图像强度的简单阈值转换导致高的假阳性和假阴性率,如图11B所示。典型的图像模糊,例如通过高斯核卷积,在某种程度上对伪影的补救,尽管存在潜在的“移动”损伤边界的代价。形态灰度重构唯一地截短图像中的高峰,倾向于保存粘连的区域的原始边界,正如由内核大小确定的。我们使用形态灰度重构来增强DWI损伤与背景之间的反差。与常规模糊过滤器相比,形态灰度重构保留边缘,参见图11C和11D,因此适合于自动DWI损伤分割。
根据例3的算法
自动生成的CSF蒙片删除DWI图中的脑室。CSF蒙片通过使用群集分析获得。k-均值群集的输入图像,是由结构图(BZERO)和由弥散加权图像(表观弥散系数,ADC)获得的图像组合的图像。向PWI图共同注册MNI(蒙特利尔神经学研究所)模板,自动生成左和右半球蒙片。通过CBV>0阈值转换生成整个脑部蒙片,通过形态细化删除整个脑部蒙片中稀疏连接的小岛。DWI体积中劣质切片中的非脑组织(例如眼睛、骨骼、空气等)对噪声成分有贡献,并增大假阳性率。因此,我们通过分析每一切片中的强度变化,删除劣质切片。
TTP损伤蒙片使用根据例1的算法自动生成。ADC图设定阈值为550mm2/sec,且删除TTP蒙片外的ADC损伤,参见图12A和12B。此外,删除每一切片上尺寸小于0.25mL的小岛。
使用形态灰度重构过滤DWI图。对于每一个存在ADC损伤的切片,使用相应的ADC蒙片确定切片特定的DWI损伤阈值,参见图12C。特别地,通过寻找ADC损伤蒙片的边缘,确定切片特定的阈值;将该阈值应用在灰度形态处理过的DWI图上;并求出包含在ADC蒙片边缘中的DWI强度的平均值。这样做的一个优势可在于获得DWI阈值的抗差估计。ADC损伤和DWI损伤不必空间匹配。为了确定DWI损伤,每一个DWI切片以切片特定的阈值进行阈值转换,且通过连通区域标记算法(CCL,在例1中描述)找到连通区域。删除没有与TTP损伤连通的DWI损伤区域。错配定义为:错配=TTParea-(DWIarea∩TTParea),因此,未包含在TTP损伤中的DWI损伤,在错配确定中是不相干的而被忽略。最后,DWI损伤中的孔被填充,边界被平滑。根据例3的分割算法,由
Figure BDA00003530197200211
公司在
Figure BDA00003530197200212
R2010a中实现。
可以理解,根据例3的算法可以用相应的方法代替。
患者和图像获取
与例1相同。
使用根据例3的算法勾画的损伤轮廓与人工勾画的DWI损伤轮廓的比较
具有丰富的临床工作经验的神经系放射学家,使用内部开发的软件半自动勾画DWI损伤轮廓。软件使用标准的阈值转换方法,其粗略地识别DWI损伤,随后神经系放射学家可以人工调整损伤边界,神经系放射学家看不到任何强度数值范围和其他脑图。基于人工勾画的DWI和TTP损伤蒙片,人工半暗带蒙片被确定;同样地,基于自动生成的DWI和TTP损伤蒙片,自动半暗带蒙片被确定。对于半暗带,敏感性和特异性在同侧半球中确定;对于DWI损伤,敏感性和特异性在TTP损伤中确定。
结果
与STM勾画的DWI损伤相比,由根据例3的算法勾画的DWI损伤,与人工勾画的DWI损伤有更好的空间和体积一致性。根据例3的算法,总的斯皮尔曼相关为R2=0.79,参见图13A。人工勾画轮廓与根据例3的算法之间的中值差是-1.4mL(±2SD:20mL)。使用根据例3的算法确定的半暗带体积,与人工勾画半暗带体积有良好的相关性R2=0.86,参见图14A。人工勾画轮廓和根据例3的算法之间的半暗带体积的中值差是-0.07mL(±2SD:73.5mL),中值灵敏度是71%(25th-75th四分位数:59-83%)和中值特异性是95%(25th-75th四分位数:92–98%)。人工勾画轮廓与STM之间的半暗带体积的中值差是-41.1mL(±2SD:116.2mL),中值灵敏度是64%(25th-75th四分位数:49-77%)和中值特异性是92%(25th-75th四分位数:88–95%)。最后,Dice系数分析表明根据例3的算法估算的半暗带(中值=0.67,25th-75th四分位数:0.51–0.75),与使用STM估算的半暗带相比(中值=0.46,25th-75th四分位数:0.27–0.60),更相似于人工勾画的半暗带(p<0.001)。
图11:预处理DWI图,图11A表示原始DWI图,图11B表示使用STM分割的原始DWI图,这在大脑的中间部分产生DWI损伤和假阳性区域。使用模糊过滤器,例如通过高斯核卷积,在某种程度上对伪影的补救,但也模糊了损伤边界(图11C),然而形态灰度重构增强了损伤和带有清晰边界(图11D)的背景之间的反差。
图12:确定DWI阈值。
图12A:ADC图在550mm2/sec进行阈值转换。
12B)ADC损伤,其不包含在TTP损伤中,被删除。图12C:ADC蒙片被用于形态灰度重构的DWI图,以确定切片特定的阈值。
图12D显示了通过根据例3的算法自动分割的DWI损伤,叠加在原始DWI图上。
图13:在TTP损伤内,通过根据例3的算法、STM和人工勾画轮廓估算的DWI损伤的体积比较。根据例3的算法的总的斯皮尔曼相关是R2=0.79,STM的斯皮尔曼相关是R2=0.62。子图表明:图13A:DWI损伤的体积比较,第一个轴表示人工确定的蒙片体积(刻度线跨度为0-120mL),第二个轴表示根据例3的算法确定的体积(刻度线跨度为0-120mL);图13B:DWI损伤的体积比较,第一个轴表示人工确定的蒙片体积(刻度线跨度为0-120mL),第二个轴表示根据STM方法确定的体积(刻度线跨度为0-120mL);图13C:人工确定的体积与根据例3的算法确定的体积相比的Bland-Altman图,第一个轴表示平均损伤体积(刻度线跨度为0-150mL),第二个轴表示差值(刻度线跨度为-120-80mL),图13D:人工确定的体积与根据STM方法确定的体积相比的Bland-Altman图,第一个轴表示平均损伤体积(刻度线跨度为0-150mL),第二个轴表示差值(刻度线跨度为-120-80mL)。在图12C-D中,实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。
图14:在最佳阈值时,通过根据例1和例3的算法、STM和人工勾画轮廓估算的半暗带的体积比较。图14A:当前神经系放射学家的最佳阈值估计为>2.8秒(R2=0.86);图14B:与人工勾画轮廓(R2=0.62)相比STM高估TTP损伤;图14C)人工与算法相比的Bland-Altman图以及图14D)人工与STM相比的Bland-Altman图。实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。图14A:半暗带的体积比较,第一个轴表示人工确定的蒙片体积(刻度线跨度为0-500mL),第二个轴表示根据例1的算法确定的体积(刻度线跨度为0-500mL);图14B:半暗带的体积比较,第一个轴表示人工确定的蒙片体积(刻度线跨度为0-500mL),第二个轴表示根据STM方法确定的体积(刻度线跨度为0-500mL);图14C:人工确定的半暗带体积与根据例1的算法确定的体积相比的Bland-Altman图,第一个轴表示平均损伤体积(刻度线跨度为0-500mL),第二个轴表示差值(刻度线跨度为-250-150mL);图14D:人工确定的半暗带体积与根据STM方法确定的体积相比的Bland-Altman图,第一个轴表示平均损伤体积(刻度线跨度为0-500mL),第二个轴表示差值(刻度线跨度为-250-150mL)。在图14C-D中,实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。
图15:分类性能,图15A:通过根据例1和例3的算法和STM进行半暗带分类的灵敏度和特异性图。对于根据例1和例3的算法,中值灵敏度是71%(25th-75th四分位数:59-83%),平均特异性是95%(25th-75th四分位数:92–98%)。STM的中值灵敏度是64%(25th-75th四分位数:49-77%),STM的中值特异性是92%(25th-75th四分位数:88–95%);图15B:与通过STM(中值=0.46)p<0.001估算的半暗带相比,在根据例1和例3的算法(中值=0.67)估算的半暗带中,Dice系数明显更高。
图16:通过a)专家,b)根据例3的算法,以及c)STM勾画DWI损伤的实例。我们看到通过根据例3的算法勾画的DWI损伤,比人工勾画的DWI损伤更一致。然而,通过STM估算的DWI损伤更分散,并包含假阳性损伤。
图17:通过a)专家,b)根据例1和例3的算法,以及c)STM勾画半暗带的实例。我们看到通过根据例1和例3的算法勾画的半暗带,与人工勾画的半暗带更一致,并显示粘连的半暗带。然而,由STM估算的半暗带更分散,并包含大面积的假阳性半暗带。
例4-包括TTP数据的图像的第一个可供选择的分割
材料与方法
形态灰度重构
如例1所示。
连通区域标记
如例1所示。
水平集
CCL生成的损伤边界,与人工勾画的认知边界相比,可不平滑。我们通过搜索平滑曲线C,间接估算认知的损伤边界,同时极小化该曲线内外的TTP变化。我们使用由Chan-Vese建议的Mumford-Shah能量极小化问题的水平集公式。详细资料可以在参考文献中找到,“无边缘的活动轮廓”,Chan TF,Vese LA著,IEEE图像处理汇刊,第10卷,第2期,第266-277页,2001,在此被包含引用。例4的水平集方法与例2的水平集方法相似。
根据例4的算法
如例1所示,除了步骤:
-平滑损伤边界,
被替换为
-应用水平集算法,水平集算法带由CCL算法创建的TTP损伤蒙片启动。
从后面的步骤产生的蒙片,可以被认为是TTP损伤蒙片,已经过在TTP>4秒的软阈值转换。
根据例4的TTP分割算法,在
Figure BDA00003530197200241
R2010a中,由
Figure BDA00003530197200242
实现,但本领域技术人员可以理解根据例4的分割算法也可以用其他编程语言实现。
可以理解,根据例4的算法可以用相应的方法代替。
患者和图像获取
我们使用一组168急性局部缺血卒中患者(♀=70)来评估APS的性能。患者准许进入参与多中心研究I-KNOW的医院。患者的中值年龄是70岁,从症状出现到最初的MRI扫描的中值时间是152分钟(表V)。
表V显示了淤血患者特征。范围表示第一和第三四分位数。
中值急性DWI损伤是58.2mL,且中值NIHSS是10。I-KNOW协议定义了七个卒中子类,其中心源性脑梗塞是主要的卒中子类,随后是未确定的子类和带有明显颈动脉狭窄的大血管疾病。标准动态磁敏感对比MRI在不同型号的扫描器上完成(GE Signa Excite1.5、GE Signa Excite3、GE Signa HDx1.5、Siemens TrioTim3、Siemens Avanto1.5、Siemens Sonata1.5、Philips Gyroscan NT1.5以及Philips Intera1.5)。PWI序列(TE30-50ms,TR1500ms,FOV24cm,矩阵128x128,切片厚度5mm),在钆对比剂(0.1mmol/kg)以5ml/s的速率静脉注射,随后30ml生理盐水以同样速率注射之后获得。获得的信号强度转换为时量关系曲线,平滑伽马变量拟合过程应用于每一个时量关系曲线。
依据例1的算法与人工TTP损伤勾画的轮廓的比较
四个具有丰富的临床工作经验的评估人(一个神经系放射学家和三个放射科医师),使用内部开发的软件人工勾画TTP损伤轮廓。评估人看不到任何强度数值范围和其他脑图。未使用自动预处理,例如阈值转换。自动确定的TTP损伤体积,与人工勾画的TTP损伤体积相比,具有一个、两个、三个和四个评估人间一致性。最后,TTP损伤在>4秒时也使用STM估算,且与人工勾画TTP损伤对比。在这项研究中,我们在TTP>4秒时估算TTP损伤,这个特定的阈值已在几项研究中评价为识别具有梗塞风险的组织,例如Sobesky,J.,et al.,“哪个达峰时间阈值最好地识别半暗带流?灌注-加权磁共振成像与正电子放射断层造影术在急性缺血性卒中上的比较”,卒中,2004.35(12):第2843-7页,在此被包含引用。为量化算法/STM与人工勾画损伤轮廓之间的几何和体积相似度,我们使用Dice系数(定义:两个蒙片覆盖蒙片的组合集的交集的两倍),Dice系数是0和1之间的一个数。一表示完全重叠,0(零)表示无重叠。敏感性和特异性在同侧半球上确定。
结果
体积比较
人工和根据例4的算法自动勾画TTP损伤之间的最佳体积相关,在于三个评评估人之间的一致性,且斯皮尔曼相关是R2=0.89,参见图18A。与人工勾画轮廓(数据未显示)相比,在一个和两个评估人之间一致时,根据例4算法的自动方法低估了TTP损伤体积,而在四个评估人之间一致时,算法高估了TTP损伤体积。在最佳体积相关时,人工勾画轮廓与算法之间的平均差为-0.26mL(±2SD:33mL)参见图18C。
TTP图阈值在TTP>4秒,与三个评估人之间一致显示出最佳相关。斯皮尔曼相关为R2=0.73,参见图18B。人工勾画轮廓与STM之间的平均差为-281mL(±2SD:46.64mL)参见图18D。
图18显示了通过根据例4的算法、在TTP>4秒时的STM和人工勾画轮廓估算的TTP损伤的体积比较。图18A显示了根据例4的算法与人工勾画轮廓之间的体积相关(R2=0.89),第一个轴(刻度线跨度为0-400)表示人工确定蒙片体积[mL],第二个轴(刻度线跨度为0-400)表示通过根据例4的算法确定的损伤体积[mL];图18B显示了STM与人工勾画轮廓之间的体积相关(R2=0.73),第一个轴(刻度线跨度为0-400)表示人工确定蒙片体积[mL],第二个轴(刻度线跨度为0-400)表示由STM方法确定的损伤体积[mL]。
图18C表示人工勾画轮廓与根据例4的算法相比的Bland-Altman图,图18D表示人工勾画轮廓与STM相比的Bland-Altman图,实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。第一个轴(刻度线跨度为0-500)表示平均损伤体积[mL],第二个轴(刻度线跨度为-250-400)表示人工确定的体积与根据例4的算法和STM方法确定的体积之间的差值[mL],分别在图18C和图18D中。
图20显示了通过根据例4的算法和STM产生的分割的实例。在上面的行中,TTP图像和随后行中表明:a)人工蒙片;b)依据例4的算法获得的蒙片;以及c)从STM在TTP图上STM覆盖(白)获得的蒙片。两个非人工图在TTP>4秒时获得。通常,与阈值转换方法相比,根据例4的算法估算边界为平滑曲线。此外,根据例4的算法能够很大程度避免假阳性损伤(图20)。
分类性能
在三个评估人之间一致时,根据例4的算法的中值灵敏度(中值:79%和25th-75th四分位数:40-85%)没有显著(p=0.22)高于STM的灵敏度(中值:72%和25th-75th四分位数:40–85)。相比之下,算法的中值特异性(95%,25th-75th四分位数:90-98%)明显更高(p<0.001),与STM95%(25th-75th四分位数:87-96%)相比,参见附录2中的表IV。此外,Dice系数分析表明通过算法估算的TTP损伤蒙片(中值Dice系数:0.74,25th-75th四分位数:0.52-0.82),明显与人工勾画TTP损伤蒙片具有更高(p<0.001)相似度,与通过STM生成的蒙片(中值Dice系数:0.54,25th-75th四分位数:0.24-0.70)相比,参见图19。
表IV描述了根据例4的算法和TTP>4秒时STM的灵敏度和特异性,用百分比表示。括号中给出了第一和第三四分位数,也进一步给出了秩和检验的统计数据。
图19显示了根据例4的算法与STM的Dice系数的箱形图。算法的中值Dice系数是0.74(0.52-0.82),STM的中值Dice系数是0.54(0.24-0.70)。
例5-通过DWI获得的图像的第一个供选择的分割
材料与方法
形态灰度重构
如例3所示。
水平集
在DWI图像中,损伤与正常组织之间的图像对比是明显的,图像对比通过形态灰度重构被进一步增强。然而,简单阈值转换存在缺乏已确定的适用于不同个体的DWI阈值,和由于DWI图中内在噪声产生的“假的”损伤的问题。水平集算法通过搜索平滑曲线C,间接估算DWI损伤边界,同时极小化该曲线内外的DWI强度变化。我们使用由Chan-Vese建议的Mumford-Shah能量极小化问题的水平集公式。详细资料可以在参考文献中找到,“无边缘的活动轮廓”,Chan TF,Vese LA著,IEEE图像处理汇刊,第10卷,第2期,第266-277页,2001,在此被包含引用。例5的水平集方法与例2的水平集方法相似。
根据例5的算法
自动生成的CSF蒙片删除DWI图中的脑室。CSF蒙片使用群集分析获得。k-均值群集的输入图像,是由结构图(BZERO)和弥散加权成像(表观弥散系数,ADC)组合的图像。向PWI图共同注册MNI(蒙特利尔神经学研究所)模板,自动生成左和右半球蒙片。通过CBV>0阈值转换生成整个脑部蒙片,通过形态细化删除整个脑部蒙片中稀疏连接的小岛。DWI体积中劣质切片中的非脑部组织(例如眼睛、骨骼、空气等)对噪声成分有贡献,并增大假阳性率。因此,我们通过分析每一切片中的强度变化,删除劣质切片。TTP损伤蒙片使用例4中描述的算法自动生成。
水平集算法必须带有一个最初的蒙片启动,最初的蒙片在550mm2/sec对ADC图进行阈值转换得到。使用内核尺寸对应于25mL的形态二值重构,删除每个切片上ADC蒙片中小于0.25mL的小岛,所述形态二值重构类似于形态灰度重构。
使用形态灰度重构过滤DWI图,水平集算法带有形态灰度重构处理过的DWI图像和ADC蒙片启动。
错配定义为TTP损伤减去DWI和TTP损伤的交集,错配=TTParea-(DWIarea∩TTParea)。自动半暗带分割算法在Matlab R2010a中实现(马萨诸塞州,纳蒂克,MathWorks公司)
可以理解,根据例5的算法可以用相应的方法代替。
患者和图像获取
我们使用一组168位急性局部缺血卒中患者(♀=70)来评估APS(自动半暗带分割)的性能。患者准许进入参与多中心研究I-KNOW的医院。患者的中值年龄是70岁,从症状出现到最初MRI扫描的中值时间是152分钟(表V)。中值急性DWI损伤是2.4mL,中值NIHSS是10。I-KNOW协议定义了七个卒中子类,其中心源性脑梗塞是主要的卒中子类,随后是未确定的子类和带有明显颈动脉狭窄的大血管疾病(表V)。MRI在不同型号的扫描器上完成(GE Signa Excite 1.5、GESigna Excite3、GE Signa HDx1.5、Siemens TrioTim3、Siemens Avanto1.5、Siemens Sonata 1.5、Philips Gyroscan NT 1.5以及Philips Intera 1.5)。在b-值=0和b-值=1.000sec/mm2时获得的平面回波DWI。
使用根据例5的算法勾画的损伤轮廓与人工勾画的DWI损伤轮廓的比较
四个具有丰富的临床工作经验评估人(一个神经系放射学家和三个放射科医师),使用内部开发的软件人工勾画DWI损伤轮廓,且评估人看不到任何强度数值范围和其他脑图。未使用自动预处理,例如阈值转换。自动确定的DWI损伤体积,与人工勾画的DWI损伤体积相比,具有一个、两个、三个和四个评估人间一致性。基于人工勾画的DWI和TTP损伤蒙片,人工半暗带蒙片被确定;同样地,基于自动生成的DWI和TTP损伤蒙片,自动半暗带蒙片被确定。对于半暗带,在同侧半球和Dice系数中,确定灵敏度和特异度。
结果
体积比较
根据例5的算法与人工勾画DWI损伤轮廓之间的最佳体积相关,表示在三个评估人间一致。斯皮尔曼相关是R2=0.61,人工与方法之间的平均差是0.23mL(±2SD:25.6mL),参见图22A和图22C。相比之下,在STM损伤体积和人工勾画轮廓的蒙片体积之间的体积相关不好。斯皮尔曼相关是R2=0.33,人工和自动方法之间的平均差是-124.4mL(±2SD:84.5mL),参见图22B和图22C。
通过根据例5的算法确定的半暗带体积,与人工勾画半暗带体积有良好的相关性R2=0.84,参见图23。人工勾画轮廓与算法之间半暗带体积的中值差为-2.1mL(±2SD:64.7mL)。
STM半暗带蒙片,与人工勾画的半暗带蒙片相比,显示相关是R2=0.62,平均差是-15.6mL(±2SD:93.6mL)。
分类性能
根据例5的算法(敏感度:中值=75%,25th-75th四分位数:63-85%;DC:中值=0.70,25th-75th四分位数:0.57-0.78)胜过STM(灵敏度:中值=52%,25th-75th四分位数:31-70%;DC:中值=0.44,25th-75th四分位数:0.26–0.60),根据灵敏度(p<0.001)且Dice系数(p<0.001)。然而,两个方法之间没有观察到特异性(p=0.22)的显著差异。
图2显示初始蒙片(子图A)到包围DWI损伤的最终蒙片(子图D)的演化,子图A、B、C和D分别对应于0、80、160和200次迭代。
图22:在TTP损伤内,通过根据例5的算法、STM和人工勾画轮廓估算的DWI损伤的体积比较。与人工勾画轮廓相比,根据例5的算法的斯皮尔曼相关是R2=0.61。子图表示:图23A:DWI损伤的体积比较,第一个轴表示人工确定的蒙片体积(刻度线跨度为0-150mL),第二个轴表示根据例5的算法确定的体积(刻度线跨度为0-150mL);图22B:DWI损伤的体积比较,第一个轴表示人工确定的蒙片体积(刻度线跨度为0-250mL),第二个轴表示根据STM方法确定的体积(刻度线跨度为0-250mL);图22C:人工确定的体积与根据例5的算法确定的体积相比的Bland-Altman图,第一个轴表示平均损伤体积(刻度线跨度为0-500mL),第二个轴表示差值(刻度线跨度为-250至150mL);图22D:人工确定的体积与根据STM方法确定的体积相比的Bland-Altman图,第一个轴表示平均损伤体积(刻度线跨度为0至500mL),第二个轴表示差值(刻度线跨度为-250至150mL)。在图22C-D中,实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。平均差为0.23mL(±2SD:25.6mL),通过STM勾画轮廓产生不利的结果,斯皮尔曼相关为R2=0.33,平均差为-124.4mL(±2SD:84.5mL)。
图23:由根据例4和例5的算法、STM和人工勾画轮廓估算的半暗带的体积比较。图23C)人工与算法相比的Bland-Altman图以及图23D)人工与STM相比的Bland-Altman图。实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。图23A:半暗带的体积比较,第一个轴表示人工确定的蒙片体积(刻度线跨度为0-500mL),第二个轴表示根据例4的算法确定的体积(刻度线跨度为0-500mL);图23B:半暗带的体积比较,第一个轴表示人工确定的蒙片体积(刻度线跨度为0-500mL),第二个轴表示根据STM方法确定的半暗带体积(刻度线跨度为0-500mL);图23C:人工确定的体积与根据例4的算法确定的体积相比的Bland-Altman图,第一个轴表示平均损伤体积(刻度线跨度为0-500mL),第二个轴表示差值(刻度线跨度为-250-150mL);图23D:人工确定的体积与根据STM方法确定的体积相比的Bland-Altman图,第一个轴表示平均损伤体积(刻度线跨度为0-500mL),第二个轴表示差值(刻度线跨度为-250-150mL)。在图23C-D中,实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。算法与人工勾画轮廓相比,斯皮尔曼相关是R2=0.84,平均差为-2.1mL(±2SD:64.7mL)。通过STM勾画半暗带轮廓产生不利的结果。
图24A-B:分类性能,依据灵敏度(p<0.05)和Dice系数(DC)(p<0.001)测量,算法与STM之间存在明显差异,然而,特异性(p=0.26)无明显差异观察到。图24A:通过根据例4和例5的算法和STM,进行半暗带分类的灵敏度和特异性图;图24B:Dice系数(DC)。
图25:通过a)专家(三个评估人之间一致),b)根据例4和例5的算法,以及c)STM勾画半暗带轮廓的实例。
例6-包括TTP数据的图像的第二个可供选择的分割
PWI分割
各种PWI度量(如去卷积曲线达峰时间(Tmax),平均通过时间(MTT)和达峰时间(TTP))被建议,以识别局部缺血半暗带,它们落入两个类别:从组织浓度曲线(例如TTP)获得的摘要图和通过使用动脉输入函数(例如Tmax和MTT)去卷积生成的图。去卷积的图假设是生理学上可解释的,而摘要图不必反映生理事件。相比之下,与去卷积图相比,由于摘要图的性质更简单,摘要图噪声更少,因而更适合于分割目的。在这个实例中,我们使用TTP的规格化的版本作为PWI度量。TTP>4秒已在几项研究中评价为识别低灌注组织,然而在4秒将TTP图进行阈值转换产生很多分散的假阳性损伤,如脑室、眼睛和敏感性伪影。在这个实例中,我们提供了一种算法,该算法在TTP图上用认知边界,识别粘连的低灌注损伤。
方法
形态灰度重构
如例1所示。
形态重构
形态重构是二值图像转换技术,参见参考文献“图像分析中的形态灰度重构:应用程序和高效算法”,文森特,L.,IEEE图像处理汇刊,1993.2(2):第176-201页,通过引用完整地包含于本文中,尤其对部分引用。在这个上下文中,形态重构用于删除假阳性TTP损伤,这项技术可以在概念上被认为是种子点的迭代扩张(参见图26,第4行包括标记为J、K、L的图像),种子点由通过TTP图的阈值转换获得的最初蒙片约束。未包含种子点的体素群集,不保留在重构的蒙片中。
图26给出了TTP损伤分割算法中步骤的示意图。对应于非脑组织删除的行(第1行,具有图像A、B、C),生成侧性蒙片(第2行,具有图像D、E、F),种子点检测(第3行,具有图像G、H、I),形态重构(第4行,具有图像J、K、L)以及水平集(第5行,具有图像M、N、O)。在最初在原始TTP图(A)上,通过应用整个脑蒙片(B)删除非脑组织,以便获得删除非脑组织的图像(C、D)。随后损伤的侧性(E)通过分析图像强度值来确定,以便获得只保留相关侧性的图像(F、G)。为在同侧半球上检测种子点,TTP图像(G)通过形态灰度重构过滤,种子点定义为具有最大强度(I)的体素或多个体素。最初的TTP损伤蒙片,由阈值转换原始TTP图(J)创建。重构由最初TTP损伤蒙片约束的种子点,保留TTP损伤蒙片并排除假阳性损伤(K)。最后,通过应用水平集算法(N),调整TTP损伤蒙片边界(L、M),以便获得损伤(O)。
水平集
自动生成的TTP损伤可不必与人工勾画的损伤具有相同的范围,此外边界与人工勾画的认知边界相比可以是不平滑的。我们通过搜索平滑曲线C,间接估算认知的损伤边界,同时极小化该曲线内外的TTP变化变化。我们使用由Chan-Vese建议的Mumford-Shah能量极小化问题的水平集公式,比较参考文献“一种无边缘的活动轮廓模型”,Chan,T.和L.Vese著,计算机视觉中的标度空间原理,1999.1682:第141-151页,通过引用完整地包含于本文中。
算法
自动生成大脑蒙片,由整个大脑蒙片和CSF蒙片组成,删除非脑组织例如脑室、眼睛等。整个大脑蒙片通过CBV图阈值转换(CBV>0)获得,未连通的或稀疏连接的区域例如眼睛和易感性伪影,通过应用形态算法被去掉;正如thinning和H-break,参见参考文献“使用Matlab的数字图像处理”,Gonzalez,R.C.,R.E.Woods和S.L.Eddins著,2004,通过引用完整地包含于本文中,尤其引用第9章。CSF蒙片在高对比度CSF-组织图上,通过体素强度的直方图分离获得。CSF-组织图通过组合结构图(T2)和弥散-加权图(表观弥散系数,ADC),(参见图26,第1行包括标记为A、B、C的图)。向PWI图共同注册MNI(蒙特利尔神经学研究所)模板,自动生成左和右半球蒙片(参见图26,第2行,包括标记为D、E、F的图)。非脑组织(例如眼睛、骨骼、空气等)和静脉(横窦和乙状窦)在PWI体积中的劣质切片中,对噪声成分有贡献,并增大假阳性率。因此,我们通过分析每一切片中的强度变化,删除劣质切片。
在自动检测损伤侧性和种子点之前,TTP图由形态灰度重构过滤。损伤侧性的检测基于在过滤的TTP图上任何一个半球上提取的三个特征。对于每一个特征,“获胜的”半球被指定1,最终该半球具有至少两个正数结果被认为同侧半球(参见图26,第2行包括标记为D、E、F的图像)。特征是:1)左和右半球的平均强度值;2)在TTP图上最大强度体素的侧性以及3)具有强度<600*10-6mm2/sec的ADC体素的数量。种子点定义为在同侧半球上具有最大强度体素(或多个体素);参见图26,第3行包括标记为G、H、I的图像。
为了平均对侧半球内的TTP延迟,TTP图被规格化。最初的TTP损伤蒙片通过同侧TTP图在4秒的阈值转换找到,使用最初的TTP损伤作为约束,种子点在三维中重构。最终,水平集算法带重构的TTP损伤蒙片启动(参见26,第5行包括标记为M、N、O的图像)。最终的蒙片可以被认为是TTP>4秒时的软阈值转换的TTP损伤蒙片。
自动TTP损伤分割算法在Matlab R2010a中实现(马萨诸塞州,纳蒂克,MathWorks公司)。
患者和图像获取
如例4所示。
自动和人工TTP损伤勾画轮廓的比较
四个在神经系放射学上具有丰富的临床工作经验的评估人(一个神经系放射学家和三个放射科医师),使用内部开发的软件人工勾画TTP损伤轮廓。评估人看不到任何强度数值范围和其他脑图。自动确定的DWI损伤体积,与人工勾画的DWI损伤共识图相比,具有一个、两个、三个和四个评估人之间的一致性。最后,在>4秒时也使用STM估算TTP损伤,且与人工勾画TTP损伤对比。在这个实例中,我们在TTP>4秒时估算TTP损伤,这个特定的阈值已在几项研究中评价为识别具有梗塞风险的E组织的最佳阈值。为量化APS/STM和人工勾画损伤轮廓之间的几何和体积相似度,我们使用体积相关和Dice系数(DC)。DC(定义:两个蒙片覆盖蒙片的组合集的两次交集),DC是一个0-1之间的数,其中1表示完全重叠,0表示无重叠。分类性能根据灵敏度和特异性评估,且在在同侧半球上确定。
结果
体积比较
图27显示了由算法、在TTP>4秒的STM和人工勾画轮廓估算的TTP损伤的体积对比。A)算法和人工勾画轮廓之间的体积相关(R2=0.92)(x轴表示人工蒙片体积[mL],y轴表示算法蒙片体积[mL]);B)STM和人工勾画轮廓之间的体积相关(R2=0.79)(x轴表示人工蒙片体积[mL],y轴表示STM蒙片体积[mL]);C)人工勾画轮廓与算法相比的Bland-Altman图(x轴表示平均损伤体积[mL],y轴表示差值)以及D)人工勾画轮廓与STM相比的Bland-Altman图(x轴表示平均损伤体积[mL],y轴表示差值)。实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。
人工和自动勾画的TTP损伤之间的最佳体积相关、dice系数(DC)和灵敏度/特异性,具有三个评估人之间一致性。APS的斯皮尔曼相关是R2=0.91,参见图27A。与人工勾画轮廓(数据未显示)相比,在一个和两个评估人之间一致时,APS低估了真实TTP损伤体积,而在四个评估人之间一致时,APS高估了真实TTP损伤体积。人工勾画轮廓和APS之间的平均差是-2.3mL,±2SD是63.6mL(p<0.001),参见图27C。
STM的斯皮尔曼相关是R2=0.78,参见图27B。人工勾画轮廓和STM之间的平均差是-22.6mL,且±2SD是85.8mL,参见图27D。
几何比较
APS在与人工勾画TTP损伤蒙片的几何一致性方面胜过STM。APS的中值DC是0.74,25th-75th内四分位数范围是0.50-0.81(p<0.001)。相比之下,STM的中值DC是0.53,25th-75th内四分位数范围是0.21-0.68,参见图28。
中值处理时间是每患者3.4秒(1.9-4.9),只依赖于TTP损伤的体积。
图28显示了APS和STM的Dice系数的箱形图。APS的中值DC是0.74(0.50-0.81),STM的中值DC是0.53(0.21-0.68)。
分类性能
在三个评估人之间一致时,APS在分类上明显胜过STM。APS的中值灵敏度是77%,25th-75th内四分位数范围(IQR)是58-87%(p<0.01),中值特异性是97%和25th-75th IQR是95-100%(p<0.05)。STM的中值灵敏度是65%,25th-75thIQR是38-80%;中值特异性是95%,25th-75th IQR是92-98%。
表VI显示了在同侧半球上确定的自动方法的灵敏度和特异性,用百分比表示。括号中给出了第一和第三四分位数,也进一步给出了秩和检验的统计数据。
讨论
本发明的一个可的优势在于:提供一种时间有效和全自动的分割算法,用以在相关TTP图上勾画TTP损伤轮廓,使用连通区域标记算法和Mumford-Shah分割的水平集实现。APS算法显示与人工勾画的损伤具有良好的一致性,在体积相关、DC和灵敏度/特异性方面明显胜过阈值TTP损伤蒙片。
图29显示通过APS(C行)和STM(D行)自动生成的TTP损伤蒙片的实例,以及由专家(B行)人工勾画的TTP损伤蒙片覆盖在TTP图上的实例。A)行显示了没有蒙片的图,三列对应于三种不同的图像。
合适的排除非脑组织的大脑蒙片,对于自动方法和人工勾画之间获得良好的一致性是有益的。使用好的脑部蒙片的效果,显示在图29中第1列。STM损伤蒙片显示了与人工勾画在识别低灌注区域上,具有良好的一致性;然而,由于缺乏合适的大脑蒙片,与APS相比假阳性率更高。大体上,来自脑室的噪声限制了APS或其他自动方法的最佳性能。在TTP图上,CSF发出高强度信号,会被分割作为低灌注区域,参见图29,第2列,STM蒙片(第4行)。对于APS我们开发了快速基于强度的脑室分割算法,该算法极大地删除CSF,为避免高估TTP损伤是有益的。基于模板注册的脑室蒙片,是基于强度的方法的可供选择的替代方法,但以时间效率为代价。
识别小的TTP损伤(<3mL)是一种挑战,因为噪声成分,例如眼睛,扩展了小损伤的体积,图29,第1列。已经建议删除小于某个阈值的区域,尽管该方法消除一些噪声成分,同时也删除了低灌注体素,因此增加了假阴性率。使用CCL的噪声抵制特性和水平集的一个优势,可在于我们克服识别小的TTP损伤的忧虑。
其中我们引入Dice系数,作为对估计自动分割方法性能的附加的方法。以前,体积相关和灵敏度/特异性已被用作评估方法,然而这些方法存在不同的缺陷。体积相关系数给出了估算的蒙片体积与人工获得的蒙片体积之间的关联的好的总结,然而,它不能反应蒙片之间的几何一致性。例如,在图29和第2列所示的情况中,对于三种方法体积非常相似(人工:16.5mL,APS:15mL和STM17mL),尽管STM的损伤体积包含许多假阳性体素。
在小的损伤中的分类度量灵敏度对假阴性敏感且可被向下降,而特异性被人工提高到刚好确定的区域,通常同侧半球例如图29和第1列(特异性APS和STM:0.99–0.99)。对于评估自动分割过程的性能,DC不是根本的方法,因为它对小的损伤很敏感,例如图29和第1列(用于APS和STM的DC:0.70和0.19),由于这个原因,我们建议至少报告体积相关(关联测量)、灵敏度/特异性(分类测量)和DC(一致性测量),在相关的研究中由于可再现性和标准参考,直到我们有一个方法在自动分割算法中覆盖性能的不同方面。
在这个实例中我们使用TTP图像,作为PWI度量标准,但考虑到不同的PWI度量标准可被用于最佳地识别具有梗塞风险的组织,考虑到扩展在这个实例中建议的方法到例如Tmax,MTT和第一瞬间图。
在患者之间识别具有梗塞风险的组织,生理上真正的TTP阈值已在几项研究中评估。TTP>4秒已经在急性缺血性卒中的MRI-PET研究经过验证的,且在临床实践中可以接受。
总之,根据例6的算法能够检测围绕定影剂-灌注区域的平滑边界,显示出与专家勾画轮廓有良好的一致性。
例7-通过DWI获得的图像的第二个供选择的分割
在这个实例中,提供了一种算法,该算法使用并发信息形成ADC图,在DWI图像上自动分割DWI损伤。结合分割DWI损伤的自动算法与分割PWI损伤的自动算法,提供一种能够呈现全自动的半暗带分割工具的方法,例如用于临床使用的工具。
形态灰度重构
如例3所示。
形态重构
形态重构是二值图像转换技术[10],在这个上下文中形态重构用于在DWI图像中消除假阳性DWI损伤。这项技术可以在概念上被认为是种子点的迭代扩张(参见图30,第5行3060),种子点由通过DWI图的阈值转换获得的最初蒙片限制。未包含在种子点中的体素群集不保留在重构的蒙片中。
算法
自动生成的脑部蒙片,由整个脑部蒙片和CSF蒙片组成,消除非脑组织例如脑室、颅骨等,参见图30第1行3052。整个大脑蒙片通过CBV图阈值转换(CBV>0)获得,未连接的或稀疏连接的组件例如眼睛和易感性人工制品通过应用形态算法被去掉;正如thinning和H-break参见参考文献“使用Matlab的数字图像处理”,Gonzalez,R.C.,R.E.Woods和S.L.艾汀斯,2004,在此完整地被引用。CSF蒙片在高对比度CSF-组织图上,通过体素强度的直方图分离获得,通过组合结构图(T2)和ADC获得。向PWI图共同注册MNI(蒙特利尔神经学研究所)模板,自动生成左和右半球蒙片。DWI体积中劣质切片中的非脑组织(例如眼睛、骨骼、空气等)对噪声成分有贡献,并增大假阳性率。因此,我们使用例6中描述的、自动生成的TTP损伤蒙片来确定DWI感兴趣的切片。最初ADC和DWI蒙片,在600*10-6mm2/sec由阈值转换ADC图而创建,DWI图像在对侧平均强度+两个标准偏差,参见图30和第2-3行3045-3056。乘最初蒙片保留DWI损伤的核心,去掉假阳性损伤(标记图像),参见图30和第4行3058。通过重构由最初DWI损伤蒙片约束的标记图像产生最终DWI损伤蒙片,参见图30和第5行3060。PWI-DWI错配定义为TTP损伤减去DWI和TTP损伤的交集。自动半暗带分割算法在MatlabR2010a中实现(马萨诸塞州,纳蒂克,MathWorks公司)
图30给出了DWI损伤分割算法中步骤的示意图。最初在原始DWI图(A)上,应用整个脑蒙片(B)删除非脑组织,以获得图像(C)。随后在600*10-6mm2/sec,通过阈值转换ADC图(D),创建ADC蒙片(E)。图像(F)表示用ADC蒙片(E)覆盖的图像(D)。在对侧平均强度+两个标准偏差,阈值转换经形态灰度重构的DWI图像(G),生成DWI蒙片(H)。图像(I)表示用ADC蒙片(H)覆盖的图像(G)。标记图像(L)通过ADC蒙片(J)和DWI蒙片(K)相乘构建,生成蒙片表明大脑中具有低ADC强度和高DWI强度的区域。重构由最初DWI蒙片(N)约束的标记图像(M),产生确定的DWI损伤蒙片(O)。(R)DWI损伤蒙片(Q)在原始DWI图像(P,A)上覆盖。
患者和图像获取
如例5所示。
自动和人工DWI损伤和半暗带勾画轮廓的比较
四个在神经系放射学上具有丰富的临床工作经验的评估人(一个神经系放射学家和三个放射科医师),在DWI图像上使用内部开发的软件人工勾画DWI损伤轮廓。评估人看不到任何强度数值范围和脑图。自动确定的DWI损伤体积与人工勾画的在三个评估人之间一致的DWI损伤共识图相比。为量化自动方法和人工勾画损伤轮廓之间的几何和体积相似度,我们使用体积相关和Dice系数(DC)。DC(定义:两个蒙片覆盖蒙片的组合集的两次交集),DC是一个0-1之间的数,其中1表示两个相同形状之间完全重叠,0表示无重叠。DWI损伤的体积和几何对比,在人工勾画的PWI损伤中完成。分类性能根据灵敏度和特异性评估,且在在同侧半球上确定。
结果
测定体积比较
图31显示由自动方法和人工勾画轮廓估算的DWI损伤的体积比较。
图32显示由自动方法和人工勾画轮廓估算的半暗带的体积比较。
对于图31-32,各自的子图A-D表示:A)APS蒙片和人工勾画轮廓之间的体积相关(x轴表示人工蒙片体积[mL],y轴表示APS蒙片体积[mL]);B)人工勾画轮廓与APS相比的Bland-Altman图(x轴表示平均损伤体积[mL],y轴表示差值);C)STM与人工勾画轮廓之间的体积相关(x轴表示人工蒙片体积[mL],y轴表示STM蒙片体积[mL])以及D)人工勾画轮廓与STM相比的Bland-Altman图(x轴表示平均损伤体积[mL],y轴表示差值)。实线表示两个比较的方法之间的平均值,而虚线表示平均差±2SD。
关于自动方法和人工勾画轮廓估算的DWI损伤的体积比较(比较,图31),对于APS与人工勾画轮廓相比,斯皮尔曼相关是R2=0.95。人工与APS勾画轮廓之间的平均差为-1.2mL(2SD:10.3mL)表明由APS轻微高估,参见图31和上面的行(子图A-B)。对于STM技术R2=0.80,人工和STM勾画轮廓之间的平均差是-11.6mL(±2SD:37.1mL)表明由STM高估,参见图31和下面的行(子图C-D)。
关于自动方法和人工勾画轮廓估算的半暗带的体积比较(比较,图32),对于APS与人工勾画轮廓相比,斯皮尔曼相关是R2=0.88(参见图32和上面的行(子图A-B))。人工与APS勾画轮廓之间的平均差为-0.2mL(±2SD:65.5mL),表明由APS较小地高估。对于STM技术R2=0.67(参见图32和较低的行(子图C-D)),人工和STM勾画轮廓之间的平均差是-31.6mL(±2SD:108.2mL),表明由STM明显高估。
Dice系数和分类性能
DWI蒙片生成的APS的中值DC是0.65(25th-75th四分位数:0.44-0.79),中值灵敏度是0.83(25th-75th四分位数:0.53-0.95)和中值特异性是1(25th-75th四分位数:100–1)。
DWI蒙片生成的STM的中值DC是0.24(25th-75th四分位数:0-0.51),中值灵敏度是0.96(25th-75th四分位数:0-1)和中值特异性是1(25th-75th四分位数:1–1)。
图33显示了半暗带分类性能评估小结。对于APS半暗带蒙片中值Dice系数(DC)是0.68(25th-75th四分位数:0.45-0.78),中值灵敏度(Sens.)是0.78(25th-75th四分位数:0.58-0.90)和中值特异性(Spec.)是0.96(25th-75th四分位数:0.92–0.98)。
对于APS半暗带蒙片中值Dice系数(DC)是0.37(25th-75th四分位数:0.06-0.59),中值灵敏度是0.54(25th-75th四分位数:0.23-0.77)和中值特异性是0.95(25th-75th四分位数:0.90–0.98)。每个患者的中值处理时间是21.4秒(8.9-33.9),主要依赖于半暗带体积和被评估的切片的数量。
讨论
根据例7的DWI损伤分割算法的一个优势在于:从ADC图和DWI图中提取并发信息,用于DWI损伤的自动分割,可产生更好的分割。只有在ADC图上呈现高信号的损伤,和在DWI图上呈现高信号的损伤,被分割为DWI损伤。通过这种方法,算法避免分割“老的”损伤,其通常在DWI图和白质结构上呈现高信号,白质结构典型地在ADC图上呈现高信号。
综上所述,本发明涉及一种方法,该方法用于估算生物组织中的半暗带的半暗带尺寸的度量,其中分别通过灌注加权成像(PWI)和弥散加权成像(DWI),获得的第一和第二图像被分析;且其中第一图像的分析包括水平集方法的应用,第二图像的分析包括灰度形态运算的应用。在本发明进一步的实施例中,连通组件标记算法可以应用在第一和第二图像中的任何一个上。本发明进一步涉及***、一个计算机程序产品和相应方法的使用。
尽管结合详细说明的实施例描述本发明,无论如何不应解释为限制到提供的实例。本发明的范围由附随的权利要求书说明。在权利要求的上下文中,术语“包含”或“包括”不排除其他的要素或步骤。而且,提到引用例如“一个”不应解释为排除多个。权利要求中关于附图中表明的部件的标号的使用,也不应解释为对本发明范围的限制。此外,不同权利要求中提到的每个特征可方便地组合,并且在不同权利要求中提到这些特征,不排除特征的组合是不可能的和无益的。
附件1
Figure BDA00003530197200391
表I
Figure BDA00003530197200392
表II
Figure BDA00003530197200393
表III
附件2
表IV
Figure BDA00003530197200402
表V
Figure BDA00003530197200403
表VI

Claims (21)

1.一种估算生物组织中半暗带的半暗带尺寸的度量的方法,所述方法包括以下步骤:
—接收第一图像,所述第一图像包括表示空间分辨灌注参数的数据;
—接收第二图像,所述第二图像包括表示空间分辨弥散参数的数据;
—在所述第一图像上识别第一区域,所述第一区域对应于灌注损伤;
—在所述第二图像上识别第二区域,所述第二区域对应于弥散损伤;
—基于第一和第二区域估算半暗带尺寸的度量;
其中所述方法进一步包括:
—在所述灌注加权图像上应用第一图像处理方法,所述第一图像处理方法包括水平集方法;
—在所述弥散加权图像上应用第二图像处理方法,所述第二图像处理方法包括灰度形态运算。
2.根据权利要求1所述的方法,所述方法进一步包括在选自所述第一图像和所述第二图像的任一个图像上应用连通组件标记。
3.根据权利要求1所述的方法,所述方法进一步包括在选自所述第一图像和所述第二图像的任一个图像上应用形态二值重构。
4.根据权利要求1所述的方法,其中所述第二图像处理方法进一步包括:
—接收由所述第一图像处理方法确定的参数;
其中识别所述第二区域的步骤是基于所述参数的。
5.根据权利要求1所述的方法,其中所述第一图像处理方法包括:
—在初级第一图像上识别初级蒙片;
且其中在所述第一图像上识别第一区域的步骤接收以下作为输入:
—次级第一图像,
—所述初级蒙片;
其中所述初级和所述次级第一图像是不相似的,由于所述次级第一图像已基于动脉输入函数被处理。
6.根据权利要求1所述的方法,其中所述第一图像处理方法进一步包括应用灰度形态运算。
7.根据权利要求1所述的方法,其中所述第二图像处理方法进一步包括应用水平集方法。
8.根据权利要求1所述的方法,其中所述第一图像处理方法的水平集方法是Mumford-Shah分割的水平集实现。
9.根据权利要求1所述的方法,所述方法进一步包括在选自所述第一图像和所述第二图像的任一个图像上应用二值形态运算。
10.根据权利要求9所述的方法,其中所述二值形态运算包括选自开运算、闭运算、填充的任一项图像处理技术。
11.根据权利要求1所述的方法,所述方法进一步包括在选自所述第一图像和所述第二图像的任一个图像上识别种子点。
12.根据权利要求11所述的方法,其中识别所述种子点的步骤包括在所述图像中识别初级点,所述初级点对应于所述图像中相对于所述图像中的其他点具有更高或更低的数据值的一个点。
13.根据权利要求11或12中任一项所述的方法,其中识别所述种子点的步骤包括损伤侧性检测。
14.根据权利要求13中任一项所述的方法,其中识别所述损伤侧性的步骤包括从多个不同成像模式中为每个半球提取一个或更多个特征。
15.根据权利要求1所述的方法,其中所述方法进一步包括:
—接收表观弥散系数图像,
且其中所述第一和/或第二图像处理方法进一步包括:
—在所述表观弥散系数图像中,识别高于或低于一个阈值的表观弥散系数区域。
16.一种用于估算受测对象的半暗带的半暗带尺寸的度量的***,所述***包括:
—用于接收第一图像的输入,所述第一图像包括表示空间分辨灌注参数的数据;
—用于接收第二图像的输入,所述第二图像包括表示空间分辨弥散参数的数据;
—设置为在所述第一图像上识别第一区域的图像处理单元,所述第一区域对应于灌注损伤;
—设置为在所述第二图像上识别第二区域的图像处理单元,所述第二区域对应于弥散损伤;以及
—基于所述第一和第二区域估算半暗带尺寸的度量;
其中所述方法进一步包括:
—设置为在所述灌注加权图像上应用第一图像处理方法的图像处理单元,所述第一图像处理方法包括水平集方法;以及
—设置为在所述弥散加权图像上应用第二图像处理方法的图像处理单元,所述第二图像处理方法包括灰度形态运算。
17.根据权利要求16所述的用于估算受测对象的半暗带的半暗带尺寸的度量的***,所述***进一步包括成像单元。
18.一种能够实现权利要求1所述的方法的计算机程序产品。
19.根据权利要求1所述的方法为估测生物组织中半暗带的半暗带尺寸的度量的使用,其中所述生物组织是脑组织。
20.根据权利要求1所述的方法,其中所述第一图像处理方法包括:
—在初级第一图像上识别初级蒙片;
—在所述第一图像上识别所述第一区域以作为次级第一图像上的次级蒙片,其中所述识别是基于所述次级第一图像和所述初级蒙片的;
其中所述初级第一图像没有为动脉输入函数而校正,而所述次级第一图像
为动脉输入函数而校正。
21.根据权利要求20所述的方法,其中在所述初级第一图像上识别初级蒙片的步骤包括在所述初级第一图像上应用水平集方法。
CN201180065488.0A 2010-12-17 2011-12-16 勾画组织损伤的方法 Active CN103339653B (zh)

Applications Claiming Priority (5)

Application Number Priority Date Filing Date Title
DKPA201070558 2010-12-17
DKPA201070558 2010-12-17
DKPA201170070 2011-02-07
DKPA201170070 2011-02-07
PCT/DK2011/050493 WO2012079593A1 (en) 2010-12-17 2011-12-16 Method for delineation of tissue lesions

Publications (2)

Publication Number Publication Date
CN103339653A true CN103339653A (zh) 2013-10-02
CN103339653B CN103339653B (zh) 2016-11-23

Family

ID=45390019

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201180065488.0A Active CN103339653B (zh) 2010-12-17 2011-12-16 勾画组织损伤的方法

Country Status (8)

Country Link
US (1) US9036878B2 (zh)
EP (1) EP2652704A1 (zh)
JP (1) JP6034299B2 (zh)
KR (1) KR20140033332A (zh)
CN (1) CN103339653B (zh)
AU (1) AU2011344876B2 (zh)
CA (1) CA2821395A1 (zh)
WO (1) WO2012079593A1 (zh)

Cited By (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107301645A (zh) * 2017-06-20 2017-10-27 上海联影医疗科技有限公司 一种数据处理方法和装置
CN107610115A (zh) * 2017-09-18 2018-01-19 彭建伟 一种直观显示缺血半暗带的图像处理方法
CN108122221A (zh) * 2016-11-29 2018-06-05 中国科学院深圳先进技术研究院 弥散加权成像图像中脑缺血区域的分割方法及装置
CN109726753A (zh) * 2018-12-25 2019-05-07 脑玺(上海)智能科技有限公司 基于时间信号曲线的灌注动态影像的分割方法及***
CN110910342A (zh) * 2018-09-12 2020-03-24 西门子医疗有限公司 通过使用深度学习来分析骨骼创伤
CN111373438A (zh) * 2017-10-17 2020-07-03 透视诊断有限公司 用于对器官成像的方法和装置
CN111407277A (zh) * 2020-03-06 2020-07-14 中国科学院武汉物理与数学研究所 一种急性缺血性脑卒中磁共振灌注-弥散影像配准方法
CN111986254A (zh) * 2020-08-21 2020-11-24 四川大学华西医院 一种靶区轮廓的分析方法、装置、存储介质及电子设备
WO2020237525A1 (zh) * 2019-05-29 2020-12-03 中国科学院重庆绿色智能技术研究院 基于相关系数的中风灌注成像病变区域检测***及方法
CN115797307A (zh) * 2022-12-13 2023-03-14 北京大学第三医院(北京大学第三临床医学院) 一种骨骼冠状位平衡参数检测***

Families Citing this family (17)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP2973405A4 (en) * 2013-03-15 2016-12-07 Seno Medical Instr Inc SYSTEM AND METHOD FOR SUPPORTING A DIAGNOSTIC VECTOR CLASSIFICATION
JP5946800B2 (ja) * 2013-07-22 2016-07-06 株式会社日立製作所 磁気共鳴イメージング装置、画像処理装置、画像処理方法、及び画像処理プログラム
EP3074949A2 (en) * 2013-11-27 2016-10-05 Universidad Politecnica De Madrid Method and system for determining the prognosis of a patient suffering from pulmonary embolism
JP6554717B2 (ja) * 2014-04-24 2019-08-07 ジーイー・メディカル・システムズ・グローバル・テクノロジー・カンパニー・エルエルシー 磁気共鳴装置およびプログラム
US9619882B2 (en) * 2014-06-19 2017-04-11 Alexander Sheung Lai Wong Correlated diffusion imaging system and method for identification of biological tissue of interest
EP2988270A1 (en) 2014-08-22 2016-02-24 Siemens Aktiengesellschaft Method for volume evaluation of penumbra mismatch in acute ischemic stroke and system therefor
KR20160118037A (ko) * 2015-04-01 2016-10-11 고려대학교 산학협력단 의료 영상으로부터 병변의 위치를 자동으로 감지하는 장치 및 그 방법
US9659368B2 (en) * 2015-05-15 2017-05-23 Beth Israel Deaconess Medical Center, Inc. System and method for enhancing functional medical images
KR101806841B1 (ko) * 2015-10-26 2017-12-11 주식회사 인피니트헬스케어 영상 기반 뇌졸중 병변 진단 시스템 및 방법
CN105654085A (zh) * 2015-12-31 2016-06-08 杭州晨鹰军泰科技有限公司 一种基于图像技术的弹孔识别方法
EP3547252A4 (en) * 2016-12-28 2019-12-04 Shanghai United Imaging Healthcare Co., Ltd. SYSTEM AND METHOD FOR PROCESSING MULTI-MODAL IMAGES
US10115197B1 (en) * 2017-06-06 2018-10-30 Imam Abdulrahman Bin Faisal University Apparatus and method for lesions segmentation
EP3493151A1 (en) * 2017-11-29 2019-06-05 Koninklijke Philips N.V. Combination of temporally resolved angiographic images with a spatially resolved angiographic image
JP7149516B2 (ja) * 2018-08-24 2022-10-07 富士通株式会社 検査支援プログラム、検査支援方法および検査支援装置
CN110362489B (zh) * 2019-07-15 2023-10-03 上海仪电(集团)有限公司中央研究院 一种测试用例生成方法、生成装置及软件自动测试***
CN110859622B (zh) * 2019-11-18 2023-07-04 东软医疗***股份有限公司 成像方法、装置及核磁***
US11436737B2 (en) * 2020-01-06 2022-09-06 Ricoh Company, Ltd. Brain segmentation using clustering and templates

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20080039895A1 (en) * 2006-04-11 2008-02-14 Northstar Neuroscience, Inc. Electromagnetic signal delivery for tissue affected by neuronal dysfunction, degradation, damage, and/or necrosis, and associated systems and methods
CN101596109A (zh) * 2009-06-12 2009-12-09 深圳先进技术研究院 获取脑部特征参数的方法、***及溶栓指征生成***与方法

Family Cites Families (14)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6994841B1 (en) * 1990-03-09 2006-02-07 Kucharczyk And Moseley Partners Relating to magnetic resonance imaging
DE102004005122A1 (de) * 2004-02-02 2005-08-18 Siemens Ag Verfahren und Vorrichtung zur verbesserten Auswertung von Magnetresonanzbildern
US7672540B2 (en) 2005-07-13 2010-03-02 Siemens Medical Solutions USA, Inc, Nonrigid registration of cardiac perfusion MR images using adaptive local template matching
US7702153B2 (en) 2005-10-07 2010-04-20 Siemens Medical Solutions Usa, Inc. Systems and methods for segmenting object of interest from medical image
WO2007046172A1 (ja) * 2005-10-20 2007-04-26 Niigata University 磁気共鳴画像処理方法および磁気共鳴画像処理装置
EP1952340B1 (en) 2005-11-21 2012-10-24 Agency for Science, Technology and Research Superimposing brain atlas images and brain images with delineation of infarct and penumbra for stroke diagnosis
WO2008000973A2 (fr) * 2006-06-29 2008-01-03 Centre National De La Recherche Scientifique - Cnrs Méthode d'estimation du potentiel de croissance des infarctus cérébraux
US8125223B2 (en) 2006-10-03 2012-02-28 Singapore Agency For Science, Technology And Research Act Segmenting infarct in diffusion-weighted imaging volume
US8064664B2 (en) * 2006-10-18 2011-11-22 Eigen, Inc. Alignment method for registering medical images
US8355552B2 (en) 2007-06-20 2013-01-15 The Trustees Of Columbia University In The City Of New York Automated determination of lymph nodes in scanned images
WO2009105530A2 (en) 2008-02-19 2009-08-27 The Trustees Of The University Of Pennsylvania System and method for automated segmentation, characterization, and classification of possibly malignant lesions and stratification of malignant tumors
JP2009247534A (ja) * 2008-04-04 2009-10-29 Ge Medical Systems Global Technology Co Llc 画像処理装置,磁気共鳴イメージング装置,および,画像処理方法
CN101332088B (zh) * 2008-07-21 2011-05-11 深圳先进技术研究院 一种获取脑部特征参数的***及溶栓决策***
KR101574376B1 (ko) * 2009-01-20 2015-12-03 케어스트림 헬스 인코포레이티드 우식을 탐지하는 방법 및 장치

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20080039895A1 (en) * 2006-04-11 2008-02-14 Northstar Neuroscience, Inc. Electromagnetic signal delivery for tissue affected by neuronal dysfunction, degradation, damage, and/or necrosis, and associated systems and methods
CN101596109A (zh) * 2009-06-12 2009-12-09 深圳先进技术研究院 获取脑部特征参数的方法、***及溶栓指征生成***与方法

Non-Patent Citations (7)

* Cited by examiner, † Cited by third party
Title
A.BARDERA等: "A framework to assist acute stroke diagnosis", 《VISION, MODELING, AND VISUALIZATION(VMV2005)》 *
A.BARDERA等: "A framework to assist acute stroke diagnosis", 《VISION, MODELING, AND VISUALIZATION(VMV2005)》, 18 November 2005 (2005-11-18), pages 1 - 8 *
NEIMARK M A等: "Sequential Acquisition and Processing of Perfusion and Diffusion MRI Data for a Porcine Stroke Model", 《ENGINEERING IN MEDICINE AND BIOLOGY SOCIETY, 2005. IEEE-EMBS 2005. 27T H ANNUAL INTERNATIONAL CONFERENCE OF THE SHANGHAI, CHINA》 *
SHEENA XIN LIU等: "Asymmetry Analysis in Rodent Cerebral Ischemia Models", 《ACADEMIC RADIOLOGY》 *
SUSHMITA DATTA等: "Segmentation of gadolinium-enhanced lesions on MRI in multiple sclerosis", 《JOURNAL OF MAGNETIC RESONANCE IMAGING》 *
TOMASZ HACHAJ等: "Automatic Detection and Lesion Description in Cerebral Blood Flow and Cerebral Blood Volume Perfusion Maps", 《JOURNAL OF SIGNAL PROCESSING SYSTEMS ; FOR SIGNAL, IMAGE, AND VIDEO TECHNOLOGY (FORMERLY THE JOURNAL OF VLSI SIGNAL PROCESSING SYSTEMS FOR SIGNAL, IMAGE, AND VIDEO TECHNOLOGY)》 *
刘国红等: "急性脑缺氧再灌注DWI及PWI的实验研究", 《中国医学影像学杂志》 *

Cited By (14)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108122221A (zh) * 2016-11-29 2018-06-05 中国科学院深圳先进技术研究院 弥散加权成像图像中脑缺血区域的分割方法及装置
CN108122221B (zh) * 2016-11-29 2020-11-17 中国科学院深圳先进技术研究院 弥散加权成像图像中脑缺血区域的分割方法及装置
CN107301645A (zh) * 2017-06-20 2017-10-27 上海联影医疗科技有限公司 一种数据处理方法和装置
CN107610115A (zh) * 2017-09-18 2018-01-19 彭建伟 一种直观显示缺血半暗带的图像处理方法
CN111373438A (zh) * 2017-10-17 2020-07-03 透视诊断有限公司 用于对器官成像的方法和装置
CN110910342A (zh) * 2018-09-12 2020-03-24 西门子医疗有限公司 通过使用深度学习来分析骨骼创伤
CN110910342B (zh) * 2018-09-12 2023-11-14 西门子医疗有限公司 通过使用深度学习来分析骨骼创伤
CN109726753A (zh) * 2018-12-25 2019-05-07 脑玺(上海)智能科技有限公司 基于时间信号曲线的灌注动态影像的分割方法及***
WO2020237525A1 (zh) * 2019-05-29 2020-12-03 中国科学院重庆绿色智能技术研究院 基于相关系数的中风灌注成像病变区域检测***及方法
CN111407277B (zh) * 2020-03-06 2022-04-26 中国科学院武汉物理与数学研究所 一种急性缺血性脑卒中磁共振灌注-弥散影像配准方法
CN111407277A (zh) * 2020-03-06 2020-07-14 中国科学院武汉物理与数学研究所 一种急性缺血性脑卒中磁共振灌注-弥散影像配准方法
CN111986254A (zh) * 2020-08-21 2020-11-24 四川大学华西医院 一种靶区轮廓的分析方法、装置、存储介质及电子设备
CN115797307A (zh) * 2022-12-13 2023-03-14 北京大学第三医院(北京大学第三临床医学院) 一种骨骼冠状位平衡参数检测***
CN115797307B (zh) * 2022-12-13 2023-08-08 北京大学第三医院(北京大学第三临床医学院) 一种骨骼冠状位平衡参数检测***

Also Published As

Publication number Publication date
EP2652704A1 (en) 2013-10-23
US9036878B2 (en) 2015-05-19
JP2014501126A (ja) 2014-01-20
AU2011344876B2 (en) 2017-05-25
AU2011344876A1 (en) 2013-08-01
WO2012079593A1 (en) 2012-06-21
JP6034299B2 (ja) 2016-11-30
US20130266197A1 (en) 2013-10-10
KR20140033332A (ko) 2014-03-18
CA2821395A1 (en) 2012-06-21
CN103339653B (zh) 2016-11-23

Similar Documents

Publication Publication Date Title
CN103339653A (zh) 勾画组织损伤的方法
Kumar et al. U-segnet: fully convolutional neural network based automated brain tissue segmentation tool
Prados et al. Spinal cord grey matter segmentation challenge
García-Lorenzo et al. Review of automatic segmentation methods of multiple sclerosis white matter lesions on conventional magnetic resonance imaging
Joliot et al. AICHA: An atlas of intrinsic connectivity of homotopic areas
Roy et al. Robust skull stripping using multiple MR image contrasts insensitive to pathology
Yoo et al. Deep learning of joint myelin and T1w MRI features in normal-appearing brain tissue to distinguish between multiple sclerosis patients and healthy controls
Devi et al. Neonatal brain MRI segmentation: A review
Wang et al. Principles and methods for automatic and semi-automatic tissue segmentation in MRI data
US6950544B2 (en) Automated measurement of anatomical structures in medical imaging
CN101111865A (zh) 用于在心脏图像中分割左心室的***和方法
US20150262372A1 (en) System and method for annotating images by propagating information
US20070167788A1 (en) Brain tissue classification
Ulas et al. DeepASL: Kinetic model incorporated loss for denoising arterial spin labeled MRI via deep residual learning
Akselrod-Ballin et al. Atlas guided identification of brain structures by combining 3D segmentation and SVM classification
Coupé et al. Nonlocal patch-based label fusion for hippocampus segmentation
US9730615B2 (en) Automated surface-based anatomical analysis based on atlas-based segmentation of medical imaging
Glaister et al. Thalamus segmentation using multi-modal feature classification: Validation and pilot study of an age-matched cohort
Qiu et al. 3D MR ventricle segmentation in pre-term infants with post-hemorrhagic ventricle dilatation (PHVD) using multi-phase geodesic level-sets
Lyu et al. Spectral-based automatic labeling and refining of human cortical sulcal curves using expert-provided examples
CN112233086A (zh) 基于脑区功能连接的fMRI数据分类识别方法及装置
Asman et al. Robust GM/WM segmentation of the spinal cord with iterative non-local statistical fusion
Studholme Dense feature deformation morphometry: Incorporating DTI data into conventional MRI morphometry
Gloger et al. Subject-Specific prior shape knowledge in feature-oriented probability maps for fully automatized liver segmentation in MR volume data
Glass et al. Hybrid artificial neural network segmentation of precise and accurate inversion recovery (PAIR) images from normal human brain☆

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant