類風(fēng)濕關(guān)節(jié)炎滑膜組織線粒體功能相關(guān)差異表達(dá)基因的分析

巴燕娜，袁 晴，賀 彤，穆 楠，吳振彪，袁志棟,3※

(1.第四軍醫(yī)大學(xué)西京醫(yī)院臨床免疫科，西安 710032； 2.第四軍醫(yī)大學(xué)基礎(chǔ)醫(yī)學(xué)院生理與病理生理學(xué)教研室，西安 710032；3.贛南醫(yī)學(xué)院基礎(chǔ)醫(yī)學(xué)院，江西 贛州 341000)

類風(fēng)濕關(guān)節(jié)炎(rheumatoid arthritis，RA)是一種以關(guān)節(jié)滑膜炎性病變?yōu)橹鞯南到y(tǒng)性自身免疫性疾病，主要表現(xiàn)為慢性、對(duì)稱性和進(jìn)行性的多關(guān)節(jié)炎[1]。該病影響了近1%的人群，女性發(fā)生率是男性的2～3倍[2]。滑膜組織中成纖維細(xì)胞樣滑膜細(xì)胞在RA中顯示出異常的行為，被認(rèn)為與關(guān)節(jié)損害有關(guān)[3-4]，這種細(xì)胞對(duì)凋亡的抵抗性是RA的一個(gè)主要特征[5]。因此，系統(tǒng)分析RA滑膜組織細(xì)胞的基因表達(dá)水平，不但可加深對(duì)RA病理性過(guò)程及臨床癥狀的理解，而且能更好地理解藥物治療的作用機(jī)制[6]。

線粒體是真核細(xì)胞中的一種重要細(xì)胞器，擁有獨(dú)立的基因組DNA，共37個(gè)基因。這些基因和核基因組編碼的定位于線粒體內(nèi)蛋白的基因?qū)€粒體的氧化磷酸化、能量轉(zhuǎn)換等功能活動(dòng)是必需的。RA患者高度失調(diào)的微血管結(jié)構(gòu)造成了缺氧微環(huán)境，異常的滑膜細(xì)胞代謝和線粒體功能障礙隨之而來(lái)，增加活性氧類，誘導(dǎo)炎癥[7]。另外，在細(xì)胞應(yīng)激或損傷時(shí)，線粒體成分可被釋放到胞質(zhì)或細(xì)胞外，可作為危險(xiǎn)信號(hào)發(fā)揮作用[8]。因而，研究線粒體相關(guān)基因的異常表達(dá)對(duì)探討RA的發(fā)病機(jī)制非常關(guān)鍵。

轉(zhuǎn)錄組芯片分析已被用來(lái)研究RA病理發(fā)生機(jī)制、治療效果和生物標(biāo)志物[9-13]。雖然線粒體與炎癥、免疫密切相關(guān)，但線粒體在RA病理發(fā)生機(jī)制中的功能一直未受到研究者的重視。為探索這些問(wèn)題，本研究通過(guò)系統(tǒng)分析RA與正常人滑膜組織的基因表達(dá)芯片數(shù)據(jù)，以期確定RA病理發(fā)生與線粒體之間的關(guān)系及鑒定關(guān)鍵的影響滑膜組織病變的線粒體相關(guān)基因。

1 材料與方法

1.1滑膜活檢組織的公共芯片數(shù)據(jù) RA滑膜活檢組織的基因表達(dá)數(shù)據(jù)下載于美國(guó)國(guó)家生物技術(shù)信息中心的GEO數(shù)據(jù)庫(kù)的芯片數(shù)據(jù)，ID為GSE77298。該數(shù)據(jù)集包含16例末期RA患者和7例無(wú)關(guān)節(jié)疾病的受試者的滑膜活檢組織，提取總RNA后在Affymetrix芯片上檢測(cè)，然后獲得滑膜組織的基因表達(dá)數(shù)據(jù)[10]。但由于該芯片設(shè)計(jì)上未考慮線粒體基因組編碼的基因，因此只能檢測(cè)線粒體相關(guān)基因而不能檢測(cè)線粒體編碼基因的表達(dá)變化。

1.2分析方法

1.2.1基因表達(dá)芯片數(shù)據(jù)的分析 下載芯片數(shù)據(jù)后，用R軟件安裝simpleaffy、affyPLM和limma等程序包后對(duì)芯片基因表達(dá)數(shù)據(jù)依次處理，利用主成分分析、熱圖聚類分析評(píng)估樣本間的總體相似度，剔除錯(cuò)誤分組的樣本，鑒定差異表達(dá)基因，繪制差異表達(dá)基因的火山圖和熱圖。

