姚萌 王奎 束長(zhǎng)龍 杜立新 李海濤 丁明月 劉榮梅

（1. 東北農(nóng)業(yè)大學(xué)生命科學(xué)學(xué)院，哈爾濱 150030；2. 中國(guó)農(nóng)業(yè)科學(xué)院植物保護(hù)研究所植物病蟲(chóng)害生物學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，北京 100193；3. 河北農(nóng)林科學(xué)院植物保護(hù)研究所，保定 071000）

作為目前研究最多、應(yīng)用最廣的生防微生物，蘇云金芽胞桿菌（Bacillus thuringiensis，Bt）廣泛分布于土壤、水體、塵埃、病死昆蟲(chóng)、植被等基質(zhì)中，其最大的特點(diǎn)是在形成芽胞的同時(shí)生成伴胞晶體，這種伴胞晶體由不同的殺蟲(chóng)晶體蛋白組成（Insecticidal crystal proteins，ICPs），包括 Cry和Cyt蛋白［1］，對(duì)多種昆蟲(chóng)具有特異的高效殺滅活性，如鱗翅目、鞘翅目、雙翅目、半翅目，及線蟲(chóng)、螨類(lèi)［2-7］等，并且對(duì)非靶標(biāo)天敵安全、對(duì)人畜無(wú)害、不污染環(huán)境，因此Bt廣泛用于防治農(nóng)、林、醫(yī)、倉(cāng)貯等領(lǐng)域的害蟲(chóng)［8］。

目前國(guó)際上Bt商品化產(chǎn)品約有數(shù)百種，其中用于防治鱗翅目害蟲(chóng)的Bt產(chǎn)品主要源自庫(kù)斯塔克亞種（B. thuringiensissubsp.kurstaki，血清型H3a，3b，3c，簡(jiǎn)稱(chēng)Btk）菌株，如HD1菌株（Bt殺蟲(chóng)劑Dipel的生產(chǎn)菌株）。該菌株已經(jīng)在全世界應(yīng)用幾十年，并且是目前用于評(píng)估鱗翅目害蟲(chóng)特異Bt殺蟲(chóng)劑毒力效價(jià)的標(biāo)準(zhǔn)參考菌株［9］。在國(guó)家有關(guān)項(xiàng)目的大力支持下，我國(guó)研究機(jī)構(gòu)也分離出一系列對(duì)鱗翅目高毒力的Bt菌株。例如，河北農(nóng)林科學(xué)院植物保護(hù)研究所分離到的鲇澤亞種（B. thuringiensissubsp.aizawai，血清型7，簡(jiǎn)稱(chēng)Bta）菌株G03，同樣對(duì)多種鱗翅目害蟲(chóng)有效，中國(guó)農(nóng)業(yè)科學(xué)院植物保護(hù)研究所在此基礎(chǔ)上構(gòu)建出了同時(shí)對(duì)鱗翅目和鞘翅目（葉甲科）害蟲(chóng)均具有高毒性的工程菌株G033A［10-11］，并且目前已經(jīng)獲得了農(nóng)藥登記證（PD20171726），商標(biāo)為“禁衛(wèi)軍”。

為了進(jìn)一步分析G03等鱗翅目害蟲(chóng)高效菌株的特點(diǎn)，本研究對(duì)G03菌株進(jìn)行了基因組測(cè)序，并與已公布序列的Bt生產(chǎn)菌株HD1以及Btk模式菌株HD73進(jìn)行比較基因組的研究，相關(guān)研究結(jié)果對(duì)G03應(yīng)用以及今后Bt菌株遺傳改良均具有重要的指導(dǎo)意義。此外，為了分析不同殺蟲(chóng)活性特性Bt菌株之間的進(jìn)化差異，本研究選擇了兩株蠐螬有效的菌株（Bt185、HBF18），以及對(duì)鱗翅目害蟲(chóng)有效但是殺蟲(chóng)活性相對(duì)較差的菌株HD12進(jìn)行了基因組進(jìn)化分析。

1 材料與方法

1.1 供試菌株

G03菌株是實(shí)驗(yàn)室保存菌株。從NCBI數(shù)據(jù)庫(kù)（https：//www.ncbi.nlm.nih.gov/）下載菌株 HD1、HD73，Bt185、HBF18以及HD12基因組數(shù)據(jù)用于后續(xù)分析。

1.2 方法

