用生物信息學(xué)分析預(yù)測絕經(jīng)后骨質(zhì)疏松癥核心基因與互作miRNA的研究

柴毅 譚峰 樊巧玲

南京中醫(yī)藥大學(xué)，江蘇 南京 210046

骨質(zhì)疏松癥(osteoporosis，OP)是最常見的骨骼疾病，以骨量低，骨組織結(jié)構(gòu)破壞，最終導(dǎo)致骨脆性增加，骨強度下降及骨折風(fēng)險增加，易發(fā)生骨折為特征的全身性骨病[1]。在美國，每年約有150萬人次發(fā)生骨折，絕大多數(shù)都發(fā)生于絕經(jīng)后婦女[2]。絕經(jīng)后骨質(zhì)疏松癥(postmenopausal osteoporosis，PMOP)是與年齡相關(guān)的衰老性疾病，多發(fā)生于絕經(jīng)2年以上，70歲以下的婦女。PMOP的病變是一個隱性的過程，全世界約有50%的絕經(jīng)后婦女受到影響，被認為是老年人發(fā)病率最高的疾病之一[3]。PMOP已成為全球經(jīng)濟的負擔(dān)之一，積極開展對PMOP的預(yù)防和治療是公共衛(wèi)生的重要任務(wù)。

OP的診斷基于全面的病史記錄、體格檢查、骨密度測定、影像學(xué)檢查和必要的生化測定，OP的診斷主要基于雙能X線吸收檢測法(dual-energy X-ray absorptiometry，DXA)骨密度測量結(jié)果與脆性骨折[4]。然而DXA檢測也存在一些問題，例如DXA檢測成本較高，不同設(shè)備的DXA檢測結(jié)果存在差異[5-6]。挖掘與OP或PMOP病理相關(guān)的核心基因，特別是高敏感性與高特異性的基因是預(yù)防和治療該疾病性價比較高的途徑之一。早期的研究發(fā)現(xiàn)了一些與OP有關(guān)的基因。例如，細胞周期蛋白E1(Cyclin E1，CCNE1)是細胞周期的調(diào)控因子。CCNE1參與了骨代謝過程，CCNE1在PMOP B細胞中的表達呈下調(diào)趨勢[7]。又如細絲蛋白α(filamin A alpha，F(xiàn)LNA)是參與破骨細胞生成過程的關(guān)鍵因子，高表達的FLNA可促進破骨細胞生成[8]。微小核糖核酸(micro-ribonucleic acid，miRNA)是一類非編碼蛋白的小RNA，miRNA通過抑制特定靶點的mRNA調(diào)控基因表達[9]。一些miRNA亦可作為OP或PMOP的敏感標志物或治療藥物的靶點，因而針對miRNA在OP或PMOP中的機制挖掘是很有必要的[10-11]。目前，治療OP或PMOP的相關(guān)機制的探索很少深入到非編碼RNA的作用層面。此外，利用生物信息學(xué)研究方法在OP或PMOP方面有關(guān)核心基因挖掘的研究數(shù)量相對較少，且現(xiàn)有關(guān)于鑒定OP或PMOP核心基因的生物信息學(xué)研究對與核心基因相關(guān)互作的miRNA預(yù)測的報道更為稀缺。

本研究通過對基因芯片GSE57273進行生物信息學(xué)分析篩選核心基因，并預(yù)測與這些核心基因相互作用的miRNA，為PMOP建立新的科學(xué)假說以及后續(xù)更深入的研究提供依據(jù)，并為PMOP診斷以及治療藥物的研發(fā)提供較為可靠的路徑和作用靶點。

1 材料與方法

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

