干旱脅迫下枇杷葉片的轉(zhuǎn)錄組分析

王 丹，龔榮高

(四川農(nóng)業(yè)大學 園藝學院，四川 成都 611130)

枇杷(EriobotryajaponicaLindl.)是亞熱帶果樹中比較重要的一種，喜溫暖氣候和肥水濕潤的土壤，在我國有著豐富的種植資源和優(yōu)良品種[1]。枇杷對干旱有一定的耐受性，但在生長過程中易發(fā)生花期、幼果期和果實成熟期旱害，導致枇杷的坐果率和產(chǎn)量嚴重降低[2-3]。而干旱是不可避免的，它是影響植物生長發(fā)育最重要的非生物因素之一[4]。關(guān)于枇杷干旱的研究在生理生化方面較多[5]，羅華建等[2]用盆栽的方法證明：干旱抑制枇杷的生長和光合作用；也有研究表明，輕度的干旱脅迫有助于誘導枇杷花芽分化，使果實成熟期提前，對枇杷的生長有促進作用[6-7]。但由于枇杷基因組背景信息缺乏，尚未見有關(guān)枇杷葉片轉(zhuǎn)錄組測序等方面的研究報道，其抗逆機理的深入研究及抗旱分子標記開發(fā)等發(fā)展相對滯后。

轉(zhuǎn)錄組測序(RNA sequencing)利用第2代高通量測序技術(shù)進行cDNA測序，能全面快速地獲得某一狀態(tài)下的全部轉(zhuǎn)錄本信息，可以用于新轉(zhuǎn)錄本預(yù)測[8]。對于缺乏基因組信息的物種而言，采用轉(zhuǎn)錄組測序技術(shù)可獲得大量的轉(zhuǎn)錄本信息，從中發(fā)掘重要的功能基因，是揭示植物優(yōu)良性狀的重要研究手段[9]。轉(zhuǎn)錄組技術(shù)還可以研究不同器官、不同環(huán)境脅迫下基因表達的差異，過去10年里，上千個基因和數(shù)十條代謝和信號通路被證明與植物的耐旱性有關(guān)，擬南芥[10]、水稻[11]和短柄草[12]中的許多基因在轉(zhuǎn)錄水平都被證實參與干旱的調(diào)控。而參與植物光合作用、ABA合成、信號轉(zhuǎn)導、滲透物質(zhì)的生物合成和轉(zhuǎn)錄因子(WRKY、bZIP家族等)調(diào)控等的基因在干旱脅迫下表達異常，表明植物對水分刺激應(yīng)答的過程受眾多分子和細胞通路的調(diào)控[13]。但在其他植物中，類似的研究非常有限，植物對干旱耐性的分子機制在很大程度上還是未知的。

本研究利用Illumina HiSeq高通量測序技術(shù)對干旱前后的枇杷葉片轉(zhuǎn)錄組進行研究，結(jié)合生物信息學方法對所獲得的差異表達基因(DEGs，Differently expressed genes)進行功能的注釋、分類和代謝途徑分析，在轉(zhuǎn)錄水平上研究干旱脅迫下枇杷關(guān)鍵基因的表達，挖掘與干旱相關(guān)的基因，為枇杷抗旱相關(guān)基因的克隆以及功能分析和抗旱育種等研究提供理論基礎(chǔ)。

1 材料和方法

1.1 試驗材料

試驗材料為四川龍泉驛區(qū)長勢良好的3年生大五星枇杷嫁接苗，定植于口徑40 cm、高30 cm的花盆中，每盆種植1株。所用土壤為菜園土和沙壤土適量混合，并加入10%腐熟的有機肥和0.1%的復合肥，后期處理過程中不再施肥。待馴化后選取高度、長勢一致的材料進行干旱脅迫控制性試驗。試驗于2015年6-9月在四川農(nóng)業(yè)大學農(nóng)場進行。試驗開始時將材料移至大棚，每盆充分灌溉，使土壤含水量一致，停止灌溉后采用自然耗水進行干旱脅迫處理。

1.2 試驗方法

1.2.1 試驗處理 以土壤含水量(占田間持水量的百分數(shù))設(shè)置正常供水(對照，田間持水量的61%～65%)、輕度干旱脅迫(田間持水量的50%～54%)、中度干旱脅迫(田間持水量的39%～43%)和重度干旱脅迫(田間持水量的28%～32%)4個水分處理。將輕度、中度和重度3個脅迫的葉片取樣后混合和對照葉片液氮速凍后于-70 ℃貯存，用于RNA的提取。

