豬源致病性大腸桿菌基因組比較與毒力因子分析

吳莉丹 冉雪琴 牛熙 黃世會 李升 王嘉福

（1. 貴州大學生命科學學院 農業生物工程研究院 山地植物資源保護與種質創新教育部重點實驗室，貴陽 550025；2. 貴州大學動物科學學院，高原山地動物遺傳育種與繁殖教育部重點實驗室，貴陽 550025）

大腸桿菌（Escherichia coli）是腸桿菌科埃希氏菌屬的革蘭氏陰性菌，是動物腸道最主要的兼性厭氧菌，通常對宿主無害。由于大腸桿菌基因組具有高度的可塑性，通過遺傳變異可能獲得和丟失一些基因［1］，導致特定的血清型對動物致病，豬尤其易感［2］。豬大腸桿菌病在全球各地普遍存在，各年齡段均可發病，尤其是新生幼豬最易感，常給豬場帶來極大的經濟損失。豬源致病性大腸桿菌主要的毒力因子包括黏附系統、鐵吸收和轉運系統、分泌系統和毒素等［3］，多種毒力因子相互協調作用，決定疾病的發生和發展進程。豬源致病性大腸桿菌首先通過黏附系統和鐵吸收轉運系統，黏附并定植于豬的小腸腔內表面，大量繁殖并產生毒素，再通過分泌系統將毒素分泌到細菌外，毒素被宿主吸收，損傷胃腸道血管和皮下組織，導致黃白痢、水腫病［4］。豬源致病性大腸桿菌的毒力不同，臨床癥狀輕重不一［5］。

毒力強弱與毒力因子的數量、結構和功能密切相關。采用比較基因組學方法，比較多殺性巴氏桿菌（Pasteurella multocida, PM）的強毒株DY120818與4 株弱毒株，表明毒力的差異可能與整合結合元件ICEPmu2 有關［6］；豬源腸外致病性大腸桿菌強毒菌株PCN033 和弱毒菌株PCN061 的基因序列分析提示，強毒菌株基因組中特有水平轉移DNAs、特有編碼黏附素、脂多糖等毒力因子基因，且擁有較多的前噬菌體序列［7］。為了系統解析本地流行的致豬腹瀉大腸桿菌的致病機制，本研究以不同毒力的6 株豬源致病性大腸桿菌為研究材料，通過全基因組測序，挖掘基因組中致病相關的毒力因子基因，分析毒力因子基因發生的變異，為大腸桿菌病的診斷和防治提供理論依據。

1 材料與方法

1.1 材料

大腸桿菌弱毒菌株P32、P111 和中等毒菌株P555、P211，由本實驗室從貴州某豬場腹瀉仔豬糞樣中分離獲得［8］；強毒菌株S10670（菌種編號：CICC 10670）、E24190（ 菌 種 編 號：CICC 24190）購自中國工業微生物菌種保藏管理中心。

1.2 方法

1.2.1 序列測定 應用LB 培養基培養各株細菌，離心收集菌體，利用天根DNA 試劑盒提取細菌基因組DNA，構建DNA 文庫，應用Illumina Hiseq 2500 測序平臺（華大基因）進行基因組測序（雙端，2×150bp）。此外，從NCBI 數據庫下載E. coli. str.K?12 的基因組序列作為參考。

1.2.2 基因組序列組裝與基因注釋 利用軟件FastQC［9］評估原始測序數據的質量，以軟件FastP［10］對原始測序數據中的接頭和低質量reads 進行過濾；通過Spades［11］軟件對過濾后的數據進行組裝，GapClose［12］軟件對未知堿基進行填充；采用SeqKit［13］軟件過濾掉組裝序列中小于500 bp 的小片段，Quast［14］軟件對組裝序列進行評價，通過多次優化，得到6 株大腸桿菌的基因組組裝序列。應用Prokka［15］軟件對組裝的基因組序列進行注釋，得到編碼基因和非編碼RNAs，使用PhiSpy［16］軟件預測組裝基因組中的前噬菌體區域，通過在線網 站CRISPRFinder（https://crispr.i2bc.paris?saclay.fr/Server/）分析組裝基因組中的CRISPR 序列。

