蔡 慧，羅佳偉

(1.四川大學數(shù)學學院，四川成都 610065；2.四川大學臨床醫(yī)學院，四川成都 610041)

0 引言

隨著生物醫(yī)學和醫(yī)療保健的大數(shù)據(jù)增長以及電子病歷記錄( EHR )數(shù)據(jù)的大量數(shù)字化，先進醫(yī)療設備及醫(yī)療智能化等原因，使得電子病歷數(shù)據(jù)呈指數(shù)級增長趨勢，對醫(yī)學大數(shù)據(jù)的準確分析將有助于早期疾病檢測.傳統(tǒng)的分析方法如回歸分析等確實十分有效，但隨著數(shù)據(jù)變得十分復雜，分析準確性會降低.為了捕捉潛在復雜信息，提高準確率，越來越多人選擇使用神經(jīng)網(wǎng)絡預測模型， RNN的變體如LSTM和GRU預測效果絕佳，但模型需要特別苛刻的條件，如數(shù)據(jù)必須完整無缺失、時間步長固定、且模型本身不捕獲時間信息等.在大數(shù)據(jù)時代，我們越來越追求模型的效率和實用性，不僅希望模型的準確率盡可能高，模型收斂越快越好，同時也希望模型能夠直接處理帶有缺失的數(shù)據(jù).

本文首先是希望能夠找到一份最真實的電子病歷數(shù)據(jù)，它滿足一般電子病歷的幾乎所有特征，數(shù)據(jù)規(guī)模大變量足夠多，有自然的缺失，觀察的次數(shù)盡可能多.其目的是希望選出的模型能夠處理足夠復雜的數(shù)據(jù)，滿足模型的通用特點.其次是希望在眾多的RNN推廣的模型中，或按照已有的推廣模型演變的變體中找到一個絕佳的模型，使得它能夠方便快捷準確的處理復雜且?guī)в腥笔У拇髷?shù)據(jù).

1 相關工作

循環(huán)神經(jīng)網(wǎng)絡(簡稱RNNs)較早出現(xiàn)在1982年[1]，后經(jīng)過一系列改進工作[2]，應用極為廣泛. RNN具有一些不錯的屬性，例如強大的預測性能以及捕獲長期時間依存關系和變長觀測值的能力[3].通過使用循環(huán)連接，信息可以在這些網(wǎng)絡中循環(huán)任意長時間并實現(xiàn)參數(shù)共享等.但與此同時也存在一些限制[4]：

(1)RNNs存在梯度消失或梯度爆炸的問題；

(2)RNNs有時甚至不收斂.

長短時間記憶單元(簡稱LSTM)[5]，是RNN的變體.它通過將門控函數(shù)納入其狀態(tài)動態(tài)中來巧妙地解決梯度消失等問題.在每個時間步中，LSTM維護一個隱藏的向量ht和一個負責控制狀態(tài)更新和輸出的記憶向量ct[6].與普通RNN相比，LSTM的收斂性更好；能夠學習長期依賴關系，可以從多變量，變長的時間序列觀測記錄中捕獲遠距離依賴性.

2015年Lipton等人[5]使用長短期記憶(LSTM)構建了一個診斷模型，該模型可以有效捕獲具有不同長度和遠距離依賴性時間序列觀測值，而無法捕獲采樣間隔的不規(guī)則性.階段性長短時間記憶單元模型Phased LSTM(P-LSTM)[7]是LSTM的變體，是在LSTM的基礎上添加了一個時間門kt，該門由具有頻率范圍的參數(shù)化振蕩控制，該頻率范圍僅在周期的一小部分期間產(chǎn)生存儲器單元的更新.Phased LSTM網(wǎng)絡實現(xiàn)了比常規(guī)LSTM更快的收斂.另外，在LSTM的基礎上，門控循環(huán)單元(GRU)[8]將LSTM中的內(nèi)部復雜的門控機制稍微做了一定的改變，簡化了通道，簡單地說是把兩條通道合并成了一條通道，在諸多RNN的變體中比較受歡迎.

