羽衣甘藍紅色系葉片形成相關microRNAs的鑒定及分析

張新鵬 徐宗大 劉林馨 夏甜甜 李曉艷*

(1.山東建筑大學 風景園林科學研究中心，濟南 250101； 2.山東農業大學 林學院，山東 泰安 271018； 3.山東建筑大學 建筑城規學院，濟南 250101)

羽衣甘藍(

Brassica

oleracea

var.

acephala

)為十字花科蕓薹屬食用甘藍的園藝變種，為二年生草本觀葉植物。羽衣甘藍葉片，特別是其心葉的葉色表型艷麗多彩，株形如牡丹盛開的花朵，業界對其有“葉牡丹”的美譽”。同時羽衣甘藍耐寒性極強，觀賞期較長，對我國北方地區景觀單調和少花的冬季而言，該植物對園林景觀綠化和美化的意義重大。羽衣甘藍心葉具體有純白、乳白、淡黃、淺粉、粉紅、玫紅和紫紅等類型，其中當屬紅色系葉片觀賞色彩最為突出，而且其花青苷含量相對最高，對人具有很好的保健功能。因此，探索和揭示紅色系葉片形成的調控機理，對于培育色彩突出且營養價值含量高的羽衣甘藍新品種具有重要意義。在許多園林植物中，葉色是其非常重要的觀賞性狀之一。對色彩鮮艷的葉片(如紅色和紫色等)而言，在生理層面，花青苷往往是呈色的關鍵花色素。先前許多研究表明花青苷合成結構基因的差異表達通常是影響植物花青苷的差異積累最直接的因素，而這些結構基因的表達模式和強度往往會受到其上游轉錄因子在不同時空上的調控。另外，在轉錄后水平上，microRNA(miRNA)也是植物中一類非常重要的基因表達調控因子，它是一類長度約21～24 nt的非編碼單鏈小分子RNA。近年來大量研究發現，許多miRNAs(如miR156、miR828和miR858等)可通過靶向某些與花青苷生物合成相關的轉錄因子(如

MYB

和

SPL

等)并抑制其表達，進而影響下游花青苷合成關鍵結構基因的表達水平，并最終調控植物組織中花青苷的合成和積累。

近年來，在許多植物中高通量測序技術已被廣泛應用于鑒定花青苷生物合成相關的重要miRNAs，如美人蕉、芍藥、蝴蝶蘭、牡丹和茶樹等。然而，對于羽衣甘藍而言，目前尚無有關葉色形成相關miRNA的報道。為了探索羽衣甘藍紅色系葉片中在轉錄后水平miRNAs對花青苷的調控機制，本研究以白色系葉片作為對比試材，以紅色系的羽衣甘藍葉片作為主要研究試材，基于高通量測序分析，以期篩選到紅色系葉片與花青苷合成相關的關鍵miRNAs，為今后在miRNA-mRNA層面系統探索羽衣甘藍葉色的分子調控機制提供一定的理論基礎。

1 材料與方法

1.1 試驗材料

紅色系和白色系的羽衣甘藍品種‘紫鴿’‘白鴿’材料于2020年種植于山東建筑大學風景園林實訓基地室外試驗田的同一區域，日常培養和外界環境條件均保持一致。在經冬季低溫葉片轉色后的觀賞期，采集這2個品種的心部(紫紅色和白色)葉片，并去掉略帶綠色的邊緣區域，分別用編號Bo-P表示紫紅色葉片，Bo-W表示白色葉片(圖1)。采集的上述組織材料，經液氮速凍后保存于-80 ℃超低溫冰箱。

圖1 研究用紫紅色(Bo-P)和白色(Bo-W)葉片Fig.1 Purple-red (Bo-P) and white (Bo-W) leaves used in this study

1.2 RNA提取和小RNA測序

將紫紅色和白色葉片組織樣品在液氮冷凍下進行研磨后，使用Trizol法提取總RNA，每種顏色葉片進行3次生物學重復。檢測合格的RNA樣品再進行進一步的文庫構建和小RNA測序。本研究構建了6個小RNA(sRNA)文庫，紫紅色葉片的文庫分別命名為Bo-P1、Bo-P2和Bo-P3，白色葉片的分別命名為Bo-W1、Bo-W2和Bo-W3。文庫構建和高通量測序均由生工生物工程(上海)股份有限公司完成。

1.3 miRNA的鑒定和靶基因預測

