低溫脅迫下不同甘蔗品種的轉(zhuǎn)錄組比較分析

唐仕云 楊麗濤 李楊瑞，

（1. 廣西農(nóng)科院甘蔗研究所，南寧 530007；2. 亞熱帶農(nóng)業(yè)生物資源保護與利用國家重點實驗室，南寧 530005）

甘蔗是我國重要的糖料作物和新型能源作物，甘蔗對低溫脅迫較為敏感，低溫寒害常給甘蔗生產(chǎn)帶來了巨大的經(jīng)濟損失［1］。選育高產(chǎn)優(yōu)質(zhì)抗寒性強的品種，可提高甘蔗抵御寒害的能力，降低蔗糖產(chǎn)業(yè)生產(chǎn)中的經(jīng)濟損失。目前甘蔗雜交育種是獲得優(yōu)良品種的主要途徑［2］，但甘蔗常規(guī)育種技術所需時間長，盲目性大。隨著現(xiàn)代生物技術的飛速發(fā)展，從分子水平改良甘蔗新品種成為了可能［3］。由于甘蔗是高度復雜的多倍體，遺傳背景十分復雜，目前甘蔗全基因組信息匱乏，現(xiàn)代生物技術輔助甘蔗常規(guī)育種的能力和幅度受限，甘蔗現(xiàn)代分子育種急需加強。轉(zhuǎn)錄組是一個生物體中所有基因表達后產(chǎn)物的綜合［4］，包括轉(zhuǎn)錄本的數(shù)量、特定發(fā)育階段的表達動態(tài)、轉(zhuǎn)錄后的修飾以及非編碼RNA的調(diào)控表達情況等，它可以從整體水平上研究基因的功能和結(jié)構(gòu)，挖掘優(yōu)良性狀相關基因。轉(zhuǎn)錄組測序技術，是通過深度測序后進行轉(zhuǎn)錄組分析的技術［5-6］，轉(zhuǎn)錄組測序能夠在單核苷酸水平上對任意物種進行檢測［7］，它可依賴或不依賴于物種的參考基因組，因而為甘蔗等非模式生物的轉(zhuǎn)錄組研究提供了新的手段［8］。本研究在前期的抗寒形態(tài)變化、抗寒光合生理等的研究基礎之上［9］，選取抗寒性具有一定差異的桂糖08-1180和ROC22為材料，采用新一代高通量Illumina轉(zhuǎn)錄組測序技術，對2個不同抗性甘蔗基因型在低溫脅迫和常溫處理下的轉(zhuǎn)錄組進行分析，研究和比較2個甘蔗品種低溫脅迫下基因表達的差異，分析與甘蔗抗寒相關的基因，為甘蔗抗寒基因挖掘、抗寒分子育種奠定基礎。

1 材料與方法

1.1 材料

以桂糖08-1180和ROC22號為材料，每個品種種植4桶，每桶種植4個蔗芽，淋足水分，待苗長出后，每桶留甘蔗4株。

1.2 方法

1.2.1 實驗處理 當甘蔗長至10片葉齡時，進行實驗處理。實驗設4個處理：桂糖08-1180常溫（G1180_CW）、桂糖08-1180低溫（G1180_DW）、ROC22常溫（ROC22_CW）、ROC22低溫（ROC22_DW），每處理設2個生物學重復。常溫處理為甘蔗在自然條件下的溫度處理，日溫度變化幅度在23.0℃至35.0℃之間。低溫處理為甘蔗在人工氣候箱中進行4℃低溫處理，常溫和低溫的處理時間均為36 h。處理完畢后，立即剪取+1葉葉片，進行液氮速凍處理。

1.2.2 RNA提取及樣品檢測 采用RNAprep pure多糖多酚植物總RNA提取試劑盒（天根）進行提取，瓊脂糖凝膠電泳檢測RNA的降解及污染程度，Nanodrop檢 測RNA純 度，Qubit Fluorometer定 量RNA濃度，Agilent 2100檢測RNA的完整性。

1.2.3 cDNA文庫構(gòu)建及測序 用帶Oligo（dT）的磁珠富集mRNA，之后加入fragmentation buffer 溶液將mRNA打斷，以打斷的mRNA為模板，用六堿基隨機引物合成一鏈cDNA，然后加入緩沖液、dNTPs和DNA polymerase Ⅰ合成二鏈cDNA，得到cDNA文庫。庫檢合格后，由北京諾禾致源生物公司進行Illumina HiSeq測序。

