洪森榮，朱盈盈，李紫瑩，胡明艷，歐陽克蕙

（1上饒師范學院生命科學學院，江西上饒 334001；2上饒農(nóng)業(yè)技術(shù)創(chuàng)新研究院，江西上饒 334001；3上饒市藥食同源植物資源保護與利用重點實驗室，江西上饒 334001；4上饒市薯芋類作物種質(zhì)保存與利用重點實驗室，江西上饒 334001；5江西農(nóng)業(yè)大學動物科學技術(shù)學院，南昌 330045）

0 引言

金花菜（Medicago polymorpha L.又稱南苜蓿）為江浙地區(qū)廣泛推廣種植的苜蓿屬一年生豆科植物，一年可刈割多次，鮮草產(chǎn)量高，營養(yǎng)豐富[1]，藥用價值、綠肥價值、飼用價值和菜用價值高，抗逆性強，固氮效率高，可作為無公害蔬菜生產(chǎn)[2]，廣泛用于草田輪作[3]。紫花苜蓿(Medicago sativa L.)是優(yōu)良的多年生豆科植物，被譽為牧草之王，產(chǎn)草量高，營養(yǎng)豐富，適口性好[4]，是國內(nèi)外人工飼草基地廣泛種植的重要草種[5]。土壤鹽漬化屬于主要非生物脅迫之一，目前在世界范圍內(nèi)日趨嚴重[6-7]。金花菜主要分布于中國長江流域以南等沿海各省區(qū)，沿海地區(qū)相對嚴重的土壤鹽漬化嚴重影響了金花菜的產(chǎn)量和品質(zhì)[8]。中國北方紫花苜蓿主栽區(qū)鹽堿地也分布較廣，嚴重制約了紫花苜蓿產(chǎn)量的穩(wěn)定提高[9]。因此，研究金花菜和紫花苜蓿的耐鹽性具有重要的現(xiàn)實意義。研究表明，轉(zhuǎn)錄組學是了解鹽脅迫下基因差異表達的較好工具[10-12]。目前，利用轉(zhuǎn)錄組技術(shù)進行鹽脅迫相關(guān)基因篩選和挖掘已經(jīng)在紫花苜蓿中進行了研究，馬進和鄭鋼[13]對鹽脅迫下紫花苜蓿的根系轉(zhuǎn)錄組分析已有報道，篩選了膽汁酸:Na+共轉(zhuǎn)運蛋白、晚期胚胎發(fā)生富集蛋白、谷胱甘肽-s-轉(zhuǎn)移酶基因和轉(zhuǎn)錄因子等與紫花苜蓿鹽脅迫相關(guān)的應(yīng)答關(guān)鍵基因；李紅等[14]對蘇打鹽堿脅迫下的紫花苜蓿進行了轉(zhuǎn)錄組學分析，發(fā)現(xiàn)蘇打鹽堿脅迫下GH3、MYB、HSF等等轉(zhuǎn)錄因子上調(diào)，而EREBP轉(zhuǎn)錄因子下調(diào)。這些研究為揭示紫花苜蓿耐鹽分子機制奠定了基礎(chǔ)，而關(guān)于鹽脅迫下金花菜和紫花苜蓿的轉(zhuǎn)錄組分析及其耐鹽基因的篩選尚未見報道。本研究將以鹽脅迫下的金花菜和紫花苜蓿試管苗為研究對象，究利用RNA-Seq技術(shù)對其的轉(zhuǎn)錄組進行分析，并篩選其耐鹽基因，為更深入地了解金花菜和紫花苜蓿的耐鹽分子機制及選育金花菜和紫花苜蓿耐鹽新品種提供理論依據(jù)。

1 材料和方法

1.1 材料

金花菜種子（編號JHC，由江西農(nóng)業(yè)大學動物科學技術(shù)學院歐陽克惠教授提供）和紫花苜蓿種子（編號ZHMX，南方型紫花苜蓿‘Millennium’種子，北京克勞沃草業(yè)技術(shù)開發(fā)中心提供）。

1.2 方法