對小RNA測序得到的raw reads進行數據過濾，得到高質量的clean reads。對得到的clean reads與甘藍基因組(https:∥www.ncbi.nlm.nih.gov/genome/?term=Brassica+oleracea)、miRbase(https:∥www.mirbase.org/)和Rfam(http:∥rfam.xfam.org/)等已知的sRNA數據庫進行比對，對sRNA分類注釋。其中，通過與miRbase 數據庫(版本號22.1)中登錄的已知植物的miRNA成熟體序列進行比對，來鑒定羽衣甘藍葉片中的保守miRNAs(known miRNAs)。同時，基于miRNA的前體能夠形成發夾二級結構，用miRdeep2軟件再對未注釋的clean reads進行miRNA的預測，鑒定新miRNAs(novel miRNAs)。

以甘藍基因組序列作為參考序列，用TargetFinder軟件對羽衣甘藍葉片中鑒定出的保守和新的miRNA進行靶基因的預測。

1.4 miRNA的差異表達分析

用R語言包DESeq2進行miRNA的差異表達分析。本研究以白色葉片樣本作為對照，采用

P

value≤0.05和|log(Fold change)|≥1這2個閾值作為篩選紫紅色和白色葉片這2組樣本中差異表達miRNA的條件。同時，對差異表達miRNA預測的靶基因進行GO和KEGG pathway富集分析。

1.5 候選miRNA的熒光定量(qRT-PCR)分析

利用TaKaRa公司的Mir-XmiRNA First-Strand Synthesis Kit反轉錄試劑盒對miRNA進行反轉錄，得到的cDNA用于候選miRNAs的進一步定量。再利用TaKaRa 公司的TB GreenPremix Ex

Taq

II試劑盒，通過qRT-PCR法對候選miRNAs(miR828和miR858-3)在紫紅色和白色葉片中的差異表達情況進行驗證分析。其中，內參基因為

U

6，每個樣品設置3次重復。以白色葉片的基因表達量為對照，通過2相對定量法計算miRNA的相對表達水平。引物合成由生工生物工程(上海)股份有限公司完成，引物的具體序列信息見表1。

表1 候選miRNA熒光定量PCR所用的引物序列
Table 1 Primer sequences of the candidate miRNAs used for qRT-PCR

miRNA名稱miRNA name上游引物(5'-3')Forward primer (5'-3')下游引物(5'-3')Reverse primer (5'-3')miR828ACACTCCAGCTGGGTCTTGCTTAAAT-GAGTGGTGTCGTGGAGTCGmiR858-3ACACTCCAGCTGGGTTCGTTGTCTGTTCTGGTGTCGTGGAGTCGU6TTGGAACGATACAGAGAAGATTAGCATTGGACCATTTCTCGATTTGTG

2 結果與分析

2.1 測序數據質量評估

本研究利用高通量測序法對紫紅色和白色葉片的6個小RNA文庫(Bo-P1、Bo-P2、Bo-P3、Bo-W1、Bo-W2和Bo-W3)進行了測序，從紫紅色和白色葉片文庫獲得的原始數據量分別為61 388 537和52 363 761 條raw reads。通過去除包含poly-A、帶接頭和低質量等reads，從紫紅色葉片文庫中最終獲得了35 183 047條Q20為97.42%的clean reads，從白色葉片文庫中獲得了33 023 141條Q20為97.49% 的clean reads(表2)。

對于sRNA的長度分布，在紫紅色和白色葉片這2個文庫中，大多數sRNA的長度分布區間均為21～24 nt。在Bo-P和Bo-W文庫中，長度為24 nt的sRNA所占比例均為最高，分別為16.4%和17.6%，其他主要長度所占比例的排序依次為21、22和23 nt(圖2)。

2.2 羽衣甘藍葉片中miRNA的鑒定及靶基因的預測

小RNA(sRNA)的組成種類比較多，通過將紫紅色和白色葉片這2個文庫得到的clean reads分別與甘藍基因組和不同的sRNA數據庫進行比對，可以將miRNA、rRNA、tRNA、snRNA和snoRNA等分類注釋出來。結果顯示，在Bo-P(圖3(a))和Bo-W(圖3(b))這2個文庫中分別有2.48%和1.99% 的clean reads比對到miRNA上。經進一步數據分析，在這6個sRNA文庫中共鑒定出1 027個保守miRNAs(known miRNAs)。在Bo-P和Bo-W文庫中分別鑒定出799和835個known miRNAs，其中這2個文庫中的大多數known miRNAs(607個)都是相同的，另外Bo-P和Bo-W文庫分別各有192和228個不同的known miRNAs。初步推測紫紅色或白色葉片中這些不同的known miRNAs可能與羽衣甘藍葉色的形成有一定的關聯。

