關(guān)于吉西他濱耐藥在浸潤性膀胱癌上的機制研究

江詩睿，吳張松，張璐，吳松，任力杰(通信作者*)

(1.安徽醫(yī)科大學(xué)深圳二院臨床學(xué)院，廣東 深圳 518035；2.深圳大學(xué)泌尿外科研究所，廣東 深圳 518116; 3.深圳大學(xué)第三附屬醫(yī)院，廣東 深圳 518001；4.深圳大學(xué)第一附屬醫(yī)院，廣東 深圳 518035)

0 引言

膀胱癌是全球最流行的癌癥之一。依據(jù)腫瘤惡性程度，膀胱癌可分為非肌層浸潤性膀胱癌和肌層浸潤性膀胱癌。大約有 25%的患者被診斷為肌層浸潤性膀胱癌，并伴有或以后發(fā)展為該疾病的轉(zhuǎn)移形式。全身化療常用于肌層浸潤性和轉(zhuǎn)移性膀胱癌中，以分別減少疾病的遠(yuǎn)處傳播和作為姑息治療手段[1-2]。近來，各種基因共表達(dá)分析方法越來越廣泛地應(yīng)用于生物信息學(xué)?；蚬脖磉_(dá)分析是指一組表達(dá)相關(guān)且相互調(diào)節(jié)的基因，使用大量基因表達(dá)數(shù)據(jù)構(gòu)建基因間的相關(guān)性，從而挖掘基因功能的一類分析方法。在很多情況下，有著相似行為/變化的物質(zhì)，會存在著一定的聯(lián)系。在不同條件下，可以根據(jù)相似的基因表達(dá)譜，數(shù)學(xué)特征和生物學(xué)意義將基因分為不同的模塊。同一模塊中的基因被認(rèn)為具有相似的功能。在它們當(dāng)中，WGCNA是一種系統(tǒng)的基因聚類方法，用于闡明疾病發(fā)生和發(fā)展的潛在分子機制。共表達(dá)網(wǎng)絡(luò)的構(gòu)建有助于篩選生物標(biāo)志物和治療靶標(biāo)。這些方法已成功地用于各種研究中，例如頭頸鱗狀細(xì)胞癌的神經(jīng)周圍浸潤的遺傳分析。WGCNA 揭示了每個模塊中高度連接的基因的分層組織，在模塊化網(wǎng)絡(luò)中可以發(fā)現(xiàn)參與其中的關(guān)鍵hub 基因[3]。

1 資料與方法

1.1 數(shù)據(jù)來源

