青稞響應燕麥鐮孢菌侵染的轉錄組分析

漆永紅，李雪萍，曹素芳，李敏權

(1.甘肅農業大學 草業學院/草業生態系統教育部重點實驗室/甘肅省草業工程實驗室/中-美草地畜牧業 可持續發展研究中心，甘肅 蘭州 730070； 2.甘肅省農業科學院 植物保護研究所， 甘肅 蘭州 730070； 3.甘肅省農業科學院 林果花卉研究所，甘肅 蘭州 730070)

青稞(Hordeumvulgarevar.nudum)是栽培大麥的一種變種，廣泛用于食品、飼料和釀造等[1-2]。青藏高原既是世界大麥起源中心之一，又是我國重要的裸大麥(青稞)種植區[3]。李雪萍等[4-5]對青稞根際土壤真菌和細菌等微生物的數量、土壤酶活性與根腐病的關系進行了研究。由燕麥鐮孢菌(Fusariumavenaceum)等引起的青稞莖基腐病是青稞苗期和成株期的主要病害之一[6]。國內一些學者對其他作物病害的轉錄組學開展了研究，柳鳳等[7]報道了杧果與杧果鐮孢菌(F.mangiferae)在轉錄水平互作過程中，鈣信號傳導、SA 信號途徑和絲裂原活化蛋白信號傳導途徑相關基因表達下調，造成了下游植物抗病基因RPM1 的表達受到抑制。葉文武等[8]首次報道通過轉錄組學方法揭示了特定的功能基因PsMPK1 沉默后大豆疫霉(Phytophthorasojae)基因組的轉錄變化，鑒定了可能與PsMPK1信號途徑相關的基因和基因家族。邱化榮等[9]基于轉錄組分析了蘋果水楊酸特異響應基因MdWRKY40的啟動子鑒定，發現SA處理影響的基因主要參與了苯丙烷類、類黃酮的生物合成，植物病原菌互作及植物激素信號傳導途徑。汪平等[10]對煙草疫霉(P.nicotianae)侵染前后劍麻葉片轉錄組學進行了研究。徐金青等[11]對青稞轉錄組SSR位點及其基因功能進行了分析，通過基因功能注釋發現青稞轉錄組中含SSR的序列主要與生物的基礎代謝相關。關于青稞與燕麥鐮孢菌在轉錄水平上的互作機制國內外報道較少。以接種燕麥鐮孢菌后抗感青稞的莖基部為材料，利用Illumina HiSeq Xten平臺對其轉錄組信息進行分析，通過對獲得的轉錄本序列進行基因功能注釋和分類，在轉錄水平上初步揭示青稞與燕麥鐮孢菌互作的分子機制，為抗病新基因的挖掘以及青稞抗病品種的選育提供理論依據。

1 材料和方法

1.1 供試材料

供試菌株：燕麥鐮孢菌[6]，經單孢分離、純化和鑒定，以常規方法保存于PDA 斜面上，接種前進行活化，來自甘肅省農業科學院植物保護研究所菌種保存室。

供試青稞品種：感病品種甘青2號，來自甘南州農業科學研究院；抗病品種NQK-01-03，來自青海省農林科學院植物保護研究所。

轉錄組測序樣品編號：抗病品種NQK-01-03，對照CK分別為T01、T02、T03，20 d分別為T04、T05、T06，40 d分別為T07、T08、T09；感病品種甘青2號，對照CK分別為T10、T11、T12，20 d分別為T13、T14、T15，40 d分別為T16、T17、T18。

1.2 試驗設計

將表面消毒的甘青2號和NQK-01-03青稞品種種植于無菌蛭石中，每盆100株。在人工氣候箱每天光照培養8 h，溫度25℃，相對濕度60%；黑暗培養16 h，溫度20℃，相對濕度80%。待青稞幼苗長至5～6葉期，在其根部接種孢子濃度為1×106的燕麥鐮孢菌50 mL，以無菌水接種作為空白對照。接菌前采集兩個青稞品種的莖基部作為轉錄測序的對照(CK)，在接菌20、40 d分別采集發病的兩個青稞品種的莖基部，將莖基部清水洗凈后，在液氮中速凍，于-80℃保存備用。

1.3 RNA提取

參照Lexogen公司SPLIT RNA試劑盒(北京賽百盛基因技術有限公司)的步驟提取青稞莖基部RNA，利用瓊脂糖凝膠電泳和Agilent 2100 檢測RNA質量及濃度。

