基于ICEEMD-PE的腦血氧降噪方法研究

吳 凱 張 歡 劉 燕 戴亞康

（1.中國科學(xué)院大學(xué) 北京 101407）（2.中國科學(xué)院蘇州生物醫(yī)學(xué)工程技術(shù)研究所 蘇州 215163）（3.南京理工大學(xué)電子工程與光電技術(shù)學(xué)院 南京 210094）

1 引言

人體腦血氧在臨床治療和診斷中起著重要的作用，例如腦外傷［1］和心腦血管疾病的治療和檢測(cè)等。腦血氧信號(hào)的檢測(cè)是通過組織中血紅蛋白等對(duì)近紅外光的吸收，得到光密度的變化情況，再結(jié)合修正的朗伯-比爾定律計(jì)算得出組織內(nèi)的氧合血紅蛋白（oxygenated hemoglobin，oxy-Hb）和脫氧血紅蛋白（deoxy-genated hemoglobin，deoxy-Hb）的濃度變化即人體腦血氧信號(hào)［2］。

然而，由于采集的腦血氧信號(hào)十分微弱，因此在測(cè)量時(shí)非常容易受到外界的干擾。這些干擾主要有被測(cè)對(duì)象運(yùn)動(dòng)和呼吸等活動(dòng)而產(chǎn)生的低頻干擾，基線漂移［3］；采樣過程中儀器或環(huán)境的磁場(chǎng)等產(chǎn)生的高頻干擾［4］，高頻噪聲。以上的干擾嚴(yán)重影響了信號(hào)的精確度，致使結(jié)果出現(xiàn)較大的誤差。

目前針對(duì)血氧的去噪方法主要有最小均方誤差（LMS）自適應(yīng)算法，小波變換和經(jīng)驗(yàn)?zāi)B(tài)分解（Empirical Mode Decomposition，EMD）［5］等。其中，LMS自適應(yīng)算法雖然結(jié)構(gòu)較為簡(jiǎn)單、計(jì)算量小且穩(wěn)定性能好，但固定步長(zhǎng)導(dǎo)致跟蹤速率、收斂速率和失調(diào)噪聲之間相互矛盾［6］；小波變換可以對(duì)空間和頻率進(jìn)行多尺度細(xì)分分析，但要事先選定小波基函數(shù)、分解層次等［7］；EMD則根據(jù)信號(hào)本身特點(diǎn)，自適應(yīng)產(chǎn)生基函數(shù)將信號(hào)分解成一系列本征模態(tài)函數(shù)（Intrinsic Mode Function，IMF）［8］，解決了小波變換對(duì)基函數(shù)的依賴問題，但依然存在模態(tài)混疊現(xiàn)象［9］。針對(duì)EMD的不足，Huang提出的總體經(jīng)驗(yàn)?zāi)B(tài)分解（Ensemble Empirical Mode Decomposition，EEMD）［10］在信號(hào)中添加不同的高斯白噪聲，在多次EMD分解后將得到IMF進(jìn)行平均，有效解決了問題。然而，在實(shí)驗(yàn)有限次的情況下，添加的高斯白噪聲沒有完全抵消，重構(gòu)分量中存在殘留噪聲［11］。因此，提出了一種基于EEMD的具有自適應(yīng)白噪聲的完全總體經(jīng)驗(yàn)?zāi)B(tài)分解（Complete Ensemble Em?pirical Mode Decomposition with Adaptive Noise，CEEMDAN）［12］將殘留噪聲基本去除，但在分解的前期仍然存在著虛假模態(tài)問題，使得信號(hào)出現(xiàn)延遲現(xiàn)象，影響了信號(hào)的重構(gòu)。

為了解決上述問題，一種新的方法——改進(jìn)的總體完備經(jīng)驗(yàn)?zāi)B(tài)分解（Improved complete ensem?ble EMD，ICEEMD）［13］被提出，利用各模態(tài)當(dāng)前的局部均值重新定義真實(shí)模態(tài)，成功解決了CEEM?DAN面臨的問題。然而，無法對(duì)ICEEMD在分解信號(hào)得到多個(gè)IMF進(jìn)一步判斷，往往根據(jù)經(jīng)驗(yàn)進(jìn)行信號(hào)重構(gòu)，存在極大的誤差。因此，本文引入一種衡量空間復(fù)雜度的算法，排列熵。

排列熵（Permutation Entropy，PE）［14］對(duì)時(shí)間非常敏感，可獲得較高的分表率，而且計(jì)算方法簡(jiǎn)單，抗噪聲能力強(qiáng)。因此，本文將ICEEMD和排列熵相結(jié)合的方法去除信號(hào)中存在的基線漂移和高頻噪聲。

