非線性格蘭杰因果關(guān)系在睡眠生理信號分析中的應(yīng)用

杜 朋 戴加飛 李 錦 王 俊 侯鳳貞

(1．南京郵電大學(xué)圖像處理與圖像通信江蘇省重點實驗室, 南京， 210003； 2．南京軍區(qū)南京總醫(yī)院神經(jīng)內(nèi)科, 南京，210002； 3．陜西師范大學(xué)物理學(xué)與信息技術(shù)學(xué)院, 西安， 710062； 4．中國藥科大學(xué)理學(xué)院, 南京，210009)

非線性格蘭杰因果關(guān)系在睡眠生理信號分析中的應(yīng)用

杜 朋1戴加飛2李 錦3王 俊1侯鳳貞4

(1．南京郵電大學(xué)圖像處理與圖像通信江蘇省重點實驗室, 南京， 210003； 2．南京軍區(qū)南京總醫(yī)院神經(jīng)內(nèi)科, 南京，210002； 3．陜西師范大學(xué)物理學(xué)與信息技術(shù)學(xué)院, 西安， 710062； 4．中國藥科大學(xué)理學(xué)院, 南京，210009)

基于非線性格蘭杰因果關(guān)系分析睡眠生理信號。分別使用多項式核函數(shù)、高斯核函數(shù)和Sigmoid核函數(shù)將低維空間數(shù)據(jù)映射到高維特征空間,在高維特征空間使用非線性格蘭杰因果方法來分析睡眠生理信號。研究結(jié)果表明，腦電信號對心電信號的影響比心電信號對腦電信號的影響更為顯著，腦電信號對血壓信號的影響比血壓信號對腦電信號的影響更為顯著，血壓對心電信號的影響比心電信號對血壓信號的影響更為顯著,而且睡眠期樣本信號間的格蘭杰因果關(guān)系更為顯著。仿真結(jié)果驗證了睡眠期信號更能客觀地反映生理信號的因果關(guān)系。

非線性格蘭杰因果關(guān)系； 睡眠信號；核函數(shù)

引 言

睡眠是一個以不動狀態(tài)為特征的基本生理過程，在典型的睡眠狀態(tài)里，肢體感覺降低。睡眠普遍存在于所有的動物中，盡管它的確切功能尚未完全理解[1, 2]，但在許多重要軀體、認知和心理過程中都有重要的作用，尤其在維持正常的大腦和身體平衡方面[3]。生理信號的研究表明從清醒期到睡眠期的過渡是平穩(wěn)的過程。腦電圖可以較為直觀地反映大腦組織的腦電活動和功能狀態(tài)，因而在睡眠信號研究過程中，腦電(Electroencephalogram, EEG)信號一直是重要的工具[4]。睡眠過程受自主神經(jīng)系統(tǒng)的調(diào)節(jié)，心電(Electrocardlogram, ECG)信號變化在反映自主神經(jīng)對心臟活動的調(diào)節(jié)及評價自主神經(jīng)功能方面，比其他生理參數(shù)有著更好的特異度和靈敏度[5]。血壓與睡眠有較為密切的關(guān)系[6]，對高血壓患者在基礎(chǔ)治療情況下進行睡眠干預(yù)，往往可以取得較好的效果。睡眠期腦電由于主觀思維的干擾較少，從而更能客觀反映生理活動狀態(tài)。本文的研究點便是心電信號、腦電信號和血壓信號這3種信號間的因果關(guān)系，以及睡眠對這3種信號間因果關(guān)系的影響。

