黃藝格，楊雨丹，賀知明

（電子科技大學(xué)廣東電子信息工程研究院，廣東東莞523808）

0 引言

胎兒心電信號(hào)監(jiān)測(cè)是獲取妊娠和分娩期間胎兒狀況信息的常規(guī)方法[1]。胎兒心電的特征如心率、波形和動(dòng)態(tài)行為在確定胎兒的生命、胎兒發(fā)育、胎兒成熟度、胎兒窘迫或先天性心臟病的存在等方面都有相當(dāng)重要的意義。

盲源分離（blind source separation，BSS）是一種不知道原始信號(hào)，但是卻能通過(guò)獲取混合信號(hào)恢復(fù)出原始信號(hào)的方法[2]。其中，獨(dú)立成分分析（independent component analysis，ICA）是盲源分離的一個(gè)主要方法，分離出的源信號(hào)數(shù)目在分離之前是確定已知的，由此可以將信號(hào)分解成幾個(gè)獨(dú)立的分量。另外，小波包變換是一種多尺度的時(shí)頻分析工具，具有多分辨率的特點(diǎn)[3]。因此可以將小波包變換與ICA 相結(jié)合用來(lái)提取胎兒心電信號(hào)。

本文首先利用小波包變換實(shí)現(xiàn)在峰值最大點(diǎn)進(jìn)行小波包分解以降低母體腹部混合信號(hào)重疊度，接著利用精度高、速度快的快速獨(dú)立成分分析（fast independent component analysis，fastICA）對(duì)處理之后的母體腹部混合信號(hào)進(jìn)行分解得到純凈的胎兒心電信號(hào)。實(shí)驗(yàn)結(jié)果表明，該方法較好地處理了常規(guī)方法因母體心電信號(hào)和胎兒心電信號(hào)重疊導(dǎo)致的提取算法精度不高的問(wèn)題。

1 胎兒心電信號(hào)波形介紹

生物醫(yī)學(xué)信號(hào)是從人體器官中獲得的電信號(hào)，它代表一個(gè)物理變量，通常被認(rèn)為是時(shí)間的函數(shù)，并且可以從它的振幅、頻率和相位來(lái)描述[4]。胎兒心電信號(hào)是一種生物醫(yī)學(xué)信號(hào)，它提供了胎兒心率的電子表示，胎兒心率包含了胎兒狀況的重要信息，可以提前預(yù)知胎兒的健康狀況以及是否出現(xiàn)呼吸窘迫等情況。但是胎兒心電信號(hào)是一個(gè)相對(duì)微弱的信號(hào)（小于20%的母體心電信號(hào)），并且經(jīng)常包含50 Hz 的工頻干擾、基線漂移、肌電干擾等。有時(shí)候母體心電信號(hào)和胎兒心電信號(hào)可能是緊密重疊的。胎兒心電信號(hào)與成人心電圖相似，包含了相同的基本波形，包括P波、QRS 復(fù)合波、T 波和U 波[5]。

2 基于小波分析和fastICA 的提取算法

2.1 提取算法原理

小波包變換和fastICA 是目前信號(hào)提取最常用也是性能較好的2 種算法，本文中將2 種算法結(jié)合起來(lái)，利用基于峰值最大化的小波包分解和基于負(fù)熵最大化的fastICA 提取胎兒心電信號(hào)。

小波包分解是一種時(shí)頻分析方法，具有多分辨率分析、表示局部信號(hào)特征等特點(diǎn)。離散小波包變換被定義為

其中，d表示小波包變換系數(shù)，j表示分辨率，k表示小波基函數(shù)ψ 的位移，f表示源信號(hào)，2-j和k2-j表示將分辨率進(jìn)行冪數(shù)級(jí)離散化，f（2-j，k2-j）為源信號(hào)的離散化表達(dá)式。小波包分解是一種把第j級(jí)的系數(shù)dj，k分解為第j-1 級(jí)的2 個(gè)小波包變換系數(shù)aj-1，k和dj-1，k的方法，以此類(lèi)推，直到分解到指定的多層小波包級(jí)數(shù)。

