基于基因芯片整合分析篩選慢性牙周炎相關基因和轉錄因子

曾小麗 李生嬌 單錚男 殷俊豪 姜吉蕊 鄭章龍 李家

1.上海牙組織修復與再生工程技術研究中心，同濟大學附屬口腔醫院·同濟大學口腔醫學院口腔頜面外科，上海200072；2.上海牙組織修復與再生工程技術研究中心，同濟大學附屬口腔醫院·同濟大學口腔醫學院修復科，上海200072

牙周炎是一種常見的口腔科疾病，慢性牙周炎的持續進展將導致牙齒的松動、脫落乃至口腔咬合功能的完全喪失[1-3]，作為被多種因素混雜影響的慢性炎癥性疾病，牙周炎的發生、發展過程中有大量基因及基因產物參與其中，因此，尋找牙周炎發病過程中的關鍵基因及相關因子是目前研究的一大熱點。

隨著高通量技術的發展，通過芯片或測序等技術研究疾病與對照組間的基因表達譜的差異成為篩選致病基因的有力工具。基因表達量的變化往往預示著生理過程的改變或病理過程的發生，與此同時，基因編碼的產物所涉及的分子功能、生物過程、信號通路則可能與疾病的發生發展相關。由于慢性牙周炎的發生發展是多種基因及其產物的共同作用效果，單個基因或蛋白的表達變化往往不足以完全反映疾病的發展過程。此外，由于樣本數目小以及芯片檢測平臺的差異，過去不同研究者獲得結果并不完全一致，甚至可能互相矛盾。因此，本實驗基于大樣本量、相同檢測平臺等檢索要求，從基因表達綜合（gene expression omnibus，GEO） 數據庫中篩選并獲得符合標準的兩個慢性牙周炎相關數據集GSE10334 和GSE16134 并深入分析其中差異表達基因的功能、相關通路以及核心模塊、基因，為未來的研究指明方向[4-5]。

1 材料和方法

1.1 基因芯片數據獲取

在GEO 數據庫（https://www.ncbi.nlm.nih.gov/geo/） 中以 “牙周炎（periodontitis）” 為關鍵詞檢索相關芯片，并篩選基于同一檢測平臺且樣本量大于100 的基因表達數據集，獲取了2 組符合標準的數據集：GSE10334 和GSE16134。兩組數據集均采用美國Affymetrix 公司提供的Human Genome U133 Plus 2.0 Array 基因芯片進行檢測，平臺號為GPL570。兩組芯片數據的樣本采集于中重度牙周炎患者炎性區域牙齦組織（試驗組） 和健康牙齦組織（對照組），兩組使用相同的采樣方法，并具有一致的樣本納入標準和臨床分型指標[4-5]。其中，GSE10334 的芯片數據來源于90 位患者的247 例樣本（183例病變，64例健康）；GSE16134的芯片數據來源于120 位患者的310 例樣本（241 例病變，69例健康）[4-5]。

1.2 篩選差異表達基因

使用GEO2R 在線分析工具（http://www.ncbi.nlm.nih.gov/geo/geo2r/） 進行試驗組與對照組之間差異表達基因（differentially expressed gene，DEG）的篩選，篩選閾值設定為調整P值（adjustedPvalue，adj.P） ＜0.05，差異倍數（fold change，FC） ＞2。使用R 語言的ggplot2 軟件包繪制火山圖和熱圖。

1.3 進行功能富集分析

對兩組數據集的差異基因取交集后進行后續分析，以增加結果可信度。使用g: Profiler （https://biit.cs.ut.ee/gprofiler/gost） 分別對上、下調基因進行功能富集分析，功能富集分析包括對基因的基因本體論（gene ontology，GO） 分析和通路分析。GO 分析涵蓋了分子功能、生物過程和細胞組分3個下級類目。通路分析中的生物通路信息來源于京都基因與基因組百科全書（Kyoto encyclopedia of genes and genomes，KEGG）、Reactome和WIKI共3個數據庫。adj.P＜0.05 視為差異具有統計學意義。

1.4 建立蛋白互作網絡

為了探索差異表達基因之間的相互作用關系，筆者通過在線工具STRING 數據庫（https://stringdb.org/） 建立蛋白相互作用（protein-protein interaction，PPI） 網絡，然后使用Cytoscape （www.cytoscape.org/） 軟件將該網絡進行可視化。

