用于分析表示生物组织的三维体积的图像数据的方法与流程

文档序号:11142391阅读:436来源:国知局
用于分析表示生物组织的三维体积的图像数据的方法与制造工艺

本发明涉及一种用于分析表示生物组织的三维体积的图像数据的方法,该图像数据包括表示第一时间点的体积的第一图像和表示不同于第一时间点的第二时间点的体积的第二图像。本发明还涉及用于分析表示生物组织的三维体积的图像数据的计算机程序。



背景技术:

发达国家中年龄超过50岁的人中,老年性黄斑变性(AMD)、尤其是新生血管AMD(nAMD)是失明的主要原因。血管透过性的增加导致视网膜内或视网膜下的异常流体聚集,该异常流体聚集在它涉及黄斑的中心时导致视觉机能障碍。这导致迅速恶化的视敏度,色素上皮疤痕,永久丧失视力或失明。

然而,抗血管生成剂、包括雷珠单抗(Ranibizumab,商品名为Novartis,巴塞尔,瑞士)的玻璃体内注射已经显示出显著改善的nAMD的过程。为了减少玻璃体内注射的负担并优化风险/利益曲线,nAMD特征的演变(progression)能通过光学相干断层扫描(OCT)进行非侵入性监测。突出的nAMD特征包括视网膜结构厚度的增加。当在视觉上比较不同时间获取的视网膜的相同区域的两个OCT图像时,这种增加可被识别,其中时间间隔是几天到几个月。

例如,用雷珠单抗治疗的患者通常经受每月一次OCT检查。如果观察到nAMD特征的显著增长,那么治疗决策表示为:病人在当天、一个月后或两个月后(治疗期)接受雷珠单抗注射。如果nAMD特征还没有完全退去,再治疗可以被指示为一个月后。否则,该患者当天不接受注射,但通常指示维持注射(维持期)。

OCT采集以及OCT采集的随后的分析通常是由技术人员、包括眼科医生进行。这意味着,所监控的患者被要求在每次采集OCT时要访问医院的医疗实践或专门单位。这使患者承受相当大的负担。而且,OCT采集的频率(例如1个月)已经是一方面nAMD发展的密切监测和另一方面患者负担的成本及负担之间的一种折衷。

原则上,OCT采集的自动分析可缓解这些问题。近年来,许多算法已被设计来自动分析视网膜的OCT采集,例如Abràmoff MD,Garvin M.,Sonka M.,Retinal的Retinal Imaging and Image Analysis(视网膜成像和图像分析),IEEE Rev Biomed Eng.2010;3:169-208。这些算法的大多数聚焦在视网膜的解剖结构的片段,如视网膜层。然而,特别是如果分析基于通过廉价的小型OCT设备所获得的稀疏的OCT信息,那么这些算法中没有一个被证明在自动提供成像的生物组织的形态学变化的信息方面、特别是层的厚度的变化方面是令人满意的。因此,仍然没有提供对于决定nAMD的进一步治疗有用的信息的自动方法。



技术实现要素:

本发明的目的是提供属于前述技术领域的用于分析图像数据的方法,该方法允许可靠地提供成像的生物组织的形态变化的信息,特别是图象数据稀疏的情况下。

通过权利要求1的特征详述本发明的方案。根据本发明,该方法包括以下步骤:

a)识别具有沿预定方向的主延伸部的第一图像所表示的体积的第一子体积和具有沿预定方向的主延伸部的第二图像所表示的体积的第二子体积,其中第一子体积和第二子体积表示三维体积的相同区域;

b)由在预定方向上的第一子体积的一系列体素(voxel)的数据图像与在预定方向上的第二子体积的一系列体素的图像数据创建距离矩阵;和

c)分析距离矩阵以获得沿预定方向的由第一子体积和第二子体积表示的生物组织的层的增长概率的至少一个局部测量值(measure)。

因此,用于分析表示生物组织的三维体积的图像数据的本发明的计算机程序包括在计算机运行时适于执行上述步骤的计算机程序代码。

