基于SVM和有監督描述子學習算法的腦MR圖像顱骨分割方法?

黃勇其 史文博 周志勇 龐樹茂 佟寶同 趙凌霄 戴亞康

（1.中國科學院蘇州生物醫學工程技術研究所蘇州215163）（2.中國科學院大學北京100049）（3.北京師范大學北京100875）（4.南方醫科大學廣州510515）

基于SVM和有監督描述子學習算法的腦MR圖像顱骨分割方法?

黃勇其1，2史文博3周志勇1龐樹茂4佟寶同1趙凌霄1戴亞康1

（1.中國科學院蘇州生物醫學工程技術研究所蘇州215163）（2.中國科學院大學北京100049）（3.北京師范大學北京100875）（4.南方醫科大學廣州510515）

神經電流源定位研究首先要解決EEG、MEG正問題的計算。在求解MEG和EEG正問題的過程中，為了精確地計算傳導矩陣，常常需要對腦組織進行分層建模。在腦MR圖像中，雖然軟組織能被清晰地成像，但顱骨卻由于缺少氫而呈現低灰度值，從而很難自動分割出顱骨。因此如何從腦MR圖像中準確、自動分割出顱骨是解決MEG、EEG正問題的關鍵。為解決上述問題，提出一種基于支持向量機的自動腦MR圖像顱骨分割方法，提取病人MR圖像的全局特征和局部特征進行訓練，并結合有監督描述子學習算法SDL，將得到的特征矩陣進行壓縮，去掉冗余的特征，得到一個緊湊的特征描述，最終利用SVM從腦MR圖像中自動識別出骨骼。實驗結果表明，采用支持向量機結合有監督描述子學習算法的分割方法與僅使用支持向量機和僅使用數學形態學方法相比，分割精度進一步提升，Dice分割精度分別為0.832，0.798，0.482，從而解決了從腦MR圖像自動分割顱骨的任務，并為解決EEG和MEG正問題的研究奠定基礎。

顱骨分割；支持向量機；有監督描述子學習算法；特征提取；特征壓縮

Class NumberTP391

1 引言

隨著MEG和EEG技術的快速發展，腦磁和腦電信號在研究腦神經生理活動中扮演著越來越重要的角色。在EEG和MEG神經電流源定位的正問題計算中，以往通常利用簡單的球模型［1］來模擬大腦結構，進而求解傳導矩陣。然而，人腦結構異常復雜，顱骨凹凸性各異，難以從簡單的球模型中求得傳導矩陣的解，進而定位出神經電流源的位置。因此若能從人腦MR圖像中提取出顱骨結構，就能針對具體的個人構建出一個真實的頭模型［2］，從而準確解決神經電流源定位的問題。

真實的顱骨模型為解決EEG正問題提供極大的便利，到目前為止，有不少有關大腦分割的研究，例如文獻［3］中作者提出一種基于貝葉斯分層混合模型的方法成功實現大腦FMRI圖像的自動分割；文獻［4］中作者基于最大后驗概率對腦MR圖像進行三維分割，成功分割出腦白質與灰質。然而僅有少數的顱骨分割算法被提出。例如Heinonen等［5］用閾值和區域生長方法從MR圖像中分割出骨骼，該方法對后腦等形狀較規則的區域能形成較好的分割，但對特定的骨骼區域，如眼廓，因為部分容積效應的原因并不能很好地分割出骨骼；Rifai等［6］應用一種可形變模型從MR圖像中分割出骨骼，該方法在仿真MRI圖像中取得了較好的分割結果，但在真實MRI圖像中，可形變模型會帶來不確定邊界的情況，例如分割的結果里面混合了皮膚，肌肉，眼睛，內耳的組織；Chu和Takaya［7］用閾值和高斯拉普拉斯算子對連續橫斷層的頭皮-顱骨的邊界和顱骨-大腦的邊界進行檢測，該方法對于大腦顱骨外邊界能達到較好的檢測，對部分頭皮-顱骨的邊界也能較好區分，但對于結構較復雜區域，如眼窩，內耳區域等不能明顯進行區分；B.Dogdas等［8］運用數學形態學開操作閉操作等方法將顱骨分割出來，該方法能較好地分割出大腦區域，對顱骨外層也能取得較好分割結果，但數學形態學開閉操作的方法受限于需要人為地設定一個閾值，魯棒性并不高。

