許樹紅　董曉強　陶　然　高　雪　高　倩　虞明星　王　彤△

【提　要】　目的　基于LASSO-Cox模型探索交叉驗證(cross validation)、pcvl法(penalized cross-validated log-likelihood)、EBIC準則(extended bayesian information criterion)、平穩(wěn)選擇(stability selection)四種方法在控制FDR(false discovery rate)方面的表現(xiàn)及其變量選擇效果。方法　通過模擬研究評價各方法在不同刪失比例、自變量間不同相關程度以及回歸系數(shù)的不同稀疏水平下的FDR和PSR(positive select rate)，并從GEO上下載DLBCL數(shù)據(jù)進行基因與預后間的關聯(lián)分析。結果　模擬結果表明，在不同刪失比例、自變量相關程度和稀疏水平的情況下，平穩(wěn)選擇法控制FDR的能力都優(yōu)于其他方法且其變量選擇效能也較高。EBIC準則在相關程度低、自變量較稀疏時表現(xiàn)較好，當樣本量較小時結果較保守。pcvl法雖然不容易漏掉有效應的變量，但其FDR仍較高。實例結果顯示，EBIC準則只選出1個基因，平穩(wěn)選擇法選出的基因中大部分有統(tǒng)計學意義且與其他方法的結果重合度高。結論　在基于LASSO-Cox模型的高維數(shù)據(jù)生存分析中平穩(wěn)選擇法能較好地控制FDR且其變量選擇效能也較高。

LASSO(least absolute shrinkage and selection operator)作為一種常用的懲罰類方法，在GWAS研究、測序數(shù)據(jù)分析中應用廣泛且發(fā)展迅速，其中心思想是將模型系數(shù)的絕對值函數(shù)作為懲罰項對模型中變量的系數(shù)進行壓縮，使得一些弱效應變量的系數(shù)變小，甚至被壓縮為0，以提供一個稀疏解[1]。LASSO類方法通過直接估計不為零的回歸系數(shù)來實現(xiàn)變量選擇，它不同于經(jīng)典的逐步回歸，并沒有采用假設檢驗的步驟，并且其變量選擇和參數(shù)估計的結果在很大程度上受調整參數(shù)的影響，然而調整參數(shù)的傳統(tǒng)選擇方法如交叉驗證法和信息準則法，在選擇調整參數(shù)時主要考慮模型的擬合以及模型的復雜程度，亦未考慮Ⅰ類錯誤估計和控制問題。目前已有研究者指出LASSO類方法存在較高的假陽性問題[2-3]，因此，如何在控制或降低假陽性的前提下，選擇合適的調整參數(shù)從而選出正確的模型尤為重要。

目前已有一些研究在多重校正以及LASSO類方法的基礎上，發(fā)展了一些控制Ⅰ類錯誤同時篩選出正確變量的方法[4]。FDR作為測序數(shù)據(jù)多重比較領域中控制Ⅰ類錯誤的常用指標，由Benjamini& Hochberg在1995年提出，指拒絕的假設檢驗結果中Ⅰ類錯誤所占比例的期望，在整體上控制Ⅰ類錯誤的同時能篩出更多有意義的基因[5]。本文將基于LASSO-Cox模型，通過模擬和實例研究比較交叉驗證法、pcvl法、EBIC準則和平穩(wěn)選擇法四種調整參數(shù)選擇方法在控制FDR方面的表現(xiàn)以及在高維數(shù)據(jù)生存分析應用中的優(yōu)劣。

原理與方法

LASSO-Cox模型的最大化目標函數(shù)為：

(1)

其中p表示協(xié)變量的維數(shù)，β=(β1,β2,…,βp)T為各協(xié)變量對應的回歸系數(shù)，l(β)為Cox模型的對數(shù)偏似然函數(shù)，λ>0為調整參數(shù)。調整參數(shù)λ的大小影響著模型的復雜程度和收斂速度，λ值過大可能會使最終的模型中協(xié)變量個數(shù)過少，導致重要的變量被遺漏，反之，λ值若過小可能會使最終的模型中協(xié)變量個數(shù)過多，使得模型包含很多噪聲變量，導致結果中FDR過高。本課題組前期研究中已對基于LASSO的Ⅰ類錯誤控制方法的基本原理進行了詳細闡述[4]，本文對交叉驗證法(cross validation,CV)、pcvl法(penalized cross-validated log-likelihood)、EBIC準則(extended bayesian information criterion)、平穩(wěn)選擇法(stability selection)的原理進行簡單介紹。

1.交叉驗證

目前最常用的交叉驗證是K折交叉驗證(K-fold cross-validation)，K為整數(shù)，1≤K≤n，通常取5或者10。該方法利用訓練集來建立模型并求得回歸參數(shù)估計值，然后用該回歸參數(shù)估計值來預測驗證集。對于對數(shù)偏似然函數(shù)，其目標函數(shù)為：

(2)

2.pcvl

Ternès[6]針對傳統(tǒng)交叉驗證法容易出現(xiàn)過度擬合這一問題提出pcvl法，公式如下：

pcvl(λ)=cvl(λ)-pen(λ)