1.4 轉錄組Illumina HiSeq Xten測序

試驗由北京百邁客生物科技有限公司完成。

1.5 數據處理

經Illumina HiSeq Xten測序獲得原始序列，去除以下類別的低質量序列：去除含有接頭的Reads 、去除低質量的Reads(包括去除N的比例大于10%的Reads；去除質量值Q≤10的堿基數占整條Read的50%以上的Reads)，經過上述一系列的質量控制之后得到高質量的Clean Data。

1.6 生物信息學分析

將下機數據進行過濾得到Clean Data，與指定的參考基因組(大麥和燕麥鐮孢菌)進行序列比對，得到的Mapped Data進行插入片段長度檢驗、隨機性檢驗等文庫質量評估。根據基因在不同樣品或不同樣品組中的表達量進行差異表達基因Unigene的注釋、差異表達基因KEGG注釋和GO分析。

1.7 新基因功能注釋

基于所選參考基因組序列，使用Cufflinks軟件對Mapped Reads進行拼接，并與原有的基因組注釋信息進行比較，尋找原來未被注釋的轉錄區，發掘該物種的新轉錄本和新基因，從而補充和完善原有的基因組注釋信息。過濾掉編碼的肽鏈過短(少于50個氨基酸殘基)或只包含單個外顯子的序列。使用BLAST[12]軟件將發掘的新基因與NR[13]，Swiss-Prot[14]，GO[15]，COG[16]，KOG[17]，Pfam[18]，KEGG[19]數據庫進行序列比對，使用KOBAS 2.0[20]得到新基因的KEGG Orthology結果，預測完新基因的氨基酸序列之后使用HMMER軟件[21]與Pfam數據庫比對，獲得新基因的注釋信息。

2 結果與分析

2.1 RNA質量檢測

由試劑盒提取的RNA，經電泳及凝膠成像觀察，發現18S rRNA 和28S rRNA條帶清晰且完整，無拖帶(圖1)。經檢測，D260/280比值均在1.8～2.2，表明所提取的RNA無降解，符合測序和建庫需求。

圖1 部分樣品RNA電泳圖Fig.1 Some samples RNA electrophoretogram 注：M：RNA MarkerⅠ，1～5部分樣品RNA，CK為空白對照

2.2 Illumina HiSeq Xten測序結果

抗病青稞NQK-01-03和感病青稞甘青2號在接菌前、接菌20 d 和40 d，經測序質量控制，共得到112.97 Gb Clean Data，各樣品Q30堿基百分比均高于88.00%，GC含量均高于53.00%(表1)。

表1 Illumina HiSeq Xten測序

2.3 DEG Set差異表達基因

抗病青稞NQK-01-03和感病青稞甘青2號接菌前、接菌20 d 和40 d共獲得DEGs 差異表達基因4 280個，其中上調基因2 037個，下調基因2 243個(表2)。接菌前、接菌20 d 和40 d抗感青稞品種之間

差異表達基因數目分別為663，2 880和737個，上調基因數目分別為349，1 274和414個，下調基因數目分別為314、1 606和323個，表明接菌后抗感青稞品種之間上調和下調差異表達基因數目均增加，其中接菌20 d上調和下調差異表達基因數最多。

表2 DEG Set差異表達基因數目

2.4 差異表達基因Unigene注釋

使用BLAST 軟件與COG、GO、KEGG、KOG、NCBI-NR、Pfam、Swiss-Prot 和eggNOG數據庫比對，共計28 869條Unigene 獲得注釋信息(表3)。抗病青稞NQK-01-03和感病青稞甘青2號之間在接菌前、接菌20 d 和40 d NCBI-NR數據庫中注釋到的Unigene 最多，達到3 981條，占全部注釋基因的13.79%；其次是Pfam、GO和Swiss-Prot，分別為2 946，2 931和2 672條，占全部注釋基因的10.20%、10.15%和9.26%；而eggNOG數據庫注釋到的Unigene 最少，為447條，僅占全部注釋基因的1.55%。與接菌前和接菌40 d相比，抗病青稞NQK-01-03和感病青稞甘青2號之間在接菌20 d各數據庫注釋到的Unigene差異表達基因數量最多，達到2 728條。

表3 注釋的差異表達基因Unigene數量