1.2.1 金花菜和紫花苜蓿種子的消毒滅菌和萌發(fā) 將金花菜種子和紫花苜蓿種子用自來水沖洗干凈，濾干后帶入無菌室，在超凈工作臺上用75%酒精消毒1 min，無菌水沖洗3遍，再用0.1%氯化汞消毒15 min，無菌水沖洗3遍，置于無菌培養(yǎng)皿內(nèi)無菌濾紙片上，待濾紙片吸去種子表面水后開始接入MS固體培養(yǎng)基（3%蔗糖+7.5%瓊脂粉）。放入恒溫恒濕培養(yǎng)箱（溫度：25±1℃；光照強度：1500～2000 Lx；光照時間：14 h/d；濕度：70%～80%）進行培養(yǎng)。10天左右，金花菜種子和紫花苜蓿種子開始陸續(xù)萌發(fā)，至30天可獲得金花菜和紫花苜蓿試管苗。

1.2.2 金花菜和紫花苜蓿試管苗的鹽脅迫處理 將金花菜和紫花苜蓿試管苗接入MS+50 mmoL/L NaCl+3%蔗糖+7.5%瓊脂粉固體培養(yǎng)基上進行鹽脅迫處理30天。培養(yǎng)條件：溫度25±1℃；光照強度1500～2000 Lx；光照時間14 h/d；濕度70%～80%。

1.2.3 總RNA提取 鹽脅迫(50 mmoL/L NaCl)處理30天金花菜和紫花苜蓿試管苗的總RNA按TRNzol Reagent試劑盒（北京天根）試劑盒說明提取。所提取的RNA經(jīng)電泳檢測合格后送至上海美吉生物醫(yī)藥科技有限公司進行轉(zhuǎn)錄組測序分析。

1.2.4 轉(zhuǎn)錄組測序與分析 轉(zhuǎn)錄組測序進行3次生物學重復(fù)。帶有polyA尾的mRNA通過Oligo(dT)磁珠富集，得到的mRNA在NEB Fragmentation Buffer中用二價陽離子隨機打斷，合成cDNA第一條鏈，RNA鏈用RNaseH降解，合成cDNA第二條鏈。雙鏈cDNA純化后經(jīng)過末端修復(fù)、加A尾、連接測序接頭、篩選cDNA(250～300 bp)、PCR擴增、PCR產(chǎn)物純化獲得文庫。測序采用Illumina平臺。

