AD-rfMRI圖像的基于改進(jìn)EDC的獨(dú)立分量數(shù)目估計

鄭偉，王軒，姚紀(jì)智，劉帥奇，張曉丹，馬澤鵬

(1.河北大學(xué) 電子信息工程學(xué)院，河北 保定 071002；2.河北省數(shù)字醫(yī)療工程重點(diǎn)實(shí)驗(yàn)室, 河北 保定 071002； 3.河北省機(jī)器視覺工程技術(shù)研究中心, 河北 保定 071002；4.河北大學(xué)附屬醫(yī)院CT-MRI診斷科， 河北 保定 071000)

功能性磁共振成像(functional magnetic resonance imaging，fMRI)和彌散張量成像等成像方式作為一種應(yīng)用于活體的、非侵入式、無損傷的功能成像方式，時間和空間分辨率相對較高，是目前實(shí)現(xiàn)阿爾茨海默病(AD)早期診斷的重要成像方式[1-2]. 而基于fMRI的功能連接可以反映腦區(qū)之間是否相關(guān)以及相關(guān)性的大小，從而反映AD患者早期及病程發(fā)展過程中腦區(qū)之間相關(guān)性的變化規(guī)律，因此具有更好的研究前景.靜息態(tài)功能磁共振成像(resting state fMRI，rfMRI)不需要額外刺激，干擾因素較少，操作簡單而且容易實(shí)現(xiàn)，可重復(fù)性較高. 2017年，Zhang等[1]掃描了96名受試者的rfMRI圖像，采用獨(dú)立分量分析將丘腦分成10個獨(dú)立分量后，通過聚類將獨(dú)立分量劃分成不同的功能簇，每個簇對應(yīng)唯一的一個時間序列，通過計算簇時間序列之間的相關(guān)系數(shù)研究丘腦的功能連接性，結(jié)果表明，獨(dú)立分量分析與聚類結(jié)合對腦區(qū)的功能連接分析優(yōu)于基于種子區(qū)的分析方法. 2018年，劉怡秋等[3]將半高全寬值為6 mm的高斯核用于對AD患者的rfMRI圖像進(jìn)行空間平滑，計算自動解剖標(biāo)簽?zāi)０?AAL)中116個腦區(qū)時間序列之間的相關(guān)性，并通過單因素方差分析，偽發(fā)現(xiàn)率(false discovery rate,FDR)校正以及雙樣本t檢驗(yàn)比較AD病程發(fā)展過程中受試者的功能連接變化. 鄭偉等[4]等針對不同年齡段女性患者靜息態(tài)功能磁共振成像，采用獨(dú)立分量分析及聚類等方法進(jìn)行盲源分離和穩(wěn)定性分析，通過統(tǒng)計檢驗(yàn)及組間對比研究女性患者各腦區(qū)之間的連接關(guān)系及其強(qiáng)弱，進(jìn)而判斷功能連接差異，為女性AD患者臨床前期診斷提供參考. 2019年，馬文潔等[5]對rfMRI圖像進(jìn)行時間層校正、頭動校正、空間歸一化以及去線性趨勢、濾波后，基于種子點(diǎn)選擇默認(rèn)模式網(wǎng)絡(luò)與背側(cè)注意網(wǎng)絡(luò)，并將相關(guān)系數(shù)用于計算功能網(wǎng)絡(luò)時間序列之間的相關(guān)性，從而研究AD患者感興趣區(qū)域的功能連接異常，為探尋AD的發(fā)病機(jī)制和發(fā)病規(guī)律提供支持. 基于負(fù)熵的快速獨(dú)立分量分析(fast independent component analysis algorithm based on negative entropy，NE-FastICA)算法收斂速度快，分離精度高，在腦功能連接分析中普遍使用， 但存在分離之前獨(dú)立分量數(shù)目的過估計問題. 本文提出對實(shí)現(xiàn)獨(dú)立分量數(shù)目估計的有效檢測準(zhǔn)則(effective detection criteria，EDC)進(jìn)行改進(jìn)，將對數(shù)函數(shù)用于懲罰函數(shù)，并與黃金分割法確定參數(shù)γ相結(jié)合，本文稱其為OIEDC2，改進(jìn)后解決了分量數(shù)目的過估計問題，給出了合理的估計數(shù)目.

1 受試者fMRI分組情況

