低溫脅迫下巴西橡膠樹三倍體全長轉(zhuǎn)錄組測序分析

毛常麗，李 玲，張鳳良，李小琴，楊 湉，吳 裕

（云南省熱帶作物科學研究所，云南景洪666100）

低溫寒害一直是限制我國橡膠樹種植面積擴大和越冬安全的主要因子，抗寒性狀也一直是我國橡膠樹良種培育的主要選擇指標之一。根據(jù)《NY/T607-2002橡膠樹選育種技術(shù)規(guī)程》中的育種程序，選出一個優(yōu)良無性系一般需要25～30年。云南省熱帶作物科學研究所培育的三倍體無性系“云研77-2”和“云研77-4”具有抗寒、高產(chǎn)特性，在云南省累積推廣種植160 000 hm2。這兩個品種從1974雜交到2000年良種審定經(jīng)歷了26年（后經(jīng)細胞學鑒定為三倍體）。一般認為，樹木三倍體具有良好的速生、抗逆和豐產(chǎn)性，實際上只能說從三倍體中更有機會獲得性狀突出的個體。橡膠樹自然三倍體的發(fā)生率特別低，要想通過傳統(tǒng)的育種程序選出優(yōu)良三倍體就更不容易。

轉(zhuǎn)錄組學是從RNA水平研究基因表達的情況，是研究細胞表型和功能的重要手段。轉(zhuǎn)錄組測序技術(shù)（RNA-Seq）通過深度測序?qū)悠分械霓D(zhuǎn)錄因子進行分析，能夠反映出轉(zhuǎn)錄水平上的基因表達情況，其已被廣泛應(yīng)用于植物種質(zhì)評價、性狀篩選等方面。篤斯越橘（Vaccinium uliginosum）、厚樸（Magnolia officinalis）、枇杷（Eriobotrya japonica）、櫻桃（Prunus avium）［1-4］等木本物種的轉(zhuǎn)錄組學抗逆性研究為橡膠樹抗逆性轉(zhuǎn)錄組學的研究提供了新的思路。

目前，橡膠樹轉(zhuǎn)錄組學主要基于第二代高通量測序平臺的RNA-Seq技術(shù)對病害、產(chǎn)膠、生殖發(fā)育等［5-7］方面進行研究，且用于測序的材料為二倍體，而對三倍體材料進行的轉(zhuǎn)錄組測序還未見報道?；谔镩g鑒定［8］、生理指標測定［9-11］等結(jié)果已確定三倍體品種云研77-2和云研77-4抗寒性比其母本GT1強，但抗寒性增強是來自于母本的2n配子中抗寒效應(yīng)基因的加倍或是其他原因則還未見報道。如果能探清其抗寒機制，將有助于三倍體抗寒品種的培育。本研究以云研77-4為材料作全長轉(zhuǎn)錄組測序及分析，旨在得到完整的轉(zhuǎn)錄組序列，為橡膠樹三倍體分子輔助選擇及抗寒三倍體育種奠定理論基礎(chǔ)。

1 材料和方法

1.1 實驗材料

材料選擇為橡膠樹三倍體無性系“云研77-4”，其親本為GT1×PR107；苗木以GT1開放授粉種子苗為砧木芽接培育；選擇長到第3蓬葉且頂蓬葉穩(wěn)定，整株無病蟲害，長勢一致的苗木，準備實施低溫脅迫處理。

1.2 方法

在4℃人工氣候箱下處理2 h、6 h和20 h，對照為不處理（0 h），每處理設(shè)3次重復(fù)。處理完后剪取第2個葉蓬中下部復(fù)葉的中間小葉，并迅速用液氮凍干后轉(zhuǎn)移至-80℃冰箱保存?zhèn)溆谩Ｏ葘Ω鳂悠贩謩e提取RNA，再將不同處理的RNA混合，對混合樣進行測序。樣品提取、文庫構(gòu)建及質(zhì)控、轉(zhuǎn)錄組測序及分析委托北京百邁客生物科技有限公司進行。參考已報道的橡膠樹全基因組序列對本研究的轉(zhuǎn)錄組序列進行注釋。

