基于TCGA數(shù)據(jù)庫(kù)膽管癌自噬-臨床預(yù)后模型的構(gòu)建和評(píng)估

鄒文強(qiáng)，符廣華，楊艷梅，張新宇

（1.哈爾濱醫(yī)科大學(xué)附屬第二醫(yī)院 普通外科，2.哈爾濱醫(yī)科大學(xué) 腫瘤防治研究所，黑龍江 哈爾濱 150086）

膽管癌（cholangiocarcinoma，CCA）是起源于膽管上皮細(xì)胞的惡性腫瘤，根據(jù)病變的解剖部位可以分為肝內(nèi)CCA、肝周CCA和遠(yuǎn)端CCA[1]。CCA的發(fā)病率約占消化道腫瘤的3%，我國(guó)及東南亞地區(qū)是最高發(fā)的地區(qū)，而且近幾十年發(fā)病率呈上升趨勢(shì)[2]。由于CCA患者早期發(fā)病癥狀隱匿且惡性程度較高，一旦發(fā)現(xiàn)，多數(shù)腫瘤已進(jìn)入晚期。研究表明根治性手術(shù)是最有效的治療方式，但由于淋巴結(jié)轉(zhuǎn)移早，易侵犯血管及神經(jīng)，術(shù)后遠(yuǎn)期預(yù)后較差[3]。CCA患者的預(yù)后評(píng)估方法主要是TNM分期和組織學(xué)診斷，但無(wú)法實(shí)現(xiàn)準(zhǔn)確的個(gè)體化生存預(yù)測(cè)[4-6]；且臨床常規(guī)的腫瘤標(biāo)志物往往敏感度和特異性不高，因此，尋找可靠的預(yù)后標(biāo)志物來(lái)優(yōu)化預(yù)后評(píng)估系統(tǒng)是必不可少的。

近年隨著對(duì)腫瘤的深入研究，發(fā)現(xiàn)自噬在腫瘤的發(fā)生發(fā)展中發(fā)揮著重要作用[7-9]。自噬是一種通過(guò)溶酶體自我分解細(xì)胞內(nèi)的老化及受損的蛋白質(zhì)和細(xì)胞器，實(shí)現(xiàn)細(xì)胞內(nèi)穩(wěn)態(tài)的途徑[10]。細(xì)胞自噬在正常生理狀態(tài)下可以促進(jìn)機(jī)體抵御疾病，但當(dāng)自噬異常時(shí)又會(huì)導(dǎo)致腫瘤、自身免疫性及神經(jīng)退行性等疾病[11]，在腫瘤中的作用取決于組織類型、腫瘤分期和致癌突變類型等因素[12]。自噬對(duì)包括CCA在內(nèi)的多種腫瘤，具有雙重調(diào)節(jié)作用。自噬通過(guò)消除受損的線粒體和減少ROS介導(dǎo)的染色體損傷，在腫瘤初期發(fā)揮抑癌功能。自噬通過(guò)抑制p53 的誘導(dǎo)和維持線粒體的代謝功能，在腫瘤進(jìn)展期促進(jìn)腫瘤細(xì)胞存活[13-15]。Sasaki等[16]的研究表明，CCA組織樣本中自噬相關(guān)蛋白（LC3、beclin-1）表達(dá)量顯著增加。CCA進(jìn)展及轉(zhuǎn)移的重要驅(qū)動(dòng)因素之一是上皮細(xì)胞向間質(zhì)細(xì)胞的轉(zhuǎn)變（EMT），Chen等[17]發(fā)現(xiàn)調(diào)節(jié)自噬的信號(hào)通路（PI3K/AKT、p53和beclin1）同時(shí)對(duì)EMT有顯著調(diào)節(jié)的作用，促使CCA細(xì)胞存活。這些研究證實(shí)自噬和CCA存在密切關(guān)系，并表明ATGs作為CCA預(yù)后標(biāo)志物的巨大潛力。既往研究表明依據(jù)自噬基因構(gòu)建的預(yù)后風(fēng)險(xiǎn)模型，可以有效評(píng)估膀胱癌、乳腺癌及胃癌等患者的預(yù)后情況[18-20]。

