Study-基于分水岭分割进行肺癌诊断

案例背景

  医学CT图像处理主要是研究医学图像中的器官和组织之间的关系，进行病理性分析。因此，借助计算机及图像处理技术对CT图像中医生所关注的区域进行精确的分割和定位是医学图像处理的关键步骤，在临床诊断中对于协助医生进行病理研判具有重要意义。

  分水岭分割是一种强有力的图像分割方法，可以有效地提取图像中我们所关注的区域。在灰度图像中使用分水岭方法可以将图像分割成不同的区域，每个区域就可能对应一个我们所关注的对象。对于这些图像的子区域可以进行进一步的处理。 除此之外，分水岭方法还可以提取目标的轮廓等特征。

  本案例通过对标记分水岭图像分割方法的实验进行改进，提出了一种简捷、高效的肺 部CT图像实质分割方法。

理论基础

  分水岭算法以数学形态学图像处理为基础，属于基于区域的图像分割算法。该算法最初的思想来源于地形学，假设将一幅图像看作 - -个立体的地形表面，且图像中每个像素点的海拔高度由该点的灰度值表示，则图像中的每个局部极小值和它的影响区域被称为集水盆地;图像的边缘灰度值变化幅度较大，对应于集水盆地的边界，即地形学中的山脊;图像的边缘位置扫描到灰度极小值点的渐变过程具有坡度的特点，因此被称为山坡口。

1.模拟浸水的过程

  假设地形图具有高低起伏的特征，对其建立地理模型并模拟浸水的过程。随着水从各个盆地的最低点向上漫溢，水面也逐步上升，相邻盆地的水将汇合于其边缘处。如果在地形图的两个集水盆汇合处修筑足够高的大坝来阻止水面的汇合，则随着水面的上升，各个盆地会完全被水淹没，但边缘处的大坝不会被淹没。因此，大坝就将地理模型分割成了不同的盆地区域，这些区域边缘建起的大坝被称为分水岭，被分割开的各个区域被称为聚水盆地，这就是模拟浸水的过程。

2.模拟降水的过程

  如果将图像视作地形图并建立地理模型，则当上空落下一滴雨珠时，雨珠降落到山体表面并顺山坡向下流，直到汇聚到相同的局部最低点。在地形图中，雨珠在山坡上经过的路线就是一个连通分支，通往局部最低点的所有连通分支就形成了-一个聚水盆地，那么山坡就被称为分水岭，这就是模拟降水的过程。

3.过度分割问题

  分水岭变换的目标是求出梯度图像的“分水岭线”，传统的差分梯度算法对近邻像素做差分运，容易受到噪声和量化误差等因素的影响，往往会在灰度均匀的区域内部产生过多的局部梯度“谷底”，这些在分水岭变换中就对应于“集水盆地”。因此，传统梯度算法最终将导致出现“过分割(Over Segmentation)”的现象，即一个灰度均匀的区域可能被过度分成多个子区域，以致产生大量的虚假边缘，从而无法确认哪些是真正边缘，这对算法的准确性造成了一定的不利影响，这就是过度分割问题。

4. 标记分水岭分割算法

  直接应用分水岭分割算法的主要缺点是会产生过分割现象，即分割出大量的细小区域，而这些区域对于图像分析可以说是毫无意义的。图像噪声等因素往往会导致图像中出现很多杂乱的低洼区域，而平滑滤波能减少局部最小点的数目，所以分割前先对图像进行平滑是避免过分割现象的有效方法之一-。此外，对分割后的图像按照某种准则进行相邻区域的合并也是一种过分割解决方法。

  基于标记(Marker) 的分水岭分割算法能够有效控制过分割现象的发生，该算法的标记包括内部标记( Internal Marker) 和外部标记( External Marker)。 其基本思想是通过引入标记来修正梯度图，使得局部最小值仅出现在标记的位置，并设置阅值h来对像素值进行过滤，删除最小值深度小于阈值h的局部区域。

  标记分水岭算法中的一一个标记对应于图像的-一个连 通成分，其内部标记与我们感兴趣的某个目标相关，外部标记与背景相关。标记的选取--般包括预处理和定义选取准则两部分，其中，选取准则可以是灰度值、连通性、大小、形状、纹理等特征。选取内部标记之后，就能以其为基础对低洼进行分割，将分割区域对应的分水线作为外部标记，之后对每:个分割出来的区域利用其他分割技术( 如二值化分割)将目标从背景中分离出来。

  首先，假设将内部标记的选取准则定义为满足以下条件。

(1)区域周围由更高的“海拔”点组成。

(2)区域内的点可以组成一个连通分量。

(3)区域内连通分量的点具有相同或相近的灰度值。

