基于生物信息學篩選活動性結核差異基因和驗證*

馬慶慶，龔亞東，劉木波，韓曉靜，陳云華，唐 竹

(貴州航天醫院中心實驗室，貴州 遵義 563000)

結核病是由結核分枝桿菌感染所致，是在世界大范圍流行并嚴重威脅人類健康的傳染病。世界衛生組織(WHO)估計，2018年約有1 000萬新增結核病患者，死亡人數高達120萬人，在致死性傳染性疾病中排名第一[1]。流行病學研究表明，全球約有17億人感染結核分枝桿菌，但其中絕大部分人處于潛伏感染狀態而不發病，5%～10%的感染者最終發展為活動性結核(PTB)，在結核病發展過程中機體免疫功能發揮了關鍵作用[2-3]，但具體分子機制至今尚未明確。因此，尋找與疾病相關的信號通路有助于闡明疾病發生、發展機制。

隨著高通量測序技術和基因組學的發展，以基因芯片為代表的基因表達分析數據呈指數增長趨勢，通過基因芯片篩選出差異基因可在短時間內獲得大量與疾病相關的基因信息，在研究疾病發生、發展的相關基因和信號通路等方面發揮著重要作用[4-5]。本研究下載并分析了美國國立生物技術信息中心的基因表達數據庫(GEO)中的2個芯片數據集(GSE29536和GSE42834)[6-8]。篩選得到PTB和健康者差異基因，再對差異基因進行基因本體(GO)功能注釋和京都基因與基因組百科全書(KEGG)信號通路富集分析，構建了蛋白質-蛋白質相互作用(PPI)網絡，對PPI網絡進行核心模塊篩選和信號通路富集分析，后續參照文獻[9-10]進行逆轉錄聚合酶鏈反應(RT-PCR)驗證。探究了與PTB相關的基因和信號通路的功能變化，旨在為明確PTB發生的分子機制提供重要依據。

1 資料與方法