2 結(jié)果與分析

2.1 數(shù)據(jù)產(chǎn)出及統(tǒng)計

云研77-4冷處理后的樣品經(jīng)PacBio三代轉(zhuǎn)錄組測序，共獲得41.85 Gb clean data的數(shù)據(jù)，按照full passes＞=3且序列準確性大于0.9的篩選條件，從原始數(shù)據(jù)中提取到環(huán)形一致序列（Circular Consensus Sequencing，CCS）504 032條，平均每條序列的長度為2 157 bp。對CCS序列進行長度分布分析，得到其長度分布范圍主要為1 800～2 200 bp（圖1）。過濾掉CCS序列中不包含5’引物、3’引物及polyA尾的非全長序列，共得到全長非嵌合序 列（Full-length Non-Chimeric sequence,FLCN）413 878條，占總CCS序列的82.11%。

圖1 橡膠樹三倍體三代轉(zhuǎn)錄組全長測序CCS長度分布

2.2 FLCN的一致序列（consensus isoform）分析

使用SMRTLink軟件中的IsoSeq模塊對FLNC序列進行相似性序列聚類，得到一致序列157 638條，從其中篩選出準確度大于99%的高質(zhì)量轉(zhuǎn)錄本154 484條，平均每條一致序列長度為2 181 bp，其consensus isoform序列長度分布如圖2。用GMAP（Genomic Mapping and Alignment Program）將校正后的一致序列與橡膠樹參考基因組進行序列比對，用cDNA_Cupcake軟件將比對后的序列進行去冗余處理，過濾identity小于0.9，coverage小于0.85的序列，共得到非冗余轉(zhuǎn)錄本序列96 437條。對非冗余的轉(zhuǎn)錄本進行可變剪接分析，檢測到基因位點23 307個，其中新基因位點4 288個，新發(fā)現(xiàn)轉(zhuǎn)錄本70 055個。為了探究更多的轉(zhuǎn)錄組信息，同時對三倍體橡膠樹中的融合轉(zhuǎn)錄本進行了分析，共得到融合轉(zhuǎn)錄本1 983個。

圖2 consensus isoform長度分布

2.3 新轉(zhuǎn)錄本分析

2.3.1 SSR分析

從新發(fā)現(xiàn)的轉(zhuǎn)錄本中篩選500 bp以上的轉(zhuǎn)錄本，利用MISA軟件進行SSR分析，結(jié)果見表1。對檢測到的SSR進行類型劃分，劃分為完美單堿基重復(fù)（p1）、完美雙堿基重復(fù)（p2）、完美三堿基重復(fù)（p3）、完美四堿基重復(fù)（p4）、完美五堿基重復(fù)（p5）、完美六堿基重復(fù)（p6）和混合SSR（c，即包含至少兩個完美SSR，且之間距離小于100 bp），并對其進行密度統(tǒng)計，統(tǒng)計結(jié)果見圖3。

由表1、圖3可知，橡膠樹新轉(zhuǎn)錄本中存在6種類型的SSR，每種重復(fù)類型的數(shù)量差異較大，主要以單堿基重復(fù)為主，且在所有的堿基重復(fù)中，以完美單堿基重復(fù)類型最多。

表1 橡膠樹新轉(zhuǎn)錄本SSR分析

圖3 橡膠樹三代轉(zhuǎn)錄組新轉(zhuǎn)錄本SSR類型分布圖

2.3.2 新基因編碼區(qū)序列預(yù)測

對可變剪接分析中得到的新轉(zhuǎn)錄本使用TransDecoder軟件［12］預(yù) 測 其編 碼 區(qū) 序 列（Coding Sequence，CDS）及其對應(yīng)的氨基酸序列，共獲得開放閱讀框（Open Reading Frame,ORF）66 182條，其中完整ORF57 516條。在獲得的ORF中，0～500 aa的有54 475條，占82.31%，500～1 000 aa的有10 537條，占15.91%，1 000～2 300 aa的有1 167條，占1.76%（圖4）。

圖4 橡膠樹轉(zhuǎn)錄組預(yù)測的CDS編碼蛋白長度分布

