基于隨機(jī)森林模型的動(dòng)態(tài)動(dòng)脈血管硬化指數(shù)估計(jì)方法*

張?？?，程云章，張?zhí)煲?/p>

(上海理工大學(xué) 上海介入醫(yī)療器械工程技術(shù)研究中心，上海 200093)

引言

根據(jù)世界衛(wèi)生組織發(fā)布的《2020世界衛(wèi)生統(tǒng)計(jì)報(bào)告》，僅在2016年，全球因非傳染性疾病(noncommunicable diseases，NCDs)導(dǎo)致的死亡就達(dá)到了4 100萬(wàn)，約占當(dāng)年全球總死亡人數(shù)的71%。心血管疾病(cardiovascular diseases，CVDs)，作為NCDs中主要疾病之一，每年造成的死亡數(shù)量約為1 700萬(wàn)[1]。動(dòng)脈粥樣硬化心血管疾病(atherosclerotic cardiovascular disease)和動(dòng)脈纖維性硬化性心血管疾病(arterioscleroticcardiovasculardisease)是心血管疾病中造成諸多臟器障礙甚至病變的重要原因之一。在動(dòng)脈粥樣硬化早期，其臨床癥狀不明顯。若能在早期及時(shí)發(fā)現(xiàn)并對(duì)患者進(jìn)行有效地治療，可顯著減緩病情惡化的程度[1-2]。

傳統(tǒng)的大動(dòng)脈彈性功能檢測(cè)的方法包括：動(dòng)脈脈搏波傳導(dǎo)速度 (pulse wave velocity, PWV)和超聲檢測(cè)內(nèi)膜-中層厚度 (intima-media thickness, IMT)。動(dòng)脈的脈搏波傳導(dǎo)速度能間接反映大動(dòng)脈的擴(kuò)張性。高分辨率超聲既可以通過(guò)觀察動(dòng)脈某一橫斷面壓力和動(dòng)脈內(nèi)徑的變化關(guān)系反映動(dòng)脈彈性的改變，又可以檢測(cè)頸動(dòng)脈內(nèi)膜-中層的厚度，是早期診斷動(dòng)脈粥樣硬化的重要依據(jù)。但是，以上兩種檢測(cè)手段都需要昂貴的檢測(cè)設(shè)備及專業(yè)的操作人員，在大規(guī)模的臨床實(shí)踐和流行病學(xué)研究中難以廣泛推廣。因此，研究更為簡(jiǎn)單有效的檢測(cè)手段具有現(xiàn)實(shí)意義[4]。

2006年，李燕等[3-4]首次提出了AASI這個(gè)能反映動(dòng)脈血管彈性的指數(shù)，并與傳統(tǒng)的反映動(dòng)脈血管彈性的指標(biāo)進(jìn)行了相關(guān)性分析，包括了脈搏波傳導(dǎo)速度、中心增強(qiáng)指數(shù)、外圍增強(qiáng)指數(shù)。之后，開(kāi)展了一些關(guān)于AASI的研究。Palmas等[5]比較了AASI和24 h脈壓尿液白蛋白排泄的相關(guān)性。Daolan等[6]比較AASI以及脈壓對(duì)心臟性死亡以及中風(fēng)的預(yù)測(cè)表現(xiàn)。結(jié)果顯示，AASI作為對(duì)中風(fēng)的預(yù)測(cè)因子優(yōu)于脈壓差，在心臟性死亡的預(yù)測(cè)上，脈壓要優(yōu)于AASI，而在正常血壓區(qū)間內(nèi)，兩者都無(wú)法作為獨(dú)立的心臟性死亡或者中風(fēng)的預(yù)測(cè)因子。Eniksen等[7]發(fā)現(xiàn)AASI對(duì)于中年人群的腎小球?yàn)V過(guò)率 (glomerular filtration rate, GFR)下降是一個(gè)重要的獨(dú)立因素。Obata等[8]研究了腦血流自動(dòng)調(diào)節(jié)下限(lower limit of cerebral autoregulation，LLA)與AASI之間的關(guān)系，發(fā)現(xiàn)LLA與AASI關(guān)系顯著，AASI可以作為術(shù)中預(yù)測(cè)LLA閾值的重要工具。

