番茄青枯雷爾氏菌致病力分化的全基因組分析

1、番茄青枯雷爾氏菌致病力分化的全基因組分析Whole-Genome Analysis of Tomato Ralstonia Solanacearum Pathogenicity Differentiation摘要青枯雷爾氏菌（Ralstonia solanacearum）是一種生長在土壤中的植物病原細菌，可在300多種單雙子葉植物中引起毀滅性的細菌性青枯病，對植物的危害極大。青枯雷爾氏菌菌株存在著致病力分化的現(xiàn)象，無致病力菌株可侵染寄主植物卻不引起寄主發(fā)病。目前國際上，尚未見無致病力青枯雷爾氏菌基因組測序的相關(guān)報道，通過本項目研究，獲得青枯雷爾氏菌無致病力菌株和強致病力菌株的基因組，為青枯雷爾

2、氏菌致病力分化研究奠定基礎(chǔ)。本文分別對青枯雷爾氏菌強致病力菌株FJAT-91和無致病力菌株FJAT-1458的全基因組進行Solexa測序，并且利用Velvet和ABySS進行de novo拼接，獲得368個FJAT-91框架序列（Scaffold），707個FJAT-1458框架序列；FJAT-91和FJAT-1458的基因組規(guī)模分別為5.6 Mb和5.7 Mb；GC含量均為66%左右。數(shù)據(jù)量和覆蓋率，共檢測到SNP。通過Prodigal預(yù)測基因共獲得6556個FJAT-91基因，7060個FJAT-1458基因，其中部分基因未達到全長；基因密度均為79%左右。通過IMG和Swiss-pro

3、t對預(yù)測基因進行同源注釋。將新測序的FJAT-91和FJAT-1458，同已測序完成的青枯雷爾氏菌GMI1000、PSI07、CFBP2957、CMR15、MolK2和IPO1609進行比較基因組學(xué)分析。青枯雷爾氏菌泛基因組共有？個基因（unigene）組成；其中？個基因在8個菌株中都存在，構(gòu)成了青枯雷爾氏菌的核心基因組，占了泛基因組的？%；此外，在無致病力菌株FJAT-1458中共有？個特異獲得基因、？個特異丟失基因，在強致病力菌株FJAT-91中，共有？特異獲得基因，？個特異丟失基因。ANI進化分歧組成成分的差異關(guān)鍵詞：青枯雷爾氏菌，全基因組Solexa測序，de novo拼接，致病力分化

4、1 前言青枯雷爾氏菌（Ralstonia solanacearum）引起的細菌性青枯?。˙acterial wilt）是危害最大、分布最廣且造成損失最為嚴重的植物病害之一（Tsuchiya et al., 2005）。青枯雷爾氏菌主要分布在熱帶、亞熱帶和溫帶地區(qū)，環(huán)境適應(yīng)能力強；其寄主中有木本植物和草本植物，一年生植物和多年生植物，雙子葉植物和單子葉植物（Janse, 2005），對茄科植物如番茄、辣椒、茄子、煙草、馬鈴薯的危害尤為嚴重，成為生產(chǎn)的主要障礙之一（Yang et al., 2006），該菌亦可危害花生、香蕉、甘薯、生姜、桉樹等植物（Sood, 2005），寄主范圍達到44個科的3

5、00多種植物（Alvarez, 2005）。因其侵染性強，寄主范圍大，地域分布廣，在水中和土壤環(huán)境中的生存能力強，青枯雷爾氏菌被認為是世界上危害最大植物病原菌之一（Genin, 2010），美國在“911”之后一度將其列在十大農(nóng)業(yè)恐怖微生物之首，多年來一直是相關(guān)學(xué)科的研究熱點。青枯雷爾氏菌表現(xiàn)出高度的種下分化的多態(tài)性，可按生理小種（Carmeille et al., 2006）、生化型（Stefani et al., 2005）、血清型、基因型、致病型等指標(biāo)進行種下分類（Thurston, 2001）。青枯雷爾氏菌的致病力分化極為顯著的，自然環(huán)境中就存在強致病力、無致病力以及介于其中的過渡菌株

6、。強致病力菌株侵入寄主植物并導(dǎo)致寄主發(fā)病，無致病力菌株可侵染寄主植物卻不引起寄主發(fā)?。˙alabel et al., 2005）?，F(xiàn)有研究表明，無致病力菌株無法自然恢復(fù)為強致病力菌株（Robertson et al., 2004），而強致病力菌株卻會自然喪失致病性，弱化形成無致病力菌株，致病力喪失的原因可能是營養(yǎng)條件和極端氣候引起的基因突變（Kanda et al., 2003）。青枯雷爾氏菌存在的復(fù)雜多態(tài)性是其內(nèi)部基因組的復(fù)雜性的體現(xiàn)，青枯雷爾氏菌能自然發(fā)生基因重組（ HYPERLINK /pubmed?term= Salanoubat et al., 2002），也存在從環(huán)境中獲得遺傳物質(zhì)

