基于加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析篩選結(jié)直腸癌潛在的預(yù)后生物標(biāo)志物

李梁珊 趙 虎 王詩雯

（復(fù)旦大學(xué)附屬華東醫(yī)院檢驗科，上海 200040）

結(jié)直腸癌（colorectal cancer，CRC）是消化系統(tǒng)惡性腫瘤，有約25%的CRC患者在確診時就已發(fā)生了轉(zhuǎn)移，預(yù)后不佳[1-2]。CRC的主要治療方法包括手術(shù)、放療、化療、免疫治療和靶向治療，惡性增殖、侵襲轉(zhuǎn)移和化療耐藥是CRC復(fù)發(fā)和預(yù)后不良的主要原因[3-4]。由于大多數(shù)CRC患者在確診時已處于中晚期階段，且治療敏感性普遍較低，因此對CRC的早期診斷和精準(zhǔn)治療對改善患者預(yù)后至關(guān)重要。目前，CRC的發(fā)病機(jī)制尚未完全闡明，用于CRC早期診斷和預(yù)后評估的生物標(biāo)志物、治療的分子靶標(biāo)和靶向藥物仍不足，因此迫切需要深入闡明CRC的發(fā)病機(jī)制，并篩選有效的CRC生物標(biāo)志物和分子治療靶點(diǎn)。加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析（weighted gene co-expression network analysis，WGCNA）是一種系統(tǒng)的生物學(xué)分析方法，旨在通過構(gòu)建基于基因表達(dá)矩陣的網(wǎng)絡(luò)，探索基因與臨床表型之間的關(guān)系，在篩選核心基因上展現(xiàn)出巨大的優(yōu)勢[5-6]。本研究擬采用差異表達(dá)基因篩選聯(lián)合WGCNA分析，挖掘與CRC發(fā)生、發(fā)展相關(guān)的核心基因，以期為CRC的診斷和預(yù)后提供潛在的分子標(biāo)志物和治療靶點(diǎn)。

1 材料和方法

1.1 數(shù)據(jù)收集

從基因表達(dá)綜合（Gene Expression Omnibus，GEO）數(shù)據(jù)庫（https：//www.ncbi.nlm.nih.gov/geo/）下載數(shù)據(jù)集GSE33113（芯片平臺GPL570）、GSE44076（芯片平臺GPL13667）、GSE110224（芯片平臺GPL570）和GSE17536（芯片平臺GPL570）的CEL格式原始表達(dá)譜芯片數(shù)據(jù)，分別包括90例CRC樣本和6例癌旁組織樣本、98例CRC樣本和98例癌旁組織樣本、17例CRC樣本和17例癌旁組織樣本、177例CRC樣本。

1.2 DEGs篩選

使用GSE33113數(shù)據(jù)集，采用R軟件中的affy程序包讀取芯片數(shù)據(jù)，通過RMA算法進(jìn)行標(biāo)準(zhǔn)化，通過impute.knn函數(shù)補(bǔ)充缺失值。采用R軟件中的limma程序包獲取差異表達(dá)基因，篩選條件為|log2FC|>1且P<0.05。結(jié)果由R軟件中的ggplot2程序包呈現(xiàn)。

1.3 差異表達(dá)基因富集分析

采用DAVID 6.8（https：//david-d.ncifcrf.gov/tools.jsp）在線工具對差異表達(dá)基因進(jìn)行基因本體論（Gene Ontology，GO）和京都基因與基因組數(shù)據(jù)庫（Kyoto Encyclopedia of Genes and Genomes，KEGG）通路富集分析。

1.4 WGCNA分析

采用R軟件中的WGCNA程序包對GSE33113數(shù)據(jù)集中標(biāo)準(zhǔn)差居前25%的基因構(gòu)建共表達(dá)網(wǎng)絡(luò)。剔除離群樣本，選擇合適的軟閾值（β值）構(gòu)建無尺度網(wǎng)絡(luò)。根據(jù)β值獲得鄰接矩陣，轉(zhuǎn)換為拓?fù)渲丿B矩陣（topological overlap matrix，TOM），利用TOM計算基因間的相異度（1-TOM），設(shè)定每個模塊最低基因數(shù)為30，采用動態(tài)剪切樹法將表達(dá)相似的基因分配到同一模塊中，最后將切割高度設(shè)為0.25，合并表達(dá)相似的模塊。計算每個模塊的模塊特征基因（module eigengene，ME）與臨床表型的相關(guān)系數(shù)，將相關(guān)系數(shù)較大且P<0.05的模塊定義為核心模塊。計算核心模塊中每個基因的基因顯著性（gene significance，GS）和基因模塊身份（module membership，MM）值。篩選核心模塊中MM >0.8且GS>0.2的基因為關(guān)鍵模塊基因。