1.2.2 總RNA的提取與測序 采用TRIzol法提取總RNA，用Agilent 2100 Technologies檢測樣品RNA的完整性和濃度，用帶有Oligo(dT)的磁珠富集mRNA，再加入緩沖液將其隨機打斷成片段，以mRNA片段為模板用隨機引物合成cDNA第一鏈，再加入緩沖液、dNTPs、RNase H和DNA Polymerase Ⅰ合成cDNA第二鏈，將洗脫純化后的雙鏈cDNA進行末端修復、加堿基A、加測序接頭處理，然后經(jīng)瓊脂糖凝膠電泳回收目的大小片段最后再進行PCR擴增，完成整個文庫制備后用Illumina HiSeq進行測序，測序策略為PE100。

1.2.3 測序數(shù)據(jù)的組裝和功能注釋 采用Trinity(20140717版)[14]將除去接頭序列、兩端低質(zhì)量序列和低度復雜序列的高質(zhì)量序列數(shù)據(jù)(Clean reads)組裝出全長轉(zhuǎn)錄本。對組裝出來的轉(zhuǎn)錄本利用TransDecoder(20140717版)鑒定編碼區(qū)域，預(yù)測開放閱讀框ORF(Opening reading fragment)。針對預(yù)測出來的ORF和Contig序列用Trinotate(20131110版)對其進行功能注釋，包括與同源性搜索(NCBI-Blast)、蛋白質(zhì)結(jié)構(gòu)域鑒定(HMMER/PFAM)、蛋白質(zhì)信號預(yù)測(SingalP/TmHMM)、Uniprot(Universal Protein：Swiss-Prot、TrEMBL和PIR-PSD)、Blast2GO軟件進行的GO(Gene Ontology)、COG(Cluster of Orthologous Groups)和KEGG(Kyoto Encyclopedia of Genes and Genomes)進行序列比對，獲得基因的功能注釋信息。

1.2.4 差異基因篩選及富集分析 基因表達量的計算使用標準化轉(zhuǎn)錄本數(shù)(RPKM)法，選取|log2Ratio|≥1和q<0.05的基因作為DEGs。應(yīng)用超幾何檢驗，以q<0.05 為閾值，篩選出在DEGs中顯著富集的條目和通路，分析DEGs行使的主要生物學功能和代謝通路。

1.2.5 實時熒光定量qRT-PCR 利用TRIzol試劑提取RNA后，采用Aidlab公司反轉(zhuǎn)錄試劑盒(TUREscript 1st Stand cDNA SYNTHESIS Kit)進行反轉(zhuǎn)錄操作，選取測序結(jié)果中表達數(shù)量較多且在各條目和通路中顯著富集的8個差異基因，依據(jù)基因的核酸序列設(shè)計熒光定量引物，采用20 μL反應(yīng)體系，對熒光定量PCR的引物濃度、退火溫度進行優(yōu)化，得到最佳的熒光定量PCR，每個樣品3個重復，引物見表1。采用2-ΔΔCT法計算各基因相對表達量。

表1 熒光定量PCR引物

2 結(jié)果與分析

2.1 測序結(jié)果和注釋數(shù)據(jù)

枇杷葉片通過新一代高通量測序技術(shù)進行轉(zhuǎn)錄組測序，共獲得114 021 518個原始reads，過濾掉接頭污染、低質(zhì)量、含N比例大于50%的序列后，得到90 517 644個干凈reads。堿基Q30為92.5%，表明測序的數(shù)據(jù)量和質(zhì)量都比較高。利用Trinity將干凈reads進行從頭組裝形成轉(zhuǎn)錄本，共發(fā)現(xiàn)88 530條轉(zhuǎn)錄本，平均長度740.64 bp，N50長度為1 189 bp，其中最長的轉(zhuǎn)錄本長度為10 189 bp，最短的長度為201 bp。結(jié)果表明，長度在200～400 bp的轉(zhuǎn)錄本數(shù)目最多，為41 141條，占所有轉(zhuǎn)錄本的46.47%(表2)。