7、的水平基因轉(zhuǎn)移（Guidot et al., 2009），模式菌株GMI1000中16%的基因可能來自于水平轉(zhuǎn)移（Nakamura et al., 2004），這些都是青枯雷爾氏菌的菌種復(fù)雜性的原因（Coupat et al., 2008）。青枯雷爾氏菌的種內(nèi)差異甚至讓人懷疑這些菌株是否還應(yīng)該屬于一個物種。2004年，強致病力的標(biāo)準(zhǔn)菌株GMI1000的基因組全序列被測序并發(fā)表（Nakamura et al., 2004），迄今共有6個公開的強致病力菌株基因組序列（ HYPERLINK /pubmed?term= Remenant et al., 2010），其中包含4個完成圖，2個草圖。通過對

8、這些基因組的分析表明，青枯雷爾氏菌基因組呈現(xiàn)鑲嵌結(jié)構(gòu)（Genin, 2004），一些持家基因及毒力基因構(gòu)成大致的骨干，但不同的生態(tài)表型或寄主特異性則可能具有不同的基因組模塊（Poueymiro et al., 2009）。然而，目前尚無一個完整的無致病力菌株基因組完成圖數(shù)據(jù)，數(shù)據(jù)的缺乏使得青枯雷爾氏菌復(fù)雜的致病力分化機制難以分析。2 材料與方法2.1 菌株的獲得通過常規(guī)組織分離方法，在番茄青枯病發(fā)生田塊，通過大量采集與分離獲得大量青枯雷爾氏菌菌株，結(jié)合TTC培養(yǎng)基鑒定和致病性檢測，獲得親緣性關(guān)系極高的番茄無致病力菌株FJAT-1458和強致病力菌株FJAT-91。青枯雷爾氏菌無致病力菌株和強致

9、病力菌株的菌落形態(tài)以及致病力分化的檢測結(jié)果見圖a。通過青枯雷爾氏菌其菌落形態(tài)來判定其致病性強弱，其中強致病性的青枯雷爾氏菌菌落形態(tài)表現(xiàn)為：流動性強，中間為粉紅色，并且白邊比較寬，表面比較濕潤。弱致病性或無致病性的青枯雷爾氏菌菌落形態(tài)表現(xiàn)為：無流動性，中間為暗紅色，并且白邊比較窄，表面比較干燥。在相同接菌量的情況下，強致病性青枯雷爾氏菌采取剪葉方式在第4天就開始發(fā)?。徊⑶以诘?天已經(jīng)完全發(fā)?。ㄋ劳雎蔬_到100%）。而弱致病性青枯雷爾氏菌在這段期間剪葉均不發(fā)病。2.2 高通量測序本文采用Solexa-illumina新一代測序技術(shù)分別對青枯雷爾氏菌FJAT-91和FJAT-1458的全基因組進行s

10、hotgun測序，其中前者是強致病力菌株，而后者是無致病力菌株，測序平臺為GA IIx。Solexa測序讀段（short reads）的讀長為54 bp，分別構(gòu)建了兩種雙末端測序文庫（Library），插入序列長度（Insert length）分別為300 bp和2500 bp。其中300 bp的插入文庫采用Paired-End測序，兩端讀段的方向相向（inward）；2500 bp文庫采用Mate-Pair測序，兩端讀段的方向相反（outward）。圖1. 生物信息學(xué)分析的技術(shù)路線adebfcadebfc圖2. 青枯雷爾氏菌致病力分化注：青枯雷爾氏菌強致病力菌株的菌落形態(tài)（a：普通相機；b：

11、體視顯微鏡），青枯雷爾氏菌無致病力菌株的菌落形態(tài)（d：普通相機；e：體視顯微鏡），致病力分化的檢測結(jié)果（c：強致病力菌株；f：無致病力菌株）2.3 de novo拼接2.3.1 拼接前的質(zhì)量控制在對高通量讀段進行拼接之前，需要進行一系列質(zhì)量控制，包括：格式轉(zhuǎn)換，質(zhì)量評估，質(zhì)量修剪和長度過濾。Solexa高通量測序生成原始數(shù)據(jù)為一套qseq文件，需要根據(jù)標(biāo)簽信息轉(zhuǎn)換為fastq格式，并根據(jù)illumina的自動篩選結(jié)果輸出高質(zhì)量的讀段。qseq轉(zhuǎn)換為fastq文件的perl腳本采用了多線程編程技術(shù)，可根據(jù)標(biāo)簽信息，并行式地從各個qseq文件中提取目標(biāo)樣本的讀段，并統(tǒng)計目標(biāo)樣本測序的總數(shù)據(jù)量。讀段

