楊?；?王靜嫻 趙 芃 李業(yè)棉 陳方堯

西安交通大學(xué)醫(yī)學(xué)部公共衛(wèi)生學(xué)院流行病與衛(wèi)生統(tǒng)計(jì)學(xué)系衛(wèi)生統(tǒng)計(jì)學(xué)教研室，710061 陜西 西安

隊(duì)列研究是流行病學(xué)研究中一種重要的研究類型。在隊(duì)列研究中，研究人員常會對大量對象進(jìn)行長期隨訪，來收集與研究目的相關(guān)的信息、記錄特定終點(diǎn)事件的發(fā)生情況及發(fā)生時(shí)間[1]。然而，在研究中通常只有一部分個(gè)體在隨訪中發(fā)生終點(diǎn)事件，另一部分個(gè)體在隨訪期內(nèi)并未觀測到終點(diǎn)事件的發(fā)生，這樣的數(shù)據(jù)稱為縱向刪失數(shù)據(jù)(censored data)[2]。而基于刪失數(shù)據(jù)的預(yù)測模型，一直是研究人員感興趣的熱點(diǎn)。本文探討基于super learner算法的集成學(xué)習(xí)在縱向刪失數(shù)據(jù)預(yù)測建模中的應(yīng)用，為刪失數(shù)據(jù)預(yù)測提供參考。

1 數(shù)據(jù)預(yù)測建模中的集成學(xué)習(xí)

目前，基于機(jī)器學(xué)習(xí)(machine learning)的統(tǒng)計(jì)分析技術(shù)已廣泛應(yīng)用于生物醫(yī)學(xué)數(shù)據(jù)的預(yù)測建模分析中[3]，它是進(jìn)行自動或半自動信息(知識)獲取的有力工具，且利用機(jī)器學(xué)習(xí)算法及相關(guān)技術(shù)對高維數(shù)據(jù)進(jìn)行處理，可對傳統(tǒng)的統(tǒng)計(jì)分析起到輔助作用。目前，在生物醫(yī)學(xué)高維數(shù)據(jù)預(yù)測建模分析中，應(yīng)用較多的機(jī)器學(xué)習(xí)算法包括回歸學(xué)習(xí)(regression learning)、正則化算法(regularization)、決策樹算法(decision tree)、貝葉斯方法(bayesian method)、核方法(kernel method)、無監(jiān)督學(xué)習(xí)(unsupervised learning)和有監(jiān)督學(xué)習(xí)(supervised learning)等。

機(jī)器學(xué)習(xí)算法的應(yīng)用為生物醫(yī)學(xué)數(shù)據(jù)的處理和分析提供了新的、有利的工具，然而機(jī)器學(xué)習(xí)算法在實(shí)踐中存在以下2個(gè)方面的問題。

(1)不同算法在處理同一個(gè)數(shù)據(jù)時(shí)得到的結(jié)果不完全一樣，即不同算法的預(yù)測結(jié)果之間存在差異；這時(shí)就需要進(jìn)行算法或模型的選擇(selection)，但基于不同的選擇標(biāo)準(zhǔn)，依舊會導(dǎo)致選出的算法或模型不同，進(jìn)而影響預(yù)測結(jié)果的穩(wěn)健性。

(2)即使是同一個(gè)算法，如人工神經(jīng)網(wǎng)絡(luò)，在初值設(shè)置不同的情況下，最終的預(yù)測結(jié)果也不盡相同；同時(shí)，機(jī)器學(xué)習(xí)算法還容易導(dǎo)致過擬合的問題，使得預(yù)測結(jié)果的外推性變差。

為了解決上述問題，研究人員提出了集成學(xué)習(xí)(ensemble learning)方法[4]。相關(guān)研究指出，當(dāng)樣本量足夠大的時(shí)候，集成學(xué)習(xí)可以整合多個(gè)機(jī)器學(xué)習(xí)算法，先獨(dú)立地基于同一個(gè)樣本進(jìn)行訓(xùn)練，然后再把各個(gè)算法所得的結(jié)果整合起來，最終得到整體的預(yù)測結(jié)果，這改善了預(yù)測結(jié)果的穩(wěn)健性和外推性[4-5]。一個(gè)較為典型的集成學(xué)習(xí)方法的實(shí)例，就是隨機(jī)森林(random forest)，它通過將不同的決策樹模型進(jìn)行整合，構(gòu)成一個(gè)“森林”，最終實(shí)現(xiàn)分類預(yù)測的目的，其預(yù)測準(zhǔn)確度較決策樹模型有較大的提升[6]。

