基于生物信息學(xué)分析探討類風(fēng)濕關(guān)節(jié)炎與骨關(guān)節(jié)炎之間的關(guān)系①

章曉云 李華南 陳 鋒 柴 源 甘 斌 陳丁鵬 （江西中醫(yī)藥大學(xué)，南昌 330004）

類風(fēng)濕關(guān)節(jié)炎（rheumatoid arthritis，RA）是一種自身免疫性疾病，可引起慢性炎癥表現(xiàn)、肥胖、心血管疾病等，影響全球0.3%～1.0%的人口，在中國每年至少500 萬人被診斷出RA［1-2］。骨關(guān)節(jié)炎（osteoarthritis，OA）是一系列能引起軟骨損傷和變性的退行性關(guān)節(jié)炎的總稱，影響了全球約2.4億人口，無論是從生活質(zhì)量還是疾病負(fù)擔(dān)方面，都已成為全球范圍內(nèi)嚴(yán)重的公共衛(wèi)生問題［3-5］。兩者是目前最常見的關(guān)節(jié)炎類型，病因均未明確，暫無有效的根治方法，給個人和社會帶來了巨大的經(jīng)濟(jì)壓力。

盡管兩種疾病的病因?qū)W不同，但均以滑膜炎癥反應(yīng)、關(guān)節(jié)軟骨破壞、骨質(zhì)流失等為特征，最終導(dǎo)致關(guān)節(jié)進(jìn)行性退變，臨床上許多患者同時患有這兩種疾病。研究表明RA 是OA 發(fā)病的重要危險因素之一，RA 患者罹患OA 的風(fēng)險比正常人高2.75 倍［4］。隨著基礎(chǔ)研究的不斷深入，已經(jīng)發(fā)現(xiàn)兩種疾病之間存在諸多共同發(fā)病機制和受調(diào)控因子，在蛋白質(zhì)水平和組織炎癥評分水平上具有很強的連續(xù)性［6-8］。因此，本研究通過GEO數(shù)據(jù)庫查找RA和OA血清基因芯片表達(dá)譜進(jìn)行生物信息學(xué)分析，篩選出差異miRNA 并進(jìn)行靶基因預(yù)測，最后交集篩選出核心mRNA，為確定RA 和OA 之間的共同發(fā)病機理創(chuàng)造條件，并為RA 和OA 治療藥物的研發(fā)提供有效的作用靶點。

1 資料與方法

1.1 數(shù)據(jù)收集 以“Rheumatoid arthritis”“Osteoar?thritis”“miRNA”為關(guān)鍵詞，在GEO 公共數(shù)據(jù)庫（https：//www.ncbi.nlm.nih.gov/geo/）［9］檢索與RA 及OA 相關(guān)的芯片，獲得編號為GSE115885 的RA 血清基因表達(dá)芯片，共包含15 例RA 患者和5 例健康對照，其所處平臺為GPL25134；獲得編號為GSE105027 的OA 血清基因表達(dá)芯片，共含12 例OA患者和12例健康對照，其所處平臺為GPL21575。

1.2 差異miRNA 的數(shù)據(jù)分析 利用Perl（V5.30.0.1）對數(shù)據(jù)進(jìn)行重注釋，并利用R 語言校正、分類，limma 包對miRNA 進(jìn)行差異分析，以P<0.05 和|logFC|≥0.5 作為過濾條件，篩選差異miRNA（differentially expressed miRNAs，DEmiRNAs）。最后運用Venny 平臺（http：//bioinfogp. cnb. csic. es/tools/venny/）將所得RA 和OA 的DEmiRNAs 進(jìn)行映射取交集，得到交集DEmiRNAs。

1.3 miRNA 轉(zhuǎn)錄因子預(yù)測及功能富集分析 Fun?Rich（V3.1.3）［10］是一個用于基因和蛋白質(zhì)功能富集和相互作用網(wǎng)絡(luò)分析的軟件。利用FunRich 軟件將篩選出的交集miRNA 進(jìn)行轉(zhuǎn)錄因子預(yù)測及功能富集分析，并根據(jù)P值顯示前10位進(jìn)行可視化。