12、的質(zhì)量值采用了illumina 1.3+標(biāo)準(zhǔn)，即PHRED值0至62分別由第64位ASCII碼字符（）至第126位字符（）表示。本文首先利用SolexaQA（Cox et al., 2010）對全部讀段的測序質(zhì)量進行評估，評測結(jié)果可為質(zhì)量修剪提供依據(jù)。本文在質(zhì)量修剪過程中，首先統(tǒng)一剪除3 末端最后一個堿基，使所有讀段只為53 bp；繼而采用Modified-Mott trimming算法（為CLC Genomics Workbench所采用質(zhì)量修剪算法）進行Perl編程，對所有讀段進行質(zhì)量修剪和過濾，過程如下：首先通過公式2，將PHRED質(zhì)量值換算成單堿基測序錯誤的概率（Pe），接著對每個堿基

13、都計算了一個新值Val，令Val=Limit-Pe，本文將Limit設(shè)為0.05。再設(shè)置一個新值sum，初始為0；從讀段的末端開始掃描序列，每掃描到一個堿基，就把堿基對應(yīng)的Val值累加，賦予sum，即sum=sum+Val；若碰到sum0的堿基開始到出現(xiàn)sum最大值的堿基結(jié)束，這段序列就是質(zhì)量修剪后（trimmed）保留下來的高質(zhì)量（High quality）讀段，而低質(zhì)量的序列末端已被切除。若掃描整條序列后，sum仍為0，則舍棄整個讀段。經(jīng)過質(zhì)量修剪之后，對含有超過2個N（Base Calling模糊不清的堿基）的讀段進行修剪，只取最長的一段少于等于2個N的讀段。最后根據(jù)長度閥值進行過濾，本

14、文只保留讀段長度大于閥值35 bp的讀段。若PE配對的兩條讀段中僅其中一條被過濾，則將另一條被保留的讀段另存至一個SE（Single End）文件中，使之不影響PE讀段的對齊。公式1公式22.3.2 ABySS和Velvet拼接在拼接之前，還需要進行以下處理：Mate-Pair文庫（MP Library），即本文的2.5 kb文庫，雙末端配對讀段的方向為outward，需要將forward和reverse讀段都進行反向互補，使兩讀段方向轉(zhuǎn)為inward。velvet和ABySS都只識別inward的PE讀段，轉(zhuǎn)換之后會使拼接質(zhì)量提高。本文先分別利用ABySS（Simpson et al., 2

15、009）和Velvet（Zerbino , Birney, 2008）單獨各自對青枯雷爾氏菌FJAT-91和FJAT-1458進行拼接。對Kmer值從21至上限35（讀段長度過濾的閥值）之間的每個奇數(shù)都進行調(diào)試，確認達到最佳拼接效果時的Kmer值為35 bp；另外，Velvet的參數(shù)min_pair_count和ABySS的參數(shù)n意義是一樣的，都設(shè)為20，表示構(gòu)成一個框架序列Scaffold至少需要20對讀段；ABySS不需要再設(shè)其他參數(shù)，Velvet需要參數(shù)exp_cov，表示期望的覆蓋度（Coverage），可設(shè)為auto，還需要指定插入文庫的長度及標(biāo)準(zhǔn)差SD，SD的設(shè)定只需要取一個盡量大

16、的值，不需要調(diào)試，本文設(shè)為20%，即300 bp的SD為60 bp，2500 bp的SD為500 bp。本文進一步將ABySS第一次拼接（Kmer=35）得到的框架序列（Scaffold）拆分為重疊群（Contig），Contig長度長于Velvet長序列拼接的長度限制（32700 bp）時，將ABySS的Contig拆分為三分之二互有重疊的兩段序列。最后將ABySS的Contig序列作為長序列（long sequence）和所有短序列（short read）一起通過Velvet進行二次拼接，Velvet的其他參數(shù)不變（Kmer=35）。2.3.3 讀段定位將讀段和拼接后的框架圖（Scaffo

17、ld）進行比對，讀段的定位（Mapping）信息可用來計算讀段的覆蓋率（Coverage）、SNP頻率和測序錯誤率（Error）。本文利用大規(guī)模高通量讀段快速聯(lián)配作圖工具Bowtie（Langmead et al., 2009）分別將青枯雷爾氏菌FJAT-91和FJAT-1458質(zhì)量修剪后的讀段，與Velvet-ABySS拼接得到的框架圖進行比對。首先，選擇序列長度大于300 bp的Scaffold，并按長度排序以及重命令，如RS91va_0001。然后，通過bowtie-build建立框架圖的ebwt索引。最后，選擇Bowtie比對模式為“-v 2”，即半全局聯(lián)配模式（全局的讀段-vs-局部