1.1 super learner算法

2006年Van der Laan等提出了一種基于損失函數(shù)的集成學(xué)習(xí)算法super learner[7]。該算法含有一個(gè)算法箱，其中包含了回歸學(xué)習(xí)、決策樹等多種機(jī)器學(xué)習(xí)算法。該算法將多種適用于不同情況的算法集合起來，并賦予不同的權(quán)重，在交叉驗(yàn)證(cross-validation)的框架下，找到使得損失函數(shù)最小的一套權(quán)重，以及在該權(quán)重下的算法組合?；静襟E如下。

Step 1：根據(jù)研究目的和數(shù)據(jù)情況，選擇并定義算法箱內(nèi)的a種系數(shù)估計(jì)方法；在完整數(shù)據(jù)集中，使用所選擇的算法，分別單獨(dú)基于完整數(shù)據(jù)進(jìn)行模型參數(shù)估計(jì)，得到a個(gè)系數(shù)估計(jì)的模型。

Step 2：將數(shù)據(jù)分為訓(xùn)練樣本和測試樣本后，進(jìn)行K層交叉驗(yàn)證。

Step 3：對訓(xùn)練集中的每一層數(shù)據(jù)使用算法箱中的全部估計(jì)方法，得到多套系數(shù)估計(jì)組合下的模型Ψa(X)。

Step 5：調(diào)整算法權(quán)重αka，得到使交叉驗(yàn)證中所得的預(yù)測結(jié)果最優(yōu)，即損失函數(shù)最小的算法權(quán)重：

需要注意的是，選擇不同的損失函數(shù)，得到的算法權(quán)重估計(jì)不全相同。通常情況下，super learner算法使用的損失函數(shù)包括絕對誤差損失函數(shù)、平方誤差損失函數(shù)、適用于二分類因變量的負(fù)對數(shù)損失函數(shù)3類。

Step 6：將調(diào)整后得到的最佳算法權(quán)重與Step1中得到的每一種算法下的模型參數(shù)估計(jì)進(jìn)行組合，得到最終的模型。

基于super learner算法進(jìn)行預(yù)測建模，優(yōu)點(diǎn)有： (1) 算法選擇靈活：可以自行選擇算法來擬合整體模型，且super learner算法會通過權(quán)重來刪除對整體預(yù)測能力無貢獻(xiàn)的算法(權(quán)重為0)；(2) 建模過程基于交叉驗(yàn)證方法，可有效提高模型預(yù)測的準(zhǔn)確度，同時(shí)降低模型過度擬合的風(fēng)險(xiǎn)；(3) super learner算法得到的模型，可以看作是不同算法下模型的加權(quán)組合，可以快速獲得原本需要花費(fèi)大量實(shí)驗(yàn)才能得出的模型權(quán)重，這使集合過程變得更加高效[7-8]。

1.2 基于super learner算法的刪失數(shù)據(jù)模型

super learner算法在數(shù)據(jù)刪失率較高或含有時(shí)間相依型協(xié)變量(time-dependent variable)時(shí)，可以更好地利用數(shù)據(jù)進(jìn)行建模，且預(yù)測能力更高。在使用super learner算法進(jìn)行刪失數(shù)據(jù)建模時(shí)與常規(guī)的cox比例風(fēng)險(xiǎn)模型不同，該算法建模時(shí)將刪失時(shí)間等距地劃分成時(shí)間段，作為應(yīng)變量納入模型；同時(shí)，針對終點(diǎn)事件發(fā)生的時(shí)間點(diǎn)建立一個(gè)指示變量，如果時(shí)間段被劃分的足夠精細(xì)，則包含有每個(gè)時(shí)間點(diǎn)指標(biāo)的參數(shù)統(tǒng)計(jì)模型，如logistic回歸將近似于cox比例風(fēng)險(xiǎn)模型[9]。

