干旱脅迫下偏關苜蓿保護酶及轉錄組差異性分析

衛 凱，朱慧森，岑慧芳，程映杰

（山西農業大學 草業學院，山西 太谷 030801）

紫花苜蓿（Medicago sativa）是黃土高原地區的當家牧草，營養物質豐富，粗蛋白、維生素和礦物質含量豐富，氨基酸的組成比較齊全，適口性好[1]。其中，偏關苜蓿是山西省農業科學院畜牧獸醫研究所和偏關縣畜牧局選育登記的地方品種，具有抗寒、抗旱性強以及營養價值高等特點[2]。山西處于山丘地帶，屬中緯度東亞季風氣候，離海洋相對較遠，旱災時有發生[3]。這使得苜蓿產量通常很低，因而，了解紫花苜?？购档纳砗头肿舆z傳機制對于紫花苜蓿的可持續生產至關重要。

干旱脅迫下，植物細胞中的酶促防御系統可以通過維持活性氧等物質的產生和清除的平衡來緩解細胞受到的損傷[4]，植物保護酶系統主要包括超氧化物歧化酶（SOD）、過氧化氫酶（CAT）、過氧化物酶（POD）和抗壞血酸過氧化物酶（APX）等[5]。保護酶活性的主要研究方法在20世紀70年代得到發明和建立[6]。張新蘭[7]測定了不同品種苜蓿葉片在干旱脅迫下抗氧化酶活性的動態變化，結果表明，抗旱性強的品種在干旱脅迫下草產量更高，葉片CAT、POD和SOD活性較高。張翠梅等[8]研究了PEG模擬干旱脅迫下紫花苜蓿葉片的膜脂過氧化程度，結果表明，抗旱性強的苜蓿品種葉片的膜脂過氧化程度更低。李紅等[9]采用PEG對12種不同抗旱性的紫花苜蓿進行干旱脅迫，結果顯示，不同抗旱性品種的葉片保護酶活性不同。因此，SOD、POD和CAT等保護酶的活性大小可以作為苜?？购敌栽u價的指標。

植物在干旱條件下，能夠誘發全部組織水平的應答，不僅涉及形態學水平的變化調整，還包括細胞水平的應答[10]。轉錄組學（Transcriptome）研究可以反映在特定時間和特定環境下植物與正常條件下的差異基因表達情況以及所在通路，有助于了解植物適應逆境脅迫的機制，從而有助于發現抵抗和緩解逆境脅迫的相關基因[11]。植物感受干旱信號后，通過信號轉導激發一系列生理生化變化，在分子水平調控特定基因表達的過程尤為關鍵[12]，最早是1997年VELCULESCU等[13]提出在轉錄水平上進行的基因表達差異分析實際上就是進行轉錄組研究。MA等[14]對紫花苜蓿Dryland品種在干旱脅迫下的轉錄組學分析，篩選出1 690個上調基因和3 827個下調基因，富集分析發現，蔗糖合酶、胺氧化酶和脂酰輔酶A還原酶相關基因在紫花苜蓿響應干旱脅迫的過程中具有重要研究價值。劉佳月[15]在干旱脅迫下紫花苜蓿中草5號的研究發現，差異表達基因在苯丙烷生物合成和谷胱甘肽代謝通路顯著富集，其數量在與代謝途徑相關的方面注釋最多，在次生代謝物的生物合成上也有相對較多的注釋數量。對甘農三號品種抗旱性的分析發現，紫花苜?？梢酝ㄟ^調控與脂質代謝、氨基酸代謝、碳水化合物代謝和次級代謝等通路相關的基因表達來響應干旱脅迫，增強抗旱能力[16]。

基于前期開展的偏關苜蓿種子萌發期的耐受性、顯微結構的響應特征以及關鍵酶的活性等研究工作[17-20]，發現偏關苜蓿表現出的抗旱潛力還有待挖掘，尤其在分子水平調控特定基因表達的過程[4]。本試驗以偏關苜蓿為研究材料，分析其在干旱脅迫和正常條件下保護酶及轉錄組的差異性，挖掘與抗旱相關的基因，為培育抗旱性強、性狀優良的苜蓿新種質提供理論依據。