在分析生理信號的過程中，不僅需要分析各個信號的同步狀態(tài)，也要研究信號間的因果關(guān)系。通過分析因果關(guān)系可以確定子網(wǎng)絡(luò)之間相互影響的方向和大小，從而降低復(fù)雜網(wǎng)絡(luò)分析的難度。目前分析信號因果關(guān)系的主流方法包括信息論中基于條件轉(zhuǎn)移熵的分析方法，然而當(dāng)參數(shù)較多時，轉(zhuǎn)移熵的計算將會非常復(fù)雜。另一個主流的方法便是由格蘭杰提出的，格蘭杰因果關(guān)系認為，如果引入第二個時間序列,可以使得第一個時間序列預(yù)期估值誤差方差減少，則認為第二個時間序列是第一個時間序列的因[7]。經(jīng)過西蒙斯的發(fā)展，格蘭杰因果檢驗作為一種計量方法已經(jīng)被經(jīng)濟學(xué)家普遍接受并廣泛使用。格蘭杰最初提出這個概念時只是對變量均值之間因果關(guān)系的一種考察，廣泛使用的方法也只驗證了變量之間的線性關(guān)系，而在非線性序列領(lǐng)域，傳統(tǒng)的格蘭杰因果檢驗不再適用。在較長時間的實驗研究中都假設(shè)所檢驗的序列是平穩(wěn)的，而實際情況下時間序列往往是非平穩(wěn)的，這時線性的分析方法往往導(dǎo)致信息的丟失，如果使用線性格蘭杰因果關(guān)系來分析數(shù)據(jù)，很可能得到錯誤的結(jié)論。為了促進睡眠生理信號的進一步研究，本文采用非線性格蘭杰因果關(guān)系指數(shù)來對睡眠生理信號進行分析，討論睡眠腦電信號和心電信號以及血壓信號之間的因果關(guān)系。

為了將線性格蘭杰因果關(guān)系擴展到非線性特征空間，本文引入核函數(shù)方法。通過核函數(shù)特征映射，將低維輸入空間中的線性不可分數(shù)據(jù)映射成高維特征空間(可再生Hilbert空間)[8]中的線性可分數(shù)據(jù)。特征空間的內(nèi)積可以通過某個核函數(shù)來表示，核函數(shù)的本質(zhì)是對應(yīng)于高維空間的內(nèi)積，從而生成高維空間的特征映射,而高維空間的數(shù)據(jù)處理可借用核方法的該對應(yīng)關(guān)系使用非線性映射來實現(xiàn)[9]

1 基本原理

1.1 線性格蘭杰因果關(guān)系模型

線性格蘭杰因果系數(shù)為

(1)

1.2 非線性格蘭杰因果關(guān)系模型

計算得到格蘭杰因果系數(shù)為

(2)

2 睡眠生理信號的非線性格蘭杰因果分析

2.1 實驗數(shù)據(jù)

本節(jié)所使用的睡眠數(shù)據(jù)均來自PhysioBank的MIT-BIH Polysomnographic Database。該庫中記錄的是多參數(shù)睡眠數(shù)據(jù)，其中包括1導(dǎo)EEG信號、1導(dǎo)ECG信號、血壓信號等多導(dǎo)睡眠信號，記錄長度為6 h，數(shù)據(jù)采樣頻率為250 Hz。本文對受試者Slp04,Slp14,Slp32,Slp48,Slp41,Slp45,Slp59,Slp60,Slp61和Slp66這10個樣本進行了測試，采用各自幾組多導(dǎo)睡眠數(shù)據(jù)中的1導(dǎo)EEG信號，提取其中的清醒期和NREM睡眠I期的信號用于實驗。為進一步驗證格蘭杰因果系數(shù)的普遍性，使用了由南京軍區(qū)總醫(yī)院從臨床診斷中采集得到的正常人腦電信號。該庫中記錄了多參數(shù)腦電數(shù)據(jù)(包括16個導(dǎo)聯(lián)信號)這些信號的采樣周期為512 Hz且記錄長度均大于1 min。本文分別從正常腦電數(shù)據(jù)中取出10組數(shù)據(jù)，取其Fp1導(dǎo)聯(lián)信號和ECG導(dǎo)聯(lián)信號，經(jīng)實驗驗證，在16個導(dǎo)聯(lián)信號中，第一導(dǎo)聯(lián)Fp1信號區(qū)分度較好，所以下面以Fp1導(dǎo)聯(lián)數(shù)據(jù)做分析對象。

2.2 分析方法

