CN105719278A - 一种基于统计形变模型的器官辅助定位分割方法 - Google Patents

Abstract

本发明公布一种基于统计形变模型的器官辅助定位分割方法，其包括如下步骤：(1)获取生物体CT数据；(2)划分高、低对比度器官，在CT数据中选取训练样本提取相应的统计先验信息；(3)建立低对比度器官的统计形变模型；(4)基于高、低对比度器官之间的相关性，借助高对比度器官辅助定位低对比度器官的初始位置；(5)在辅助定位的基础上，沿标记点法线方向进行器官的搜索分割。本发明采用的器官分割方法可快速自动地找寻器官的初始化位置，利用统计形变模型融合组织器官的位置、形状等先验信息，快速系统地完成器官的分割，大大提高了图像分割的效率，是一种有效的器官分割方法。

Description

一种基于统计形变模型的器官辅助定位分割方法

技术领域

本发明属于医学影像领域，涉及一种基于统计形变模型的器官辅助定位分割方法。

背景技术

医学图像分割在现代医学研究、临床诊断、病理分析和治疗方面有着极为重要的作用，它是生物医学图像分析的先决条件，也是人类对于组织器官结构、功能认识的必要步骤。但由于人体器官的多样性和复杂性，分割问题一直作为研究人员努力探究的热点问题。

图像分割中常用的方法包括：边缘检测法，阈值分割法，区域生长法、活动轮廓模型法等。但这些方法都有它们各自的局限性：边缘检测法是根据引起图像灰度变化的物理过程来描述图像中灰度变化的过程，然而边缘定位与噪声滤除是两个相互矛盾的部分，很难同时得到很好的解决；阈值分割法适用于目标和背景占据不同灰度级范围的图像，因此会受到分割区域灰度对比度的限制；区域生长法可能由于噪声和灰度不均一导致空洞和过分割；活动轮廓模型是一种自顶向下定位图像特征的机制，需要确定初始轮廓位置，因此分割结果很大程度上受到初始轮廓的影响。在医学图像处理过程中成像设备分辨率的局限和生物组织蠕动产生的伪影和噪声，都会导致一些器官组织边缘迷糊、分辨率低，这就给器官分割造成了相当大的困难，并且单一的分割方法也都难以取得令人满意的效果。

对此，近年来，基于先验知识的统计模型分割方法已被广泛应用于各种与图像相关的领域中。上世纪九十年代Cootes等人提出了ActiveShapeModel，后又针对全局统计模型进行了改进，提出多变量的混合模型用以替代初始的单一模型，避免了单个静态模型的局限性，具有更高科学性，参见：CootesTF,HillA,TaylorCJ,etal.Theuseofactiveshapemodelsforlocatingstructuresinmedicalimages，InformationProcessinginMedicalImagingSpringerBerlinHeidelberg，1993:33-47.此外，研究者们从模型和算法角度均对统计形变模型在医学图像分割上的应用进行了深入探索研究。

但就发展现况来看，仍存在一些缺陷：首先是先验模型的初始化问题，统计形变模型的分割方法中，器官的初始化定位不易，若初始定位偏离真实形状越远越不易分割，使得后续分割难度加大；其次是分割的自动化程度仍有待提高，需在保障分割精度的同时，更加高效的完成器官分割。本发明主要针对以上两个问题，利用统计形变模型更高效、准确地分割组织器官。

发明内容

为了解决器官图像分割工作量大、复杂度高、效率低的问题，本发明提出了一种基于统计形变模型的器官辅助定位分割方法。为了达到能够系统分割器官组织的目的，本发明采用了基于统计形变模型的辅助定位分割方法，借由高对比度器官图像来辅助分割低对比度器官图像，从而系统地完成器官图像的分割。

为了实现上述目的，本发明采用的技术方案是：

一种基于统计形变模型的器官辅助定位分割方法，包括如下步骤：

(1)获取生物体CT数据；

(2)划分高、低对比度器官图像，在CT数据中选取训练样本提取相应的统计先验信息；

(3)建立低对比度器官图像的统计形变模型；

(4)基于高、低对比度器官之间的相关性，借助高对比度器官辅助定位低对比度器官的初始位置；

