短短芽胞桿菌B011基因組中芽胞形成調(diào)控基因spo0E的鑒定與表達(dá)分析

周向平，王運(yùn)生，陳 武，劉天波，王凱歌，劉 峰，周 立，袁志輝

（1.湖南農(nóng)業(yè)大學(xué)植物保護(hù)學(xué)院，長沙 410128；2.湖南省煙草公司永州市公司，永州 425100；3.湖南省煙草科學(xué)研究所，長沙 410004；4.湖南科技學(xué)院化學(xué)與生物工程學(xué)院，永州 425199；5.湖南省銀杏工程技術(shù)研究中心，永州 425199）

形成芽胞是產(chǎn)芽胞類細(xì)菌在長期的自然進(jìn)化過程中形成的抵抗逆境的能力體現(xiàn)。當(dāng)芽胞類細(xì)菌感知到逆境脅迫時(shí)，如細(xì)胞密度過高、營養(yǎng)物質(zhì)匱乏等，細(xì)胞首先一般可通過三種途徑來減輕環(huán)境壓力，一是激活運(yùn)動(dòng)裝置，如鞭毛蛋白，通過趨化性尋找新的食物；二是分泌抗生素、抗菌肽和其他化學(xué)物質(zhì)，清除或驅(qū)除同一生境中的競爭性微生物；三是分泌大量水解酶類，降解胞外蛋白質(zhì)和糖類[1]。對枯草芽胞桿菌Bacillussubtilis168 的研究顯示，如果上述3 種途徑?jīng)]有減輕環(huán)境壓力，50%～70%的細(xì)胞將形成芽胞[2]。芽胞形成屬于群體感應(yīng)行為，即細(xì)胞向胞外分泌信號分子，細(xì)胞能感應(yīng)這些信號分子的濃度（與細(xì)胞密度成正比），當(dāng)達(dá)到一定閾值后，通過信號轉(zhuǎn)導(dǎo)通路激活或抑制相關(guān)基因的表達(dá)，使細(xì)胞進(jìn)入不同發(fā)育途徑[3]。

spo0A是調(diào)控芽胞形成的開關(guān)基因[4]，調(diào)控因子Spo0A 磷酸化為Spo0A～P，Spo0A～P 是Spo0A 的活性形態(tài)，是芽胞形成的關(guān)鍵調(diào)控因子[5]，一方面，Spo0A～P 可以正調(diào)控spoIIE基因，正式啟動(dòng)芽胞形成途徑[6]，另一方面Spo0A～P 負(fù)調(diào)控abrB基因[7]，AbrB 是一個(gè)全局性轉(zhuǎn)錄負(fù)調(diào)控因子，同時(shí)Spo0A～P 還可以正調(diào)控abbA基因，而AbbA 可抑制AbrB 活性[8]，AbrB 可負(fù)調(diào)控群體感應(yīng)關(guān)鍵基因comK[9]。細(xì)胞內(nèi)Spo0A～P 的濃度水平?jīng)Q定了細(xì)胞的發(fā)育途徑，Spo0E 是一類磷酸酶，可將活性態(tài)Spo0A 即Spo0A～P，轉(zhuǎn)為非活性態(tài)Spo0A，從而起到負(fù)調(diào)控胞內(nèi)Spo0A～P 的濃度[10,11]，同時(shí)spo0E受AbrB 負(fù)調(diào)控，Spo0A-AbrB-Spo0E 形成一個(gè)調(diào)控回路[12]。Spo0E在芽胞桿菌基因組中多以基因家族（SMFP,spo0Emultigene family of phosphatases）存在，基因上游一般存在多個(gè)啟動(dòng)子（promoter）區(qū)域，可受不同轉(zhuǎn)錄因子調(diào)控[13]。

短短芽胞桿菌BrevibacillusbrevisB011 菌株是本實(shí)驗(yàn)室從煙草根系篩選到的內(nèi)生細(xì)菌，其對煙草青枯病菌、煙草黑脛病菌、辣椒白絹病菌等均有較明顯的抑制效果，通過形態(tài)鑒定、生理生化試驗(yàn)、16S rDNA測序及比對分析，確定該生防菌株為短短芽胞桿菌[14]。本文在對短短芽胞桿菌B011 菌株進(jìn)行全基因組測序的基礎(chǔ)上，對不同生長階段的B011 菌株進(jìn)行轉(zhuǎn)錄組測序，通過生物信息學(xué)工具預(yù)測得到了9 個(gè)spo0E同源基因，并通過與其他101 個(gè)短芽胞桿菌Brevibacillus菌株全基因組進(jìn)行比較分析，結(jié)合轉(zhuǎn)錄組數(shù)據(jù)，推測spo0E同源基因功能，為進(jìn)一步驗(yàn)證基因功能和B011 菌株的遺傳改造提供候選基因。

