基于ESPCA模型的結(jié)腸癌差異基因篩選研究*

哈爾濱醫(yī)科大學(xué)衛(wèi)生統(tǒng)計(jì)學(xué)教研室(150081) 李 稱 張秋菊 孫 琳 叢雨欣 趙 敏 劉美娜

【提 要】 目的 驗(yàn)證ESPCA模型降維及變量選擇的效果優(yōu)于SPCA；利用ESPCA模型分析結(jié)腸癌基因數(shù)據(jù)，尋找不同主成分的通路信息，篩選結(jié)腸癌的差異基因。方法 模擬實(shí)驗(yàn)比較ESPCA模型與SPCA的降維及變量選擇能力，通過靈敏度、特異度、準(zhǔn)確度評(píng)價(jià)其變量篩選效果；收集TCGA結(jié)腸癌基因數(shù)據(jù)，利用ESPCA模型對(duì)結(jié)腸癌基因數(shù)據(jù)進(jìn)行降維及變量篩選。結(jié)果 模擬實(shí)驗(yàn)獲得的ESPCA模型具有較高的靈敏度、特異度和準(zhǔn)確度，降維及變量選擇效果優(yōu)于SPCA。結(jié)腸癌數(shù)據(jù)分析結(jié)果顯示：降維后ESPCA模型第一主成分中基因與18條GO-BP通路及7條KEGG通路有關(guān)，ESPCA模型第二主成分中基因與19條GO-BP通路有關(guān)；根據(jù)度中心性及中介中心性對(duì)主成分基因排序篩選出6個(gè)差異基因：MYC、CD44、PKM、FBL、PSMA7、RPS2，經(jīng)GSE137327數(shù)據(jù)集驗(yàn)證具有較高的AUC值。結(jié)論 ESPCA模型在降維過程中既考慮數(shù)據(jù)本身信息，又考慮了生物網(wǎng)絡(luò)信息，具有良好的降維及變量選擇效果。結(jié)腸癌數(shù)據(jù)分析中ESPCA第一主成分中基因參與癌癥相關(guān)的通路；ESPCA第二主成分中基因參與的通路與免疫相關(guān)；通路中基因參與癌癥的發(fā)生發(fā)展過程、可能通過免疫反應(yīng)相關(guān)通路調(diào)節(jié)結(jié)腸癌的發(fā)生發(fā)展過程；獲得的6個(gè)結(jié)腸癌差異基因可為研究疾病發(fā)生機(jī)制和鑒別診斷提供依據(jù)。

隨著組學(xué)技術(shù)發(fā)展，高通量組學(xué)數(shù)據(jù)為尋找生物標(biāo)志物、疾病鑒別診斷提供新機(jī)遇?！案呔S、小樣本”為組學(xué)數(shù)據(jù)的重要特征，成為影響組學(xué)研究結(jié)果準(zhǔn)確性的重要因素。為解決這一問題，許多數(shù)據(jù)降維、變量選擇方法被應(yīng)用于組學(xué)數(shù)據(jù)分析。目前，稀疏主成分分析(sparse principal component analysis，SPCA)廣泛用于高維組學(xué)數(shù)據(jù)降維，如LASSO正則化[1]、加權(quán)稀疏主成分分析[2]。SPCA模型在基因篩選中，假定所有基因有同等概率被選至各主成分中，沒有考慮變量間的權(quán)重、忽略基因在網(wǎng)絡(luò)中的相互作用。邊組稀疏主成分分析[3](edge-group sparse PCA，ESPCA)以基因調(diào)控網(wǎng)絡(luò)中相互作用關(guān)系為先驗(yàn)信息對(duì)高維組學(xué)數(shù)據(jù)進(jìn)行分析；基因調(diào)控網(wǎng)絡(luò)[4]可以反應(yīng)出基因間相互作用關(guān)系，將調(diào)控網(wǎng)絡(luò)引入到稀疏主成分分析中，獲得的主成分更加接近真實(shí)的生物模式，但該方法尚未應(yīng)用在癌癥基因數(shù)據(jù)的變量篩選中。

