2個四倍體栽培棉種冷脅迫響應基因的鑒定與比較分析

祁博文，郭夢蘭，盧合均，梅歡，趙汀，方磊

（浙江大學農業(yè)與生物技術學院，杭州 310058）

棉花是天然可紡織纖維的主要來源，是一種重要的經濟作物。陸地棉（Gossypium hirsutum）和海島棉（G.barbadense）都是四倍體栽培棉種，由同一祖先馴化而來。陸地棉具有產量高、對逆境耐受能力好的特點，但是纖維品質一般；海島棉的棉鈴小，相對產量低，卻能提供更細長強韌的優(yōu)質纖維。陸地棉是最為廣泛的栽培品種，占我國棉花種植面積的95%，而長絨的海島棉只占1%[1]。棉花主要生長在熱帶和亞熱帶地區(qū)，是一種對低溫敏感的作物[2]。隨著生產需求的不斷增加，棉花的種植區(qū)域漸漸轉移到了氣溫較寒冷的高緯度地區(qū)[3]，更易受到冷害的影響。

冷脅迫對植物形態(tài)的破壞主要表現為植株明顯矮化[4]，葉片失水萎蔫、表面出現黃斑或壞死，嚴重時可導致死亡[5]。植物遭受冷脅迫后，細胞膜的穩(wěn)定結構被破壞，流動性受阻，滲透壓發(fā)生改變，從而引起生長發(fā)育的失衡[6]。冷脅迫可以誘發(fā)植物細胞內產生大量活性氧（reactive oxygen species,ROS），使得植物內DNA、脂類等物質過氧化，破壞代謝平衡，因此，多數植物內部存在一套活性氧清除系統(tǒng)[7]。冷脅迫還會改變植物內源激素水平，從而引發(fā)一系列基因的表達形成反應機制，如脫落酸（abscisic acid,ABA）積累可以增加膜透性，誘導冷脅迫響應基因的表達[8]；生長素內源水平是植物耐冷性的另一考察指標，通常認為其內源含量與耐冷性呈負相關，這也許與生長素的極性運輸有關[9]。

目前，已有許多針對棉花冷脅迫響應相關基因的報道。ZHENG 等的研究表明，冷脅迫可以使棉花纖維長度顯著縮短，且可以誘導蘋果酸代謝、可溶糖代謝、纖維素合成等途徑基因的差異表達[10]。王俊娟通過對陸地棉抗冷品種的轉錄組測序，成功克隆到一個編碼脫水素蛋白的基因GhDHNA1，該基因在抗冷棉花品種和冷敏感品種中有顯著表達差異，且表達水平與耐冷性呈正相關[11]。LU等從棉花中鑒定到一個新的G蛋白受體基因，將其轉入擬南芥后發(fā)現轉基因株系在冷脅迫下表現出更高的萌發(fā)率[12]。HU 等對陸地棉TM?1 逆境處理材料進行轉錄組測序，共鑒定出8 725 個與耐冷性相關的基因，其中ABA信號轉導途徑相關基因，如ABA受體、SnRK蛋白和PP2C蛋白磷酸化激酶等在冷脅迫后顯著上調表達[13]。目前的研究更多地是對單個棉花品種冷脅迫響應的基因進行鑒定，而對造成棉種間耐冷性差異的分子機制的研究較少。已有的研究結果顯示，陸地棉和海島棉具有冷害耐受性差異[13]，而其背后的分子機制依然未知。本研究通過對陸地棉TM?1 和海島棉Hai7124 在冷脅迫處理后不同階段的葉片進行轉錄組測序，旨在揭示棉花受冷脅迫誘導的基因，并比較二者冷脅迫響應基因的分子功能和代謝途徑差異，以解析TM?1和Hai7124冷脅迫響應差異的分子基礎，為棉花耐冷和優(yōu)質新品種的培育提供基因資源。

1 材料與方法

1.1 棉花材料

2 種四倍體栽培棉花材料：陸地棉遺傳標準系TM?1和海島棉Hai7124。

1.2 棉花的冷脅迫處理與表型觀察

將黏土和蛭石按照體積比2∶1 混合，分裝于小號盆缽（直徑10 cm）中。用水濕潤后，選取長勢均勻、飽滿的棉花種子，在每個盆缽中播種3～4 粒。種子生長條件為溫度22 ℃，14 h 光照/10 h 黑暗，相對濕度70%。待幼苗生長均勻并萌發(fā)至兩葉一心時，分為2 組：一組為冷脅迫處理組，另一組為對照（CK）組。冷脅迫處理條件為生長于恒溫培養(yǎng)箱中，溫度4 ℃，14 h光照/10 h黑暗，相對濕度70%。對照組條件為生長于恒溫培養(yǎng)箱中，溫度28 ℃，其他條件與冷脅迫處理一致。