1.2.2線粒體相關(guān)基因功能注釋與蛋白網(wǎng)絡(luò)分析 本研究定義核基因組編碼但其蛋白產(chǎn)物定位于線粒體，進(jìn)而與線粒體結(jié)構(gòu)、功能相關(guān)的基因?yàn)榫€粒體功能相關(guān)基因。線粒體功能相關(guān)基因信息來(lái)源于MitoMiner(MitoMiner 4.0v2018JUN)數(shù)據(jù)庫(kù)[14]，該數(shù)據(jù)庫(kù)收錄了1 626個(gè)基因，包括人的1 184個(gè)已知其蛋白產(chǎn)物定位于線粒體的基因和442個(gè)預(yù)測(cè)其蛋白產(chǎn)物定位于線粒體的基因。將其與鑒定的全轉(zhuǎn)錄組差異表達(dá)基因列表取交集，即提取顯著性差異表達(dá)的線粒體功能相關(guān)基因。差異表達(dá)的線粒體相關(guān)基因的功能注釋與蛋白質(zhì)網(wǎng)絡(luò)分析采用在線分析工具M(jìn)etascape[15]，而Metascape應(yīng)用了圖論聚類算法——分子復(fù)合物檢測(cè)MCODE(Molecular Complex Detection)[16]。

2 結(jié) 果

2.1滑膜組織差異表達(dá)的基因 首先對(duì)芯片數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化，并根據(jù)主成分分析和熱圖聚類分析結(jié)果評(píng)估樣本間的總體相似度，發(fā)現(xiàn)RA1、RA2、RA9和RA15樣本數(shù)據(jù)與總體的差異較大，進(jìn)而剔除這4例樣本，并分析基因差異表達(dá)，獲得1 986個(gè)在RA滑膜組織中顯著差異表達(dá)的基因(校正后的P<0.05且表達(dá)值變化2倍以上)，其中表達(dá)水平2倍以上變化的上調(diào)基因有1 490個(gè)，下調(diào)的基因有496個(gè)，見(jiàn)圖1、圖2。這些差異表達(dá)基因中的IGLC1、SPP1、MMP1、JCHAIN和CHI3L2為表達(dá)水平上調(diào)倍數(shù)變化6以上的基因；PLIN1和ADIPOQ(脂聯(lián)素)是表達(dá)水平下調(diào)倍數(shù)變化6以上的基因。

Up：上調(diào)(紅色)；Down：下調(diào)(綠色)；logFC：RA患者與正常人標(biāo)準(zhǔn)化后的基因表達(dá)值比值，即以2為底倍數(shù)變化的對(duì)數(shù)；-log10(adj.P.Val)：以10為底校正后P值的對(duì)數(shù)的相反數(shù)；NotSig：無(wú)顯著變化

圖1 基因表達(dá)的火山圖

HC1-7：正常對(duì)照；RA：類風(fēng)濕關(guān)節(jié)炎患者；藍(lán)色：基因表達(dá)下調(diào)；紅色：基因表達(dá)上調(diào)

圖2 差異表達(dá)基因的層級(jí)聚類熱圖

2.2差異表達(dá)的線粒體功能相關(guān)基因 本研究結(jié)果顯示MitoMiner 4.0數(shù)據(jù)庫(kù)收錄的1 355個(gè)線粒體功能相關(guān)基因在芯片表達(dá)數(shù)據(jù)中有表達(dá)，但由于芯片設(shè)計(jì)的局限性，線粒體基因組編碼的基因未檢測(cè)到表達(dá)。在前面鑒定的1 986個(gè)顯著差異表達(dá)基因中發(fā)現(xiàn)有169個(gè)線粒體功能相關(guān)基因，占目前所知的蛋白產(chǎn)物定位于線粒體的基因總數(shù)量的10.39%，其中表達(dá)上調(diào)的基因數(shù)為131個(gè)，下調(diào)的基因數(shù)為38個(gè)。

2.3差異表達(dá)線粒體相關(guān)基因的功能注釋 為了解169個(gè)差異表達(dá)線粒體相關(guān)基因的功能，對(duì)這些基因進(jìn)行注釋后獲得了以條形圖形式列出了P值排列前20的富集通路，離散色階代表不同的統(tǒng)計(jì)學(xué)差異。P值最顯著的通路是核苷酸代謝過(guò)程通路(GO:0009117)，其次是蛋白質(zhì)在內(nèi)質(zhì)網(wǎng)中的蛋白質(zhì)加工(KEGG通路：HSA04141)和線粒體結(jié)構(gòu)組織通路(GO:0007005)。另外，還包括凋亡信號(hào)的調(diào)控(GO:2001233)、內(nèi)質(zhì)網(wǎng)應(yīng)激反應(yīng)(GO:0034976)、氧化應(yīng)激反應(yīng)(GO:0006979)和活性氧類代謝過(guò)程(GO:0072593)通路等(圖3A)。差異表達(dá)線粒體相關(guān)基因的功能通路網(wǎng)絡(luò)中顏色代表了聚類簇的ID，聚類簇ID與圖3A富集的通路一致，每個(gè)節(jié)點(diǎn)代表了一個(gè)富集的術(shù)語(yǔ)，但功能通路網(wǎng)絡(luò)聚類簇更體現(xiàn)了富集的術(shù)語(yǔ)間的功能相關(guān)性(圖3B)。這些異常表達(dá)線粒體功能相關(guān)基因的功能富集通路體現(xiàn)了RA相關(guān)功能通路出現(xiàn)了擾動(dòng)。

