遗传算法实现图像分割(MATLAB)

 本文是对于Omar Banimelhem and Yahya Ahmed Yahya 发表论文《Multi-Thresholding Image Segmentation Using Genetic Algorithm》的翻译。

用遗传算法对图像进行多阈值分割(Multi-Thresholding Image Segmentation Using Genetic Algorithm)

摘要:图像分割是计算机视觉和图像处理中的一个重要问题。将一幅数字图像分成多个区域或者集合。在多个问题和应用中,图像分割变得越来越重要,从而使得研究者提出并改进图像分割过程的算法。图像分割的方法有很多。本文使用基于遗传算法的阈值技术来寻找不同目标和背景间的最优阈值。

关键词:图像分割;阈值;遗传算法。

1 引言

图像分割是区分一幅图像中的目标和背景的过程。对于许多依赖于计算机视觉的应用,如医学成像、卫星图像中的目标定位、机器视觉、指纹和人脸识别、农业成像等,图像分割是一个重要的预处理任务。图像分割阶段的准确性会对图像处理的后续阶段产生很大的影响。许多研究者对于图像分割问题已研究了多年;然而,由于图像具有如不同形式的直方图的特点,所以图像分割问题仍是一个公开的研究问题,需要进一步的研究。

最近,提出了许多技术,包括基于图形的算法边缘检测的算法进和基于阈值的算法[1-5],用于图像分割。由于阈值技术简单,在过去几年中得到了很大的关注。

另一方面,由于遗传算法能够为许多实际应用提供近似最优解,所以最近用其解决图像分割问题。一般来说,遗传算法(Genetic Algorithm,GA)是模拟自然选择的生物进化过程一种软计算模型(a soft computational model)[6]。

本文采用基于遗传算法的图像阈值方法,通过寻找阈值,将阈值问题转换为优化问题。该方法的目标是最大化目标间的类内方差和最小化目标像素之中的背景像素间的类间方差[7]。在每个目标中,像素的灰度间的同质性会影响图像分割质量。如果目标像素灰度变化,结果将有可能变差[7]。

论文的其余部分组织如下。第2节详细描述了图像阈值方法。第3节讨论了遗传算法的概念。第4节讨论所提方法。第5节给出仿真结果。最后,第6节对本文进行总结。

2 图像阈值

数字图像可以视为二维矩阵或者二元函数。它由称为像素的离散点组成。在彩色图像中,每个像素有三个值:红、绿和蓝。每个值在0到L-1的范围中,其中L是精度的数量。另一方面,灰度图像由像素组成,其中每个像素只有0到L-1之间的一个值,称之为灰度。对于许多图像处理问题,处理灰度图像比处理彩色图像更简单、有效。因此,在使用图像处理算法前,会将彩色图像转换为灰度图像。最常用的灰度为256(即每个像素值在0到255之间)。

图像阈值是一种用于灰度图像的图像分割方法。该想法是为了找到一个阈值,如果像素低于阈值,就认为是背景,否则认为是目标的一部分。例如,图1-a中的图像有一个目标和背景。图像分割的结果如图1-b所示。

图1:原始图像(a)和分割结果(b)

基于阈值的方法分为单级阈值(single-level thresholding)和多级阈值(multi-level thresholding)。多阈值方法通常是通过寻找将多个目标分开的多个阈值作为图像阈值。例如,图2-a中的图像有三个目标,图像分割的结果如图2-b所示。一般而言,分割一幅有n个目标和背景的图像,可以使用n个阈值。

图2:一幅图像(a)有3个目标和图像(b)是分割结果

对于找到分开目标的最优阈值,处理图像的统计比处理图像本身要简单。图像直方图统计是图像中每个灰度出现的频率的一维表示。通过计算具有相同灰度的像素的数目得到。图1中的图像的图像直方图如图3所示。

图3:图1的图像直方图

将彩色图像转换为灰度图像,然后用图像的一维直方图简化阈值分割过程,提高计算效率,可以用于许多实时应用。

从图像直方图中找到阈值的方法有很多。在单阈值的情况下,我们可以尝试0到L-1之间的所以值,然后我们选择得到最好的分割结果的值作为阈值。