18、的染色體），允許有2個錯配，bowtie尚不支持Gapped比對，比對結(jié)果中不會含有InDel；其他參數(shù)有-phred64-quals，即本文的讀段質(zhì)量標(biāo)準(zhǔn)為illumina 1.3+；比對結(jié)果保存為SAM格式（參數(shù)-S）；此外，采用Paired-End模式（參數(shù)-q -1 m1 -2 m2），不能達到理想的比對效果，因而本文未采用該模式。2.3.4 計算覆蓋率（深度和廣度）Bowtie的比對結(jié)果保存為SAM格式，可利用Perl語言讀取SAM文件計算定位的讀段（Mapped read）總數(shù)量和堿基總數(shù)。利用SAMTools（Li et al., 2009）處理SAM文件（文本文件），過濾未定位

19、的讀段，并轉(zhuǎn)換為BAM格式（壓縮的二進制SAM文件），對BAM文件排序并建立索引，最后將多個基于同一個框架圖的BAM文件合而為一。本文利用BEDTools（Quinlan , Hall, 2010）的genomeCoverageBed程序計算per-base coverage，即單堿基覆蓋率，需要參數(shù)“-d”。根據(jù)per base coverage通過Perl編程計算所有定位讀段對整個框架圖（含Gap）的總覆蓋率，包括測序深度（Depth）和廣度（Breadth）。測序深度為所有定位讀段覆蓋框架圖的倍數(shù)，廣度為所有定位讀段覆蓋框架圖的完整程度，在本文研究中，框架圖中未覆蓋讀段的部分即為Gap區(qū)

20、域。通過掃描per base coverage，計算各個框架序列（Scaffold）的覆蓋度（深度Depth和廣度Breadth），同時計算各個重疊群（Contig，不含Gap）的測序深度?？梢愿鶕?jù)覆蓋率的高低評估框架圖的拼接質(zhì)量，濾去覆蓋率過低（本文設(shè)為30）的Scaffold；而覆蓋率數(shù)倍于平均覆蓋率的序列，則可能是重復(fù)序列。2.3.5 檢測SNP和測序錯誤讀段的Bowtie比對結(jié)果已轉(zhuǎn)換為SAMTools處理之后的BAM格式，利用SAMTools的mpileup生成pileup文件，在此文件中包括了每個匹配位點的定位詳情，pileup格式說明詳見SAMTools網(wǎng)站上的說明：HYPERL

21、INK /pileup.shtml/pileup.shtml通過perl語言編程來處理pileup文件。首先從pileup文件中列出所有差異位點（讀段與參照序列之間有錯配）。其次，對這些差異位點進行SNP過濾，SNP或者Error的區(qū)別只是頻率不同。某個位點發(fā)生多次測序錯誤的概率會隨著錯配次數(shù)的增多而急劇降低，錯配次數(shù)越多，說明該位點存在SNP的概率越高。比如說，參照序列為堿基“A”，讀段錯配堿基為“T”和“G”，A至T的次數(shù)和頻率都很高，則有存在SNP的可能性就高，而A至G的次數(shù)僅一次，則完全不能排除是測序錯誤的可能性。因此，進行SNP過濾需要設(shè)定一系列閥值，本文設(shè)定的閥值如下：差異位點的p

22、er base coverage在30以上，500以下，錯配次數(shù)（絕對數(shù)量）在5次以上（包括5次），SNP頻率大于等于2%。最后，對不符合SNP篩選條件的錯配當(dāng)作測序錯誤，錯配數(shù)進行累積，從而計算測序錯誤率。測序的錯誤率為所有測序錯誤堿基數(shù)除以所有定位讀段的總堿基數(shù)。另外，框架圖中少數(shù)幾個位點需要進行SNP校正。本文進行讀段定位的參照序列是新拼接的框架圖序列，來自于這些定位讀段，在計算差異位點的錯配比率時，部分SNP頻率高于參照堿基的頻率，并且在所有SNP中頻率最高，將該SNP和參照堿基交換，進行校正。對框架圖進行SNP校正之后，再計算SNP頻率和錯誤率。2.3.6 框架圖拼接后的質(zhì)控統(tǒng)計青枯

23、雷爾氏菌FJAT-91和FJAT-1458基因組拼接之后，需要對其拼接結(jié)果（框架圖）進行評估，從而判斷結(jié)果是否理想。首先，通過框架序列（Scaffold）的長度閥值和最低覆蓋率對框架圖進行質(zhì)量控制，本文設(shè)定的長度閥值為300 bp，最低覆蓋率為30。然后，對質(zhì)控之后的框架圖進行評估，計算以下幾個數(shù)值：Scaffold的數(shù)量，長度總和（基因組規(guī)模），Scaffold長度的最大值和平均值，Gap的比例，Scaffold的N50、N90序列數(shù)量和N50、N90長度值，Contig的數(shù)量，Contig的長度總和，Contig長度的最大值和平均值，Contig的N50、N90序列數(shù)量和N50、N90長度