在冷脅迫施加后的不同時間點（0、1、3、6、12、24 h），分別取冷脅迫和對照組的5 株植株，進行快速拍照和表型觀察，每個時間點的植株在表型觀察后及時丟棄，不放回培養(yǎng)箱。

1.3 棉花的取樣與RNA 提取

在冷處理開始前，分別從處理組和對照組選取生長較為一致的健壯幼苗作為0 h的材料。冷處理開始后，在1、3、6、12 和24 h 時分別從處理組和對照組取3 株長勢一致且苗期棉花的第2 片真葉（即為3個生物學重復），立即放置在－80 ℃液氮中，用于RNA 的提取，以進行轉錄組測序。測序采用Illumina HiSeq 4000測序平臺（美國Illumina公司），生成150 bp 的雙端測序讀長?？俁NA 的提取使用RK16?50T 型植物RNA 快速提取試劑盒（南京鐘鼎生物技術有限公司），具體操作按使用說明書進行。

1.4 原始數據的處理和冷脅迫下受誘導的基因的鑒定

得到fastq格式的原始數據后，使用軟件fastp[14]過濾以獲得高質量讀長。使用HISAT×2（v2.1）[15]將高質量讀長比對到TM?1 參考基因組（v2.1）或Hai7124 參考基因組（v1.1）（http://cotton.zju.edu.cn/download）上，并基于基因轉移格式（gene transfer format, GTF）的注釋文件和DESeq2 軟件包鑒定差異表達基因（differentially expressed genes, DEGs），接著使用featureCounts 2.0.0[16]對基因表達進行定量，計算每百萬讀長中來自某一基因每千堿基長度的片段數目（expected number of fragments per kilobase of transcript sequence per million base pairs sequenced,FPKM），以評估基因的表達水平。差異基因的篩選標準為取偽發(fā)現率（false discovery rate,FDR）矯正后的P＜0.05 且差異表達倍數≥2（上調）或≤0.05（下調）的基因，上調和下調的DEGs均認為是受誘導的響應基因，并以FPKM=1 進行過濾，即FPKM＞1 的基因視為表達，FPKM＜1 的基因視為不表達?；虮倔w（gene ontology,GO）富集分析和京都基因與基因組百科全書（Kyoto Encyclopedia of Genes and Genomes, KEGG）富 集 分 析[17]使 用TBtools 1.046軟件[18]。轉錄組原始數據在本課題組之前的研究中已上傳至NCBI 數據庫（https://www.ncbi.nlm.nih.gov/）[13]，數據編號為PRJNA490626。

1.5 實時熒光定量反轉錄聚合酶鏈反應（real time quantitative reverse polymerase chain reaction,RT-qPCR）驗證

以Histone 3作為內參基因，引物為Y8991（F，5′?CGGTGGTGTGAAGAAGCCTCAT?3′；R，5′?AA TTTCACGAACAAGCCTCTGGAA?3′）。每個樣品設置3 個技術重復。使用Primer?Blast 工具（http://www.ncbi.nlm.nih.gov/tools/primer?blast）設計引物。引物序列如表1 所示。RT?qPCR 使用的是HiScript?ⅡQ RT SuperMix 試劑盒（南京諾唯贊生物科技股份有限公司），反應在不含RNase 的PCR反應管中進行：向反應管中加入RNA 模板（500 ng）、gDNA wiper 混合物（4 μL），并用無RNA 酶的雙蒸餾水補至16 μL，于42 ℃條件下反應2 min 后，加入HiScript?ⅡqRT SuperMix Ⅱ（4 μL），然后經過25 ℃引物結合10 min，50 ℃延伸30 min，85 ℃酶失活5 min 的程序，以此來合成cDNA 的第一條鏈。使用引物Y347 對得到的cDNA 進行PCR 擴增和電泳檢測，以驗證cDNA樣品的質量。

表1 用于定量檢測、基因克隆和載體構建的引物序列Table 1 Primer sequences for quantitative detection,gene cloning and vector construction

實時熒光定量PCR 使用的熒光染料是SYBR Green Master Mix（南京諾唯贊生物科技股份有限公司），在StepOnePlus 型PCR 儀（美國Thermo FisherScientific 公司）上進行反應，方法為2 步法PCR 循環(huán)，基因的相對表達水平通過2－ΔΔCT方法計算得到。