1.4 靶基因預(yù)測 通過miRTarBase 數(shù)據(jù)庫（http：//mirtarbase. mbc. nctu. edu. tw/php/index. php）、Target Scan 數(shù)據(jù)庫（http：//www. targetscan. org/vert_72/）和miRDB數(shù)據(jù)庫（http：//mirdb.org/）分別預(yù)測DEmiRNA的靶基因，使用Venny 平臺取3 個數(shù)據(jù)庫均能預(yù)測到的靶基因，生成交集基因，并與步驟1.2所得的交集DEmiRNAs 構(gòu)成關(guān)系網(wǎng)絡(luò)。整合得到miRNAmRNA 調(diào)控網(wǎng)絡(luò)，并導(dǎo)入Cytoscape（V3.7.2）繪制網(wǎng)絡(luò)圖。

1.5 蛋白互作網(wǎng)絡(luò)構(gòu)建 為進(jìn)一步探究miRNA 所調(diào)控靶基因的作用機制，將靶基因?qū)隨tring 數(shù)據(jù)庫（http：//string-db. org/cgi/input. pl）［11］，限定研究物種為“Homo Sapiens”獲得蛋白互作關(guān)系，設(shè)置連接評分（Combined score）>0.95，再導(dǎo)出蛋白互作關(guān)系的數(shù)據(jù)文件，選取文件中的node1、node2 和Com?bined score 三列導(dǎo)入Cytoscape 中，最后利用Cyto?scape 軟件的“NetworkAnalyzer”工具進(jìn)行可視化處理，構(gòu)建蛋白互作（PPI）網(wǎng)絡(luò)圖。

1.6 基因富集分析 利用DAVID 數(shù)據(jù)庫（https：//david. ncifcrf. gov/）［12］對mRNA 進(jìn)行GO 功能富集分析，研究靶基因的主要生物功能；對潛在靶點進(jìn)行KEGG 信號通路富集分析，研究靶基因的主要信號通路，P<0.05 代表富集結(jié)果顯著。最后利用R 語言繪制GO富集分析圖及KEGG富集分析氣泡圖。

2 結(jié)果

2.1 DEmiRNA 的鑒定 利用R語言對芯片進(jìn)行重注釋和數(shù)據(jù)校正后，再對其進(jìn)行差異分析。結(jié)果顯示，與健康對照相比，RA 患者的血清中共存在41個明顯改變的miRNAs，其中上調(diào)31個，下調(diào)10個；OA患者的血清中共存在265個明顯改變的miRNAs，其中上調(diào)124 個，下調(diào)141 個；對二者取交集后共得到4個交集DEmiRNAs，詳見表1。

表1 交集DEmiRNAsTab.1 Intersection DEmiRNAs

2.2 miRNA 調(diào)控的轉(zhuǎn)錄因子及富集分析結(jié)果 通過FunRich 軟件對4 個交集DEmiRNA 進(jìn)行轉(zhuǎn)錄因子預(yù)測，共得到197個轉(zhuǎn)錄因子，根據(jù)其所調(diào)控的基因數(shù)與P值，選取富集程度顯著的10 個轉(zhuǎn)錄因子進(jìn)行展示，見圖1。差異最顯著的10 個轉(zhuǎn)錄因子分別為LHX3、SP4、NFIC、VSX2、HOXA7、TCF3、MYC、HOXB4、ETS1、SP1，其基本信息見表2；對4 個交集DEmiRNAs 進(jìn)行功能富集分析，其生物過程主要涉及肽代謝、轉(zhuǎn)錄、翻譯等；細(xì)胞成分主要涉及受體復(fù)合體、Ⅰ型膠原、Ⅸ型膠原等；分子功能主要涉及細(xì)胞外基質(zhì)結(jié)構(gòu)成分、蛋白質(zhì)酪氨酸/絲氨酸/蘇氨酸磷酸酶活性、輔助轉(zhuǎn)運蛋白活性等，見圖2。

圖2 交集DEmiRNAs的功能富集分析Fig.2 Functional enrichment analysis of intersection DEmiRNAs

表2 關(guān)鍵轉(zhuǎn)錄因子的基本信息Tab.2 Basic information on key transcription factors