實際應用中的多元時間序列數(shù)據(jù)，特別是電子病歷的記錄存在各種缺失值.處理缺失的方法，例如平滑，插值[9]和樣條曲線[10]，無法捕獲變量相關性，也可能無法捕獲復雜的模式來執(zhí)行插補.最近，出現(xiàn)了一個名為GRU-D的模型[3]，該模型通過使用衰減率將標記觀察值缺失與否的向量和觀察時間間隔應用于GRU的更新狀態(tài)中，并希望找到內(nèi)部的規(guī)律.該模型不僅能夠填補缺失，還能利用缺失的模式來改善預測結果.在本文中，這種填補缺失的方法被運用到了每一個LSTM的推廣模型中.

2 神經(jīng)網(wǎng)絡模型

在這個部分中，主要描述了四種LSTM變體的預測模型，它們分別是GRU-D、LSTM-D、P'-LSTM-D和P'-GRU-D，其具體的介紹如下。

2.1 定義[11]

假設一個具有N個樣本的數(shù)據(jù)(本文中即N個病人)，每一個樣本都具有D個變量，每個樣本的每一個變量具有Tn個時間點(n∈1,2,3,….,N)，不同的樣本可以具有不同的時間長度. n代表被采樣的第n個樣本，其中，第n個樣本的數(shù)據(jù)形式可以描述如下：

Xn=(x1,x2,…,xTn)T∈RTn×D、Mn=(m1,m2,…,mTn)T∈RTn×D、Sn=(s1,s2,…,sTn)∈RTn

其中，s=(s1,…,sTn) 表示第 n 個樣本的所有時間點.

為了方便理解，下面列舉了一個具有3個樣本數(shù)(P1,P2,P3)，兩個變量( x1,x2 )和七個時間點(s1,s2,…,s7)的數(shù)據(jù)例子如圖1所示：

其中(a)表示st和xt；(b)表示st和mt；(c)表示st和δt；這三個部分具有相通的數(shù)據(jù)結構，并且P1,P2,P3可以具有不同的時間點，這里直觀的表現(xiàn)是表格可以不一樣長.需要注意的是，真實帶有缺失的數(shù)據(jù)xt在輸入模型之前，需要把缺失的地方用 0 填補占位，但最終要利用mt識別出真正缺失的數(shù)據(jù)并做填補處理.

2.2 GRU-D和LSTM-D神經(jīng)網(wǎng)絡預測模型

在這節(jié)中，本文層層遞進介紹從RNN到LSTM再到GRU，最后GRU推出GRU-D，另附上本文推出的LSTM-D神經(jīng)網(wǎng)絡預測模型.

2.2.1 LSTM神經(jīng)網(wǎng)絡預測模型 首先，從簡單的循環(huán)神經(jīng)網(wǎng)絡(RNN)開始，它及其它的變體的網(wǎng)絡結構和操作始終保持完全共享. 換句話說，所有RNN都具有重復單元鏈的形式，只是在不同的時間步長更新權重.在最初的循環(huán)神經(jīng)網(wǎng)絡中，重復模塊的結構都很簡單，例如使用一個tanh層(如圖2(a)).

ht=tanh(Wxxt+Whht-1+b)

(1)

而標準的長短時間記憶網(wǎng)絡(LSTM)[5]，這是一種特殊的RNN，它使用門控機制和專用記憶單元，通過利用長時間依賴性，可以更好地解決梯度消失或梯度爆炸問題[13].在這項工作中，我們使用具有窺孔連接的版本[14]，并且使用最受歡迎的LSTM架構之一[12](如圖3(b))，這種結構復雜許多.它的主要進步是引入了兩個信息通道，上：存儲單元狀態(tài)通道；下：隱藏單元狀態(tài)通道，以及使用的各種門控機制.其定義的更新公式如下：

圖3 GRU-D、LSTM-D神經(jīng)網(wǎng)絡模塊結構圖Fig.3 Structural Diagram of GRU-D,LSTM-D Neural Network Module

ft=σ(Wxfxt+Whfht-1+bf)

(2)

it=σ(Wxixt+Whiht-1+bi)

(3)