2 結果與分析

2.1 TM-1和Hai7124在冷脅迫下的響應表型觀察

在冷處理前，陸地棉TM?1葉片生長狀態(tài)良好。冷脅迫處理1 h后，未產生明顯的表型變化；3 h后，TM?1 的子葉表面較為平整，真葉葉尖出現輕微下垂；6 h 后葉片逐漸失去光澤，子葉在萎蔫的同時出現輕微卷縮，真葉下垂程度加劇，葉片出現向背面卷縮；12 h后整個植株萎蔫明顯，整個子葉從葉柄處完全下垂，真葉葉面的卷縮程度進一步加劇，但真葉仍能勉強保持大部分葉面平整；24 h后，整個植株萎蔫，真葉葉面不再平整，皺縮嚴重，幼嫩主莖開始出現較為嚴重的倒伏現象（圖1）。

海島棉Hai7124在4 ℃冷脅迫下的表型變化趨勢和TM?1 一致，但更快地出現了表型變化。在冷脅迫下處理3 h后，相比TM?1，Hai7124真葉萎蔫下垂的程度較為明顯；6 h后二者的表型差異進一步加劇，且Hai7124的子葉卷縮更加明顯，子葉下垂程度更加嚴重；12 h 后Hai7124 幼嫩的主莖出現輕微倒伏，真葉皺縮嚴重，但TM?1的幼嫩主莖仍能維持直立狀態(tài)，子葉仍能保持較為平整的葉面；24 h后整個Hai7124 植株的子葉葉柄下垂程度、葉片皺縮程度均明顯超過TM?1。通過對冷脅迫處理不同時間的表型比較發(fā)現，TM?1比Hai7124對冷脅迫的耐受力更強，冷脅迫表型更晚出現（圖1）。

圖1 TM-1和Hai7124棉花幼苗在冷脅迫下的表型比較Fig.1 Comparisons of phenotypes of TM-1 and Hai7124 cotton seedlings under 4 ℃treatment

采用Illumina HiSeq 高通量測序技術對冷脅迫處理和對照組的TM?1 和Hai7124 的幼苗進行轉錄組測序，結果顯示：平均每個樣本產生54.69百萬個讀長，轉錄組比對率在93.9%～97.8%之間。Q20值（堿基錯誤識別概率小于1%）均在98.5%以上，Q30值（堿基錯誤識別概率小于0.1%）均在96%以上。使用DESeq2軟件包比對處理組和對照組的測序數據，得到TM?1和Hai7124在冷脅迫處理后各個時間點的DEGs（圖2）。根據不同處理時間進行統(tǒng)計，得到TM?1 在冷脅迫處理后的1、3、6、12、24 h 分別有745、853、1 349、3 804、4 693 個DEGs，TM?1 響應冷脅迫的基因數隨處理時間延長呈現上升趨勢；Hai7124 在冷脅迫處理后的1、3、6、12、24 h 分別有113、3 936、3 085、5 674、4 301個DEGs，說明Hai7124在冷脅迫處理初始的反應比TM?1緩慢，而在3 h時有大量基因受誘導，說明3 h時Hai7124已經對冷脅迫做出強烈的反應，有較多基因差異表達。

圖2 冷脅迫下各個時間段TM-1和Hai7124的差異表達基因數Fig.2 Numbers of differentially expressed genes (DEGs) of TM-1 and Hai7124 at different stages under cold stress

一個基因如果在1、3、6、12、24 h中的任意一個時間點受到誘導，該基因即被視為總體上受誘導的基因。以此統(tǒng)計了TM?1 和Hai7124 中受誘導的基因數目。從圖3A可知：冷脅迫處理下TM?1中共有8 713 個基因受誘導，63 649 個基因不受誘導；Hai7124中共有12 461個基因受誘導，62 065個基因不受誘導。通過對TM?1 和Hai7124 之間的直系同源基因比較發(fā)現，兩者之間有3 884 個共同受誘導的基因，TM?1 擁有1 518 個特異受誘導的基因，而Hai7124擁有3 577個特異受誘導的基因（一些基因并沒有對應的直系同源基因）（圖3B）。

圖3 受誘導基因數Fig.3 Numbers of induced genes

對TM?1中1 518個特異受誘導的基因進行GO富集分析，結果顯示這些基因多富集在信號受體活性、光合作用、蛋白質代謝過程、細胞蛋白修飾過程、對光刺激反應、前體代謝產物和能量產生、脂質代謝過程等功能中（圖4A）；KEGG 富集分析顯示，這些基因主要集中在光合作用?天線蛋白、晝夜節(jié)律——植物、光合作用、信號轉導、植物激素信號轉導、環(huán)境信息處理、半胱氨酸和蛋氨酸的代謝、磷酸肌醇代謝、三羧酸循環(huán)等代謝途徑中（圖4B）。