然而,对于多阈值分割,尝试所有可能的组合需要L×(L-1)×…×(L-t+1)次,其中t是阈值的个数。很明显,尝试错有可能的阈值的朴素方法(naïve method)对于更大的t而言,并不是有效的,也不是可行的。事实上,朴素方法的计算复杂度表明,要使用更有效的算法来寻找t个阈值。

通过测量阈值处理的性能(goodness)来获知一组阈值是否是有效的阈值也是很重要的。在自然图像中,可以看到相同的物体(objects)像素相似,不同的物体像素不相似,这表明可以考虑使用类内方差和类间方差进行测量;即属于一个物体的像素间的灰度上的方差和不同物体中的像素间的方差。

朴素方法的计算复杂度和多阈值情况中的性能测量(goodness measure)的存在促使我们使用一个有效的搜索算法。在下一节,将描述如何用遗传算法实现它。

3 遗传算法

遗传算法(GA)是模拟基因的自然地元启发式算法。遗传算法在不存在确定性方法或者如果确定性方法计算复杂,GA是一个基于种群(population)的算法(即每次迭代产生多个解)。每次迭代的解的个数称为种群大小。每个解用染色体表示,每一条染色体都是由基因组成。

对于一个种群大小为n的遗传算法,是从n个随机解开始。然后,选择最优解杂交产生新的解。将产生的最优解添加到下一个迭代中,而舍弃不好的结果。当算法对解进行了迭代时,可以实现将这些解在一定程度上改进到收敛至接近最优解。

当使用遗传算法时,有许多因素需要考虑。第一个因素是染色体和基因的表示,因为不好的表示可能会减慢收敛过程。遗传算法的另一个重要因素是从旧的解中产生新解的机制。最常用的机制是交叉和突变。第三个因素是如何找到一个适应度函数(即评价解的方法)用于接受或者拒绝解,和如何选择最优的解来杂交。

一般来说,遗传算法有四个阶段:种群初始化、适应度评估、复制(reproduction)和终止。

初始化是创建初始随机解的过程,可以通过将基因设置为随机值来完成。在初始化的过程中,创造n条染色体作为第一代解。初始化后,用适应度函数评价每条染色体的适应度(解的性能(solution goodness))。

复制过程有四个步骤:选择、交叉、突变和接受解。

•在选择步骤中,选择当前种群中的最适解来复制得到新的解。然而,更少的适应度将有机会被选择。可用通过多种机制实现选择步骤。一种流行的方法是rollet wheel。该方法的工作原理如下:假设一个rollet wheel被分为n部分,每个部分表示一个染色体。每个部分的面积是根据染色体的适应度来映射。然后,在0到360之间随机选择一个数字,并将其映射到rollet。这给最是染色体更多的机会,而不太给不太合适的那些提供可能的机会。在实践中,这可以通过评估每一个染色体的适应度来实现,然后将总数归一化到1。这个过程将适应度测量转换为概率值,然后发现累积分布,并生成一个0到1之间的随机数。选择累积中的对应于随机数的染色体。这一选择在两条染色体上执行,每次复制两条新的染色体。

•选择染色体后,通过选择染色体上的一个随机的点和该点后交换基因来完成交叉操作,如图4所示。

图4:交叉操作

•交叉可能陷入局部最优。为了克服该问题,对随机选择的一个基因进行突变操作,改变其值,从而打破平局。基因的一个常用的表示是比特(bits),每个基因都用一个bit来表示。在这种情况下,通过在染色体上随机翻转一个bit来完成变异。

•在交叉和突变之后,复制得到了两条新的染色体。最后一步是在新的种群中接受这两条染色体。通常,如果新的染色体比它们的父代好,就会被接受。

终止是遗传算法的最后一步。通常,当满足某个条件时,遗传算法的迭代就停止了。最常用停止条件是迭代次数。当满足一个预定义的迭代次数时,遗传算法就停止了。一般的遗传算法如图5所示[6]。

图5:一般的遗传算法[6]的流程图

4 所提算法