1.5 關(guān)鍵基因篩選

取關(guān)鍵模塊基因與差異表達(dá)基因的交集作為候選關(guān)鍵基因，采用GEPIA數(shù)據(jù)庫進(jìn)行檢索，將在CRC組織與癌旁組織中表達(dá)有顯著性差異、且與疾病預(yù)后相關(guān)的基因作為關(guān)鍵基因。

1.6 關(guān)鍵基因驗證

在GSE33113數(shù)據(jù)集、GSE44076數(shù)據(jù)集和GSE110224數(shù)據(jù)集中篩選交集的差異表達(dá)基因。采用GEPIA數(shù)據(jù)庫驗證關(guān)鍵基因在CRC組織中的表達(dá)情況及其在預(yù)后評估中的價值。采用R軟件pROC程序包繪制受試者工作特征（receiver operating characteristic，ROC）曲線，通過曲線下面積（area under curve，AUC）判斷從GSE23878數(shù)據(jù)集中獲得的關(guān)鍵基因診斷CRC的效能。

1.7 關(guān)鍵基因集富集分析（gene set enrichment analysis，GSEA）

根據(jù)GSE17536數(shù)據(jù)集EXO1表達(dá)量的中位數(shù)將樣本分為高表達(dá)組和低表達(dá)組。使用MSigDB數(shù)據(jù)庫的c2.cp.kegg.v7.5.1.symbols.gmt [Curated]數(shù)據(jù)集和GSEA 4.2.3軟件，篩選NOM p-val<0.05、FDR q-val<0.05的富集基因集和信號通路。

1.8 不同CRC細(xì)胞中關(guān)鍵基因表達(dá)量的檢測

人正常腸上皮細(xì)胞系HIEC和CRC細(xì)胞系HCT116、HCT15均購自中國科學(xué)院細(xì)胞庫。常規(guī)培養(yǎng)細(xì)胞。采用RNA提取試劑盒（江蘇康為世紀(jì)生物科技股份有限公司）提取細(xì)胞總RNA。采用SimpliAmp PCR熱循環(huán)儀（美國ThermoFisher Scientific公司）和逆轉(zhuǎn)錄試劑盒（南京諾唯贊公司）將RNA逆轉(zhuǎn)錄為cDNA。以β-actin為內(nèi)參，采用實時熒光定量聚合酶鏈反應(yīng)（real-time fluorescence quantitative polymerase chain reaction，RTqPCR）檢測關(guān)鍵基因相對表達(dá)量，試劑盒購自南京供維贊公司，檢測儀器為QuantStudio 5 qPCR儀（美國ThermoFisher Scientific公司）。反應(yīng)條件：95 ℃ 10 min；95 ℃ 15 s，60 ℃ 30 s，共40個循環(huán)；95 ℃ 15 s，60 ℃1 min，95 ℃ 15 s。引物序列：EXO1上游引物為5'-TGAGGAAGTATAAAGGGCAGGT-3'，下游引物為5'-AGTTTTTCAGCACAAGCAATAGC-3'；β-actin上游引物為5'-TGACGTGGACATCCGCAAAG-3'，下游引物為5'-CTGGAAGGTGGACAGCGAGG-3'。引物由鉑尚生物技術(shù)（上海）有限公司合成。

采用RIPA裂解液（上海碧云天生物技術(shù)有限公司）提取細(xì)胞總蛋白，采用二喹啉甲酸（bicinchoninic acid，BCA）法測定蛋白濃度。取一定量總蛋白進(jìn)行十二烷基硫酸鈉-聚丙烯酰胺凝膠電泳，轉(zhuǎn)印至聚偏氟乙烯膜上，室溫條件下用5%脫脂奶粉封閉1 h，用1×TBST緩沖液清洗聚偏氟乙烯膜，加入相應(yīng)一抗（EXO1抗體購自美國Proteintech公司，β-actin抗體購自杭州華安生物技術(shù)有限公司），4 ℃孵育過夜。1×TBST緩沖液洗膜3次，每次10 min，加入二抗（美國Cell Signaling Technology公司），室溫孵育1 h，1×TBST緩沖液洗膜3次，每次10 min，采用ECL發(fā)光液（上海圣爾生物科技有限公司）在Tanon 5200化學(xué)發(fā)光成像儀（上海天能公司）上顯影。

