張梅 崔超 馬千里 干宗良? 王俊?
近年來,心腦血管等身體內部器官疾病的發病率不斷地增加,為了對身體內部器官的健康狀況進行準確的預報和診斷,需要行之有效的方法來采集和分析人體生物電信號.目前生物電信號采集和記錄技術有了很大進步.相比較于生物電信號采集技術的突飛猛進,生物電信號的分析技術的進展[1-4]卻不盡如人意.目前常用的生物電分析方法主要有時域分析[5-11]、頻域分析[6-9,12]、功率譜分析[13,14]、神經網絡[12,15-19]分析方法等,但大多數是停留在理論方面,沒有廣泛地應用于實踐.
本文提出了一種新的生物電信號分析方法——基于符號化部分互信息的多參數生物電信號的分析方法,這種方法用于分析生物電信號之間的耦合信息.分析兩個或兩個以上系統的耦合關系十分有益,以分析生物電信號為例,由于目前尚不能完全掌握這些系統,所以我們的分析局限于已記錄的生物電信號上.部分互信息可以反映這些信號間的相關性,這對于疾病的檢測和診斷具有重要意義.
如果直接對生物電信號原始序列進行部分互信息分析[20],由于噪聲等因素的影響,往往不能準確地對信號序列進行分析.而且和其他基于熵的算法一樣,基于部分互信息的耦合分析方法也存在不足,要求數據量較大,而且時間序列具有很好的平穩性,這些都給實際應用帶來困難.
針對以上問題,本文提出了使用符號化的部分互信息來對生物電信號進行耦合分析.符號化[21,22]是動力學系統研究的一個重要手段,其主要的意義是對噪聲的影響不敏感,得到的結論較為嚴格,而其應用的關鍵在于如何劃分符號區域,使得處理后的信號不會丟失時間序列的動力學特性.
對于一個離散隨機變量X,有概率分布{px},其信息熵定義如下:

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

I(X,Y)反映了隨機變量X和Y間的相關性,其中聯合熵H(X,Y)由XY的聯合分布{pxy}計算得到.
部分互信息I(X,Y|Z)是XY相交但不包括Z的部分,其定義式如下:

I(X,Y|Z)表示在已知Z的情況下Y提供的關于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的函數時.隨機變量間的互信息越大,其耦合程度越大.
為了克服噪聲等因素的影響,本文提出了使用符號化的部分互信息來對生物電信號進行耦合分析.原始序列的符號化定義如下[23,24]:

其中μ1是原始序列中大于等于零的取樣信號的平均值,μ2是小于零的平均值,取a=0.05.符號化主要的意義是對噪聲的影響不敏感,得到的結論較為嚴格,而其應用的關鍵在于如何劃分符號區域,使得處理后的信號不會丟失時間序列的動力學特性.當a取0.04到0.07時都可得到原信號的大尺度信息,如果a取小于0.04或大于0.07時,會得不到較為合理的統計特性.這是因為在符號化的過程中,如果a的值過大或過小,會導致細節信息的丟失.故取a=0.05,既能去掉原信號噪聲影響,又能較好地捕捉信號中的動態信息.
為了取得更好的分析效果,本文將原始序列符號化與部分互信息分析法結合,即先對多參數的生物電原始信號腦電X,心電Y,肌電Z進行符號化及編碼處理,得到符號化序列˙X,˙Y,˙Z,再計算其部分互信息來獲得多參量生物電信號的耦合信息,符號化的部分互信息定義如下:

本文使用的睡眠數據來自PhysioBank的MITBIH Polysomnographic Database.該庫中的記錄是多參數睡眠數據,包括1導EEG(electroencephalogram)信號,1導ECG(electrocardiosignal)信號,1導EOG(eyectro-oculogram)信號,1導EMG(electromyographic signal)信號等多導睡眠信號,記錄長度為6 h,數據采樣頻率250 Hz,每份記錄的數據都附帶著以30 s為一個分期的注釋信息,本文的實驗結論判定依此注釋為參考.
本文所用數據,分別采用了15組受試者睡眠和清醒期的多參數生物電信號數據中的1導EEG(C3-O1)腦電信號、1導ECG心電信號、1導EMG肌電信號,提取其中的清醒期和NREM睡眠I期的若干組信號,分別記為樣本Sleep和Weak.
根據本文提出的算法,我們首先根據(4)式對生物電原始信號各參量EEG,ECG,EMG進行相應的符號化及編碼處理,再根據(5)式計算生物電信號樣本Sleep和Weak各參量間的部分互信息I(EEG,ECG|EMG),從而分別獲得睡眠期和清醒期生物電信號的耦合信息.最后使用spss統計分析軟件對計算結果進行了假設檢驗,驗證了該算法的有效性.以生物電信號的耦合信息作為參數,我們可以判斷該生物機制是處在活躍的還是消極的狀態,這對于一些身體內部器官健康狀況的診斷具有重要意義.
3.3.1 研究多參量生物電信號耦合程度與數據長度L的關系
對受訓者樣本15組Sleep和15組Wake中的每組生物電信號序列分別取數據長度L=500,1000,1500,2000,2500,3000,并對每個時間序列進行符號化及編碼處理,分別計算每個個體的部分互信息,并分別對睡眠期部分互信息和清醒期部分互信息進行平均,得到部分互信息與數據長度L的關系如圖2所示.