1.2.1 基因組制備與測(cè)序 按照張彥蕊等方法提?。?2］G03菌株的基因組DNA并檢測(cè)基因組DNA的質(zhì)量，合格后送交公司測(cè)序。利用TruSeq DNA無(wú)PCR文庫(kù)制備試劑盒構(gòu)建測(cè)序文庫(kù)，并交由華大基因公司利用Illumina HiSeq 2500 測(cè)序平臺(tái)測(cè)序。

1.2.2 測(cè)序數(shù)據(jù)的組裝、注釋 對(duì)測(cè)序獲得的原始Reads進(jìn)行質(zhì)量控制，刪除低質(zhì)量的Reads，利用IDBA_UD軟件包［13］對(duì)剩余高質(zhì)量的pair-end clean reads數(shù)據(jù)進(jìn)行基因組序列組裝。并利用NCBI 原核生物基因組注釋流程PGAP（Prokaryotic Genome Annotation Pipeline，https：//www.ncbi.nlm.nih.gov/genome/annotation_prok/）對(duì)組裝基因組數(shù)據(jù)進(jìn)行基因預(yù)測(cè)與注釋。

利用泛基因組分析PGAP軟件包［14］進(jìn)行基因家族聚類(lèi)。利用維恩圖在線繪制網(wǎng)站（http：//bioinfogp.cnb.csic.es/tools/venny/）進(jìn)行菌株間共享基因及特有基因的統(tǒng)計(jì)分析。利用Blast軟件包［15］、TBtools軟件［16］進(jìn)行Go功能富集。并利用Blast軟件包［15］進(jìn)行殺蟲(chóng)基因及功能基因的注釋。利用MEGAX軟件包［17］和CVTree軟件包［18］分別進(jìn)行蛋白質(zhì)系統(tǒng)發(fā)生分析和全基因組系統(tǒng)發(fā)生分析。

2 結(jié)果

2.1 基因組測(cè)序、組裝與注釋

利用IDBA_UD軟件包對(duì)G03菌株進(jìn)行基因組組裝。結(jié)果顯示，基因組大小為6.38 M，組裝contig數(shù)量為395，GC含量為34.71%，Contig N50為45.68 kb，Contig L50為41。相關(guān)組裝結(jié)果已經(jīng)提交到NCBI GenBank數(shù)據(jù)庫(kù)，基因組登錄號(hào)為NHNQ00000000.1。進(jìn)一步利用NCBI 原核生物基因組注釋流程PGAP對(duì)基因組進(jìn)行了注釋?zhuān)Y(jié)果顯示G03菌株包含6916個(gè)基因，其中蛋白質(zhì)編碼基因（Coding sequence，CDS）6 848 個(gè)。

2.2 CVTree系統(tǒng)進(jìn)化分析

CVTree（Composition vector tree，簡(jiǎn)稱(chēng) CVTree）是郝柏林研究組于2003年提出的一種全新的用于推測(cè)物種之間親緣關(guān)系與分類(lèi)的研究方法［18］。這種方法基于全基因組，不需要進(jìn)行序列比對(duì)，根據(jù)全基因組序列中不同核苷酸堿基或氨基酸殘基短串的出現(xiàn)頻率來(lái)構(gòu)建系統(tǒng)發(fā)育樹(shù)。迄今為止，CVTree已經(jīng)廣泛用于分析多種不同物種及數(shù)據(jù)，包括古生菌［19-20］、原核生物［18，21］、真菌［22］、病毒［23］、葉綠體序列［24］、tRNA 序列［25］、癌癥序列［26］等。

我們利用CVTree對(duì)包括G03、HD1、HD73菌株在內(nèi)的7株Bt菌株的全基因組構(gòu)建組分矢量并作出Neighbor-joining系統(tǒng)發(fā)生樹(shù)（氨基酸短串K值取6）。結(jié)果（圖1）發(fā)現(xiàn)，Btk菌株HD1和HD73聚在同一個(gè)分支，菌株G03與Btk菌株HD1和HD73所在分支相聚較近，屬于同一大分支。而同樣對(duì)鱗翅目害蟲(chóng)有效的菌株HD12則相距較遠(yuǎn)。此外，對(duì)鞘翅目害蟲(chóng)有效的菌株Bt185與HBF18與G03所在分支也有較大的距離。

圖1 菌株CVTree系統(tǒng)進(jìn)化分析樹(shù)圖

2.3 比較基因組分析

