干旱脅迫下割手密根系轉錄組差異表達分析

劉洪博，劉新龍，蘇火生，陸鑫，徐超華，毛鈞，林秀琴，李純佳，李旭娟，字秋艷

（云南省農業(yè)科學院甘蔗研究所/云南省甘蔗遺傳改良重點實驗室，云南開遠661699）

干旱脅迫下割手密根系轉錄組差異表達分析

劉洪博，劉新龍，蘇火生，陸鑫，徐超華，毛鈞，林秀琴，李純佳，李旭娟，字秋艷

（云南省農業(yè)科學院甘蔗研究所/云南省甘蔗遺傳改良重點實驗室，云南開遠661699）

【目的】探明云南割手密 82-114受干旱脅迫的內在分子機制，挖掘與抗旱密切相關的功能基因，提高割手密抗旱基因的研究與利用。【方法】以干旱脅迫24、48和72 h的云南割手密82-114的根系進行Illumina HiSeqTM 4000高通量轉錄組測序分析，將組裝得到的Unigene分別與Swiss-Prot、Nr、KOG、Pfam和KEGG數(shù)據(jù)庫進行比對，利用RPKM來衡量各樣本間的基因表達豐度，以FDR≤0.05和|log2 fold change|≥1來評估樣本間的差異表達基因，通過Gene Ontology（GO）數(shù)據(jù)庫、KEGG pathway數(shù)據(jù)庫對不同干旱脅迫樣本差異表達基因的功能和參與的調控路徑進行分析?！窘Y果】未脅迫處理的CK與脅迫24、48和72 h后的樣本轉錄組分別檢測到134 724、130 368、133 564和131 321個表達基因，與CK相比，各脅迫樣本顯著上調（或下調）的差異表達基因分別有3 061（1 302）、2 304（2 841）和3 236（2 525）個；以FDR值= 0為篩選條件，3個樣本的極顯著差異表達基因主要參與轉錄激活、水分運輸、DNA結合、ATP結合及細胞膜、跨膜運輸和防御反應等有關代謝活動。GO富集分析表明：3個脅迫處理的樣本在生物學途徑中富集最多的類別并不完全相同，其中，脅迫24 h與48 h的樣本富集最多的類別是與DNA依賴型轉錄和轉錄調節(jié)子相關的基因，脅迫72 h樣本富集最多的類別是與翻譯相關的基因；而在細胞組件和分子功能中，3個脅迫樣本富集最多的基因類別均是與內在的膜、核和ATP結合相關的基因。KEGG富集分析表明：脅迫24 h的樣本有2 248個差異表達基因參與了107條代謝途徑，脅迫48 h的樣本有2 114個差異表達基因參與了130條代謝途徑，脅迫72 h的樣本有2 392個差異表達基因參與了144條代謝途徑；經KEGG顯著性富集分析，篩選到鞘糖脂生物合成途徑、MAP激酶信號途徑和ABC轉運蛋白等與非生物脅迫或逆境相關?！窘Y論】獲得3個干旱脅迫樣本與CK樣本間差異表達基因的變化及其功能信息，發(fā)現(xiàn)參與云南割手密82-114干旱脅迫響應相關的細胞外信號調節(jié)激酶、磷脂酶D、β-氨基己糖苷酶和ATP結合盒B亞族（MDR/TAP）-1等，獲得在整個干旱脅迫時期穩(wěn)定上調表達基因70個，穩(wěn)定下調表達基因11個，通過同源基因序列的比對分析，闡明這些基因在各自的代謝通路中被強烈誘導或抑制表達，顯示與干旱脅迫存在密切關系。

割手密；干旱；轉錄組；差異表達基因

0 引言

【研究意義】割手密（Saccharums pontaneum L.）又稱甘蔗細莖野生種，是禾本科蜀黍族甘蔗亞族甘蔗屬多年生草本植物[1]，成熟莖微甜，因此又稱“甜根子草”，具有發(fā)達的地下橫走莖，宿根性好、分蘗能力強、固土力強，一般生長于田埂、坡地和干旱沙地中，是甘蔗品種改良最重要的野生親本和甘蔗抗旱、耐瘠等重要基因源[2-3]。在現(xiàn)代甘蔗品種的基因組中，割手密基因組約占10%—25%[4-5]，對甘蔗產量、抗性育種具有重要作用。由于近年來全球氣候變化的影響，農作物干旱問題日益突出，對山坡地種植的甘蔗作物尤為明顯[6]。割手密發(fā)達的地下橫走莖，葉片窄小的形態(tài)特征，且適宜在熱帶、亞熱帶的沙壤環(huán)境中生存，對干旱具有極強的耐受性。因此，甘蔗育種工作者十分重視割手密的抗旱、強分蘗等優(yōu)良特性對甘蔗遺傳改良的作用?！厩叭搜芯窟M展】由于割手密資源分布范圍廣，遺傳多樣性豐富[7]，在資源利用上具有一定的盲目性。蘇火生等[8]、齊永文等[9]通過割手密初級核心種質庫的構建，成功篩選出代表不同生態(tài)類型和基因型的核心種質，為割手密后續(xù)的研究利用提供了參考；同時國內外學者對割手密不同無性系的抗旱水平進行了較多的研究，認為在旱坡地、干旱的沙壤土中生長的割手密具有較好的抗旱和耐旱能力[10-12]；這類土壤環(huán)境下的土壤持水量差，短時間的高溫天氣即可導致土壤和植物缺水，因此割手密在生長旺盛時期的抗旱性可能具有短時效應；同時割手密在干旱脅迫下的電導率較低，細胞破壞程度小，抗旱能力強于甘蔗品種[13]。與其他作物的野生親本相比，割手密抗旱相關基因的研究還很少見報道。姚艷麗等[14]通過同源基因比較的方法，克隆了割手密 2個與抗旱相關的DREB轉錄因子，通過原核表達分析，證明這兩個基因是功能基因；國外學者PARK等[15]通過轉錄組學研究，證明割手密TUS05-05的SspNIP2在冷處理30 min后表達量明顯增加，同時通過轉SspNIP2植株的水分脅迫測試表明，轉基因植株比正常植株的蒸騰速率低，證明了該基因與抗寒和抗旱性狀具有一定的關系?！颈狙芯壳腥朦c】植物的抗旱機理是一個非常復雜的過程，僅通過部分抗旱相關基因的克隆與功能研究，不能更好地詮釋抗旱機理，割手密作為甘蔗遺傳改良重要的野生親本，其本身受到野外環(huán)境的脅迫，抗旱性與甘蔗相比具有一定的遺傳優(yōu)勢。然而割手密不同干旱脅迫時期的轉錄組數(shù)據(jù)及其表達譜等方面的研究還未見報道，且與割手密抗旱相關功能基因的研究也非常有限；目前，中國甘蔗主產區(qū)主要分布在廣西和云南兩省，且主產區(qū)大部分為濕熱氣候生態(tài)型。因此，非常有必要利用該生態(tài)類型的割手密資源進行干旱脅迫下的轉錄組分析研究，它不僅可以探尋到適合該區(qū)域的抗旱類型，從育種的角度分析，還可以使選育出的品種具有本土生態(tài)基因型，更適合于當?shù)胤N植。云南割手密82-114屬于濕熱氣候生態(tài)型資源，株高、莖徑等性狀表現(xiàn)優(yōu)異，是割手密種質中利用較為成功的種質；通過前期的抗旱性篩選評價，云南割手密82-114的葉片持綠性與其他割手密材料相比有較好的表現(xiàn)，水分澆灌后恢復快，耐短時高溫及干旱環(huán)境，適合中國甘蔗種植區(qū)域的氣候類型；且利用該種質雜交獲得的F2種質中有5份材料抗旱性超過中國主栽品種ROC22，3份材料抗性中等[16]，說明該種質在抗旱性育種方面具有一定的潛力?！緮M解決的關鍵問題】本研究以內陸濕熱生態(tài)型割手密云南82-114為試驗材料，通過?lluminaHiSeqTM4000高通量測序平臺，比較該種質在拔節(jié)期干旱脅迫24、48和72 h的根部轉錄組水平差異，分析可能參與抗旱或耐旱調控的相關基因，并進一步揭示割手密受干旱脅迫的響應機理，為割手密在甘蔗抗旱育種研究中的應用奠定分子基礎。

