一種基于水平集的前列腺磁共振圖像分割方法與流程

本發明涉及一種基于水平集的前列腺磁共振圖像分割方法，屬于磁共振圖像分割技術領域。

背景技術：

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

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

技術實現要素：

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

上述目的主要通過以下方案實現：

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

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

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

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

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

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

其中，p是能量密度函數，

能量密度函數構造為：

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

式(2)中函數R_p(φ)的加托導數為：

其中，函數d_p定義為：

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

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

H_ε的導數δ_ε為：

輪廓驅動能量函數ε_drive(φ)的加托導數為：

求解梯度流方程的穩態解，

其中，是函數ε(φ)的加托導數；

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

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

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

步驟二：外輪廓分割

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

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

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

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

其中，

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

其中，c₀為正常數；

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

然后，確定外輪廓邊界約束函數：

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

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

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

步驟三：內部區域分割

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

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

其中，c₀為正常數；

然后，確定內輪廓邊界約束函數：

采用全向邊界梯度作為邊界指示器來描述前列腺中央腺圖像的邊界特征，假定I為前列腺圖像，I_i,j為I的某一元素，設定為中心元素，其相鄰的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元素與中心元素的差值，定義如下對應的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的全向邊界梯度函數定義為：

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)中的邊界約束函數為：

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

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

本發明的有益效果為：

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

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

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

附圖說明

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

圖1為本發明方法的流程圖；

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

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

圖4為本發明的人機交互界面。

具體實施方式

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

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

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

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

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

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

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

其中，p是能量密度函數，

為了避免邊界效應，能量密度函數構造為：

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

式(2)中函數R_p(φ)的加托導數為：

其中函數d_p定義為：

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

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

H_ε的導數δ_ε為：

輪廓驅動能量函數ε_drive(φ)的加托導數為：

為了求取能量函數ε(φ)的最小值，常規方法就是求解梯度流方程的穩態解：

其中，是函數ε(φ)的加托導數；

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

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

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

步驟二：外輪廓分割

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

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

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

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

其中，

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

其中，c₀為正常數。

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

然后，確定外輪廓邊界約束函數：

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

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

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

步驟三：內部區域分割

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

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

其中，c₀為正常數；

然后，確定內輪廓邊界約束函數：

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

假定I為前列腺圖像，I_i,j為I的某一元素，設定為中心元素。其相鄰的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元素與中心元素的差值，定義如下對應的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的全向邊界梯度函數定義為：

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)中的邊界約束函數為：

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

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

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