張梅 崔超 馬千里 干宗良? 王俊?

1 引言

近年來,心腦血管等身體內(nèi)部器官疾病的發(fā)病率不斷地增加,為了對身體內(nèi)部器官的健康狀況進行準確的預(yù)報和診斷,需要行之有效的方法來采集和分析人體生物電信號.目前生物電信號采集和記錄技術(shù)有了很大進步.相比較于生物電信號采集技術(shù)的突飛猛進,生物電信號的分析技術(shù)的進展[1-4]卻不盡如人意.目前常用的生物電分析方法主要有時域分析[5-11]、頻域分析[6-9,12]、功率譜分析[13,14]、神經(jīng)網(wǎng)絡(luò)[12,15-19]分析方法等,但大多數(shù)是停留在理論方面,沒有廣泛地應(yīng)用于實踐.

本文提出了一種新的生物電信號分析方法——基于符號化部分互信息的多參數(shù)生物電信號的分析方法,這種方法用于分析生物電信號之間的耦合信息.分析兩個或兩個以上系統(tǒng)的耦合關(guān)系十分有益,以分析生物電信號為例,由于目前尚不能完全掌握這些系統(tǒng),所以我們的分析局限于已記錄的生物電信號上.部分互信息可以反映這些信號間的相關(guān)性,這對于疾病的檢測和診斷具有重要意義.

如果直接對生物電信號原始序列進行部分互信息分析[20],由于噪聲等因素的影響,往往不能準確地對信號序列進行分析.而且和其他基于熵的算法一樣,基于部分互信息的耦合分析方法也存在不足,要求數(shù)據(jù)量較大,而且時間序列具有很好的平穩(wěn)性,這些都給實際應(yīng)用帶來困難.

針對以上問題,本文提出了使用符號化的部分互信息來對生物電信號進行耦合分析.符號化[21,22]是動力學系統(tǒng)研究的一個重要手段,其主要的意義是對噪聲的影響不敏感,得到的結(jié)論較為嚴格,而其應(yīng)用的關(guān)鍵在于如何劃分符號區(qū)域,使得處理后的信號不會丟失時間序列的動力學特性.

2 基本原理

對于一個離散隨機變量X,有概率分布{px},其信息熵定義如下：

兩個隨機變量X和Y的互信息定義為

I(X,Y)反映了隨機變量X和Y間的相關(guān)性,其中聯(lián)合熵H(X,Y)由XY的聯(lián)合分布{pxy}計算得到.

部分互信息I(X,Y|Z)是XY相交但不包括Z的部分,其定義式如下：

I(X,Y|Z)表示在已知Z的情況下Y提供的關(guān)于X的平均信息量,其物理意義如圖1所示.

圖1 部分互信息I(X,Y|Z)(間斷條紋部分)

對于相同的Z部分信息熵是對稱的,即I(X,Y|Z)=I(Y,X|Z),同樣有0≤I(X,Y|Z),當且僅當XY相互獨立時取零,特別是當X或Y是Z的函數(shù)時.隨機變量間的互信息越大,其耦合程度越大.

為了克服噪聲等因素的影響,本文提出了使用符號化的部分互信息來對生物電信號進行耦合分析.原始序列的符號化定義如下[23,24]：

其中μ1是原始序列中大于等于零的取樣信號的平均值,μ2是小于零的平均值,取a=0.05.符號化主要的意義是對噪聲的影響不敏感,得到的結(jié)論較為嚴格,而其應(yīng)用的關(guān)鍵在于如何劃分符號區(qū)域,使得處理后的信號不會丟失時間序列的動力學特性.當a取0.04到0.07時都可得到原信號的大尺度信息,如果a取小于0.04或大于0.07時,會得不到較為合理的統(tǒng)計特性.這是因為在符號化的過程中,如果a的值過大或過小,會導(dǎo)致細節(jié)信息的丟失.故取a=0.05,既能去掉原信號噪聲影響,又能較好地捕捉信號中的動態(tài)信息.