(4)

ot=σ(Wxoxt+Whoht-1+bo)

(5)

(6)

ht=0t⊙tanh(ct)

(7)

其中ft,it,ot分別代表在時間步長t上的遺忘門，輸入門和輸出門.ct是細胞記憶狀態(tài)向量，而xt和ht分別代表輸入特征向量和隱藏層輸出向量.Whi,Whf,Wxi,Wxf,Who和Wxo是權重參數(shù).bi,bf,bc, 和bo是偏差參數(shù)，sigmoid(σ) 和tanh是激活函數(shù).

2.2.2 GRU預測模型 門控循環(huán)單元(GRU)[8]，如圖2(c)所示.和圖2(b)相比，它只有一個輸出，合并了一些通道，對LSTM中的門控單元做了一定的修改整合，它的參數(shù)比長短時間記憶模型要少.加入偏差項后其定義的更新公式如下：

zt=σ(Wxzxt+Whzht-1+bz)

(8)

rt=σ(Wxrxt+Whrht-1+br)

(9)

(10)

(11)

其中zt和rt分別叫做更新門和重置門，均通過上一個隱藏狀態(tài)ht-1和這一時刻的輸入xt以及學習所有模型的參數(shù)(W，U，b)，并結合激活函數(shù)σ來獲取(見公式8和9).

2.2.3 GRU-D神經(jīng)網(wǎng)絡預測模型 2018年車正平(Zhengping Che)等人[3]在GRU的基礎上拓展了模型，使得模型能夠直接處理存在缺失的數(shù)據(jù), 文中將這種模型簡稱為GRU-D(其中重復模塊單元部分如圖3(a)所示)，下面具體來介紹這個模型：

某些領域的數(shù)據(jù)，尤其是在醫(yī)療保健領域，序列中缺失存在兩個重要屬性：(1)如果變量的最后一次觀察是在很久以前發(fā)生的，則該變量的值往往接近某個固定值；(2)如果輸入變量已經(jīng)缺失了一段時間，那么輸入變量的影響將隨著時間的流逝而逐漸消失.為了從根本上解決序列中缺失的問題，他們在神經(jīng)網(wǎng)絡結構中引入衰減率，將缺失表達的兩種模式(觀察值缺失與否的標記向量和時間間隔向量)嵌入到GRU的網(wǎng)絡結構中，并使用反向傳播訓練所有模型組件，希望捕獲到以上兩種特性，達到填補缺失的目的.

需要注意的是，電子病歷數(shù)據(jù)中的每個輸入變量在醫(yī)療應用中都有其自身的含義和重要性，衰減率應因變量而異.GRU-D模型根據(jù)數(shù)據(jù)本身存在的缺失狀況從訓練數(shù)據(jù)中學習衰減率，而不是先行固定.為了使每個衰減率保持在0到1的合理范圍內(nèi)，選擇用指數(shù)負整流器，輸入變量xt衰減率γ定義如下：

γxt=exp{-max(0,Wγxδt+bγx)}

(12)

(13)

有時，輸入衰減可能無法完全捕獲丟失的模式，因為并非所有丟失信息都可以用衰減的輸入值表示. 為了從缺失中獲取更多的信息，他們在GRU-D中還添加了一個隱藏狀態(tài)衰減率γh，這是有對GRU隱藏狀態(tài)衰減的作用，它在計算新的隱藏狀態(tài)ht之前，通過衰減先前的隱藏狀態(tài)ht-1來實現(xiàn).其中隱藏狀態(tài)衰減率γh的計算和輸入衰減率γx類似，如下：

γht=exp{-max(0,Wγhδt+bγh)}

(14)

衰減隱藏狀態(tài)的中間變量公式為：

(15)

另外，將標志向量直接輸入到GRU模型中，得到的GRU-D更新后的公式如下：

(16)

(17)

(18)

(19)

