基于自適應區域增長的fMRI腦功能激活區檢測

李 敏 曾衛明

(上海海事大學信息工程學院 上海 201306)

基于自適應區域增長的fMRI腦功能激活區檢測

李 敏 曾衛明

(上海海事大學信息工程學院 上海 201306)

區域增長作為一種有效的功能磁共振(fMRI)分析方法，由于受到諸如噪聲、生長準則等因素的影響，限制了它在腦功能激活區檢測方面的應用與發展。為了克服這些問題，提出一種自適應的區域增長方法，用于fMRI腦功能激活區的檢測。該方法首先利用主成分分析(PCA)對預處理后的圖像進行降噪；然后通過裂分合并與模板匹配相結合的方法自動獲取初始生長點；最后使用典型相關系數和皮爾森相關系數作為生長準則進行區域增長。通過模擬數據實驗和真實數據實驗驗證了該方法的有效性。實驗結果表明，相比于其他的fMRI分析方法(如ICA、SPM和基于裂分合并的區域增長方法)，自適應區域增長方法能夠獲得更加準確有效的結果。而且該方法還能應用到靜息態數據的分析中，進一步證明了該方法的可行性與有效性。綜上，自適應區域增長方法能夠拓寬區域增長在功能磁共振數據分析中的應用。

區域增長 裂分合并 皮爾森相關系數 典型相關系數 功能磁共振成像

0 引 言

區域增長是一種傳統的圖像分割方法，其基本思想是將具有相似性質的像素合并構成區域[1]。在真實的功能磁共振數據中，神經元的活動往往是以團簇的形式出現的，也就是說，真正的激活區更有可能是一些空間相連的體素的集合，而不是單個的體素。2003年，區域增長方法就已經運用在功能磁共振數據分析中，并且在此基礎上又提出了一種基于裂分合并的區域增長方法[2]。

與經典的激活區檢測方法如一般線性模型(GLM)和聚類分析相比，區域增長方法充分利用了數據的空間域信息[3]，因此得到的激活區的連續性更好，同時也不需要線性假設和分類個數的先驗條件[1]。但仍然有許多重要因素制約著區域增長方法在功能磁共振數據中的應用，比如：僅依賴皮爾森相關的生長準則對噪聲十分敏感，影響了最終結果的準確性與可靠性；由于靜息態數據沒有已知任務刺激的緣故，至今仍未使用區域增長方法對其進行分析。

本文在改進和完善傳統區域增長方法的基礎上，提出了一種新的自適應區域增長方法，并用它來進行功能磁共振數據方面的分析。該方法首先對預處理后的數據進行PCA降噪，然后通過裂分合并與模板匹配相結合的方式自動得到初始生長點。在隨后的區域增長過程中，選擇典型相關與皮爾森相關相結合的方法作為生長準則，以便能夠從fMRI中檢測出更多有意義的腦激活區。

1 方 法

1.1 PCA 降噪

由于噪聲的孤立性、離散性及大小不一，使其對激活區檢測有很強的敏感性，從而噪聲成為功能磁共振數據分析中最嚴重的影響因素之一[4]。因此，為了提高區域增長算法的精確性與可靠性，去噪是必不可少的預處理之一。

主成分分析(PCA)最早由Karl Pearson提出，該方法的主要思想是從數據集中提取主要信息來表示原數據，這些主要信息是一些正交的變量稱為主成分[5-6]。因此，如果對功能磁共振數據進行主成分提取，那么作為次要成分的噪聲在數據重構時就可以被有效地剔除，從而實現了去噪的效果。

主成分分析主要依靠特征分解，本文中主要分為以下四個步驟：

(1) 計算原始數據的協方差矩陣；

(2) 提取協方差矩陣的特征向量；

(3) 將特征向量的轉置矩陣乘以原始數據得到壓縮矩陣；

(4) 將特征向量乘以壓縮矩陣得到重構后的數據。

1.2 自適應初始生長點選擇

自適應初始生長點選擇分為兩個部分：裂分合并選擇初始生長區域以及模板匹配選擇初始生長點。首先，裂分合并方法先大致獲得一些初始區域，然后通過一個已知的模板來進一步從中選出初始生長點。

1) 裂分合并初始生長區域選擇

裂分合并是最早被用作圖像分割的一種方法，該方法把整幅圖像分割成一系列不相鄰的區域，主要步驟[7]如下：

(1) 令R表示整幅圖像區域并選擇一個謂詞P作為生長準則。

