基于TCGA 數(shù)據(jù)庫(kù)的肺腺癌自噬基因預(yù)后模型建立與應(yīng)用

陳彬彬 葉威 黃約諾 方濤 林雪嬌 董曉亭

肺腺癌（lung adenocarcinoma，LUAD）是肺癌中常見(jiàn)的病理類型，建立準(zhǔn)確預(yù)測(cè)LUAD 患者預(yù)后的工具有利于指導(dǎo)臨床診治[1]。自噬（autophagy）是溶酶體將細(xì)胞內(nèi)受損、變性或衰老的蛋白質(zhì)以及細(xì)胞器進(jìn)行消化降解的過(guò)程[2]。自噬在腫瘤發(fā)生、發(fā)展中發(fā)揮著雙刃劍的作用：在正常情況下，細(xì)胞自噬防止致癌的蛋白質(zhì)和細(xì)胞器等積累，維持細(xì)胞的穩(wěn)定；在腫瘤細(xì)胞中，細(xì)胞自噬可以降解非功能細(xì)胞器提供營(yíng)養(yǎng)給癌細(xì)胞促進(jìn)生長(zhǎng)[2]。ATGs 的下調(diào)能夠通過(guò)抑制Hedgehog 信號(hào)通路增強(qiáng)抗LUAD 的療效[3]。本研究利用出公共數(shù)據(jù)庫(kù)篩選出生存相關(guān)的自噬基因（autophagy gene，ATG），建立LUAD-ATG 預(yù)后模型并加以驗(yàn)證其應(yīng)用價(jià)值。

1 材料與方法

1.1 數(shù)據(jù)下載整理 從癌癥基因數(shù)據(jù)庫(kù)（The cancer Genome Atlas，TCGA）中選擇LUAD 患者的RNA 表達(dá)數(shù)據(jù)，提取與患者相關(guān)的臨床數(shù)據(jù)。從人類自噬數(shù)據(jù)庫(kù)中（Human Autophagy Database http：//www.autophagy.lu）查找自噬相關(guān)基因。本文的所有下載數(shù)據(jù)截至2021 年1 月23 日。

1.2 自噬差異基因表達(dá) 從TCGA 數(shù)據(jù)庫(kù)中篩選出LUAD 組織和正常組織樣本中ATGs 的表達(dá)數(shù)據(jù)。通過(guò)R（v4.0.3）軟件運(yùn)用Wilcox 檢驗(yàn)篩選LUAD 癌組織樣品和正常樣品中差異表達(dá)的ATGs。篩選標(biāo)準(zhǔn)為：｜logFC｜＞1，F(xiàn)DR＜0.05。為了進(jìn)一步探索潛在生物學(xué)功能，調(diào)用中“org.Hs.eg.db”“enrichplot”和“gplot2”包進(jìn)行GO 功能注釋分析，P＜0.05 被認(rèn)為是有意義的生物學(xué)通路。

1.3 預(yù)后相關(guān)的自噬基因 采用Survival 包先后運(yùn)用單、多因素COX 回歸篩選出影響患者生存期的ATGs（P 值＜0.05）。此外，運(yùn)用多因素Cox 回歸計(jì)算出預(yù)后相關(guān)ATGs 的相關(guān)系數(shù)（coefficient，coef）并計(jì)算每個(gè)樣本的風(fēng)險(xiǎn)值（risk score，RS），從而建立預(yù)后模型。

RS=a1*coef1+a2*coef2......an*coefn

a 代表預(yù)后相關(guān)TAGs 的表達(dá)量，coef 代表預(yù)后相關(guān)TAGs 的相關(guān)系數(shù)。

