基于改进蚁群算法的医学图像边界提取

陈 超

(内江师范学院四川省数据恢复重点实验室 四川 内江 641112)

0 引 言

医学图像边界信息检测是医学图像分析的重要环节。一般的医学图像边缘是通过寻找最佳阈值来分割图像，而边界是一个重要的信息点。因为边界往往是正常组织和病变组织相邻的位置，实际检测病变组织的边界不准确或者不精确，对于诊断病情没有实际辅助意义。一些国内学者[1-2]使用权重因子改进概率转移函数与混沌系统动态更新信息素矩阵检测出目标的真实边缘，结果表明识别边缘效果好但没有做其他目标对比实现，无法说明其泛化能力。詹永丰等[3]用边缘梯度与二维高斯滤波器对CT图像进行滤波处理实现CT 图像边界的完整分割，但分割出来的目标完整性不理想。张小莉等[4]改进蚁群算法(ACO)，结合聚类动态更新信息素浓度和参数变量把实验医学图像精确分割，但是没有针对噪声影响下分割是否有效进行分析。许梅梅等[5]提出自适应边缘算法优化医学图像边缘信息检测出目标的真实边缘，但有时丢失了部分真实的边缘。郑延斌等[6]基于蚁群算法的虚拟植物优化研究，针对植物的分割取得优秀效果。朱立军等[7]提出一种改进蚁群算法的医学图像提取。张耀楠等[8]利用ACO结合网格分割目标边缘。吴学礼等[9]应用改进ACO实现无人机避险对比性实验，取得了良好效果，但没有对时间复杂度做定量分析。黄志华[10]利用自适应ACO求解旅行商最优路径问题(TSP)。文献[11-13]也提出从改进蚁群算法的参数实现越野路径规划或者引入混沌系统预测控制的无人机群避碰对机动性的影响或者多阈值蚁群算法实现各种图像目标的分割，取得了显著效果，但是实时性或者分割的有效性没有得到明显的提高。本文改进蚁群算法的信息素重要程度参数和挥发系数等，实时地根据聚类前后期不同，分阶段地动态修改相关参数，再结合原图和罗伯特边缘最近的变化速度动态地调整蚁群算法的信息素，指导更多的蚂蚁向最优路径聚集。通过五组医学图像的实验结果做图像分割客观分析后[14-17]，结果表明改进后的蚁群算法结合图像的边缘变化等信息可以很好地实现医学图像的有效分割。

1 算法设计

1.1 基本蚁群算法

ACO主要是模拟自然界蚂蚁的群体智能，其中每个蚂蚁的信息在蚁群内部不断更新并互相学习交流,使整体不断优化，算法通过外激素的浓度对最优解进行筛选。问题的解在蚁群中通过不断迭代逐步被优化。对于各种复杂优化问题,蚁群算法通过定义解的分量来解决。开始时,解的分量对应每只蚂蚁,备选解以重复增加解的分量的方式来产生。当处于选择点上时,决定哪个解的分量要加到它的当前部分解，主要取决于信息量的浓度以及禁忌表的情况，在构造完一个完整的解后,解的分量会被蚂蚁涂上信息素,指导接下来整体和个体的活动。ACO中主要包含两个环节：每个蚂蚁的状态转移规则；信息素更新规则。

1) 状态转移规则:

人工蚂蚁每次从当前位置i，在八邻域范围内按照如下转移概率公式选择下一步要到达的位置j：

(1)

式中：Ω表示第i只蚂蚁允许选择可能位置的集合，在图像中一般是八邻域范围内；τij(t)是第t次迭代，医学图像某个像素点具体的位置的信息素浓度；α和β为常数；ηij(t)为期望程度。

2) 信息素更新规则：

(2)

(3)

式中：ρ为每只蚂蚁信息素的保留因子,取(0,1)范围内的常数；Q为常数；Ck(t)是第k只蚂蚁在t次迭代寻优过程中所经过路径的实际长度。

1.2 改进蚁群算法

1.2.1改进状态转移规则

在迭代中，每只人工蚂蚁都会在当前像素点位置的八邻域范围内以当前的转移概率选择最优坐标点。针对不同迭代次数NC,2个参数会对应取不同的值,以应对不同时期算法不同的特点,大致分为四个阶段进行调整，两个参数α、β动态调整具体如下：

(4)

(5)

式中:MAX为最大迭代次数；NC为当前迭代次数。这些参数的设置都是经过大量实验给定的。α′和β′分别为本文改进后的重要程度和启发因子。则改进后的转移概率公式为：

(6)

1.2.2改进信息素更新规则

医学图像的每只蚂蚁信息素值初始化为：

(7)

其中：

gk=|f(i+1,j)-f(i,j)|+|f(i，j+1)-f(i,j)|

(8)