1 材料和方法

1.1 試驗材料

供試材料為籽粒飽滿、外形完整、大小均勻的偏關苜蓿種子，由山西農業大學草業學院牧草種質資源庫提供。

1.2 試驗方法

試驗在山西農業大學草業學院實驗室中進行。將偏關苜蓿種子播種于營養缽，置于14 h光照、10 h黑暗的培養箱中，晝夜溫度分別為25、20?C，用Hoagland營養液每日定時澆灌。30 d后挑選生長情況一致的偏關苜蓿，在20%PEG溶液（用Hoagland營養液配制）中進行72 h的脅迫處理。未進行干旱處理的材料為對照組，脅迫72 h的為干旱組，采集偏關苜蓿葉片，用于SOD、POD、CAT、APX活性的測定及RNA提取。

1.3 保護酶活性的測定

SOD活性測定參照文獻[21]進行，采用NBT光還原法測定其在560 nm波長處的吸光值，酶活性以抑制3 mL NBT反應液光化學反應的50%所需酶量為1個活性單位；POD活性測定采用愈創木酚顯色法[21]，在分光光度計上測量470 nm下5 min時的吸光值，以每分鐘每克樣品鮮質量吸光值的變化為1.0表示1個活性單位；CAT活性測定參照高俊鳳[22]中的H2O2分解法，以1 min內A240減少0.1的酶量記為1個酶活性單位，共計時4 min；APX活性采用高俊鳳[22]的方法測定，在20℃下測定10～30 s內A290的變化，計算單位時間內抗壞血酸減少量及酶活性。

1.4 RNA提取及高通量測序

試驗采用TRIzol法提取RNA。檢測合格的樣品采用Illumina Hiseq TM 2500進行高通量測序，將得到的原始數據進行過濾，再利用TRINITY[23]軟件進行組裝，采用CD-HIT[24]方法選擇相似類中最長的作為唯一基因，得到基因的序列集合（Unigene）。

1.5 基因功能注釋

將Unigene序列與NR（NCBI non-redundant protein sequences）、NT（NCBInucleotidesequences）、COG（Clusters of Orthologous Groups of proteins）、Swiss-Prot、KEGG（Kyoto Encyclopedia of Genes and Genomes）、GO（Gene Ontology）六大數據庫進行基因功能注釋[25]。

1.6 差異表達基因分析

基因表達量的計算使用RPKM法[26]（Reads per kilobase million）。差異基因的篩選條件為：FDR＜0.05和|log2FC≥1|。GO分 析 的 方 法 是GOseq[27]。利用KEGG數據庫進行通路（Pathway）分析。

1.7 數據分析

試驗數據采用Microsoft Excel 2010和SAS 9.0統計分析軟件進行處理。

2 結果與分析

2.1 干旱脅迫下偏關苜蓿SOD、POD、CAT和APX活性變化

從表1可以看出，偏關苜蓿葉片的干旱組SOD、POD、CAT和APX活性相對于對照組均顯著上升（P＜0.05），其中，干旱組的POD活性是對照組的4.87倍。

表1 對照組與干旱組葉片中保護酶活性比較Tab.1 Comparison of protective enzyme activities in leaves of the control and drought groups

2.2 測序數據統計

測序結果如表2所示，對照組和干旱組的Q20百分率分別達到98.09%和98.05%，GC百分率平均值分別為51.56%和51.08%，這些數據說明測序質量偏好。

表2 測序數據輸出質量情況Tab.2 The output quality of sequencing data

2.3 轉錄組的拼裝

過濾后的測序數據經拼裝后最終得到79 794條Unigene，總長度為81201 250 bp，平均長度為1 018bp，N50為1 599 bp（表3）。采用組裝序列的長度分布來衡量轉錄組組裝質量，可知序列滿足分析的要求。

表3 Unigene的長度及數量統計Tab.3 Unigene length and quantity statistics

2.4 基因功能注釋

由表4可知，成功注釋到56 911條Unigenes，注釋到6個數據庫的數量分別是53 963（94.82%）、50 745（89.17%）、36 588（64.29%）、35 212（61.87%）、23 866（41.94%）、40 252個（70.73%），NR數據庫的注釋百分比最高。