(2) 對于每一個區域Ri，如果P(Ri)=false，就將每個區域拆分為8個相鄰的子區域。

(3) 將P(Ri∪Rj)=true的任意兩個相鄰區域Ri和Rj進行合并。

(4) 重復步驟(2)-(3)直到無法再進行裂分合并。

fMRI數據F={fxyzt}MNOT是一個時空的四維信號，其中(x,y,z)是像素的三維笛卡兒坐標，t是時間指標，M×N×O表示一次掃描中的像素的個數，T表示所有的掃描個數，fxyzt是像素圖像(x,y,z)在t時刻的灰度值[1]。

在裂分合并階段，整個fMRI數據被看作一個四維的塊，每個塊都有一個度量標準PSM。如果一個區域R不滿足度量標準就將該區域裂分成8個同樣大小的子區域[R1,R2,…,R8]，子區域Ri重復這一步驟直到滿足度量標準，或者達到最小區域面積[2,7]。如果相鄰子區域合并后滿足度量標準，就將它們合并成一個區域[8]。

由于靜息態fMRI數據沒有任務態fMRI數據中已知的任務刺激，因此兩種狀態下度量標準PSM的選擇是不同的。在任務狀態下，PSM的定義依據已知任務刺激和區域的平均時間序列的相關性，具體定義為：

(1)

其中，fxyz:是已知的任務刺激，MR是區域的平均時間序列，V和C分別表示方差和協方差，Tsm是閾值。

在靜息狀態下，PSM的定義依據該區域每個體素點時間序列的標準差，具體定義如下：

(2)

其中，fxyz:是區域R中每個體素點的時間序列，MR是區域R的平均時間序列，n是區域R中體素點的數量，Tsm是閾值。

2) 模板匹配

上述裂分合并方法能夠選擇出一些大致的初始生長區域。但對于靜息態fMRI數據，這些初始生長區域只是一些相關性較高的區域，仍然無法確定具體某個功能網絡的初始生長點，因此，通過模板匹配，即將落在模板區域內的初始生長區域作為最終的初始生長點[9]。例如：在靜息態fMRI默認網絡檢測的實驗中，使用默認網絡的模板匹配裂分合并的結果，即將屬于模板范圍內的裂分合并后的初始生長區域，作為最終初始生長點。

1.3 自適應區域增長

皮爾森相關系數是衡量兩個變量相關性強弱的重要標準，但由于它忽略了功能磁共振數據中存在的空間域信息的影響，因此皮爾森相關系數對噪聲很敏感，往往會造成最終檢測的激活區中出現偽激活區[10-11]。本文采用將典型相關分析與皮爾森相關分析相結合的方法作為生長準則。其中，皮爾森相關用于檢驗體素點之間的相關性，典型相關則充分利用空間信息進行檢測，兩種方式相結合可以有效地規避噪聲敏感的問題。

1) 皮爾森相關系數

Pearson相關系數是經典的比較兩個變量之間相關性的度量標準，具體定義如下：

(3)

其中，fxyz:是每個種子區域R鄰域中的體素點的時間序列，MR是區域R的平均時間序列，V和C分別表示方差和協方差。

2) 典型相關系數

典型相關分析是由Hotelling于1936年最先提出的，該方法的主要思想是從兩組數據中找到一組基向量，該組向量能夠最大程度地反應原始數據之間的相關性。典型相關分析常用于分析多維變量間相關性問題，在fMRI中，它充分利用了血液動力學的子空間模型和空間域信息進行腦區域一致性檢測[10,12]。典型相關系數的定義如下：

(4)

其中，f是種子區域R鄰域中的體素點的時間序列，MR是種子區域R的平均時間序列，Wf和WMR分別是f和mR的特征向量，Cff和CMRMR分別是f和MR方差，CfMR是協方差矩陣。

在本文中，生長準則的定義如下：

PRG=P{corr(fxyz:,MR)>TRG&maxρ>TRG}

(5)

其中，TRG是閾值。式(5)表示只有當鄰域中的體素點同時滿足式(3)和式(4)時，該體素點才能合并入生長區域。

2 實 驗

本文中的預處理過程和算法都是在Matlab平臺上實現的(Matlab, 2013a,Math-WorksInc.,Sherborn,MA,USA)。任務態和靜息態的fMRI數據的預處理都是由軟件包SPM12完成的，預處理過程包括：層間時間校正、頭動校正、空間標準化和平滑。

2.1 混合模擬實驗

