哈爾濱醫(yī)科大學(xué)衛(wèi)生統(tǒng)計學(xué)教研室(150081) 叢雨欣 張秋菊 田 偉 張 奇 李 稱 趙 敏 劉美娜

【提 要】 目的 構(gòu)建lncRNA-mRNA共表達網(wǎng)絡(luò)，探索與胃癌分期相關(guān)的lncRNA及其調(diào)控關(guān)系，為研究胃癌的進展機制及尋找胃癌治療的潛在靶點提供依據(jù)。方法 利用TCGA數(shù)據(jù)庫，收集胃癌RNA-seq數(shù)據(jù)及臨床信息數(shù)據(jù)；采用表達數(shù)量性狀位點(eQTL)分析與加權(quán)基因共表達網(wǎng)絡(luò)分析(WGCNA)相結(jié)合的方法，構(gòu)建lncRNA-mRNA共表達網(wǎng)絡(luò)模塊，結(jié)合臨床信息篩選與胃癌分期相關(guān)的模塊；采用Kruskal-Wallis秩和檢驗篩選模塊內(nèi)胃癌不同分期差異表達的lncRNA。結(jié)果 獲得286例胃癌組織和30例癌旁對照組織樣本的RNA-seq數(shù)據(jù)；eQTL分析得到5118對順式作用和1 953 109對反式作用lncRNA-mRNA；2 999個lncRNA和3 884個mRNA納入WGCNA，產(chǎn)生25個共表達模塊，其中與胃癌分期高度相關(guān)的模塊有3個；模塊midnightblue、orange內(nèi)18個樞紐lncRNA中有14個lncRNA在胃癌不同分期差異表達。結(jié)論 本研究篩選出14個與胃癌分期相關(guān)的lncRNA，這14個lncRNA可能通過調(diào)控mRNA的表達影響胃癌的進展，分析其對應(yīng)的網(wǎng)絡(luò)調(diào)控關(guān)系，為研究lncRNA-mRNA調(diào)控機制及探索胃癌治療靶點提供了參考和方向。

lncRNA(long non-coding RNA)是一類轉(zhuǎn)錄長度超過200nt的非編碼序列，通常由RNA聚合酶II轉(zhuǎn)錄形成。由于不編碼蛋白質(zhì)，lncRNA起初被認(rèn)為是沒有功能的轉(zhuǎn)錄垃圾。近年來發(fā)現(xiàn)，lncRNA在表觀遺傳學(xué)、轉(zhuǎn)錄及轉(zhuǎn)錄后水平等方面調(diào)控基因的表達，如干擾鄰近基因的表達、作為共因子調(diào)節(jié)轉(zhuǎn)錄因子的活性、在轉(zhuǎn)錄后水平可與mRNA形成雙鏈復(fù)合物調(diào)控基因表達水平。更值得關(guān)注的是，lncRNA在細胞異常調(diào)節(jié)、誘導(dǎo)腫瘤發(fā)生過程中起關(guān)鍵作用[1]。多項研究表明，lncRNA 與胃癌的發(fā)生、發(fā)展、轉(zhuǎn)移和預(yù)后密切相關(guān)[2-3]。因此，探索與胃癌分期相關(guān)的lncRNA，明確其對靶基因的調(diào)控關(guān)系對胃癌機制研究及治療具有重要意義。

加權(quán)基因共表達網(wǎng)絡(luò)分析(weighted gene co-expression network analysis，WGCNA)作為一種高通量數(shù)據(jù)挖掘算法，已廣泛應(yīng)用于篩選疾病相關(guān)基因的研究領(lǐng)域，如利用WGCNA構(gòu)建與胃癌進展相關(guān)的lncRNA網(wǎng)絡(luò)[4]。本研究在此基礎(chǔ)上將表達數(shù)量性狀位點(expression quantitative trait loci，eQTL)分析與WGCNA相結(jié)合，收集胃癌RNA-seq數(shù)據(jù)及臨床信息數(shù)據(jù)，利用eQTL篩選出的相關(guān)程度高的lncRNA-mRNA共表達對，對應(yīng)的lncRNA和mRNA納入WGCNA，構(gòu)建共表達網(wǎng)絡(luò)模塊，獲得與胃癌分期相關(guān)的lncRNA并分析其調(diào)控關(guān)系，為胃癌進展及治療靶點研究奠定基礎(chǔ)。

數(shù)據(jù)與方法

1.數(shù)據(jù)來源

在TCGA(the cancer genome atlas)數(shù)據(jù)庫中，納入286例胃癌組織和30例癌旁對照組織樣本的RNA-seq reads count數(shù)據(jù)。收集匹配的臨床數(shù)據(jù)，包括年齡、性別、腫瘤TNM分期等信息。