1.2.3 進化類型分析 參照文獻［17］合成chuA、yjaA和TspE4.C2基因的特異性引物（上海生工）。采用煮沸法提取受試菌株的基因組DNA 為模板，進行三重PCR 擴增。PCR 反應體系（20 μL）：2×TaqPCR MasterMix 10 μL；上、下游引物（10 μmol/L）各0.5μL；DNA 模板1 μL；ddH2O 為6 μL。PCR 反應條件：94℃ 5 min；94℃ 45 s，退火60℃ 30 s，72℃18 s，30 次循環；72℃ 5 min；4℃保存。取3 μL 擴增產物經1.5%瓊脂糖凝膠電泳檢測，據檢測結果鑒定大腸桿菌的系統進化群。用Abricate 軟件將基因組組裝序列比對到數據庫ECOH［18］（https://github.com/katholt/srst2/tree/master/data），據基因序列相似性鑒定菌株的血清型。利用mlst［19］軟件對組裝的基因組實施多位點序列分型（multilocus sequence typing，MLST）。

1.2.4 毒力因子基因的篩選與驗證 將測序菌株的基因集，應用Abricate 軟件與VFDB［20］數據庫（http://www.mgc.ac.cn/VFs/main.htm）的核心序列集進行比對（BLAST 參數：e?value 10-5，最小覆蓋度80%，最小相似度80%），得到毒力因子候選基因。從基因組注釋文件中提取各菌株分泌系統基因簇，繪制分泌系統的基因簇分布圖，比較不同菌株各類分泌系統的異同。根據基因組與毒力因子數據庫的比對結果，從每種毒力因子基因類型中挑選代表性的毒力因子基因（共23 個）進行PCR 驗證（表1），基因組DNA 提取、PCR 反應和瓊脂糖凝膠電泳同1.2.3。

表1 毒力因子基因檢測引物Table 1 Primers for the detection of virulence factor genes

1.2.5 具有高影響變異的毒力因子基因的篩選與驗證 將質控后各菌株的測序數據通過軟件BWA 比對到參考基因組E. coli. str.K?12 上，應用SAMtools和GATK 軟件檢測菌株基因組中的單核苷酸多態性（SNP）和短序列插入缺失突變（InDel），以SnpEff［21］軟件對SNPs 和InDels 影響的基因進行注釋，從中篩選出具有高影響變異的毒力因子基因。在毒力因子基因保守區域設計引物進行PCR 擴增（表1）。擴增產物進行雙向測序（上海生工），用DNAMAN［22］軟件對基因序列進行拼接，根據測序圖譜校正拼接序列，將各基因的拼接序列與參考基因組進行比對，分析基因中發生的堿基變異。

2 結果

2.1 基因組序列組裝與注釋

經拼接（表2），毒力不同的6 株大腸桿菌最終獲得108-191 條scaffolds；N50 為91-161 kb，確保了基因組組裝的質量。菌株的基因組全長為4.62-5.3 Mb，對參考基因的覆蓋率約為90%，其中菌株S10670 對參考基因組的覆蓋率最低（87.45%），P32的覆蓋率最高（91.98%）。

表2 菌株基因組組裝質量分析Table 2 Quality analysis of the assembly genomes in E. coli strains

經Prokka 注釋（表3），從6 株大腸桿菌中鑒定出4 308-5 040 個CDS 和3 364-3 557 個編碼基因，普遍高于參考基因組的數量。參考菌株有88個tRNA 編碼基因，檢測菌株的tRNA 基因數目為83-98 個，預測到菌株的rRNA 基因顯著少于參考菌株。菌株S10670 的前噬菌體數目最多（3 個），菌株E24190 和P555 基因組中沒有預測到前噬菌體序列。P211 的CRISPR 序列最多，其余菌株基因組中可能存在的CRISPR 序列數量高于參考菌株。

表3 基因組注釋分析Table 3 Genome annotation

2.2 進化類型分析

通過三重PCR 檢測和序列相似性分析，菌株S10670 屬于系統進化D 群，其余菌株屬于A 群。各菌株的血清型不同（表4），基于7 個管家基因的等位基因序列比對得到各菌株的多位點序列類型（sequence typing, ST）和克隆復合體類型（clonal complex, Cplx），E24190、P555、P32 為ST10 Cplx，P211 是一種新的ST 型。