1 材料與方法

1.1 試驗材料的培育及脅迫處理

以來自于國家甘蔗種質資源圃內的濕熱型割手密云南82-114為材料，于2015年4月10日種植于云南省農業(yè)科學院甘蔗研究所抗旱溫室，土壤以沙和紅土按1﹕1的比例混合（此土壤基質保水能力弱，易于水分脅迫），栽種于12個直徑為40 cm且含有15 cm厚的珍珠巖桶中?；谇捌陬A試驗結果，植株在這種栽培條件下停水24 h后即可出現(xiàn)根部缺水癥狀。試驗在割手密拔節(jié)期前按正常澆水處理，待拔節(jié)期開始第15天早晨對桶栽材料進行正常灌溉，下午14時取3桶對照（CK）材料的根系，此后其余材料停止?jié)菜⒃?4 h（植株表現(xiàn)為葉片萎蔫，早晨植株無吐水現(xiàn)象，標記為A）、48 h（植株表現(xiàn)為中度干旱，葉片萎蔫下垂，下部葉鞘脫水，標記為B）、72 h（植株表現(xiàn)為重度干旱，植株下部葉片、葉鞘枯黃，土壤極度干旱，標記為C）后分別取一次根系組織，每次取3桶，混合取樣，并將根系樣品放入液氮中迅速冷凍，置于-70℃超低溫冰箱中保存?zhèn)溆谩?/p>

1.2 割手密根總RNA的提取及cDNA的合成

將采取的樣品用干冰盒包裝好，送杭州聯(lián)川生物技術有限公司提取RNA和測序。利用Nanodrop檢測RNA質量，并經瓊脂糖凝膠檢測RNA完整性。樣品總RNA提取后，經Oligo（dT）的磁珠富集純化mRNA。經抽提的mRNA被隨機打斷成短片段后，以片段化的mRNA為模板，合成cDNA，并通過純化、粘性末端修復為平末端，3′末端加堿基A并加接頭，片段選擇，最后進行PCR擴增，經文庫質量檢測合格后，安排上機測序。高通量測序用?lluminaHiSeqTM4000進行，整個流程由杭州聯(lián)川生物技術有限公司完成。

1.3 測序數(shù)據(jù)質量的評估和 RAW數(shù)據(jù)的處理、Trinity組裝

由于高通量測序錯誤率會隨著測序序列的讀長的增加而升高[16]，需要對原始數(shù)據(jù)進行質量評估，分析不同長度片段測序的錯誤率。將由測序儀造成的錯誤率小于0.1%的堿基所占比例（Q30）來衡量測序的質量[17]。將測序儀測得的初始數(shù)據(jù)（raw reads）以FASTQ文件進行存儲，并將Raw reads去除adapter，去除Raw reads中的雜質數(shù)據(jù)，得到clean reads，計算Raw Data和過濾后的Valid Data的比率。利用短reads組裝軟件Trinity進行序列的組裝[18]，將組裝好的Contig進行聚類，用reads驗證每個Component 的reads支持情況，并通過Butterfly合并在de Bruijn圖中有連續(xù)節(jié)點的線性路徑，剔除read支持少的路徑。通過統(tǒng)計每個樣品由raw reads 得到的測序堿基數(shù)、clean reads占Raw reads的比例、clean reads比對到contig的比例以及測序飽和度分析，來評估測序的質量[19]。

1.4 割手密干旱脅迫下差異基因的篩選

利用組裝好的基因作庫，用序列相似性比對的方法求各基因在各樣本中的表達豐度，基因的表達量利用RPKM[20]（Reads Per Kilobase of exon model per Million mapped reads）值來衡量，RPKM值計算方法如下：

