長雙歧桿菌W13 的全基因組序列分析

許婷婷，司靜，閆夢娜，逯彬，劉翔，王艷萍，耿偉濤

（天津科技大學(xué) 食品科學(xué)與工程學(xué)院，天津 300457）

哺乳動物的胃腸道內(nèi)定植著復(fù)雜的微生物群落，這些微生物在調(diào)控宿主營養(yǎng)吸收、腸道功能及免疫應(yīng)答等方面具有重要作用[1]。其中，益生菌作為胃腸道中的有益細菌，有助于改善腸道微生態(tài)平衡、促進人體健康[2-3]。據(jù)統(tǒng)計，目前人體中所含益生菌種類超過400 種，其中主要有4 類：乳酸菌、雙歧桿菌、酵母菌和芽孢桿菌[4]。

雙歧桿菌是較早發(fā)現(xiàn)的一種重要的益生菌。目前鑒定的絕大多數(shù)雙歧桿菌來源于人類和其他哺乳動物的胃腸道，其在胃腸道中的數(shù)量取決于個體的年齡和飲食方式[5]。Arboleya 等[6]發(fā)現(xiàn)腹瀉、絞痛、過敏、壞死性小腸結(jié)腸炎、肥胖等疾病患者的胃腸道中雙歧桿菌豐度顯著降低，這表明雙歧桿菌對腸道菌群穩(wěn)態(tài)具有重要貢獻[7]。

隨著高通量測序技術(shù)的發(fā)展，對雙歧桿菌的研究進入了基因組學(xué)時代[8]。了解菌株的基因組基本信息，可以更進一步認識基因與蛋白質(zhì)、代謝功能、個體行為之間的關(guān)系。比較基因組是在高通量測序的基礎(chǔ)上發(fā)展而來，通過比較不同物種、菌株或同一物種的不同菌株的基因組的差異和相似性，能夠了解基因功能之間的關(guān)系[9]。

本課題組在前期試驗中從嬰兒糞便中篩選分離得到了一株具有抗氧化特性的長雙歧桿菌。長雙歧桿菌作為人體腸道內(nèi)存在的分布最廣、豐度最高的一種雙歧桿菌，且具有一定的益生特性，近年來引起了廣泛的關(guān)注。本文通過生物信息學(xué)技術(shù)研究長雙歧桿菌W13 的基因組特征，挖掘菌株的功能特性基因，以期為長雙歧桿菌長亞種資源的開發(fā)與應(yīng)用提供參考。

1 材料與方法

1.1 材料

長雙歧桿菌W13：天津科技大學(xué)發(fā)酵食品與益生菌資源開發(fā)實驗室前期自嬰兒糞便中分離得到。菌株W13 的基因組序列已提交至美國國立生物技術(shù)信息中心（national center for biotechnology information，NCBI）的GenBank 數(shù)據(jù)庫，登錄號為CP096771.1。

1.2 方法

1.2.1 菌株W13 的基因組特性

自NCBI 下載菌株W13 的全基因組序列，利用CGVIEW（https://proksee.ca/projects/new）預(yù)測基因組的基因結(jié)構(gòu)，并進行CRISPR 序列的預(yù)測。

1.2.2 生物信息學(xué)分析

使用本地Blast 工具將預(yù)測得到的編碼基因蛋白序列與直系同源蛋白分組比對（evolutionary genealogy of genes：non-supervised orthologous groups，eggNOG）數(shù)據(jù)庫、基因組百科全書（kyoto encyclopedia of genes and genomes，KEGG）進行比對，獲得編碼基因的功能信息；通過碳水化合物活性酶（carbohydrate-active enzymes，CAZy）數(shù)據(jù)庫的比對分析基因組中的碳水化合物活性酶。通過使用子系統(tǒng)技術(shù)快速注釋（rapid annotation using subsystem technology，RAST）數(shù)據(jù)庫、抗性基因數(shù)據(jù)庫（comprehensive antibiotic resistance database，CARD）完成菌株W13 基因組中耐藥基因的預(yù)測。使用IslandViewer4 （https://www.pathogenomics.sfu.ca/islandviewer/upload/）對菌株W13 中存在的基因組島進行預(yù)測[10]。利用PHASTER（https://phaster.ca/）在線網(wǎng)站分析鑒定和注釋微生物基因組內(nèi)的噬菌體序列[11-13]。