本研究擬通過(guò)癌癥基因組圖譜數(shù)據(jù)庫(kù)（The Cancer Genome Atlas，TCGA）分析與CCA預(yù)后密切相關(guān)的ATGs，同時(shí)構(gòu)建自噬基因預(yù)后風(fēng)險(xiǎn)模型，觀察其對(duì)CCA患者的預(yù)后價(jià)值，從而為臨床診斷和指導(dǎo)治療CCA提供個(gè)性化的模型。

1 材料和方法

1.1 數(shù)據(jù)下載及整理

在TCGA中獲取CCA的轉(zhuǎn)錄組數(shù)據(jù)和臨床數(shù)據(jù)，包括35 例CCA組織樣本和9 例正常組織樣本。使用Strawberry Perl軟件將各個(gè)獨(dú)立樣本的基因表達(dá)量匯總在一個(gè)文件中，并將基因的“ID”轉(zhuǎn)換為“Symbol”。自噬相關(guān)基因（ATGs）從人類自噬數(shù)據(jù)庫(kù)（Human Autophagy Database，HADb）中下載，包括232個(gè)ATGs的信息。在TCGA數(shù)據(jù)庫(kù)已經(jīng)下載的CCA樣本信息中，利用Perl腳本將232個(gè)ATGs的基因表達(dá)量提取出來(lái)。

1.2 差異表達(dá)的自噬相關(guān)基因（DEATGs）的識(shí)別

通過(guò)R語(yǔ)言（R 4.0.2）使用“Willcox.test”來(lái)篩選DEATGs，篩選條件為：FDR＜0.05，|logFC|＞1。通過(guò)R語(yǔ)言繪制DEATGs的火山圖和箱線圖。

1.3 DEATGs的富集分析

基因本體功能（gene ontology，GO）富集分析和京都基因與基因組百科全書通路（Kyoto Encyclopedia of Genes and Genomes pathway，KEGG pathway）分析，可以明確DEATGs的主要生物學(xué)屬性。調(diào)用R語(yǔ)言中“org.Hs.eg.db”“enrichplot”和“gplot2”包，進(jìn)行GO富集分析和KEGG通路分析；使用R語(yǔ)言中“GOplot”包，繪制GO氣泡圖和KEGG圈圖。

1.4 預(yù)后相關(guān)ATGs的識(shí)別和基于ATGs的預(yù)后風(fēng)險(xiǎn)評(píng)分模型構(gòu)建

使用Strawberry Perl軟件將DEATGs表達(dá)量和生存狀態(tài)及時(shí)間，匯總在一個(gè)文件中。調(diào)用R軟件的“survival”包，篩選標(biāo)準(zhǔn)為P＜0.05，對(duì)DEATGs進(jìn)行單因素Cox回歸分析，篩選出與CCA預(yù)后相關(guān)的自噬基因（ATGs）。再進(jìn)行多因素Cox回歸分析，刪除相關(guān)性較高的基因來(lái)優(yōu)化模型，最后得到用于構(gòu)建預(yù)后風(fēng)險(xiǎn)評(píng)分模型的風(fēng)險(xiǎn)基因。使用R軟件獲取風(fēng)險(xiǎn)基因的評(píng)分系數(shù)（Coef）。模型的計(jì)算公式為：患者風(fēng)險(xiǎn)評(píng)分=（基因1的表達(dá)量×coef1+基因2的表達(dá)量×coef2+…）。按照患者風(fēng)險(xiǎn)評(píng)分，以風(fēng)險(xiǎn)值中位數(shù)為界限，將患者分為高風(fēng)險(xiǎn)和低風(fēng)險(xiǎn)組。

1.5 預(yù)后模型的驗(yàn)證

使用R 4.0.2軟件進(jìn)行生存分析、獨(dú)立預(yù)后分析和臨床相關(guān)性分析，繪制風(fēng)險(xiǎn)曲線和多指標(biāo)ROC曲線。

2 結(jié)果

2.1 自噬相關(guān)基因的差異分析