統(tǒng)計CK、A、B和C 4個樣本的基因表達豐度，從整體水平觀察樣本間的表達情況和基因表達水平的離散程度；以P-value（假設認為2個不同樣本間的所有基因都不存在差異的情況下出現(xiàn)的極端情況概率）＜0.01為非常顯著，利用錯誤發(fā)現(xiàn)率（false discovery rate，F(xiàn)DR）≤0.001為閾值[21]（FDR值越低，表明基因表達差異越顯著），篩選出樣本間的差異表達基因（differentially expressed genes，DEGs），并從差異倍數(shù)（log2foldchange≥2）和顯著性水平（P-value）進行評估，以未進行水分脅迫處理的CK為對照，對A、B、C 3個脅迫樣本表達基因的上下調關系進行統(tǒng)計。

1.5 差異表達基因的功能注釋與分析

Gene Ontology功能顯著性富集分析首先把所有得到的差異表達基因向Gene Ontology數(shù)據(jù)庫（http://www. geneontology.org）的各個 term映射，計算每個 term的基因數(shù)目，然后應用超幾何檢驗，找出與整個基因組背景相比，在DEGs中顯著富集的GO條目，其計算公式為：

其中，N為所有基因總數(shù)；n為N中DEGs的數(shù)目；M為所有基因注釋為某特定GO term的基因數(shù)目；m為注釋為某特定GO term的DEGs數(shù)目。以P-value≤0.05為閾值，滿足此條件的GO term定義為在DEGs中顯著富集的GO term。通過GO功能顯著性富集分析確定DEGs行使的主要生物學功能。同時將DEGs序列在 KEGG數(shù)據(jù)庫（kyoto encyclopedia of gene and genomes，京都基因與基因組百科全書）中進行注釋[22]，其顯著性富集分析的P-value值計算和閾值同GO富集分析相同。通過Pathway顯著性富集確定DEGs參與的最主要生化代謝途徑和信號轉導途徑。

2 結果

2.1 轉錄組測序質量檢測和RAW數(shù)據(jù)整理分析

通過對照樣品CK和A、B、C脅迫組樣品的轉錄組測序，分別獲得50 408 133、46 705 498、43 087 971和49 777 187條Raw reads，按照測序儀造成的錯誤率小于0.1%的堿基來衡量測序質量，如圖1所示，CK、A、B和C樣品的質量檢測，測序儀讀取1—99個堿基質量的平均值均高于Q30，測得序列的GC含量占47%，其他樣本的Q30和GC含量均符合上述質量要求，說明轉錄組測序質量較高，可以進行下一步的數(shù)據(jù)整理和分析。

經過原始Raw數(shù)據(jù)的整理和篩選分析，CK、A、B和C樣本的有效數(shù)據(jù)量均有所下降，與Raw數(shù)據(jù)相比，篩選后的有效數(shù)據(jù)量分別占98%、98.5%、97.8%和98.7%（表1）。經過denovo的拼接結果顯示，Transcript的總量為346 744，而unigene的總量為152 225，其中CK樣本總量為134 724個，A樣本總量為130 368個，B樣本總量為133 564個，C樣本總量為131 321個。轉錄本的長度分布主要集中在200 bp左右，并隨著長度的增加片段數(shù)量迅速下降，在片段長度達到 1 900 bp時降低到5 426個，此后大于2 000 bp的片段總計又增加到45 307個，與unigene長度的分布相比，整體分布較一致。但 unigene的頻數(shù)有所下降，說明存在一定量的可變剪接或單倍型序列。

2.2 樣本基因的注釋和差異表達基因的篩選分析

將獲得的 unigene基因分別在 Swiss-Prot、Nr、Pfam、KEGG和KOG中注釋，結果云南割手密82-114在對照和不同脅迫處理組中共包含152 225個基因，其中70 355個基因得到注釋，81 870個基因未得到注釋（表2），在注釋的基因中，來自Nr數(shù)據(jù)庫注釋的基因最多，達到70 186個。

利用RPKM來衡量各樣本間的基因表達豐度，將各樣本中RPKM值用最大值、上四分位數(shù)、中值、下四分位數(shù)和最小值來統(tǒng)計，以 RPKM盒型圖來表示（圖3），樣本A和C的表達豐度高于CK和B樣本，但B樣本的中位數(shù)水平整體高于其他樣本；從脅迫樣本間基因表達的離散程度來看，A和C樣本的離散程度較高，B樣本的離散程度相對較低。

表1 數(shù)據(jù)產量統(tǒng)計Table 1 The statistics of data yield

圖1 CK、A、B和C樣本測序錯誤率檢測Fig. 1 The sequencing error rate detection of CK, A, B and C samples

表2 割手密unigene基因的BLAST注釋結果Table 2 The unigene annotation result of Saccharum spontaneum for BLAST

以FDR≤0.05和|log2fold change|≥1來評估樣本間的差異基因，與CK相比，A、B、C 3個樣本在分別脅迫24、48和72 h后，上調和下調基因均明顯增加，其中，C樣本的上調基因為78 644個，是3個樣本中上調表達基因數(shù)目最多的，其次為B樣本，A樣本最少；下調表達最多的樣本為A，最少的樣本為C（表3）。

以各樣本中前50個DEGs為篩選對象，F(xiàn)DR值=0為篩選條件，共篩選出79個DEGs，其中上調表達占39個，下調表達占40個。3個脅迫樣本間的極顯著上調和極顯著下調表達基因數(shù)不同，其中A樣本包含2個上調，10個下調基因；B樣本包含23個上調和18下調基因；C樣本包含13上調和12下調基因。從A、B和C 3個樣本的前50個差異基因來看，大部分涉及轉錄激活、水分運輸、DNA結合、ATP結合和蛋白酶結合相關基因，以及與細胞膜、跨膜運輸、蔗糖代謝和防御反應相關的基因。從不同的文獻報道來看，這些基因的激活和過量表達與干旱脅迫有著十分密切的關系[23]。

圖2 樣本表達水平RPKM盒形圖Fig. 2 The RPKM box map of expression levels for samples

表3 差異基因的比較分析Table 3 The comparative analysis of differential genes

2.3 差異表達基因的Gene Ontology富集分析

