徐州醫(yī)科大學公共衛(wèi)生學院流行病與衛(wèi)生統計學系(221004) 楊家驥 余星皓 曾 平 黃水平

【提 要】 目的 將整合通路信息的sparse group LASSO方法與近年來發(fā)表的表型預測方法進行比較，通過模擬各種復雜疾病可能的遺傳結構，比較各方法的預測能力，期望通過TCGA數據找到高效和穩(wěn)健的統計方法。方法 本研究利用SGL方法整合基因途徑信息和基因表達數據，并與傳統模型(LASSO、Enet、GSSLASSO)進行比較。通過乳腺癌真實基因型數據模擬表型數據：考慮不同分組(分組k=50，200，300，328)和不同遺傳度對模型的影響(遺傳度h2=0.3，0.5，0.8)。采用相關系數R評價幾種模型的預測能力，進一步通過結直腸癌(CRC)、胰腺癌(PAAD)、乳腺癌(BRCA)三個真實數據比較各方法表型預測的準確性。結果 模擬結果表明，隨著遺傳度的增高，各方法的預測準確性也逐漸增高。整合通路信息的SGL方法和GSSLASSO方法比傳統的LASSO和Enet方法有著更高的預測精度。而兩種整合通路信息的方法中，SGL方法有著更好的預測能力和穩(wěn)定性。在50，200，300分組情況下，GSSLASSO預測效果和LASSO以及Enet相近，但是在考慮通路信息的328分組下，GSSLASSO表現出了較好的預測效果。實例數據分析CRC，PAAD數據中，SGL方法具有最優(yōu)的預測精度，其次是GSSLASSO，LASSO和Enet方法預測效果最差。結論 整合通路信息的預測方法預測效果明顯優(yōu)于一般模型，而無論是在模擬數據還是實例數據中SGL的方法具有最優(yōu)的預測精度。

生物技術的飛速發(fā)展產生了大量高通量測序數據，這不僅僅為研究遺傳因素與復雜疾病、特征之間的關系提供了廣闊的空間，并且可以通過整合一系列組學信息，進一步促進了復雜表型遺傳風險預測和評估的發(fā)展[1-6]。與傳統預測模型不同，遺傳預測模型研究會帶來數據高維度(變量的數目p要遠遠大于樣本量n)的問題，這就使得傳統的分類和預測方法預測精度下降，計算負擔加重[7]。針對高維數據，研究者通常會利用正則化的方法來提高統計模型的預測準確性和可解釋性(例如LASSO、Elastic net)，通過增加一個l1或者l2懲罰項對一部分模型系數進行壓縮，以達到變量選擇的目的。這些方法被廣泛應用于大規(guī)模分子數據的疾病預測和診斷中[8-10]。

在遺傳預測方面，近幾年研究者們提出了許多利用分組信息進行預測的方法，絕大多數是對LASSO方法進行改進的模型選擇方法，但這些正則化的方法對分組進行懲罰不可避免地會導致遺傳信息的丟失。例如Yuan和Lin提出的group LASSO的方法，該方法首先將所有變量分組，然后在目標函數中懲罰每一組的l2范數，這樣就可以將一整個組剔除[11]。2010年Friedman提出了一種稀疏分組LASSO(sparse group LASSO，SGL)的方法[12]，這種方法對其分組和組內變量均進行正則化，以達到變量選擇和模型選擇雙重目的[13]。另外，研究者們也發(fā)展了多種利用外部分組信息的高維數據分析方法。Tang等人在group LASSO基礎上提出了一種分組的穗和板套索廣義線性模型(group spike-and-slab LASSO，GSSLASSO)，該方法發(fā)現在模型擬合中納入KEGG通路信息，可以有效地提高預測的準確性[14]。本研究將SGL模型應用于連續(xù)型的高維遺傳數據中，進一步整合KEGG通路信息，利用模擬研究和真實數據分析與整合分組信息的GSSLASSO模型以及不考慮分組信息的LASSO方法[15]、Enet模型[16]比較，評價其預測精度及穩(wěn)定性。

方法與材料

1.方法