1 材料與方法

1.1 短芽胞桿菌基因組序列分析

短短芽胞桿菌B011 菌株為本實(shí)驗(yàn)室從土壤環(huán)境中分離獲得的菌株，經(jīng)鑒定其為短短芽胞桿菌[14]。對其進(jìn)行了全基因組測序，并提交序列至NCBI，GenBank 登錄號為CP041767。從NCBI RefSeq 數(shù)據(jù)庫[15]（https://www.ncbi.nlm.nih.gov/）下載所有短芽胞桿菌已經(jīng)公布全基因組序列的菌株核苷酸序列和蛋白序列，統(tǒng)計(jì)各個(gè)菌株基因組的GC 含量、蛋白編碼基因個(gè)數(shù)、基因組大小等信息。

1.2 短短芽胞桿菌基因家族分析

用OrthoMCL[16]構(gòu)建基因家族，選取已公布完成圖的4個(gè)短短芽胞桿菌菌株（B011、X23、NBRC 100599、DZQ7）用來做基因家族分析，OrthoMCL 首先通過blastp 對所有蛋白序列進(jìn)行兩兩比對，將任何兩兩比對e-value 小于1E-5 的序列保留下來，然后根據(jù)score 值進(jìn)行MCL 聚類，聚在同一個(gè)cluster 的基因被認(rèn)為是同一個(gè)基因家族。

1.3 Spo0E 家族基因預(yù)測

從Pfam 數(shù)據(jù)庫[17]（https://pfam.xfam.org）下載SpoOE-like (PF09388)模型文件，然后用hmmsearch[18]搜索上述Brevibacillus蛋白庫，取E 值小于0.05 的蛋白為候選基因，然后經(jīng)BLASTp[19]進(jìn)一步驗(yàn)證NR 數(shù)據(jù)庫中相似序列是否含有Spo0E功能描述。

1.4 多序列比較和構(gòu)建進(jìn)化樹

用MAFFT[20]進(jìn)行多序列比對，參數(shù)選用自動(dòng)模式（-auto），然后用FastTree[21]構(gòu)建系統(tǒng)進(jìn)化樹，參數(shù)用默認(rèn)參數(shù)，進(jìn)化樹上傳到在線分析平臺iTOL[22]進(jìn)行進(jìn)一步分析與顯示作圖。

1.5 基序（motif）分析

基序采用MEME Suite 軟件包[23]進(jìn)行分析，選取輸出前5 個(gè)基序。MEME Suite 是集合眾多預(yù)測和注釋motif 工具的在線網(wǎng)站，其中MEME 算法是基于最大期望值（EM）算法來識別motif。保守基序用seqLogo[24]作圖。

1.6 B011 菌株中spo0E 類似基因的表達(dá)情況分析

B011 菌株在LB 培養(yǎng)液中，37 ℃、200 r/min 條件下分別振蕩培養(yǎng)12、18、30、36 h 后，送交測序公司進(jìn)行轉(zhuǎn)錄組測序。芽胞形成調(diào)控相關(guān)基因相對表達(dá)量采用RNA-Seq 的TPM（Transcripts Per Million）值進(jìn)行表征，TPM 與RPKM/FPKM 類似，是將原始reads 除去基因長度的影響后標(biāo)準(zhǔn)化到每百萬條reads，相當(dāng)于重新標(biāo)準(zhǔn)化的文庫，保證每個(gè)樣本中所有TPM 的總和是相同的[25]。RNA-Seq 分析流程可簡單描述為：首先將質(zhì)控合格的reads 回帖到參數(shù)基因組，然后用TPMCalculator[26]計(jì)算每個(gè)基因的TPM 值。

2 結(jié)果與分析

2.1 短芽胞桿菌基因組序列及特征