应注意:生物组织的三维体积可用二维图片(即,二维网格的像素)来表示,示出仅一片组织。在这种情况下,“体素”实际上是像素并被认为表示围绕一片组织的相应位置的三维体积(例如,立方体)。但是,本发明的方法特别适合由三维网格的体素构成的图像数据。

特别地,预定的方向与三维图像数据的轴相一致。这有利于图像数据的处理。第一时间点和第二时间点的典型的时间间隔是12小时-150天,特别是5-90天。

应注意:不要求第一子体积和第二子体积表示被分析的生物组织的完全相同的区域,而是要求存在被两个子体积所涵盖的区域。它是将在随后的步骤中被研究的区域(或其中的一部分)。

距离矩阵包含在用于成对的体素(第一子体积的一个体素和第二子体积的一个体素)的图像数据方面的距离测量值(或成本函数(cost function))。在最简单的情况下,距离矩阵是二维的,或它可具有如下所示的更高的维度。应注意:通常不必要建立全部的矩阵,而是相对于预定方向具有非常不同的位置的第一子体积和第二子体积的体素不需要彼此比较,因为它们不对应于生物组织的相同部分,即使是在被比较的图像之间有重大变化的情况。

上述步骤的结果是在第一子体积和第二子体积所表示的生物组织内的层的增长概率的局部测量值,其可以是第一图像中已经存在的层或仅第二图像中的已经被检测的新的层。在本文中,“局部”是指彼此垂直且垂直于预定方向的两个轴主要限定的平面上的特定位置。

它表明,基于距离矩阵的该方法比传统的图像识别技术更具有鲁棒性(robust)并明显更简单,传统的图像识别技术试图识别生物组织的特定结构,以用于随后确定这些结构的尺寸。基于距离矩阵,组织的层的厚度变化可被自动且可靠地识别。

特别地,图像数据是光学相干断层成像(OCT)的图像数据。在OCT扫描,第一子体积和第二子体积的情况下,每者均对应于一定位置处的(超声术语中)A型扫描。预定方向对应于深度轴,即进行A型扫描所沿的轴。在下文中,深度轴用y表示,而表示对应的A型扫描的位置的垂直于y轴的两个轴分别用x和z表示。特别地,3D OCT图像数据由光谱域OCT获得。

通常,本发明的方法可应用于由诸如血管造影术、计算机断层扫描、磁共振成像(MRI)等其他技术获得的三维图像数据。

本发明的方法特别适合人类或哺乳动物的视网膜组织的研究。特别地,该方法允许视网膜、尤其是黄斑区的各层的增长的自动研究。相应地,处理过的数据有利于nAMD特征的检测,并提供对于决定nAMD的进一步治疗有用的信息。因此,该方法可用作决策支持工具,以便帮助临床医师评估玻璃体内注射的需要。

本发明的方法的目的不在于检测所有nAMD特征,即视网膜囊肿、视网膜下液,弥漫性视网膜水肿、视网膜色素上皮脱离和视网膜组织。相反,其检测3D OCT体积中的局部视网膜的增长。相比于基于图像识别的方法,该检测对图像质量的改变不敏感,因而更具有鲁棒性(robust)。

本发明不限于视网膜组织的研究。利用本发明的方法可研究其他生物组织,特别是这些组织的层的厚度随时间的变化(增加或减少)受关注的情况。一个示例是在MRI图像中观察到的阿尔茨海默氏病的演变。

优选地,分析距离矩阵的步骤包括获得对应于将距离矩阵的预定起点和预定终点相连接的路径的翘曲函数,距离矩阵给出的距离的总和沿该路径最小。通过采用寻路算法、例如已知的Dijkstra算法可确定最短距离。

该方法受到根据例如语音识别或诸如地震学领域中的时间顺序分析而已知的“动态时间规整”方法的启发。翘曲函数对应于第一子体积的体素和第二子体积的体素之间的最佳匹配。一般说来,第一子体积的每个体素与第二子体积的一个体素配对。如果第一子体积中的一个体素与第二子体积中的多个体素配对,那么这可能意味着在那个体素中记载的组织已经在两个采集之间增长或流体已经出现在体素的正下面或正上面。