2.3.3 長鏈非編碼RNA（Longnon-coding RNA，lncRNA）預(yù)測及轉(zhuǎn)錄因子分析

通過對轉(zhuǎn)錄本中不編碼蛋白的lncRNA進行編碼潛能篩選，從而判定該轉(zhuǎn)錄本是否為lncRNA。利用CPC［13］、CNCI［14］、pfam、CPAT［15］四種分析方法對新發(fā)現(xiàn)的轉(zhuǎn)錄本進行l(wèi)ncRNA預(yù)測，分別預(yù)測到5 598、7 164、15 207、17 488個lncRNA，其中4種方法均預(yù)測到的lncRNA共3 575個（圖5）。同時，根據(jù)lncRNA在參考基因組注釋信息上的位置對lncRNA進行分類繪圖，結(jié)果見圖6。由圖6可知，本研究材料中的lncRNA主要以lincRNA為主，占53.19%，其次為sense_lncRNA，占27.30%。

圖5 四種篩選方法維恩圖

圖6 lncRNA位置分類

LncRNA的存在為下一步研究橡膠樹表觀遺傳、轉(zhuǎn)錄水平及轉(zhuǎn)錄后水平調(diào)控基因的表達提供了新的思路。

轉(zhuǎn)錄因子是指能夠結(jié)合在某基因上游特異核苷酸序列上的蛋白質(zhì)，這些蛋白質(zhì)可以調(diào)控RNA聚合酶與DNA模板的結(jié)合，從而調(diào)控基因的轉(zhuǎn)錄。使用iTAK軟件［16］對三倍體橡膠樹三代測序轉(zhuǎn)錄因子進行預(yù)測，共預(yù)測到轉(zhuǎn)錄因子12 064個，且預(yù)測到的轉(zhuǎn)錄因子包括了在其他物種中有報道的RLK-Pelle_DLSV、C3H、NAC、bHLH、MYBrelated及RLK-Pelle_LRR-VIII-1等 類 型。C3H、bHLH、MYB-related等轉(zhuǎn)錄因子在植物抗逆、生長發(fā)育及有效成分合成方面均有報道，這為今后開展橡膠樹相關(guān)研究提供了重要信息。同時，為了更直觀地呈現(xiàn)全長轉(zhuǎn)錄組測序所得到的注釋序列在全基因組上的分布密度及全長轉(zhuǎn)錄組測序的完整性，我們將注釋到的轉(zhuǎn)錄本與全基因組進行了Circos可視化繪圖（圖7）。從圖7可知，本研究材料全長轉(zhuǎn)錄組序列中的基因密度、轉(zhuǎn)錄本密度和全基因組中的基因密度、轉(zhuǎn)錄本密度相當，lncRNA在各條染色體上均有分布。

圖7 基因組水平Cicros可視化圖

2.4 新轉(zhuǎn)錄本功能注釋

為探清三倍體橡膠樹新轉(zhuǎn)錄本的功能，用BLAST軟件（version 2.2.26）［17］將得到的新轉(zhuǎn)錄本與NR［18］、Swissprot［19］、GO［20］、COG［21］、KOG［22］、Pfam［23］、KEGG［24］數(shù)據(jù)庫比對，進行新轉(zhuǎn)錄本功能注釋，共注釋到功能基因66 684個。7個數(shù)據(jù)庫注釋到的轉(zhuǎn)錄本數(shù)量見表2。

從表2可以看出，7個功能數(shù)據(jù)庫中注釋到的轉(zhuǎn)錄本長度大部分在1 000 bp以上，這也進一步說明全長轉(zhuǎn)錄組測序能測得較長的序列。綜合7個功能數(shù)據(jù)庫中注釋到的轉(zhuǎn)錄本，其功能分類主要涉及細胞組分、分子功能及生物學過程三大方面，這為抗寒性相關(guān)基因的篩選奠定了基礎(chǔ)。

表2 7個數(shù)據(jù)庫注釋到的轉(zhuǎn)錄本數(shù)量

3 討論