1.5 核心模塊與核心基因確立

使用Cytoscape 軟件中的MCODE 插件并基于默認閾值（K-core=2，Node Score Cutoff=0.2，Degree Cutoff=2，Maximum Depth=100） 篩選PPI 網絡中的重要功能模塊。使用cytoHubba 插件進行PPI 網絡中關鍵節點的篩選。cytoHubba 插件使用12 種拓撲分析策略進行重要節點的預測和評估，在6 種以上算法中得分均位于前10 的節點被視作核心基因。為進一步了解PPI網絡中重要子網絡的功能，筆者對MCODE核心模塊中的基因進行了功能富集分析，并繪制通絡分析氣泡圖。

1.6 預測轉錄因子

使用Cytoscape 軟件中的iRegulon 插件進行轉錄因子的預測。該插件使用標準化富集分數（normalized enrichment score，NES） 來評估預測結果的可信度，NES 值越大，說明可信度越高，筆者將NES＞4.5的轉錄因子將用于調控網絡的建立。

2 結果

2.1 篩選出的差異表達基因

基于adj.P＜0.05 且FC 絕對值＞2 的篩選標準，在GSE10334 中篩選出221 個DEG，其中150 個上調，71 個下調（圖1A）；在GSE16134 中篩選出297 個DEG，其中212 個上調， 85 個下調（圖1C）；兩組DEG 在試驗組與對照組的表達量變化及其基因、樣本聚類結果見圖1B、D。

圖1 試驗組與對照組間差異基因的篩選Fig 1 Identification of differentially expressed genes between test group and control group

2.2 差異表達基因的功能富集分析結果

在兩組數據集中，196 個DEG 表達趨勢一致，包括139 個上調基因和57 個下調基因（圖2A、D）。GO 分析結果顯示：上調基因共富集于10 項分子功能、99 項生物過程以及20 項細胞組分，下調基因共富集于1項分子功能、11項生物過程以及7項細胞組分。圖2B、E展示了上、下調DEG adj.P前10 的GO 條目。在分子功能方面，上調基因參與趨化因子受體結合、免疫球蛋白結合，而下調基因參與皮膚表皮結構形成。在生物過程方面，上調基因涉及免疫應答、刺激應答、白細胞活化等，而下調基因涉及皮膚及表皮的發育和角化。在細胞組分方面，上調基因主要位于胞外間隙、細胞表面、胞外體、質膜外側，而下調基因主要位于超分子纖維、超分子聚合物、聚合細胞骨架纖維、透明角化顆粒等處。

通路分析結果顯示：上調基因富集于26 條信號通路，其中包括12 條KEGG 通路、12 條Reactome 通路和2 條WIKI 通路；而下調基因富集于2條Reactome 通路。圖2C、F 分別展示了上調基因adj.P前15 的通路條目和下調基因所富集的2 條通路條目。上調基因參與細胞因子-細胞因子受體相互作用、白細胞經內皮細胞遷移、趨化因子信號通路以及人補體系統、造血干細胞分化等免疫反應相關通路，而下調基因涉及角質化包膜形成、角化等2條通路。

圖2 差異基因功能富集分析結果Fig 2 Functional enrichment analysis of differentially expressed genes

2.3 PPI網絡建立及核心模塊、核心基因篩選結果

由交集基因構成的PPI 網絡包含132 個節點和519 個節點作用關系（圖3）。圖4A、B 顯示：PPI網絡的第一模塊涉及病毒蛋白與細胞因子及細胞因子受體的相互作用、細胞因子-細胞因子-受體相互作用和趨化因子信號通路等通路。圖4C、D 顯示：第二模塊參與B細胞受體信號通路。圖4E、F顯示：第三模塊與白細胞介素介導的信號通路、JAK-STAT信號通路等通路有關。

圖3 構建PPI網絡并篩選核心基因Fig 3 Construction of PPI network and identification of hub genes

圖4 PPI網絡核心模塊篩選與通路分析結果Fig 4 Identification and pathway analysis of key modules in PPI network