24、值。其中N50、N90值的概念如下：將各個序列按長度大小排序，從大至小逐一掃描各個序列的長度值，進行累積，當(dāng)該累積值第一次超過所有序列總長的50%時，此時掃描到的序列，其長度值即為N50值，而此時已掃描的序列數(shù)量即為N50序列數(shù)量。N90值亦同理。2.4 基因預(yù)測與注釋本文利用原核生物基因預(yù)測軟件Glimmer3、GeneMarkS和Prodigal分別在青枯雷爾氏菌FJAT-91和FJAT-1458基因組框架圖的重疊群序列（Contig）上預(yù)測編碼基因序列。因本文測得的框架圖尚不完整，含有較多的Gap，基因預(yù)測需包括未達全長的基因片段（Fragment）?；蚬δ茏⑨?，將FJAT-91和FJ

25、AT-1458的基因分別與Swiss-Prot和trEMBL數(shù)據(jù)庫進行BLASTP同源比對，選擇最相似的比對的結(jié)果，并且要求序列比對相似性在30%以上，比對雙方中較短的那條序列的全長，聯(lián)配長度需占其60%以上，期望值E-value在0.001以上。通過Swiss-prot的GO關(guān)聯(lián)（GOA，GO Association），分別對FJAT-91和FJAT-1458中具有GO注釋的基因利用WEGO進行GO功能分類。蛋白質(zhì)結(jié)構(gòu)域分析，利用Hmmer3的hmmscan程序?qū)JAT-91和FJAT-1458所有預(yù)測基因的蛋白質(zhì)序列與本地化PFAM數(shù)據(jù)庫進行比對，從而鑒定這些基因的蛋白質(zhì)結(jié)構(gòu)域。非編碼R

26、NA的預(yù)測包括核糖體RNA（rRNA）和轉(zhuǎn)運RNA（tRNA）。rRNA可通過青枯雷爾氏菌GMI1000中rRNA序列進行同源預(yù)測，tRNA利用rRNAScan-SE和aragorn進行預(yù)測。重復(fù)序列注釋，通過RepeatModeler（調(diào)用RECON和RepeatScout）分析FJAT-91和FJAT-1458的框架圖序列，發(fā)現(xiàn)重復(fù)序列，構(gòu)建物種特異的重復(fù)序列數(shù)據(jù)庫，聯(lián)同RepBase，利用RepeatMasker檢測框架圖中的重復(fù)序列，通過和已有的TE數(shù)據(jù)庫進行比對，對重復(fù)序列進行注釋。2.5 比較基因組學(xué)研究青枯雷爾氏菌CMR15、PSI07、CFBP2957、GMI1000和Po82

27、全基因組完成圖（Completed Genome），青枯雷爾氏菌IPO1609、UW551和MolK2基因組草圖（Genome Draft），以及青枯雷爾氏菌近緣物種Ralstonia pickettii 12J全基因組完成圖，均可從NCBI或者Genoscope RalstoniaScope網(wǎng)站下載數(shù)據(jù)，后者的網(wǎng)址如下： HYPERLINK s.fr/agc/microscope/about/collabprojects.php?P_id=67 s.fr/agc/microscope/about/collabprojects.php?P_id=67將新測序的青枯雷爾氏菌強致病力菌株FJAT-

28、91和無致病力菌株FJAT-1458聯(lián)同以上9個菌株的基因組一起進行比較基因組學(xué)研究。2.5.1 雙向最佳比對在全基因組水平上，兩個物種之間的蛋白質(zhì)序列進行雙向最佳比對（BBH，Bidirectional Best Hit），可以鑒定得到這兩個物種之間的直系同源基因（Orthologous Genes）。利用BLAST+程序包中的BLASTP程序，將以上每個菌株全基因組中的所有蛋白質(zhì)序列分別和其他10個菌株進行BLASTP，只選擇最佳匹配項，并要求序列比對相似性在30%以上，比對雙方中較短的那條序列的全長，聯(lián)配長度需占其60%以上，期望值E-value在0.001以上。通過Perl編程，各個菌

29、株兩兩之間的最佳比對結(jié)果進行聚類2.5.2 全基因組比對（megablast）和基因組島鑒定2.5.3 計算AAI并構(gòu)建物種進化樹通過蛋白質(zhì)序列雙向最佳比對（BBH）鑒定各個物種兩兩之間的直系同源基因（Orthologous Genes），從而計算這些直系同源基因的平均氨基酸一致性（AAI，Average Amino acid Identity）。各菌株之間的分歧程度通過100%減去AAI值來表示，以此生成距離矩陣，并保存為Phylip所能讀取的格式。之后，利用phylip的neighbor程序，采用鄰接法（NB）構(gòu)建物種進化樹。3 結(jié)果與分析3.1 框架圖拼接結(jié)果分析3.1.1 數(shù)據(jù)量和測序