為了取得更好的分析效果,本文將原始序列符號化與部分互信息分析法結(jié)合,即先對多參數(shù)的生物電原始信號腦電X,心電Y,肌電Z進行符號化及編碼處理,得到符號化序列˙X,˙Y,˙Z,再計算其部分互信息來獲得多參量生物電信號的耦合信息,符號化的部分互信息定義如下：

3 基于符號化部分互信息的生物電信號耦合分析

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

本文使用的睡眠數(shù)據(jù)來自PhysioBank的MITBIH Polysomnographic Database.該庫中的記錄是多參數(shù)睡眠數(shù)據(jù),包括1導(dǎo)EEG(electroencephalogram)信號,1導(dǎo)ECG(electrocardiosignal)信號,1導(dǎo)EOG(eyectro-oculogram)信號,1導(dǎo)EMG(electromyographic signal)信號等多導(dǎo)睡眠信號,記錄長度為6 h,數(shù)據(jù)采樣頻率250 Hz,每份記錄的數(shù)據(jù)都附帶著以30 s為一個分期的注釋信息,本文的實驗結(jié)論判定依此注釋為參考.

本文所用數(shù)據(jù),分別采用了15組受試者睡眠和清醒期的多參數(shù)生物電信號數(shù)據(jù)中的1導(dǎo)EEG(C3-O1)腦電信號、1導(dǎo)ECG心電信號、1導(dǎo)EMG肌電信號,提取其中的清醒期和NREM睡眠I期的若干組信號,分別記為樣本Sleep和Weak.

3.2 實驗方法

根據(jù)本文提出的算法,我們首先根據(jù)(4)式對生物電原始信號各參量EEG,ECG,EMG進行相應(yīng)的符號化及編碼處理,再根據(jù)(5)式計算生物電信號樣本Sleep和Weak各參量間的部分互信息I(EEG,ECG|EMG),從而分別獲得睡眠期和清醒期生物電信號的耦合信息.最后使用spss統(tǒng)計分析軟件對計算結(jié)果進行了假設(shè)檢驗,驗證了該算法的有效性.以生物電信號的耦合信息作為參數(shù),我們可以判斷該生物機制是處在活躍的還是消極的狀態(tài),這對于一些身體內(nèi)部器官健康狀況的診斷具有重要意義.

3.3 實驗結(jié)果及分析

3.3.1 研究多參量生物電信號耦合程度與數(shù)據(jù)長度L的關(guān)系

對受訓(xùn)者樣本15組Sleep和15組Wake中的每組生物電信號序列分別取數(shù)據(jù)長度L=500,1000,1500,2000,2500,3000,并對每個時間序列進行符號化及編碼處理,分別計算每個個體的部分互信息,并分別對睡眠期部分互信息和清醒期部分互信息進行平均,得到部分互信息與數(shù)據(jù)長度L的關(guān)系如圖2所示.

圖2 部分互信息熵與數(shù)據(jù)長度L的關(guān)系

分析圖2睡眠期與清醒期的部分互信息曲線可知：

1)在數(shù)據(jù)長度相等的情況下,清醒期的生物電信號耦合程度要高于睡眠期;

2)當數(shù)據(jù)長度由0增至2000的過程中,睡眠期與清醒期的耦合程度分別呈遞增趨勢,數(shù)據(jù)長度大于2000時,耦合程度趨于平穩(wěn),基本不變.

從數(shù)據(jù)的精準性考慮,數(shù)據(jù)長度越大,其統(tǒng)計概率分布越接近實際分布,相應(yīng)的準確性越高.但從計算量以及計算速度來講,數(shù)據(jù)長度越小,那么算法速度也越快,綜合我們的實驗結(jié)果,取數(shù)據(jù)長度L=2000時,既能兼顧處理速度,又可以保證實驗精度.

3.3.2 研究多參量生物電信號耦合程度與編碼長度N的關(guān)系