1.2.4 轉(zhuǎn)錄組原始數(shù)據(jù)的讀取及拼接 對Illumina高通量測序獲得的原始圖像數(shù)據(jù)，經(jīng)CASAVA堿基識別分析后，轉(zhuǎn)化為原始測序序列，稱之為原始讀子（Raw reads）。對原始讀子去除接頭、去除N（N表示無法確定堿基的信息）和去除低質(zhì)量reads，獲得clean reads。由于目前甘蔗的全基因組測序信息尚未見公布，本研究在獲得clean reads后，采用Trinity軟件進行de novo從頭拼接。對拼接得到的轉(zhuǎn)錄本序列，取最長的轉(zhuǎn)錄本作為該基因的Unigene，對所有轉(zhuǎn)錄本及Unigene進行統(tǒng)計，并以此為基礎進行后續(xù)的生物信息分析。

1.2.5 基因功能注釋 使用Blast方法，在Nr、Nt、KOG/COG、Swiss-prot、Pfam、GO、KEGG數(shù)據(jù)庫對所有的Unigene進行注釋。若在這些數(shù)據(jù)庫中未能比對上的Unigene，采用estscan（3.0.3）軟件預測其ORF。

1.2.6 基因表達水平的分析 采用RSEM軟件，以Trinity拼接獲得的轉(zhuǎn)錄組作為參考序列（ref），將每個樣品的clean reads往參考序列（ref）上做映射（mapping），并對其進行FPKM轉(zhuǎn)換，分析各基因的表達水平。本研究采用DESeq進行基因差異表達分析，設定的篩選閾值為：padj＜0.05。

1.2.7 基因GO富集和KEGG富集分析 采用GOseq進行GO富集分析。使用KOBAS（2.0）進行KEGG Pathway富集分析。

2 結(jié)果

2.1 測序數(shù)據(jù)及質(zhì)量情況

從表1可以看出，本實驗8個樣品共獲得394 907 890條raw reads，過濾后共獲得382 684 356條clean reads，測序產(chǎn)生的clean reads的堿基數(shù)合計為57.41 G。堿基錯誤率均在0.02%的較低水平，Q20值和Q30值均較高，GC含量在53.28%-55.16%之間。

表1 8個樣本的數(shù)據(jù)產(chǎn)出情況匯總表

2.2 轉(zhuǎn)錄本拼接

從表2可以看出，對8個樣品拼接后，共獲得360 404條轉(zhuǎn)錄本和183 515個Unigenes，平均長度分別為1 073 bp和741 bp，N50長度分別為1 903 bp和1 284 bp。

2.3 基因功能注釋

從表3可以看出，在NR、NT、KO、Swiss-Prot、PFAM、GO和KOG數(shù)據(jù)庫中，對所有unigenes進行比對，有110 021個unigenes能在7個數(shù)據(jù)庫中的至少1個數(shù)據(jù)庫里顯著匹配，注釋成功率為59.95%。10 825個unigenes在7個數(shù)據(jù)庫中均能匹配，占總unigenes的5.89%。

表2 轉(zhuǎn)錄本數(shù)量及長度

表3 基因注釋成功率統(tǒng)計表

2.4 基因差異表達分析

從圖1可以看出，低溫脅迫下桂糖08-1180有8 965個上調(diào)表達基因，7 180個下調(diào)表達基因。ROC22有10 898個上調(diào)表達基因，9 419個下調(diào)表達基因。

本研究進一步對2個品種的特有差異表達基因及共有差異表達基因進行分析，結(jié)果如圖2，可以看出，有13 113個差異表達基因為2個品種共有。7 204個差異表達基因只在ROC22品種中出現(xiàn)，為ROC22特有。而3 032個差異表達基因只在桂糖08-1180品種中出現(xiàn)，為桂糖08-1180特有。

2.5 桂糖08-1180 和ROC22差異表達基因的異同分析

2.5.1 桂糖08-1180 和ROC22共有差異表達基因GO和KEGG富集分析 對桂糖08-1180和ROC22共有的差異表達基因進行GO分析，結(jié)果如圖3，可以看出，共有差異表達基因富集在生物過程大類中的2個小類，分別是非生物刺激反應和低溫反應。

共有差異表達基因還富集在細胞成分大類中的20個小類，其中細胞組分、細胞、細胞部件中富集的差異基因較多。與光合作用相關的葉綠體、質(zhì)體等富集的基因較多。

對共有的差異表達基因進行KEGG Pathway分析，共有的差異表達基因富集在前20條代謝通路的結(jié)果如圖4，可以看出，富集差異基因數(shù)最多的途徑是核糖體，光合作用天線蛋白、光合作用等與光合作用相關的pathway中富集的差異基因也較多，說明低溫對光合作用的代謝途徑產(chǎn)生了較大的影響。

圖1 基因差異表達的火山圖分析

