儲(chǔ)昭陽， 朱向明， 魯柯兵， 龔儒杰， 馬平川， 徐凱慧

（1.皖南醫(yī)學(xué)院研究生學(xué)院， 安徽 蕪湖241002； 2.皖南醫(yī)學(xué)院弋磯山醫(yī)院超聲醫(yī)學(xué)科； 3.臺(tái)州市黃巖區(qū)婦幼保健院超聲醫(yī)學(xué)科）

膀胱癌（bladder cancer， BC） 2018 年全球新發(fā)患者近54.9 萬人， 相關(guān)死亡20 萬人［1］， 患者的5 年存活率持平［2］。 開發(fā)新的診斷和治療BC 的方法至關(guān)重要， 而識(shí)別新的與BC 診斷和預(yù)后相關(guān)的生物標(biāo)志物也同樣重要。 RNA 結(jié)合蛋白（RNA binding proteins， RBPs） 是一類與多種類型的RNA相互作用的蛋白質(zhì)， 包括rRNAs、 ncRNAs、 snRNAs、miRNAs、 mRNAs、 tRNAs 和snoRNAs。 迄今為止，在人類基因組中已經(jīng)發(fā)現(xiàn)了超過1 500 個(gè)RBPs 基因［3］。 這些RBPs 在維持細(xì)胞的生理平衡方面發(fā)揮著重要作用， 尤其是在發(fā)育過程和應(yīng)激反應(yīng)中。RBPs 可以以一種結(jié)構(gòu)或序列依賴的方式與目標(biāo)RNA 結(jié)合， 形成核糖核酸蛋白復(fù)合物， 在轉(zhuǎn)錄后水平調(diào)控mRNA 的穩(wěn)定性、 RNA 加工、 剪接、 定位、 輸出和翻譯［4］。 在過去的幾十年里， 許多研究揭示了RBPs 在腫瘤中表達(dá)異常， 影響了mRNA向蛋白的轉(zhuǎn)化， 并參與腫瘤的發(fā)生［5-7］。 其中， 只有少數(shù)RBPs 被深入研究并發(fā)現(xiàn)在人類癌癥中發(fā)揮關(guān)鍵作用［8-12］。 RBPs 的系統(tǒng)功能研究將有助于我們?nèi)媪私馄湓谀[瘤中的作用。 因此， 本研究從癌癥基因組圖譜 （Cancer Genome Atlas， TCGA）數(shù)據(jù)庫下載了BC 的RNA 測序和臨床病理數(shù)據(jù)。隨后通過高通量生物信息學(xué)分析鑒定了BC 組織和正常膀胱組織間異常表達(dá)的RBPs， 并系統(tǒng)地探索了它們的潛在功能和分子機(jī)制。 本研究確定了一些可能影響B(tài)C 預(yù)后的RBPs， 促進(jìn)了對(duì)BC 進(jìn)展的分子機(jī)制的理解。 這些RBPs 可能為診斷和預(yù)后提供潛在的生物標(biāo)志物。

1 材料與方法

1.1 數(shù)據(jù)處理 從TCGA （https： ／ ／portal.gdc.cancer.gov／） 數(shù)據(jù)庫下載19 例正常膀胱組織和411 例BC 組織的RNA 測序數(shù)據(jù)集以及相應(yīng)的臨床數(shù)據(jù)， 其中數(shù)據(jù)的篩選為包含組織的病理分期以及生存期。 為了鑒定正常膀胱組織和BC 組織之間的差異表達(dá)基因， 所有原始數(shù)據(jù)使用LIMMA（ http： ／ ／www.bioconductor.org／packages／release／bioc／html／limma.html） 軟 件 包 進(jìn) 行 預(yù) 處 理， 根據(jù)| logFC （fold change）| ≥1 和偽發(fā)現(xiàn)率（FDR）＜0.05 標(biāo)準(zhǔn)， 進(jìn)行Wilcox 檢驗(yàn)篩選出RBPs 相關(guān)的差異基因。

1.2 KEGG 途徑與GO 富集分析 通過GO 富集和京都基因與基因組百科全書（KEGG） 途徑分析，綜合檢測這些差異表達(dá)的RBPs 的生物學(xué)功能。GO 分析術(shù)語包括細(xì)胞成分 （CC）、 分子功能（MF） 和生物過程（BP）。 所有富集分析均利用org.Hs.eg.db、 DOSE、 clusterProfiler、 enrichplot 軟件包完成， P 和FDR 值均＜0.05 為差異有統(tǒng)計(jì)學(xué)意義。

1.3 PPI 網(wǎng)絡(luò)的構(gòu)建和模塊篩選 將差異表達(dá)的RBPs 提交到STRING 數(shù)據(jù)庫（http： ／ ／www.string-db.org／）［13］， 以識(shí)別蛋白質(zhì)-蛋白質(zhì)相互作用信息。 使用Cytoscape 3.8.0 軟件進(jìn)一步構(gòu)建PPI 網(wǎng)絡(luò)并進(jìn)行可視化。 使用聚類構(gòu)建功能模塊（MCODE）插件得到PPI 網(wǎng)絡(luò)中的重要模塊和基因， MCODE得分和節(jié)點(diǎn)數(shù)均大于5。 P＜0.05 為差異有統(tǒng)計(jì)學(xué)意義。

1.4 BC 患者RBPs 篩選及預(yù)后生存模型構(gòu)建 使用Survival 軟件包對(duì)TCGA 數(shù)據(jù)集差異表達(dá)的RBPs進(jìn)行單因素Cox 回歸分析。 采用對(duì)數(shù)秩檢驗(yàn)進(jìn)一步篩選有意義的候選基因。 隨后， 在初步篩選出有意義的候選基因的基礎(chǔ)上， 構(gòu)建了多因素Cox風(fēng)險(xiǎn)回歸模型， 并計(jì)算風(fēng)險(xiǎn)分?jǐn)?shù)來評(píng)估患者的預(yù)后。 每個(gè)樣本的風(fēng)險(xiǎn)評(píng)分公式： Risk score ＝β1×Exp1＋β2×Exp2 ＋βi×Expi。 其中β 表示系數(shù)值，Exp 表示基因表達(dá)水平。 根據(jù)中位風(fēng)險(xiǎn)評(píng)分生存分析， 將BC 患者分為低危組和高危組。 采用對(duì)數(shù)秩檢驗(yàn)比較2 個(gè)亞組之間的總生存期（OS） 差異。此外， 使用Survival ROC 軟件包進(jìn)行ROC 曲線分析以評(píng)估上述模型的預(yù)測能力。 將TCGA 數(shù)據(jù)隨機(jī)分為2 組， 一組作為train 組分析， 另一組作為test 組進(jìn)行模型驗(yàn)證， 以證實(shí)該預(yù)測模型的預(yù)測能力。 最后， 使用RMS 軟件包繪制列線圖， 預(yù)測患者的生存期。

