通過(guò)加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析(WGCNA)研究食管腺癌的關(guān)鍵基因

姚 歌，王新帥，張晨曦，李新陽(yáng)，王子明，王鈺鯤，張衛(wèi)國(guó)

食管癌是發(fā)生在食管上皮組織的惡性腫瘤，是全世界范圍內(nèi)癌癥相關(guān)死亡的主要原因之一，占所有惡性腫瘤的2%，死亡率高達(dá)90%，發(fā)生率和死亡率在惡性腫瘤中分別居第六位和第四位[1]。食管癌主要有兩種組織學(xué)亞型：鱗狀細(xì)胞癌和腺癌，都是高度侵襲性腫瘤。最新的流行病學(xué)數(shù)據(jù)表明，全球79%的食管鱗癌發(fā)生在東南部和中亞，而46%的食管腺癌發(fā)生于北歐、西歐、北美洲和大洋洲[2]。在過(guò)去的幾年中，發(fā)達(dá)國(guó)家腺癌發(fā)病率有所上升，而鱗癌的發(fā)病率下降。同時(shí)食管癌仍然是世界范圍內(nèi)最難治愈的惡性腫瘤之一，由于診斷時(shí)多為晚期，只有不到50%的食管癌患者有手術(shù)切除的機(jī)會(huì)。近年來(lái)隨著手術(shù)、化療和放療技術(shù)的發(fā)展，食管癌的臨床治療有了很大進(jìn)展，但5 a總生存率仍然很低。分子靶向治療是在基因分子水平上，針對(duì)已經(jīng)明確的致癌位點(diǎn)來(lái)設(shè)計(jì)相應(yīng)的治療藥物，藥物進(jìn)入體內(nèi)會(huì)特異地結(jié)合致癌位點(diǎn)發(fā)生作用，使腫瘤細(xì)胞特異性壞死，較少對(duì)正常組織影響。隨著腫瘤分子靶向治療在肺癌、結(jié)直腸癌等疾病的成功應(yīng)用，加上靶向藥物本身選擇性高、抗腫瘤活性強(qiáng)優(yōu)點(diǎn)，分子靶向治療已成為人們?cè)噲D解決食管癌治療問(wèn)題的新途徑。但目前食管癌可應(yīng)用的靶點(diǎn)較少，因此尋找新的分子靶點(diǎn)顯得尤為重要。加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析(weighted gene co-expression network analysis,WGCNA)是鑒定表型性狀的基因模塊和關(guān)鍵基因的理想方法，已有研究采用此方法對(duì)乳腺癌、肺癌、結(jié)直腸癌等惡性腫瘤進(jìn)行分析，但尚未有應(yīng)用WGCNA對(duì)食管腺癌進(jìn)行分析的研究。因此，本研究采用WGCNA算法識(shí)別與食管腺癌發(fā)展高度相關(guān)的基因模塊，并從中篩選出樞紐基因，為進(jìn)一步發(fā)現(xiàn)新的生物標(biāo)志及食管腺癌診斷和治療的潛在靶點(diǎn)提供重要的理論依據(jù)?，F(xiàn)報(bào)道如下。

1 材料與方法

1.1 數(shù)據(jù)來(lái)源及數(shù)據(jù)預(yù)處理基因表達(dá)譜數(shù)據(jù)(GSE92396)及其臨床信息來(lái)源于美國(guó)國(guó)家生物信息中心下的GEO數(shù)據(jù)庫(kù)，平臺(tái)編號(hào)為GPL6244，共納入22例樣本的基因表達(dá)數(shù)據(jù)用于分析，包含9例正常食管組織樣本、12例食管腺癌組織及1例介于兩者之間的組織樣本基因表達(dá)譜數(shù)據(jù)。在進(jìn)行WGCNA前，對(duì)下載的矩陣信息進(jìn)行了預(yù)處理：利用平臺(tái)信息將探針名轉(zhuǎn)化為基因名，將對(duì)應(yīng)多個(gè)探針的基因表達(dá)量作log2處理及標(biāo)準(zhǔn)化處理后作為該基因的最終表達(dá)量。最終獲得行名為樣本名，列名為基因名的矩陣用于后續(xù)共表達(dá)網(wǎng)絡(luò)的構(gòu)建。

1.2 數(shù)據(jù)篩選在R語(yǔ)言下運(yùn)行WGCNA包，進(jìn)行共表達(dá)網(wǎng)絡(luò)的構(gòu)建，為獲取有效的共表達(dá)網(wǎng)絡(luò)，計(jì)算各基因在所有樣本中表達(dá)量的方差，取方差前25%的基因用于共表達(dá)網(wǎng)絡(luò)的構(gòu)建。對(duì)樣本進(jìn)行聚類分析，檢測(cè)并去除離群異常值。