2.2.4 LSTM-D神經(jīng)網(wǎng)絡預測模型 同理，本文把缺失模式下的數(shù)據(jù)應用到LSTM中去，可以得到LSTM-D(其中重復模塊單元部分如圖3(b))模型，對xt的缺失填補方法與GRU-D一樣，也是需要連同神經(jīng)網(wǎng)絡內(nèi)部所有參數(shù)一起訓練學習得到兩個衰減率γxt和γht，且每個變量在每個時刻的衰減率可能不同，具體詳情可以參考2.2.3小節(jié)的描述，這里就不再過多贅述了.LSTM-D最終的更新方程如下：

(20)

(21)

(22)

(23)

(24)

ht=0t⊙tanh(ct)

(25)

對比公式2-7可以發(fā)現(xiàn)，以上公式和原始LSTM預測模型的公式變化可參照2.2.3小節(jié)中GRU-D與GRU的變化，這里不再做過多的介紹.

2.3 P'-LSTM-D和P'-GRU-D預測神經(jīng)網(wǎng)絡模型

在這節(jié)中，我們主要從P-LSTM模型中替換時間門函數(shù)加入缺失模式后推導出P'-LSTM-D模型，在衍生出P'-GRU-D模型.

P-LSTM(階段性LSTM)神經(jīng)網(wǎng)絡模型[7]， 即在LSTM基礎上添加了一個新的門，叫“時間門kt”(如圖4(a)所示)，這個門由具有一定頻率范圍的參數(shù)化振蕩控制，該頻率范圍僅在周期的一小部分期間產(chǎn)生存儲記憶信息的更新，簡而言之，時間門就像一個篩選器，控制數(shù)據(jù)信息進入模型訓練的比例，只有當這個時間門是開放的階段，數(shù)據(jù)才參與訓練并且存儲記憶信息.這個篩選的比例不是人為固定，而是通過分別學習kt中定義的幾個參數(shù)來確定.原模型中的kt定義如下：

圖4 P'-LSTM-D、P'-GRU-D神經(jīng)網(wǎng)絡模塊結構圖Fig.4 Structural Diagram of P'-LSTM-D、P'-GRU-D Neural Network Module

(26)

其中，第一個參數(shù)τ控制振蕩的實時周期；第二個參數(shù)ron控制開放時間的比例；第三個參數(shù)s控制振蕩的相移.φt是一個輔助變量，代表節(jié)奏周期內(nèi)的相位.kt分成三個階段，0到1、1到0和最后的不通過階段，這里的α通常設置為很小的數(shù)，比如0.001等.

本文所用的時間門kt，更換為另一個更為簡單的函數(shù)來代替，原理和思想和原文中的門一樣，同樣有三個參數(shù)：

kt=max {0,Asin(ωt+c)}

(27)

其中，這里A、ω、c需要在訓練過程中不斷學習，需要注意的是，這里需要限制 |A|1 ,其中不同變量的參數(shù)可能不同，不同時間的參數(shù)也可能不同，數(shù)據(jù)只在開放的階段才參與訓練.再結合缺失處理，得到P'-LSTM-D 神經(jīng)網(wǎng)絡模型的更新公式如下，具體還可以參考[11].

(28)

(29)

(30)

(31)

(32)

(33)

(34)

同理，也可以得出類似的 P'-GRU-D 的更新公式：

(35)

(36)

(37)

(38)

這里公式35-37與2.2.3 小節(jié)中GRU-D神經(jīng)網(wǎng)絡預測模型中公式16-19表達含義相同，唯一不同的是加入了P'-LSTM-D定義中的時間門kt，其作用與之前描述相同.

2.4 損失函數(shù)

(39)

w代表在輸出時t時刻的權重參數(shù)，by是偏差參數(shù)，hnt是前幾節(jié)中隱藏狀態(tài)的輸出ht，這里的n僅注明是表示第n個病人的ht. 目標是想要使得下面的損失函數(shù)到達最小：

(40)

2.5 模型評價指標

在本文中，使用5個評價模型預測效果的指標，二分類問題中引入混淆矩陣，可參考維基百科上及參考文獻[15]的解釋，它們分別是：召回率(recall)；精度(precision)；準確度(accuracy)；F1分數(shù)(F1-measure)；AUC (AUROC).

