一株紅霉素耐藥百日咳鮑特菌的完整基因組分析

于丹 栗東芳 袁林 馮欣 史偉 姚開虎★

百日咳鮑特菌（Bordetella pertussis，B.pertussis）又稱為百日咳桿菌，是百日咳的病原菌。20世紀40年代開始推廣接種百日咳疫苗后，百日咳發病率明顯下降。但20世紀末，一些高疫苗接種率國家陸續報道百日咳發病率再度升高，稱為百日咳再現（pertussis re-emerge），引起全球廣泛關注［1］。百日咳再現涉及多方面的原因，其中百日咳鮑特菌的適應性變化是重要原因之一，包括對疫苗免疫選擇壓力和抗生素選擇壓力的適應性變化［1］。從不同國家報道數據看，盡管1994年第一例紅霉素耐藥百日咳鮑特菌出現在美國，但至今只有中國內地報道了紅霉素耐藥百日咳鮑特菌的廣泛流行［2］。前期，本課題組研究結果提示紅霉素耐藥百日咳鮑特菌大多具有相同的抗原基因型，均為ptxA1/ptxP1/prn1型，與之前的研究結果相一致［3］，提示耐藥菌株具有相似遺傳背景。為進一步明確耐藥菌株的基因組特征，本研究對紅霉素耐藥株P2013109 進行了完整基因組測序，并與我國疫苗株CS 等已公開的完整基因組數據比較，為疫苗研究等提供參考。

1 材料與方法

1.1 菌株

菌株P2013109 為本實驗室從百日咳患兒鼻咽拭子標本中分離。菌株的特點：藥物敏感性檢測顯示該菌株紅霉素最低抑菌濃度（Minimum inhibitory concentration，MIC）＞256 mg/L，紅霉素KB 紙片（15 μg/片）擴散法檢測藥物敏感性未見抑菌環，前期研究確定其疫苗相關抗原基因型為ptxA1/ptxC1/ptxP1/prn1/fim2-1/fim3-1/tcfA2［3］。

1.2 試劑/儀器

硅膠模型TM 基因組DNA 提取試劑盒（北京賽百盛公司），碳瓊脂培養基（CM0119，英國OXOID 公司），培養箱，4℃冰箱，-80℃冰箱，離心機。

1.3 細菌培養以及基因組DNA 提取

凍存菌株標本接種于含10%羊血的碳瓊脂培養基。培養基置于溫度為35℃的培養箱中孵育3天后觀察菌落特征和純度，轉種一次。刮取瓊脂上的菌落，采用試劑盒提取細菌基因組DNA。

1.4 基因組測序和組裝

將P2013109 的基因組DNA 委托上海美吉生物醫藥公司分別利用二代Illumina Hiseq 2000 平臺和三代Pacbio 平臺進行全基因組測序。采用SOAPdenovo v2.04 軟件對二代測序數據進行初步組裝；利用blasR 將三代測序數據比對到初步組裝結果上，并對其進行矯正與糾錯；然后利用糾正后的三代測序數據進行組裝，使用的軟件為Celera Assembler8.0。最后再次利用二代測序數據進行校驗，同時進行局部內洞填充和堿基校正。

1.5 基因組注釋

對組裝得到的基因組序列，利用Glimmer 3.02軟件進行開放閱讀框（open reading frame，ORF）預測，并與數據庫比對獲得蛋白功能注釋。利用Barrnap 0.4.2 和tRNAscan-SE v1.3.1 軟件對基因組中包含的rRNA 和tRNA 進行預測。

1.6 比較基因組分析

從NCBI 下載已完成全基因組測序的15 株百日咳桿菌及其注釋信息（見下表1），包括中國百日咳疫苗株CS［4］、該菌種最早完成測序的Tohama I菌株［5］，以及歐美發達國家已經發表的當前流行的ptxP3型菌株的基因組序列［6-7］。利用MUMmer3.23 軟件將其他基因組序列與參考基因組序列Tohama I 進行比對，得到每個基因組的單核苷酸多態性（SNPs）信息，去除落在重復序列區域的SNP 后，將在所有樣本中都出現的SNP 位點來構建系統發育樹，采用TreeBeST 軟件中的最大似然法，迭代1 000 次。用blast v2.5.0 比對軟件將其他菌株比對到菌株P2013109 上，得到每株菌在P2013109 的比對上的區域，然后將所有菌的比對結果用circos 軟件做圖。

2 結果

2.1 基因組組裝與基因組預測