2 super learner算法的R語言實(shí)現(xiàn)

super learner算法在R語言中的實(shí)現(xiàn)，可以通過R包“SuperLearner”(https://CRAN.R-project.org/package=SuperLearner)進(jìn)行。該軟件包的當(dāng)前版本為2.0-26，發(fā)布時(shí)間為2019年12月。應(yīng)用該軟件包實(shí)現(xiàn)基于super learner算法預(yù)測建模的基本語法和參數(shù)如下：

2.1 預(yù)測模型的建立

使用R軟件包“SuperLearner”中的“SuperLearner”函數(shù)進(jìn)行建模時(shí)，函數(shù)的基本參數(shù)和格式如下：

#模型的建立>Fit<- SuperLearner(Y,#應(yīng)變量,可以是二分類變量或連續(xù)型變量>X,#自變量矩陣>family = ,#鏈接函數(shù),目前可選gaussian()或binomial()>SL.library=,#算法參數(shù)>method=<"METHORD">#定義損失函數(shù)>)

其中，參數(shù)Y為應(yīng)變量，且必須是數(shù)值型向量；X是自變量的集合，通常是一個(gè)數(shù)據(jù)集。

函數(shù)參數(shù)“method”用來指定損失函數(shù)類型，可以選擇的損失函數(shù)包括基于lawson-hanson算法的最小二乘損失函數(shù)(method.NNLS，默認(rèn)損失函數(shù))、基于goldfarb和idnani雙重算法的最小二乘損失函數(shù)(method.NNLS2)、負(fù)對數(shù)二項(xiàng)似然損失函數(shù)(method.NNloglik)和1-AUC損失函數(shù)(method.AUC)、基于二次規(guī)劃算法來計(jì)算權(quán)重組合的損失函數(shù)(method.CC_LS)和通過最小化負(fù)二項(xiàng)式對數(shù)似然確定最佳權(quán)重組合的損失函數(shù)(method.CC_nloglik)。

參數(shù)“SL.library”為算法參數(shù)，用來指定該super learner模型中希望組合的算法。super learner算法的算法箱中包含多個(gè)可選算法，如隨機(jī)森林算法等，使用相應(yīng)的算法時(shí)需要R已經(jīng)安裝了相應(yīng)的R程序包，否則在調(diào)用算法時(shí)會報(bào)錯(cuò)。同時(shí)，該函數(shù)也允許使用者調(diào)用自己設(shè)定的算法，或?qū)浖幸讯x算法的超參數(shù)(hyperparameter)進(jìn)行自定義。相關(guān)研究認(rèn)為，在使用super learner算法進(jìn)行建模分析時(shí)，根據(jù)數(shù)據(jù)情況和分析目的，對算法超參數(shù)進(jìn)行適當(dāng)?shù)恼{(diào)整，會對建模結(jié)果帶來有利的影響[10]。

需要注意的是，基于同一個(gè)數(shù)據(jù)集，選擇不同的算法組合或損失函數(shù)，建立得到的模型很可能是不同的。需要特別說明的是，選擇1-AUC損失函數(shù)時(shí)，應(yīng)變量Y需要是二分類變量。

2.2 模型評價(jià)

模型評價(jià)采用該軟件包中的另一個(gè)函數(shù)“CV.SuperLearner”來實(shí)現(xiàn)，其基本格式和參數(shù)的含義如下：

#模型評價(jià)>Fit<- CV.SuperL-earner(Y,#應(yīng)變量,可以是二分類變量或連續(xù)型變量>X,#自變量矩陣>V=,#交叉驗(yàn)證次數(shù),默認(rèn)為10次>family=,#鏈接函數(shù),目前可選gaussian()或binomial()>SL.library=, #算法參數(shù)>method=<"METHORD">#定義損失函數(shù)>)

其中，參數(shù)“V”為交叉驗(yàn)證的折數(shù)，不進(jìn)行指定時(shí)，默認(rèn)為10折。交叉驗(yàn)證的折數(shù)不宜過高，否則會將數(shù)據(jù)過度分割。該函數(shù)的其余參數(shù)，含義與2.1部分中相同。模型驗(yàn)證的結(jié)果，可以采用“summary()”函數(shù)和“plot()”函數(shù)進(jìn)行輸出。

3 實(shí)例分析

3.1 數(shù)據(jù)來源

實(shí)例數(shù)據(jù)獲取自TCGA數(shù)據(jù)庫(https://cancergenome.nih.gov/)[11]，研究所使用的數(shù)據(jù)集包含了結(jié)腸腺癌患者的生存結(jié)局及其microRNA測序結(jié)果。該數(shù)據(jù)刪失率為88.8%，屬于高刪失數(shù)據(jù)。參考2019年Chen等[12]的分析結(jié)果，選取其中15個(gè)在單因素分析中有統(tǒng)計(jì)學(xué)意義的microRNA納入分析，采用super learner算法建立預(yù)測模型。

3.2 數(shù)據(jù)預(yù)處理

在建模前首先需對刪失時(shí)間進(jìn)行處理。應(yīng)用super learner算法對刪失數(shù)據(jù)建模前，必須將數(shù)據(jù)從常規(guī)的結(jié)構(gòu)轉(zhuǎn)換為計(jì)數(shù)過程(counting process)格式，即需要對時(shí)間變量進(jìn)行拆分，使得在給定的時(shí)間間隔下，每個(gè)小的時(shí)間區(qū)間內(nèi)最多僅有1個(gè)終點(diǎn)事件發(fā)生[9]。

在較早期版本的R里，“SuperLearner”包中的“createDiscrete()”函數(shù)可以完成該轉(zhuǎn)換，但更新后的版本中，該函數(shù)已經(jīng)不可使用?？梢钥紤]使用“survival”包中的“survSplit()”函數(shù)達(dá)到相同的效果，亦可基于R語言自行編程進(jìn)行數(shù)據(jù)格式的轉(zhuǎn)換處理。

建模時(shí)，將新產(chǎn)生的區(qū)間間隔的起點(diǎn)或終點(diǎn)作為自變量，不再將原有的時(shí)間變量納入分析，模型的應(yīng)變量為小區(qū)間內(nèi)事件的發(fā)生與否。

3.3 建模分析

經(jīng)過預(yù)處理后，首先將數(shù)據(jù)集按照7∶3的比例分割為訓(xùn)練集和測試集：

# 將數(shù)據(jù)集(s)分割為訓(xùn)練集(70%)和驗(yàn)證集(30%)

>num <- round(nrow(s)*0.7)

>index <- sample(nrow(s),size = num)

>s_train <-s[index,] #訓(xùn)練集

>s_test <-s[-index,] #驗(yàn)證集

#基于訓(xùn)練樣本建立模型

#指定結(jié)局變量和自變量集合

y_train <- as.vector(s_stran$y) #定義結(jié)局變量，應(yīng)為數(shù)值型向量

x_train <-s_train[,-1] #定義自變量集合，應(yīng)為數(shù)據(jù)集

#模型建立

>model <- SuperLearner(y_train,x_train, family=binomial(),

> SL.library=list("SL.ranger","SL.ksvm","SL.ipredbagg","SL.bayesglm"))

訓(xùn)練結(jié)果見圖1：

圖1 基于super learner算法建模時(shí)的算法權(quán)重和對應(yīng)的風(fēng)險(xiǎn)值

以上輸出結(jié)果中，Risk為每一種算法下?lián)p失函數(shù)的值，默認(rèn)Coef為每個(gè)算法最終獲得的權(quán)重。由該結(jié)果可知，所選的4個(gè)算法中，“SL.ranger”權(quán)重為0.622 598 9，“SL.ipredbagg”算法權(quán)重為0.377 401 1，二者之和為1；而其余2個(gè)算法權(quán)重為0，被判定為無效算法。以上結(jié)果相當(dāng)于，當(dāng)前模型是基于“SL.ranger”和“SL.ipredbagg”2個(gè)算法得到的。