Gene Ontology（簡稱GO）是一個國際標準化的基因功能分類體系，能夠全面描繪生物體中基因和基因產物的屬性。根據(jù)云南割手密82-114在干旱脅迫下篩選的不同DEGs，研究DEGs在GO中的分布狀況，闡明試驗中樣本處理差異在基因功能上的體現(xiàn)。以CK作為對照，A、B、C 3個脅迫處理的樣本經過GO功能顯著性富集分析，在參與的生物過程、細胞組件和分子功能中，A和C樣本均有47個小類別，而B樣本則缺少一個代謝通路（圖3）。與CK對比，在生物學途徑中，A樣本共有1 480個DEGs被注釋，其中富集最多的是與 DNA依賴型轉錄（transcription DNA-dependent）相關和轉錄調節(jié)子（regulation of transcription）相關的基因；B樣本有2 499個DEGs被注釋，其中富集最多的類別與A相同；C樣本有3 165個DEGs被注釋，其中富集最多的類別是與翻譯（translation）相關的基因。在細胞組件中，A樣本共有1 517個DEGs被注釋，其中與內在的膜（integral to membrane）和核（nucleus）相關的基因富集最多；B樣本共有2 629個DEGs被注釋，富集最多的基因類別與A樣本相同；C樣本共有3 326個DEGs被注釋，其中富集最多的基因類別同樣與A相同。在分子功能中，A樣本有1 760個DEGs被注釋；B樣本有3 093個DEGs被注釋；C樣本有3 825個DEGs被注釋；而3個樣本中富集最多的類別均是與ATP結合相關的基因。從圖3中可以看出，A、B和C樣本的差異基因數(shù)小于GO注釋的差異基因數(shù)目，說明差異基因中存在多拷貝現(xiàn)象，符合異源多倍體植物含有多個單倍型基因的特點。

2.4 差異表達基因的KEGG富集分析

通過 KEGG富集分析，確定云南割手密 82-114不同干旱脅迫時期，DEGs參與的最主要生化代謝途徑和信號途徑。找出與整個基因組背景相比，在DEGs中顯著性富集的Pathway。結果表明，A樣本有2 248個DEGs注釋到KEGG數(shù)據(jù)庫中107個分類代謝途徑中，B樣本有2 114個DEGs注釋到KEGG數(shù)據(jù)庫中130個分類代謝途徑中，C樣本有2 392個差異表達基因注釋到KEGG通路中144個分類代謝途徑中。以CK作為對照，經KEGG顯著性富集分析，如表4所示，A樣本上調的代謝通路有2個，下調的代謝通路有3個；B樣本上調的代謝通路有1個，下調的代謝通路有6個；C樣本下調的代謝通路有2個，上調的代謝通路有3個。

在上調的代謝途徑中，MAP激酶（促分裂原活化蛋白激酶）信號途徑在A樣本和C樣本中同時出現(xiàn)，該途徑可將不同的細胞膜感受器與細胞應答聯(lián)系起來，響應各種生物和非生物脅迫，在植物激素信號傳遞、細胞發(fā)育和分裂過程中起著重要的作用[24]。近年來，隨著功能獲得型突變體的獲得，MAP激酶級聯(lián)途徑在干旱脅迫下信號傳導中的功能和作用得到了廣泛的驗證。利用獲得的MAP激酶序列進行蛋白相似性比對和進化分析，所獲得的 MAP激酶與水稻的OsMAPK5有 84%的同源性，與栗（Setariaitalical）的MAPK5有88%的同源性，相似性最高，在蛋白進化分析中距離最近；而在水稻中的研究結果表明，水稻OsMAPK5的過表達對非生物脅迫的抵抗力增加，對非生物脅迫具有正調控作用[25]，在鹽、干旱、37℃處理都能誘導水稻OsMAPK5表達[26]，與割手密云南82-114所獲得的上調控MAP代謝途徑相一致，說明A和C樣本所獲得的MAP激酶有可能參與非生物脅迫的正調控。

同時在KEGG富集中，3份脅迫樣本中僅有A樣本含有醚脂類代謝途徑，且表現(xiàn)為上調，參與該途徑調控的磷脂酶D（PLD）不僅是生物膜的骨架成分，而且能夠通過水解后產生的產物來參與多種生理過程以及環(huán)境刺激引起的細胞反應，在PLD活性受反義抑制后，擬南芥葉片對ABA的敏感性下降，因此會阻礙氣孔關閉，導致水分缺失，這種現(xiàn)象在擬南芥 PLD缺失突變株[27]和更蘇（Craterostigma plantagineum）植物的PLD活性檢測[28]中已經得到證實；本研究在A樣本中篩選到的差異表達基因PLD與玉米PLD1有94%的同源性，且在干旱脅迫中處于上調表達，表明PLD在割手密云南82-114干旱脅迫進程中同樣發(fā)揮著重要調控作用。

在A與B樣本中同時存在下調的代謝途徑有鞘糖脂生物合成途徑、糖胺聚糖降解途徑，這些途徑都涉及β-氨基己糖苷酶（β-HEX），且都是起到下調作用；鞘脂是細胞質膜和葉泡膜結構的重要成分，在植物膜穩(wěn)定性、細胞信號和逆境反應中發(fā)揮重要的作用[29-32]，如圖4所示，鞘糖脂生物合成的2個途徑都受到β-氨基己糖苷酶和α-半乳糖苷酶的調節(jié)，且在A和B樣本中都處于下調表達，在擬南芥的β-氨基己糖苷酶的研究中，已經發(fā)現(xiàn)多個亞型，其中HEXO2和HEXO3蛋白主要位于質膜上[33]，說明水分脅迫后鞘脂在膜的穩(wěn)定性和細胞信號傳導的作用受到影響。在圖4中，β-氨基己糖苷酶下調導致催化形成的神經酰胺類物質減少，而神經酰胺具有啟動細胞的能力，促進細胞的新陳代謝，促使角質蛋白有規(guī)律再生的功能，說明受干旱脅迫的A和B樣本在根系細胞的新陳代謝、能力的恢復中產生了負面影響，增加上述基因的表達活性可對延緩細胞衰老，增加根系再生具有重要作用。

