CN101971213A - 图像中解剖结构分割与建模的方法及系统 - Google Patents

Abstract

本发明提供了一种三维大脑扫描图像(如MR或CT)中一个或多个脑室的分割方法。该图像被配准到大脑模型，为一个或多个脑室中每一个的脑室模型。基于脑室模型确定各自的相关区域。首先对相关区域实施区域增长过程得到对象区域，然后基于解剖学知识对其进行修剪。根据分割结构，在3D空间中建立一个或多个对象的3D表面模型。利用用户在3D空间中选择的表示丢失细节特征的修正点，编辑和改善3D表面。将被选点附近的3D表面区域平滑的向修正点扭曲，将修正路径与3D表面的其余部分结合得到精确的解剖结构模型。

Description

图像中解剖结构分割与建模的方法及系统

技术领域

本发明涉及一种图像中解剖结构分割的方法及系统，以及一种为分割结构建立3D表面模型的方法及系统。具体应用的例子是医学图像中大脑脑室系统的分割与建模，例如MR图像和CT图像。

背景技术

如图1所示，人类大脑脑室系统由4个互相联络的室组成，即左侧脑室、右侧脑室、第三脑室和第四脑室。脑室内含有脑脊液(CSF)，周围还有白质(WM)和灰质(GM)。位于大脑中的两个侧脑室较大，呈C型，大约围绕在基底节的背向。每一个侧脑室通过前角、下角和后角分别伸向额叶、颞叶和枕叶。侧脑室均经室间孔与第三脑室(建立在中脑的中央)相通，第三脑室经脑导水管(位于中脑内)与第四脑室(建立在后脑内)相通。图中的缩写定义为：AC：前连合(Anterior Commissure)；BC：基底池(Basal Cistern)(脚间池)；CC：胼胝体(Corpus Callosum)；CP：大脑脚(Cerebral Peduncle)；CQ：四叠体(Corpora Quadrigemina)；HP：脑下垂体(Hypophysis)(脑垂腺)；ICV：大脑内静脉(Internal Cerebral Vein)(在横裂内)；IS：下丘脑漏斗(Infundibular Stalk)；LT：综板(Lamina Terminalis)；LV：侧脑室(Lateral Ventricles)；MI：中间块(Massa Intermedia)(中连合)；MO：延髓(Medulla Oblongata)；OC：视交叉(Optic Chiasma)；PC：后连合(Posterior Commissure)；PG：松果体(Pineal Gland)；SP：透明隔(Septum Lucidum)；TC：灰结节(Tuber Cinereum)；TF：横裂(Transverse Fissure)(位于胼胝体下的蛛网膜下隙)；V3：第三脑室；V4：第四脑室。

MR成像使得非侵入性地获得人脑的3D图像成为可能。由于CSF容量和脑室形状的变化通常与几种内在和外在的病理有关，因此在MR图像中对脑室系统的分割和量化尤为重要。

由于脑室手动分割方法消耗时间，具有主观性和不可重复生产性(或者不可重复性)，对于脑室分割提出了许多自动方法。一般而言，根据是否采用3D脑室模型，脑室分割方法可归类为基于模型方法或非基于模型方法。

非基于模型方法，如强度界定[17]和区域增长[12，13，19]适用于脑室系统形状和大小变化的情况。但是，由于这些方法不采用脑室形状的先验知识，可能会产生从脑室区域到非脑室区域的“泄露”。而且，由于图像的不均匀性或者图像中噪声和部分容积伪影(volume artifacts)的存在，这些方法可能会遗漏一些脑室区域。由于第三脑室的精确边界依赖于它们的形状和拓扑约束以及与周围物体的关系，所以当采用这些非基于模型方法，第三脑室的准确分割尤其具有挑战性。

相反，基于模型方法，如图谱扭曲[4]或几何和参数模型变形[3，6，8]，采用显式或隐式模型作为脑室先验知识的形状模型。当模型与研究对象形状之间的差别很小时，这些方法对噪声具有鲁棒性，并能够实现精确分割。但是，由于脑室形状和大小具有很大差异，要想实现适应于每一差异的模型变形，设计一个合理的能量或相似性函数非常困难。而且，这些方法中不可避免地存在局部最小化问题，将导致错误的分割。

一般来说，图像中解剖结构分割主要存在两个难点。第一，由于部分容积效应，导致被研究结构(如脑室系统)及其周围组织(如灰质)之间出现过渡区域。如果这些过渡区域被完全排除，该结构可能被欠分割或分裂成几个不连接的部分。第二，被研究结构及其周围组织之间的一些界线很弱而在图像中不易被检测。因此，一些对象区域可能“泄露”(即连接)至其他非对象区域。目前，没有现成的方法可以检测出过渡区域，并且同时防止对象区域“泄露”至非对象区域。

发明内容

本发明目的在于提供一种图像中结构分割与建立3D表面模型的方法及系统。

特别地，本发明提供了一种由大脑扫描数据组成的三维大脑扫描图像中一个或多个脑室的分割方法，该方法包括步骤：

(a)根据该图像，配准一大脑模型，该大脑模型包括一个或多个脑室中每一个脑室的一个或多个各自的脑室模型；从而形成大脑模型中的位置与大脑扫描图像中各自位置的对应关系；

(b)根据所述对应关系，基于一个或多个脑室模型在图像中确定一个或多个各自的相关区域；

(c)利用大脑扫描数据，在一个或多个相关区域上执行区域增长，用以形成表示各自脑室的各自容量；以及

(d)利用各自的容量，分割大脑扫描图像。

本发明可以进一步包括建立分割解剖结构的表面模型，以及编辑该表面模型以精确地描述分割时丢失的特征和细节的步骤。