LASSO是把一個懲罰項加到回歸系數絕對值之和上，使其滿足總和小于等于一個常數的約束條件，它通過構造一個罰函數得到一個較為精簡的模型，使得一些系數被壓縮，使殘差平方和最小化，從而能夠產生某些嚴格等于0 的回歸系數，最終得到一個解釋力較強的模型。使用LASSO的原因主要有兩個：一是為了提高模型的預測精度，通過將一系列回歸系數設置為0，使得預測值的方差減少，因此可以提高整體的預測精度；二是為了滿足模型的可解釋性，通過變量選擇的方法找出影響較大的變量。Enet是一種LASSO與嶺回歸組合后的回歸分析[17-18]，即將嶺回歸引入的l2正則項與LASSO回歸引入的l1正則項組合，通過構造罰函數，使得殘差平方最小，最終得到合理的模型。一方面達到了嶺回歸對重要特征選擇的目的，另一方面又像LASSO回歸那樣，刪除了對因變量影響較小的特征，取得了很好的效果。GSSLASSO回歸是Tang等人[14]提出的一種納入了外部注釋信息的模型。該模型在系數上引入一個新的先驗分布，即混合尖峰和平板雙指數先驗。其通過自適應的調整收縮量來提高系數估計和預測的準確性，根據不同的基因表達與表型數據，產生合適的收縮系數，去除與疾病不相關的基因，同時保留系數較大的基因，并把期望最大化步驟整合到循環(huán)坐標下降算法中，可以很好地識別重要的預測因子并從大量的候選分組中構建有效的預測模型。稀疏組LASSO(SGL)在分組LASSO的基礎上增加一個l1懲罰項，既考慮到組內系數的稀疏性，也考慮到分組的稀疏性來決定選擇重要的分組；類似于彈性網方法，參數α通常設置為0～1，用于在組LASSO(α=0)和LASSO(α=1)之間建立聯系。

另外，在基因數據中同一個基因可能屬于不同的通路，因此采用復制變量的方法對屬于不同通路的基因進行復制，以達到更好的預測效果。

2.數據來源及質量控制

(1)模擬數據

(2)實例數據