表4 KEGG代謝途徑顯著性富集的代謝通路Table 4 The significant enrichment in metabolic pathways from KEGG

3 討論

本研究應用 RNA-seq技術對云南割手密 82-114進行干旱脅迫前后的轉錄組比較分析，經 GO和KEGG的顯著性富集分析，篩選的差異基因類型與前人在甘蔗[34]、小麥[35]、楸樹[36]等分析的結果基本一致，主要表現(xiàn)為植物膜結合和信號傳導相關的基因最多。由于甘蔗和割手密的基因組還沒完成測序，本研究從不同脅迫時間段對割手密進行轉錄組測序分析，探討割手密隨著干旱脅迫時間的變化，其差異表達基因、代謝通路是否有一定的差異，尋找與干旱脅迫相關的抗性基因；同時3個脅迫組分別測序，對增加數(shù)據(jù)分析的可靠性具有一定好處。通過GO富集分析，3個脅迫樣本在參與的生物學途徑并不完全相同，其中A和B樣本富集的DNA依賴型轉錄相關和轉錄調節(jié)子相關的基因，從材料的脅迫程度來看，屬于對干旱脅迫的生理應答反應，但C樣本富集最多的類別則是與翻譯相關的基因，這可能與部分細胞膜、各細胞器受到破壞有一定的關系。但這些變化是否隨著干旱脅迫時間的延長代謝通路發(fā)生改變，目前還并不清楚。

值得肯定的是，3個脅迫樣本在GO注釋的分子功能和細胞組件中富集最多的類別是相同的，通過 3個脅迫樣本的交集分析，如圖5所示，在上調差異表達基因中，A、B、C 3個樣本共有的差異表達基因有70個，主要包含氧化壓力應答、碳水化合物代謝過程、電子載體活性、脂類轉運、細胞壁組織和膜等相關的基因；在下調的差異表達基因中，A、B、C 3個樣本共有的差異表達基因僅有 11個，其中以基因序列comp117244_c0_seq7參與的代謝路徑最多，主要包括過氧化氫應答、高溫應答、病毒應答、高光強度應答、液泡膜、ATP結合、細胞壁和泛素蛋白連接酶結合相關代謝過程；其次為comp114714_c0_seq3，GO功能注釋顯示，在擬南芥中參與鹽脅迫應答、蛋白結合、光形態(tài)發(fā)生、細胞核和鋅離子結合相關的代謝過程。而這些富集的差異基因與上述KEGG富集分析結果和甘蔗干旱脅迫下轉錄組變化的結果是一致的[24]。說明割手密與甘蔗在干旱脅迫后的調控機制方面存在一定的共性。從3個樣本在GO富集分析結果來看，A、B、C 3個樣本差異基因數(shù)量都是不同的，割手密隨著干旱脅迫時間的變化，其差異表達基因和代謝通路同樣發(fā)生變化，在不同的時間段啟動的干旱脅迫應答機制不同，說明在研究植物干旱脅迫轉錄組分析時，不同脅迫時間上調和下調的表達基因不同，響應機理也不完全相同，C樣本由于處于極度的干旱水平，部分葉片已經枯死，因此，隨著細胞內毒素物質的增多，對細胞的代謝產生嚴重影響，很多基因的表達處于一種紊亂的狀態(tài)，抗旱功能基因的差異表達難以準確判斷。通過3個脅迫樣本的交集分析，70個上調基因和11個下調基因很可能對割手密抗旱性起著重要的調節(jié)作用，其所參與的代謝路徑值得進一步的研究。

圖4 鞘糖脂生物合成途徑-環(huán)球系列Fig. 4 Glycosphingolipid biosynthesis-globo series

圖5 A、B、C 3個脅迫樣本的交集分析Fig. 5 The intersection analysis of A, B, and C samples under drought stress

在差異表達基因所參與的代謝通路中，除了A、B、C 3個樣本都有交叉的代謝通路外， B和C樣本中也富集了一些不同代謝通路的差異基因，如B樣本中獨有的脂氧合酶（lipoxygenase），C樣本中的甘氨酸羥甲基轉移酶（Glycine hydroxymethyl-transferase）和乙酰乳酸合成酶（Acetolactate synthase）。B樣本富集的脂氧合酶（LOX）是植物十八碳酸代謝途徑中的關鍵酶，該酶對植物的生長發(fā)育及其對環(huán)境脅迫的反應具有重要的影響，可降解細胞膜，與植物的衰老功能相聯(lián)系[37]；在擬南芥[38]和水稻[39]水分脅迫的試驗中，定量表達分析證明了水分脅迫樣本LOX表達活性上升，且水稻的LOX3啟動子中也證明存在一些逆境誘導元件。但本研究中割手密水分脅迫處理 48 h后LOX處于下調表達，與已報道的結果不相符合，需要進一步研究驗證。在C樣本中富集的甘氨酸羥甲基轉移酶在干旱脅迫中的報道并不多，但在水稻干旱脅迫的蛋白組學研究中篩選到該蛋白，在植物干旱脅迫后的具體功能尚不清楚。

ABC轉運蛋白廣泛存在于植物的細胞質膜、葉泡、線粒體和過氧化物酶中，其蛋白數(shù)目已超過 100多種[40]，隨著轉運底物的差異在植物生長素的外源毒素的解毒、氣孔調節(jié)等生理活動中發(fā)揮作用，是植物器官生長、環(huán)境脅迫等重要的轉運蛋白。經過序列的比對分析，該基因編碼的蛋白屬于ABC轉運蛋白 B亞族MDR/TAP-1類蛋白成員，在擬南芥中該基因與除草劑交叉抗性有關，參與擬南芥體內有毒物質的外排過程和次生代謝物質的跨膜運輸[41]。本研究在KEGG富集的通路中，B和C樣本都含有ABC轉運蛋白，不同的是，B樣本在脅迫48 h后處于上調表達，C樣本在脅迫72 h后處于下調表達，說明不同脅迫時間該蛋白的表達水平不同，B樣本的干旱脅迫較C樣本輕，可能處于對干旱反應的生理調節(jié)階段，對植株干旱處于正向調節(jié)，而C樣本的ABC轉運蛋白處于下調水平，可能細胞的代謝處于一種極度紊亂水平，影響了該基因的正常表達，使該基因不能發(fā)揮正常的功能，但這些猜測還需要后續(xù)的進一步研究驗證。