對于系統(tǒng)發(fā)育樹的構(gòu)建，先使用軟件BPGA 預(yù)測核心基因，再使用軟件MEGA 11.0.13 構(gòu)建系統(tǒng)發(fā)育進化樹。根據(jù)核心基因建樹結(jié)果，選取與菌株W13 親緣關(guān)系較近的7 株菌，利用EZBioCloud 網(wǎng)站（https://www.ezbiocloud.net/tools/ani） 比較菌株W13 與7 株菌之間的平均核苷酸一致性（average nucleotide identity，ANI）[14]，并繪制熱圖。長雙歧桿菌長亞種的泛基因組和核心基因組預(yù)測使用PanGP 軟件完成[15]，通過BRIG 軟件比較分析基因組。

2 結(jié)果與分析

2.1 菌株W13 基因組特性

2.1.1 菌株W13 基因組基本特征

對菌株W13 的基因組序列進行分析，其基本特征見表1。

表1 W13 基因組基本特征Table 1 Basic characteristics of Bifidobacterium longum W13 genome

由表1 可知，長雙歧桿菌W13 的基因組序列全長為2 335 916 bp。此外，長雙歧桿菌W13 的GC 含量（GC 含量指DNA 中鳥嘌呤和胞嘧啶所占比例的比率）為60.95%。對菌株W13 的基因組基本特征做可視化處理，繪制菌株W13 的染色體基因組圈圖，見圖1。

圖1 長雙歧桿菌W13 的基因組圈圖Fig.1 Genomic cycle of Bifidobacterium longum W13

圖1 展示了預(yù)測的GC 含量、編碼序列區(qū)（coding sequence，CDS）功能、轉(zhuǎn)運RNA（transfer RNA，tRNA）、核糖體RNA（ribosomal RNA，rRNA）、信使RNA（transfermessenger RNA，tmRNA）、成簇規(guī)律間隔的短回文重復(fù)序列（clustered regularly interspaced short palindromic repeats，CRISPR）。該菌株的基因組由一個不包含質(zhì)粒的完整環(huán)狀染色體構(gòu)成，依據(jù)CGVIEW 預(yù)測菌株W13 基因組中的編碼基因，共包含1 961 個編碼基因，占整個基因組的86.76%。CRISPR 序列的具體排布方式見圖2。

圖2 菌株W13 中的CRISPR 基因簇Fig.2 CRISPR gene cluster in strain W13

由圖2 可知，在菌株W13 中，依據(jù)其Cas 效應(yīng)蛋白可將菌株含有的內(nèi)源CRISPR 系統(tǒng)分類為Ⅰ-C 型。CRISPR 序列廣泛分布在原核生物中，是大多數(shù)細菌及古細菌中一種不斷進化適應(yīng)的免疫防御機制，其分布具有菌株特異性，可用于進行菌株的分型。

2.1.2 菌株W13 的基因組島預(yù)測

基因組島是微生物基因組中可水平轉(zhuǎn)移的基因簇，基因組島的轉(zhuǎn)移是微生物基因交換的一種方式，可以提高微生物的多樣性和對環(huán)境的適應(yīng)性。

使用IslandViewer4 基因組島在線預(yù)測網(wǎng)站得到的菌株W13 基因組島預(yù)測結(jié)果見圖3。

圖3 菌株W13 的基因組島預(yù)測結(jié)果Fig.3 Predicted genomic islands in strain W13

從圖3 可以看出，菌株W13 的染色體上共含有33 個基因組島。在這33 個基因組島中，島4 是最長的基因組島，長度為25 714 bp；島13 是最短的基因組島，長度為4 321 bp。

其中，島23 屬于碳水化合物活性酶GH2 家族；島16 屬于碳水化合物活性酶GH25 家族；島2、島26 均屬于碳水化合物活性酶GH43 家族；島4、5、6 均屬于碳水化合物活性酶GT2 家族。碳水化合物活性酶相關(guān)基因占總基因組島的21%，這表明菌株W13 具有較高的糖合成能力。

此外，島7 與CRISPR 序列有關(guān)，包含Cas 效應(yīng)蛋白，這與2.1.1 的分析對應(yīng)。

2.2 菌株W13 的基因功能注釋

2.2.1 eggNOG 數(shù)據(jù)庫注釋

eggNOG 數(shù)據(jù)庫收集了全面的物種和大量的蛋白序列數(shù)據(jù)，并且能進行同源基因分類以及功能注釋，是目前最先進、最完善的數(shù)據(jù)庫。

通過分析，通過eggNOG 數(shù)據(jù)庫注釋，預(yù)測菌株W13 的蛋白質(zhì)功能見圖4。

圖4 eggNOG 數(shù)據(jù)庫注釋分類統(tǒng)計圖Fig.4 eggNOG database annotation classification statistics