從TCGA數(shù)據(jù)庫(kù)中下載了35 個(gè)腫瘤組織樣本和9 個(gè)正常組織樣本的轉(zhuǎn)錄組數(shù)據(jù)和臨床數(shù)據(jù)。在HADb數(shù)據(jù)庫(kù)中提取了232 個(gè)自噬相關(guān)基因，并提取自噬相關(guān)基因表達(dá)量。按照FDR＜0.05和log2FC＞1 的篩選標(biāo)準(zhǔn)，過(guò)濾得到49 個(gè)自噬相關(guān)差異基因（DEATGs），其中上調(diào)的基因有48個(gè)，下調(diào)的基因有1個(gè)（圖1）。

圖1 49個(gè)差異自噬基因的箱線圖

2.2 差異自噬基因的富集分析

對(duì)49個(gè)差異表達(dá)的自噬相關(guān)基因（DEATGs）進(jìn)行GO富集分析和KEGG通路分析。經(jīng)過(guò)GO富集分析，在生物學(xué)過(guò)程方面，基因主要富集在自噬、利用自噬機(jī)制的過(guò)程、宏觀自噬、自噬調(diào)控、自噬信號(hào)通路的調(diào)控、細(xì)胞生長(zhǎng)及內(nèi)在凋亡信號(hào)通路等；在細(xì)胞成分方面，基因主要富集在晚期內(nèi)涵體、細(xì)胞的焦點(diǎn)粘連、核被膜及液泡膜上；在分子功能方面上，基因主要富集在細(xì)胞黏附分子結(jié)合、泛素樣蛋白連接酶結(jié)合、鈣黏素結(jié)合及蛋白絲氨酸/蘇氨酸激酶活性等（圖2）。KEGG通路分析如圖3，這些基因主要與自噬-動(dòng)物、志賀氏癥、人乳頭瘤病毒感染、細(xì)胞凋亡、PI3K-AKT信號(hào)通路、鉑類耐藥、肝細(xì)胞癌、慢性粒細(xì)胞白血病、胰腺癌等通路顯著相關(guān)。

圖2 GO富集分析的氣泡圖

圖3 KEGG通路分析的圈圖和熱圖

2.3 預(yù)后相關(guān)自噬基因的篩選

將差異自噬基因的表達(dá)量與生存時(shí)間和生存狀態(tài)合并，按照P＜0.05的篩選標(biāo)準(zhǔn)，通過(guò)單因素Cox回歸分析，過(guò)濾得到5個(gè)與預(yù)后相關(guān)的差異自噬基因（RHEB、PPP1R15A、ATG101、BNIP3、NRG1），見圖4。

圖4 單因素Cox分析森林圖

2.4 自噬基因的預(yù)后模型的構(gòu)建

將這5個(gè)預(yù)后相關(guān)基因通過(guò)多因素Cox回歸分析進(jìn)行優(yōu)化，刪除相關(guān)性較高基因，結(jié)果這5個(gè)基因都與患者的預(yù)后顯著相關(guān)，可以用來(lái)構(gòu)建預(yù)后模型（見表1）。自噬基因的預(yù)后風(fēng)險(xiǎn)模型公式為（NRG1表達(dá)量×1.1207+BNIP3表達(dá)量×1.4002+ATG101表達(dá)量×1.3457-PPP1R15A表達(dá)量×1.1613-RHEB表達(dá)量×1.3279）。

表1 多因素Cox分析CCA預(yù)后影響因素

2.5 自噬預(yù)后模型的驗(yàn)證