本研究于2019年10月從TCGA數(shù)據(jù)庫中下載膀胱癌表達(dá)譜數(shù)據(jù)，該數(shù)據(jù)包括412浸潤性膀胱癌樣本[4]。下載以RPKM呈現(xiàn)形式原始表達(dá)譜數(shù)據(jù)，用于篩選差異表達(dá)的基因及后續(xù)WGCNA分型分析。另外，對應(yīng)臨床相關(guān)資料，包括無病生存(disease free progression，DFS)數(shù)據(jù)可從網(wǎng)頁下載(http：//www.cbioportal.org)。浸潤性膀胱癌細(xì)胞系的 IC50 值來自腫瘤藥物敏感性基因組學(xué)數(shù)據(jù)庫(Genomics of Drug Sensitivity in Cancer，DSC)。

1.2 納入和排除標(biāo)準(zhǔn)

根據(jù)下載的412浸潤性膀胱癌的臨床及用藥數(shù)據(jù)，有83例接受了吉西他濱治療，而239例患者未使用該藥治療。根據(jù)患者中有無進展生存數(shù)據(jù)記錄的，排除使用吉西他濱之前就發(fā)生腫瘤進展的病例(N=15)。去除吉西他濱治療期間接受放射治療的病例(N=3)。排除診斷浸潤性膀胱癌后一年內(nèi)未接受吉西他濱治療的病例(N=5)。最終，有60例符合入組條件。最后按照DFS的中位數(shù)，分別把接受吉西他濱標(biāo)準(zhǔn)化治療組(N=60)及未接受吉西他濱治療組(N=239)分成高低DFS組進行差異分析。

1.3 GO功能富集和KEGG信號通路分析

利用R軟件的“clusterprofiler”程序包GO及KEGG功能富集，包括基因執(zhí)行的分子功能(Molecular Function，MF)，基因所處的細(xì)胞組分(Cellular Component，CC)，基因以及參與的生物學(xué)過程(Biological Process，BP)這三部分[5]。

1.4 加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析(WGCNA)

利用 R 軟件的“WGCNA” 程序包構(gòu)建加權(quán)基因共表達(dá)網(wǎng)絡(luò)[6]。首先，通過方差分析確定的 TOP 20%的變異基因(4101個基因)作為候選基因，用作后續(xù) WGCNA 分析。繪制患者的層次聚類樹以查找和刪除異常值。然后，根據(jù)標(biāo)準(zhǔn)無標(biāo)度網(wǎng)絡(luò)選擇軟閾值β=7。最后，根據(jù)每個基因之間的相關(guān)性，進行平均連鎖層次聚類，將變異基因分為不同的模塊。

1.5 對吉西他濱有特殊反應(yīng)的基因功能模塊的選擇

根據(jù)每個模塊中基因表達(dá)水平的高低，將患者分為高表達(dá)組和低表達(dá)組。根據(jù)無進展生存時間繪制生存曲線。選擇具有顯著意義的基因功能模塊(P＜0.05)作為吉西他濱特殊反應(yīng)模塊。

1.6 共表達(dá)網(wǎng)絡(luò)的可視化和 hub 基因的定義

使用Cytoscape 3.7.2軟件對2.6部分篩選目標(biāo)功能模塊行基因共表達(dá)網(wǎng)絡(luò)可視化[7]。hub基因是由Cytoscape的插件CytoHubba計算的，它提供了四種類型的基于本地的方法和七種類型的基于全局的方法[8]。每個模塊中的高度連接的基因是通過四種基于本地的方法來計算的，包括Degree，Maximum Neighborhood Component (MNC)，Density of Maximum Neighborhood Component (DMNC)，and Maximal Clique Centrality (MCC)。用以上四種方法的交集獲得的基因被鑒定為hub基因或關(guān)鍵控制基因，用韋恩圖可視化。

1.7 基于 mRNA 表達(dá)的預(yù)后標(biāo)記的構(gòu)建和準(zhǔn)確性驗證

對所選的 24 個 hub 基因進行Lasso回歸分析和單變量Cox 回歸分析。由此，我們確定了三個與DFS顯著相關(guān)的基因(P＜0.05)。使用以下公式計算簽名的風(fēng)險評分：其中 Coefi是系數(shù)，xi 是每個選定基因的 z 分?jǐn)?shù)轉(zhuǎn)換的相對表達(dá)值。利用 R 軟件中的“survival”和“survminer”程序包探索最佳風(fēng)險分值并繪制 Kaplan-Meier 生存曲線。重要的是使用“survminer”R軟件包的“surv_cutpoint”功能來將患者分為高風(fēng)險和低風(fēng)險組的最佳臨界值，并使用R包“FactoMineR”和“factoextra”行主成分分析來驗證分組結(jié)果。此外，使用Cox回歸分析評估風(fēng)險評分與DFS之間的關(guān)聯(lián)，其中年齡，階段和性別為協(xié)變量。

1.8 基因功能驗證和富集分析

從GDSC數(shù)據(jù)庫獲得基因表達(dá)數(shù)據(jù)和藥物相關(guān)信息，并探索吉西他濱的IC50與浸潤性膀胱癌細(xì)胞系中風(fēng)險基因表達(dá)的相關(guān)性[9]。進行GSEA分析來研究與浸潤性膀胱癌不同亞組相關(guān)的功能[10]。

1.9 數(shù)據(jù)分析

所有統(tǒng)計分析使用GraphPadPrism8.0和使用R軟件(version 3.5，www.r-project.org)進行。所有數(shù)據(jù)均表示為(±SD)。本研究中使用了 Student’st-test 和 Pearson’s χ2test來分析數(shù)據(jù)，當(dāng)P＜0.05具有顯著性差異。

2 結(jié)果

2.1 分組情況

在生存分析中，死亡或腫瘤轉(zhuǎn)移以及任何原因引起的復(fù)發(fā)均被視為終末事件。在無進展生存數(shù)據(jù)記錄的患者中，經(jīng)過吉西他濱治療的患者共83例，未使用該藥的患者為239例。排除使用吉西他濱之前發(fā)生腫瘤進展的病例(N=15)，去除吉西他濱治療期間接受放射治療的病例(N=3)，排除診斷后一年內(nèi)未接受吉西他濱治療的病例(N=5例)。最終，篩選出60例符合條件患者。根據(jù)DFS的中位數(shù)，將60例吉西他濱治療的患者分為兩組，在較長生存期組中有30例患者，在較短生存期組中有30例患者。其中，無病生存期相對較長的組中，有8人的生存期超過5年，有3人的生存期超過10年(圖1A)。根據(jù)中位無病生存期，將未經(jīng)吉西他濱治療的239例患者分為兩組，較長生存期組120例，較短生存期組119例。其中，無病生存期相對較長的組中，有28人的生存期超過5年，有2人的生存期超過10年(圖1B)。

圖1 吉西他濱給藥功效的差異具有其獨特的差異表達(dá)基因