1.9 統(tǒng)計學(xué)方法

采用GraphPad Prism 8.4.3軟件進(jìn)行統(tǒng)計分析和作圖。多組間比較采用單因素方差分析。以P<0.05為差異有統(tǒng)計學(xué)意義。

2 結(jié)果

2.1 差異表達(dá)基因篩選

GSE33113數(shù)據(jù)集共篩選出差異表達(dá)基因1 211個，其中表達(dá)上調(diào)505個、表達(dá)下調(diào)706個。見圖1。

圖1 GSE33113數(shù)據(jù)集DEGs火山圖

2.2 差異表達(dá)基因的GO和KEGG富集分析

差異表達(dá)基因的GO富集分析結(jié)果顯示，差異表達(dá)基因主要分布在細(xì)胞外間隙、胞外區(qū)、細(xì)胞外外泌體和細(xì)胞外基質(zhì)，涉及趨化因子介導(dǎo)的信號通路、炎癥應(yīng)答、細(xì)胞分裂、調(diào)控細(xì)胞增殖和調(diào)控細(xì)胞生長等生物過程，影響趨化因子活性、細(xì)胞外基質(zhì)結(jié)合和CXCR趨化因子受體結(jié)合等分子功能。KEGG富集分析結(jié)果顯示，差異表達(dá)基因主要富集在細(xì)胞因子與細(xì)胞因子受體的相互作用、細(xì)胞周期、過氧化物酶體增殖物激活受體（peroxisome proliferator-activated receptor，PPAR）信號通路和代謝途徑等。見圖2。

圖2 差異表達(dá)基因涉及的生物學(xué)功能

2.3 WGCNA分析

R軟件WGCNA程序包分析結(jié)果顯示，GSE33113數(shù)據(jù)集剔除1個離群樣本，樣本的聚類特征熱圖見圖3（a）。選擇軟閾值（β值）=6（r2=0.87），將基因劃分為20個模塊，在對表達(dá)相似的基因模塊進(jìn)行合并后，得到了18個模塊，其中與CRC臨床特征相關(guān)最顯著的模塊為模塊1（r=-0.81，P<0.05）和模塊2（r=0.55，P<0.05），分別有24和96個關(guān)鍵模塊基因。見圖3。

圖3 WGCNA的構(gòu)建

2.4 關(guān)鍵基因的篩選

取關(guān)鍵模塊基因與差異表達(dá)基因的交集，分別獲得24和62個候選關(guān)鍵基因。將這些基因在GEPIA數(shù)據(jù)庫上進(jìn)行檢索，發(fā)現(xiàn)模塊1中的AQP8和模塊2中的PBK、EXO1、CCNB1、DEPDC1B和KPNA2表達(dá)在CRC組織與癌旁組織之間差異有統(tǒng)計學(xué)意義（P<0.01），且與疾病預(yù)后相關(guān)（P<0.05），因此確定這6個基因為最終的關(guān)鍵基因。見圖4。

圖4 癌組織和癌旁組織6個關(guān)鍵基因相對表達(dá)量的比較

2.5 關(guān)鍵基因驗證

檢索文獻(xiàn)后發(fā)現(xiàn)CRC中EXO1的相關(guān)研究較少，因此選擇EXO1作進(jìn)一步驗證。EXO1處于GSE33113數(shù)據(jù)集、GSE44076數(shù)據(jù)集和GSE110224數(shù)據(jù)集的交集中。GEPIA數(shù)據(jù)庫分析結(jié)果顯示，CRC患者癌組織EXO1 mRNA表達(dá)水平顯著高于癌旁組織（P<0.05），且與總生存期（overall survival，OS）顯著相關(guān)（P=0.022），但與無病生存期（disease free survival，DFS）和TNM分期無關(guān)（P>0.05）。ROC曲線分析結(jié)果顯示，EXO1診斷CRC的AUC為0.913。見圖5、圖6。

圖5 關(guān)鍵基因表達(dá)CRC與OS患者的關(guān)系

圖6 關(guān)鍵基因EXO1驗證

2.6 關(guān)鍵基因EXO1的GSEA分析

GSEA分析結(jié)果顯示，在EXO1基因高表達(dá)的CRC樣本中，EXO1主要富集于細(xì)胞周期和DNA復(fù)制等基因集上。見圖7。