表4 Unigene在各大數據庫中的功能注釋情況Tab.4 Statistics of Unigene functional annotation accor ding to the major databases

利用NR數據庫進行同源物種比對，結果發現（圖1），與偏關苜蓿同源基因相似度最高的是高粱，相似率達38.60%；其次為玉米（Zea mays），相似率為26.10%；日本水稻（Oryza sativa）的相似率為14.30%；二穗短柄草（Brachypodium distachyon）的相似率為6.90%。

對所有的Unigenes進行SSR分析，共識別出7 868個SSR。其中，二堿基重復2 570個，三堿基重復4 931個，四堿基重復85個，五堿基重復169個，六堿基重復113個，各堿基的分布頻率如表5所示。

表5 SSR分布頻率Tab.5 SSR distr ibution fr equency

2.5 基因差異性表達分析

進行差異基因篩選，獲得5 867個差異表達基因，上調基因占總數的46.28%，下調基因占總數的53.72%。對基因的表達量制作火山圖進行分析，如圖2所示。

2.6 差異表達基因的GO分析

對處理組和干旱組之間的差異表達基因進行GO功能分析,結果如圖3所示，圖中從左往右依次為生物附著、生物調節、細胞組分組織或生物合成、細胞過程、發育過程、定位系統的建立、生長、免疫系統過程、定位、代謝過程、多元有機體過程、多細胞生物過程、生物過程負調控、生物過程正調控、生物過程的調節、繁殖、繁殖過程、應急響應、節律過程、信號傳導、單生物過程、細胞、細胞連接、細胞部分、細胞外基質、細胞外基質部分、胞外器、胞外區部分、高分子復合物、膜、膜部分、膜包圍腔、擬核、細胞器、細胞器部分、共質體、病毒體、病毒體部分、抗氧化活性、結合、催化活性、電子載體活性、酶調節活性、金屬伴侶活性、分子轉導活性、核酸結合轉錄因子、營養貯藏、蛋白質結合轉錄因子、受體活性、結構分子、轉運活性，注釋到51個類別，3個部分中數目最多的類別分別是代謝過程（2 145個）、細胞（2 785個）、催化活性（2 035個）；進行GO富集分析，選擇富集度高的前30個如圖4所示，3個部分中富集最多的類別分別是作用于糖基鍵水解酶活性、染色體和DNA包裝。

對照組與干旱組蛋白激酶、O-糖基化合物水解酶、過氧化物酶、作用于糖基鍵的水解酶、木葡聚糖內糖基轉移酶、抗氧化酶和蛋白磷酸酶活性等相關基因功能有明顯差異?；蚬δ茏⑨尡砻?，代謝過程、催化活性、金屬伴侶活性、抗氧化活性和酶調節活性等與干旱脅迫聯系緊密的生物學過程和分子功能具有明顯的響應變化。

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

偏關苜蓿共有35 215條Unigenes可在KEGG中得到注釋，占Unigene總數的61.88%，注釋到代謝途徑的數目最多，為9 106個。KEGG富集分析表明，有3 471個差異表達基因可注釋到122條代謝途徑。

統計分析可知，122條代謝途徑中最顯著的4個是次生代謝產物的生物合成、植物病原體相互作用、苯丙烷的生物合成和其他聚糖降解，最顯著的前20個通路相關的差異表達基因占總量的59.69%。由表6和圖5可知，信號轉導、物質合成與降解、代謝途徑和次生代謝等通路相關基因的數量最多。差異表達基因在植物激素信號轉導和次生代謝物的生物合成富集程度相對較高。

表6 差異表達基因富集程度排名前20的Pathway條目Tab.6 Top 20 enrichment pathway in differentially expressed genes

從次生代謝物的生物合成、苯丙烷生物合成和植物信號轉導等代謝通路中挑選出部分與MYC2轉錄因子、苯丙氨酸酶、類胡蘿卜素結合蛋白、光合作用蛋白、淀粉酶、蔗糖合酶有關的基因作為偏關苜蓿響應干旱脅迫的候選基因，主要有Unigene12274_All、Unigene10105_All、Unigene35365_All、CL7158.Contig1_All、CL 11234.Contig1_All、CL 12445.Conti g1_All、Unigene12274_All、CL4041.Contig1_All。