所提算法用遗传算法来寻找近似最优阈值。我们的目标是最大化类内方差并最小化类间方差。所提算法使用编码为bits的一个向量的然染色体。每个向量有L*n bits,其中L是log(灰度的个数),n是所用阈值的个数。染色体重大每个L bits表示一个阈值。染色体结构如图6所示。

图6:染色体结构

首先,找到一个图像的直方图。直返图提供的信息将用于评价适应度。然后,随机产生K个初始染色体。之后,遗传算法迭代固定的迭代次数,最终选出最优染色体作为解。图所提算法的伪代码如图7所示。

 

所提图像分割算法

输入:图像Im、种群大小、交叉率、变异率、迭代次数、阈值个数

输出:分割图像

图7:所提图像分割算法

适应度函数测量分割的性能(goodness)。通常,在自然图像中,不同物体的像素之间的方差大,而同一物体中像素间的方差小。这影响了适应度函数的设计。在所提算法中,我们使用如下的函数F,作为阈值i的适应度[7]:

其中

物体间的方差SBetween objects如下:

其中

是用阈值thsld i所分割的像素的平均,表示如下:

分割i的概率表示如下:

图像的全局平均是:

mg = mi × pi

分割i的的方差是:

图像的全局方差是

5 实现及结果

该算法在Matlab中实现。初始种群的大小是10条染色体。交叉率和变异率分别是0.95和0.05。算法的迭代次数是500次。它运行在2.27GHz核i3的处理器上。

该算法对于分割图像有高质量结果。一组输入图像和分割结果的样本如图8-11所示。

图8:(a)摄影师的原始图像(b)一个分割后的结果(c)直方图图像和在83 levels的阈值

   

图9:(a)米粒的原始图像(b)一个分割后的结果(c)直方图图像和在130 levels的阈值

 

图10:(a)渔夫的原始图像(b)两个分割后的结果(c)直方图图像和在111与168 levels的阈值

图11:(a)升力体(liftingbody)的原始图像(b)两个分割后的结果(c)直方图图像和在119与191 levels的阈值

四幅图像的仿真时间和阈值如表1所示。

 

表1:四幅图像的仿真时间和阈值

6 总结

本文中,我们提出了使用遗传算法的多阈值图像分割的一个新方法。由表示基因的bits的向量构造染色体,每个向量用t levels建模,每个level用Log L bits表示。在未来的工作中会研发其他的适应度函数。

%**************************************************************************主函数
function main()
clear all
close all
clc
global chrom oldpop fitness lchrom  popsize cross_rate mutation_rate yuzhisum
global maxgen  m n fit gen yuzhi A B C oldpop1 popsize1 b b1 fitness1 yuzhi1
A=imread('house.jpg');     %读入要进行分割的图像
%A=imresize(A,0.4);
B=rgb2gray(A);         %灰度化
C=B;
C=imresize(B,0.1);     %将读入的图像缩小
lchrom=8;              %染色体长度
popsize=10;            %种群大小
cross_rate=0.7;        %交叉概率
mutation_rate=0.4;     %变异概率
maxgen=150;            %最大代数
[m,n]=size(C);
initpop;    %初始种群
for gen=1:maxgen
    generation;  %遗传操作
end
findresult; %图象分割结果
%输出进化各曲线
figure;
gen=1:maxgen;
plot(gen,fit(1,gen)); 
title('最佳适应度值进化曲线');

%**************************************************************************初始化种群
%初始化种群%
function initpop()
global lchrom oldpop popsize chrom C
imshow(C);
for i=1:popsize
    chrom=rand(1,lchrom);
    for j=1:lchrom
        if chrom(1,j)<0.5
            chrom(1,j)=0;
       else 
           chrom(1,j)=1;
        end
    end
    oldpop(i,1:lchrom)=chrom;    
end

%**************************************************************************产生新一代个体
%产生新一代个体%
function generation()
fitness_order;                  %计算适应度值及排序
select;                         %选择操作
crossover;                      %交叉
mutation;                       %变异

%**************************************************************************计算适度值并且排序
%计算适度值并且排序%
function fitness_order()
global lchrom oldpop fitness popsize chrom fit gen C m n  fitness1 yuzhisum
global lowsum higsum u1 u2 yuzhi gen oldpop1 popsize1 b1 b yuzhi1 
if popsize>=5
    popsize=ceil(popsize-0.03*gen);