(5)在步骤(4)辅助定位的基础上，沿标记点法线方向进行器官的搜索分割。

2.根据权利要求1所述的基于统计形变模型的器官辅助定位分割方法，其特征在于：所述的步骤(1)，使用Micro-CT成像系统采集生物体多个角度的投影数据，利用滤波反投影方法对投影数据进行重建，获取生物体三维CT体数据；

所述的步骤(2)，根据X射线对生物体不同组织器官图像穿透能力的差异，将器官图像分为高、低对比度器官图像；选取训练样本数据，对其中高、低对比度器官图像分别进行分割，通过有限元剖分获取相应的三维点云数据；分别将高、低对比度器官的三维点云数据进行相似性变换对齐配准；配准后求取两者的均值模型；

所述的步骤(3)，计算步骤(2)中所得低对比度器官图像配准结果的协方差矩阵，并进行特征分解，获得对应的特征值和特征向量，对特征值的绝对值降序排列，选取前t个主成分，将主成分特征向量单位化处理，此时，低对比度器官的统计形变模型可被描述为其均值模型与特征向量和形变参数的线性组合；

所述的步骤(4)，选取待测样本，对其中高对比度器官阈值分割，有限元剖分后获取三维点云数据；求解待测样本中高对比度器官与步骤(2)求得的相应均值模型之间的位置形状变换；利用高对比度器官之间的变换辅助定位待测样本中低对比度器官的初始位置；

所述的步骤(5)，在待测数据中获取初始位置后，利用低对比度器官的统计形变模型搜索查找边界点，使其逐渐往目标区域逼近，搜索的区域范围选定为每个边界点沿当前轮廓的法线方向(即与该点相邻的两点连线的垂直方向上)，直至搜寻到最佳匹配的变形模型，完成最终分割。

2.根据权利要求1所述的基于统计形变模型的医学图像分割方法，其特征在于：所述步骤(2)的具体方法包括以下步骤：

(2a)根据X射线对生物体不同组织器官穿透能力的差异，将器官分为高、低对比度器官，组织密度较大X射线不易透过的器官为高对比度器官，组织密度较小X射线穿透性较大的器官为低对比度器官；

(2b)选取N个小鼠的三维体数据作为训练样本，对其中高、低对比度器官分别进行分割，前者采用自动阈值分割，后者采用手动分割，并将分割结果进行有限元离散剖分，获取相应的三维点云数据；

(2c)对步骤(2b)所得的高、低对比度器官三维点云数据在形状和位置上进行归一化，分别选取基准模板将所有样本进行旋转、平移等相似性变换，得到两者对齐配准结果；

(2d)将步骤(2c)配准对齐后的结果求和平均，分别计算高对比度器官均值模型和低对比度器官均值模型

${\stackrel{\‾}{X}}^H=\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}{X_i}^H,{\stackrel{\‾}{X}}^L=\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}{X_i}^L;$

其中，为第i个高对比度器官对齐配准后的样本，为第i个低对比度器官对齐配准后的样本，N为训练样本数，均为K行乘3列的矩阵，K为配准后单个样本三维点云中点的个数。

3.根据权利要求1所述的基于统计形变模型的医学图像分割方法，其特征在于：所述步骤(3)的具体方法包括以下步骤：

(3a)计算低对比度器官配准结果的协方差矩阵S^L：

$S^L=\frac{1}{N-1}\underset{i=1}{\overset{N}{\Σ}}({X_i}^L-{\stackrel{\‾}{X}}^L){({X_i}^L-{\stackrel{\‾}{X}}^L)}^T$

S^L为3K行乘3K列的矩阵，K为单个样本三维点云中点的个数。其中N为训练样本数，为低对比度器官均值模型，为低对比度器官对齐配准后的样本，T代表矩阵的转置运算；

(3b)对协方差矩阵S^L进行特征分解，求解对应的特征值和特征向量，将特征值按绝对值降序排列，选取前t(1≤t≤3K)个特征值，以确保能够反应样本变化的主要模式，将主成分对应的特征向量单位化处理；