組裝出的轉(zhuǎn)錄本共預(yù)測到41 748條ORF，占所有轉(zhuǎn)錄本總數(shù)的47.16%，ORF平均長度為891.54 bp，N50長度為1 137 bp，長度400～600 bp所占的比例最高。將Contigs和ORF序列注釋到Uniprot、GO和KEGG數(shù)據(jù)庫中，對注釋到每個庫中的數(shù)目進行統(tǒng)計，BlastX數(shù)據(jù)庫注釋的最多，為22 039條，占所有轉(zhuǎn)錄本的24.89%；BlastP數(shù)據(jù)庫注釋的數(shù)量為16 352條；KEGG數(shù)據(jù)庫注釋的數(shù)量為5 484條，占6.19%；注釋到GO數(shù)據(jù)庫的有20 912條。

表2 枇杷葉片轉(zhuǎn)錄組Transcript和ORF長度統(tǒng)計

2.2 COG功能分類

COG分類結(jié)果表明，13 285條轉(zhuǎn)錄本具有功能注釋信息且被注釋到24個功能類別。涵蓋基因最多的是一般功能預(yù)測(18.14%)，其次是翻譯后修飾、蛋白質(zhì)折疊和分子伴侶(10.65%)，翻譯、核糖體結(jié)構(gòu)和生物合成(10.64%)，復制、重組和修復(9.73%)，碳水化合物運輸與代謝(7.19%)，信號轉(zhuǎn)導機制(6.86%)和能量產(chǎn)生與轉(zhuǎn)化(6.57%)。所有類別中，細胞運動(0.11%)和核酸結(jié)構(gòu)(0.01%)的數(shù)量最少(圖1)。

A.加工與修飾RNA ；B.染色體結(jié)構(gòu)與變化；C.能量產(chǎn)生與轉(zhuǎn)化；D.細胞周期調(diào)控與分裂、染色體重排；E.氨基酸運輸與代謝；F.核苷酸運輸與代謝；G.碳水化合物運輸與代謝；H.輔酶運輸與代謝；I.脂質(zhì)運輸與代謝；J.翻譯、核糖體結(jié)構(gòu)與生物合成；K.轉(zhuǎn)錄；L.復制、重組與修復；M.細胞壁、膜生物合成Cell；N.細胞運動；O.翻譯后修飾、蛋白質(zhì)折疊和分子伴侶；P.無機離子運輸與代謝；Q.次生代謝物合成、運輸和代謝；R.一般功能預(yù)測；S.未知功能；T.信號轉(zhuǎn)導機制；U.細胞內(nèi)轉(zhuǎn)運、分泌和囊泡運輸；V.防御機制；Y.核酸結(jié)構(gòu)；Z.細胞骨架。

2.3 差異基因的GO分類和富集分析

通過計算基因的RPKM值，將兩樣品間RPKM的log2倍數(shù)絕對值大于1和q<0.05的基因作為差異基因，共篩選出25 197個差異表達的基因，其中上調(diào)基因13 482個，下調(diào)基因11 715個。54個GO分類注釋結(jié)果(圖2)表明：生物學過程分類中，注釋數(shù)量最多的是細胞進程(上調(diào)基因為8.23%，下調(diào)基因為8.84%)和代謝途徑(上調(diào)基因7.40%，下調(diào)基因8.0%)。在細胞組分中，細胞和細胞部分擁有的基因數(shù)最多(上調(diào)基因均為10.08%，下調(diào)基因均為11.94%)。分子功能中結(jié)合作用(上調(diào)基因分別為7.55%，下調(diào)基因為9.51%)的數(shù)量最多。

通過顯著富集計算，生物學過程中有345條GO條目顯著富集，細胞組分有87條，分子功能中有173條。其中分子功能中富集差異基因最多的前3條GO條目分別是催化活性、轉(zhuǎn)移酶活性和水解酶活性(表3)。

表3 分子功能中DEGs富集最多的前10 類GO注釋

2.4 代謝通路分析

KEGG分析結(jié)果表明，1 618個DEGs在291個標準代謝通路中得到注釋。對KEGG中每個代謝通路用超幾何檢驗進行富集分析，找出DEGs顯著性富集的30條代謝通路(表4)，主要有：代謝途徑、次生代謝物的二次代謝、不同生境微生物代謝作用和核糖體。這些途徑為研究枇杷干旱脅迫的具體過程、功能和途徑提供了珍貴的資源。其中，與激素代謝相關(guān)的3條代謝途徑：雙萜類、油菜素類固醇和類胡蘿卜素生物合成途徑中的差異基因呈現(xiàn)出較為一致的表達。雙萜類生物合成途徑的13個基因中有6個貝殼杉烯酸氧化酶和4個貝殼杉烯氧化酶基因均下調(diào)表達。類胡蘿卜素生物合成途徑共14個DEGs，其中有2個9-順式環(huán)氧類胡蘿卜素雙加氧酶基因和10個上調(diào)表達的脫落酸(ABA)-8′羥化酶基因。油菜素內(nèi)固醇生物合成中DEGs共14個，12個均表達下調(diào)，其中有9個基因為CYP734A1。