E(f)为当前待检测医学图像的所有像素的期望值。信息素浓度的调整因子ρ改进为ρ′，具体如下：

(1-gk)

(9)

式中：学习因子cc1=0.85；rr1=0.75。

改进后的信息素更新公式为：

(10)

式中：m是所有像素的综述；n是在聚类时小于指定半径的像素总个数，则：

(11)

本文提出动态自适应的方法来控制：

(12)

式中：ε=0.001。算法终止，反之继续运行寻找最优阈值。最后使用改进ACO边缘与聚类算法将样本到聚类中心的欧式加权距离:

(13)

式中：cjk为聚类半径；m为蚂蚁数。

2 基于改进ACO医学图像边界提取步骤

(1) 初始化蚁群算法各参数，如α、β、Q等参数；

(2) 求出图像直方图的峰值点作为初始聚类中心；

(3) 根据信息素和梯度信息由式(6)计算蚂蚁选择最优路径；

(4) 根据式(13)算欧氏加权距离；

(5) 使用式(10)更新路径上的信息素浓度；

(6) 计算该类更新后的聚类中心；

(7) 当聚类中心变化量小于阈值时，停止计算，并输出聚类中心及聚类个数；

(8) 由第(7)步得到的聚类数目和中心作为FCM的原始聚类个数和中心位置；

(9) 根据最新的信息素和欧式距离更新聚类中心；

(10) 计算最优分割阈值或达到最大迭代次数时，停止计算，否则转到(3)继续迭代。

3 实 验

3.1 实验平台

实验环境是Window XP系统、MATLAB 7.0、英特尔奔腾CPU 2.60.00 GHz(2 CPUs)、内存8 GB，64位操作系统。

3.2 核心代码及结果分析

针对三大经典边缘算子以及最近两年国内外成果(文献[1]和文献[12])做了5组对比性实验,所有实验数据都是256×256的医学图像，相关的核心代码如下：

function T=AntCC(x)

%初始化相关参数

alpha=1;

beta=1;

maxant=numel(x);

r=20;

//聚类半径

rho=0.9;

//挥发系数

ranta=0.9;

//隶属度

c=105;

//食物源中心(初始聚类中心)

new_m=0;

//类集合中所有像素灰度值总和

new_n=0;

//类集合中所有像素个数

x=double(x);

j=1;

for i=1:maxant

distance(i)=sqrt((x(i)-c)^2);

//蚂蚁i到食物源c的距离

if distance(i)<=r;

distance3(j)=distance(i);

//小于指定欧式距离的保留，其余的赋值为0

j=j+1;

end

if distance(i)<=r;

ph(i)=1;

else

ph(i)=0;

end

if distance(i)==0;

//判断为一类

similar(i)=1;

else

similar(i)=r/sqrt((x(i)-c)^2);

//不为一类，那么就对相似度赋值为r/sqrt((x(i)-c)^2)，

//越不相似的，相似度越低

end

if distance(i)<=r;

//判为一类

new_n=new_n+1;

newph(new_n)=ph(i)*(1-rho^new_n)/(1-rho);

%rho=0.9(改进)

newph1(new_n)=1;

//信息浓度,更新为1

end

if distance(i)<=r;

similar1(i)=similar(i);

end

end

similar2=similar1(find(similar1));

allsum=newph1.*similar2;

for n=2:new_n

allsum(n)=allsum(n-1)+allsum(n);

end

for n=1:new_n

newallsum(n)=allsum(n)/n;

end

cc1=1.95;

rr1=0.75;

for n=1:new_n

if distance3(n)==0;

newp(n)=1;

//是同一类，必然要选

else

//不是同一类，就根据当前概率来选择

newp(n)=((newph1(n)^alpha)*

(similar2(n)^beta))/newallsum(n);

V= cc1*rr1*(1-arccos(-2*n/new_n+1));

V=0.729*V;

newp(n)=newp(n)+V;

end

end

for i=1:maxant

if distance(i)<=r;

//差异太大的设置为0

new_m(i)=x(i);

end

end

new_m1=new_m(find(new_m));

//只保留非0

//元素，两个蚁群算法都是使用这个来计算new_m1

//改进蚁群算法求类间分割阈值

new_m2=0;

new_n2=0;

for n=1:new_n

//目标类

if newp(n)>ranta;

//大于隶属度0.9,对该类像素值求和(改进)

new_m2=new_m2+new_m1(n);

//类集合中所有像素灰度值总和

new_n2=new_n2+1;

//统计一共有多少个像素

end

end

//改进蚁群算法求类间分割阈值

//计算目标所有像素的灰度值的平均值

T=new_m2/new_n2;

在本次实验中一共五张医学图像，如图1-图5所示。其中每张图中包含原图、对应的直方图、原图的灰度化图像、roberts算子、sobel算子、prewitt算子、文献[1]、文献[12]、本文改进算法分割的效果图。图1为CT医学图像，从分割效果可以看出改进后的蚁群算法可以很好地检测出图像的边缘信息，包含的信息量较文献[1]和文献[12]增加很多，对于诊断有很大的辅助作用。

