骨質疏松hub基因篩選驗證及互作miRNA預測

張亦朦，張宇新△，田發明

骨質疏松癥是以骨量下降和骨組織微結構退化為特征的全身性骨骼疾病。成骨、破骨代謝平衡是維持骨穩態的關鍵因素，隨著年齡增長，間充質干細胞向成骨分化減少，引起骨質疏松［1］。早期診斷和干預是防止骨質疏松向椎體或髖關節骨折進展的關鍵。因而尋找敏感、特異的生物學標志物，明確其對間充質干細胞分化、代謝的調控機制已成為老年性骨質疏松癥診療亟待解決的問題。微小RNA（miRNA）是一類小片段的單鏈非編碼RNA［2］。已有研究證實miRNA廣泛參與調節細胞增殖、分化和代謝等過程，與骨質丟失等骨骼疾病密切相關［3］。但以往研究主要集中于絕經后骨質疏松癥的miRNA表達譜分析和骨質疏松與骨關節炎等疾病之間的分子機制。針對老年性骨質疏松癥診斷標志物研究較少。此外，miRNA參與細胞分化調控過程的信號通路較復雜，具體調節機制尚未明確。本研究利用生物信息分析技術篩選與老年性骨質疏松癥相關的hub基因，并研討核心miRNA參與細胞分化調控的分子機制，為老年性骨質疏松癥診斷提供潛在的生物學標志物，為臨床治療提供理論依據。

1 資料與方法

1.1 數據資料 從GEO數據庫（http://www.ncbi.nlm.nih.gov/geo）中下載GLP570平臺的GSE35956基因芯片。該數據集包含5例健康者［男1例，女4例，平均年齡（57.6±8.5）歲］和5例骨質疏松癥患者［均為女性，平均年齡（86.2±5.3）歲］的骨髓間充質干細胞樣本。通過股骨頭置換術獲得股骨頭及骨髓，從中分離出骨髓間充質干細胞。

1.2 實驗動物 選取SPF級3～4月齡（青年組）和18～20月齡（老年組）的C57BL/6小鼠各6只，2組均為雌性5只，雄性1只。實驗小鼠購自北京維通利華實驗動物技術有限公司，動物生產許可證號：SYXK（京）：2020-0004，動物飼養于華北理工大學實驗中心動物房。

1.3 主要試劑與儀器 α-MEM培養基購自武漢普諾賽生命科技有限公司；10%胎牛血清購自上海李記生物科技有限公司；胰蛋白酶-EDTA購自北京索萊寶科技有限公司；青霉素-鏈霉素溶液、逆轉錄試劑盒和熒光定量PCR（qPCR）試劑盒購自蘭杰柯科技有限公司；ARKTIK96定性PCR儀和PIKOREAL96實時熒光定量PCR儀購自美國Thermo Scientific公司。

1.4 差異表達基因篩選及富集分析 利用R語言（版本3.6.3）的Limma包對GEO數據庫下載的GSE35956基因芯片進行差異表達基因（differentially expressed genes，Degs）分析，差異表達的基因篩選條件為P＜0.05和|log2fold change|＞1.5，通過ggplot2包對差異表達的結果進行可視化處理，利用complexheatmap包繪制差異表達前20位的基因熱圖。利用DAVID 6.8（Database for Annotation，Visualization and Integrated Discovery 6.8）數據庫對得到的Degs進行基因本體（gene ontology，GO）富集和KEGG通路分析［4-5］，GO功能富集和KEGG通路分析均以P＜0.05作為篩選標準。

1.5 篩選hub基因 利用String數據庫（https://string-db.org/）［6］和Cytoscape軟件對差異基因進行網絡拓撲分析，并繪制PPI網絡圖，通過MCODE插件，以MCODE scores＞4，Max depth=100，node score cutoff=0.3和k-score=2為標準確定核心模塊。通過cytohubba內的5種算法（degree、bottleneck、radiality centrality、stress centrality、closeness centrality）來篩選hub基因［7］，并對篩選出的核心模塊及hub基因進行富集分析。使用Targetscan數據庫對得到的hub基因進行miRNA互作分析，預測與hub基因相關的miRNA及其之間的關系［8］，并通過funrich數據庫（3.1.3版本）對獲得的miRNA進行GO和KEGG富集分析。

1.6 組織細胞培養 小鼠頸椎脫臼處死后，取雙側股骨、脛骨，眼科剪剪開兩端，用α-MEM培養液（含10%胎牛血清、1%青霉素-鏈霉素）沖洗髓腔，移液槍輕柔吹打至均勻。1 000 r/min離心5 min，棄上清液后加入培養液重懸細胞，細胞懸液鋪于培養瓶中，置于37℃、5%CO2加濕培養箱內培養，待細胞生長至培養瓶壁80%～90%時，用0.25%胰蛋白酶消化，制成細胞懸液后傳代培養，培養至第1代細胞用于后續實驗。

