【图像分割】基于matlab粒子群算法优化模拟退火算法图像分割【含Matlab源码 2020期】

原创

佛怒火莲 2022-08-07 00:49:03 博主文章分类：matlab完整代码 ©著作权

文章标签 算法 matlab 模拟退火算法方差粒子群算法 文章分类 数据结构与算法人工智能

©著作权归作者所有：来自51CTO博客作者佛怒火莲的原创作品，请联系作者获取转载授权，否则将追究法律责任

一、粒子群算法优化模拟退火算法图像分割简介（具体理论见参考文献）

1基于模拟退火思想的粒子群算法
1.1 基本PSO算法
首先, 粒子群算法是由Eberhan博士和Kennedy博士最先提出的全局优化进化算法。该算法源于对鸟群捕食行为的灵感, 其基本思想是通过群体中个体之间的协作和共享来寻求最优解。

粒子群算法中, 每个优化问题的解看作是搜索空间中的一只鸟, 称之为“粒子"。鸟被抽象为没有质量和体积的微粒 (点) , 并将其延伸到N维空间, 粒子i在N维空间里的位置表示为矢量Xi= (xi1, xi2, …, xiN) , 飞行速度表示为矢量Vi= (vi1, vi2, …, viN) , 每个粒子都有一个由目标函数决定的适应值 (FitnessValue) , 并且知道自己到目前为止发现的最好位置pbest和现在的位置Xi, 这个可以看作是粒子自己的飞行经验。除此之外, 每个粒子还知道到目前为止整个群体中所有粒子发现的最好位置gbest (gbest是pbest中的最好值) 。这个可以看作是粒子同伴的经验。粒子就是通过自己的经验和同伴中最好的经验来决定下一步的运动。

对于第k次迭代, PSO中的每一个粒子是按照式 (1) 和式 (2) 进行变化的。

【图像分割】基于matlab粒子群算法优化模拟退火算法图像分割【含Matlab源码 2020期】_模拟退火算法

粒子通过不断的学习更新, 最终飞向解空间中最优解所在位置, 搜索过程结束。最后输出的g是全局最优解。在更新过程中, 粒子每一维的速率被限制到vmin-vmax之内。

1.2 模拟退火思想的粒子群算法 (SAPSO)
基本粒子群优化算法中, 虽然粒子速度作了限制, 不会变化太大, 但位置更新时未作限制, 就有可能新的位置会变得很坏, 引起收敛速度缓慢, 所以要对更新的位置作限制。限制的思路有两种方法:一种类似于速度的限制方法, 给每一维变量限制一个范围;另一种思路采用模拟退火算法思想, 模拟退火算法用于优化问题的出发点是基于物理中固体物质的退火过程与一般优化算法的相似性。算法的基本思想是从一给定解开始的, 从邻域中随机产生另一个解, 接受准则允许目标函数在有限范围内变坏, 它由一控制参数t决定, 其作用类似于物理过程中的温度T, 对于控制参数t的每一取值, 算法持续进行“产生新解—判断—接受或舍弃”的迭代过程, 对应着固体在某一恒定温度下趋于热平衡的过程。经过大量的解变换后, 可以求得给定控制参数t值时优化问题的相对最优解。然后减小控制参数t的值, 重复执行上述迭代过程。当控制参数逐渐减少并趋于零时, 系统也越来越趋于平衡状态, 最后系统状态对应于优化问题的整体最优解, 该过程也称冷却过程[8]。

方法1: 在基本粒子群算法的基础上, 在更新当前位置时, 把它限制在xmin-xmax内。

采用模拟退火算法思想有3种方法改进, (在此以求目标函数的极小值为例, 即适应值越小越好) 分析如下:

方法2:按式 (2) , 计算新的位置, 然后计算两个位置所引起的适应值的变化量ΔE;若ΔE≤0, 接受新值, 否则若exp (-ΔE/T) >rand (0, 1) 也接受新值。其中rand (0, 1) 表示0-1之间的随机数, T是起始温度, α为退火系数。否则就拒绝。

方法3: 接受准则允许目标函数在有限范围内变坏, 并不直接按概率取舍, 而是按ΔE<e, e为按允许目标函数变坏的范围。

方法SAPSO:把方法1与方法3结合在一起, 现提出一种基于模拟退火思想微粒群算法 (SAPSO) , 既限定位置的范围, 又设定允许目标函数变坏的范围e, 接受准则允许目标函数在有限范围内变坏, 即若ΔE<e, 就接受新值, 否则就拒绝。把SAPSO方法结合最大类间方差法应用于图像分割, 这充分发挥了模拟退火、PSO、最大类间方差法三种算法的优势, 非线性快速且稳定地求解最大类间方差及对应的灰度阈值, 在取得较好的图像分割效果的同时, 显著地提高了计算速度。