根據(jù)模型計(jì)算每個(gè)患者風(fēng)險(xiǎn)值，18 例屬于高風(fēng)險(xiǎn)組，18 例為低風(fēng)險(xiǎn)組。為了驗(yàn)證自噬預(yù)后模型在預(yù)測(cè)臨床預(yù)后方面的可靠性，進(jìn)行以下分析。Kaplan-Meier生存分析中P=8.52e-05（P＜0.05），因此高、低風(fēng)險(xiǎn)組的患者在生存率上有差異的可能性較大（圖5），高風(fēng)險(xiǎn)組3 年生存率為18.7%，低風(fēng)險(xiǎn)組3年生存率為70.3%。圖6A為風(fēng)險(xiǎn)曲線圖，從左到右患者的風(fēng)險(xiǎn)值逐漸增大；圖6B為生存狀態(tài)圖，描述了風(fēng)險(xiǎn)值與生存時(shí)間/生存狀態(tài)的關(guān)系，隨著風(fēng)險(xiǎn)值增大死亡患者增多，隨著風(fēng)險(xiǎn)值減小死亡患者減少；圖6C為風(fēng)險(xiǎn)熱圖，描述了風(fēng)險(xiǎn)值與基因表達(dá)量的關(guān)系，隨著風(fēng)險(xiǎn)值的增加，4個(gè)基因的表達(dá)量增加（RHEB、ATG101、BNIP3、NRG1），表明這4個(gè)基因?yàn)楦唢L(fēng)險(xiǎn)基因；隨著風(fēng)險(xiǎn)值的增加，1個(gè)基因的表達(dá)量減少（PPP1R15A），表明這1個(gè)基因?yàn)榈惋L(fēng)險(xiǎn)基因。在單因素獨(dú)立預(yù)后分析中，風(fēng)險(xiǎn)值與生存時(shí)間及生存狀態(tài)顯著相關(guān)（HR1.382，95%CI1.160～1.645，P＜0.001）（圖7A），在多因素獨(dú)立預(yù)后分析中風(fēng)險(xiǎn)值與生存時(shí)間及生存狀態(tài)也顯著相關(guān)（HR1.427，95%CI1.161～1.756，P＜0.001）（圖7B），表明風(fēng)險(xiǎn)值可以作為獨(dú)立預(yù)后因子，進(jìn)一步驗(yàn)證了Cox模型在預(yù)測(cè)臨床預(yù)后方面的可靠性。在多指標(biāo)ROC曲線中，風(fēng)險(xiǎn)值的AUC值最大（圖8），表明風(fēng)險(xiǎn)得分有更好的預(yù)測(cè)臨床預(yù)后的能力。

圖5 生存曲線分析

圖6 風(fēng)險(xiǎn)值與生存時(shí)間/生存狀態(tài)、基因表達(dá)量的關(guān)系

圖7 單因素（A）和多因素（B）分析患者預(yù)后影響因素

3 討論

CCA作為臨床常見的惡性腫瘤，預(yù)后極差，可能與肝血管的侵犯、遠(yuǎn)處轉(zhuǎn)移、較早淋巴結(jié)的轉(zhuǎn)移及術(shù)后早期復(fù)發(fā)有關(guān)，因此，尋找可靠的早期診斷標(biāo)志物來(lái)完善預(yù)后評(píng)估系統(tǒng)，對(duì)于CCA患者預(yù)后有極大幫助。本研究在TCGA數(shù)據(jù)庫(kù)中篩選出49個(gè)差異表達(dá)的自噬相關(guān)基因DEATGs，且通過(guò)富集分析證實(shí)自噬參與CCA發(fā)生發(fā)展及耐藥的過(guò)程，本研究應(yīng)用單因素及多因素Cox回歸分析最終篩選出5個(gè)與CCA預(yù)后顯著相關(guān)的自噬基因（RHEB、PPP1R15A、ATG101、BNIP3、NRG1）。

RHEB基因?qū)儆赗AS家族的一員，編碼的蛋白屬于GTP結(jié)合蛋白，在與GTP結(jié)合的情況下Rheb被激活，當(dāng)GTP被水解時(shí)Rheb受到抑制。在TSC1-TSC2-TBC1D7/RHEB/mTORC1經(jīng)典的信號(hào)通路調(diào)節(jié)中，扮演重要角色。在參與調(diào)節(jié)自噬的發(fā)生發(fā)展、調(diào)控細(xì)胞分化及抑制腫瘤發(fā)面發(fā)揮重要作用[21]。在結(jié)腸癌組織中，p27KIP1受到RHEB的調(diào)控，RHEB/p27 激活mTORC1 信號(hào)通路來(lái)增強(qiáng)自噬，促進(jìn)腫瘤細(xì)胞的存活，加速結(jié)腸癌的進(jìn)展[22]。He等[23]的研究表明，F(xiàn)ADD基因在乳腺癌中顯著上調(diào)，F(xiàn)ADD基因敲減降低了乳腺癌細(xì)胞系中Rheb的表達(dá)，抑制了mTORC1的活性，抑制腫瘤細(xì)胞的存活。Chen等[24]在宮頸癌組織樣本中的研究表明，circMYLK與miR-1301-3p結(jié)合可以促宮頸癌細(xì)胞生長(zhǎng)，其機(jī)制是RHEB/mTORC1信號(hào)通路的激活。

