低溫脅迫下苦瓜葉片轉錄組差異基因分析及生理響應特征

杜文麗 陳中釤 許端祥 徐同偉 高 山 溫慶放

(1福州市蔬菜科學研究所，福建 福州 350111；2福建省農業科學院作物所，福建 福州 350013)

苦瓜(Momordica charantiaL.)又名涼瓜，起源于東印度，現主要種植于熱帶、亞熱帶和溫帶區域，加勒比海和南美地區也有栽培。近年來，人們逐步認識到苦瓜的營養價值及保健功能。我國南方地區生產的苦瓜除供本地和北運外，還遠銷東南亞地區，己成為南方部分地區出口創匯的主要蔬菜之一[1-2]。苦瓜屬喜溫蔬菜，在實際生產中，苦瓜苗期常常遭受早春寒流或秋季氣溫驟降的侵襲，10℃以下的冷害導致生長緩慢甚至形成僵苗，嚴重制約苦瓜種植范圍的擴大和上市時節的調節，生產上亟需選育耐低溫且綜合性狀好的新品種[3]。

與植物耐低溫相關的研究逐漸由細胞學、生理生化轉向分子生物學，且植物耐低溫分子機制的研究已取得了長足進步，現已有100 余種轉基因耐低溫植物，共涉及50 多個物種[4-5]。目前，低溫脅迫涉及一系列復雜的信號傳導網絡，誘導植物體內相關基因的表達，但對于植物將低溫脅迫信號傳遞到細胞核內激活轉錄因子，以及激活下游基因表達使機體出現耐低溫特性的機制仍不清楚。前人研究發現CBF1 通過結合CRT/DRE(c-repeat/dehydration responsive element)順式作用元件，調控耐低溫特性基因表達來參與植物中低溫脅迫下的信號傳導，這種通過轉錄因子所調控的信號路在學術界已被廣泛驗證[6-7]。有研究學者先后在擬南芥中獲得CBF2、CBF3 和CBF4 基因轉入擬南芥中，發現擬南芥植株耐低溫脅迫的能力得到顯著提升，且導致許多與低溫脅迫相關的生理生化指標發生變化，如脯氨酸和可溶性糖含量增加等[8]。周麗霞等[9]發現WRKY7、WRKY1 等基因對油棕的耐低溫性發揮調控作用。李佳[10]研究發現馬鈴薯在4℃處理后 St2RMYB80 基因表達量下調，St4RMYB1 和St3RMYB2 基因表達量上調。

第二代測序技術已廣泛應用于轉錄組測序(RNAseq)領域，通過比對植物生物和非生物脅迫轉錄組基因reads 數獲得植物體大量基因信息，能夠在分子水平上了解基因表達量的變化，揭示特定的生物學過程和分子機制[11-13]。Chen 等[13]通過Illumina nextseq 2500 測序平臺對低溫脅迫下的甘蔗葉片進行研究，發現2 583 個基因上調表達，3 302 個基因下調表達。黃超[11]利用Illumina Hiseq-2000 測序平臺對4℃低溫處理0.8 和32 h 的甜瓜葉片進行測序，篩選出19個包括編碼WRKY、b HLH、NAC、AP2/ERF等轉錄因子在內可能與耐低溫具有潛在相關的候選基因。此外，轉錄組測序技術已應用于大白菜、花椰菜、黃瓜等蔬菜作物研究中[14-17]。目前，關于苦瓜在低溫脅迫方面的研究還限于表型、生理生化指標鑒定和遺傳規律等方面[18-19]，通過RNA-seq 揭示苦瓜的耐低溫機制方面研究鮮有報道，在一定程度上約束了苦瓜育種進程。

本研究以福州市蔬菜科學研究所選育的較耐低溫苦瓜自交系43 為研究對象，依托百邁客Illumina HiSeq2500 平臺對26℃(對照)和模擬福建早春夜溫8℃氣候[20]處理12 h 的兩個樣品葉片進行轉錄組測序分析，再通過生理指標分析和分子手段對苦瓜響應低溫脅迫的分子機制進行探討，同時進一步挖掘苦瓜耐低溫相關的候選基因，以期為苦瓜耐低溫基因工程提供一定的理論支撐。

1 材料與方法

1.1 試驗材料

