湯亞琪，蘇文文，李 瑩，王保通，李 強

(西北農(nóng)林科技大學植物保護學院/旱區(qū)作物逆境生物學國家重點實驗室，陜西楊凌 712100)

小麥條銹病是世界范圍內(nèi)小麥主產(chǎn)區(qū)最重要的病害之一[1]，其病原菌條形柄銹菌小麥?；?PucciniastriiformisWest. f. sp.triticiEriks & Henn，Pst)為專性寄生菌[2]，生理小種致病性變異迅速且頻繁，容易造成主栽小麥品種抗條銹的頻繁“喪失”[3]，導致小麥條銹病的爆發(fā)流行，嚴重威脅著全球小麥的安全生產(chǎn)[4]。

中4是以普通小麥(TriticumaestivumL.)(2n=6x=42，AABBDD)與中間偃麥草(Thinopyrumintermedium)(2n=6x=42，E1E1E2E2XX)經(jīng)有性雜交、選育獲得的異源八倍體小偃麥(AABBDDXX)[5]。自1983年列為中國小麥條銹菌鑒別寄主至今，對我國現(xiàn)有小麥條銹菌流行小種均表現(xiàn)為抗病[6]，其抗性穩(wěn)定且持久，具有較高的研究和利用價值。但由于中4是異源八倍體小偃麥，與六倍體普通小麥經(jīng)常規(guī)雜交后獲得穩(wěn)定遺傳后代的難度較大，導致其抗性研究嚴重滯后[7]。藺瑞明等[7]將中4分別與中國春phb1突變體和條銹病感病品種銘賢169雜交，對其雜交后代進行遺傳分析，發(fā)現(xiàn)其對條銹菌小種CYR31和CYR32的抗病性由1對顯性基因控制。曹世勤等[8]利用基因推導法發(fā)現(xiàn)中4含有未知抗條銹病基因。甘肅省天水市農(nóng)業(yè)科學研究所以中4為親本，先后選育出中梁22號[9]、中93444[10]、中93447[11]等多個具有較好豐產(chǎn)性和抗病性的品種(系)材料。但是其抗條銹病機制目前尚未見報道。

本研究利用轉(zhuǎn)錄組測序，篩選中4接種條銹菌后的差異表達基因，以期為進一步解析中4抗條銹性分子機制提供科學依據(jù)。

1 材料與方法

1.1 供試材料

供試小麥品種為我國小麥條銹菌鑒別寄主中4和條銹病感病品種銘賢169，小麥條銹菌為我國小麥條銹菌流行小種CYR32，均由西北農(nóng)林科技大學植物保護學院小麥病原真菌監(jiān)測及抗病遺傳實驗室提供。

1.2 樣品采集

挑選籽粒發(fā)育良好、顏色形狀一致的中4種子，催芽后播種在7 cm × 7 cm的小花盆內(nèi)。待小麥幼苗長至第一葉完全展開時，接種條銹菌CYR32，接種0、12、24和48 h后，分別采集接種葉片樣品，以未接種葉片(0 h)為對照，每個時間點采集3個生物學重復，用于轉(zhuǎn)錄組數(shù)據(jù)質(zhì)量測定。接種48 h后的材料在人工氣候箱(溫度15～17 ℃、濕度90%、16 h光照/8 h黑暗、光照強度10 000 lx)繼續(xù)培養(yǎng)16 d，用于觀察表型，以銘賢169作為感病對照。

1.3 RNA提取及文庫構建

樣品總RNA的提取按照Plant RNA Midi Kit(OMEGA，美國)試劑盒說明書進行，獲取總RNA后，檢測RNA樣品的純度、濃度和完整性，用帶有Oligo(dT)的磁珠富集真核生物mRNA，加入Fragmentation buffer將mRNA隨機打斷，以First Strand cDNA Synthesis Kit(OMEGA，美國)合成第一條cDNA鏈，然后加入緩沖液、dNTPs、RNase H和DNA polymerase I合成第二條cDNA鏈，利用AMPure XP beads純化雙鏈cDNA，純化后的雙鏈cDNA須先進行末端修復、加A尾并連接測序接頭，再用AMPure XP beads進行片段大小選擇，最后通過PCR富集得到 cDNA文庫。

1.4 轉(zhuǎn)錄組測序及功能注釋