我們對(duì)CVTree分析中親緣關(guān)系較近的Bta菌株G03以及Btk菌株HD1、HD73這3株Bt菌株進(jìn)行了比較基因組的分析，首先利用PGAP軟件包對(duì)其進(jìn)行基因家族聚類(lèi)，并利用維恩圖在線繪制網(wǎng)站（http：//bioinfogp.cnb.csic.es/tools/venny/）對(duì)這 3株菌株的共享基因及特有基因家族進(jìn)行繪圖統(tǒng)計(jì)分析，如圖2所示。PGAP分析共鑒定出了6 860個(gè)基因家族，其中核心基因家族共有 4 720個(gè)（全部菌株均含有的基因家族），占全部基因家族的68.8%。而G03、HD1和HD73菌株分別有540、433和327個(gè)特有基因（僅在一株菌株中存在的基因）。G03與HD1有5122基因家族是共有的，這大于HD1與HD73共有基因家族（5101）。

2.4 菌株功能分類(lèi)比較分析

我們利用Blast程序包分別將G03和HD1菌株的全部基因、共有基因以及特有基因序列與uniprot數(shù)據(jù)庫(kù)（https：//www.uniprot.org/）進(jìn)行序列比對(duì)，并利用TBtools軟件進(jìn)行Go功能注釋以及l(fā)evel 2水平上的功能解析。圖3統(tǒng)計(jì)了菌株G03基因組Go功能注釋在level 2水平上的統(tǒng)計(jì)情況，紅色是菌株G03基因組特有基因的統(tǒng)計(jì)結(jié)果，藍(lán)色部分是G03和HD1菌株都有的基因的統(tǒng)計(jì)結(jié)果。

圖2 3株Bt菌株共享及特有基因的維恩圖分析

結(jié)果分析發(fā)現(xiàn)生物學(xué)途徑（biological process）相關(guān)功能部分，G03特有基因在cellular process和metabolic process兩個(gè)方面數(shù)量最多；在細(xì)胞學(xué)組件（cellular component）相關(guān)功能部分，G03特有基因在cell和cell part兩個(gè)方面數(shù)量最多；在分子功能（molecular function）方面，G03特有基因在binding和catalytic activity兩個(gè)方面數(shù)量最多。其中生物學(xué)途徑與分子功能相關(guān)的Go富集反映了G03在生物學(xué)水平與分子水平上的特點(diǎn)。

cellular process是指在細(xì)胞水平上發(fā)生的過(guò)程相關(guān)功能富集。例如，細(xì)胞對(duì)環(huán)境的應(yīng)答過(guò)程，在該方面的功能有較多Go富集顯示菌株G03在環(huán)境適應(yīng)等細(xì)胞過(guò)程方面有較多的特點(diǎn)。metabolic process是指生物體轉(zhuǎn)化化學(xué)物質(zhì)的化學(xué)反應(yīng)和途徑相關(guān)功能富集，包括合成代謝和分解代謝，這方面較多Go富集說(shuō)明G03菌株在物質(zhì)分解利用或化合物合成方面有較特別的能力。

Binding和catalytic activity分別指結(jié)合與催化方面的功能富集，這兩個(gè)功能是緊密相關(guān)的，G03特有基因在這兩個(gè)方面有較多的富集，顯示其在分子識(shí)別催化方面有特別的功能。

2.5 G03菌株殺蟲(chóng)基因分析

我們對(duì)Bta及Btk菌株包含的殺蟲(chóng)基因種類(lèi)以及數(shù)量做了比較分析，HD73菌株僅有cry1Ac一個(gè)殺蟲(chóng)基因，菌株G03與HD1均含有cry1Aa、cry1Ac、cry1Ia、cry2Ab、vip3Aa這5個(gè)高毒力殺蟲(chóng)基因［3，27-29］。已有研究表明 Cry1A、Cry1I、Cry2A、Vip3A這幾類(lèi)殺蟲(chóng)蛋白具有不同的受體，并且殺蟲(chóng)蛋白之間還有協(xié)同增效的作用，這可能是G03與HD1具有高活性的原因。此外，與HD1相比，G03菌株還含有cry1Ca、cry1Da、cry9Ea基因。

圖3 G03菌株Go功能注釋

2.6 G03、HD1菌株特有基因分析