表4 系統進化群、血清型及MLST 分型Table 4 Phylogenetic groups, serotypes and MLST typing

2.3 毒力因子基因的篩選與驗證

6 株大腸桿菌的毒力因子基因相比（表5），菌株之間毒力因子基因數量為112-280 個，基因數量相差1 倍以上，強毒菌株的毒素基因數量明顯多于中等毒菌株和弱毒菌株。除強毒菌株S10670 外，其余菌株的鞭毛相關基因占總毒力因子基因的一半以上。菌株P211 的鐵攝取和生物膜相關基因較多。各菌株的分泌系統基因數量差異較大，弱毒菌株也可以擁有豐富的分泌系統基因。

表5 菌株基因組中毒力因子分析Table 5 Virulence factors in genomes of strains

6 株致病菌分泌系統的差異集中于T2SS 和T6SS（圖1）。強毒菌株中預測出全部的T2SS 基因簇（11-13 個基因），參考基因組和弱毒菌株P32 中的T2SS基因簇不完整，2 個中等毒菌株中未檢測到T2SS 基因。大腸桿菌三型分泌系統（E. colitype III secretion systems, ETT）分為ETT1 和ETT2，其中ETT1 型只在強毒菌株S10670 基因組中分布；菌株S10670、P32 和P211 含有ETT2 基因，其中富含假基因，其余菌株不含ETT2 基因。菌株S10670 和P555 中預測出全部的T6SS 基因簇（19-20 個基因），E24190和P32 的T6SS 基因較少，只有3-6 個，其余菌株未檢測到T6SS 基因。

圖1 六株菌的分泌系統基因簇比較Fig. 1 Comparison of the secretion system gene cluster from six strains

從各類毒力因子基因中挑選代表性的23 個毒力因子基因，采用PCR 技術進行驗證（圖2，表6）。從強毒菌株S10670 基因組中檢測到T2SS 的外膜管道蛋白基因etpD、ETT1 的緊密黏附素基因eae、ETT2 中的轉錄因子eivF和脂蛋白前體基因eprK、T3SS 的效應物espX1和espO1?2基因、T6SS 溶血素共溶蛋白家族基因aec16、血紅素吸收系統基因chuA、T5SS 的ehaA和ehaB基因、志賀樣毒素基因stx2A和stxB、侵襲素EaeH。強毒菌株E24190 基因組中檢測到T2SS 的GspD基因和六聚體ATP 酶基因GspE，其中GspD與etpD基因同源；含有T3SS 的效應物基因espX1，T6SS 的溶血素共溶蛋白家族基因hcp1、自分泌黏附素基因aatB、耐熱腸毒素和不耐熱腸毒素基因、腸聚集黏附性耐熱毒素基因east1，溶血素基因hlyA、侵襲素基因EaeH。

圖2 23 個毒力因子基因PCR 擴增分析Fig. 2 PCR amplification of 23 virulence factor genes

表6 23 個毒力因子基因驗證匯總分析Table 6 Summary of 23 validated virulence factor genes

中等毒菌株P555 基因組中檢測出espX1、aec16、hcp1、EaeH和ehaB基因；P211 中含有eprK、espX1、鐵載體內膜轉運系統基因sitA、ehaB基因和CFA/I 菌毛基因cfaB。弱毒菌株P32 基因組中含有GspE、eprK、hcp1、EaeH和ehaB基因；而P111基因組中只有espX1和ehaB兩個毒力因子基因。上述PCR 檢測結果與基因組測序和注釋分析的毒力因子基因一致。

2.4 毒力因子基因變異分析與驗證