PPP1R15A基因又稱為GADD34，在應(yīng)激生長(zhǎng)停止條件和DNA損傷劑治療后其轉(zhuǎn)錄水平增加。PPP1R15A編碼的蛋白質(zhì)可以調(diào)節(jié)電離輻射后的凋亡，在促進(jìn)細(xì)胞死亡和未折疊蛋白反應(yīng)中起著重要的作用。Holczer等[25]在肝癌細(xì)胞中的研究發(fā)現(xiàn)，GADD34的表達(dá)被抑制后，導(dǎo)致mTOR的快速激活和自噬水平下調(diào)，隨后凋亡細(xì)胞死亡。GADD34表達(dá)的減少顯著抑制腫瘤，并導(dǎo)致MDSCs和T細(xì)胞的積累減少，抑制GADD34減少了MDSCs分泌的血管上皮生長(zhǎng)因子α和轉(zhuǎn)化生長(zhǎng)因子β[26]。

ATG101基因位于染色體12q13.13，編碼的蛋白是ULK1 自噬復(fù)合物的重要組成部分，與ATG13結(jié)合維持ULK1復(fù)合物的穩(wěn)定性，進(jìn)而調(diào)控自噬[27]。研究表明，PTCH1是已知的腫瘤抑制因子，抑制腫瘤細(xì)胞增殖和誘導(dǎo)細(xì)胞凋亡，減少PTCH1的表達(dá)量會(huì)增加ATG101介導(dǎo)的自噬通量來(lái)促進(jìn)癌細(xì)胞增殖[28]。

BNIP3基因的全稱為BCL2相互作用蛋白3基因，該基因編碼一種線粒體蛋白，含有BH3 結(jié)構(gòu)域，可以充當(dāng)促凋亡因子。細(xì)胞內(nèi)存在很多蛋白可以誘導(dǎo)線粒體自噬和凋亡，但BNIP3能將線粒體自噬和凋亡串聯(lián)在一起，功能更加特殊[29]。m6A去甲基化酶（FTO）在乳腺癌中顯著高表達(dá)，BNIP3是其下游靶點(diǎn)，促進(jìn)腫瘤細(xì)胞增殖和轉(zhuǎn)移[30]。Li等[31]在胰腺癌的研究中發(fā)現(xiàn)，BNIP3甲基化抑制了線粒體介導(dǎo)的胰腺癌細(xì)胞凋亡的誘導(dǎo)。

NRG1基因編碼一種膜糖蛋白，介導(dǎo)多種細(xì)胞信號(hào)通路，該基因有多種異構(gòu)體，該基因的失調(diào)與癌癥等疾病密切相關(guān)。關(guān)于前列腺癌中一項(xiàng)研究發(fā)現(xiàn)，腫瘤微環(huán)境下產(chǎn)生的NRG1，可以促進(jìn)腫瘤細(xì)胞的耐藥性。研究表明，阻斷NRG1會(huì)減少其介導(dǎo)的HER3的激活，抑制NRG1信號(hào)通路的誘導(dǎo)，起到抑制胰腺腫瘤細(xì)胞（PC）和腫瘤相關(guān)成纖維細(xì)胞（CAFS）相關(guān)的促腫瘤細(xì)胞增殖和轉(zhuǎn)移的作用[32]。Trombetta 等[33]研究發(fā)現(xiàn)，肺腺癌中NRG1的高表達(dá)會(huì)通過(guò)PI3K-AKT和MAPK信號(hào)通路異常激活ErbB2/ErbB3，使肺腺癌細(xì)胞侵襲性增強(qiáng)。

本研究基于上述5 個(gè)自噬基因構(gòu)建預(yù)后風(fēng)險(xiǎn)模型，發(fā)現(xiàn)風(fēng)險(xiǎn)模型的AUC值（0.906）要明顯優(yōu)于其他指標(biāo)，因此可以作為CCA的獨(dú)立預(yù)后指標(biāo)。當(dāng)然本研究還存在一定的局限性，譬如本研究屬于純數(shù)據(jù)庫(kù)分析，缺少相關(guān)的功能實(shí)驗(yàn)進(jìn)行驗(yàn)證。