圖像的特征壓縮能大大縮短顱骨分割的時間，目前比較流行的特征壓縮方法有：Golub和Van Loan［9］提出的奇異值分解（SVD）方法，該方法用特征矩陣的方式取代特征向量的表示方法，在降維的同時保留了數據的關鍵信息，實現了原始數據的低維近似描述。Dhillon和Modha［10］，Srebro和Jaakko?la［11］，Castelli等［12］，都對數據壓縮進行了廣泛的研究，在提取相似度最高的數據和移除噪聲方面取得了較好的結果，然而隨著數據量的增大，它們的算法均不同程度地受到時間復雜度和空間復雜度的限制。Achilioptas和McSherry［13］采用隨機采樣的方法來加速奇異值分解以得到原始數據的低維近似描述，然而該研究的效率很大程度上取決于數據矩陣的能譜結構。另一個應用廣泛的數據壓縮方法是壓縮感知［14］，壓縮感知理論利用少量的線性測量值就能感知信號的初始結構，然后通過求解最優化問題實現信號的精確重構。盡管壓縮感知得到了較廣泛的應用，但仍有一些需要解決的問題，例如在稀疏表示方面，如何引入流形結構，在重建算法方面，如何解決重建時間過長的問題等。

為了從腦MR圖像中分割出顱骨，解決EEG、MEG源成像問題，本文提出一種基于支持向量機結合有監督描述子學習算法SDL的自動腦圖像顱骨分割算法，以圖像的局部和全局信息作為輸入，提取MR圖像的概率圖譜、鄰域強度、統計矩等有代表性的多個特征值，并利用SDL算法對得到的特征矩陣進行特征選擇，去除冗余的特征值，最后將得到的經壓縮后的特征矩陣輸入到支持向量機進行訓練，達到自動分割腦MR顱骨的目的，并與僅用支持向量機的方法和數學形態學的方法進行對比，證明了本方法的有效性。

2 方法

本文分割方法主要包括三部分內容，首先將圖像數據進行配準并提取多個關鍵特征，再利用有監督描述子學習算法（SDL）進行特征壓縮，最后利用支持向量機SVM訓練得到的特征矩陣。

2.1 圖像特征提取

機器學習算法最重要的步驟之一就是特征的選取，特別地，隨著數據量的增加，特征空間的擴大，大小合適的特征數量就顯得尤為重要。這樣的特征矩陣不僅能減少SVM訓練時間，而且也不會出現過擬合的現象。

本文選取的特征主要包括如下三類特征：

1）概率圖譜：該特征為全局特征，選取數據集中的某個病人的CT圖像作為參考CT，并將其他病人的CT圖像以剛性配準的方法配準到參考病人的CT圖像中，然后以300H單元的閾值分割出骨骼。再從分割出的所有病人的CT顱骨圖像中求得平均值，得到的結果即為概率圖譜。結果如圖1所示。

算法的全部數據都是基于參考MR坐標系的，因此需要將得到的概率圖譜配準到參考MR坐標系中。配準流程和訓練流程如圖2所示，每個病人的訓練CT均采用剛性配準的方法配準至相應的MR中，并最終配準至參考MR。待分割目標MR也要配準到參考MR中，之后再利用本文算法進行分割。

圖18組CT圖像分別分割顱骨后求得的平均值，圖中為第34層，取值范圍［0，1］

圖2來自不同病人數據的多模態的配準概述

2）鄰域強度：以體素為中心，提取每個坐標軸的7個體素的灰度值作為局部特征（總共19個灰度值）。在建立體素模型時，不同大小的體素區域其三維模型差別非常明顯，當體素區域過小時，不能準確表達三維圖像的表面形態，而當體素區域過大時，則會占用大量的物理空間，造成浪費，因此選擇合適的體素區域非常重要。本文選擇圖3所示體素區域用以獲取鄰域強度特征值。

圖3用于計算鄰域強度的體素區域

3）統計矩：為了反映圖像的局部紋理特征，選取統計矩μn作為第三個特征值，計算如下：