相应地,图像数据表示的生物组织的各个层的厚度的变化将导致翘曲函数表示的路径的性状改变。因此翘曲函数的分析产生对于获得沿预定方向的生物组织的层的增长概率的局部测量值有用的信息。

在研究视网膜组织的情况下,沿预定方向的起点特别对应于视网膜的内界膜,终点将对应于脉络膜。

有利地,分析距离矩阵的步骤包括将该路径与基准路径比较的另一个子步骤,增长概率的局部测量值对应于路径和基准路径之间的总体差异。特别地,基准路径被选择为对应于在第一子体积的一系列体素的被比较图像数据和第二子体积的一系列的体素的图像数据之间没有变化的情况。在最简单的情况下,使用二维距离矩阵,基准路径是将起点与终点连接的直线,即基准路径显示没有翘曲。如果矩阵的列表示基准(第一)子体积的y位置,并且矩阵的行表示当前(第二)子体积的y位置,那么对角线从左上角(起点)到右下角(终点)延伸。具有在基准路径(对角线)下方的路径的翘曲函数表示检测到组织的层的增长。因此,在翘曲函数限定的路径和基准路径之间的区域是生物组织的层的厚度增大的测量值。

用于比较路径和基准路径之间的总体差异的不同的测量值可被应用。最优选的是与由路径和基准路径界定的区域的大小成比例的测量值。在二维距离矩阵的情况下,该区域是以任意单元测量的简单普通几何区域,在更高维度的距离矩阵的情况下必须使用超曲面的区域的测量值。局部测量值是被确定的区域乘以预定的比例常数。

特别根据距离矩阵的距离值的子集,其他局部测量值可包括加权。此外,除了路径和基准路径之间的区域,可使用其他标准,如路径和基准路径之间的长度差异或路径长度。

优选地,获得增长图,该增长图根据垂直于预定方向的平面中的位置由增长概率的局部测量值构成。在由(真正)三维网格体素构成三维图像数据的情况下,增长图将是二维的。这意味着,对于(x,z)平面中的很多点来说,如上所述地确定局部测量值,即,对于所有这些点来说,沿预定方向的增长概率的测量值在二维增长图中获得并被收集。如果图像数据是二维的,增长图将是一维的。增长图对于进一步分析来说是非常有用的基础。

特别地,沿预定方向由第一图像和第二图像表示的生物组织的层的增长概率的全局测量通过将机器学习算法应用于增长图而获得。利用机器学习算法允许系统地包括在各图像数据分析领域的专家的知识。机器学习算法可直接应用于增长图,即在图中储存的局部测量值,或者增长图的值可被预处理,例如,在应用机器学习算法之前进一步整理。

在优选的实施例中,该方法包括以下步骤:从增长图提取全局增长特征向量,并通过将支持向量机器算法应用于全局增长特征向量以便对全局增长概率分类。与提供主要沿给定(x,z)坐标处开始的直线的增长概率的局部测量值相反,全局增长特征向量涉及在研究的生物组织的整个体积内的增长概率。提取全局增长特征向量是机器学习算法的预处理步骤,并且机器学习算法应用于全局增长特征向量,以便将全局增长概率分类,即全局增长特征向量通过支持向量分类(SVC)模型被映射到最终的全局增长概率。可能的分类可简单地是"无增长"或"增长"(二进制),或者可定义中间类。进一步的分析和/或治疗方面的决定可基于分类由临床医师做出。

优选地,全局增长特征向量包括局部测量值的至少一个平均值和局部测量值的空间自相关的至少一个测量值。从广义上理解术语“平均”,因为归一化不被严格要求,在本文中“平均’也可简单地是局部测量值的总和、局部测量值的平方和或类似。空间自相关的合适的测量值例如是莫兰指数(Moran's I)或吉尔里指数(Geary's C)。后者能将分散或随机的局部增长(这很可能是由于噪声)与nAMD特征出现或增长时通常观察到的群生局部增长进行区分。

