基于生物信息學分析肝細胞癌中l(wèi)ncRNA-miRNA-mRNA的調(diào)控網(wǎng)絡(luò)

張 政，盧鴻健，李明慧，高文博，田 宇，王嘉祺，張榮花，章廣玲

(1.華北理工大學臨床醫(yī)學院/河北省醫(yī)工融合精準醫(yī)療重點實驗室，河北 唐山 063000；2.華北理工大學基礎(chǔ)醫(yī)學院/河北省慢性疾病基礎(chǔ)醫(yī)學重點實驗室，河北 唐山 063210)

肝細胞癌（hepatocellular carcinoma，HCC）是世界上常見的惡性腫瘤之一，占原發(fā)性肝癌的75%～85%。GLOBOCAN 的最新統(tǒng)計數(shù)據(jù)顯示[1]，2020 年全球HCC 新發(fā)病例達91 萬，死亡病例83 萬，分別居腫瘤新發(fā)病率第6 位、死亡率第3 位；我國人口占世界的18%，但新發(fā)病例占全球的45%。HCC 的發(fā)生發(fā)展過程復(fù)雜，是一種異質(zhì)性較高的癌癥，大多數(shù)發(fā)現(xiàn)時已是中晚期，錯過了手術(shù)的最佳時間且預(yù)后較差，其5 年總生存率只有12%[2]。因此，尋找HCC 發(fā)展過程中有效的靶點和信號通路并進行相應(yīng)干預(yù)，對延緩HCC 的進展或者逆轉(zhuǎn)具有重要意義。研究表明[3]，非編碼RNA（non-coding RNA，ncRNA）包括長鏈非編碼RNA（long non-coding RNA，lncRNA）和微小RNA（microRNA，miRNA）等，其可能在HCC 的發(fā)生發(fā)展中發(fā)揮了重要的作用。miRNAs 分子由18～25個bp 組成的非編碼RNA，能夠通過堿基互補配對與靶基因的3' 端非編碼區(qū)（3'untranslated region，3'UTR）特異性結(jié)合，沉默或降解靶mRNAs 從而影響其翻譯功能。lncRNAs 分子是一種長度超過200 bp 的非編碼核苷酸，含有miRNAs 的結(jié)合位點，可以作為miRNAs 的“海綿”吸附miRNAs，從而競爭性抑制miRNAs 對靶mRNAs 的作用[4]。本研究擬通過生物信息學方法探究癌癥基因組圖譜（TCGA）數(shù)據(jù)庫中HCC 組織中差異表達的基因，并構(gòu)建lncRNAmiRNA-mRNA 調(diào)控網(wǎng)絡(luò)，為靶向診斷和治療HCC提供科學依據(jù)。

1 資料與方法

1.1 HCC 組織中差異表達基因的篩選 利用R 軟件DESeq2 包分析TCGA 數(shù)據(jù)庫中374 例HCC 組織與50 例正常肝組織中l(wèi)ncRNAs、miRNAs、mRNAs 的表達譜。使用R 軟件ggplot2 包對表達譜數(shù)據(jù)繪制火山圖，以|Log2FC|＞1（fold change 差異倍數(shù)大于2 倍）且P＜0.05 作為篩選標準。

1.2 miRNAs 轉(zhuǎn)錄因子的預(yù)測 利用FunRich 軟件預(yù)測差異表達miRNAs 的轉(zhuǎn)錄因子，根據(jù)其所調(diào)控的基因數(shù)和值進行排序，選取富集程度最顯著的前10 個轉(zhuǎn)錄因子進行展示。

1.3 基因功能（GO）和京都基因與基因組百科全書（KEGG）富集分析 利用R 軟件clusterProfiler 包分別分析miRNAs、mRNAs 各基因之間潛在的影響，并以調(diào)整后的P＜0.05 為標準，進行GO 和KEGG分析，分別預(yù)測其主要生物功能和參與的主要信號通路。

1.4 miRNA-mRNA 調(diào)控網(wǎng)絡(luò) 使用FunRich 軟件預(yù)測miRNA 的候選靶基因，結(jié)合HCC 組織中差異表達的miRNAs 和mRNAs 的FC 值，篩選具有相互作用的miRNA-mRNA 關(guān)系。

1.5 PPI 網(wǎng)絡(luò)構(gòu)建 將篩選出的基因?qū)隨TRING 數(shù)據(jù)庫，分析蛋白之間的相互作用（protein-protein interaction，PPI），并將分析結(jié)果導(dǎo)入Cytoscape 中，構(gòu)建可視化PPI 蛋白互作網(wǎng)絡(luò)圖。

