CN100583149C - 医学图像中器官的保形分段 - Google Patents

Abstract

公开了一种用于对图像中的对象的轮廓进行分段的方法和设备。它特别地应用于医学图像中描绘的身体器官或部分的分段。含有至少一个对象的输入图像包括至少二维的像素数据组。从输入图像得到一个经边缘检测的图像。在所述经边缘检测的图像中选择标记，并被指定以相应的固定电位值。生成了一个电位图以作为对在具有根据经边缘检测的图像的电阻率的电阻性媒体中静电问题的解，其中各标记规定了相应的固定电位值的电极。根据生成的电位图通过例如阈值来估计对象轮廓。

Description

医学图像中器官的保形分段

本发明涉及图像的分段。更特定地，本发明提出了一种用于对特别是医学图像的数字图像中描绘的不同的和分立的对象的一大类边界进行分段的鲁棒的过程。

这样的分段技术处理一个数字图像以便检测、分类、和枚举其中所描绘的分立的对象。它包括在感兴趣的区域内确定对象的轮廓，即，它们的轮廓或边界，这对于对象的识别和定位以及分析它们的形状、形式、尺寸、运动等等是有用的。

分段是在图像处理中最具有挑战性的任务之一，因为在某种程度上它需要对图像进行语义理解。它也表示了一种困难的问题，因为数字图像通常缺乏足够的信息来把分段问题的可能的解决方案限制在包括了正确的解的解决方案的一个小的集合中。

图像分段在医学图像领域找到普遍的应用，尤其是计算X射线层析摄影术(CT)图像、X射线图像、超声(US)图像、磁共振(MR)图像等等。特别希望对在这样的医学图像中可以看见的各种组织对象(诸如心脏、前列腺、肾脏、肝脏、胰腺等等)进行定位。一个组织对象相对于周围事物的位置可被用来规划和执行医学手术，诸如外科手术、用于癌症的放射治疗等等。

图像分段是在数字形式的医学图像上进行的。诸如人的身体的一部分那样的目标的数字图像是一个包括一系列数据单元的数据组，每个数据单元具有相应于目标特性的数字数据值。该特性可以由图像传感器在有规律的时间段内在图像传感器的视场上测量。它也可以按照像素网格根据投影数据进行计算。数据值所对应的特性可以是黑白照片的光的强度、彩色图像的分开的RGB分量、X射线衰减系数、MR中的氢含量等。

在已知的分段技术中，流域变换对交互的和自动的分段提供了方便的和有力的工具。这种经典的“形态学的分段”可以认为是淹没被看成为地形表面的梯度图像。首先，构建梯度图像。其次，对于感兴趣的对象，自动地或交互地规定里面的标记。最后，构建与标记相联系的流域，对地形表面逐渐地淹没(引入时间参数)。与感兴趣区域中标记的形状和布置无关的结果可通过对蕴含表面地形和完整的标记组这两者的全局最小化而得到。

在“Multiscale Morphological Segmentation based onWatershed，Flooding and Eikonal PDE(基于流域、淹没和镜像PDE的多尺度形态学的分段)”，作者F.Meyer和P.Maragos，Scale-Space’99，LNCS 1682，pp.351-362中，根据用于对被称为镜象方程(Eikonal Equation)的偏微分方程(PDE)进行求解的有效方法，提出流域变换的一个变例。它允许计算被确定的自适应超度量距离，以使得由一条边缘所分开的各点被认为比没有被该边缘分开的点分隔得更远。

当把这样的技术用于医学器官分段时，在部分的对象边界具有弱的对比度时会遇到限制。当处理要被分段的多个对象时这些缺点是很明显。该技术可导致丢失边缘所引起的泄漏的弱边缘，或一个“泄漏”的分段，因为某些对象没有被看作为闭合的对象。图1a显示一个具有由它们相应的边缘(经边缘检测的图像)表示的两个同心的特征的图像的例子，即，一个具有闭合的边缘的星形外部特征10，和一个完全被包含在星形外部特征10里面的圆形的特征20。由于缺乏对比度，不能检测到圆形特征20的小的边缘部分22。被应用于图1a的图像的上述的流域技术导致图1b的结果，在其中只确定了外部特征10。的确，圆形特征20丢失，因为对于外部特征10构建的流域将通过小的边缘部分22“泄漏”到圆形特征20中。