圖1 交集DEmiRNAs的轉(zhuǎn)錄因子預(yù)測Fig.1 Transcription factor prediction of intersection DEmiRNAs

2.3 “miRNA-基因”調(diào)控網(wǎng)絡(luò)預(yù)測 本研究運用miRTarBase、TargetScan、miRDB數(shù)據(jù)庫預(yù)測上述4個DEmiRNAs 的靶基因，通過Venny 平臺篩選靶基因在3 個數(shù)據(jù)庫的交集，共得到433 個靶基因，見圖3。最后將得到的靶基因與交集DEmiRNAs構(gòu)成關(guān)系網(wǎng)絡(luò)導(dǎo)入Cytoscape 軟件構(gòu)建miRNA-mRNA 調(diào)控網(wǎng)絡(luò)圖，見圖4。

圖3 交集DEmiRNAs的靶基因韋恩圖Fig.3 Target gene Venn diagram of intersection DEmiRNAs

圖4 交集miRNA-mRNA調(diào)控網(wǎng)絡(luò)Fig.4 Intersection miRNA-mRNA regulatory network

2.4 PPI網(wǎng)絡(luò) 借助String數(shù)據(jù)庫和Cytoscape軟件構(gòu)建PPI網(wǎng)絡(luò)，見圖5。圖中共涉及40 個節(jié)點、43 條邊，其中節(jié)點代表蛋白，邊則代表蛋白間的互作關(guān)系，節(jié)點越大、顏色越深則度值越大；邊越粗，則蛋白間的關(guān)系越緊密。度值排名前5 的蛋白質(zhì)為CDC5L、HNRNPU、GATAD2A、CHD4、ACTB。這些度值較大的蛋白質(zhì)在整個網(wǎng)絡(luò)中發(fā)揮關(guān)鍵作用，可能是RA 跟OA 發(fā)生發(fā)展的關(guān)鍵基因，其基本信息見表3。

表3 關(guān)鍵靶點的基本信息Tab.3 Basic information on key targets

圖5 蛋白互作網(wǎng)絡(luò)Fig.5 Protein interaction network

2.5 GO 功能與KEGG 信號通路富集分析 將靶基因?qū)隓AVID 數(shù)據(jù)庫，分析獲得疾病基因和蛋白信息所富集的通路，并將基因信息進(jìn)行可視化處理，見圖6。結(jié)果顯示，靶基因的生物過程主要涉及mRNA代謝過程的調(diào)控、細(xì)胞周期過程的負(fù)調(diào)控、有絲分裂細(xì)胞周期等；細(xì)胞成分主要涉及轉(zhuǎn)錄因子復(fù)合物、轉(zhuǎn)移酶復(fù)合物、蛋白激酶復(fù)合物等；分子功能主要涉及泛素蛋白轉(zhuǎn)移酶活性、蛋白激酶調(diào)節(jié)劑活性、DNA 解旋酶活性等；信號通路主要涉及調(diào)控干細(xì)胞多能性的信號通路，Hippo信號通路、FoxO 信號通路、Apelin信號通路、p53信號通路等。

圖6 GO和KEGG信號通路可視化結(jié)果Fig.6 Visualization results of GO and KEGG signaling pathway

3 討論

近年來，miRNA 的研究引起了人們的極大興趣，已經(jīng)發(fā)現(xiàn)了2 000 多種miRNAs，它們在真核基因表達(dá)調(diào)控中具有廣泛作用，參與細(xì)胞分化，增殖及凋亡等多種生物過程，人類基因組中約有1/3 的基因處于miRNA的控制之下［13-14］。研究表明miRNA通過調(diào)控自身抗體與炎癥細(xì)胞因子釋放、軟骨細(xì)胞自噬凋亡、調(diào)節(jié)免疫等機制在RA 與OA 疾病的發(fā)病過程中發(fā)揮至關(guān)重要的作用［13-15］。隨著以基因芯片為代表的高通量測序技術(shù)的發(fā)展，目前可在短時間獲取疾病的相關(guān)信息，進(jìn)而從基因?qū)用孢M(jìn)行更全面的生物標(biāo)志物篩選。本研究通過生物信息學(xué)分析了解RA 與OA 之間的相關(guān)性，了解兩病的相關(guān)發(fā)病機理，為藥物研發(fā)提供一定的參考資料。