首先讀取各組原始睡眠腦電信號，由于原始信號中包括腦電信號，血壓信號以及心電信號3組信號，因此需要對這3組信號分別進行兩兩處理。從3組信號中取出2組，分別記為時間序列{x}和{y}，由于實驗計算得到因果系數(shù)的計算結(jié)果在m=2時出現(xiàn)峰值，這里取m=2時間序列，對序列分別進行均值、歸一化以及m次平移處理，組成矩陣X和Y，參考1.2小節(jié)中的方法，分別代入3類不同的核函數(shù)并取不同的核函數(shù)參數(shù)，即可求得相應(yīng)的格蘭杰因果系數(shù)矩陣。

核函數(shù)系數(shù)較為合適時可以得到明顯的區(qū)分度，對于多項式核函數(shù)參數(shù),當(dāng)p=3時格蘭杰因果系數(shù)計算結(jié)果較為顯著;對于高斯核函數(shù)，在σ=3到σ=8的區(qū)間中計算結(jié)果相差較小，因而采用了σ=5的計算結(jié)果;對于Sigmoid核函數(shù)，需要考慮到Sigmoid核函數(shù)是有條件的正定函數(shù)，因此需要滿足正定條件，本文取α=1/9，β=0可以滿足正定函數(shù)的條件，計算結(jié)果也較為顯著。

2.3 數(shù)值計算結(jié)果與分析

實驗1驗證耦Henon網(wǎng)絡(luò)函數(shù)之間的格蘭杰因果關(guān)系

由兩個耦合函數(shù)生成的Henon網(wǎng)絡(luò)表達式為

分析{x}序列和{y}序列，從表達式的形式中可以看出{x}序列是{y}序列的因。取m=2，利用3種核函數(shù)，選取不同的參數(shù)來計算Henon在1 000次循環(huán)中生成序列的格蘭杰因果系數(shù),結(jié)果分別如表1,2和3所示。

表1 多項式核函數(shù)計算得到的Henon網(wǎng)絡(luò)x→y方向格蘭杰因果系數(shù)

表2 高斯核函數(shù)計算得到的Henon網(wǎng)絡(luò)x→y方向格蘭杰因果系數(shù)

表3 Sigmoid核函數(shù)計算得到的Henon網(wǎng)絡(luò)x→y方向格蘭杰因果系數(shù)

從表1可以看出x→y方向的格蘭杰因果系數(shù)總是非零,而且在p=2時取得最大值，這說明了Heanon網(wǎng)絡(luò)間的因果影響是平方級數(shù)的，這也和表達式相一致。從表2和表3可以看出，x→y方向的格蘭杰因果系數(shù)總是非零，參數(shù)σ和β在一定范圍內(nèi)的變動對結(jié)果影響較小。這表明理想情況下使用3種核函數(shù)計算格蘭杰因果系數(shù)均可以得到正確的結(jié)論，對比計算結(jié)果可以發(fā)現(xiàn)當(dāng)使用合適的參數(shù)時，多項式核函數(shù)計算的因果系數(shù)更大，而且最大值對應(yīng)的參數(shù)與表達式的階數(shù)相一致。同時，高斯核函數(shù)與Sigmoid核函數(shù)對參數(shù)的敏感度較小，因而計算結(jié)果的均方差值較小。

實驗2驗證樣本不同時期腦電和心電信號的格蘭杰因果關(guān)系

通過分別計算清醒期和睡眠期不同樣本腦電和心電信號的格蘭杰因果系數(shù)，可以得到腦電和心電信號格蘭杰因果關(guān)系圖，使用3種核函數(shù)的結(jié)果分別如圖1和圖2所示。

圖1 清醒期腦電信號與心電信號格蘭杰因果系數(shù)棒圖 圖2 睡眠期腦電信號與心電信號格蘭杰因果系數(shù)棒圖

Fig.1 Granger causality index error bar chart of awake sample's EEG and ECG signals Fig. 2 Granger causality index error bar chart of sleep sample's EEG and ECG signals