從圖4 可以看出，在菌株W13 的基因組中，除占比最多的假定蛋白外，關(guān)于碳水化合物的運輸和代謝（G）、氨基酸的運輸和代謝（E）的相關(guān)基因最多，分別占比10.30%、9.95%，表明菌株可能具有較高的糖代謝和氨基酸代謝能力；其次為K（轉(zhuǎn)錄）、L（復(fù)制、重組和修復(fù)）、J（翻譯、核糖體結(jié)構(gòu)和生物發(fā)生），分別占比9.36%、9.07%、8.37%。

2.2.2 KEGG 數(shù)據(jù)庫注釋

KEGG 數(shù)據(jù)庫是收集了生物的基因組、通路和化合物信息的綜合性的數(shù)據(jù)庫，該數(shù)據(jù)庫主要將基因組信息分為三大類，分別為環(huán)境信息處理（environmental information processing）、新陳代謝（metabolism）、遺傳信息處理（genetic information processing），通過KEGG 數(shù)據(jù)庫注釋，菌株W13 中關(guān)于環(huán)境信息處理、新陳代謝、遺傳信息處理中分別涉及基因數(shù)目為151、630、90 個。將這三大類中涉及較多的基因統(tǒng)計，得到的注釋信息見圖5。

圖5 KEGG 數(shù)據(jù)庫注釋分類統(tǒng)計圖Fig.5 KEGG database annotation classification statistics

從圖5 可以看出，環(huán)境信息處理類占比最多的為腺苷三磷酸結(jié)合盒轉(zhuǎn)運蛋白（ATP-binding cassette transporter proteins）。該家族蛋白是已知最大的蛋白質(zhì)家族之一，廣泛存在于細菌、古菌和真核生物中。特別地，群體感應(yīng)（quorum sensing，QS）相關(guān)基因占比較多，QS 是細菌根據(jù)細胞密度變化進行基因表達調(diào)控的一種生理行為，可通過微生物的信息交流調(diào)控腸道屏障功能與營養(yǎng)素代謝，進而維持機體對營養(yǎng)素的吸收和腸道健康，促進腸道穩(wěn)態(tài)。新陳代謝類數(shù)量最多的為嘌呤代謝途徑相關(guān)基因，很多細菌在人體嘌呤代謝中起著重要作用，人體內(nèi)三分之二的尿酸是由腎臟排出的，其余的主要由腸道排出，尿酸酶由尿酸轉(zhuǎn)化為尿囊素和尿素，廣泛存在于乳酸菌中。因此菌株W13 可能具有延緩腎臟疾病進程的潛力。遺傳信息處理類數(shù)量最多的為核糖體相關(guān)基因，主要負責(zé)基因的翻譯。

2.2.3 RAST 數(shù)據(jù)庫注釋

RAST 數(shù)據(jù)庫可用于注釋細菌和古菌基因組，該服務(wù)器通過識別蛋白質(zhì)編碼基因、rRNA 和tRNA 基因，來預(yù)測不同基因的功能，使用這些信息重建代謝網(wǎng)絡(luò)[16]。

在菌株W13 的基因組中，RAST 數(shù)據(jù)庫注釋到了28%的基因，預(yù)測結(jié)果見圖6。

圖6 RAST 數(shù)據(jù)庫注釋Fig.6 RAST database annotation

從圖6 可以看出，在菌株W13 的基因組功能預(yù)測結(jié)果中，占比最高的三大類分別為蛋白質(zhì)代謝（protein metabolism）、氨基酸及其衍生物（amino acids and derivatives）、碳水化合物（carbohydrates），分別占比21%、19.92%、16.09%。特別地，在脅迫應(yīng)答調(diào)控基因類下注釋到了谷胱甘肽氧化還原循環(huán)相關(guān)酶，谷胱甘肽是一種三肽，在微生物的大多數(shù)活細胞中以高濃度存在，谷胱甘肽在眾多細胞功能中起著關(guān)鍵作用，包括清除自由基、參與氧化還原反應(yīng)、參與脫氧核糖核苷酸的形成等，這表明菌株W13 可能具有抗氧化潛力，印證了本課題組前期對菌株W13 抗氧化能力的研究[17]。同時，還預(yù)測到了膽鹽水解酶相關(guān)基因，這表明菌株W13 可能具有耐膽鹽的能力。

2.2.4 CAZy 數(shù)據(jù)庫注釋

