活動性肺結核患者與潛伏感染者差異表達基因的生物信息學分析

楊秉芬 翟斐 安紅娟 曹志紅 劉艷華 王若 程小星

結核病作為世界性傳染病已嚴重威脅人類健康，世界衛生組織(WHO)估計，2017年約有1000萬新增結核病患者，死亡人數高達130萬人，在致死性傳染性疾病中排名第一[1]。流行病學研究表明，全球約有17億人感染結核分枝桿菌，但其中絕大部分人處于潛伏感染狀態而并不發病，約有百分之五到十的感染者最終會發展為活動性結核，在結核病的發病和發展過程中，人體的機體免疫功能發揮了關鍵作用[2]，但是目前并不是十分清楚具體的分子機制。本文通過轉錄組測序，獲取活動性肺結核患者和潛伏感染者PBMCs差異表達基因，并對這些基因進行生物信息學分析，期望找到活動性肺結核患者與潛伏感染者免疫功能差異的重要分子，為結核病的發病機制研究尋找新的突破點。

資料與方法

一、病例收集

患者樣本均來自解放軍總醫院第八醫學中心結核病研究所住院治療的患者，通過影像學檢查與痰菌培養確診為活動性肺結核患者，并排除其他免疫性疾病或使用免疫抑制劑，共6例，男性2例，女性4例，平均年齡(36.0±4.34)歲。潛伏感染者樣本選自總醫院第八醫學中心體檢中心的體檢人員，無發熱、咳嗽和咳痰等癥狀，影像學檢查肺部無病灶，通過干擾素釋放實驗，篩選出潛伏感染者，共6例，男性2例，女性4例，平均年齡(30.0±2.59)歲。

二、分離外周血單個核細胞(PBMCs)

采集EDTA抗凝外周血2 mL，按照1∶1的比例使用1640基礎培養基 稀釋，然后使用等比Ficoll(GE Biosciences, USA)，1000 g，18 min密度梯度離心，中間細胞層為單個核細胞，收集后使用1640基礎培養基洗兩遍后備用。

三、干擾素釋放實驗

將PBMCs懸液細胞濃度調整為2.5×106/mL，預包被γ-干擾素抗體的ELISPOT板每孔加入100μL細胞懸液，用結核分枝桿菌特異混合抗原多肽A和B進行刺激，設置陰性對照孔(不加任何抗原)和陽性對照空(PHA)，設置孵箱溫度37 ℃、5%CO2培養20 h。用去離子水裂解細胞，扣干水后再用洗滌液洗5遍，拍干后加又入親和素標記的γ-干擾素抗體，放入培養箱37 ℃孵育60 min。洗去未結合的抗體，然后加入HRP-標記的鏈霉素，37℃孵育60 min。洗去未結合的鏈霉素，拍干，加入酶底物顯色，室溫避光顯色，觀察斑點形成，用去離子水終止反應，晾干后用免疫斑點計數儀(Cellular Technology Ltd, USA)進行圖像和斑點計數。如果陰性對照孔斑點數為0～5點，用測試孔的斑點數減去陰性對照孔斑點數≥6；或者當陰性對照孔斑點數≥6點時，測試孔斑點數≥2倍陰性對照孔斑點數，以上兩種情況為陽性反應，可判斷為結核菌感染。

四、RNA提取

將細胞用無血清培養基AIM Ⅴ進行重懸，加入到96孔細胞培養板中，按終濃度10 μg/mL加入結核分枝桿菌H37Rv裂解液，37℃、5%CO2培養箱培養過夜。使用TRIzolTMReagent(Invitrogen，USA)進行RNA提取：每個樣本加入1 mL TRIzol，0.2 mL氯仿，劇烈振蕩15 s，室溫放置3 min；12 000 g 4 ℃離心15 min，將上層無色水相分離，并加入異丙醇0.5 mL，上下顛倒混勻，室溫內靜置10 min；12 000 g 4℃離心10 min，棄上清，使用1 mL無核酸酶水稀釋的75%乙醇清洗RNA沉淀；7500 g 4℃離心10 min，棄上清，RNA沉淀置于冰上晾干(5～10 min)，用20μL無核酶水溶解RNA沉淀。使用Agilent 2100對RNA樣本進行濃度、純度和完整性檢測，挑選合格的樣本，并根據測序樣本總量要求混合樣本，每個混合樣本包含2例女性和1例男性。

五、RNA測序及數據分析

