東南大學醫學院附屬江陰醫院 影像科,江蘇 江陰 214400
MRI成像由于良好的空間分辨率、較高的組織對比度和非侵入性等特性,被廣泛應用于圖像引導介入、外科手術引導、放療計劃等方面[1-3]。MRI圖像處理是診斷腦部疾病的基礎,其中MR腦部圖像分割作為圖像處理的第一步,能輔助識別腦白質、灰質和腦脊液形態,提取病變區域,達到定量評估癡呆、精神分裂癥、阿爾茲海默等神經性疾病,但MRI腦部圖像分割一直面臨著諸多挑戰。
MR圖像分割算法主要分為手工分割、半自動分割技術和全自動分割技術[4-6]。臨床醫師常用手工分割方法,當病灶數量較少且對比度較大時,主要操作簡單、分割精度相對較高,但分割結果主觀差異大,耗時較多。半自動分割技術主要包括圖像閾值法、區域增長法、聚類分析法、基于特定理論法等,其中閾值法簡單高效,借助圖像灰度直方圖選取分割閾值,將圖像劃分為背景和目標區域,但對灰度不均圖像經常無法確定閾值,且對噪聲敏感,分割精度很難得到保證。基于區域分割算法是根據相似性準則將像素或區域聚合為更大區域的過程,能有效克服空間不連續的缺點,但對相似性準則要求較高,且容易過度分割。聚類分析將圖像空間像素用對應特征空間點表示,然后根據距離準則對特征空間進行分割,適用于存在不確定性和模糊性的圖像,對初始聚類中心設置極為敏感。全自動分割技術不必借助腦圖譜等分割參考圖像,往往需要借助圖像紋理、位置、方向等經驗知識對圖像數據進行訓練,然后根據概率論構建自動分割模型。
隱馬爾可夫隨機場是一種有效的圖像分割模型[7-8],主要基于最大后驗概率準則模擬目標函數最小值。而共軛梯度算法(Conjugate Gradient,CG)是一種最實用的參數優化算法,能快速準確計算目標函數最小值。基于此,本研究嘗試將隱馬爾可夫隨機場和共軛梯度算法應用于腦部MRI圖像分割。
馬爾科夫隨機場是一個由無向圖表示的概率分布模型[9],圖中每個結點表示一個或者一組變量,結點之間的邊表示兩個變量之間的依賴關系。隱馬爾可夫隨機場則是一種統計模型,用來描述一個隱含未知參數的馬爾可夫過程[10]。假定S={s1,s2,….sM}是圖像節點的集合,圖像由M個像素構成,每個節點s有一個鄰域集合Vs(S),鄰域系統V(S)符合公式(1)描述的特質。

則第r個鄰域系統X可由公式(2)定義,式中d(s,t)表示像素s與t之間的歐氏距離,僅由像素位置決定。

定義c為S的一個子集,其中所有元素互為鄰域,對于一個非獨立點集合,滿足公式(3),其中第p個子集Cp包含p個點。

采用y=(y1,y2,…,yM)表示待分割圖像的像素點,x=(x1,x2,…,xM)表示已分割完成圖像的像素點,M為像素點數量。yi和xi分別為第i個像素點灰度值,y和x可看成馬爾可夫隨機場家族Y的觀測數據和X的對應真實分割結果,其中馬爾可夫變量灰度空間Ey={0,1,…,255},離散空間Ex={1,…,K},K為圖像區域分類數量。分割圖像過程即可描述為在y發生時尋找最優X中x發生的概率,隱馬爾可夫隨機場可以通過最大化條件概率模擬此過程,見公式(4)。

根據貝葉斯定理[11-12],可將公式(4)轉換成公式(5),P[Y=y]是常數。

由于所有像素點均相互獨立,且根據中心極限定理,相同類別數據在大樣本下總是趨于高斯分布,因此所有像素點灰度值滿足正態分布,可得公式(6)和(7)。

由于所有像素點相互獨立且同步,依據Hammerskey-Clifford定理,該隱馬爾可夫隨機場等價于Gibbs分布[13],可得公式(8),其中T為控制參數。

此處U(x)是由Potts模型定義的條件似然能量函數[14],見公式(9),其中β是常數。δ是克羅內克函數,見公式(10)。

將公式(7)和(8)帶入公式(5),可得公式(11),μxs,σxs分別是類別xs的均值和標準差,當β>0時,優先被分割均勻度較高區域,且區域大小由β控制,A是常數。

求概率函數P[X=x|Y=y]最大值等價于求函數Ψ(x,y)的最小值,如公式(12)。

直接準確計算分割結果x*是不可能的,因此采用最優化算法計算最佳估計值是x? 必要的。
設 μ=(μ1,μ2,…,μK)和 σ=(σ1,σ2,…,σK)分 別 為 分 割 圖 像x=(x1,x2,…,xM)中K類區域像素灰度均值和標準差,見公式(13),則計算函數Ψ(x,y)的最小值可轉化為求Ψ(μ)的最小值。可根據ys歸類到最接近的μj來計算x,如果最接近ys的均值是μj,則xs=j,因此尋找μ*即可得到x*,μ集合為μ=[0,…255]K,見公式 (14)。