30、深度經(jīng)過Solexa高通量測序，青枯雷爾氏菌強致病力菌株FJAT-91共產(chǎn)生了1.3 Gb數(shù)據(jù)量（測得堿基總數(shù)），無致病力菌株FJAT-1458共產(chǎn)生了1.1 Gb數(shù)據(jù)量。經(jīng)過illumina自動篩選之后，F(xiàn)JAT-91的有效讀數(shù)，即高質(zhì)量讀段的數(shù)量，為2109萬個，堿基總數(shù)為1.1 Gb，占了原始數(shù)據(jù)的87.6%，測序深度為203（已知青枯雷爾氏菌基因組平均規(guī)模為5.6 Mb）；而FJAT-1458的有效讀數(shù)為1751萬個，堿基總數(shù)為0.9 Gb，是原始數(shù)據(jù)的87.5%（詳見表1），測序深度為169。如表1所示，讀段經(jīng)過質(zhì)量修剪之后，青枯雷爾氏菌FJAT-91和FJAT-1458實際用于拼接

31、的測序深度分別為156和131，各自修剪了23.2%和22.5%的原始序列，讀段的平均長度分別為49.8 bp和50.0 bp，原始讀段長度為54 bp。在de novo拼接之后，將讀段定位至框架圖上，未能聯(lián)配至框架圖上的讀段，可以認為這些讀段沒有在拼接中用到。讀段定位結(jié)果顯示，實際用于拼接的讀段覆蓋率深度在FJAT-91和FJAT-1458中分別有129和103，分別下降了17%和21%，用于拼接的讀段分別僅占原始數(shù)據(jù)的56%和53%。表1. 青枯雷爾氏菌FJAT-91和FJAT-1458讀段數(shù)量和堿基總量FJAT-91（virulent）FJAT-1458（avirulent）Raw da

32、taread length54 bp54 bpreads number2407275020001978total bases1299928500 bp1080106812 bpgenome sizeavg. 5.6 Mbavg. 5.6 Mbdepth 1)232193Illumina filteredread length54 bp54 bpreads number2109410617507820total bases1139081724 bp945422280 bpgenome sizeavg. 5.6 Mbavg. 5.6 Mbdepth203169Quality trimmedread

33、 length3553 bp, mean=49.8 bp3553 bp, mean=50 bpreads number1750103814740010total bases871237877 bp736276730 bpgenome sizeavg. 5.6 Mbavg. 5.6 Mbdepth156131Mapping to scaffoldread length3553 bp, mean=49.9 bp3553 bp, mean=49.9 bpreads number1454335111935129total bases725271690 bp595253263 bpscaffold si

34、ze 2)5620164 bp5734490 bpdepth 3)129103breadth 4)88.36%88.74%注：1) 青枯雷爾氏菌基因組規(guī)模平均值約為5.6 Mb，測序深度為所有讀段可覆蓋基因組的倍數(shù)，測序深度=堿基總數(shù) / 5.6 Mb；2) 大于300 bp的所有框架序列的總和，包括Gap；3) 測序深度=堿基總數(shù) / Scaffold總長。4) 廣度Breadth即所有讀段覆蓋參照基因組序列的程度，在本文中，未覆蓋的區(qū)域即Gap部分，測得部分的重疊群序列，實際的測序深度=depth / breadth。圖3. 高通量測序的覆蓋深度3.1.2 de novo拼接基因組框架圖經(jīng)

35、過調(diào)試，最佳Kmer值設(shè)為35 bp，分別采用了三種拼接方法：單獨利用ABySS進行拼接，單獨利用Velvet進行拼接，結(jié)合Velvet和ABySS進行拼接，在同樣的參數(shù)設(shè)置下，Scaffold長度在300 bp以上的拼接結(jié)果，三者的拼接效果表2所示。表2. 拼接結(jié)果的質(zhì)控統(tǒng)計FJAT-91FJAT-1458VelvetABySSV-A 1)VelvetABySSV-Aall scf 2) number126161641022 4)231361641804total length578398761659125700755615908861659125878500scaffold length

36、= 300 bpscf number590268336811222894707scf length569820757751615620164600038655203115734490max scf length5884254107460309728900831134329864avg. scf length9657.982152.515272.185347.941907.58111.02N rate (Gap)14.23%4.49%11.61%16.82%0.44%11.22%scf N50 number18273105242736scf N50 length59152505611205524

37、052304430933scf N90 number1621558784601792241scf N90 length36888041032719867823619ctg 3) number377729522655457931702978ctg length488738055160294967576499094154960965091150ctg max length168583496634691170923113426925ctg avg. length1293.981868.571871.031089.961733.781709.59ctg N50 number64831230278045