1.7 qPCR驗證 收集老年組和青年組小鼠骨髓間充質干細胞，PBS沖洗2遍，加入1 mL Trizol試劑提取總RNA，取1μL測定其濃度，反轉錄生成cDNA后進行qPCR反應。引物序列見表1。反應體系為20μL，包含SYBR Green qPCR Mix 10μL、上下游引物（10μmol）各1μL、cDNA模板2μL和RNase free H2O 6μL。反應條件：95℃預變性30 s；95℃變性10 s，60℃退火延伸30 s，40個循環。每個樣本重復3次。采用2-ΔΔCt法計算目的基因經內參GAPDH標化后的相對表達量。

1.8 統計學方法 采用Graphpad Prism 9軟件進行數據分析，符合正態分布的計量資料以均數±標準差（±s）表示，組間比較采用t檢驗，P＜0.05為差異有統計學意義。

2 結果

2.1 Degs篩選及富集分析 在骨質疏松癥患者骨髓間充質干細胞樣本中共篩選出1 081個Degs，其中982個表達上調，99個表達下調，見圖1。GO富集分析結果顯示，Degs參與了質膜黏附分子介導的同型細胞黏附、神經嵴細胞遷移、骨化等生物過程，Degs在鈣離子結合、神經肽受體和主要組織相容性復合體（major histocompatibility complex、MHC）Ⅱ類受體活性及細胞因子活性等分子功能中富集程度較高；Degs主要富集到細胞外空間、質膜等細胞組分。KEGG通路分析表明，Degs主要集中在神經活性配體-受體相互作用、Janus激酶-信號轉導和轉錄激活因子（janus kinase-signal transducer and activator of transcription，JAK-STAT）信號通路和谷氨酸能突觸等通路中，這些通路跟骨質疏松癥密切相關，見圖2。

Tab.1 Primer sequence for qPCR表1 qPCR引物序列

2.2 蛋白互作網絡構建及hub基因篩選 將Degs導入String進行分析并利用Cytoscape構建PPI網絡圖，得到586個節點和1 262條蛋白-蛋白互作線網絡。依據蛋白互作分析結果，通過5種算法篩選出8個hub基因：瘦素（Leptin，LEP）、淋巴細胞特異性蛋白酪氨酸激酶（lymphocyte-specific protein tyrosine kinase，LCK）、WT1轉 錄 因 子（WT1 transcription factor，WT1）、SRY-Box轉 錄 因 子10（SRY-box transcription factor 10，SOX10）、整合素β2（integrin subunit beta 2，ITGB2）、白蛋白（albumin，ALB）、膽囊收縮素（cholecystokinin，CCK）和骨形態發生蛋白7（bone morphogenetic protein 7，BMP7）。利 用Targetscan數據庫鑒定出6個hub基因與105個預測miRNA存在互作關系（以2為閾值），其中，LEP、LCK、WT1、ITGB2、ALB和BMP7與105個miRNA存在調控關系，見圖3。

2.3 預測miRNA的富集分析結果6個hub基因相關的miRNA主要富集到信號轉導、細胞通訊等生物過程，在轉錄因子活性、蛋白絲氨酸/蘇氨酸激酶活性和細胞質分子功能和細胞組分中富集程度高。預測miRNA參與了蛋白多聚糖Syndecan介導的信號轉導、腫瘤壞死因子相關的凋亡誘導配體（TNFrelated apoptosis-inducing ligand、TRAIL）信號通路和酪氨酸激酶受體（receptor tyrosine kinases，erbb）信號網絡等代謝通路調控，見圖4。

Fig.1 The volcano plot and expression profiles of Degs圖1差異基因火山圖及表達水平分析

2.4 小鼠骨髓間充質干細胞原代提取和培養結果 小鼠原代細胞換液后，鏡下可見均勻鋪滿瓶底的細胞，4 d后可見細胞呈多級生長，鏡下可見散在高折射細胞為分裂期細胞，7 d后觀察約40%～50%細胞呈梭形多級生長，可見部分細胞聚集呈簇狀生長。10 d左右細胞融合率可達90%，可進行傳代培養，見圖5。

Fig.2 GO and KEGG annotations of Degs圖2差異基因的GO和KEGG富集分析

Fig.3 Key genes-miRNA interaction network diagram of differential genes圖3核心基因與miRNA互作網絡圖

Fig.4 GO and KEGG annotations of predicted miRNA圖4預測miRNA的GO和KEGG富集分析

Fig.5 The cultivation of bone marrow mesenchymal stem cells(×40)圖5小鼠骨髓間充質干細胞的培養（×40）

2.5 2組小鼠hub基因的表達水平比較qPCR結果顯示，與青年組相比，老年組骨髓間充質干細胞的BMP7、CCK、ITGB2、SOX10基 因 上 調 表 達（P＜0.01），見表2。

Tab.2 Comparison of expression levels of BMP7,CCK,ITGB2 and SOX10 genes between the two groups of mice表2 2組小鼠BMP7、CCK、ITGB2、SOX10基因表達水平比較 （n=3，±s）

Tab.2 Comparison of expression levels of BMP7,CCK,ITGB2 and SOX10 genes between the two groups of mice表2 2組小鼠BMP7、CCK、ITGB2、SOX10基因表達水平比較 （n=3，±s）

