基于平滑最小方差無失真響應的一致性同步算法研究

顧廣華 崔 冬* 王 娟 齊順愛 李小俚

1(燕山大學信息科學與工程學院，河北 秦皇島 066004)2(河北省信息傳輸與信號處理重點實驗室，河北 秦皇島 066004)3(北京師范大學認知神經(jīng)科學和學習國家重點實驗室，北京 100875)

引言

一致性分析廣泛應用于探討腦電信號的功能連接性上，在知覺、行動和認知等方面的神經(jīng)同步方面起著重要作用[1]。電極對間的一致性值越大，說明兩路信號相互依賴的程度越強，反映大腦相應位點具有較高的功能連接性[2]。

糖尿病是輕度認知功能障礙(mild cognitive impairment, MCI)的危險因素之一[3]。MCI是介于正常衰老和癡呆的中間狀態(tài)，分為遺忘型輕度認知障礙(amnestic mild cognitive impairment, aMCI)和非遺忘型輕度認知障礙(non-amnestic mild cognitive impairment, nMCI)，其中aMCI最易發(fā)展成為癡呆，每年約有10%～15% 的患者發(fā)展成癡呆[4]，對aMCI盡早進行干預可預防和延緩發(fā)展成癡呆的情況。

一致性可應用于分析靜息狀態(tài)下大腦功能連接的情況，也可以應用于檢查執(zhí)行認知任務時大腦功能連接的改變。Locatelli對比分析10例癡呆患者和10 例正常老人，發(fā)現(xiàn)癡呆患者在顳-頂區(qū)α頻段一致性下降，顳-頂-枕區(qū)兩半球間θ頻段和β頻段一致性下降，而在額葉和枕部區(qū)域間δ頻段的一致性升高[5]。Babiloni等對大腦內(nèi)部的電極對(F3、Fz、F4、C3、Cz、C4、P3、Pz、P4)進行一致性研究，發(fā)現(xiàn)AD患者在大腦兩半球間一致性下降，而在前額-頂區(qū)一致性升高[6]。Jelles等對15 例AD患者和17 例年齡匹配的正常人進行一致性研究，發(fā)現(xiàn)在α頻段和β頻段患者的一致性低于正常組[7]。Brunovsky發(fā)現(xiàn)，癡呆患者的α頻段一致性下降，δ頻段的一致性增加，且這些特征和癡呆嚴重程度有關[8]。Stam發(fā)現(xiàn)，AD患者在β頻段一致性下降并且和MMSE得分呈相關性，在α頻段AD患者的一致性值低于正常對照組[9]。Jelic等發(fā)現(xiàn)，AD患者在顳-頂區(qū)的一致性值與記憶障礙患者和正常老年人相比下降[10]。

本研究基于最小方差無失真響應(minimum variance distortionless response, MVDR)一致性算法[11-12]，提出結合CB核平滑功率譜的平滑最小方差無失真響應算法(smoothing minimum variance distortionless response, SMVDR)。采用線性模型，在窄帶信號和寬帶信號情況下，對MVDR和SMVDR算法的性能進行仿真分析，包括抗噪性能、信號幅度。仿真結果顯示，SMVDR算法在準確性及抗噪性能方面明顯優(yōu)于MVDR算法。利用SMVDR算法，分析糖尿病患者輕度認知障礙腦電數(shù)據(jù)，對aMCI組和對照組在不同頻段、不同通道組間進行獨立樣本t檢驗統(tǒng)計分析，找出有顯著特征的通道組及頻帶，并將一致性結果與認知量表進行相關性分析，探討能夠區(qū)分糖尿病人中具有輕度認知障礙和認知正?；颊叩奶卣?。

1 方法

1.1 MVDR一致性算法

一致性分析可以反映兩個隨機過程在各頻率成分上的同步關系，基于譜分析描述大腦在半球間以及大腦半球內(nèi)不同區(qū)域、不同節(jié)律上電活動的同步情況[13-15]。MVDR算法是一種非參數(shù)譜的估計算法，主要是基于濾波器組的分解，令信號通過一組受限的帶通濾波器，其中每個濾波器都固定一個分析頻率，從而在每個頻帶估計信號的頻譜。

設x(n)是零均值平穩(wěn)隨機過程，表示為

x(n)=[x(n)x(n-1)…x(n-L+1)]T

(1)

式中，上標T表示向量的轉置。

令x(n)是K個濾波器的輸入，每個濾波器的長度為L[11-12]，有

(2)

濾波器gk的輸出yk(n)的功率為

(3)