其中zv表示區域V中的體素v的灰度值，區域V用一個邊長為3體素的立方體區域表示，體素區域的大小選擇同鄰域強度，本文選取26-鄰域的體素區域用以計算統計矩，如圖4所示。m為該區域的平均灰度值，其中n≥2，n取2時為二階矩，即方差，度量的是圖像邊界偏離體素區域均值的程度，三階矩描述的是圖像紋理灰度的起伏分布，四階矩描述紋理灰度的反差，五階矩和更高矩不容易和圖像的直方圖連續起來。本文取二階矩作為第三個圖像特征值。

2.2 有監督描述子學習算法

有監督描述子學習算法［15］（SDL）是一種特征表示學習法，該方法能去除不相關和冗余的特征值，并將高維的特征轉換為較低維的特征表示，最終得到一個緊湊且高區別度的特征表示。

對于給定的特征向量{X1，X2，…，XL}，和對應的標簽值{Y1，Y2，…，Yn}，L表示特征向量的個數。目標是找到變換矩陣W∈RM×m，V∈RN×n，使得向量Di∈Rm×n經過變換WDiVT近似于特征向量Xi。即相當于求解如下最優問題：

為了提高近似特征值的精度，引入如下有監督的流形正則化項：

其中Sij為一核函數，如下式所示：

結合式（1）和式（2）可得：

β∈(0，∞)為一常量。

為了求得上述結果，采取迭代優化的方法，經優化Di=WTXiV，結合矩陣的秩，式（3）轉化為求解如下最大值：

對于任意給定一個向量V，可根據下式求得最優的向量W：

其中

同樣的，給定一個向量W，可根據下式求得最優的向量V：

其中

通過單一值分解方法（SVD）可求得式（5）和式（7）的解。

本文首先利用有監督流形正則化算法（SMR）［16］生成所得特征矩陣的低維近似矩陣描述，再利用有監督描述子學習算法訓練得到的低維特征，除去冗余和不相關的特征值，最后得到一個緊湊的區別度高的特征描述。

3 實驗結果及分析

3.1 實驗準備

在實驗中，本文共選取8組臨床MR圖像進行實驗，并在Matlab R2013a的編程環境下進行實驗，采用LIBSVM［17］工具箱進行編程訓練。然后將SDL特征壓縮后得到的分割結果和特征壓縮前得到的分割結果進行對比，并與數學形態學方法、線性判別式分析法（LDA）、主成分分析法（PCA）分割結果進行對比，分別計算其Matthews相關系數（MCC），Dice相關系數（Dice），真陽性率（Pt），假陽性率（Pf）。

評價參數簡介：

1）Matthews相關系數（MCC）。定義如下：

其中P11，P10，P01，P00由如下混淆矩陣定義給出：

由式（9）可知，MCC取值范圍在［-1，+1］，值越高，代表兩幅圖像的相似程度越高。

2）Dice相關系數（D）。定義如下：

其中，SMR表示用本文方法從MR圖像分割得到的顱骨圖，SCT表示用閾值分割方法從腦CT圖像中分割得到的顱骨圖。

3）分割誤差（Pt和Pf）定義為：。其中，Nvp為標準圖像（閾值分割結果）中標記為分割目標的像素個數和；Nuvp為標準圖像中被標記為非分割目標的像素個數總和；TN為本文分割結果中正確分割的像素數和，FN為本文分割結果中錯誤分割的像素數和。

3.2 實驗結果分析

特征的有效性決定了分割的精確度，為驗證所提取特征的有效性，對本文提出的三類特征進行不同組合并利用SVM進行訓練，得到如圖5所示的分割準確率與不同特征組合的關系圖。其中1代表概率圖譜值，2代表鄰域強度值，3代表統計矩，1+2表示概率圖譜與領域強度值的組合，以此類推。從圖中可以看出，概率圖譜特征值在識別顱骨區域和非顱骨區域上具有很強的代表性，其與其他特征值的不同組合均顯示很高的分割準確率。而僅使用領域強度或統計矩或兩者的組合作為特征，識別率并不高，不能達到區分顱骨和非顱骨的目的。本文選擇三種特征的組合確保了分割的準確性。