由此,隱馬爾科夫隨機場模型將圖像分割過程轉化為尋求目標函數Ψ(μ)最小值,本研究采用共軛梯度算法進行求解。共軛梯度算法是一種最小化類的迭代算法,主要思想是[15-17]:選擇一個優化方向后,采用本次選擇的步長更新完此方向的誤差,隨后在優化更新過程中不再朝此方向更新。由于每次將一個方向優化到極小,后面優化過程不再影響之前優化方向上的極小值,因此理論上對N維問題求極小只用對N個方向都求出極小即可,但需要每次優化方向共軛正交。
采用共軛梯度算法求解Ψ(μ)最小值時,需對Ψ(μ)進行一階求導,此處一階導數由中心差逼近形式表達,見公式(15),則此一階導數的充分逼近取決于ε值,經測試ε取0.01時效果最佳。

圖像分割準確性采用Dice相似性系數和特異性(Sensitivity)評估[18-19],表示分割結果與真實結果的接近程度,取值范圍均為[0,1],評價指標越接近于1表明分割結果越精確,見公式(16)和(17)。其中A代表分割結果,B代表真實結果,TP、TN、FP和FN分別代表真陽性、假陽性、假陽性和假陰性。

為了驗證本研究提出分割算法的可行性和準確性,選用兩組不同分辨率的腦部MR-T1加權圖像進行仿真實驗。①選自IBSR腦圖庫臨床實例MRI圖像30幅[20],圖像大小為256×256×63 個體素,體素大小1 mm×3 mm×1 mm,算法中常數β設置為1,溫度T為10,初始點μ0=(1,5,140,190);② 選自BrainWeb腦圖庫仿真圖像20幅[21],圖像大小為181×217×181 個體素,體素大小1 mm×1 mm×1 mm,且在圖像中加入(0,0)、(3%,20%)和(5%,20%)的噪聲。算法中參數設置如下:常數β為1,Brain Web-3控制參數T分別為10、4和1,初始點μ0均為(1,45,110,150)。并將本研究提出算法與經典馬爾可夫隨機場算法和三種改進的馬爾可夫隨機場算法進行比較,所有的仿真實驗均在MATLAB平臺上實現。
本研究提出算法能完整清晰地識別出IBSR庫腦部MRI臨床圖像的白質、灰質、腦脊液和背景部分(圖1),表明本研究提出算法的可行性。經典MRF算法所得平均Dice相似性系數和特異性值均最低,表面該算法分割效果最差。其次,兩種改進MRF算法所得平均Dice相似性系數和特異性值稍高于MRF算法,且大致相等,表明改進的MRF算法分割性能有所提升。基于本研究提出算法分割所得灰質、白質、腦脊液和總體平均Dice相似性系數分別達到0.895、0.855、0.832和0.848,特異性值分別達到0.902、0.912、0.923和0.912,且Dice相似性系數在各功能分區均有所提升:灰質(10.41%~11.12%)、白質(3.14%~3.39%)、腦脊液(4.26%~6.26%)和總體平均(5.87%~6.67%);特異性值提升如下:灰質(1.69%~7.89%)、白質(2.93%~6.54%)、腦脊液(5.01%~9.75%)和總體平均(3.17%~8.16%),表明本研究提出算法的有效性和優越性。定量分析結果,見表1。

圖1 IBSR庫腦部MRI圖像分割結果

表1 不同算法分割IBSR庫腦部MR圖像所得Dice相似性系數
由圖2可知,本研究提出算法可將BrainWeb庫不同噪聲水平MRI腦部合成圖像準確分割成灰質、白質、腦脊液和背景部分,表明本研究提出算法的準確性與抗噪性。定量分析結果如表2所示,在噪聲水平為0時,MRF改進算法所得腦區各部分Dice相似性系數和特異性值均低于本研究提出的HMRF-CG算法,表明本研究提出算法的優越性。噪聲水平從(0,0)增加到(5%,20%)時,MRF改進算法所得平均Dice相似性系數有所升高,從0.705上升到0.918,平均特異性值呈下降趨勢,從0.864下降到0.816,主要是該算法誤將噪聲誤劃分為圖像目標區域像素所致,導致分割特異性下降。基于本研究提出算法不同噪聲水平所得Dice相似性系數特異性值均微弱下降,但分割精度仍大于0.9,且均高于MRF改進算法所得對應值。不同噪聲水平下不同腦分區Dice相似性系數下降幅度為:灰質(39.17%~0.66%)、白質(48.43%~0.11%)、腦脊液(27.96%~3.47%)和總體平均值(38.15%~1.42%),特異性值提升幅度為:灰質(14.20%~15.86%)、白質(14.13%~14.27%)、腦脊液(10.20%~12.35%)和總體平均值(12.85%~14.22%),表明本研究提出算法的抗噪性和特異性均較強。

圖2 BrainWeb庫不同噪聲水平MRI腦部圖像分割結果
本文聯合使用隱馬爾可夫隨機場和共軛梯度算法分割MRI腦部圖像,其中隱馬爾可夫隨機場為圖像分割過程建模,共軛梯度算法提供目標函數最優化解法。定性與定量分析顯示,基于本研究提出算法能獲得清晰準確地劃分腦部各功能分區,臨床實例和包含噪聲MRI腦部圖像分割所得Dice相似性系數和特異性值均高于其他算法,表明本研究提出算法的可行性、抗噪性和魯棒性,適用于臨床MRI腦部圖像分割。

表2 Brainweb庫不同噪聲水平腦部圖像分割所得Dice系數