CN105405136A - 一种基于粒子群算法的自适应脊柱ct图像分割方法 - Google Patents

Abstract

本发明涉及一种基于粒子群算法的自适应脊柱CT图像分割方法，该方法包括如下步骤：(1)输入原始图像；(2)初始化第一代种群；(3)计算个体最优和全局最优；(4)产生新的个体；(5)计算新的个体最优和全局最优；(6)判断是否满足最大迭代次数，如果满足则转到步骤(7)，否则，返回步骤(4)；(7)进行初步的图像分割，作为初始的分割结果；(8)采取拓扑学操作，对脊柱进一步准确的分割，输出图像。本发明优化了搜索粒度，由粗到细控制算法收敛的速度，大大加快了算法的收敛速度，提高了分割精度，并利用待分割区域的先验知识，采用拓扑学对脊柱进一步进行准确的分割，大大提升了图像分割的效率和质量。

Description

一种基于粒子群算法的自适应脊柱CT图像分割方法

技术领域

本发明涉及一般的图像处理，具体涉及图像分析，尤其是涉及通过图像属性的分析来分割图像。

背景技术

随着计算机科学技术的不断发展，图像处理技术也得到了一定的发展，逐步形成了自己的科学理论体系。图像分割是图像处理的重要技术，在实际应用和理论研究中受到越来越多人的更加关注。

图像分割技术是把图像划分为若干个有意义的区域，把感兴趣的目标提取出来，以便后续工作的有效进行，是图像处理到图像分析的关键步骤。现有的图像分割方法可分以下几类:基于区域的分割方法，基于边缘的分割方法，基于聚类的方法。在基于聚类的方法中，是将聚类的理论应道到图像分割上，将图像分割问题，转化为一个求解目标函数的最优值问题。现有的聚类算法解决图像分割问题主要包括，蚁群算法，神经网络算法，遗传算法和粒子群算法，其中粒子群优化算法的目标是解决现实世界中的各种最优化问题，该特点正好和图像分割思想符合，因此粒子群优化算法广泛被学者应用到图像分割中去。但是粒子群算法虽然实现简单，参数少，但是容易陷入局部最优，参数的选取的合适与否直接影响到图像分割的速度和精度。

章慧等人在“ImageThresholdSegmentationMethodBasedonImprovedParticleSwarmOptimization”(ComputerScience39(2012)289-301)中提出一种基于混沌粒子群算法的图像阈值分割方法，该方法引入混沌思想，对全局最优值进行混动扰动，因此增强了算法的全局搜索能力，能够在一定程度上避免陷入局部最优的状态，但是算法的分割精度和时间复杂度度提高了。

发明内容

本发明所要解决的技术问题的提供一种基于粒子群算法的自适应脊柱CT图像分割方法，该方法可改善算法的自适应性，提高求解速度和精度，加速算法收敛。

本发明解决上述问题的技术方案如下：

一种基于粒子群算法的自适应脊柱CT图像分割方法，该方法包括如下步骤：

(1)输入原始图像，计算出图像的最小灰度值和最大灰度值；

(2)初始化第一代种群：在图像的最小灰度值至最大灰度值之间分别随机产生一个整数来初始化第一代种群任一个体的每一维，再以同样的方法初始化第一代种群下一个体的每一维，依次完成第一代种群每一个体的初始化；

(3)按以下方法计算个体最优值和全局最优值：

3a)利用OSTU法得到图像的类间距方差，将类间距方差作为适应度函数；

3b)将种群个体带入到适应度函数计算出种群个体的适应度函数值，然后，先挑选种群个体中适应度最大的作为个体得到个体最优，再从个体最优中挑选适应度最大的个体作为全局最优；

(4)按以下方法产生新的个体：

4a)先按下述公式公式(Ⅰ)计算出k是自适应权重：

k＝w·(1-(t/n)^1/r)(Ⅰ)，

式(Ⅰ)中，t表示当前的迭代次数，n表示最大迭代次数，r是大于等于1的常数，w为惯性因子；

4b)再将自适应权重按下述公式(Ⅱ)产生新一代的个体：

x_id(t+1)＝k·x_id(t)+c₁·rand₁()·(p_id-x_id(t))+c₂·rand₂()·(p_gd-x_id(t))(Ⅱ)，

式(Ⅱ)中，c₁和c₂为基本粒子群算法中的加速常数，rand₁()和rand₂()表示两个在[0，1]范围里变化的随机值，p_id是个体最优值，p_gd是全局最优值，x_id(t)为第t代的第i个个体的第d维，x_id(t+1)为第t+1代的第i个个体的第d维；