(3c)根据步骤(3b)所得，低对比度器官统计形变模型中的任意形状X^L均可描述为其均值模型与特征向量P^L和形变参数b^L的线性组合，即为

$X^L\≈{\stackrel{\‾}{X}}^L+P^L{b}^L;$

为了使形状变化限制在合理的范围内，形变参数b^L需满足：

$-3\sqrt{\mathrm{\&lambda;}}<b^L<3\sqrt{\mathrm{\&lambda;}}$

其中，λ为第一主成分特征值。

4.根据权利要求1所述的基于统计形变模型的医学图像分割方法，其特征在于：所述步骤(4)的具体方法包括以下步骤：

(4a)选取M个小鼠的三维CT体数据作为待分割的测试样本，将待测样本中的高对比度器官进行自动阈值分割，并进行有限元离散剖分，获取相应的三维点云数据；

(4b)步骤(4a)所得待测样本的高对比度器官三维点云形状为Test^H，计算Test^H与高对比度器官均值形状之间的变换，即Trans^H，公式如下：

(4c)基于高、低对比度器官之间的相关性，将步骤(4b)中计算所得的高对比度器官之间的变换Trans^H近似为低对比度器官在待测样本中相对位置形状的变换Trans^L，即Trans^L≈Trans^H，结合低对比度器官均值形状模型求解低对比度器官初始位置Test^L，计算公式为：

本发明与现有技术相比具有如下优点：

第一，本发明采用的是一种基于高、低对比度器官位置形变的相关性，利用前者辅助定位后者初始位置，建立统计形变模型进行器官图像分割的方法，可快速自动地对器官进行初始化定位，有效地统计不同器官的空间与几何等先验形状信息，并通过不同的形变参数来表示同类器官的任意形变过程。

第二，本发明解决了系统分割器官图像工作量大的问题，只需预先建立一个器官形状统计模型，然后在辅助定位的基础上进行形变即可有效地完成器官图像分割，极大地减少了逐张、逐个器官图像分割的工作量，同时也避免了疲劳因素带来的分割错误。

附图说明

图1为基于统计形变模型的器官图像辅助定位分割方法流程图。

图2为采集的CT断层数据三维重建后的图像。

图3为训练样本的高对比度器官图像阈值分割结果及三维点云显示图。

图4为训练样本的低对比度器官图像手动分割结果及三维点云显示图。

图5为对训练样本高对比度器官图像对齐配准的结果及其均值模型图。

图6为对训练样本低对比度器官图像对齐配准的结果及其均值模型图。

图7为采用本发明方法分割低对比度器官图像的分割结果图。

具体实施方式

下面结合附图对本发明作进一步的详细描述，应指出的是，所描述的实施例仅旨在便于对本发明的理解，而对其不起任何限定作用。

下面结合附图对本发明作进一步说明：

步骤1：获取小鼠CT断层数据：

将已注射造影剂的实验小鼠固定在Micro-CT成像系统的成像台上，调整X射线管、旋转台以及X射线平板探测器位置，使得三者中心在一条直线上，对小鼠进行360度照射扫描，采集投影数据，利用滤波反投影方法对投影数据三维重建，得到小鼠CT断层数据。

用于实验的小鼠CT体数据如附图2所示。用Micro-CT成像系统360度扫描拍摄获取了小鼠的CT数据，利用3DMed把扫描的CT数据进行三维体重建，得到raw格式数据。图2(a)为三维展示图，图2(b)、(c)、(d)分别为x-y平面、y-z平面、x-z平面显示的二维切片图。

步骤2：划分高、低对比度器官图像，选取训练样本提取相应的统计先验信息：

(2a)根据X射线穿透组织器官的差异性将小鼠器官图像分为高、低对比度器官图像，实验以小鼠脊柱图像作为高对比度器官图像，肾脏图像作为低对比度器官图像；

(2b)选取N个小鼠的三维体数据作为训练样本，对其高对比度器官(脊柱)图像进行自动阈值分割，对低对比度器官(肾脏)图像进行手动分割，将分割结果进行有限元离散剖分，分别获取相应的高、低对比度器官三维点云数据；

(2c)对步骤(2b)中得到的训练样本的高、低对比度器官三维点云数据分别进行包括旋转、平移在内的相似性变换，对其在形状和位置上归一化。