對6 株致病性大腸桿菌的基因組進行SNP 和InDel 變異分析（表7），從中檢測到豐富的SNP 位點，菌株S10670 擁有的SNP 最多。除菌株P111 外，其他菌株基因組中的缺失突變比插入突變數量多。進一步分析變異對基因的影響，得出受變異影響大的毒力因子基因只有7 個，分別是T3SS 效應因子基因espY1、espL4和espX4，Fim 鞭毛基因FimD和FimI，轉錄因子基因rpoS和尿囊素酶基因allB。采用PCR 和Sanger 測序（表8），證實菌株S10670 和P555 的espY1基因發生了一個移碼突變（圖3?A）；菌株E24190、P211、P32、P111 的espL4基因都缺失了612 bp（圖3?B），菌株P111 的espL4基因中還有一段超過900 bp 的插入序列，菌株P555 的espL4基因中檢測到一段不同于參考基因的序列；菌株P555 的espX4基因也檢測到一段超過900 bp 插入序列，S10670 的espX4基因還有兩個移碼突變（圖3?C）；E24190 中FimD基因檢測到一個無義突變（圖3?D）；菌株P32 的FimI有一個移碼突變（圖3?E）；菌株S10670 中allB基因含有一個無義突變（圖3?F）。E24190 的rpoS基因發生771 bp 的插入變異。這些變異將導致相關基因編碼的蛋白失活。

圖3 毒力因子基因的高影響變異位點Fig. 3 High-impact variant loci in virulence factor genes

表7 各菌株的基因組變異位點Table 7 Variation sites in the genome of each strain

表8 七個具有高影響變異的毒力因子基因Table 8 Seven virulence factor genes with high impact variants

3 討論

豬源致病性大腸桿菌是一種常見的傳染性病原體，可引起嚴重腹瀉等疾病［23］。為了有效地預防和控制豬大腸桿菌病，本研究通過基因組測序、組裝，得到6 株不同毒力的豬源致病性大腸桿菌的基因組信息。相比之下，強毒菌株S10670 和E24190 的基因組較大（5.3 Mb 和4.96 Mb），基因組最小的是弱毒菌株P111，僅有4.62 Mb，與參考基因組的大小相似。參考基因組E. coli. str.K?12 屬于非致病的腸道共生菌，全長4.64 Mb，有4 288 個編碼基因，與致病性大腸桿菌（O157:H7，菌株EDL933）相比少了1 387 個基因［24］。菌株S10670、E24190、P555、P211、P32、P111 的CDS 數目依次為5 040、4 732、4 531、4 483、4 837、4 308，編碼基因從3 519 個降至3 364 個。總體上，隨著菌株的毒力從高到低變化，基因組長度、CDS 和編碼基因數量表現出逐漸降低的趨勢。

大腸桿菌的致病力與種系發育分群有一定關系，B2 和D 群是主要的腸道外致病性大腸桿菌，A和B1 群多為腸道共生型大腸桿菌，四川地區豬源致病性大腸桿菌多屬于A 和B1 群［25］，新疆地區的多為B1 和A 群［26］，本研究檢測的菌株中，本地分離的4 株菌均屬于A 群，提示養豬業應關注腸道共生型大腸桿菌對大腸桿菌病的影響。大腸桿菌的某些特殊血清型具有致病性［27］。本研究中，從患病豬糞樣中分離的菌株血清型為H4、O4: H45、O9: H4、O26: H30，從云南分離得到101 株豬源大腸桿菌中優勢血清型為O4［28］，國內山東地區檢測出血清型O9［29］，本研究檢測菌株的血清型有一定的多樣性。此外，本研究檢測到3 個菌株屬于ST10 克隆復合群，新發現一種ST 型。山東360 株豬源大腸桿菌中，優勢的復合群為ST10 克隆復合群［30］。表明ST10 克隆復合群可能是豬源大腸桿菌的優勢ST 型。

大腸桿菌攜帶的毒力因子是決定其致病性的關鍵因素［31］。在感染初期，黏附素可介導細菌定植于宿主的腸道上皮細胞［32］，不同致病性大腸桿菌的菌毛類型不同。本研究的強毒菌株與弱毒菌株相比擁有更多的黏附系統基因，強毒菌株擁有毒素相關的菌毛黏附素和一些毒性較強的非菌毛黏附素，以協助大腸桿菌入侵宿主細胞［33］，如CFA/I 菌毛基因僅在P211 基因組中存在。

細菌的鐵吸收和轉運系統能幫助細菌從宿主獲取鐵元素，本研究檢測的6 株大腸桿菌中，均含有鐵吸收和轉運系統中的腸桿菌素基因，鐵載體內膜轉運系統僅在P211 基因組中存在，血紅素吸收系統僅從S10670 中檢測到。構建血紅素吸收系統基因chuT基因缺失突變株，證明細菌的致病性不受chuT基因的影響［34］。推測鐵吸收和轉運系統可能不是大腸桿菌毒力的決定性因素。