38、3431ctg N50 length209540283842182729463001ctg N90 number232617351502279918831681ctg N90 length628730823509741809注：1) V-A = Velvet-ABySS; 2) scf = Scaffold, 3) ctg = Contig; 4) 加粗的數(shù)字表示最佳Velvet-ABySS拼接的總體特征顯示，de novo拼接共獲得1022個FJAT-91的Scaffold，長度在300 bp以上的Scaffold有368個，所有Scaffold總長達到5.6 Mb，其中Gap（堿基未測得，

39、標(biāo)記為“N”）序列的比例為11.6%，所有Contig序列總長為4.97 Mb，包括2655個Contig；在FJAT-1458中，共有1804個Scaffold，長度在300 bp以上的Scaffold有707個，Scaffold總長為5.7 Mb，Gap序列比例為11.2%，Contig序列總長為5.09 Mb，包括2978個Contig。Scaffold和Contig的N50和N90統(tǒng)計結(jié)果（見表2），在FJAT-91中，Scaffold的N50序列數(shù)量為10，N50序列長度值為112 kb，Contig的N50序列數(shù)量為302，N50序列長度值為3.8 kb；在FJAT-1458中，S

40、caffold的N50序列數(shù)量為36，N50值為31 kb，Contig的N50數(shù)量為431，N50值為3 kb。N90的統(tǒng)計結(jié)果詳見表2。此外，Velvet拼接獲得的Scaffold，數(shù)量和長度均好于ABySS的拼接結(jié)果；而ABySS在Contig拼接上的效果要好于Velvet，Gap比例較低；從表2的結(jié)果可見，Velvet-ABySS的Scaffold拼接效果超過了Velvet，Contig拼接效果也與ABySS較為接近，達到了Velvet和ABySS優(yōu)缺點互補的拼接效果。3.1.3 讀段定位結(jié)果分析讀段分別定位至FJAT-91和FJAT-1458的基因組框架圖上，總體統(tǒng)計結(jié)果見表1。其中

41、各個Contig的覆蓋深度分布如圖4所示，在FJAT-91（左）和FJAT-1458（右）中，多數(shù)Contig的覆蓋深度在50至200以內(nèi)，F(xiàn)JAT-91中各個Contig的覆蓋深度比FJAT-1458略高。從圖4亦可知，F(xiàn)JAT-91和FJAT-1458的GC含量基本上在60%至75%之間。此外，圖4中FJAT-91和FJAT-1458的散點分布集中，并無明顯分區(qū)現(xiàn)象，表明測序中基本沒有污染。圖4. 覆蓋深度和GC含量散點圖（左：FJAT-91，右：FJAT-1458）不滿足SNP過濾條件的錯配屬于測序錯誤，從而評估Solexa測序的錯誤率。在FJAT-91中共有1413371個單堿基測序錯

42、誤，測序錯誤率為0.16%；在FJAT-1458中共有1548491個單堿基測序錯誤，測序錯誤率為0.21%。在青枯雷爾氏菌FJAT-91中共檢測到22,144個SNP，其基因組框架圖大小約為5.6 Mb，平均全基因組的SNP密度為3.96/kb；在FJAT-1458中檢測到19,017個SNP，SNP密度為3.34/kb。本文未檢測單核苷酸的插入和缺失。SNP在全基因組的平均發(fā)生頻率，在FJAT-91中為10.32%，在FJAT-1458中為10.68%，頻率相近。青枯雷爾氏菌FJAT-91和FJAT-1458的SNP中，屬于轉(zhuǎn)換（C/T、G/A）的SNP分別只有1700和1400多個，所占

43、比例較少，而顛換四種類型C/G、A/T、C/A和G/T中，C/A和G/T兩者總共約占了全部SNP的70%，F(xiàn)JAT-91整體的轉(zhuǎn)換與顛換比為15.81，F(xiàn)JAT-1458的轉(zhuǎn)換與顛換比為15.64。但是，轉(zhuǎn)換的SNP發(fā)生頻率分別為18%和16%，大于顛換的SNP突變頻率。表3. 轉(zhuǎn)換和顛換的SNP分布數(shù)量和頻率TransitionsTransversionAllC/TG/AC/GA/TC/AT/GFJAT-91SNP Number1707154329052208320745122144Frequency18.28%18.66%11.37%10.94%7.76%7.91%10.32%FJAT-1

44、458SNP Number1437142426102856667659419017Frequency16.30%16.08%11.66%8.21%8.95%8.81%10.68%3.2 基因組特征概述青枯雷爾氏菌CMR15、PSI07、CFBP2957、GMI1000和Po82已獲得全基因完成圖，菌株IPO1609、UW551和MolK2只得到基因組草圖，共有八個青枯雷爾氏菌基因組已知；本文新測序了兩個基因組FJAT-91和FJAT-1458，前者為強致病力菌株，后者為無致病力菌株；分析這些菌株的基因組特征，包括基因組規(guī)模，GC含量，編碼基因數(shù)量，基因密度，基因序列平均長度，rRNA及tRNA