圖7 EXO1高表達(dá)GSEA富集分析結(jié)果

2.7 CRC細(xì)胞中EXO1的表達(dá)

CRC細(xì)胞系HCT116和HCT15中的EXO1 mRNA和蛋白表達(dá)水平均顯著高于人正常腸上皮細(xì)胞系HIEC（P<0.05）。見圖8。

圖8 CRC細(xì)胞中EXO1的表達(dá)情況

3 討論

CRC的發(fā)病機(jī)制復(fù)雜，診斷和治療難度較大，因此迫切需要尋找有意義的診斷和預(yù)后標(biāo)志物，以提高CRC患者的早期診治水平[7]。有研究發(fā)現(xiàn)，CRC正常組織樣本和腫瘤組織樣本之間的差異表達(dá)基因可能對CRC具有診斷和預(yù)后評估價值[8]。

本研究利用R語言對GSE33113數(shù)據(jù)集進(jìn)行分析，篩選出1 211個差異表達(dá)基因，其中表達(dá)上調(diào)505個、表達(dá)下調(diào)706個。KEGG通路富集分析結(jié)果顯示，差異表達(dá)基因主要參與細(xì)胞因子與細(xì)胞因子受體的相互作用、細(xì)胞周期、PPAR信號通路和代謝途徑等。細(xì)胞周期的正常調(diào)控對細(xì)胞生長至關(guān)重要，包括檢測和修復(fù)DNA損傷，阻止不受控制的細(xì)胞分裂等[9]。當(dāng)細(xì)胞周期發(fā)生紊亂時，細(xì)胞生長失去控制，進(jìn)而導(dǎo)致腫瘤的發(fā)生。有研究發(fā)現(xiàn)，PPAR-γ過表達(dá)可抑制上皮-間質(zhì)轉(zhuǎn)化過程，有助于CRC患者獲得更好的臨床結(jié)局[10]。隨后，PPAR信號通路的腫瘤抑制功能再次被證實[11]。本研究篩選到的差異表達(dá)基因可能通過影響細(xì)胞因子與細(xì)胞因子受體的相互作用、細(xì)胞周期、PPAR信號通路和代謝途徑等參與CRC的發(fā)生、發(fā)展。

本研究篩選出6個在CRC組織與癌旁組織之間表達(dá)有顯著差異，且與CRC預(yù)后相關(guān)的關(guān)鍵基因（AQP8、PBK、EXO1、CCNB1、DEPDC1B和KPNA2）。因CRC中EXO1的相關(guān)研究較少，因此選擇EXO1作進(jìn)一步分析。EXO1屬于核酸外切酶Rad2家族，具有5'雙鏈DNA核酸外切酶和5'瓣狀核酸內(nèi)切酶活性，可參與DNA錯配修復(fù)、DNA雙鏈斷裂修復(fù)、DNA復(fù)制、端粒維持、細(xì)胞周期調(diào)控和核苷酸切除修復(fù)等生物進(jìn)程[12-14]。已有研究證實，EXO1可通過泛素介導(dǎo)的蛋白酶體途徑快速降解，以響應(yīng)DNA損傷[15]。EXO1已被報道在乳腺癌、肝細(xì)胞肝癌和肺腺癌中高表達(dá)[14，16-18]，是一種潛在的生物標(biāo)志物。本研究結(jié)果顯示，EXO1表達(dá)與CRC患者OS顯著相關(guān)；ROC曲線分析結(jié)果顯示，EXO1對CRC具有較好的診斷價值；GSEA分析結(jié)果提示EXO1可能參與細(xì)胞周期和DNA復(fù)制過程，進(jìn)而調(diào)控CRC的發(fā)生、發(fā)展。RT-qPCR和免疫印跡法結(jié)果顯示，CRC細(xì)胞系HCT116和HCT15中EXO1表達(dá)顯著上調(diào)，說明EXO1可能是CRC潛在的生物標(biāo)志物，或可用于CRC患者的預(yù)后評估。

綜上所述，本研究采用WGCNA等生物信息學(xué)分析方法篩選出在CRC中高表達(dá)，且與疾病預(yù)后相關(guān)的6個關(guān)鍵基因（AQP8、PBK、EXO1、CCNB1、DEPDC1B和KPNA2），其中EXO1或可作為CRC潛在的預(yù)后評估生物標(biāo)志物。本研究結(jié)果可為進(jìn)一步探索相關(guān)基因在CRC發(fā)生、發(fā)展中的作用機(jī)制提供參考。