基于GEO和TCGA數(shù)據(jù)庫食管鱗狀細(xì)胞癌預(yù)后模型構(gòu)建與驗(yàn)證

管沛文， 李 雪， 尤崇革

蘭州大學(xué)第二醫(yī)院檢驗(yàn)醫(yī)學(xué)中心，甘肅 蘭州 730000

食管癌(esophageal carcinoma，ESCA)是世界上第七大常見癌癥，死亡率居世界第六位[1]。食管癌包括兩種主要病理類型，其中食管鱗狀細(xì)胞癌(esophageal squamous cell carcinoma，ESCC)占食管癌病例的90%。ESCC在發(fā)展中國家發(fā)病率較高，且男性發(fā)病率為女性的2～4倍[2-3]。由于ESCA的5年總體生存率低于13%[4-6]，因此，尋找新的預(yù)后標(biāo)志物至關(guān)重要。本研究將基因表達(dá)數(shù)據(jù)庫(Gene Expression Omnibus，GEO)中獲得的差異表達(dá)基因整合到癌癥基因組圖譜(The Cancer Genome Atlas，TCGA)數(shù)據(jù)庫中，對ESCC進(jìn)行數(shù)據(jù)挖掘分析。我們在差異基因的基礎(chǔ)上進(jìn)行了功能富集分析、蛋白質(zhì)-蛋白質(zhì)相互作用分析、生存分析、Cox回歸分析、基因集富集分析等，并發(fā)現(xiàn)了全新的生物標(biāo)志物，為ESCC分子機(jī)制提供新的見解。

1 資料與方法