CAZy 數(shù)據(jù)庫是收錄碳水化合物活性酶的數(shù)據(jù)庫。該數(shù)據(jù)庫將蛋白質(zhì)按功能類別分為六大類：糖苷水解酶（glycoside hydrolases，GH）、糖基轉(zhuǎn)移酶（glycosyl transferases，GT）、多糖裂解酶（polysaccharide lyases，PL）、碳水化合物酯酶（carbohydrate esterases，CE）、輔助氧化還原酶（auxiliary activities，AA）和非催化的結(jié)合碳水化合物的功能域（carbohydrate-binding modules，CBM）。菌株W13 的注釋結(jié)果見圖7。

圖7 碳水化合物酶分布比例圖Fig.7 Distribution ratio of carbohydrate enzymes

在菌株W13 中共注釋到80 個GH 家族，其中占比較多的是GH2、GH3、GH13 和GH43，數(shù)量分別為5、6、13 和8。GH2 分布在哺乳動物組織、植物和微生物中，能夠產(chǎn)生具有轉(zhuǎn)糖基活力的β-半乳糖苷酶，這是以母乳或奶粉為主要營養(yǎng)來源的嬰幼兒體內(nèi)最重要的消化酶，GH2 家族的β-半乳糖苷酶包括LacZ 型和異二聚體LacLM 型酶，雙歧桿菌菌株多含有LacZ 型酶[18]。在食品工業(yè)中，β-半乳糖苷酶因其水解乳糖和合成低聚半乳糖的能力而被利用；GH3 家族包含β-葡萄糖苷酶，該酶在微生物中尤為普遍，可使乳糖水解為葡萄糖和半乳糖。GH3 和GH43 家族含具有β-木糖苷酶活性的糖苷水解酶，該酶是重要的木糖分解酶，用于許多生物技術(shù)過程，可應(yīng)用于海洋產(chǎn)品加工；GH13 為淀粉水解酶家族，是雙歧桿菌屬糖生物組中最豐富的GH 家族之一[19-21]。

此外，在菌株W13 中共注釋到32 個GT 家族，其中占比較多的糖基轉(zhuǎn)移酶家族是GT2 和GT4，數(shù)量分別為17、9。多種合成型糖基轉(zhuǎn)移酶屬于GT2 家族，例如纖維素、幾丁質(zhì)、透明質(zhì)酸或葡聚糖合酶。GT4 具有蔗糖合成酶、蔗糖磷酸合成酶等活性[22]，綜上說明W13有較高的糖合成能力，具備合成胞外多糖的潛力。

2.2.5 抗生素抗性基因預(yù)測

將菌株W13 基因的蛋白序列上傳到CARD，使用工具軟件RGI 進行比對并預(yù)測潛在的抗性基因。預(yù)測結(jié)果見表2。

表2 CARD 預(yù)測抗性基因Table 2 Resistance genes predicted by CARD

根據(jù)表2 的比對發(fā)現(xiàn)，菌株W13 中含有利福平（rifampicin，RFP）抗性基因，與青春雙歧桿菌rpoB 突變體的序列相似性達到了92.65%。

2.2.6 前噬菌體預(yù)測

一部分微生物的橫向轉(zhuǎn)移基因是通過噬菌體的侵染而獲得的。在PHASTER 中共預(yù)測出2 條前噬菌體片段，見表3。

表3 前噬菌體預(yù)測基本情況Table 3 Basic information of prophage prediction

由表3 可知，菌株W13 中所包含的2 條前噬菌體長度均小于30 kbp（1 kbp=1 000 bp），且不完整，因此認為這2 條噬菌體不具有活性。將預(yù)測到的前噬菌體片段對比基因組島分析，發(fā)現(xiàn)基因組島中包含部分與前噬菌體相關(guān)的基因片段，但不具有噬菌體活性。因此認為菌株W13 中并無可進行水平轉(zhuǎn)移的致病島。

2.3 比較基因組學(xué)分析

2.3.1 基于核心基因的系統(tǒng)發(fā)育進化樹分析

利用BPGA 軟件將菌株W13 的蛋白序列與8 株長雙歧桿菌豬亞種、8 株長雙歧桿菌嬰兒亞種及8 株長雙歧桿菌長亞種的蛋白序列分析比對，結(jié)果見圖8。

圖8 基于核心基因建立的菌株W13 系統(tǒng)發(fā)育樹Fig.8 Phylogenetic tree of strain W13 based on core gene

圖8 的比對結(jié)果顯示，菌株W13 的核心基因數(shù)為990，必需基因數(shù)為857。利用核心基因建立系統(tǒng)發(fā)育樹，由圖8 可知，菌株W13 與模式菌株Bifidobacterium longum subsp. longum JCM 1217 進化遺傳距離最近，推測菌株W13 為長雙歧桿菌長亞種。Bifidobacterium longum subsp.infants CECE 7210 與長亞種的進化遺傳距離較近，這一結(jié)果與Blanco 等[24]的研究結(jié)果一致，是由于CECE 7210 被錯誤地分類到了嬰兒亞種中。