從GEO公共數(shù)據(jù)庫(https://www.ncbi.nlm.nih.gov/geo/)下載基因表達芯片GSE57273。該芯片包含3組PMOP藥物干預(yù)前后樣本。其所處的平臺(Platforms)為GPL4133Agilent-014850 Whole Human Genome Microarray 4x44 K G4112F(Feature Number version)。另下載Series Matrix File以便后續(xù)使用。

1.2 差異基因的數(shù)據(jù)分析

使用GEO2R(https://www.ncbi.nlm.nih.gov/geo/geo2r/)和Morpheus(https://software.broadinstitute.org/morpheus/)在線分析軟件進行差異基因(differentially expressed genes，DEGs)的甄別和篩選。GEO2R是GEO數(shù)據(jù)庫自帶的公共在線分析工具，它可將GEO數(shù)據(jù)進行復(fù)雜的R語言分析，從而呈現(xiàn)出每個基因的計算結(jié)果。經(jīng)原始數(shù)據(jù)進行成組t檢驗統(tǒng)計學(xué)分析，以adj.P<0.01和|logFC|≥3作為DEGs的篩選條件[12]。

1.3 基因功能富集分析

GO(gene ontology)是常用的分析方法，其主要功能是注釋基因或其產(chǎn)物并識別高通量基因組或轉(zhuǎn)錄組數(shù)據(jù)的特征生物學(xué)特性。GO按照生物途徑(biological process，BP)、分子功能(molecular function，MF)、細胞定位(cellular component，CC)對基因進行注釋和分類。此外，KEGG(Kyoto Encyclopedia of Genes and Genomes)數(shù)據(jù)庫可供查詢通路信息和信號通路檢索等。KEGG通路分析是另一種常用的基因功能富集分析方法。本研究應(yīng)用DAVID(Database for Annotation, Visualization and Integrated Discovery)進行在線分析提供所需的GO和KEGG生物功能富集數(shù)據(jù)。本研究使用的DAVID數(shù)據(jù)庫版本為6.8，地址為https://david.ncifcrf.gov/，由美國國立變態(tài)反應(yīng)與傳染病研究所提供研究服務(wù)。使用Fisher Exact或EASE Score統(tǒng)計方法，GO各項以P<0.05且FDR<0.05為篩選條件，KEGG各項以P<0.05為篩選條件[13]。

1.4 蛋白互作網(wǎng)絡(luò)與模型分析

STRING(Search Tool for the Retrieval of Interacting Genes)是一款可以用來呈現(xiàn)與評估蛋白互作(protein-protein interaction，PPI)的在線分析工具。STRING中的所有數(shù)據(jù)和下載文件都可以在“Creative Commons BY 4.0”許可下免費獲取。本研究將所篩選的所有DEGs植入STRING(版本10.5，https://string-db.org/)分析工具試探它們之間潛在的聯(lián)系。置信度(confidence score)≥0.4，互作最大值(maximum number of interactors)=0設(shè)為篩選條件[14]。此后，把STRING的計算結(jié)果導(dǎo)入Cytoscape(版本3.6.0)進行MCODE(Molecular Complex Detection)分析以挖掘PPI中連接最為緊密的集簇。本研究使用的MCODE版本為1.5.1，設(shè)置參數(shù)為degree=2，node score=0.2，k-core=2，max. depth=100[15]。

1.5 miRNA-基因調(diào)控網(wǎng)絡(luò)的預(yù)測

CyTargetLinker可以擴展生物調(diào)控互作網(wǎng)絡(luò)(regulatory interaction networks，RegINs)，由荷蘭系統(tǒng)生物學(xué)聯(lián)合會提供支持。它涵蓋了miRNA—靶點、轉(zhuǎn)錄因子—靶點和藥物—靶點之間的互作關(guān)系。本研究下載了人類物種基因數(shù)據(jù)集(https://projects.bigcat.unimaas.nl/cytargetlinker/regins/)。選用該數(shù)據(jù)集中基于實驗驗證的miRTarBase 4.4數(shù)據(jù)庫(含20 942個RegINs)，基于預(yù)測功能的TargetScan 6.2數(shù)據(jù)庫(含511 040個RegINs)和MicroCosm 5數(shù)據(jù)庫(含541 039個RegINs)預(yù)測核心基因與miRNA之間的調(diào)控關(guān)系。

2 結(jié)果

2.1 DEGs的鑒定

本研究選用基因表達芯片GSE57273，經(jīng)GEO2R初步分析共獲得32 996個DEGs，隨后由Morpheus分析并經(jīng)條件篩選，最終獲得841個DEGs，其中包含826個下調(diào)基因和15個上調(diào)基因。

2.2 GO功能與KEGG信號通路富集分析

根據(jù)GO的分析結(jié)果，本研究以P<0.05，F(xiàn)DR<0.05為篩選條件，并按照計數(shù)值從大到小排列，在BP、CC與MF類別中各選取前3項列為表1?？梢钥闯?，在生物學(xué)過程中，這些DEGs主要參與了基因表達，細胞大分子生物合成和RNA代謝過程；在細胞定位中，這些DEGs富集于核漿、細胞質(zhì)基質(zhì)以及粘附連接；從分子功能上看，這些DEGs具有使有機環(huán)狀化合物結(jié)合、雜環(huán)化合物結(jié)合和核酸結(jié)合的作用。根據(jù)KEGG分析結(jié)果，本研究以P<0.05為篩選條件，按照計數(shù)值從大到小將DEGs富集的信號通路列為表2。結(jié)果顯示，這些DEGs主要富集于癌癥信號通路、病毒致癌通路、粘附斑、rap1信號通路和內(nèi)質(zhì)網(wǎng)蛋白加工通路。

表1 與PMOP相關(guān)的DEGs的GO富集分析

表2 與PMOP相關(guān)的DEGs的KEGG富集分析

圖2 蛋白質(zhì)互作網(wǎng)絡(luò)的前3個集簇模塊

2.3 核心基因的鑒定與篩選的集簇模塊

通過STRING的PPI構(gòu)建，經(jīng)由Cytoscape對網(wǎng)絡(luò)的計算工具得出所有DEGs的連接度(degree)(圖1)。degree值表示網(wǎng)絡(luò)中某一基因與周圍基因的關(guān)系數(shù)量，因此degree越大代表與它相互作用關(guān)系的基因數(shù)量就越多。本研究按degree從高至低進行排序，以排名前10位的高degree基因定為核心基因，它們分別是HSP90AA1(degree=75)、EP300(degree=55)、SMARCA2(degree=44)、RANBP2(degree=41)、ASH1L(degree=36)、EIF4E(degree=35)、PTEN(degree=31)、CNOT6L(degree=30)、RPL7(degree=29)、KRAS(degree=29)。此外，MCODE分析發(fā)現(xiàn)17組集簇，共包含523個節(jié)點(node)與2 026條連線(edge)。本研究以score為依據(jù)展示前3組集簇并分析各個集簇所富集的通路(圖2)。其中集簇模塊A富集核糖體(hsa03010：Ribosome)與mRNA監(jiān)視通路(hsa03015：mRNA surveillance pathway)(P<0.05)，集簇模塊B富集泛素介導(dǎo)的蛋白質(zhì)水解4通路(hsa04120：Ubiquitin mediated proteolysis 4)(P<0.05)，而集簇模塊C未鑒定出具有統(tǒng)計學(xué)意義的通路(P>0.05)(表3)。

圖1 DEGs的PPI網(wǎng)絡(luò)

2.4 miRNA—核心基因調(diào)控網(wǎng)絡(luò)模塊的鑒定

本研究應(yīng)用CyTargetLinker預(yù)測可以與上述篩選出的10個核心基因相互作用的miRNA。結(jié)果顯示，在MicroCosm數(shù)據(jù)庫中有258個預(yù)測的miRNA靶點互作關(guān)系，在TargetScan數(shù)據(jù)庫中有1 171個預(yù)測的miRNA靶點互作關(guān)系，總共有875個節(jié)點與1 429條連線。另外，閾值(threshold)可對結(jié)果顯示的可視化調(diào)控網(wǎng)絡(luò)進行支持數(shù)據(jù)庫的疊加篩選，通過調(diào)設(shè)閾值可控制調(diào)控網(wǎng)絡(luò)的顯示結(jié)果，其設(shè)置范圍一般為1～3[16]。本研究將閾值設(shè)為2，結(jié)果顯示共有37個miRNA與7個靶基因存在互作關(guān)系。這些基因與預(yù)測的miRNA如表3所示。

3 討論

本研究通過對GEO數(shù)據(jù)庫中的基因芯片GSE57273進行生物信息學(xué)分析獲得DEGs，并分析了有關(guān)這些基因富集的生物過程、細胞定位、分子功能和信號通路為機制研究提供了理論依據(jù)與研究方向，通過挖掘與核心基因互作的miRNA為PMOP的研究提供新思路。

GO與KEGG分析有助于更深入地認識并篩選出DEGs的功能和作用。由于受到數(shù)據(jù)呈現(xiàn)空間的限制，本研究未能把這些DEGs所富集的生物學(xué)過程與信號通路全部列出。本研究由KEGG篩選出的富集通路以癌癥通路為首。從參與的DEGs來看，該通路僅包含了少部分核心基因，由于該通路涉及了大部分非核心基因，這可能導(dǎo)致它們富集的通路與PMOP不同。另外，就僅涉及的核心基因而言，現(xiàn)有關(guān)于這些基因功能的研究也存在局限性，大部分集中于腫瘤領(lǐng)域的研究，而本研究結(jié)果為這些核心基因參與PMOP提供了一定依據(jù)，有助于拓展核心基因的功能。

表3 由CyTargetLinker擴展網(wǎng)絡(luò)分析預(yù)測的與7個核心基因互作的miRNA

本研究羅列的核心基因主要富集于RNA代謝過程、基因表達過程、PI3K-Akt信號通路和癌癥相關(guān)的信號通路等。HSP90AA1編碼的蛋白質(zhì)是一種功能類似于同型二聚體的誘導(dǎo)型分子。HSP90AA1參與細胞的生長發(fā)育過程，有研究證實使用HSP90AA1抑制劑可迅速導(dǎo)致細胞死亡，表明HSP90AA1在細胞活動中發(fā)揮重要作用[17-18]。因而，驗證HSP90AA1是否參與PMOP的骨偶聯(lián)的調(diào)控過程值得深入探究。EP300編碼與腺病毒E1A相關(guān)的細胞p300轉(zhuǎn)錄共激活蛋白。它與組蛋白乙酰轉(zhuǎn)移酶的功能類似，可以通過染色質(zhì)重塑調(diào)節(jié)轉(zhuǎn)錄并且在細胞增殖和分化過程發(fā)揮重要作用。EP300在骨髓中高表達，并參與了細胞成骨分化和骨量減少的調(diào)控過程。EP300在髓核細胞中可受到骨形態(tài)發(fā)生蛋白2和骨形態(tài)發(fā)生蛋白7的調(diào)控，為EP300與PMOP的關(guān)聯(lián)提供理論支撐[19-20]。由SMARCA2基因編碼的蛋白屬于SWI/SNF家族蛋白，并且這種蛋白與果蠅的brahma蛋白高度相似。該家族蛋白具有解旋酶和ATP酶活性，其通過改變基因周圍的染色質(zhì)結(jié)構(gòu)來發(fā)揮調(diào)節(jié)基因轉(zhuǎn)錄的功能。SMARCA2在卵巢、大腦等處高表達，而這類蛋白的減少會影響間充質(zhì)干細胞的成骨分化和成脂分化平衡，進而導(dǎo)致OP[21]。RANBP2編碼與核孔復(fù)合物免疫定位有關(guān)的RAN結(jié)合蛋白，主要在睪丸、甲狀腺、骨髓和大腦等處高表達。RAN是與核膜相關(guān)的RAS超家族的小GTP結(jié)合蛋白，它通過與蛋白質(zhì)的相互作用來調(diào)控多種細胞功能。RANBP2細胞定位為核孔復(fù)合體，并且在神經(jīng)視網(wǎng)膜中的表達非常豐富[22]。ASH1L編碼轉(zhuǎn)錄激活因子的三空腔結(jié)構(gòu)蛋白質(zhì)，它同樣在睪丸、甲狀腺、骨髓和大腦等處高表達。現(xiàn)有研究證實ASH1L可以參與細胞分化以及骨髓造血的調(diào)控[23-24]。EIF4E基因編碼的蛋白質(zhì)是真核生物翻譯起始因子4F復(fù)合物的組成部分。研究發(fā)現(xiàn)，EIF4E在病理狀態(tài)下的骨髓間充質(zhì)干細胞(bone marrow stromal cells，BMSCs)中表達水平會發(fā)生變化[25]。當然，EIF4E作為一種原癌基因，其表達和激活與轉(zhuǎn)化和腫瘤發(fā)生密切相關(guān)。PTEN編碼的蛋白質(zhì)是磷脂酰肌醇-3，4，5-三磷酸3-磷酸酶，該基因被認為是在大量癌癥中高頻率突變的腫瘤抑制因子。近年來有研究已經(jīng)發(fā)現(xiàn)PTEN/PI3K/AKT信號通路可能是調(diào)控BMSCs增殖和分化的途徑之一[26]。CNOT6L可能參與包括細胞增殖的各種細胞活動，有研究證實CNOT6L可通過p53介導(dǎo)的信號通路途徑來調(diào)節(jié)細胞周期阻滯和衰老[27-28]。RPL7編碼一種由60S亞基組成的核糖體蛋白，該蛋白屬于核糖體蛋白L30P家族，主要在卵巢與骨髓高表達。KRAS編碼一種屬于小GTP酶超家族成員的蛋白。有研究證實KRAS可以與骨形態(tài)發(fā)生蛋白4靶向性結(jié)合[29]。

此外，本研究還擴展了與核心基因相互作用的miRNA。與HSP90AA1相互調(diào)控的miRNA參與了類固醇生物合成過程。與EP300相互作用的miRNA參與D-谷氨酰胺和D-谷氨酸代謝和細胞周期調(diào)控過程。與SMARCA2靶向調(diào)控的miRNA富集于調(diào)節(jié)干細胞多能性的信號通路、肌動蛋白細胞骨架的調(diào)控與甲狀腺激素信號通路等。與RANBP2互作的miRNA功能富集于蛋白的內(nèi)質(zhì)網(wǎng)加工過程和雌激素信號通路等。hsa-miR-137、hsa-miR-219-5p、hsa-miR-548 d-3p可調(diào)控ASH1L，參與了細胞氮化合物代謝過程。hsa-miR-134、hsa-miR-935可調(diào)控KRAS，目前尚未發(fā)現(xiàn)它們所富集的信號通路，因此仍需進一步研究探索。與EIF4E互作miRNA的功能富集于細胞周期調(diào)控和糖胺聚糖的生物合成等。

本研究預(yù)測的37個miRNA與核心基因為研究PMOP提供了新的角度。然而，本研究所得到的結(jié)論尚需要更多的研究驗證這些核心基因和miRNA在PMOP病理進展中的特異性與可靠性。隨著未來針對這些核心基因及miRNA不斷地深入研究，它們或可應(yīng)用于PMOP病理的篩查，或為PMOP藥物研發(fā)提供作用靶點，或有助于對PMOP合并腫瘤患者進行基因診斷。此外，上述核心基因及miRNA或可作為遺傳致病因子的篩查點并為此定制個性化的防治措施。由于PMOP機制復(fù)雜，影響因素眾多，探尋并梳理具有邏輯關(guān)系和實際意義的基因及互作集簇，找出后續(xù)研究工作的切入點是需要進一步思考的問題。