為了降低冗余特征對分類器的影響并獲取合適的特征維數，首先利用SDL方法對得到的21維特征矩陣進行不同程度的壓縮，隨后對另兩種壓縮方法PCA和LDA進行實驗，得到支持向量機的分類準確率隨壓縮特征維數的變化關系圖。從圖6可看出，在本文方法中，SDL比LDA和PCA取得的分割準確率要高，并且可看出當特征維數壓縮到8維至14維可以取得較高的分割準確率。

圖5不同特征組合與分割準確率的關系圖

圖6SVM分類準確率與特征數目的關系

隨后分別與特征壓縮前的SVM分割結果和運用數學形態學方法得到的分割結果進行對比，得到五種分割方法的數據對比，本文選擇將特征矩陣壓縮至9維進行對比實驗。在以上Matthews相關系數（MCC）、Dice相關系數（D）、分割誤差（Pt和Pf）中，MCC、Dice和Pt的數值越高，表示分割結果越好，兩幅圖像的相似程度越高，而Pf數值越低表示分割結果越好，相似程度越高。

圖7為分割結果對比圖，均選取了同一層面進行對比分析，圖中灰色區域標識分割結果與CT閾值分割結果的差異之處，從圖中可以看出，特征壓縮前，分割的誤差主要集中在大腦兩側區域和后腦區域，經LDA和PCA特征壓縮后的分割結果改善不明顯，而經過SDL特征壓縮后再進行分割，雖然后腦區域仍存在些許誤差，但較好地解決了大腦兩側的分割問題，進一步提高了分割精度。而從數學形態學的分割結果看，其分割結果明顯劣于特征壓縮前的分割結果和特征壓縮后的分割結果。從下圖中也可以明顯地看出數學形態學分割結果與CT閾值分割存在較大的誤差。在大腦的各個部位均不能很好的達到分割目的。

圖7分割結果對比圖，圖中均為第72層圖像

為了定量分析五種分割方法的區別，利用直方圖的形式將五種分割方法的評價參數進行對比分析。如下列直方圖所示，實驗中SDL特征壓縮后MCC相關系數平均值達到0.788，Dice系數均值為0.832，均高于特征壓縮前的0.741、0.798，LDA方法的0.69、0.74，PCA方法的0.72、0.73，數學形態學方法的0.478和0.482，Pt和Pf進行對比也顯示了SDL特征壓縮后的分割優勢。

圖8五種分割方法MCC對比

圖8為五種分割方法的MCC結果對比直方圖，從圖中可以明顯看出經過SDL特征壓縮后取得的分割結果要高于特征壓縮前的結果（即原始數據的分割結果），同時高于LDA和PCA方法的結果，且遠高于數學形態學方法分割結果，其中LDA和PCA方法分割精度相似。類似于MCC結果，五種分割方法Dice系數直方圖如圖9所示，從圖中可知數學形態學分割結果較差，Dice均值只有0.482不足0.50，LDA和PCA分割結果次之，分別為0.74和0.73，壓縮前的分割結果Dice均值為0.798，SDL特征壓縮后的分割結果效果最好，Dice均值達到了0.832。

圖9五種分割方法Dice對比

圖10為五種分割方法的Pt對比直方圖，Pt代表五種方法的正確分割顱骨數與標準圖像中標記為顱骨的數目之比，數值越高代表分割精度越好，兩幅圖像相似度越高。從圖中可看出采用支持向量機SVM的方法取得的分割效果要優于數學形態學的分割方法。

圖10五種分割方法Pt對比

同圖10，圖11為五種分割方法的Pf對比直方圖，Pf代表三種方法的錯誤分割顱骨數與標準圖像中標記為非顱骨的數目之比，數值越高代表分割精度越低，兩幅圖像相似度越低。如圖所示，數學形態學方法的Pf值要高于其他四種方法，而經SDL特征壓縮后的Pf值最小，代表分割結果最好。

圖11五種分割方法Pf對比

綜上所述，本文提出的算法成功地從腦MR圖像中自動分割出顱骨，并取得了較高的精確度，較好地吻合了真實的大腦頭模型，為解決MEG/EEG的源成像問題提供了方法指導。