本发明的目的是提供鲁棒的避免流域变换的上述的缺点的分段技术。本发明的另一个目的是提供具有高精度边界确定的和能够处理具有不同的对比度的多个边界的分段技术。

因此，本发明提供按照权利要求1的分段设备和包括这样的设备的医学成像系统。本发明还提供按照权利要求14的方法和按照权利要求15的计算机程序产品。

该分段问题被变换成等价的稳态静电流问题，其中标记的固定电位被用作为边界条件。本发明利用所谓的保形映射，即一种归一化的图像数据自适应拉普拉斯偏微分方程(PDE)的解。它避免当一些边缘部分丢失时的由边缘引起的泄漏，并导致平滑的内插形状，从而在边缘对比度局部减弱或不存在时使边缘闭合。所以，可以达到高精度的边界确定，以及本发明的一个或多个实施例的分段技术可以处理具有不同对比度的多个边界。

当结合附图考虑此后的说明时，将明白本发明的其它特性和优点，其中：

图1a是图像中两个同心特征的示意图，

图1b是使用已知的流域变换的图1a的特征的最终分段的示意图，

图2是按照本发明的分段方法的流程图，

图3是使用按照本发明的分段方法根据图1a的特征生成的电位映射的示意图，

图4是使用按照本发明的分段方法时图1b的特征的最终分段的示意图，

图5是按照本发明的教导编程的通用计算机的示意图。

本发明处理图像上的对象的分段。虽然本发明是以软件实施方案显示的，但它也可用例如在计算机系统的图形卡中的硬件部件来实现。

现在参照附图(更特别地是图2)，图上显示按照本发明的示例性分段方法的流程图。总的方案包括在步骤100初始地获取一个包含要被分段的对象的数字医学2D或3D图像。如需要的话，这个步骤可包括使用文件变换器部件来把图像从一个文件格式变换成另一个必须的格式。最终得到的输入图像此后被称为M，而p是在图像内的像素标号。下面，图像和它的数据将用同一个名称表示，因此，M既是指输入图像，也是指像素p的输入数据。此后的方法也可以容易地推广到维数大于3D的图像。

在第二步骤110，输入图像M被进行处理以便提供一个经边缘检测的图像ED(p)。经边缘检测的图像ED(p)是通过使用技术上已知的方法创建的，其诸如局部变量方法或梯度强度方法。初始输入数据M经受到边缘检测，这样，它的边缘被加以检测从而确定边缘强度数据ED(p)，以便区分一个对象的边缘区域与其它区域。替换地，输入图像M可以首先由适当的技术进行边缘鲜明度和特征增强，以产生具有增强鲜明度的图像。

在经边缘检测的图像中的像素值ED(p)说明ROI中的特征。它们表示特征凸出量，后者可以是像素强度值，或与图像M中特征强度有关的任何适当的数据。在像素强度中的局部梯度是在经边缘检测的图像中要被使用的最典型的量。

经边缘检测的图像ED(p)包括一定的数目的不同的对象。在此后的说明中，除非另外提出，将以两个对象作为说明。本领域技术人员可以很容易地把它推广到大的数目。

经边缘检测的图像ED(p)可被修改成把在感兴趣区域(ROI)以外(即，不太可能是器官轮廓的地方)的边缘强度数据设置为零。

为了藉助于求解静电问题把对象分段，对经边缘检测的图像ED(p)或ROI的每个像素p规定一个电位函数Φ。按照本发明的分段方法使用从分段问题到静电问题的转换，它的解是电位函数Φ(p)。

在步骤120，在经边缘检测的图像ED(p)上选择标记。固定的电位值被定义和被附着到不同的标记上。这些电位值被用作为静电问题的固定的边界条件。