(5)按以下方法产生新的个体最优和全局最优：

将新个体带入适应度函数，计算新的种群新个体的适应度函数值；当新个体的适应度值比个体最优大，那么则用新个体代替个体最优，否则个体最优不变；当新个体的适应度值比全局最优的适应度函数值大则用个体最优代替全局最优，否则全局最优不变；

(6)判断是否满足最大迭代次数，如果满足，得到全局最优值，否则返回步骤(4)；

(7)以全局最优的每一维数据作为阈值，对得到的图像矩阵进行分割；

(8)采取拓扑学操作，对脊柱进一步准确的分割，输出图像。

上述方案中，步骤3a)所述的OSTU法是：设图像被{t_{_1},t_{_2},…,t_{_M-1}}个阈值分割为M个部分，则最佳分割阈值{t^* _{_1},t^* _{_2},…,t^* _{_M-1}}需要满足如下条件：

{t^* _{_1},t^* _{_2},…,t^* _{_M-1}}＝argmax{σ²(t₁,t₂,…,t_M-1)}(Ⅲ)

其中，σ²(t₁,t₂,…,t_M-1)为图像灰度值的方差值， C_k为分割图像，k＝1,2,...,M，PR_i＝h(i)/N为各阈值的计算概率，为图像的总像素数，其中，L为所述图像的灰度值，h(i)为灰度值为i的像素数；

分割图像C₁对应阈值0,…,t_{_1}，分割图像C₂对应阈值t_{_1}+1,…,t_{_2}，...，分割图像C_M对应阈值t_{_M-1},…,L，0≤t_{_1}≤t_{_2}≤…≤t_{_M-1}≤L。

本发明与现有的技术相比具有以下优点：

本发明基于以下思路：假设空间上有一个起始点和一个目标点，在算法的刚开始阶段，粒子距离目标较远，可以加大搜索幅度，大步接近目标点，随着迭代次数的增加，粒子距离目标越来越近时，粒子搜索幅度变小，搜索的幅度更加精细，局部搜索能力变强，小幅度的逼近目标，最终到达目标点。也就是说，在给定相同的迭代次数的情况下，本发明和其他粒子群算法的区别在于，对于每一次迭代，本发明给分配给每次迭代的搜索幅度是不同的，本发明搜索幅度随着迭代次数的增加自适应的减少，而其他粒子群算法，每次迭代的搜索幅度都是不变的。

粒子群算法的速度的更新幅度会直接影响到算法的搜索性能，而公式(Ⅰ)中的自适应权重k关系到搜索的幅度大小，因此本发明通过改进权重k，通过控制权重的大小来控制搜索幅度，进而将当前迭代次数和总迭代次数考虑进去，自适应改变搜索幅度，在算法刚开始时，搜索幅度较大，搜索面积较广，可以较快的接近全局最优，经过不断的迭代，算法搜索面积减小，搜索精度提高，最终更加逼近最优解，从而具有更快的收敛速率，大大提高了图像分割精度。

另外，脊柱的分割容易受到肾脏，肋骨和一些病灶点的干扰，本发明利用一系列拓扑学的操作，通过剔除非连通区域，以及通过检测脊柱的宽度变化去除肾脏，肋骨的影响，以达到对脊柱的精确分割。

附图说明

图1是本发明的流程图。

图2是公式(Ⅰ)的函数图像。

图3是脊柱的原始CT图像。

图4是初步分割图，其中，(a)图是本发明所述方法的初步分割图，(b)图为传统方法的初步分割图。

图5是本发明所述方法的最终分割的脊柱区域图像。

图6是本发明所述方法的最终分割的脊柱区域映射到原始图像的边界图。

具体实施方式

下面结合附图对本发明做进一步的描述。

参照图1，本发明的具体实现步骤如下：

步骤1，输入图3所示的原始图像，求出图像的最小灰度值和最大灰度值。

步骤2，初始化第一代种群；最小灰度值至最大灰度值之间分别随机产生一个整数来初始化种群的个体的每一维；本发明实施例中将维数取为1。

步骤3，计算个体最优值和全局最优。

设图像被{t_{_1},t_{_2},…,t_{_M-1}}个阈值分割为M个部分，则最佳分割阈值{t^* _{_1},t^* _{_2},…,t^* _{_M-1}}需要满足如下条件：

{t^* _{_1},t^* _{_2},…,t^* _{_M-1}}＝argmax{σ²(t₁,t₂,…,t_M-1)}(Ⅲ)