1.2.5 差異表達基因的篩選和富集分析 差異表達基因的篩選標準為差異倍數(shù)(FC)≥2且錯誤發(fā)現(xiàn)率(FDR)＜0.01。基于描述基因功能的綜合性數(shù)據(jù)庫(http://www.geneontology.org/)進行GO(Gene Ontology)富集分析；基于有關(guān)Pathway主要公共數(shù)據(jù)庫進行KEGG(Kyoto Encyclopedia of Genes and Genomes)富集分析。

1.2.6 耐鹽基因的篩選及其qRT-PCR驗證 在1.2.5的基礎(chǔ)上，篩選與鹽脅迫相關(guān)的差異基因，以GAPDH為內(nèi)參基因，設(shè)計基因引物（表1），進行熒光定量PCR驗證。PCR反應(yīng)程序：預(yù)變性95℃30 s，變性95℃10 s，退火60℃ 30 s，40個循環(huán)。使用2-△△Ct法計算基因表達水平。實驗重復(fù)3次。所有數(shù)據(jù)表示為平均值，并使用SPSS19.0軟件進行統(tǒng)計分析，應(yīng)用單因素方差分析(One-way ANOVA)檢驗與鹽脅迫相關(guān)的差異表達基因在金花(JHC)和紫花苜蓿(ZHMX)試管苗中表達的差異顯著性(P＜0.05)。

表1 與鹽脅迫相關(guān)的差異表達基因的引物

2 結(jié)果與分析

2.1 測序數(shù)據(jù)統(tǒng)計

本次測序共獲得15.53 Gb Clean Data，金花菜組和紫花苜蓿組Clean Data均達到7.64 Gb以上，Q30堿基百分比在92.87%以上。

2.2 轉(zhuǎn)錄本長度分布

金花菜和紫花苜蓿組裝得到的unigene個數(shù)為68093，Transcript個數(shù)為 112419，N50 平均長度為1270 bp。

2.3 與組裝結(jié)果比對統(tǒng)計

金花菜組和紫花苜蓿組本次分析比對率從79.12%到83.16%不等。

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

金花菜組和紫花苜蓿組表達基因數(shù)量為66139，表達轉(zhuǎn)錄本數(shù)量為109244，所有基因數(shù)量為68093，所有轉(zhuǎn)錄本數(shù)量為112419。

2.5 轉(zhuǎn)錄因子家族分析

金花菜組和紫花苜蓿組轉(zhuǎn)錄因子家族主要有LOB、NAC、BBR-BPC、ZF-HD、CAMTA、Nin-like、CPP、EIL、C2C2、bZIP、HSF、B3_superfamily、GRAS、TCP、NF-X1、MADS、WRKY、E2F/DP、SBP、GRF、SRS、AP2/ERF、C2H2、bHLH、LBD(AS2/LOB)、C3H、NF-Y、MYB_superfamily、Whirly、FAR1、GeBP、BES1等。

2.6 表達量分布

金花菜組和紫花苜蓿組的表達量分布見圖1。從圖1可知，金花菜組和紫花苜蓿組表達量FPKM的對數(shù)值在-2～4.4之間，表達量密度在0～0.6之間。

圖1 金花菜組和紫花苜蓿組unigene/transcript的表達量密度分布情況

2.7 樣本間Venn分析

金花菜組unigene/transcript數(shù)目為41748，紫花苜蓿組unigene/transcript數(shù)目為43592，兩者共有的unigene/transcript數(shù)為20886。

2.8 樣本間相關(guān)性熱圖

金花菜組與紫花苜蓿組樣本間的相關(guān)性為0.634。

2.9 表達量差異統(tǒng)計

金花菜組和紫花苜蓿組差異基因數(shù)為26722，金花菜組有15850個基因下調(diào)，有10872個基因上調(diào)。

2.10 差異基因火山圖

金花菜組和紫花苜蓿組差異基因散點圖和火山圖見圖2。由圖2可知，金花菜組和紫花苜蓿組差異基因表達量排名前十的基因有腺苷激酶2(adenosine kinase 2)、伸長因子 1-γ(elongation factor 1-gamma)、晚期胚胎發(fā)生蛋白2(late embryogenesis abundant protein 2)、氯通道蛋白CLC-c亞型X1(chloride channel protein CLC-c isoform X1)、非特征蛋白LOC11444183亞型X2(uncharacterized protein LOC11444183 isoform X2)、非特征蛋白LOC25487653(uncharacterized protein LOC2 5487653)、硫氧還蛋白(thioredoxin)、葉綠體可能的GTP二磷酸激酶RSH2(probable GTP diphosphokinase RSH2,chloroplastic）)、假 定 的 跨 膜 蛋 白 (transmem brane protein,putative)和假定的蛋白絲氨酸/蘇氨酸激酶CK1-CK1-Pl家族(putative protein-serine/threonine kinase CK1-CK1-Pl family)。

圖2 金花菜組和紫花苜蓿組差異基因散點圖和火山圖

2.11 差異基因COG分類統(tǒng)計

金花菜組和紫花苜蓿組差異基因牽涉的主要功能有RNA加工和修飾、染色質(zhì)結(jié)構(gòu)和動力學、能源生產(chǎn)和轉(zhuǎn)換、細胞周期控制、細胞分裂、染色體分配、氨基酸轉(zhuǎn)運與代謝、核苷酸轉(zhuǎn)運與代謝、碳水化合物運輸和代謝、輔酶運輸和代謝、脂質(zhì)轉(zhuǎn)運與代謝、翻譯、核糖體結(jié)構(gòu)與生物發(fā)生、轉(zhuǎn)錄本、復(fù)制、重組和修復(fù)、細胞壁/膜/包膜生物發(fā)生、翻譯后修飾、蛋白質(zhì)周轉(zhuǎn)、伴侶、無機離子轉(zhuǎn)運與代謝、次生代謝物生物合成、運輸和分解代謝、信號轉(zhuǎn)導(dǎo)機制、細胞內(nèi)運輸、分泌和囊泡運輸、防御機制、核結(jié)構(gòu)、細胞骨架等。

2.12 差異基因GO富集分析

金花菜組和紫花苜蓿組差異基因GO富集分析結(jié)果如圖3所示。從圖3可知，生物過程(BP)富集顯著排名前五主要包括細胞殺傷cell killing、免疫系統(tǒng)過程immune system process、碳利用carbon utilization、生物調(diào)節(jié)biological regulation、代謝過程metabolic process；細胞成分(CC)富集顯著排名前五主要包括腔上包膜membrane-enclosed lumen、含蛋白質(zhì)復(fù)合物protein-containing complex、細胞 cell、膜部件 membrane part、胞外區(qū)部分extracellular region part；分子功能(MF)富集顯著排名前五主要包括翻譯調(diào)節(jié)器活動translation regulator activity、轉(zhuǎn)錄調(diào)節(jié)活性transcription regulator activity、結(jié)構(gòu)分子活性 structural molecule activity、蛋白質(zhì)折疊伴侶protein folding chaperone、分子載體活性molecular carrier activity。

圖3 金花菜組和紫花苜蓿組差異基因GO富集分析圖

2.13 差異基因KEGG富集分析

金花菜組和紫花苜蓿組差異基因KEGG富集分析結(jié)果如圖4所示。從圖4可知，金花菜組和紫花苜蓿組差異基因KEGG富集顯著的前十的Pathways有核糖體Ribosome、光合作用天線蛋白Photosynthesisantenna proteins、糖 酵 解/糖 異 生 Glycolysis/Gluconeogenesis、光合作用Photosynthesis、苯丙烷生物合成Phenylpropanoid biosynthesis、丙酮酸代謝Pyruvate metabolism、纈氨酸、亮氨酸和異亮氨酸降解Valine,leucine and isoleucine degradation、乙醛酸和二羧酸代謝Glyoxylate and dicarboxylate metabolism、精氨酸和脯氨酸代謝Arginine and proline metabolism、光合生物的固碳作用Carbon fixation in photosynthetic organisms。

圖4 金花菜組和紫花苜蓿組差異基因KEGG富集分析圖

2.14 CDS長度分布

金花菜組和紫花苜蓿組差異基因CDS長度為200～200的unigene/transcript數(shù)目為5853，長度為201～400的unigene/transcript數(shù)目為22301，長度為401～600的 unigene/transcript數(shù)目為 8085，長度為 601～800的unigene/transcript數(shù)目為 4761，長度為 801～1000 的unigene/transcript數(shù)目為 3362，長度為 1001～1200 的unigene/transcript數(shù)目為 2735，長度為 1201～1400 的unigene/transcript數(shù)目為 2025，長度為 1401～1600 的unigene/transcript數(shù)目為 1723，長度為 1601～1800 的unigene/transcript數(shù)目為1081，長度＞1800的unigene/transcript數(shù)目為3375，CDS的總數(shù)為55301。

2.15 耐鹽基因的篩選及其qRT-PCR檢測

選取鹽脅迫條件下金花菜組和紫花苜蓿組基因表達差異性最大的10個差異基因，并對其進行qRT-PCR驗證（表2）。由表2可知，qRT-PCR測定的10個DEGs在金花菜組和紫花苜蓿組中的表達水平趨勢與轉(zhuǎn)錄組測序基因的表達趨勢是一致，表明本研究中得到的DEGs是可信的。從表2還可知，與紫花苜蓿組相比，金花菜組的WRKY轉(zhuǎn)錄因子(WRKY transcription factor)、蛋白TIFY 8(protein TIFY 8)、轉(zhuǎn)錄因子MYB13(transcription factor MYB13)、核運輸因子2A(nuclear transport factor 2A)、低溫鹽響應(yīng)蛋白(low temperature and salt responsive protein)等基因下調(diào)。而堿性亮氨酸拉鏈43(basic leucine zipper 43)、NAC轉(zhuǎn)錄因子47(NAC transcription factor 47)、ABC轉(zhuǎn)運A家族成員7(ABC transporter A family member 7)、晚期胚胎發(fā)生豐富蛋白2(late embryogenesis abundant protein 2)、乙烯反應(yīng)轉(zhuǎn)錄因子ERF110(ethylene-responsive transcrip tionfactor ERF110)等基因上調(diào)。

表2 與耐鹽相關(guān)差異基因的qRT-PCR驗證

3 討論

植物的耐鹽脅迫涉及多基因參與表達且調(diào)控復(fù)雜，一般通過改變基因的表達來調(diào)節(jié)代謝與信號轉(zhuǎn)導(dǎo)[15]。熱依麥阿依·阿布都艾尼等通過鹽脅迫下棉花根系的轉(zhuǎn)錄組分析篩選出其耐鹽基因，結(jié)果表明‘新陸中69號’差異表達基因共1557個，上調(diào)基因891個，下調(diào)基因666個，富集到苯丙烷類合成代謝途徑的差異基因最多[16]。在本試驗中，鹽脅迫下金花菜組和紫花苜蓿組2個樣本的差異基因數(shù)為26722，金花菜組有15850個基因下調(diào)，有10872個基因上調(diào)，富集到核糖體、糖酵解/糖異生、丙酮酸代謝和苯丙烷生物合成等代謝途徑的差異基因較多。這說明在不同植物抵御鹽脅迫過程中，起到調(diào)控工作的代謝途徑也不盡相同。

WRKY轉(zhuǎn)錄因子主要參與植物的生物和非生物脅迫調(diào)控[17]。大量研究表明，WRKY轉(zhuǎn)錄因子能響應(yīng)調(diào)控植物鹽脅迫，緩解鹽脅迫傷害[18]。有研究表明，紫花苜蓿MsWRKY11可通過調(diào)節(jié)細胞K+/Na+內(nèi)穩(wěn)態(tài)，保護過氧化物酶活性，降低Na+對細胞膜結(jié)構(gòu)的破壞，從而提高紫花苜蓿的耐鹽能力[19]。在本試驗中，金花菜組的WRKY轉(zhuǎn)錄因子下調(diào)，而紫花苜蓿組的WRKY轉(zhuǎn)錄因子上調(diào)，表明鹽脅迫下WRKY轉(zhuǎn)錄因子在紫花苜蓿組起到重要調(diào)控作用。這與鹽脅迫下南方型紫花苜蓿‘Millennium’WRKY轉(zhuǎn)錄因子表達上調(diào)[20-22]的研究結(jié)果一致。

TIFY基因家族是植物特異性基因家族，主要調(diào)控植物生長發(fā)育、響應(yīng)逆境脅迫、協(xié)同信號轉(zhuǎn)導(dǎo)。研究表明，紫花苜蓿MsTIFY基因?qū)}脅迫有明顯的響應(yīng)，在調(diào)控鹽脅迫發(fā)揮著關(guān)鍵性作用[23]。在本試驗中，金花菜組的蛋白TIFY 8基因下調(diào)，而紫花苜蓿組的蛋白TIFY 8基因上調(diào)，表明鹽脅迫下蛋白TIFY 8基因在紫花苜蓿組起到重要調(diào)控作用。這與鹽脅迫下新疆大葉紫花苜蓿MsTIFY6、MsTIFY32、MsTIFY46、MsTIFY48、MsTIFY65五條基因表達上調(diào)[24]的研究結(jié)果一致。

