劉廣鑫，董晏君，趙麗娟，鄧益琴，程長洪，馬紅玲，江建軍，馮 娟，郭志勛，林 蠡

1. 中國水產科學研究院南海水產研究所/農業(yè)農村部南海漁業(yè)資源開發(fā)利用重點實驗室，廣東 廣州 510300

2. 肇慶市水產技術推廣中心，廣東 肇慶 526040

3. 仲愷農業(yè)工程學院 動物科技學院/廣東省水環(huán)境與水產品安全工程技術研究中心/廣州市水產病害與水禽養(yǎng)殖重點實驗室，廣東 廣州 510225

高密度的養(yǎng)殖模式中，池塘中殘余的飼料蛋白和水生動物的排泄物等均含有大量的含氮化合物。這些累積的含氮化合物在養(yǎng)殖水體中被慢慢分解，并產生氨基酸，氨基酸再被微生物脫氨產生氨氮，最終養(yǎng)殖水體的氮源會以分子氨 (NH3-N)、離子銨、亞硝態(tài)氮、硝態(tài)氮、氮氣 (N2) 及有機氮化物 (R-NH2) 的形式以動態(tài)平衡的方式存在[1-3]。而NH3和之間可相互轉化：，當水體 pH值和溫度升高時，水體中NH3的比例增加。NH3和對經濟動物組織、器官有直接毒害作用，是造成機體免疫力下降并導致疾病爆發(fā)的關鍵因素之一[4-5]。通過綠色、高效的生物方式降低水體中無機氮的濃度達到水質凈化的目的，是提高水生動物成活率的主要方法之一。

芽孢桿菌 (Bacillusspp.) 作為異養(yǎng)微生物，是一種廣泛存在于自然界中，并能夠高效利用無機氮實現(xiàn)自身生長、增殖的革蘭氏陽性菌，也是水產養(yǎng)殖領域中研究最多、應用最廣泛的菌種之一[6]。隨著研究的深入，芽孢桿菌同化無機氮的代謝途徑逐漸被挖掘，同化代謝過程中涉及的多種特異性功能蛋白的作用被發(fā)現(xiàn)，例如特異性感應蛋白、轉運蛋白、氧化還原蛋白、同化蛋白等，但其作用機理尚不清晰[7]。

全基因組測序 (Complete genome sequence) 是指將生物的基因組進行隨機打斷、測序和拼接，從而獲得完整基因組序列的技術，通過序列分析和功能基因注釋 (包括NR、KEGG、eggNOG、GO和CARD等)，以便從基因水平深入了解菌株表型的機理[8-9]。目前，隨著生物技術和信息學的有機結合，全基因組測序技術已經發(fā)展成為一種高效、準確的信息獲得手段，在生物研究領域中應用廣泛，例如腸道菌群分析、微生態(tài)群系組成、特異性表型分析、功能基因的預測及核心基因的挖掘等[10-13]。

利用芽孢桿菌降低水體中無機氮的濃度具有綠色、持續(xù)、無二次污染等特點，被學者和養(yǎng)殖戶普遍認可[14]。盡管芽孢桿菌同化無機氮的效果得到了很好驗證，但其代謝途徑和特異性功能蛋白等仍有待深入解析和驗證。貝萊斯芽孢桿菌 (B. velezensis) LG37篩選于華中農業(yè)大學教學實習基地養(yǎng)殖池塘中，前期研究發(fā)現(xiàn)LG37可通過無機氮作為唯一氮源的最小培養(yǎng)基進行培養(yǎng)，并具有高效同化無機氮的能力[7]?；诖?，本研究結合三代PacBio RS II和二代 Illumina HiSeq 2000測序技術，對 LG37進行全基因組測序分析，并對獲得的數(shù)據進行NR、KEGG、eggNOG、GO、CARD數(shù)據庫注釋、分析，旨在深入揭示LG37在無機氮代謝過程中的關鍵功能基因和代謝通路信息，為其在無機氮代謝研究和水產養(yǎng)殖中的應用提供參考信息和科學依據。

1 材料與方法

1.1 實驗材料

1.1.1 菌株

貝萊斯芽孢桿菌LG37為實驗室于華中農業(yè)大學教學實習基地水產養(yǎng)殖池塘中篩選獲得。

1.1.2 相關試劑和儀器設備

LB (Luria-Bertani) 培養(yǎng)基；細菌基因組 DNA提取試劑盒，天根生化科技 (北京) 有限公司。Pac-Bio RS II (Pacific Biosciences, Menlo Park, CA, USA)。

1.2 實驗方法

1.2.1 LG37 基因組 DNA 提取