在模擬實驗中，我們在靜息態數據中選取了7×7×7大小的立方體添加刺激，刺激強度范圍為0.02至0.05，步長為0.005。模擬實驗的結果如圖1所示，實驗的結果比較如圖2所示。

其中，(1)至(3)分別代表ICA,SPM,ARGM在同一切片的實驗結果，從左至右分別表示刺激強度范圍0.02至0.05，步長為0.005。

圖1 模擬實驗結果

圖2 模擬數據實驗結果比較

圖1展示了ARGM能夠成功檢測出添加的任務刺激，從圖2可以看出，ARGM方法可以減小偽激活率，并且與ICA和SPM方法相比也能夠獲得較好的結果。

2.2 任務態實驗

為了說明自適應區域增長方法的有效性，我們在真實fMRI視覺試驗的數據上對該方法進行了驗證。BOLD/EPI圖像來源于Siemens2T掃描系統，TR=2秒，圖像尺寸為64×64×64，像素大小為3×3×3毫米。刺激序列以休息開始并且刺激和休息交替出現，每個組塊20秒，總共選擇了70個掃描體數據，完整刺激序列為[(0,1)×3,0]。在刺激狀態時，將藍/黃棋盤保持旋轉頻率7Hz，作為被試者的視覺刺激，在休息狀態時，要求被試者集中在屏幕的中心十字上。

利用以上真實實驗數據對本文提出的自適應區域增長方法(ARGM)進行驗證，實驗結果證明了該方法的可行性。在本文的討論部分，我們將會對ARGM與ICA、SPM及裂分合并區域增長方法進行比較。實驗相關參數,分別為0.9、0.95，實驗結果如圖3所示。

圖3 任務態自適應區域增長結果

2.3 靜息態實驗

本文對靜息態閉眼試驗所得到的fMRI數據上對ARGM進行了驗證。BOLD/EPI圖像來源于Siemens2T掃描系統，TR=2秒，圖像尺寸為64×64×64，像素大小為3×3×3毫米，總共選擇了150個掃描體數據。

利用以上fMRI靜息態數據，ARGM對腦功能網絡：DMN、DPN2、WMN、AUN以及VIN進行檢測，實驗中參數Tsm、TRG分別為10、0.95，實驗結果如圖4所示。

圖4 靜息態自適應區域增長結果

(1)至(5)依次為DMN、DPN2、WMN、AUN和VIN。

靜息態實驗結果表明，ARGM可以成功地檢測出一些腦功能網絡，如DMN、DPN2、WMN、AUN以及VIN等。也就是說，盡管腦信號強度在靜息態下比較微弱，但通過自適應區域增長的方法仍然可以有效地進行激活區檢測，因此，可以利用區域增長方法對靜息態下的功能磁共振數據進行分析。

3 討 論

為了驗證在任務態下自適應區域增長算法(ARGM)的可靠性，我們將其與經典的激活區檢測方法ICA、SPM和裂分合并區域增長算法(SMRG)進行比較，其中，ICA的order參數為40，SPM使用的是SPM12，實驗結果如圖5所示。

圖5 實驗結果比較

圖5(1)-(3)曲線圖分別顯示了ICA、SPM、SMRG與ARGM的結果對比。其中橫軸表示體素點與已知任務刺激的相關性，縱軸表示體素點的數量，方形線表示ARGM相比其他方法多檢測到的體素點的數量，菱形線表示其他方法相比ARGM多檢測到的體素點的數量。由于ICA、SPM、SMRG檢測出的激活區體素點的數量各不相同，因此ARGM相比每種方法多檢測出的體素點的數量也不相同，三種方法相比ARGM多檢測出的數量也不相同，所以圖5(1)-(3)中縱坐標的值不相同。

與ICA、SPM兩種方法相比，ARGM能夠檢測出更多與已知任務刺激相關性高的體素點，(圖5的3個子圖中ARGM多檢測出的體素點主要集中于相關性0.6以上)，比如ARGM比ICA多檢測出相關性0.7至0.8的體素點612個，比SPM多檢測出相關性0.8至0.9的體素點407個。相比之下，ICA相比ARGM多檢測出的體素點主要集中在0.4至0.5之間，SPM相比ARGM多檢測出的主要集中在0.2至0.5之間。此外，從圖5(3)中可以看出，ARGM相比SMRG檢測出了更多相關性0.6以上的體素點，大大改進了現有的區域增長方法。

4 結 語

