原發(fā)性骨質(zhì)疏松癥相關(guān)核心基因的生物信息學(xué)分析

楊洋 薛建華 虞俊波 顧琪琪 張亞峰 劉佳佳*

1. 南通大學(xué)附屬醫(yī)院創(chuàng)傷中心，江蘇 南通 226001

2. 南通大學(xué)附屬醫(yī)院骨科，江蘇 南通 226001

原發(fā)性骨質(zhì)疏松癥(osteoporosis，OP)是一種以骨量減少、骨的微結(jié)構(gòu)退化，骨骼脆性增加和骨折發(fā)生風(fēng)險增加為特征的全身性骨代謝疾病[1]。OP常無明顯的臨床癥狀，一般在骨折發(fā)生后才被發(fā)現(xiàn)，主要累及脊柱、髖部、腕部和肱骨近端等部位，故OP又被稱為“無聲的浩劫”。因此，及時地預(yù)測OP對防治骨質(zhì)疏松骨折尤為重要[2]。

骨代謝的平衡由成骨細(xì)胞的骨形成作用和破骨細(xì)胞的骨吸收作用共同維持，兩者之間的動態(tài)平衡共同維持骨骼微環(huán)境的穩(wěn)態(tài)[3]。既往研究表明，過度激活的破骨細(xì)胞以及功能降低的成骨細(xì)胞是OP的重要病因。針對Wnt、NF-κB、MAPK等通路設(shè)計開發(fā)藥物，重塑骨代謝平衡，在OP的治療中具有重要意義[4-5]。然而，現(xiàn)有抗OP藥物如雙膦酸鹽、特立帕肽等，仍具有非典型性骨折、下頜骨壞死、費用昂貴等特點[6]。探究OP發(fā)病機制，尋找新的藥物靶點，在OP的臨床治療中具有重要意義[7-8]。

因此，為進一步探索OP的分子機制，本研究擬通過分析基因表達(dá)數(shù)據(jù)庫(gene expression omnibus，GEO)提供的OP相關(guān)基因芯片數(shù)據(jù)，篩選差異表達(dá)基因，并對這些基因進行核心基因和信號通路的生物信息學(xué)分析，為進一步在分子水平研究OP的發(fā)病機制提供理論依據(jù)。

1 材料和方法

1.1 一般資料