1.3 構(gòu)建基因共表達(dá)網(wǎng)絡(luò)通過(guò)選擇一個(gè)合適的加權(quán)系數(shù)β(軟閾值)使構(gòu)建網(wǎng)絡(luò)中的基因之間的連接服從無(wú)尺度網(wǎng)絡(luò)分布，并利用基因之間的相關(guān)系數(shù)構(gòu)建分層聚類樹(shù)，聚類樹(shù)的不同分支代表不同的基因模塊，不同顏色代表不同的模塊。然后基于基因的加權(quán)相關(guān)系數(shù)，將基因按照表達(dá)模式進(jìn)行分類，將模式相似的基因歸為一個(gè)模塊，從而將基因通過(guò)基因表達(dá)模式歸類成不同的模塊來(lái)進(jìn)行下一步分析。最后應(yīng)用此系數(shù)將相關(guān)矩陣轉(zhuǎn)換為鄰接矩陣，進(jìn)一步轉(zhuǎn)化為拓?fù)渲丿B矩陣(Topological Overlap Matrix,TOM)，并隨機(jī)選取400個(gè)基因做TOM熱圖來(lái)證明模塊之間的高度獨(dú)立性以及每個(gè)模塊中基因表達(dá)的相對(duì)獨(dú)立性。

1.4 模塊與臨床特征關(guān)聯(lián)分析通過(guò)WGCNA算法計(jì)算各模塊基因和樣本構(gòu)成的矩陣(Module Eigengene,ME)與臨床性狀的皮爾森相關(guān)系數(shù)及其P值，通過(guò)皮爾森相關(guān)系數(shù)衡量不同模塊與臨床性狀之間的關(guān)系，并取相關(guān)系數(shù)最高的模塊用于后續(xù)分析。然后利用基因顯著性(Gene Significance,GS)和模塊顯著性(Module Significance,MS)計(jì)算與樣本種類相關(guān)模塊的表達(dá)模式。GS為每個(gè)基因的表達(dá)量與某一臨床信息的皮爾森相關(guān)系數(shù)，MS為模塊中所有基因的GS取絕對(duì)值均值。為了篩選目的模塊中的樞紐基因，計(jì)算目的模塊中每一基因與該模塊ME的相關(guān)系數(shù)，即模塊隸屬度(module membership,MM)，用以衡量目的模塊內(nèi)的某一基因從屬于該模塊的程度。

1.5 樞紐基因的尋找和分析通過(guò)對(duì)各個(gè)模塊和樣本臨床信息進(jìn)行關(guān)聯(lián)分析后做出相應(yīng)熱圖，選取相關(guān)性最高的模塊，將該模塊中的基因信息導(dǎo)入到Cytoscape(3.7.1版本)軟件中，使用Cytoscape將模塊中的基因信息可視化為網(wǎng)絡(luò)，并通過(guò)連接度對(duì)模塊中基因進(jìn)行排序，篩選出前7個(gè)基因用于進(jìn)一步分析。最后把篩選出的基因用GEPIA數(shù)據(jù)庫(kù)做生存分析，觀察其對(duì)生存預(yù)后的影響。

1.6 統(tǒng)計(jì)學(xué)處理本研究使用R統(tǒng)計(jì)軟件(3.6版)和WGCNA包進(jìn)行統(tǒng)計(jì)計(jì)算。本文采用的食管組織的相關(guān)臨床特征和各個(gè)共表達(dá)模塊的ME間相關(guān)系數(shù)的計(jì)算是基于R語(yǔ)言平臺(tái)Rstuio下實(shí)現(xiàn)。WGCNA即加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析，可被應(yīng)用于尋找功能相似的基因。WGCNA使用了一種軟閾值(soft thresholding)定義一個(gè)權(quán)重值，判斷每一組基因連接的可能性。在此概念下，形成加權(quán)共表達(dá)網(wǎng)絡(luò)。WGCNA使用一個(gè)以生物學(xué)意義作為基準(zhǔn)的評(píng)判標(biāo)準(zhǔn)——無(wú)尺度拓?fù)鋵W(xué)準(zhǔn)則。

2 結(jié)果

2.1 數(shù)據(jù)來(lái)源及數(shù)據(jù)預(yù)處理從GEO數(shù)據(jù)庫(kù)中得到行名為22個(gè)樣本，列名為33 297個(gè)核酸探針的矩陣。利用平臺(tái)文件GPL6244將探針名轉(zhuǎn)換為基因名，并對(duì)其基因表達(dá)量作log2處理及標(biāo)準(zhǔn)化處理，最終得到包含個(gè)基因的樣本基因表達(dá)譜矩陣，作為共表達(dá)網(wǎng)絡(luò)構(gòu)建的輸入文件。