2.4差異表達(dá)的線粒體相關(guān)基因的蛋白相互作用網(wǎng)絡(luò) 利用Metascape生成了差異表達(dá)的線粒體相關(guān)基因的蛋白質(zhì)互作網(wǎng)絡(luò)(圖4A)。Metascape檢測(cè)的蛋白-蛋白相互作用網(wǎng)絡(luò)的密集連接區(qū)，它們代表了分子復(fù)合物。每一個(gè)MCODE成分分別進(jìn)行了通路和過(guò)程富集分析，這些MCODE成分中根據(jù)P值打分最高的術(shù)語(yǔ)是KEGG通路HSA04141[內(nèi)質(zhì)網(wǎng)蛋白質(zhì)加工，log10(P)=-17.1]、GO:0009117[核苷酸代謝過(guò)程，log10(P)=-14.1]和GO:0006753[核苷磷酸鹽代謝過(guò)程，log10(P)=-14.0]。蛋白質(zhì)互作網(wǎng)絡(luò)中4種不同的顏色表示識(shí)別到的4個(gè)MCODE成分(圖4A、圖4B)，分別是MCODE_1：Reactome 通路R-HSA-1799339[靶向膜的SRP依賴的共翻譯蛋白，log10(P)=-11.1]、Reactome通路R-HSA-72766[翻譯，log10(P)=-8.6]和HSA04141[內(nèi)質(zhì)網(wǎng)蛋白質(zhì)加工，log10(P)=-4.2]；MCODE_2：Reactome通路R-HSA-446203 [天冬酰胺N-連接糖基化，log10(P)=-9.5]、KEGG通路HSA04141[內(nèi)質(zhì)網(wǎng)蛋白質(zhì)加工，log10(P)=-8.8]和KEGG通路HSA00510[N-聚糖生物合成，log10(P)=-6.2]；MCODE_3：Reactome通路R-HSA-382551[小分子轉(zhuǎn)運(yùn)，log10(P)=-3.1]和MCODE_4：GO:0050878[體液水平的調(diào)控，log10(P)=-3.8]。

4A：線粒體相關(guān)基因的蛋白質(zhì)互作網(wǎng)絡(luò)圖； 4B：蛋白質(zhì)網(wǎng)絡(luò)中密集連接的蛋白質(zhì)組

圖4 差異表達(dá)的線粒體相關(guān)基因編碼的蛋白質(zhì)互作用網(wǎng)絡(luò)

3 討 論