非高斯函數(shù)是描述信號(hào)偏離高斯分布的一個(gè)度量，用于測(cè)量在BSS 中混合信號(hào)之間的獨(dú)立程度。非高斯性越強(qiáng)，信號(hào)的獨(dú)立程度越高[6]。根據(jù)小波包分解的性質(zhì)，只有少數(shù)小波包系數(shù)是重要的，因?yàn)橹挥泻苌俚狞c(diǎn)的幅值比0 大，這些幅值比0 大的點(diǎn)的非高斯性比源信號(hào)更強(qiáng)。因此，本文選擇了小波包分解處理混合信號(hào)。

在BSS 問(wèn)題中，峰度是第四階累積量，根據(jù)隨機(jī)變量的峰值k（x）的大小，信號(hào)y可以分為以下3類(lèi)：（1）如果kurt（y）＜0，稱(chēng)y為亞高斯信號(hào)；（2）如果kurt（y）＞0，稱(chēng)y為超高斯信號(hào)；（3）如果kurt（y）=0，稱(chēng)y為高斯信號(hào)。

本文利用小波包分解和峰度來(lái)確定最好的小波包分解節(jié)點(diǎn)，這個(gè)節(jié)點(diǎn)的系數(shù)是低重疊度源信號(hào)的近似表示[3]。每個(gè)節(jié)點(diǎn)的小波包系數(shù)的峰度可以通過(guò)下式計(jì)算得到：

其中，c表示每個(gè)節(jié)點(diǎn)的小波包系數(shù)，E表示求解平均值。

將混合信號(hào)用小波包變換分解之后，去除了小波包系數(shù)小于零的信號(hào)點(diǎn)，得到非高斯性比原混合信號(hào)高得多的信號(hào)，為之后fastICA 提取減少迭代次數(shù)。

FastICA 的工作原理是最大化非高斯性，即源信號(hào)的獨(dú)立性，母體腹部混合信號(hào)來(lái)源于母體心電信號(hào)、胎兒心電信號(hào)和環(huán)境的干擾信息[7]。FastICA 利用采用高階統(tǒng)計(jì)數(shù)據(jù)的固定點(diǎn)算法來(lái)恢復(fù)獨(dú)立成分。當(dāng)觀測(cè)信號(hào)數(shù)量少時(shí)，fastICA 對(duì)于分離信號(hào)非常有效。理論上講，母體心電信號(hào)、胎兒心電信號(hào)和復(fù)雜的噪聲相互之間是獨(dú)立的。母體心電信號(hào)和胎兒心電信號(hào)的信號(hào)形態(tài)相似并且具有幾乎相同的光譜，唯一的區(qū)別是胎兒的心率幾乎是母親的1.5～2 倍[8]。

假設(shè)已經(jīng)觀察到的信號(hào)x（t）=[x1（t），x2（t），…，xn（t）]T是由一些未知的、相互獨(dú)立的信號(hào)s（t）=[s1（t），s2（t），…，sn（t）]T混合而成的，混合公式如下：

其中，A=[a1，a2，…，an]，是混合矩陣；n（t）表示噪聲信號(hào)。FastICA 算法步驟就是求出一個(gè)n×n的分離矩陣W，W=A-1，再求出y（t）=Wx（t），求得的y（t）=[y1（t），y2（t），…，yn（t）]T就是對(duì)s（t）的估計(jì)，其中某一路數(shù)據(jù)即為胎兒信號(hào)[9]?？偟膩?lái)說(shuō)，ICA 可以被建模為求解一個(gè)最優(yōu)化問(wèn)題，它可以對(duì)一個(gè)輸出為獨(dú)立信號(hào)的函數(shù)進(jìn)行建模。

本文的FastICA 基于負(fù)熵最大化，負(fù)熵近似由以下目標(biāo)函數(shù)估計(jì)：

其中，W 是分離矩陣；V是標(biāo)準(zhǔn)高斯隨機(jī)變量，其均值為0，方差為1；函數(shù)G（V）為一個(gè)非二次函數(shù)，表達(dá)式如下：

當(dāng)負(fù)熵達(dá)到最大時(shí)，就表明信號(hào)是可以分離的，所以本文算法求得W 需要滿(mǎn)足負(fù)熵最大和估計(jì)信號(hào)y的均值為0 且方差為1，也就是E{G（WTx）2}=1，這樣混合信號(hào)中的幾個(gè)獨(dú)立分量便可以估計(jì)出來(lái)。

本文中使用的是單位化的fastICA 算法，首先需要注意的是，上文中提到的JG（W）的最大值是在E{G（WTx）}的某些最優(yōu)值中得到的[10]。被E{G（WTx）2}=||W||2=1 約束的E{G（WTx）}的最優(yōu)解是在公式（7）成立時(shí)得到的。