使用統(tǒng)計產(chǎn)品與服務(wù)解決方案軟件(Statistical product and service solutions, SPSS)對清醒期的樣本均值進行樣本獨立t檢驗，即對兩樣本均數(shù)差別的顯著性進行檢驗。多項式核函數(shù)，高斯核函數(shù)和Sigmoid核函數(shù)對應(yīng)格蘭杰因果系數(shù)的顯著性水平p值分別為0.018，0.023和0.033，三者的顯著性水平均小于0.05，因此可以通過置信度為0.95的t檢驗。同樣使用SPSS對睡眠期樣本均值進行樣本獨立t檢驗。多項式核函數(shù)高斯核函數(shù)和Sigmoid對應(yīng)格蘭杰因果系數(shù)的顯著性水平更p值為0.018，0.006和0.008，三者的顯著性水平均小于0.05，因此可以通過置信度為0.95的t檢驗。

對比清醒期和睡眠期的因果關(guān)系顯著性水平，可以發(fā)現(xiàn)睡眠期的生理信號因果關(guān)系顯著性水平更為明顯。腦電信號對于心電信號有著較為明顯的因果關(guān)系，即腦電信號是心電信號的因，而心電信號不是腦電信號的因。

實驗3驗證樣本不同時期血壓和腦電信號的格蘭杰因果關(guān)系

通過分別計算清醒期和睡眠期不同樣本腦電和血壓信號的格蘭杰因果系數(shù)，可以得到腦電和血壓信號因果關(guān)系圖，使用3種核函數(shù)的結(jié)果分別如圖3和圖4所示。

圖3 清醒期腦電信號與血壓信號格蘭杰因果系數(shù)棒圖 圖4 睡眠期腦電信號與血壓信號格蘭杰因果系數(shù)棒圖

Fig.3 Granger causality index error bar chart of awake sample′s blood pressure and EEG signals Fig.4 Granger gausality index error bar chart of sleep sample′s blood pressure and EEG signals

使用SPSS對清醒期樣本均值進行樣本獨立t檢驗，多項式核函數(shù)，高斯核函數(shù)和Sigmoid核函數(shù)對應(yīng)格蘭杰因果系數(shù)的顯著性水平p值分別為0.001，0.34和0.030，多項式核函數(shù)和Sigmoid核函數(shù)對應(yīng)的顯著性水平均小于0.05，因此可以通過置信度為0.95的t檢驗，而高斯核函數(shù)則無法通過t檢驗，即無法利用高斯核函數(shù)來確定腦電血壓信號的因果關(guān)系；使用SPSS對睡眠期樣本均值進行樣本獨立t檢驗，多項式核函數(shù)高斯核函數(shù)和Sigmoid核函數(shù)對應(yīng)格蘭杰因果系數(shù)的顯著性水平p值分別為0.006，0.35和0.045，其中多項式核函數(shù)和高斯核函數(shù)對應(yīng)的顯著性水平均小于0.05，因此可以通過置信度為0.95的t檢驗，Sigmoid核函數(shù)則無法通過t檢驗，即無法利用Sigmoid核函數(shù)來確定腦電血壓信號的因果關(guān)系。

對比清醒期和睡眠期的樣本腦電血壓信號可以發(fā)現(xiàn)，由多項式核函數(shù)計算得到的格蘭杰因果關(guān)系均具有較高顯著性水平，即腦電信號是血壓信號的因，血壓信號不是腦電信號的因。

實驗4驗證樣本不同時期血壓和心電信號的格蘭杰因果關(guān)系

分別計算清醒期和睡眠期不同樣本血壓和心電信號的格蘭杰因果系數(shù)，可以得到血壓和心電信號因果關(guān)系圖，使用3種核函數(shù)的結(jié)果分別如圖5和圖6所示。

圖5 清醒期血壓信號與心電信號格蘭杰因果系數(shù)棒圖 圖6 睡眠期血壓信號與心電信號格蘭杰因果系數(shù)棒圖

Fig.5 Granger causality index error bar chart of awake sample′s blood pressure and ECG signals Fig.6 Granger causality index error bar chart of sleep sample′s blood pressure and EEG signals