4 結語

本文提出了一種基于支持向量機和有監督描述子學習算法的自動腦MR圖像顱骨分割方法，創新之處在于提出多特征值以區別顱骨區域和非顱骨區域，在確保特征值的有效性基礎上將支持向量機SVM方法與有監督描述子學習算法結合在一起，去除了冗余特征對分類器的影響，進一步提高了分割精確度。首先提取待分割MR圖像的全局特征（概率圖譜值）、局部特征（鄰域強度和統計矩），構建出特征矩陣，并利用有監督描述子學習算法進行特征壓縮，去除冗余和不相關的特征值，得到降維后的特征矩陣，在本文中特征壓縮前得到的特征矩陣是21維的，經有監督描述子學習算法壓縮后的特征矩陣維度為9，極大地降低了特征矩陣的維度，最后利用支持向量機SVM對21維的特征矩陣和9維的特征矩陣進行訓練，最終實現自動分割腦MR圖像顱骨的目的。之后分別將SDL特征壓縮后的分割結果與特征壓縮前的分割結果和LDA、PCA、數學形態學的分割結果進行對比，分別求取評價參數Matthews相關系數（MCC）、Dice相關系數（D）和分割誤差（Pt和Pf）進行對比，并利用柱狀圖進行數據對比分析，從柱狀圖中可以明顯看出采用支持向量機的訓練方法進行腦MR顱骨分割，其結果要優于數學形態學方法的分割結果。而經有監督描述子學習算法進行特征壓縮后取得的分割準確率也要高于特征壓縮前的分割準確率，這充分證明了本文分割算法的有效性，并且體現了SVM結合有監督描述子學習算法的優勢。另外，由于本文方法依賴于圖像的配準技術，需要將病人的MR圖像與參考MR進行配準，也需要將病人的CT圖像配準到參考CT圖像，造成程序存在時效性問題，因此下一步工作可集中于研究高效的配準算法，進一步提高本文算法的時效性。

［1］Potthast R，Kühn L.On the convergence of the finite inte?gration technique for the anisotropic boundary value prob?lem of magnetic tomography［J］.Mathematical methods inthe applied sciences，2003，26（9）：739-57.

［2］Cuffin BN.A method for localizing EEG sources in realis?tic head models［J］.Biomedical Engineering，IEEE Trans?actions on，1995，42（1）：68-71.

［3］劉金山，陳鎮坤，黃來華.基于貝葉斯分層混合模型的大腦FMRI圖像分割［J］.數理統計與管理，2015，34（4）：603-611. LIO Jinshan，CHEN Zhenkun，HUANG Laihua.Sege?mentation of Brain FMRI Based on Bayesian Hierarchical Mixture Model［J］.Journal of Applied Statistics and Man?agement，2015，34（4）：603-611.

［4］劉凡，高上凱，高小榕.基于最大后驗概率的大腦磁共振圖像的三維分割［J］.清華大學學報：自然科學版，2002（S1）：43-48. LIU Fan，GAO Shangkai，GAO Xiaorong.3-D Segmenta?tion of Brain Magnetic Resonance Images Based on Maxi?mum a Posterior［J］.Journal of Tsinghua University（Sci& Tech），2002（S1）：43-48.

［5］Heinonen T，Eskola H，Dastidar P，Laarne P，Malmivuo J.Segmentation of T1 MR scans for reconstruction of resis?tive head models［J］.Computer methods and programs in biomedicine，1997，54（3）：173-81.

［6］Rifai H，Bloch I，Hutchinson S，Wiart J，Garnero L.Seg?mentation of the skull in MRI volumes using deformable model and taking the partial volume effect into account［J］.Medical image analysis，2000，4（3）：219-33.

［7］Chu C，Takaya K，editors.3-dimensional rendering of MR（magnetic resonance）images［C］//WESCANEX 93‘Com?munications，Computers and Power in the Modern Envi?ronment'Conference Proceedings，IEEE，1993：IEEE.

［8］Dogdas B，Shattuck DW，Leahy RM.Segmentation of skull and scalp in 3-D human MRI using mathematical morphol?ogy［J］.Human brain mapping，2005，26（4）：273-85.