盡管學(xué)者們對(duì)AASI做了大量的研究，但大部分的研究都停留在臨床數(shù)據(jù)統(tǒng)計(jì)與對(duì)AASI表征生理信息的驗(yàn)證上。AASI作為一個(gè)獲取難度相對(duì)較低的能反映心血管健康狀況的生理指標(biāo)，對(duì)于AASI的實(shí)際應(yīng)用成果卻不多。因此，本研究基于動(dòng)態(tài)動(dòng)脈硬化指數(shù) (ambulatory arterial stiffness index, AASI)，提出了一種新的測(cè)量患者動(dòng)脈硬化程度的方法。作為一種動(dòng)脈硬化的早篩方法，該方法能有效地降低測(cè)量成本，并對(duì)后續(xù)的臨床診斷和檢測(cè)具有一定的參考意義。與傳統(tǒng)的血管彈性探測(cè)方法相比，基于機(jī)器學(xué)習(xí)得到的AASI估計(jì)值無(wú)需超聲診斷裝置或者侵入式探測(cè)，降低了重要生理信息獲取的難度，對(duì)心血管疾病的早篩有重要意義。

1 材料與方法

基于機(jī)器學(xué)習(xí)的動(dòng)態(tài)動(dòng)脈血管彈性指數(shù)的估計(jì)流程圖，見(jiàn)圖1。

圖1 AASI估計(jì)系統(tǒng)流程圖

1.1 動(dòng)態(tài)動(dòng)脈硬化指數(shù)(AASI)

早在20世紀(jì)初，一些研究已經(jīng)注意到收縮壓和舒張壓的變化可以反映動(dòng)脈血管的彈性功能。當(dāng)一個(gè)人的動(dòng)脈血管健康有彈性時(shí)，收縮壓升高，舒張壓也相應(yīng)地升高；當(dāng)血管彈性降低后，收縮壓增高，而舒張壓增高不明顯、甚至降低。這說(shuō)明收縮壓和舒張壓兩者之間的動(dòng)態(tài)變化關(guān)系在一定程度上反映動(dòng)脈的彈性功能。收縮壓和舒張壓的對(duì)應(yīng)數(shù)值，可以分析舒張壓和收縮壓這之間的回歸關(guān)系。李燕等把一段時(shí)間內(nèi)采集得到的收縮壓與舒張壓做回歸分析，其中收縮壓作為自變量，舒張壓作為因變量。AASI的計(jì)算方式為1減去回歸斜率b。動(dòng)脈彈性好，則AASI趨近于0；動(dòng)脈彈性僵硬，則AASI趨近于1[3-4]。AASI的數(shù)學(xué)表達(dá)式見(jiàn)式⑴、式⑵：

DBP=a+b×SBP

(1)

AASI=1-b

(2)

其中，DBP代表一個(gè)心動(dòng)周期的舒張壓，SBP代表一個(gè)心動(dòng)周期的收縮壓。a、b是回歸方程的回歸系數(shù)。

1.2 數(shù)據(jù)源

本研究所用公共數(shù)據(jù)庫(kù)來(lái)自于UCI機(jī)器學(xué)習(xí)資料檔案庫(kù)(UCI Machine Learning Repository)中的無(wú)袖帶式血壓估計(jì)數(shù)據(jù)集。該數(shù)據(jù)集實(shí)質(zhì)上是經(jīng)過(guò)去噪處理后的MIMICII(multi-parameter intelligent monitoring intensive careII)數(shù)據(jù)集[9-10]。數(shù)據(jù)集包含125 Hz的從指尖獲取的光電容積脈搏波(PPG)信號(hào)、125 Hz的有創(chuàng)動(dòng)脈血壓(ABP)信號(hào)、125 Hz的心電信號(hào)(ECG)。其中，由于MIMICII數(shù)據(jù)庫(kù)自身帶有相當(dāng)可觀的危重患者的數(shù)據(jù)，血壓處于極高值或極低值的病例被移除(SBP≥180 mmHg,DBP≥130 mmHg)∩(SBP≤80 mmHg,DBP≤60 mmHg)。此研究所使用的數(shù)據(jù)集基本信息見(jiàn)表1。在本研究中，80%的樣本(共計(jì)2 036例患者數(shù)據(jù))被用于模型的訓(xùn)練；20%的樣本(共計(jì)510例患者)被用于模型的測(cè)試。

表1 所用數(shù)據(jù)集的基礎(chǔ)信息

1.3 特征提取