式中，E{·}為數(shù)學期望，上標H表示向量或矩陣的轉置共軛，協(xié)方差矩陣Rxx為

Rxx=E{x(n)x(n)H}

(4)

(5)

在頻率ωk處x(n)的頻譜為

(6)

對于MVDR譜估計,每個濾波器gk的設計相當于最小化下面的函數(shù)，有

(7)

則有

(8)

因此，MVDR譜為

(9)

對于兩個零均值平穩(wěn)隨機信號x1(n)和x2(n),兩個濾波器可表示為

(10)

在頻率ωk處的自譜Sxixi(ωk) (i=1,2)和交叉譜Sx1x2(ωk)為

(i,j=1,2)

(11)

式中，上標*表示復共軛。

最后，兩個信號x1(n)和x2(n)之間的MVDR一致性函數(shù)定義為

(12)

1.2 SMVDR一致性算法

本研究提出新的計算一致性的算法，通過CB核平滑的最小方差無失真地響應算法MVDR頻譜。CB核廣泛應用在時頻分析中，克服了高斯核的缺點，是一種緊致支撐核[16]。CB核通過控制核的參數(shù)，可以有效地抑制交叉項，提高時頻分辨率，恢復由于截斷引起的信息丟失，同時提高算法的準確性[17]。通過計算CB核和瞬時自相關函數(shù)的卷積來產(chǎn)生廣義瞬時自相關函數(shù)，可以減少有效信息的損失，使能量集中，提高準確性。CB核的定義[18]如下：

(13)

式中：N為采樣點，n=1,2,…,N；M代表窗長，通常取奇數(shù)值，m=1,2,…,M。

核的寬度通過參數(shù)M/B控制，峰值通過參數(shù)C調整，以上參數(shù)使CB核具有較高的分辨率和足夠的交叉相抑制。

對于計算SMVDR一致性，首先把通過MVDR方法計算出來的自譜的平方值和交叉頻譜的平方值作為信號，再計算這3種平方譜的瞬時頻率自相關函數(shù)，有

(14)

式中，m是頻率延遲。

給定該CB核的窗長M，計算核與瞬時頻率自相關函數(shù)的卷積，即進行核平滑操作，定義廣義瞬時頻率自相關函數(shù)如下：

(15)

式中，?表示卷積操作。

令頻率延遲m為某一常數(shù)a時，得到新的SMVDR一致性計算方法為

這種改進的計算一致性的SMVDR方法 也是關于頻率的函數(shù)，并且取值范圍同樣在0～1之間,可以評估兩信號間的相關性程度，適用于雙通道腦電信號的同步分析。

1.3 算法仿真分析

選取線性仿真模型，包含兩個正弦信號x1(n)和x2(n)，表達形式[11]如下：

(17)

式中：fi是正弦信號的頻率，其中i=0,1,…,Nf-1；fs為采樣頻率；φi是隨機相位；r1(n)和r2(n)是兩個獨立的均值為0、方差為1的高斯白噪聲。

為獲得窄帶信號，設定5個一致性頻率：f0=10,f1=15,f2=25,f3=30,f4=35 Hz。理論上，在這5個頻率點處一致性值應該接近1，而在其他頻率處一致性值應該接近0。

濾波器個數(shù)K設為1 500，濾波器長度L設為150。這種參數(shù)設置可以保證一致性估計值較為準確，并且計算復雜度不會太高[19]。CB核中參數(shù)的選取，C=1.1，M/B=1.2，可以保證分辨率較高，并且能較好地抑制交叉項[20]。平滑操作中，M=3，a=1。

1.3.1抗噪性能

為了比較兩種一致性算法的抗噪性能，在線性模型產(chǎn)生的窄帶信號中加高斯白噪聲，其信噪比SNR分別為-10、-5、5和20 dB。噪聲對MVDR和SMVDR算法的影響如圖1所示?？梢钥闯?，隨著信噪比的增加，兩種一致性算法性能更好，并且在5個一致性頻率處的峰值也隨著信噪比的增大而增大，如圖1 (b)～(d) 所示。SMVDR算法在一致性頻率處具有更窄的峰值，而在其他頻率處更接近0，因此該算法在抗噪性能方面要優(yōu)于MVDR算法。

圖1 噪聲對窄帶信號一致性算法的影響。(a) SNR=-10 dB；(b) SNR=-5 dB；(c) SNR=5 dB；(d) SNR=20 dBFig.1 Effects of the noise impact on coherence methods in narrow-band signals. (a) SNR=-10 dB; (b) SNR=-5 dB; (c) SNR=5 dB; (d) SNR=20 dB