步骤(2c)采用的配准算法是迭代最近点法(Iterativeclosestpoint,ICP)。其主要是找到基准模型和待配准模型之间的旋转R与平移T变换，使得两模型满足度量准则下的最优匹配。设基准点集为P{P_i,i＝1,2,...,N_p}(N_p为基准点集个数，P_i为基准点集P中第i个点)，待配准点集为Q{Q_i,i＝1,2,...,N_q}(N_q为待配准点集个数，Q_i为基准点集Q中第i个点)，它们之间的变换可以表示为P＝RQ+T。为了使它们配准起来，需要找到最优的旋转R与平移T，使得同源点间距离最小，满足目标表达式：

D_min＝min∑[P-(RQ_i+T)]²

其中，D_min为配准后两个点集的最近距离，min为最小距离运算。

分别选取高、低对比度器官的三维点云模板，按以上步骤对所有训练样本进行配准，使两者均在形状和位置上达到归一化。

(2d)分别计算步骤(2c)中高、低对比度器官配准结果(N个高对比度器官对齐配准后的所有样本和N个低对比度器官对齐配准后的所有样本的均值模型：

${\stackrel{\&OverBar;}{X}}^H=\frac{1}{N}\underset{i=1}{\overset{N}{\&Sigma;}}{X_i}^H,{\stackrel{\&OverBar;}{X}}^L=\frac{1}{N}\underset{i=1}{\overset{N}{\&Sigma;}}{X_i}^L$

其中，均为K行乘3列的矩阵，K为配准后单个样本三维点云中点的个数。

对训练样本的高对比度器官图像进行阈值分割的结果如附图3。在具体实施例中高对比度器官图像选取的是脊柱图像，图3(a)为对训练样本的脊柱器官图像进行阈值分割后的三维展示图，将分割好的脊柱在Amira中剖分，获取三维点云数据，在matlab中展示，如图3(b)。

对训练样本的低对比度器官图像进行手动分割的结果如附图4。在具体实施例中低对比度器官选取的是肾脏图像，图4(a)为对训练样本的肾脏器官图像进行手动分割后的三维展示图，将分割好的脊柱图像在Amira中剖分，获取三维点云数据，在matlab中展示，如图4(b)。

对训练样本的高、低对比度器官三维点云数据对齐配准的结果分别如图5(a)和图6(a)所示。在得到脊柱、肾脏的三维点云数据后选择各自的配准模板，求解出相应的旋转R与平移T变换，把所有样本对齐配准到一起，然后计算各自的均值模型。图5(b)和图6(b)分别为脊柱和肾脏的均值模型。

步骤3：建立低对比度器官的统计形变模型：

计算低对比度器官图像配准后的协方差矩阵，进行主成分分析建立低对比度器官图像的统计形变模型：

(3a)计算低对比度器官形状向量的协方差矩阵S^L：

$S^L=\frac{1}{N-1}\underset{i=1}{\overset{N}{\&Sigma;}}({X_i}^L-{\stackrel{\&OverBar;}{X}}^L){({X_i}^L-{\stackrel{\&OverBar;}{X}}^L)}^T$

S^L为3K行乘3K列的矩阵，K为单个样本三维点云中点的个数。其中N为训练样本数，为低对比度器官图像均值模型，为低对比度器官图像对齐配准后的样本，T代表矩阵的转置运算；

(3b)对协方差矩阵S^L进行特征分解，求出对应的特征值λ_i(i＝1,2,…,3N)和特征向量P_i(i＝1,2,…,3N)。(3N为特征值和特征向量的个数)将特征值按绝对值从大到小排列，取前t(1≤t≤3K)个特征值所对应的特征向量来构建主成分标准正交基，用来近似表示任意的形状向量。t的值由下式而定：

$\frac{\underset{i=1}{\overset{t}{\&Sigma;}}{\mathrm{\&lambda;}}_i}{\underset{j=1}{\overset{3N}{\&Sigma;}}{\mathrm{\&lambda;}}_j}\&GreaterEqual;\mathrm{\&eta;}$

其中，η表示压缩后的模型能够解释原来模型中所有变化的比例，一般取80％～90％；

(3c)将主成分特征向量做单位化处理，记为P＝[p₁,p₂,…,p_t]，最终的低对比度器官图像统计形状模型表示为：

$X^L={\stackrel{\&OverBar;}{X}}^L+P^L{b}^L,$

其中X^L为任意形状，为步骤(2d)中所求低对比度器官图像均值模型，P^L为特征向量，b^L为形变参数，且为了保障形变在合理的范围内，形变参数需要控制在一定范围内：

$-3\sqrt{\mathrm{\&lambda;}}<b^L<3\sqrt{\mathrm{\&lambda;}}$

步骤4：基于高、低对比度器官图像之间相关性，借助高对比度器官辅助定位低对比度器官图像的初始位置

(4a)选取M个小鼠体数据作为待分割的测试样本，仅对其高对比度器官(脊柱)图像进行自动阈值分割，并进行有限元剖分，获取三维点云数据Test^H；

(4b)求解待分割样本中高对比度器官图像的Test^H与步骤(2d)中得到的高对比度器官图像均值模型之间的位置形状变换Trans^H，该变换是包括旋转、平移、缩放的相似性变换，计算公式为：

(4c)根据生物组织高、低对比度器官之间的相关性，将步骤(4b)中求解的Trans^H作为低对比度器官待测样本中相对位置形状的变换Trans^L，计算获取低对比度器官的初始位置形状Test^L。

因为Trans^H≈Trans^L，将其作用于肾脏的均值模型上，求出待测样本中低对比度器官初始位置Test^L，即：

步骤5：在步骤(4)初始定位的基础上，利用步骤(3)建立的低对比度器官图像统计形变模型Model^L沿标记点法线方向搜索查找最佳边界点，获取待分割样本中低对比度器官图像的最佳匹配模型，完成器官图像的分割。

采用本发明方法分割器官图像的分割结果如附图7。待分割测试样本的脊柱图像分割完成后，首先求解出待测样本的脊柱点云数据与脊柱均值模型之间相匹配的变换Trans^H，Trans^H≈Trans^L，作用于肾脏图像均值模型，找到其在待测样本中的初始化位置，利用建立的统计形变模型获取与待测样本匹配的最佳模型。图7(a)为测试样本的肾脏分割结果三维点云图，为了更好的对比分割的准确性，图7(b)把实际的肾脏位置与分割结果放在一起作对比，红色的点云数据是本发明方法得到的分割结果，绿色的点云数据是测试样本实际的肾脏点云数据，图7(c)是更直观地把两个点云放置在三维切片中显示的结果。图7(d)和(e)是从二维x-y平面上验证分割的切片图，其中红色的离散点是本发明的分割结果，绿色的离散点是实际肾脏的轮廓点。

以上所述，仅为本发明中的一个具体实例，但本发明的保护范围并不局限于此，任何熟悉该技术的人在本发明所揭露的技术范围内，可理解和想到的变换或替换，都应涵盖在本发明的包含范围之内，因此，本发明的保护范围应该以权利要求书的保护范围为准。

Claims (5)

1.一种基于统计形变模型的器官辅助定位分割方法，包括如下步骤：

(1)获取生物体CT数据；

(2)划分高、低对比度器官图像，在CT数据中选取训练样本提取相应的统计先验信息；

(3)建立低对比度器官图像的统计形变模型；

(4)基于高、低对比度器官之间的相关性，借助高对比度器官辅助定位低对比度器官的初始位置；

(5)在步骤(4)辅助定位的基础上，沿标记点法线方向进行器官的搜索分割。

2.根据权利要求1所述的基于统计形变模型的器官辅助定位分割方法，其特征在于：所述的步骤(1)，使用Micro-CT成像系统采集生物体多个角度的投影数据，利用滤波反投影方法对投影数据进行重建，获取生物体三维CT体数据；

所述的步骤(2)，根据X射线对生物体不同组织器官图像穿透能力的差异，将器官图像分为高、低对比度器官图像；选取训练样本数据，对其中高、低对比度器官图像分别进行分割，通过有限元剖分获取相应的三维点云数据；分别将高、低对比度器官的三维点云数据进行相似性变换对齐配准；配准后求取两者的均值模型；

所述的步骤(3)，计算步骤(2)中所得低对比度器官图像配准结果的协方差矩阵，并进行特征分解，获得对应的特征值和特征向量，对特征值的绝对值降序排列，选取前t个主成分，将主成分特征向量单位化处理，此时，低对比度器官的统计形变模型可被描述为其均值模型与特征向量和形变参数的线性组合；

所述的步骤(4)，选取待测样本，对其中高对比度器官阈值分割，有限元剖分后获取三维点云数据；求解待测样本中高对比度器官与步骤(2)求得的相应均值模型之间的位置形状变换；利用高对比度器官之间的变换辅助定位待测样本中低对比度器官的初始位置；

所述的步骤(5)，在待测数据中获取初始位置后，利用低对比度器官的统计形变模型搜索查找边界点，使其逐渐往目标区域逼近，搜索的区域范围选定为每个边界点沿当前轮廓的法线方向(即与该点相邻的两点连线的垂直方向上)，直至搜寻到最佳匹配的变形模型，完成最终分割。

3.根据权利要求1所述的基于统计形变模型的医学图像分割方法，其特征在于：所述步骤(2)的具体方法包括以下步骤：

(2a)根据X射线对生物体不同组织器官穿透能力的差异，将器官分为高、低对比度器官，组织密度较大X射线不易透过的器官为高对比度器官，组织密度较小X射线穿透性较大的器官为低对比度器官；

(2b)选取N个小鼠的三维体数据作为训练样本，对其中高、低对比度器官分别进行分割，前者采用自动阈值分割，后者采用手动分割，并将分割结果进行有限元离散剖分，获取相应的三维点云数据；

(2c)对步骤(2b)所得的高、低对比度器官三维点云数据在形状和位置上进行归一化，分别选取基准模板将所有样本进行旋转、平移等相似性变换，得到两者对齐配准结果；

(2d)将步骤(2c)配准对齐后的结果求和平均，分别计算高对比度器官均值模型和低对比度器官均值模型

其中，为第i个高对比度器官对齐配准后的样本，为第i个低对比度器官对齐配准后的样本，N为训练样本数，均为K行乘3列的矩阵，K为配准后单个样本三维点云中点的个数。

4.根据权利要求1所述的基于统计形变模型的医学图像分割方法，其特征在于：所述步骤(3)的具体方法包括以下步骤：

(3a)计算低对比度器官配准结果的协方差矩阵S^L：

S^L为3K行乘3K列的矩阵，K为单个样本三维点云中点的个数。其中N为训练样本数，为低对比度器官均值模型，为低对比度器官对齐配准后的样本，T代表矩阵的转置运算；

(3b)对协方差矩阵S^L进行特征分解，求解对应的特征值和特征向量，将特征值按绝对值降序排列，选取前t(1≤t≤3K)个特征值，以确保能够反应样本变化的主要模式，将主成分对应的特征向量单位化处理；

(3c)根据步骤(3b)所得，低对比度器官统计形变模型中的任意形状X^L均可描述为其均值模型与特征向量P^L和形变参数b^L的线性组合，即为

为了使形状变化限制在合理的范围内，形变参数b^L需满足：

其中，λ为第一主成分特征值。

5.根据权利要求1所述的基于统计形变模型的医学图像分割方法，其特征在于：所述步骤(4)的具体方法包括以下步骤：

(4a)选取M个小鼠的三维CT体数据作为待分割的测试样本，将待测样本中的高对比度器官进行自动阈值分割，并进行有限元离散剖分，获取相应的三维点云数据；

(4b)步骤(4a)所得待测样本的高对比度器官三维点云形状为Test^H，计算Test^H与高对比度器官均值形状之间的变换，即Trans^H，公式如下：

(4c)基于高、低对比度器官之间的相关性，将步骤(4b)中计算所得的高对比度器官之间的变换Trans^H近似为低对比度器官在待测样本中相对位置形状的变换Trans^L，即Trans^L≈Trans^H，结合低对比度器官均值形状模型求解低对比度器官初始位置Test^L，计算公式为：