目前，結(jié)腸癌的發(fā)病率逐年遞增，基因治療[5]正被廣泛用于腫瘤等疾病的治療中，但其潛在的基因通路及結(jié)腸癌差異基因仍不明確。本研究通過模擬比較ESPCA模型與SPCA在降維、變量選擇的效果，并利用ESPCA模型對(duì)結(jié)腸癌基因數(shù)據(jù)集進(jìn)行降維及變量選擇，獲得結(jié)腸癌相關(guān)通路及差異基因，為研究結(jié)腸癌的相關(guān)機(jī)制及鑒別診斷與治療提供依據(jù)。

分析方法

1.ESPCA 模型

ESPCA模型[3]是一種以網(wǎng)絡(luò)交互作用為導(dǎo)向處理高維組學(xué)數(shù)據(jù)的方法，將稀疏主成分與先驗(yàn)網(wǎng)絡(luò)相結(jié)合，通過選擇具有交互作用的組結(jié)構(gòu)實(shí)現(xiàn)變量選擇，算法步驟為：

(1)假設(shè)g為組結(jié)構(gòu)，g={e1，…，eM}表示基因調(diào)控網(wǎng)絡(luò)中所有邊組的集合。當(dāng)g為非重疊的組結(jié)構(gòu)時(shí)，組稀疏懲罰為L(zhǎng)1范數(shù)或L0范數(shù)。然而實(shí)際應(yīng)用中需要考慮重疊組結(jié)構(gòu)[6]的情況，基因相互作用網(wǎng)絡(luò)中兩個(gè)相連基因被視為一個(gè)邊組。顯然，邊組為重疊組結(jié)構(gòu)，因此邊組稀疏懲罰(ES懲罰)如下：

‖u‖ES=minimize|g′|

?g′?g，support(u)?V(g′)

(1)

其中，g′為g的子集，V(g′)為g′的頂點(diǎn)集，|g′|表示g′的元素個(gè)數(shù)，support(u)為u中非零元素指數(shù)集。

(2)由于ES懲罰的應(yīng)用，ESPCA載荷非零元素的選擇基于g中重要基因相互作用的邊。最終形成帶有ES懲罰的主成分分析即ESPCA，其表達(dá)式如下：

(2)

其中，X為基因表達(dá)矩陣，u為主成分載荷，v為主成分，k為選擇邊的個(gè)數(shù)，‖·‖2為歐氏距離。

(3)模型參數(shù)估計(jì)的核心問題在于用固定的v、z=Xv，求解k邊組稀疏投影問題。該問題為NP-難問題，將通過貪婪算法[7]求解下式的最優(yōu)解：

(3)

(4)

其中，g(i)為包含基因i的邊組指數(shù)集；supp(z，k)表示z的絕對(duì)值中較大的k個(gè)元素的指數(shù)集。

(4)采用交替迭代策略直到收斂得到最優(yōu)解：

(5)應(yīng)用Gram-Schmidt正交化方法，分別求出ESPCA多個(gè)主成分和主成分載荷。

(6)

2.模型參數(shù)選擇

邊個(gè)數(shù)k和主成分?jǐn)?shù)l是ESPCA模型最重要的參數(shù)，這兩個(gè)參數(shù)的確定是開放性的問題，常以最小化貝葉斯信息準(zhǔn)則(bayes′ information criterion，BIC)實(shí)現(xiàn)模型參數(shù)估計(jì)為可行的策略之一。根據(jù)已有文獻(xiàn)[3]，本研究將k=150，l=2設(shè)為默認(rèn)值。

模擬實(shí)驗(yàn)

1.模擬參數(shù)設(shè)置

(1)生成兩個(gè)主成分載荷u1、u2，長(zhǎng)度為200的列向量，u1中前100個(gè)變量的載荷隨機(jī)產(chǎn)生，后100個(gè)變量載荷均為0；u2反之；

(2)生成兩個(gè)主成分v1、v2，表示隨機(jī)產(chǎn)生的長(zhǎng)度為100的列向量；

(4)相互作用網(wǎng)絡(luò)g由g1，g2，g3組成；g1，g2分別為前100和后100個(gè)變量以概率p=0.3組成的網(wǎng)絡(luò)；g3為前、后100個(gè)變量之間以p=0.04組成的網(wǎng)絡(luò)；

(5)k設(shè)定為300，為了便于比較，提取前兩個(gè)主成分，且SPCA每個(gè)主成分提取變量個(gè)數(shù)與ESPCA的相同。

2.模擬評(píng)價(jià)指標(biāo)

(1)靈敏度：實(shí)際為陽(yáng)性(系數(shù)非零)的變量中，判斷為陽(yáng)性的比例。