截止2020 年8 月1 日，在NCBI 的RefSeq 數(shù)據(jù)庫中已收錄并公布全基因組序列的短芽胞桿菌菌株有101 個(gè)，加上本實(shí)驗(yàn)室的B011 菌株，一共包含了102 個(gè)菌株的全基因組序列。包含有20 個(gè)種，其中側(cè)孢短芽胞桿菌B.laterosporus23 個(gè)，短短芽胞桿菌13 個(gè)，波茨坦短芽胞桿菌B.borstelensis7 個(gè)，土壤短芽胞桿菌B.agri5 個(gè)，有26 個(gè)菌株沒有種分類信息。短芽胞桿菌基因組大小范圍為4.00～7.13 Mb，GC 含量范圍為40%～8.5%。從圖1 可以看出，不同種之間基因組大小和GC 含量差異明顯，如側(cè)孢短芽胞桿菌基因組GC 含量偏低（40%～41.3%），基因組大小也偏小（4.53～5.80 Mb），而短短芽胞桿菌基因組GC含量偏高，為46.8%～47.5%，基因組大小也偏大，為6.02～6.73 Mb，在基因組大小-GC 含量圖上能明顯分開（圖1）。而短短芽胞桿菌，抗生素短芽胞桿菌B.antibioticus, 美麗短芽胞桿菌B.formosus, 強(qiáng)音短芽胞桿菌B.fortis和波塔短芽胞桿菌B.porteri聚在一起，這與平均核苷酸一致性（ANI, Average Nucleotide Identity）和基于全基因組聚類結(jié)果一致（結(jié)果未公布），表明這幾個(gè)種親緣關(guān)系較近。

圖1 102 個(gè)Brevibacillus 菌株的聚類分布圖Fig.1 The cluster distribution of 102 Brevibacillus genomes using the genome size as X-axis and GC content as Y-axis

2.2 短短芽胞桿菌基因家族分析

比較了4 個(gè)短短芽胞桿菌基因組，用OrthoMCL 對4 個(gè)短短芽胞桿菌的編碼基因進(jìn)行了聚類，構(gòu)建基因家族。首先對4 個(gè)菌株的蛋白序列用BLAST 進(jìn)行兩兩比對，然后用MCL 進(jìn)行聚類，共構(gòu)建了5933 個(gè)基因家族，其中4749 個(gè)基因家族在4 個(gè)菌株都存在，可視為短短芽胞桿菌的核心基因組，其中4718 個(gè)基因家族為單拷貝基因家族，每個(gè)基因組只含一份（表1，圖2）。

表1 4 個(gè)短短芽胞桿菌基因家族概況Table 1 The gene families of the 4 B.brevis with complete genomic sequences

圖2 短短芽胞桿菌基因家族比較分析Fig.2 OrthoMCL gene families of 4 Brevicillus brevis

2.3 Spo0E 基因預(yù)測

從102 個(gè)短芽胞桿菌菌株基因組中一共預(yù)測得到837 個(gè)spo0E基因，平均每個(gè)菌株含8.2 個(gè)spo0E類似基因，最少的僅有3 個(gè)，最多的有23 個(gè)。不同種之間spo0E同源基因數(shù)量具有一定差異，側(cè)孢短芽胞桿菌平均含有13 個(gè)，變化范圍為8～23，而短短芽胞桿菌平均含有8 個(gè)spo0E基因，變化范圍為6～10 個(gè)（圖3）。Spo0E 蛋白因子的平均長度為68 個(gè)氨基酸，最小的24 個(gè)氨基酸，最長的為162 個(gè)氨基酸，但不同種之間長度差異不顯著。菌株B011 預(yù)測含有9 個(gè)spo0E類似基因，編碼肽鏈的長度范圍為52～101個(gè)氨基酸。這9 個(gè)基因在基因組上為分散分布，多數(shù)與旁邊基因并不形成操縱子結(jié)構(gòu)。

圖3 不同Brevibacillus 種spo0E 類似基因個(gè)數(shù)Fig.3 The number of spo0E-like genes in different species of Brevibacillus

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