圖2 部分互信息熵與數據長度L的關系
分析圖2睡眠期與清醒期的部分互信息曲線可知:
1)在數據長度相等的情況下,清醒期的生物電信號耦合程度要高于睡眠期;
2)當數據長度由0增至2000的過程中,睡眠期與清醒期的耦合程度分別呈遞增趨勢,數據長度大于2000時,耦合程度趨于平穩,基本不變.
從數據的精準性考慮,數據長度越大,其統計概率分布越接近實際分布,相應的準確性越高.但從計算量以及計算速度來講,數據長度越小,那么算法速度也越快,綜合我們的實驗結果,取數據長度L=2000時,既能兼顧處理速度,又可以保證實驗精度.
3.3.2 研究多參量生物電信號耦合程度與編碼長度N的關系
對受訓者樣本15組Sleep和15組Weak中的每組生物電信號序列分別取編碼長度N=1,2,3,4,5,6,并對每個時間序列進行符號化及編碼處理,分別計算每個個體的部分互信息,并分別對睡眠期部分互信息和清醒期部分互信息進行平均,得到生物電信號耦合程度與編碼長度N的關系如圖3所示.

圖3 部分互信息熵與編碼長度N的關系
分析圖3部分互信息與編碼長度N的關系可知:
1)在編碼長度N從1到6的遞增過程中,生物電信號的部分互信息都成增大趨勢;
2)清醒期的耦合程度明顯高于睡眠期的,且在N=6時清醒期與睡眠期耦合程度差異最顯著.
考慮到實驗效果的明顯性,N應該取大些;考慮到算法的復雜度會影響實驗處理速度以及在臨床應用上的實時性,N應該取小一些.綜合考慮上述因素,N取6.
3.3.3 研究在數據長度L=2000,編碼長度N=8條件下青老年能量耗散的差異
對受訓者樣本15組Sleep和15組Weak中的每組生物電信號序列分別取編碼長度N=6及數據長度L=2000,計算每個個體的部分互信息,得到生物電耦合程度與睡眠期清醒期的關系如圖4所示.
分析圖4睡眠期與清醒期的生物電信號耦合程度對比可知:
1)清醒期生物電信號耦合程度明顯高于睡眠期的,生物電信號的部分互信息熵可以作為衡量一個過程是否處于積極有序狀態的參數,耦合程度越高,說明該物理過程越有序;
2)上圖中0—14號生物電信號個體間的生物電信號耦合程度有很大差異,不同個體處于不同的年齡段和具有不同的體重引起的.
3.3.4 統計分析與假設檢驗
為進一步驗證本文計算結果的準確性及算法的有效性,使用spss統計分析軟件對計算結果進行了假設檢驗,主要方法如下.
對清醒期與睡眠期生物電信號耦合程度的差異顯著性進行假設檢驗,分別將睡眠期和清醒期的部分互信息值記作樣本S和W,使用spss對兩組部分互信息數據S,W進行獨立樣本T檢驗,結果如表1所示.

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

表1 清醒期與睡眠期生物電信號耦合程度的差異顯著性假設檢驗
本文關注的是準確地對清醒期與睡眠期生物電信號進行判斷和辨別,根據表1分析清醒期與睡眠期生物電信號耦合程度的差異顯著性,假設清醒期與睡眠期的均值相等:
1)由表1可知,在假設方差相等的情況下,自由度(d f)為28,查找t值表可知理論t值t(d f)0.05=t(28)=2.048,而樣本t值為3.291,大于t(28),且Sig(雙側)=0.003<0.05,假設不成立,所以清醒期與睡眠期生物電信號耦合程度差異顯著;
2)由表1可知,在假設方差不相等的情況下,自由度(d f)為23,查找t值表可知理論t值t(d f)0.05=t(23)=2.069,而樣本t值為3.291,大于t(23),且Sig(雙側)=0.003<0.05,假設不成立,所以清醒期與睡眠期生物電信號耦合程度差異顯著,證明該算法可以有效地對清醒期與睡眠期生物電信號進行區分.
以上實驗結果證明,使用符號化的部分互信息算法可以很好地對多參數生物電信號的耦合程度進行分析,以生物電信號各參量之間的耦合程度作為參數,能夠有效地辨別該生理機制是否是活躍的,這對于身體內部器官健康狀況的診斷具有重要意義.
1)提出了一種新的部分互信息的計算方法——符號化部分互信息,經驗證該算法能夠很好地分析離散隨機變量間的耦合信息.
2)使用符號化的部分互信息來計算生物電時間序列間的耦合程度,實現了利用統計學來量化研究生命現象.
3)基于該算法分別對睡眠期和清醒期的EEG,ECG,EMG多參量生物電信號進行分析.分析結果表明,清醒期的生物電信號的耦合程度顯著高于睡眠期的.并進行了假設驗證,證明睡眠期和清醒期的生物電耦合信息具有顯著差異,表明耦合程度可以作為衡量一個物理過程是否處于積極狀態以及睡眠分期的參數.
[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)[孟慶芳,周衛東,陳月輝,彭玉華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]