1.5 驗(yàn)證預(yù)后模型RBPs 的預(yù)測價(jià)值 分別對(duì)train 組及test 組進(jìn)行單因素和多因素獨(dú)立預(yù)后分析。 單因素獨(dú)立預(yù)后分析是將每個(gè)因素跟生存時(shí)間和生存狀態(tài)進(jìn)行比較， 觀察與患者生存是否相關(guān)； 多因素獨(dú)立預(yù)后分析是將這里面的因素一起輸入， 跟生存時(shí)間和生存狀態(tài)進(jìn)行比較， 并考慮這些因素之間的一個(gè)相關(guān)性。 如果train 組及test組分析結(jié)果P 值均＜0.05， 則說明模型可以獨(dú)立于其他的臨床現(xiàn)狀作為獨(dú)立的預(yù)后因子。 之后分別探討預(yù)后模型中的各個(gè)RBPs 在BC 中的預(yù)后價(jià)值， 使用Kaplan-Meier Plotter Server （https： ／ ／kmplot.com／analysis／） 在線工具繪圖儀來驗(yàn)證預(yù)后模型中每個(gè)RBPs 與OS 間的關(guān)系。

1.6 統(tǒng)計(jì)學(xué)方法 使用Perl 5.30.2 進(jìn)行數(shù)據(jù)表達(dá)矩陣的處理， 應(yīng)用R4.0.0 軟件進(jìn)行統(tǒng)計(jì)學(xué)分析及相應(yīng)圖形繪制， 采用Wilcox 檢驗(yàn)篩選差異基因，采用ggplot 包進(jìn)行圖形繪制， Survival 包進(jìn)行單因素和多因素Cox 比例回歸模型篩選， 并建立多基因預(yù)后模型。 使用Survival ROC 包計(jì)算ROC 曲線評(píng)價(jià)模型的有效性， 并計(jì)算曲線下面積。 P＜0.05為差異有統(tǒng)計(jì)學(xué)意義。

