花魔芋響應(yīng)軟腐病菌侵染初期的轉(zhuǎn)錄組分析

何 斐馮明月王菁軒邵永春王 甜王 瑩盧謝敏方姝穎

(安康學(xué)院 現(xiàn)代農(nóng)業(yè)與生物科技學(xué)院,陜西安康 725000)

魔芋系天南星科魔芋屬(Amorphophallus)多年生草本植物,它是目前為止自然界唯一能大量提供葡甘聚糖的重要經(jīng)濟(jì)作物[1]。這種天然高分子優(yōu)質(zhì)膳食纖維-葡甘聚糖,具有凝膠性、增稠性、持水性和成膜性等獨(dú)特的特性,被廣泛應(yīng)用于生物醫(yī)藥、食品和化工等領(lǐng)域[2]。由于其較高的經(jīng)濟(jì)價(jià)值,近年來魔芋已發(fā)展成為云南、貴州、四川、陜南及湖北等山區(qū)農(nóng)民脫貧致富的支柱產(chǎn)業(yè)[3]?；в笫侵袊闹髟云贩N,其產(chǎn)量大、葡甘聚糖含量高,但抗病性較差[4]。其中,土傳和種傳細(xì)菌性軟腐病是花魔芋生產(chǎn)過程中危害最大且極難防控的病害,在植株生長期和采后貯藏期均可發(fā)病,可危害花魔芋葉片、葉柄及地下塊莖,輕則造成減產(chǎn)20%～30%,重則導(dǎo)致80%以上減產(chǎn)甚至絕收[4],嚴(yán)重制約魔芋產(chǎn)業(yè)的可持續(xù)發(fā)展。隨著科技的不斷進(jìn)步,利用抗病基因進(jìn)行分子育種已成為解決作物抗病性的有效途徑,這一途徑的關(guān)鍵是對(duì)花魔芋新型內(nèi)源抗病基因的挖掘和篩選以及進(jìn)一步對(duì)其抗病防御機(jī)制的解析。

目前,對(duì)花魔芋軟腐病相關(guān)的研究雖然已有一定報(bào)道,但是主要集中在病原菌鑒定與檢測[3]、病害防治[5]及傳統(tǒng)的抗病品種選育[6]等方面。Wu等[7]研究發(fā)現(xiàn),芽孢桿菌產(chǎn)生的?；呓z氨酸內(nèi)酯酶(Aii A)通過干擾其群體感應(yīng)系統(tǒng)的信號(hào)分子N-?；呓z氨酸內(nèi)酯(AHL),從而對(duì)花魔芋軟腐病菌產(chǎn)生抑制作用。雷珍珍等[8]從抗軟腐病花魔芋植株基因組中分離獲得一條NBSLRR 類抗病基因同源序列。上述結(jié)果為花魔芋對(duì)軟腐病菌侵染響應(yīng)分子機(jī)制的研究提供了理論基礎(chǔ),但這些研究僅局限于病原菌群體感應(yīng)的信號(hào)分子或花魔芋中部分抗病基因的分析,而植物對(duì)病原菌侵染的抗性反應(yīng)涉及多種途徑和基因的表達(dá)調(diào)控,有必要系統(tǒng)鑒定花魔芋中參與軟腐病脅迫響應(yīng)的多種基因,以揭示花魔芋抗病的分子機(jī)制及挖掘重要的抗病基因與調(diào)控功能。