對受訓(xùn)者樣本15組Sleep和15組Weak中的每組生物電信號序列分別取編碼長度N=1,2,3,4,5,6,并對每個時間序列進行符號化及編碼處理,分別計算每個個體的部分互信息,并分別對睡眠期部分互信息和清醒期部分互信息進行平均,得到生物電信號耦合程度與編碼長度N的關(guān)系如圖3所示.

圖3 部分互信息熵與編碼長度N的關(guān)系

分析圖3部分互信息與編碼長度N的關(guān)系可知：

1)在編碼長度N從1到6的遞增過程中,生物電信號的部分互信息都成增大趨勢;

2)清醒期的耦合程度明顯高于睡眠期的,且在N=6時清醒期與睡眠期耦合程度差異最顯著.

考慮到實驗效果的明顯性,N應(yīng)該取大些;考慮到算法的復(fù)雜度會影響實驗處理速度以及在臨床應(yīng)用上的實時性,N應(yīng)該取小一些.綜合考慮上述因素,N取6.

3.3.3 研究在數(shù)據(jù)長度L=2000,編碼長度N=8條件下青老年能量耗散的差異

對受訓(xùn)者樣本15組Sleep和15組Weak中的每組生物電信號序列分別取編碼長度N=6及數(shù)據(jù)長度L=2000,計算每個個體的部分互信息,得到生物電耦合程度與睡眠期清醒期的關(guān)系如圖4所示.

分析圖4睡眠期與清醒期的生物電信號耦合程度對比可知：

1)清醒期生物電信號耦合程度明顯高于睡眠期的,生物電信號的部分互信息熵可以作為衡量一個過程是否處于積極有序狀態(tài)的參數(shù),耦合程度越高,說明該物理過程越有序;

2)上圖中0—14號生物電信號個體間的生物電信號耦合程度有很大差異,不同個體處于不同的年齡段和具有不同的體重引起的.

3.3.4 統(tǒng)計分析與假設(shè)檢驗

為進一步驗證本文計算結(jié)果的準確性及算法的有效性,使用spss統(tǒng)計分析軟件對計算結(jié)果進行了假設(shè)檢驗,主要方法如下.

對清醒期與睡眠期生物電信號耦合程度的差異顯著性進行假設(shè)檢驗,分別將睡眠期和清醒期的部分互信息值記作樣本S和W,使用spss對兩組部分互信息數(shù)據(jù)S,W進行獨立樣本T檢驗,結(jié)果如表1所示.

圖4 睡眠期與清醒期的生物電信號耦合程度對比

表1 清醒期與睡眠期生物電信號耦合程度的差異顯著性假設(shè)檢驗

本文關(guān)注的是準確地對清醒期與睡眠期生物電信號進行判斷和辨別,根據(jù)表1分析清醒期與睡眠期生物電信號耦合程度的差異顯著性,假設(shè)清醒期與睡眠期的均值相等：

1)由表1可知,在假設(shè)方差相等的情況下,自由度(d f)為28,查找t值表可知理論t值t(d f)0.05=t(28)=2.048,而樣本t值為3.291,大于t(28),且Sig(雙側(cè))=0.003＜0.05,假設(shè)不成立,所以清醒期與睡眠期生物電信號耦合程度差異顯著;

2)由表1可知,在假設(shè)方差不相等的情況下,自由度(d f)為23,查找t值表可知理論t值t(d f)0.05=t(23)=2.069,而樣本t值為3.291,大于t(23),且Sig(雙側(cè))=0.003＜0.05,假設(shè)不成立,所以清醒期與睡眠期生物電信號耦合程度差異顯著,證明該算法可以有效地對清醒期與睡眠期生物電信號進行區(qū)分.

以上實驗結(jié)果證明,使用符號化的部分互信息算法可以很好地對多參數(shù)生物電信號的耦合程度進行分析,以生物電信號各參量之間的耦合程度作為參數(shù),能夠有效地辨別該生理機制是否是活躍的,這對于身體內(nèi)部器官健康狀況的診斷具有重要意義.