1.3.2信號幅度

為了觀察信號不同頻率間相對功率對一致性算法的影響，將最初選定的5個一致性頻率分成兩組，設定兩個信號在前3個頻率點f0=10、f1=15和f2=25 Hz處以及后兩個頻率點f3=30、f4=35Hz處的信號幅度之比分別為10∶1和5∶1，一致性對比的仿真圖形如圖2所示。

如圖2(a) 所示，當兩組信號的幅度對比是10∶1時，后兩個頻率點f3=30 Hz和f4=35 Hz處的一致性值不明顯；如圖2(b) 所示，當信號幅度增大時，后兩個頻率點處的一致性值有所提高?？梢?，兩種一致性算法均受信號幅度的影響。因此，在對實際腦電信號進行數(shù)據(jù)分析時，為減小信號不同頻段振幅的影響，先對實際的腦電信號進行分頻段的濾波處理，再對每一頻段的信號用SMVDR算法進行一致性計算，使分析的結果更準確。

圖2 信號幅度對一致性算法的影響。(a) 幅度對比10∶1；(b) 幅度對比5∶1Fig.2 Effects of the amplitude on coherence methods. (a) Amplitude contrast is 10∶1; (b) Amplitude contrast is 5∶1

1.4 腦電信號一致性分析

1.4.1糖尿病患者資料

本次研究的實驗對象選取在人民解放軍火箭軍總醫(yī)院就診的31 名糖尿病患者，均符合aMCI的納入及排除標準。所有受試者均被告知研究內(nèi)容并同意參加實驗，經(jīng)北京師范大學倫理委員會批準。

記錄患者的姓名、年齡、性別、受教育程度、MOCA得分情況。依據(jù)MOCA得分對患者進行分組：MOCA分數(shù)<26為aMCI組，18 例；MOCA分數(shù)≥26為對照組，13例。運用SPSS統(tǒng)計軟件，對兩組患者的基本人口特征和神經(jīng)心理學量表進行獨立樣本t檢驗，結果用“均值±標準差”的形式表示，如表1所示?？梢钥闯觯瑑山M患者在人口統(tǒng)計特征中不存在顯著差異，即年齡、受教育程度方面是匹配的，并且MOCA分數(shù)達到顯著水平。

表1人口學特征和神經(jīng)心理學檢查統(tǒng)計結果

Tab.1Theresultsofdemographiccharacteristicsandneuropsychologicalexamination

因素aMCI組對照組P值年齡/歲70.44±7.2671.92±6.700.568受教育程度13.06±4.5613.38±2.730.818MOCA21.94±1.2626.92±1.190.000???

1.4.2腦電信號采集與預處理

采用EGI的GES300 型128 通道腦電記錄儀，試者處于閉眼靜息狀態(tài)進行靜息EEG數(shù)據(jù)采集，雙側乳突為參考電極，采樣頻率500 Hz，時長5 min，通過0～200 Hz帶通濾波，阻抗<20 kΩ。采用基于小波增強的獨立分量分析算法[21]對信號進行預處理，去除眼電等偽跡。

1.4.3腦電信號數(shù)據(jù)分析

本研究選擇19 路電極并劃分在5 個大腦區(qū)域，分別為前額(frontal, F)、頂區(qū)(parietal, P)、枕區(qū)(occipital, O)、左顳(left temporal, LT)、右顳(right temporal, RT)，如圖3 所示。

圖3 腦電信號19電極安置區(qū)域Fig.3 The location of nineteen EEG electrodes

為了研究和比較aMCI組和對照組在大腦不同節(jié)律、不同腦區(qū)的一致性特征，將所選的雙通道電極對分為3個組，分別是大腦兩側電極對、大腦一側短距離電極對和大腦一側長距離電極對，具體選取情況如表2所示。

表2腦區(qū)間電極對的選取情況

Tab.2Theselectionofelectrodepairsbetweenbrainregions

區(qū)域分組通道組大腦兩側Fp1Fp2,F3F4,C3C4,P3P4,O1O2,F7F8,T3T4,T5T6,F7T4,T3T6,F8T5,P4O1,F7T6,F8T3,Fp1F4,F8C3,PzT4,F4O1,F8O1,T3O2,F8Pz,C4F3,FzPz,F4T5大腦一側短距離P3O1,P4O2,F3C3,F4C4,C4P4,C3P3,F7T3,F8T4,T3T5,T4T6,F8F4大腦一側長距離F3T3,F4T4,F3T5,F4T6,F3P3,F4P4,F3O1,F4O2,T3O1,T4O2,T3P3,T4P4,T6C4,T6P6,T4C4,F8P4,Fp1O1,F8C4