对于经边缘检测的图像ED(p)的每个被识别的对象，需要内部标记Min和外部标记Mout。这些内部和外部标记可被看作为与到每个对象相联系的内部和外部边界。对于每个对象，内部边界可仅限于在这个对象里面的任何地方的一个点，而外部边界可以是图像的外部边沿，或处在这个对象与相邻的对象之间的空间。内部和外部标记可以由用户(例如，通过鼠标或类似的指点设备)输入到计算机系统，或被自动设置。有关提取形状的现有的知识可以通过给出预期的对象轮廓的内部和外部标记大体形状而被引入到按照本发明的方法，以使算法更加鲁棒。

对于每个对象，里面的标记Min作为具有电位Vin的电极对待，类似地，外面的标记Mout作为具有不同的电位Vout的电极对待。电位Vin和Vout的数值是任意的，例如可以是-100伏和+100伏。

在步骤130，使用ED(p)来生成电阻率图，即，电阻率分布(或，等效地，电导率分布)被映射到经边缘检测的图像ED(p)上：每个像素p根据它在经边缘检测的图像中的数值接收一个电阻率数值ρ(p)。典型地，电阻率数值可以等于经边缘检测的图像中的像素强度值，在假设的对象边缘上比在介入的介质中具有更高的数值。电介质的电阻率ρ(p)通常被取为边缘强度值ED(p)的递增函数，它绝不会取小于正的最小值的数值。

在优选实施例中，电阻性介质电导率γ(p)＝1/ρ(p)(在2D中每单位平方，或在3D中每单位立方)被取为其中A是被选择为在边缘检测图像ED(p)中预期噪声电平的量级的参数，换句话说，在像素p处的电阻率ρ(p)正比于A²+ED(p)²。

在步骤140，求解2D静电问题。给定2D电阻率分布和在标记位置处的电位值Vin、Vout后，可以使用已知的数值分析方法来求解稳态静电流方程，以得出电位分布。

在流域变换技术中，“流问题”是通过淹没被看作为地形表面的梯度图像而被解决的。要求解的偏微分方程(PDE)是依赖于时间的，因为淹没是逐步盖满地形表面的。这里，按照本发明的方法的概念在于求解一个静态稳电流问题，它在其电阻率与边缘检测图像数据相联系的电阻性介质上被规定。这种介质包括对于设置的边界条件先前作为电极而被规定的标记。步骤140为边缘检测图像ED(p)的每个像素生成电位函数Φ或保形电位映射，以作为对于上述的静电流问题的解。

令J是在电阻性媒体内的电流密度向量。J(p)，Φ(p)和γ(p)＝1/ρ(p)是像素p的函数，如前所述，p指的是在初始图像M或经边缘检测的图像ED(p)中的位置。用

表示梯度向量算符，这三个量通过以下关系式相联系：

$J\left(p\right)=\mathrm{\γ}\left(p\right)\·\▿\mathrm{\Φ}\left(p\right)---\left(1\right)$

由于假设在感兴趣区域(ROI)中没有源与接收器，向量J(p)满足零散度(divergence)条件：

$\▿\·J\left(p\right)=0---\left(2\right)$

其中

代表散度算符，必须求解通过把来自公式(1)的J(p)的表达式代入公式(2)而得到的偏微分方程(PDE)。

$\▿\·(\mathrm{\γ}\left(p\right)\·\▿\mathrm{\Φ}\left(p\right))=0---\left(3\right)$

这是一个包括了在ROI内的电位函数Φ(p)的一个未知数的PDE，它受到在内部和外部标记上规定的边界条件的限制。

在ROI中没有包括电流源或接收器(sink)。寻求电阻性介质的电位分布可导致想要的保形电位图。这个电位图产生在对象边界处的快速改变，而在边缘之间缓慢变化。当部分边缘丢失时，按照本发明的方法允许在丢失的边缘部分的邻居中平滑插入电位。

可以使用求解电问题的多网格算法。例如，可以实施有限差分法。也可以使用其它现有技术的方法。

