牙周病的表觀遺傳學分析：生物信息學研究

王田田 劉香瓊 李思敏

牙周炎是一種主要由菌斑生物膜導致的多因素疾病，該病可導致牙齒支持組織的進行性破壞[1]。遺傳和表觀遺傳學因素可通過調(diào)節(jié)牙周致病菌引起的宿主免疫反應，影響牙周炎的進展[2]。為闡明牙周炎的病因機制，確定牙周炎發(fā)病過程中遺傳和表觀遺傳因素的改變顯得非常必要。

非編碼RNA，作為表觀遺傳學的重要調(diào)控因素，有兩種主要類型，微小RNA（miroRNAs，miRNAs）和長鏈非編碼RNA（long noncoding RNAs，lncRNAs）。長鏈非編碼RNA，微小RNA，可與信使RNA，三者之間相互作用，形成調(diào)節(jié)性競爭性內(nèi)源性RNA網(wǎng)絡[3]。微小RNA含有22個堿基對，能夠負向調(diào)控信使RNA的表達[4]。最近使用芯片分析的研究[5，6]表明，牙周炎患者的炎癥牙齦組織和健康牙齦組織之間在微小RNA的表達譜上存在很大差異。長鏈非編碼RNA，作為非編碼RNA的一種，包含200個以上的堿基對，可以通過正向或者負向的方式來調(diào)控編碼蛋白質(zhì)的信使RNA的表達[7]。研究表明長鏈非編碼RNA在慢性牙周炎患者的炎癥牙齦組織中差異性表達[8]。

應用RNA芯片分析方法研究各種RNA在牙周炎中的表達，所得到的結果存在很大的異質(zhì)性[9]。因此，有必要通過生物信息學的技術，將這些既往研究中的所有的芯片數(shù)據(jù)集進行整合，從而預測出重要的與牙周炎發(fā)病風險有關的生物學標記物。近些年來，已有生物信息學研究[10]確定了與牙周炎發(fā)病過程相關的重要的基因、蛋白質(zhì)，以及生物學通路。然而，基于競爭性內(nèi)源性RNA理論，對牙周炎發(fā)病過程中競爭性內(nèi)源性RNA網(wǎng)絡的研究尚無報道。

本研究采用生物信息學方法，基于牙周炎的高通量RNA測序及芯片數(shù)據(jù)，通過整合長鏈非編碼RNA，微小RNA以及信使RNA的表達譜數(shù)據(jù)，旨在構建參與牙周炎發(fā)病機制的競爭性內(nèi)源性RNA（ceRNA）網(wǎng)絡。

資料和方法

數(shù)據(jù)采集與差異表達分析：從基因表達數(shù)據(jù)庫（GEO）中下載與牙周炎樣品匹配的基因表達數(shù)據(jù)和miRNA表達譜數(shù)據(jù)。將炎癥性牙齦樣品作為實驗組，健康牙齦樣品作為對照組。使用R語言的limma包中的Linear模型，對牙周炎樣本和牙周健康樣本之間差異表達的mRNA和miRNA進行差異表達分析。P值＜0.05，F(xiàn)DR＜0.05，|log2FC（倍數(shù)變化）|≥0.58的mRNA和miRNA分別被認為是顯著差異表達的基因（DEG），以及顯著差異表達的miRNA（DE-miRNA）。

1.RNA互作對數(shù)據(jù)的提取

從TarBase v7.0，miRTarBase Release 6.0，和miRecords三個數(shù)據(jù)集，提取已獲得實驗驗證的324507對miRNA-mRNA互作對。從miRcode和starBase database v2.0數(shù)據(jù)集，提取出1756對lncRNA-miRNA互作對，整合后得到DE-miRNA的mRNA靶標。

2.差異表達基因的共表達分析

為構建共表達網(wǎng)絡和確定潛在的風險基因，利用R語言進行加權基因共表達網(wǎng)絡分析。首先，獲得P值小于0.05的顯著差異表達基因的表達譜，并對基因之間的相互關系的鄰接矩陣進行計算。基于公式（1），對基因互作對間的連接系數(shù)進行計算。

χi，χj：基因i和基因j的表達載體；cor：兩個基因載體的皮爾遜相關系數(shù)；β：皮爾遜相關系數(shù)的指數(shù)變換，將皮爾遜相關系數(shù)進一步覆蓋到連接系數(shù)αij。

基于網(wǎng)絡的拓撲性質(zhì)，通過加權基因共表達網(wǎng)絡分析，獲得共表達網(wǎng)絡的模塊。該網(wǎng)絡分析了兩種互作的基因以及與其互作的其他基因之間的關系。根據(jù)以下公式（2），這兩種互作的基因將連接系數(shù)αij覆蓋到加權系數(shù)wij上?；诩訖嘞禂?shù)，獲得了顯著的基因模塊。