2.5 差異表達基因KEGG通路注釋

將DEGs與KEGG數據庫進行比對，找到DEGs 中抗病青稞NQK-01-03和感病青稞甘青2號之間在接菌前、接菌20 d和40 d中差異表達基因的KEGG通路60條，其中顯著性Q-value≤0.05為閾值的代謝途徑有5條，分別為苯丙氨酸代謝、異喹啉類生物堿的生物合成、酪氨酸代謝、苯丙基類生物合成和光合作用天線蛋白(表4) 。接菌前，抗感青稞品種間差異表達基因的KEGG通路富集20條，但Q-value≤0.05的沒有；接菌20 d，抗感青稞品種間差異表達基因的KEGG通路富集在苯丙氨酸代謝、異喹啉類生物堿的生物合成和酪氨酸代謝3種途徑；接菌40 d，兩者之間的KEGG通路富集在苯丙基類生物合成和光合作用天線蛋白兩種途徑，表明在接種燕麥鐮孢菌后抗感青稞品種間差異表達基因的KEGG通路富集發生了變化。抗感青稞之間差異表達基因KEGG通路參與的生理和信號傳導方式有3類，分別參與代謝、合成和光合過程，其中與代謝有關的2條，包括酪氨酸代謝和苯丙氨酸代謝。與合成有關的2條，包括苯丙基類生物合成和異喹啉生物堿合成。與光合有關的1條，光合作用天線蛋白。

表4 差異表達基因KEGG通路富集

2.6 差異表達基因GO分析

抗病青稞NQK-01-03和感病青稞甘青2號之間在接菌前、接菌20 d 和40 d差異表達基因GO的平均基因數為1.3、9.3、97和977個，在0.1%、1%、10%和大于10%所占的比率分別為4.76%、3.24%、3.37%和3.39%(表5)。接菌前和接菌40 d，抗感青稞之間差異表達基因GO的數目及比率接近；而接菌20 d，抗感青稞之間差異表達基因GO的數目及比率最高。

抗病青稞NQK-01-03和感病青稞甘青2號之間GO注釋系統包含生物學過程、細胞組分和分子功能3個主要分支(圖2)。生物學過程包含多細胞生物過程、多有機體過程、免疫系統過程、細胞成分組織或生物發生、晝夜節湊過程等20個條目，“生物學過程”富集的GO分析主要包括代謝過程、 細胞進程、單生物過程、刺激反應和生物調節。細胞組分包含細胞外基質部分、核狀體、高分子復合物、膜封閉腔、膠原三聚體等17個條目，“細胞組分”富集的GO分析主要包括細胞部、細胞、細胞器、膜和細胞器部分。分子功能包含酶調節活性、抗氧化劑活性、蛋白結合轉錄因子活性、鳥苷核苷酸交換因子活性、通道調節活性、通道調節活性、翻譯調節活性、金屬血紅素活性和蛋白標記等17個條目，“分子功能”富集的GO分析主要包括催化活性、轉運體活性、結構分子活性、核酸結合轉錄因子活性和電子載體。

表5 差異表達基因GO的數目及比率

2.7 新基因功能注釋

共發現新基因數目3 095個，各數據庫注釋的新基因數量和所占的比率表明COG、KEGG、KOG和Swiss-Prot注釋的新基因數目分別為287、407、704和806個，所占比例9.27%～26.04%；Pfam和GO注釋的新基因數目分別為1 021和1 194個，所占比例分別為32.99%和38.58%；eggNOG和NCBI-NR注釋的新基因數目最多分別達到1 833和2 396個，所占比例分別為59.22%和77.42%(表6)。

表6 新基因功能注釋

3 討論

植物與病原菌互作的研究，早期主要采用顯微技術觀察病原菌侵染后兩者超微結構的變化、病原菌侵染對寄主生理代謝的影響和植物的誘導抗性等。隨著高通量測序技術的發展，從轉錄水平上探索植物與病原菌互作的分子機制已成為該項研究的主要方法[22]。轉錄組學是從RNA水平上研究細胞整個基因組的表達情況，能夠全面獲得植物與病原菌互作的分子信息，篩選出兩者之間DEGs，研究病原菌在不同環境下基因表達變化，進而摸索出植物對病原菌的影響等[23]。轉錄組測序技術適用于缺乏基因組序列數據的非模式生物基因表達圖譜的構建，被廣泛應用于解析生物性狀形成的分子機制、發掘基因資源、植物與病原的互作等方面[7]。Li等[24]利用Solexa第2代測序技術研究了水稻和稻瘟病菌(Magnaportheoryzae)互作早期轉錄組的測定和分析。本試驗首次以接種燕麥鐮孢菌后的抗感青稞品種莖基部為材料，利用Illumina HiSeq Xten平臺對其轉錄組信息進行分析，對獲得的轉錄本序列進行基因功能注釋和分類，在轉錄水平上初步揭示了青稞與燕麥鐮孢菌互作的分子機制。