步骤(c)可以包括生成形为连通区域的容量，在步骤(d)之前具有基于具体涉及到脑室的解剖知识修剪所述容量的步骤。

本发明可以表现为执行所述方法的计算机系统。该计算机系统可集成获取图像的设备。本发明还可以表现为计算机程序产品，例如记录有形计算机语言的介质，其包括供计算机系统操作用以执行所述方法各步骤的程序指令。

附图说明

仅参考以下附图的例子，详细说明本发明实施例，其中：

图1(a)-(c)示出了人类大脑脑室系统的一个例子；

图2示出了本发明一实施例中系统的主要流程图，包括方法步骤202和204；

图3示出了分割图像中脑室系统的方法202的流程图；

图4示出了根据方法202的分割输出结果，生成精确的脑室结构的3D表面模型的方法204的流程图；

图5示出了本发明实施例提供的采用修正点修正表面模型的过程；

图6示出了利用方法202分割数据集IBSR-18中左侧脑室得到的结果；

图7示出了利用方法202分割数据集IBSR-18中第三脑室得到的结果；

图8示出了利用方法202分割数据集IBSR-18中第四脑室得到的结果；

图9示出了利用方法202从数据集BIL-20中4个不同的脑容量集中分割出的4个脑室结构。

具体实施方式

参考图2，示出了本发明一实施例的方法200的步骤，该方法用以生成脑室的3D表面模型。

方法200的输入是一容量图像。在步骤202中，容量图像中的脑室被分割。在步骤204中，建立每一个脑室的3D表面模型，并且编辑该3D表面模型，以提升其精确性。需要注意的是，在其他实施例中，步骤204可以不在步骤202之后。此外，步骤204的方法并非仅应用于方法200中，还可以有其他可能的应用，在3D建模的领域中，可以单独执行或合并执行。

步骤202：在容量图像中分割脑室

参考图3，示出了本发明一实施例的方法202的步骤，该方法用以生成表明脑室系统的容量图像。

方法202的输入是一容量图像。在步骤302中，根据多个(如10个)可自动识别的脑室标定点，将图像重组到标准的塔莱拉什(Talairach)空间，该图像的标准脑室模型被扭曲。在步骤304中，采用变形的脑室模型为每一个脑室指定相关区域。在步骤308、310和312中，侧脑室、第三和第四脑室被分割。在步骤306a、306b和306c中执行迟滞性界定(hysteric thresholding)(即，具有滞后作用的界定)，用以生成包括部分脑室的CSF连通区域，该CSF区域还包括极小的非脑室区域。

步骤302：重组图像

给定一容量图像I，通常采用塔莱拉什变换[9]将该图像I重组到标准塔莱拉什空间[14]，这样根据解剖学知识能够处理或理解该图像。但是，当塔莱拉什标定点没能自动定位时，不能自动执行塔莱拉什变换。

因此，在实施例中，采用基于脑皮层轮廓的配准方法重组图像。大脑的皮层轮廓接近于其皮层表面的凸壳。利用形态分析[11]自动提取出图像中皮层轮廓S₁，通过插入[8]3D TT(Talairach-Tournoux)大脑图谱的2D数字电子版，生成3D塔莱拉什空间中的皮层轮廓S₂，采用3D TT大脑图谱[8]中的脑室系统作为标准脑室容量模型。