4 結論

通過轉錄組分析，獲得云南割手密82-114在干旱脅迫下的基因表達譜數(shù)據(jù)，闡明了參與干旱調控的信號傳遞、脂類代謝、多糖代謝、水分運輸、DNA結合等相關基因的表達差異；通過富集的差異基因同源比對分析，證明了ATP結合盒B亞族(MDR/TAP) -1、脂氧合酶、β-氨基己糖苷酶、磷脂酶D和細胞外信號調節(jié)激酶等都與干旱脅迫存在一定的相關性，并在其他作物上得到證實；獲得在整個干旱脅迫時期穩(wěn)定上調表達基因70個，穩(wěn)定下調表達基因11個，闡明了這些基因在各自的代謝通路中被強烈誘導或抑制表達，顯示與干旱脅迫存在密切關系。

[1] 郭本兆. 中國植物志. 北京: 科學出版社, 1987, 10(2): 40. GUO B Z. Flora of China. Beijing: Science Press, 1987, 10(2): 40. (in Chinese)

[2] 陳義強, 鄧祖湖, 郭春芳, 陳如凱, 張木清. 甘蔗常用親本及其衍生品種的抗旱性評價. 中國農業(yè)科學, 2007, 40(6): 1108-1117.

CHEN Y Q, DENG Z H, GUO C F, CHEN R K, ZHANG M Q. Drought resistant evaluations of commonly used parents and their derived varieties. Scientia Agricultura Sinica, 2007, 40(6): 1108-1117. (in Chinese)

[3] BREMER G. Problems in breeding and cytology of sugar cane. Euphytica, 1961, 12(1): 178-188.

[4] GLASZMANN J C, FAUTRET A, NOYER J L, FELDMANN P, LANAUD C. Biochemical genetic markers in sugarcane. Theoretical and Applied Genetics, 1989, 78(4): 537-543.

[5] D’HONT A, LU Y H, FELDMANN P, GLASZMANN J C. Cytoplasmic diversity in sugarcane revealed by heterologous probes. Sugar Cane, 1993, 1: 12-15.

[6] 鄧軍, 李孝翠, 張躍彬. 云南甘蔗生產現(xiàn)狀與發(fā)展對策. 甘蔗糖業(yè), 2011(4): 77-80.

DENG J, LI X C, ZHANG Y B. The current situation and development countermeasures of sugarcane production in Yunnan. Sugarcane and Canesugar, 2011(4): 77-80. (in Chinese)

[7] 劉洪博, 應雄美, 毛鈞, 陸鑫, 蘇火生, 馬麗, 林秀琴, 蔡青. 11份割手密遺傳多樣性的 SSR分析. 植物遺傳資源學報, 2013, 14(3): 542-546.

LIU H B, YING X M, MAO J, LU X, SU H S, MA L, LIN X Q, CAI Q. Genetic diversity analysis for 11 clones of Saccharum spontaneum L.. Journal of Plant Genetic Resource, 2013, 14(3): 542-546. (inChinese)

[8] 蘇火生, 劉新龍, 毛鈞, 應雄美, 陸鑫, 馬麗, 范源洪, 蔡青. 割手密初級核心種質取樣策略研究. 湖南農業(yè)大學學報(自然科學版), 2011, 37(3): 253-259.

SU H S, LIU X L, MAO J, YING X M, LU X, MA L, FAN Y H, CAI Q. Sampling strategy of pre-core collection for Saccharum spontaneum. Journal of Hunan Agricultural University(Natural Science Edition), 2011, 37(3): 253-259. (in Chinese)

[9] 齊永文, 樊麗娜, 羅青文, 王勤南, 陳勇生, 黃忠興, 劉睿, 劉少謀, 鄧海華, 李奇?zhèn)? 甘蔗細莖野生種核心種質構建. 作物學報, 2013, 39(4): 649-656.

QI Y W, FAN L N, LUO Q W, WANG Q N, CHEN Y S, HUANG Z X, LIU R, LIU S M, DENG H H, LI Q W. Establishment of Saccharum spontaneum L. core collections. Acta Agronomica Sinica, 2013, 39(4): 649-656. (in Chinese)

[10] 陳能武, 楊榮仲, 吳才文, 黃久凱, 王貴華. 四川割手密(S. spontaneum)資源的雜交利用潛力研究. 甘蔗, 1996, 3(4): 1-7.

CHEN N W, YANG R Z, WU C W, HUANG J K, WANG G H. Studies on breeding potential of S. spontaneum of Sichuan province. Sugar Cane, 1996, 3(4): 1-7. (in Chinese)

[11] 許文花, 楊清輝. 甘蔗割手密無性系抗旱性鑒定. 亞熱帶農業(yè)研究, 2005, 1(2): 22-26.

XU W H, YANG Q H. Drought resistance of S. spontaneum L. clones. Subtropical Agriculture Research, 2005, 1(2): 22-26. (in Chinese)

[12] AMINATUN M, TARYONO, ENDANG S, PAUL H, SISMINDAR I. Tolerance of accessions of glagah (Saccharum spontaneum) to drought stress and their accumulation of proline. American Journal of Agricultural & Biological Science, 2013, 8(1): 1-11.

[13] NEUFELD H S, GRANTZ D A, MEINZER F C, GOLDSTEIN G, CRISOSTO G M, CRISOSTO C. Genotypic variability in vulnerability of leaf xylem to cavitation in water-stressed and well-irrigated sugarcane. Plant Physiology, 1992, 100(2): 1020-1028.

[14] 姚艷麗, 徐磊, 胡小文, 邢淑蓮, 劉洋. 割手密DREB2類轉錄因子的克隆與比較分析. 分子植物育種, 2016, 14(9): 2261-2267.

