趙振利， 張靖曼， 鄭秋莉， 楊海波， 范國強(qiáng)

(河南農(nóng)業(yè)大學(xué)泡桐研究所，河南 鄭州 450000)

泡桐是中國重要的速生用材、農(nóng)田防護(hù)和園林綠化樹種，在緩解木材短缺、保障糧食安全、改善生態(tài)環(huán)境等方面發(fā)揮重要作用[1-2]。然而，植原體侵染后可引起泡桐發(fā)生叢枝病，泡桐叢枝病(PaWB)是一種危害極大的泡桐病害，患病后導(dǎo)致泡桐幼樹死亡、大樹生長緩慢，嚴(yán)重影響泡桐產(chǎn)業(yè)的發(fā)展。因此，研究PaWB的發(fā)生機(jī)制是解決泡桐這一病害的關(guān)鍵。表觀遺傳作為一種重要的修飾對基因表達(dá)進(jìn)行調(diào)控進(jìn)而影響植物的生長發(fā)育，而DNA甲基化是重要的表觀遺傳調(diào)控機(jī)制之一，其在調(diào)節(jié)植物生長發(fā)育、應(yīng)對生物或非生物脅迫、維持基因組穩(wěn)定等方面發(fā)揮著重要作用[3-5]。目前，關(guān)于DNA甲基化的研究大多在草本植物中進(jìn)行的，如在鐵皮石斛[6]、懷地黃[7]、四倍體長春花[8]等草本植物中探究DAN甲基化在植物體生長發(fā)育過程中的調(diào)控作用及對逆境響應(yīng)等方面的研究。但對木本植物DNA甲基化研究甚少，阻滯了木本植物表觀遺傳學(xué)研究的進(jìn)展。本研究利用全甲基化測序技術(shù)(WGBS)對白花泡桐健康苗(PF)和叢枝病苗(PFI)進(jìn)行高通量測序，系統(tǒng)分析PaWB發(fā)生過程中甲基化的變化情況，并通過DNA甲基化和轉(zhuǎn)錄組關(guān)聯(lián)分析，篩選出與PaWB發(fā)生相關(guān)的關(guān)鍵調(diào)控基因，以期為進(jìn)一步闡明PaWB發(fā)生分子機(jī)理提供科學(xué)依據(jù)。

1 材料與方法

1.1 試驗材料

試驗材料來自河南農(nóng)業(yè)大學(xué)泡桐研究所培養(yǎng)的PF和PFI。將植物材料放在1/2 MS培養(yǎng)基(蔗糖20 g·L-1，瓊脂粉4.0 g·L-1)中培養(yǎng)，每種材料的外植體葉片培養(yǎng)20瓶，每瓶接種3個外植體，待幼苗長至高度為2 cm時轉(zhuǎn)入1/2 MS生根培養(yǎng)基(蔗糖20 g·L-1，瓊脂粉3.0 g·L-1)中，置于溫度(25±2) ℃、光照強(qiáng)度1 300 lx、光照時間16 h·d-1的條件下培養(yǎng)30 d后，剪取組培苗頂端2～4片嫩葉，用液氮冷凍后置于-80 ℃超低溫冰箱內(nèi)保存?zhèn)溆谩?/p>

1.2 試驗方法

1.2.1 DNA提取、建庫及測序 將上述生長發(fā)育狀態(tài)一致的葉片混合，10株為一個樣本，進(jìn)行葉片DNA提取，DNA的提取按照張延召等[9]的方法進(jìn)行。DNA樣品經(jīng)深圳華大基因公司質(zhì)檢合格后進(jìn)行亞硫酸氫鹽(BS)文庫構(gòu)建，文庫質(zhì)量檢測后通過Illumina HiSeq 2000(Illumina，San Diego，CA，USA)測序儀進(jìn)行雙末端測序，具體方法參照文獻(xiàn)[10]。