在某些事例中，按照本发明的方法可以导致电位图的解，它们显示在围绕选择的内部标记各点处电位会快速改变，这会改变分段的寻求。一种方法必须被设计成使得在ROI的边缘点处而不是在内部标记处发生快速电位改变。小的内部标记等价于小的内部电极，并且可能加剧在这些标记附近这样的不想要的快速电位改变。

在按照本发明的为避免在标记附近的不想要的快速电位改变(产生不想要的伪像)的方法的优选实施例中，实施迭代算法，以便得到更精确的电位图。从实际的观点看来，在第五步骤140后，即，在电位图生成后，在附加步骤145使用生成的电位图Φ(p)来更新电导率γ(p)。在本发明的优选实施例中，修改的数值γ^m(p)是：

${\mathrm{\γ}}^m\left(p\right)=B\·\mathrm{\γ}{\left(p\right)}^2\·|\▿\mathrm{\Φ}\left(p\right)|---\left(4\right)$

其中γ(p)是在用来生成电位图Φ(p)的以前的步骤中使用的电导率，以及B是任意常数。

这相应于电导率γ(p)的重新归一化，并且是在ROI上执行的。它允许局部弛豫操作，以避免某些标记的电位吸引。

相同的标记步骤12被保持，在规定了在具有更新的电阻率γ^m(p)的电阻性介质中静电问题后，生成了修改的电位图以作为对于新的静电问题的解。

为了收敛到稳定的电位图，可能需要使用按照(4)的更新的电导率进行几次迭代。

除了在边界电极里面或在边界电极外面以外，这个电位图的解没有显示最小值或最大值。图3显示最终得到的收敛的电位图，它可以根据图1a的特征生成，不同的灰度代表各电位值。来自图1a的初始边缘检测图像的两个对象的形状都出现在电位图上。丢失的边缘部分22相应于平滑的电位变化的区域。中心的深色圆23相应于其中实行了弛豫操作的区域。

在图2上由步骤146(选择非线性任选项)表示的替换实施例中，方程(2)的线性静电PDE可以由非线性PDE代替，其中电导率γ(p)本身是电位图Φ(p)的函数。这个任选项或者紧接在步骤130，即，电导率图生成后被选择，或者在步骤140的电位图生成后，即，在迭代算法期间被生成。

在优选实施例中，γ(p)由非线性函数γ^nl(p)代替：

${\mathrm{\γ}}^{\mathrm{nl}}\left(p\right)=\frac{\mathrm{\γ}\left(p\right)}{{(1+{|\▿\mathrm{\Φ}\left(p\right)|}^2)}^{1/2}}---\left(5\right)$

当它代入到公式(2)时产生：

$\▿\·(\frac{\mathrm{\γ}\left(p\right)}{{(1+{|\▿\mathrm{\Φ}\left(p\right)|}^2)}^{1/2}}\·\▿\mathrm{\Φ}\left(p\right))=0---\left(6\right)$

这是一个非线性PDE，它包括一个作为处在ROI内的电位函数Φ(p)的一个未知数，受到在内部和外部标记上规定的边界条件的支配。

在这个步骤后面的是在分段方法中的第五步骤140，其中生成作为(6)的非线性静电问题解的电位图。

这个步骤是任选的。它可以紧接在电阻率图生成后实施以代替线性PDE(2)，并且不用以前描述的步骤145的电导率更新，以便也避免电极快速电位改变的伪像。

它也可以与步骤145的电导更新相关联，即，在求解第一线性PDE(2)而得出电位图Φ(p)后，电导率按照(4)被更新为γ^m(p)，然后要求解的新的静电问题是非线性PDE(6)。与前面描述相同的迭代算法然后结合非线性PDE(6)的使用以代替线性PDE(2)的使用而被实施。

在步骤150，由于在以前的步骤中计算的保形电位图Φ(p)，生成输入图像M的分段的轮廓。

