本發(fā)明涉及一種基于水平集的前列腺磁共振圖像分割方法，屬于磁共振圖像分割技術(shù)領(lǐng)域。

背景技術(shù)：

隨著人口增長和生活習(xí)慣的改變，前列腺癌的發(fā)病率和死亡率近年來呈明顯上升趨勢。臨床經(jīng)驗表明，前列腺癌若能盡早被發(fā)現(xiàn)，得到及時治療，有很高的存活率，因此，對于前列腺癌診斷和治療的相關(guān)研究具有重要意義。醫(yī)學(xué)影像作為前列腺癌診斷和治療的重要手段之一，發(fā)揮著越來越重要的作用。磁共振圖像因其具有對軟組織分辨率高，可多參數(shù)成像，能對任意斷層進行掃描的特點，被認為是目前前列腺癌診斷和輔助治療的最佳醫(yī)學(xué)影像。

圖像分割是基于圖像早期診斷和治療的基礎(chǔ)，是首要解決的關(guān)鍵問題。目前，前列腺磁共振圖像的分割研究大多集中在前列腺外輪廓分割，其分割方法主要有圖論、形變模型、特定理論以及混合的圖像分割方法等四大類。對基于磁共振圖像的前列腺內(nèi)外輪廓全分割研究近年來才開始展開。2011年，法國的Makni等人采用C-均值的方法最早實現(xiàn)了基于磁共振圖像的前列腺內(nèi)外輪廓全分割。2012年，荷蘭的Litjens等人通過模式識別的方法利用解剖、灰度值和紋理特征對前列腺內(nèi)部體素進行分類，實現(xiàn)了內(nèi)外輪廓全分割。這兩種方法都是先通過手動的方式實現(xiàn)外輪廓分割，并將其作為內(nèi)部分割的初始化，使前列腺的全分割費時費力。2013年，美國的Toth等人利用多重耦合水平集活動輪廓模型的方法進行了內(nèi)外輪廓全分割，該方法的分割效果更多依賴于所分割圖像的質(zhì)量，而且分割一張圖像需要訓(xùn)練大量圖集，消耗時間較長。2014年，加拿大的Qiu等人利用優(yōu)化連續(xù)最大流模型的方法進行了前列腺的內(nèi)外全分割，改善了分割效果并提高了分割效率。由于前列腺內(nèi)部區(qū)域磁共振圖像存在噪聲，灰度顯示不均勻，區(qū)域邊界模糊不清，使基于磁共振圖像的前列腺內(nèi)外輪廓全分割研究一直是個挑戰(zhàn)。

技術(shù)實現(xiàn)要素：

針對上述問題，本發(fā)明要解決的技術(shù)問題是提供一種基于水平集的前列腺磁共振圖像分割方法，在構(gòu)建統(tǒng)一水平集能量函數(shù)的基礎(chǔ)上，首先基于縱向弛豫時間圖像，將前列腺從其周圍組織中分割出來，即實現(xiàn)前列腺的外輪廓分割，其次基于橫向弛豫時間圖像，將外輪廓作為約束條件來實現(xiàn)前列腺的內(nèi)部區(qū)域分割，進而實現(xiàn)前列腺的內(nèi)外輪廓全分割。

上述目的主要通過以下方案實現(xiàn)：

本發(fā)明的一種基于水平集的前列腺磁共振圖像分割方法，其特征在于：所述方法的具體實現(xiàn)過程為：

步驟一：定義水平集演化方程

在Ω域內(nèi)定義一個水平集函數(shù)能量函數(shù)ε(φ)定義為：

ε(φ)＝μR_p(φ)+αε_drive(φ) (1)

其中，R_p(φ)是水平集的距離調(diào)整項，ε_drive(φ)是輪廓驅(qū)動能量項，μ＞0，α＜0，都為常數(shù)；

水平集的距離調(diào)整項R_p(φ)定義為：

其中，p是能量密度函數(shù)，

能量密度函數(shù)構(gòu)造為：

能量密度函數(shù)p(s)具有兩個極值點，分別是s＝0和s＝1，其一階導(dǎo)數(shù)和二階導(dǎo)數(shù)為：

式(2)中函數(shù)R_p(φ)的加托導(dǎo)數(shù)為：

其中，函數(shù)d_p定義為：

輪廓驅(qū)動能量項ε_drive(φ)定義為：

其中，g是邊界約束函數(shù)，H是單位階躍函數(shù)，通常將單位階躍函數(shù)H近似地用函數(shù)H_ε來代替，且定義為：

H_ε的導(dǎo)數(shù)δ_ε為：