植物染病后,基因轉(zhuǎn)錄和翻譯水平通常發(fā)生明顯改變,而基因表達(dá)水平的變化則早于病癥的出現(xiàn)。李玲[9]在研究2個(gè)稻瘟病菌生理小種侵染時(shí),發(fā)現(xiàn)在3個(gè)水稻品種上接種不同病菌后24 h和48 h分別有12 456個(gè)和4 851個(gè)基因發(fā)生差異表達(dá),說明在染病早期草本植物的基因表達(dá)已經(jīng)明顯改變。此外,軟腐病對(duì)花魔芋造成的危害嚴(yán)重且不可逆,增強(qiáng)植株自身染病初期的抗性十分必要。轉(zhuǎn)錄組測序技術(shù)(RNA-seq)作為一種高效快捷、準(zhǔn)確度高、成本低的研究手段,被廣泛應(yīng)用于模式植物擬南芥[10]和玉米[11]、非模式植物甘蔗[12]和花生[13]等多種植物抗病機(jī)理研究。本研究通過轉(zhuǎn)錄組測序技術(shù)探究花魔芋感染軟腐病菌初期基因表達(dá)的整體變化情況,通過對(duì)差異表達(dá)基因進(jìn)行GO 功能分類和KEGG 富集性分析,篩選花魔芋響應(yīng)細(xì)菌性軟腐病菌侵染的相關(guān)基因、抗病相關(guān)的信號(hào)傳導(dǎo)途徑和生物活性物代謝通路,以期為深入研究魔芋抗病分子機(jī)制提供理論參考,也為魔芋抗病分子育種提供候選基因。

1 材料與方法

1.1 材料

植物材料為1 a生的安康市本地花魔芋盆栽苗?；в筌浉【展z桿菌(Pectobacterium chrysanthemi)CZS-B6菌株由中國典型培養(yǎng)物保藏中心(China Center for Type Culture Collection,簡稱CCTCC,也稱武漢大學(xué)保藏中心)提供。

1.2 方法

1.2.1 病原菌接種與取樣 花魔芋軟腐病菌的培養(yǎng)和接種體的制備參考徐煒等[14]的方法。參照何斐等[15]的方法進(jìn)行盆栽花魔芋的種植與管理,于盆栽花魔芋生長至苗高約10 cm 左右,采用澆灌法[14]將稀釋至108CFU/m L 的接種體菌懸液澆灌到健康的盆栽花魔芋植株根際,每株澆灌200 m L,以澆灌無菌水的花魔芋作為對(duì)照。根據(jù)預(yù)試驗(yàn)效果,分別選取澆灌處理第5天的健康和軟腐病菌侵染初期的花魔芋球莖,于液氮中迅速冷凍后至-80℃保存。每6個(gè)生物學(xué)重復(fù)當(dāng)作1個(gè)處理,每處理重復(fù)2次。

1.2.2 RNA 提取、cDNA 文庫構(gòu)建和轉(zhuǎn)錄組測序 將凍存樣品在液氮中研磨,用RNApure Plant Kit(含DNaseⅠ)試劑盒(康為世紀(jì)生物科技有限公司,北京)提取總RNA,分別采用Nanodrop 2000 和Agilent 2100 生物分析儀(Agilent Technology,美國)檢測RNA 樣品的濃度和質(zhì)量。經(jīng)檢測合格后送至北京組學(xué)生物科技有限公司構(gòu)建cDNA 文庫,采用實(shí)時(shí)定量PCR 方法庫檢合格后,使用Illumina HiSeq 2500 測序平臺(tái)(Illumina Inc.,美國)進(jìn)行轉(zhuǎn)錄組測序(RNAseq)。

1.2.3 測序數(shù)據(jù)處理 花魔芋無參考基因組,采用Trinity軟件[16]對(duì)原始測序數(shù)據(jù)(Raw data)過濾后獲得高質(zhì)量數(shù)據(jù)(Clean data),對(duì)其進(jìn)行組裝拼接成轉(zhuǎn)錄本。去冗余后獲得unigenes,使用COG、GO[17]、KEGG[18]、KOG[19]、Pfam[20]、Swiss-Prot[21]和Nr數(shù)據(jù)庫進(jìn)行基因功能注釋。

差異表達(dá)基因篩選:用FPKM 值(每千個(gè)堿基的轉(zhuǎn)錄每百萬映射讀取的片段數(shù))[22]計(jì)算基因表達(dá)量。差異表達(dá)基因篩選標(biāo)準(zhǔn):花魔芋健康株與軟腐病株基因表達(dá)水平差異倍數(shù)FC(Fold change)≥2,且錯(cuò)誤發(fā)現(xiàn)率FDR(False discovery rate)<0.01[23]。

