類風濕關節炎差異表達基因篩選、生物信息學分析及其潛在治療藥物的預測

孫勝南，廖苑君，藍樹金，李讓，趙小蕾，覃繼恒，饒紹奇

1廣東醫科大學醫學系統生物學研究所，廣東東莞 523808；2廣東醫科大學公共衛生學院

類風濕關節炎(RA)是一種慢性炎癥性關節自身免疫性疾病，病理特征為關節滑膜細胞增生，新生血管形成血管翳，炎癥細胞浸潤，軟骨和骨結構破壞[1,2]。RA患者存在持續、反復的滑膜炎癥，以及嚴重的關節損傷，關節炎嚴重者可致殘，使患者勞動力下降，嚴重影響生活質量[3]。文獻[4]報道，每年每10萬人中有5～50人罹患RA，其中最常發生于婦女和老年人。文獻[5]報道，RA患病率和發病率呈逐年上升趨勢，RA與其他風濕病是一個日益嚴重的公共衛生問題。雖然抗風濕藥物甲氨蝶呤等能減弱或緩解RA的病情[6]，但患者治愈率仍不高。因此，RA早期階段的檢測非常必要[7]，然而現有的RA生物標志物預測價值有限，迫切需要進一步深入研究RA發生發展的分子機制，從而實現RA的早發現、診斷、治療，以提高治療效果，改善患者預后。基因芯片是基于堿基互補配對原則，利用化學的合成等技術，出現的高通量測序方法。至今為止，已經建成了很多的基因序列、表達量公共數據庫，但這些基因表達量數據在不同疾病模型及生物學發育過程中所代表的生物學意義卻知道的少之又少。大量資料表明，疾病的發生是受特異基因表達量的變化影響，這些基因變化的表達水平與人類正常的表達水平有明顯的不同[8]，利用統計學、生物學、計算機等不同領域的幫助，提取不同生理、病理條件下基因表達量所代表的意義，研究不同基因所在通路，從而發現關鍵性的致病基因，這對RA的早期發現和診斷有巨大的意義。本研究通過篩選RA差異表達基因，并對差異表達基因進行生物信息學分析，且預測潛在的RA治療藥物，為RA的分子發生機制、臨床診斷和藥物治療提供參考。

1 資料與方法