(2)特異度：實(shí)際為陰性(系數(shù)為零)的變量中，判斷為陰性的比例。

(3)準(zhǔn)確度：正確判斷為陽(yáng)性和陰性變量的比例。

3.模擬實(shí)驗(yàn)結(jié)果

模擬實(shí)驗(yàn)獲得兩個(gè)主成分，每一個(gè)主成分在不同的噪聲水平下，ESPCA靈敏度、特異度、準(zhǔn)確度的效果優(yōu)于SPCA；隨著噪聲水平γ的增大，靈敏度、特異度、準(zhǔn)確度呈下降趨勢(shì)；由于PC1大于PC2的權(quán)重值，故PC1有更高的靈敏度、特異度、準(zhǔn)確度，見圖1。

圖1 在不同噪聲水平下SPCA、ESPCA的靈敏度、特異度、準(zhǔn)確度

結(jié)腸癌差異基因篩選

下載TCGA(TCGA，http：//cancergenome.nih.gov/)結(jié)腸癌RNA-seq表達(dá)數(shù)據(jù)：包括60488個(gè)基因，512例樣本(471個(gè)為癌組織，41個(gè)為癌旁組織)。

1.數(shù)據(jù)預(yù)處理

基因篩選：(1)根據(jù)“80原則”剔除缺失變量；(2)相同探針保留變異較大的基因；(3)剔除t檢驗(yàn)P>0.05的變量，最終剩余7837個(gè)基因。邊的獲得：利用7837個(gè)基因在STRING網(wǎng)站中獲得基因網(wǎng)絡(luò)相互作用的邊，共87493條。

2.ESPCA模型降維及富集分析

ESPCA模型應(yīng)用于預(yù)處理后結(jié)腸癌數(shù)據(jù)集，提取前兩個(gè)主成分：ESPCA1和ESPCA2，繪制結(jié)腸癌得分圖，結(jié)果見圖2，圖中兩個(gè)主成分可以很好地區(qū)分癌組織與癌旁組織。

圖2 基于ESPCA前兩個(gè)主成分結(jié)腸癌的得分圖

ESPCA1中有95個(gè)基因變量，ESPCA2中有92個(gè)基因變量，將兩個(gè)主成分中的基因分別進(jìn)行富集分析發(fā)現(xiàn)：ESPCA1中基因與18條GO-BP通路和7條KEGG通路有關(guān)，其中細(xì)胞分裂、PI3K-Akt信號(hào)、EMC-受體交互、癌癥中心的碳代謝等通路與癌癥的發(fā)生發(fā)展相關(guān)，如 PKM，SLC1A5，SLC2A1，MYC等基因參與KEGG中的癌癥中心碳代謝通路，結(jié)果見表1；ESPCA2的基因與19條GO-BP通路有關(guān)，富集分析發(fā)現(xiàn)通路與免疫過程相關(guān)，如：補(bǔ)體激活、免疫反應(yīng)調(diào)節(jié)、B細(xì)胞受體信號(hào)通路、免疫應(yīng)答等，結(jié)果見表2。

表1 ESPCA第一主成分富集分析結(jié)果

表2 ESPCA第二主成分富集分析結(jié)果

3.篩選結(jié)腸癌的差異基因

主成分中的基因做PPI網(wǎng)絡(luò)圖(圖3)，根據(jù)度中心性和中介中心性排序，選擇前6個(gè)基因，分別為MYC、CD44、PKM、FBL、PSMA7、RPS2，其AUC值分別為0.97、0.99、0.98、0.96、0.97、0.96，見表3。

圖3 基因相互作用網(wǎng)絡(luò)

表3 6個(gè)基因排序及AUC值

下載GEO結(jié)腸癌數(shù)據(jù)集GSE137327共18例樣本，其中9個(gè)為對(duì)照組，9個(gè)為病例組進(jìn)行外部驗(yàn)證，6個(gè)基因的AUC值分別為0.72、0.69、0.74、0.77、0.68、0.74，結(jié)果見圖4。

圖4 篩選6個(gè)基因的ROC曲線

討 論

本研究模擬實(shí)驗(yàn)的驗(yàn)證結(jié)果顯示，在不同的噪聲水平下，ESPCA的靈敏度、特異度、準(zhǔn)確度均優(yōu)于SPCA，在邊的信息下提取的樣本信息更加接近真實(shí)模式。ESPCA模型不僅能用于基因組學(xué)數(shù)據(jù)，也可用于對(duì)蛋白質(zhì)組學(xué)、代謝組學(xué)數(shù)據(jù)的降維及變量選擇。

ESPCA模型對(duì)結(jié)腸癌數(shù)據(jù)進(jìn)行降維分析，獲得兩個(gè)主成分ESPCA1和ESPCA2，基于兩個(gè)主成分可以很好地區(qū)分癌組織和癌旁組織。富集分析發(fā)現(xiàn)ESPCA1中的基因參與PI3K-Akt信號(hào)、EMC-受體交互、癌癥中心的碳代謝等與癌癥相關(guān)的通路。研究表明：PI3K-Akt信號(hào)通路[8]是癌癥中影響生存進(jìn)展的主要通路，Akt[9]的失調(diào)會(huì)導(dǎo)致癌癥、糖尿病、心血管和神經(jīng)系統(tǒng)等疾??；EMC[10]對(duì)細(xì)胞粘附、增殖、凋亡起重要作用，影響腫瘤的進(jìn)展過程；與正常細(xì)胞相比，癌癥細(xì)胞的中心代謝通路[11]存在明顯差異，癌癥細(xì)胞即使在正常的氧氣濃度下也能將大部分葡萄糖轉(zhuǎn)換為乳酸，這種差異影響癌癥的發(fā)生發(fā)展。富集分析發(fā)現(xiàn)ESPCA2中的基因參與補(bǔ)體激活、免疫反應(yīng)調(diào)節(jié)、B細(xì)胞受體信號(hào)、免疫應(yīng)答等免疫相關(guān)通路。有研究發(fā)現(xiàn)，免疫浸潤(rùn)有關(guān)的mRNA與I～I(xiàn)II期的結(jié)腸癌的診斷及預(yù)后有關(guān)；腫瘤浸潤(rùn)淋巴細(xì)胞是三陰性乳腺癌[12]重要的預(yù)后因素。目前尚未發(fā)現(xiàn)免疫通路影響結(jié)腸癌的發(fā)生，提示基因通過調(diào)控免疫反應(yīng)通路可能影響結(jié)腸癌的發(fā)生發(fā)展，免疫反應(yīng)機(jī)制可能是發(fā)生癌癥的潛在機(jī)制之一。

篩選出的6個(gè)基因在結(jié)腸癌數(shù)據(jù)集及GSE137327驗(yàn)證集中有較高的AUC值。其中，PKM、MYC、CD44和PSMA7已被證實(shí)與結(jié)腸癌有關(guān)，PKM[13]基因通過調(diào)控STAT3相關(guān)的信號(hào)通路促進(jìn)結(jié)腸癌細(xì)胞的黏附和遷移；MYC[14]基因?yàn)榻Y(jié)腸癌的關(guān)鍵基因，并作為L(zhǎng)EF1的轉(zhuǎn)錄調(diào)節(jié)因子，通過激活LEF1的表達(dá)來調(diào)節(jié)結(jié)腸癌細(xì)胞的增殖過程；CD44[15]基因已被多個(gè)研究證實(shí)是結(jié)腸癌和胃癌干細(xì)胞的生物標(biāo)志物；PSMA7[16]基因在結(jié)腸癌組織中過表達(dá)，能誘導(dǎo)HT-29細(xì)胞的凋亡。目前尚未發(fā)現(xiàn)FBL、RPS2與結(jié)腸癌發(fā)生發(fā)展有關(guān)，本研究為進(jìn)一步深入了解結(jié)腸癌發(fā)生發(fā)展過程及基因治療提供了參考和依據(jù)。

本文雖然擴(kuò)展了現(xiàn)有文獻(xiàn)的參數(shù)設(shè)置，但仍需探索更多的參數(shù)組合下模型的效果；GSE137327驗(yàn)證集中癌組織和癌旁組織的比例與TCGA結(jié)腸癌數(shù)據(jù)集中的比例相差較大，結(jié)果可能會(huì)有一定的偏差；只考慮了網(wǎng)絡(luò)中兩個(gè)相互作用的基因作為一個(gè)組結(jié)構(gòu)，未來可以考慮將通路中具有類似功能的多個(gè)基因作為一個(gè)組結(jié)構(gòu)進(jìn)行ESPCA模型分析。