2 算法原理與分析

2.1 ICEEMD算法

記原始信號(hào)為 x（n），算子 Ek（.）和 M（.），Ek（.）表示EMD分解得到的第k個(gè)模態(tài)；M（.）則滿足M（x（n））=x（n）-E1（x（n）），ICEEMD 分解得到的得到的模態(tài)分量個(gè)數(shù)。其分解步驟如下：

1）對(duì)x（n）添加高斯白噪聲生成新的信號(hào)，wi表示第i組添加的高斯白噪聲，如下式：

其中，i=1，2，…，I，βk使得每個(gè)分解階段選擇合適的信噪比。

重復(fù)上述步驟便可得到原始信號(hào)所有的模態(tài)分量。

2.2 排列熵算法

將ICEEMD分解得到的模態(tài)分量IMFk(t)，k=1，2…，K表示模態(tài)分量個(gè)數(shù)，t=1，2，…，n，進(jìn)行空間重構(gòu)，得到矩陣如下：

其中N=n-(m-1)τ，m為嵌入維數(shù)，τ為延遲時(shí)間。矩陣中每行都可以看成是一個(gè)重構(gòu)分量，根據(jù)數(shù)值大小，將其按升序進(jìn)行排序，得到新的重構(gòu)分量，如下式表示：

其中，每個(gè)元素滿足如下關(guān)系：

式中j=2，3，…，m，因此 IMFk重構(gòu)所組成的矩陣中，每行都含有一組符號(hào)序列S(q)=(i1。i2。…。im)，式中q=1，2，…，N，且N ≤m！種不同排列方式，不同的序列出現(xiàn)的概率為P，得到IMFk分量的排列熵為

為了方便，將Hp進(jìn)行歸一化處理，即：

根據(jù)上式計(jì)算可得出各IMF分量的排列熵，HpE越大表明該分量越隨機(jī)；反之，則越規(guī)則。選取閾值剔除高頻噪聲和基線漂移，得到有用分量IMFm～I(xiàn)MFn，將其重構(gòu)產(chǎn)生去噪后的信號(hào)x'，表達(dá)式為

2.3 性能對(duì)比實(shí)驗(yàn)

為了驗(yàn)證算法的優(yōu)越性，在這里設(shè)計(jì)了一個(gè)對(duì)比實(shí)驗(yàn)。設(shè)原始信號(hào)S分別由信號(hào)S1、S2所示組成，如下式：

將原始信號(hào)分別進(jìn)行EMD分解、EEMD分解、CEEMDAN分解和ICEEMD分解，然后觀察其IMF分量。

圖1 原始信號(hào)S

模擬信號(hào)的分解實(shí)驗(yàn)結(jié)果如圖2所示，從圖（a）可以看出，EMD分解后得到的IMF分量中IMF1、IMF2和IMF3明顯是由多種不同時(shí)間尺度的信號(hào)組成，即存在模態(tài)混跌現(xiàn)象；從圖（b）可以看出，EEMD分解后的IMF分量可以很好地將兩種不同頻率、不同幅值的信號(hào)分解出來，分別為IMF1和IMF2，但是存在IMF3～I(xiàn)MF9的噪聲殘留；從圖（c）可以看出，CEEMDAN分解后的IMF分量IMF1和IMF3之間明顯存在虛假模態(tài)IMF2；從圖（d）可以看出，ICEEMD可以將兩個(gè)不同的模擬信號(hào)很完美地分解出來。因此，ICEEMD可以將非線性非平穩(wěn)的時(shí)頻信號(hào)很好的分解為一系列不同時(shí)間尺度的IMF分量。

圖2 分解得到IMF分量

3 實(shí)驗(yàn)與結(jié)果

3.1 實(shí)驗(yàn)步驟

利用自制的設(shè)備對(duì)腦血氧信號(hào)進(jìn)行采樣，即氧合血紅蛋白和脫氧血紅蛋白的濃度變化信號(hào)，之后按本文提出的算法進(jìn)行處理，驗(yàn)證方法的效果。圖3為設(shè)備的電路結(jié)構(gòu)框圖。200組的高斯白噪聲，標(biāo)準(zhǔn)差為0.2，最大的迭代次數(shù)為500次。

3）將通過ICEEMD分解而產(chǎn)生的兩組IMF分量引入排列熵計(jì)算，本文選擇的m=7，τ=1。

