基于變分模態(tài)分解的脈率變異性信號(hào)提取方法

楊文龍，李 鑫，丑麗娟,2，楊海萍，丑永新*

(1.常熟理工學(xué)院 電氣與自動(dòng)化工程學(xué)院， 江蘇 蘇州 215500；2.東北石油大學(xué) 計(jì)算機(jī)與信息技術(shù)學(xué)院，黑龍江 大慶 163711)

0 引言

人體心臟的收縮與舒張驅(qū)動(dòng)血液在心血管系統(tǒng)周而復(fù)始地循環(huán)，引起人體皮膚淺表處(如手指、腕部、耳垂等[1])脈管的搏動(dòng)，這種搏動(dòng)稱為脈搏。脈搏信號(hào)中含有人體心血管系統(tǒng)豐富的生理和病理信息，診脈是中醫(yī)進(jìn)行疾病診斷的重要依據(jù)之一。從脈搏信號(hào)中提取的脈率[2]，是評(píng)價(jià)人體健康狀態(tài)的重要指標(biāo)之一。人體心臟搏動(dòng)規(guī)律隨著人體狀態(tài)的變化而不斷調(diào)整，因此，脈率也在不斷發(fā)生變化。脈率之間的這種變化稱為脈率變異性(Pulse Rate Variability,PRV)，PRV是評(píng)價(jià)人體心血管系統(tǒng)和神經(jīng)系統(tǒng)的重要依據(jù)。研究表明，從PRV中提取參數(shù)可用于心源性猝死、冠心病、心力衰竭、高血壓及呼吸暫停綜合癥等疾病的診斷[3-6]，具有十分重要的臨床應(yīng)用價(jià)值。

目前，已有多種PRV信號(hào)提取方法，主要可分為3種：① 基于脈搏波時(shí)域征的提取方法，這類方法通過(guò)設(shè)置幅度、時(shí)間或者差分閾值來(lái)提取脈搏波波峰或者波谷等特征[7]，進(jìn)而計(jì)算脈搏信號(hào)的主波間期。這類方法的速度快，準(zhǔn)確性高，但是抗噪性能差。② 基于脈搏信號(hào)頻譜的PRV估計(jì)法[8]，從脈搏信號(hào)的頻譜成分中估計(jì)PRV，抗噪性很強(qiáng)，但是誤差太大。③ 從脈搏信號(hào)的分解分量中提取。這類方法先通過(guò)傅里葉級(jí)數(shù)[9]、小波分解[10]、EMD分解[11]、形態(tài)學(xué)濾波[12-13]等方法將脈搏信號(hào)分解成不同頻段的分量，再?gòu)牡皖l分量中提取PRV信號(hào)。這種方法的準(zhǔn)確性和抗噪性很高，但是分量的準(zhǔn)確提取選擇決定著PRV信號(hào)的提取精度，傅里葉級(jí)數(shù)的窗寬、小波分解的頻段、形態(tài)學(xué)濾波法結(jié)構(gòu)元素的長(zhǎng)度等因素選擇主觀，使其在實(shí)際應(yīng)用中受到限制。EMD算法可實(shí)現(xiàn)脈搏波的自適應(yīng)分解，然而其存在模態(tài)混疊。

近年來(lái)，有學(xué)者結(jié)合維納濾波、Hilbert變換和頻率混合變分問(wèn)題的求解提出了變分模態(tài)分解(Variational Mode Decomposition，VMD)算法，在震動(dòng)[14]、氣象[15]、電力[16]、生物醫(yī)學(xué)等信號(hào)處理方面表現(xiàn)出優(yōu)越性能。因此，本文以VMD算法為核心，提出一種從脈搏信號(hào)中提取PRV信號(hào)的新算法。首先，通過(guò)VMD算法對(duì)脈搏信號(hào)分解，從低頻分量中獲取脈搏信號(hào)主波定位的參考信號(hào)；然后，以此為依據(jù)，對(duì)脈搏信號(hào)的波谷進(jìn)行定位，計(jì)算PRV信號(hào)。設(shè)計(jì)實(shí)驗(yàn)，通過(guò)實(shí)測(cè)脈搏信號(hào)對(duì)所提出算法的準(zhǔn)確性進(jìn)行評(píng)估。

1 算法原理

1.1 VMD分解算法