wij：代表與基因i相互作用的基因，與基因j相互作用的基因兩者之間的重復。

3.共表達網(wǎng)絡中風險基因模塊的確定

為了篩選出特定的基因模塊，構建了評分分類器。對每個模塊進行差異表達分析，分析了在牙周炎標本和牙周健康標本之間差異表達的每個基因的表達水平，P值被設定為整合表達評分。此外，我們還計算了模塊中每個基因的度（degree）與其-log10（P值）之間的關系。之后，可以篩選出明顯富集的基因模塊。

4.牙周炎ceRNA網(wǎng)絡的構建和分析

從TarBase，miRTarBase和miRecords三個數(shù)據(jù)庫中獲取到差異表達的9498對mRNA-miRNA相互作用關系對。從miRcode和starBase v2.0兩個數(shù)據(jù)庫中下載了1726對DEmiRNA-lncRNA的互作對。通 過 將DEG-mRNA-miRNA互 作 對 和DEmiRNA-lncRNA互作對整合到一起，構建了競爭性內(nèi)源性RNA（ceRNA）網(wǎng)絡。

5.分析結果預測精確度的評估

為了評估本研究中分析方法的預測精確度，我們應用了含有核函數(shù)的支持向量機（SVM）。通過計算受試者工作特征（ROC）曲線（AUC）下的面積，預測模型的診斷能力得以評估。通過使用R包中的ROCR.61，整合假陽性率和敏感度的數(shù)據(jù)，繪制出ROC曲線。此外，為了可視化地顯示分類器的性能，我們進行了無監(jiān)督的分層聚類分析。通過將完整的聯(lián)結和皮爾遜相關性作為距離測量，實現(xiàn)表達譜的分層聚類。

結 果

1.差異表達的mRNA，miRNAs和lncRNAs的確定

通過進行差異表達分析，我們確定了704個差異表達基因和47個差異表達的miRNA。通過使用TarBase，miRTarBase和miRecords這三個數(shù)據(jù)庫，我們獲得了3282對差異表達的miRNA與其靶基因的互作對，這些互作對是由37個差異表達的miRNA和2650個靶基因構成的。

2.共表達分析

通過差異表達分析，獲得P值＜0.05的基因并對其進行共表達分析。為了建立無標度網(wǎng)絡，我們選擇了相鄰矩陣（β）的權重，此外我們還對不同β值（范圍：1-20）的基因的相鄰矩陣的相關系數(shù)和平均值進行了計算。當軟權重為12時，無標度R∧2趨于穩(wěn)定。

3.顯著模塊的確定

在核查已構建網(wǎng)絡的參數(shù)后，通過對加權共表達網(wǎng)絡進行建模，將基因分為不同模塊。從加權關聯(lián)網(wǎng)路中，檢索出11個具有不同剪切靜態(tài)顏色的基因模塊。隨機抽取出100個基因之后，對所有基因模塊的基因相關性進行計算。每個基因的度k和-log10（p）之間的相關性被加權，加權結果顯示這四個模塊（綠色，紫色，藍綠色，黃色）中的基因與牙周炎的發(fā)病機制顯著相關（P＜0.05，圖1b-e）。

4.ceRNA網(wǎng)絡的構建

在所有11個模塊（圖1a）中，綠色和藍綠色模塊的顯著性最高，因此被用于構建ceRNA網(wǎng)絡（圖2a，b）。綠色模塊包含2種DE-miRNA（miR-125a-3p和miR-142-3p）和三種DEGs（HSPA4L，EVI2B和PANK3）（圖2a）。藍綠色模塊含有3種DE-miRNAs（miR-142-3p，miR-125a-3p和miR-200a），和三種DEGs（YOD1，ITGAL，和CTNNBIP1）（圖2b）。在ceRNAs網(wǎng)絡中，3種lncRNA（TUG1，F(xiàn)GD5-AS1，MALAT1）能夠同時直接靶向兩種微小RNA，包括miR-125a-3p和miR-142-3p。沒有一種lncRNA能夠同時直接靶向miR-125a-3p和miR-200a。

圖1 共表達網(wǎng)絡中的模塊

圖2 綠色模塊（a）和藍綠色模塊（b）的競爭性內(nèi)源性網(wǎng)絡

5.使用兩種額外的獨立隊列進行進一步的基因驗證