4 結(jié)束語

1)提出了一種新的部分互信息的計算方法——符號化部分互信息,經(jīng)驗證該算法能夠很好地分析離散隨機變量間的耦合信息.

2)使用符號化的部分互信息來計算生物電時間序列間的耦合程度,實現(xiàn)了利用統(tǒng)計學來量化研究生命現(xiàn)象.

3)基于該算法分別對睡眠期和清醒期的EEG,ECG,EMG多參量生物電信號進行分析.分析結(jié)果表明,清醒期的生物電信號的耦合程度顯著高于睡眠期的.并進行了假設(shè)驗證,證明睡眠期和清醒期的生物電耦合信息具有顯著差異,表明耦合程度可以作為衡量一個物理過程是否處于積極狀態(tài)以及睡眠分期的參數(shù).

[1]Fang X L,Jiang Z L 2007 Acta Phys.Sin.56 7330(in Chinese)[方小玲,姜宗來2007物理學報56 7330]

[2]Meng Q F,Zhou W D,Chen Y H,Peng Y H 2010 Acta Phys.Sin.59 123(in Chinese)[孟慶芳,周衛(wèi)東,陳月輝,彭玉華2010物理學報59 123]

[3]Ma Q L,Bian C H,Wang J 2010 Acta Phys.Sin.59 4480(in Chinese)[馬千里,卞春華,王俊2010物理學報59 4480]

[4]Bian H R,Wang J,Han C X,Deng B,Wei X L,Che Y Q 2011 Acta Phys.Sin.60 118701(in Chinese)[邊洪瑞,王江,韓春曉,鄧斌,魏熙樂,車艷秋2011物理學報60 118701]

[5]Wang J,Ma Q L 2008 Chin.Phys.B 17 4424

[6]Hsu W Y 2012 Int.J.Neural Syst.22 51

[7]Nevado-Holgado A J,Marten F,Richardson M P,Terry J R 2012 Neuro Image 59 2374

[8]Petrantonakis P C,Hadjileontiadis L J 2012 IEEE Trans.Signal Proces.60 2604

[9]Thatcher R W 2012 Dev.Neuropsychol.37 476

[10]Wang J,Yu Z F 2012 Chin.Phys.B 21 018702

[11]Wang J,Zhao D Q 2012 Chin.Phys.B 21 028703

[12]Musselman M,Djurdjanovic D 2012 Exp.Syst.Appl.39 11413

[13]Shao S Y,Shen K Q,Yu K,Wilder-Smith E P V,Li X P 2012 Clin.Neurophysiol.123 2042

[14]Tarokh L,van Reen E,Acebo C,Le Bourgeois M,Seifer R,Fallone G,Carskadon M A 2012 Alcohol.-Clin.Exp.Res.36 1530

[15]Wang R F,Zhang J H,Zhang Y,Wang X Y 2012 Biomed.Signal Proc.Control 7 490

[16]Orhan U,Hekim M,Ozer M 2012 J.Med.Syst.36 2219

[17]Acharya U R,Sree S V,Alvin A P C,Suri J S 2012 Exp.Syst.Appl.39 9072

[18]Siuly S,Li Y 2012 IEEE Trans.Neur.Syst.Reh.Eng.20 526

[19]Acharya U R,Molinari F,Sree S V,Chattopadhyay S,Ng K H,Suri J S 2012 Biomed.Signal Proc.Control.7 401

[20]Frenzel S,Pompe B 2007 Phys.Rev.Lett.99 204101

[21]Staniek M,Lehnertz K 2008 Phys.Rev.Lett.100 158101

[22]Li J,Ning X B 2006 Phys.Rev.E 73 052902

[23]Wessel N,Ziehmann C,Kurths J,Meyerfeldt U,Schirdewan A,Voss A 2000 Phys.Rev.E 61 733

[24]Shen W,Wang J 2011 Acta Phys.Sin.60 118702(in Chinese)[沈韡,王俊2011物理學報60 118702]