3.4 模型損失函數(shù)值

以下是模型交叉驗(yàn)證的結(jié)果。

#模型交叉驗(yàn)證，默認(rèn)為10折

>CV.model <- CV.SuperLearner(y_train, x_strain, family=binomial(),

> SL.library=list("SL.ranger","SL.ksvm","SL.ipredbagg","SL.bayesglm"))

以下2個(gè)函數(shù)分別輸出不同算法下，模型損失函數(shù)的量化(圖2)和圖形化結(jié)果(圖3)。

#輸出交叉驗(yàn)證結(jié)果

>summary(cv.model)

圖2 不同算法模型損失函數(shù)值(Risk)的量化結(jié)果

#將結(jié)果以統(tǒng)計(jì)圖的形式體現(xiàn)

>plot(cv.model)

圖3 不同算法下模型損失函數(shù)值(Risk)及其95%CI的圖形化結(jié)果

3.5 基于測試集的模型評價(jià)

模型建立后，可以基于拆分好的測試集，對模型的預(yù)測能力進(jìn)行檢驗(yàn)，以下過程同時(shí)需要加載R包“dplyr”和“caret”，方法如下：

#基于測試集進(jìn)行預(yù)測

#指定結(jié)局變量和自變量集合

y_test <- as.vector(s_test$y) #定義結(jié)局變量，應(yīng)為數(shù)值型向量