供試材料為耐低溫苦瓜自交系43，由福州市蔬菜科學研究所苦瓜課題組選育。選取飽滿健壯的苦瓜種子采用溫湯浸種法催芽，待種子露白后轉移至營養缽中(椰糠＋草木灰)中，隨后放置于26℃培養箱(立德泰勒科學儀器有限公司，型號:LT-ACC)中培養，光強8 000 lx，相對濕度70%)。約2 周后苦瓜苗生長至四葉一心時，將培養箱溫度調至8℃，選取生長良好、發育狀態一致的苦瓜苗進行低溫處理，處理前澆透水。觀察每個時期的植株形態，其他培養條件同上[20]。處理0(對照)、4、8、12 h 后剪取幼苗真葉測定相關生理指標。剪取處理后剪取12 h 幼苗真葉用錫箔紙包好，以0 h 為對照，立即放置于液氮中速凍后送北京百邁客生物科技有限公司進行RNA-seq。

1.2 生理指標測定方法

硫代巴比妥酸(thiobarbituric acid，TAB)法用于測定丙二醛(malondialdehyde，MDA)含量。酸性茚三酮法用于測定游離脯氨酸(proline，Pro)含量。電解質滲出率(rate of electrolyte efflux，REC)的測定參考文獻[21]。

1.3 RNA-seq 方法

總RNA 按照提取后經消化處理后及時對RNA 樣品的純度、濃度和完整性進行檢測，保證使用高質量的RNA 進行轉錄組測序。依托Illumina HiSeq2500 高通量測序平臺進行上機測試工作，PE150 為測序讀設置值。利用FASTX 軟件對獲得的數據進行質量控制，即去除Raw Data 中的接頭和低質量序列后得到Clean Data，再使用Trinity 軟件進行組裝最終獲得Unigene序列。利用BLAST 程序(設置E-value= 10-5)使用NR、COG、KOG 及KEGG 等數據庫對Unigene 進行注釋及和分類。基因表達衡量指標采用FPKM 法[22]，使用DESeq[23]對樣品間差異表達進行分析，以樣品間表達量符合錯誤發現率數≥2 且≤0.5 作為篩選標準獲得差異表達基因。最后對差異表達的基因進行GO 功能注釋、COG 功能注釋和富集分析。

1.4 實時熒光定量PCR 檢測

參照植物RNA 提取試劑盒(RNAprep Pure Plant Kit，天根生化科技有限公司)獲得苦瓜葉片的總RNA，然后參照[Prime ScriptTM RT reagent Kit，寶生物工程(大連)有限公司]試劑盒說明反轉錄合成cDNA。用Primer 5.0 軟件設計與低溫逆境相關的苦瓜差異表達基因引物，以McCYP(GeneBank 登錄號:HQ171897)作為內參基因，候選基因的引物序列見表1，引物由北京六合華大基因科技股份有限公司合成。實時熒光定量PCR(quantitative real-tirne PCR，qRTPCR)采用ABI 7500 Fast Real-Time PCR System(美國)進行PCR 擴增，每次試驗設3 次重復。反應體系參照[Super Real Pre Mix Plus (SYBR Green)，天根生化科技有限公司)]試劑盒說明書進行，數據采用2-ΔΔCt方法計算[24]。

1.5 數據處理

采用SPSS 16 軟件對試驗數據進行統計，利用Excel2016.1.0 軟件作圖。

2 結果與分析

2.1 低溫脅迫下苦瓜的表型與生理生化指標變化

苦瓜經低溫脅迫4 h 后，苦瓜的真葉就開始出現萎蔫失水的現象，這是植物受到低溫脅迫的典型反應；低溫脅迫8 h 后，葉柄軟化、向下彎曲，未成熟莖出現更嚴重的萎蔫現象(圖1-A)。但在低溫脅迫12 h 后，又出現明顯的恢復現象，表明苦瓜43 號自交系具有良好的低溫耐受能力，適于作為本研究的材料。

圖1 低溫脅迫前后苦瓜表型變化和生理生化變化Fig.1 Phenotypic and physiological changes in cold-stressed bitter gourd