圖4 TM-1中特異受誘導基因的GO和KEGG富集分析Fig.4 GO and KEGG enrichment analysis of specifically induced genes in TM-1

對Hai7124 中3 577 個特異受誘導的基因進行GO富集分析，結果顯示，這些基因多富集在水解酶活性、催化活性、葉綠體、質體、類囊體、脂質代謝過程等功能中（圖5A）；KEGG 富集分析顯示，這些基因主要集中在各種代謝過程中，如脂質代謝，甘氨酸、絲氨酸和蘇氨酸代謝，鞘脂代謝，β?丙氨酸代謝，組氨酸代謝，精氨酸和脯氨酸代謝，賴氨酸降解等（圖5B）。

圖5 Hai7124中特異受誘導基因的GO和KEGG富集分析Fig.5 GO and KEGG enrichment analysis of specifically induced genes in Hai7124

2.2 冷脅迫下脂肪酸代謝途徑種間差異比較分析

脂肪酸是植物細胞的必需組分之一，是植物激素等許多重要生物分子的中間合成體，也是細胞能量儲存的基礎。脂肪酸的含量和組成能夠影響細胞膜的流動性，從而影響植物對低溫環(huán)境的耐受能力。根據現有研究報道，我們從TM?1 和Hai7124里鑒定出了脂肪酸合成途徑關鍵家族中的乙酰輔酶A羧化酶（acetyl?CoA carboxylase,ACC）、脂肪酸合成酶（fatty acid synthase,FAS）、十八烷?ACP去飽和酶（stearoyl?ACP desaturase, SAD）、脂肪酸去飽和酶（fatty acid desaturase, FAD）、長鏈酰基輔酶A合成酶（long chain acyl?CoA synthetase, LACS）、酰基?ACP硫酯酶（acyl?ACP thioesterase,FAT）的基因并繪制了棉花脂肪酸代謝途徑的網絡圖（圖6～7）。在脂肪酸合成途徑關鍵基因家族中，在TM?1 中共鑒定到91 個基因，其中有17 個受冷脅迫誘導；Hai7124中共鑒定到84個基因，其中有24個受到冷脅迫誘導。在TM?1中，我們發(fā)現16個ACC家族基因中有2個基因受冷脅迫誘導；26個FAD家族基因中有10 個受冷脅迫誘導；6 個FAT 家族基因中有2個受冷脅迫誘導；24 個LACS 家族基因中有3 個受冷脅迫誘導。共有11 個基因在整個冷脅迫過程中呈現上調表達趨勢，也有一些基因展現出波動變化的表達模式（圖6）。

圖6 TM-1在冷脅迫下脂肪酸代謝通路圖和基因表達熱圖Fig.6 Diagram of fatty acid metabolism pathway and heat maps of expressed genes in TM-1 under cold stress

相比之下，在Hai7124中，我們發(fā)現了4個ACC家族基因、11個FAD家族基因、2個FAT家族基因、4個LACS家族基因以及3個SAD家族基因受冷脅迫誘導。在冷脅迫下，只有一小部分基因在某些時間段上調表達（如FAD家族基因），而有10個基因在整個冷處理過程中呈下調表達趨勢，也有一些基因呈波動表達模式（圖7）?？梢酝茰y，在冷脅迫下，TM?1中有更多的脂肪酸代謝途徑關鍵基因上調表達；而Hai7124 中則相反，包含更多下調表達的基因。這些受誘導的脂肪酸代謝途徑基因的差異表達模式可能造成了TM?1和Hai7124對冷脅迫的耐受差異。

圖7 Hai7124在冷脅迫下脂肪酸代謝通路圖和基因表達熱圖Fig.7 Diagram of fatty acid metabolism pathway and heat maps of expressed genes in Hai7124 under cold stress

從TM?1 和Hai7124 脂肪酸代謝通路中挑選了5 個冷脅迫響應基因進行RT?qPCR 驗證，結果（圖8）表明，GH_A13G2126（GhFAD6?1）、GH_D06G0458（GhLACS9）、GB_D13G2209（GbFAD6）、GB_D03G0951（GbLACS4）表 達 整 體 上 調，GH_D13G2107（GhSFD6?2）表達先上調后下調。

圖8 5個脂肪酸代謝途徑基因差異表達的RT-qPCR驗證Fig.8 RT-qPCR validation of five differentially expressed genes involved in fatty acid metabolism pathway