2 結(jié)果

2.1 篩選得到差異表達(dá)的BC 組織中RBPs 將TCGA 得到的數(shù)據(jù)通過perl 和R 語言處理后， 處理數(shù)據(jù)共納入差異表達(dá)的RBPs 有1 495 個(gè)［3］， 其中109 個(gè)符合本研究篩選標(biāo)準(zhǔn)， 其中包括52 個(gè)上調(diào)的RBPs 和57 個(gè)下調(diào)的RBPs。 這些差異表達(dá)的RBPs 的分布用熱圖和火山圖來表示。

2.2 差異表達(dá)的BC 組織中RBPs 的GO 和KEGG通路富集分析 為研究所鑒定的RBPs 的功能和機(jī)制， 將這些不同表達(dá)的RBPs 分為2 組： 表達(dá)上調(diào)組和表達(dá)下調(diào)組。 然后， 將這些差異表達(dá)的RBPs進(jìn)行富集分析。 結(jié)果表明， 上調(diào)差異的RBPs 的GO 富集在催化活性、 核糖核酸作用活性、 核苷酰轉(zhuǎn)移 酶 活 性、 mRNA 3′ 端- 非 編 碼 區(qū) 結(jié) 合 區(qū)（UTR） 等分子功能中； 就細(xì)胞組成而言， 上調(diào)差異的RBPs 的GO 富集在細(xì)胞質(zhì)核糖核蛋白顆粒、核糖核蛋白顆粒； 生物過程的富集主要包括大分子甲基化、 甲基化、 RNA 分解過程等（圖1A）；下調(diào)差異的RBPs 的GO 富集分析中分子功能主要包含mRNA 3′端-UTR 結(jié)合區(qū)、 mRNA 3′端-UTR富含AU 的結(jié)合區(qū)、 富含AU 的元素結(jié)合區(qū)等； 細(xì)胞組成主要包含細(xì)胞質(zhì)核糖核蛋白顆粒、 核糖核蛋白顆粒、 P-顆粒體等結(jié)構(gòu)； 生物過程主要包括調(diào)節(jié)mRNA 代謝過程、 調(diào)節(jié)RNA 剪接、 調(diào)節(jié)翻譯等（圖1B）。 上調(diào)差異的RBPs 相關(guān)基因的KEGG通路分析結(jié)果表明， DEGs 主要參與癌癥中的miR?NAs、 mRNA 監(jiān)測途 徑、 DNA 復(fù)制等信號(hào)通路（圖1C）； 下調(diào)差異的RBPs 相關(guān)基因的KEGG 通路分析結(jié)果表明， DEGs 相關(guān)基因主要參與孕酮介導(dǎo)的卵母細(xì)胞成熟、 卵母細(xì)胞減數(shù)分裂、 mRNA監(jiān)測途徑等信號(hào)通路（圖1D）。

圖1 差異表達(dá)RBPs 的GO 富集分析（A、 B） 和KEGG 途徑（C、 D）

2.3 蛋白工作網(wǎng)絡(luò)及核心模塊 Cytoscape 得出88個(gè)RBPs （45 個(gè)上調(diào)， 43 個(gè)下調(diào)） 對(duì)應(yīng)的166 個(gè)蛋白互作的關(guān)系。 使用MCODE 模塊構(gòu)建出7 個(gè)子網(wǎng)絡(luò)， 并得出16 個(gè)基因是與BC 相關(guān)的關(guān)鍵RBPs（這些RBPs 在DNA 烷基化、 DNA 甲基化或去甲基化、 翻譯的負(fù)調(diào)控等方面都有顯著的富集）。