堿性亮氨酸拉鏈(basic leucine zipper，bZIP)蛋白是最大且功能最具多樣性的轉(zhuǎn)錄因子家族之一，主要參與植物生長發(fā)育以及脅迫響應(yīng)等生理過程，特別是在非生物脅迫的抗逆反應(yīng)中起重要作用，例如鹽脅迫等[25]。研究表明，紫花苜蓿‘中苜1號’MsbZIP1基因?qū)}脅迫處理有顯著響應(yīng)，推測該基因可能參與調(diào)控紫花苜蓿鹽脅迫[26]。在本試驗中，金花菜組的bZIP43基因上調(diào)，而紫花苜蓿組的bZIP43基因下調(diào)，表明鹽脅迫下bZIP43基因在金花菜組起到重要調(diào)控作用。這與紫花苜蓿鹽誘導(dǎo)HD-Zip類轉(zhuǎn)錄因子MsHB2表達量顯著升高從而負調(diào)控鹽脅迫[27]的研究結(jié)果一致，鹽脅迫下南方型紫花苜?！甅illennium’Ms bZIP等轉(zhuǎn)錄因子被篩選為與鹽脅迫應(yīng)答相關(guān)候選轉(zhuǎn)錄因子[21]。

MYB類轉(zhuǎn)錄因子是真核生物中分布最廣泛的家族之一，廣泛參與包括高鹽、干旱、極端溫度在內(nèi)的多種逆境脅迫應(yīng)答[28]。有研究表明，鹽脅迫下紫花苜蓿大量出現(xiàn)鹽脅迫應(yīng)答的轉(zhuǎn)錄因子候選基因如Ms MYB等[20]。在本試驗中，金花菜組的轉(zhuǎn)錄因子MYB13下調(diào)，而紫花苜蓿組的轉(zhuǎn)錄因子MYB13上調(diào)，表明鹽脅迫下轉(zhuǎn)錄因子MYB13在紫花苜蓿組起到重要調(diào)控作用。

