爱华网[摘要]：图像分割是一种重要的数字图像处理技术。本文首先介绍了图像分割技术，其次总结了目前图像分割技术中所用到的阈值、边缘检测、区域提取等方法以及分水岭算法。针对各种阈值分割算法，本文在最后做了详细的研究，并给出了相应的MATLAB程序源代码。
[关键词]：数字图像处理；图像分割；阈值；算法
ResearchonMATLAB-BasedThresholdImageSegmentationAlgorithm
Abstract：Imagesegmentationisaveryimportantdigitalimageprocessingtechnique.Inthisthesis,thedigitalimageprocessingtechniquewillbeillustratedinthefirstplace.Besides,itsummarizesthemethodofimagesegmentation,suchasthreshold,edgedetection,regionextraction,andwatershedalgorithm.Lastbutnotleast,itmakesastudyofvariousthresholdimagesegmentationalgorithms.Meanwhile,thecodeofMATLABwillbeshowed.
Keywords：DigitalImageProcessingTechnique;ImageSegmentation;Threshold;Algorithm
引言
随着科学技术的发展，人们对信息处理和信息交流的要求越来越高，最典型的表现是手机及相关数码产品的更新换代速度加快。图像信息具有直观、形象、易懂和信息量大等特点，因此它是人们日常生活、生产中接触最多的信息种类之一。正是因为图像所带给人们的直观信息的特点，使得数字图像处理技术随着计算机技术、多媒体技术的发展取得了长足的进步。
在图像分析中，通常需要将关心的目标从图像中提取出来，这种从图像中将某个区域与其他部分进行分离并提取出来的处理，就是图像分割。因为图像分割实际上就是区分图像中的“前景目标”和“背景”，所以通常又称之为图像的二值化处理。图像分割是图像处理中的一项关键技术，自20世纪70年代起一直受到人们的高度重视，至今已提出上千种分割算法，但因尚无通用的分割理论，现提出的分割算法大都是针对具体问题的，并没有一种适合所有图像的通用分割算法。另外，还没有制定出选择适用分割算法的标准，这给图像分割技术的应用带来许多实际问题。最近几年又出现了许多新思路、新方法或改进算法。总的来说，图像分割是图像识别和图像分析的基本前提步骤，图像分割的质量好坏直接影响后续图像处理的效果，甚至决定成败。因此，图像分割在数字图像处理技术中占有非常重要的地位。
1.图像分割概述
1.1.图像分割的定义
图像分割时指将一副图像分解为若干互不交叠的、有意义的、具有相同性质的区域。好的图像分割应具备以下特征：
⑴分割出来的各个区域对某种特性（例如灰度和纹理）而言具有相似性，区域内部是连通的且没有过多小孔
⑵相似区域对分割所依据的性质有明显的差异
⑶区域边界是明确的
图像分割更形式化的定义如下：令表示图像，表示具有相同性质的谓词，图像分割把分解成个区域，，满足：
⑴，Ф，，，
⑵，，
⑶，，，
条件⑴表明分割区域要覆盖整个图像且各区域互不重叠，条件⑵表明每一个区域都具有相同的性质，条件⑶表明相邻的两个区域性质相异不能合并为一个区域。
实际的图像处理和分析都是面向某种具体应用的，所以上述条件中的各种关系也要视具体情况而定。目前，还没有一种通用的方法可以很好的兼顾这些约束条件，也没有一种通用的方法可以完成不同的图像分割任务。原因在于实际的图像时千差万别的，还有一个重要原因在于图像数据质量的下降，包括图像在获取和传输过程引入的各种噪声以及光照不均与等因素。到目前为止，对图像分割的好坏和评价还没有统一的准则。因此，图像分割使图像分析和计算机视觉中的经典难题。至今，提出的分割算法已有上千种，每年还有不少新算法出现。这些算法的实现方式各不相同，然而大都基于图像在像素级的两个性质：不连续性和相似性。属于同一目标的区域一般具有相似性，而不同区域在边界出现不连续性。
1.2.图像分割的方法
图像分割方法依照分割时所依据的图像特性不同，大致可以分为四大类。第一类是阈值方法，这种方法是根据图像灰度值得分布特性确定某个阈值来进行图像分割的。第二类为边缘检测方法，这类方法是通过检测出封闭某个区域的边界来进行图像分割的。通俗地讲，这类方法实际上是沿着闭合的边缘线将其包围的区域剪切出来。第三类方法是区域提取方法，这类方法的特点是根据特定区域与其他背景区域特性上的不同来进行图像分割的。另外，还有一种基于形态学的分水岭算法。
1.2.1.阈值方法
所谓阈值分割方法就是确定某个阈值，根据图像中每个像素的灰度值大小或小于该阈值，来进行图像分割。阈值方法的数学模型如下：
设原图像为，经过分割处理后的图像为，为二值图像，则有：
根据上式可知，阈值方法的核心就是阈值的确定方法。
1.2.2.边缘检测方法
边缘检测技术对于处理数字图像非常重要，因为边缘时所要提取目标和背景的边界线，提取出边缘才能将目标和背景区分出来。在图像中，边界表明一个特征区域的始终和另一个特征区域的开始，边界所分开区域的内部特征或属性是一致的，而不同区域内部的特征或属性是不同的，边缘检测正是利用物体和背景在某种图像特性上的差异来实现的，这些差异包括灰度、颜色或纹理特征。边缘检测实际上就是检测图像特性发生变化的位置。
最简单的边缘检测方法是边缘检测算子，它利用相邻区域的像素值不连续的性质，采用一阶或二阶导数来检测边缘点。近年来还提出了基于曲面拟合的方法，基于边界曲线拟合的方法。
1.2.3.区域提取方法
数字图像的像素分布在一定程度上可以反映图像内容的特征，利用像素值的分布特征进行图像分类成了目前对图像进行分类的常用手段。图像分类的性能主要取决于图像的特征提取，数字图像的特征提取主要从颜色、纹理、形状等几个方面提取图像的特征。
1.2.4.分水岭算法
分水岭算法是一种借鉴了形态学的分割算法，在该方法中，将一副图像看成是一个拓扑地图其中灰度值对应的高度图。高灰度值对应着山峰，低灰度对应着山谷。水总是朝着低的地方流动，直到某一局部低洼处才停下来，这个低洼处被称为吸水盆地。最终所有的水会分聚在不同的吸水盆地，吸水盆地之间的山脊称为分水岭。水从分水岭流下时，它朝不同可能的吸水盆地流去的可能性是相等的。将这种想法用于图像分割，就是要在灰度图像中找出不同的吸水盆地和分水岭，由这些不同的吸水盆地和分水岭组成的区域即为要分割的目标。MATLAB图像处理工具箱中的watershed函数可以用于显示分水岭算法。
2.基于MATLAB的图像阈值分割算法的研究
2.1.人工阈值选择法
阈值分割最简单的方法就是人工选择法。基于灰度阈值的分割方法，其关键是如何合理的选择阈值。人工选择方法是通过人眼的观察，应用人对图像的知识，在分析图像直方图的基础上，人工选择出合理的阈值。也可以在人工选择出阈值后，根据分割的效果，不断地进行交互操作，从而选择出最佳的阈值。例如，分析图2.1，做出其灰度直方图（见图2.2），这里选择0.3为阈值，则可得到较为理想的分割结果（见图2.3）。
图2.1国际标准测试图片图2.3阈值分割后的图片
图2.2图2.1的灰度直方图
显然，灰度直方图的峰谷阈值方法是一种有效且非常简单的阈值方法，但是该方法有一个局限性，就是要求图像的灰度直方图必须具有双峰型。
2.2.自动阈值选择算法
虽然人工选择方法可以选取出令人满意的阈值，但是在无人介入的情况下自动选区阈值是大部分应用的要求，自动阈值选择方法通常使用灰度直方图来分析图像中灰度值的分布，结合特定的应用领域知识来选取最合适的阈值。
2.2.1.-参数法
-参数法是针对预先已知图像中目标物所占比例的情况下，所采用的一种简单且有效的方法。-参数法的设计思想是，选择一个值，使前景目标物所占的比例为，背景所占比例为。
根据上面的原理：-参数法的具体步骤如下：
①首先获得理想状态下的目标物所占画面的比例
②尝试性地给定一个阈值
③计算在下判定的目标物的像素点数，此时目标物所占比例为
④当足够接近时，此时得到最佳分割阈值
2.2.2.迭代法
迭代法的的设计思想是，开始时选择一个阈值作为初始估计值，然后按某种策略不断的改进这一估计值，直到满足给定的准则为止。在迭代过程中，关键之处在于选择什么样的阈值改进策略。好的改进策略应该具备两个特征：一是能够快速收敛，二是在每一个迭代过程中，新产生的阈值优于上一次的阈值。下面介绍一种迭代法：
①选择图像灰度的中值作为初始阈值
②利用阈值把图像分割为两个区域和，用下式计算区域和的灰度均值和
，为第类中的像素个数
③计算出和后，用下式计算出新的阈值
④重复②和③，直到和的差小于某个特定的值
2.2.3.均匀性度量法
均匀性度量法的设计思想是，假设当图像被分为目标物和背景两个类别时，属于同一类别内的像素值分布应该具有均匀性。在这里采用方差来度量像素间的均匀性。设原图像为，结果图像为。通过图像分割将原图像分为和（即背景与目标）两类，则算法步骤如下：
①给定一个初始阈值
②分别计算两类中的方差和
其中，为第类中的像素个数
③分别计算两类在图像中的分布概率和
④选择最佳的阈值，使得图像按照该阈值分为和两类后，满足
2.2.4.类间最大距离法
类间最大距离法德设计思想是，在某个适当的阈值下，图像分割后的前景目标与背景两个类之间的差异最大为最佳阈值。在这里两个类别（目标与背景）的差异，用两个类别中心与阈值之间的距离差来度量。具体步骤如下：
①给定一个初始阈值，将图像分为和两类
②计算两类的灰度均值和
③计算两类的相对距离值
④选择最佳的阈值，使得图像按照该阈值分为和两类后，满足
2.2.5.Otsu法
Otsu是一种使类间方差最大的自动确定阈值的方法，该方法具有简单、处理速度快的特点，是一种常用的阈值选取方法。MATLAB中的graythresh函数求取阈值采用的就是Otsu法。具体步骤如下：
①给定一个初始阈值，将图像分为和两类
②计算两类的灰度均值和，以及图像的总体灰度均值
③计算两类的概率和
，
④计算类间方差
⑤选择最佳的阈值，使得图像按照该阈值分为和两类后，满足
2.2.6.最大熵法
熵是信息论中对不确定性的度量，是对数据中所包含信息量大小的度量。熵取最大值时，就表明获得的信息量为最大。最大熵法的设计思想是，选择适当的阈值将图像分为两类，两类的平均熵之和为最大时，可以从图像中获得最大信息量，以此来确定最佳阈值。具体步骤如下：
①求出图像中的所有像素的分布概率，图像灰度的分布范围为
②给定一个初始阈值，将图像分为和两类
③分别计算两个类的平均相对熵
④选择最佳的阈值，使得图像按照该阈值分为和两类后，满足
2.2.7.最大类间、类内方差比法
从统计意义上讲，方差是表征数据分布不均衡性的统计量，要通过阈值对两类问题进行分割，显然，适当的阈值使得两类数据间的方差越大越好，表明该阈值的确将两类不同的区域分开了，同时希望属于同一类的方差越小越好，表明同一类区域有一定的相似性。因此可以采用类内，类间方差比作为选择阈值的评价参数。具体步骤如下：
①求出图像中的所有像素的分布概率，图像灰度的分布范围为
②给定一个初始阈值，将图像分为和两类
③计算两类的方差和，和灰度均值和，以及图像的总体灰度均值
④计算两类的概率和
，
⑤计算类间方差和类内方差
⑥选择最佳的阈值，使得图像按照该阈值分为和两类后，满足
2.2.8.聚类方法
所谓聚类算法，是采用模式识别中的聚类思想，以类内保持最大相似性以及类内保持最大距离为目标，通过迭代优化获得最佳的图像分割阈值。具体步骤如下：
①给定一个初始的代表两个类别和的中心灰度值和，两类的方差和
②进行分类处理，对某个像素，若，则，否则，
③对和中的所有像素分别重新计算其中心点的值和，中心点的值取每个类的均值，同时计算两类的方差和
④如果，则返回②，否则将该分类结果作为最终的结果输出。其中和为新分类的两类像素分布概率
2.2.9.局部阈值方法
前面给出的阈值算法，对于较为简单的图像（即目标与背景比较容易区分），上面的方法简单有效。但是对于较为复杂的图像，则往往会产生一些问题。图2.4给出了国际标准测试图像（图2.1）的局部阈值分割结果。
图2.3局部阈值分割后的图片
2.3.各种阈值算法的评价
对于目前的各种阈值分割算法，尚无一种固定的分割评价标准。分割的好与坏取决于实际需求和技术人员的主观因素。
参考文献
[1]朱虹.数字图像处理基础[M].北京：科学出版社，2009.125～146.
[2]秦襄培.MATLAB图像处理与界面编程宝典[M].北京：电子工业出版社，2009.441～477.
附录（MATLAB程序源代码）
%迭代法
clc;
clearall;
f=imread('man.bmp');
f=rgb2gray(f);
f=im2double(f);
T=0.5*(min(f(:))+max(f(:)));
done=false;
while~done
g=f>=T;
Tn=0.5*(mean(f(g))+mean(f(~g)));
done=abs(T-Tn)<0.1;
T=Tn;
end
T
r=im2bw(f,T);
figure,imshow(f),title('原图');
figure,imshow(r),title('迭代法');
%均匀性度量法
clc;
clearall;
I=imread('man.bmp');
I=rgb2gray(I);
I=double(I);
[m,n]=size(I);
Smin=-1;
forT=0:255
sum1=0;num1=0;
sum2=0;num2=0;
fori=1:m
forj=1:n
ifI(i,j)>=T
sum2=sum2+I(i,j);
num2=num2+1;
else
sum1=sum1+I(i,j);
num1=num1+1;
end
end
end
ave1=sum1/num1;
ave2=sum2/num2;
d1=-1;
d2=-1;
fori=1:m
forj=1:n
ifI(i,j)>=T
d=(I(i,j)-ave2)^2;
ifd2==-1
d2=d;
else
d2=d2+d;
end
else
d=(I(i,j)-ave1)^2;
ifd1==-1
d1=d;
else
【】d1=d1+d;
end
end
end
end
p1=num1/(m*n);
p2=num2/(m*n);
S=p1*d1+p2*d2;
if(Smin==-1)
Smin=S;
else
if(S
Smin=S;
Th=T;
end
end
end
Th
figure,imshow(uint8(I)),title('原图');
fori=1:m
forj=1:n
ifI(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I),title('均匀性度量法');
%类间最大距离法
clc;
clearall;
I=imread('man.bmp');
I=rgb2gray(I);
I=double(I);
[m,n]=size(I);
Smax=0;
forT=0:255
sum1=0;num1=0;
sum2=0;num2=0;
fori=1:m
forj=1:n
ifI(i,j)>=T
sum2=sum2+I(i,j);
num2=num2+1;
else
sum1=sum1+I(i,j);
num1=num1+1;
end
end
end
ave1=sum1/num1;
ave2=sum2/num2;
S=((ave2-T)*(T-ave1))/(ave2-ave1)^2;
if(S>Smax)
Smax=S;
Th=T;
end
end
Th
figure,imshow(uint8(I)),title('原图');
fori=1:m
forj=1:n
ifI(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I),title('类间最大距离法');
%Otsu法
clc;
clearall;
I=imread('man.bmp');
I=rgb2gray(I);
I=double(I);
figure,imshow(uint8(I)),title('原图');
[m,n]=size(I);
Th=Otsu(I);
Th
fori=1:m
forj=1:n
ifI(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I),title('Otsu法');
functionT=Otsu(I)
[m,n]=size(I);
I=double(I);
count=zeros(256,1);
pcount=zeros(256,1);
fori=1:m
forj=1:n
pixel=I(i,j);
count(pixel+1)=count(pixel+1)+1;
end
end
dw=0;
fori=0:255
pcount(i+1)=count(i+1)/(m*n);
dw=dw+i*pcount(i+1);
end
Th=0;
Thbest=0;
dfc=0;
dfcmax=0;
while(Th>=0&Th<=255)
dp1=0;
dw1=0;
fori=0:Th
dp1=dp1+pcount(i+1);
dw1=dw1+i*pcount(i+1);
end
ifdp1>0
dw1=dw1/dp1;
end
dp2=0;
dw2=0;
fori=Th+1:255
dp2=dp2+pcount(i+1);
dw2=dw2+i*pcount(i+1);
end
ifdp2>0
dw2=dw2/dp2;
end
dfc=dp1*(dw1-dw)^2+dp2*(dw2-dw)^2;
ifdfc>=dfcmax
dfcmax=dfc;
Thbest=Th;
end
Th=Th+1;
end
T=Thbest;
%最大熵法
clc;
clearall;
I=imread('man.bmp');
I=rgb2gray(I);
figure,imshow(I);title('原始图像');
h=imhist(I);
h1=h;
len=length(h);
[m,n]=size(I);
h1=h1/(m*n);
fori=1:(len-1)
ifh(i)~=0
P1=sum(h1(1:i));
P2=sum(h1((i+1):len));
else
continue;
end
H1(i)=-(sum(P1.*log(P1)));
H2(i)=-(sum(P2.*log(P2)));
H(i)=H1(i)+H2(i);
end
m1=max(H);
Th=find(H==m1);
Th
fori=1:m
forj=1:n
ifI(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I);title('最大熵法');
%最大类内类间方差比法
clc;
clearall;
I=imread('man.bmp');
I=rgb2gray(I);
I=double(I);
[m,n]=size(I);
Smax=-1;
forT=0:255
sum1=0;num1=0;
sum2=0;num2=0;
fori=1:m
forj=1:n
ifI(i,j)>=T
sum2=sum2+I(i,j);
num2=num2+1;
else
sum1=sum1+I(i,j);
num1=num1+1;
end
end
end
ave1=sum1/num1;
ave2=sum2/num2;
ave=(sum1+sum2)/(m*n);
d1=-1;
d2=-1;
fori=1:m
forj=1:n
ifI(i,j)>=T
d=(I(i,j)-ave2)^2;
ifd2==-1
d2=d;
else
d2=d2+d;
end
else
d=(I(i,j)-ave1)^2;
ifd1==-1
d1=d;
else
d1=d1+d;
end
end
end
end
p1=num1/(m*n);
p2=num2/(m*n);
S1=p1*(ave1-ave)^2+p2*(ave2-ave)^2;
S2=p1*d1+p2*d2;
S=S1/S2;
ifS>Smax
Smax=S;
Th=T;
end
end
Th
figure,imshow(uint8(I)),title('原图');
fori=1:m
forj=1:n
ifI(i,j)>=Th
I(i,j)=255;
else
I(i,j)=0;
end
end
end
figure,imshow(I),title('最大类内类间方差比法');
%局部阈值
clc;
clearall;
I=imread('man.bmp');
I=rgb2gray(I);
I=im2double(I);
figure,imshow(I),title('原图');
se=strel('disk',10);
ft=imtophat(I,se);
Th=graythresh(ft);
Th
G=im2bw(ft,Th);
figure,imshow(G),title('局部阈值');