end
if gen==75     %当进化到末期的时候调整种群规模和交叉、变异概率
    cross_rate=0.3;        %交叉概率
    mutation_rate=0.3;     %变异概率
end
%如果不是第一代则将上一代操作后的种群根据此代的种群规模装入此代种群中
if gen>1   
    t=oldpop;
    j=popsize1;
    for i=1:popsize
        if j>=1
            oldpop(i,:)=t(j,:);
        end
        j=j-1;
    end
end
%计算适度值并排序
for i=1:popsize
    lowsum=0;
    higsum=0;
    lownum=0;
    hignum=0;
    chrom=oldpop(i,:);
    c=0;
    for j=1:lchrom
        c=c+chrom(1,j)*(2^(lchrom-j));
    end
    b(1,i)=c*255/(2^lchrom-1);  %转化到灰度值        
    for x=1:m
        for y=1:n
            if C(x,y)<=b(1,i)
                lowsum=lowsum+double(C(x,y));%统计低于阈值的灰度值的总和
                lownum=lownum+1; %统计低于阈值的灰度值的像素的总个数
            else
                higsum=higsum+double(C(x,y));%统计高于阈值的灰度值的总和
                hignum=hignum+1; %统计高于阈值的灰度值的像素的总个数
            end
        end
    end
    if lownum~=0
        u1=lowsum/lownum; %u1、u2为对应于两类的平均灰度值
    else
        u1=0;
    end
    if hignum~=0
        u2=higsum/hignum;
    else
        u2=0;
    end   
    fitness(1,i)=lownum*hignum*(u1-u2)^2; %计算适度值
end
if gen==1 %如果为第一代,从小往大排序
    for i=1:popsize
        j=i+1;
        while j<=popsize
            if fitness(1,i)>fitness(1,j)
                tempf=fitness(1,i);
                tempc=oldpop(i,:);
                tempb=b(1,i);
                b(1,i)=b(1,j);
                b(1,j)=tempb;
                fitness(1,i)=fitness(1,j);
                oldpop(i,:)=oldpop(j,:);
                fitness(1,j)=tempf;
                oldpop(j,:)=tempc;
            end
            j=j+1;
        end
    end
    for i=1:popsize
        fitness1(1,i)=fitness(1,i);
        b1(1,i)=b(1,i);
        oldpop1(i,:)=oldpop(i,:);
    end
    popsize1=popsize;
else %大于一代时进行如下从小到大排序
    for i=1:popsize
        j=i+1;
        while j<=popsize
            if fitness(1,i)>fitness(1,j)
                tempf=fitness(1,i);
                tempc=oldpop(i,:);
                tempb=b(1,i);
                b(1,i)=b(1,j);
                b(1,j)=tempb;
                fitness(1,i)=fitness(1,j);
                oldpop(i,:)=oldpop(j,:);
                fitness(1,j)=tempf;
                oldpop(j,:)=tempc;
            end
            j=j+1;
        end
    end
end 
%下边对上一代群体进行排序
for i=1:popsize1
    j=i+1;
    while j<=popsize1
        if fitness1(1,i)>fitness1(1,j)
            tempf=fitness1(1,i);
            tempc=oldpop1(i,:);
            tempb=b1(1,i);
            b1(1,i)=b1(1,j);
            b1(1,j)=tempb;
            fitness1(1,i)=fitness1(1,j);
            oldpop1(i,:)=oldpop1(j,:);
            fitness1(1,j)=tempf;
            oldpop1(j,:)=tempc;
        end
        j=j+1;
    end
end
%下边统计每一代中的最佳阈值和最佳适应度值
if gen==1
    fit(1,gen)=fitness(1,popsize);
    yuzhi(1,gen)=b(1,popsize);
    yuzhisum=0;
else
    if fitness(1,popsize)>fitness1(1,popsize1)
        yuzhi(1,gen)=b(1,popsize); %每一代中的最佳阈值
        fit(1,gen)=fitness(1,popsize);%每一代中的最佳适应度
    else
        yuzhi(1,gen)=b1(1,popsize1); 
        fit(1,gen)=fitness1(1,popsize1);
    end