已经证实:通过机器学习算法,平均值以及空间自相关的测量值的组合实现对全局增长概率的可靠的分类。在优选的实施例中,全局增长特征向量由以下四个量构成:

a)局部增长分数的总和;

b)局部增长分数的平方和;

c)Moran's I;以及

d)Geary's C。

全局增长特征向量可以包括其他的或额外的量,只要他们被选择为提供有关研究的生物组织的层的增长概率的信息。

在优选的实施例中,第一子体积具有沿垂直于预定方向的至少一个方向的延伸部,该延伸部大于沿第二子体积的至少一个方向的延伸部。相应地,距离矩阵是至少三维的,其包括与第二子体积中沿预定方向的位置相关的第一维度和第一子体积中沿垂直于预定方向的至少一个方向以及沿预定方向的位置相关的至少两个其他维度。更优选地,距离矩阵是四维的,其包括第二子体积中沿预定方向的位置相关的第一维度和第一子体积中沿垂直于预定方向的两个方向以及沿预定方向的位置相关的三个其他维度。

如果确定了翘曲函数,相应的路径将具有至少一个额外的自由度。使用动态时间规整产生的算法的该3D扩展版允许补偿例如由一定的跟踪不准确引起的比较图像之间的小的偏移或倾斜或者相对于成形装置的生物组织的取向变化引起的轻微扭曲。在超声术语中,在位置(x,z)处的目标A型扫描(第二子体积)不仅仅与基准图像(第一子体积)中的位置(x,z)处的A型扫描匹配,还与多个邻近的A型扫描匹配,例如那些属于第一次提到的(x,z)A型扫描中心的薄圆柱体积。

这意味着这两个体积仍然需要大体记录。例如,在不同时期获取的视网膜的OCT图像的情况下,这能通过现在的设备实现。作为示例,采用Heidelberg的基于SLO的跟踪设备,获得两个体积,这两个体积沿垂直于眼睛的视轴的方向被充分地记录。通过翘曲算法的3D扩展的小残余记录误差被考虑。

利用该技术,这两个体积不必在深度方向、即沿y轴记录。通过将内界膜(ILM)分成两个体积,例如通过使用基于图切割的算法可获得沿该方向的记录。然后,所有的A型扫描沿y轴位移以便将ILM设置为y=0的坐标。

在最简单的情况下,第一体积和第二体积仅沿预定方向延伸,即他们对应于一系列的单个体素。这就要求在两个图像记录方面有一定的精度。而且,为避免由于小的残留记录误差引起结果的恶化,沿x和z方向的分辨率应足够精细,以避免相邻的(x,z)位置之间的层几何形状的突变。

有利地,有关第一子体积和第二子体积的图像数据的强度在创建距离矩阵之前被归一化。如果构成图像数据的体素的强度(亮度值)用于创建距离矩阵,那么这就是特别重要的。归一化确保通过本发明被彼此比较的两个图像部分在共同的基础上不被采集之间的照明和噪声变化影响。根据成像技术和分析的生物组织,归一化可被一次性应用到整个图像或单独应用到用来创建距离矩阵的每个子体积的图像数据。合适的归一化法是应用中值滤波器。

如果量不用于创建对图像强度的恒定变换(constant shift)不敏感的距离矩阵或者如果被提供的图像数据已经在图像采集和预处理的构架(framework)中被归一化,那么归一化可以是不必要的。

在优选的实施例中,距离矩阵基于沿预定方向的第一子体积的一系列体素的图像数据的强度的局部梯度和第二子体积的一系列体素的图像数据的强度的局部梯度而创建,即用于建立距离矩阵的距离测量值基于这些局部梯度的差异。优选在前述的归一化之后,梯度滤波增强例如主视网膜特征的研究的生物组织的特征的边界(在nAMD特征的视网膜层、边界等之间的界面)。因此,局部梯度构成用于形成鲁棒性距离矩阵的基础,该距离矩阵很大程度上不受全局范围上影响图像的影响因素影响。可行测量值是沿y轴的局部强度梯度的绝对差异。