表2 6個小RNA文庫測序數據
Table 2 Sequencing data of small RNAs (sRNAs) in the six libraries

樣品名稱Sample name過濾后的readsClean readsQ20值/%Q20 valueGC含量/%GC content平均長度/bpAverage lengthBo-P111 610 46097.5544.2222.30Bo-P211 980 71797.1943.9922.63Bo-P311 591 87097.5144.0822.64Bo-W111 443 05897.4343.5723.00Bo-W210 014 61697.4843.6122.97Bo-W311 565 46797.5743.7822.95Bo-P總計 Bo-P summary35 183 04797.4244.1022.52Bo-W總計 Bo-W summary33 023 14197.4943.6522.97

圖2 sRNA在紫紅色葉片(a)和白色葉片(b)中的長度分布Fig.2 Length distribution of sRNA in Bo-P leaf (a) and Bo-W leaf (b)

對于那些與已知sRNA數據庫未比對上的clean reads，對其進行了新miRNA(novel miRNA)的預測。結果發現，在這6個sRNA文庫中共鑒定出100個novel miRNAs。其中，在Bo-P和Bo-W文庫中分別鑒定出100和99個novel miRNAs。另外，對于在紫紅色和白色葉片的6個sRNA文庫中鑒定出的known miRNAs和novel miRNAs，通過分析發現它們分別有892和86個miRNAs存在靶基因，其相應的靶基因個數為8 428和1 923個。

括號內的數字和百分號數字分別代表該項的reads數和其所占clean reads的比例。 The numbers and percentages in brackets represent the number of reads of this item and its proportion in clean reads, respectively.圖3 sRNA在紫紅色葉片(a)和白色葉片(b)中的分類注釋Fig.3 Classification of sRNA annotations in Bo-P leaf (a) and Bo-W leaf (b)

2.3 不同顏色葉片差異表達的miRNA分析

為了解影響羽衣甘藍紅色系葉片形成的關鍵miRNA，對紫紅色和白色葉片中鑒定出的miRNAs進行了差異表達分析。基于閾值

P

value≤0.05和|log(Fold change)|≥1，在這2種顏色葉片中共鑒定出110個(屬于33個家族)差異表達的miRNAs(DEMs)，其中32個在紫紅色葉片中表達上調，78個在紫紅色葉片中表達下調(圖4)。

2.4 差異表達miRNA靶基因的功能富集分析

為了更好地明確紫紅色和白色葉片中這些靶基因行使的主要生物學功能，對其進行了GO和KEGG功能富集分析。GO功能顯著富集(矯正后的

P

value≤0.05)分析結果(圖5)顯示，對于生物學過程類別而言，富集主要集中于轉錄(Transcription)和轉錄調控(Regulation of transcription)這2個子類別。對于細胞成分類別而言，富集主要集中于CCAAT結合位點復合體(CCAAT-binding factor complex)和細胞核(Nucleus)這2個子類別。在分子功能這1大類中，富集主要集中于轉錄因子活性及序列特異性DNA結合(Transcription factor activity, sequence-specific DNA binding)和DNA結合(DNA binding)這2個子類別。

對差異表達miRNA的靶基因也進行了KEGG pathway富集分析，結果顯示這些差異表達miRNA的靶基因主要在植物MAPK信號通路(MAPK signaling pathway-plant)、植物激素信號轉導(Plant hormone signal transduction)和類黃酮合成(Flavonoid biosynthesis)等通路富集(圖6)。

圖4 差異表達miRNAs的火山圖Fig.4 Volcanic map of differentially expressed miRNAs

圖5 差異表達miRNAs的靶基因的GO功能富集Fig.5 GO functional enrichment of target genes of differentially expressed miRNAs

圖6 差異表達miRNA的靶基因的KEGG pathway富集Fig.6 KEGG pathway enrichment of target genes of differentially expressed miRNAs (DEMs)

2.5 羽衣甘藍紅色系葉片形成相關的差異表達miRNA

在這110個差異表達miRNAs中，基于靶基因的功能注釋和前人花青苷合成相關miRNA的研究，篩選出了2個可能與羽衣甘藍紅色系葉片形成相關的重要miRNAs(miR828和miR858-3)。其中，miR828的靶基因為XM_013735358.1 [R2R3-MYB: Transcription factor MYB114-like]。miR858-3的靶基因包括XM_013758669.1 [Transcription factor MYB12-like]和XM_013763976.1 [R2R3-MYB: Transcription factor MYB86-like]。小RNA測序結果顯示，與白色葉片相比，miR828和miR858-3在紫紅色葉片中的表達均明顯上調(表3)。

