文章題目：Utilizing PacBio Iso-Seq for Novel Transcript and Gene Discovery of Abiotic Stress Responses in Oryza sativa L.

中文題目：利用Pacbio Iso-Seq測(cè)序技術(shù)發(fā)現(xiàn)水稻非生物脅迫下的新轉(zhuǎn)錄本和新基因

發(fā)表期刊：International Journal of Molecular Sciences

發(fā)表時(shí)間：2020年10月31日

影響因子：4.556

研究背景

全球氣候變化導(dǎo)致高溫、干旱和夜間高溫等非生物脅迫條件的嚴(yán)重程度和頻率增加，這些都造成了作物產(chǎn)量的降低。隨著世界人口的增長(zhǎng)，植物育種專家面臨著開(kāi)發(fā)高產(chǎn)、優(yōu)質(zhì)、減少環(huán)境污染的新品種的艱巨任務(wù)。水稻是世界上一半以上人口的主要卡路里來(lái)源，特別是對(duì)亞洲最貧窮的人來(lái)說(shuō)。世界各地的基因庫(kù)中保存著23萬(wàn)多份水稻及其野生近緣種的廣泛自然遺傳多樣性的種質(zhì)資源，是一種無(wú)價(jià)的可用于作物改良的重要基因庫(kù)。

雖然世界上近80%的水稻種植是基于indica（秈稻亞種）品種，但目前的標(biāo)準(zhǔn)基因組及其注釋來(lái)自粳稻亞種Nipponbare。由于缺乏適當(dāng)?shù)幕蚪M，不同水稻亞種的研究大多都基于Nipponbare基因組。例如，在3000水稻基因組計(jì)劃中將測(cè)序序列比對(duì)到Nipponbare基因組上，丟棄了不能比對(duì)到該參考基因組的序列。這可能會(huì)導(dǎo)致非粳稻亞種特有的遺傳信息的丟失。另外，最近已經(jīng)對(duì)其它水稻亞種的栽培品種的基因組進(jìn)行了測(cè)序，例如indica（Shuhui498，Zhenshan 97, Minghui 63）、aus（Kasalath，N22），但其完整性和注釋程度仍存在差異。值得注意的是aus亞種是抗病、耐磷酸鹽缺失、耐澇、耐厭氧發(fā)育和抗旱等潛在性狀的寶貴基因來(lái)源。例如，在aus品種基因組中發(fā)現(xiàn)了耐磷酸鹽缺失相關(guān)基因OsPSTOL1、耐澇相關(guān)基因OsSNORKEL1/2和OsSUB1A。值得注意的是，這些基因在粳稻的Nipponbare亞種基因組序列中是不存在的。

在過(guò)去的幾年里，RNA測(cè)序（特別是基于illumina的短序列RNA-seq）已經(jīng)成為分析轉(zhuǎn)錄組的有力工具，用來(lái)識(shí)別在非脅迫控制和各種環(huán)境脅迫條件下差異表達(dá)的基因。然而，需要基于參考基因組或轉(zhuǎn)錄組序列對(duì)RNA-seq數(shù)據(jù)進(jìn)行比對(duì)和注釋來(lái)確定轉(zhuǎn)錄水平。在水稻中，參考基因組決定了可以鑒別的差異表達(dá)基因和轉(zhuǎn)錄本亞型。顯然，參考基因組/轉(zhuǎn)錄組中沒(méi)有的基因的表達(dá)信息在分析過(guò)程中會(huì)丟失。這在研究耐脅迫的外來(lái)品種、陸地品種或野生稻種時(shí)尤其相關(guān)，因?yàn)樗鼈兛赡芎袇⒖计贩NNipponbare不存在的耐受基因。這將嚴(yán)重限制識(shí)別支持作物改良計(jì)劃的新候選基因的可能性。
解決這一問(wèn)題的一個(gè)顯而易見(jiàn)的辦法是對(duì)所需的基因組進(jìn)行測(cè)序、組裝和注釋。但是，這種方法比較昂貴和耗時(shí)。在這篇文章里，我們探索了一種更有針對(duì)性的RNA-Seq序列方法來(lái)測(cè)序和重建了三個(gè)不同亞種的水稻品種的部分轉(zhuǎn)錄本作為參考，Pacific Bioscience（PacBio）屬于提供高通量全長(zhǎng)轉(zhuǎn)錄本序列的新一代測(cè)序方法。該方法已成功應(yīng)用于對(duì)現(xiàn)有植物轉(zhuǎn)錄本及注釋的探索和擴(kuò)展，如高粱、小麥、甘蔗、野生棉花、不同的穗型草、苜蓿等。