用三角网格表示轮廓S₁和S₂，分别在顶点标记标号为Q₁和Q₂。采用迭代最近点(ICP)的方法[2]配准点集Q₁至Q₂，获得一线性变换，用于将图像重组到塔莱拉什空间。采用标准放射科协议(http://www.grahamwideman.com/gw/brain/orientation/orientterms.htm)定义塔莱拉什空间的坐标系统(xyz)，其原点定位于3D TT图谱的前连合，x轴从对象的右边延伸到左边，y轴从对象的后面延伸到前面，z轴从对象的下面延伸到上面。

步骤304：指定相关区域

在实施例中，为每一个部分脑室指定相关区域，首先要在图像和3D TT图谱中识别出10个脑室标定点[7]。采用基于模型的半全局方法自动识别图像中的10个脑室标定点，作为本领域知识，采用医学图像理解环境(MIUE)工具[7，8]在3D TT图谱中交互式的指定这些标定点。

举例来说，在每个侧脑室中有4个标定点，它们是每个侧脑室的最后侧点、最优点、额叶前外侧极点以及后部中线交点。该标定点还包括第三脑室中的前部极点以及第四脑室中的后部优越点。

基于图像中的10个脑室标定点和TT大脑图谱，标准的脑室模型被配准到图像中。由于自动检测到的标定点的位置可能并不精确[7，10]，可采用薄板样条逼近方法[10]获得配准(或扭曲)函数。

扭曲(或变形)的脑室容量模型被分成四个子容量：V₁(左侧脑室)，V₂(右侧脑室)，V₃(第三脑室)和V₄(第四脑室和导水管)。根据公式(1)扩张相应的扭曲子容量V_i得到每一部分脑室对应的相关区域Ω_i。

${\mathrm{\Ω}}_i=\left\{\begin{array}{cc}\left\{p\right|s(V_i,p)\≤d_0\}-V_0& (i=\mathrm{1,2})\\ \left\{p\right|s(V_i,p)\≤d_0\}& (i=\mathrm{3,4})\end{array}\right.---\left(1\right)$

在公式(1)中，相关区域Ω₁至Ω₄分别用于分割左侧脑室、右侧脑室、第三和第四脑室。s(V_i，p)表示体素p(p＝(p_x，p_y，p_z)∈R³)至容量V_i边界的有符号的最小欧式距离，s(V_i，p)为正值表示体素p位于容量V_i之外，s(V_i，p)为负值表示体素p位于容量V_i之内。举例来说，设d₀为6mm，这样每个区域正好足够大到包括三种类型的脑组织：灰质，白质以及CSF，还包括相关的部分脑室。这使得在该区域内能随后估计出相关部分脑室的界限。另外，V₀表示中间矢状片。举例来说，根据公式(2)设V₀的厚度为8mm，从Ω₁和Ω₂中排除V₀，以防止两个侧脑室“泄露”到跨半球的CSF，或者两个侧脑室“泄露”到彼此内。

V₀＝{p|-4≤x≤4}                    (2)

步骤306a，306b和306c：执行迟滞性界定

尽管几种方法[5，15，16，21]适用于从脑容量中分割出CSF区域，通常提取的CSF区域不仅包括脑室区域，还包括一大部分非脑室区域。从大量互相连通的非脑室区域中分割出脑室区域很困难。因此，这些方法可能不能定位脑室CSF和非脑室组织之间的过渡区域，从而导致欠分割。尽管现有方法[20]适用于提取过渡区域，这些方法要么基于梯度，要么基于局部熵。因此，它们很有可能提取一大部分非脑室CSF区域作为过渡区域。

在实施例中，步骤304中指定的相关区域Ω₁(Ω₁至Ω₄)被用作引导步骤306a，306b和306c采集CSF连通区域X，该CSF连通区域包括与其对应的部分脑室。在步骤306a，306b和306c中，根据以下子步骤，采用迟滞性界定采集区域Ω_i对应的区域X：

步骤1：在每一个区域Ω_i中分别计算部分脑室的两对强度界限。

举例来说，根据以下步骤执行上述步骤1。

首先，根据体素的强度，利用模糊C均值方法[1]将区域中图像的所有体素归类为5组。这5组代表三种类型的组织(GM，WM和CSF)以及两个过渡区域CSF_GM(位于CSF和GM之间)和GM_WM(位于GM和WM之间)。

然后，将一个组k的强度g记为u_k(g)，将每一组的中心强度记为c_k，(k＝1，2，...，5)。为了不失一般性，假设c₁＜c₂＜...c₅，计算出两个隶属函数u_k和u_k+1的交点g_k，如u_k(g_k)＝u_k+1(g_k)，其中k＝1，2，...，4。k组的最低界限t_kL和最高界限t_kH被分别设为g_k-1和g_k，g₀和g₅被分别设为可能的强度最小值和最大值。

根据本领域知识，挑选出与CSF和CSF_GM所对应的两组。举例来说，在T1-MR图像中，具有强度界限[t_1L，t_1H]的第一组被选为CSF，具有强度界限[t_2L，t_2H]的第二组被选为CSF_GM。包括CSF组的界限作为较窄界限[T_L1，T_H1]中，包括CSF和CSF_GM组的界限作为较宽界限[T_L2，T_H2]。也就是说，T_L1＝t_1L，T_H1＝t_1H，T_L2＝min{t_1L，t_2L}，T_H2＝max{t_1H，t_2H}。

步骤2：根据较窄界限[T_L1，T_H1]，为每一Ω_i采集部分脑室对应的核心区域K。

举例来说，根据以下步骤执行上述步骤2。

首先，根据最低界限T_L1和最高界限T_H1对图像I进行二值化，得到CSF组{p|T_L1≤I(p)≤T_H1}。然后，根据6个相邻体素的连通性，从CSF组中提取出最大连通区域K。由于区域Ω_i是通过扩张变形的部分脑室而生成的，它与图像中对应的部分脑室大致匹配，自然的，区域K至少包括区域Ω中相关侧脑室的主要部分。也就是说，根据区域Ω₃得到的区域K包括第三脑室的主要部分，当根据每一个其他区域获得区域K时，区域K是左侧脑室、右侧脑室或第四脑室的主要部分。区域K被记为相关部分脑室的核心区域。

步骤3：根据较宽界限[T_L2，T_H2]，区域K可自适应的扩张到包括过渡区域。

在一个实施例中，采用基于区块的边界区域增长方法执行上述步骤3，可自适应的扩张区域K使之包括部分脑室的过渡区域，同时防止区域K“泄露”到非脑室区域。

如果容量内一边界体素p的26个最接近的相邻体素中至少一个q满足：q∈Ω-K且T_L2≤I(q)≤T_H2，那么该边界体素p被认为是激活体素。根据26个相邻体素的连通性，K中的边界激活体素被分组形成边界区块集

其中n表示区块的数量。一个区块中的所有体素与26个相邻体素连接，两个不同的区块

和(i≠j)不相连接。

分别对每一个区块

实施区域增长。最初，被设为

根据公式(3)由

重复地生成

在公式(3)中，N₂₆(p)表示体素p的26个相邻体素。

${\∂}_{i,k+1}=\underset{p\∈{\∂}_{i,k}}{\∪}\left\{q\right|q\∈N_{26}\left(p\right),T_{L2}\≤I\left(q\right)\≤T_{H2},q\∈\mathrm{\Ω}\}-\left(K{\∪}_{j=0}^k{\∂}_{k,j}\right)$

(k＝0，1，2，...)            (3)

继续执行由

生成

的过程，直至k＝k_i时为空，或者直至

中体素数量高于

中体素数量的两倍，即由于过渡区域的大小被期望是与

的规模相同，该终止条件

是为了防止过渡区域“泄露”到非脑室区域。

在该过程的最后，由区块

获得新的容量

最后，根据公式(4)，该新的扩张容量V_i和核心区域K组合在一起，生成CSF连通区域X。

$X={\∪}_{i=0}^n\left({\∪}_{n=0}^{k_i}{\∂}_{i,n}\right)\∪K---\left(4\right)$

在分别位于步骤308，310和312的步骤306a，306b和306c中，根据关于以下脑室系统的形状、强度及解剖学的本领域知识，进一步修剪CSF连通区域。

步骤308：侧脑室分割

为了分割两个侧脑室，分别对区域Ω₁和Ω₂实施迟滞性界定，以获得两个容量X₁和X₂，它们是左侧脑室和右侧脑室的主要部分。为了检测到中间矢状片V₀中侧脑室的可能的其余部分，首先分别定位在区域V₀中的X₁和X₂的边界区块

和

利用基于区块的边界区域增长方法自适应的扩张区域V₀中和

因而得到两个新的容量X ₁和X ₂，它们包括左侧脑室和右侧脑室的其余部分。

步骤308a：侧脑室分离

当在两个侧脑室之间的透明隔足够大时(例如，在矢状方向上至少为一个体素的厚度。在测试数据集中，这大约在30％的测试对象中发生)，X ₁和X ₂分离(即，容量X ₁和X ₂的交集X ₁₂＝X ₁∩X ₂为空)，X₁∪X ₁和X₂∪X ₂分别被认为是左侧脑室和右侧脑室。当透明隔非常薄时，通过非空的重叠区域X ₁₂，两个区域X ₁和X ₂可能连接在一起，那么有必要根据以下步骤分离左侧和右侧脑室：首先，从X ₁和X ₂中去除X ₁₂，以获得两个区域X′₁＝X ₁-X ₁₂和X′₂＝X ₂-X ₁₂。然后，对于X ₁₂中的每一个体素p，如果它到区域X′₁的边界的距离小于到区域X′₂的边界的距离，即s(X′₁，p)＜s(X′₂，p)，体素p被分配到X′₁，否则如果s(X′₁，p)＞s(X′₂，p)，体素p被分配到X′₂。如果s(X′₁，p)＝s(X′₂，p)，p被认为是透明隔中的体素。最后，并集X₁∪X′₁和X₂∪X′₂被认为是左侧和右侧脑室的分割，即X₁和X₂被分别更新为X₁∪X′₁和X₂∪X′₂。

步骤310：第三脑室分割

首先，对区域Ω₃实施迟滞性界定，以获得CSF连通容量X₃。

然后，将识别出的属于左侧或右侧脑室部分的体素从X₃中去除。也就是说，X₃被更新为X₃-(X₁+X₂)。最后，从X₃中去除其他外脑室体素。

步骤310a：基于投影的非脑室区域修剪

举例来说，采用一种基于投影的修剪方法从X₃中去除非脑室体素。由于第三脑室是位于大脑中间的狭窄的裂隙，沿着矢状方向(从左至右)，容量X₃中包括的非脑室部分远宽于第三脑室。以下为基于投影的修剪方法的步骤。

步骤1：根据公式(6)，通过将容量X₃投影到的中间矢状平面x＝0上，生成二维图像f(y，z)。

f(y，z)＝#{p|p_y＝y，p_z＝z，p∈X₃}                        (6)

在公式(6)中，#表示一个集基数，即平面x＝0中的点(y，z)，f(y，z)表示在点(y，z)的投影线上容量X₃中体素的数量。

步骤2：采用模糊C均值方法将所有非零值{f(y，z)≠0}归类为两组，获得一自适应阈值h，使得一组中f(y，z)小于h，另一组中大于h。

步骤3：对于每一体素p∈X₃，如果f(p_y，p_z)＞h，从X₃中去除该体素。

步骤310b：标定点引导非脑室区域修剪

在采用基于投影的外脑室修剪方法之后，在第三脑室的前下部分X₃可能仍包含小部分狭窄的非脑室区域。举例来说，采用标定点引导非脑室区域修剪方法去除这个非脑室区域。在标定点引导非脑室区域修剪方法中，在第三脑室前极之前的所有体素被去除。采用基于模型的方法[7]在图像中识别出标定点(第三脑室的前极)。

步骤310c：基于形状的非脑室区域修剪

在第三脑室的上部，X₃包括由横裂和ICV组成的狭窄的C型区域。而且，从PC(或PG)到下后部，X₃包括一条或多条小而狭窄的路径“泄露”到基底池。举例来说，基于该“泄露”的条状特征，采用基于形状的非脑室区域修剪方法去除这些“泄露”。首先，通过从后到前分组冠状切片上的连通区域，定位需要去除的所有候选分量。然后，从这些候选分量中，识别和去除条状的“泄漏”。举例来说，在基于形状的非脑室区域修剪方法中执行以下子步骤。

步骤1：在初始步骤中，将候选泄漏分量设为

将临时分量设为

在步骤1中，定位容量X₃中最后部的冠状切片y₀＝min{y|p(x，y，z)∈X₃}，

被置为空，而

被置为{{C₀}|C₀∈S₀}，其中C₀表示索引为y₀的冠状切片上X₃中所有8个相邻体素连通区域S₀中的一个，{C₀}为由区域C₀组成的候选泄露分量。

步骤2：通过跟踪中的每一分量生成

定位所有的候选泄漏分量。

对于每个分量

如果在索引y₀+k+1的冠状切片上具有8个相邻体素连通区域C_k+1，并且C_k+1连接至C_k，这意味着，至少存在一个体素p_k+1∈C_k+1为另一体素p_k∈C_k的26个相邻体素之一，区域C_k+1加上分量L_k组成了一个新的分量L_k+1＝{C₀，C₁，...，C_k，C_k+1}。

如果C_k+1中体素与C_k中体素的面积比率大于给定阈值r(例如，r设为3)，

要加上L_k+1。否则，为了进一步增长，

要加上L_k+1。如果

非空，重复步骤2由

生成

步骤3：在步骤3中，去除X₃上部的C型泄露分量。

如果满足以下三个条件，候选分量

被识别为由横裂和ICV组成的C型泄露分量。

(1)C_k+1为一分支区域，即，存在另一个连通区域C′_k∈S(k)且C′_k≠C_k。

(2)若P_k，P_k+1，P′_k分别是区域C_k，C_k+1 and C′_k的质心，∠P_kP_k+1P′_k的度数小于30°，以及

(3)每个区域C_i∈L_k+1(i＝0，1，...k+1)在索引为y₀+i的冠状切片中X₃的所有体素的上部，即max{z|p(x，y₀+i，z)∈C_i}＞max{z|p(x，y₀+i，z)∈X₃-C_i}。

如果L_k+1＝{C₀，C₁，...，C_k，C_k+1}被识别为C型泄漏分量，从X₃中去除C₀，C₁，...，C_k，从

中去除L_k+1。

步骤4：在步骤4中，去除X₃后部的条状泄露分量。

对于每个候选分量如果它位于区域Ω₃的质心G(x，y，z)的后部，即y₀+k+1＜G_y，那么它被识别为泄露分量。从X₃中去除C₀，C₁，...，C_k，从中去除L_k+1。

最终的X₃区域就是第三脑室的分割结果。

步骤312：第四脑室和导水管(adequate)分割

由于第四脑室与导水管之间没有明确的界限，它们被同时分割。对区域Ω₄实施迟滞性界定，得到容量X₄。在导水管和第四脑室的连合处，由于导水管的后壁(即四叠体)变得非常薄，可能从图像中识别不出来，X₄会从第四脑室“泄露”到小脑周围的基底池。同时，由于导水管是连接第三脑室和第四脑室的唯一条狭窄路径，部分导水管或全部的导水管可能不包括在X₄中。

步骤312a：对第四脑室进行基于形状的修剪

为了去除“泄露”，计算容量X₄的每个轴切片中索引为坐标z的体素f(z)的个数。当f(z)达到它的最大值时，切片z_max被定位。对于f(z)＞0的切片，根据公式(7)计算在上部(或背部)方向上从切片z_max到后续切片的相对增长率。

q(z)＝[f(z+1)-f(z)]/f(z)                    (7)

由于在“泄露”开始之处，f(z)快速增长，因此从z_max到腹侧方向(记为切片z_leak轴向)上，在q(z)达到正的最大值的地方，定位第一个泄露切片。如果q(z)的最大值不是正的，表示X₄没有“泄露”到基底池。在这种情况下，z_leak被置为V₄中体素的z坐标的最大值。

由于导水管向前倾斜连接至第三脑室，y_leak表示在泄露切片z_leak上容量X₄的最后部，在背向上从z_leak向前的y坐标小于y_leak的所有体素将作为“泄露”并从容量X₄中去除。

从切片z_max向下(或背部的)的方向，要求f(z)不再增长。因此如果存在一个切片z_min，如f(z_min+1)＞f(z_min)，z坐标小于z_min的所有体素将从X₄中修剪掉。

最终的区域X₄就是第四脑室的分割结果。

为了找到导水管，将从区域Ω₄中得到的较宽界限记为[T_L2，T_H2]，将体素p(x，y，z)的定向邻居记为N_z+(p)＝{(p_x+i，p_y+j，p_z+k)|i，j＝-1，0，1，k＝0，1}，将切片z_leak中X₄的所有体素记为S₀，根据公式(8)，通过定向区域增长由S_n生成S_n+1。

$S_{n+1}=\underset{p\∈S_n}{\∪}\left\{q\right|q\∈Z_{z+}\left(p\right),T_{L2}\≤I\left(q\right)\≤T_{H2}\},(n=\mathrm{0,1,2},...)---\left(8\right)$

重复地由S_n生成S_n+1，直至S_n+1为空或者S_n+1中体素的个数大于S_n中体素的个数(即，#(S_n+1)＞#(S₀))。将S₁∪S₂...∪S_n作为导水管容量。如果当#(S_n+1)＞#(S₀)时，由S_n生成S_n+1的重复过程终止，这是因为检测到的导水管到达第三脑室。另一方面，如果当S_n+1为空时，由S_n生成S_n+1的重复过程终止，这是因为检测到的导水管由于部分容积效应未能到达第三脑室。在大多数情况下，这个过程是由于#(S_n+1)＞#(S₀)而停止。

步骤204：建立并改善脑室表面模型

图4示出了方法204的流程图，这个方法用于由方法202的分割输出结果生成精确地脑室结构3D表面模型。在步骤402中，为脑室结构建立了一个表面模型，在步骤404中，通过局部正弦扭曲改善了表面模型的细节。

步骤202产生了脑室的分割结果。然而，如果切片之间的距离太大和/或研究图像的质量较差，一些细节可能仍会遗失或是不精确。在这种情况下，为了描述图像切片之间丢失了细微特征或者图像质量影响了细微特征的解剖结构，几何表面模型更为灵活和平滑。为了建立精确的脑室表面模型，可采用著名的移动立方体方法[22]，根据步骤202的脑室容量输出结果生成初始表面模型。表现为三角网格的初始表面模型被简化，用以减少计算时间并提高后续处理的效率。

为了加强精确度，实施例提供的系统支持用户采用局部正弦扭曲方法修改表面模型。基于本领域知识，用户能够通过在3D模型空间中设置修正点来指明丢失的细微特征。局部正弦变形(LSD)函数扭曲限定区域平滑地接近修正点，用以恢复丢失的细微特征，且不失解剖结构的连续性(如图5所示)。

图5示出了本发明一实施例提供的采用修正点修正表面模型的过程。

假设用户在模型M附近设置了一个修正点A，表示该模型中丢失了细节特征。该模型表现为多角网格，对于网格上每个顶点V，从A到V的距离记为d(A，V)。A与该模型之间的距离为d(A，M)＝min(d(A，V)|V∈M)。给定半径R＞d(A，M)，(系统中R可调)，建立有限的顶点集P＝{p₁，p₂，...p_k |d(A，p_i)＜R}。对于集合P中的每个点p_i，计算对应的点q_i，如下：q_i位于A到p_i{i＝1，2，...k}的线上，根据公式(9)，通过LSD函数计算A到q_i的距离：

$d(A,q_i)=\mathrm{Sin}\left(\frac{\mathrm{\πd}(A,p_i)}{2R}\right)$ (i＝1，2，...k)              (9)

将每个p_i替换为按照上式计算的q_i{i＝1，2，...k}，表面模型的局部区域被扭曲到修正点A，因此恢复了丢失的解剖细微特征。

步骤204中表面模型加强过程是一交互过程，并且能够迭代执行得到满意的输出结果。

步骤202的实验结果

图6示出了采用方法202以及本发明一个实施例提供的侧脑室分割方法，分割数据集IBSR-18(IBSR-18-02图，切片144a，57c，120s)中的一左侧脑室得到的结果。轮廓线1202-1218表示为左侧脑室自动确定的相关区域。轮廓线1220-1230表示扩张到中间矢状片的左侧脑室模型的相关区域。轮廓线1232-1242表示通过较窄界限获得的区域，轮廓线1231-1242周围的轮廓线表示通过较宽界限获得的其他区域。

图7示出了采用本发明一个实施例提供的方法202，分割数据集IBSR-18(IBSR-18-02图，切片142a，60c，128s)中的一第三脑室得到的结果。从左至右的这四列示出了轴向、冠状、矢状和3D视图。第一行示出了为第三脑室自动确定的相关区域。第二行示出了在ROI中迟滞性界定的结果。第三行示出了通过基于投影的修剪获得的结果。第四行示出了在对前部进行标定点引导修剪和对其他外脑室部分进行基于形状的修剪之后的最终结果。

图8示出了采用本发明一个实施例提供的方法202，分割数据集BIL-20(BIL-Ja03，切片44a，102c，129s)中的一第四脑室得到的结果。从左至右的这四列示出了轴向、冠状、矢状和3D视图。第一行示出了为第四脑室自动确定的相关区域，第二行示出了在ROI中执行迟滞性界定获得的结果，第三行示出了在“泄露”去除之后的最终结果。

图9示出了采用本发明一个实施例提供的方法202，从数据集BIL-20中的4个不同的大脑容量集中分割的4个脑室结构。第一至第四行示出了一异常的成人大脑(具有脑肿瘤)、一正常的成人大脑、一孩童大脑以及一老年人大脑的容量图像。从左至右的这四列分别示出了轴向、冠状、矢状原始扫描视图以及提取的脑室系统的3D视图。

实施例的优点

在步骤202中，采用容量变形模型作为本领域知识，自动的确定被研究结构分割的相关区域，如实施例的例子中的脑室结构。对于精确的分割，适当的ROI是关键。如果ROI太小，它可能不包括被研究结构。另一方面，如果ROI太大，它可能包括太多不相关信息，导致错误的分割。在步骤202中，首先通过3D点基于标定点的扭曲方法，将模型进行变形，使其大致匹配图像中与之对应的结构，然后通过扩张(或扩大)变形模型，确定ROI。ROI采用被研究结构的先验形状，所以ROI中不相关信息量被最小化。因此，步骤202可抗噪声和抗形状和大小的较大差异。

此外，在步骤202给定的相关区域中，区域增长过程采用迟滞性界定方法。在迟滞性界定方法中，采用两对强度界限，即较窄界限和较宽界限。较宽界限的范围包括较窄界限的范围。该对较窄界限用于采集不包括过渡区域的核心部分，而该对较宽界限用于采集该结构的过渡区域。当检测到“泄露”，区域增长过程停止。当最小化“泄露”时，步骤202的区域增长过程能够检测到过渡区域。由于对于正确分割来说，过渡区域检测能力很关键，而且“泄露”最小化大大减轻了区域修剪过程的负担，所以这是一个优势。

另外，步骤202的区域修剪过程采用了多种知识为基础的策略，例如基于投影、标定点引导、及基于形状的修剪，这对第三脑室的正确分割也很关键。

又，步骤202优于现有技术的方法，如[19]。[19]提供的方法依赖于精确识别出AC，PC和MSP，因此如果AC，PC和MSP提供的位置的精确度不高(AC和PC位置的误差需要小于3mm)，该方法可能会失效。另外，在[19]提供的方法中，只有一对界限用于相关区域，因此该方法不能处理部分容积问题，这会导致组成部分的一些部分不连通。特别地，方法[19]所采用的ROI的形状为三角形，不同于脑室的形状。因此，ROI包括大量非脑室组织，导致[19]中存在潜在的分割误差和“泄露”。相反地，在步骤202中，采用10个脑室标定点对脑室模型进行扭曲，以匹配图像中与之对应的脑室结构。由于采用薄板样条逼近方法[10]获得扭曲函数，变形模型进一步扩张到厚度为6mm，步骤202更能接受较大的标定点识别误差(在IBSR-18([7]中表2所示)中可大到3.4mm)。尽管第三脑室的前极识别错误会影响步骤202中第三脑室分割的精确性，但由于第三脑室的前极只用于修剪第三脑室的后部，该后部相对地占整个第三脑室的一小部分，所以上述影响很小且是局部的。而且，步骤202中迟滞性界定采用两对较宽的和较窄的界限，用以自适应的开发相关区域，该迟滞性界定的应用保证ROI包括过渡区域，并且同时最小化非脑室区域。又由于步骤202所采用的ROI来源于大脑图谱中的脑室形状，ROI的形状与目标结构的形状非常接近，因此大大减少了潜在的分割误差和“泄露”。

参考文献：

[1]J.C.Bezdek，基于模糊目标函数算法的模式识别，Plenum出版社.1981，纽约。

[2]P.Besl and N.McKay，“一种三维形状配准方法”，IEEE模式分析与机器智能汇刊，14(2)，pp.239-256，1992。

[3]C.Baillard，P.Hellier，C.Barillot，2000，“利用水平集与密集配准的3D大脑结构分割”，IEEE生物医学图像分析中的数学方法专题会议录，94-101。

[4]M.Holden，J.A.Schnable，D.L.G.Hill，2001，“利用非刚性配准的大脑脑室容量微小变化定量方法”，医学影像计算与计算机辅助介入大会，49-56。

[5]K.Held，E.R.Kops，B.J.Krause，(1997)，脑磁共振图像的马尔可夫随机场分割，IEEE医学影像汇刊16，876-886。

[6]M.R.Kaus，S.K.Warfield，A.Nabavi，P.M.Black，F.A.Jolesz，R.Kikinis，2001，脑肿瘤MR图像的自动分割”，放射学杂志218，586-591。

[7]J.Liu，W.Gao，S.Huang，and W.L.Nowinski，“神经图像中3D点自动标定定位的基于模型的半球分割方法”，TMI-2008(in press)。

[8]J.Liu，W.L.Nowinski，(2006)，“形状插值人类大脑立体定向图谱的混合方法”，《神经信息学》，vol.4，no.2，pp.177-198。

[9]W.L.Nowinski，G.Qian，K.N.B.Prakash，Q.Hu，A.Aziz 2006，“磁共振神经图像的快速塔莱拉什变换”，计算机辅助X线断层摄影术杂志，7-8月；30(4)：629-41。

[10]K.Rohr，基于界标图像分析，克吕韦尔学术出版集团，2001.

[11]S.M.Smith，“快速鲁棒自动大脑提取”，人脑图像描记，17(3)：143-155，2002年11月。

[12]H.G.Schnack，P.H.E.Hulshoff，W.F.C.Baare，M.A.Viergever，R.S.Kahn，2001，“人类大脑的MR图像中脑室系统的自动分割”，神经影像第14期，95-104。

[13]M.Sonka，S.K.Tadikonda，S.M.Collins，1996，“MR大脑图像的阐述”，IEEE医学影像汇刊.15，443-452。

[14]J.Talairach，P.Tournoux(1988)，人类大脑的共面立体定向图谱，斯图加特，蒂梅出版社。

[15]H.Tang，E.X.，Wu，Q.，Ma，Y.，Gallagher，D.，Perera，G.M.，and Zhuang，T.(2000).“采用多分辨率边缘检测和区域选择的MRI大脑图像分割”，计算机化医学影像与图学24，349-357。

[16]W.Wells，W.Grimson，R.Kikinis，and F.Jolesz(1996)“MRI数据的自适应分割”.IEEE医学影像汇刊15，429-442。

[17]A.J.Worth，N.Makris，M.R.Patti，J.M.Goodman，E.A.Hoge，V.S.Caviness，D.N.Kennedy，1998，“利用自动驱动直方图的MR大脑图像中侧脑室与尾状核的精确分割”，IEEE医学影像汇刊17，303-310。

[18]Y.Wang，L.H.Staib，1998，利用统计学形状模型的边界查找，“关于计算机视觉与模式识别的IEEE会议公告”，IEEE出版社，华盛顿，DC，pp.338-345。

[19]Y.Xia，Q.Hu，A.Aziz and W.L.Nowinski，“快速自动提取MR神经图像中人脑脑室系统的知识导向算法”，神经影像，21，2004，pp.269-282。

[20]C.Yan，N.Sangt，and T.Zhang，2003，“局部基于熵的过渡区域提取与界定”，模式识别快报，24(2003)，pp.2935-2941。

[21]Y.Zhang，M.Brady，and S.Smith(2001)，“通过隐马尔可夫随机场模型与最大期望值算法的大脑MR图像分割”，IEEE医学图像汇刊20，45-57。

[22]Lorensen WE and Cline HE，“移动立方体：高分辨率3D表面构建方法”，计算机图形学，1987，21(4)，163-169。

[23]P.Lindstrom and G.Turk，“快速且存储高效的多边形简化方法”，Proc.Visualization’98，IEEE计算机协会出版社，1998年10月，pp.279-286。

Claims (19)

1.一种由大脑扫描数据组成的三维大脑扫描图像中一个或多个脑室的分割方法，所述方法包括步骤：

(a)根据所述图像配准一大脑模型，所述大脑模型包括一个或多个脑室中每一个脑室的一个或多个各自的脑室模型；从而形成大脑模型中的位置与大脑扫描图像中各自位置的对应关系；

(b)根据所述对应关系，基于一个或多个脑室模型在图像中确定一个或多个各自的相关区域；

(c)利用大脑扫描数据，在一个或多个相关区域上执行区域增长，用以形成表示各自脑室的各自容量；以及

(d)利用各自的容量，分割所述大脑扫描图像。

2.根据权利要求1所述的方法，其特征在于，步骤(a)包括子步骤：

(i)将图像重组到大脑模型的坐标系统；

(ii)在重组图像中识别与大脑模型中标定点对应的标定点；以及

(iii)基于识别出的标定点，配准图像的大脑模型。

3.根据权利要求2所述的方法，其特征在于，步骤(i)包括子步骤：

(iv)确定图像中的第一皮层轮廓；

(v)确定大脑模型中的第二皮层轮廓；

(vi)配准第一皮层轮廓中的点与第二皮层轮廓中的点，以获得线性变换；以及

(vii)利用该线性变换，将图像重组到大脑模型的坐标系统。

4.根据权利要求3所述的方法，其特征在于，在步骤(vi)中采用迭代最近点方法配准所述点。

5.根据权利要求2-4任一项所述的方法，其特征在于，步骤(iii)包括子步骤：

(viii)利用薄板样条逼近方法获得配准函数；

(ix)利用所述配准函数，基于识别出的标定点配准图像的大脑模型。

6.根据前述任一项权利要求所述的方法，其特征在于，对于一个或多个相关区域中的每一个，步骤(c)包括子步骤：

(x)分别计算一对较窄强度界限和一对较宽强度界限；

(xi)根据该对较窄强度界限，确定相关区域的核心区域；

(xii)根据该对较宽强度界限，扩张核心区域使之包括相关区域周围的过渡区域，以形成作为连通区域的所述容量。

7.根据权利要求6所述的方法，其特征在于，步骤(x)包括子步骤：

(xiii)根据强度，对相关区域内的体素分组；

(xiv)基于每一组与相邻组的交点，为每一组计算一对强度界限；

(xv)确定包括相关区域强度的组的强度界限对为该对较窄强度界限；以及

(xvi)确定包括相关区域强度的组与包括过渡区域强度的组的强度界限对组合的上限和下限，作为该对较宽强度界限。

8.根据权利要求6或7所述的方法，其特征在于，步骤(ii)包括子步骤：

(xvii)根据该对较窄强度界限对图像进行二值化，得到一组；

(xviii)根据6-相邻体素连通性方法，从该组中提取出最大连通区域作为该核心区域。

9.根据权利要求6-8任一项所述的方法，其特征在于，步骤(iii)包括子步骤：

(xix)确定核心区域的边界激活体素；

将边界激活体素分组到边界路径将边界激活体将边界激活体素分组到边界路径；

(xxi)对每个边界路径实施区域增长，以获得扩张的边界路径；

(xxii)扩张核心区域使之包括扩张的边界路径，以形成连通区域。

10.根据权利要求6-8任一项所述的方法，其特征在于，在步骤(d)之前包括一修剪步骤，具体为通过对相应的脑室的处理步骤，从连通区域中去除体素。

11.根据权利要求10所述的方法，其特征在于，所述相关区域中的一个与侧脑室对应，该相关区域的修剪步骤包括以下步骤：

(xxiii)根据一中间矢状片，定位连通区域的共同边界体素集；

(xxiv)利用共同边界体素集的26-相邻体素最大连通区域，生成第一新区域作为起始点；

(xxv)重复地由前一个新区域生成后一个新区域，直至后一个新区域为空；以及

(xxvi)确定连通区域为共同边界体素集与所述新区域之和。

12.根据权利要求10或11所述的方法，其特征在于，所述相关区域中的一个与第三脑室对应，所述修剪步骤包括子步骤：

(xxvii)将连通区域投影到中间矢状面上，获得一投影图像，其中该投影图像中每一个像素表示连通区域上沿着投影线到该像素的体素的个数；

(xxviii)获得投影图像中将像素分成两组的分组阈值；

(xxix)从连通区域中去除值高于阈值的像素对应的体素。

13.根据权利要求10-12任一项所述的方法，其特征在于，所述修剪步骤包括子步骤：

(xxx)识别图像中的标定点；以及

(xxxi)从连通区域中去除定位位置相对于标定点的体素。

14.根据权利要求10-13任一项所述的方法，其特征在于，所述修剪步骤包括子步骤：

(xxxii)在图像的每一个切片上重复地定位8-相邻体素连通区域，上述像素属于所述连通区域，直至下一个切片中体素与当前切片中体素的面积比率大于预设阈值；以及

(xxxiii)确定所述8-相邻体素连通区域的组合为C型泄漏分量；以及

(xxxiv)从连通区域中去除该C型泄漏分量。

15.根据权利要求10-13任一项所述的方法，其特征在于，所述修剪步骤包括子步骤：

(xxxv)识别在连通区域具有最大个数像素的图像中的第一切片；

(xxxvi)从第一切片开始，计算下一个切片在连通区域的像素个数的增量；

(xxxvii)识别在连通区域具有最大像素个数增量的切片为泄漏切片；以及

从从连通区域中去除从位于泄露切片之外的连通区域中去除体素。

16.根据权利要求10-13任一项所述的方法，其特征在于，所述修剪步骤包括子步骤：

(xxxix)识别第一切片，所述第一切片在连通区域的像素个数高于前一个切片；以及

从从连通区域中从位于第一切片之外的连通区域中去除体素。

17.一种图像中结构的3D表面模型的建立方法，所述方法包括步骤：

根据权利要求1的方法分割所述结构；

在3D空间中建立分割结构的3D表面模型，其中所述3D表面模型包括多个脑室；以及

重复地编辑所述3D表面模型：

(xxxxi)在3D表面模型遗漏特征之处，在3D表面模型上设置修正点；

(xxxxii)计算修正点到3D表面模型中每一个顶点的距离；

(xxxxiii)对于3D表面模型中的每一个顶点，如果修正点到3D表面模型中该顶点的距离小于预设阈值，计算该顶点的对应点，该对应点位于修正点到顶点的线上，修正点到对应点的距离为sin(πd(A，p_i)/2R)，其中d(A，p_i)为修正点到顶点p_i的距离，R为预设阈值。

18.一种包括用于执行权利要求1-17任一项所述方法的处理器的计算机系统。

19.一种计算机程序产品，为计算机可读，包括供计算机系统中的处理器操作的指令，使得处理器执行根据权利要求1-17任一项所述的方法。

Images