線粒體在細(xì)胞內(nèi)不僅參與物質(zhì)的代謝、能量轉(zhuǎn)化和鈣離子存儲(chǔ)等生命活動(dòng)，還在病理狀態(tài)下由損壞的線粒體導(dǎo)致先天性免疫系統(tǒng)異常激活而引起自身性免疫疾病[17]。然而，在RA中的線粒體功能一直未引起研究者的重視。為了探究線粒體功能相關(guān)的基因在RA中的表達(dá)變化和功能，本研究分析了GEO數(shù)據(jù)庫(kù)基因表達(dá)芯片數(shù)據(jù)GSE77298，獲得了1 986 個(gè)在RA滑膜組織中顯著差異表達(dá)的基因，這些基因與DisGeNET數(shù)據(jù)庫(kù)[18]收錄的1 832個(gè)RA相關(guān)基因的交集為322個(gè)，與Liu等[19]報(bào)道的RA患者中顯著差異表達(dá)上調(diào)的228個(gè)(實(shí)際是227個(gè))基因和下調(diào)的36個(gè)基因相比，相一致的上調(diào)基因110個(gè)，相一致的下調(diào)基因10個(gè)；另一項(xiàng)研究報(bào)道在滑膜組織中RA組比正常組高表達(dá)的基因1 018個(gè)，下調(diào)的基因893個(gè)[20]，與本研究結(jié)果相同的分別為686個(gè)和257個(gè)。另外，在RA滑膜組織中顯著差異表達(dá)、倍數(shù)變化6以上的上調(diào)基因SPP1、MMP1在以前的研究中[21-23]也得到了驗(yàn)證。而ADIPOQ是倍數(shù)變化6以上的表達(dá)下調(diào)基因，Harshan等[20]的研究結(jié)果顯示脂聯(lián)素在RA患者滑膜組織中表達(dá)水平較正常人顯著降低，本研究結(jié)果與上述文獻(xiàn)一致。在這1 986個(gè)顯著差異表達(dá)基因中，本研究鑒定了169個(gè)與正常對(duì)照存在統(tǒng)計(jì)學(xué)差異表達(dá)的線粒體功能相關(guān)基因。而這些線粒體功能相關(guān)基因部分也得到了DisGeNET數(shù)據(jù)庫(kù)[18]及RA相關(guān)研究報(bào)道[19-20,24-28]的數(shù)據(jù)支持與驗(yàn)證。169個(gè)線粒體相關(guān)基因的25個(gè)與DisGeNET數(shù)據(jù)庫(kù)中的一致，其中FHL1、G0S2和SLIT3為表達(dá)下調(diào)的3個(gè)基因，其他的22個(gè)為上調(diào)基因。SLIT3能夠抑制RA中Robo3誘導(dǎo)的滑膜成纖維細(xì)胞的入侵[24]，因此SLIT3的表達(dá)下調(diào)可能促進(jìn)RA的病情惡化。另外，G0S2可作為生物標(biāo)志物用來(lái)預(yù)測(cè)響應(yīng)抗腫瘤壞死因子治療RA患者[27]。而表達(dá)上調(diào)倍數(shù)變化最大的為TYMS，其遺傳多態(tài)性與RA患者的甲氨蝶呤臨床治療的效果相關(guān)[25-26]。

本研究發(fā)現(xiàn)的部分差異表達(dá)基因不僅與已往研究報(bào)道[18-28]一致，而且通過(guò)功能通路注釋和互作網(wǎng)絡(luò)分析，結(jié)果顯示GO:0009117(核苷酸代謝過(guò)程)和KEGG通路HSA04141(內(nèi)質(zhì)網(wǎng)蛋白質(zhì)加工)不僅分別是差異表達(dá)的線粒體相關(guān)基因功能富集通路中P值最顯著的通路之一，也分別是蛋白互作網(wǎng)絡(luò)MCODE成分中根據(jù)P值得分最高的生物學(xué)過(guò)程術(shù)語(yǔ)之一。將蛋白質(zhì)互作數(shù)據(jù)與通路和生物學(xué)過(guò)程富集數(shù)據(jù)相結(jié)合顯示了差異表達(dá)的線粒體功能相關(guān)基因的蛋白質(zhì)互作網(wǎng)絡(luò)可能在RA中存在異常。而有報(bào)道發(fā)現(xiàn)青蒿素能通過(guò)抑制PDK1(GO:0009117共有34個(gè)差異表達(dá)的線粒體功能相關(guān)基因之一)誘導(dǎo)的蛋白激酶B和核糖體S6激酶家族成員2磷酸化激活來(lái)抑制RA成纖維樣滑膜細(xì)胞的遷移和入侵[28]。

另外，富集的通路還包括線粒體結(jié)構(gòu)組織通路(GO:0007005)、凋亡信號(hào)的調(diào)控(GO:2001233)、內(nèi)質(zhì)網(wǎng)應(yīng)激反應(yīng)(GO:0034976)、氧化應(yīng)激反應(yīng)(GO:0006979)和活性氧類代謝過(guò)程(GO:0072593)通路，以及蛋白質(zhì)互作網(wǎng)絡(luò)中識(shí)別到的MCODE成分R-HSA-1799339 (靶向膜的SRP依賴的共翻譯蛋白)等。以上功能與通路的富集注釋結(jié)果以及已有相關(guān)研究報(bào)道[18-28]的驗(yàn)證顯示了這些線粒體功能相關(guān)基因的異常表達(dá)可能與RA滑膜組織的病理變化有關(guān)。

需要指出的是，本研究分析的是芯片數(shù)據(jù)，由于該芯片探針未覆蓋線粒體基因組編碼的基因，所以分析結(jié)果僅鑒定了核基因組編碼的線粒體功能相關(guān)基因，并且這些基因均為編碼蛋白基因。因此，后續(xù)更完整的研究應(yīng)囊括1 626個(gè)定位線粒體內(nèi)的相關(guān)編碼基因、線粒體基因組編碼的基因，以及非編碼RNA基因，以期能夠全面系統(tǒng)地揭示線粒體對(duì)RA的發(fā)生、發(fā)展及滑膜組織病變所起的作用，找到有效的RA疾病診斷分子標(biāo)志物和治療的藥物靶標(biāo)。