基于高通量芯片對急性髓系白血病競爭內源性RNA網絡的構建和分析

潘玉卿,孫運艷,李軼勛,王麗英,斯南卓瑪,陳 睿,張 曦,杜 艷

（1.昆明醫科大學第一附屬醫院檢驗科 云南省實驗診斷研究所 云南省檢驗醫學重點實驗室，云南 昆明650032；2.云南省腫瘤醫院 昆明醫科大學第三附屬醫院血液科，云南 昆明650118；3.云南省腫瘤醫院 昆明醫科大學第三附屬醫院檢驗科，云南 昆明650118）

急性髓系白血病（acute myeloid leukemia，AML）是一組造血干/祖細胞惡性克隆性血液病，也是成人中最常見的白血病類型。臨床上，因白血病細胞無控性增殖引起骨髓正常造血受到抑制和骨髓外組織器官浸潤，從而引發貧血、感染、出血和肝脾腫大，甚至彌漫性血管內凝血等癥狀，危害性極大［1］。統計［2］顯示：1990—2017年，全世界白血病發病總數中AML占比由18.0%升高到23.1%。但由于治療方案的改進，AML患者的預后明顯改善。其中，小兒AML的預后明顯優于成人，且長期生存率已接近70.0%［3-4］。對AML遺傳學異常及預后的研究［5］表明：伴t（8；21）、t（15；17）和inv（16）/t（16；16）/del（16）的AML患者預后良好，但-5/del（5q）、-7/del（7q）、累及3q和11q23等復雜核型患者的預后不佳已成為國際共識［6］。因此，尋找新的AML分子標記物并全面揭示其發病機制對于患者的診斷和治療尤為重要。ceRNA網絡是由非編碼RNA（noncoding RNAs，ncRNAs）所組成并對基因發揮重要的調控機制［7］，網絡中大量ncRNAs在癌癥發病中發揮重要作用，其中長鏈非編碼RNA（long ncRNAs，lncRNAs）能夠作為強大的“海綿結構”對下游的miRNA進行招募和吸附，通過進一步抑制miRNA對mRNA的結合，從而間接抑制相關網絡中靶基因的轉錄及轉錄后調控［8］。研究［9］顯示：ceRNA能夠通過改變癌基因或抑癌基因的表達，在促進白血病細胞惡性轉化的過程中發揮重要作用。目前國內外對于AML的治療主要以放、化療及造血干細胞移植為主，但關于ceRNAs在AML發病機制中作用的報道相對較少。因此，本研究根據GEO數據庫中AML患者ncRNAs的高通量測序數據（GSE96535），采用生物信息學方法構建lncRNA-miRNA-mRNA的ceRNAs調控網絡，探討與該疾病相關的生物標記物，為揭示ceRNAs在AML發病中的作用機制提供理論基礎。

1 材料與方法

1.1 AML轉錄組數據的獲取以“acute myeloid leukemia”和“ncRNAs”為關鍵詞從GEO數據庫（http：//www.ncbi.nlm.nih.gov/geo）中篩選符合要求的數據集。篩選條件：①表達矩陣為ncRNAs；②樣本來源為AML患者的骨髓；③正常骨髓樣本為對照。選擇GPL 21242平臺的GSE96535數據集Agilent-060228 Human LncRNA microarray，選擇其中的AML患者（6例）和健康對照者（6名）的測序數據進行研究。

1.2 數據集中差異表達RNA的篩選采用SeqMap軟件［10］將AML患者和正常對照者表達矩陣中探針序列對應的lncRNA及mRNA進行重新注釋，并與GENCODE數據庫（Release 34，GRCh38.p13）中的信息進行比對，采用R語言的“edgeR”函數將得到的lncRNA及mRNA進行差異分析，并采用“pheatmap”函數繪制差異表達較大的前40個lncRNA及mRNA的熱點圖。差異表達的條件：|log2fold change（FC）|≥0.5和錯誤發現率（false discovery rate，FDR）＜0.05。