YAO Y L, XU L, HU X W, XING S L, LIU Y. The analysis of cloning and comparing of DREB2 group transcription factor genes from Saccharum spontaneum L.. Molecular Plant Breeding, 2016, 14(9): 2261-2267. (in Chinese)

[15] PARK J W, BENATTI T R, MARCONI T, YU Q, SOLISGRACIA N, MORA V, SILVA J A. Cold responsive gene expression profiling of sugarcane and Saccharum spontaneum with functional analysis of a cold inducible Saccharum Homolog of NOD26-like intrinsic protein to salt and water stress. PLoS ONE, 2015, 10(5): 51-56.

[16] 桃聯(lián)安, 楊李和, 經艷芬, 段慧芬, 董立華, 安汝東, 周清明, 朱建榮. 云南割手密血緣 F2代抗旱性隸屬函數(shù)法綜合評價. 西南農業(yè)學報, 2011, 24(5): 1676-1680.

TAO L A, YANG L H, JING Y F, DUAN H F, DONG L H, AN R D, ZHOU Q M, ZHU J R. Comprehensive evaluation of drought resistance for descendant F2of Saccharum spotaneum in yunnan by subjection function. Southwest China Journal of Agricultural Sciences, 2011, 24(5): 1676-1680. (in Chinese)

[17] ERLICH Y, MITRA P P, DELABASTIDE M, HANNON G J. Alta-Cyclic: A self-optimizing base caller for next-generation sequencing. Nature Methods, 2008, 5(8): 679-682.

[18] GRABHERR M G, HAAS B J, YASSOUR M, LEVIN J Z, THOMPSON D A, AMIT I, ADICONIS X, FAN L, RAYCHOWDHURY R, ZENG Q D, CHEN Z H, MAUCELI E, HACOHEN N, GNIRKE A, RHIND N, PALMA F, BIRREN B W, NUSBAUM C, LINDBLADTOH K, FRIEDMAN N, REGEV A. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nature Biotechnology, 2011, 29(7): 644-652.

[19] WANG Z, GERSTEIN M, SNYDER M. RNA-Seq: A revolutionary tool for Transcriptomics. Nature Reviews Genetics, 2009, 10: 57-63.

[20] MORTAZAVI A, WILLIAMS B A, MCCUE K, SCHAEFFER L, WORL B. Mapping and quantifying mammalian transcriptomes by RNA-Seq. Nature Methods, 2008, 5(7): 621-628.

[21] AUDIC S, CLAVERIE J M. The significance of digital gene expression profiles. Genome Research, 1997, 7(10): 986-995.

[22] KANEHISA M, ARAKI M, GOTO S, HATTORI M, HIRAKAWA M, ITOH M, KATAYAMA T, KAWASHIMA S, OKUDA S, TOKIMATSU T, YAMANISHI Y. KEGG for linking genomes to life and the environment. Nucleic Acids Research, 2008, 36(Suppl_1): 480-484.

[23] 程志遠, 李穆, 石瑩,何平, 何麗蓮, 李富生. 利用高通量測序分析甘蔗的抗旱響應機制. 分子植物育種, 2015(9): 2018-2028.

CHENG Z Y, LI M, SHI Y, HE P, HE L L, LI F S. Research of drought-response mechanism in sugarcane by high throughput sequencing-based digital gene expression profiling. Molecular Plant Breeding, 2015(9): 2018-2028. (in Chinese)

[24] 張騰國, 劉玉冰, 夏小慧. 植物 MAP激酶級聯(lián)途徑研究進展. 西北植物學報, 2008, 28(8): 1704-1714.

ZHANG T G, LIU Y B, XIA X H. RESEARCH ADVANCES ABOUT MAPK PATHWAYS IN PLANTS. Acta Botanica Boreali- Occidentalia Sinica, 2008, 28(8): 1704-1714. (in Chinese)

[25] XIONG L, YANG Y. Disease resistance and abiotic stress tolerance in rice are inversely modulated by an abscisic acid-inducible mitogenactivated protein kinase. The Plant Cell, 2003, 15(3): 745-759.

[26] AGRAWAL G K, AGRAWAL S J, IWAHASHI H, RAKWAL R. Novel rice MAP kinases OsMSRMK3 and OsWJUMK1 involved in encountering diverse environmental stresses and developmental regulation. Biochemical & Biophysical Research Communications, 2003, 300(3): 775-783.

[27] SANG Y M, ZHENG S Q, LI W Q, HUANG B R, WANG X M. Regulation of plant water loss by manipulating the expression of phospholipase D. The Plant Journal, 2001, 28 (2): 135-144.

[28] FRANK W, MUNNIK T, KERKMANN K, SALAMINI F, BARTELS D. Water deficit triggers phospholipase D activity in the resurrectionplant Craterostigma plantagineum. The Plant Cell, 2000, 12(1): 111-124.

[29] MUKHTAR N, IQBAL K, ANIS I, MALIK A. Sphingolipids from Conyza canadensis. Phytochemistry, 2002, 61(8):1005-1008.

[30] DAUM G. Lipid metabolism and membrane biogenesis. Springer, 2004, 6: 21-45.

[31] NG C K, CARR K, MCAINSH M R, POWELL B, HETHERINGTON A M. Drought-induced guard cell signal transduction involves sphingosine-1-phosphate. Nature, 2001, 410(6828): 596-599.

[32] NAPIER J A, MICHAELSON L V, DUNN T M. A new class of lipid desaturase central to sphingolipid biosynthesis and signalling. Trends in Plant Science, 2002, 7(11): 475-478.

[33] STRASSER R, BONDILI J S, SCHOBERER J, SVOBODA B, LIEBMINGER E, GL?SSL J, ALTMANN F, STEINKELLNER H, MACH L. Enzymatic properties and subcellular localization of Arabidopsis beta-N-acetyl- hexosaminidases. Plant Physiology, 2007, 145(1): 5-16.

[34] LEMBKE C G, NISHIYAMA M Y, SATO P M, ANDRADE R F D, SOUZA G M. Identification of sense and antisense transcripts regulated by drought in sugarcane. Plant Molecular Biology, 2012, 79(4/5): 461-477.