cytoHubba 插件的分析結果顯示：滿足篩選標準的基因為CXC 基序趨化因子配體（C-X-C motif chemokine ligand，CXCL） 8、CXCL1、CXC 基序趨化因子受體（C-X-C motif chemokine receptor，CXCR） 4、選擇素（selectin，SEL） L， CD19 分子（CD19 molecule，CD19） 和IKAROS 家族鋅指（IKAROS family zinc finger，IKZF） 1。在PPI 網絡中，CXCL1、CXCL8、CXCR4 和SELL 位于第一模塊，IKZF1位于第三模塊。

2.4 轉錄因子預測結果與調控網絡的建立

依據iRegulon預測結果，NES＞4.5的轉錄因子為干擾素調節因子（interferon regulatory factor，IRF） 4、組蛋白乙酰轉移酶（E1A binding protein，EP） P300、易洛魁同源盒蛋白（iroquois homeobox，IRX） 4 和TATA 序列結合蛋白（TATA-box binding protein，TBP），分別調控60、18、24 和36個DEG（圖5），其中IRF4屬于上調DEG且位于第二核心模塊，參與B細胞受體信號通路（圖4 C、D）。

圖5 轉錄因子及其調控網絡Fig 5 Transcription factors and their regulatory networks

3 討論

牙周炎不僅是成年人失牙的主要原因，也是多種疾病如糖尿病、動脈粥樣硬化、慢性阻塞性肺疾病等系統性疾病的潛在風險因素[3,6]。鑒定與牙周炎進展相關的相關分子及信號通路，有利于了解牙周炎的發病機制，為牙周炎的早期診治及伴發疾病的預防提供理論基礎。本研究基于GEO數據庫公開的牙周炎基因表達譜，通過系統性生物信息學分析，為牙周炎的機制研究提供理論依據和研究方向。由于樣本量小會在一定程度上降低試驗結果的可信度，因此本研究篩取大樣本數據集GSE10334 和GSE16134 （分別含有247 例和310 例樣本） 進行生物信息學分析。聯合分析兩組芯片數據集，在兩組中均差異表達的基因共196個，其中139 個為上調基因，57 個為下調基因。功能富集分析結果顯示，上調基因參與免疫系統，病毒蛋白與細胞因子及細胞因子受體的相互作用，細胞因子-細胞因子受體相互作用，白細胞跨內膜遷移和趨化因子受體結合趨化因子等免疫反應相關途徑，由此可見免疫應答相關途徑在牙周炎發病過程中發揮重要作用。而下調基因涉及角質化包膜形成、角化等途徑，說明上皮屏障形成能力降低，提示牙周炎的發生與牙齦上皮完整性喪失相關。由差異基因構建的PPI 網絡中，MCODE 得分前三的模塊被視作核心模塊，分別富集于趨化因子信號通路、B細胞受體信號通路以及白細胞介素介導的信號通路等途徑，由此推測，這些分子事件在牙周炎發病過程中活躍。CXCL8、CXCL1、CXCR4、SELL、CD19和IKZF1為核心基因，位于PPI網絡的中心位置，因此最有可能成為牙周炎的治療靶點。

CXCL8 和CXCL1 屬于趨化因子CXC 亞家族，通過募集白細胞來增強或延長炎癥反應。多項研究[7-9]發現，相較于健康牙齦，炎癥牙齦組織切片中CXCL8 和CXCL1 的mRNA 水平和蛋白水平均顯著升高。這些結果提示：炎性環境下牙周細胞可能通過分泌更多的趨化因子以對抗牙周病原體的感染，而這種過度的免疫應答正是牙周組織破壞的主要原因。此外Souto 等[1]研究發現，CXCL8水平與牙周炎患者的探診出血、探診深度和臨床附著喪失等臨床檢測指標正相關，其結果提示CXCL8 表達量與牙周炎嚴重程度相關，有望成為牙周炎進展標志物。

CXCR4 是CXC 趨化因子受體家族成員，Hajishengallis 等[10]發現，牙齦卟啉單胞菌可通過結合巨噬細胞表面CXCR4，逃避宿主的免疫清除。口服CXCR4 的拮抗劑AMD3100，可抑制牙齦卟啉單胞菌的定殖，進而預防或減輕小鼠實驗性牙周炎及其引起的牙槽骨吸收[11]。因此推測，CXCR4是牙齦卟啉單胞菌在牙周炎致病作用的一個關鍵位點。