數據全部來源于加利福尼亞大學基因組瀏覽器UCSC Xena(https：//xenabrowser.net/)，下載其中腫瘤基因圖譜數據庫(TCGA)，一共包括三份癌癥數據集，即乳腺癌(BRCA)、結直腸癌(CRC)、胰腺癌(PAAD)。數據包括癌癥患者的臨床數據和RNAseq基因表達水平數據。對于每種癌癥，首先合并從原發(fā)癌組織測量的臨床數據和基因表達水平的數據；然后移除了零表達值超過50%的基因并對剩余的基因表達數據進行標準化。

乳腺癌：原始數據包括1247例患者的臨床數據和1218例患者的20530基因表達數據，對兩份數據進行合并，刪除重復的患者和男性患者，同時刪除零表達值超過50%的基因，最終獲得1083例患者的17675個基因表達數據。

結直腸癌：原始數據包括736例患者的臨床數據和434例患者的20530基因表達數據，對兩份數據進行合并，刪除重復的患者，同時刪除零表達值超過50%的基因，最終獲得275例患者的17493個基因表達數據

胰腺癌：原始數據包括196例患者的臨床數據和183例患者的20530基因表達數據，對兩份數據進行合并，刪除重復的患者，同時刪除零表達值超過50%的基因，最終獲得178例患者的18009個基因表達數據。

表1 TCGA數據集中的每種癌癥的樣本大小和基因數量

3.統計分析

文中所有分析均使用R 3.5.2軟件，LASSO、Enet使用glmnet(version 2.0-16)軟件包，通過100折交叉驗證選擇最優(yōu)懲罰參數，Enet設置α為0.05；GSSLASSO使用BhGLM(version 1.1.0)軟件包，設置s1為1，通過十折交叉驗證從s0=0.01×m，m=0.1，1，2，…，9中選擇最優(yōu)s0；SGL使用SGL(version 1.2)軟件包?；騅EGG通路注釋使用clusterProfiler軟件包。本研究采用100次五折交叉驗證評估模型的預測精度，每次交叉驗證隨機抽取80%的數據集作為訓練集，剩余20%作為測試集；在訓練集數據中擬合預測模型，并在測試集中對連續(xù)表型進行預測，通過相關系數R評估預測性能。

結 果

1.模擬研究

圖1是遺傳度分別為0.3、0.5、0.8，三種模型與SGL的預測能力比較，預測性能由相關系數R衡量；每種方案重復100次。圖1可以看出，與其他方法相比(LASSO、Enet、GSSLASSO)，在不同的遺傳度背景中SGL表現能力最佳。在分組為50的不同遺傳力下，四種方法預測能力相近。而在按照KEGG分組的情況下，利用分組信息的方法明顯優(yōu)于傳統模型，SGL預測能力最優(yōu)。在50、200、300分組的情況下，LASSO、Enet和GSSLASSO預測能力相差不大。

圖1 三種模型與SGL的預測能力比較

2.真實數據結果

根據先前的研究，發(fā)病年齡可能是更常見于遺傳起源的重要癌癥指標，本研究使用乳腺癌初始病理診斷時的年齡(即發(fā)病年齡)作為表型。首先利用KEGG通路信息對基因進行分組，接著應用SGL和GSSLASSO兩種方法，如果同一個基因出現在不同的通路中，就把它復制到各通路中；LASSO和Enet方法不進行分組，并刪除重復的基因。

圖2是使用來自TCGA數據集的三種表型對四種模型的預測性能進行比較，預測能力通過相關系數R來衡量；每種方法重復100次。從圖2可以看出在兩個真實數據集中，利用外部信息的SGL和GSSLASSO方法預測精度高于LASSO和Enet。在結直腸癌數據中，SGL方法預測能力最強，其次是GSSLASSO，SGL的預測精度最高，相關系數在0.25左右；其次是GSSLASSO，相關系數在0.22左右；LASSO、Enet預測能力最差，相關系數低于0.2；在胰腺癌數據中SGL方法預測能力最強，其次是GSSLASSO，SGL的相關系數在0.3左右，GSSLASSO的相關系數在0.2左右，LASSO和Enet的相關系數在0.16左右。

圖2 三種表型對四種模型的預測性能進行比較

討 論

本研究比較的四種方法均是稀疏模型，其中LASSO和Enet已被廣泛應用于遺傳學研究，在遺傳統計學的多個方面都有較好的應用價值。目前，國內外關于遺傳風險預測模型的主要研究大多數只考慮了變量的稀疏性，只是把單獨的遺傳位點納入模型，沒有考慮到位點之間可能存在的相關性和復雜結構。在實際情況中，經常會遇到具有分組結構的變量，如果忽視其中的分組結構，統計效能就會大打折扣。大量研究表明無論是SNP數據還是基因表達數據都具有復雜的遺傳結構，與常見人類疾病相關的單個遺傳變異不會直接導致疾病，而是作用于中間的分子表型或與其他遺傳位點共同作用，進而引起高階疾病特征的變化，忽略基因運作的分子網絡和功能結構以及這些網絡和結構變化如何導致疾病特征的變化會使得模型的預測精度不佳[19-21]。

SGL和GSSLASSO方法作為分組稀疏模型，它們考慮到了組間的稀疏性和組內的稀疏性。從模擬結果和真實數據結果我們可以看出，對組內基因進行懲罰可以提高預測精度。一般來說，模型的預測能力取決于數據的真實結構，模型假設與數據結構的吻合度越高，預測準確性越高。從模擬數據的分析結果可以看到各模型在不同遺傳度和不同分組設置下預測能力的表現。本研究共有12種模擬方案，在9種方案中，SGL的預測能力要明顯優(yōu)于其他方法，在其余方案中，四種模型預測能力相近。當分組數設置為200組時，SGL模型明顯優(yōu)于其他方法，當分組數設置為50/300時，幾種模型的預測精度接近，SGL模型要略好；當分組設置為328組，并且數據通過基因進行分組時，SGL模型有較強的預測能力和穩(wěn)定性。真實數據的結果也可以看出，三種整合外部注釋信息的SGL和GSSLASSO方法要優(yōu)于LASSO和Enet方法。

本研究從模擬和實際出發(fā)，探索整合外部注釋信息能否提高對連續(xù)型表型的預測準確性，通過對比分析，整合外部注釋信息可以顯著提高對遺傳表型的預測，并且SGL方法有較高的預測精度。