委托華大基因對活動性肺結核患者及潛伏感染者PBMCs進行轉錄組測序，在得到測序數據raw reads后先去除原始序列中帶有的adaptor序列、N的比例大于5%和低質量的reads，獲得clean reads，再使用SOAPaligner/SOAP2工具來做參考序列比對。通過統計基因在不同位置的reads數獲得reads的分布隨機性。基因表達量的計算使用RPKM法(reads per kilobase transcriptome per million mapped reads)，根據基因的RPKM值計算該基因在不同樣本間的差異表達倍數。對差異檢驗的P值作多重假設檢驗校正，通過控制FDR(False Discovery Rate)來決定p-value的域值。在本次分析中，差異表達基因定義為FDR≤0.001，且倍數差異在2倍以上的基因，共獲得差異表達基因98個，其中上調基因67個，下調基因31個。

六、GO及KEGG分析

使用DAVID 6.8(the Database for Annotation, Visualization and Integrated Discovery, https://david.ncifcrf.gov/)對差異表達基因進行GO分析：細胞組分(cell complement)、分子功能(molecular function)和生物進程(biology process)[3]。使用Catespace插件GoClue分析差異表達基因的KEGG通路富集[4]。

七、蛋白相互作用網絡分析

使用STRING(Search Tool for the Retrieval of Interacting Genes, http://string db.org/)數據庫對差異表達基因進行蛋白相互作用分析，篩選出綜合得分(combined scores)大于或等于0.4的相互作用蛋白用于構建蛋白相互作用網絡[5]，使用軟件Cytoscape 3.7.1(www.cytoscape.org)將相互作用網絡繪制出來[6]。用Cytoscape插件MCODE(Molecular Complex Detection, http://apps.cytoscape.org/apps/mcode)對蛋白相互作用網絡進行蛋白復合物分析[7]，設定degree cut-off為3，node score cut-off為0.2，k-score為2，maximum depth為100。

結 果

一、差異表達基因的聚類分析

根據華大基因測序結果，對活動性肺結核患者和潛伏感染者差異表達基因使用MeV4.9.0(Multi Experiment Viewer 4.9.0)軟件進行聚類分析[8]，結果如圖1所示，67個基因在活動性肺結核患者上調顯著，31個基因下調顯著。

圖1 差異表達基因的聚類分析

二、差異表達基因的GO分析和KEGG分析

使用DAVID 6.8對差異表達基因進行GO分析，根據富集分數——-log10(p-value)進行富集條目排序，篩選出富集程度最高的10個條目，結果如圖2所示藍色條形表示“生物進程”富集，在“免疫反應”、“天然免疫反應”和“中性粒細胞趨化”3個詞條富集顯著；綠色條形表示“細胞組分”富集，在“細胞外空間”、“胞外區”和“外泌體”3個詞條富集顯著；紅色條形表示“分子功能”富集，在“絲氨酸型內肽酶活性”、“糖結合”和“RAGE受體結合”3個詞條富集顯著。使用Catespace插件GoClue分析差異表達基因的KEGG通路富集如圖3所示，在紫色圓點所示“自然殺傷細胞介導的細胞毒作用”、“抗原處理和遞呈”，綠色圓點所示的“IL17信號通路”，青色圓點所示的“鏈球菌感染”4個信號通路顯著富集。其中細胞外細胞因子γ-干擾素(interferon-gama, IFNG)，膜蛋白整合素亞單位αM(intergrin subunit alpha M, ITGAM)，胞內蛋白基質金屬蛋白酶1(matrix metallopeptidase 1, MMP1)和轉錄因子FOS連接多個信號通路，是重要的分子，其功能有待進一步研究。

圖2 差異表達基因的GO分析 藍色條形表示“生物進程”富集，綠色條形表示“細胞組分”富集，紅色條形表示“分子功能”富集

三、差異表達基因的蛋白相互作用網絡分析

使用STRING數據庫對差異表達基因進行蛋白相互作用分析，篩選出綜合得分大于或等于0.4的相互作用蛋白，用于構建蛋白相互作用網絡，使用Cytoscape 3.7.1軟件繪制相互作用網絡，用NetworkAnalyzer分析網絡拓撲學[9]，結果如圖4所示：STRING分析共有95個節點(dots)、423對連接(edges)構成網絡，平均聚集程度(degree)為8.91。節點的顏色表示基因變化：紅色是上調表達基因，其中防御素α1(defensin alpha1, DEFA1)、葉酸受體γ(folate receptor gama, FOLR3)、肽聚糖識別蛋白1(peptidoglycan recognition protein 1, PGLYRP1)和觸珠蛋白(haptoglobin, HP)上調倍數(fold change)的大于32，綠色是下調表達基因，其中KIR3DL1、KIR3DS1和KIR2DL4下調較明顯，這些基因與NK細胞殺傷作用相關。節點的大小表示節點聚集程度：其中中性粒細胞表達彈性蛋白酶(elastase, neutrophil expresse，ELANE)(degree=32)，ITGAM(degree=31)，基質金屬蛋白酶9(matrix metallopeptidase，MMP9)(degree=30)和 精氨酸酶1(arginase 1，ARG1)(degree=26)聚集程度最高。連接線的粗細表示相互作用的綜合得分。

圖3 差異表達基因的KEGG分析 不同顏色表示不同的KEGG信號通路富集，大的圓點是信號通路，由大到小表示富集程度的變化，最小的圓點是代表性基因

圖4 差異表達基因的蛋白相互作用網絡分析 圓點的顏色表示基因表達差異，紅色表示基因表達上調，綠色表示基因表達下調；圓點的大小表示基因在分子網絡中的聚集程度；連接線的粗細表示相互作用的綜合得分

圖5 差異表達基因的蛋白復合物分析 使用Cytoscape插件MCODE對蛋白相互作用網絡進行蛋白復合物分析獲得A、B、C三個復合物。

四、差異表達基因的蛋白復合物分析

使用Cytoscape插件MCODE對蛋白相互作用網絡進行蛋白復合物分析，設定degree cut-off為3，node score cut-off為0.2，k-score為2，maximum depth為100，結果如圖5所示，共有3個蛋白復合物，其中復合物1由上調表達基因PGLYRP1、FOLR3和HP等組成，復合物2上調表達基因CD177、MS4A3和CEACAM1等組成，復合物3由下調表達基因KIR3DL1、KIR3DL2和KIR3DL3等組成，復合物3組成基因與NK細胞介導的殺傷作用相關，說明活動性肺結核患者NK細胞可能比潛伏感染者有顯著功能異常，提示我們可以研究NK細胞的功能在結核病發病中的作用。

討 論

機體的免疫功能在結核病的發病和發展中起了關鍵的作用[2]，但是具體的分子機制目前并不是十分清楚。生物信息學的發展，特別是一些數據庫的開放為科研人員提供了平臺，科研人員可以結合自己的關注點進行分析，生成的可視化的結果能更直觀地找到新的研究切入點。本文目的旨在通過轉錄組測序，獲取活動性肺結核患者和潛伏感染者PBMCs差異表達基因，并對這些基因進行生物信息學分析，結果顯示GO分析生物進程富集的詞條是“免疫反應”、“天然免疫反應”和“中性粒細胞趨化”。KEGG分析結果顯示富集詞條是“自然殺傷細胞介導的細胞毒作用”、“抗原處理和遞呈”和“IL17信號通路”等信號通路，其中IFNG、ITGAM、MMP1和FOS基因連接多個信號通路。構建蛋白相互作用網絡，篩選到聚集程度高的分子ELANE、ITGAM、MMP9和ARG1。本文測序使用的PBMCs在體外經過結核抗原刺激，活動性結核病患者和潛伏感染者對結核抗原的反應性得到放大，以上的結果提示結核病患者PBMCs在天然免疫如“自然殺傷細胞”和“抗原遞呈”與潛伏感染者存在顯著差異。篩選到重要的節點蛋白，其中IFN-γ在抗結核免疫具有重要的作用，IFN-γ和TNF-α釋放實驗可用于鑒別潛伏感染和活動性肺結核，IFN-γ促進結核菌感染的巨噬細胞形成脂滴，加強宿主的保護性免疫反應[10-11]，MMP1和MMP9在結核病的肺損傷中具有重要的作用，可作為結核病診斷的標志物，機制研究發現miR-223可以下調MMP1和MMP9的表達[12-14]。在寄生蟲和結核菌共感染模型，ARG1表達上調引起炎癥增強和結核病情加重，Ⅰ型干擾素可以通過誘導NOS2和抑制ARG1起到抗結核的作用，又有研究表明抑制肉芽腫中缺氧巨噬細胞ARG1表達可導致肉芽腫病理變化增強和荷菌量增加[15-17]，這些基因的差異與文獻報道相符，說明轉錄組測序和信息學分析結果的可靠。而我們發現新的關鍵分子ELANE、ITGAM和FOS，其在結核病免疫中的作用知之甚少，為我們下一步研究提供依據。而蛋白復合物分析結果顯示由上調基因組成的復合物1和復合物2、由下調基因組成的復合物3，值得注意的是模塊3中KIR3DL1、KIR3DL2和KIR3DL3等基因與NK細胞介導的殺傷作用相關，最近多中心研究表明NK細胞亞群在活動性結核病患者和潛伏感染者存在顯著差異，并且能預示疾病的進展和抗結核治療的效果[18]，我們分析獲得的功能分子KIR3DL1、KIR3DL2和KIR3DL3為NK細胞在結核病發病中的作用機制研究提供依據。