VMD分解算法的目的是將原始信號(hào)分解為多個(gè)本征模態(tài)函數(shù)(Intrinsic Mode Function，IMF)的和，每個(gè)IMF分量為特定帶寬的調(diào)幅-調(diào)頻信號(hào)。分解后第k個(gè)IMF分量的定義如下：

imfk(t)=ak(t)cos(φk(t))，

(1)

式中，ak(t)為第k個(gè)IMF分量的包絡(luò)，φk(t)為其相位，則相位對(duì)應(yīng)的角頻率為ωk(t)=φ′k(t)。

對(duì)于給定的脈搏信號(hào)ppg(t)，可以分解成K個(gè)IMF分量之和，每個(gè)IMF分量對(duì)應(yīng)的約束變分模型為：

(2)

式中，{imfk}={imf1,imf2,…,imfk}為VMD分解得到所有IMF的集合，{ωk}={ω1,ω2,…,ωk}為各IMF的中心角頻率集合。

對(duì)于式(2)所示問(wèn)題的最優(yōu)解問(wèn)題，引入拉格朗日函數(shù)后，可轉(zhuǎn)化為非約束變分問(wèn)題，即

(3)

(4)

式中，ppg為時(shí)域表達(dá)式，PPG為其在頻域內(nèi)的表達(dá)式。式(4)經(jīng)過(guò)Parseval傅里葉等距變換后，可以得到：

(5)

式中，Λn+1(ω)的更新方式為：

(6)

同理可得：

(7)

在此基礎(chǔ)上，VMD分解分為以下幾步：

② 根據(jù)式(5)～式(7)更新imfk、ωk和λ；

③ 重復(fù)前兩步，直到滿足以下約束條件：

(8)

1.2 基于VMD分解的PRV信號(hào)提取原理

脈搏信號(hào)經(jīng)過(guò)VMD分解后，可得到一系列IMF分量，不同的IMF分量中包含脈搏信號(hào)不同的頻段。其中，imf1中含有脈搏信號(hào)的基波分量，基波分量的主要成分為脈搏信號(hào)的主波，即圖1中的虛線信號(hào)。該信號(hào)不含脈搏信號(hào)的高頻噪聲、重搏波和潮波等波谷定位的干擾信息，相比直接從脈搏信號(hào)中定位波谷的難度明顯降低。同時(shí)，imf1能自適應(yīng)地隨著脈搏信號(hào)趨勢(shì)改變而變化，其兩個(gè)峰值(圖1中黑色圓圈)之間的最小值正好對(duì)應(yīng)脈搏信號(hào)的波谷。因此，imf1分量可作為PRV信號(hào)檢測(cè)的重要依據(jù)。

圖1 從imf1分量中定位波谷的示意圖Fig.1 Process of locating troughs fromimf1 component

依據(jù)圖1所示的原理，提出圖2所示的PRV信號(hào)檢測(cè)流程。

圖2 算法流程Fig.2 Process of proposed method

首先，對(duì)脈搏信號(hào)進(jìn)行VMD分解得到imf1分量；對(duì)imf1分量進(jìn)行20點(diǎn)平滑濾波，去除殘余高頻噪聲的干擾；然后，計(jì)算imf1分量的極大值，存儲(chǔ)極大值位置,以相鄰極大值位置為參考，在原始脈搏信號(hào)中計(jì)算最小值位置即為脈搏信號(hào)波谷位置；對(duì)波谷位置進(jìn)行一階差分，計(jì)算PRV信號(hào)，PRV信號(hào)計(jì)算過(guò)程見(jiàn)式(9);最后，對(duì)PRV信號(hào)的提取結(jié)果進(jìn)行評(píng)估。

(9)

式中，fs為脈搏信號(hào)采樣頻率，PPI(j)為由波谷位置一階差分得到的脈搏信號(hào)間期，PRV(j)為第j個(gè)PRV信號(hào)樣本點(diǎn)。

2 實(shí)驗(yàn)數(shù)據(jù)

采用自主研制的可穿戴式生理信號(hào)采集系統(tǒng)獲取實(shí)驗(yàn)數(shù)據(jù)，采集系統(tǒng)如圖3所示。系統(tǒng)由上位機(jī)和下位機(jī)組成，上下位機(jī)通過(guò)Zigbee組網(wǎng)進(jìn)行無(wú)線傳輸。下位機(jī)分為兩個(gè)模塊，分別為放置在胸部的心電采集模塊和放置在腕部的脈搏信號(hào)采集模塊。上位機(jī)實(shí)時(shí)接收下位機(jī)數(shù)據(jù)，對(duì)其進(jìn)行存儲(chǔ)與處理。