替代性地,距离测量值基于其他量特别是强度本身,优选是在归一化之后。可行测量值是归一化的强度的绝对差异。原则上,它也能借助结合对应的距离矩阵而通过加权求和来结合两个或更多个不同的距离测量值。

根据下面的详细说明和全部权利要求,其他有利的实施方式和特征的组合变得明显。

附图说明

用于说明实施例的附图示出︰

图1由两个不同时间点的图像数据所表示的生物组织的体积的示意图,该体积包括待用本发明的方法分析的子体积;

图2a)强度梯度的基于体素的距离;b)对应于基于体素的距离的距离最短的路径;c)距离最短的路径和参考路径之间的区域;

图3增长图的示意图,包括根据垂直于深度轴线的平面上的位置的增长概率的局部测量值;

图4由两个不同的时间点的图像数据所表示的生物组织的体积的示意图,该体积包括用本发明的方法的变型分析的子体积;以及

图5根据本发明的方法的流程图。

在附图中,相同的部件被赋予相同的附图标记。

具体实施方式

图1是由两个不同的时间点的图像数据表示的生物组织的体积的示意图,该体积包括待用本发明的方法分析的子体积。在下文中,参照从人的视网膜黄斑区域取得的由光谱域OCT获得的3D OCT图像描述该方法。然而,该方法的基本步骤可以应用到生物组织等的其它三维图像。根据本发明的方法的流程图在图5中给出。

首先,获得实质上表示在不同时间点、例如30天的时间间隔获取的给定组织的相同体积的两个图像10、20(步骤101)。图像由布置为规则网格的多个体素组成,沿网格的方向分别用x、y和z表示。两个图像10、20沿垂直于眼睛的视轴的方向,即x-z平面被大体记录,这可通过用现在的OCT采集设备、例如借助那些具有眼底SLO跟踪特征的设备(由Heidelberg engineering公司开发)实现。Y轴是深度轴并平行于眼睛的视轴延伸。

在给定的示例中,这两个体积由49x49x496(即,较稀疏)的体素构成。应该注意:这两个体积不必在深度方向、即沿y轴记录。因此,内界膜(ILM)使用基于图割的算法被分成两个体积,该基于图割的算法在Dufour PA,Ceklic L,Abdillahi H等人的Graph-based multi-surface segmentation of OCT data using trained hard and soft constraints(利用受过训练的硬性和软性约束的OCT数据的基于图的多面分割),IEEE Trans Med Imaging.2013;32:531-43中描述。其他已知的算法是适用的。然后,所有的A型扫描沿y轴位移以便将ILM设置为y=0的坐标。

在将来自连续采集的两个A型扫描对齐时,由于采集之间的照明和噪声变化,体素强度的直接比较是不可靠的。作为第一预处理步骤102,单位半径的中值滤波器因而被用于去噪。作为第二个预处理步骤103,梯度滤波器用于增强主视网膜特征的边界(视网膜层之间的界面、nAMD特征的边界等)。

3D图像是每个A型扫描的比较A型扫描,即获得由第一图像10表示的体积中在坐标(x,z)处的A型扫描11(步骤104)和获得由第二图像20表示的体积中在坐标(x,z)处的A型扫描21(步骤105)。两个A型扫描11、21是图像10、20所表示的体积的子体积并由具有坐标(x,y,z)的线性的一系列体素构成,其中x和z是常数(彼此比较的子体积是相同的),即子体积的A型扫描11、21的主延伸部是沿着y方向的。A型扫描11、21然后互相比较。

基于滤波的图像,对于每一对在各自(x,z)坐标处的A型扫描来说,相强度梯度基于体素的距离被确定(步骤106),即基准A-扫描11中的每个体素与目标A型扫描21中的每个体素相比较,除非它们的y坐标太远。由此产生的距离矩阵30在图2的a)中示出,其中,对于基准A型扫描11来说轴沿y坐标(y1轴12)延伸并对于目标A型扫描21来说沿y坐标(y2轴22)延伸。在矩阵角部处黑色三角形区域31、32是由于|y1-y2|超过一定阈值的区域的排除。在距离矩阵30中的亮线33对应于视网膜特征的增强边界。