取樣材料

樣品取自三個(gè)水稻品種的10個(gè)不同亞種的不同組織部位（表1）：

分析結(jié)果

1.重構(gòu)轉(zhuǎn)錄本

使用PacBio Sequel I測(cè)序平臺(tái)對(duì)每個(gè)品種進(jìn)行SMRT測(cè)序，得到15.49~24.51GB的轉(zhuǎn)錄本數(shù)據(jù)。用IsoSeq3軟件對(duì)原始測(cè)序數(shù)據(jù)進(jìn)行ccs和lima處理，每個(gè)品種SMRT cell分別得到460340~736747條全長(zhǎng)非嵌合序列（full-length non-chimeric reads簡(jiǎn)稱FLNC，包含3 ‘ 引物、5 ‘ 引物以及polyA尾）。全長(zhǎng)非嵌合經(jīng)過(guò)IsoSeq3聚類和polish分別得到37951~54684高質(zhì)量轉(zhuǎn)錄本（HQ）以及1233 ~2170低質(zhì)量轉(zhuǎn)錄本（LQ）。先將HQ與NCBI核苷酸數(shù)據(jù)庫(kù)進(jìn)行blastn比對(duì)（E<=1e-10），再將上一步為比對(duì)上的轉(zhuǎn)錄本序列與NCBI蛋白數(shù)據(jù)庫(kù)進(jìn)行blastx比對(duì)（E<=1e-10），去除未比對(duì)上兩個(gè)數(shù)據(jù)庫(kù)的轉(zhuǎn)錄本序列，最終得到37535~54594條HQ用于后續(xù)分析（表2）。

Pacbio RSII平臺(tái)聲明使用RNA-seq二代測(cè)序數(shù)據(jù)對(duì)轉(zhuǎn)錄本數(shù)據(jù)進(jìn)行矯正，可以得到更多的HQ序列，因?yàn)長(zhǎng)Q序列中含有大量的插入和缺失。然而與RSII相比，PacBio Sequel I測(cè)序平臺(tái)的測(cè)序結(jié)果更好。為了驗(yàn)證這一結(jié)果，我們用minimap2將未矯正的HQ比對(duì)到相應(yīng)亞種的基因組，結(jié)果表明缺失的比例很小，所以進(jìn)一步的分析中沒(méi)有包含LQ序列。

2.轉(zhuǎn)錄本去冗余

在文庫(kù)準(zhǔn)備過(guò)程中，會(huì)產(chǎn)生5 ‘ RNA降解產(chǎn)物，并進(jìn)行測(cè)序。這些降解產(chǎn)物具有相同的外顯子結(jié)構(gòu)，但缺乏5 ‘序列信息，因此產(chǎn)生與技術(shù)偏差或生物學(xué)背景無(wú)關(guān)的冗余異構(gòu)體。我們測(cè)試了三種不同的去冗余方法，包括cogent、cDNA cupcake和TAMA，其中cDNA cupcake和TAMA需要基于參考基因組，而cogent不需要基于參考基因組。cogent基于pacbio全長(zhǎng)轉(zhuǎn)錄本序列重構(gòu)一個(gè)參考基因組，然后將相同的序列比對(duì)到重建的基因組，基于比對(duì)結(jié)果利用cDNA cupcake算法對(duì)轉(zhuǎn)錄本去冗余。cDNA cupcake和TAMA直接用minimap2軟件和各自的亞種參考基因組進(jìn)行比對(duì)。基于這三種方法，只有很少的轉(zhuǎn)錄本不能回比到基因組上（表3）?？偟膩?lái)說(shuō)，這三種去冗余方法均能顯著減少異構(gòu)體的數(shù)量，分別為47.6% （cDNA cupcake，Nipponbare）和68.3%（cogent，Dular）。

