基于改進蟻群算法的醫學圖像邊界提取

陳 超

(內江師范學院四川省數據恢復重點實驗室 四川 內江 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個像素的時候，既能保證醫學圖像檢測的速度又能保證醫學圖像檢測的效果。因本文算法搜索的目標是灰度小、梯度大的像素，易造成漏判和誤判，未來還需進一步研究和解決。