對測序得到的百日咳桿菌P2013109 二代和三代數據進行質控和組裝后，得到一個完整的基因組序列，總長度4 126 010 bp，GC 含量為67.72%，預測P2013109 的基因共有4 085 個，基因的總長度為3 571 551 bp，占全基因組的86.6%，平均基因長度為874 bp。基因區GC 含量為68.4%，基因間區GC 含量為62.8%。共找到61 個tRNA 和3 組rRNA。P2013109 基因組序列已經提交到NCBI 的Genbank 數據庫，其登錄號為CP038790。與其他15 株菌的基因組信息見下表1。

表1 百日咳桿菌P2013109 與其他15 個菌株基因組的基本特征Table 1 The genome characteristics of B.pertussis P2013109 and the others

2.2 COG 聚類和KEGG 代謝通路分析

通過與String 數據庫比對，在P2013109 基因組上，共有3 388 個CDS 獲得COG 和NOG 功能注釋，根據注釋結果可以將這些功能歸為4 大類，23小類。如表2 所示。在這些注釋結果中，與代謝相關的功能最多，主要以氨基酸、無機離子、脂類和碳水化合物的轉運和代謝為主；其次是細胞過程與信號傳導類的功能，主要以細胞壁/膜/包體的生物合成、翻譯后修飾等為主；在信息儲存和處理類的功能中，以轉錄和翻譯為主。

通過與KEGG 數據庫比對，在P2013109 基因組中共找到代謝通路相關的基因2 160 個。共注釋上41 個類型，歸屬6 個大類，即細胞過程、環境信息處理、遺傳信息處理、人類疾病，代謝和有機體系統相關的功能。

2.4 比較基因組分析

2.4.1 抗原基因型和耐藥基因型確定

比對百日咳鮑特菌P2013109 抗原基因ptxA、ptxP、fim3和prn的核苷酸序列，得到P2013109 的基因組的抗原基因型結果，與前期抗原基因型測序結果一致。其他15 株菌通過搜集文獻和比對結果相結合的方式得到這些抗原基因型的分析結果，見表1 所示。細菌染色體23s rRNA 基因有3個拷貝，每個拷貝中均具有A2047G 變異。

2.4.2 全基因組比對

將其他所有菌株與P2013109 比對結果如下圖1 所示。與P2013109 相比，其他基因組主要有3 個比較大的缺失區域，R1 約4kb，編碼轉座酶；R2 是菌株H321 和I979 缺失的一個大約27 kb 的區域，是編碼百日咳毒素和IV 型分泌系統的區域；R3 區域主要包含R3-1 和R3-2 兩個子區域，除了CS 和Tohama I 菌株外，其他菌株均缺失。R3-1 主要包含BP1948-1951、BP1953 和BP1954 等6 個基因，前4個基因涉及支鏈氨基酸轉運底物結合蛋白、轉運系統滲透酶、APC 轉運體ATP 結合蛋白；BP1953 和BP1954，分別編碼氧化還原酶和單氧酶。R3-2 區域也包含4 個基因，BP1959 是IS1663 的轉座酶基因，BP1961 編碼一種細胞色素flavocytochrome，bfrI和BP1965 分別與鐵載體受體和跨膜蛋白有關。

2.4.3 系統發育分析

將這16 株菌與廣東報道的一株百日咳全基因序列一起做系統發育分析，共得到326 個SNPs。基于這些SNPs 構建系統發育樹如圖2 所示，由圖可知，這16 株菌可以分成3 個型別，TohamaⅠ和CS 這兩個疫苗株作為I 型，分別于1951年和1952年分離，與其他2 個型別有76 個SNP 的差異（23.3%），此外，TohamaⅠ與其他所有菌株比較呈現54 個SNPs（16.5%）。P2013109 自成一個分支，并且與其他所有菌株相比，呈現25 個SNPs（7.7%）。其他在2000年之后采集自國外的菌株自成一個型別，與另外兩個型別有21 個SNP 的差異（6.4%），其中這14 株菌又可以分為兩個亞型，主要區別是Ⅲa 亞型的fim3基因是2 型，而Ⅲb 型菌株是1 型。

分支旁的數字為自展1 000 次的置信值；標尺代表7%的序列差異。系統發育樹分析時納入了廣東分離株BP_GD2017，見參考文獻［8］。

圖1 P2013109 與其他15 株菌的全基因組比較結果Figure 1 The whole genome comparison analysis of B.pertussis P2013109 and the other 15 genomes