致病性大腸桿菌的外毒素主要包括志賀毒素（Stx）和腸毒素等。本實驗室前期從貴州本地豬場的患病仔豬中分離到志賀樣毒素和腸毒素基因［35-36］。本研究中僅在S10670 基因組中檢測到志賀樣毒素基因，E24190 中含有耐熱腸毒素、不耐熱腸毒素基因以及腸聚集黏附性耐熱毒素east1基因。此外，細菌釋放的溶血素（hly）可進入宿主的血液循環，在靶細胞膜上形成孔道，使大部分哺乳動物的紅細胞和免疫細胞破裂［37］，本研究僅從菌株S10670 基因組中檢測到溶血素基因hlyA，中等毒菌株和弱毒菌株中都沒有檢測到溶血素基因。大腸桿菌通過志賀樣毒素、溶血素等改變宿主腸道的完整性，進而引發炎癥等，最終導致豬的腹瀉［38］。本研究結果提示強毒力菌株可能通過多種毒素的綜合作用，破壞宿主的腸道屏障致病。弱毒菌株含有的毒素基因較少，致病力較弱。

細菌利用不同的分泌系統可以將蛋白運送到細菌外。迄今從革蘭氏陰性菌中發現了6 種分泌系統（即I-VI 型）［39］。T1SS 主要運輸抗生素和毒素等小分子［40］，與ABC 轉運蛋白非常相似［41］，溶血素和自分泌黏附素均屬于T1SS，自分泌黏附素是由aatPABCD基因簇編碼的外膜蛋白，能促進致病性大腸桿菌對宿主細胞的黏附，破壞宿主的自我保護屏障［42］。從107 株鴨源致病性大腸桿菌中，檢測出28.94%的菌株含有aatA基因［43］。本研究中這兩種T1SS 基因只存在于強毒菌株中。T2SS 不僅具有豐富的蛋白酶分泌活性，還能介導多種毒素向細菌外傳遞［44］。本研究發現強毒菌株中的T2SS 基因數量多，弱毒菌株的T2SS 基因較少。T3SS 可以將細菌的物質注入宿主細胞內，干擾宿主的免疫功能。大腸桿菌有兩種T3SS，分為ETT1 和ETT2，ETT1 由LEE毒力島編碼，擁有LEE 毒力島的致病性大腸桿菌可損傷動物的腸黏膜細胞。ETT2 在腸致病性大腸桿菌的侵襲、胞內存活和毒力等方面發揮重要作用［45］，如ETT2 結構基因epaPQR簇參與致病性大腸桿菌的鞭毛形成并調控相關基因的轉錄，進而影響細菌的定殖能力［46］。本研究從強毒菌株S10670 基因組中檢測到幾乎全部的T3SS 基因，中等毒菌株P211 和弱毒菌株P32 中的ETT2 基因減少，其余菌株中未發現T3SS 基因。大腸桿菌的T5SS 是一個自主轉運系統，主要由膜外功能區和跨膜區、信號肽及自分泌重復序列3 部分組成［47］，主要與黏附和生物被膜的形成有關。本研究中所有菌株都檢測到T5SS 的ehaB基因，但基因ehaA僅在S10670 基因組中分布。T6SS 與細菌的致病性密切相關，25%的革蘭氏陰性菌擁有T6SS 基因，其中大多數為致病菌［48］。本研究中只有菌株S10670 和P555 檢測到豐富的T6SS 基因20 個和19 個，菌株E24190、P32 中的T6SS 基因較少，僅有6 個和3 個。

4 結論

豬源致病性大腸桿菌基因組長度呈多態性，強毒菌株的基因組長度、CDS 和編碼基因的數量高于弱毒菌株，并且擁有更多的毒力因子基因，尤其毒素基因；且強毒菌株的毒力因子基因具有較多的高影響變異。結果表明，豬源致病性大腸桿菌毒力的強弱與菌株之間毒力因子基因的數量和基因結構相關。