研究號 GSE35958 基因表達(dá)譜芯片數(shù)據(jù)來源于美國國立生物技術(shù)信息中心的公共基因表達(dá)數(shù)據(jù)庫(gene expression omnibus ，GEO)(https://www.ncbi.nlm.nih.gov/geo/)。該數(shù)據(jù)集由Benisch 等[9]于2012 年發(fā)表，由HGU133 Plus 2.0 ( Affymetrix公司,美國)注釋，芯片平臺為 GPL570。該數(shù)據(jù)集包含2組來自不同人群的人股骨頭骨髓間充質(zhì)干細(xì)胞(HMSC)，共計9例樣本。其中4例樣本的供體為行全髖關(guān)節(jié)置換術(shù)的非OP老年患者，另5例樣本的供體為股骨頸低能量骨折伴有OP的老年患者。

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

利用R軟件對芯片原始數(shù)據(jù)進行注釋和過濾，使用Bioconductor (http://www.bioconductor.org)提供的RMA算法對各原始芯片數(shù)據(jù)進行背景校正及歸一化等預(yù)處理，得到標(biāo)準(zhǔn)化的數(shù)據(jù)。選取校正后的P<0.05且∣logFC∣≥2為閾值，獲得該數(shù)據(jù)集的差異表達(dá)基因(differentially expressed genes，DEGs)。

1.3 GO功能和KEGG信號通路富集分析

本研究采用Metascape(http://metascape.org/gp/)在線分析工具進行差異基因富集分析。Metascape整合多個數(shù)據(jù)庫資源(GO、KEGG、UniProt和DrugBank等)，本研究使用的是基因本體( gene ontology，GO)及京都基因和基因組百科全書( kyoto encyclopedia of genes and genomes，KEGG) 兩個數(shù)據(jù)庫，分別對差異基因進行GO功能富集以分析其參與的生物學(xué)進程和KEGG 信號通路富集分析。

1.4 PPI與核心模塊分析

通過STRING(search tool for the retrieval of interacting genes /proteins) 11.0在線工具對差異基因進行蛋白質(zhì)-蛋白質(zhì)相互作用( protein-protein interaction，PPI) 網(wǎng)絡(luò)分析。將STRING的計算結(jié)果導(dǎo)入Cytoscape 3.6.1軟件，利用插件MCODE和Cytohubba分別尋找蛋白相互作用網(wǎng)絡(luò)中與OP相關(guān)的核心基因和核心模塊。

2 結(jié)果

2.1 差異基因的獲取

本研究選用基因表達(dá)譜芯片GSE35958，通過R軟件采用RMA算法對數(shù)據(jù)進行預(yù)處理后，根據(jù)篩選標(biāo)準(zhǔn)篩選出569個差異基因，其中上調(diào)7個，下調(diào)562個。通過ImageGP(http://www.ehbio.com/ImageGP/)在線工具繪制OP樣本與非OP樣本差異表達(dá)基因熱圖。熱圖的顏色反映基因表達(dá)量的水平，白色代表升高，藍(lán)色代表降低，顏色的亮度代表基因升高及降低的程度，見圖 1。

2.2 基因富集分析

將篩選出的差異基因?qū)隡etascape后，分別勾選GO中的生物學(xué)過程(biological processes)、細(xì)胞定位(cellular components)、分子功能(molecular functions)和KEGG中的信號通路(pathway)提交分析?？梢钥闯?，差異基因功能主要富集在mRNA前體內(nèi)含子結(jié)合、核體、組蛋白修飾、mRNA 3’端加工和調(diào)控結(jié)合等(圖2)，通路主要富集在葉酸生物合成、粘合連接和致病性大腸桿菌感染等(圖3)。然而，這3條信號通路中并不包含本研究篩選出的核心基因?；诖?，本研究繼續(xù)對10個核心基因再次進行KEEG通路富集分析，得到一條核心基因顯著富集的信號通路，即泛素介導(dǎo)蛋白水解信號通路(hsa04120)(圖4)。

2.3 PPI構(gòu)建及分析

利用STRING在線工具對569個差異表達(dá)基因進行 PPI構(gòu)建，設(shè)定最低相互作用分?jǐn)?shù)為 0.8，得到了由 517個節(jié)點、363條連線構(gòu)成的相互作用關(guān)系圖(圖5)。將STRING 的計算結(jié)果導(dǎo)入Cytoscape 3.6.1軟件后，采用MCC算法得到前10個核心基因，按值的大小由高到低分別是TCEB1、CUL2、KBTBD6、KBTBD7、ASB8、KLHL42、ASB5、FBXO11、ANAPC10、CDC23，且這10個基因全部為下調(diào)的差異表達(dá)基因。此外，用MCODE插件對該蛋白相互作用網(wǎng)絡(luò)進行分析，得到了由上述10個核心基因構(gòu)成的核心模塊(圖6)。

3 討論

骨質(zhì)疏松癥是一種年齡相關(guān)的慢性骨骼代謝疾病，主要以骨量減少、骨微結(jié)構(gòu)改變和骨脆性增加為特點。目前，OP的發(fā)病機制尚未完全闡明，但目前已知OP的發(fā)生具有一定的遺傳因素，與特定基因的表達(dá)有關(guān)[10-11]。

本研究分析了2組來自不同人群的人股骨頭骨髓間充質(zhì)干細(xì)胞(HMSC)，共計9例樣本，通過GEO數(shù)據(jù)庫中的基因表達(dá)譜芯片GSE35958獲得差異基因，并分析差異基因GO功能富集和KEGG通路富集，發(fā)現(xiàn)泛素介導(dǎo)蛋白水解信號通路顯著激活，為OP發(fā)病機制的研究提供了理論依據(jù)和研究方向。在本研究中，KEGG 最終篩選出的富集通路為泛素介導(dǎo)蛋白水解信號通路，該通路包含CUL2、CDC23、ANAPC10 3個核心基因。而在骨代謝過程中，破骨細(xì)胞首先粘附在骨組織，然后通過酸化和蛋白水解完成骨吸收[12]，提示該通路可能與OP有關(guān)聯(lián)。另外，就分析得到的核心基因而言，現(xiàn)有研究主要集中在腫瘤領(lǐng)域，而本研究為這些基因參與OP的發(fā)生提供了一定依據(jù)，有助于拓展核心基因的功能。

CUL2編碼的蛋白質(zhì)為Cullin蛋白家族成員，包括Elongin B和C，Rbx1和各種底物識別受體的支架蛋白，形成E3泛素連接酶[13-14]。有學(xué)者發(fā)現(xiàn)，E3泛素連接酶可以與HIFβ結(jié)合，促進HIFα的泛素化和蛋白酶降解[15-16]，而HIFα又可與HIF1β結(jié)合形成HIF轉(zhuǎn)錄因子，調(diào)節(jié)基因表達(dá)[17]。HIF下游靶基因中包括血管內(nèi)皮生長因子A(VEGFA)，可促進細(xì)胞生長和微血管的增殖[15]。以上結(jié)果提示CUL2編碼的蛋白質(zhì)可能會通過調(diào)控HIF表達(dá)以影響骨血管再生和成骨細(xì)胞活性，從而影響OP的發(fā)生發(fā)展。

TCEB1編碼的蛋白質(zhì)Elongin C為轉(zhuǎn)錄因子B(SIII)復(fù)合物的一個亞基，SIII復(fù)合物由Elongin A / A2、B和C組成。其中，Elongin B和C是調(diào)節(jié)亞基，其支架蛋白為CUL2編碼的蛋白質(zhì)。因此，TCEB1可能與上述的HIF調(diào)控下游靶基因VEGFA有關(guān)，具體結(jié)果有待進一步探討。

CDC23編碼的蛋白質(zhì)與釀酒酵母Cdc23有很強的相似性，后者為細(xì)胞周期G2/M轉(zhuǎn)換所必需的蛋白質(zhì)。該蛋白是后期促進因子(APC)的組成部分，APC由8個蛋白亞基組成，在真核細(xì)胞中高度保守，可催化細(xì)胞周期蛋白的形成。有研究發(fā)現(xiàn)，CDC23在甲狀腺癌組織中過度表達(dá)，而在正常和增生的甲狀腺組織中缺失，并且CDC23為細(xì)胞周期和細(xì)胞生長的關(guān)鍵調(diào)節(jié)因子[18]。因此，驗證 CDC23是否參與OP的骨偶聯(lián)調(diào)控值得深入研究。

ANAPC10編碼的蛋白質(zhì)為后期促進因子(APC)的一個亞基，由185個氨基酸組成，具有一個保守的核心，是細(xì)胞周期調(diào)控的泛素蛋白連接酶，在細(xì)胞從有絲分裂到G1的進程中起著至關(guān)重要的作用[19-20]。由于APC10與SP1蛋白存在功能相關(guān)，而SP1蛋白在乳腺癌、甲狀腺癌和胰腺癌等癌組織中與正常組織相比呈現(xiàn)高表達(dá)[21-22]。因此，APC10的失調(diào)可能與腫瘤的發(fā)生有關(guān)。但是，ANAPC10與OP是否存在相關(guān)性，目前尚無證據(jù)支撐，可作為研究的一個方向。

綜上所述，本研究分析得到的10個核心基因和一條信號通路為研究OP 的發(fā)生機制和藥物治療提供了新的思路，或為OP藥物研發(fā)提供作用靶點。然而，本研究所得到的結(jié)論尚需通過基礎(chǔ)實驗來作進一步驗證。