基于植物中430個(gè)高度保守的同源蛋白利用BUSCO軟件對(duì)TAMA算法去冗余前后的HQ進(jìn)行完整性評(píng)估（圖一），由于取樣不完全，缺失了54%~27%的重要蛋白，其中Nipponbare參考基因組（IRGSP）只缺失了6種。cDNA cupcake和TAMA的結(jié)果相似，而對(duì)于cogent，超過(guò)50%的蛋白缺失，最有可能的原因是轉(zhuǎn)錄本沒(méi)有回比到重建的基因組。

去冗余后轉(zhuǎn)錄本長(zhǎng)度中值都有所增長(zhǎng)，長(zhǎng)度分布和轉(zhuǎn)錄本長(zhǎng)度中值與Nipponbare參考基因組相似。統(tǒng)計(jì)了去冗余后10個(gè)品種每個(gè)基因相應(yīng)的轉(zhuǎn)錄本數(shù)量，其中基因只有一個(gè)轉(zhuǎn)錄本的比例，TAMA最高達(dá)到了75%，cDNA cupcake在60%左右，cogent只有50%。同時(shí)計(jì)算了Nipponbare參考基因組每個(gè)基因?qū)?yīng)的轉(zhuǎn)錄本數(shù)量進(jìn)行比較，該參考基因中基因只有一個(gè)轉(zhuǎn)錄本的比例達(dá)到了85%（圖2）。

來(lái)自同一亞種的不同品種親緣關(guān)系更近，我們使用系統(tǒng)發(fā)育樹(shù)評(píng)估亞種之間的遺傳距離。利用去冗余后的轉(zhuǎn)錄本序列基于Nipponbare參考基因組識(shí)別SNPs，使用SNPhylo繪制進(jìn)化樹(shù)（圖3）。SNPhylo提取高質(zhì)量并且具有代表性的SNPs進(jìn)行后續(xù)分析，cDNA cupcake算法大約30000個(gè)SNPs，cogent算法大約23200個(gè)SNPs，TAMA算法大約16000個(gè)SNPs。三種方法中，同一亞種的不同品種聚類在了一起，cDNA cupcake算法和cogent的聚類結(jié)果更相似。三種方法均能將aus與另外兩個(gè)亞種區(qū)分開(kāi)，但cogent和TAMA對(duì)indica和japonica種間的區(qū)分不如cDNA cupcake明顯。

3.評(píng)估重構(gòu)的轉(zhuǎn)錄本

基于TAMA算法得到的HQ進(jìn)行轉(zhuǎn)錄本的評(píng)估。由于TAMA只對(duì)比對(duì)到參考基因組上的轉(zhuǎn)錄本進(jìn)行去冗余，我們用cogent對(duì)沒(méi)比對(duì)上參考基因組的轉(zhuǎn)錄本進(jìn)行去冗余。合并結(jié)果后，每個(gè)品種最終得到10511（Dular）~15011（IR64）個(gè)基因，14255（Dular）~20803（Moroberekan）個(gè)轉(zhuǎn)錄本（表4）。與Nipponbare參考基因組相比，大約三分之一的基因位點(diǎn)和大約一半的轉(zhuǎn)錄模型被重建。每個(gè)品種每個(gè)基因的平均轉(zhuǎn)錄數(shù)約為1.4～1.5，略高于參考基因組的1.2。中位轉(zhuǎn)錄本長(zhǎng)度為986 bp（Dular）~1394 bp（Nipponbare），與Nipponbare參考值1385 bp相似。平均GC含量在50.87%（Dular）~52.76%（IR64），與Nipponbare參考值51.24%相似。利用gffcompare軟件與Nipponbare參考基因組進(jìn)行比較識(shí)別新基因與轉(zhuǎn)錄本。

4.功能注釋

為了深入了解重建轉(zhuǎn)錄本的生物學(xué)信息，我們進(jìn)行了功能注釋。使用TransDecoder軟件預(yù)測(cè)開(kāi)放閱讀框（ORFs），包括blast和PFAM，結(jié)果表明大約有60%~70%的完整ORFs（包括啟動(dòng)和終止密碼子）。此外還發(fā)現(xiàn)了26%~38%的5 ‘ ORF、很少比例的3 ‘ ORF和中間ORF（既沒(méi)有起始密碼子也沒(méi)有終止密碼子）（圖4）。