[35] MOHAMMADI M, KAV N N V, DEYHOLOS M K. Transcriptional profiling of hexaploid wheat (Triticum aestivum L.) roots identifies novel, dehydration-responsive genes. Plant Cell & Environment, 2007, 30(5): 630-645.

[36] 麻文俊. 楸樹優(yōu)良無性系 2-8苗期生理變化與基因表達對干旱脅迫的響應[D]. 北京: 中國林業(yè)科學研究院, 2013.

MA W J. Physiological changes and gene expression response of Catalpa bungei superior clone 2-8 seedlings to drought stress[D]. Beijing: Chinese Academy of Forestry, 2013. (in Chinese)

[37] MELAN M A, DONG X, ENDARA M E, DAVIS K R, AUSUBEL F M, PETERMAN T K. An Arabidopsis thaliana lipoxygenase gene can be induced by pathogens, abscisic acid, and methyl jasmonate. Plant Physiology, 1993, 101(2): 441-450.

[38] GIGON A, MATOS A R, LAFFRAY D, ZUILY-FODIL Y, PHAM-THI A T. Effect of drought stress on lipid metabolism in the leaves of Arabidopsis thaliana (ecotype Columbia). Annals of Botany, 2004, 94(3): 345-351.

[39] 劉南南, 江玲, 張文偉, 劉玲瓏, 翟虎渠, 萬建民. 水稻種胚 LOX3基因在逆境脅迫中的作用. 中國水稻科學, 2008, 22(1):8-14.

LIU N N, JIANG L, ZHANG W W, LIU L L, ZHAI H Q, WAN J M. Role of embryo LOX3 gene under adversity stress in rice (Oryza sativa). Chinese Journal of Rice Science, 2008, 22(1): 8-14. (in Chinese)

[40] PONTE-SUCRE A. Availability and applications of ATP-binding cassette (ABC) transporter blockers. Applied Microbiology & Biotechnology, 2007, 76(2): 279-286.

[41] DUDLER R, HERTIG C. Structure of an mdr-like gene from Arabidopsis thaliana evolutionary implications. Journal of Biological Chemistry, 1992, 267(9): 5882-5888.

（責任編輯 李莉）

Transcriptome Difference Analysis of Saccharum spontaneum Roots in Response to Drought Stress

LIU HongBo, LIU XinLong, SU HuoSheng, LU Xin, XU ChaoHua, MAO Jun, LIN XiuQin, LI ChunJia, LI XuJuan, ZI QiuYan
(Sugarcane Research Institute, Yunnan Academy of Agricultural Sciences/Yunnan Key Laboratory of Sugarcane Genetic Improvement, Kaiyuan 661699, Yunnan)

【Objective】The objective of this study is to understand the inner molecular mechanisms of Saccharum spontaneum clone named Yunnan 82-114 in response to water stress, and mine these genes closely related to drought tolerance, improve the utilization efficiency of Saccharum spontaneum in sugarcane drought-resistant breeding program. 【Method】lluminaHiSeqTM 4000,a high-through transcriptome sequencing technology, was applied to obtain the transcriptome differential expression data of Yunnan82-114 roots under drought stress treatment for 24, 48, and 72h. These assembled unigenes above were compared with database of Swiss-Prot, Nr, KOG, Pfam, and KEGG, respectively, then the abundance of gene expression among different samples were screened according to transcriptome data by using RPKM method, and the differentially expressed genes among the treated samples were estimated by referring to the standard of FDR≤0.05 & |log2 fold change|≥1. Function and pathway of those different expression genes were also investigated using Gene Ontology（GO）database and KEGG pathway database.【Result】Total 134 724, 130 368, 133 564, and 131 321 expressing genes were got, respectively, from the untreated control and three treated samples(24 h, 48 h and 72 h). Compared the control with the treated samples, about 3 061 (1 302), 2 304 (2 841) and 3 236 (2525) genes whose expression was significantly up-regulated (or down-regulated) were detected, respectively. When the FDR value was set to be 0 as a selection criterion, some genes performing extreme significantly different expression between control and treated samples were acquired, they mainly involved in transcriptional activation, water transport, DNA binding, ATP binding, membrane/ transmembrane transport, and defense response. GO enrichment analysis showed that some differences in the biological approach categories existed among the treated samples, for the samples treated for 24 h and 48h, the most enriched category was DNA-dependent transcription and transcription factor, then to the samples treated for 72 h, the enriched category was the translation-related genes. At the aspects of cellular components and molecular functions, the enriched category was inner membrane, nucleus, and ATP binding related genes in the three treated samples. KEGG enrichment analysis illustrated that 2 248 differentially expressed genes involved in 107 pathways in 24 h-treated sample, 2 114 in 130 in 48h-treated sample, and 2 392 in 144 in 72 h-treated sample. Finally, some genes related to abiotic stress and resistance were screened via KEGG significance enrichment analysis, which include glycosylsphingolipid biosynthesis, MAP kinase signal transduction, and ABC translocator.【Conclusion】Expression pattern and function information of differentially expressed genes between three water-stressed samples and a control were analyzed. Some genes like extracellular signal-regulated kinase, phospholipase D, hexosaminidase, and ATP-binding cassette B-1 were found to involve in response to water stress for Yunnan82-114. About 70 significant up-regulated expression genes at transcription level were obtained, but there were 11 genes exhibiting down-regulated expression, this implied that these genes induced or suppressed fiercely by drought stress have obvious relationships with drought tolerance of Yunnan82-114.

Saccharum spontaneum; drought; transcriptomes; differentially expressed genes

2016-11-24；接受日期：2017-01-23

國家自然科學基金（31540084）、云南省應用基礎研究計劃（2013FZ153）、云南省中青年學術技術帶頭人后備人才（2014HB038）

聯(lián)系方式：劉洪博，E-mail：liuhongbo1982@126..com。劉新龍，E-mail：lxlgood868@163.com。劉洪博和劉新龍為同等貢獻作者。通信作者蘇火生，

E-mail：shs304@163.com