特征提取是機(jī)器學(xué)習(xí)中最重要的部分之一。根據(jù)Sharifi等所做的研究[11-15]，本研究一共選用了13個(gè)特征。由于源信號(hào)PPG和ECG是時(shí)域信號(hào)而AASI則相對(duì)是一個(gè)常值，因此,本研究中所有的時(shí)域信號(hào)皆以在整個(gè)時(shí)軸上的均值呈現(xiàn)，以此使輸入信號(hào)與被估計(jì)值得以匹配。整個(gè)特征提取過(guò)程在Matlab R2020a平臺(tái)上實(shí)現(xiàn)。表2是所用全部特征名,圖2—圖4是對(duì)信號(hào)的提取特征過(guò)程。

注：A和B作為周期內(nèi)的最大值和最小值，被定義為是一個(gè)心動(dòng)周期內(nèi)的收縮壓以及舒張壓。

表2 特征名及其定義

注：C和E是一個(gè)心動(dòng)周期內(nèi)的PPG最大值和PPG最小值； D是PPG信號(hào)這個(gè)周期的均值；F和G是FDPPG和SDPPG的均值。圖3 PPG信號(hào)及其衍生的FDPPG和SDPPGFig.3 PPG, FDPPG and SDPPG

1.4 隨機(jī)森林回歸

隨機(jī)森林(random forest ，RF)回歸是由決策樹(shù)學(xué)習(xí)衍生而來(lái)的。即使特征與目標(biāo)呈非線性關(guān)系，RF也能有效地估計(jì)輸出。RF算法采用bootstrap分類器來(lái)估計(jì)響應(yīng)變量。它將每一組決策樹(shù)模型都與數(shù)據(jù)集進(jìn)行匹配。每棵樹(shù)都使用訓(xùn)練數(shù)據(jù)不同的bootstrap樣本訓(xùn)練，以此作為自己的回歸函數(shù)，由回歸給出的最終輸出與各個(gè)決策樹(shù)輸出的平均值對(duì)應(yīng)，即為RF糾正決策樹(shù)算法常見(jiàn)的過(guò)擬合現(xiàn)象的過(guò)程。而未進(jìn)入bootstrap分類器的樣本被稱為out-of-bag (OOB)樣本,它們被用來(lái)測(cè)試決策樹(shù)的準(zhǔn)確性和評(píng)估整體模型的性能[16-20]。

注：H和I是ECG信號(hào)的最大值以及最小值。M和N是

整個(gè)回歸過(guò)程在Matlab平臺(tái)上使用Treebagger函數(shù)完成。其形式見(jiàn)式(3)：

B=Treebagger(NumTrees,X,y)

(3)

其中，AASI作為響應(yīng)y，不同形式的特征組合作為輸入X，樹(shù)的數(shù)量(NumTrees)設(shè)置為100。

1.5 遺傳算法

表2中的13個(gè)特征的組合將被作為RF模型的輸入。但若使用遍歷思路解決該特征篩選問(wèn)題，將會(huì)有213-1輪RF模型計(jì)算。考慮到整個(gè)數(shù)據(jù)集大小，該思路是不可取的。因此，引入遺傳算法(geneticalgorithm,GA)進(jìn)行特征篩選，以此找到最優(yōu)特征子集組合[21-25]。

在本研究中，RF模型的平均絕對(duì)誤差(MAE)作為GA的適應(yīng)度函數(shù)。同時(shí)，所有的13個(gè)特征均以“0”、“1”的形式表明其被選擇與否(“0”代表該特征不被選擇?！?”代表該特征被選擇)，它們共同組成13位二進(jìn)制數(shù)字作為其特征子集組合編碼，以此量化特征篩選問(wèn)題,見(jiàn)圖5。

圖5 特征子集編碼過(guò)程

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

本研究所用衡量方法的指標(biāo)是均方根誤差(root mean squared error, RMSE)、絕對(duì)平均誤差(mean absolute error, MAE)、皮爾森相關(guān)系數(shù)(pearson correlation coefficient,r)和組內(nèi)相關(guān)系數(shù)(Intraclass correlation coefficient, ICC)。實(shí)驗(yàn)結(jié)果，見(jiàn)表3。