45、數(shù)量，分析結(jié)果詳見表4。表4. 基因組概述StrainPhylotypeGenomesize (bp)GCcontentCDSnumberGenedensityCDS avg.length (bp)rRNA operonstRNAnumberCMR15III5,593,04166.8%514986.7%904.2359PSI07IV5,606,28866.3%524786.8%808.6149CFBP2957IIA5,683,40266.4%531085.2%972.9156GMI1000I5,810,92266.9%563589.2%955.5457IPO1609IIB5,523,29266

46、.7%520380.7%990.8129UW551IIB5,952,91864.2%530181.3%913.31Molk2IIB5,862,10166.7%543883.8%903.3134Po82IIB5,430,26366.6%508091.3%975.63FJAT-91I5,620,16466.0%655687.6%656.0148FJAT-1458I5,734,49066.0%706086.3%622.4254Mean5,681,688 66.3% 5295 1)85.6%928.0 2)1.8注：1) CDS數(shù)量的均值不包括FJAT-91和FJAT-1458; 2) CDS平均長度

47、不包括FJAT-91和FJAT-1458.圖5. 十個青枯雷爾氏菌基因組中的CDS數(shù)量FJAT-91FJAT-91圖6. FJAT-91基因組注釋注：最外層黑白條碼表示10 kb掃描窗口的基因密度，白色表示基因密度在60%以下，顏色越深表示基因密度越高；中間層表示所有預(yù)測基因的位置，其中紅色表示基因不完整（Partial，F(xiàn)ragment），即有缺少起始密碼子或終止密碼子，藍色表示基因達到全長（Complete，F(xiàn)ull Length）；最內(nèi)層，將具有基因名（Symbol Name）的基因按其位置標(biāo)識在圖中。新測序的青枯雷爾氏菌FJAT-91和FJAT-1458，與GMI1000同屬于進化I型

48、；框架圖（Scaffold序列長度大于300 bp）的堿基總數(shù)分別為5.6 Mb和5.7 Mb，與青枯雷爾氏菌基因組平均規(guī)模較為接近。統(tǒng)計重疊群序列（不含Gap）的GC含量，青枯雷爾氏菌FJAT-91和FJAT-1458的GC含量均為66%。10個青枯雷爾氏菌的GC含量普遍在一個較高的水平，均值為66.3%。已知的8個青枯雷爾氏菌各有5000多個基因，其中GMI1000最多，具有5635個基因，Po82最少，具有5080個基因，均值為5295個。新測序的FJAT-91和FJAT-1458分別預(yù)測得到6556和7060個基因，其中有較多的基因片段（Fragment），比例超過50%（如圖5所示）

49、，這些基因不完整（Partial），缺少起始密碼子或終止密碼子，未達全長。FJAT-91中全長基因共有3010個，基因片段共有3546個；FJAT-1458中全長基因共有3075個，基因片段共有3985個。因為含有較多基因片段，F(xiàn)JAT-91和FJAT-1458的CDS平均長度分別僅為656 bp和622 bp，而其他8個菌株的CDS平均長度的均值為928 bp，相差較多；FJAT-91全長基因的CDS平均長度為801 bp，F(xiàn)JAT-1458全長基因的CDS平均長度為766 bp，也相差了100多bp，說明這些全長基因中也有一部分基因未達到全長。原核生物具有較高的基因密度，在青枯雷爾氏菌中，

50、基因密度一般可在80%以上，甚至超過90%。統(tǒng)計FJAT-91和FJAT-1458重疊群序列的基因密度，分別達到87.6%和86.3%，基因密度水平基本一致。經(jīng)過rRNA同源比對和tRNA預(yù)測，在FJAT-91中共發(fā)現(xiàn)1個rRNA，48個tRNA，在FJAT-1458中共發(fā)現(xiàn)2個rRNA，54個tRNA。FJAT-1458FJAT-1458圖7. FJAT-1458基因組注釋注：最外層黑白條碼表示10 kb掃描窗口的基因密度，白色表示基因密度在60%以下，顏色越深表示基因密度越高；中間層表示所有預(yù)測基因的位置，其中紅色表示基因不完整（Partial，F(xiàn)ragment），即有缺少起始密碼子或終止

51、密碼子，藍色表示基因達到全長（Complete，F(xiàn)ull Length）；最內(nèi)層，將具有基因名（Symbol Name）的基因按其位置標(biāo)識在圖中。3.3 青枯雷爾氏菌泛基因組Fig. 3 pan-genome, Core genome and Specific genes of eight R. solanacearum strains. The pan-genome of R. solanacearum and special genes of every strain were identified. 11888 unigenes constituted the pan-genome of R. solanacearum. The pan-genome was composed of the core-genome (genes presented in all strains), the dispensable genome (genes pre