利用基因家族分析的方法，對(duì)菌株G03與 HD1進(jìn)行比較分析，結(jié)果顯示G03有597個(gè)特有基因家族，而HD1有814個(gè)特有基因家族。其中，G03菌株特有的基因家族中有232個(gè)可以注釋到功能，而HD1菌株特有的基因家族中有404個(gè)可以注釋到功能。進(jìn)一步，對(duì)這些注釋到功能的基因家族進(jìn)行分析，發(fā)現(xiàn)噬菌體與轉(zhuǎn)座酶相關(guān)的基因家族最為豐富。HD1菌株特有的基因家族中有34個(gè)基因家族與噬菌體相關(guān)，而G03菌株特有的基因家族中只有 10個(gè)與噬菌體有關(guān)。

此外，特有基因家族中，HD1菌株轉(zhuǎn)座酶基因也比較多，有26個(gè)基因家族屬于轉(zhuǎn)座酶，而G03菌株特有的基因家族中只有8個(gè)。轉(zhuǎn)座酶是轉(zhuǎn)座子執(zhí)行轉(zhuǎn)座功能的酶，通常由轉(zhuǎn)座子編碼。本研究進(jìn)一步比較了兩個(gè)菌株含有轉(zhuǎn)座酶的數(shù)量與多樣性。結(jié)果顯示，G03菌株僅僅含有33個(gè)轉(zhuǎn)座酶或片段，而HD1多達(dá)216個(gè)轉(zhuǎn)座酶或片段。這些轉(zhuǎn)座酶按照序列一致性（按照序列一致性50%設(shè)定閾值）分為58個(gè)家族，其中16個(gè)基因家族是兩個(gè)菌株共有的。這些基因家族中，最大的基因家族含有45個(gè)轉(zhuǎn)座酶基因，其中有4個(gè)來(lái)源于菌株G03，41個(gè)來(lái)源于菌株HD1。該家族轉(zhuǎn)座酶含有TNP1 DDE 結(jié)構(gòu)域，是 IS4、IS421、IS5377、IS427、IS402、IS1355 以及IS5等多種插入序列（Insert sequence，IS）的轉(zhuǎn)座酶。進(jìn)一步用UPGMA對(duì)這些含有TNP1 DDE 結(jié)構(gòu)域的轉(zhuǎn)座酶進(jìn)行系統(tǒng)進(jìn)化分析（去除了4個(gè)轉(zhuǎn)座酶片段），結(jié)果（圖 4）顯示這些轉(zhuǎn)座酶分為3個(gè)類(lèi)群，其中類(lèi)群1是兩個(gè)菌株共有的；類(lèi)群2 有6個(gè)基因，是HD1特有的；類(lèi)群3 只有1個(gè)基因，是G03特有的。

噬菌體與插入序列、轉(zhuǎn)座子是微生物獲得外源基因的重要途徑，同時(shí)也是微生物基因組變異的重要因素。上述分析結(jié)果表明G03含有較少的噬菌體與插入序列、轉(zhuǎn)座子相關(guān)的基因，顯示菌株G03與HD1相比，應(yīng)該具有較低的變異概率。

3 討論

目前市場(chǎng)上常用的Bt生產(chǎn)菌株種類(lèi)相對(duì)較少，主要以Btk菌株HD1為主，相關(guān)研究報(bào)道較多。Bta菌株G03是我國(guó)科學(xué)家自主分離并且已經(jīng)用于害蟲(chóng)防控的菌株，已通過(guò)基因工程改造拓展了其殺蟲(chóng)譜，在鱗翅目和鞘翅目（葉甲科）害蟲(chóng)防治中發(fā)揮重要作用，然而該菌株遺傳學(xué)、基因組相關(guān)研究鮮見(jiàn)報(bào)道。本研究利用Illumina測(cè)序技術(shù)進(jìn)行了菌株G03基因組測(cè)序，并進(jìn)行了基因功能注釋。在此基礎(chǔ)上與NCBI數(shù)據(jù)庫(kù)中已經(jīng)公開(kāi)的Btk生產(chǎn)菌株HD1、Btk模式菌株HD73、兩株蠐螬有效的菌株（Bt185、HBF18）以及含有多種殺蟲(chóng)基因類(lèi)型的菌株HD12進(jìn)行了比較基因組研究，旨在從基因組角度分析G03等鱗翅目害蟲(chóng)高效菌株的特點(diǎn)、比較與Btk菌株的異同，為進(jìn)一步對(duì)Bt菌株進(jìn)行改良、推進(jìn)我國(guó)自主分離研發(fā)的Bt菌株及其殺蟲(chóng)劑產(chǎn)品的應(yīng)用奠定基礎(chǔ)。

圖4 含有TNP1 DDE結(jié)構(gòu)域轉(zhuǎn)座酶進(jìn)化關(guān)系分析。箭頭標(biāo)出的是G03轉(zhuǎn)座酶