［9］Golub GH，Van Loan CF.Matrix computations［M］.JHU Press，2012.

［10］Dhillon IS，Modha DS.Concept decompositions for large sparse text data using clustering［J］.Machine learning，2001，42（1-2）：143-75.

［11］Srebro N，Jaakkola T，editors.Weighted low-rank ap?proximations［M］.ICML，2003.

［12］Thomasian A，Castelli V，Li C-S，editors.Clustering and singular value decomposition for approximate index?ing in high dimensional spaces［C］//Proceedings of the seventh international conference on Information and knowledge management，1998：ACM.

［13］Achlioptas D，McSherry F.Fast computation of low-rank matrix approximations［J］.Journal of the ACM（JACM），2007，54（2）：9.

［14］Shi G，Liu D，Gao D，Liu Z，Lin J，Wang L.Advances in theory and application of compressed sensing［J］.Acta Electronica Sinica，2009，37（5）：1070-81.

［15］Zhen X，Wang Z，Yu M，Li S，editors.Supervised de?scriptor learning for multi-output regression［C］//Pro?ceedings of the IEEE Conference on Computer Vision and Pattern Recognition，2015.

［16］Ye J.Generalized low rank approximations of matrices［J］.Machine Learning，2005，61（1-3）：167-91.

［17］Chang C-C，Lin C-J.LIBSVM：a library for support vec?tor machines［J］.ACM Transactions on Intelligent Sys?tems and Technology（TIST），2011，2（3）：27.

Automated Segmentation Based on Support Vector Machine and Supervised Descriptor Learning from Brain MR Image

HUANG Yongqi1，2SHI Wenbo3ZHOU Zhiyong1PANG Shumao4TONG Baotong1ZHAO Lingxiao1DAI Yakang1
（1.Suzhou Institute of Biomedical Engineering and Technology of Chinese Academy of Sciences，Suzhou215263）（2.University of Chinese Academy of Sciences，Beijing100049）（3.Beijing Normal University，Beijing100875）（4.Southern Medical University，Guangzhou510515）

A solution of EEG/MEG forward problem is essential and important in stereotactic neurosurgery applications.It is necessary to build a multi-layer brain model to distinguish different tissues for MEG/EEG forward problem.Although soft tissues can be clearly seen in MR images，but the intensity of skull is so low because of a lack of hydrogen in skull that can't be segmented auto?matically and accurately from MR image.Extracting skull form MR image automatically end up to be a key problem when calculating the MEG/EEG forward problem.In order to solve the above problem，a support vector machine（SVM）is proposed based segmentation algorithm using global features and local features of MR image.Moreover，the supervised descriptor learning（SDL）algorithm is com?bined that can transform the feature matrix into a compact one，and finally the skull from brain MR image is extrated by training on multi-modal images from the same patient whose CTs and MRs are available.Compared to the algorithm based on SVM only and math?ematical morphology based algorithm，the proposed method shows a considerable improvement on segmentation accuracy.The pro?posed method achieves an accuracy with Dice coefficient 0.832 compared with the other two methods 0.798 and 0.482.The proposed hybrid algorithm extract the skull successfully，so that the EEG，MEG source imaging problem can be solved easily in future work.

kull segmentation，support vector machine（SVM），supervised descriptor learning（SDL），feature extraction，feature compression

TP391

10.3969/j.issn.1672-9722.2017.07.034

2017年1月13日，

2017年2月28日

中國科學院百人計劃項目；國家自然科學基金（編號：61301042）；國家863計劃（編號：2015AA020514）；國家自然科學基金青年基金項目（編號：61501452）；江蘇省博士后基金項目（編號：1501089C）資助。

黃勇其，男，碩士研究生，研究方向：智能醫學圖像處理。史文博，男，研究方向：圖像處理。周志勇，男，博士，研究方向：醫學圖像處理與分析算法。龐樹茂，男，碩士研究生，研究方向：醫學圖像分析，醫學大數據。佟寶同，男，碩士，研究方向：醫學圖像處理軟件開發。趙凌霄，男，博士，研究方向：醫學影像分析處理與可視化算法研究和軟件開發。戴亞康，男，研究方向：醫學圖像處理算法研究和軟件開發。