2基于模拟退火思想的PSO算法用于图像分割
2.1 最大类间方差分割方法
最大类间方差阈值分割法的基本思路是将直方图在某一阈值处分割成两组, 计算两组的方差信息, 因为方差是灰度分布均匀性的一种度量, 方差值越大, 说明构成图像的两部分差别越大, 当部分目标错分为背景或部分背景错分为目标都会导致两部分差别变小, 因此使类间方差最大的分割意味着错分概率最小, 此时的阈值即为最佳阈值。

最大类间方差法以式 (3) :

【图像分割】基于matlab粒子群算法优化模拟退火算法图像分割【含Matlab源码 2020期】_算法_02

2.2采用SAPSO算法优化最大类间方差分割法的具体实现

从上面讨论可以看出, 最大类间方差法的计算量很大, 严重影响了阈值的选取效率, 难以满足图像分割的实时处理要求和大幅图像处理的需要。由于最大类间方差法实际上也是一种寻优方法, 所以文中使用前文提出的SAPSO算法优化最大类间方差法的寻优过程。

2.2.1 SAPSO算法流程
下面是SAPSO算法的基本步骤:

初始化:设定微粒群规模M, 每个微粒代表 (0～L-1) 区间内的某个灰度值, 其中L为需要处理图像的灰度级。随机对微粒群各微粒的初始位置和速度进行初始设定, 预先设定算法的运行参数。
根据式 (1) 和式 (2) , 由当前位置和速度产生各个微粒的新的位置。并且计算微粒群中每个微粒的当前最好位置pbest和全局最好位置gbest。

While (不满足迭代终止条件) do

选择式 (3) 的F (k) 作为微粒群算法的适应函数, 计算每个粒子新位置的适应值。
对各个微粒, 将其适应值与所经历过的最好位置pbest的适应值进行比较, 若较好, 则将其作为当前的最好位置;根据各个微粒的个体极值pbest找到全局最好位置gbest。对适应值进行比较, 若较好, 则将其作为当前的全局最好位置。进行一步迭代。
根据式 (1) 进化每个微粒的Vi, 并把它限制在vmax内。
根据式 (2) 进化每个微粒的Xi, 并把它限制在0～L-1内。
计算两个位置所引起的适应值的变化量ΔE;若ΔE>-e, e为允许目标函数变坏范围, 接受新值;否则就拒绝, xk+1仍为xk。

End

将gbest对应的阈值K作为分割阈值门限, 进行图像分割。利用灰度图像直方图可以实现灰度图像的二值化, 一般而言, 如果一幅图像用f (i, j) 表示, 图像的二值化阈值处理方式如下:

f (i, j) =1;f (i, j) ≥K
f (i, j) =0;f (i, j) <K

通常, 用f (i, j) =1的部分表示图像, 用f (i, j) =0的部分表示背景。

2.2.2 SAPSO算法的主要参数的选择
虽然粒子群算法发展迅速并取得了可观的研究成果, 但其理论基础仍相对薄弱, 尤其是算法基本模型中的参数设置和优化问题还缺乏成熟的理论论证和研究。目前比较通用的参数设置方法是经验法和实验法。

粒子群规模的选择, A EI-Gallad 在其研究论文中利用统计方法分析了微粒群规模、速度极限和迭代次数对算法性能的影响, 但是这三个参数与具体问题的特性密切相关, 文中算法经过实验得到粒子数为5比较合适。

惯性权重w、加速常数c1、c2和最大速度vmax共同维护微粒对全局和局部搜索能力的平衡。这四个算法参数的设置与具体问题密切相关, 目前比较常用的方法是针对具体求解问题, 利用充分的实验来确定。

加速常数c1、c2, 低的值允许微粒在被拉回之前可以在目标区域外徘徊, 而高的值则导致微粒突然地冲向或越过目标区域, 一般都取2, 不过在文献中也有其它的取值, 但一般c1=c2并且范围在0～4之间[11], J Kennedy和R Eberhart曾在论文中建议设置c1=c2=2以保证式 (1) 中的随机乘积均值等于1, 文中算法中统一取c1=c2=2。

惯性权重因子w使微粒保持运动惯性, 使微粒有扩展搜索空间的趋势, 有能力探索新的区域。在此取值与经典粒子群算法一致, w=1。

速度限制最大值vmax决定当前位置与最好位置之间的区域的分辨率 (或精度) 。如果vmax太高, 微粒可能会飞过好解;如果vmax太小, 微粒不能在局部好区间之外进行足够的探索, 导致陷入局部优值。

该限制有3个目的:

(1) 防止计算溢出;

(2) 实现人工学习和态度转变;