保形电位图Φ(p)可以服从动态阈值算法(即，通过与位置有关的阈值的定义)。对于这个求阈值操作的一个解决方案包含为每个对象的每个内部像素p设置一个邻居关系，并且在这个邻居关系中根据ED(p)来确定哪个电位值Φ(p)相应于边缘强度数据的最高的和值(这可以通过适当的直方图技术完成)。当只有一个边界边缘要检测时，被这样确定了阈值的这个图可提供所需的对象分段。

当有多个器官边界需要检测时，这个算法需要一个在电位Φ递增的间隔内执行动态求阈值步骤。

也可以使用其它求阈值方法，例如在以前提到的F.Meyer和P.Maragos的论文中给出的用于流域变换的方法。

替换地，不同的器官边界可以按顺序地规定。例如，具有一个内部标记的分段方法的第一应用允许规定里面的轮廓。然后，这个最里面的轮廓可被作为新的内部标记以用于分段方法的第二应用，等等。

图4给出在求阈值操作已被应用于到图3所示的电位图以后的分段轮廓。来自图1的两个对象即星形对象10和圆形对象20都被识别。这时当图1a的丢失的边缘部分22完成时，圆形对象现在就是一个闭合的对象。

关于分段的图像的显示，在本领域中已知的技术可被用来观看结果。例如，在显示归类为某个灰度级的对象里面的所有像素时，同时把被归类为在所述对象外面的所有的像素设置成与以前的灰度级别极为不同的另一个灰度级别，从而使轮廓变得很明显。

按照本发明的方法的整体计算复杂性是低的。提出了使用分立电路近似的分级结构的静电问题的从粗糙到精细的非常有效的解决工具。对于收敛需要有限次数的迭代(等价于对于执行局部弛豫操作的每个内部点的约10次过程)。

本发明还提供用于把图像中的对象分段的设备，它包括获取装置，用来接收包含至少一个对象的输入图像，这个图像包括至少两维的像素数据组。按照本发明的设备还包括处理装置，它用来实施按照本发明的方法，以便根据输入图像提供边缘检测图像，以及生成作为对一个从分段问题变换成的静电问题的解的保形电位图。该设备还包括选择装置，用来选择在经边缘检测的图像上的至少一个标记，以及把标记与固定的静电位数值相关联。

本发明可以使用按照本申请的教导而编程的传统通用数字计算机或微处理器来方便地实施。该设备有利地是医学成像系统的一部分。

图5是按照本发明的计算机系统300的框图。计算机系统300可包括CPU(中央处理单元)310、存储器320、输入设备330、输入/输出传输通道340、和显示设备350。其它设备，如附加盘驱动器、附加存储器、网络连接...等等也可以包括在内，但没有被表示。

存储器320包括一个源文件，该文件包含具有要被分段的对象的输入图像M。存储器320还可包括计算机程序产品，它是要由CPU 310来执行的。这个程序包括执行以上描述的方法的指令。输入设备被用来接收来自用户的指令，以便例如选择ROI和/或在图像内的标记位置，和/或运行或不运行该方法的实施例的不同的阶段。输入/输出通道可被用来接收被存储在存储器320中的输入图像M，以及发送分段的图像(输出图像)到其它设备。显示设备可被用来显示输出图像，输出图像包括来自输入图像的最终得到的经分段的对象。

Claims (14)

1.用于把图像中的对象分段的设备，包括：

-输入装置，用于接收包含至少一个对象(10，20)的输入图像，所述图像包含至少二维的像素数据组；

-用于从所述输入图像生成经边缘检测的图像的装置；

-用于选择在所述经边缘检测的图像中的标记和用于给所述标记指定相应的固定电位值的装置；

-计算装置，用于生成电位图以作为对在其电阻率取决于所述经边缘检测的图像的电阻性介质中的静电流问题的解，其中的标记以相应的固定电位值来规定电极；以及

-用于根据所述电位图估计对象轮廓的装置。

2.如在权利要求1中要求的设备，还包括用于使用生成的电位图来更新所述电阻性介质的电阻率的装置，其中该计算装置被安排来生成一个另外的电位图以作为对在具有更新后的电阻率的电阻性介质中新的静电流问题的解。

3.如在权利要求2中要求的设备，其中在像素p处更新的电阻率正比于

