李 敏 曾衛(wèi)明
(上海海事大學(xué)信息工程學(xué)院 上海 201306)
基于自適應(yīng)區(qū)域增長的fMRI腦功能激活區(qū)檢測
李 敏 曾衛(wèi)明
(上海海事大學(xué)信息工程學(xué)院 上海 201306)
區(qū)域增長作為一種有效的功能磁共振(fMRI)分析方法,由于受到諸如噪聲、生長準(zhǔn)則等因素的影響,限制了它在腦功能激活區(qū)檢測方面的應(yīng)用與發(fā)展。為了克服這些問題,提出一種自適應(yīng)的區(qū)域增長方法,用于fMRI腦功能激活區(qū)的檢測。該方法首先利用主成分分析(PCA)對預(yù)處理后的圖像進(jìn)行降噪;然后通過裂分合并與模板匹配相結(jié)合的方法自動(dòng)獲取初始生長點(diǎn);最后使用典型相關(guān)系數(shù)和皮爾森相關(guān)系數(shù)作為生長準(zhǔn)則進(jìn)行區(qū)域增長。通過模擬數(shù)據(jù)實(shí)驗(yàn)和真實(shí)數(shù)據(jù)實(shí)驗(yàn)驗(yàn)證了該方法的有效性。實(shí)驗(yàn)結(jié)果表明,相比于其他的fMRI分析方法(如ICA、SPM和基于裂分合并的區(qū)域增長方法),自適應(yīng)區(qū)域增長方法能夠獲得更加準(zhǔn)確有效的結(jié)果。而且該方法還能應(yīng)用到靜息態(tài)數(shù)據(jù)的分析中,進(jìn)一步證明了該方法的可行性與有效性。綜上,自適應(yīng)區(qū)域增長方法能夠拓寬區(qū)域增長在功能磁共振數(shù)據(jù)分析中的應(yīng)用。
區(qū)域增長 裂分合并 皮爾森相關(guān)系數(shù) 典型相關(guān)系數(shù) 功能磁共振成像
區(qū)域增長是一種傳統(tǒng)的圖像分割方法,其基本思想是將具有相似性質(zhì)的像素合并構(gòu)成區(qū)域[1]。在真實(shí)的功能磁共振數(shù)據(jù)中,神經(jīng)元的活動(dòng)往往是以團(tuán)簇的形式出現(xiàn)的,也就是說,真正的激活區(qū)更有可能是一些空間相連的體素的集合,而不是單個(gè)的體素。2003年,區(qū)域增長方法就已經(jīng)運(yùn)用在功能磁共振數(shù)據(jù)分析中,并且在此基礎(chǔ)上又提出了一種基于裂分合并的區(qū)域增長方法[2]。
與經(jīng)典的激活區(qū)檢測方法如一般線性模型(GLM)和聚類分析相比,區(qū)域增長方法充分利用了數(shù)據(jù)的空間域信息[3],因此得到的激活區(qū)的連續(xù)性更好,同時(shí)也不需要線性假設(shè)和分類個(gè)數(shù)的先驗(yàn)條件[1]。但仍然有許多重要因素制約著區(qū)域增長方法在功能磁共振數(shù)據(jù)中的應(yīng)用,比如:僅依賴皮爾森相關(guān)的生長準(zhǔn)則對噪聲十分敏感,影響了最終結(jié)果的準(zhǔn)確性與可靠性;由于靜息態(tài)數(shù)據(jù)沒有已知任務(wù)刺激的緣故,至今仍未使用區(qū)域增長方法對其進(jìn)行分析。
本文在改進(jìn)和完善傳統(tǒng)區(qū)域增長方法的基礎(chǔ)上,提出了一種新的自適應(yīng)區(qū)域增長方法,并用它來進(jìn)行功能磁共振數(shù)據(jù)方面的分析。該方法首先對預(yù)處理后的數(shù)據(jù)進(jìn)行PCA降噪,然后通過裂分合并與模板匹配相結(jié)合的方式自動(dòng)得到初始生長點(diǎn)。在隨后的區(qū)域增長過程中,選擇典型相關(guān)與皮爾森相關(guān)相結(jié)合的方法作為生長準(zhǔn)則,以便能夠從fMRI中檢測出更多有意義的腦激活區(qū)。
1.1 PCA 降噪
由于噪聲的孤立性、離散性及大小不一,使其對激活區(qū)檢測有很強(qiáng)的敏感性,從而噪聲成為功能磁共振數(shù)據(jù)分析中最嚴(yán)重的影響因素之一[4]。因此,為了提高區(qū)域增長算法的精確性與可靠性,去噪是必不可少的預(yù)處理之一。
主成分分析(PCA)最早由Karl Pearson提出,該方法的主要思想是從數(shù)據(jù)集中提取主要信息來表示原數(shù)據(jù),這些主要信息是一些正交的變量稱為主成分[5-6]。因此,如果對功能磁共振數(shù)據(jù)進(jìn)行主成分提取,那么作為次要成分的噪聲在數(shù)據(jù)重構(gòu)時(shí)就可以被有效地剔除,從而實(shí)現(xiàn)了去噪的效果。
主成分分析主要依靠特征分解,本文中主要分為以下四個(gè)步驟:
(1) 計(jì)算原始數(shù)據(jù)的協(xié)方差矩陣;
(2) 提取協(xié)方差矩陣的特征向量;
(3) 將特征向量的轉(zhuǎn)置矩陣乘以原始數(shù)據(jù)得到壓縮矩陣;
(4) 將特征向量乘以壓縮矩陣得到重構(gòu)后的數(shù)據(jù)。
1.2 自適應(yīng)初始生長點(diǎn)選擇
自適應(yīng)初始生長點(diǎn)選擇分為兩個(gè)部分:裂分合并選擇初始生長區(qū)域以及模板匹配選擇初始生長點(diǎn)。首先,裂分合并方法先大致獲得一些初始區(qū)域,然后通過一個(gè)已知的模板來進(jìn)一步從中選出初始生長點(diǎn)。
1) 裂分合并初始生長區(qū)域選擇
裂分合并是最早被用作圖像分割的一種方法,該方法把整幅圖像分割成一系列不相鄰的區(qū)域,主要步驟[7]如下:
(1) 令R表示整幅圖像區(qū)域并選擇一個(gè)謂詞P作為生長準(zhǔn)則。
(2) 對于每一個(gè)區(qū)域Ri,如果P(Ri)=false,就將每個(gè)區(qū)域拆分為8個(gè)相鄰的子區(qū)域。
(3) 將P(Ri∪Rj)=true的任意兩個(gè)相鄰區(qū)域Ri和Rj進(jìn)行合并。
(4) 重復(fù)步驟(2)-(3)直到無法再進(jìn)行裂分合并。
fMRI數(shù)據(jù)F={fxyzt}MNOT是一個(gè)時(shí)空的四維信號,其中(x,y,z)是像素的三維笛卡兒坐標(biāo),t是時(shí)間指標(biāo),M×N×O表示一次掃描中的像素的個(gè)數(shù),T表示所有的掃描個(gè)數(shù),fxyzt是像素圖像(x,y,z)在t時(shí)刻的灰度值[1]。
在裂分合并階段,整個(gè)fMRI數(shù)據(jù)被看作一個(gè)四維的塊,每個(gè)塊都有一個(gè)度量標(biāo)準(zhǔn)PSM。如果一個(gè)區(qū)域R不滿足度量標(biāo)準(zhǔn)就將該區(qū)域裂分成8個(gè)同樣大小的子區(qū)域[R1,R2,…,R8],子區(qū)域Ri重復(fù)這一步驟直到滿足度量標(biāo)準(zhǔn),或者達(dá)到最小區(qū)域面積[2,7]。如果相鄰子區(qū)域合并后滿足度量標(biāo)準(zhǔn),就將它們合并成一個(gè)區(qū)域[8]。
由于靜息態(tài)fMRI數(shù)據(jù)沒有任務(wù)態(tài)fMRI數(shù)據(jù)中已知的任務(wù)刺激,因此兩種狀態(tài)下度量標(biāo)準(zhǔn)PSM的選擇是不同的。在任務(wù)狀態(tài)下,PSM的定義依據(jù)已知任務(wù)刺激和區(qū)域的平均時(shí)間序列的相關(guān)性,具體定義為:
(1)
其中,fxyz:是已知的任務(wù)刺激,MR是區(qū)域的平均時(shí)間序列,V和C分別表示方差和協(xié)方差,Tsm是閾值。
在靜息狀態(tài)下,PSM的定義依據(jù)該區(qū)域每個(gè)體素點(diǎn)時(shí)間序列的標(biāo)準(zhǔn)差,具體定義如下:
(2)
其中,fxyz:是區(qū)域R中每個(gè)體素點(diǎn)的時(shí)間序列,MR是區(qū)域R的平均時(shí)間序列,n是區(qū)域R中體素點(diǎn)的數(shù)量,Tsm是閾值。
2) 模板匹配
上述裂分合并方法能夠選擇出一些大致的初始生長區(qū)域。但對于靜息態(tài)fMRI數(shù)據(jù),這些初始生長區(qū)域只是一些相關(guān)性較高的區(qū)域,仍然無法確定具體某個(gè)功能網(wǎng)絡(luò)的初始生長點(diǎn),因此,通過模板匹配,即將落在模板區(qū)域內(nèi)的初始生長區(qū)域作為最終的初始生長點(diǎn)[9]。例如:在靜息態(tài)fMRI默認(rèn)網(wǎng)絡(luò)檢測的實(shí)驗(yàn)中,使用默認(rèn)網(wǎng)絡(luò)的模板匹配裂分合并的結(jié)果,即將屬于模板范圍內(nèi)的裂分合并后的初始生長區(qū)域,作為最終初始生長點(diǎn)。
1.3 自適應(yīng)區(qū)域增長
皮爾森相關(guān)系數(shù)是衡量兩個(gè)變量相關(guān)性強(qiáng)弱的重要標(biāo)準(zhǔn),但由于它忽略了功能磁共振數(shù)據(jù)中存在的空間域信息的影響,因此皮爾森相關(guān)系數(shù)對噪聲很敏感,往往會(huì)造成最終檢測的激活區(qū)中出現(xiàn)偽激活區(qū)[10-11]。本文采用將典型相關(guān)分析與皮爾森相關(guān)分析相結(jié)合的方法作為生長準(zhǔn)則。其中,皮爾森相關(guān)用于檢驗(yàn)體素點(diǎn)之間的相關(guān)性,典型相關(guān)則充分利用空間信息進(jìn)行檢測,兩種方式相結(jié)合可以有效地規(guī)避噪聲敏感的問題。
1) 皮爾森相關(guān)系數(shù)
Pearson相關(guān)系數(shù)是經(jīng)典的比較兩個(gè)變量之間相關(guān)性的度量標(biāo)準(zhǔn),具體定義如下:
(3)
其中,fxyz:是每個(gè)種子區(qū)域R鄰域中的體素點(diǎn)的時(shí)間序列,MR是區(qū)域R的平均時(shí)間序列,V和C分別表示方差和協(xié)方差。
2) 典型相關(guān)系數(shù)
典型相關(guān)分析是由Hotelling于1936年最先提出的,該方法的主要思想是從兩組數(shù)據(jù)中找到一組基向量,該組向量能夠最大程度地反應(yīng)原始數(shù)據(jù)之間的相關(guān)性。典型相關(guān)分析常用于分析多維變量間相關(guān)性問題,在fMRI中,它充分利用了血液動(dòng)力學(xué)的子空間模型和空間域信息進(jìn)行腦區(qū)域一致性檢測[10,12]。典型相關(guān)系數(shù)的定義如下:
(4)
其中,f是種子區(qū)域R鄰域中的體素點(diǎn)的時(shí)間序列,MR是種子區(qū)域R的平均時(shí)間序列,Wf和WMR分別是f和mR的特征向量,Cff和CMRMR分別是f和MR方差,CfMR是協(xié)方差矩陣。
在本文中,生長準(zhǔn)則的定義如下:
PRG=P{corr(fxyz:,MR)>TRG&maxρ>TRG}
(5)
其中,TRG是閾值。式(5)表示只有當(dāng)鄰域中的體素點(diǎn)同時(shí)滿足式(3)和式(4)時(shí),該體素點(diǎn)才能合并入生長區(qū)域。
本文中的預(yù)處理過程和算法都是在Matlab平臺上實(shí)現(xiàn)的(Matlab, 2013a,Math-WorksInc.,Sherborn,MA,USA)。任務(wù)態(tài)和靜息態(tài)的fMRI數(shù)據(jù)的預(yù)處理都是由軟件包SPM12完成的,預(yù)處理過程包括:層間時(shí)間校正、頭動(dòng)校正、空間標(biāo)準(zhǔn)化和平滑。
2.1 混合模擬實(shí)驗(yàn)
在模擬實(shí)驗(yàn)中,我們在靜息態(tài)數(shù)據(jù)中選取了7×7×7大小的立方體添加刺激,刺激強(qiáng)度范圍為0.02至0.05,步長為0.005。模擬實(shí)驗(yàn)的結(jié)果如圖1所示,實(shí)驗(yàn)的結(jié)果比較如圖2所示。
其中,(1)至(3)分別代表ICA,SPM,ARGM在同一切片的實(shí)驗(yàn)結(jié)果,從左至右分別表示刺激強(qiáng)度范圍0.02至0.05,步長為0.005。
圖1 模擬實(shí)驗(yàn)結(jié)果
圖2 模擬數(shù)據(jù)實(shí)驗(yàn)結(jié)果比較
圖1展示了ARGM能夠成功檢測出添加的任務(wù)刺激,從圖2可以看出,ARGM方法可以減小偽激活率,并且與ICA和SPM方法相比也能夠獲得較好的結(jié)果。
2.2 任務(wù)態(tài)實(shí)驗(yàn)
為了說明自適應(yīng)區(qū)域增長方法的有效性,我們在真實(shí)fMRI視覺試驗(yàn)的數(shù)據(jù)上對該方法進(jìn)行了驗(yàn)證。BOLD/EPI圖像來源于Siemens2T掃描系統(tǒng),TR=2秒,圖像尺寸為64×64×64,像素大小為3×3×3毫米。刺激序列以休息開始并且刺激和休息交替出現(xiàn),每個(gè)組塊20秒,總共選擇了70個(gè)掃描體數(shù)據(jù),完整刺激序列為[(0,1)×3,0]。在刺激狀態(tài)時(shí),將藍(lán)/黃棋盤保持旋轉(zhuǎn)頻率7Hz,作為被試者的視覺刺激,在休息狀態(tài)時(shí),要求被試者集中在屏幕的中心十字上。
利用以上真實(shí)實(shí)驗(yàn)數(shù)據(jù)對本文提出的自適應(yīng)區(qū)域增長方法(ARGM)進(jìn)行驗(yàn)證,實(shí)驗(yàn)結(jié)果證明了該方法的可行性。在本文的討論部分,我們將會(huì)對ARGM與ICA、SPM及裂分合并區(qū)域增長方法進(jìn)行比較。實(shí)驗(yàn)相關(guān)參數(shù),分別為0.9、0.95,實(shí)驗(yàn)結(jié)果如圖3所示。
圖3 任務(wù)態(tài)自適應(yīng)區(qū)域增長結(jié)果
2.3 靜息態(tài)實(shí)驗(yàn)
本文對靜息態(tài)閉眼試驗(yàn)所得到的fMRI數(shù)據(jù)上對ARGM進(jìn)行了驗(yàn)證。BOLD/EPI圖像來源于Siemens2T掃描系統(tǒng),TR=2秒,圖像尺寸為64×64×64,像素大小為3×3×3毫米,總共選擇了150個(gè)掃描體數(shù)據(jù)。
利用以上fMRI靜息態(tài)數(shù)據(jù),ARGM對腦功能網(wǎng)絡(luò):DMN、DPN2、WMN、AUN以及VIN進(jìn)行檢測,實(shí)驗(yàn)中參數(shù)Tsm、TRG分別為10、0.95,實(shí)驗(yàn)結(jié)果如圖4所示。
圖4 靜息態(tài)自適應(yīng)區(qū)域增長結(jié)果
(1)至(5)依次為DMN、DPN2、WMN、AUN和VIN。
靜息態(tài)實(shí)驗(yàn)結(jié)果表明,ARGM可以成功地檢測出一些腦功能網(wǎng)絡(luò),如DMN、DPN2、WMN、AUN以及VIN等。也就是說,盡管腦信號強(qiáng)度在靜息態(tài)下比較微弱,但通過自適應(yīng)區(qū)域增長的方法仍然可以有效地進(jìn)行激活區(qū)檢測,因此,可以利用區(qū)域增長方法對靜息態(tài)下的功能磁共振數(shù)據(jù)進(jìn)行分析。
為了驗(yàn)證在任務(wù)態(tài)下自適應(yīng)區(qū)域增長算法(ARGM)的可靠性,我們將其與經(jīng)典的激活區(qū)檢測方法ICA、SPM和裂分合并區(qū)域增長算法(SMRG)進(jìn)行比較,其中,ICA的order參數(shù)為40,SPM使用的是SPM12,實(shí)驗(yàn)結(jié)果如圖5所示。
圖5 實(shí)驗(yàn)結(jié)果比較
圖5(1)-(3)曲線圖分別顯示了ICA、SPM、SMRG與ARGM的結(jié)果對比。其中橫軸表示體素點(diǎn)與已知任務(wù)刺激的相關(guān)性,縱軸表示體素點(diǎn)的數(shù)量,方形線表示ARGM相比其他方法多檢測到的體素點(diǎn)的數(shù)量,菱形線表示其他方法相比ARGM多檢測到的體素點(diǎn)的數(shù)量。由于ICA、SPM、SMRG檢測出的激活區(qū)體素點(diǎn)的數(shù)量各不相同,因此ARGM相比每種方法多檢測出的體素點(diǎn)的數(shù)量也不相同,三種方法相比ARGM多檢測出的數(shù)量也不相同,所以圖5(1)-(3)中縱坐標(biāo)的值不相同。
與ICA、SPM兩種方法相比,ARGM能夠檢測出更多與已知任務(wù)刺激相關(guān)性高的體素點(diǎn),(圖5的3個(gè)子圖中ARGM多檢測出的體素點(diǎn)主要集中于相關(guān)性0.6以上),比如ARGM比ICA多檢測出相關(guān)性0.7至0.8的體素點(diǎn)612個(gè),比SPM多檢測出相關(guān)性0.8至0.9的體素點(diǎn)407個(gè)。相比之下,ICA相比ARGM多檢測出的體素點(diǎn)主要集中在0.4至0.5之間,SPM相比ARGM多檢測出的主要集中在0.2至0.5之間。此外,從圖5(3)中可以看出,ARGM相比SMRG檢測出了更多相關(guān)性0.6以上的體素點(diǎn),大大改進(jìn)了現(xiàn)有的區(qū)域增長方法。
在本文中,我們提出了一種自適應(yīng)的區(qū)域增長功能磁共振激活區(qū)檢測方法ARGM。該方法主要對噪聲、初始生長點(diǎn)的選擇以及生長準(zhǔn)則三個(gè)方面做了改進(jìn),不僅降低了噪聲對激活區(qū)檢測的影響,還將區(qū)域增長方法應(yīng)用到靜息態(tài)數(shù)據(jù)分析中,而且檢測出更多更準(zhǔn)確的激活區(qū)。在模擬數(shù)據(jù)實(shí)驗(yàn)中,自適應(yīng)區(qū)域增長方法成功檢測出了激活區(qū);在任務(wù)態(tài)數(shù)據(jù)實(shí)驗(yàn)中,與幾種經(jīng)典的fMRI分析方法相比,自適應(yīng)區(qū)域增長方法能夠得到更多真實(shí)的激活區(qū);在靜息態(tài)數(shù)據(jù)中,自適應(yīng)區(qū)域增長方法也成功地檢測出了腦功能網(wǎng)絡(luò),說明該方法在腦信號較弱的靜息態(tài)下也仍然適用。因此,自適應(yīng)區(qū)域增長方法ARGM是一種可靠且有效的fMRI激活區(qū)檢測方法。
[1] 呂英立.腦功能核磁共振激活區(qū)檢測算法研究[D].北京:中國科學(xué)院自動(dòng)化研究所,2004.
[2]LuYL,JiangTZ,ZangYF.ASplit-Merge-BasedRegion-GrowingMethodforfMRIActivationDetection[J].HumanBrainMapping,2004,22(4):271-279.
[3] 肖飛,綦星光.圖像分割方法綜述[J].可編程控制器與工廠自動(dòng)化,2009,31(11):77-79.
[4]LiJX,JiangSP.AdaptiveThresholdImageDenoisingAlgorithmBasedonPrincipalComponentAnalysis[J].InfraredTechnology,2014,36(4):312-315.
[5] 芮挺,王金巖,沈春林,等.基于PCA的圖像小波去噪方法[J].小型微型計(jì)算機(jī)系統(tǒng),2006,27(1):158-161.
[6] 李俊秀,姜三平.基于主成分分析的圖像自適應(yīng)閾值去噪算法[J].紅外技術(shù),2014,36(4):311-315.
[7]RafaelCGonzalez,RechardEWoods.數(shù)字圖像處理[M].2版.北京電子工業(yè)出版社,2011:498-499.
[8] 胡亞斌,劉建平.基于區(qū)域增長和改進(jìn)分水嶺算法的胃腺癌細(xì)胞圖像分割[J].計(jì)算機(jī)與數(shù)字工程,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] 孫權(quán)森,曾生根,王平安,等.典型相關(guān)分析的理論及其在特征融合中的應(yīng)用[J].計(jì)算機(jī)學(xué)報(bào),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。上??莆攸c(diǎn)項(xiàng)目(14590501700)。李敏,碩士生,主研領(lǐng)域:區(qū)域增長,功能磁共振成像。曾衛(wèi)明,教授。
TP3
A
10.3969/j.issn.1000-386x.2017.03.030