本研究首次利用基因組測(cè)序技術(shù)，從基因組水平上全面的分析了G03殺蟲(chóng)相關(guān)的基因。與PCR殺蟲(chóng)基因鑒定技術(shù)相比，基因組測(cè)序技術(shù)避免了PCR鑒定技術(shù)由于引物設(shè)計(jì)不全面、PCR擴(kuò)增條件不合適等原因，導(dǎo)致基因鑒定的誤差。通過(guò)基因組分析發(fā)現(xiàn)，菌株G03含有8個(gè)殺蟲(chóng)基因（cry1Aa、cry1Ac、cry1Ca、cry1Da、cry1Ia、cry2Ab、cry9Ea、vip3Aa），其中cry1Ca、cry1Da是菌株HD1不具備的。目前，表達(dá)Cry1A、Cry2A和Vip3A等蛋白的轉(zhuǎn)基因的抗蟲(chóng)作物已經(jīng)在多個(gè)國(guó)家商業(yè)化種植（International Service for the Acquisition of Agri-biotech Applications，www.isaaa.org/），使用表達(dá)與上述蛋白無(wú)交互抗性的Bt殺蟲(chóng)劑在害蟲(chóng)抗性治理方面具有重要的意義。利用小菜蛾抗Cry1Ac種群研究表明，Cry1Ca與Cry1Ac蛋白沒(méi)有交互抗性［30］，因此，菌株G03比HD1在害蟲(chóng)的Bt殺蟲(chóng)劑抗性風(fēng)險(xiǎn)控制方面更具優(yōu)勢(shì)。

Bt菌株中與殺蟲(chóng)活性相關(guān)的免疫抑制因子A（Immune Inhibitor A，InhA）、 雙 效 菌 素 ZwA（Zwittermicin A，ZwA）、酰基高絲氨酸內(nèi)酯水解酶（AHLs Hydrolase，aiiA）以及幾丁質(zhì)酶（Chitinase）對(duì)Bt菌株殺蟲(chóng)活性及環(huán)境適應(yīng)性都有重要意義。InhA是Bt分泌的一種金屬蛋白酶，能夠降解昆蟲(chóng)產(chǎn)生的抗菌肽從而逃避宿主的免疫系統(tǒng)［31-32］。aiiA可以水解細(xì)菌的群體感應(yīng)相關(guān)的信號(hào)分子?；呓z氨酸內(nèi)酯（Acylated Homoserine Lactones，AHLs），對(duì)多種細(xì)菌有抑制作用，因此有助于提高Bt在昆蟲(chóng)腸道內(nèi)競(jìng)爭(zhēng)優(yōu)勢(shì)［33-34］。Chitinase是一種可溶性的胞外蛋白質(zhì)類(lèi)殺蟲(chóng)活性物質(zhì)，可以幫助Bt菌株降解昆蟲(chóng)腸道圍食膜中的幾丁質(zhì)成分，使其能通過(guò)穿孔的腸道進(jìn)入血腔引起昆蟲(chóng)敗血癥，對(duì)Bt殺蟲(chóng)蛋白具有增效作用［35-36］。前期研究發(fā)現(xiàn)，菌株G03與HD1都含inhA、aiiA、chitinase基因以及合成ZwA的基因簇，并且aiiA、chitinase基因進(jìn)行系統(tǒng)發(fā)生分析發(fā)現(xiàn)鱗翅目害蟲(chóng)高效的菌株aiiA、chitinase基因序列比較接近，說(shuō)明G03菌株在感染寄主、適應(yīng)環(huán)境方面與HD1類(lèi)似［37］。

此外，本研究還發(fā)現(xiàn)與HD1相比，G03含有較少的噬菌體與轉(zhuǎn)座子基因，說(shuō)明G03基因組應(yīng)該具有更好的遺傳穩(wěn)定性。作為生產(chǎn)菌株，遺傳穩(wěn)定性對(duì)產(chǎn)品發(fā)酵生產(chǎn)尤為重要，基因組比較分析結(jié)果顯示G03菌株在該方面更有優(yōu)勢(shì)。

綜上所述，本研究對(duì)生產(chǎn)菌株G03進(jìn)行的基因組測(cè)序、注釋與比較基因組分析，不僅可以從殺蟲(chóng)基因角度揭示生產(chǎn)菌株高活性的原因，還可以為進(jìn)一步菌株改良提供基因組數(shù)據(jù)支撐。