表3種挑選了四種經(jīng)典回歸模型進(jìn)行AASI的估計(jì)，全特征輸入代表了未經(jīng)篩選的特征子集組合。其中，脊回歸損失函數(shù)的正則項(xiàng)系數(shù)設(shè)定為α=5。由表3可知，隨機(jī)森林模型在估計(jì)AASI問(wèn)題時(shí)誤差最小(MAE、RMSE)、估計(jì)值與真實(shí)值的線性相關(guān)程度最高(r)。為了進(jìn)一步提高本研究方法的性能，引入GA對(duì)綜合性能表現(xiàn)最好的隨機(jī)森林模型進(jìn)行特征子集篩選。

2.1 遺傳算法特征篩選

引入GA進(jìn)行特征篩選，GA的迭代過(guò)程見(jiàn)圖6。

表3 四種回歸模型估計(jì)結(jié)果

圖6 尋找最優(yōu)特征子集的GA過(guò)程

由GA獲得的最優(yōu)特征子集見(jiàn)表4，其估計(jì)結(jié)果表現(xiàn)見(jiàn)表5。

表4 由GA獲得的最優(yōu)特征子集

2.2 在血壓以及心率維度下AASI估計(jì)分布情況

為了分析進(jìn)行過(guò)特征篩選后的RF模型對(duì)于不同患者的準(zhǔn)確度，本研究從血壓、心率兩個(gè)維度對(duì)所有樣本進(jìn)行了誤差統(tǒng)計(jì)，見(jiàn)圖7。

圖7的橫坐標(biāo)軸代表每個(gè)樣本的平均動(dòng)脈血壓，縱坐標(biāo)軸代表每個(gè)樣本的心率，色塊代表絕對(duì)誤差大小：由淺入深代表絕對(duì)誤差由大至小。一個(gè)黑點(diǎn)代表一個(gè)病例，為了使統(tǒng)計(jì)結(jié)果方便可視，且能展現(xiàn)出總體樣本的趨勢(shì)，在此引入了局部加權(quán)回歸使圖像完整且平滑。

表5 RF模型下經(jīng)過(guò)GA特征篩選與未經(jīng)過(guò)GA篩選的結(jié)果對(duì)比

3 討論與總結(jié)

由表3可知，隨機(jī)森林模型在估計(jì)AASI時(shí)，其綜合性能顯著優(yōu)于表中的其他三種模型。由表5可知，當(dāng)引入GA完成特征篩選后，RF模型估計(jì)AASI的能力提高，誤差明顯減小。

在本研究中，RF模型的最優(yōu)特征子集為Max of PPG、Min of PPG、FDPIR、HR和RR。

2.2節(jié)分析了在最優(yōu)特征子集和隨機(jī)森林模型下，AASI的絕對(duì)誤差在血壓和心率維度上的分布情況。由圖7可知，AASI估計(jì)值的絕對(duì)誤差在樣本平均動(dòng)脈血壓突破了100 mmHg以后，與其呈正相關(guān)，在圖7表現(xiàn)為色塊在橫軸上由深入淺地變化，而另一維度心率與估計(jì)AASI誤差的關(guān)系在圖7上尚不顯著。產(chǎn)生此現(xiàn)象的原因可能為：(1)訓(xùn)練回歸模型的樣本中，高血壓的樣本數(shù)要低于正常血壓的樣本數(shù)，因此，模型對(duì)高血壓樣本產(chǎn)生的誤差大于正常血壓區(qū)間的樣本；(2)AASI對(duì)高血壓患者的表現(xiàn)要弱于對(duì)正常血壓患者。根據(jù)李燕等[4]的研究，AASI與頸動(dòng)脈-股動(dòng)脈PWV的相關(guān)系數(shù)在面對(duì)正常血壓患者時(shí)為0.51，而面對(duì)高血壓患者時(shí)僅為0.28[26]。

圖7 在最優(yōu)特征子集下AASI誤差分布圖

為對(duì)抗心血管疾病早發(fā)的態(tài)勢(shì)，本研究基于動(dòng)態(tài)動(dòng)脈硬化指數(shù)與隨機(jī)森林模型，完成了動(dòng)態(tài)動(dòng)脈硬化指數(shù)的初步估計(jì)，并且優(yōu)化了在估計(jì)過(guò)程中的特征選擇問(wèn)題，同時(shí)分析了動(dòng)態(tài)動(dòng)脈硬化指數(shù)的估計(jì)絕對(duì)誤差在血壓和心率維度上的分布情況。但是模型的可靠性仍需要多中心的外部資料進(jìn)行驗(yàn)證，并在推廣過(guò)程中不斷優(yōu)化和提高其性能。