從837 個(gè)spo0E類似基因的系統(tǒng)進(jìn)化樹上可以看出，這些基因明顯聚成6 個(gè)大的分枝，B011 菌株的9個(gè)基因分散在不同的5 個(gè)分枝上，且每個(gè)分枝上均分布有不同的種，說明這些基因在種分化之前就已經(jīng)存在了。與其他種不同的是，側(cè)孢短芽胞桿菌趨向于分布在單獨(dú)的分枝上，特別是如圖4 所示的分枝4 絕大部分為側(cè)孢短芽胞桿菌成員，說明該分枝上的基因可能為側(cè)孢短芽胞桿菌特異基因。

圖4 短桿芽胞桿菌spo0E 類似基因系統(tǒng)進(jìn)化樹Fig.4 The phylogenetic tree of 837 spo0E-like genes from 102 Brevibacillus species

2.5 短桿芽胞桿菌Spo0E 編碼蛋白基序特征

經(jīng)MEME 的基序分析，發(fā)現(xiàn)短短桿芽胞桿菌的837 個(gè)spo0E基因編碼蛋白序列中有825（98.6%）含有基序1，基序1 當(dāng)中最保守的序列為SQQLD，這與已報(bào)道的芽胞桿菌Spo0E 的SQELD 基序高度保守，特別是基序中的第8 位S 和第12 位D 在所有序列中高度保守，推測可能為其活性位點(diǎn)。菌株B011 的9個(gè)基因中，有5 個(gè)為SQQLD，2 個(gè)為SQELD，1 個(gè)為SVQLD，1 個(gè)為SIQID（圖5）?；? 位于基序1的前面，但缺乏完全保守的位點(diǎn)，在673 個(gè)基因中存在該基序。這兩個(gè)基序均存在于B011 菌株的9 個(gè)基因中。

圖5 短桿芽胞桿菌Spo0E 蛋白因子基序特征及B011 菌株的9 個(gè)Spo0E 類似蛋白序列多序列比對Fig.5 The main two motifs in the Spo0E-like proteins of Brevibacillus and multiple sequences alignment of that of strain B011

2.6 B011 菌株中spo0E 類似基因的表達(dá)情況

B011 菌株的9 個(gè)spo0E同源基因除了FO446_11425在5 個(gè)轉(zhuǎn)錄組樣本中沒有檢測到外，其他8 個(gè)基因均檢測到有表達(dá)（圖6），其中基因FO446_19435的表達(dá)量非常高，在發(fā)酵到18 h 達(dá)到最高峰，18 h為其對數(shù)生長期，其TPM 值達(dá)到了5013，隨后下降。基因FO446_04050表達(dá)模式與FO446_19435類似，均在18 h 達(dá)到最大值，表達(dá)量也較高，18 h 時(shí)TPM 值達(dá)到了1003。而FO446_23370、FO446_27395兩個(gè)基因表達(dá)模式相似，但在18 h 時(shí)表達(dá)量與30、36 h 時(shí)的測量值相當(dāng)，屬于中等表達(dá)豐度基因，另外4 個(gè)基因的表達(dá)量較低，其中基因FO446_01695在12 h 表達(dá)量相對更高，隨后下降。

圖6 菌株B011 菌株中9 個(gè)spo0E 類似基因表達(dá)熱圖Fig.6 The gene expression heatmap of 9 spo0E-like genes in strain B011

3 討論

短短芽胞桿菌是能夠產(chǎn)生多種植物病原菌拮抗物質(zhì)的重要生防菌，其活性物質(zhì)鑒定和基因組測序分析已有多個(gè)報(bào)道[27,28]。從已報(bào)道的具有完成圖的短短芽胞桿菌基因組來看，B011 菌株的基因組大小是最小的，蛋白編碼基因也是最少的，預(yù)測有5663 個(gè)蛋白編碼基因，而本實(shí)驗(yàn)室篩選到的另一個(gè)生防菌X23 菌株同為短短芽胞桿菌[29]，基因組大小為6643473，編碼6367 個(gè)蛋白，比B011 菌株基因組大了約500000，蛋白編碼基因也相應(yīng)多了約700 個(gè)基因。B011 菌株是從植物組織內(nèi)分離得到的植物內(nèi)生菌，而X23 菌株是從土壤中分離得到的土棲菌，植物內(nèi)生菌由于所處的生存環(huán)境相對于土壤來說要穩(wěn)定和簡單，對于植物內(nèi)生菌來說，某些基因可能對其生存作用不大，這部分基因在進(jìn)化的過程中就有可能丟掉或者退化，這也與植物內(nèi)生菌或者病原微生物的基因組相對于土棲菌來說，其基因組相對來說要小一些的現(xiàn)象是一致的。