其中，σ²(t₁,t₂,…,t_M-1)为图像灰度值的方差值， C_k为分割图像，k＝1,2,...,M，PR_i＝h(i)/N为各阈值的计算概率，为图像的总像素数，其中，L为所述图像的灰度值，h(i)为灰度值为i的像素数；

分割图像C₁对应阈值0,…,t_{_1}，分割图像C₂对应阈值t_{_1}+1,…,t_{_2}，...，分割图像C_M对应阈值t_{_M-1},…,L，0≤t_{_1}≤t_{_2}≤…≤t_{_M-1}≤L。

最大类间方差法由以下公式实现：

这里M取值为2，则最大类间方差法公式变为：

σ²＝ω₁·(u₁-u)²+ω₂·(u₂-u)²

其中，σ²表示类间距方差，ω₁表示第1类灰度值的概率，ω₂表示第2类灰度值的概率，u₁表示第1类灰度值的平均值，u₂表示第2类灰度值的平均值，u是整幅图像的平均灰度值。

其次，将类间距方差作为适应度函数；

再次，将种群个体带入到适应度函数计算出种群个体的适应度函数值；

然后，挑选种群个体中适应度最大的个体得到个体最优；

最后，挑选个体最优中适应度最大的个体得到全局最优；

步骤4，产生新的个体。

首先，建立粒子群的迭代公式，粒子按照下述公式进行更新：

k＝w·(1-(t/n)^1/r)(Ⅰ)

x_id(t+1)＝k·x_id(t)+c₁·rand₁()·(p_id-x_id(t))+c₂·rand₂()·(p_gd-x_id(t))(Ⅱ)

其中，自适应权重k可以给粒子提供一个自适应的搜索幅度，从而增强算法前期的全局搜索能力和后期的局部搜索能力；

上式(Ⅰ)中，t表示当前的迭代次数，本例中取t＝1,2...,360，n表示最大迭代次数，r是大于等于1的常数，w为惯性因子。r取值大小关系到粒子的移动速率的变化大小，r越小，k值变化范围就越大，那么粒子位置变化幅度越大，算法收敛速度相对越快，但是精度降低，r越大，k值变化范围就越小，那么粒子位置变化幅度越小，算法搜索精度较高，收敛速度较慢。根据实验表明取值3可以得到最好的效果。

上式(Ⅱ)中，c₁和c₂为基本粒子群算法中的加速常数，rand₁()和rand₂()表示两个在[0，1]范围里变化的随机值，p_id是个体最优值，p_gd是全局最优值，x_id(t)为第t代的第i个个体的第d维，x_id(t+1)为第t+1代的第i个个体的第d维；

步骤5，根据适应度值越大越好的原则，即计算每个粒子的的适应度值，从而得出个体最优值，全局最优值；具体的说，首先，根据步骤3建立的适应度函数，计算每个粒子的适应度值。那么，挑选种群个体中适应度最大的个体得到个体最优；挑选个体最优中适应度最大的个体得到全局最优；

然后根据上述粒子群迭代公式更新种群中各粒子的位置，产生新一代粒子。

步骤6，判断是否满足最大迭代次数360，如果满足，则得到最终的全局最优p_g，根据阈值对图像进行分割，输出图像矩阵。否则返回步骤(4)。

步骤7，以p_g的每一维的值，作为阈值，对读入的图像进行初步分割，得到分割后的图像。

步骤8，对步骤7的产生的图像，采用拓扑学操作，进一步精确分割。

首先，对于与脊柱不相连的热点、肩胛骨，通过连通性限制，剔除与脊柱不连通的区域；

其次，对于和脊柱直接相连的肋骨或者未被完全去除的肾脏，通过检测脊柱宽度的剧烈变化进行去除，如果脊柱宽度变化大于预先设定值width的话，我们认为有肋骨、肾脏等与脊柱相连造成的，于是对脊柱边界进行修正。

最后，输出图像。

本发明的效果可以通过以下仿真实验进一步说明：

1.仿真条件

本发明在Windows7，CPU-I52410m2.3GHZ，Matlab2010b。实验参数设置，类别取值为2，种群数取值为60，迭代次数取值为360。

2.实验内容

图2是公式k＝w·(1-(t/n)^1/r)的函数图像，横坐标是迭代次数，纵坐标是函数值k，其中，w取值为2，迭代次数n取值为60，r取值3。由图像可以看出，随着迭代次数的不断增加，权重越来越小，这样开始时，速度变化较大，搜索幅度也较大，全局搜索能力较强，到算法的后期，距离目标越来越近时，权重越来越小，搜索幅度变小，局部搜索能力较强，这样算法可以快速，精确的逼近目标。