圖2 百日咳桿菌P2013109 與其他百日咳桿菌的的全基因組進化樹Figure 2 Genome-wide phylogenetic tree of Genome-wide phylogenetic tree of Bordetella pertussis P2013109 and the other genomes（including the BP_GD2017 isolated from Gangdong）

表2 百日咳桿菌P2013109 基因組的COG 功能分類Table 2 Gene distribution based on COG and NOG functional classification of B.pertussis P2013109

3 討論

本研究顯示百日咳鮑特菌P2013109 株基因組4，126，010 bp，與已經報道的百日咳鮑特菌全基因組［9-10］大小相當。該紅霉素耐藥株表型明確，全基因組序列也表明，細菌染色體23S rRNA 基因有3個拷貝，3 個拷貝中均具有A2047G 變異。23S rRNA 是核糖體RNA，核糖體是細胞內蛋白質合成的場所，紅霉素可以通過與位于23S rRNA 區域的肽鏈形成的部位（肽酰轉移酶區域）結合進而抑制蛋白質合成。23S rRNA 突變是目前百日咳鮑特菌耐藥的主要機制［11］，基因型與表型相互應證。與P2013109 株全基因組比較，其他菌株主要有3 個缺失區，尤其R1 區涉及轉座酶。這些基因區在耐藥變異過程中是否發揮作用，需要深入研究。基于完整基因組SNPs 構建系統發育樹發現，納入分析的菌株明確分為3 支。分離于20世紀50年代的百日咳流行菌株CS 株［4］和Tohama I 株形成獨立的進化分支，與當前流行菌株在進化分支上存在較大差異。

近期，廣東李柏生等［8］對廣東分離到的1 株百日咳鮑特菌BP-GD2017 進行全基因組分析報道，該菌株的基因組草圖序列數據納入基因組比較分析發現，該菌株處于國外菌株分支Ⅲ中（見圖3），它們的抗原基因型也完全相同，為ptxA1/ptxP3。P2013109 獨屬一支，抗原基因型為ptxA1/ptxP1。課題組前期檢測了99 株百日咳鮑特菌的疫苗相關基因型，ptxA1/ptxP1占91.9%（91/99）［3］。這些數據明確提示長期疫苗使用之后，百日咳致病菌株百日咳毒素相關基因型發生了明顯變化，與既往研究報告［12-13］相符。提示近年來流行的百日咳菌株的基因改變除了疫苗壓力造成相關抗原性發生適應進化，也有可能在進化分支上出現改變。本課題組前期研究表明，臨床分離菌株耐藥株抗原基因型均為ptxA1/ptxP1，而基因型ptxA1/ptxP3的分離株全部對紅霉素敏感，這提示國內外臨床分離菌株耐藥性的差異，與流行菌株的遺傳背景不同相關。

當前百日咳疫苗生產用菌株基本上都是分離自上個世紀。研究已經發現，百日咳疫苗長期使用以后，臨床分離菌株的抗原基因型與疫苗株明顯不同［12-13］。尤其在無細胞疫苗使用國家，一些毒力抗原缺失的菌株，如缺失PRN、FIM 和PTX 的菌株逐漸被發現，甚至出現流行的情況，菌株疫苗抗原的變異和缺失等適應性變化可能是百日咳再現的重要因素［14］。但課題組在前期研究中沒有發現國內百日咳鮑特菌分離株有抗原缺失，與疫苗株比較主要是抗原基因變異［3］。本研究從基因組序列證明P2013109 有抗原基因變異，無缺失，同時與疫苗株比較遺傳進化分支明顯不同。近期，Xu 等［15］分析了71 株中國分離的ptxP1菌株（其中包括46 株紅霉素耐藥株）的基因組重測序結果，基于2 744 個SNPs 位點構建遺傳進化樹，結果發現耐藥株只出現在3 個遺傳分支中，這些耐藥株只分離自中國內地，包括香港、臺灣在內的其他地區和國家公開的菌株基因組，均沒有處于同一遺傳支的分離株。

本研究第一次呈現了耐紅霉素百日咳鮑特菌ptxP1型菌株的完整基因組序列，序列分析明確顯示其在遺傳學關系上與疫苗株和國外當前主要流行的ptxP3型菌株遺傳關系較遠，屬于特有的一個遺傳分支。已有基因組重測序研究提示還需要對其他耐藥菌株進行完整基因組的構建和研究，才能反映紅霉素耐藥鮑特菌基因組的全貌，為國內將來新型百日咳疫苗的研發提供基礎信息。

致謝：廣東省疾病預防控制中心李柏生提供其課題組發表的1 株百日咳鮑特菌分離株的全基因組序列數據。