由圖1-B 可知，在經過低溫脅迫后MDA 含量先升高后降低，在低溫脅迫8 和12 h 后，比0 h 分別增加了96.80%和60.34%，峰值約為對照的2 倍，說明此時苦瓜的細胞膜受到損傷導致MDA 含量增加，但在低溫脅迫12 h 后損傷有所恢復。由圖1-C 可知，隨著低溫脅迫時間延長，Pro 含量整體呈現先升高后降低的趨勢，在低溫脅迫4、8 和12 h 時，Pro 含量相比于對照分別提高了114.47%、125.20%和87.07%，這與MDA 含量的變化情況相一致。由圖1-D 可知，REC 隨脅迫時間的延長而持續增加并保持在較高水平，在低溫脅迫4、8 和12 h 時，比對照分別增加了7.55、9.38 和12.47個百分點。

2.2 測序數據整體評估及差異表達基因的獲得

為深入探索苦瓜在轉錄水平上對低溫脅迫的應答，剪取8℃低溫條件下處理0 和12 h 的苦瓜葉片，采用Illumiana HiSeq2500 平臺進行轉錄組測序，共獲得1.77 Gb Clean Data，Q30 值均大于91.00%(表2)，GC含量均在48%的水平。Clean Data 經過濾、拼接和組裝后獲得的Unigene 數量為89 757 條，N50 長度為1 239 bp，平均長度達到700.75 bp (表3)。Unigene序列數量最多的分布在200 ～300 bp 之間，在500 ～1 000 bp 之間次之(42.64%)，Unigene 長度達到1 000 bp 以上的占18.2%以上，共計16 329 條。以上結果表明RNA-seq 結果符合下一步生物信息學分析的要求。

表2 測序數據結果統計Table 2 Statistic of evaluating the sequencing data

表3 測序數據結果統計Table 3 Statistic of evaluating the sequencing data

將得到的Unigene 序列比對結果經NR、KEGG等多個數據平臺得到對應的注釋信息(表4)。在預設的E 值范圍最終有48 662 個注釋信息的Unigene，占全部Unigenes 的54.2%。比對得到的Unigene 長度大于1 000 nt 的達到14 964 個，占所有得到注釋Unigene 的30.5%。但仍有41 095 個Unigenes 在數據庫中沒有得到注釋。通過對差異表達分析結果進行統計，在低溫脅迫處理后共獲得1 285 個差異表達基因，其中916 個基因上調、369個基因下調。

表4 測序數據結果統計Table 4 Statistic of evaluating the sequencing data

2.3 差異表達基因GO 分類

對差異表達基因進行GO 聚類分析聚類，結果有1 495 條Unigenes 映射到GO 不同功能節點上，顯著富集主要包括在細胞組分、分子功能、生物過程中的50個亞類(圖2)。其中最大的主類群是參與生物學過程的差異表達基因，共有749 個；歸入細胞組分和分子功能類的分別為347 和399 個；細胞部分、細胞器和膜結構在細胞組分中是排名前三的亞類，分別為123、62 和46 個差異表達基因。在分子功能主類群，大多數的差異表達基因歸屬于催化活性(177 個)和結合活性(149個)。在生物學過程類群，由191 條Unigene 組成的代謝進程是數量最多的亞類，其次是細胞過程(174 個)和單一生物過程(134 個)。生物過程分組中共富集到14 個GO 條目，差異基因所占比例較大的有代謝過程、細胞過程、單個有機體過程、應激反應、生物調節。在細胞組分分組中共富集到8 個GO 條目，差異表達基因所占比例較高的有細胞部分細胞、膜、細胞器、膜結構。在分子功能分組中共富集11 個GO 條目，催化活性、結合活性占較大比例。這些大量積累的基因表明，苦瓜在遭受低溫脅迫過程中進行著復雜的新陳代謝和酶促反應。

2.4 差異表達基因的COG 分類

圖2 苦瓜轉錄組GO 功能分類Fig.2 Gene ontology classification of Momordica charantia L.transcriptome

將轉錄組獲得的差異表達基因與COG 數據庫進行比對，對其所編碼的蛋白質進行直系同源分類，進一步預測基因功能并進行分類統計(表5)。共得到了239 個COG 功能注釋的信息，按其功能劃分為24 類，其中有3 類功能類型在本研究中無Unigene 比對。其中，僅一般功能預測(16.74%)和轉錄(10.04%)類所占的比例最大，隨后依次是次生代謝物合成，運輸和代謝(9.62%)及重復，重組和修飾(8.79%)。

2.5 轉錄組差異表達基因代謝途徑分析