圖2 加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析(WGCNA)

2.2 吉西他濱共表達(dá)網(wǎng)絡(luò)的構(gòu)建

通過R軟件的差異分析，我們發(fā)現(xiàn)了數(shù)量有限的差異表達(dá)基因。富集的KEGG途徑與免疫細(xì)胞，造血細(xì)胞和細(xì)胞因子有關(guān)。這些途徑對此問題都有一定的解釋，但是這些結(jié)果仍然相對分散。為了更好地解釋該問題，我們使用了WGCNA來構(gòu)建對吉西他濱有特殊反應(yīng)的基因共表達(dá)網(wǎng)絡(luò)。先前選擇的接受吉西他濱治療的患者的分層聚類樹(N=60)顯示，所有樣本均被聚集成一個大的聚類，而沒有需要刪除的離群值。樣本之間最變異的基因的前20%通過平均連鎖層次聚類被分組為基因功能模塊。如圖2C所示，總共鑒定出15個基因功能模塊。未分配的基因被分類為灰色模塊。在所有模塊中，最大的模塊是綠色模塊，其中包含338個基因。最小的模塊是深藍(lán)色模塊，僅包含42個基因。不同模塊基因表達(dá)的KM生存分析表明，黃色模塊(P=0.0084)和綠色模塊(P=0.032)與DFS的相關(guān)性最高。我們通過Cytoscape插件CytoHubba識別了24個hub基因(分別在綠色模塊中選擇了10個，在黃色模塊中選擇了14個)，其特征是對網(wǎng)絡(luò)中的節(jié)點進行排名。對于黃色模塊，GO和KEGG分析顯示主要途徑顯著富集了粘著斑，ECM受體相互作用，跨膜受體蛋白，細(xì)胞外結(jié)構(gòu)組織，膠原生物合成過程等。hub 基因是 COL3A1，HTRA1，THY1，ADAM12，COL5A2，ITGA11，OLFML2B，SP ARC，F(xiàn)BN1，ANGPTL2，CDH11，COL1A1，CTSK，PXDN。該模塊包含差異表達(dá)的基因 PCOLCE 和 EFEMP2，它們在圖3D 中以黑色標(biāo)記。對于綠色模塊，GO 和 KEGG 分析顯示主要途徑顯著富集于血管平滑肌收縮，肌動蛋白細(xì)胞骨架的調(diào)節(jié)，粘著斑，組織遷移，平滑肌收縮。hub 基因是 CNN1，ALDH1B1，ACTG2，KANK2，ACTA2，MYLK，TAGLN，CSRP1，TPM1，F(xiàn)LNA。該模塊包含差異表達(dá)的基因 FLNC，其在圖 4F 中以黑色標(biāo)記。

圖3 對吉西他濱的特殊應(yīng)答基因共表達(dá)網(wǎng)絡(luò)的構(gòu)建

2.3 吉西他濱在浸潤性膀胱癌治療的風(fēng)險評估的預(yù)測價值

為了更好地了解尿路上皮癌中與吉西他濱耐藥相關(guān)的潛在調(diào)節(jié)因子，我們鑒定了 24個 hub 基因。接下來，我們對24 個 hub 基因進行了單變量 Cox 回歸分析。結(jié)果表明，浸潤性膀胱癌患者中 KANK2，ITCA11，F(xiàn)BN1，PXDN，COL5A2和 THY1 的高表達(dá)具有較差的 DFS(圖 4A)。此外，我們還將 LASSO Cox 回歸算法應(yīng)用于 24 個 hub 基因，以便準(zhǔn)確識別涉及臨床結(jié)果的調(diào)節(jié)基因。根據(jù)最低標(biāo)準(zhǔn)選擇了五個基因 (ALDH1B1，COL5A2，F(xiàn)LNA，ITGA11 和 KANK2)(圖 4B，C)。最終，我們選擇 KANK2，ITCA11 和 COL5A2 進行多元Cox回歸分析，并使用從多元 Cox 回歸算法獲得的系數(shù)來計算風(fēng)險評分，以建立浸潤性膀胱癌患者吉西他濱反應(yīng)的預(yù)測模型(圖4D)，并指出 KANK2 可以作為不良 DFS 的獨立危險因素。

2.4 預(yù)后風(fēng)險評分模型與浸潤性膀胱癌臨床病理表現(xiàn)出顯著相關(guān)性