2.分析方法

(1)預(yù)處理方法

RNA-seq數(shù)據(jù)的預(yù)處理：包括提取樣本信息、構(gòu)建基因表達矩陣、將探針名轉(zhuǎn)化為基因名，最終獲得行名為樣本名、列名為基因名的矩陣；利用edge R進行差異分析，設(shè)置log2 fold change界值為1，F(xiàn)DR界值為0.05。

(2)表達數(shù)量性狀位點分析

表達數(shù)量性狀位點分析由數(shù)量性狀位點分析(quantitative trait loci，QTL)發(fā)展而來，QTL指定位控制數(shù)量性狀的基因在基因組中的位置[5]，eQTL是將每個基因的表達水平作為數(shù)量性狀進行分析[6]。通過全基因組mRNA表達量測序得到特定組織樣本的基因表達量，以全部SNP為自變量、每種mRNA表達量為因變量進行線性回歸，得到每一個SNP位點和mRNA表達量之間的關(guān)系。由于eQTL分析中表型數(shù)據(jù)不限于離散型數(shù)據(jù)，因此eQTL可作為lncRNA-mRNA共表達分析方法[7]。本研究將SNP數(shù)據(jù)替換為lncRNA表達量數(shù)據(jù)，分析mRNA與lncRNA表達量的相關(guān)性，獲得mRNA-lncRNA共表達對。

eQTL分析的優(yōu)點是可以區(qū)分順式作用關(guān)系和反式作用關(guān)系，有利于具體分析lncRNA對mRNA的調(diào)控作用。順式調(diào)控[8]是指lncRNA對染色體臨近位置(小于1Mb)的mRNA的表達調(diào)控，位于編碼蛋白上下游的lncRNA可能與啟動子或共表達基因的其他順式作用元件有交集，從而在轉(zhuǎn)錄或者轉(zhuǎn)錄后水平對基因的表達進行調(diào)控。反式調(diào)控是指對不同染色體或者染色體遠端位置的mRNA進行表達調(diào)控，調(diào)控關(guān)系不受空間距離的限定，因此trans-eQTL數(shù)量龐大，且可能存在假陽性相關(guān)?；诖?，本研究選取eQTL分析獲得的全部順式作用共表達對和前5%具有顯著性的反式作用共表達進行WGCNA分析，構(gòu)建lncRNA與mRNA共表達網(wǎng)絡(luò)。

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

加權(quán)基因共表達網(wǎng)絡(luò)分析是一種在高通量基因表達數(shù)據(jù)中，利用系統(tǒng)生物學(xué)思想，計算基因表達相關(guān)性，構(gòu)建基因共表達模塊，進而發(fā)現(xiàn)具有生物學(xué)意義模塊的高通量數(shù)據(jù)挖掘算法[9]。該算法首先計算各基因間Pearson相關(guān)系數(shù)，構(gòu)建加權(quán)鄰接矩陣：

aij=power(Sij，β)≡|Sij|β

其中aij代表基因i與基因j的鄰接系數(shù)；Sij代表基因i與基因j的Pearson相關(guān)系數(shù)。β為軟閾值，軟閾值的確定應(yīng)滿足使共表達網(wǎng)絡(luò)服從無尺度網(wǎng)絡(luò)分布，即出現(xiàn)連接度為k的節(jié)點的對數(shù)lgk與該節(jié)點出現(xiàn)的概率的對數(shù)lg[P(k)]呈負相關(guān)，且R2應(yīng)>0.8。

基因間的鄰接矩陣轉(zhuǎn)換為拓撲矩陣(topological matrix，TOM)，將某個基因與所有直接、間接相關(guān)的基因間的關(guān)系聯(lián)系起來，基于拓撲重疊性計算基因與基因間的相異度，根據(jù)相異度利用動態(tài)剪切樹法劃分基因的功能模塊，完成共表達網(wǎng)絡(luò)的基礎(chǔ)構(gòu)建。繼而計算模塊特征值(即該模塊的第一主成分)與臨床表型信息的相關(guān)系數(shù)，篩選出與胃癌分期相關(guān)的模塊作為關(guān)鍵模塊，分析關(guān)鍵模塊內(nèi)的lncRNA及其共表達的mRNA。根據(jù)網(wǎng)絡(luò)模塊節(jié)點中度的分析，確定模塊中的樞紐lncRNA，分析其在胃癌分期是否差異表達。由于目前已知功能的lncRNA微乎其微，這種分析策略有效縮小候選lncRNA的范圍，有助于尋找與胃癌分期相關(guān)的lncRNA。

(4)統(tǒng)計分析