本文將65～80歲的89名女性AD患者的rfMRI圖像(來源于ADNI數(shù)據(jù)庫https://ida.loni.usc.edu/login.jsp?project=ADNI)作為研究對象，按照疾病的發(fā)展階段分為4組，對照組受試者(subjects of normal control，fsub_CN)17人、早期輕度認(rèn)知障礙受試者(subjects of early mild cognitive impairment，fsub_EMCI)23人、晚期輕度認(rèn)知障礙受試者(subjects of late mild cognitive impairment，fsub_LMCI)19人和AD受試者(subjects of AD，fsub_AD)30人的rfMRI圖像，受試者情況如表1所示. rfMRI的相關(guān)參數(shù)如表2所示.

表1 受試者基本信息

表2 圖像基本參數(shù)

2 基于EDC的獨(dú)立分量數(shù)目估計

基于信息理論準(zhǔn)則的估計方法主要包括貝葉斯信息準(zhǔn)則(bayesian information criterion，BIC)[6]、赤池信息準(zhǔn)則(akaike’s information criterion，AIC)[7]、最小描述長度(minimum description length，MDL)[8]、kullback-leibler信息論準(zhǔn)則(kullback-leibler information criterion,KIC)[9]和EDC[10]，這些信息理論準(zhǔn)則通過最小化目標(biāo)函數(shù)實(shí)現(xiàn)對獨(dú)立分量數(shù)目的自適應(yīng)估計，EDC估計如式(1)所示.

(1)

其中，-2L(x(n)|θk)可衡量EDC模型對數(shù)據(jù)的擬合程度，η(θk,N)CN為懲罰函數(shù)項(xiàng)，γ∈[0.1,1]為調(diào)整懲罰項(xiàng)的參數(shù)，可以調(diào)整估計性能，懲罰項(xiàng)中CN是可變的.2種不同懲罰函數(shù)下的EDC分別用EDC1和EDC2表示，如式(2)和式(3)所示.

EDC1(k)=-2L(x(n)|θk)+η(θk,N)(NlgN)γ,γ∈[0.1,1]

，

(2)

EDC2(k)=-2L(x(n)|θk)+η(θk,N)Nγ,γ∈[0.1,1]

.

(3)

信息理論準(zhǔn)則用于估計獨(dú)立分量數(shù)目時，若對rfMRI估計的獨(dú)立分量數(shù)目在20～40，表明估計結(jié)果較合理[10]. 將AIC、KIC、BIC、MDL、EDC1以及EDC2用于估計4組受試者的獨(dú)立分量數(shù)目，結(jié)果如圖1所示.

圖1 不同準(zhǔn)則下對4組受試者的估計結(jié)果對比Fig.1 Comparison of estimation results of four groups of subjects based on different criterion

由圖1可知，AIC、KIC、BIC和MDL的估計性能接近，且在fsub_EMCI、fsub_LMCI和fsub_AD中均存在嚴(yán)重的過估計問題；EDC1和EDC2對每組受試者估計的獨(dú)立分量數(shù)目相對較少，雖然EDC2存在過估計問題，但估計性能優(yōu)于EDC1.

EDC2的過估計問題是由于其懲罰函數(shù)和γ值的非最優(yōu)性導(dǎo)致的，本文提出將對數(shù)函數(shù)引入懲罰函數(shù)，并結(jié)合黃金分割法確定γ最優(yōu)值的EDC2改進(jìn)策略.

3 改進(jìn)懲罰函數(shù)的EDC2及實(shí)驗(yàn)對比分析

3.1 基于對數(shù)懲罰函數(shù)的改進(jìn)EDC2

式(1)中的懲罰函數(shù)項(xiàng)η(θk,N)CN的CN項(xiàng)是可變的，滿足式(4)所示的條件.

(4)

其中，N為rfMRI的混合信號數(shù)，即體素數(shù).

由式(3)可知，EDC2中CN=Nγ對γ求一階導(dǎo)數(shù)的結(jié)果如式(5)所示.

CN=NγlnN,γ∈[0.1,1]

，

(5)

其中，N較大時，隨著γ值增大，CN的變化率較快，EDC2估計的獨(dú)立分量數(shù)減小的較快，容易造成低估計；γ值減小時，EDC2估計的獨(dú)立分量數(shù)增加的較快，容易造成過估計.

由于對數(shù)函數(shù)在自變量較大時變化趨勢較為平緩，本文將對數(shù)函數(shù)引入懲罰函數(shù)項(xiàng)，用于改進(jìn)EDC2，稱為改進(jìn)有效檢測準(zhǔn)則(improvement of effective detection criteria，IEDC2)，如式(6)所示.

CIEDC2=[loga(N)]γ,γ∈[0.1,1]

，

(6)

其中,a為對數(shù)函數(shù)的底數(shù)，可以取2、e和10.

改進(jìn)的懲罰函數(shù)項(xiàng)需滿足式(4)，下面將給出證明：

首先計算關(guān)于式(6)的2個極限，CIEDC2與N比值的極限如式(7)所示.

(7)

函數(shù)y=loga(N)是一個單調(diào)遞增函數(shù)，相同N下，γ=1時CIEDC2達(dá)到最大值，由于N→∞時，N達(dá)到∞的速度比loga(N)快，式(7)中CIEDC2取最大值時的極限為

(8)

當(dāng)γ=0.1時CIEDC2達(dá)到最小值，N→∞時，N達(dá)到∞的速度比[loga(N)]0.1快，式(7)中CIEDC2取最小值時的極限如式(9)所示.

(9)

CIEDC2與CIEDC2max、CIEDC2min之間滿足CIEDC2min

(10)

CIEDC2與lglgN比值的極限如式(11)所示.

(11)

對數(shù)函數(shù)換底公式如式(12)所示.

(12)

將式(12)代入式(11)，結(jié)果如式(13)所示.

(13)

若令lgN=u，則式(13)可轉(zhuǎn)化為式(14).

(14)

其中(lga)γ為與u無關(guān)的常數(shù)，則式(14)可寫成式(15).

(15)

當(dāng)u→∞時，指數(shù)函數(shù)達(dá)到∞的速度比對數(shù)函數(shù)快，因此式(11)的推導(dǎo)結(jié)果為

(16)

由式(10)、式(16)可知，式(6)作為懲罰函數(shù)時滿足式(4)的條件，因此將對數(shù)函數(shù)引入EDC2以改進(jìn)懲罰函數(shù)項(xiàng)，是合理的.

式(6)中，底數(shù)a不同時，對數(shù)函數(shù)變化趨勢不同，a取2時記為f1=log2x,a取e時記為f2=lnx,a取10時記為f3=lgx，3種常用對數(shù)函數(shù)變化趨勢對比如圖2所示.

圖2 3種函數(shù)變化趨勢對比Fig.2 Comparison of three logarithmic functions

由圖2可知，底數(shù)越大，函數(shù)值越小，3個函數(shù)中f3的值最小，反映為自變量較大時f3的變化趨勢較平穩(wěn)，本文將f3用于IEDC2懲罰函數(shù)得式(17)，則IEDC2估計可表示為式(18).

CIEDC2=η(θk,N)(lgN)，γ∈[0.1,1].

(17)

IEDC2(k)=-2L(x(n)|θk)+η(θk,N)(lgN)，γ∈[0.1,1]

.

(18)

3.2 IEDC2實(shí)驗(yàn)結(jié)果及對比分析

為了衡量IEDC2改進(jìn)后的估計性能，以fsub_EMCI、fsub_LMCI 2組受試者為例，將EDC2、IEDC2在不同γ值下的估計結(jié)果進(jìn)行對比.

EDC2的估計結(jié)果如圖3a所示，圖3a中EDC2對2組受試者的估計結(jié)果存在過估計問題；γ值越小，估計的獨(dú)立分量數(shù)目越多，γ值越大，估計的獨(dú)立分量數(shù)目越少，不同γ參數(shù)下，EDC2對fsub_EMCI估計的獨(dú)立分量數(shù)目有3個值為20～40，對fsub_LMCI估計的獨(dú)立分量數(shù)目有4個值為20～40，選擇的獨(dú)立分量數(shù)目不同會導(dǎo)致盲源分離結(jié)果也不同.

IEDC2的估計結(jié)果如圖3b所示，2組受試者的估計結(jié)果雖然在參數(shù)γ較小時仍存在過估計問題，但與圖3(a)相比，相同γ值下，IEDC2估計的獨(dú)立分量數(shù)目明顯下降，表明IEDC2改善了過估計問題；不同參數(shù)γ下，IEDC2對fsub_EMCI、fsub_LMCI估計的獨(dú)立分量數(shù)目均只有1個值在20～40，表明IEDC2避免了因獨(dú)立分量數(shù)目不同導(dǎo)致的盲源分離結(jié)果出現(xiàn)差異，提高了估計結(jié)果準(zhǔn)確性.

a.EDC2;b.IEDC2.

4 基于黃金分割法的γ值確定

EDC2的估計性能除了受到不同懲罰函數(shù)的影響外，還與參數(shù)γ有關(guān)，γ越大EDC2估計的獨(dú)立分量數(shù)越少，越容易造成低估計；γ越小，EDC2估計的獨(dú)立分量數(shù)越多，越容易造成過估計. 將EDC2中的CN記為y1=Nγ=271 633γ，IEDC2中的CN記為y2=[lg(271 633)]γ，則y1、y2隨γ值變化的趨勢如圖4所示.

圖4 y1和y2隨γ值的變化趨勢對比Fig.4 Comparison of the change trend of y1 and y2 with γ values

圖4中橫軸為參數(shù)γ，左邊縱軸對應(yīng)y1=271 633γ的函數(shù)值，右邊縱軸對應(yīng)y2=[lg(271 633)]γ的函數(shù)值，對比可知，體素數(shù)N一定時，相同γ值下，y1函數(shù)值的變化率快，表現(xiàn)為γ值增加時，EDC2估計的獨(dú)立分量數(shù)減小較快，容易造成低估計；γ值減小時，EDC2估計的獨(dú)立分量數(shù)增加較快，容易造成過估計.

本文針對不同γ值下EDC2對4組受試者的估計結(jié)果存在過估計問題提出改進(jìn)，在IEDC2懲罰函數(shù)的基礎(chǔ)上，利用黃金分割法確定γ的最優(yōu)值，改善IEDC2的估計性能.將黃金分割法用于確定最優(yōu)γ值優(yōu)化改進(jìn)有效檢測準(zhǔn)則，本文簡稱為OIEDC2.

黃金分割法每次將搜索區(qū)間縮小0.618倍，根據(jù)精度要求以0.618n的縮減速率尋找最優(yōu)點(diǎn).γ值的范圍是確定的記為[a,b]，黃金分割法首先確定2個初始實(shí)驗(yàn)點(diǎn)r1和r2，其中實(shí)驗(yàn)點(diǎn)r1由式(19)確定.實(shí)驗(yàn)點(diǎn)r2由式(20)確定.

r1=a+0.618×(b-a).

(19)

r2=a+b-r1=b-0.618×(b-a)

.

(20)

黃金分割法第一次確定2個實(shí)驗(yàn)點(diǎn)，之后通過對比、迭代依次確定下一個實(shí)驗(yàn)點(diǎn)，同時更新變量區(qū)間，縮小變量范圍，最終通過評價函數(shù)確定最優(yōu)點(diǎn).黃金分割法通過不斷更新γ的區(qū)間最終確定參數(shù)γ的最優(yōu)值，直至IEDC2(γ)滿足式(21).將滿足式(21)的γ稱為參數(shù)γ的最優(yōu)值，記為γbest.將γbest用于式(18)可得OIEDC2估計如式(22)所示.

20≤IEDC2(γ)≤40

.

(21)

OIEDC2(k)=2L(x(n)|θk)+η(ηk,N)[lgN]γbest.

(22)

5 實(shí)驗(yàn)結(jié)果及對比分析

為了衡量OIEDC2改進(jìn)后的估計性能，將EDC2、IEDC2和OIEDC2分別用于估計4組受試者的獨(dú)立分量數(shù)目，結(jié)果如表3所示.

表3 EDC2、IEDC2和OIEDC2的估計結(jié)果

由表3可知，EDC2對4組受試者的估計結(jié)果均存在不同程度的過估計問題；IEDC2對fsub_CN估計的獨(dú)立分量數(shù)在20～40，估計結(jié)果較合理，對其余3組受試者雖然存在不同程度的過估計問題，但對fsub_LMCI估計的獨(dú)立分量數(shù)與40較接近，估計結(jié)果比EDC2好，對fsub_EMCI和fsub_AD估計的獨(dú)立分量數(shù)比EDC2低，表明與EDC2相比，IEDC2的估計性能有所改善；OIEDC2對4組受試者估計的獨(dú)立分量數(shù)目均處于20～40，估計結(jié)果都比較合理，因此本文改進(jìn)后的OIEDC2避免了過估計問題，實(shí)驗(yàn)結(jié)果比改進(jìn)前的EDC2要好.

6 結(jié)論

本文對基于EDC2的獨(dú)立分量數(shù)目估計方法進(jìn)行改進(jìn).將對數(shù)函數(shù)引入EDC2的懲罰函數(shù)項(xiàng)，稱為IEDC2，將黃金分割法引入IEDC2確定最優(yōu)參數(shù)γ，稱為OIEDC2，并將結(jié)果進(jìn)行對比分析，衡量EDC2改進(jìn)前后的估計性能.實(shí)驗(yàn)結(jié)果表明，本文方法有效解決了4組受試者獨(dú)立分量數(shù)目的過估計問題，估計結(jié)果更合理.為后續(xù)的基于獨(dú)立分量的盲源分離，獨(dú)立分量選擇，感興趣區(qū)域功能連接計算，實(shí)現(xiàn)4組受試者的功能連接分析，并通過統(tǒng)計檢驗(yàn)獲取隨疾病發(fā)展進(jìn)程感興趣區(qū)域功能連接的變化規(guī)律，進(jìn)而為尋找導(dǎo)致功能連接變化的內(nèi)外部因素以及性別與AD發(fā)病的關(guān)系，即為AD發(fā)病機(jī)制的研究奠定基礎(chǔ).