為了研究三基因風(fēng)險因子的預(yù)后模型的有效性，我們根據(jù)風(fēng)險評分中位數(shù)將浸潤性膀胱癌患者分為高風(fēng)險(吉西他濱耐藥)和低風(fēng)險(吉西他濱敏感)組，三基因預(yù)測模型風(fēng)險評分分布如圖5A 所示。此外，計算了風(fēng)險評分與浸潤性膀胱癌狀況之間的關(guān)系(圖5B)。此外，結(jié)果表明，高風(fēng)險組浸潤性膀胱癌患者的生存率較差(圖 5C)。為了判斷我們的分類是否正確，我們將通過 PCA 分析這兩個子類，結(jié)果顯示高風(fēng)險亞組可以聚集在一起，低風(fēng)險亞組也可以聚集在一起(圖5D)，并且在 3D 模型中得到驗證。接下來，風(fēng)險評分模型與其他臨床數(shù)據(jù)結(jié)合，我們進行了單因素和多因素 Cox 回歸分析，以確定風(fēng)險特征是否為獨立的預(yù)后指標(biāo)。根據(jù)證據(jù)，危險評分與 DFS 顯著相關(guān) (圖5E，F(xiàn))。

2.5 浸潤性膀胱癌中吉西他濱耐藥的三基因風(fēng)險分子特征和潛在的信號通路

我們首先顯示了高風(fēng)險和低風(fēng)險組(圖6A)中三種基因風(fēng)險分子(KANK2，ITCA11和COL5A2)的表達(dá)水平。與DFS較短的組相比，DFS較短的組患者通常更容易表達(dá)的KANK2(圖6B)。為了研究 KANK2 對膀胱癌細(xì)胞化療反應(yīng)的影響，我們利用高等級膀胱癌細(xì)胞系癌癥項目數(shù)據(jù)集(GDSC)中的藥物敏感性基因組學(xué)，分析 KANK2 表達(dá)水平與對化療敏感性分析。有趣的是，如圖6C 所示，KANK2 水平顯示出與吉西他濱的 IC50呈正相關(guān)的趨勢，但是P值不顯著。此外，為揭示吉西他濱在浸潤性膀胱癌中低應(yīng)答的潛在機制，GSEA 揭示了幾種眾所周知的信號通路可能與吉西他濱耐藥相關(guān)，包括上皮間充質(zhì)轉(zhuǎn)化 (EMT)(NES=1.95，F(xiàn)DRq＜0.001)，干細(xì)胞分化 (NES=1.47，F(xiàn)DRq=0.021)，hedgehog(Hh)信號傳導(dǎo) (NES=1.57，F(xiàn)DRq=0.031)和Wnt信號通路(NES=1.67，F(xiàn)DRq=0.006)與高危亞組顯著相關(guān)(圖6D)。所有，這些發(fā)現(xiàn)為吉西他濱在浸潤性膀胱癌中的耐藥性提供了潛在的靶點，需要將來進一步研究。

3 討論

圖5 預(yù)后風(fēng)險評分顯示浸潤性膀胱癌與臨床病理特征密切相關(guān)

我們根據(jù)來自TCGA數(shù)據(jù)庫的基因表達(dá)譜數(shù)據(jù)，通過WGCNA鑒定了與吉西他濱耐藥相關(guān)的兩個不同模塊，即綠色和黃色模塊。通過利用Lasso回歸分析和Cox回歸分析，我們還從綠色和黃色模塊中選擇了三個hub基因(例如，COL5A2，ITGA11和KANK2)得出了預(yù)后風(fēng)險特征，從而將浸潤性膀胱癌患者分為高風(fēng)險和低風(fēng)險亞組。此外，高風(fēng)險和低風(fēng)險亞組不僅與患者疾病預(yù)后相關(guān)，而且與惡性膀胱癌中吉西他濱耐藥相關(guān)的關(guān)鍵信號通路密切相關(guān)。我們當(dāng)前工作的可以開發(fā)專家決策支持系統(tǒng)，該系統(tǒng)可以特別預(yù)測患者對吉西他濱的敏感性，并指導(dǎo)臨床醫(yī)生選擇吉西他濱進行化學(xué)治療的最佳人選。此外，我們推斷靶向EMT，Hh和Wnt信號通路可以用來改善患者的吉西他濱治療結(jié)局，這顯然是未來研究的方法。這項研究是從臨床到分子水平全面分析吉西他濱在浸潤性膀胱癌患者中療效差異的因素。這是從組學(xué)角度首次驗證影響吉西他濱化療藥物治療浸潤性膀胱癌療效的各種因素。但是，考慮到所涉及患者的數(shù)量有限以及缺乏相應(yīng)實驗數(shù)據(jù)來支持我們的結(jié)論，未來需要進一步研究以證明其可行性。更精確地來講，這些治療方案可以通過首先建立患者來源的 3D 類器官和/或在異種移植模型中進行驗證，顯示對個體化給藥的持續(xù)監(jiān)測及療效驗證。

圖6 浸潤性膀胱癌中吉西他濱耐藥性涉及的三基因風(fēng)險特征和潛在的信號通路