SELL 編碼淋巴細胞歸巢受體L-選擇素，該因子屬于細胞表面黏附分子，調控白細胞在血管內的黏附與遷移。L-選擇素主要表達于中性粒細胞，在炎性疾病中，L-選擇素結合其相應配體，調控中性粒細胞活化與聚集，并使其從血管中遷移至炎癥區域。中性粒細胞是固有免疫中的重要組成細胞，在免疫防御中起先驅作用，但近年來的研究[12]發現，過度活化的中性粒細胞可通過產生大量活性氧和基質降解酶造成牙周組織的破壞。因此，SELL 可能通過促進中性粒細胞遷移和活化參與牙周炎病理發生。

CD19 是特異性表達于B 細胞表面的免疫球蛋白基因超家族成員[10-11,13-14]，為B細胞抗原受體復合物的共受體，可降低激活下游信號通路和觸發B細胞對抗原的反應的閾值[14]。因此，CD19 在試驗組的表達量升高，提示在炎性牙齦中B 細胞數目和免疫應答活性增加。有研究[15]表明，B 細胞是牙周炎骨缺損區分泌RANKL的主要細胞之一。使用抗B 細胞藥物如利妥昔單抗則可以緩解牙周炎的多種癥狀[2]。因此，以CD19 為靶點拮抗B 細胞的藥物研發，可能成為藥物治療牙周炎的一個研究方向。

IKZF1 編碼淋巴系特異性轉錄因子IKAROS，IKAROS 是造血細胞分化的關鍵調節因子，調控T淋巴細胞和B 淋巴細胞抗原受體介導的增殖反應以及相關效應細胞的功能[16]，T 細胞與B 細胞參與牙周炎結締組織區域的炎癥浸潤及骨喪失[2,15]。

IRF4 為iRegulon 預測結果中可信度最高的轉錄因子， 調控包括核心基因CXCR4、 CD19 和IKZF1在內的60個DEG的轉錄，IRF4 mRNA與蛋白質在炎性牙齦中的表達已在過去的文獻[7]報道中得到支持。雖然IRF4 在牙周炎發病過程中的具體作用尚不明確，但IRF4 參與系統性紅斑狼瘡、類風濕性關節炎等自身免疫性疾病的發生[17]，提示其與機體免疫亢進相關。

在本研究中，為在多樣本高通量數據中快速高效找到顯著上調或下調的基因，筆者使用了嚴格的閾值篩選方法，但值得注意的是，這種方法可能會遺漏部分具有重要生物學意義但表達量變化并不顯著的基因。例如，白細胞介素-1β 和腫瘤壞死因子-α 作為牙周炎疾病中重要的炎性因子和骨吸收介質，可激活下游核因子κB 信號通路參與牙周炎相關骨丟失[18]。然而由于不滿足FC＞2 的閾值標準，這些牙周炎相關基因在本實驗中不作為“差異表達基因” 納入后續分析，其在牙周炎發病過程中的作用也未得到體現。因此本研究尚存在一定局限性，相關分析結果需要通過后續的分子生物學實驗和動物實驗進一步驗證。

綜上所述，本研究基于大樣本芯片數據集研究，通過生物信息學分析篩選、鑒定與牙周炎發病相關的易感基因、轉錄因子及相關信號通路。共發現196 個牙周炎相關差異基因，包括6 個核心基因（CXCL8、CXCL1、CXCR4、SELL、CD19和IKZF1） 和1 個關鍵轉錄因子（IRF4），這些分子通過細胞因子-細胞因子受體相互作用、趨化因子受體結合趨化因子等多條免疫相關途徑參與牙周炎的病理發生。本研究首次報道了IKZF1、IRF4與中重度牙周炎患者牙齦組織炎性病變相關的研究。病毒蛋白與細胞因子及細胞因子受體的相互作用、白細胞跨內膜遷移、B細胞受體信號通路等途徑是本研究新發現的與牙周炎發病相關的信號通路。針對以上分子及通路展開深入研究，將有助于牙周炎發病機制的闡明和治療靶標的確立。

利益沖突聲明：作者聲明本文無利益沖突。