基于二次三维区域生长的腹部器官分割方法

1.一种基于二次三维区域生长的腹部器官分割方法，包括如下步骤：
(1)输入一套完整的DICOM格式腹部CT切片，读取每幅切片中的像素数据，得到一套切片图像数据；
(2)根据腹部CT成像特点，截取每幅切片图像数据中的身体区域，得到一套身体区域切片，并将这些身体区域切片按切片序号从小到大存成三维数据体D1；
(3)对三维数据体D1中的每幅身体区域切片，用二维高斯滤波器去除其噪声，得到去噪三维数据体D2；
(4)利用Canny算法提取去噪三维数据体D2中每幅切片的边缘，得到三维边缘数据体Y；
(5)根据腹部CT图像中感兴趣器官的解剖位置及其灰度值分布，确定三维区域生长的初始种子点，结合三维边缘数据体Y和初始种子点的灰度信息，确定三维区域生长的生长规则，根据已确定的初始种子点和生长规则，对去噪三维数据体D2进行第一次三维区域生长，得到生长后的三维数据体D3：
所述的确定三维区域生长的初始种子点，按如下步骤进行：
(5a)从去噪三维数据体D2中选取连续10幅包含感兴趣器官面积较大的切片，将这些切片作为种子切片；
(5b)根据感兴趣器官在腹部的解剖位置，在选取的种子切片中确定种子区域的范围，并根据感兴趣器官的灰度值分布，选取待定初始种子点；
(5c)对待定初始种子点进行三维形态学腐蚀，并提取每幅种子切片中的最大连通区域作为初始种子点；
所述的确定三维区域生长的生长规则，是将满足||g(x)-v||＜cσ的去噪三维数据体中的像素x标记为生长点，但对与三维边缘数据体Y对应的去噪三维数据体D2中的那部分像素不做标记，其中g(x)表示像素x的灰度值，v表示初始种子点的均值，σ表示初始种子点的标准差，c为控制常数；
(6)采用三维形态学开运算，提取生长后的三维数据体D3的三维形态学边缘D4；
(7)结合三维形态学边缘D4，重新确定三维区域生长的生长规则，利用第一次三维区域生长的初始种子点，对去噪三维数据体D2重新进行三维区域生长，得到第二次三维区域生长后的三维数据体D5：
所述的重新确定三维区域生长的生长规则，是将满足||g(x)-v||＜cσ的去噪三维数据体中的像素x标记为生长点，但排除满足该条件的以下三种情况的像素：
7a)与三维边缘数据体Y对应的去噪三维数据体D2中的那部分像素；
7b)与三维形态学边缘D4对应的去噪三维数据体D2中的那部分像素；
7c)计算第一次三维区域生长后的三维数据体D3中的每幅切片上各个连通区域的面积，与面积小于阈值T的连通区域对应的去噪三维数据体D2中的那部分像素；
(8)采用三维形态学膨胀方法，对第二次三维区域生长后的三维数据体D5进行平滑处理，得到感兴趣器官的最终分割结果。