轉(zhuǎn)錄組是研究動植物生命過程的有效方法之一，隨著測序手段的不斷更新，近幾年發(fā)展起來的三代全長轉(zhuǎn)錄組測序因其能讀取更長的序列而被廣泛應(yīng)用于植物的生長發(fā)育、非生物脅迫等領(lǐng)域。筆者所在單位云南省熱帶作物科學研究所在國際上首次獲得了達到染色體級別的高質(zhì)量巴西橡膠樹優(yōu)良品種GT1的參考基因組序列［25］這為本研究的轉(zhuǎn)錄組注釋、組裝提供了重要參考。

橡膠樹轉(zhuǎn)錄組測序已經(jīng)應(yīng)用于自根幼態(tài)無性系體胚發(fā)育［26］、橡膠樹叢枝病機理［5］及苗期生長優(yōu)勢［27］等方面，并從中篩選了相關(guān)的調(diào)控基因，為橡膠樹分子輔助育種奠定了重要的理論基礎(chǔ)，但是現(xiàn)有的報道均基于Illumina HiSeq的二代轉(zhuǎn)錄組測序技術(shù)，且測序的樣品均為二倍體。本研究通過PacBio三代轉(zhuǎn)錄組測序技術(shù)對三倍體橡膠樹抗寒無性系進行測序，共產(chǎn)生41.85 Gb clean data的數(shù)據(jù)，分析得到高質(zhì)量轉(zhuǎn)錄本154 484條，平均每條序列長度為2 181 bp，數(shù)據(jù)得量、轉(zhuǎn)錄本總數(shù)以及單條轉(zhuǎn)錄本序列長度均大于先前報道的橡膠樹二代轉(zhuǎn)錄組數(shù)據(jù)［28-29］。將得到的數(shù)據(jù)分析后得到新轉(zhuǎn)錄本70 055條。對94 084條轉(zhuǎn)錄本進行分析，發(fā)現(xiàn)了潛在的80 432個SSR分子標記，高于Wirulda等［30］用454高 通 量 測 序 技 術(shù) 獲 得 的24 000個SSRs，這些標記可廣泛應(yīng)用于橡膠樹種質(zhì)資源評價、分子輔助育種、抗逆性鑒定等。利用CPC、CNCI、pfam、CPAT四種方法預(yù)測到共有的lncRNA 3 575個，同時對新轉(zhuǎn)錄本進行了NR、Swissprot、GO、COG、KOG、Pfam、KEGG功能注釋，共注釋到功能基因66 684個，高于前人注釋到的37 432個功能基因［28］。轉(zhuǎn)錄因子可與基因5’端上游特定序列專一性結(jié)合，從而保證目的基因在特定的時間與空間進行表達。本研究中共預(yù)測到轉(zhuǎn)錄因子12 064個，其中包含了與抗逆性相關(guān)的MYB、C3H轉(zhuǎn)錄因子，下一步將從中篩選出與三倍體橡膠樹抗寒性相關(guān)的轉(zhuǎn)錄因子，為研究三倍體抗寒性提供參考。從全長轉(zhuǎn)錄組數(shù)據(jù)和全基因組數(shù)據(jù)Circos可視化圖可以看出，全長轉(zhuǎn)錄組數(shù)據(jù)在全基因組測序得到的數(shù)據(jù)覆蓋率相當，說明全長轉(zhuǎn)錄組序列可靠，可為下一步二代轉(zhuǎn)錄組序列的拼接、組裝及注釋提供參考。

4 結(jié)論

本研究基于三代測序技術(shù)建立了三倍體橡膠樹的轉(zhuǎn)錄組數(shù)據(jù)庫，并對序列進行分析，為后續(xù)的橡膠樹分子輔助育種、二代轉(zhuǎn)錄組數(shù)據(jù)組裝，特別是三倍體橡膠樹的轉(zhuǎn)錄組分析提供重要的基礎(chǔ)，下一步將繼續(xù)從lncRNA、轉(zhuǎn)錄因子等方面分析三倍體橡膠樹抗寒性標記及強抗寒性的分子機制，為三倍體選育種提供有用標記及參考。