2.3.2 平均核苷酸一致性分析

基于核心基因組，選定了7 株與菌株W13 親緣關(guān)系較近的菌株，在NCBI 數(shù)據(jù)庫中下載的這7 株長雙歧桿菌的基因組數(shù)據(jù)見表4。

表4 ANI 分析選用的7 株菌基因組信息Table 4 Genome information of 7 strains selected for ANI analysis

將菌株W13 與表4 中的7 株菌進行了ANI 值的計算，對比結(jié)果見圖9。

圖9 菌株W13 與親緣關(guān)系較近的8 株菌的ANI 值Fig.9 ANI value of strain W13 and 8 strains with close relationship

由圖9 可知，菌株W13 與5 株長雙歧桿菌長亞種的ANI 值高于98%，而與2 株長雙歧桿菌嬰兒亞種的ANI 值均在95%以下，與系統(tǒng)發(fā)育樹的結(jié)果具有一致性。這說明菌株W13 屬于長雙歧桿菌長亞種。

2.3.3 泛基因組、核心基因組分析

使用PanGP 軟件預(yù)測了長雙歧桿菌長亞種的泛基因組和核心基因組隨基因組數(shù)目的變化，結(jié)果見圖10。

圖10 Bifidobacterium longum subsp. longum 的泛基因組特征曲線Fig.10 Pan-genome characteristic curves of Bifidobacterium longum subsp. longum

由圖10 可以看出，隨著菌株數(shù)量的增加，泛基因組的擬合曲線呈上升趨勢，這說明長雙歧桿菌長亞種的泛基因組是開放性的，這可能是由于長雙歧桿菌長亞種在自然界中分布廣泛，與外界各種遺傳物質(zhì)發(fā)生交換，說明該物種有較高的遺傳多樣性。相反，其核心基因的數(shù)量隨菌株數(shù)量增加逐漸趨于穩(wěn)定，說明該物種的穩(wěn)定性較高。

2.3.4 長雙歧桿菌長亞種的基因組比較

利用BRIG 軟件對菌株W13 與8 株長雙歧桿菌長亞種的全基因組進行了比較，這8 株菌的基因組信息見表5。

表5 BRIG 軟件比較基因組使用的菌株信息Table 5 Strain genome information compared by BRIG software

以菌株W13 作為參考基因組，菌株W13 與從NCBI 數(shù)據(jù)庫下載的8 株長雙歧桿菌長亞種的基因組進行比較，結(jié)果見圖11。

圖11 比較基因組圖Fig.11 Comparative genome map

由圖11 可以看出，菌株W13 的染色體與其它菌株相比有其特定的功能基因。圖中空白區(qū)域代表在參考基因組中存在而在其它基因組中不存在的基因，即區(qū)域1～7 是菌株W13 所特有的。經(jīng)過分析，區(qū)域2、4、7 為假定蛋白，區(qū)域1、5 與復(fù)制、重組和修復(fù)區(qū)域；細胞壁/膜/包膜生物發(fā)生有關(guān)；區(qū)域3 為CRSPR 序列區(qū)域；區(qū)域6 與輔酶轉(zhuǎn)運與代謝相關(guān)，表明W13 可能具有更強的抗氧化、抗衰老、抗炎等能力。

3 結(jié)論

本研究從基因組特性出發(fā)，對長雙歧桿菌W13 進行了系統(tǒng)化的分析。通過對菌株W13 進行全基因組測序，在KEGG、eggNOG、RAST 數(shù)據(jù)庫中系統(tǒng)分析了其功能基因與代謝通路，發(fā)現(xiàn)該菌株可能具有抗氧化、耐膽鹽的能力，W13 含有較多的群體感應(yīng)基因，可促進腸道穩(wěn)態(tài)。除此之外，W13 含有內(nèi)源Ⅰ-C 型CRISPR 系統(tǒng)，可用于進行菌株的分型，未來可應(yīng)用于進行內(nèi)源CRISPR 系統(tǒng)基因編輯。通過分析碳水化合物活性酶基因，預(yù)測該菌株具有合成胞外多糖的潛力，能夠高產(chǎn)胞外多糖。通過基因組學(xué)分析，認為長雙歧桿菌長亞種在自然界中分布廣泛，有較高的遺傳多樣性和穩(wěn)定性。綜上，認為長雙歧桿菌W13 是一株益生功能強、具有商業(yè)化價值的菌株，為長雙歧桿菌的開發(fā)與應(yīng)用提供了參考。