圖2 干旱脅迫下枇杷葉片DEGs的GO分類

表4 DEGs富集的30條KEGG通路

2.5 差異基因及其實時熒光定量驗證

通過功能分析，除了部分未知功能的基因外，可能與干旱脅迫響應(yīng)相關(guān)的基因主要有信號轉(zhuǎn)導相關(guān)基因：促分裂原活化蛋白激酶(MAPK)、富含亮氨酸重復序列類受體蛋白激酶(LRR-RLKs)、G-type-S受體蛋白、ERF、MYB、WRKY轉(zhuǎn)錄因子和熱激蛋白等；抗氧化酶基因：過氧化氫酶、谷胱甘肽過氧化物酶、過氧化物酶；滲透調(diào)節(jié)物質(zhì)合成相關(guān)基因：海藻糖-6-磷酸合酶、乙醛脫氫酶、谷氨酸脫氫酶；激素合成相關(guān)基因：細胞色素P450家族、貝殼杉烯酸氧化酶、貝殼杉烯氧化酶、黃烷酮醇-4-還原酶、三角狀五肽重復(PPR)蛋白等。

選擇了差異表達數(shù)量豐富且在GO分類和KEGG代謝通路中都顯著富集的8個差異基因(表1)進行實時熒光定量分析，驗證這8個基因在葉片中的表達情況(圖3)，其中的6個基因c1047(過氧化物酶-上調(diào))、c65100(蛋白激酶byr2-上調(diào))、c56685(絲氨酸蘇氨酸蛋白激酶-上調(diào))、c68887(脫落酸8′-羥化酶-上調(diào))、c22815(細胞色素P450 734A1-下調(diào))和c68870(吲哚-3-乙酰乙酸合酶-上調(diào))的表達結(jié)果與轉(zhuǎn)錄組測序一致，這也證明了測序方法的正確性。

圖3 八個差異表達基因的定量PCR表達

3 討論

抗旱性一個復雜的多基因特征[15]，高通量測序能夠為植物基因表達的全面分析提供合理的數(shù)據(jù)資源，尤其適合沒有參考基因組信息的物種[16]。本研究采用高通量測序共獲得枇杷葉片轉(zhuǎn)錄本88 530條，預(yù)測出41 748條ORF序列，將測序結(jié)果與Uniprot、GO和KEGG等數(shù)據(jù)庫進行比對。通過篩選，得到25 197個DEGs(P<0.01)，然而有大量的DEGs未得到注釋且未被報道過，這很大程度上豐富了枇杷的基因信息資源庫[17-18]。注釋結(jié)果表明，大部分差異表達的基因參與了如催化酶、轉(zhuǎn)移酶和水解酶活性的調(diào)節(jié)，大分子物質(zhì)的代謝，滲透調(diào)節(jié)物質(zhì)合成，信號轉(zhuǎn)導和激素水平的調(diào)節(jié)。

受體蛋白激酶(RLKs)是一組膜蛋白，盡管其功能和詳細的分子機制還不清楚，但RLK基因轉(zhuǎn)錄變化表明他們可能在干旱中起到重要的作用[19]。有研究表明，蛋白激酶信號基因在干旱時被誘導[20]。水稻受體蛋白激酶4在脅迫中被發(fā)現(xiàn)上調(diào)表達[11]。本研究中的絲氨酸/蘇氨酸蛋白激酶、LRR-RLKs等大部分基因均不同程度的下調(diào)，MAPK則表現(xiàn)為上調(diào)，說明這幾類蛋白激酶對干旱脅迫應(yīng)答調(diào)控起不同的作用。蔡國華[21]的研究結(jié)果表明，在干旱條件下，過表達玉米MAPKK在擬南芥中通過ROS和依賴ABA途徑增強了對干旱脅迫的抗性，另有研究表明，一些RLKs，如擬南芥RPK1在干旱脅迫時與ABA和BR(油菜素內(nèi)固醇)信號通路有聯(lián)系，而不同受體蛋白類型在信號通路中的作用還有待進一步探討[22]。