為了獨立驗證我們預測出來的這6種基因的診斷能力，基于GPL6244這一平臺，我們在GSE33774這一數(shù)據(jù)集（含有7個牙周炎樣品和8個牙周健康對照樣本）中，對這些基因的鑒別性能進行了檢測。基于這六種基因的表達譜，我們進行了SVM分析。分析結果（圖3）顯示，HSPA4L基因（被hsa-miR-125a-3p所靶向）的ROC為0.9643，是這6個基因ROC中最高的。這表明，HSPA4L是在牙周炎樣本與牙周健康對照樣本中最為顯著性差異表達的基因。

圖3 ROC曲線（a），預測出的六種風險基因的熱圖（b）

討 論

本研究圍繞三種類型的RNA（信使RNA，微小RNA，長鏈非編碼RNA），對牙周炎的遺傳和表觀遺傳學機制進行剖析。六種信使RNA（HSPA4L，CTNNBIP1，ITGAL，YOD1，PANK3，和EVI2B）被 預測為牙周炎的生物標記物。HSPA4L，熱沖擊蛋白家族的一員，在牙周炎發(fā)病機制中發(fā)揮雙重作用：一方面可以維持牙齦健康，防止受到口腔微生物的破壞[11]；另一方面，HSPA在某種條件下會引發(fā)牙齦炎，加重后向牙周炎轉化，進而增加牙槽骨的損失[12]。CTNNBIP1，又名β-鏈互作蛋白1，可抑制核內(nèi)β-鏈蛋白的轉錄能力[13]。核內(nèi)β-鏈蛋白的激活能夠通過抑制非經(jīng)典的Wnt通路，促進牙周韌帶干細胞的增殖，并導致其成骨分化的減弱[14]。ITGAL，整合素亞基αL，可編碼整合素淋巴細胞功能相關抗 原-1（LFA-1/CD11a）的 整 合 素αL鏈[15]。LFA-1/CD11a是一種重要的免疫調(diào)節(jié)細胞因子，可通過介導T細胞的召集和保留，調(diào)節(jié)牙周炎的免疫過程[16]。YOD1，泛蛋白酵素，在牙周炎的樣本中表達下調(diào)，可被miR-30b所靶向。MiR-30家族能夠負向地調(diào)控骨髓間充質(zhì)干細胞的骨形態(tài)蛋白-2引起的成骨分化[17]。PANK3，泛酰酸激酶3，已被證明可被miR-142-3p所靶向；MiR-142-3p能夠調(diào)節(jié)脂多糖引起的牙周炎癥[9]。EVI2B，同向性病毒結合位點2B，已被證明在多種微生物感染的大鼠骨模型中表達量增加[18]。

三種微小RNA（miR-125a-3p，miR-200a，和miR-142-3p）被預測可作為牙周炎的風險生物標記物。miR-125a-3p來源于pre-miR-125a的3'末端，在牙周炎樣本中表達上調(diào)。MiR-125a的上調(diào)可通過調(diào)節(jié)TRAF6/NFATc1/miR-125a調(diào)節(jié)反饋回路，抑制破骨細胞的生成，進而抑制牙周炎的進展[19]。miRNA-200a通過減少牙槽骨骨髓基質(zhì)細胞（BMSCs）的成骨分化，從而抑制牙周骨再生[20]。MiR-142-3p是介導破骨細胞死亡的RANKL的獨立誘導物，可通過調(diào)節(jié)Wnt通路促進成骨細胞的分化[21]。

三種長鏈非編碼RNA（MALAT1，TUG1，和FGD5-AS1）被預測可作為牙周炎的風險生物標記物。MALAT1，肺腺癌轉移相關轉錄子1，能夠與轉錄因子NF-κB相互作用，調(diào)節(jié)脂多糖引發(fā)的炎癥反應[22]。NF-κB信號通路可以調(diào)節(jié)RANK配體引起的破骨細胞生成，進而引起牙周炎的骨吸收[23]。TUG1，牛磺酸上調(diào)基因1，能夠通過促進轉化生長因子TGF-β，從而抑制細胞的增殖[24]。TGF-β信號通路在牙周炎的發(fā)病過程中，參與了牙齦成纖維細胞與牙齦上皮細胞之間的相互作用[25]。FGD5已被證明能夠在牙周炎癥進展過程中，通過調(diào)節(jié)血管內(nèi)皮生長因子VEGF的促血管生成作用，來促進血管網(wǎng)絡的延伸[26]。

結論：六種mRNA（HSPA4L，PANK3，YOD1，CTNNBIP1，EVI2B，和ITGAL），三 種miRNA（miR-125a-3p，miR-200a，和miR-142-3p），三種lncRNA（MALAT1，TUG1，和FGD5-AS1）參與了牙周炎的競爭性內(nèi)源性RNA網(wǎng)絡。這些結果有助于解釋牙周炎的發(fā)病機制，為未來研究提供了潛在的治療靶標。