利用LB固體培養(yǎng)基對LG37進行活化培養(yǎng)，挑取長勢良好的單菌落接種于LB 液體培養(yǎng)基，培養(yǎng)至對數(shù)期時，通過離心的方式對培養(yǎng)菌液進行集菌及DNA提取，DNA 提取的具體操作流程根據天根生化試劑盒說明書進行。利用NanoDrop 2000/2000c (Thermo Scientific, Wilmington, DE, USA) 對提取、純化的LG37基因組質量進行分析。其DNA(OD260/280) 值在 1.8～2.0，質量濃度≥20 ng·μL?1可用于后續(xù)的建庫測序。

1.2.2 系統(tǒng)進化樹分析

將擴增獲得的 LG37 16S rDNA 序列與 NCBI中的芽孢桿菌16S rDNA序列運用DNAMAN軟件進行比對分析。應用MEGA 6.0軟件對LG37和下載的16S rDNA進行遺傳進化分析。其中序列比對采用Clustal W 多重比較法進行，基于鄰接法 (Neighbor-joining) 構建系統(tǒng)進化樹，確定菌株的進化分類[15]。

1.2.3 文庫構建和庫檢

采用三代 PacBio RS II結合二代 Illumina HiSeq 2000測序的方法對LG37進行全基因組測序，其完成單位為上海歐易生物醫(yī)學科技有限公司[16]。經質檢合格的 DNA 樣品采用 PacBio RS II single-molecule real-time (SMRT) 測序技術 (Pacific Biosciences, Menlo Park, CA, USA) 對其進行全基因組測序。然后通過分級基因組裝配(The hierarchical genome-assembly process, HGAP) 對隨機打斷的測序讀段進行裝配，以期獲得超長的測序讀段。接下來使用二代測序技術Illumina HiSeq 2000進行測序[17]，并利用 Bowtie2 (v2.3.0) 將 Illumina 基因片段與組裝的SMRT-contigs進行比對[18]，從而獲得完整的環(huán)狀全基因組序列信息。用Glimmer 3.02進行基因預測和注釋[19]，用tRNAscan-SE進行tRNA genes的預測和篩選[20]，用RNAmmer進行rRNA 的預測分析[21]。

1.2.4 本地 Blast+功能基因比對