然后根據(jù)RS 的中位數(shù)將所有患樣本分為高風(fēng)險(xiǎn)組和低風(fēng)險(xiǎn)組，運(yùn)用Kaplan-Meier 生存曲線評(píng)價(jià)高、低風(fēng)險(xiǎn)組的生存差異，P＜0.05 被認(rèn)為生存差異有統(tǒng)計(jì)學(xué)意義。進(jìn)一步判斷獨(dú)立預(yù)后因素，通過(guò)Survival 包對(duì)年齡、性別、TNM 分期、臨床分期和RS 多個(gè)因素進(jìn)行單因素COX 分析。TMN 分期能夠?qū)δ[瘤進(jìn)行準(zhǔn)確的分期，指導(dǎo)臨床腫瘤治療，其中T 是指照腫瘤大小，N 是指淋巴狀態(tài)，M 是指遠(yuǎn)處轉(zhuǎn)移狀態(tài)。為了排除混雜因素的干擾，對(duì)上述多個(gè)影響因素進(jìn)行多因素COX 分析。最后采用受試者工作特征（receiver operator characteristic，ROC）曲線評(píng)估預(yù)測(cè)模型的準(zhǔn)確性。其中曲線下面積（area under the curve，AUC），AUC＞0.06 被認(rèn)為是可以接受的預(yù)測(cè)。

1.4 臨床相關(guān)性分析 運(yùn)用R 軟件中beeswarm 包繪采用Wilcox 檢驗(yàn)來(lái)比較預(yù)后相關(guān)的ATGs 與腫瘤大小、淋巴結(jié)狀態(tài)、遠(yuǎn)處轉(zhuǎn)移、臨床分期、性別、年齡之間的相關(guān)性。P＜0.05 被認(rèn)為差異有統(tǒng)計(jì)學(xué)意義。

2 結(jié)果

2.1 自噬相關(guān)差異基因的篩選 TCGA 數(shù)據(jù)庫(kù)中共下載了535 例LUAD 組織樣本和59 例癌旁組織樣本的mRNA 表達(dá)數(shù)據(jù)及對(duì)應(yīng)的臨床信息。采用Wilcox 檢驗(yàn)所得P 值的對(duì)數(shù)-lg（pvalue）為縱坐標(biāo)，以差異倍數(shù)的對(duì)數(shù)值log2（fold change，F(xiàn)C）為橫坐標(biāo)繪制火山圖。以｜logFC｜＞1，F(xiàn)DR＜0.05 為標(biāo)準(zhǔn)篩選出30 個(gè)差異表達(dá)ATGs，其中HSPA5、IFNG、P4HB、ATIC、PTK6、ERO1A 等18 個(gè)基因上調(diào)，DLC1、CCL2、DAPK2、PRKCQ、FOS 等12 個(gè)基因下調(diào)，進(jìn)一步通過(guò)箱線圖展示ATG 在腫瘤組織和癌旁組織的差異表達(dá)情況（P＜0.05），見(jiàn)圖1。

圖1 肺腺癌和癌旁樣本的差異表達(dá)的自噬基因

為了揭示差異TAGs 的潛在功能，用clusterProfiler 包對(duì)差異表達(dá)的ATGs 進(jìn)行功能分析，分別展示了分子功能（molecular function）、生物過(guò)程（biological process）、細(xì)胞組成（cellular component）三個(gè)部分GO分類條目（P＜0.05）。圖2 顯示，主要與內(nèi)源性凋亡、未折疊蛋白反應(yīng)、ErbB-2 信號(hào)通路等相關(guān)。

圖2 差異表達(dá)的自噬基因的GO 富集分析

2.2 ATG 預(yù)后模型的建立 通過(guò)單因素Cox 回歸分析以P ＜0.01 為閾值篩選出ERN1、ATG16L2、VEGFA、CCR2、ATG2A、DDIT3、CFLAR、IRGM 等共有28 個(gè)可能影響患者總體生存期的ATGs，其中HR＞1表示高風(fēng)險(xiǎn)基因，HR＜1 表示低風(fēng)險(xiǎn)因素，P＜0.05 表示差異有統(tǒng)計(jì)學(xué)意義（見(jiàn)圖3）。進(jìn)一步采用多因素Cox 分析優(yōu)化，篩選出RAC1、SQSTM1、CD46、NRG3、IKBKB、VMP1、WIPI1、FKBP1B、IKBKE 在內(nèi)的9 個(gè)預(yù)后相關(guān)ATGs。其中CD46、NRG3、IKBKB、VMP1、IKBKE、WIPI1 被認(rèn)為是危險(xiǎn)因素，高表達(dá)與較差生存期相關(guān)；高表達(dá)的RAC1、SQSTM1、FKBP1B 被認(rèn)為是保護(hù)因素，過(guò)度表達(dá)能夠延長(zhǎng)生存期。見(jiàn)表1。

圖3 肺腺癌差異表達(dá)的自噬基因單因素COX 回歸分析注：Hazard ration（HR）為風(fēng)險(xiǎn)比；pvalue 為P值

表1 9 個(gè)自噬基因的多因素COX 回歸分析

根據(jù)風(fēng)險(xiǎn)值公式計(jì)算的各個(gè)樣本RS 值并將LUAD 患者分為高風(fēng)險(xiǎn)組、低風(fēng)險(xiǎn)組。為了進(jìn)一步驗(yàn)證結(jié)果準(zhǔn)確性，圖4A 運(yùn)用Kaplan-Meier 生存曲線進(jìn)行分析發(fā)現(xiàn)低RS 組具有更長(zhǎng)的生存期（P＜0.05）。此外，通過(guò)將所有的患者跟進(jìn)RS 進(jìn)行排名，以觀察其生存率。隨著RS 越高，死亡患者數(shù)量越多（見(jiàn)圖4B-C）。圖4D 顯示，LUAO 自噬相關(guān)模型中ATGs 在高風(fēng)險(xiǎn)組、低風(fēng)險(xiǎn)組的表達(dá)量。

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

2.3 生存分析及ROC 曲線 納入Cox 回歸分析的自變量中的連續(xù)變量為RS，年齡、TMN 分期、病理分期等為分類變量。單因素Cox 回歸模型分析結(jié)果顯示，RS、腫瘤原發(fā)灶情況是影響癌總體生存期的危險(xiǎn)因素（P＜0.05），見(jiàn)圖5A。多因素Cox 回歸模型分析結(jié)果顯示，RS、腫瘤原發(fā)灶情況是LUAD 患者預(yù)后的獨(dú)立影響因素（P＜0.05），見(jiàn)圖5B。

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

ROC 曲線分析影響因素與5 年的生存率之間的關(guān)系，RS 的ROC 曲線的AUC 值最大，為0.659，認(rèn)為RS 可作為L(zhǎng)UAD 獨(dú)立預(yù)后指標(biāo)。見(jiàn)圖6。

圖6 LUAD 預(yù)后預(yù)測(cè)模型的ROC 曲線分析

2.4 不同臨床特征下自噬相關(guān)基因表達(dá)水平情況進(jìn)一步分析9 個(gè)預(yù)后相關(guān)ATGs 與臨床特征的相關(guān)性（見(jiàn)圖7-8），IKBKB 與臨床分期、腫瘤原發(fā)灶情況、腫瘤淋巴結(jié)情況顯著相關(guān)（P＜0.05）；NRG3 在腫瘤淋巴結(jié)情況顯著相關(guān)（P＜0.05）；FKBPIB 與性別具有相關(guān)性（P＜0.05），在女性中FKBPIB 的表達(dá)量更多。

圖7 CD46、FKBP1B、IKBKB、IKBKE、NRG3 與臨床特征之間的相關(guān)性

圖8 RAC1、SQSTM1、VMP1、WIPI1 與臨床特征之間的相關(guān)性

3 討論

自噬既是一種促細(xì)胞生存機(jī)制，也是一種程序性細(xì)胞死亡機(jī)制。越來(lái)越多證據(jù)證明自噬不僅能夠?yàn)榘┘?xì)胞提供能量，而且能夠從各個(gè)途徑促進(jìn)癌癥侵襲及轉(zhuǎn)移[4]。Wang等[5]對(duì)393 例接受了放射治療的北美非小細(xì)胞肺癌患者進(jìn)行研究，通過(guò)對(duì)ATG2B、ATG10、ATG12 ATG16L2 在內(nèi)的4 個(gè)ATGs 中的9個(gè)潛在功能單核苷酸多態(tài)性進(jìn)行基因分型，發(fā)現(xiàn)ATG16L2 rs10898880 CC 變異基因型患者的無(wú)局部復(fù)發(fā)生存率、無(wú)局部區(qū)域進(jìn)展生存率和總生存率較好，表明該位點(diǎn)的多態(tài)性參與了抗腫瘤治療和肺部炎癥的過(guò)程，也提示自噬以某種方式參加腫瘤的病理變化。NRG3 是神經(jīng)調(diào)節(jié)蛋白家族之一（其它成員包括NRG1、NRG2 等），是ErbB 家族受體酪氨酸激酶的配體[6]。近年研究顯示，NRG-1 與非小細(xì)胞肺癌、食管癌、甲狀腺癌等多種腫瘤的發(fā)生、發(fā)展密切相關(guān)，其作用過(guò)程可能涉及腫瘤相關(guān)的信號(hào)通路分子改變[7-9]。但是關(guān)于NRG3 尚未針對(duì)癌癥進(jìn)行過(guò)研究。IKBKB 是組成IKK 復(fù)合體的重要催化亞基之一，其在IKK 復(fù)合體激活NF-κB 過(guò)程中發(fā)揮重要作用[10]。而NF-κB 異常激活存在腫瘤細(xì)胞和腫瘤微環(huán)境中，被認(rèn)為在乳腺癌、肝癌、肺癌、前列腺癌在內(nèi)的多種腫瘤發(fā)生、發(fā)展發(fā)揮著重要的作用[11-14]。此外，IKBKB的高表達(dá)與乳腺癌耐順鉑耐藥性相關(guān)，可能通過(guò)抑制FOXO3 而導(dǎo)致耐藥性的產(chǎn)生[15]。Qi等[16]通過(guò)A549細(xì)胞中轉(zhuǎn)染pcDNA3.1/IKBKB 質(zhì)粒過(guò)表達(dá)了IKBKB基因，發(fā)現(xiàn)IKBKB 高表達(dá)的A549 細(xì)胞耐藥性明顯增加，明顯抑制細(xì)胞凋亡。這可能與激活NF-κB 信號(hào)通路有關(guān)，導(dǎo)致耐藥性的產(chǎn)生。