通過Illumina高通量測序獲得12個cDNA文庫，用FastQC對cDNA文庫的堿基質(zhì)量值(Q30)、堿基含量分布(GC含量)和重復相關性(r2)等指標進行檢測，以判斷文庫數(shù)據(jù)質(zhì)量并評估其是否滿足后續(xù)分析。用HISAT2對reads進行高效比對，用StringTie進行從頭組裝，基于中國春全基因組(https://wheat-urgi.versailles.inra.fr/Seq-Repository)進行有參組裝。

將完成組裝的轉(zhuǎn)錄本數(shù)據(jù)用非冗余(non-redundant，NR)數(shù)據(jù)庫和基因本體(gene ontology，GO)數(shù)據(jù)庫進行比對注釋，用真核生物同源蛋白簇(eukaryotic orthologous groups，KOG)、蛋白質(zhì)直系同源簇(cluster of orthologous groups of proteins，COG)以及基因的進化譜系：非監(jiān)督直系群(evolutionary genealogy of genes:non-supervised orthologous groups，eggNOG)等數(shù)據(jù)庫進行基因富集聚類，用京都基因與基因組百科全書(kyoto encyclopedia of genes and genomes，KEGG)數(shù)據(jù)庫進行代謝通路的富集分析。

1.5 實時熒光定量PCR

選取轉(zhuǎn)錄組測序得到的差異表達基因，根據(jù)其全長cDNA序列，使用Primer 5設計特異性引物(表1)，以小麥TaEF為內(nèi)參基因。提取中4各個樣品的總RNA并反轉(zhuǎn)錄為cDNA，按照ChamQTMSYBR○RqPCR Master Mix(Vazyme，中國)說明書的反應體系和條件，使用QuantStudio 5(ABI，美國)進行qRT-PCR，按照2-△△Ct計算基因相對表達量，重復3次。

表1 熒光定量PCR所用引物Table 1 Primers used for fluorescence quantitative PCR

2 結果與分析

2.1 中4抗條銹菌鑒定結果

從圖1可以看出，中4葉片只有少許壞死斑且不連片，表現(xiàn)為高抗，而銘賢169產(chǎn)生大量夏孢子堆。

2.2 轉(zhuǎn)錄組數(shù)據(jù)質(zhì)量

取接種條銹菌0、12、24和48 h后的中4葉片樣品，每個處理3個重復，共計12個樣品進行轉(zhuǎn)錄組測序分析，共獲得130.51 Gb Clean Data，各樣品Clean Data均達到8.42 Gb；堿基質(zhì)量值Q30的堿基百分比為92.83%～94.14%，表明堿基識別可靠，準確度高；G和C、A和T的含量為58.33%～60.83%。分別將各樣品的Clean Reads與中國春全基因組IWGSC RefSeq v1.0進行序列比對，比對效率為60.29%～69.03%，由于中4為異源八倍體小偃麥(AABBDDXX)[6]，其X染色體組與普通小麥親緣關系較遠，所以轉(zhuǎn)錄組的比對效率沒有出現(xiàn)高度匹配?；诒葘Y果進行基因表達量分析、功能注釋、富集分析以及新基因的發(fā)掘[8]，發(fā)掘新基因20 788個，其中 15 069個得到功能注釋。關于樣品相關性，r2值越接近1，說明兩個重復樣品相關性越強且樣品重復性良好，測序結果顯示，中4接種條銹菌CYR32 0 h后，3個樣品的r2值為0.960～ 0.983；接種條銹菌CYR32 12 h后，3個樣品的r2值為0.958～0.975；接種條銹菌CYR32 24 h后，3個樣品的r2值為0.972～ 0.982；而接種條銹菌CYR32 48 h后，3個樣品的r2值為0.824～0.932，表明樣品間整體相關性較強，重復性良好。

2.3 差異表達基因數(shù)目及染色體定位

以Fold Change≥2且FDR<0.01為篩選標準[12]，對不同時間點差異表達基因的分布進行調(diào)查，發(fā)現(xiàn)0 h VS 12 h(表示兩個侵染時間點差異表達基因的比較，下同)、0 h VS 24 h、0 h VS 48 h分別篩選出7 816、10 412、9 563個上調(diào)表達差異基因和7 209、10 159、9 612個下調(diào)表達差異基因；12 h VS 24 h、12 h VS 48 h、24 h VS 48 h分別篩選出937、1 150、57個上調(diào)表達差異基因和 2 866、 1 643、229個下調(diào)表達差異基因。結果表明，以0 h為對照，24 h獲得的差異表達基因最多，48 h次之，12 h最少；以12 h為對照，24 h獲得的差異表達基因多于48 h；48 h較24 h新增加了286個差異表達基因。

為獲得潛在的抗病相關基因，本研究重點關注上調(diào)表達差異基因，從圖2可以看出，共獲得13 674個上調(diào)表達差異基因，12 h與24 h、48 h共有的上調(diào)表達差異基因分別有6 241、5 611個，24 h與48 h共有的上調(diào)表達差異基因7 416個，0～48 h篩選到持續(xù)上調(diào)表達的差異基因5 151個，這些基因可能與中4的的抗條銹性有關。對5 151個上調(diào)差異基因進行染色體定位，其中 4 482個獲得定位信息，其余669個無法確定其染色體位置，為新基因。進一步分析發(fā)現(xiàn)，差異表達基因在B染色體組上分布最多，為1 581個；D染色體組次之，為1 523個；A染色體組最少，為 1 378個。上調(diào)表達差異基因在2B(7.12%)、5B (6.38%)、3B(6.34%)、2D(5.85%)、7D(5.42%)染色體上所占比例較高。

2.4 共同上調(diào)表達差異基因的分類注釋及富集

2.4.1 共同上調(diào)表達差異基因的GO功能注釋

中4轉(zhuǎn)錄組測序數(shù)據(jù)中有58 059個轉(zhuǎn)錄本在GO數(shù)據(jù)庫中被注釋，5 151個共同上調(diào)表達差異基因中3 411個基因獲得GO注釋。獲得注釋的轉(zhuǎn)錄本可分為分子功能、細胞組分和生物學過程3大類，46個亞類(圖3)。分子功能包含13個亞類，其中蛋白結合、催化活性通路富集到的上調(diào)表達差異基因數(shù)目較多，此外，在轉(zhuǎn)運活性、受體活性、酶調(diào)節(jié)活性、抗氧化活性、蛋白結合轉(zhuǎn)錄因子活性等通路也富集到上調(diào)表達差異基因；細胞組分包含13個亞類，其中細胞部件、細胞、細胞器等通路富集到的上調(diào)表達差異基因數(shù)目較多，在大分子復合物、細胞連接、超分子復合物等通路也富集到上調(diào)表達差異基因；生物學過程包含20個亞類，其中細胞過程、代謝過程、單有機體過程等通路富集到的上調(diào)表達差異基因數(shù)目較多，此外，在生殖過程、多有機體過程、信號轉(zhuǎn)導等通路也富集到上調(diào)表達差異基因。

2.4.2 共同上調(diào)表達差異基因的功能富集分析

利用COG、KOG、eggNOG數(shù)據(jù)庫對5 151個上調(diào)表達差異基因分別進行直系同源分類，以期獲得更準確的分類及功能注釋。COG數(shù)據(jù)庫中有2 638個上調(diào)表達差異基因獲得注釋分類，最大類別為R類(僅一般功能預測)，包含555個基因，占21.04%；其次是K類(轉(zhuǎn)錄)，包含306個基因，占11.60%(圖4A)；KOG數(shù)據(jù)庫中有 2 486個上調(diào)表達差異基因獲得注釋分類，最大類別為R類，包含436個基因，占17.54%，其次是O類(翻譯后修飾，蛋白質(zhì)轉(zhuǎn)換，分子伴侶)，包含271個基因，占10.90%(圖4B)；eggNOG數(shù)據(jù)庫中有3 906個上調(diào)表達差異基因獲得注釋分類，最大類別為S類(未知功能)，包含859個基因，占21.99%，其次是R類，包含666個基因，占 17.05%(圖4C)。此外，COG、KOG、eggNOG三個數(shù)據(jù)庫分別注釋到與防御相關的V類(防御機制)上調(diào)表達差異基因，分別為47、28、25個。雖然eggNOG數(shù)據(jù)庫獲得的注釋基因最多，但是注釋到與防御相關的基因最少，在未知功能和一般功能預測富集到的基因最多，為1 525個，占 39.04%；KOG數(shù)據(jù)庫在所占比例為 5.83%～ 5.39%處出現(xiàn)分類富集，且數(shù)據(jù)庫25個分類中均注釋到上調(diào)表達差異基因；COG和KOG數(shù)據(jù)庫雖然在最大類別注釋相同，但是其余類別差異很大。三個數(shù)據(jù)庫注釋分類結果有較大差異，因此，應綜合分析三個數(shù)據(jù)庫彌補單一數(shù)據(jù)庫的不足。

2.4.3 共同上調(diào)表達差異基因的KEGG通路富集分析

5 151個差異表達基因在KEGG數(shù)據(jù)庫中注釋到431個上調(diào)表達差異基因，主要分為5大類20個小類，共涉及109個代謝通路。其中，細胞過程、環(huán)境信息處理、遺傳信息加工、代謝和有機系統(tǒng)5大類注釋到的上調(diào)表達差異基因數(shù)目分別為39、24、146、203和19。有機系統(tǒng)類別注釋到14個參與植物與病原菌互作通路的差異表達基因。將顯著性設為qvalue<0.05，對上調(diào)差異基因進行富集分析，發(fā)現(xiàn)差異基因在剪接體、苯丙氨酸代謝、過氧化物酶體、糖和氨基酸的代謝等通路顯著富集 (圖5)。

2.5 R基因和新基因的發(fā)掘

基于所選參考基因組序列的比對結果，共注釋到包括7種類型的33個R(resistance)基因，其中，12個與富含半胱氨酸類受體蛋白激酶有關，10個與抗病性蛋白有關(表2)。此外，還篩選到β-半乳糖苷酶、纖維素合酶、過氧化物酶等病程相關蛋白。

表2 中4抗條銹性上調(diào)差異基因中R基因的功能分析Table 2 Functional analysis of R gene in up-regulated differential genes of Zhong 4 resistance to stripe rust

將篩選到的669個上調(diào)表達的差異新基因在7種數(shù)據(jù)庫注釋分析，結果顯示，有335個獲得GO注釋，123個獲得COG注釋，79個獲得KOG注釋，511個獲得eggNOG注釋；其中，eggNOG數(shù)據(jù)庫最大的兩類注釋分別為功能未知和一般功能預測，占48.57%，KOG數(shù)據(jù)庫中細胞周期控制、細胞分裂和染色體分裂通路中獲得注釋信息的基因比例較高。151個上調(diào)表達差異新基因獲得KEGG注釋，且其中85個差異新基因共涉及89個信號通路；562個獲得NR注釋。其中，苯丙素的生物合成和苯丙氨酸代謝通路分別富集到7個上調(diào)表達差異基因，谷胱甘肽代謝通路分別富集到6個，內(nèi)質(zhì)網(wǎng)蛋白加工和碳代謝通路分別富集到5個，植物與病原菌互作、剪接體、糖酵解和糖異生以及淀粉和蔗糖代謝通路分別富集到4個，說明顯著富集的差異新基因主要與糖的合成與代謝、氨基酸的代謝以及植物與病原菌互作有關。上調(diào)表達差異新基因還注釋到R基因、WRKY 和NAC轉(zhuǎn)錄因子等應答非生物脅迫的基因。

2.6 qRT-PCR對RNA-seq數(shù)據(jù)的驗證

為驗證轉(zhuǎn)錄組測序RNA-seq數(shù)據(jù)的準確度，在差異基因中隨機選取與NB-ARC、鉀離子轉(zhuǎn)運蛋白和熱激蛋白相關的3個基因進行qRT-PCR，其中，轉(zhuǎn)錄本數(shù)據(jù)經(jīng)NCBI Blast確定其特異序列位置，用Premier 6設計與內(nèi)參引物TaEF-1α 退火溫度值相近的特異引物。以0 h為對照，計算基因相對表達水平與轉(zhuǎn)錄組數(shù)據(jù)表達量并進行對比。結果(表3)表明，轉(zhuǎn)錄組測序RNA-seq數(shù)據(jù)與qRT-PCR數(shù)據(jù)整體趨勢一致，轉(zhuǎn)錄組測序數(shù)據(jù)準確可信。

表3 熒光定量PCR結果與轉(zhuǎn)錄組測序結果對比Table 3 Comparison of real-time PCR results with RNA-seq results

3 討 論

條銹病的流行嚴重威脅小麥的產(chǎn)量和品質(zhì)，培育抗病品種是防治該病最有效、經(jīng)濟和環(huán)保的措施[1-4]。中4自育成至今，對我國現(xiàn)有條銹菌小種均表現(xiàn)為抗病，然而其抗條銹病機制尚不明確。

本研究利用轉(zhuǎn)錄組測序，分析了條銹菌CYR32侵染中412、24和48 h后的差異表達基因，獲得 5 151個持續(xù)上調(diào)表達差異基因。經(jīng)多個數(shù)據(jù)庫注釋及富集分析，篩選得到R基因、PR基因、蛋白酶體、WRKY和NAC轉(zhuǎn)錄因子等與調(diào)控小麥條銹菌有關的基因。其中33個上調(diào)表達差異基因參與植物與病原菌互作過程；大量上調(diào)表達差異新基因獲得注釋，在植物與病原菌互作過程、糖的合成與代謝以及氨基酸的代謝通路顯著富集。

Chen等[13]通過轉(zhuǎn)錄組學分析小種?；钥剐浴⑷诳剐院透邷爻芍昶诳剐?，發(fā)現(xiàn)高溫成株期抗性中基因多樣性最高，但是基因誘導表達緩慢，解釋了其抗性持久性的原因。Coram等[14]分析Yr39介導的與條銹菌親和與不親和的兩個F7重組自交系的轉(zhuǎn)錄本，發(fā)現(xiàn)條銹菌侵染48 h后時，轉(zhuǎn)錄本數(shù)量均達到峰值；對轉(zhuǎn)錄組數(shù)據(jù)進行注釋后發(fā)現(xiàn)，50%以上的轉(zhuǎn)錄本參與防御和信號轉(zhuǎn)導過程，例如R基因同源物、PR蛋白、苯丙氨酸的生物合成和蛋白激酶等。Coram等[15]分析Yr5介導的與條銹菌親和與不親和的相關轉(zhuǎn)錄本，發(fā)現(xiàn)在侵染24 h時轉(zhuǎn)錄本積累量出現(xiàn)峰值，注釋后發(fā)現(xiàn)被誘導的轉(zhuǎn)錄本均參與基礎防御相關過程，Yr5快速響應的信號通路和防御反應有R基因參與，包括蛋白激酶信號轉(zhuǎn)導和活性氧導致的過敏性壞死反應。本研究結果發(fā)現(xiàn)，條銹菌侵染24 h時，差異表達基因數(shù)目出現(xiàn)峰值，表明小麥在應對條銹菌的誘導下，大量與寄主抗性及相關抗性基因轉(zhuǎn)錄表達來調(diào)控免疫反應。Hao等[16]研究發(fā)現(xiàn)，在條銹菌侵染24 h時，KEGG顯著富集的通路有苯丙氨酸的生物合成和代謝，該途徑可以產(chǎn)生水楊酸、類黃酮和木質(zhì)素等化合物，并且苯丙氨酸代謝途徑在條銹菌侵染初期參與植物防御反應；在條銹菌侵染120 h時，KEGG富集到與淀粉、蔗糖和碳水化合物代謝相關的基因，推測小麥光合作用的強度影響其對條銹菌的抗性。本研究KEGG通路分析發(fā)現(xiàn)，上調(diào)表達差異基因顯著富集在蛋白質(zhì)加工、苯丙氨酸代謝、過氧化物酶體、糖和氨基酸代謝等通路，與前人研究結果一致。Garnica等[17]研究發(fā)現(xiàn)，ABC轉(zhuǎn)運蛋白基因在孢子萌發(fā)中差異表達，并在附著期高度富集；MFS轉(zhuǎn)運蛋白優(yōu)先在孢子萌發(fā)期表達，主要影響其生長過程；己糖轉(zhuǎn)運蛋白和氨基酸轉(zhuǎn)運體在吸器中高度表達，條銹菌的生活方式對轉(zhuǎn)運體類型起決定作用。葡萄糖轉(zhuǎn)運體在發(fā)芽夏孢子中高表達，說明條銹菌侵染小麥后消耗寄主碳水化合物吸收營養(yǎng)。

PR蛋白在系統(tǒng)獲得抗性和植物防御中尤為重要。研究表明，PR1蛋白能夠抑制小麥條銹菌的生長[18]。另外，TaWRKY49負調(diào)控小偃6號對小麥條銹菌的抗性，并誘導水楊酸、茉莉酸和活性氧響應相關基因的表達，同時抑制乙烯相關基因的表達；TaWRKY45和aWRKY62正調(diào)控小偃6號對小麥條銹菌的抗性，且水楊酸、乙烯、茉莉酸甲酯、過氧化氫和脫落酸均能夠誘導TaWRKY62基因表達[19-20]。NAC基因家族編碼的蛋白質(zhì)是目前發(fā)現(xiàn)的最大的植物特異性轉(zhuǎn)錄因子之一[21-23]，在植物應答非生物和生物脅迫以及調(diào)控植物發(fā)育等方面發(fā)揮著重要的作用，TaNAC4[21]、TaNAC21/22[22]和TaNAC30[23]負調(diào)控小麥對條銹菌的抗性，且茉莉酸甲酯、乙烯和脫落酸能夠誘導TaNAC4基因的表達，而水楊酸對TaNAC4基因的表達無影響[21]；TaNAC21/22作為tae-miR164的靶標在小麥對條銹菌的抗性方面發(fā)揮著重要作用[22]；沉默TaNAC30可顯著增加H2O2的積累量，說明TaNAC30可增加小麥對條銹菌的抗性[23]。

本研究分析了條銹菌CYR32侵染中4的轉(zhuǎn)錄組數(shù)據(jù)，通過基因組序列比對獲得7種類型R基因并得到大量的新基因，新基因中也預測到R基因和抗條銹性相關基因，表明小麥品種中4中蘊含豐富的抗性基因。