使用SPSS對樣本均值進行樣本獨立t檢驗，多項式核函數(shù),高斯核函數(shù)和Sigmoid核函數(shù)對應(yīng)格蘭杰因果系數(shù)的顯著性水平p值為0.001，0.002和0.098。由于多項式核函數(shù)和高斯核函數(shù)對應(yīng)的顯著性水平均小于0.05，因此可以通過置信度為0.95的t檢驗，而Sigmoid核函數(shù)則無法通過t檢驗，即無法利用Sigmoid核函數(shù)來確定心電信號和血壓信號的因果關(guān)系。使用SPSS對樣本均值進行樣本獨立t檢驗，多項式核函數(shù)對應(yīng)格蘭杰因果系數(shù)的顯著性水平p值分別為0.002，0.003和0.008，三者的顯著性水平均小于0.05，因此可以通過置信度為0.95的t檢驗。

對比清醒期和睡眠期的因果關(guān)系系數(shù)，可以發(fā)現(xiàn)睡眠期樣本的格蘭杰因果關(guān)系系數(shù)顯著性水平更高。血壓信號對心電信號因果影響顯著性較高，即心電信號是血壓信號的因，血壓信號不是心電信號的因。

實驗5腦電采集樣本腦電和心電信號的格蘭杰因果關(guān)系

通過計算南京軍區(qū)總醫(yī)院臨床采樣正常樣本的腦電信號和心電信號的格蘭杰因果系數(shù)可以得到圖7。

圖7 3種核函數(shù)計算得到的正常樣本心電信號和腦電信號格蘭杰因果系數(shù)棒圖Fig.7 Granger causality index error bar chart of normal sample′s EEG and ECG signals

從圖7可知，使用SPSS對樣本均值進行樣本獨立t檢驗，多項式核函數(shù)高斯函數(shù)和Sigmoid核函數(shù)對應(yīng)格蘭杰因果系數(shù)的顯著性水平p值分別為0.001，0.001和0.002，三者的顯著性水平均小于0.05，因此可以通過置信度為0.95的t檢驗。

從實驗結(jié)果可以看出，正常樣本腦電信號對心電信號的格蘭杰影響較為顯著，心電信號對腦電信號的影響很小。比較3種核函數(shù)得到的計算結(jié)果可以發(fā)現(xiàn)，使用合適的參數(shù)時多項式核函數(shù)計算的結(jié)果顯著性水平較高，高斯核函數(shù)和Sigmoid核函數(shù)計算結(jié)果均方差較小。

3 結(jié)束語

基于非線性格蘭杰因果關(guān)系的睡眠腦電分析，本文利用核函數(shù)方法將線性格蘭杰因果關(guān)系擴展到非線性空間，從而可以更為準確地分析生理信號。通過對實驗結(jié)果的分析發(fā)現(xiàn)，心電信號、腦電信號和血壓信號之間存在著較為顯著的因果關(guān)系，腦電對心電信號的影響比心電對腦電信號的影響更為顯著，腦電對血壓信號的影響比血壓對腦電信號的影響更為顯著，血壓對心電信號的影響比心電對血壓信號的影響更為顯著。即腦電信號是心電信號的因，心電信號是血壓信號的因，反之則無法確定。此外，相對于樣本清醒期的數(shù)據(jù)分析結(jié)果，樣本睡眠期的這種顯著程度更為強烈，從而驗證了睡眠期的生理信號可以更為清晰地反映生理活動狀態(tài)。

[1] Fortune M H. State transitions between wake and sleep, and within the ultradian cycle, with focus on the link to neuronal activity[J]. Sleep Medicine, 2004, 8(6): 473-485.

[2] Saper C, Scammell T. The sleep switch: Hypothalamic control of sleep and wakefulness[J]. Trends Neurosci, 2001, 24(12): 726-731.

[3] Urrila A S, Paunio T, Palom E, et al. Sleep in adolescent depression: Physiological perspectives[J]. Acta Physiol, 2015, 213:758-777.

[4] 張莉, 何傳紅, 何為, 等. 基于典型相關(guān)分析與低通濾波的肌電偽跡去除 [J] .數(shù)據(jù)采集與處理, 2010, 25(2): 255-258.

Zhang Li, He Chuanhong, He Wei, et al. Method for removing EMG artifacts based on CCA and low-pass filtering [J]. Journal of Data Acquisition and Processing, 2010, 25(2): 255-258.

[5] Lee A H, Mojica J, Thomas R J, et al. Dissociation of EEG and autonomic sleep physiology[J]. Sleep Medicine, 2015, 16(2): 308-309.