在真實腦電數(shù)據(jù)分析中，算法參數(shù)設置與仿真相同，即濾波器長度L和濾波器個數(shù)K分別取150和1 500，CB核中參數(shù)M和a分別取3 和1。實際腦電信號噪聲干擾較多，且頻率成分在各頻段分布不均勻，SMVDR算法具有較好的抗噪性能，但估計出的一致性值會受到各頻率振幅的影響。因此，在用SMVDR算法進行腦電信號數(shù)據(jù)分析前，首先對實際腦電信號進行濾波。本研究利用Gabor小波變換，將EEG信號分為Delta(1～4 Hz)、Theta(4～8 Hz)，alpha(8～13 Hz)和beta(13～30 Hz)4個頻段。在每個頻段上，各電極對之間的一致性值通過下式進行計算，有

(18)

式中，U和L分別為指定頻段的上限頻率和下限頻率。

數(shù)據(jù)長度截取為5 s一段，計算出每對電極的一致性值，剔除超出均值±3倍標準方差的異常值，求得的均值即為每對信號的一致性值，可用于后續(xù)的統(tǒng)計分析。

利用SPSS 20.0統(tǒng)計軟件對數(shù)據(jù)進行統(tǒng)計分析，兩組間的比較采用獨立樣本t檢驗，P<0.05認為差異有統(tǒng)計學意義。采用皮爾森線性相關分析所有患者一致性值與MOCA測試得分的相關性，利用Bonferroni方法對顯著相關值進行了校正，以排除由于概率的原因出現(xiàn)的顯著相關，探究腦電信號一致性值與認知功能的內(nèi)在聯(lián)系。

2 結果

2.1 aMCI組與對照組統(tǒng)計分析

統(tǒng)計分析aMCI組和對照組在各個頻段中有顯著特征的通道組的情況，結果如圖4所示?？梢钥闯?，delta頻段aMCI組的一致性值要普遍低于對照組；theta頻段aMCI組在大腦兩側的右額-左枕(F4O1)以及右顳-左枕(F8O1)通道組間的一致性值顯著高于對照組，且在F4O1處具有較強的顯著性；alpha頻段aMCI組的一致性值都要普遍低于對照組，且在T4O2處具有較強的顯著性；beta頻段aMCI組的一致性值都要普遍高于對照組，且通道組都集中在大腦兩側，如T5T6、F8T5、F8O1和F4O1。

圖4 各頻段aMCI組與對照組的一致性統(tǒng)計分析結果(每列上為腦地形圖，下為一致性對比分析；腦地形圖中，?表示aMCI組一致性值小于對照組，表示aMCI組一致性值大于對照組；一致性對比分析圖中，*表示P<0.05，**表示P<0.01)。(a) delta頻段；(b) theta頻段；(c) alpha頻段； (d)beta頻段Fig.4 The coherence statistic analysis results of aMCI group and normal control group in each frequency band(The top line of each row is brain map, the bottom line is the graph of coherence comparison analysis between the two groups; in brain map, ? donates that coherence values of aMCI is lower than that of NC, ★ donates that coherence values of aMCI is lager than that of NC; in coherence comparison analysis graph, * donates P<0.05, ** donates P<0.01). (a) Delta frequency band; (b)Theta frequency band; (c) Alpha frequency band; (d) Beta frequency band

2.2 一致性與認知功能的相關性分析

腦電信號一致性值與MOCA得分的相關性分析結果如圖5所示?？梢钥闯?，在F7F8和T4T6的δ頻段一致性均與MOCA得分存在顯著的正相關；在θ頻段，F(xiàn)4O1和F8O1通道組的一致性均與MOCA得分存在顯著的負相關；在F7F8、F3O1、T4O2、F7T3、F8F4和FP1O1的α頻段，一致性均與MOCA得分存在顯著的正相關；在左右顳間T5T6、F8T5，右額-左枕F4O1和右顳-左枕F8O1的β頻段，一致性均與MOCA得分呈顯著的負相關。

3 討論

本研究提出基于CB核平滑的SMVDR一致性算法。在算法的仿真分析中，利用線性模型從抗噪性能、信號振幅影響進行兩種方法的對比。從仿真結果中可以看出，SMVDR算法比MVDR算法具有較好的抗噪性及準確性。