1.1數據來源 2021年3月在GEO數據庫(https://www.ncbi.nlm.nih.gov/geo/)中輸入關鍵詞“Active /Pulmonary Tuberculosis,blood”進行表達譜數據集搜索，篩選條件：(1)包含PTB患者和健康者樣本；(2)包含外周血樣本。經過篩選得到GSE29536和GSE42834芯片數據集，包括49例PTB患者和224例健康對照者，樣本均來源于人全血。基于GPL10558(Illumina HumanHT-12 V4.0 expression beadchip)和GPL6102(Illumina human-6 v2.0 expression beadchip)平臺，表達數據均為Expression profiling by array，種屬為Homo sapiens。

1.2方法

1.2.1微陣列數據挖掘 在數據來源中選擇2個微陣列數據集(GSE29536，GSE42834)用于鑒定差異表達基因。將P<0.01作為標準，并用R語言heatmap，ggpubr，ggthems等R包分別對2個微陣列數據集進行可視化聚類分析及火山圖繪制。

1.2.2差異基因篩選 通過GEO2R篩選差異基因，限定條件為|log2 fold change(log2FC)|>1 且P<0.01。在http://bioinformatics.psb.ugent.be/web tools/Venn/網站進行差異基因篩選。

1.2.3GO功能富集與KEGG通路分析 GO主要是用于注釋基因和蛋白質功能。包括3個生物學部分，即生物過程、細胞組件和分子功能。使用David網站(https://david.ncifcrf.gov/home.jsp)將篩選好的166個差異共表達基因進行GO功能分析與KEGG通路富集。并使用R語言中的ggplot2包進行氣泡圖可視化，應用Funrich tool軟件分析通路富集相關性。P<0.01、FDR<0.05 設定為顯著性基因富集臨界值。

1.2.4蛋白互作網絡 通過STRING 數據庫(https://string-db.org/)構建已確定的差異基因進行PPI網絡分析，然后使用 Cytoscape軟件中的cytohHubba，ClueGO等插件對網絡進行可視化和通路富集分析，篩選關鍵(hub)基因。

1.2.5RT-PCR驗證關鍵(hub)基因在PTB患者血液中的表達 驗證血液樣本系本院結核病科收治的45例PTB患者(PTB組)和15例健康者(健康對照組)。使用Paxgene Blood RNA Extraction Kit(Qiagen) 提取和純化全血中總RNA，經反轉錄后獲得cDNA；以GAPDH作為內參基因，采用FastStart Universal SYBR Green Master(Roche)試劑盒，利用Primer 3在線引物設計工具(http://primer3.ut.ee)設計各基因qPCR引物，qPCR檢測PTB組和健康對照組研究對象hub基因表達。根據2-ΔΔCT值計算hub基因在PTB患者中的相對表達量。

2 結 果

2.1差異基因篩選 GSE29536共有461個差異基因，上調基因229個，下調基因232個；GSE42834共有343個差異基因，上調基因282個，下調基因61個，見圖1A、B。上、下調基因與在線分析結果一致。見圖1C、D。2個芯片數據集差異共表達基因166個，包括114個上調基因和52個下調基因。見圖1E。

A.GSE29536 差異基因火山圖;B.GSE42834 差異基因火山圖，紅色代表上調，藍色代表下調；C.GSE29536聚類圖；D.GSE42834聚類圖；E.差異共表達基因維恩圖。

2.2GO與KEGG通路富集分析 差異基因主要富集在防御病毒、免疫反應、先天免疫反應、對病毒的反應等生物過程，見圖2A。介導蛋白質結合、腺苷三磷酸結合、相同的蛋白質結合、鈣離子結合、蛋白激酶結合等分子功能，見圖2B。富集在細胞質、胞質溶膠、質膜、細胞外區域、質膜的組成部分、線粒體等細胞組分。見圖2C。差異基因主要涉及甲型流行性感冒(流感)、單純皰疹感染、麻疹、丙型肝炎、胞質DNA傳感途徑、Toll樣受體信號通路、視黃酸誘導基因I(RIG-I)樣受體信號通路等過程。見圖2D。免疫系統信號通路與結核病患者相關性最高(44.6%)，差異有統計學意義(P=0.001)。

A.生物過程圖；B.細胞組件圖；C.分子功能圖；D.KEGG通路富集圖。

2.3蛋白互作網絡分析 信號傳導及轉錄活化因子1(STAT1)、GBP1、TRIM5、CXCL10、TLR5、CCR7、CD3E、MAPK14、CD28、CCR1等在結核病患者中異常表達。10個hub基因分別為STAT1、CXCL10、干擾素(IFN)調節因子7(IRF7)、ISG15、IFIH1、IFIT1、IFIT3、GBP1、OAS1和OAS2。見圖3C。細胞因子生成、先天免疫、病毒防御等通路及生物過程明顯交互。見圖3D。主要功能涉及調節IFN、細胞因子、免疫反應等方面。見表1。

A.蛋白互作網絡圖；B.hub基因；C.ClueGo富集分析通路餅狀圖；D.生物過程交互網絡圖。

表1 hub基因功能注釋

2.4RT-PCR驗證 STAT1、ISG15、OAS1表達明顯上調，差異均有統計學意義(P<0.01)。見圖4。

aP<0.01。

3 討 論

結核分枝桿菌能通過多種途徑操縱巨噬細胞活化并在其胞內建立適宜的增殖環境以防止被消除[11-12]；同時，采取一系列機制逃逸宿主細胞的免疫應答機制[13-14]，與其易發生耐藥[15]及宿主免疫系統的相互適應性改變等密切相關。基于基因表達數據庫挖掘相關差異表達基因和信號通路對結核病的診治具有至關重要的作用[16]。

為進一步了解PTB發生、發展的分子機制，本研究選取GEO 數據庫中與PTB密切相關的2張芯片數據集(GSE29536和GSE42834)參照文獻[17-18]的方法進行了生物信息學分析，樣本包括49例PTB患者和224例健康對照者(均為人全血樣本)，結果顯示，GSE29536數據集共有461個差異基因，其中上調229個，下調232個；GSE42834數據集共有343個差異基因，其中上調282個，下調61個。采用GEO2R在線分析上、下調基因結果得到差異表達基因166個，進一步說明PTB的發病機制是極其復雜的，可能是由眾多基因和(或)蛋白質分子等相互作用的結果[19-20]。

本研究通過GO富集分析發現，差異基因主要富集于以下幾個方面:(1)生物過程，差異基因顯著富集于防御病毒、免疫反應、先天免疫反應、對病毒的反應等;(2)細胞單位，差異基因顯著富集于細胞質、胞質溶膠、質膜、細胞外區域、質膜的組成部分、線粒體等;(3)分子功能，差異基因顯著富集于蛋白質結合、腺苷三磷酸結合、相同的蛋白質結合、鈣離子結合、蛋白激酶結合等。KEGG 信號通路主要富集在甲型流感、單純皰疹感染、麻疹、丙型肝炎、胞質DNA傳感途徑、Toll樣受體信號通路、RIG-I樣受體信號通路等。

由于傳統篩選表達差異基因的方法會掩蓋一些具有重要生物學意義但表達上調倍數較低而被掩蓋的基因，本研究繼續運用FunRich軟件分析了全部差異基因相關的信號通路，結果顯示，免疫系統信號通路與結核病患者相關性最高。通過文獻發現，本研究結果大部分可以在目前已發表的文獻中找到證據[21- 22]。結核分枝桿菌侵入機體后免疫調節功能的變化仍有待于進一步研究，而現有研究證明，固有免疫細胞如自然殺傷細胞、吞噬細胞等首先反應，啟動免疫反應機制，隨之啟動適應性免疫，共同完成免疫防御。機體抗結核感染的過程是一個精細而復雜的免疫反應過程，其相互影響，共同維持機體的正常免疫防御機制[23]。

本研究借助STRING在線數據庫和Cytoscape軟件對差異基因進行了蛋白質相互作用網格分析，篩選得到10個關鍵基因，即STAT1、IRF7、CXCL10、ISG15、IFIH1、IFIT1、OAS1、GBP1、IFIT3和OAS2。

本研究驗證結果證實，STAT1、ISG15、OAS1 mRNA在PTB患者外周血中表達明顯上調。有研究表明，在結核病感染的早期STAT1可通過磷酸化促進下游凋亡因子激活轉錄。同時，STAT1在促進巨噬細胞極化到M1極化巨噬細胞方面也非常重要，與M2極化巨噬細胞比較，STAT1可更有效地清除結核分枝桿菌[24]。STAT1可結合特定的含磷酸酪氨酸的肽段，當STAT被磷酸化時會聚集成同源的二聚體，從而參與IFN-γ引發的信號通路。當STAT進入細胞核時IFN誘導的早期基因表達被與啟動子結合的IFN-γ激活序列激活[25]。因此認為，STAT1可能在結核病免疫防御過程中具有至關重要的作用。IRF、STAT均具有共有的靶基因。有研究表明，IRF7是誘導IFN表達的最重要的調節因子之一[26]。IFN尤其是IFN-γ是與結核的發生、發展、轉歸、診治均密切相關。目前，已被用于結核病的診斷及免疫治療等方面。另外，有研究探索了CXCL10(IP-10)作為IFN-γ的替代標記物的可行性，結果證實CXCL10的診斷準確性與IFN-γ相當，且CXCL10在幼兒和細胞免疫低下的人類免疫缺陷病毒感染者中可能更為準確[27]。ISG15是IFN誘導產生的一種類泛素蛋白分子，其在抗病毒反應中具有重要作用，DOS SANTOS等[28]研究證實,ISG15與PTB的炎性反應和疾病嚴重程度密切相關，表明其可能作為結核病的生物標志物。

綜上所述，本研究通過對結核病芯片表達譜數據的挖掘分析，獲得結核病患者基因表達差異基因，并初步對差異基因進行了功能富集和探討，構建出蛋白質相互作用網絡圖，并篩選出關鍵差異基因。結合RT-PCR驗證結果也支持了部分數據分析結果，但更深入的機制仍值得探究。本研究篩選到的基因為今后更進一步研究以上生物學過程在結核病發病機制中的作用提供了新的可能和方向。