图1 第一张医学图像边界提取效果图

图2也是CT图像，只是清晰度不同，左边图像灰度值较低，右边较高。使用改进后的蚁群算法分割后，可以看出边界明显较前面两种算法更加完整，尤其是第二排第二幅图像拥有的信息量提高了。

图2 第二张医学图像边界提取效果图

图3是病人的头部MRI图像，对于后脑勺部分的纹理比较复杂，灰度信息差别不大，使用改进后的算法分割后可以从分割的边界信息辅助诊断，靠近嘴巴和后脑勺部分的边界得到了清晰的呈现。

图3 第三张医学图像边界提取效果图

图4是胸部到盆腔的整段X影像。在胸部靠近脊椎部分受噪声污染导致图像边界的清晰度受到严重影响，使用改进后的蚁群算法能检测出真实边界。

图4 第四张医学图像边界提取效果图

图5是一个病患的TCD颅脑超声脑电图,脑血管的边界信息可以很好地辅助医生确诊疾病，检测出的边界清晰度和完整性得到了提高，尤其是脑颅内部的血管边界得到了清晰的展示。

图5 第五张医学图像边界提取效果图

3.3 结果评价

分别对改进ACO结合图像边缘信息分割后做客观评价,主要的评价指标如下[14-18]：

(1) 概率随机指数(PRI)计算了一对像素的分数，这些像素的标签在计算出的分割和真值之间是一致的，在多个真值分割中求平均值，以解释人类感知的尺度变化。PRI越大检测效果越好。

(2) 信息变化(VOI)度量将两个分段之间的距离定义为一个分段给定另一个分段的平均条件熵，从而粗略地测量一个分段中不能由另一个解释的随机性数量。在原图上能够检测出真实而精准的边界效果为最佳，VOI越大，检测边缘信息越好。

(3) 全局一致性错误(GCE)衡量一个分段被视为另一个分段的细化程度。 GCE越大，检测边缘信息越好。

表1-表5给出了文献[1]、文献[12]、本文算法的五幅医学图像的定量分析。

表1 第一幅医学边缘提取参数比较统计

表2 第二幅医学边缘提取参数比较统计

表3 第三幅医学边缘提取参数比较统计

表4 第四幅医学边缘提取参数比较统计

表5 第五幅医学边缘提取参数比较统计

改进ACO提取五幅医学图像边缘阈值如表6所示。

表6 改进ACO提取五幅医学图像边缘阈值

计算机算法复杂性的衡量标准是该算法完成相应功能时所需要消耗的计算机资源的多少，主要包括时间和空间两个部分，即时间复杂性和空间复杂性两个部分构成[18-19]。

(1) 时间复杂度比较。蚁群算法的时间复杂度分析如表7所示。

表7 蚁群算法的时间复杂度分析

综上所述时间复杂度为：

T(n)=O(Ncn2m2)

(14)

式中：Nc为最大循环数；n为蚂蚁数；m为图像高度和宽度。

改进后的算法在检测目标时，没有增加额外时间开销，且有一定时间的减少，在普通的计算机上即可完成实验。在上述数据集检测时没有出现数量级的增加，对于有些医学图像分割时间有略微的增加。分割运动目标的时间、整个过程所有的实验结果都直接写入磁盘需要的时间等，对比文献[1]和文献[12]算法的时间见表8。

表8 改进ACO_FCM提取五张图像边界所需时间 s

(2) 空间复杂度比较。改进后蚁群算法中的数据主要在于问题的描述和寻找最优分割阈值数据，没有增加额外的内存空间，依然使用二维数组，空间复杂度为[18-19]:

S(n)=O(n2)+O(nm)

(15)

(3) 鲁棒性分析。在以上仿真实验中，所有数据在实现目标检测的过程中没有出现不稳定现象，也没有出现过拟合现象，只有在检测后期目标位置不精准，导致检测出其他的目标，因此改进ACO_FCM算法有较好的鲁棒性。

4 结 语

本文提出了一种新的基于蚁群算法的医学图像提取算法。该算法不是通过人工预设阈值来判断医学图像，而是先通过人工蚂蚁的遍历来增强目标(医学图像边缘)的信息强度，然后使用聚类算法自动获得最佳医学图像分割阈值。通过定性和定量分析，不管是检测出的边缘展现还是细节分析都验证了本文改进后的蚁群算法具有良好的检测边界的作用。通过实验验证，当人工蚂蚁间距取 2个像素的时候，既能保证医学图像检测的速度又能保证医学图像检测的效果。因本文算法搜索的目标是灰度小、梯度大的像素，易造成漏判和误判，未来还需进一步研究和解决。