1.1 ESCC數(shù)據(jù)收集與處理ESCC相關(guān)基因組表達(dá)數(shù)據(jù)集從GEO數(shù)據(jù)庫(https://www.ncbi.nlm.nih.gov/)中下載。TCGA-ESCA中ESCC的轉(zhuǎn)錄組數(shù)據(jù)和臨床數(shù)據(jù)從Genomic Data Commons(GDC)(https://portal.gdc.cancer.gov/)中下載。

鑒于TCGA中只有一個(gè)正常樣本，我們基于GEO的數(shù)據(jù)確定差異表達(dá)基因。為了監(jiān)測在不相關(guān)輸入的零假設(shè)下更好地挖掘差異基因，我們選擇“Robust Rank Aggreg”包處理數(shù)據(jù)集[7]。另外，R中的軟件包“normalizeBetweenArrays”和“l(fā)imma”(版本3.6.2，http：//www.r-project.org/)用于標(biāo)準(zhǔn)化矩陣文件中樣本并篩選差異基因；“impute”包用于處理缺失數(shù)據(jù)。最終，“Robust Rank Aggreg”包用以整合5個(gè)數(shù)據(jù)集中的共表達(dá)差異基因。

1.2 KEGG信號通路富集和GO功能注釋分析為了進(jìn)一步探索差異基因，DAVID(https://david.ncifcrf.gov/)(版本6.8)用于進(jìn)行“京都基因與基因組百科全書”(Kyoto Encyclopedia of Genes and Genomes，KEGG)和“基因本體論”(Gene Ontology，GO)分析；STRING(https://string-db.org/)在線數(shù)據(jù)庫和Cytoscape軟件用于篩選前15個(gè)關(guān)鍵基因。

1.3 差異基因的預(yù)后分析我們使用在線工具Kaplan-Meier plotter(https://kmplot.com/analysis/)探索差異基因在ESCC中總生存期(overall survival，OS)的差異，并通過單變量Cox回歸分析進(jìn)一步驗(yàn)證預(yù)后價(jià)值，最終篩選出有預(yù)后價(jià)值的關(guān)鍵基因。

1.4 實(shí)時(shí)熒光定量PCR(qRT-PCR)ESCC組織和正常對照(各6例)取自蘭州大學(xué)第二醫(yī)院萃英中心樣本庫。使用TRIzol試劑提取總RNA，然后用FastKing gDNA DispellingRTSuperMix(天根，中國北京)進(jìn)行逆轉(zhuǎn)錄。所有qRT-PCR反應(yīng)均使用羅氏LightCycler96 PCR儀進(jìn)行?；虻南鄬Ρ磉_(dá)通過2-ΔΔCt計(jì)算。引物和Tm值如表1所示。

表1 引物序列

1.5 預(yù)后模型的構(gòu)建與驗(yàn)證多元Cox回歸算法用于開發(fā)預(yù)后模型。風(fēng)險(xiǎn)評分計(jì)算如下：

(其中Coefi為系數(shù)，xi為基因表達(dá)量)。并將TCGA數(shù)據(jù)分為建模組和驗(yàn)證組，分別確定模型的臨床應(yīng)用價(jià)值。此外，R軟件中“survival”包用于展示兩個(gè)風(fēng)險(xiǎn)亞組OS的差異。最終，ROC曲線、靈敏度、特異度、一致性分析中C-統(tǒng)計(jì)量以及基于患者死亡結(jié)局的靈敏度和特異度用于確定模型的預(yù)測能力。

1.6 基因集富集分析基因集富集分析用于研究Cox回歸中兩個(gè)風(fēng)險(xiǎn)亞群的相關(guān)通路。我們選擇C2(C2.cp.kegg.v7.1.symbols.gmt)作為參考基因集進(jìn)行基因集富集分析。|NES|>1.5、adj.P<0.05和FDR<0.25為差異有統(tǒng)計(jì)學(xué)意義。

1.7 統(tǒng)計(jì)學(xué)分析所有統(tǒng)計(jì)分析均使用R軟件(版本3.6.2)進(jìn)行。單變量Cox分析用以篩選預(yù)后相關(guān)的關(guān)鍵基因。中位風(fēng)險(xiǎn)評分用作臨界值，將患者分為高風(fēng)險(xiǎn)組和低風(fēng)險(xiǎn)組，Kaplan-Meier方法用于分析高風(fēng)險(xiǎn)與低風(fēng)險(xiǎn)之間的OS差異組。通過卡方檢驗(yàn)比較不同亞組的臨床特征。采用多元Cox回歸分析建立ESCC的預(yù)后模型。

2 結(jié)果

2.1 數(shù)據(jù)的下載與處理鑒于所需樣本為ESCC的組織，微陣列數(shù)據(jù)集GSE20347[8]、GSE23400[9]、GSE26886[10]、GSE45168和GSE77861[11]被納入研究中(見表2)。

表2 ESCC相關(guān)基因組表達(dá)數(shù)據(jù)集

每個(gè)芯片的差異基因如圖1A所示。ESCC組織中共有42個(gè)高表達(dá)基因和56個(gè)低表達(dá)基因(adj.P<0.05和logFC≥1)，前40個(gè)差異基因表達(dá)情況如圖1B所示。

圖1 差異表達(dá)基因 A：火山圖；B：熱圖展示了RRA方法從5個(gè)獨(dú)立的基因芯片中選擇的20個(gè)高表達(dá)基因和20個(gè)低表達(dá)基因，方框中標(biāo)記基因的FC值

2.2 GO和KEGG富集分析結(jié)果GO和KEGG富集分析結(jié)果如表3～4所示，GO富集分析中，上調(diào)基因主要富集在細(xì)胞外區(qū)域、細(xì)胞外空間等；下調(diào)基因則主要集中在細(xì)胞外外泌體、絲氨酸型肽酶活性等。KEGG富集分析中，上調(diào)基因主要富集于細(xì)胞外基質(zhì)-受體結(jié)合等，而下調(diào)基因中無有意義的KEGG通路。最終篩選出GABRP、MAGEA11、SERPINH1、WDR72、INHBA、IFI6、TGFBI、CEACAM6、MAGEA6、HMGA2、VCAN、ACPP、CRABP2、KIAA1429和ALOX12。

表3 上調(diào)基因GO和KEGG富集分析主要結(jié)果展示

表4 下調(diào)基因GO和KEGG富集分析主要結(jié)果展示

2.3 關(guān)鍵基因的生存分析Kaplan-Meier Plotter結(jié)果表明，IFI6、APOBEA3B、TGFBI、CEACAM6、ACPP和GABRP的高表達(dá)表明ESCC患者的預(yù)后較差(P<0.05)，而INHBA、WDR72、HMGA2、MAGEA11、MAGEA6、VCAN、SERPINH1、CRABP2和ALOX12高表達(dá)則表明ESCC患者預(yù)后良好(P<0.05)。另外，通過單變量Cox回歸再次驗(yàn)證食管癌中上述基因的預(yù)后價(jià)值發(fā)現(xiàn)，ACPP為危險(xiǎn)基因(HR>1)，而MAGEA11、MAGEA6、HMGA2、VCAN和ALOX12是保護(hù)基因(HR<1)(見表5)。如上所述，我們選擇了MAGEA11、MAGEA6、HMGA2、VCAN、ALOX12和ACPP作為我們的研究重點(diǎn)。

表5 關(guān)鍵基因的單變量Cox回歸結(jié)果

2.4 qRT-PCR為了進(jìn)一步驗(yàn)證數(shù)據(jù)挖掘結(jié)果，本研究通過qRT-PCR驗(yàn)證了基因在ESCC組織與正常組織中的表達(dá)水平。雖然樣本量有限，但除MAGEA11外，其余5個(gè)候選基因均在腫瘤組織中呈顯著高表達(dá)，差異均有統(tǒng)計(jì)學(xué)意義(P<0.05)(見圖2)。

圖2 qRT-PCR驗(yàn)證預(yù)后相關(guān)基因在ESCC組織及正常組織中的差異表達(dá)

2.5 預(yù)后模型的建立與驗(yàn)證多元Cox回歸分析上述所選基因的預(yù)后價(jià)值并構(gòu)建模型，風(fēng)險(xiǎn)評分公式如下：Risk Score=-0.0479×(VCAN)-0.1605×(ALOX12)+0.1318×(ACPP)。此外，基于“高風(fēng)險(xiǎn)”和“低風(fēng)險(xiǎn)”的生存曲線(見圖3A、4A)，我們發(fā)現(xiàn)建模組和驗(yàn)證組中，高風(fēng)險(xiǎn)亞組的患者生存率較低風(fēng)險(xiǎn)亞組的患者低(P<0.01)。

如熱圖所示(見圖3C、4C)，ACPP在高風(fēng)險(xiǎn)組中水平較高，VCAN、ALOX12則在低風(fēng)險(xiǎn)組中水平較高。建模組ROC曲線下方面積為0.873，當(dāng)風(fēng)險(xiǎn)值為0.502時(shí)，預(yù)測患者預(yù)后不良靈敏度為84.6%，特異度為82.4%(見圖3B)；驗(yàn)證組ROC曲線下方面積為0.897，當(dāng)風(fēng)險(xiǎn)值為0.644時(shí)，預(yù)測患者不良的靈敏度為92.3%，特異度為82.4%(見圖4B)。且建模組中該模型預(yù)測死亡結(jié)局的靈敏度為70.59%，特異度為65.21%；驗(yàn)證組中靈敏度為94.12%，特異度為74.07%。另外，列線圖和一致性分析中C-統(tǒng)計(jì)量(實(shí)驗(yàn)組為0.804，驗(yàn)證組為0.834)表明該模型預(yù)測效能良好(見圖3D、4D)。

圖3 預(yù)后模型在實(shí)驗(yàn)組中的效能 A：兩亞組的生存曲線；B：1年期ROC曲線；C：3基因表達(dá)量熱圖；D：預(yù)測TCGA隊(duì)列1、3和5年OS的列線圖

圖4 預(yù)后模型在驗(yàn)證組中的效能 A：兩亞組的生存曲線；B：1年期ROC曲線；C：3基因表達(dá)量熱圖；D：預(yù)測TCGA隊(duì)列1、3和5年OS的列線圖

2.6 基因集富集分析為了深入了解所選調(diào)控因子對癌癥預(yù)后的影響機(jī)制，本研究基于“高風(fēng)險(xiǎn)”和“低風(fēng)險(xiǎn)”亞組進(jìn)行基因集富集分析。富集結(jié)果如圖5所示。高風(fēng)險(xiǎn)亞組主要與ABC轉(zhuǎn)運(yùn)蛋白、丁酸甲酯代謝等通路相關(guān)；而低風(fēng)險(xiǎn)亞組主要與基底轉(zhuǎn)錄因子等相關(guān)。

注：A：高風(fēng)險(xiǎn)亞組的基因集富集分析結(jié)果；B:低風(fēng)險(xiǎn)亞組的基因集富集分析結(jié)果。

3 討論

ESCC是全球侵襲性最強(qiáng)的癌癥之一，也是中國第二大最常見的癌癥[12]。由于ESCC初次診斷時(shí)常伴隨癌細(xì)胞轉(zhuǎn)移，多模式療法也未能提高ESCC的存活率[4，13]。因此，ESCC預(yù)后相關(guān)標(biāo)志物的研究至關(guān)重要。本研究依據(jù)公開數(shù)據(jù)庫集中確定了ESCC和正常食管樣品之間的98個(gè)DEG。并且，通過Kaplan-Meier Plotter分析和Cox回歸分析發(fā)現(xiàn),MAGEA11、MAGEA6、HMGA2、VCAN、ALOX12和ACPP對ESCC的預(yù)后有統(tǒng)計(jì)學(xué)意義。癌癥/睪丸抗原家族中的黑色素瘤抗原-A6(MAGEA6)通過操縱自噬和由AMP激活的蛋白激酶信號傳導(dǎo)來調(diào)節(jié)癌癥的進(jìn)展[14-15]。本研究中，MAGEA6和MAGEA11在ESCC組織中高表達(dá)且對預(yù)后具有重要意義，這與其他研究結(jié)果一致[16]。HMGA2屬于非組蛋白染色體高遷移率族蛋白家族，可促進(jìn)各種哺乳動(dòng)物啟動(dòng)子的形成[17]。研究表明，HMGA2通過與microRNA[18]和非編碼RNA[19]相互作用來調(diào)控腫瘤的發(fā)生和癌變過程。

ALOX12來自哺乳動(dòng)物的脂氧合酶家族[20]，ALOX12及其代謝物12S-羥基二十碳四烯酸影響腫瘤的轉(zhuǎn)化和進(jìn)展，有助于克服癌癥的化學(xué)耐藥性[21]。VCAN是硫酸軟骨素基質(zhì)蛋白聚糖[22]，它通過TGF-β[23-24]、TNF信號通路[25]和基質(zhì)金屬蛋白酶[23，25]來調(diào)節(jié)癌癥的進(jìn)展。在我們的研究中，ACPP被發(fā)現(xiàn)是ESCC預(yù)后的危險(xiǎn)基因，鑒于有研究發(fā)現(xiàn)ACPP在口腔鱗狀細(xì)胞癌可促進(jìn)癌細(xì)胞轉(zhuǎn)移[26]，提示ACPP可能通過促進(jìn)ESCC癌細(xì)胞的轉(zhuǎn)移而影響預(yù)后。基于上述三種基因建立的預(yù)后模型經(jīng)ROC曲線下面積與列線圖的C統(tǒng)計(jì)量驗(yàn)證，模型預(yù)測能力較好(AUC>0.8且C-index>0.8)。

GSEA結(jié)果表明，高風(fēng)險(xiǎn)組主要富集在藥物代謝細(xì)胞色素P450，視黃醇代謝和類固醇激素生物合成途徑。研究發(fā)現(xiàn)[27]，參與視黃醇代謝的細(xì)胞結(jié)合蛋白改變有助于促進(jìn)乳腺癌的易感性，并與前列腺癌、卵巢腺癌和成膠質(zhì)細(xì)胞瘤的不良預(yù)后相關(guān)。另外，固醇激素的合成途徑[28]和藥物代謝細(xì)胞色素P450[29]已被證實(shí)與多種癌癥有關(guān)。如今，尚未有文獻(xiàn)報(bào)道以上途徑對ESCC預(yù)后的影響。綜上所述，本研究開發(fā)的預(yù)測模型可有效地預(yù)測ESCC患者的預(yù)后，但相關(guān)的作用機(jī)制與臨床適用性仍需進(jìn)一步研究。