＊＊P＜0.01。

組別青年組老年組t BMP7 1.00±0.00 4.24±0.31 23.373＊＊CCK 1.00±0.00 3.93±0.38 17.340＊＊ITGB2 1.00±0.00 2.08±0.25 9.650＊＊SOX10 1.00±0.00 2.92±0.46 9.257＊＊

3 討論

據世界衛生組織統計，骨質疏松癥已成為常見疾病之一，患病率隨著年齡的增長而增加。目前，關于骨質疏松癥基因組學、表觀遺傳學等分子機制研究已成熱點。miRNA是細胞和組織衰老過程的調節因子，其表達水平與心血管疾病、阿爾茨海默病等老年疾病密切相關［9］。本研究在健康者、骨質疏松癥患者及不同月齡小鼠中均鑒定到差異miRNA，表明miRNA在老年性骨質疏松癥的診療和預后領域具有潛在價值。

GO和KEGG富集分析可以深入解析關鍵基因的調控功能和代謝通路。Rachner等［10］發現骨質疏松癥患者miRNA靶基因主要參與蛋白域特異性結合、酶結合、細胞連接、蛋白質定位和蛋白質轉運等過程。本研究發現，Degs主要參與黏附分子介導的同型細胞黏附、鈣離子結合等功能，定位于細胞外空間、質膜等位置；黏附分子介導的信號可以影響細胞的生長、分化及細胞功能。同時，KEGG通路分析結果表明，叉頭框蛋白O（fork-head box protein O，FOXO）、絲裂原活化蛋白激酶（mitogen-activated protein kinase，MAPK）等信號通路同樣受到靶基因調控，推測JAK-STAT、MAPK等信號通路失調可能是骨細胞分化失衡的重要原因，可能在骨老化過程中發揮作用。hub基因可能由于改變表達水平來影響間充質干細胞之間黏附信號，進而影響間充質干細胞的增殖、向成骨分化和凋亡。

LEP是一種由白色脂肪細胞分泌的蛋白質，在循環中可以和大腦中瘦素受體結合。目前LEP既可直接作用于骨髓間充質干細胞，促進其向成骨細胞分化［11］，還可以通過影響中樞神經系統下的神經肽Y和Y2受體調節骨穩態［12］。Bai等［13］通過KEGG通路分析發現LEP在神經活性配體-受體相互作用、JAK-STAT信號通路等通路中均有富集，其中JAKSTAT信號通路可以通過調節核因子-κB受體激活因 子 配 體（receptor activator of NF-κB ligand，RANKL）破骨細胞的分化成熟而影響骨穩態。LCK編碼的蛋白質屬于蛋白酪氨酸激酶Scr家族，該家族的蛋白具有細胞轉化和細胞內的信號傳遞的功能［14-15］。非受體酪氨酸激酶可以通過影響RANK來調節破骨細胞的分化，還可以與整合素αvβ3結合來影響細胞骨架的生成［16］。另有研究表明，ITGB2參與細胞表面黏附通路，是介導成骨細胞和破骨細胞前體接觸及黏附的重要分子［17］，主要調節骨髓基質細胞和破骨前體細胞的黏附［18］。相關研究顯示ITGB2基因缺失的小鼠體現出成骨缺陷［19］。本研究中qPCR驗證結果也支持上述觀點。

BMP7編碼蛋白質的轉化生長因子-β超家族的分泌配體，通過影響Runt相關的轉錄因子2（runt related transcription factor2，Runx2）來調節干細胞的成骨分化［20］。有研究發現，BMP7可以增加成骨活性，促進骨折愈合，并且在骨折模型骨折端的骨組織中發現BMP7的表達量增加［21］。本研究發現，老年組BMP7 mRNA表達量高于青年組，推測可能由于反饋性調節導致老年組中BMP7表達升高。SOX10屬于SOX轉錄因子家族之一，與胚胎發育、神經嵴分化和周圍神經系統發育相關，具有誘導脂肪來源的干細胞向軟骨細胞方向分化的功能［22］，并且可通過促進神經嵴細胞成骨分化來參與骨損傷后的修復［23］。ALB是代謝、血漿滲透壓的編碼蛋白［24］。WT1是一種抑癌基因，在骨髓異常增生綜合征患者中呈高表達，WT1可能通過調控細胞周期來調控髓系細胞的分化和發育［25］。目前還沒有直接證據證明ALB、WT1基因與骨質疏松有直接關系，相關分子機制還需進一步研究。

綜上所述，本研究通過GEO芯片對骨質疏松癥的各種基因進行篩選，找出6個hub基因，并預測相關的miRNA，為減低骨質疏松癥的功能障礙的研究提供了全新的視角。本研究篩選出的hub基因和miRNA有可能成為骨質疏松癥的早期診斷標志物和重要靶點，為骨質疏松癥的診斷和治療提供新的目標點和研究方向。本研究結果對骨質疏松癥后續的試驗研究提供了理論依據。