輪廓驅(qū)動能量函數(shù)ε_drive(φ)的加托導(dǎo)數(shù)為：

求解梯度流方程的穩(wěn)態(tài)解，

其中，是函數(shù)ε(φ)的加托導(dǎo)數(shù)；

將式(6)和式(11)代入(12)中，可以得到能量函數(shù)ε(φ)的梯度流表達式為：

式(13)所示的偏微分方程就是基于式(1)的前列腺內(nèi)外輪廓分割的水平集演化方程；

瞬態(tài)偏導(dǎo)數(shù)可以近似采用正向有限差分方程進行求解，時變函數(shù)φ(x,y,t)的離散形式用來表示，則水平集演化方程可以離散為如下所示的有限差分方程：

步驟二：外輪廓分割

讀取原始的縱向弛豫時間圖像，選擇外分割初始化方法-變形橢圓法：

基本橢圓參數(shù)方程如式(15)所示：

其中，a_x是x方向的半軸長，a_y是y方向的半軸長；

沿著y軸通過轉(zhuǎn)換基本的橢圓方程獲得變形橢圓的參數(shù)方程ψ(x_d,y_d)，如式(16)所示：

其中，

將已確定的前列腺變形橢圓所確定的區(qū)域設(shè)為S_e，則初始水平集函數(shù)為：

其中，c₀為正常數(shù)；

在式(16)和式(17)中，(x_c,y_c)是變形橢圓的中心坐標，t_y∈[-1,1]是描述橢圓上部沿著y軸方向線性變尖的參數(shù)，b_y∈[-1,0)∪(0,1]是描述橢圓下部沿著y軸方向內(nèi)凹彎曲的參數(shù)，調(diào)整式(16)和式(17)相應(yīng)的參數(shù)，使得可變形橢圓最大限度逼近前列腺的外輪廓形狀；

然后，確定外輪廓邊界約束函數(shù)：

在縱向弛豫時間圖像中，假定I為前列腺圖像，定義圖像I的邊界指示器為：

其中，G_σ是方差為σ的高斯核，將式(19)作為前列腺外輪廓分割的邊界約束函數(shù)，并給定參數(shù)值。

最后對水平集演化方程(14)進行迭代求解，實現(xiàn)前列腺的外輪廓分割；

步驟三：內(nèi)部區(qū)域分割

讀取原始橫向弛豫時間圖像，選擇內(nèi)分割初始化方法-多線段擬合法：

在中央腺內(nèi)依次選取N個點，使這N個點首尾相連形成一封閉區(qū)域，設(shè)為S_N，則初始水平集函數(shù)為：

其中，c₀為正常數(shù)；

然后，確定內(nèi)輪廓邊界約束函數(shù)：

采用全向邊界梯度作為邊界指示器來描述前列腺中央腺圖像的邊界特征，假定I為前列腺圖像，I_i,j為I的某一元素，設(shè)定為中心元素，其相鄰的8元素分別為I_i-1,j-1，I_i-1,j，I_i-1,j+1，I_i,j-1，I_i,j+1，I_i+1,j-1，I_i+1,j，I_i+1j+1，為求取這8元素與中心元素的差值，定義如下對應(yīng)的8個卷積模板，

中心元素與相鄰8元素的差值計算為：

Dif_lu＝conv2(I,Temp_lu,'same') (29)

Dif_u＝conv2(I,Temp_u,'same') (30)

Dif_ru＝conv2(I,Temp_ru,'same') (31)

Dif_l＝conv2(I,Temp_l,'same') (32)

Dif_r＝conv2(I,Temp_r,'same') 33)

Dif_ld＝conv2(I,Temp_ld,'same') (34)

Dif_d＝conv2(I,Temp_d,'same') (35)

Dif_rd＝conv2(I,Temp_rd,'same') (36)

conv2是卷積運算符，圖像I的全向邊界梯度函數(shù)定義為：

Grad_I＝[Grad_I_x Grad_I_y Grad_I_xy- Grad_I_xy+] (37)

其中，各項分別定義為：

圖像I的全向邊界梯度模定義為：

|Grad_I|＝sqrt(Grad_I_x²+Grad_I_y²+Grad_I_xy-²+Grad_I_xy+²) (42)

式(12)中的邊界約束函數(shù)為：

式(43)稱為前列腺內(nèi)輪廓分割的邊界約束函數(shù)，并給定參數(shù)值；

最后對水平集演化方程(14)進行迭代求解，獲得前列腺內(nèi)部中央腺的輪廓；將第二步得到的外輪廓與第三步所得到的中央腺輪廓進行區(qū)域相減，便得到前列腺外周帶區(qū)域，進而實現(xiàn)了前列腺的全面分割。