對(duì)差異表達(dá)的ATGs 進(jìn)行GO 分析發(fā)現(xiàn)，較多的基因富集在凋亡信號(hào)通路、未折疊蛋白反應(yīng)，這與以前的研究結(jié)論是一致的：腫瘤細(xì)胞通過(guò)引起的錯(cuò)誤折疊的蛋白質(zhì)和受損的細(xì)胞器來(lái)消除細(xì)胞器的壓力[17]。ErbB-2 信號(hào)通路在腫瘤發(fā)生中發(fā)揮著重要的作用，ERBB2 種族突變存在于肺癌患者中[18]。本研究發(fā)現(xiàn)，ErbB-2 信號(hào)通路上富集到多個(gè)ATG，側(cè)面說(shuō)明ATGs 也可通過(guò)這種經(jīng)典的癌癥通路影響腫瘤進(jìn)展。

進(jìn)一步使用單、多COX 回歸分析后確定了9 個(gè)與生存相關(guān)的ATGs（RAC1、SQSTM1、CD46、NRG3、IKBKB、VMP1、WIPI1、FKBP1B、IKBKE）并構(gòu)建預(yù)后風(fēng)險(xiǎn)模型。根據(jù)風(fēng)險(xiǎn)評(píng)分中位值分為高、低風(fēng)險(xiǎn)組，不同的RS 分組與生存相關(guān)，擁有不同臨床特征（P＜0.05）。通過(guò)繪制RS 分布、Kaplan-Meier 分析證實(shí)，相較于高RS 組，低RS 組患者的預(yù)后更好。ROC曲線具有簡(jiǎn)單、直觀的特點(diǎn)，用來(lái)評(píng)價(jià)模型的準(zhǔn)確性。根據(jù)ROC 曲線，9 個(gè)ATGs 建立的ATGs 預(yù)后模型AUC 為0.659，進(jìn)一步證實(shí)預(yù)后模型的穩(wěn)定性和可靠性。此外，對(duì)RS 和其他因素進(jìn)行了單、多因素Cox 回歸分析，證實(shí)RS 可作為L(zhǎng)UAD 患者的獨(dú)立預(yù)后預(yù)測(cè)因子。

近年出現(xiàn)了許多不同的統(tǒng)計(jì)學(xué)方法構(gòu)建的腫瘤預(yù)后模型，包括乳腺癌、肺癌、結(jié)直腸癌等[19-21]，有利于臨床診治，推進(jìn)精準(zhǔn)醫(yī)學(xué)發(fā)展。本研究構(gòu)建的ATG-LUAD 預(yù)后模型仍存下述局限：（1）本模型處于初步建立階段，仍需要更多臨床數(shù)據(jù)驗(yàn)證模型的可靠性；（2）本模型以TCGA 數(shù)據(jù)庫(kù)為基礎(chǔ)構(gòu)建，缺少亞洲人群的轉(zhuǎn)錄組數(shù)據(jù)，在預(yù)測(cè)國(guó)人生存期方面的準(zhǔn)確性較差。