图4(a)是本发明中自适应粒子群算法中的初步分割图，可以看出，脊椎区域基本能够分割出来，产生的噪声较少，4(b)传统粒子群算法的初步图像分割图，脊椎区域分割较为粗糙，噪声点较多。图5和图6是本发明的最终分割结果图，可以看出脊椎区域被准确的分割出来，脊椎附近的噪声被全部消除，边缘分割结果较好。

表1为本发明在不同种群数量以及迭代次数的情况下，将本发明和基于粒子群的图像分割算法做对比，获得的适应度函数值。

从表1可以看出，本发明与基于粒子群的图像分割算法相比，在相同的种群数量和迭代次数的情况下，本发明的适应度值优于传统的粒子群算法，因此本发明具有更有效更快速的全局搜索能力。

Claims (2)

1.一种基于粒子群算法的自适应脊柱CT图像分割方法，该方法包括如下步骤：

(1)输入原始图像，计算出图像的最小灰度值和最大灰度值；

(2)初始化第一代种群：在图像的最小灰度值至最大灰度值之间分别随机产生一个整数来初始化第一代种群任一个体的每一维，再以同样的方法初始化第一代种群下一个体的每一维，依次完成第一代种群每一个体的初始化；

(3)按以下方法计算个体最优值和全局最优值：

3a)利用OSTU法得到图像的类间距方差，将类间距方差作为适应度函数；

3b)将种群个体带入到适应度函数计算出种群个体的适应度函数值，然后，先挑选种群个体中适应度最大的作为个体得到个体最优，再从个体最优中挑选适应度最大的个体作为全局最优；

(4)按以下方法产生新的个体：

4a)先按下述公式公式(Ⅰ)计算出k是自适应权重：

k＝w·(1-(t/n)^1/r)(Ⅰ)，

式(Ⅰ)中，t表示当前的迭代次数，n表示最大迭代次数，r是大于等于1的常数，w为惯性因子；

4b)再将自适应权重按下述公式(Ⅱ)产生新一代的个体：

x_id(t+1)＝k·x_id(t)+c₁·rand₁()·(p_id-x_id(t))+c₂·rand₂()·(p_gd-x_id(t))(Ⅱ)，

式(Ⅱ)中，c₁和c₂为基本粒子群算法中的加速常数，rand₁()和rand₂()表示两个在[0，1]范围里变化的随机值，p_id是个体最优值，p_gd是全局最优值，x_id(t)为第t代的第i个个体的第d维，x_id(t+1)为第t+1代的第i个个体的第d维；

(5)按以下方法产生新的个体最优和全局最优：

将新个体带入适应度函数，计算新的种群新个体的适应度函数值；当新个体的适应度值比个体最优大，那么则用新个体代替个体最优，否则个体最优不变；当新个体的适应度值比全局最优的适应度函数值大则用个体最优代替全局最优，否则全局最优不变；

(6)判断是否满足最大迭代次数，如果满足，得到全局最优值，否则返回步骤(4)；

(7)以全局最优的每一维数据作为阈值，对得到的图像矩阵进行分割；

(8)采取拓扑学操作，对脊柱进一步准确的分割，输出图像。

2.根据权利要求1所述的一种基于粒子群算法的自适应脊柱CT图像分割方法，其特征在于步骤3a)所述的OSTU法是：设图像被{t_{_1},t_{_2},…,t_{_M-1}}个阈值分割为M个部分，则最佳分割阈值{t^* _{_1},t^* _{_2},…,t^* _{_M-1}}需要满足如下条件：

{t^* _{_1},t^* _{_2},…,t^* _{_M-1}}＝argmax{σ²(t₁,t₂,…,t_M-1)}(Ⅲ)

其中，σ²(t₁,t₂,…,t_M-1)为图像灰度值的方差值， C_k为分割图像，k＝1,2,...,M，PR_i＝h(i)/N为各阈值的计算概率，为图像的总像素数，其中，L为所述图像的灰度值，h(i)为灰度值为i的像素数；

分割图像C₁对应阈值0,…,t_{_1}，分割图像C₂对应阈值t_{_1}+1,…,t_{_2}，...，分割图像C_M对应阈值t_{_M-1},…,L，0≤t_{_1}≤t_{_2}≤…≤t_{_M-1}≤L。