圖2 所有差異基因維恩圖

2.5.2 桂糖08-1180和ROC22特有差異表達基因GO富集分析 對桂糖08-1180和ROC22低溫脅迫下特有的差異表達基因進行GO富集分析，結(jié)果如圖5和圖6，可以看出，桂糖08-1180的特有差異表達基因在生物過程大類的17個小類得到富集，其中DNA代謝過程、DNA整合、DNA重組中富集的差異基因較多；細胞成分大類的4個小類也得到富集，其中，以RNA為目標的RNA聚合酶復合物中富集的差異基因較多；分子功能大類的9個小類得到富集，其中，在ADP結(jié)合、天冬氨酸型內(nèi)肽酶活性、天冬氨酸型肽酶活性中富集的差異基因較多。

ROC22的特有差異表達基因富集在生物過程大類的14個小類中，其中細胞代謝過程、生物合成過程、有機物生物合成過程富集的差異基因較多；細胞成分大類的8個小類也得到富集，其中細胞內(nèi)成分、細胞、細胞組分富集的差異基因較多；分子功能大類的8個小類得到富集，但富集的差異表達基因均較少。

2.5.3 桂糖08-1180和ROC22特有差異表達基因KEGG富集分析 從圖7可以看出，桂糖08-1180特有差異表達基因富集的前20條KEGG通路中，植物激素信號轉(zhuǎn)導中富集的差異基因數(shù)最多，淀粉與蔗糖的代謝、苯丙素生物合成富集的基因較多。而ROC22特有差異表達基因富集的前20條KEGG通路中，植物激素信號轉(zhuǎn)導、甘油磷酯代謝、嘧啶代謝、半胱氨酸和蛋氨酸代謝等富集的差異基因數(shù)較多。

3 討論

RNA-seq是研究生物體轉(zhuǎn)錄組信息的有效方法，測定生物體在不同處理條件下的表達信息，可以幫助我們發(fā)掘相關的功能基因［5，7］。轉(zhuǎn)錄組學是功能基因組學研究的一個重要內(nèi)容，利用轉(zhuǎn)錄組數(shù)據(jù)庫研究非模式植物的冷調(diào)節(jié)基因?qū)⑹腔蚩寺〉囊粋€重要手段，也是特殊抗性植物基因發(fā)掘的重要方法［10］。目前轉(zhuǎn)錄組測序技術在甘蔗的抗旱［11-12］、抗黑穗?。?3］、分蘗習性方面［14］應用均有報道，且均篩選到了大量具有相應功能的差異表達基因。本研究對低溫脅迫下2個不同甘蔗品種（系）進行轉(zhuǎn)錄組測序及基因差異表達分析，獲得了大量的具有一定功能的差異表達基因，這有助于從整體上了解不同甘蔗品種對低溫脅迫的不同反應機制，可為甘蔗品種的抗寒性鑒定提供理論依據(jù)。

圖3 桂糖08-1180和ROC22共有差異表達基因GO富集分析

圖4 桂糖08-1180和ROC22共有的差異表達基因的KEGG富集結(jié)果

圖5 桂糖08-1180特有差異表達基因的GO富集分析

低溫首先發(fā)生在細胞膜及膜系統(tǒng)相關系統(tǒng)中［16-17］，細胞膜最先感知低溫脅迫，通過信號傳導，啟動冷響應基因和轉(zhuǎn)錄因子，以緩解低溫脅迫的破壞［18］，植物在低溫脅迫下，通過膜運輸和感知系統(tǒng)穩(wěn)定膜結(jié)構(gòu)，維持細胞內(nèi)環(huán)境的穩(wěn)定［19］。但這種緩解能力是有限度的，一旦超出植株的忍受范圍，體內(nèi)的保護系統(tǒng)受到破壞，低溫引起膜系統(tǒng)的透性和活性變化，從而使光合作用減弱，呼吸速率大起大落，同時細胞間冰晶對細胞膜和細胞壁產(chǎn)生機械損傷，水分平衡受到破壞，最終導致植物死亡［20］。Li等［21］選用抗寒性強的甘蔗品種GT28和抗寒性弱的品種YL6，研究了低溫脅迫下甘蔗微管與葉綠體膜超微結(jié)構(gòu)的變化，認為抗寒性強的品種GT28的微管和葉綠體膜結(jié)構(gòu)穩(wěn)定且正常，而抗寒性較弱的YL6微管數(shù)顯著減少，葉綠體膜模糊降解，GT28具有更高的光合能力。在本研究中，2個品種中共有差異表達基因均在膜系統(tǒng)功能小類中得到了顯著的富集。在與光合作用相關的葉綠體、質(zhì)體、膜系統(tǒng)和細胞器等相關的GO條目和類胡蘿卜素生物合成、光合作用天線蛋白、光合作用、吡啉與葉綠素代謝等KEGG代謝通路上，富集到了大量的差異基因。說明低溫脅迫對桂糖08-1180和ROC22的光合作用相關代謝途徑產(chǎn)生了重要的影響。筆者曾以桂糖08-1180和ROC22等9個育種材料，進行低溫脅迫和常溫處理，測定了各品種的形態(tài)指數(shù)、光合及葉綠素熒光等指標，結(jié)果表明隨著低溫脅迫處理時間的延長，各品種的冷害指數(shù)不斷增大，葉片葉綠素含量均顯著降低，葉片凈光合速率、氣孔導度在常溫與低溫處理之間差異顯著，F(xiàn)v/Fm、ΦPSⅡ均顯著下降，且其變化大小因品種不同表現(xiàn)不一，桂糖08-1180的冷害指數(shù)及光合特性參數(shù)顯示其抗寒性均優(yōu)于 ROC22［9］。