2.4 預(yù)后RBPs 篩選及生存模型構(gòu)建 通過篩選得到11 個(gè)與預(yù)后相關(guān)的RBPs 相關(guān)基因（圖2A），通過COX 回歸降維分析構(gòu)建預(yù)后模型（圖2B），由2 個(gè)高風(fēng)險(xiǎn)的RBPs 基因（TRIM71、 DARS2） 及3 個(gè)低風(fēng)險(xiǎn)的RBPs 基因（APOBEC3H、 PABPC3、MTG1） 組成； 模型在train 組及test 組中ROC 曲線下面積分別為0.701、 0.664 （圖2C、 D）， 反映了模型預(yù)測BC 患者的生存具有良好的診斷效能，train 組及test 組風(fēng)險(xiǎn)曲線（圖2E、 F） 也驗(yàn)證了模型預(yù)測BC 患者的生存的準(zhǔn)確性。 繪制train 組及test 組的生存曲線（圖2G、 H）， 可以看出在train 組及test 組中， 高低風(fēng)險(xiǎn)2 組生存率的差異有統(tǒng)計(jì)學(xué)意義， 說明得到的模型可以預(yù)測BC 患者的生存。

2.5 列線圖 根據(jù)5 個(gè)RBPs 繪制出列線圖（圖3）， 可以根據(jù)列線圖去預(yù)測患者的生存期。

2.6 獨(dú)立預(yù)后分析 將train 組及test 組通過多元回歸分析結(jié)果繪制森林圖（圖4）， 得出預(yù)后模型P 值均＜0.05， 說明模型可以獨(dú)立于其他的臨床現(xiàn)狀作為獨(dú)立的預(yù)后因子。

2.7 模型RBPs 的預(yù)測價(jià)值的驗(yàn)證 Kaplan-Meier Plotter Server 鑒定了預(yù)后模型中的5 個(gè)（TRIM71、DARS2、 APOBEC3H、 PABPC3 和MTG1） RBPs 與BC 患者的總生存期顯著相關(guān)（圖5）。

圖2 A： 單變量Cox 回歸分析中與BC 預(yù)后相關(guān)的RBPs； B： 多變量Cox 回歸分析中預(yù)后相關(guān)的RBPs； C、 D： 基于風(fēng)險(xiǎn)評(píng)分預(yù)測OS 的ROC 曲線； E、 F： 風(fēng)險(xiǎn)評(píng)分分布、 表達(dá)熱圖和生存狀態(tài)； G、 H： 低風(fēng)險(xiǎn)和高風(fēng)險(xiǎn)亞組的生存曲線

圖3 TCGA 隊(duì)列中預(yù)測BC 患者1 年、 2 年和3年OS 的列線圖

3 討論

RBPs 失調(diào)已被報(bào)道在各種惡性腫瘤［6-13］。 然而， 只有一部分RBPs 被深入研究， 部分證實(shí)它們參與了癌癥［9-13］的發(fā)生和發(fā)展。 本研究基于TCGA的BC 數(shù)據(jù)確定了腫瘤組織和正常膀胱組織之間的109 個(gè)不同表達(dá)的RBPs， 系統(tǒng)分析了相關(guān)的生物學(xué)途徑， 構(gòu)建了這些RBPs 的PPI 網(wǎng)絡(luò)， 對(duì)異常表達(dá)RBPs 進(jìn)行單變量Cox 回歸分析、 生存分析、 多元逐步Cox 回歸分析和ROC 分析， 進(jìn)一步探討其生物學(xué)功能和臨床意義。 本研究構(gòu)建了一個(gè)基于5個(gè)預(yù)后相關(guān)的RBPs 基因預(yù)測BC 預(yù)后的風(fēng)險(xiǎn)模型。 這些發(fā)現(xiàn)可能有助于開發(fā)新的BC 患者診斷和預(yù)后的生物標(biāo)志物。