其中，I 為單位矩陣。為了簡(jiǎn)化這個(gè)矩陣的求逆，我們決定近似第一項(xiàng)。由于數(shù)據(jù)已被白化，即E{xxT}=1，所以E{xxTG'（WTx）}≈E{xxT}E{G'（WTx）}=E{G'（WTx）}I，可以用來(lái)得出近似值。這樣導(dǎo)數(shù)矩陣就變成了對(duì)角矩陣，大大簡(jiǎn)化了求逆過(guò)程。同時(shí)用W 的當(dāng)前值取代W0來(lái)求出β 的近似值。因此，我們得到了以下近似的牛頓迭代公式：

其中，W*表示W(wǎng) 的更新值，β=E{WTxG（WTx）}，并且標(biāo)準(zhǔn)化提高了穩(wěn)定性。這個(gè)算法可以通過(guò)在第一個(gè)等式的兩邊同時(shí)乘β-E{G'（WTx）}來(lái)進(jìn)一步簡(jiǎn)化。因此得到了fastICA 的迭代公式：

在實(shí)踐中，fastICA 的期望值將被它們的估計(jì)值所取代。理想情況下，所有可用的數(shù)據(jù)都應(yīng)該被使用，但這會(huì)導(dǎo)致迭代時(shí)間太長(zhǎng)。所以一般用一部分樣本的均值進(jìn)行估計(jì)，其數(shù)目多少可能對(duì)最后估計(jì)的準(zhǔn)確性有相當(dāng)大的影響。因此，在每次迭代時(shí)應(yīng)該分別選擇取樣點(diǎn)，如果收斂不盡如人意，可以選擇增加樣本數(shù)量。

2.2 本文算法與fastICA 算法的比較

（1）本文算法與現(xiàn)有的fastICA 算法相比，前期增加了小波包變換，在ICA 之前就減少了不相關(guān)信號(hào)點(diǎn)的數(shù)量，使得fastICA 提取的時(shí)候樣本數(shù)大大減少，從而具有非?？斓氖諗克俣取?/p>

（2）小波包變換和fastICA 結(jié)合算法與fastICA算法相比，不需要選擇步長(zhǎng)參數(shù)。改進(jìn)算法可以使用任何一個(gè)非線性函數(shù)G得到非高斯分布的獨(dú)立分量，這與fastICA 算法形成了對(duì)比[11]。

2.3 算法步驟

在BSS 方法中，源信號(hào)的混合被假定為固定的和線性的[12]，但是胎兒心電信號(hào)和其他干擾因素總是不固定和線性可分的，因此ICA 并不能把源信號(hào)完全分開(kāi)。在這種情況下，本文算法是基于fastICA的基礎(chǔ)上提出的，信號(hào)處理步驟如下：

（1）輸入采集的原始信號(hào)。原始信號(hào)為包含母體心電信號(hào)和胎兒心電信號(hào)以及各種噪聲的母體腹部混合信號(hào)，將其定義為原始混合心電信號(hào)X。

（2）對(duì)輸入的原始混合心電信號(hào)X進(jìn)行中心化處理，使其均值為0，再對(duì)信號(hào)進(jìn)行白化處理，得到處理后信號(hào)Z。

（3）對(duì)信號(hào)Z進(jìn)行小波包分解處理，得到混合心電數(shù)據(jù)矩陣A，即小波包分解后的混合心電信號(hào)（混合母體和胎兒的心電信號(hào)）。

（4）采用基于負(fù)熵最大化的fastICA 對(duì)混合心電信號(hào)A 進(jìn)行母體與胎兒心電信號(hào)分離處理，獲取含有噪聲的胎兒心電信號(hào)。

（5）對(duì)含有噪聲的胎兒心電信號(hào)進(jìn)行小波去噪處理，去除其中的噪聲和基線漂移，輸出清晰的胎兒心電信號(hào)。

3 采樣與結(jié)果分析

本文使用的臨床實(shí)驗(yàn)信號(hào)來(lái)自DaISy（Database for the Identification of Systems）數(shù)據(jù)庫(kù)。圖1 為觀測(cè)到的胎兒心電信號(hào)與母體心電信號(hào)的混合信號(hào)，通過(guò)觀察可以看到該混合信號(hào)中包含母體心電信號(hào)、胎兒心電信號(hào)以及各種噪聲，為了更方便地觀測(cè)各路信號(hào)，使用MATLAB 將各通道信號(hào)分離出來(lái)進(jìn)行分析。

圖2 為DaISy 數(shù)據(jù)庫(kù)中某位孕婦的8 個(gè)導(dǎo)聯(lián)的母體體表心電信號(hào)，（a）～（e）是母體腹部混合信號(hào)，（f）～（h）是母體胸部心電信號(hào)。5 個(gè)母體腹部混合信號(hào)中，第1 通道可以明顯地觀察到胎兒心電，可見(jiàn)腹部信號(hào)是母親心電和胎兒心電的混合；第4 通道中噪聲多，幾乎不含胎兒心電；第5 通道中的胎兒心電成分很少，不能用來(lái)分離胎兒心電。

圖2 8 通道母體心電信號(hào)（源信號(hào)）

本文采用第1、3、6 通道的信號(hào)進(jìn)行fastICA 處理，首先將這3 個(gè)通道的信號(hào)混合，在小波包處理階段，選擇db4 多貝西小波函數(shù)對(duì)混合信號(hào)進(jìn)行兩層小波包分解，從而減少重疊。經(jīng)過(guò)小波包分解，利用峰度極大化進(jìn)行最優(yōu)節(jié)點(diǎn)選擇，得到2 條清晰的直線，即源信號(hào)的非高斯表示，減少了重疊度。圖3（a）為混合心電信號(hào)的聯(lián)合分布圖，可以看出圖上母體心電信號(hào)、胎兒心電信號(hào)重疊在一起；圖3（b）為經(jīng)過(guò)小波包分解處理之后的混合心電信號(hào)的聯(lián)合分布圖，可以明顯地看出信號(hào)的重疊程度減少了。

通過(guò)fastICA 解混后，可以觀察到由第1、3、6 通道的觀測(cè)信號(hào)分離出來(lái)了4 路信號(hào)，如圖4 所示，由于胎兒心電較母體心電頻率快、幅值低，可以判斷出（a）為母體信號(hào)，（b）為包含了少量胎兒心電的母體心電，且（b）中含有的母體信號(hào)也僅為（a）信號(hào)的1/10 左右，（c）為分離出的較清晰的胎兒心電信號(hào)，（d）為噪聲信號(hào)。

圖3 母體腹部混合信號(hào)小波包分解前后的小波包系數(shù)的聯(lián)合分布圖

圖4 解混后的信號(hào)

為方便觀察，本文特意從圖4 中截取了胎兒的心電圖，如圖5 所示，其中仍然包含了少量母體心電信號(hào)、噪聲和基線漂移，對(duì)于后期觀察會(huì)造成一定影響，并且大大降低了信噪比。

由于解混之后的胎兒心電信號(hào)包含太多噪聲，我們利用小波包分解重構(gòu)得到去噪后的心電信號(hào)，如圖6 所示。

圖5 本文算法獲得的胎兒心電

圖6 小波去噪之后的胎兒心電

4 結(jié)語(yǔ)

本文提出了一種利用峰度最大化小波包分解結(jié)合fastICA 算法進(jìn)行胎兒心電分離。在小波包分解之后混合心電信號(hào)的重疊度大大減弱，之后對(duì)信號(hào)進(jìn)行fastICA 處理，最終獲取對(duì)源信號(hào)的估計(jì)。臨床實(shí)驗(yàn)結(jié)果表明該算法擁有更好的性能。在進(jìn)行fastICA分離信號(hào)之后，通過(guò)小波去噪、對(duì)胎兒心電圖信號(hào)基線漂移和噪聲的噪聲抑制來(lái)重建信號(hào)，得到清晰的胎兒心電信號(hào)。實(shí)驗(yàn)結(jié)果表明，該算法可以將胎兒心電信號(hào)與母體心電信號(hào)分離開(kāi)來(lái)。在基于負(fù)熵的fastICA 中即使只用到負(fù)熵一般形式的近似，其收斂速度也是二次化的，算法更簡(jiǎn)單可靠，收斂性更好。但是該算法需要使用3 導(dǎo)聯(lián)母親體表信號(hào)才能進(jìn)行信號(hào)分離，增加了獲取數(shù)據(jù)的復(fù)雜性，如果應(yīng)用到臨床醫(yī)學(xué)中實(shí)時(shí)性還有待改進(jìn)。