(3)

“體育是人的類存在的主體象征”[11]，是人類主體本性需求下的產(chǎn)物，是本體體驗與身體經(jīng)驗的結合，體育不僅是肉體強蠻的塑造之法，同時也是主體自由的解放之徑。體育行為以對生命的養(yǎng)護為最終的實踐目的，以主體的自由為最高的價值旨歸，見證人類個體自我造就的每一時刻。以銅為鏡，可以正衣冠；以史為鏡，可以知興替；以人為鏡，可以明得失。體育相當于反觀人自身的明鏡，在體育競賽中反觀自我的品性，在體育訓練中反觀肉體的健康，在體育游戲中反觀主體的自由。在體育的實踐行為當中認識自我，發(fā)揮自身的優(yōu)點，摒棄自身的不足，避免自身的異化。

(4)

圖1　cvl(λ)和pcvl(λ)以及模型中非零自變量個數(shù)隨λ變化趨勢圖(n=100,p=1000,q=6,L=3時的模擬結果)

3.EBIC準則

(5)

當γ=0時，EBIC準則等于傳統(tǒng)的BIC準則；當γ=1時，EBIC準則與mBIC準則類似。γ的取值大小影響著EBIC準則的選擇一致性。Luo和Chen[8]基于線性模型對EBIC準則的選擇一致性進行了深入的探討。

4.平穩(wěn)選擇法

(6)

模擬試驗

1.模擬試驗設置

參考Luo、Song和Fan 等的模擬情形[13-15]，模擬試驗中自變量個數(shù)p=1000，樣本量n=(100,120,140,160,180,200)，模擬次數(shù)B=100。

模擬方案一：探索生存數(shù)據(jù)刪失比例的不同和自變量間的相關程度不同對調整參數(shù)選擇方法的影響，具體設置如下：

(7)

(8)

(3)刪失時間是由均數(shù)為Uexp(-xβ)的指數(shù)分布產(chǎn)生，其中U是[1,L]上的均勻分布，不同的L值對應的刪失比例不同，設置L=(2,3,4,5)。

模擬方案二：研究自變量的不同稀疏情況對調整參數(shù)選擇方法的影響，具體設置如下：

(1)自變量服從均數(shù)為0的正態(tài)分布，自變量間相關結構與模擬方案一的區(qū)組相關結構相同。

(3)刪失時間的指數(shù)分布中L=3。

2.模擬試驗分析與評價

利用R 3.3.2軟件中glmnet package建立LASSO-Cox模型，分別采用CV法、pcvl法、EBIC準則、平穩(wěn)選擇法進行調整參數(shù)的選擇和變量篩選。利用cv.glmnet函數(shù)進行10折交叉驗證分析。EBICγ1準則中γ取值為1，EBICγ2準則中γ取值略大于1-lnn/2lnp(在其基礎上向上取整)。pcvl法由于還沒有R軟件包可實現(xiàn)，本次研究主要參考Ternès，Rotolo 和 Michiels提供的R程序[6]。利用lol package進行平穩(wěn)選擇，平穩(wěn)選擇的重復抽樣次數(shù)設為100，截斷點πthr取0.6。

評價指標主要為錯誤發(fā)現(xiàn)率(FDR)和PSR(positive select rate)。FDR作為變量選擇的假陽性指標。PSR作為變量選擇效能指標，用來評價模型發(fā)現(xiàn)確實存在的有效應變量的能力。FDR和PSR定義如下：

FDR=FP/(TP+FP)

(9)

PSR=TP/(TP+FN)

(10)

其中FP(false positive)表示在模擬試驗的真實模型中回歸系數(shù)為零，但估計的回歸系數(shù)是非零的個數(shù)；TP(true positive)表示在模擬試驗的真實模型中系數(shù)是非零，估計的結果也是非零的個數(shù)；FN(false negative)表示在模擬試驗的真實模型中回歸系數(shù)為非零，但估計的結果為零的個數(shù)。FDR越接近0、PSR越接近1，說明該調整參數(shù)選擇方法越好。

3.模擬試驗結果

(1)如圖2和圖3，在樣本量、刪失比例、相關系數(shù)均一定的情況下，各方法的FDR從低到高依次為：平穩(wěn)選擇法≤EBICγ1準則

(2)如圖4和圖5，在樣本量、自變量間相關系數(shù)和稀疏水平均一定時，各方法的FDR相比較結果以及PSR相比較結果與方案一的結果一致。隨著真實非零自變量個數(shù)的增多即稀疏水平的降低，平穩(wěn)選擇法的FDR基本保持不變，pcvl法的FDR略升高且與CV法之間的差距逐漸減小，EBIC準則對稀疏水平的變化比較敏感，其結果波動較大；在樣本量較小且真實非零自變量個數(shù)增多的情況下，各方法的PSR值均出現(xiàn)不同程度的降低，平穩(wěn)選擇法的PSR值始終高于EBIC準則。在樣本量較大的情況下，隨著真實非零自變量個數(shù)的增多，CV法、pcvl法和平穩(wěn)選擇法的PSR值基本保持不變且接近于1。