2.8 抗旱相關轉錄因子的預測

對材料中與轉錄調控相關的差異基因進行了功能注釋，共得到136個轉錄因子（Transcription factors，TFs）基因，被分為18個基因家族，其中，數目較多的是bHLH（29個）、WRKY（22個）、ERFs（18個）和MYB（11個）。篩選出表達差異倍數2以上的轉錄因子基因（表7），這些轉錄因子在植物響應干旱脅迫中可能發揮重要的作用。

表7 可能與抗旱相關的轉錄因子Tab.7 Possibly dr ought-r esistance-related transcr iption factors

3 討論

3.1 干旱脅迫對偏關苜蓿保護酶的影響

干旱脅迫下植物體內產生過量的活性氧，活性氧的氧化作用會對植物造成損傷，植物清除活性氧的能力強弱是抗旱性評價的重要生理指標[28]。SOD是一種廣泛存在于植物內的酶，可以清除活性氧自由基從而使植物免受傷害。OLGA等[29]在干旱脅迫下研究不同抗旱性玉米品種的SOD活性差異，結果表明，正常條件下不同抗性品種的SOD活性無明顯差異，但在干旱脅迫下抗旱性強的品種SOD活性顯著高于抗旱性差的品種，可能的原因是抗旱性強的品種保護酶活性強，SOD相關基因表達量更高。本研究結果表明，偏關苜蓿葉片SOD活性顯著高于未脅迫時，與韓剛等[30]研究結果一致，說明干旱脅迫下SOD活性提高有利于偏關苜蓿清除植物體內產生的過量活性氧，增強偏關苜蓿的干旱脅迫耐受性。

POD是酶促反應系統的保護酶，通過防御活性氧來緩解干旱脅迫對植物的傷害。植物感受干旱信號后會調控相關基因表達，多數POD的合成屬于誘導表達型，與植物信號轉導有關[31]。干旱處理后小麥幼苗根和葉中的POD活性顯著增強[32]，與本研究一致，偏關苜蓿在干旱脅迫下POD活性顯著升高，發揮其清除活性氧的生理功能。

CAT對H2O2具有較高的親和力，可以清除植物體內過多的H2O2，將H2O2分解為對細胞無害的水和分子氧。研究表明，干旱脅迫下植物體內CAT活性會發生改變，高粱CAT活性在干旱脅迫下顯著增加，證明CAT是高粱重要的保護性酶[33]。也有研究發現，小麥抗旱性與干旱脅迫下CAT基因表達量密切相關，CAT等保護酶活性較高的品種抗旱性更強[34]，與本試驗結果相一致，干旱脅迫下偏關苜蓿葉片CAT活性顯著高于未脅迫時，說明偏關苜?？梢栽鰪奀AT活性來響應干旱脅迫。

APX是植物體內清除H2O2的關鍵酶。研究表明，甘薯在干旱脅迫下APX活性增強，H2O2含量減少，從而削弱其對光合的抑制，保護葉綠體[35]。轉APX棗樹在干旱條件下，APX活性高于非轉基因植株，表明轉基因棗樹耐旱能力強于非轉基因植物[36]。本研究結果表明，偏關苜蓿葉片APX活性顯著高于未脅迫時，APX通過有效清除活性氧在增強偏關苜蓿的耐旱能力中起一定作用。

3.2 干旱脅迫下偏關苜蓿轉錄組的差異性