表3 2個可能參與羽衣甘藍紅色系葉片形成的候選差異表達miRNAs的表達水平
Table 3 Expression levels of two candidate differentially expressed miRNAs possibly involved in formation of red-color leaf in var.

miRNA名稱miRNA namelog2(Bo-P/Bo-W)靶基因Target gene注釋AnnotationmiR8286.06XM_013735358.1Transcription factor MYB114-likemiR858-31.85XM_013758669.1Transcription factor MYB12-likeXM_013763976.1Transcription factor MYB86-like

本研究進一步采用qRT-PCR的方法對這2個候選重要miRNAs的表達模式進行了驗證分析。結果顯示，它們均與小RNA測序的結果呈現出一致的趨勢，其中高通量測序和qRT-PCR這2次結果顯示miR828在紫紅色葉片中的表達量分別為白色葉片中的66.7倍和44.4倍(圖7)。

不同字母代表在0.05水平上存在顯著性差異。 Different lowercase letters indicate significant differences at 0.05 level. 圖7 2個候選差異表達miRNAs的定量分析Fig.7 qRT-PCR analysis of 2 candidate differentially expressed miRNAs

3 討論與結論

羽衣甘藍紅色系內葉的色彩表型極為突出，該色系的葉片富含花青素且可食用，這使其在一定程度上兼具觀賞、營養和商業價值。因此，探索其紅色系葉片形成的機理具有重要意義。對于羽衣甘藍紅色系在生理層面的形成機理，前人研究發現花青苷，特別是矢車菊素類糖苷，是紅色系葉片呈色的關鍵花色素。羽衣甘藍紅色系葉片在分子層面的形成機理，目前多數研究主要集中在mRNA(結構基因和轉錄因子)層面。對花青苷合成相關結構基因而言，研究發現

BoDFR

應該是羽衣甘藍紅色系葉片形成的關鍵結構基因。任杰通過對

DFR

啟動子區域序列分析發現，該基因上存在花青苷合成相關轉錄因子MYB的結合位點。對花青苷合成相關轉錄因子而言，研究發現MYB類轉錄因子(如BoMYB2、BoMYB1、BoPAP1和BoMYBL2.1)可能對羽衣甘藍紅色系葉片中的花青苷積累起著重要的調控作用。但在轉錄后水平上，特別是miRNA層面目前尚未見報道，本研究首次報道了羽衣甘藍紅色系葉片中與葉色調控相關的miRNAs。

近年來在植物花青苷生物合成領域(包括葉色、花色和果色等)，miRNAs作為在轉錄后層面重要的基因調控者，引起了眾多研究者的關注。本研究分別構建了‘紫鴿’和‘白鴿’羽衣甘藍品種紅色和白色系葉片的小RNA文庫，并進行了高通量測序和相關分析。在這2種顏色葉片的文庫中，均顯示21～24 nt長度的sRNAs占多數，而24 nt的sRNA比例最高。這一結果與枸杞、甘薯和茶樹等的研究結果一致，而與大麗花、芍藥和牡丹等的報道不同，這表明在不同植物中sRNA的生物合成途徑可能存在一些差異。

本研究在紫紅色和白色葉片文庫中共鑒定了110個顯著差異表達的miRNAs，同時對其進行了靶基因的預測及功能分類。基于前人花青苷合成相關miRNAs的研究和對相應靶基因的進一步分析，本研究在Bo-P和Bo-W 文庫中共篩選到2個與花青苷合成相關的DEMs(miR828和miR858-3)，并發現它們的靶基因中均有R2R3-MYB轉錄因子(XM_013735358.1和XM_013763976.1)。對MYB類轉錄因子家族而言，許多研究已表明其參與調控花青苷生物合成的亞族多為R2R3-MYB轉錄因子。

對于miR828和miR858這2個家族的一些miRNAs，前人許多研究也已表明它們均可以通過調控其靶基因R2R3-MYB轉錄因子的表達，進而最終影響花青苷的生物合成和積累。本研究的qRT-PCR分析結果顯示，miR828和miR858-3的表達趨勢與高通量測序結果基本一致。同時這2個結果均顯示，miR828在紫紅色和白色葉片中的差異高達44倍以上。

因此，本研究認為miR828和miR858-3(特別是miR828)應該是轉錄后水平上調控羽衣甘藍紅色系葉片形成的關鍵miRNAs。該結果將有助于人們進一步理解miRNAs在羽衣甘藍葉色形成中潛在的功能角色，以期為今后深入揭示羽衣甘藍紅色系葉片形成的分子機理提供一定程度的理論參考。