圖2　不同刪失比例和相關情況下模型變量選擇的FDR(100次模擬的平均值)

實例分析

1.數(shù)據(jù)來源及整理

圖3　不同刪失比例和相關情況下模型變量選擇的PSR(100次模擬的平均值)

2.實例分析結果

(1)這些基因之間存在比較高的相關關系，每個變量與其他變量間的最大相關系數(shù)的平均值是0.810。

(2)上述方法選出與預后有關的基因數(shù)分別有51、28、13、1個，其中EBICr1準則和EBICr2準則結果一致。分別以各方法選出的基因為協(xié)變量，擬合Cox比例風險模型，得到的模型擬合結果如表1所示：平穩(wěn)選擇法選出的13個基因中有11個基因具有統(tǒng)計學意義，且與CV法相同的基因有12個，與pcvl法相同的基因有10個；pcvl法選出的28個基因中12個有統(tǒng)計學意義，與CV法相同的基因有26個，與EBIC準則相同的基因有1個。CV法選出的51個基因中只有16個基因具有統(tǒng)計學意義。其中MMP12、CXCL2、P2RY12、MAL基因已有文獻報道與DLBCL 的預后和疾病進展有關[17-21]。值得注意的是有文獻報道CD163基因與DLBCL預后有關[22]，但只有pcvl法選出了CD163。實例結果表明CV法選出變量很多且大部分無統(tǒng)計學意義。pcvl法一定程度上減少了假陽性，且不容易漏掉有效應的自變量，但選出的變量仍比較多。CV法、pcvl法選出的有統(tǒng)計學意義的基因大部分被平穩(wěn)選擇法選出，且占后者結果的絕大部分，說明平穩(wěn)選擇法選出的自變量假陽性很低，且與其他方法結果重合度高。EBIC準則結果則過于保守。

圖4　自變量不同稀疏水平下模型變量選擇的FDR(100次模擬的平均值)

討　　論

模擬研究和實例數(shù)據(jù)分析結果顯示傳統(tǒng)的CV法假陽性問題確實比較嚴重，原因是CV法在選擇調整參數(shù)時只考慮了模擬的擬合優(yōu)度未考慮模型的稀疏性。模擬結果顯示pcvl法既能將PSR保持在較高水平又能一定程度地減少FDR，同時其結果受自變量間相關的影響較小，但與EBIC準則、平穩(wěn)選擇法相比其FDR仍然比較高，而且當真實有效應的自變量較多時pcvl法控制FDR的能力略微下降，可能是因為pcvl法通過在CV法的目標函數(shù)中加入懲罰項達到擬合優(yōu)度與稀疏程度的折中，其結果一定程度上依賴于CV法。

EBIC準則γ參數(shù)取值不同時結果不相同，當γ=1時其控制FDR的能力與平穩(wěn)選擇法相差不大，但其PSR在小樣本時比較低，提示樣本量很小時可能會遺漏重要的變量；與γ=1的結果相比，γ取值略大于1-lnn/2lnp時的PSR較高但其FDR也比較高；EBIC準則結果受自變量間相關和稀疏水平影響較大；實例研究結果中EBIC準則只選出1個基因，遺漏掉的重要變量比較多。Luo等[15]將EBIC準則應用于SIS-Adaptive LASSO的調整參數(shù)的選擇，在其實例分析中γ取值越小篩選的基因數(shù)越多，當0.3≤γ≤0.7時只篩選出兩個基因，說明該方法確實比較保守。

平穩(wěn)選擇法將FDR控制在低水平的同時仍使PSR保持在可接受水平。該方法降低了模型變量選擇結果對調整參數(shù)值的依賴性，因此在自變量間不同相關程度和不同稀疏水平的情況下表現(xiàn)均比較穩(wěn)定，雖然在樣本量小且有意義變量多的情況下效能有所降低，但仍高于EBIC準則。

圖5　自變量不同稀疏水平下模型變量選擇的PSR(100次模擬的平均值)

本研究主要基于LASSO-Cox模型進行四種方法變量選擇的比較和評價，然而這些方法不僅可以用于LASSO調整參數(shù)的選擇，還可以應用于其他懲罰類方法如adaptive LASSO、SCAD(smoothly clipped absolute deviation)、MCP(minimax concave penalty)等，這部分內容將在今后的研究工作中進一步探討。在利用變量選擇方法篩選變量時，除了考慮其統(tǒng)計學意義，還需要結合問題本身的專業(yè)知識和實際意義，尤其是面對目前日益龐大的生物學數(shù)據(jù)，利用變量選擇方法篩選出很有可能與癌癥患者生存有關聯(lián)的候選基因后，其具體的生物學意義仍需要進一步實驗加以驗證。

表1　 各方法篩選出的基因

注：粗體字表示平穩(wěn)選擇法、pcvl法、CV法選出的相同基因；斜體字表示pcvl法、CV法選出的相同基因；斜體字加粗表示平穩(wěn)選擇法、CV法選出的相同基因；字體背景灰色表示EBIC準則和pcvl法選出的相同基因；***表示P值小于0.001，**表示P值小于0.01，*表示P值小于0.05