end

%**************************************************************************精英选择
%精英选择%
function select()
global fitness popsize oldpop temp popsize1 oldpop1 gen b b1 fitness1
%统计前一个群体中适应值比当前群体适应值大的个数
s=popsize1+1;
for j=popsize1:-1:1
    if fitness(1,popsize)<fitness1(1,j)
        s=j;
    end
end
for i=1:popsize
    temp(i,:)=oldpop(i,:);
end
if s~=popsize1+1
    if gen<50  %小于50代用上一代中用适应度值大于当前代的个体随机代替当前代中的个体
        for i=s:popsize1
            p=rand;
            j=floor(p*popsize+1);
            temp(j,:)=oldpop1(i,:);
            b(1,j)=b1(1,i);
            fitness(1,j)=fitness1(1,i);
        end
    else
        if gen<100  %P~100代用上一代中用适应度值大于当前代的个体代替当前代中的最差个体
            j=1;
            for i=s:popsize1
                temp(j,:)=oldpop1(i,:);
                b(1,j)=b1(1,i);
                fitness(1,j)=fitness1(1,i);
                j=j+1;
            end
        else %大于100代用上一代中的优秀的一半代替当前代中的最差的一半,加快寻优
            j=popsize1;
            for i=1:floor(popsize/2)
                temp(i,:)=oldpop1(j,:);
                b(1,i)=b1(1,j);
                fitness(1,i)=fitness1(1,j);
                j=j-1;
            end
        end
    end
end
%将当前代的各项数据保存
for i=1:popsize
    b1(1,i)=b(1,i); 
end    
for i=1:popsize
    fitness1(1,i)=fitness(1,i); 
end
for i=1:popsize
    oldpop1(i,:)=temp(i,:);
end
popsize1=popsize;

%**************************************************************************交叉
%交叉%
function crossover()
global temp popsize cross_rate lchrom
j=1;
for i=1:popsize
    p=rand;
    if p<cross_rate
        parent(j,:)=temp(i,:);
        a(1,j)=i;
        j=j+1;
    end
end
j=j-1;
if rem(j,2)~=0
    j=j-1;
end
if j>=2
    for k=1:2:j
        cutpoint=round(rand*(lchrom-1));
        f=k;
        for i=1:cutpoint
            temp(a(1,f),i)=parent(f,i);
            temp(a(1,f+1),i)=parent(f+1,i);
        end
        for i=(cutpoint+1):lchrom
            temp(a(1,f),i)=parent(f+1,i);
            temp(a(1,f+1),i)=parent(f,i);
        end
    end
end

%**************************************************************************变异
%变异%
function mutation()
global popsize lchrom mutation_rate temp newpop oldpop
sum=lchrom*popsize;    %总基因个数
mutnum=round(mutation_rate*sum);    %发生变异的基因数目
for i=1:mutnum
    s=rem((round(rand*(sum-1))),lchrom)+1; %确定所在基因的位数
    t=ceil((round(rand*(sum-1)))/lchrom); %确定变异的是哪个基因
    if t<1
        t=1;
    end
    if t>popsize
        t=popsize;
    end
    if s>lchrom
        s=lchrom;
    end
    if temp(t,s)==1
        temp(t,s)=0;
    else
        temp(t,s)=1;
    end
end
for i=1:popsize
    oldpop(i,:)=temp(i,:);
end

%**************************************************************************查看结果
%查看结果%
function findresult()
global maxgen yuzhi m n C B A 
result=floor(yuzhi(1,maxgen)) %result为最佳阈值
C=B;
%C=imresize(B,0.3);
imshow(A);
title('原始裂缝图像')
figure;
subplot(1,2,1)
imshow(C);
title('原始裂缝灰度图')
[m,n]=size(C);
%用所找到的阈值分割图象
for i=1:m
    for j=1:n
        if C(i,j)<=result
            C(i,j)=0;
        else
            C(i,j)=255;
        end
    end
end
subplot(1,2,2)
imshow(C);
title('阈值分割后的裂缝图');

 

已标记关键词 清除标记
©️2020 CSDN 皮肤主题: 技术黑板 设计师:CSDN官方博客 返回首页