NAC轉(zhuǎn)錄因子是植物特有轉(zhuǎn)錄調(diào)控因子，主要參與植物次生生長、激素信號轉(zhuǎn)導(dǎo)和逆境御反應(yīng)[29]。鹽脅迫下的表達分析顯示，紫花苜蓿Ms NAC1基因表達量呈上調(diào)趨勢，且根中的表達量上調(diào)幅度高于葉片，說明Ms NAC1基因參與調(diào)控鹽脅迫的生理響應(yīng)[30]。在本試驗中，金花菜組的NAC轉(zhuǎn)錄因子47上調(diào)，而紫花苜蓿組的NAC轉(zhuǎn)錄因子47下調(diào)，表明鹽脅迫下NAC轉(zhuǎn)錄因子47在金花菜組起到重要調(diào)控作用。這與鹽脅迫下紫花苜蓿通過轉(zhuǎn)錄組測序分析發(fā)掘10個鹽誘導(dǎo)表達NAC轉(zhuǎn)錄因子[31]的研究結(jié)果一致。

ABC(ATP-binding cassette)轉(zhuǎn)運蛋白因其蛋白質(zhì)中含有的核苷酸結(jié)合域而得名，是植物器官生長、營養(yǎng)物質(zhì)運輸、抵抗病原體、脅迫反應(yīng)等過程中至關(guān)重要的轉(zhuǎn)運蛋白[32]。核運輸因子(nuclear transport factor，NTF)是細胞核轉(zhuǎn)運蛋白家族IV亞組的成員，表達局限于根部，并受到逆境脅迫的刺激，NTF可通過影響植物中的ROS含量來提高植物的抗逆性[33]。低溫鹽響應(yīng)蛋白(Low temperature and salt responsive protein，LTSR)也與植物調(diào)控鹽脅迫密切相關(guān)，可提高植物的抗逆性[34]。在本試驗中，金花菜組的ABC轉(zhuǎn)運A家族成員7基因上調(diào)、NTF2A和LTSR下調(diào)，而紫花苜蓿組的ABC轉(zhuǎn)運A家族成員7基因下調(diào)、NTF2A和LTSR上調(diào)，表明鹽脅迫下ABC轉(zhuǎn)運A家族成員7基因在金花菜組起到重要調(diào)控作用，而鹽脅迫下NTF2A和LTSR在紫花苜蓿組起到重要調(diào)控作用。