差異表達(dá)基因GO 注釋和KEGG 富集分析:利用Blast2GO(https://www.blast2go.org/)進(jìn)行差異表達(dá)基因的GO 注釋[24]。通過P值(P<0.05)的Benjamini-Hochberg矯正方法確定差異表達(dá)基因中顯著富集的GO 分類,使用Cytoscape(http://www.cytoscape.org/)注釋通路。利用KEGG 數(shù)據(jù)庫(http://www.genome.jp /kegg/),以P<0.05且Q≤0.05作為顯著富集標(biāo)準(zhǔn),進(jìn)行差異表達(dá)基因的通路富集分析[25]。

2 結(jié)果與分析

2.1 測序數(shù)據(jù)質(zhì)量分析

利用Illumina測序平臺(tái)對(duì)健康和軟腐病菌侵染初期的花魔芋球莖轉(zhuǎn)錄組進(jìn)行測序,經(jīng)過對(duì)測序原始數(shù)據(jù)質(zhì)量控制,在健康和發(fā)病魔芋中分別獲得19 236 352和21 116 175條clean reads(NCBI SRA 數(shù)據(jù)庫登錄號(hào):PRJNA738708,各樣本序列登錄號(hào) SAMN19759221 - SAMN 19759224),GC 含量分別為51.32%和51.36%,Q30值均高于93%。各樣品的clean reads與拼接得到的轉(zhuǎn)錄本參考序列的比對(duì)效率達(dá)到75%以上,數(shù)據(jù)顯示轉(zhuǎn)錄組測序質(zhì)量良好,可進(jìn)行下一步分析(表1)。

表1 花魔芋球莖轉(zhuǎn)錄組測序數(shù)據(jù)Table 1 Transcriptome data assembly for corm of A.konjac

2.2 Unigene功能注釋

將unigene序列與COG、GO、KEGG、KOG、Pfam、Swissprot和Nr數(shù)據(jù)庫比對(duì)和注釋,共獲得57 525條有注釋信息的unigenes,占全部unigenes的51.6%。其中COG 數(shù)據(jù)庫注釋信息最少,有18 357條;Nr數(shù)據(jù)庫注釋信息最多,高達(dá)53 611條,占總注釋unigenes的93.2%(表2)。

表2 Unigene注釋分析Table 2 Unigene annotation based on different databases

2.3 差異表達(dá)基因分析

轉(zhuǎn)錄組分析(圖1)表明,軟腐病菌侵染初期的花魔芋球莖中共出現(xiàn)差異表達(dá)基因3 222個(gè)。與健株相比,病株中有2 660 個(gè)基因表達(dá)上調(diào),562個(gè)基因表達(dá)下調(diào),上調(diào)基因比下調(diào)表達(dá)的基因多373.3%。差異基因表達(dá)量上調(diào)、下調(diào)最大的前10個(gè)基因中(表3),多數(shù)基因編碼特定功能的蛋白質(zhì)或酶,如上調(diào)表達(dá)基因TRINITY_DN674_c0_g3 編碼Mi AMP1 抗菌蛋白,下調(diào)表達(dá)基因TRINITY_DN3252_c0_g3 編碼熱激蛋白,TRINITY_DN11231_c0_g1編碼胰蛋白酶與蛋白酶抑制劑、胰凝乳蛋白酶抑制劑。研究表明這些蛋白質(zhì)、酶與抵御病菌侵染密切相關(guān)[26-27]。

表3 花魔芋在軟腐病菌脅迫下表達(dá)量上升和下降最大的前10個(gè)基因Table 3 Top 10 up-and down-regulated genes differentially expressed in A.konjac under P.chrysanthemi stress

圖1 軟腐病菌侵染初期花魔芋球莖中的差異表達(dá)基因分析Fig.1 Differential expression analysis of genes in corm of A.konjac with early response to P.chrysanthemi infection

2.4 差異表達(dá)基因功能分類和富集分析

對(duì)差異表達(dá)基因進(jìn)行GO 功能分類(圖2),共分為3大類53個(gè)功能組。在生物學(xué)過程、細(xì)胞組分和分子功能3大類中,2個(gè)差異表達(dá)最多的功能組分別是新陳代謝過程(706個(gè)差異表達(dá)基因)和細(xì)胞進(jìn)程(620個(gè))、細(xì)胞組分(627個(gè))和細(xì)胞(626個(gè))、催化活性(644個(gè))和結(jié)合(589個(gè))。

圖2 基于轉(zhuǎn)錄組測序數(shù)據(jù)的差異表達(dá)基因GO 功能分類Fig.2 GO functional classifications of differentially expressed genes(DEGs)based on transcriptome sequencing data

通過KEGG 對(duì)花魔芋各通路進(jìn)行注釋分析(圖3),發(fā)現(xiàn)1 148個(gè)差異表達(dá)基因在細(xì)胞過程、環(huán)境信息處理、基因信息處理、代謝和有機(jī)系統(tǒng)5大類的信號(hào)轉(zhuǎn)導(dǎo)及植物-病原菌互作等通路中富集。在29 條通路中,與代謝有關(guān)的通路最多(17條),占比58.6%。

圖3 差異表達(dá)基因的KEGG 通路分析Fig.3 KEGG pathway enrichment analysis of differentially expressed genes

在顯著富集的7條通路(表4)中,有3條與生物合成相關(guān),2條與代謝有關(guān),剩余的2條分別是糖酵解/糖異生和氧化磷酸化通路。富集程度最大的2 條通路分別是類黃酮生物合成[ko00941]和苯丙烷生物合成通路[ko00940]。在對(duì)花魔芋健株與軟腐病株差異表達(dá)基因分類和途徑分析的基礎(chǔ)上,推測軟腐病菌侵染花魔芋誘導(dǎo)的代謝途徑中的差異表達(dá)可能與苯丙烷和類黃酮生物合成等抗性有關(guān),對(duì)防御病菌侵染發(fā)揮重要作用[28]。

表4 差異表達(dá)基因中顯著富集的通路Table 4 Significantly enriched pathways of differentially expressed genes(DEGs)

2.5 趨勢表達(dá)基因通路分析

為解析花魔芋抵御軟腐病菌侵染的關(guān)鍵通路和基因,逐條分析差異表達(dá)基因映射在KEGG 數(shù)據(jù)庫中的參考通路,結(jié)果發(fā)現(xiàn)一些通路的差異基因存在特定的表達(dá)變化趨勢(表5)。其中,硒化合物代謝和維生素B6代謝通路中,所有差異基因相比健株均呈現(xiàn)表達(dá)量上升的趨勢。在氨基糖和核苷酸糖代謝,甘氨酸、絲氨酸和蘇氨酸代謝,類黃酮生物合成等14條通路中,60%以上的差異基因相比健株表達(dá)量上調(diào)。但甾類激素生物合成(57.14%)和類胡蘿卜素生物合成(66.67%)通路中的大部分差異表達(dá)基因相比健株呈現(xiàn)下調(diào)趨勢。

表5 趨勢表達(dá)基因KEGG 通路富集分析Table 5 Enrichment analysis of down-and up-regulated genes based on KEGG database

3 討論與結(jié)論

本研究對(duì)健康和軟腐病菌侵染初期的花魔芋球莖進(jìn)行轉(zhuǎn)錄組測序分析,獲得111 404條unigenes和3 222個(gè)差異表達(dá)基因,其中2 660個(gè)基因(82.6%)上調(diào)表達(dá),562個(gè)基因(17.4%)下調(diào)表達(dá)(病株/健株)。劉夢洋[29]研究發(fā)現(xiàn),大白菜抗軟腐病突變體受軟腐病菌侵染初期12 h 時(shí),412 個(gè)基因(81.9%)表達(dá)量上調(diào),91 個(gè)基因(18.1%)表達(dá)量下調(diào)。于秋佟[30]通過轉(zhuǎn)錄組測序發(fā)現(xiàn),軟腐病菌侵染24 h的擬南芥有579個(gè)基因(89.0%)表達(dá)量上調(diào),72個(gè)基因(11.1%)表達(dá)量下調(diào)。這些研究與本結(jié)果類似,即上調(diào)表達(dá)基因數(shù)量明顯多于下調(diào)表達(dá)基因的數(shù)量,說明更多的基因被誘導(dǎo)表達(dá)來抵御軟腐病菌的入侵和發(fā)展。本試驗(yàn)結(jié)果為后續(xù)鑒定花魔芋抵御軟腐病菌侵染初期發(fā)揮關(guān)鍵作用的基因縮小了篩選范圍,也為魔芋抗病育種奠定基礎(chǔ)。

在分析軟腐病菌脅迫下花魔芋中差異表達(dá)變化較大的基因時(shí)(表3),鑒定到部分響應(yīng)軟腐病菌脅迫的差異表達(dá)基因。例如,在上調(diào)表達(dá)的基因中TRINITY_DN674_c0_g3編碼MiAMP1抗菌蛋白,已有研究表明其與植物抗病有關(guān),Mi-AMP1抗菌蛋白重組畢赤酵母能誘導(dǎo)梨對(duì)擴(kuò)展青霉產(chǎn)生抗性[31]。TRINITY_DN 1022_c0_g3編碼類枯草桿菌蛋白酶家族,有研究顯示類枯草桿菌蛋白酶在病原物侵染后被特異性誘導(dǎo)并參與宿主植物細(xì)胞程序性死亡,同時(shí)在病原物識(shí)別和信號(hào)轉(zhuǎn)導(dǎo)級(jí)聯(lián)反應(yīng)中發(fā)揮重要作用[32]。值得關(guān)注的是,筆者發(fā)現(xiàn)也存在與抗病相關(guān)的下調(diào)表達(dá)基因,如TRINITY_DN3252_c0_g3編碼熱激蛋白,已有報(bào)道顯示眾多熱激蛋白參與調(diào)控植物抗逆和防衛(wèi)反應(yīng)[33]。TRINITY_DN11231_c0_g1 編碼胰蛋白酶和蛋白酶抑制劑,多種植物的蛋白酶抑制劑具有抗病蟲功能[34]。據(jù)此,推測花魔芋可能根據(jù)病菌脅迫的特點(diǎn),通過上調(diào)表達(dá)具有針對(duì)性的抗性基因,同時(shí)下調(diào)表達(dá)廣譜抗性基因來響應(yīng)軟腐病菌的脅迫。

本研究在分析KEGG 代謝通路時(shí),發(fā)現(xiàn)花魔芋中16條通路的差異表達(dá)基因在軟腐病菌脅迫下全部或大部分上調(diào)表達(dá)(表5)。經(jīng)逐條分析,可按照功能將其分成兩類:第1類與宿主植物抗逆有關(guān);第2類與細(xì)胞合成原料相關(guān)。在第1類通路中,硒化合物和維生素B6代謝通路中差異表達(dá)基因全部上調(diào)。已有研究證明,硒是植物生長發(fā)育不可缺少的微量元素,其通過維持植物體內(nèi)較高的抗氧化酶活性來誘導(dǎo)宿主產(chǎn)生抗病性[35];維生素B6作為多種代謝酶的重要輔助因子,參與植物氨基酸、乙烯合成等多種基本代謝活動(dòng),亦能通過調(diào)整細(xì)胞內(nèi)抗氧化能力從而調(diào)控番茄植株對(duì)灰霉病的抗性[36]。黃酮類化合物是植物自身抗病、抗逆的主要成分[37],在本研究中,涉及黃酮類化合物合成的類黃酮生物合成通路中有91.0%差異基因上調(diào)表達(dá)。N-聚糖生物合成通路合成的寡糖是一類植物抗逆誘導(dǎo)劑[38],其中木寡糖能誘導(dǎo)擬南芥對(duì)半活體病原細(xì)菌丁香假單胞菌番茄致病變種、活體病原煙草花葉病毒和死體病原真菌核盤菌產(chǎn)生抗性[39]。在本研究中N-聚糖生物合成通路83.3%差異基因上調(diào)表達(dá)。而鏈霉素生物合成通路合成的鏈霉素主要用于防治由革蘭氏陰性細(xì)菌引起的煙草野火病和青枯病、大白菜軟腐病、柑橘潰瘍病、水稻白葉枯病和條斑病等植物病害[40],在本研究中鏈霉素生物合成通路60%差異基因上調(diào)表達(dá),推測該通路合成的鏈霉素可能在抵御軟腐病菌侵入過程中發(fā)揮重要作用。第1大類參與宿主植物抗逆相關(guān)的信號(hào)轉(zhuǎn)導(dǎo)和物質(zhì)運(yùn)輸?shù)闹饕愤€包括氧化磷酸化、MAPK 信號(hào)通路、氨基糖和核苷酸糖代謝通路、植物-病原菌互作和植物激素信號(hào)轉(zhuǎn)導(dǎo)[41],這些通路75%以上的差異基因相比于健株表達(dá)量上調(diào)。