其中，AUC是ROC(受試者工作特征曲)曲線下的面積，常常用來評價二分類問題的優(yōu)劣，本文主要看這個綜合指標.

3 真實數(shù)據(jù)描述

3.1 數(shù)據(jù)概述

在這項工作中，患有充血性心力衰竭(CHF)的患者的數(shù)據(jù)來自由美國麻省理工學院計算生理學研究室以及貝斯以色列迪康醫(yī)學中心和飛利浦醫(yī)療共同構建的重癥監(jiān)護醫(yī)學數(shù)據(jù)庫(MIMIC-III，v1.4)[16]. MIMIC-III包含2001年至2012年在波士頓貝斯以色列迪康醫(yī)學中心接受診治的40 000多名患者的電子病歷記錄數(shù)據(jù).MIMIC-III包含高分辨率數(shù)據(jù)，包括帶有時間戳的生命體征，實驗室值 ，治療指標等.其中，貝斯以色列迪康醫(yī)學中心和麻省理工學院的機構審查委員會批準了MIMIC數(shù)據(jù)庫的創(chuàng)建和使用.

圖5為從數(shù)據(jù)庫中初步篩選病人的過程圖，需要注意的是，本文僅納入患有充血性心力衰竭且進入ICU的成年病人，最終納入了10 315名病人數(shù)據(jù).

其中，表1為描述病人的基本情況:

3.2 臨床終點預測

數(shù)據(jù)庫中針對充血性心力衰竭病人所記錄的變量超過2 000個，本文做了一個排序處理，剔除數(shù)據(jù)規(guī)模不大的變量，先保留每個變量至少大于10萬條數(shù)據(jù)的標準，這樣做的目的是當轉換成模型能夠識別的數(shù)據(jù)框架時，缺失較大的變量能夠被過濾掉，這樣填補的數(shù)據(jù)比較少一點，從而保證數(shù)據(jù)的質量盡可能真實.

數(shù)據(jù)庫中單獨標記的變量通常對病人十分重要，針對充血性心力衰竭病人，MIMIC III數(shù)據(jù)庫中記錄的標記變量中，記錄特別完善的四個變量分別是：中心靜脈壓(CVP)；動脈收縮壓(ABPs)；呼吸頻率(RR)；心率(HR).這四個變量不僅完整的記錄著真實值，還記錄著該變量是處于正?；蚴钱惓５臓顟B(tài).由于前兩個變量都與血壓有關，后兩個變量均與呼吸有關，本文最終鎖定對充血性心力衰竭病人影響較為敏感的兩個結局，分別是呼吸頻率(RR)和中心靜脈壓(CVP)作為臨床終點.根據(jù)實際數(shù)據(jù)情況納入與臨床終點相關性最大的變量如下：

(1)針對呼吸頻率(RR)是否正常作為臨床終點參與訓練的自變量有26個：動脈血壓均值、中心靜脈壓、GCS評分、心率、吸入氣溫度、呼吸暫停時間間隔、無創(chuàng)血壓均值、肺動脈收縮壓、(呼吸機)吸氣壓力、先前體重、血氧飽和度、體溫(C)、血糖、血紅蛋白、凝血激活酶時間、血小板計數(shù)、白血細胞、肺動脈壓警報、動脈舒張壓、肺動脈舒張壓、平均血紅蛋白、紅細胞分布寬度、FiO2、氣道壓力、動脈收縮壓、無創(chuàng)血壓收縮壓;

(2)針對中心靜脈壓(CVP)正常與否作為臨床終點參與訓練的自變量有31個: 無創(chuàng)血壓收縮壓、動脈收縮壓、氣道壓力、FiO2、氧分壓、紅細胞分布寬度、平均血紅蛋白、尿素氮、肺動脈舒張壓、無創(chuàng)血壓舒張壓、動脈舒張壓、肺動脈壓警報、白血細胞、血小板計數(shù)、凝血激活酶時間、血細胞比容、肌酐、動脈CO2、潮氣量、體溫(C)、呼吸頻率、先前體重、(呼吸機)吸氣壓力、肺動脈收縮壓、無創(chuàng)血壓均值、分鐘量(每分鐘肺呼出氣體量)、吸入氣溫度、心率、GCS評分、動脈血壓均值、呼吸暫停時間間隔.