(3) 决定问题空间搜索的粒度。在此取vmax=4。

粒子长度, 由优化问题决定, 就是问题解的长度, 粒子的范围, 文中粒子位置代表阈值, 对于灰度级为L的图像, 阈值范围在0～ (L-1) , 故位置限制最小值xmin=0, 位置限制最大值xmax=L-1。

模拟退火算法所用到的参数, 初始温度T和退温方式对算法的优化性能具有一定影响, 在此采用如下经验公式[16]:T=-fpg/ln (0.2) , Tk+1=α*Tk, 其中, fpg为初始种群中最佳微粒的目标值, α为退火系数。

二、部分源代码

%% Cleaning the Stage
 clc;
 clear;
 close all;
 warning(‘off’);%% Reading Image
 MainOrg=imread(‘tst.jpg’);
 Gray=rgb2gray(MainOrg);
 InpMat= double(MainOrg);%% Basics
 [s1,s2,s3]=size(InpMat);
 R = InpMat(:,:,1);
 G = InpMat(:,:,2);
 B = InpMat(:,:,3);
 X1 = (R-min(R(😃))/(max(R(😃)-min(R(😃));
 X2 = (G-min(G(😃))/(max(G(😃)-min(G(😃));
 X3 = (B-min(B(😃))/(max(B(😃)-min(B(😃));
 X = [X1(😃 X2(😃 X3(😃];%% Cluster Numbers
 clusteres = 7;%% Cost Function and Parameters
 % Cost Function
 CostFunction=@(m) CLuCosPSOSA(m, X, clusteres);
 % Decision Variables
 VarSize=[clusteres size(X,2)];
 % Number of Decision Variables
 nVar=prod(VarSize);
 % Lower Bound of Variables
 VarMin= repmat(min(X),1,clusteres);
 % Upper Bound of Variables
 VarMax= repmat(max(X),1,clusteres);%% PSO-SA Clustering Option and Run
 % PSO-SA Options
 % Iterations (more value means: slower runtime but, better result)
 Itr=50;
 % SA solver + PSO body
 SA_opts = optimoptions(‘simulannealbnd’,‘display’,‘iter’,‘MaxTime’,Itr,‘PlotFcn’,@pswplotbestf);
 options.SwarmSize = 250;
 % PSO-SA Run
 disp(['SA-PSO Segmentation Is Started … ']);
 [centers, Error] = particleswarm(CostFunction, nVar,VarMin,VarMax,SA_opts);
 disp(['SA-PSO Segmentation Is Ended. ']);%% Calculate Distance Matrix
 % Create the Cluster Center
 g=reshape(centers,3,clusteres)';
 % Create a Distance Matrix
 d = pdist2(X, g);
 % Assign Clusters and Find Closest Distances
 [dmin, ind] = min(d, [], 2);
 % Sum of Cluster Distance
 WCD = sum(dmin);
 % Fitness Function of Centers Sum
 z=WCD;
 % Final Segmented Image
 SA_Segmented=reshape(ind,s1,s2);
 PSOSAuint=uint8(SA_Segmented);
 ColorSeg = labeloverlay(Gray,PSOSAuint);
 %
 medgray = medfilt2(SA_Segmented,[5 5]);
 %
 redChannel = ColorSeg(:,:,1); % Red channel
 greenChannel = ColorSeg(:,:,2); % Green channel
 blueChannel = ColorSeg(:,:,3); % Blue channel
 medcolor1 = medfilt2(redChannel,[4 6]);
 medcolor2 = medfilt2(greenChannel,[4 6]);
 medcolor3 = medfilt2(blueChannel,[4 6]);
 medrgb = cat(3, medcolor1, medcolor2, medcolor3);%% Plot PSO-SA Segmented Result
 disp(['Error Is: ’ num2str(Error)]);
 figure(‘units’,‘normalized’,‘outerposition’,[0 0 1 1])
 subplot(2,3,1)
 subimage(MainOrg);title(‘原始’);
 subplot(2,3,2)
 subimage(Gray);title(‘灰色’);
 subplot(2,3,3)
 imshow(SA_Segmented,[]);
 title(['粒子群算法优化模拟退火算法灰色分段, 集群 = ’ num2str(clusteres)]);
 subplot(2,3,4)
 imshow(ColorSeg,[]);
 title(['粒子群算法优化模拟退火算法彩色分段, 集群 = ’ num2str(clusteres)]);
 subplot(2,3,5)
 imshow(medgray,[]);
 title(['粒子群算法优化模拟退火算法灰色中值滤波 ']);
 subplot(2,3,6)
 imshow(medrgb,[]);
 title([‘粒子群算法优化模拟退火算法彩色中值滤波’]);
 % That’s it, GoodBye 😐