表達(dá)的功能通路富集分析顯示， 不同RBPs 在翻譯的調(diào)控、 RNA 加工、 RNA 代謝中具有較高富集程度。 已有研究證明翻譯的調(diào)控、 RNA 加工、RNA 代謝與多種人類疾病的發(fā)生和發(fā)展有關(guān)［14-17］。 RNA 穩(wěn)定性的轉(zhuǎn)錄后調(diào)控是基因表達(dá)過程中的一個(gè)重要環(huán)節(jié)。 RBPs 可以與RNA 相互作用形成核糖核酸蛋白復(fù)合物， 從而增加靶mRNA的穩(wěn)定性， 促進(jìn)基因表達(dá)， 在各種疾病的進(jìn)展中發(fā)揮關(guān)鍵作用。 此外， 我們還建立了這些不同表達(dá)的RBPs 的蛋白-蛋白相互作用網(wǎng)絡(luò)， 接下來利用TCGA 隊(duì)列進(jìn)行分析， 通過多元逐步COX 回歸分析， 建立預(yù)測BC 預(yù)后的風(fēng)險(xiǎn)模型。 ROC 曲線分析顯示， 這5 個(gè)RBPs 具有較好的診斷能力， 可以選擇預(yù)后較差的BC 患者。 在預(yù)后模型的5 個(gè)RBPs中， 許多已被證明在腫瘤的發(fā)展和進(jìn)展中發(fā)揮重要作用［18-21］。 APOBEC3H 結(jié)構(gòu)揭示了一種與雙鏈RNA 相互作用的異常機(jī)制［18］； TRIM71 通過降解突變體p53 抑制卵巢腫瘤發(fā)生［19］， 并且TRIM71通過調(diào)節(jié)Lin28B-let-7-HMGA2 信號(hào)傳導(dǎo)抑制腫瘤發(fā)生［20］； PABPC3 在濾泡性甲狀腺癌細(xì)胞系外顯子組中被認(rèn)為是反復(fù)突變的癌癥驅(qū)動(dòng)基因［21］；乙型肝炎病毒（HBV） 對(duì)DARS2 的上調(diào)通過miR-30e-5p、 MAPK、 NFAT5 途徑促進(jìn)肝癌發(fā)生［22］；人類心肌細(xì)胞以及發(fā)育中的斑馬魚的MTG1 沉默揭示了早期的心血管損傷［23］。 然而， 這5 種RBPs的分子機(jī)制對(duì)BC 的發(fā)生仍知之甚少， 進(jìn)一步探索其可能的機(jī)制可能是有價(jià)值的。 隨后， 本研究建立了一個(gè)列線圖來幫助更直觀地預(yù)測BC 患者1年、 2 年和3 年的OS， 并使用Kaplan-Meier Plotter檢測了5 個(gè)RBPs 編碼基因的預(yù)后價(jià)值， 結(jié)果與TCGA 隊(duì)列的預(yù)后分析結(jié)果基本一致。 提示5 個(gè)RBPs 的預(yù)后模型對(duì)BC 患者治療方案的調(diào)整具有一定的參考價(jià)值。 這可能有利于臨床醫(yī)師為BC 的治療制定臨床治療決策。

圖4 train 組（A） 和test 組（B） 單因素獨(dú)立分析（左側(cè)） 和多因素獨(dú)立分析（右側(cè)） 結(jié)果

圖5 Kaplan-Meier Plotter Server 對(duì)BC 中預(yù)后模型RBPs 預(yù)后價(jià)值的驗(yàn)證

總體而言， 本研究基于5 個(gè)RBPs 構(gòu)建的預(yù)后模型， 顯著降低測序成本， 更有利于臨床應(yīng)用。此外， 5 個(gè)RBPs 預(yù)測模型在BC 患者的生存預(yù)測中有更好的表現(xiàn)。 此外， RBPs 相關(guān)基因也顯示出重要的生物學(xué)功能， 這表明它們有可能用于臨床輔助治療。 盡管如此， 本研究仍有一些局限性：（1） 預(yù)后模型僅基于TCGA 數(shù)據(jù)庫的數(shù)據(jù)， 未在臨床患者隊(duì)列等數(shù)據(jù)庫中得到驗(yàn)證； （2） 本研究是在回顧性分析的基礎(chǔ)上設(shè)計(jì)的， 未來應(yīng)進(jìn)行前瞻性研究來驗(yàn)證結(jié)果； （3） 數(shù)據(jù)集臨床信息不夠完善， 可能會(huì)降低多元逐步Cox 回歸分析的統(tǒng)計(jì)效度和信度。

綜上所述， 本研究通過一系列生物信息學(xué)分析系統(tǒng)地探討了不同表達(dá)的RBPs 在BC 中的表達(dá)及預(yù)后價(jià)值。 這些RBPs 在BC 的發(fā)生、 發(fā)展中發(fā)揮著重要作用。 通過構(gòu)建5 個(gè)RBPs 編碼基因的預(yù)后模型， 對(duì)未來BC 治療和預(yù)后有著積極指導(dǎo)意義。 目前還鮮見關(guān)于RBPs 相關(guān)的BC 預(yù)后模型的報(bào)道， 本研究結(jié)果對(duì)揭示BC 的發(fā)病機(jī)制， 開發(fā)新的診斷思路、 尋找新的治療靶點(diǎn)以及預(yù)后分子標(biāo)志物提供一定證據(jù)。