在本文中，我們提出了一種自適應的區域增長功能磁共振激活區檢測方法ARGM。該方法主要對噪聲、初始生長點的選擇以及生長準則三個方面做了改進，不僅降低了噪聲對激活區檢測的影響，還將區域增長方法應用到靜息態數據分析中，而且檢測出更多更準確的激活區。在模擬數據實驗中，自適應區域增長方法成功檢測出了激活區；在任務態數據實驗中，與幾種經典的fMRI分析方法相比，自適應區域增長方法能夠得到更多真實的激活區；在靜息態數據中，自適應區域增長方法也成功地檢測出了腦功能網絡，說明該方法在腦信號較弱的靜息態下也仍然適用。因此，自適應區域增長方法ARGM是一種可靠且有效的fMRI激活區檢測方法。

[1] 呂英立.腦功能核磁共振激活區檢測算法研究[D].北京:中國科學院自動化研究所,2004.

[2]LuYL,JiangTZ,ZangYF.ASplit-Merge-BasedRegion-GrowingMethodforfMRIActivationDetection[J].HumanBrainMapping,2004,22(4):271-279.

[3] 肖飛,綦星光.圖像分割方法綜述[J].可編程控制器與工廠自動化,2009，31(11):77-79.

[4]LiJX,JiangSP.AdaptiveThresholdImageDenoisingAlgorithmBasedonPrincipalComponentAnalysis[J].InfraredTechnology,2014,36(4):312-315.

[5] 芮挺,王金巖,沈春林,等.基于PCA的圖像小波去噪方法[J].小型微型計算機系統,2006,27(1):158-161.

[6] 李俊秀,姜三平.基于主成分分析的圖像自適應閾值去噪算法[J].紅外技術,2014,36(4):311-315.

[7]RafaelCGonzalez,RechardEWoods.數字圖像處理[M].2版.北京電子工業出版社,2011:498-499.

[8] 胡亞斌,劉建平.基于區域增長和改進分水嶺算法的胃腺癌細胞圖像分割[J].計算機與數字工程,2006,34(8):151-153.

[9]VanhaudenhuyseA,NoirhommeQ,TshibandaLJ,etal.DefaultNetworkConnectivityReflectstheLevelofConsciousnessinNon-communicativeBrainDamagedPatients[J].Brain,2010,133(1):161-171.

[10]FrimanO,CedefamnJ,LundbergP,etal.DetectionofNeuralActivityinFunctionalMRIUsingCanonicalCorrelationAnalysis[J].MagneticResonanceinMedicine,2001,45(2):323-330.

[11]DeleusF,VanHulleMM.FunctionalConnectivityAnalysisoffMRIDataBasedonRegularizedMultisetCanonicalCorrelationAnalysis[J].JournalofNeuroscienceMethods,2011,197(1):143-157.

[12] 孫權森,曾生根,王平安,等.典型相關分析的理論及其在特征融合中的應用[J].計算機學報,2005,28(9):1524-1533.

ACTIVE REGION DETECTION OF BRAIN FUNCTION BY FMRI BASEDON SELF-ADAPTIVE REGION GROWING

Li Min Zeng Weiming

(CollegeofInformationEngineering,ShanghaiMaritimeUniversity,Shanghai201306,China)

Region growing has been utilized in the analysis of functional magnetic resonance imaging (fMRI) data for many years, while some influential factors, such as the noise problem, definition of the homogeneity criterion, restricting the application and development in brain function active region detection. In order to overcome these disadvantages, an adaptive region growing method (ARGM) is proposed to detect fMRI brain function active region, where PCA was firstly used to de-noise the fMRI data as a step of preprocessing. Then the region seed was automatically selected by the split-merge method combined with a prior template. Next, an improved homogeneity criterion defined by Canonical correlation coefficient and Pearson correlation coefficient were used to judge the region growing. Compared with the typical fMRI data analysis methods such as ICA and SPM and the classical region growing method, ARGM generates a more accurate and reliable result in task-related experiment. In addition, the resting-state experiment has also demonstrated the effectiveness and usefulness of the proposed method. To conclude, the proposed method is able to broaden the application of region growing in analyzing fMRI data.

Region growing Split-merge Pearson correlation coefficient Canonical correlation coefficient Functional magnetic resonance imaging

2015-12-11。上海科委重點項目(14590501700)。李敏，碩士生，主研領域：區域增長，功能磁共振成像。曾衛明，教授。

TP3

A

10.3969/j.issn.1000-386x.2017.03.030