1.3 ceRNA網絡構建和功能富集分析通過在線數 據 庫miRTarBase（http：//mirtarbase.mbc.nctu.edu.tw）［11］、miRDB（http：//www.mirdb.org/miRDB/）［12］和TargetScan（http：//www.targetscan.org/）［13］預測差異表達的mRNA與差異表達的miRNA的互作關系。將3個數據庫的預測結果取交集后用于后續調控網絡的構建。采用miRcode（http：//www.mircode.org/）在線數據庫［14］預測差異表達的lncRNA和miRNA之間的互作關系。采用Cytoscape軟件（www.cytoscape.org/）構建和可視化lncRNA-miRNA-mRNA的ceRNA調控網絡，并根據基因表達譜和注釋對分子相互作用網絡進行可視化。為了探討AML的發生機制，本文作者采用R軟件中的ClusterProfile函數對ceRNA網絡中差異表達的基因進行基因本體論（Gene Ontology，GO）功能富集及京都基因和基因組百科全書（Kyoto Encyclopedia of Genes and Genomes，KEGG）通路富集分析，以FDR＜0.05作為cut-off值，得到GO富集中關鍵的生物學過程（biological process，BP）、細 胞 成 分 （cellular component，CC）和分子功能（molecular function，MF）及KEGG中關鍵的信號通路。

1.4 AML患者Hub基因的篩選首先將差異表達基因導入數據庫STRING（http：//string-db.org/）分析蛋白之間的相互作用網絡，挖掘核心的調控基因。再將得到的網絡信息導入Cytoscape軟件并進行可視化。采用Cytoscape的插件cytoHubba計算每個蛋白節點的連接度，本研究中前10個基因被確定為Hub基因［15］。

1.5 AML患者Hub基因的表達和生存分析基因表達譜交互分析數據庫（GEPIA，http：//gepia.cancer-pku.cn）可用于分析來自癌癥基因信息數據庫（TCGA，網址為https：//portal.gdc.cancer.gov/）和 基因組織表達數據庫（GTEx，網址為www.gtexportal.org/）中的腫瘤和正常樣本的RNA測序表達數據。本文作者采用GEPIA比較AML患者和正常對照者中前10個Hub基因表達，并繪制Hub基因的生存曲線。

1.6 統計學分析采用R語言軟件（版本為3.6.0，網址為www.r-project.org）進行統計學分析。采用Kaplan-Meier法和Log rank檢驗進行生存分析，數據以危險比（hazard ratio，HR）和95%可信區間（confidence interval，CI）表示，組間比較采用Mann-Whitney U檢驗。以P＜0.05為差異有統計學意義。

2 結 果

2.1 AML患者中差異mRNA和lnc RNA的表達情況以FDR＜0.05和|log2fold change（FC）|≥0.5為篩選標準，本研究共篩選出1 105個差異表達的mRNA，其中表達上調374個，表達下調731個；387個差異表達的lncRNA，其中表達上調203個，表達下調184個。數據集中差異較大的前40個mRNA和lncRNA見圖1。

2.2 構建lncRNA-miRNA-mRNA的ceRNA調控網絡通過miRcode、miRDB、miRTarBase和TargetScan數據庫預測與差異表達的lncRNA發生互作的miRNA以及差異表達的mRNA與miRNA之間的靶向調控關系，得到由45個lncRNA、31個miRNA、89個mRNA和636條連線所組成的ceRNA調控網絡，采用Cytoscape將其可視化。見圖2。

圖1 AML患者和正常對照者差異表達mRNA(A)和lncRNA(B)的熱點圖Fig.1 Heat maps of differentially expressed mRNA(A)and lnc RNA(B)in AML patients and normal controls

圖2 AML的mRNA-miRNA-lncRNA ceRNA調控網絡Fig.2 m RNA-miRNA-lnc RNA ceRNA network of AML

2.3 功能富集分析 GO功能富集結果顯示：在BP層面，差異基因主要涉及上皮細胞增殖、細胞因子調控（圖3A）；MF層面，主要涉及跨膜受體蛋白絲氨酸活性（圖3B）；CC層面，主要涉及細胞膜微區等（圖3C）。KEGG通路富集結果顯示：差異基因主要富集在癌癥的蛋白聚糖和錯誤轉錄調控等信號通路中（圖3D）。

2.4 Hub基因篩選本研究的Hub基因分別是GATA 3、WT 1、PPARG和小GTP酶蛋白家族成員RND3（RND3）、Fas配體（FASLG）、SOX4、干擾素γ（IFN-γ）、MYB、蛋氨酸氨基轉移酶2A（MAT 2A）及E2F7。見圖4。

2.5 ceRNA網絡的Hub基因表達量和生存分析GEPIA分析結果顯示：與正常對照者比較，ceRNA網絡中AML患者的GATA 3、MYB、PPARG、SOX4、E2F7和WT 1表達水平明顯升高（圖5）。生存分析結果顯示：E2F7低表達組AML患者的生存率明顯高于E2F7高表達組。見圖6。