本發(fā)明的有益效果為：

1、提出了基于縱向弛豫時間圖像與橫向弛豫時間圖像相結(jié)合的前列腺兩步分割法，綜合了縱向弛豫時間圖像內(nèi)外信號對比鮮明，可實現(xiàn)外輪廓分割和橫向弛豫時間圖像內(nèi)部結(jié)構(gòu)顯示清晰，外周帶與中央腺信號形成明顯對比，可實現(xiàn)內(nèi)區(qū)域的分割的優(yōu)點，克服了縱向弛豫時間圖像難以區(qū)分內(nèi)部結(jié)構(gòu)和橫向弛豫時間圖像內(nèi)部清晰的多區(qū)域信號灰度值將干擾外輪廓的提取與分割的缺點。

2、所建立的能量函數(shù)中融入了距離調(diào)整項，在演化過程中可以不斷的進行調(diào)整，引起的周圍擴散效應(yīng)可以維持期望的形狀與期望輪廓附近的距離，在不必重新初始化的同時避免了通用水平集方法由于不斷的初始化引起的數(shù)值錯誤。

3、采用了初始輪廓接近內(nèi)外輪廓的外分割初始化-變形橢圓法和內(nèi)分割初始化-多線段擬合法來分別初始化水平集函數(shù)，可以提高分割效果。

附圖說明

為了易于說明，本發(fā)明由下述的具體實施及附圖作以詳細描述。

圖1為本發(fā)明方法的流程圖；

圖2為本發(fā)明外分割初始化-變形橢圓法中不同a_x參數(shù)的輪廓示意圖；

圖3為本發(fā)明外分割初始化-變形橢圓法中不同b_y參數(shù)的輪廓示意圖；

圖4為本發(fā)明的人機交互界面。

具體實施方式

為使本發(fā)明的目的、技術(shù)方案和優(yōu)點更加清楚明了，下面通過附圖中示出的具體實施例來描述本發(fā)明。但是應(yīng)該理解，這些描述只是示例性的，而并非要限制本發(fā)明的范圍。此外，在以下說明中，省略了對公知結(jié)構(gòu)和技術(shù)的描述，以避免不必要地混淆本發(fā)明的概念。

如圖1、圖2、圖3、圖4所示，本具體實施方式采用以下技術(shù)方案：一種基于水平集的前列腺磁共振圖像分割方法，其特征在于：所述方法的具體實現(xiàn)過程為：

步驟一：定義水平集演化方程

在Ω域內(nèi)定義一個水平集函數(shù)能量函數(shù)ε(φ)定義為：

ε(φ)＝μR_p(φ)+αε_drive(φ) (1)

其中，R_p(φ)是水平集的距離調(diào)整項，ε_drive(φ)是輪廓驅(qū)動能量項，驅(qū)動水平集函數(shù)曲線向前列腺輪廓邊界運動，μ＞0，α＜0，都為常數(shù)；

水平集的距離調(diào)整項R_p(φ)定義為：

其中，p是能量密度函數(shù)，

為了避免邊界效應(yīng)，能量密度函數(shù)構(gòu)造為：

能量密度函數(shù)p(s)具有兩個極值點，分別是s＝0和s＝1，其一階導(dǎo)數(shù)和二階導(dǎo)數(shù)為：

式(2)中函數(shù)R_p(φ)的加托導(dǎo)數(shù)為：

其中函數(shù)d_p定義為：

輪廓驅(qū)動能量項ε_drive(φ)定義為：

其中，g是邊界約束函數(shù)，H是單位階躍函數(shù)，通常將單位階躍函數(shù)H近似地用函數(shù)H_ε來代替，且定義為：

H_ε的導(dǎo)數(shù)δ_ε為：

輪廓驅(qū)動能量函數(shù)ε_drive(φ)的加托導(dǎo)數(shù)為：

為了求取能量函數(shù)ε(φ)的最小值，常規(guī)方法就是求解梯度流方程的穩(wěn)態(tài)解：

其中，是函數(shù)ε(φ)的加托導(dǎo)數(shù)；

將式(6)和式(11)代入(12)中，可以得到能量函數(shù)ε(φ)的梯度流表達式為：

式(13)所示的偏微分方程就是基于式(1)的前列腺內(nèi)外輪廓分割的水平集演化方程；

瞬態(tài)偏導(dǎo)數(shù)可以近似采用正向有限差分方程進行求解，時變函數(shù)φ(x,y,t)的離散形式用來表示，這樣水平集演化方程可以離散為如下所示的有限差分方程：

步驟二：外輪廓分割