基于二次三维区域生长的腹部器官分割方法\n技术领域\n[0001] 本发明属于医学图像处理领域，涉及一种腹部器官分割方法，可用于腹部CT图像中，对感兴趣腹部器官如肝脏、脾脏和肾脏进行提取，以辅助临床医疗诊断。\n背景技术\n[0002] 腹部器官的分割具有重要的理论价值和广泛的临床应用前景。从复杂背景中提取感兴趣器官是进行三维可视化的前提和基础。更重要的是，确定感兴趣器官的病变位置和区域有助于进行外科手术、放射性治疗等。医学图像的对比度低，边缘模糊，以及不同病人器官形状的不确定性，使得长期以来医学图像的分割大部分集中在人机交互水平，处理时间长，而且处理结果易受人为因素影响。因此，图像自动分割一直是医学图像处理的研究重点和难点。\n[0003] John E.Koss等人1999年在《IEEE Transactions On Medical Imaging》杂志上提出利用Hopfield网络对图像的纹理特征进行聚类分割，分割得到的器官取决于Hopfield网络中指定的聚类数目。该方法需要多次迭代才能收敛，时间消耗大，得到的区域往往不连续。Hyunjin Park等人2003年在《IEEE Transactions On Medical Imaging》杂志上提出利用腹部器官的解剖学形态图谱对图像作配准，从而估计器官灰度值分布的高斯模型参数，并利用贝叶斯方法得到分割结果。该方法需要手工设置多个控制点，影响了方法的自动化。周永新等人2005年在《IEEE Transactions On Information Technology In Biomedical》杂志上提出了一种腹部器官的自动分割方法，该方法通过将图像和解剖学形态图谱做配准，并采用模糊连接分割来提取感兴趣器官。但解剖学图谱的建立需要用到大量的数据，也会花费大量的时间。意大利的Paola Campadelli等人在2007年IEEE会议上发表的论文《Automatic segmentation of abdominal organs from CT scans》中，提出了利用三维区域生长法进行腹部器官的自动分割。该方法能比较有效的提取出感兴趣器官，但是由于腹部器官的灰度值与其周围组织的灰度值很相近，使结果出现过分割现象；另外，该方法中各个感兴趣器官的分割是按其解剖位置从上到下依次进行的，上一个器官的分割结果直接影响后续器官的分割；而且对各个感兴趣器官分割的具体实施步骤也大不相同，缺少一个适用于所有感兴趣器官分割的统一框架，使得该方法缺少鲁棒性。\n发明内容\n[0004] 本发明的目的在于克服上述已有技术的不足，提供一种基于二次三维区域生长的腹部器官提取方法，以实现对感兴趣器官的准确分割，并有效地提高了腹部器官分割方法的鲁棒性。\n[0005] 为实现上述目的，本发明包括如下步骤：\n[0006] (1)输入一套完整的DICOM格式腹部CT切片，读取每幅切片中的像素数据，得到一套切片图像数据；\n[0007] (2)根据腹部CT成像特点，截取每幅切片图像数据中的身体区域，得到一套身体区域切片，并将这些身体区域切片按切片序号从小到大存成三维数据体D1；\n[0008] (3)对三维数据体D1中的每幅身体区域切片，用二维高斯滤波器去除其噪声，得到去噪三维数据体D2；\n[0009] (4)利用Canny算法提取去噪三维数据体D2中每幅切片的边缘，得到三维边缘数据体Y；\n[0010] (5)根据腹部CT图像中感兴趣器官的解剖位置及其灰度值分布，确定三维区域生长的初始种子点，结合三维边缘数据体Y和初始种子点的灰度信息，确定三维区域生长的生长规则，根据已确定的初始种子点和生长规则，对去噪三维数据体D2进行第一次三维区域生长，得到生长后的三维数据体D3：\n[0011] 所述的确定三维区域生长的初始种子点，按如下步骤进行：\n[0012] (5a)从去噪三维数据体D2中选取连续10幅包含感兴趣器官面积较大的切片，将这些切片作为种子切片；\n[0013] (5b)根据感兴趣器官在腹部的解剖位置，在选取的种子切片中确定种子区域的范围，并根据感兴趣器官的灰度值分布，选取待定初始种子点；\n[0014] (5c)对待定初始种子点进行三维形态学腐蚀，并提取每幅种子切片中的最大连通区域作为初始种子点；\n[0015] 所述的确定三维区域生长的生长规则，是将满足||g(x)-v||＜cσ的去噪三维数据体中的像素x标记为生长点，但对与三维边缘数据体Y对应的去噪三维数据体D2中的那部分像素不做标记，其中g(x)表示像素x的灰度值，v表示初始种子点的均值，σ表示初始种子点的标准差，c为控制常数；\n[0016] (6)采用三维形态学开运算，提取生长后的三维数据体D3的三维形态学边缘D4；\n[0017] (7)结合三维形态学边缘D4，重新确定三维区域生长的生长规则，利用第一次三维区域生长的初始种子点，对去噪三维数据体D2重新进行三维区域生长，得到第二次三维区域生长后的三维数据体D5：\n[0018] 所述的重新确定三维区域生长的生长规则，是将满足||g(x)-v||＜cσ的去噪三维数据体中的像素x标记为生长点，但排除满足该条件的以下三种情况的像素：\n[0019] 7a)与三维边缘数据体Y对应的去噪三维数据体D2中的那部分像素；\n[0020] 7b)与三维形态学边缘D4对应的去噪三维数据体D2中的那部分像素；\n[0021] 7c)计算第一次三维区域生长后的三维数据体D3中的每幅切片上各个连通区域的面积，与面积小于阈值T的连通区域对应的去噪三维数据体D2中的那部分像素；\n[0022] (8)采用三维形态学膨胀方法，对第二次三维区域生长后的三维数据体D5进行平滑处理，得到感兴趣器官的最终分割结果。\n[0023] 本发明与现有技术相比具有以下优点：\n[0024] 1.本发明由于根据腹部CT图像中感兴趣器官的解剖位置及其灰度值分布自动地确定三维区域生长的初始种子点，克服了现有方法需要人工指定初始种子点的缺点；\n[0025] 2.本发明由于结合Canny边缘进行第一次三维区域生长，并结合了第一次三维区域生长得到的分割结果的形态学边缘进行第二次三维区域生长，有效地抑制了传统三维区域生长方法出现的过分割现象；\n[0026] 3.本发明自动地提取各个感兴趣器官的初始种子点，对不同器官进行完全独立地分割，不会出现现有分割方法中提取不同器官时出现的累加错误，提高了方法的鲁棒性。\n附图说明\n[0027] 图1为本发明的主要操作过程示意图；\n[0028] 图2为腹部感兴趣器官种子区域在CT图像中的位置示意图；\n[0029] 图3为腹部感兴趣器官的灰度直方图统计结果图；\n[0030] 图4为腹部感兴趣器官的单幅分割结果示意图；\n[0031] 图5为腹部感兴趣器官的分割结果三维重建效果图。\n具体实施方式\n[0032] 参照图1，本发明的具体实施步骤如下：\n[0033] 步骤1：输入待分割图像。\n[0034] 输入一套包含完整的肝脏、脾脏和肾脏器官的DICOM格式的腹部CT切片，由于DICOM格式的图像带有大量与器官分割无关的信息，为了减少数据的存储空间，这里只读取与图像有关的信息，主要包括每幅切片的像素数据和后续处理中所需的切片序列标识号信息，从而得到一套待分割的切片图像数据。\n[0035] 步骤2：图像预处理。\n[0036] 对待分割的切片图像数据进行预处理，主要包括：截取身体区域和去除图像噪声。\n[0037] 截取身体区域：由于腹部CT图像中包含大量非身体区域的像素，为了进一步减少存储图像的空间和处理时间，这里只需截取身体的大致区域，可采用以下简单方法：任意指定一幅切片图像，通过沿不同方向进行扫描找出其最上、最下、最左、最右的非零点坐标。由最上和最下非零点纵坐标，以及最左和最右非零点横坐标确定的区域即为该幅切片上的身体区域。由于螺旋CT序列成像时速度较快，因此切片上身体区域的位置基本固定，所以其他切片上身体区域的位置都采用该幅切片上的位置。这样，便得到一套身体区域切片。由于切片的输入顺序不同会对分割结果造成影响，需要对输入的所有切片进行排序。每幅切片都自带有与实际身体位置相对应的序列标识号信息，对输入的图像按该序列号从小到大进行排序，并存成三维数据体D1。\n[0038] 去除图像噪声：医学图像通常都含有大量的噪声，为了减少图像噪声对分割结果带来的影响，需要对图像进行去噪处理。由于高斯滤波器对于去除各种类型的噪声都有很好的效果，为此采用该滤波方法去除图像噪声。利用高斯函数的可分离性，二维高斯滤波可通过一维高斯模板函数逐次沿水平方向和沿垂直方向的卷积来实现。其具体实施步骤如下：\n[0039] 2a.设计一维高斯模板：这里采用的模板函数为 模板\n窗口大小取为7σ，其中x表示像素值，取[-3σ，3σ]内的整数，m表示模板窗口中像素x的平均值，σ表示平滑尺寸，这里取σ＝2.0；\n[0040] 2b.利用一维高斯模板，对三维数据体D1中的每幅身体区域切片沿x方向进行卷积，得到第一次卷积后图像；\n[0041] 2c.利用一维高斯模板，对第一次卷积后图像，沿y方向进行卷积，得到最终滤波结果，并将其存为去噪三维数据体D2。\n[0042] 步骤3：利用Canny算法提取去噪三维数据体D2中每幅切片的边缘，得到三维边缘数据体Y。\n[0043] 这里希望能得到尽量精确的图像边缘，而Canny边缘检测算法正能做到这一点。\n该方法不容易受噪声的干扰，能够检测到真正的弱边缘。医学图像的边缘模糊，采用Canny算法能较好的提取感兴趣器官的边缘。Canny边缘检测算法的第一步是对图像做高斯滤波，而步骤2中已经对每幅切片利用高斯滤波器去除图像噪声，所以不需重新进行这一处理，只需对去噪三维数据体D2中的每幅切片图像进行Canny算法的其他后续处理，便得到三维边缘数据体Y。\n[0044] 步骤4：结合三维边缘数据体Y，对去噪三维数据体D2进行第一次三维区域生长，得到生长后的三维数据体D3。具体实施如下：\n[0045] 4a.确定初始种子点：\n[0046] (4a1)从去噪三维数据体D2中选取连续10幅包含感兴趣器官面积较大的切片，将这些切片作为种子切片。\n[0047] (4a2)根据感兴趣器官在腹部的解剖位置，在选取的种子切片中确定种子区域的范围，并根据各感兴趣器官的灰度值分布，选取待定初始种子点，具体实施过程如下：\n[0048] 由于感兴趣器官在腹部的解剖位置是相对固定的，所以在腹部CT图像中感兴趣器官的位置也是相对固定的。比如肝脏位于腹部CT图像的左上方，脾脏位于腹部CT图像的右下方，两个肾脏分别位于腹部CT图像的左、右下方。根据这些信息，在种子切片上选取种子区域，选取的感兴趣器官的种子区域如图2所示。图2中由白线画出的最大矩形所围区域是由步骤2截取的身体区域，x表示身体区域的长度，y表示身体区域的宽度，其中图(a)中由白线画出的较小矩形所围区域是肝脏的种子区域，图(b)中由白线画出的较小矩形所围区域是脾脏的种子区域，图(c)中由白线画出的两个较小矩形所围区域是两个肾脏的种子区域，左边矩形对应的是右肾脏种子区域，右边矩形对应的是左肾脏种子区域。肝脏种子区域的最左位置选在距身体区域最左位置x*3/16处，最上位置选在距身体区域最上位置y*3/16处，肝脏种子区域的长度为x/4，宽度为y/2。脾脏种子区域的最右位置选在距身体区域最右位置x/8处，最下位置选在距身体区域最下位置y/5处，脾脏种子区域的长度为x*5/16，宽度为y/3。右肾脏种子区域的最左位置选在距身体区域最左位置x/6处，左肾脏种子区域的最右位置选在距身体区域最右位置x/6处，两个肾脏的最下位置都选在距身体区域最下位置y/3处，两个肾脏种子区域的长度都为x/4，宽度都为y/3。\n[0049] 统计由上述方法得到的各个感兴趣器官种子区域内的灰度直方图，统计结果如图\n3所示，其中图(a)为肝脏的灰度直方图，图(b)为脾脏的灰度直方图，图(c)为肾脏的灰度直方图。根据灰度直方图，取种子区域内满足T1≤g(x)≤T2的像素x作为待定的初始种子点，其中g(x)表示像素x的灰度值，T1为最小阈值，T2为最大阈值。参照图3，对于肝脏，根据图(a)取T1为110，T2为210；对于脾脏，根据图(b)取T1为10，T2为170；对于肾脏，根据图(c)取T1为20，T2为200；\n[0050] (4a3)通过上述处理，待定初始种子点中可能仍含有不是待分割器官的像素，需要做进一步的处理。具体做法是：对待定的初始种子点，采用3×3×3的结构元素进行三维形态学腐蚀处理，计算腐蚀后每幅种子切片上各个连通区域的面积，并提取最大连通区域内的像素作为最终的初始种子点。\n[0051] 4b.确定生长规则：\n[0052] 将满足||g(x)-v||＜cσ的去噪三维数据体中的像素x标记为生长点，但对与三维边缘数据体Y对应的去噪三维数据体D2中的那部分边缘像素不做标记，其中g(x)表示像素x的灰度值，v表示初始种子点的均值，σ表示初始种子点的标准差，c为控制常数。\n[0053] 4c.利用已确定的初始种子点和生长规则，对去噪三维数据体D2进行三维区域生长，具体生长过程如下：\n[0054] (4c1)将初始种子点存入堆栈中；\n[0055] (4c2)从堆栈中取出一个像素，搜索其周围26邻域中的像素，将满足生长规则的像素标记为生长点，并将其存入堆栈中，这里的26邻域是指当前像素在所处切片上的8邻域，与当前像素位置相对应的上层切片上的像素及其周围8邻域，再加上与当前像素位置相对应的下层切片上的像素及其周围8邻域；\n[0056] (4c3)重复步骤(4c2)，直到堆栈为空，或者所有满足生长规则的像素都已标记为生长点，则停止生长，那些标记为生长点的像素构成生长后的三维数据体D3。\n[0057] 步骤5：采用结构元素为3×3×3的三维形态学开运算，提取第一次三维区域生长后的三维数据体D3的三维形态学边缘D4；\n[0058] 步骤6：结合三维边缘数据体Y和三维形态学边缘D4，对去噪三维数据体D2进行第二次三维区域生长，得到生长后的三维数据体D5。其具体实施过程是：以步骤4a确定的初始种子点为第二次三维区域生长的初始种子点，重新确定生长规则，按照步骤4c进行生长，得到生长后的三维数据体D5。其中重新确定的生长规则为：将满足||g(x)-v||＜cσ的去噪三维数据体中的像素x标记为生长点，但排除满足该条件的以下三种情况的像素：\n[0059] 情况1：与三维边缘数据体Y对应的去噪三维数据体D2中的那部分像素；\n[0060] 情况2：与三维形态学边缘D4对应的去噪三维数据体D2中的那部分像素；\n[0061] 情况3：计算第一次三维区域生长后的三维数据体D3中每幅切片上各个连通区域的面积，与面积小于阈值T的连通区域对应的去噪三维数据体D2中的那部分像素，其中g(x)表示像素x的灰度值，v表示初始种子点的均值，σ表示初始种子点的标准差，c为控制常数。\n[0062] 步骤7：采用结构元素为5×5×5的三维形态学膨胀，对第二次三维区域生长后的三维数据体D5进行平滑处理，得到腹部器官的最终分割结果。\n[0063] 本发明的效果通过以下仿真实验进一步说明。\n[0064] 1.实验条件及内容：\n[0065] 本实验数据来自于北京肿瘤医院64排螺旋CT扫描图像，成像层厚0.625mm，每幅切片都是512×512的16位DICOM格式图像。总共测试了3套腹部CT切片数据，每套数据大约都有350幅切片。各个感兴趣器官的分割实验都在VC++6.0的MFC环境下进行，并在Matlab环境下对最终分割结果做了三维重建，由于切片数太多，不能在Matlab环境下对一套完整的切片分割结果都进行三维重建，在此采用隔层重建的方法。\n[0066] 2.实验结果及结果分析\n[0067] 对上述的一套腹部CT切片数据，用本发明所述方法进行各个感兴趣器官的分割，实验结果如图4和图5所示。\n[0068] 图4是分割结果的单幅显示，这里只分别给出一幅肝脏、脾脏和肾脏的分割结果图。其中图(a)为肝脏分割的原图，图(b)为肝脏分割结果图，图(c)为脾脏分割的原图，图(d)为脾脏分割结果图，图(e)为肾脏分割的原图，图(f)为肾脏分割结果图。\n[0069] 图5是分割结果的三维重建显示，其中图(a)是肝脏分割结果的三维重建效果图，图(b)是脾脏分割结果的三维重建效果图，图(c)是肾脏分割结果的三维重建效果图。\n[0070] 从两种形式的分割结果图可以看出：本发明能够对腹部CT切片中的感兴趣器官进行准确地分割，有效地抑制了过分割现象。