接下来,类似于信号处理中的动态时间规整(DTW)的序列对齐策略应用到待比较的两个A型扫描11、21。基于基于体素的距离矩阵30,利用在与基于体素的距离矩阵30相关的曲线中的已知的Dijkstra算法,发现从ILM 34(左下角)到脉络膜35(右上角)的两个A型扫描11、21之间的最小失真路径40(步骤107)。产生的路径40在图2的b)中以与图2的a)相同的坐标示出。它是从ILM 34到脉络膜35的路径,其中沿该路径由距离矩阵30给出的基于体素的距离值的总和最小。

如果在第一A型扫描中的一个体素与第二A型扫描中的多个体素配对,这可能意味着,在那个体素中记载的组织已经在两个采集之间增长,或者流体已经出现在该体素的正下方或正上方。如在图2的b)中能看到的,这种情况在路径40的区域41中给出。应注意,梯度滤波器的使用迫使主视网膜层在最小扭曲路径40中被对齐。

沿距离矩阵30的对角线36延伸的路径对应于相同图像(没改变)的情况。在对角线36上方延伸的最小失真路径40表示检测到视网膜的增长。因此,在对角线36上方和最小失真路径40下方的区域42用作局部视网膜增长分数,见图2的c),其再次使用了与图2的a)和2的b)相同的坐标(步骤108)。

只要不是所有(x,z)坐标已经被分析,重复上述步骤(决定109)。

其次,基于对于每个(x,z)的局部视网膜增长分数,通过重复上述步骤,建立增长图50(步骤110)。增长图50基本上表示作为x和z坐标、即根据垂直于眼睛的视轴的平面中的位置的函数的局部视网膜增长分数s(x,z)的值。在图3中,增长图50被示意性地示出。应该注意,示出的数字只是为了说明的目的,通常,分数也不是整数,而只是指示上一步骤中确定的区域42的值。

下一步,连续体积的全局增长概率被计算(步骤111)。为此目的,全局增长特征向量从视网膜增长图50提取。该特征向量由以下组成:

其中,利用总和或平方和,第一两个特征使局部增长分数s(x,z)(由视网膜增长图给出)平均。最后两个特征是空间自相关特性特征,即Moran's I和Geary's C。后两个特征能将分散或随机的局部增长(这很可能是由于噪声)与nAMD特征出现或增长时通常获得的群生局部增长相区分。

最后,利用支撑向量分类(SVC)模型,该全局增长特征向量被映射到全局增长概率(步骤112)。SVC也对全局增长概率提供最优截止(cut-off)。因此,该算法将产生二进制输出︰“增长”或者“没有增长”。如果输出(步骤113)是“增长”,保健专业人士应进行额外的检查并最终作出治疗决定。

利用由熟练的专业眼科医生评估的OCT图像的训练数据集,SVC算法已经被训练。对SVC模型已经采用标准训练策略:通过两重交叉验证,已经发现参数值(最佳软边缘参数C、最优核函数和最优核函数参数)。然后,利用参数的最优设定,分级者已经在整个训练集中重新训练。在训练集中已有象(x,z)位置那样多的训练样本。

图4是由两个不同的时间点的图像数据表示的生物组织的体积的示意图,该体积包括由本发明的方法的变型分析的子体积。

为了处理在垂直于眼睛的视轴的平面中的小的(残余)记录误差,在第二图像220的位置(x,z)处的目标A-扫描221不仅仅与早前的OCT图像210中位置(x,z)处的A型扫描匹配,而还与那个(x,z)A型扫描中心的薄的准圆柱体积211匹配。最小失真路径的搜索是类似的,除了它以四维的基于体素的矩阵(来自目标图像220的一个维度和来自早前的图像210的三个维度)进行操作,而不是二维矩阵。

本发明并不限于上述的方法。如上所述,该方法的一些步骤可由替代算法取代。此外,该方法是适用于其它组织和源于其他成像技术的图像数据。

总之,应当注意:本发明创建了用于图像数据分析的方法,该方法允许可靠地提供被成像的生物组织的厚度变化的信息。

当前第1页1 2 3 
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1