干旱脅迫下，植物細胞通過信號轉導來調控抗旱相關基因表達，使自身結構、生理生化以及物質代謝發生一系列變化。利用高通量測序技術可以研究偏關苜蓿干旱響應基因，對測序結果的統計和評估可知本次數據質量較好，本研究共得到79 794條Unigenes，平均長度1 018 bp，所得到的Unigenes數量低于黃花苜蓿的114 347條[37]，而高于中苜一號紫花苜蓿的41 734條[15]，造成這種差異的原因可能是品種和試驗環境不同。本研究所有的Unigenes在六大功能注釋數據庫中進行了注釋，其中在NR數據庫的注釋率最高，與劉佳月[15]對紫花苜蓿的研究結果相同。Unigenes在NR數據庫中的注釋率高達99.1%，結果發現與偏關苜蓿同源率最高的是高粱，其次是玉米和日本水稻。SSR標記有多態性高、穩定性好和不受環境影響等優點，是遺傳變異研究最好的標記之一，目前已在苜蓿種質鑒定、親緣關系、DNA遺傳圖譜構建、分子標記育種和遺傳多樣性等領域廣泛應用[38]。本研究識別的7 868個SSR標記可為偏關苜蓿遺傳圖譜構建和分子育種提供依據。

本研究中GO功能分析結果可知，分子功能顯著差異在作用于糖基鍵的水解酶活性、蛋白質異二聚化活性、O-糖基化合物水解酶活性等；細胞組分顯著差異在細胞質囊、微管、DNA蛋白復合物等；生物學過程顯著差異在DNA包裝、DNA-蛋白質復合體亞基結構、DNA構象變化等。說明偏關苜蓿通過復雜的基因調控網絡來響應干旱脅迫。白雪梅[39]對紫花苜蓿的研究結果表明，蛋白激酶基因是紫花苜蓿在響應干旱脅迫中的關鍵信號傳遞體，參與磷酸化和去磷酸化，在抗旱中起關鍵作用。本研究發現，對照組與干旱組在過氧化物酶活性、抗氧化酶活性和蛋白磷酸酶活性相關基因功能有顯著差異。其中，抗氧化酶和過氧化物酶的基因差異性表達可以反映偏關苜蓿對干旱的分子響應，抵抗活性氧誘導的氧化應激損傷。另外，有研究表明，蛋白磷酸酶活性相關基因參與了草地早熟禾對干旱脅迫的響應[40]。

KEGG分析發現，本研究中偏關苜蓿差異表達基因顯著富集在植物激素信號轉導與淀粉和蔗糖代謝通路，這表明植物感受干旱信號后會將信號傳遞給相關細胞器來引起抗旱基因表達，使淀粉和蔗糖代謝相關途徑向有利于滲透調節的方向激活或抑制。一般而言,Ca2+信號傳導、植物激素信號、蛋白激酶和轉錄因子是植物體內重要的信號傳遞體[41]，在本研究得到的差異表達基因中均有注釋。淀粉和蔗糖代謝能夠影響植物的生長及其脅迫響應，這表明偏關苜蓿受到干旱脅迫后通過形成一些糖類小分子化合物來進行滲透調節[42]，其原理是干旱脅迫下淀粉和蔗糖可以降解為易溶于水的葡萄糖，同等條件下單糖溶液滲透壓大于雙糖和多糖溶液，從而避免植物細胞損失更多的水分[43]。

能夠顯著響應干旱脅迫的基因可以作為苜蓿處在逆境時分子水平的應答特征，同時具有提高苜蓿耐旱和抗旱的潛質。本研究中干旱脅迫下偏關苜蓿蔗糖合成酶基因的下調以及淀粉水解酶基因的上調表達，一方面可以為植物響應干旱脅迫提供更直接的能量，另一方面可以防止細胞內滲透壓下降造成的失水。JACOBSEN等[44]研究也發現，隨著葉片水勢下降，大麥幼苗葉片淀粉酶活性增加。WANG等[45]對茶樹在干旱脅迫下轉錄組的表達進行分析，發現差異表達基因在黃酮類生物合成和苯丙烷類生物合成途徑大量富集。本研究中差異表達基因在苯丙烷生物合成中顯著富集，其中苯丙氨酸酶基因上調表達，表明了與苯丙烷生物合成途徑相關的關鍵基因在偏關苜蓿響應干旱脅迫中起重要作用。類胡蘿卜素相關基因表達的差異也是苜蓿品種抗旱差異的重要原因之一[16]。相關上述研究的差異表達基因有助于提高耐旱性，可作為與偏關苜蓿干旱脅迫相關的潛在抗旱基因。