1.6 lncRNA-miRNA-mRNA 調(diào)控網(wǎng)絡(luò)構(gòu)建 利用Starbase 數(shù)據(jù)庫，預(yù)測miRNAs 的上游lncRNAs，結(jié)合HCC 中差異表達lncRNAs 的FC 值，篩選具有相互作用的lncRNA-miRNA 關(guān)系對?；诠蚕淼膍iRNAs 通過Cytoscape 軟件對lncRNA-miRNAmRNA 調(diào)控網(wǎng)絡(luò)進行可視化。

2 結(jié)果

2.1 差異表達lncRNAs、miRNAs、mRNAs 的篩選 利用R 軟件DESeq2 包分析374 例HCC 組織和50 例正常肝 組織中l(wèi)ncRNAs、miRNAs、mRNAs 的表達譜，根據(jù)P＜0.05 和|Log2FC|＞1 的標準篩選篩選出HCC 組織中差異表達的lncRNAs 2850 個，其中上調(diào)2475 個，下調(diào)375 個；差異表達的miRNAs 38個，其中上調(diào)31 個，下調(diào)7 個；差異表達的mRNAs 4455 個，其中上調(diào)3199 個，下調(diào)1256 個，見圖1。

圖1 差異表達基因的火山圖

2.2 轉(zhuǎn)錄因子-miRNA 調(diào)控網(wǎng)絡(luò) 利用FunRich 軟件對38 個差異表達miRNAs 進行轉(zhuǎn)錄因子預(yù)測，共得到201 個轉(zhuǎn)錄因子，根據(jù)其所調(diào)控的基因數(shù)和P值，選取富集程度顯著的前10 個轉(zhuǎn)錄因子進行展示，見圖2，其基本信息見表1，前10 個轉(zhuǎn)錄因子分別為EGR1、SP1、NKX6-1、POU2F1、LHX3、SP4、MEF2A、ARID3A、MSX2、HOXA9。

表1 差異表達前10 名轉(zhuǎn)錄因子的基本信息

圖2 差異表達miRNAs 轉(zhuǎn)錄因子的預(yù)測

2.3 miRNAs 的GO 和KEGG 富集分析 對38 個差異表達的miRNAs 進行GO 和KEGG 功能富集分析，其生物過程主要涉及血管內(nèi)皮細胞遷移的調(diào)節(jié)、血管內(nèi)皮細胞增殖調(diào)控與新生血管生成、內(nèi)皮細胞增殖的調(diào)節(jié)、血管內(nèi)皮細胞增殖與新生血管生成；分子功能主要涉及mRNA 結(jié)合參與轉(zhuǎn)錄后基因沉默、RNA 聚合酶Ⅱ復(fù)合物結(jié)合、RNA 聚合酶核心酶結(jié)合、基礎(chǔ)轉(zhuǎn)錄機制結(jié)合；KEGG 分析主要富集于miRNAs 在癌癥中的作用，見圖3。

圖3 差異表達miRNAs 的GO 和KEGG 功能富集分析

2.4 miRNA-mRNA 調(diào)控網(wǎng)絡(luò) 通過FunRich 軟件預(yù)測差異表達的38 個miRNAs 的下游靶基因，共得到1530 個靶基因，然后與差異表達的4455 個mRNAs取交集獲得286 個候選靶基因，見圖4。通過Starbase 和MiRDB 數(shù)據(jù)庫比對分析，獲得14 302 對miRNA-mRNA，與差異表達的miRNA-mRNA 匹配后，最后篩選出交叉表達負調(diào)控的12 個miRNAs 和108 個mRNAs 作為后續(xù)研究對象。

圖4 差異表達miRNAs 的下游靶基因與差異表達mRNAs 的韋恩圖

2.5 108 個mRNAs 的蛋白PPI 調(diào)控網(wǎng)絡(luò)及富集分析利用STRING 數(shù)據(jù)庫和Cytoscape 軟件構(gòu)建篩選出的108 個基因的PPI 網(wǎng)絡(luò)，見圖5。圖中節(jié)點代表mRNAs，有相互作用關(guān)系的mRNAs 之間有連線，共有44 個節(jié)點和134 條邊；節(jié)點越大、顏色越深則介數(shù)中心度越大，邊的粗細代表連接評分，邊越粗，mRNAs 間的互作關(guān)系越大。

圖5 108 個mRNAs 的蛋白互作網(wǎng)絡(luò)圖