3 討論與結論

冷脅迫對棉花生產的影響巨大，我國有5%～12%種植面積的棉花是因為受冷害而補種或重播種的[19]。因此，解析陸地棉和海島棉冷脅迫響應差異形成的分子機制，對培育高產、高品質且抗逆性好的棉花新品種極為重要。已有的研究證明，陸地棉對冷害的耐受性更好。在本研究中，經過冷脅迫處理后3 h 時，陸地棉TM?1 和海島棉Hai7124 的表型變化開始出現差異，此時Hai7124 的真葉萎蔫更明顯，到6 和12 h 時表型差異更大，Hai7124 的葉片萎蔫和主莖倒伏更加劇烈，說明TM?1 比Hai7124對冷脅迫的耐受能力更強。

目前，通過同源克隆等方法，已經在棉花中鑒定出了許多冷脅迫相關響應基因，如GhERF1[20]、GhAGP[21]、GhFAD[22?23]等。然而棉花的耐冷、耐旱等抗逆性狀大多由多基因控制，無法通過單基因的變化來提升抗逆性，因此越來越多的研究人員致力于從全基因組或全轉錄組水平上對棉花的冷脅迫響應機制進行研究。ZHANG等從全基因組水平對陸地棉TM?1 進行正向選擇基因研究發(fā)現，馴化后的TM?1中，A亞基因組在纖維發(fā)育方面受到正向選擇，而D 亞基因組則在逆境耐受方面受到正向選擇[24]。LU等通過對亞洲棉、雷蒙德氏棉和陸地棉的全基因組MATE 基因家族進行分析，分別鑒定出了70、72和128 個MATE 家族基因，將基因Gh_D06G0281在擬南芥中過表達后，顯著提高了擬南芥植株對冷脅迫的抗性[25]。這些結果表明陸地棉經過正向選擇和人為馴化而擁有了更好的環(huán)境耐受力。本研究中，我們從轉錄組水平鑒定出了冷脅迫處理后TM?1和Hai7124 中總體受誘導的基因，并通過直系同源基因的比較鑒定出了兩者之中特異受誘導的基因。TM?1 和Hai7124 在冷脅迫處理后受誘導的基因隨處理時間增加而增加，Hai7124 在3、12 h 時間點的受誘導基因數量增加明顯。冷脅迫下，TM?1 中有1 518個特異受誘導的基因參與激素信號轉導、晝夜節(jié)律等調節(jié)過程，而Hai7124中有3 577個特異受誘導的基因參與賴氨酸、精氨酸、絲氨酸、蘇氨酸等氨基酸代謝，這反映了TM?1 和Hai7124 在冷脅迫下受誘導基因的分子功能和代謝途徑的差異。

除了冷脅迫響應的分子機制以外，棉花的耐冷研究還包括耐冷性的測定以及冷脅迫下的各類生理指標等[26]。細胞膜的滲透性和流動性在一定程度上反映了植物在冷脅迫下的受損害程度，也可以從側面反映植物對冷脅迫的耐受能力。而細胞膜的滲透性、流動性與膜脂中的不飽和脂肪酸比例密切相關。提高細胞膜內不飽和脂肪酸含量的一個重要方法就是飽和脂肪酸的去飽和反應，主要通過去飽和酶催化膜內脂肪酸來實現。在本研究中，篩選出TM?1 和Hai7124 里的冷脅迫相關基因后，我們著重關注了脂肪酸相關通路中的基因受冷脅迫誘導的表達模式差異：TM?1 脂肪酸途徑相關受誘導基因中共有11 個基因在整個冷處理過程中呈上調表達趨勢；而Hai7124 中共有10 個基因呈下調表達趨勢。TM?1 和Hai7124 這些脂肪酸途徑相關基因在冷脅迫下的差異表達，可能是造成其耐冷性差異的生物學基礎。

陸地棉和海島棉經過了同一起源、獨立馴化的過程，形成了兩者在產量、品質和適應性方面的差異。本研究從陸地棉TM?1 和海島棉Hai7124 這2個栽培品種中鑒定出冷脅迫響應基因，為未來棉花耐冷機制的研究提供了基因資源；同時解析了TM?1和Hai7124冷脅迫響應基因在分子功能、代謝途徑以及表達模式上的差異，這可能與兩者適應性差異的形成密切相關，然而還需要在未來進行進一步的功能驗證。本研究為創(chuàng)制轉基因新材料并培育耐冷、高產、優(yōu)質的棉花新品種提供了理論基礎。