2.2 數(shù)據(jù)篩選本研究選取方差前25%的14 124個(gè)基因用于共表達(dá)網(wǎng)絡(luò)的構(gòu)建。通過(guò)構(gòu)建22個(gè)食管組織樣本的14 124個(gè)基因的層次聚類樹(shù)(圖1)，發(fā)現(xiàn)無(wú)明顯離群樣本。

圖1 22個(gè)食管組織樣本系統(tǒng)聚類圖

2.3 構(gòu)建基因共表達(dá)網(wǎng)絡(luò)模塊更符合無(wú)尺度特征，選取β值為4用于構(gòu)建共表達(dá)網(wǎng)絡(luò)(圖2)。通過(guò)動(dòng)態(tài)剪切樹(shù)法進(jìn)行模塊識(shí)別，最終得到45個(gè)模塊(圖3)。隨機(jī)選取400個(gè)基因做TOM熱圖(圖4)，描述基因之間的拓?fù)渲丿B矩陣，確定每個(gè)模塊獨(dú)立存在。

圖2 加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析軟閾值(β)的確定

圖3 基因聚類樹(shù)狀圖

圖4 拓?fù)渲丿B矩陣(TOM)圖

2.4 模塊與臨床特征關(guān)聯(lián)分析對(duì)各個(gè)模塊和樣本臨床信息進(jìn)行關(guān)聯(lián)分析，從模塊和性狀熱圖中可以發(fā)現(xiàn)深紅色(darkred)模塊與食管腺癌相關(guān)程度最高(圖5)。計(jì)算深紅色模塊中基因GS和MM相關(guān)系數(shù)(cor=0.77，P=5.6e-25)進(jìn)一步驗(yàn)證此結(jié)果的可信度(圖6)。結(jié)果深紅色模塊和食管腺癌的關(guān)聯(lián)程度最大，因此后續(xù)研究選取深紅色模塊篩選樞紐基因。

圖5 臨床特征與模塊特征相關(guān)性熱圖

圖6 深紅色模塊基因模塊隸屬度與基因顯著性相關(guān)的散點(diǎn)圖

圖7 矩陣網(wǎng)絡(luò)熱圖

2.5 尋找并分析樞紐基因進(jìn)一步分析發(fā)現(xiàn)深紅色模塊共有121個(gè)基因，通過(guò)WGCNA對(duì)這 121個(gè)基因進(jìn)行連通性分析后導(dǎo)入Cytoscape軟件做進(jìn)一步分析。使用Cytoscape軟件對(duì)目標(biāo)模塊的網(wǎng)絡(luò)進(jìn)行可視化(圖8)，并依據(jù)連接度選出樞紐基因(圖9)，結(jié)果顯示深紅色模塊中的樞紐基因?yàn)閂ASP、BEX4、KDELR3、POLR3B、MYO15B、FBXO27、ABCC3。最后將樞紐基因?qū)隚EPIA數(shù)據(jù)庫(kù)中做生存分析(圖10)。

圖8 深紅色模塊中基因網(wǎng)絡(luò)可視化

3 討論

食管癌是世界上最為常見(jiàn)的消化系統(tǒng)惡性腫瘤之一，中國(guó)食管癌的發(fā)病率和病死率居世界第一位[3]。食管癌惡性程度較高，且早期癥狀不典型，導(dǎo)致早期診斷困難，臨床就診患者多為中晚期，預(yù)后極差。隨著食管外科手術(shù)技術(shù)的提高和放射治療、化學(xué)治療技術(shù)的不斷創(chuàng)新，食管癌可切除率以及患者生存率有了一定提高，但化學(xué)治療藥物毒副作用較大，很多患者無(wú)法耐受，且對(duì)化學(xué)治療藥物有耐藥趨勢(shì)[4]。因此，不斷有研究嘗試食管癌生物靶向治療及免疫治療手段，探索相關(guān)分子作用的靶點(diǎn)。傳統(tǒng)的基因共表達(dá)模式的研究通常通過(guò)計(jì)算任意兩個(gè)基因之間的相關(guān)系數(shù)，設(shè)立硬閾值來(lái)衡量?jī)蓚€(gè)基因是否具有相似表達(dá)模式，缺乏說(shuō)服力。WGCNA通過(guò)在高通量數(shù)據(jù)中識(shí)別出功能相關(guān)或表達(dá)相似的基因構(gòu)成的模塊，從生物功能整體考慮基因功能及其聯(lián)系，彌補(bǔ)了傳統(tǒng)方法的缺陷。通過(guò)使用這種方法，研究者不僅可以發(fā)現(xiàn)某模塊內(nèi)基因間的聯(lián)系，還可以關(guān)注該模塊與其他模塊內(nèi)基因的相關(guān)性。此外，通過(guò)將臨床信息與模塊相關(guān)聯(lián)，還可進(jìn)一步獲得與臨床特征相關(guān)的基因，有助于建立疾病某臨床特征的研究基礎(chǔ)[5]。