利用SMVDR一致性算法分析糖尿病患者輕度認知障礙的腦電信號，發(fā)現(xiàn)在有顯著特征的通道組中，在δ頻段aMCI組的一致性值要普遍低于對照組，尤其是在大腦兩側的左右顳區(qū)(F7F8)δ頻段，與Babiloni研究結果相同[6]。近幾年對老年癡呆患者進行一致性研究，發(fā)現(xiàn)其中的某些頻帶內(nèi)一致性值會增大，而在某些頻帶中一致性值會減小。

在θ頻段，aMCI組的一致性值要高于對照組。Locatelli等的研究結果也顯示，癡呆患者在大腦兩側θ頻段一致性值增大[5]。θ頻段aMCI的一致性值增大，可能與皮層下結構的皮層傳入神經(jīng)有關[21]，也可能與補償機制有關，因為認知障礙損傷患者在記憶能力以及認知功能等方面都會減弱，使得在q頻段要對其進行補償[22]。這些原因可能導致某些通道組在θ頻段的一致性值升高。

α頻段aMCI組的一致性值都要普遍低于對照組。從近些年針對老年癡呆患者腦電信號的一致性研究結果來看，Locatelli、Knott以及Stam等得出的結論表明，癡呆患者的腦電信號一致性值要低于正常對照組，尤其是在α頻段[5,9,23]。其中，Locatelli提出癡呆患者在大腦皮層的前額-枕區(qū)以及右顳-枕區(qū)的一致性值要低于對照組，本研究F3O1和T4O2通道組的研究結果與其相符[5]。有研究顯示，2型糖尿病認知障礙患者在耦合強度方面，大腦皮層的前額-左枕以及右額-右顳的α頻段一致性值要低于對照組[24]，本研究Fp1O1和F8F4通道組與其研究結果相符。α頻段的一致性情況可以揭示皮質區(qū)域之間的功能連接水平，糖尿病患者中aMCI組的一致性值降低可能與不同腦區(qū)大腦皮質之間的連接能力受阻有關[21]。

圖5 不同頻段神經(jīng)心理學量表(MOCA)與一致性值線性分析結果(每個子圖中第一行為相關性分析腦地形圖，其余為一致性值與MOCA測試得分散點圖；腦地形圖中，?表示一致性值與MOCA得分呈顯著性正相關，表示一致性值與MOCA得分呈顯著性負相關 )。(a)delta頻段；(b)theta頻段; (c)alpha頻段；(d)beta頻段Fig.5 The results of correlation analysis results of neuropsychological scale and coherence values in each frequency band(The first line of each subgraph is brain map of correlation analysis, the others are the scatter graphs of correlation between coherence values and MOCA scores; in brain map, ? donates a significant positive correlation between coherence values and MOCA scores, donates a significant negative correlation between coherence values and MOCA scores). (a) Delta frequency band；(b) Theta frequency band；(c) Alpha frequency band；(d)Beta frequency band

β頻段aMCI組的一致性值要高于對照組。有研究表明，糖尿病認知障礙患者的腦電信號在耦合強度方面，位于大腦區(qū)域左右顳間的aMCI組一致性值要高于對照組，與本研究中左右顳間的T5T6和F8T5通道組與糖尿病認知障礙患者的腦電信號在耦合強度方面的研究結果相符合[25]。

本研究主要對新算法SMVDR方法和原有方法做仿真分析，并將其應用于實際腦電信號的分析中。仿真實驗和實際數(shù)據(jù)的分析分別得出了一些結論，但研究中尚存在一些局限性：樣本數(shù)量較少，只有31人；通道數(shù)目較少，只選取了19路腦電信號等。在以后的工作中，將對此做進一步分析。

4 結論

本研究提出了一種新的計算雙通道腦電信號的一致性算法——SMVDR方法，該方法具有較好的抗噪性能及準確性。通過對糖尿病患者輕度認知障礙的腦電信號一致性的統(tǒng)計分析，發(fā)現(xiàn)這些患者在δ頻段的左右顳間一致性下降，在θ頻段的前額-枕區(qū)一致性增加，在α頻段的右顳-枕區(qū)和前額-右顳區(qū)域一致性下降，在β頻段的左右顳間一致性增加。通過一致性值與MOCA得分的相關性分析，發(fā)現(xiàn)在特征通道下δ和α頻段的一致性與MOCA得分存在顯著的正相關，θ和β頻段的一致性與MOCA分數(shù)呈負相關。SMVDR方法以及相應節(jié)律下的不同腦區(qū)之間的一致性特征，可以輔助醫(yī)生盡早地診斷糖尿病輕度認知障礙，對預測輕度認知障礙并及時進行干預具有重要意義。