使用Trinotate和Mercator4進(jìn)行功能注釋。Mercator4是專門為植物開(kāi)發(fā)的，它使用了一種簡(jiǎn)單的層次樹(shù)結(jié)構(gòu)，被稱為“容器”，用來(lái)描述生物學(xué)概念。主要的生物過(guò)程如光合作用，都是由頂層的容器來(lái)表示的，每個(gè)子容器描述的是一個(gè)更詳細(xì)的子過(guò)程。目前本體包括27個(gè)功能類別，代表了植物中不同的生物過(guò)程。N22、IR64和Nipponbare三個(gè)品種作為各自亞種的代表與植物中所有水稻基因的分類進(jìn)行比較分析，結(jié)果顯示三個(gè)品種的注釋結(jié)果分布相似（圖5）。超過(guò)28000個(gè)水稻已知基因在Mercator庫(kù)中沒(méi)有注釋分類信息，因此三個(gè)品種有8000~10000個(gè)轉(zhuǎn)錄本沒(méi)有分類注釋到Mercator庫(kù)。

5.品種間共有和特有的轉(zhuǎn)錄本

為了鑒定品種特異性轉(zhuǎn)錄本，以N22、IR64和Nipponbare三個(gè)品種的轉(zhuǎn)錄本作為blast比對(duì)庫(kù)，其它9個(gè)品種與其進(jìn)行比對(duì)（圖6）。識(shí)別到N22特有轉(zhuǎn)錄本652個(gè)，IR64特有轉(zhuǎn)錄本2426個(gè)，Nipponbare特有轉(zhuǎn)錄本349個(gè)。

6.aus特有轉(zhuǎn)錄本的差異表達(dá)分析

aus品種N22特別抗旱和耐熱脅迫，因此我們想知道在這些條件下是否有aus特異轉(zhuǎn)錄本受到調(diào)控。以N22為研究對(duì)象，分析干旱和熱脅迫下差異表達(dá)的基因。利用從發(fā)育種子中分離的RNA進(jìn)行Illumina測(cè)序，將測(cè)序數(shù)據(jù)回比到重構(gòu)的N22轉(zhuǎn)錄本。使用DESeq2基于參數(shù)FDR<0.1和|log2FC|>=1軟件識(shí)別出56個(gè)aus特異的差異表達(dá)基因。這56個(gè)差異基因進(jìn)行blast比對(duì)，其中46%比對(duì)上擬南芥，27%沒(méi)有任何注釋信息，11%僅描述了一個(gè)PFAM域或與其它植物物種的序列同源，而在水稻中僅有5%已知同源基因。

舉個(gè)例子，在高溫和干旱雙重脅迫下顯著上調(diào)的基因B12288。它在japonica和indica中均有同源基因RAB21，這個(gè)基因受干旱的誘導(dǎo)，其編碼的蛋白屬于LEA脫氫蛋白家族。與水稻其它脫氫蛋白進(jìn)行多序列比對(duì)研究（圖7），N22實(shí)際上的基因與野生稻、O. sativa ssp. japonica中其它4種脫氫酶的親緣關(guān)系密切。序列覆蓋率為89.5%，序列同源性86.0%，其中包含脫氫酶高度保守的重復(fù)區(qū)。相比japonica蛋白，N22蛋白序列與野生稻更接近。

總結(jié)
本文主要探討了Pacbio Iso-Seq獲得的轉(zhuǎn)錄本相比于Nipponbare參考基因組是否可以用于aus等水稻亞種的下游分析。此外通過(guò)這些轉(zhuǎn)錄本，我們希望發(fā)現(xiàn)水稻非生物脅迫下新的轉(zhuǎn)錄本和基因。我們的分析表明所有品種都可以鑒定出特異的轉(zhuǎn)錄本，還確定了aus亞種特異的差異表達(dá)基因。Pacbio Iso-Seq這種方法適用于其它沒(méi)有基因組或者基因組質(zhì)量不高的物種，相比對(duì)基因組組裝，這種方法更省時(shí)便宜。