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

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

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

1 受試者fMRI分組情況

本文將65～80歲的89名女性AD患者的rfMRI圖像(來(lái)源于ADNI數(shù)據(jù)庫(kù)https://ida.loni.usc.edu/login.jsp?project=ADNI)作為研究對(duì)象，按照疾病的發(fā)展階段分為4組，對(duì)照組受試者(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ù)目估計(jì)

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

(1)

其中，-2L(x(n)|θk)可衡量EDC模型對(duì)數(shù)據(jù)的擬合程度，η(θk,N)CN為懲罰函數(shù)項(xiàng)，γ∈[0.1,1]為調(diào)整懲罰項(xiàng)的參數(shù)，可以調(diào)整估計(jì)性能，懲罰項(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)則用于估計(jì)獨(dú)立分量數(shù)目時(shí)，若對(duì)rfMRI估計(jì)的獨(dú)立分量數(shù)目在20～40，表明估計(jì)結(jié)果較合理[10]. 將AIC、KIC、BIC、MDL、EDC1以及EDC2用于估計(jì)4組受試者的獨(dú)立分量數(shù)目，結(jié)果如圖1所示.

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

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

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

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

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

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

(4)

其中，N為rfMRI的混合信號(hào)數(shù)，即體素?cái)?shù).

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

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

，

(5)

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

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

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

，

(6)

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

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

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

(7)

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

(8)

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

(9)

CIEDC2與CIEDC2max、CIEDC2min之間滿(mǎn)足CIEDC2min

(10)

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

(11)

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

(12)

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

(13)

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

(14)

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

(15)

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

(16)

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

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

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

由圖2可知，底數(shù)越大，函數(shù)值越小，3個(gè)函數(shù)中f3的值最小，反映為自變量較大時(shí)f3的變化趨勢(shì)較平穩(wěn)，本文將f3用于IEDC2懲罰函數(shù)得式(17)，則IEDC2估計(jì)可表示為式(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é)果及對(duì)比分析

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

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

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

a.EDC2;b.IEDC2.

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

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

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

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

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

黃金分割法每次將搜索區(qū)間縮小0.618倍，根據(jù)精度要求以0.618n的縮減速率尋找最優(yōu)點(diǎn).γ值的范圍是確定的記為[a,b]，黃金分割法首先確定2個(gè)初始實(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個(gè)實(shí)驗(yàn)點(diǎn)，之后通過(guò)對(duì)比、迭代依次確定下一個(gè)實(shí)驗(yàn)點(diǎn)，同時(shí)更新變量區(qū)間，縮小變量范圍，最終通過(guò)評(píng)價(jià)函數(shù)確定最優(yōu)點(diǎn).黃金分割法通過(guò)不斷更新γ的區(qū)間最終確定參數(shù)γ的最優(yōu)值，直至IEDC2(γ)滿(mǎn)足式(21).將滿(mǎn)足式(21)的γ稱(chēng)為參數(shù)γ的最優(yōu)值，記為γbest.將γbest用于式(18)可得OIEDC2估計(jì)如式(22)所示.

20≤IEDC2(γ)≤40

.

(21)

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

(22)

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

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

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

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

6 結(jié)論

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