圖6 ROC22特有差異表達基因的GO富集分析

植物抗寒性是一種受微效多基因控制的數(shù)量性狀，此類基因是一種誘發(fā)性基因，只有在特定的條件下才能表達，多個基因的共同表達才能提高植物的抗寒性［15］。冷脅迫中根據(jù)基因產(chǎn)物的不同分為2大類：一類在脅迫過程中直接發(fā)揮作用的功能基因編碼的產(chǎn)物，包括脯氨酸等調(diào)節(jié)物質(zhì)和抗凍蛋白等功能蛋白；另一類為調(diào)節(jié)基因，其基因產(chǎn)物參與調(diào)控下游基因表達與信號傳導［22］。植物的抗寒生理反應主要來自于冷誘導基因（COR）的表達，COR基因的表達受多種轉(zhuǎn)錄因子的控制，其中，CBF/DREB通路是誘發(fā)其表達的一種主要方式［23］。植物激素在冷脅迫中對控制依賴CBF基因通路途徑和不依賴CBF/DREB通路途徑中均為重要的信號因子［24］。甘蔗是一種起源于熱帶的喜溫性植物，其抗寒能力和機制不同于起源于低溫帶的抗凍性植物，不同品種之間的抗性機制也必然不同。本研究使用轉(zhuǎn)錄組測序的方法，對2個抗寒性不同的品種進行的差異基因的表達譜分析，結(jié)果表明了2個品種具有一定相同的抗寒性機制，同時，2個品種又獨具各自的抗性機制。桂糖08-1180表現(xiàn)在DNA代謝過程、DNA整合、DNA重組、ADP結(jié)合、天冬氨酸型內(nèi)肽酶活性、天冬氨酸型肽酶活性、淀粉與蔗糖的代謝、苯丙素生物合成途徑中，富集的差異基因較多。而ROC22表現(xiàn)在細胞代謝過程、生物合成過程、有機物生物合成過程、甘油磷酯代謝、嘧啶代謝、半胱氨酸和蛋氨酸代謝等途徑中，富集的差異基因數(shù)較多。但對不同品種富集的差異表達基因的功能、代謝途徑及相關生理特性，還需做進一步深入研究。

4 結(jié)論

通過轉(zhuǎn)錄組測序技術，共獲得了大小為57.41 G的clean bases，進行從頭拼接后，獲得了183 515個Unigenes?；虿町惐磉_分析表明，低溫脅迫下，桂糖08-1180有16 145個基因的表達發(fā)生了變化，ROC22有20 317個基因的表達發(fā)生了變化，13 113個差異表達基因為2個品種所共有，7 204個差異表達基因為ROC22品種特有，3 032個差異表達基因為桂糖08-1180品種特有。

圖7 桂糖08-1180和ROC22特有差異表達基因的KEGG富集結(jié)果

對桂糖08-1180和ROC22共有差異表達基因的功能富集分析表明，低溫脅迫下，ROC22與桂糖08-1180，在對寒害逆境反應、膜系統(tǒng)、光合作用方面具有許多相同的功能小類，富集了大量的差異基因，說明對低溫脅迫的應急響應、低溫脅迫對膜系統(tǒng)的傷害及光合作用的下降是2個品種的共性。

對桂糖08-1180和ROC22特有的差異表達基因功能富集分析表明，ROC22在有機化合物的合成、運輸和轉(zhuǎn)運蛋白活性相關的條目富集得較多，而桂糖08-1180在DNA整合、RNA聚合酶、ADP結(jié)合及代謝酶中富集的差異表達基因較多。