4）多次重復(fù)實(shí)驗(yàn)，統(tǒng)計(jì)出實(shí)驗(yàn)結(jié)果。根據(jù)排列熵的結(jié)算結(jié)果，選擇合適的閾值，進(jìn)行信號(hào)重構(gòu)。

3.2 實(shí)驗(yàn)結(jié)果

根據(jù)上述的實(shí)驗(yàn)步驟，得到oxy-Hb和de?oxy-Hb的原始濃度變化信息和經(jīng)過ICEEMD分解后的IMF分量信息，如圖4所示；將得到的IMF分量引入PE計(jì)算，得到的結(jié)果如表1所示。

圖3 腦血氧檢測(cè)電路結(jié)構(gòu)圖

圖4 血紅蛋白濃度變化的原始信號(hào)和經(jīng)過ICEEMD分解后得到的IMF分量

電路中，微處理器選擇美國德州儀器生產(chǎn)的CC250芯片，光源LED選擇EPITEX公司生產(chǎn)的SMT 660/910雙波長(zhǎng)LED，檢測(cè)器采用光敏管OPT101芯片。首先，微處理器輸出PWM波控制不同波長(zhǎng)的LED輪流發(fā)光射入大腦，以距離光源3cm的光敏管將出射光信號(hào)轉(zhuǎn)化為電信號(hào)；然后，經(jīng)過電流電壓轉(zhuǎn)換，低通濾波，比例放大，AD轉(zhuǎn)換等一系列處理后，由CC2540通過USB和BLE4.0兩種方式傳入上位機(jī)，并對(duì)數(shù)據(jù)進(jìn)行處理。具體的處理步驟如下：

1）利用近紅外光譜技術(shù)，以10Hz的采樣頻率采集腦部入射光和出射光的光密度信息，將結(jié)果做集成平均運(yùn)算。然后，利用修正的朗伯比爾定律計(jì)算出oxy-Hb和deoxy-Hb的濃度變化信息。

2）分別對(duì)oxy-Hb和deoxy-Hb進(jìn)行ICEEMD分解，得到一系列IMF分量。其中，加入試驗(yàn)次數(shù)為

表1 ICEEMD分解得到IMF分量的排列熵

3.3 結(jié)果分析

根據(jù)上述實(shí)驗(yàn)進(jìn)行多次實(shí)驗(yàn)，經(jīng)過統(tǒng)計(jì)計(jì)算后，這里將PE值超過0.5的IMF分量視為高頻噪聲；PE值低于0.2的IMF分量視為低頻的基線漂移。所以，選在該范圍內(nèi)的IMF分量進(jìn)行信號(hào)重構(gòu)，得到降噪后的信號(hào)。

圖5、6分別為oxy-Hb和deoxy-Hb的原始濃度變化信號(hào)及其對(duì)應(yīng)的頻譜信息和降噪之后的濃度變化信號(hào)和頻譜信息的對(duì)比。從兩張圖的（a）和（c）可以直觀地看出去噪后的信號(hào)光滑度更好，且從（b）和（d）中可以看出，高于10Hz的高頻噪聲被有效地抑制了。

圖5 oxy-Hb的濃度變化情況

圖6 deoxy-Hb的濃度變化情況

基線漂移主要是影響信號(hào)的穩(wěn)定性，但是在人體的血液中，血紅蛋白的總體濃度趨于穩(wěn)定［15～16］。分別計(jì)算出去噪前后的血紅蛋白的變化量，結(jié)果如圖7所示。從圖中可以看出，去噪之后的血紅蛋白濃度變化更為穩(wěn)定，很好地抑制了基線漂移問題。

圖7 降噪前與降噪后的血紅蛋白濃度變化對(duì)比

4 結(jié)語

本文提出了一種ICEEMD和排列熵相結(jié)合的腦血氧降噪方法，利用了ICEEMD根據(jù)信號(hào)自身從而自適應(yīng)的將信號(hào)分解為一系列不同尺度的模態(tài)分量。該方法與EMD相比，解決了模態(tài)混疊問題；與EEMD相比，解決了殘余噪聲問題；與CEEM?DAN相比，則解決了虛假模態(tài)問題。此外，利用排列熵判斷IMF分量的時(shí)間復(fù)雜度，為信號(hào)的重構(gòu)提供了依據(jù)，進(jìn)一步提高了重構(gòu)信號(hào)的精確度。通過實(shí)驗(yàn)證明，本文提出的方法可以有效抑制腦血氧信號(hào)中存在的基線漂移和高頻噪聲問題。