本研究通過GEO 數(shù)據(jù)庫查找獲取RA 血清基因表達(dá)芯片GSE115885 和OA 血清基因表達(dá)芯片GSE105027，對兩種芯片進(jìn)行生物學(xué)信息分析共獲得4 個DEmiRNAs，分別為hsa-miR-98-5p、hsa-miR-3613-3p、hsa-miR-6858-3p、hsa-miR-4281。研究結(jié)果顯示除hsa-miR-98-5p 對兩種疾病的作用變化趨勢相反，其他三種miRNA 的作用趨勢一致。相關(guān)研究表明RA 與OA 疾病發(fā)生和發(fā)展過程中涉及多種轉(zhuǎn)錄因子，主要涉及滑膜炎癥反應(yīng)，軟骨與骨組織的破壞［16］。本研究對4 個交集DEmiRNA 轉(zhuǎn)錄因子預(yù)測結(jié)果顯示與其相關(guān)的轉(zhuǎn)錄因子共有197 個，差異最顯著的為LHX3、SP4、NFIC、VSX2、HOXA7、TCF3、MYC、HOXB4、ETS1、SP1。研究表明miRNAs參與了幾乎所有已知細(xì)胞過程的轉(zhuǎn)錄后調(diào)控機制，包括發(fā)育、分化、凋亡以及對病原體感染的先天性和適應(yīng)性免疫反應(yīng)，而轉(zhuǎn)錄因子與RNA 聚合酶Ⅱ形成的轉(zhuǎn)錄起始復(fù)合體在轉(zhuǎn)錄起始過程中具有重要的作用［17-18］。本研究結(jié)果對交集DEmiRNAs 進(jìn)行功能富集分析顯示其生物過程主要涉及肽代謝、轉(zhuǎn)錄、翻譯等，與既往研究完全符合［17-18］。因此，筆者認(rèn)為這4種交集DEmiRNAs與RA 和OA 的發(fā)病過程中有著緊密的關(guān)系，可以為闡明RA 與OA 之間的關(guān)系以及后續(xù)研究提供依據(jù)。

對自身靶基因的調(diào)控是miRNAs 在疾病發(fā)病過程中發(fā)揮作用的研究重點。本研究結(jié)果顯示4 個DEmiRNAs 共調(diào)控433個靶基因，其中有359個靶基因為共同調(diào)控，其中CDC5L、HNRNPU、GATAD2A、CHD4、ACTB 排名較靠前，在整個網(wǎng)絡(luò)中起著關(guān)鍵的作用。CDC5L 參與細(xì)胞周期調(diào)控，影響炎癥及細(xì)胞凋亡過程，因此其可通過調(diào)節(jié)炎癥細(xì)胞及軟骨細(xì)胞的凋亡從而影響RA 與OA 疾病的發(fā)生與發(fā)展［19-20］。研究表明HNRNPU 參與免疫基因轉(zhuǎn)錄調(diào)節(jié)過程，對RA 發(fā)病具有關(guān)鍵作用［21-22］。HNRNPU 可與功能性基因間重復(fù)RNA 元件相互作用，增強mRNA的穩(wěn)定性，可促進(jìn)炎癥基因的表達(dá)水平，從而可減輕RA 與OA 患者的炎癥反應(yīng)［23］。CHD4是低氧誘導(dǎo)因子（hypoxia inducible factor，HIF）的共激活因子，CHD4 可增強HIF1α 的募集，從而促進(jìn)HIF 目標(biāo)基因轉(zhuǎn)錄［24］。研究表明HIF1α 介導(dǎo)了免疫炎癥細(xì)胞的激活，在RA 患者滑膜中表現(xiàn)出缺氧性并呈現(xiàn)HIF1α 上調(diào)［25］，HIF1α 也被認(rèn)為是OA 的關(guān)鍵調(diào)節(jié)因子，因此筆者認(rèn)為CHD4 可通過與HIF1 相互作用從而影響兩種疾病的發(fā)生。ACTB 參與軟骨細(xì)胞自噬與凋亡過程，從而影響RA 與OA 患者軟骨退變情況［26］。目前對于GATAD2A 的研究主要集中在癌癥疾病與精神分裂疾病，其可影響患者精神狀態(tài)［27］，尚未有關(guān)于其在RA 與OA 疾病中的研究，但RA 與OA 患者因為疾病的難治性可影響患者的心理及精神情緒，因此未來對于GATAD2A 與RA 與OA 疾病之間的關(guān)系有待進(jìn)一步研究。