所有分析通過R 3.6.1和Cytoscape 3.7.2實現(xiàn)。①edgeR包篩選差異表達基因；②篩選出的基因通過MatrixEQTL包進行eQTL分析，獲得lncRNA-mRNA共表達對；③利用WGCNA包實現(xiàn)共表達網(wǎng)絡(luò)的構(gòu)建并篩選與胃癌分期相關(guān)的關(guān)鍵模塊；④將模塊內(nèi)基因共表達網(wǎng)絡(luò)的權(quán)重信息導(dǎo)入Cytoscape 3.7.2軟件，篩選模塊樞紐基因并繪制基因共表達網(wǎng)絡(luò)圖；⑤對樞紐lncRNA進行Kruskal-Wallis秩和檢驗，分析各樞紐lncRNA在不同胃癌分期是否差異表達；⑥模塊內(nèi)基因映射至在線網(wǎng)站DAVID(http：//david-d.ncifcrf.gov/)中，進行GO和KEGG富集分析。

結(jié) 果

1.預(yù)處理結(jié)果

差異分析獲得4 767個差異表達的mRNA，其中2 466個mRNA上調(diào)，2 301個mRNA下調(diào)；3 542個差異表達的lncRNA，其中2 767個lncRNA上調(diào)，775個lncRNA下調(diào)。對差異基因進行火山圖的可視化，見圖1。

圖1 差異表達基因火山圖

2.構(gòu)建胃癌lncRNA-mRNA共表達網(wǎng)絡(luò)

以年齡、性別為協(xié)變量對差異基因進行eQTL分析，獲得5 118對順式作用共表達對和1 953 109對反式作用共表達對，取全部的順式作用共表達對和前5%具有顯著性的反式作用共表達對，去除重復(fù)基因，最終獲得2 999個lncRNA和3 884個mRNA進入下一步的WGCNA分析。

選取軟閾值β=4，通過動態(tài)剪切樹法進行模塊初步識別并合并相似模塊，設(shè)置每個基因網(wǎng)絡(luò)模塊最少的基因數(shù)目為50，模塊合并閾值為0.3，最終得到25個基因網(wǎng)絡(luò)模塊，灰色模塊是無法聚類到其他任何模塊的基因集合，見圖2。

圖2 加權(quán)基因共表達網(wǎng)絡(luò)模塊聚類圖

根據(jù)各模塊的特征向量分別計算每個模塊與胃癌分期的相關(guān)性，繪制模塊樣本性狀相關(guān)性熱圖，見圖3。結(jié)果顯示，midnightblue、orange、yellow、purple、blue、royalblue、green共七個模塊與胃癌分期的相關(guān)性有統(tǒng)計學(xué)意義(P<0.05)。

圖3 臨床信息與模塊相關(guān)性熱圖

3.篩選與胃癌分期相關(guān)lncRNA

選擇與胃癌分期相關(guān)性較強的前三個模塊，midnightblue、orange、yellow模塊，各模塊按度的大小排序前50的基因為樞紐基因，見圖4。三個模塊中分別包含10個、8個、3個樞紐lncRNA和40個、42個、47個樞紐mRNA。秩和檢驗結(jié)果獲得14個與胃癌分期相關(guān)lncRNA：midnightblue模塊有AP002954.4、AC002331.1、LINC01272、RP11-44K6.4、RP11-638I2.9、AC069363.1、LINC01094，orange模塊有RP11-443C10.1、RP11-13P5.2、AC093850.2、RP11-576I22.2、AC007750.5、RP11-95H3.1、RP11-867G23.1；其中AP002954.4在T3、T4期表達水平相較于T1期上調(diào)，AC002331.1、RP11-44K6.4在T4期表達水平相較于T1、T2期上調(diào)，AC069363.1在T4期表達水平相較于T1期上調(diào)，RP11-638I2.9在T4期表達水平相較于T2期上調(diào)，LINC01272在T2期表達水平相較T1期上調(diào)，LINC01094、RP11-443C10.1、RP11-13P5.2、AC093850.2、RP11-576I22.2、AC007750.5、RP11-95H3.1、RP11-867G23.1在T2、T3、T4表達水平相較于T1期上調(diào)，差異均具有統(tǒng)計學(xué)意義，見表1。與14個lncRNA共表達的mRNA共有59個，其中有56個mRNA秩和檢驗有統(tǒng)計學(xué)意義。

表1 樞紐lncRNA秩和檢驗結(jié)果

圖4 midnightblue(a)、orange(b)、yellow(c)模塊的共表達網(wǎng)絡(luò)圖

4.富集分析

lncRNA-mRNA共表達模塊中的59個mRNA進行GO和KEGG分析，結(jié)果見表2。GO分析結(jié)果顯示，模塊基因主要涉及到免疫應(yīng)答(immune response)、防御反應(yīng)(defense response)、趨化作用(chemotaxis)和炎性反應(yīng)(inflammatory response)等生物過程；KEGG分析結(jié)果顯示，模塊基因主要富集于細胞因子-受體相互作用通路(cytokine-cytokine receptor interaction)、Toll樣受體信號通路(toll-like receptor signaling pathway)、趨化因子信號通路(chemokine signaling pathway)、黏著斑通路(focal adhesion)等，這些通路與機體炎性反應(yīng)有關(guān)。

表2 富集分析結(jié)果

討 論

本研究首次將WGCNA與eQTL相結(jié)合應(yīng)用在癌癥數(shù)據(jù)中，篩選與胃癌分期相關(guān)的lncRNA。由于目前大部分lncRNA功能未知，通過分析lncRNA與mRNA表達相關(guān)性，利用已知功能的mRNA推導(dǎo)lncRNA的功能成為一種重要的研究策略。相較于既有的單純通過WGCNA構(gòu)建與胃癌進展相關(guān)的lncRNA網(wǎng)絡(luò)的研究[4]，本研究結(jié)合eQTL的優(yōu)勢在于，一方面篩選出相關(guān)程度較高的lncRNA與mRNA進行共表達網(wǎng)絡(luò)分析，提高lncRNA功能推斷的準(zhǔn)確性；另一方面將lncRNA與mRNA的共表達區(qū)分為順式共表達對與反式共表達對，有利于后續(xù)進一步研究生物學(xué)上lncRNA與mRNA的調(diào)控關(guān)系。本研究獲得14個與胃癌分期相關(guān)的lncRNA，總體上隨著胃癌分期的進展，lncRNA表達水平呈上升趨勢。與14個lncRNA共表達的59個mRNA中56個mRNA在不同分期表達水平差異有統(tǒng)計學(xué)意義，這說明eQTL與WGCNA相結(jié)合的方法能夠有效篩選出表達相關(guān)性高的lncRNA-mRNA，篩選出的lncRNA可能通過調(diào)控mRNA的表達在胃癌的發(fā)展進程中具有關(guān)鍵作用，可作為潛在的生物標(biāo)志物。研究發(fā)現(xiàn)，midnightblue模塊的樞紐lncRNA LINC01272在肺鱗癌與癌旁組織中顯著差異表達，且可以作為診斷肺鱗癌早期與晚期的標(biāo)志物[10]；orange模塊的樞紐lncRNA AC093850.2在肺鱗癌、肺腺癌、乳腺癌及胰腺癌中表達水平均高于癌旁組織[11]。目前仍未有研究發(fā)現(xiàn)這些樞紐lncRNA與胃癌分期有關(guān)，本研究結(jié)果為深入研究胃癌分期相關(guān)lncRNA提供了參考和依據(jù)。

lncRNA通過調(diào)控mRNA的轉(zhuǎn)錄、穩(wěn)定性等多種方式影響靶基因的表達，lncRNA-mRNA共表達網(wǎng)絡(luò)是研究lncRNA功能和調(diào)控機制的重要方式。GO分析和KEGG分析結(jié)果顯示lncRNA-mRNA共表達網(wǎng)絡(luò)中的mRNA主要參與免疫應(yīng)答、防御反應(yīng)、炎性反應(yīng)等生物過程，可以推斷這些樞紐lncRNA可能通過調(diào)控mRNA影響機體的炎性反應(yīng)、免疫應(yīng)答等，在一定程度影響胃癌的進展過程。已有多項研究發(fā)現(xiàn)，非編碼RNA可以通過影響腫瘤和免疫細胞中免疫調(diào)節(jié)分子的表達影響抗腫瘤免疫應(yīng)答、參與炎性反應(yīng)過程[12]。midnightblue模塊中的樞紐lncRNA LINC01272已被證實與炎性反應(yīng)有關(guān)，Wang Sen等的研究發(fā)現(xiàn)LINC01272在炎癥性腸病患者中的表達水平顯著高于健康對照，是炎癥性腸病的潛在診斷標(biāo)志物[13]。

綜上，利用eQTL和WGCNA相結(jié)合的方法構(gòu)建lncRNA-mRNA共表達網(wǎng)絡(luò)，探索與胃癌分期相關(guān)的lncRNA及其調(diào)控關(guān)系。篩選出的14個與胃癌分期相關(guān)的lncRNA可為研究胃癌的進展機制提供依據(jù)，其對應(yīng)的網(wǎng)絡(luò)調(diào)控關(guān)系為后續(xù)更深入研究lncRNA-mRNA調(diào)控機制及探索胃癌治療靶點提供了參考和方向。