讀取原始的縱向弛豫時間圖像，選擇外分割初始化方法-變形橢圓法：

基本橢圓參數(shù)方程如式(15)所示：

其中，a_x是x方向的半軸長，a_y是y方向的半軸長；

由于前列腺橫斷軸位各層的外輪廓形狀是沿著y軸改變的，因此沿著y軸通過轉(zhuǎn)換基本的橢圓方程，可以獲得變形橢圓的參數(shù)方程ψ(x_d,y_d)，如式(16)所示：

其中，

將已確定的前列腺變形橢圓所確定的區(qū)域設(shè)為S_e，則初始水平集函數(shù)為：

其中，c₀為正常數(shù)。

在式(16)和式(17)中，(x_c,y_c)是變形橢圓的中心坐標，t_y∈[-1,1]是描述橢圓上部沿著y軸方向線性變尖的參數(shù)，b_y∈[-1,0)∪(0,1]是描述橢圓下部沿著y軸方向內(nèi)凹彎曲的參數(shù)。通過調(diào)整式(19)和式(20)相應(yīng)的參數(shù)，使得可變形橢圓最大限度逼近前列腺的外輪廓形狀；

然后，確定外輪廓邊界約束函數(shù)：

在縱向弛豫時間圖像中，假定I為前列腺圖像，定義圖像I的邊界指示器為：

其中，G_σ是方差為σ的高斯核，式中卷積用來平滑前列腺圖像，減小噪聲的影響，將式(19)作為前列腺外輪廓分割的邊界約束函數(shù)，并給定參數(shù)值。

最后對水平集演化方程(14)進行迭代求解，實現(xiàn)前列腺的外輪廓分割；

步驟三：內(nèi)部區(qū)域分割

讀取原始橫向弛豫時間圖像，選擇內(nèi)分割初始化方法-多線段擬合法：

在中央腺內(nèi)依次選取N個點，使這N個點首尾相連形成一封閉區(qū)域，設(shè)為S_N，則初始水平集函數(shù)為：

其中，c₀為正常數(shù)；

然后，確定內(nèi)輪廓邊界約束函數(shù)：

采用全向邊界梯度作為邊界指示器來描述前列腺中央腺圖像的邊界特征，

假定I為前列腺圖像，I_i,j為I的某一元素，設(shè)定為中心元素。其相鄰的8元素分別為I_i-1,j-1，I_i-1,j，I_i-1,j+1，I_i,j-1，I_i,j+1，I_i+1,j-1，I_i+1,j，I_i+1,j+1。為求取這8元素與中心元素的差值，定義如下對應(yīng)的8個卷積模板，

中心元素與相鄰8元素的差值計算為：

Dif_lu＝conv2(I,Temp_lu,'same') (29)

Dif_u＝conv2(I,Temp_u,'same') (30)

Dif_ru＝conv2(I,Temp_ru,'same') (31)

Dif_l＝conv2(I,Temp_l,'same') (32)

Dif_r＝conv2(I,Temp_r,'same') 33)

Dif_ld＝conv2(I,Temp_ld,'same') (34)

Dif_d＝conv2(I,Temp_d,'same') (35)

Dif_rd＝conv2(I,Temp_rd,'same') (36)

conv2是卷積運算符，圖像I的全向邊界梯度函數(shù)定義為：

Grad_I＝[Grad_I_x Grad_I_y Grad_I_xy- Grad_I_xy+] (37)

其中，各項分別定義為：

圖像I的全向邊界梯度模定義為：

|Grad_I|＝sqrt(Grad_I_x²+Grad_I_y²+Grad_I_xy-²+Grad_I_xy+²) (42)

式(12)中的邊界約束函數(shù)為：

式(43)稱為前列腺內(nèi)輪廓分割的邊界約束函數(shù)，并給定參數(shù)值；

最后對水平集演化方程(14)進行迭代求解，獲得前列腺內(nèi)部中央腺的輪廓；將第二步得到的外輪廓與第三步所得到的中央腺輪廓進行區(qū)域相減，便得到前列腺外周帶區(qū)域，進而實現(xiàn)了前列腺的全面分割。

以上顯示和描述了本發(fā)明的基本原理和主要特征和本發(fā)明的優(yōu)點。本行業(yè)的技術(shù)人員應(yīng)該了解，本發(fā)明不受上述實施例的限制，上述實施例和說明書中描述的只是說明本發(fā)明的原理，在不脫離本發(fā)明精神和范圍的前提下，本發(fā)明還會有各種變化和改進，這些變化和改進都落入要求保護的本發(fā)明范圍內(nèi)。本發(fā)明要求保護范圍由所附的權(quán)利要求書及其等效物界定。