本研究采用WGCNA對(duì)GEO數(shù)據(jù)庫(kù)中的22個(gè)食管組織樣本基因表達(dá)譜進(jìn)行分析，最終獲得了45個(gè)共表達(dá)模塊。將共表達(dá)模塊與臨床信息相關(guān)聯(lián)，從臨床特征與模塊特征相關(guān)性熱圖中得到深紅色模塊中的基因與食管腺癌組織相關(guān)性系數(shù)最高(0.66，P=7e-04)，且與食管正常組織相關(guān)性系數(shù)最低(-0.59，P=0.004)，深紅色模塊中的基因與食管腺癌的形成高度相關(guān)。將深紅色模塊中121個(gè)基因間的連通性信息導(dǎo)入到Cytoscape軟件中，進(jìn)一步識(shí)別出7個(gè)可能在食管腺癌中發(fā)揮關(guān)鍵作用的樞紐基因，分別為VASP、POLR3B、BEX4、KDELR3、ABCC3、FBXO27、MYO15B。VASP作為一種肌動(dòng)蛋白結(jié)合蛋白，參與一系列依賴細(xì)胞骨架重塑和細(xì)胞極性的過(guò)程，如軸突導(dǎo)向、片長(zhǎng)軸突和絲狀軸突動(dòng)力學(xué)、血小板激活和細(xì)胞遷移，VASP與絲狀肌動(dòng)蛋白的形成有關(guān)，可能在細(xì)胞黏附和運(yùn)動(dòng)中起著廣泛的作用[6]；ABCC3是MDR相關(guān)蛋白家族中的一員,能轉(zhuǎn)運(yùn)谷胱甘肽、硫酸鹽、膽汁酸和葡糖苷酸聚合物等多種內(nèi)源性和外源性化合物，被發(fā)現(xiàn)在多個(gè)惡性腫瘤中有高表達(dá)[7]；BEX4通過(guò)解除微管動(dòng)力學(xué)和染色體完整性的調(diào)控，可抵抗凋亡和細(xì)胞死亡，導(dǎo)致非整倍體的獲得，同時(shí)增加了腫瘤的增殖潛能和生長(zhǎng)[8]；POLR3B在感應(yīng)和限制細(xì)胞內(nèi)細(xì)菌和DNA病毒感染方面起著關(guān)鍵作用，在先天免疫反應(yīng)中起核質(zhì)和胞質(zhì)DNA傳感器的作用[9]；KDELR3是編碼KDEL內(nèi)質(zhì)網(wǎng)蛋白保留受體家族成員，在酵母和動(dòng)物細(xì)胞中，通過(guò)從順式高爾基體或高爾基前腔室中不斷獲取內(nèi)質(zhì)網(wǎng)腔內(nèi)的可溶性蛋白，可以實(shí)現(xiàn)內(nèi)質(zhì)網(wǎng)腔內(nèi)可溶性蛋白的保持[10]；FBXO27為SCF (SKP1-CUL1-F-box蛋白)型E3泛素連接酶復(fù)合物的底物識(shí)別組分，能夠識(shí)別和結(jié)合變性糖蛋白[11]。MYO15B在食管癌及其他腫瘤中的作用及功能尚未見(jiàn)報(bào)道，其作用需要進(jìn)一步驗(yàn)證。

圖9 樞紐基因網(wǎng)絡(luò)可視化

P<0.05，差異有統(tǒng)計(jì)學(xué)意義。圖10 利用GEPIA數(shù)據(jù)庫(kù)做樞紐基因的生存分析

將7個(gè)樞紐基因?qū)氲紾EPIA數(shù)據(jù)庫(kù)中做生存分析，發(fā)現(xiàn)VASP、ABCC3、MYO15B基因高表達(dá)患者總體生存率低于低表達(dá)患者，且P<0.05，說(shuō)明VASP、ABCC3、MYO15B可能為食管腺癌的致癌基因。BEX4基因高表達(dá)患者總體生存率高于低表達(dá)患者，BEX4可能為食管腺癌的抑癌基因。因此這些基因可能是影響食管腺癌生存和預(yù)后的關(guān)鍵基因，可能作為食管腺癌診斷和治療的潛在靶點(diǎn)，對(duì)進(jìn)一步理解食管腺癌的分子機(jī)制可起到一定幫助。需要進(jìn)一步對(duì)這些基因的作用進(jìn)行單獨(dú)的分析和實(shí)驗(yàn)驗(yàn)證。