在芽胞形成調(diào)控基因方面，B011 菌株中含有9 個(gè)spo0E家族基因，除了FO446_11425在5 個(gè)轉(zhuǎn)錄組樣本中沒有表達(dá)外，其余8 個(gè)基因均有表達(dá)，說明這些基因可能都在行使功能。Spo0E基因在產(chǎn)芽胞細(xì)菌中一般是以多基因家族存在的，相互之間存在功能冗余，如在枯草芽胞桿菌中敲除spo0E后對其芽胞形成并沒有明顯影響，因?yàn)樵诳莶菅堪麠U菌中還存在另外兩個(gè)spo0E同源基因yisI和ynzD，具有類似功能[30]。在所有B011 的9 個(gè)spo0E同源基因中均存在芽胞桿菌SQELD 類似基序，說明SQELD 基序可能是芽胞桿菌共有基序，三維結(jié)構(gòu)顯示，基序SQELD 位于活性位點(diǎn)中心[31]，其中S 和D 高度保守，中間三個(gè)氨基酸QEL 在不同菌株間存在一定變異，如在B011 菌株中有5 個(gè)為SQQLD，但表達(dá)量最高的兩個(gè)spo0E同源基因FO446_19435和FO446_04050均為SQELD。B011 菌株中這9 個(gè)spo0E同源基因是否具有功能互補(bǔ)？在體內(nèi)又是如何進(jìn)行調(diào)控的需進(jìn)一步研究。

鑒于基因FO446_19435表達(dá)量非常高，且表達(dá)模式是先升高后下降，符合spo0E的一般表達(dá)模式，Spo0E 的生物學(xué)功能是對Spo0A～P 進(jìn)行去磷酸化，從而降低細(xì)菌體內(nèi)Spo0A～P 的含量，抑制芽胞的形成。推測基因FO446_19435為B011 菌株的主效spo0E基因，但同時(shí)基因FO446_04050可能是spo0E的候選補(bǔ)充基因，它們的表達(dá)模式類似，且表達(dá)量也較高。FO446_19435與下游基因FO446_19440屬于同一個(gè)操縱子，與FO446_19435類似，基因FO446_19440的表達(dá)量非常高，18 h 的TPM 達(dá)到了27181，但基因FO446_19440功能未知，且只在菌株短短芽胞桿菌GZDF3.1 中找到有同源基因，但奇怪的是與FO446_19440同為一個(gè)操縱子的FO446_19435并未在短短芽胞桿菌 GZDF3.1 找到同源基因。

枯草芽胞桿菌中spo0E受到多個(gè)轉(zhuǎn)錄因子調(diào)控，受AbrB 和Rok 負(fù)調(diào)控，同時(shí)還受sigma B 轉(zhuǎn)錄因子正調(diào)控[32]。在B011 菌株中，abrB(FO446_00625)基因表達(dá)隨著發(fā)酵時(shí)間是升高的，發(fā)酵后期表達(dá)量明顯上調(diào)，總體趨勢是與spo0E的表達(dá)相反，說明與枯草芽胞桿菌類似，B011 菌株的spo0E可能也受AbrB 負(fù)調(diào)控，而sigma B (FO446_02910)轉(zhuǎn)錄因子的表達(dá)量變化不大，sigma B 轉(zhuǎn)錄因子主要與脅迫有關(guān)[33]，由于本試驗(yàn)轉(zhuǎn)錄組樣本沒有脅迫因素，可能是導(dǎo)致sigmaB轉(zhuǎn)錄因子表達(dá)量變化不大的原因。B011 菌株中缺少rok基因，rok基因主要與感受態(tài)細(xì)胞形成相關(guān)[34]，而B011 菌株中同時(shí)也缺少感受態(tài)細(xì)胞形成關(guān)鍵基因comK?？紤]到B011 菌株的spo0E基因在發(fā)酵18 h 時(shí)達(dá)到最高，而abrB的表達(dá)量是一直上升的，所以推測在B011菌株中，spo0E基因除了受AbrB 調(diào)控外，可能還受到其他基因調(diào)控。