[6] Schmidt E M, Linz B, Diekelmann S, et al. Effects of an interleukin-1 receptor antagonist on human sleep, sleep-associated memory consolidation, and blood monocytes[J]. Brain Behavior and Immunity, 2015, 47:178-185.

[7] Granger C W J. Investigating causal relations by econometric models and cross-spectral meth[J]. Econometrica, 1969, 37(3): 424-438.

[8] Marinazzo D, Pellicoro M, Stramaglia S. Kernel method for nonlinear granger causality[J]. Phys Rev Lett, 2008, 100(14): 144103.

[9] Schoelkopf B, Muandet K, Fukumizu K, et al. Computing functions of random variables via reproducing kernel Hilbert space representations[J]. Statistics and Computing, 2015, 25(4): 755-766.

[10] 劉柳 ，陶大程. Lasso問題的最新算法研究 [J] .數(shù)據(jù)采集與處理, 2015, 30(1): 35-46.

Liu Liu, Tao Dacheng. Review on recent mehtod of solving Lasso problem [J]. Journal of Data Acquisition and Processing, 2015, 30(1): 35-46.

[11] Marinazzo D， Pellicoro M, Stramaglias S. Kernel-Granger causality and the analysis of dynamical networks[J]. Physical Review E, 2008, 77(2)： 056215.

[12] Camps-Vallsa G, Martin-Guerreroa J D, Rojo-Alvarez J L, et al. Fuzzy sigmoid kernel for support vector classifiers[J]. Neurocomputing, 2004, 62(1): 501-506.

ApplicationofNonlinearGrangerCausalityinAnalysisofPhysiologicalSignalsDuringSleep

Du Peng1， Dai Jiafei2， Li Jin3， Wang Jun1， Hou Fengzhen4

(1．Image Processing and Image Communications Key Lab. , Nanjing University of Posts & Telecomm, Nanjing, 210003, China； 2．Nanjing General Hospital of Nanjing Military Command, Nanjing, 210002, China； 3．College of Physics and Information Technology, Shanxi Normal University, Xi′an， 710062, China；4．School of Science, China Pharmaceutical University, Nanjing， 210009, China)

A method based on the nolinear Granger causality is used to analyze sleep physiological signal. Polynomial kernel function, Gaussian kernel function and sigmoid kernel function are used to map the linear data in low dimensional input space into high dimensional feature space in which linear Granger method can be used to analyse the biomedical signals. The analysis results show that the causal effect of electrocardlogram (ECG) signals to electroencephalogram (EEG) signals, ECG signals to blood pressure signals and blood pressure signals to ECG signals are more significant than that of the opposite direction. In addition, the results of sleep subjects have more significant difference than that of normal subjects. The simulation results validate that the sleep physiological signal reflects the causality more objectively.

nonlinear Granger causality； sleep signal； kernel function

國家自然科學(xué)基金(61271082， 61401518)資助項目；江蘇省重點研發(fā)計劃(BE2015700)資助項目；江蘇省自然科學(xué)基金(BK20141432資助項目)；南京軍區(qū)南京總醫(yī)院基金(2014019)資助項目；南京市醫(yī)療衛(wèi)生科技項目(201503009)資助項目；中國藥科大學(xué)中央高?；究蒲袠I(yè)務(wù)費專項資金(FY2014LX0039)資助項目。

2015-11-18；

2016-03-23

TP911.73

A

杜朋(1991-)，男，碩士研究生，研究方向：醫(yī)學(xué)信號處理,E-mail:dupeng_graduate@163.com。

戴加飛(1981-)，男，醫(yī)師，研究方向：腦電圖方面對癲斷及鑒別診斷。

李錦(1973-)，女，博士，副教授，研究方向：信號與信息處理、超聲醫(yī)學(xué)電子。

王俊(1973-)，通信作者，男，博士，教授，博士生導(dǎo)師，研究方向：生物醫(yī)學(xué)信息處理，E-mail:wangj@njupt.edu.cn。

侯鳳貞(1980-)，通信作者，女，博士，副教授，研究方向：生物醫(yī)學(xué)信息處理，E-mail:houfz@cpu.edu.cn。