3 討 論

目前國內外對AML的治療主要采用阿糖胞苷/蒽環霉素（“7+3”）聯合吉妥珠單抗的強化化療方案［16］，阿扎胞苷/地西他濱/維奈托克聯合低甲基化療藥物或低劑量阿糖胞苷的低強度化療方案［17］。但因化療藥物容易產生耐藥性和各種不良反應，而針對復發性難治性AML（R/R AML）的新型靶向藥對約40%的患者無法達到完全緩解，治療費用居高等問題進一步影響了患者的治療［18］。因此，探討白血病的發生機制及其特異的診療靶點對于AML的治療具有重要意義。

研究［19］表明：包括miRNAs和lncRNAs在內的ncRNA促進了AML的發生發展，并可作為AML診斷和預后的生物標志物，而lncRNAs可以作為內源性競爭性RNA，通過miRNA反應元件（MREs）與miRNA發生相互作用，從而影響miRNA誘導的下游靶基因的轉錄沉默［20］。目前，臨床上仍缺乏AML早期診斷的分子標志物，因此，本研究通過GEO數據庫中AML患者的表達矩陣，構建lncRNA-miRNA-mRNA的ceRNA網絡，旨在闡明AML的發病機制并尋找特異的治療靶點，為AML的發病機制及預后評估奠定理論基礎。本研究所構建的ceRNA由差異表達的89個mRNA、45個lncRNAs及靶向調控的31個miRNA所構成。對其中的89個mRNA進行GO和KEGG功能富集后，通過STRING數據庫得到前10個Hub基因，并采用GEPIA數據庫驗證了Hub基因的表達水平，并進一步對預后相關的E2F7基因進行了分析。

圖3 差異基因的GO富集和KEGG富集分析Fig.3 GO and KEGG enrichment analysis of differentially expressed genes

圖4 ceRNA網絡中的Hub基因Fig.4 Hub genes of ceRNA network

對ceRNA網絡中的ncRNA和Hub基因進行文獻檢索后得到部分節點促進了AML的發生發展。其中，lncRNA CRNDE在AML細胞株中高表達，下調其表達能夠抑制AML細胞增殖并促進其凋亡［21］。lncRNA LINC00239在AML來源的KG?1和HL?60細胞中異常表達，并且通過抑制阿霉素誘導的細胞凋亡從而增加AML細胞對阿霉素的耐藥性［22］。過表達lncRNA SNHG4能夠通過調控miR-10a/PTEN軸抑制AML細胞的增殖［23］。作為淋巴細胞發育至關重要的GATA 3基因［24］，EP300-ZNF384復合物的介導使急性淋巴細胞白血病獲得特異的免疫表型［25］。原癌基因MYB在AML中可以作為重要的轉錄因子促進疾病的發生發展，當使用其抑制劑莫能菌素有效降低MYB的表達后，AML細胞喪失活力，細胞分化和凋亡也受到誘導［26］。PPARG作為核受體超家族成員之一，可聯合不同的受體激動劑和化療藥物在AML治療過程中發揮重要的協同抗白血病作用［27］。lncRNA SNHG5可通過與miR-489-3p競爭性結合進而影響SOX4的表達，從而調控AML細胞的增殖和凋亡［28］。AML患者初診和鞏固治療過程中外周血轉錄因子WT 1的表達量與中、低風險AML患者預后明顯相關［29］。E2F7作為細胞周期的負性調控因子，在AML患者的原始細胞中表達升高，下調E2F7的表達能夠促進細胞分化并抑制增殖［30］。本課題組通過GEPIA數據庫統計可見：E2F7高表達患者比低表達患者預后差，提示該基因能夠作為評估AML療效的分子標記物。此外，本研究所構建的ceRNA網絡是利用生物信息學的手段對GEO數據庫的相關數據做出預測，對于基因間的交互關系和調控機制還需進一步的基礎實驗和臨床標本加以驗證。

圖5 AML癌患者中6個差異表達的基因Fig.5 Six differentially expressed genes in AML patients

圖6 Hub基因E2F7的預后分析Fig.6 Prognostic analysis of Hub gene E2F7

本研究基于GEO數據庫中AML患者的ncRNA表達芯片成功構建了lncRNA-miRNAmRNA的ceRNA網絡，其中ncRNA及靶向的Hub基因可能為AML的精準治療提供有效的理論依據和數據支持。