這里需要做一次時間點的篩選，將采取四舍五入的方式將周圍點均聚集到整點上，即我們最終呈現(xiàn)的數(shù)據(jù)都是以整點小時為單位，另外臨床終點不能有任何缺失，簡言之是需要以臨床終點數(shù)據(jù)y的時間點作為標準，抽取其余與y相對應的時間點，多余時間點需舍去.

4 結論

在能夠處理缺失的情況下，本文主要關注神經(jīng)網(wǎng)絡預測模型的效果以及模型的收斂速度這兩個重要指標.為了使模型具有可比性，本文中有意改變兩個變量：(1)模型類型；(2)時間序列長度，即模型納入的時間點總數(shù).

除此之外，本文盡可能保持其他設置參數(shù)相同，如：把大小(batch size)，即一次訓練的樣本數(shù)目固定為100；將時期(epochs)設置成100，其中一個時期代表所有訓練樣本的一個正向傳遞和一個反向傳遞；學習率(lr)固定為0.001；隱藏層固定為100等.具體的結果及結論如下：

圖6中記錄著以呼吸頻率作為臨床終點的模型表現(xiàn)效果，這里指AUC的值.從圖6的結果可以看出，LSTM-D和GRU-D的表現(xiàn)效果極佳，平均AUC能達到0.96，并且隨著時間序列長度的改變，表現(xiàn)效果改變不大，且LSTM-D和GRU-D表現(xiàn)不相上下.相比而言，P'-LSTM-D和P'-GRU-D的效果稍顯劣勢，平均AUC在0.86左右，所有病人的時間點超過60個后，由于服務器內(nèi)存限制，模型無法完成訓練，且隨著納入的時間點越多，表現(xiàn)效果越差.

為了再次驗證隨著時間序列長度的改變，LSTM-D和GRU-D的表現(xiàn)效果會不會出現(xiàn)較大的改變，故以CVP為臨床終點訓練，驗證集上的測試結果如圖7所示，隨著時間序列長度的改變，LSTM-D和GRU-D的表現(xiàn)效果也不相上下.

模型的收斂速度情況，可以從圖8結果中得知，在收斂速度這個指標上，GRU-D表現(xiàn)最佳，LSTM-D其次，而P'-LSTM-D和P'-GRU-D出人意料收斂較慢，或許是因為P'-LSTM-D和P'-GRU-D的模型參數(shù)太多.

為了考查所得到的GRU-D和LSTM-D收斂規(guī)律是不是偶然，本文再次查看以CVP為臨床終點的模型收斂速度，圖9結果亦表明，在GRU-D和LSTM-D表現(xiàn)效果相當?shù)那闆r下，GRU-D的收斂速度要快一些.故從四個模型比較結果來看，GRU-D模型更勝一籌，LSTM-D也很不錯，只是收斂速度稍微較GRU-D慢一點.

下面是以呼吸頻率為臨床終點的具體效果圖：

圖12 納入100個時間點時驗證表現(xiàn)效果圖Fig.12 Verified Performance Rendering with 100 Time Points

從以上結果可以很明顯的了解到模型的整體效果，這里僅展現(xiàn)納入時間點為10、30和100的效果圖，分別如圖10-12所示，所有的模型針對我們充血性心力衰竭病人的數(shù)據(jù)實驗上看，均出現(xiàn)GRU-D和LSTM-D神經(jīng)網(wǎng)絡預測模型效果很好，P'-LSTM-D和P'-GRU-D神經(jīng)網(wǎng)絡預測模型表現(xiàn)較差.

本文建議，數(shù)據(jù)長度不是很長的數(shù)據(jù)，可直接使用GRU-D和LSTM-D神經(jīng)網(wǎng)絡預測模型其中任何一個，它們不僅可以直接處理缺失，而且表現(xiàn)效果都極佳，在另一方面它們能夠處理不等距的數(shù)據(jù)，捕獲時間信息.