通過NCBI下載本地Blast+，即“基本局部對比搜索工具”(Basic Local Alignment Search Tool, https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/)。首先對于基因組測序獲得的蛋白注釋信息進行模型建庫，然后將NCBI (https://www.ncbi.nlm.nih.gov/)和 UniProt (https://www.uniprot.org/) 數(shù)據庫相關信息進行下載、比對、篩選和分析，獲得無機氮代謝相關候選基因，并注釋其GO分子功能[22]。

2 結果

2.1 LG37 基因組組裝

通過三代 PacBio RS II結合二代 Illumina HiSeq 2000測序平臺對LG37 進行全基因組測序分析(表1)。97 382 個平均長度 3 016 nt的測序讀段被最終拼裝成一條重疊基因組序列，長度為3 949 250 bp。通過測序比對分析，最終獲得一條 3 929 697 bp的環(huán)狀染色體基因組，GC含量為46.5%，總基因數(shù)3 967 個 (圖1)。其中蛋白編碼基因 3 854 個，編碼區(qū)域長度為 3 495 864 bp，占總長度的 89.0%，GC含量為47.3%。RNA 基因113個，其中tRNA基因86個，rRNA基因27個?；蜷g隔區(qū)域長度433 833 bp，占總長度的11.0%。該LG37菌株無載體。

表1 貝萊斯芽孢桿菌LG37 基因組特性Table 1 Genome features of B. velezensis LG37 strain

圖1 LG37 環(huán)狀基因組圖譜Fig. 1 Circular genome map of B. velezensis LG37 strain

2.2 LG37 系統(tǒng)進化樹分析

經全基因組測序結果顯示，其16S rDNA基因序列長度為 1 445 bp。通過 NCBI→Blast→Nucleotide Blast在線分析軟件比對數(shù)據庫中已收錄的16S ribosomal RNA sequences，并結合基因組測序數(shù)據匹配分析。結果表明，篩選的芽孢桿菌菌株為水源貝萊斯芽孢桿菌 (B. velezensis)，菌株命名為LG37。將LG37菌株與NCBI數(shù)據庫中進化關系較近的芽孢桿菌的16S rDNA序列進行比對，并通過MEGA 6.0軟件中的鄰接法構建系統(tǒng)發(fā)育樹 (圖2)。發(fā)現(xiàn) LG37 與B. velezensisFZB42 親緣關系最近，16S rDNA序列相識度高達100%。

圖2 貝萊斯芽孢桿菌基于16S rDNA 序列構建的系統(tǒng)進化樹Fig. 2 Neighbor-joining tree of B. velezensis LG37 based on 16S rDNA sequences

2.3 LG37 基因功能注釋

對全基因組測序獲得的蛋白編碼基因，通過搜尋 NCBI的 nr庫、KEGG (Kyoto Encyclopedia of Genes and Genome) 蛋白數(shù)據庫以及 SEED 蛋白數(shù)據庫進行基因功能注釋分析，利用CDD數(shù)據庫進行 COG (Clusters of Orthologous Groups) 分類，并通過KEGG數(shù)據庫構建代謝通路。全基因組測序共獲得3 854個蛋白編碼基因，具有明確生物學功能的蛋白編碼基因3 057個，具有KEGG的ortho-log蛋白編碼基因2 108個，具有COG分類的蛋白編碼基因2 890個。所有蛋白匹配的同源蛋白來源于 75 個物種，其中B. velezensisFZB42 比例最高(78.4%)。

2.3.1 COG 功能分類

LG37完整基因組的蛋白編碼基因的氨基酸序列數(shù)據與測序所得數(shù)據進行匹配獲得蛋白序列注釋，并根據功能進行COG class聚合分析。COG分類的蛋白編碼基因共2 890個，排名前3的分別為：General function prediction only 324個，Cell motility and secretion 312 個，Coenzyme metabolism 283 個，共占比 31.8% (圖3)。

圖3 COG功能分類Fig. 3 COG classification

2.3.2 KEGG 代謝通路分析

LG37全基因組測序獲得的氨基酸數(shù)據經過蛋白注釋后，與KEGG數(shù)據進行比對分析，將目的基因注釋蛋白與想匹配的蛋白功能注釋信息相結合，獲得目的基因的翻譯蛋白信息。具有KEGG的ortholog蛋白編碼基因2 108個，共匹配到160個KEGG代謝通路中。其中排名前6的KEGG通路匹配基因數(shù)分別為：碳水化合物代謝(Carbohydrate metabolism) 和氨基酸代謝 (Amino acid metabolism) 各 157 個，膜轉運 (Membrane transport) 104個，輔酶因子和維生素的代謝(Metabolism of cofactors and vitamins) 93 個，信號轉導 (Single transduction) 90 個，核苷酸代謝 (Nucleotide metabolism) 82 個，共計 683 個，占總數(shù)的32.4%。其中涉及到氮代謝相關的KEGG通路主要集中在GO的Cell Metabolism中，包括：氨基酸代謝，核苷酸代謝，輔酶因子和維生素的代謝，能量代謝 (Energy metabolism)，其他次生代謝產物的生物合成 (Biosynthesis of other secondary metabolites and amino acid metabolism)，詳見圖4。

圖4 LG37基因功能注釋KEGG代謝通路Fig. 4 Gene KEGG pathway classification map of LG37

LG37同化無機氮的主要代謝通路為氮代謝(ko00910)，通過基因組測序分析，共篩選得到18個基因被匹配至氮代謝通路中 (表2)。其中涉及的Pathway ID及匹配的通路基因分別參考 https://www.genome.jp/pathway/map00910和圖5。

表2 LG37 基因組氮代謝通路及其相關基因Table 2 Related genes of nitrogen metabolism pathways of LG37

圖5 氮代謝通路(ko00910)Fig. 5 Nitrogen metabolism (ko00910)

2.4 無機氮代謝相關候選基因的篩選

測序數(shù)據分析結果顯示，共有18個基因涉及到氨氮代謝通路中。然后利用本地Blast+將測序獲得的全基因組信息與NCBI和UniProt數(shù)據庫中無機氮代謝相關數(shù)據進行比對分析，額外篩選出無機氮代謝候選基因76個。其中共篩選出轉運蛋白編碼基因7個，氧化還原蛋白編碼基因13個，調節(jié)蛋白編碼基因5個，結合蛋白和轉錄調控因子編碼基因8個，同化蛋白編碼基因6個及其他無機氮代謝相關基因 (表3)。

表3 無機氮代謝候選基因Table 3 Candidate genes of inorganic nitrogen metabolism

續(xù)表3 Ａto be continued

3 討論

近年，隨著高密度養(yǎng)殖等模式的逐步完善，我國水產養(yǎng)殖單位面積產量呈逐年上升趨勢，但高產出的同時也伴隨著飼料營養(yǎng)的大量投入，其每天產生氨氮的量可根據公式初步計算[23-25]：Q=F×R×0.144 [Q為氨氮產生的總量；F為配合飼料投喂總量 (kg·d?1)；R為飼料中蛋白質所占比例；式中蛋白氮含量約16%，0.144是指90%的被經濟動物吸收的氮作為氨和尿素氮排出 ]。養(yǎng)殖水體中由于4種無機氮之間可相互轉化，且因水體中多種理化指標及微生物等多重作用因素而處于動態(tài)平衡狀態(tài)，因此，雖然水體中NH3和對水生動物有毒害作用，但不能僅降低NH3和的濃度，而應在降低NH3或濃度的同時改良水體中的微生物菌落組成和水體中總無機氮含量[2-3]，實現(xiàn)降低無機氮濃度的目的。

貝萊斯芽孢桿菌不僅可降低無機氮的濃度，而且對養(yǎng)殖水體中病原微生物的生長、繁殖起到很好的抑制效果。雷陽等[26]在對蝦養(yǎng)殖池塘中篩選獲得一株可高效同化氨氮的芽孢桿菌，后經鑒定并命名為貝萊斯芽孢桿菌FS017。Zhu等[27]篩選獲得一株貝萊斯芽孢桿菌CPA1-1，可以很好地降低養(yǎng)殖水體中氨氮和亞硝態(tài)氮的濃度。同時，貝萊斯芽孢桿菌可對副溶血弧菌 (Vibrio parahaemolyticus)、哈維氏弧菌 (V. harveyi)、嗜水氣單胞菌 (Aeromonas hydrophila) 等常見水生動物病原菌起到很好的抑制效果[28]。表明貝萊斯芽孢桿菌在水產養(yǎng)殖領域具有較高的應用價值，因此對其無機氮同化機理進行研究非常必要。

為從基因水平探討貝萊斯芽孢桿菌同化無機氮的機理，本研究對LG37的基因組進行了測序分析。結果顯示，該菌株基因組序列全長 3 929 697 bp，經過基因功能注釋，共有3 057個蛋白具有明確的生物學功能，所有蛋白匹配的同源蛋白來源于75個物種，其中貝萊斯芽孢桿菌FZB42比例最高(78.4%)[29]。通過基因注釋和本地Blast+的篩選、分析，3 057個功能蛋白編碼基因中，與NCBI和UniProt數(shù)據庫中記錄的無機氮代謝相關基因相匹配的共有94個，約占功能基因的3.1%。其中涉及多種特異性無機氮同化功能蛋白的編碼基因，包括感應蛋白 (GlnK)、轉運蛋白(MnrA、NarT和NrgA)、氧化還原蛋白 (NasABD、NarIGHK)等相應編碼基因[7]，其通路下游涉及的同化蛋白和調節(jié)蛋白等不僅參與無機氮的同化代謝，同時也參與有機氮的同化代謝。本結果表明LG37同化無機氮的代謝途徑與有機氮的代謝途徑之間，僅在通路上游的感應、轉運、氧化還原至氨的相關節(jié)點的功能蛋白上具有特異性，而對氨下游涉及的氮化合物合成代謝相關途徑中，不論氮源是無機氮還是有機氮，其代謝途徑均一致，涉及的蛋白的功能和作用也一致。

氮元素在微生物新陳代謝過程中必不可少[30]，進化過程中為適應環(huán)境因子的變化和氮源的多樣性，微生物對氮元素的同化代謝已逐漸發(fā)展成為多元、互補、交叉及全局性的補償性代謝網絡途徑，進而實現(xiàn)對不同種類氮源的利用[31-33]。微生物在代謝過程中涉及的相關功能蛋白，其在基因組上的編碼基因或基因簇往往處于相臨或相近的位置[34-35]，以此提高微生物的代謝效率，實現(xiàn)營養(yǎng)物質的高效利用，同時降低自身的能量損耗[36-37]。本研究發(fā)現(xiàn)，orf00316—orf00630和orf03597—orf03963兩個區(qū)段基因翻譯的蛋白功能主要涉及無機氮的轉運和氧化還原，并最終生成NH3或說明該區(qū)段內的基因翻譯蛋白主要行使的功能是微生物對環(huán)境中無機氮的感知、轉運、氧化還原、吸收同化等；而篩選的其他基因翻譯的蛋白在無機氮或有機氮的代謝途徑中均行使相應的功能。

綜上，貝萊斯芽孢桿菌LG37基因組序列的解析為深入研究其同化無機氮的關鍵功能基因提供了參考依據，為進一步研究其無機氮的代謝途徑提供了理論數(shù)據。后續(xù)將聚焦于orf00316—orf00630和orf03593—orf03963基因區(qū)段的研究，即LG37感知、轉運、氧化還原等無機氮同化途徑。