柳鳳等[7]采用高通量測序技術比較了杧果病健組織轉錄組信息的差異，發現了糖、核酸和蛋白質、活性氧、內源激素和酚類等多種代謝途徑，參與抗氧化、離子運輸、信號轉導和逆境響應等生物進程。董璐[25]在不同抗感品種菊花的白色銹病病原菌(Pucciniahoriana)的轉錄組和表達譜研究中，發現了與戊糖磷酸途徑、苯丙烷生物合成、苯丙氨酸的代謝及碳代謝等相關的差異表達基因。試驗研究發現，青稞接種燕麥鐮孢菌后，抗感青稞品種間差異表達基因的KEGG通路富集中在苯丙氨酸代謝、異喹啉類生物堿的生物合成、酪氨酸代謝、苯丙基類生物合成和光合作用天線蛋白。

Ward和Weber[26]采用RNA-Seq 技術對覆盆子抗感品系感染根腐病病原菌(Phytophthorarubi)后進行轉錄組測序，結果發現WRKY 轉錄因子和病程相關蛋白基因及與三羧酸循環和木質素合成途徑相關的基因表達量都上調有關。Bagnaresi等[27]利用RNA-Seq技術對抗稻瘟病(M.oryzae)水稻品種Gigante vercelli(GV)和感病品種Vialone nano(VN)進行轉錄組分析，發現GV中編碼植保素合成酶、含黃素單加氧酶、糖基水解酶以及幾丁質酶的基因特異性上升。本試驗僅獲得了抗病青稞NQK-01-03和感病青稞甘青2號接菌前、接菌后的DEGs 差異表達上調基因和下調基因數目，為了進一步驗證本轉錄產生的差異表達基因的實際功能，應當從不同的代謝途徑中隨機挑選一些基因，通過設計不同的引物進行qRT-PCR驗證，這部分工作有待于進一步完善。

試驗發現，抗病青稞NQK-01-03和感病青稞甘青2號之間在接菌20 d時的差異表達基因DEGs數目、Unigene數目和GO數目高于接菌40 d，表明這些差異表達的基因數量與接種的時間有密切的關系。通過盆栽接種來采集發病的青稞莖基部開展研究，在開展植物與病原在轉錄水平互作研究時，通過刺傷植物的莖基部接種既可以縮短侵染時間，又可能得到更好的差異表達基因效果。轉錄組測序技術對植物寄主與病原菌互作時基因的差異表達情況研究，可獲得大量可靠的與抗病有關的基因資源，為深入了解植物的抗病機制以及改良植物抗病性的分子育種奠定基礎[28]。試驗發掘了一些新基因，這些新基因的序列信息對進一步分析抗病基因結構和功能提供理論依據和奠定了基礎，有關這些新基因的功能和作用機制有待于進一步研究。

青稞與燕麥鐮孢菌的互作是一個復雜的過程，僅通過轉錄組分析兩者間的互作機制是不夠的，應綜合運用生理代謝、生物化學和分子生物學等手段進行系統分析，這將有助于闡明青稞與燕麥鐮孢菌的互作機制，建立青稞防御燕麥鐮孢菌的作用模型，為青稞莖基腐病的防治提供一種新的防治策略。

4 結論

抗病青稞NQK-01-03和感病青稞甘青2號之間經過測序質量控制，共得到112.97Gb Clean Data。DEGs基因表達注釋共獲得差異表達上調基因2 037個，下調基因2 243個。NCBI-NR數據庫中注釋到的Unigene數量最多，占全部注釋基因的13.79%。差異表達基因的KEGG通路富集Q-value≤0.05主要有5條，GO注釋系統包含3個主要分支。eggNOG和NCBI-NR注釋的新基因數目最多，所占比例分別為59.22%和77.42%。