參與脅迫信號傳遞網(wǎng)絡(luò)的轉(zhuǎn)錄因子，尤其是ERF、WRKY、MYB等家族是信號轉(zhuǎn)導通路上的關(guān)鍵調(diào)控因子，對干旱脅迫響應(yīng)和水分刺激應(yīng)答發(fā)揮了重要調(diào)節(jié)作用[23]。ERF 類轉(zhuǎn)錄因子是AP2/EREBP 轉(zhuǎn)錄因子家族的一個亞家族，通過參與乙烯、ABA、茉莉酸等多種信號轉(zhuǎn)導途徑，在非生物脅迫應(yīng)答調(diào)控中發(fā)揮重要的作用，一些ERFs在轉(zhuǎn)基因植物中能加強對干旱的耐受性，過表達的ERF類轉(zhuǎn)錄因子WXP1、WXP2可以增加葉表皮蠟質(zhì)的積累，顯著減少葉片脫水、提高植株對干旱脅迫的耐受性[24]。也有研究表明，在耐旱玉米基因型中的轉(zhuǎn)錄因子MYB33和MYB101能增強玉米的耐旱性[25]，在柑橘胚珠的研究中與信號傳遞相關(guān)的基因WRKY被驗證在脅迫發(fā)生時表達下調(diào)[26]。本試驗中ERF、WRKY、MYB是干旱脅迫處理下數(shù)量豐富、表達量變化最顯著的幾類轉(zhuǎn)錄因子，下調(diào)顯著的主要是ERF類轉(zhuǎn)錄因子，上調(diào)顯著的是MYB相關(guān)家族編碼基因，說明這些轉(zhuǎn)錄因子對干旱脅迫的應(yīng)答顯著，這與Chen等[27]的研究結(jié)果一致，原因可能是通過上調(diào)基因來實現(xiàn)葉片對水分刺激的應(yīng)答調(diào)節(jié)，通過下調(diào)基因數(shù)量的平衡來實現(xiàn)葉片對干旱脅迫的響應(yīng)[28]。

ABA的生物合成和分解代謝率決定ABA的積累以及對脅迫的反應(yīng)強度，盡管有各種各樣的途徑合成ABA，但類胡蘿卜素生物合成途徑在被子植物中是唯一已經(jīng)確定的途徑[29]。本研究中類胡蘿卜素生物合成通路中的10個ABA 8′-羥化酶上調(diào)表達，而ABA 8′-羥化酶屬于細胞色素P450的CYP707A家族，它是合成植物激素ABA的前體[30]且CYP707A8′-羥基化基因在擬南芥的研究中參與抗旱[31]、種子休眠和萌發(fā)[32]。ABA 8′-羥化酶可能參與合成類胡蘿卜素，并積累額外的ABA以響應(yīng)脅迫。ABA與赤霉素(GA)有拮抗作用，在干旱脅迫下，激素的比值會發(fā)生不同程度的改變，有研究表明，GA3/ABA值能體現(xiàn)玉米萌發(fā)期適應(yīng)脅迫的能力[33]。在本研究中，雙萜類生物合成途徑中的GA合成上游的2種關(guān)鍵酶均下調(diào)表達：貝殼杉烯酸氧化酶和貝殼杉烯氧化酶，這可能與ABA 8′-羥化酶基因的上調(diào)表達相互平衡，共同響應(yīng)干旱，但二者的具體變化情況還不明確。油菜素內(nèi)酯(BRs)是一種植物甾醇類激素，有研究發(fā)現(xiàn)，油菜素內(nèi)酯能提高擬南芥幼苗忍耐低溫的能力，還能增加植物的環(huán)境抗逆性[34]。而本研究中，油菜素內(nèi)固醇生物合成途徑中下調(diào)基因中有9個都是細胞色素P450家族的CYP734A1基因，而擬南芥中的CYP734A1能夠使體內(nèi)BR失活[35]，說明此時植物可能通過降低CYP734A1基因的表達量來調(diào)節(jié)BR的活性以應(yīng)對干旱脅迫，關(guān)于CYP734A1基因和其他未知功能基因與植物響應(yīng)干旱之間的具體關(guān)系都還有待更深入的研究。