其中ρ(p)是用来第一次生成电位图Φ(p)的所述电阻性介质的电阻率，其中Φ(p)是在所述经边缘检测的图像的像素p处规定的电位图，以及

代表梯度向量算符。

4.如在权利要求1中要求的设备，其中所述经边缘检测的图像是在包含所述对象的所述输入图像的感兴趣区域中被规定的。

5.如在权利要求1中要求的设备，其中用于选择标记的装置被安排成选择在所述对象里面的至少一个内部标记和至少一个在所述对象外面的外部标记，所述内部标记被指定第一固定的电位值，以及所述外部标记被指定第二固定的电位值。

6.如在权利要求1中要求的设备，其中在像素p处所述电阻率正比于A²+ED(p)²，其中ED(p)表示对所述像素p的经边缘检测的数据的值，以及A是在经边缘检测的图像中预期噪声电平的量级的常数。

7.如在权利要求1中要求的设备，其中所述静电流问题具有以下形式：

$\▿.(\mathrm{\γ}\left(p\right).\▿\mathrm{\Φ}\left(p\right))=0$

其中Φ(p)是在所述经边缘检测的图像的像素p处规定的所述电位图，γ(p)是所述电阻性介质的电导率，即电阻率的倒数，

代表梯度向量算符，以及

代表散度算符。

8.如权利要求7中要求的设备，其中在像素p处所述电阻率正比于A²+ED(p)²，其中ED(p)表示对所述像素p的经边缘检测的数据的值，以及A是在经边缘检测的图像中预期噪声电平的量级的常数。

9.如在权利要求1中要求的设备，其中所述静电流问题具有以下形式：

$\▿.(\frac{\mathrm{\γ}\left(p\right)}{{(1+{|\▿\mathrm{\Φ}\left(p\right)|}^2)}^{1/2}}.\▿\mathrm{\Φ}\left(p\right))=0$

其中Φ(p)是在所述经边缘检测的图像的像素p处规定的所述电位图，γ(p)是所述电阻性介质的电导率，即电阻率的倒数，

代表梯度向量算符，以及

代表散度算符。

10.如权利要求9中要求的设备，其中在像素p处所述电阻率正比于A²+ED(p)²，其中ED(p)表示对所述像素p的经边缘检测的数据的值，以及A是在经边缘检测的图像中预期噪声电平的量级的常数。

11.如在权利要求3中要求的设备，其中新的静电流问题具有以下形式：

$\▿.(\mathrm{\γ}\left(p\right).\▿\mathrm{\Φ}\left(p\right))=0\mathrm{or}\▿.(\frac{\mathrm{\γ}\left(p\right)}{{(1+{|\▿\mathrm{\Φ}\left(p\right)|}^2)}^{1/2}}.\▿\mathrm{\Φ}\left(p\right))=0$

其中γ(p)是所述电阻性介质的更新后的电导率，即所述更新后的电阻率的倒数，

代表梯度向量算符，以及代表散度算符。

12.如权利要求11中要求的设备，其中在像素p处所述电阻率正比于A²+ED(p)²，其中ED(p)表示对于所述像素p的经边缘检测的数据的值，以及A是在经边缘检测的图像中预期噪声电平的量级的常数。

13.一种医学成像系统，包括用于获取描绘身体器官的输入图像的装置，以及如在前述权利要求的任一项中要求的用于把所述输入图像中的对象(10，20)进行分段的设备。

14.一种用于把所述输入图像中的对象分段的方法，包括以下步骤：

-接收包含至少一个对象(10，20)的输入图像(100)，所述图像包含至少二维的像素数据组；

-根据所述输入图像生成经边缘检测的图像(110)；

-选择在所述经边缘检测的图像中的标记和用于给所述标记(120)指定相应的固定的电位值；

-生成电位图以作为对于在具有取决于经边缘检测的图像(130)的电阻率的电阻性介质中静电流问题的解(140)，其中的标记以相应的固定的电位值来规定电极；以及

-根据所述电位图估计对象轮廓(150)。

Images