1.1 芯片數據的篩選 截止至2019年8月，以“Rheumatoid arthritis(RA)，homo sapiens”為關鍵詞從公共Gene Expression Omnibus(GEO)數據庫(https://www.ncbi.nlm.nih.gov/GEO/)搜索目標基因芯片，經過篩選，最終選取GSE97779和GSE1919作為本次研究的目標數據基因芯片。GSE97779共有14個樣本，由9例RA患者樣本和5例健康樣本組成；GSE1919共有15例樣本，由5例RA患者樣本、5例健康樣本、5例骨關節患者樣本(不納入本研究)組成。

1.2 RA差異表達基因篩選 GEO2R(https://www.ncbi.nlm.nih.gov/geo/GEO2R/)是一個基于R語言的差異表達基因篩選工具，可以比較疾病和正常組之間的差異。通過輸入序列號、定義組別、分配組別(疾病組和健康組)和單擊“Top 250”，得到疾病組和健康組差異分析結果。在此結果中篩選差異表達基因，需同時滿足以下兩個條件:①P<0.05;②|log2FC|>1，即差異倍數(Fold change)FC>2或FC<-2。

1.3 RA差異表達基因蛋白質相互作用(PPI)網絡構建及核心基因提取 將所得差異表達基因作為初始種子基因上傳到String數據庫(http://string-db.org/)，得到基因互作關系網絡圖，以置信分數(confidence)>0.4為篩選標準，將符合條件的互作基因關系導入到Cytoscape軟件進行可視化，得到PPI網絡。運用Cytoscape的cytoHubba插件，獲取排名前十位的核心基因。

1.4 RA差異表達基因基因本體(GO)功能富集分析、京都基因與基因組百科全書(KEGG)通路富集分析 將篩選出的差異表達基因上傳至DAVID數據庫(https://david.ncifcrf.gov/)，進行GO功能富集分析和KEGG通路富集分析。GO分析包括生物學過程、細胞成分、分子功能三方面內容，主要對差異表達基因進行生物學過程、細胞成分、分子生物學功能進行分析；KEEG主要針對差異表達基因的信號通路進行富集分析。

1.5 RA潛在治療藥物篩選 藥物-基因相互作用數據庫(DGIdb，http：//dgidb.org/)是一個整合了來自組織和展示論文、數據庫和網絡資源的藥物-基因相互作用和基因可藥性信息的數據庫。DGIdb將來自30個不同來源的內容標準化，包含關于查詢藥物、基因和藥物-基因相互作用更詳細的信息。將核心基因列表上傳到DGIdb數據庫基因清單窗口框，進行尋找基因-藥物作用關系，得到PubMed id列表、相互作用類型和得分等結果數據，可直接將結果TSV格式下載到本地，根據得分篩選潛在治療藥物。

2 結果

2.1 RA差異表達基因篩選結果 按照篩選條件，得到GSE97779和GSE1919共有的100個差異表達基因，其中上調基因77個，下調基因23個(詳見圖1)。

2.2 RA差異表達基因的PPI網絡構建和核心基因篩選結果 最終得到具有81個節點、383條邊的PPI網絡(圖2A)，經過Cytoscape的cytoHubba插件篩選得到10個核心基因(圖2B)，分別是C-C基序趨化因子受體5(CCR5)、選擇素(SELL)、顆粒酶B(GZMB)、CD2分子(CD2)、白細胞介素2受體亞單位(IL2RB)、CD44分子(CD44)、白細胞介素受體7(IL7R)、白細胞介素2受體亞單位(IL2RG)、淋巴細胞胞質蛋白2(LCP2)、LCK原癌基因(LCK)。

注：Up_List1、down_List1表示GSE97779集合中的上調、下調差異表達基因，Up_List3、down_List4表示GSE1919集合中的上調、下調差異表達基因。

圖1GSE97779和GSE1919共有的100個差異表達基因

注：A為RA差異表達基因的PPI網絡，B為10個核心基因。

圖2RA差異表達基因的PPI網絡構建和關鍵基因篩選結果

2.3 RA差異表達基因GO功能富集分析、KEGG通路富集分析結果 RA差異表達基因GO功能富集分析結果以P<0.05為篩選條件，按照P值從大到小的順序排列；在生物學過程中，這些差異表達基因主要參與了免疫反應、趨化因子介導信號通路、趨化作用、細胞趨化作用、炎癥反應的正調控、炎癥反應、T細胞受體信號通路、信號轉導、白細胞遷移、細胞外基質、T細胞共刺激、調節細胞增殖、細胞間信號、細胞表面受體信號通路、細胞防御反應、調節免疫反應、T細胞分化、T細胞分化的調控、細胞黏附、適應性免疫反應和T細胞趨化性的正調控過程；由細胞組分可知，這些差異表達基因表達主要位于質膜外側面、細胞外區域、原生質膜、α-βT細胞受體復合物、細胞外空間、細胞表面、膜、免疫突觸和質膜部位；從分子功能上看，這些差異表達基因具有使趨化因子的活動、CXCR3趨化因子受體結合的作用。RA差異表達基因KEGG通路富集分析結果以P<0.05為篩選條件，按照P值從大到小的順序排列，這些差異表達基因表達通路主要富集于細胞因子—細胞因子受體相互作用通路、趨化因子信號通路、類風濕關節炎、原發性免疫缺陷、自然殺傷細胞介導的細胞毒性、造血細胞系、Jak-STAT信號通路、T細胞受體信號通路。

2.4 RA潛在治療藥物預測篩選結果 根據PPI網絡中篩選出的10個核心基因，分別為CCR5、SELL、GZMB、CD2、IL2RB、CD44、IL7R、IL2RG、LCP2、LCK。在DGIdb數據庫中尋找相關潛在藥物，最后搜索發現除LCP2外的9個基因都找到相關藥物。由基因CCR5預測到的潛在治療藥物有胍基乙酸、依諾肝素、安克里維克、阿帕維洛克、塞尼克羅克、馬拉維克、酚普羅考蒙、阿普洛克、馬來酸維里韋洛克；由基因SELL5預測到的潛在治療藥物有雙咪唑糖、里維潘塞爾、阿塞單抗；由基因GEZM5預測到的潛在治療藥物有(2S)-2-(甲胺基)-1-苯丙-1-醇；由基因CD25預測到的潛在治療藥物有阿利法西普、西普利單抗、奧曲肽、葡萄牛磺酸；由基因IL2RB5預測到的潛在治療藥物有艾地氯喹、地尼白介素、巴利昔單抗、賽尼哌、硫柳單抗；由基因CD445預測到的潛在治療藥物有是牛磺酸、干擾素GAMA-1B、甲酚酸、慶大霉素；由基因IL7R5預測到的潛在治療藥物有魯索利替尼；由基因IL2RG5預測到的潛在治療藥物有賽尼哌、巴利昔單抗、地尼白介素；由基因LCK5預測到的潛在治療藥物有達沙替尼、波納替尼、阿卡魯替尼、十二醇、伊布替尼、薩拉卡蒂尼布、鹽酸帕佐帕尼、伊洛拉塞提布、葡萄孢素。我們選取得分排名前三的藥物馬拉維克、阿利法西普、艾地氯喹作為治療RA的候選治療藥物。

3 討論

本研究通過對GEO數據庫中的基因芯片GSE97779和GSE1919進行生物信息學分析獲得差異表達基因，并分析了有關這些基因富集的生物學過程、細胞成分、分子功能和信號通路，為RA的機制研究提供了理論依據與研究方向。在我們的研究中，RA病例組和對照組之間共100個差異表達基因。功能分析顯示，差異表達基因主要富集于細胞因子活性、G 蛋白偶聯受體結合和細胞因子受體結合。KEGG通路分析顯示，差異表達基因主要富集于細胞因子-細胞因子受體相互作用、趨化因子信號通路、病毒蛋白與細胞因子、細胞因子受體相互作用。

CD44抗原是一種參與細胞間相互作用、黏附和遷移的細胞表面糖蛋白。在人類中，CD44抗原由11號染色體上的CD44基因編碼。CD44被認為粘連分子、吞噬糖蛋白-1、Hermes抗原、淋巴細胞歸巢受體，一系列研究表明CD44可能是RA的標志物，與RA的發生發展相關[9～11]。這些研究表明，本研究的篩選結果與前人研究一致。CD2是在T細胞和NK細胞表面發現的一種細胞黏附分子，參與了CD8+T細胞和NK細胞介導的人類和其他幾種哺乳動物細胞免疫的免疫調節[12]。CCR5是一種趨化因子受體，被多個參與炎癥反應的免疫細胞廣泛表達，我們推測CCR5可能通過這些炎癥途徑通路參與了RA發生。GZMB基因沉默通過調節炎癥因子、凋亡相關因子、血管生成相關因子的表達，抑制MAPK信號通路，緩解RA引起的滑膜組織增生和關節軟骨組織損傷[13]。LCP2是淋巴細胞胞質蛋白2，在PubMed數據庫中未發現LCP2與RA有直接關系。但有研究證實其與腫瘤密切相關，包括結直腸癌[14,15]、多形性膠質母細胞瘤[16,17]、彌漫性大b細胞淋巴瘤[18]。癌癥與RA可能存在共同的免疫通路。LCK是一種非受體src家族激酶，在細胞周期控制、細胞黏附、運動、增殖和分化等多種細胞過程中發揮重要作用。有研究表明，LCK是調節T細胞功能的關鍵分子之一，調節TCR信號的啟動、T細胞的發育和T細胞內穩態[19]。我們推測LCK通過調節T細胞免疫通路影響RA的發生發展。到目前為止，還沒有研究顯示有SELL與RA的發生、發展有直接的聯系。因此，對該基因的進一步研究可能會有新的發現，為RA的分子發病機制研究提供了新的方向。

功能富集分析中，100個差異表達基因主要富集于免疫、炎癥、信號通路，如白細胞遷移、細胞黏附的正向調節、受體配體活性、細胞因子-細胞因子受體相互作用通路等。目前藥物控制是RA患者的重要治療手段，尋找有效、敏感的抗RA藥物有助于改善RA患者的預后。10個核心基因共發現76種化學物質，其中選擇得分前3的化合物為候選藥物。馬拉維克是一種一級選擇性CCR5拮抗劑，僅被批準與其他抗逆轉錄病毒藥物聯合治療艾滋病病毒感染的[20]，已有研究報道，馬拉維克有助于免疫綜合征的初始預防和免疫綜合征的積極治療[21]，減輕神經性疼痛癥狀的發展[22]。但其抗RA作用及作用機理至今尚未闡明。在本研究中，我們發現了馬拉維克作為RA治療劑的潛力。關于阿利法西普和艾地氯喹，已有研究報道它們與銀屑病和新發型1型糖尿病有關[23,24]。我們的研究為RA的治療提供了新的藥物靶點，當然還需要更多的實驗研究來證實這些發現。

總之，本文從GEO數據庫獲得RA相關數據集，基于GEO2R采用limma包分析篩選出100個差異表達基因、PPI網絡構建及核心基因篩選，并對其進行GO和KEGG功能分析和潛在藥物的預測。綜上所述，我們采用大數據挖掘和計算生物學的綜合生物信息學方法篩選出了10個與RA發生有關的核心基因(CCR5、SELL、GZMB、CD2、IL2RB、CD44、IL7R、IL2RG、LCP2、LCK)和3種潛在治療藥物。本研究為RA的發病機制、早期診斷和臨床藥物治療提供有價值的參考。