將差異基因比對到KEGG 數據庫獲得基因注釋結果以進一步探索基因功能(表6)。結果顯示共有237 個差異基因歸入94 種代謝路徑。其中參與核糖體代謝途徑的差異表達基因有19 個；其次有18 個差異表達Unigene 與苯丙素生物合成相關；與苯丙氨酸代謝路徑相關的差異表達基因占第三位；此外，參與植物激素信號傳導以及淀粉和糖類代謝途徑的差異表達基因均為11 個。與脯氨酸合成相關的精氨酸和脯氨酸代謝通路與低溫脅迫相關的氧化磷酸化(oxidative phosphorylation)通路放入差異表達基因分別獲得1 個和2 個。代謝途徑包含碳水化合物、信號傳導、光合等涉及生命活動的各個階段，表明這些代謝途徑在苦瓜遭受低溫脅迫時可能發揮著重要作用。

2.6 苦瓜耐低溫相關差異表達候選基因

為了挖掘參與耐低溫脅迫信號網絡相關基因，在差異表達基因中篩選到了包括WRKY、ERF/AP2(乙烯應答)、ZAT轉錄因子等可能參與苦瓜耐低溫脅迫的10 個重要基因(表7)。10 個候選基因的序列和相似性分析見附表1。

2.7 qRT-PCR 分析

為了驗證RNA-seq 數據的準確性并進一步分析苦瓜候選基因的耐低溫機制，對挑選出的10 個候選基因利用qRT-PCR 觀測低溫脅迫下動態變化情況，結果如圖3 所示。WRKY33 的表達量在脅迫4 h 后急劇上升，8 h 時稍下調隨后又有所上升，但都始終高于對照(0 h)；MYB108 的表達量隨脅迫時間的延長呈逐步上升趨勢；ERF017 的表達量隨低溫脅迫時間的延長呈現先下降后上升的趨勢；ERF021 和ERF053 的表達量在受到低溫脅迫4 h 時急劇上升，且均在低溫脅迫處理8 h 時達到峰值，在低溫脅迫處理12 h 時稍有下降但仍高于對照；NAC25 的表達量在受到低溫脅迫后下降，低溫脅迫處理8 h 表達量上升后又下調；Ga2＋ATPase和Ga2＋ATPase2 的表達量隨著低溫脅迫時間的延長均呈現先上升后下降的趨勢；ZAT11 和ZAT12的表達量在受到低溫脅迫后先上調后下降，且峰值均出現在低溫脅迫8 h，低溫脅迫12 h 時仍高于對照。所有基因的表達變化與RNA-seq 所得的變化情況一致，說明轉錄組結果數據可靠。這些基因表達量的上調或下調也預示著它們在苦瓜耐低溫脅迫中可能發揮著重要作用。

表5 差異表達基因的COG 注釋分類情況Table 5 Classification of COG annotation about DEGs

3 討論

苦瓜經低溫脅迫后出現萎蔫失水的現象，隨后又明顯恢復，表明該苦瓜材料具有良好的低溫耐受能力。MDA 是評價植物細胞膜系統損傷程度的重要指標，本研究中苦瓜受到一定程度的冷害后MDA 含量先升高后降低，說明植物損傷后逐漸得到恢復，這與在低溫脅迫下甜瓜不同品種間的耐低溫性與MDA 含量成顯著負相關的研究結果一致[10]。Pro 是植物細胞內重要的滲透調節物質，通過提高細胞持水能力從而提高植物的耐低溫特性[1，25]。本研究中Pro 含量先上升后下降，說明在低溫脅迫前期，苦瓜體內能夠快速促進Pro含量的提升使植物細胞膜結構免受損傷。MDA 和Pro的變化情況與苦瓜的表型變化一致，說明苦瓜在低溫脅迫時可能刺激相關通路基因的表達以保護細胞膜系統不受損害。電解質參透率能反映植物在逆境條件下細胞膜透性程度[1]。本研究中苦瓜葉片的電解質參透率在低溫脅迫的過程中始終維持較高水平，表明低溫造成的膜損傷并沒有在12 h 內得到完全修復。

表6 苦瓜部分差異表達Unigene 的代謝通路Table 6 Some metabolic pathways involving DEG of Momordica charantia L

表7 差異表達中耐低溫脅迫候選基因Table 7 Low temperature candidate gene of the DGEs