圖3 脈搏信號(hào)采集過(guò)程Fig.3 Recording process of pulse signal

招募年齡為21～26歲的健康在校大學(xué)生(男生18名，女生13名)作為實(shí)驗(yàn)對(duì)象。數(shù)據(jù)實(shí)驗(yàn)前在常熟市第一人民醫(yī)院進(jìn)行體檢，確保實(shí)驗(yàn)對(duì)象的健康狀態(tài)。在常熟第一人民醫(yī)院倫理委員會(huì)的監(jiān)督之下，實(shí)驗(yàn)對(duì)象佩戴實(shí)驗(yàn)設(shè)備，處于靜坐狀態(tài)，采集各路信號(hào)4 min，采樣頻率為1 000 Hz。

3 實(shí)驗(yàn)結(jié)果及分析

對(duì)實(shí)驗(yàn)數(shù)據(jù)按照?qǐng)D2流程處理，其中，VMD分解的參數(shù)配置為：分解層數(shù)為5，τ=0.001，迭代停止誤差ε=0.1。對(duì)于如圖4所示的脈搏信號(hào)，進(jìn)行VMD分解后得到IMF分量，如圖5所示，分解后得到5個(gè)IMF分量。相比于圖4，imf1分量含有脈搏信號(hào)的主波成分，imf2、imf4和imf5中含有噪聲成分，而imf3中含有重搏波和潮波成分。

圖4 原始脈搏信號(hào)Fig.4 A raw pulse signal

(a)imf1

(b)imf2

(c)imf3

(d)imf4

(e)imf5圖5 脈搏信號(hào)的VMD分解結(jié)果Fig.5 VMD decomposition results of a raw pulse signal

因此，選用imf1分量按照?qǐng)D2所示流程計(jì)算PRV信號(hào)，結(jié)果如圖6所示。可以看出，只有第一個(gè)波谷沒(méi)有定位出來(lái)，其他波谷都可以準(zhǔn)確定位，這是由于要從imf1分量的兩個(gè)極大值之間定位波谷導(dǎo)致的。圖7為從圖2脈搏信號(hào)對(duì)應(yīng)完整信號(hào)中提取的PRV信號(hào)，可以看出所提方法可以準(zhǔn)確提取PRV信號(hào)。

圖6 脈搏波波谷定位結(jié)果Fig.6 Pulse wave and its trough locations

圖7 PRV信號(hào)提取結(jié)果Fig.7 A PRV signal extracted from pulse signal

將所提算法用于31組實(shí)驗(yàn)數(shù)據(jù)的PRV信號(hào)提取，檢出數(shù)、漏檢數(shù)、誤檢數(shù)和準(zhǔn)確率如表1所示。

表1 PRV信號(hào)提取準(zhǔn)確率Tab.1 Accuracy of PRV signal extraction

可以看出，31組數(shù)據(jù)中，25組PRV信號(hào)的提取準(zhǔn)確率為100%，其余各組的準(zhǔn)確率為97%以上，其中有4組的準(zhǔn)確率在99.5%以上。只有1組的準(zhǔn)確率為97.35%，這是由于原始脈搏信號(hào)中運(yùn)動(dòng)偽跡造成的干擾段導(dǎo)致。總體而言，所提算法的準(zhǔn)確率很高，這是因?yàn)橥ㄟ^(guò)VMD分解后，原始脈搏信號(hào)的高頻噪聲和重搏波、潮波等干擾特征分解到其他IMF分量中，對(duì)主波波谷定位的干擾降低，從而得到較高的準(zhǔn)確率。

4 結(jié)束語(yǔ)

PRV信號(hào)準(zhǔn)確提取十分重要，本研究提出一種基于VMD分解的脈率變異性提取算法，可準(zhǔn)確地從原始脈搏信號(hào)中提取PRV信號(hào)。具有抗干擾性強(qiáng)不需要設(shè)置閾值等優(yōu)點(diǎn)，可用于日常工程生活環(huán)境下PRV信號(hào)的準(zhǔn)確提取。然而，VMD分解速度較慢，耗用系統(tǒng)資源較多，嚴(yán)重限制了所提算法的應(yīng)用場(chǎng)景。因此，未來(lái)工作將對(duì)VMD分解算法進(jìn)行改進(jìn)，提高其運(yùn)算速度，進(jìn)而拓闊所提算法的應(yīng)用領(lǐng)域。