植物處于干旱脅迫時轉錄因子起不可或缺的作用，本研究中bHLH、WRKY、ERFs和MYB家族基因在偏關苜蓿響應干旱脅迫過程意義重大。通常認為，WRKY、bHLH、NAC等轉錄因子家族與干旱脅迫有關[15]。MYB轉錄因子對干旱脅迫有負調控作用[46]，干旱脅迫誘導NAC轉錄因子上調[47]。LIU等[48]采用PEG脅迫小麥1、6 h后取葉片進行轉錄組測序，差異基因主要富集到FAR1、NAC、bZIP、bHLH、AP2/ERF、WRKY、Myb-related和Myb轉錄因子家族。對杜梨干旱脅迫的轉錄組分析結果顯示，數目最大的轉錄因子家族是MYB，有72個MYB基因差異表達，其次是bHLH，有57個bHLH基因響應干旱脅迫，45個WRKY基因在干旱的誘導下差異表達，39個NAC基因、15個HSF基因在干旱條件下差異表達[49]。對不抗旱香蕉品種和抗旱香蕉品種進行干旱脅迫后取樣測序，結果表明，與香蕉干旱脅迫相關的轉錄因子大多是MYB、NAC、WRKY、bHLH、bZIP、ERF和HSFs轉錄因子家族的成員[50]。

植物復雜的信號網絡可以使其迅速地感知外界環境變化，調控基因的表達，從而形成應對干旱脅迫的多種機制。WANG等[51]研究了過表達葡萄V vbHLH 1擬南芥與非轉基因擬南芥在干旱脅迫下的差異性，結果表明，VvbHLH 1的過表達增加了植株中SOD和POD的活性，H2O2和丙二醛含量顯著下降；與脯氨酸生物合成、黃酮生物合成和抗氧化酶相關的基因上調。這說明VvbHLH 1通過提高抗氧化酶活性和增加黃酮類化合物的積累增強了擬南芥對干旱脅迫的適應能力。在本研究中也發現，黃酮和黃酮醇的生物合成相關基因和抗氧化相關基因顯著表達，并且挖掘出bHLH轉錄因子，這表明在偏關苜蓿中可能也存在類似的抗旱過程。有研究發現，在干旱脅迫下過表達甘薯IbWRKY 2的擬南芥中SOD活性較高，并且MDA和H2O2的含量較低，與脫落酸信號轉導途徑、抗氧化酶和脯氨酸生物合成有關的基因差異顯著[52]。這表明WRKY轉錄因子具有提高轉基因擬南芥耐旱性的作用，本研究中偏關苜蓿的WRKY轉錄因子具有提高苜?？购敌缘臐摿?。有研究表明，苯丙烷類代謝主要受MYB轉錄因子家族的調控[53]；在本研究中也發現，苯丙烷生物合成以及黃酮和黃酮醇的生物合成途徑相關基因顯著表達，所挖掘的偏關苜蓿bHLH及MYB轉錄因子可能通過調控保護酶活性相關基因和注釋到上述代謝通路的基因來提高抗旱性?；谄P苜蓿響應干旱脅迫的差異表達基因，根據功能注釋結果可知，偏關苜?？赡芡ㄟ^調控抗氧化活性和酶調節活性相關的基因以及bHLH等轉錄因子，使SOD、POD、CAT及APX酶活性上升，減少干旱脅迫造成的損傷。本研究中受干旱脅迫誘導而表達的基因，可在今后苜蓿耐旱育種中作進一步研究。

4 結論

本研究以偏關苜蓿為材料，以正常條件和干旱脅迫下的偏關苜蓿葉片進行保護酶活性測定和轉錄組分析，結果表明，偏關苜?？梢酝ㄟ^調節抗氧化酶活性來響應干旱脅迫，轉錄水平上總共獲得56 911條基因序列，識別出7 868個SSR位點。通過差異表達基因的分析得出，Unigene12274_All、Unigene10105_All、Unigene35365_All、CL7158.Contig1_All、CL 11234.Contig1_All、CL 12445.Con tig1_All、Unigene12274_All、CL 4041.Contig1_All以及bHLH 35、ERF25、PCF7和RF2b可作為偏關苜蓿潛在抗旱基因，為后續苜?？购捣肿訖C制研究提供理論依據。