圖3 低溫脅迫下苦瓜候選基因的相對表達量Fig.3 Relative expression of candidate gene under cold stress in Momordica charantia L.

通過對低溫脅迫條件下苦瓜葉片進行RNA-seq分析，共獲得1 285 個差異表達基因，受低溫誘導表達的基因數目明顯多于低溫抑制的基因數目。差異表達基因GO 富集到的大量基因存在于代謝、催化活性和應激反應亞類中，表明苦瓜在低溫脅迫過程中進行著復雜的新陳代謝和應激反應。差異表達基因歸類到物質代謝、植物激素信號傳導、光合作用相關蛋白及植物晝夜節律等通路，表明苦瓜可能通過改變信號傳導、環境應激及次級代謝物合成來參與低溫脅迫的響應，這與黃瓜、番茄、葡萄通過轉錄組測序技術研究低溫響應相關差異表達基因的主要通路結果類似[7，26-27]。研究表明WRKY、MYB、AP2/ERF、ZAT和NAC這幾類轉錄因子在調節植物低溫脅迫應答網絡中起著關鍵性的調控作用[28-30]。對篩選出表達量變化較為明顯的幾類轉錄因子進行qRT-PCR 驗證并對這些基因在苦瓜低溫脅迫過程進行動態跟蹤，發現WRKY33 基因表達量先上調后下降隨后又上升，MYB108 表達量持續上升、ZAT11 和ZAT12 上調表達后下調，與前人研究結果也相一致[31]。ERF021 和ERF053 上調表達，之后出現下降趨勢，ERF017 表達量隨著脅迫時間的延長先略下調后上調，這與Du 等[32]證實在冷脅迫后期大量ERF上調表達的研究結果一致。NAC25 表達量在低溫脅迫后劇烈下降，預示其可能負調控植物低溫應答。Ca2＋作為第二信使在植物感應外界逆境和信號傳遞過程中發揮重要的作用，Ca2＋ATPase的功能研究證明該家族基因響應各種逆境脅迫[33-34]，本研究Ca2＋ATPase基因表達受低溫誘導明顯，說明其可能在苦瓜應對低溫逆境時發揮重要作用。qRT-PCR 和轉錄組測序基因表達的上調和下調一致，具體的定量值有差異可能是機器或者不同批次處理植株中取樣引起的誤差，表明轉錄組測序的結果真實可信，基因表達和生理變化能夠相互印證。此外，發現脯氨酸合成、氧化磷酸化過程通路中3 個差異基因表達被激活，這與生理生化測試結果一致，可作為苦瓜耐低溫基因挖掘的補充資源。

由本研究結果推測在苦瓜耐低溫脅迫中Ca2＋作為信號分子激活相關轉錄因子的表達，被激活的轉錄因子再調控下游功能基因的表達及修飾，減弱了低溫脅迫對植物細胞膜產生的破壞，促使植物體內失衡的活性氧和激素水平快速恢復，從而提高了植物的耐低溫性。本研究篩選出的候選基因可能參與苦瓜機體耐低溫脅迫過程，為下一步開展苦瓜耐低溫相關分子機制研究提供了一定的理論支持。

4 結論

本研究以耐低溫品種自交系43 號為材料對低溫冷脅迫下苦瓜進行形態和相關生理指標測定及分析，并對26、8℃處理下的苦瓜葉片進行RNA-seq 獲得了豐富的轉錄本信息，對差異表達基因的基因功能進行注釋，了解苦瓜葉片在低溫脅迫下生理變化過程，篩選出包括WRKY、NAC、AP2/ERF等轉錄因子家族候選基因，可能參與苦瓜機體耐低溫脅迫過程調控。10 個差異表達基因的qRT-PCR 動態分析結果與RNA-seq 變化一致，說明轉錄組數據可靠，基因表達變化和生理指標的變化也一致。本研究利用轉錄組測序技術對低溫脅迫下的苦瓜材料進行測序，通過生物信息學方法對苦瓜耐低溫性的多種代謝途徑進行分析，并從分子水平闡述了苦瓜低溫脅迫下的響應機制，為進一步研究苦瓜耐低溫脅迫分子生物學研究功能奠定了一定的理論基礎。

附表1: 10 個候選基因相似性比較結果Schedule 1: Similarity comparison results of 10 candidate genes