為充分了解受DEmiRNAs 調(diào)控的mRNA 在RA與OA 中涉及的共同通路及功能，本研究進(jìn)行了GO及KEGG 分析。GO 結(jié)果顯示靶基因主要涉及mRNA代謝過程的調(diào)控、細(xì)胞周期過程的負(fù)調(diào)控、有絲分裂細(xì)胞周期等生物過程，主要涉及泛素蛋白轉(zhuǎn)移酶活性、蛋白激酶調(diào)節(jié)劑活性、DNA 解旋酶活性等分子功能，這些生物過程及分子功能均與炎癥細(xì)胞與軟骨細(xì)胞的增殖與分化，自噬與凋亡有關(guān)，從而影響RA 與OA 的炎癥變化與關(guān)節(jié)退變程度。KEGG 分析結(jié)果顯示信號通路主要涉及調(diào)控干細(xì)胞多能性的信號通路，Hippo 信號通路、FoxO 信號通路、Apelin 信號通路、p53 信號通路等。研究發(fā)現(xiàn)通過對誘導(dǎo)性多能干細(xì)胞進(jìn)行編程，可使其具有對炎癥刺激作出反應(yīng)的能力，并產(chǎn)生強大的自主調(diào)節(jié)的抗炎細(xì)胞因子，因此調(diào)控干細(xì)胞多能性的信號通路是未來研究RA 與OA 發(fā)病機制的主要途徑［28］。研究表明Hippo 信號通路及p53 信號通路均可影響細(xì)胞增殖和凋亡，抑制其表達(dá)可減少軟骨細(xì)胞凋亡以維持軟骨形態(tài)，F(xiàn)oxO 信號通路可調(diào)節(jié)軟骨細(xì)胞自噬和防御氧化應(yīng)激，從而有效改善RA 與OA 患者的病情，減緩關(guān)節(jié)退變［29-31］。Apelin 可通過自分泌/旁分泌和內(nèi)分泌途徑有效地調(diào)節(jié)炎癥，調(diào)節(jié)軟骨、滑膜、骨骼等組織和各種免疫細(xì)胞的變化，對于RA 和OA的發(fā)病機制至關(guān)重要［32］。

總而言之，RA 和OA 均存在龐大的疾病基因網(wǎng)絡(luò)，但這些網(wǎng)絡(luò)并非各自獨立，而是存在密不可分的聯(lián)系。本研究獲得RA 和OA 各自疾病基因網(wǎng)絡(luò)的交叉部分正是二者相關(guān)聯(lián)的部分，也是研究RA和OA 相互關(guān)系的突破點，同時也可為藥物同時干預(yù)兩種相互關(guān)聯(lián)的疾病提供較為可靠的路徑和作用靶點。且從側(cè)面反映出不同疾病可通過應(yīng)用同一藥物治療的特色，與中醫(yī)辨證論治中的“異病同治”思想不謀而合，以此促進(jìn)中藥現(xiàn)代化發(fā)展。此外本研究尚存在一些缺點：首先，由于篩選條件的限制，只能對主要轉(zhuǎn)錄因子、作用靶點、信號通路進(jìn)行分析，在一定程度上使得研究結(jié)果具有局限性；其次，對差異基因信息進(jìn)行整合處理依賴于生物信息學(xué)技術(shù)的發(fā)展，疾病數(shù)據(jù)庫的完整性、準(zhǔn)確性直接決定著整合后信息的可靠性。故本文雖然能夠篩選出大量靶點和通路，但仍需后續(xù)結(jié)合體內(nèi)外實驗進(jìn)一步驗證與支持，使理論更加可靠。