2.6 108 個mRNAs 的GO 與KEGG 富集分析 利用GO 功能分析篩選出的108 個mRNAs 的生物功能，生物過程主要涉及肌肉適應(yīng)、MAP 激酶活性的調(diào)節(jié)、胚胎器官發(fā)育、蛋白絲氨酸/蘇氨酸激酶活性的正調(diào)控；分子功能主要涉及MAP 激酶酪氨酸/絲氨酸/蘇氨酸磷酸酶活性、MAP 激酶磷酸酶活性、激酶調(diào)節(jié)活性、DNA 結(jié)合轉(zhuǎn)錄激活物活性和RNA 聚合酶Ⅱ特異性。KEGG 富集分析108 個mRNAs 參與的信號通路，主要涉及調(diào)節(jié)干細胞多能性的信號通路、催乳素信號通路、MAPK 信號通路、EGFR 酪氨酸激酶抑制劑耐藥性，見圖6。

圖6 108 個mRNAs 的GO 和KEGG 富集分析圖

2.7 構(gòu)建lncRNA-miRNA-mRNA 調(diào)控網(wǎng)絡(luò) 利用Starbase 數(shù)據(jù)庫，分析篩選出的12 個miRNAs 的上游lncRNAs 共221 個，結(jié)合差異表達lncRNAs 的FC 值，篩選出5 對lncRNA-miRNA，包括2 個miRNAs（miR-23A 和miR-146B）與5 個lncRNAs（GAS5、MCM3AP-AS1、LINC00173、KCNQ1OT1、SNHG7）；結(jié)合篩選出的2 個miRNAs 的與組織中差異表達mRNAs 的FC 值，共匹配42 對miRNAmRNA（包括40 個mRNAs）。通過Cytoscape 軟件對lncRNA-miRNA-mRNA 調(diào)控網(wǎng)絡(luò)進行可視化處理，見圖7。

圖7 HCC 組織中l(wèi)ncRNA-miRNA-mRNA 的調(diào)控網(wǎng)絡(luò)

3 討論

HCC 是由肝細胞惡性病變引起的癌癥，患者發(fā)現(xiàn)時大多數(shù)已是中晚期，具有發(fā)病隱匿、進展快、易耐藥、預(yù)后差以及易復(fù)發(fā)等特點，且現(xiàn)有的治療效果并不理想[5]。因此尋找開發(fā)有效的HCC 診斷和預(yù)后的生物標記物具有重要意義。越來越多的研究表明，lncRNAs 和miRNAs 等非編碼RNAs 與HCC 的發(fā)病密切相關(guān)。LncRNA n339260 通過降低miR-30e-5p的表達使TP53INP1 的表達增加，從而促進肝細胞癌的進展[6]。LncRNA ZEB1-AS1 通過調(diào)節(jié)miR-1224-5p/MAP4K4 軸，調(diào)節(jié)持續(xù)沙眼衣原體感染中線粒體介導(dǎo)的HeLa 細胞凋亡[7]。抑制lncRNA GAS6-AS2 通過調(diào)節(jié)miR-136-5p/OXSR1 軸能夠在體外和體內(nèi)水平減輕敗血癥相關(guān)的急性腎損傷[8]。然而在HCC 中尚沒有報道完整的lncRNA-miRNA-mRNA 調(diào)控網(wǎng)絡(luò)，因此深入研究lncRNAs 和miRNAs 在HCC 發(fā)生發(fā)展中的作用及調(diào)控機制具有重要意義。

已有研究表明，轉(zhuǎn)錄因子調(diào)控miRNAs 的表達，且這些轉(zhuǎn)錄因子在HCC 的發(fā)病過程中具有重要的調(diào)控作用。本研究預(yù)測了調(diào)控HCC 中差異表達miRNAs 的上游轉(zhuǎn)錄因子共201 個，根據(jù)值篩選出差異顯著的前10 個轉(zhuǎn)錄因子分別為EGR1、SP1、NKX6-1、POU2F1、LHX3、SP4、MEF2A、ARID3A、MSX2、HOXA9，其中EGR1 差異最明顯。EGR1 是一種轉(zhuǎn)錄因子，可能參與腫瘤細胞的增殖、侵襲和轉(zhuǎn)移以及腫瘤血管生成，與癌癥的發(fā)生和發(fā)展密切相關(guān)[9]。有報道表明[10]，EGR1 可能以轉(zhuǎn)錄因子方式調(diào)控hsa-miR-146b 的表達。EGR1 在miR-195 的上游啟動子區(qū)有結(jié)合位點，且調(diào)控miR-195 的表達，EGR1 可抑制miR-195/AKT3 軸并抑制胃癌細胞凋亡[11]。EGR1 還能夠通過負調(diào)控miR-491-5p 促進細胞增殖、生長和遷移，并抑制細胞凋亡，從而抑制胃癌的進程[12]。因此，研究轉(zhuǎn)錄因子對miRNAs 的調(diào)控作用可在系統(tǒng)層次上對疾病的發(fā)生提供理論依據(jù)。