晚期胚胎發(fā)生豐富蛋白(Late embryogenesis abundant，LEA)基因是生物體中廣泛存在的一類廣泛響應(yīng)逆境脅迫的功能基因，主要參與調(diào)控植物對鹽脅迫等的耐受[35]。研究表明，在逆境脅迫下，紫花苜蓿MsLEA4-4基因的表達量隨處理時間的延長呈逐漸上升趨勢，表明MsLEA4-4基因參與了紫花苜蓿的抗逆性調(diào)控[36]。在本試驗中，金花菜組的LEA2基因上調(diào)，而紫花苜蓿組的LEA2基因下調(diào)，表明鹽脅迫下LEA2基因在金花菜組起到重要調(diào)控作用。

乙烯反應(yīng)轉(zhuǎn)錄因子ERFs轉(zhuǎn)錄因子主要響應(yīng)逆境脅迫[37]。研究表明，紫花苜蓿MsERF003基因在各個組織中的表達量均有上調(diào)，在莖和莢果中尤為顯著，暗示MsERF003基因在響應(yīng)逆境脅迫時，主要在莖和胚胎發(fā)育中發(fā)揮作用[38]。在本試驗中，金花菜組的ERF110基因上調(diào)，而紫花苜蓿組的ERF110基因下調(diào)，表明鹽脅迫ERF110基因在金花菜組起到重要調(diào)控作用。