然后，对平滑滤波后的图像应用分水岭算法，并将满足条件的内部标记为所允许的局部最小值，再将分水岭变换得到的分水线结果作为外部标记。最后，内部标记对应于每个感兴趣目标的内部，外部标记对应于背景。根据这些标记结果将其分割成互不重叠的区域，每个区域包含唯一的目标和背景。

  因此，标记分水岭算法的显著特点和关键步骤就是获取标记的过程，本案例将采用基于梯度边缘检测与标记分水岭算法相结合的方法来对肺部图像进行分割。

程序实现

  分水岭算法的主要目标在于找到图像的连通区域并进行分割。在实际处理过程中，如 果直接以梯度图像作为输入，则容易受到噪声的干扰，产生多个分割区域;如果对原始图 像进行平滑滤波处理后再进行梯度计算，则容易将某些原本独立的相邻区域合成一个区 域。当然，这里的区域主要还是指图像内容变化不大或者灰度值相近的连通区域。为了易 于调用，本案例将标记分水岭分割算法封装为子函数，具体代码如下所示。

function Watershed_Fun(fileName)
rgb = imread(fileName);
if ndims(rgb) == 3
    I = rgb2gray(rgb);
else    
    I = rgb;
end
sz = size(I);
if sz(1) ~= 256
    I = imresize(I, 256/sz(1));
    rgb = imresize(rgb, 256/sz(1));
end
hy = fspecial('sobel');
hx = hy';
Iy = imfilter(double(I), hy, 'replicate');
Ix = imfilter(double(I), hx, 'replicate');
gradmag = sqrt(Ix.^2 + Iy.^2);
se = strel('disk', 3);
Io = imopen(I, se);
Ie = imerode(I, se);
Iobr = imreconstruct(Ie, I);
Ioc = imclose(Io, se);
Iobrd = imdilate(Iobr, se);
Iobrcbr = imreconstruct(imcomplement(Iobrd), imcomplement(Iobr));
Iobrcbr = imcomplement(Iobrcbr);
fgm = imregionalmax(Iobrcbr);
se2 = strel(ones(3,3));
fgm2 = imclose(fgm, se2);
fgm3 = imerode(fgm2, se2);
fgm4 = bwareaopen(fgm3, 15);
bw = im2bw(Iobrcbr, graythresh(Iobrcbr));
D = bwdist(bw);
DL = watershed(D);
bgm = DL == 0;
gradmag2 = imimposemin(gradmag, bgm | fgm4);
L = watershed(gradmag2);
Lrgb = label2rgb(L, 'jet', 'w', 'shuffle');
[pathstr, name, ext] = fileparts(fileName);
filefolder = fullfile(pwd, '实验结果', [name, '_实验截图']);
if ~exist(filefolder, 'dir')
    mkdir(filefolder);
end
h1 = figure(1);
set(h1, 'Name', '图像灰度化', 'NumberTitle', 'off');
subplot(1, 2, 1); imshow(rgb, []); title('原图像');
subplot(1, 2, 2); imshow(I, []); title('灰度图像');
fileurl = fullfile(filefolder, '1');
set(h1,'PaperPositionMode','auto');
print(h1,'-dtiff','-r200',fileurl);
h2 = figure(2);
set(h2, 'Name', '图像形态学操作', 'NumberTitle', 'off');
subplot(1, 2, 1); imshow(Iobrcbr, []); title('图像形态学操作');
subplot(1, 2, 2); imshow(bw, []); title('图像二值化');
fileurl = fullfile(filefolder, '2');
set(h2,'PaperPositionMode','auto');
print(h2,'-dtiff','-r200',fileurl);
h3 = figure(3);
set(h3, 'Name', '图像梯度显示', 'NumberTitle', 'off');
subplot(1, 2, 1); imshow(rgb, []); title('待处理图像');
subplot(1, 2, 2); imshow(gradmag, []); title('梯度图像');
fileurl = fullfile(filefolder, '3');
set(h3,'PaperPositionMode','auto');
print(h3,'-dtiff','-r200',fileurl);
h4 = figure(4); imshow(rgb, []); hold on;
himage = imshow(Lrgb);
set(h4, 'Name', '图像分水岭分割', 'NumberTitle', 'off');
set(himage, 'AlphaData', 0.3);
hold off;
fileurl = fullfile(filefolder, '4');
set(h4,'PaperPositionMode','auto');
print(h4,'-dtiff','-r200',fileurl);

通过对病变的肺部图像和正常的肺部图像分别进行标记分水岭处理，可得到分割结果 并进行彩色标记显示。

在实验过程中，为了避免传统分水岭算法所产生的过分割问题，函 数应用了对梯度图像进行形态学滤波处理的技巧，有效地提高了算法的运行效率和准确 性。实验表明，标记分水岭分割算法对肺部图像进行分割具有良好的效果，能在一定程度 上突出病变区域，起到辅助医学诊断的目的，具有一定的使用价值。