x_test <- s_test[,-1] #定義自變量集合，應(yīng)為數(shù)據(jù)集

>predictions <- predict.SuperLearner(model, newdata=x_test)

>conv.preds <- ifelse(predictions$pred>=0.5, 1, 0)

# 輸出預(yù)測結(jié)果

>cm <- confusionMatrix(factor(conv.preds), factor(y_test))

可以得到，當(dāng)前模型的預(yù)測準(zhǔn)確度為0.873 7(95%CI：0.789 7～0.933 0)?；谟?xùn)練集的模型預(yù)測結(jié)果顯示，模型的預(yù)測準(zhǔn)確度良好。此外，我們也可以估計(jì)模型的C-index[13]進(jìn)行模型預(yù)測結(jié)果的評價(jià)，C-index的估計(jì)需要調(diào)用R包“Hmisc”來估計(jì)。

我們同時(shí)使用常規(guī)的cox回歸模型基于以上數(shù)據(jù)建模，得到相對應(yīng)回歸模型的C-index為0.784，而基于super learner算法的到的預(yù)測模型，其C-index為0.883，比常規(guī)的cox回歸模型高。

4 討論

super learner算法得到的預(yù)測模型，可以看作為基于其算法箱內(nèi)的1個(gè)或多個(gè)算法得到預(yù)測模型的加權(quán)組合，預(yù)測準(zhǔn)確度高。因此，近年來，在生物醫(yī)學(xué)的預(yù)測建模分析中，super learner算法的使用越來越多[14]。

在處理刪失數(shù)據(jù)時(shí)，使用該算法及R包“SuperLearner”對刪失數(shù)據(jù)進(jìn)行預(yù)測建模與常規(guī)的cox回歸不同，具體的不同之處在于，預(yù)測得到終點(diǎn)事件的發(fā)生概率，是矯正了事件發(fā)生時(shí)間后的結(jié)果。由于該軟件包暫時(shí)不支持直接導(dǎo)出回歸系數(shù)的估計(jì)，因此目前使用該軟件包進(jìn)行基于super learner算法的集成學(xué)習(xí)建模，更適合以預(yù)測為目的的分析。同時(shí)，相關(guān)模擬研究也證明，在刪失數(shù)據(jù)的預(yù)測分析中，super learner算法的表現(xiàn)都優(yōu)于其他相關(guān)算法和模型[15]。

綜上所述，基于super learner的集成學(xué)習(xí)算法，為縱向刪失數(shù)據(jù)的預(yù)測建模提供了一個(gè)新的選擇，其預(yù)測效果良好，對于數(shù)據(jù)分析實(shí)踐來說，是良好的統(tǒng)計(jì)預(yù)測建模工具。