近年來，伴隨著高通量測序技術(shù)的發(fā)展產(chǎn)生了大量的基因數(shù)據(jù)，這為進一步研究HCC 的發(fā)生發(fā)展及預(yù)后監(jiān)測提供了有力的基礎(chǔ)。本研究首先分析TCGA 數(shù)據(jù)庫中差異表達基因的表達譜，得到的差異表達的lncRNAs 2850 個、miRNAs 38 個、mRNAs 4455 個。對38 個差異表達的miRNAs 進行GO 和KEGG 功能富集分析，顯示主要涉及血管內(nèi)皮細胞遷移和增殖、新生血管生成、內(nèi)皮細胞增殖的調(diào)節(jié)等生物過程，mRNA 結(jié)合參與轉(zhuǎn)錄后基因沉默等分子功能，且KEGG 分析顯示這些miRNAs 主要在癌癥中發(fā)揮作用。通過FunRich 軟件預(yù)測得到的miRNAs 靶基因與差異表達的mRNAs 交叉匹配，對篩選出的108 個mRNAs 進行GO 分析發(fā)現(xiàn)，主要在胚胎器官發(fā)育、MAP 激酶活性的調(diào)節(jié)等生物學過程富集，DNA 結(jié)合轉(zhuǎn)錄激活物活性和RNA 聚合酶Ⅱ特異性、激酶調(diào)節(jié)活性等分子功能富集。KEGG 分析發(fā)現(xiàn)，靶基因參與的信號通路主要涉及調(diào)節(jié)干細胞多能性的信號通路、MAPK 信號通路等。因此可以得出在HCC 中差異表達的mRNAs 其生物功能主要與細胞增殖、形態(tài)分化和細胞遷移有關(guān)。同樣地，在PPI 網(wǎng)絡(luò)的Top 10 基因中KLF4、NRG1、FOSB、MAP2K1 也參與了上述生物學功能以及相關(guān)信號通路。研究發(fā)現(xiàn)[13]，miR-146B 通過靶向KLF4 調(diào)節(jié)肝星狀細胞活化。NRG1 是鱗狀細胞癌分化的關(guān)鍵調(diào)節(jié)因子，抑制腫瘤細胞角化和終末鱗狀分化并促進腫瘤細胞增殖[14]。TGF-β1通過誘導(dǎo)FOSB 促進前列腺癌細胞遷移和侵襲[15]。CircRNA ZFR 通過上調(diào)MAP2K1 促進HCC 的增殖[16]。

本研究尋找與miRNAs 相匹配的上游lncRNAs和下游mRNAs 構(gòu)建HCC 中l(wèi)ncRNA-miRNA-mRNA 的調(diào)控網(wǎng)絡(luò)，其中包含2 個關(guān)鍵的miRNAs（miR-23A、miR-146B）以及與其匹配的5 個lncRNAs（GAS5、MCM3AP-AS1、LINC00173、KCNQ1OT1、SNHG7）和40 個mRNAs（ROBO1、LIN28B 等）。研究表明[17]，miR-23A 是HCC 中下調(diào)最顯著的miRNA，可能是改善肝癌患者索拉非尼反應(yīng)性的潛在靶點。miR-146B 通過靶定PTP1B 抑制胃癌細胞的增殖并促進其凋亡[18]。miR-146B 通過NF-κB（核因子κB）途徑靶向TLR4（Toll 樣受體4），抑制膽囊癌腫瘤的發(fā)生和轉(zhuǎn)移[19]。以上研究表明，miR-23A、miR-146B具有抑制腫瘤的作用。此外，lncRNA GAS5 靶向miR-23A 抑制CCl4 誘導(dǎo)的肝纖維化[20]。LncRNA SNHG7 通過調(diào)節(jié)miR-146B-5p/ROBO1 軸促進胰腺癌進展[21]。LncRNA SNHG17 通過抑制miR-23A-3p 來調(diào)節(jié)CXCL12 介導(dǎo)的血管生成，從而促進結(jié)直腸腺癌細胞的增殖和遷移[22]。LncRNA KCNQ1OT1通過miR-146A-5p/ACER3 軸抑制HCC 的放射敏感性并促進腫瘤發(fā)生[23]。

綜上所述，lncRNA-miRNA-mRNA 調(diào)控網(wǎng)絡(luò)的構(gòu)建，有望闡明HCC 的發(fā)生發(fā)展機制，為HCC 的臨床診斷和預(yù)后監(jiān)測提供有價值的標記物。但是lncRNA-miRNA-mRNA 調(diào)控網(wǎng)絡(luò)的功能以及各差異基因的功能，還有待通過體內(nèi)外實驗和臨床實驗進一步證實。