1.2.2 測序數(shù)據(jù)的過濾、比對分析 測序結(jié)果經(jīng)過識別堿基后轉(zhuǎn)化為原始測序序列(Raw Reads)，然后除去Raw Reads中含接頭序列、低質(zhì)量、含N比例>10%的Reads后得到高質(zhì)量序列(Clean Reads)。通過軟件BSMAP (version bsmap-2.73) 將獲得的Clean Reads比對到泡桐參考基因組上(http://paulownia.genomics.cn/page/species/index.jsp)，使用唯一比對數(shù)據(jù)得到全基因組C堿基甲基化信息并進(jìn)行后續(xù)分析。具體方法參照文獻(xiàn)[11]。

1.2.3 滑動窗口法鑒定DMR 采用滑動窗口法確定DMR[12]。具體方法如下：(1)DMR的識別標(biāo)準(zhǔn)選擇2組間測序深度≥4×(×表示測序read覆蓋基因組的深度)的CpG位點作為候選位點，并對每個位點的差異顯著甲基化進(jìn)行Fisher精確檢驗。(2)將第一個甲基化差異顯著的CpG(自動校正P值<0.05)作為DMR的初始位點，并將后續(xù)候選位點合并為候選DMR，其中候選位點滿足以下條件：① 相鄰2個候選位點之間的距離≦300 bp；② 候選DMR中的所有候選CpG位點保持相同的甲基化趨勢；③ 候選的DMR至少有5個CpG位點且其中至少3個CpG位點的Fisher精確檢驗P值<0.01。(3)對于上述候選DMR，依據(jù)平均甲基化水平再次進(jìn)行Fisher精確檢驗。(4)候選DMR的錯誤發(fā)現(xiàn)率根據(jù)P值進(jìn)行調(diào)整，最后鑒定出的Fisher檢驗錯誤發(fā)現(xiàn)率P值<0.05且2個樣本之間平均甲基化水平差異小于20%的區(qū)域即是DMR。

1.2.4 DNA甲基化與轉(zhuǎn)錄組關(guān)聯(lián)分析鑒定 DMG啟動子區(qū)域含有DMR的基因定義為DMG，將獲得的DMR與DONG等[13]前期研究的轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行關(guān)聯(lián)，找到甲基化和轉(zhuǎn)錄組共同差異表達(dá)的基因，并從中篩選出與泡桐叢枝病發(fā)生相關(guān)的基因，然后對這些共同差異表達(dá)基因進(jìn)行GO和KEGG功能分析。

2 結(jié)果與分析

2.1 白花泡桐全基因組DNA甲基化測序數(shù)據(jù)分析

使用Illumina HiSeq 2000測序平臺對樣品PF、PFI進(jìn)行測序，共得到71.63 G的Clean reads(表1)。使用軟件BSMAP(version bsmap-2.73)將PF、PFI樣品得到的Clean reads比對到白花泡桐基因組序列上，比對率分別為90.40%、80.72%，全基因組平均測序深度分別是53.15×和45.45×(表2)。質(zhì)控檢測結(jié)果顯示，測序結(jié)果可進(jìn)行后續(xù)分析。

表1 PF、PFI測序數(shù)據(jù)統(tǒng)計結(jié)果

表2 Clean reads比對到參考基因組結(jié)果統(tǒng)計Table 2 Statistics result of clean reads map to reference genome results

2.2 白花泡桐全基因組DNA甲基化C特征分析

C位點的DNA甲基化主要有3種類型(mCG、mCHG和mCHH)，這3種類型的數(shù)量及構(gòu)成比例反映了特定物種的全基因組甲基化的特征[14]。本研究發(fā)現(xiàn)，PF和PFI中分別有37.65%和36.96%的C位點發(fā)生了甲基化，且DNA甲基化變化主要有mCG、mCHH和mCHG 3種類型，其中PF的C位點和CHH的甲基化比率均高于PFI，而CG和CHG的甲基化比率則均低于PFI(表3)。對DNA甲基化類型的研究結(jié)果表明，泡桐發(fā)病前后其mCG甲基化位點在染色體上出現(xiàn)的頻率最高(均達(dá)到71%以上)，mCHG位點和mCHH位點出現(xiàn)頻率依次降低(圖1)。研究結(jié)果表明，白花泡桐被植原體侵染后，C位點發(fā)生甲基化呈現(xiàn)下降趨勢，不同類型的mC模式均發(fā)生了變化，這種變化趨勢體現(xiàn)了DNA甲基化水平變化保持相對穩(wěn)定。

表3 PF、PFI中C位點甲基化比例統(tǒng)計Table 3 Statistical analysis of methylation ratio at C site in PF,PFI

2.3 白花泡桐全基因組DNA甲基化水平分析

全基因組C位點DNA甲基化水平分析結(jié)果表明，白花泡桐在植原體侵染后其DNA甲基化水平由14.53%下降到12.47%。此外，白花泡桐發(fā)病后mCG、mCHG和mCHH水平也均呈現(xiàn)下降趨勢，分別由50.84%降至45.91%、由22.97%降至19.98%和由8.76%降至7.13%(表4)。研究結(jié)果表明，DNA甲基化主要發(fā)生在CG序列，CHH序列在mC中占比最低，從而說明PaWB的發(fā)生與mCG的類型變化關(guān)系更密切。

2.4 白花泡桐DNA甲基化位點在不同基因區(qū)域的分布

對白花泡桐發(fā)病前后基因不同區(qū)域上發(fā)生甲基化水平的分析發(fā)現(xiàn)，在各基因區(qū)域上，CG甲基化水平最高，CHG甲基化水平次之，CHH甲基化水平最低(圖2)。由圖2可以看出，在白花泡桐感染植原體后許多基因區(qū)域的CG甲基化水平都有所下降，說明在侵染過程中存在多個位點的去甲基化過程。研究結(jié)果表明，白花泡桐發(fā)病過程中內(nèi)含子區(qū)域發(fā)生甲基化水平較高，這可能是因為內(nèi)含子參與基因表達(dá)調(diào)控及防止了內(nèi)含子序列的異常轉(zhuǎn)錄所導(dǎo)致的。

2.5 白花泡桐DNA甲基化組的DMR分析

對樣品PF和PFI進(jìn)行甲基化區(qū)域的差異研究分析，共檢測到109 334個DMR(表5)。從表5可以看出，白花泡桐在感染植原體過程中hyper-DMR和hypo-DMR分別有37 408個和71 926個，且CG、CHG、CHH序列都發(fā)生低甲基化的比率更大，其中CHH序列發(fā)生低甲基化比率最大為83.4%，CHG序列發(fā)生低甲基化比率最小為52.7%。PF與PFI中DMR的比較結(jié)果表明，植原體感染后泡桐體內(nèi)甲基化水平普遍降低。

A：白花泡桐健康苗不同類型mC堿基的分布比例；B：叢枝病苗不同類型mC堿基的分布比例A：Distribution ratio of mC bases in different sequence types of PF；B：Distribution ratio of mC bases in different sequence types of PFI圖1 白花泡桐不同類型mC堿基的分布比例Fig.1 Distribution ratio of mC bases in different sequence types of Paulownia fortunei

表4 PF、PFI中C、CG、CHG 和 CHH 甲基化水平Table 4 Methylation levels of C, CG, CHG and CHH in PF,PFI %

DMR在基因不同區(qū)域的分布研究結(jié)果見圖3。從圖3可以看出，位于基因間區(qū)的CG-DMR最多，CHG-DMR和CHH-DMR傾向于發(fā)生在基因間區(qū)和轉(zhuǎn)錄終止區(qū)(TE)上，說明基因間區(qū)(intergenic)甲基化水平較高，并且TE可能參與CHG和CHH的甲基化調(diào)控過程。此外，在上游5k(upstream5k)區(qū)域，CHH-DMR的發(fā)生水平較高，這可能是CHH-DMR參與了調(diào)控基因表達(dá)造成的結(jié)果。

A：白花泡桐健康苗中不同基因區(qū)域上各類型C的甲基化水平；B：叢枝病苗中不同基因區(qū)域上各類型 C 的甲基化水平。A:Methylation levels of various types of C in different gene regions of PF; B:Methylation levels of different types of C in different gene regions of PFI.圖2 不同基因區(qū)域上各類型C的甲基化水平Fig.2 Methylation level of various types of C in different gene regions

2.6 白花泡桐DNA甲基化與轉(zhuǎn)錄組關(guān)聯(lián)分析

通過對DNA甲基化和轉(zhuǎn)錄組的差異表達(dá)基因進(jìn)行關(guān)聯(lián)分析，共得到1 693個共同差異表達(dá)的基因(表6)。CAO等[15]研究表明，在植原體侵染泡桐后多個與脅迫應(yīng)答相關(guān)的基因的表達(dá)水平和甲基化水平呈負(fù)相關(guān)。在此研究結(jié)果基礎(chǔ)上，找到了6個與泡桐叢枝病發(fā)生相關(guān)的DMG，功能分別涉及到調(diào)節(jié)植物防御反應(yīng)、光合作用、防御反應(yīng)和ABA合成等方面(表7)。其中，PAU024721.1在啟動子區(qū)域其CG和CHH甲基化水平下降，基因表達(dá)量上調(diào)；PAU026702.1存在多個位點的CHH、CHG甲基化水平下降，基因表達(dá)量上調(diào)；PAU000493.1和PAU017490.2存在多個位點的CHH甲基化水平降低，基因表達(dá)水平上調(diào)；PAU009265.1存在CG甲基化水平升高，基因表達(dá)量下調(diào)；PAU009776.1存在CHG甲基化水平升高，基因表達(dá)量下調(diào)。

表5 DMR在PFvsPFI間的變化情況Table 5 Details of DMR during PFvsPFI

圖3 不同序列類型DMR在基因組區(qū)域中所占比例Fig.3 Proportion of DMR of different sequence types in the genome region

2.7 泡桐叢枝病差異表達(dá)基因功能分析

將得到的叢枝病發(fā)生相關(guān)的基因進(jìn)行GO注釋，共有45條GO條目被顯著富集(P<0.05)(圖4)，分為分子功能(molecular function)、細(xì)胞組分(celluar component)、生物過程(biological process)3類?；虻腒EGG通路聚類分析(圖5)，結(jié)果顯示，共有122條通路富集，其中主要富集在代謝通路、次生代謝產(chǎn)物生物合成、植物病原體互作、植物激素信號轉(zhuǎn)導(dǎo)等生物學(xué)過程。上述研究結(jié)果表明，差異甲基化基因可能和植物激素、植物病原互作具有密切關(guān)聯(lián)。

表6 在PF、PTF中鑒定到的差異表達(dá)的DMG數(shù)目統(tǒng)計Table 6 Statistics of the number of differentially expressed DMGs identified in PF, TF

3 討論

DNA甲基化作為一種重要的表觀遺傳修飾方式，它能導(dǎo)致基因表達(dá)發(fā)生可遺傳的變化，在植物的生長發(fā)育過程中起重要調(diào)控作用。DNA甲基化水平在不同植物種類中變化明顯，在同物種不同品種之間也有差異。本研究發(fā)現(xiàn)，白花泡桐中3種不同類型CG、CHG和CHH的DNA甲基化水平分別為50.84%、22.97%和8.76%，與木本植物蘋果[16]、白樺[17]和毛果楊[18]的甲基化水平相近，而跟草本植物擬南芥[19]差別較大，其中，白花泡桐和蘋果的CG、CHH的甲基化水平接近，且都比白樺和毛果楊的高，白花泡桐的CHG的甲基化水平和白樺接近。本研究對PF和PFI 的DNA甲基化進(jìn)行比較發(fā)現(xiàn)PFI整體甲基化水平呈下降趨勢，獲得的109 334個DMR中CG、CHG、CHH都發(fā)生低甲基化現(xiàn)象，且主要發(fā)生在CG序列?；蚪M功能元件區(qū)域分析研究也發(fā)現(xiàn)，CDS、內(nèi)含子等元件上的DNA甲基化水平在泡桐患病前后變化不大，但在泡桐感染植原體后，許多基因結(jié)構(gòu)的CG甲基化水平都有所下降，說明在侵染過程中存在多個位點的去甲基化過程，從而影響到基因的表達(dá)情況。

表7 PaWB發(fā)生相關(guān)差異甲基化基因Table 7 Differentially expressed methylated genes related to the occurrence of PaWB

圖4 PaWB發(fā)生相關(guān)基因GO分類圖Fig.4 Histogram of GO classifications of PaWB related genes

圖5 PaWB發(fā)生相關(guān)基因KEGG通路聚類圖

病原體入侵植物時，植物激素會通過相互連接的信號通路協(xié)同工作并誘導(dǎo)復(fù)雜的防御反應(yīng)[20-21]。脫落酸(ABA)是植物逆境反應(yīng)中的一種重要激素，它通過影響胼胝質(zhì)的沉積，阻礙病原菌在植物細(xì)胞中的運動從而對病原菌侵染有一定的抗性[22-23]。9-順式環(huán)氧類胡蘿卜素雙加氧酶(NCED)是合成ABA的關(guān)鍵酶，它催化9-順式葉黃素分解為ABA的前體黃素，然后通過脫落酸醛的兩個反應(yīng)將黃素轉(zhuǎn)化為ABA[24]。HAO 等[25]發(fā)現(xiàn)，在干旱脅迫下，過度表達(dá)AtNCED3的植物具有更高的耐旱性，而抑制和破壞AtNCED3表達(dá)的植物則產(chǎn)生干旱表型。本研究發(fā)現(xiàn)，編碼NCED的基因發(fā)生去甲基化且基因表達(dá)水平升高，從而導(dǎo)致ABA含量的升高，說明在植原體侵染過程中泡桐通過調(diào)控體內(nèi)ABA信號通路從而調(diào)控植原體的侵染。

ROS是參與防御反應(yīng)的調(diào)節(jié)因子之一，但高濃度的ROS會對植物造成一定傷害，因此去除ROS的抗氧化劑細(xì)胞色素 P450對保護(hù)植物起著重要作用[26-28]。HUBBARD 等[29]研究證實了細(xì)胞色素P450在植物防御狀態(tài)時含量升高的情況。本研究也發(fā)現(xiàn)，白花泡桐患叢枝病過程中，編碼細(xì)胞色素P450基因甲基化水平降低致使基因表達(dá)量升高，這可能是由于在植原體侵染泡桐后，植物產(chǎn)生ROS參加防御反應(yīng)，又激活編碼細(xì)胞色素P450基因的表達(dá)。WRKY轉(zhuǎn)錄因子33(WRKY33)屬于第I組WRKY蛋白家族，對植物體抵抗病原體的侵染發(fā)揮重要的作用[30]。ZHENG等[31]在擬南芥中發(fā)現(xiàn)編碼AtWRKY33 基因的過表達(dá)增強(qiáng)了對壞死性真菌病原體灰葡萄孢和鏈格孢菌的抗性。在本研究中，編碼WRKY33基因由于發(fā)生甲基化水平降低導(dǎo)致其基因表達(dá)量上調(diào)，這種變化可能是為了抵抗植原體的侵染，DNA甲基化通過負(fù)調(diào)控相關(guān)基因從而提高WRKY33基因的表達(dá)以增強(qiáng)泡桐的抗性。

本研究利用WGBS技術(shù)對PaWB發(fā)生過程中全基因組DNA甲基化的變化進(jìn)行了研究，分析了DNA甲基化的水平和模式變化情況，并通過DNA甲基化和轉(zhuǎn)錄組的差異表達(dá)基因關(guān)聯(lián)分析，找到了與泡桐叢枝病發(fā)生相關(guān)的DMG，這些基因涉及植物防御反應(yīng)、光合作用、植物病原互作、植物激素信號轉(zhuǎn)導(dǎo)等生物學(xué)過程。