在第2類通路中,半胱氨酸和蛋氨酸代謝,甘氨酸、絲氨酸和蘇氨酸代謝,谷胱甘肽代謝,苯丙氨酸代謝,苯丙氨酸、酪氨酸和色氨酸的生物合成均與氨基酸的代謝或合成相關(guān),這些通路均有超過66.7%的差異基因受誘導(dǎo)上調(diào)表達(dá),這可能引起相應(yīng)的氨基酸合成加快,從而為合成抗病相關(guān)的蛋白和酶提供充足原料[42]。然而甾類激素和類胡蘿卜素生物合成通路的大部分差異基因表達(dá)下調(diào)(表5)。廣泛存在于植物體內(nèi)的類胡蘿卜素和甾類激素主要參與調(diào)控植物生長發(fā)育及生理功能[43],且類胡蘿卜素含量也與植物抗病性相關(guān)[44]。類胡蘿卜素和甾類激素生物合成通路中基因表達(dá)的下調(diào)可能引起細(xì)胞內(nèi)甾類激素和類胡蘿卜素供應(yīng)緊缺,從而影響細(xì)胞包括抗軟腐病相關(guān)基因表達(dá)等在內(nèi)的正常生理功能,這可能是染病花魔芋無法有效抵御軟腐病菌入侵的原因之一。

本研究通過RNA-seq測序得到3222個(gè)在花魔芋感染細(xì)菌性軟腐病初期出現(xiàn)差異表達(dá)的基因,篩選出與抗病相關(guān)且表達(dá)量差異大的基因。進(jìn)一步進(jìn)行GO 注釋和KEGG 通路分析發(fā)現(xiàn)2條通路的富集程度最大且合成物質(zhì)具有抗病性-苯丙烷生物合成和類黃酮生物合成通路。通過這2條通路合成的物質(zhì)可能在花魔芋抵御軟腐病菌侵入過程中發(fā)揮關(guān)鍵作用。同時(shí),研究發(fā)現(xiàn)16條通路中60%以上的差異基因上調(diào)表達(dá),而2條通路中大部分差異基因下調(diào)表達(dá)。其中,硒化合物代謝、維生素B6代謝、類黃酮生物合成、N-聚糖生物合成及鏈霉素生物合成這5條通路生成的相關(guān)物質(zhì)與花魔芋抗軟腐病菌脅迫密切相關(guān)。這些抗病相關(guān)通路和差異表達(dá)基因的分析結(jié)果為花魔芋抗病機(jī)制的研究奠定了基礎(chǔ),也為抗病相關(guān)基因的功能研究提供一定參考。