孫燕 王金剛 臧丹丹 趙恒田 劉淑華

（1.東北農(nóng)業(yè)大學(xué)園藝園林學(xué)院，哈爾濱 150030；2.中國科學(xué)院東北地理與農(nóng)業(yè)生態(tài)研究所黑土區(qū)農(nóng)業(yè)生態(tài)重點實驗室，哈爾濱 150081）

藍果忍冬（Lonicera caeruleaL.）為忍冬科忍冬屬植物，原產(chǎn)于北半球溫帶地區(qū)，現(xiàn)廣泛分布于我國各地，是我國重要的野生果樹資源［1-3］。藍果忍冬果實出汁率高，有香氣，營養(yǎng)價值豐富，具有調(diào)節(jié)血壓、減緩衰老等作用［4-5］，因其富含花青素、氨基酸、維生素、類黃酮等多種活性物質(zhì)，被譽為“飲料之王”［6］?！锌扑{1 號’是我國第一個具有自主知識產(chǎn)權(quán)和國家備案的藍果忍冬新品種，植株健壯，豐產(chǎn)性強，果實為長橢圓形，口感酸甜無苦味，是優(yōu)良的鮮果品種，適宜東北寒區(qū)廣泛種植［7］。

花青素是植物體內(nèi)一類重要的次生代謝物質(zhì)，屬于類黃酮類化合物，在賦予植物色彩時起主要作用［8-9］。在果實生長發(fā)育過程中，果皮中色素會發(fā)生一系列變化，隨著果實的成熟，葉綠素在水解酶的作用下逐漸分解消失，類胡蘿卜素略有增加；花青素在細胞質(zhì)中合成，以苯丙氨酸為直接前體，由一系列結(jié)構(gòu)基因編碼的合成酶催化形成，并經(jīng)過一系列的運輸后貯存于液泡中；漿果的大小、色澤、硬度均會隨著果實的成熟發(fā)生顯著變化［10-11］。近年來，關(guān)于藍果忍冬果實花青素的研究受到廣泛關(guān)注，Zhou 等［12］對藍果忍冬品種‘蓓蕾’果實中的花青素提取并純化，發(fā)現(xiàn)提取的花青素ABL-2 的主要成分是花青素-3-葡萄糖苷、花青素-3，5-二葡萄糖苷、花青素-3-蘆丁糖苷和芍藥-3-葡萄糖苷，且能夠顯著抑制細胞生長，誘導(dǎo)DNA 損傷，殺死腫瘤細胞，改善小鼠的生存狀態(tài)；Paereckait 等［13］利用分光光度法對11 種藍果忍冬果實中的酚類化合物和花青素含量等進行提取研究并比較不同藍果忍冬品種的抗菌活性。然而，目前有關(guān)藍果忍冬基因信息尚不明確，對藍果忍冬果實轉(zhuǎn)錄組的相關(guān)研究少有報道，嚴重限制了藍果忍冬的開發(fā)和利用。

本研究以藍果忍冬‘中科藍1 號’果實為試驗材料，使用高通量測序平臺對綠熟期果實和成熟期藍果忍冬果實進行轉(zhuǎn)錄組測序分析，通過組裝建立藍果忍冬轉(zhuǎn)錄組數(shù)據(jù)庫，探究藍果忍冬果實不同發(fā)育時期的功能基因和調(diào)控網(wǎng)絡(luò)，旨在為后續(xù)了解藍果忍冬果實發(fā)育提供理論參考。

1 材料與方法

1.1 材料

試驗材料為藍果忍冬‘中科藍1 號’品種，由中國科學(xué)院東北地理與農(nóng)業(yè)生態(tài)研究所培育，分別采取綠熟期（花后28 d，果實成長到體積不再發(fā)生變化，果皮顏色為綠色）和成熟期（花后38 d，果皮顏色全部變?yōu)樗{紫色，果實未軟化）兩個時期的果實（果皮和果肉），各取3 個生物學(xué)重復(fù)，用液氮速凍后保存在-80℃冰箱中用于后續(xù)轉(zhuǎn)錄組測序和花色苷含量的測定（圖1）。

圖1 不同發(fā)育時期的‘中科藍1 號’果實Fig.1 Fruits of ‘Zhongkelan No.1’ at different developmental stages

1.2 方法

1.2.1 不同時期的花色苷總量測定 參照高效液相色譜法（high performance liquid chromatography，HPLC）［14］，利用安捷倫液相色譜系統(tǒng)（Agilent Technologies 1260 Infinity Ⅱ Prime，US）檢測花青素種類及含量，流動相A 為含3% 甲酸的雙蒸水，流動相B 為100% 乙腈。洗脫程序為0 min 90%A，10% B；16 min 86% A，14% B；20 min 90% A，10% B；流速為1 mL/min，柱溫為30℃，檢測波長為520 nm，每個樣品設(shè)3 個重復(fù)；分別以矢車菊素-3，5-雙葡萄糖苷、飛燕草素-3-O-葡萄糖苷、飛燕草素-3-O-蕓香糖苷、矢車菊素-3-O-葡萄糖苷、矢車菊素-3-O-蕓香糖苷、天竺葵素-3-O-葡萄糖苷、芍藥素-3-O-葡萄糖苷、芍藥素-3-O-蕓香糖苷標準品繪制標準曲線，含量由峰面積計算方法得出，單位為mg/100 g FW。

1.2.2 RNA 提取 采用CTAB 法提取藍果忍冬不同發(fā)育時期的總RNA，經(jīng)PrimeScript One Step RT-PCR Kit Ver.2（TaKaRa）反轉(zhuǎn)錄為cDNA 并將cDNA 稀釋10 倍用于后續(xù)RT-qPCR 實驗。樣品總RNA 經(jīng)檢測合格后送至北京百邁客生物技術(shù)有限公司運用高通量平臺進行RNA-Seq 測序。轉(zhuǎn)錄組樣品編碼為G：“綠熟期”，B：“成熟期”。

1.2.3 轉(zhuǎn)錄組測序與分析 為了獲得高質(zhì)量的轉(zhuǎn)錄組數(shù)據(jù)，采用illumina novaseq 6000 平臺進行雙末端測序，截除掉原始序列接頭和引物序列，對低質(zhì)量的序列進行過濾處理，并利用Trinity 軟件對生成的高質(zhì)量測序數(shù)據(jù)進行組裝，得到轉(zhuǎn)錄本序列Unigene；利用BLAST 軟件在NR、GO、Swiss-Prot、COG、KEGG 和eggNOG4.5 數(shù)據(jù)庫中進行比較分析并預(yù)測Unigene 的氨基酸序列，使用HMMER軟件將Unigene 序列在Pfam 數(shù)據(jù)庫比對從而獲得Unigene 的功能信息，并使用 TBtools 軟件繪制熱圖。

1.2.4 差異表達基因的篩選和富集分析 用DESeq2軟件對藍果忍冬綠熟期和成熟期果實的轉(zhuǎn)錄組數(shù)據(jù)進行組間差異表達分析，篩選標準為FDR（False Discovery Rate）＜0.01 且FC（Fold Change）≥2，并利用Goatools 軟件對篩選出的差異表達基因進行COG 注釋分析、GO 功能富集分析和KEGG 富集分析。

1.2.5 差異表達基因的RT-qPCR 驗證 根據(jù)藍果忍冬轉(zhuǎn)錄組數(shù)據(jù)隨機篩選出10 個差異基因，在NCBI數(shù)據(jù)庫中比對并設(shè)計定量引物（表1）。選用LcACT1（MT344113）、LcTUB1（MT344114）作為內(nèi)參基因，利用Agilent AriaMx 實時熒光定量PCR 儀（Agilent Stratagene，USA）對其進行實時熒光定量PCR 反應(yīng)，程序為94℃預(yù)變性30 s，94℃變性15 s，57℃退火30 s，72℃延伸35 s，共44 個循環(huán)結(jié)束，65-95℃熔解曲線讀取，用2-ΔΔCT算法計算處理。

表1 轉(zhuǎn)錄組數(shù)據(jù)RT-qPCR 驗證引物序列Table 1 Primer sequences used in RT-qPCR validation of transcriptome data

2 結(jié)果

2.1 藍果忍冬‘中科’兩種時期的花色苷含量變化

利用高效液相色譜（high performance liquid chromatography，HPLC）檢測技術(shù)和紫外檢測器對520 nm 波長下進行掃描，檢測到的峰圖如圖2所示，可以看出，藍果忍冬‘中科’綠熟期果實未檢測出花色苷成分（圖2-a），成熟期果實中檢測到的花青素種類共有8 種（A1-A8）（圖2-b），總花色苷含量為598.55 mg/100 g FW。其中，矢車菊素-3-O-葡萄糖苷含量最高，占總花色苷的86.28%，是藍果忍冬‘中科’果實中主要的花色苷（表2），說明不同發(fā)育階段藍果忍冬‘中科’果實中花色苷類化合物組成和含量有較大差異。

表2 ‘中科藍1 號’成熟期果實中花色苷含量Table 2 Anthocyanins contents in the ripening stage fruits of ‘Zhongkelan No.1’

圖2 ‘中科藍1 號’不同發(fā)育時期果實高效液相色譜Fig.2 HPLC of ‘Zhongkelan No.1’ fruit at different developmental stages

2.2 轉(zhuǎn)錄組測序數(shù)據(jù)統(tǒng)計

選擇藍果忍冬綠熟期果實和成熟期果實各3 次生物學(xué)重復(fù)，共6 個樣品，轉(zhuǎn)錄組測序共得到40.32 Gb Clean Data，每個樣品Clean Data 均大于6.07 Gb，GC 含量為46.32%-49.38%，Q30 均值在93.89%及以上，表明測序數(shù)據(jù)較好，可用于后續(xù)分析（表3）。

表3 測序數(shù)據(jù)統(tǒng)計結(jié)果Table 3 Statistics of sequencing data

6 個藍果忍冬果實的轉(zhuǎn)錄組數(shù)據(jù)經(jīng)過質(zhì)量分析處理后共得到236 033 條轉(zhuǎn)錄本序列，平均長度為1 371.99 bp，N50 長2 186 bp；組裝后獲得81 777個單基因序列，且長度大于1 kb 的Unigene 共有17 339 個（表4），對所獲得的單基因序列與NR、Swiss-Prot、KEGG、COG、KOG、GO 和Pfam 數(shù)據(jù)庫進行功能注釋比對，注釋成功的Unigene 共有28 992 條，占全部Unigene 的35.45%（表5）。

表4 組裝結(jié)果統(tǒng)計Table 4 Assembly result statistics

表5 單基因功能注釋統(tǒng)計Table 5 Single-gene function annotation statistics

2.3 不同發(fā)育時期差異表達基因篩選結(jié)果

根據(jù)評估標準FDR（False Discovery Rate）＜0.01 且FC（Fold Change）≥2，共篩選出3 247 個差異表達基因，其中有1 642 個基因在上調(diào)表達，1 605 個基因在下調(diào)表達，從圖3可以清楚地看出，在兩個生長時期藍果忍冬果實中基因表達水平存在較大差異。

圖3 差異表達基因的火山圖Fig.3 Volcano map of differentially expressed genes

2.4 不同發(fā)育時期差異表達基因COG功能分類

利用COG 數(shù)據(jù)庫對藍果忍冬果實兩個發(fā)育時期的差異表達基因進行對比分析（圖4）發(fā)現(xiàn)，注釋到COG 通路中的差異表達基因有1 081 個，分布在23 條途徑中，占全部差異基因的37.88%，其中，碳水化合物運輸和代謝（carbohydrate transport and metabolism）是最大通路，所占數(shù)量最多，其次是一般功能預(yù)測（general function prediction only），差異基因在脂質(zhì)運輸和代謝（lipid transport and metabolism）通路和次級代謝物生物合成、運輸和分解代謝（secondary metabolites biosynthesis，transport and catabolism）通路占比也較高，說明在果實成熟發(fā)育過程中代謝物質(zhì)變化較大。

圖4 不同發(fā)育時期差異表達基因COG 注釋分類統(tǒng)計圖Fig.4 COG annotation classification of differentially expressed genes at different developmental stages

2.5 不同發(fā)育時期差異表達基因GO功能富集結(jié)果

對藍果忍冬果實兩個發(fā)育時期的差異表達基因進行對比分析，發(fā)現(xiàn)注釋到GO 通路中的差異表達基因有1 908 個，分別映射到51 個代謝通路中（圖5），在生物過程（biological process）中，主要富集在代謝過程（metabolic process）和細胞過程（cellular process），其次是單個有機體過程（single-organism process）；在細胞組分（cellula component）通路中，膜部分（membrane）占比最多、其次是細胞（cell）和細胞部分（cell part）；在分子功能（molecular function）通路中，催化活性（catalytic activity）和結(jié)合功能（binding）較其他通路顯著富集。花青素是植物次生代謝產(chǎn)物，糖苷衍生物。推測可能歸類于代謝過程（metabolic process）和催化活性（catalytic activity）的Unigenes 中，可能有參與藍果忍冬花青素生物合成積累差異相關(guān)的良好候選基因。

圖5 不同發(fā)育時期差異表達基因GO 功能富集結(jié)果Fig.5 Gene ontology annotation of differentially expressed genes at different developmental stages

2.6 不同發(fā)育時期差異表達基因KEGG功能富集結(jié)果

以q-value＜0.05 作為篩選差異表達基因在KEGG數(shù)據(jù)庫通路的標準，由圖6可知，2 854 個差異基因中共有1 034 個DEGs 注釋到KEGG 通路中，其中植物激素信號轉(zhuǎn)導(dǎo)（plant hormone signal transduction）途徑是最大通路，有55 個，占比10.07%；其次是氨基酸生物合成（biosynthesis of amino acid）和碳代謝（carbon metabolism）途徑，均為46 個，占比為8.42%；此外，淀粉和蔗糖代謝（starch and sucrose metabolism）通路與果實的發(fā)育也顯著富集，有37 個，占比6.78%。

圖6 不同發(fā)育時期差異表達基因KEGG 分類圖Fig.6 KEGG classification of differentially expressed gene at different developmental stages

花青素是類黃酮化合物，且與苯丙氨酸代謝途徑有關(guān)。對KEGG 代謝通路進行富集顯著性篩選，得出富集最顯著的20 條代謝通路（圖7），其中，苯丙氨酸，酪氨酸和色氨酸的生物合成（phenylalanine，tyrosine and tryptophan biosynthesis）、類黃酮生物合成（flavonoid biosynthesis）、半胱氨酸和蛋氨酸代謝（cysteine and methionine metabolism）這3 條代謝途徑被顯著富集，說明可能存在大量與藍果忍冬花青素生物合成積累相關(guān)的差異基因。

圖7 不同發(fā)育時期差異表達基因KEGG 功能富集結(jié)果Fig.7 KEGG analysis of differentially expressed genes at different developmental stages

2.7 藍果忍冬轉(zhuǎn)錄組測序數(shù)據(jù)的RT-qPCR驗證

為了驗證RNA-seq 的結(jié)果，對測序結(jié)果隨機挑選了10 個差異表達基因進行RT-qPCR 驗證，分析RT-qPCR 與RNA-seq 結(jié)果之間的相關(guān)性，每個基因進行3 次生物學(xué)重復(fù)。結(jié)果顯示（圖8），隨機挑選的10 個差異表達基因的表達量變化與轉(zhuǎn)錄組基因的變化完全一致，說明轉(zhuǎn)錄組的測序數(shù)據(jù)和RT-qPCR的數(shù)據(jù)之間具有良好的一致性，表明本研究RNAseq 數(shù)據(jù)具有較高的可靠性。

圖8 轉(zhuǎn)錄組測序數(shù)據(jù)的RT-qPCR 驗證Fig.8 Validation of transcriptome sequencing data by RT-qPCR

3 討論

近年來，隨著分子生物學(xué)技術(shù)的發(fā)展，新一代高通量測序技術(shù)已被廣泛應(yīng)用于生物分析和醫(yī)學(xué)研究等領(lǐng)域，因其具有成本低，測序速度快，能夠估計不同發(fā)育時期和不同組織中的整體基因的表達量等特點，是當前無參考基因組的非模式生物的經(jīng)濟有效的選擇［15-18］。目前，已有研究表明在藍果忍冬果實或種子中分離和鑒定出酚類物質(zhì)、抗氧化成分和揮發(fā)性物質(zhì)等多種生物活性成分，包括黃酮類化合物、花青素、脂肪酸、甘油三酯和總甾醇等化合物，并證實花青素-3-葡萄糖苷為其主要的抗氧化成分［19-22］。本研究基于RNA-seq 測序技術(shù)，對‘中科藍1 號’藍果忍冬兩個發(fā)育時期的果實進行轉(zhuǎn)錄組分析，序列經(jīng)過組裝之后，共獲得81 777 條平均長度為763.35 nt 的Unigene，與前人報道的‘蓓蕾’品種果實發(fā)育階段的基因數(shù)目相近［23］，比馬尾松、枇杷、荔枝等木本植物組裝的Unigene 多［24-26］，結(jié)果表明，組裝的質(zhì)量以及片段長度可以滿足轉(zhuǎn)錄組數(shù)據(jù)分析的需要。有28 992 條Unigene 得到注釋，占比35.45%，雖然比大葉女貞（50.20%）和藍莓（73.60%）的注釋率低，但與文冠果（36.15%）和紅掌（38.87%）相近［27-30］，造成注釋率低的原因可能是缺乏近緣物種的參考基因組信息，未得到注釋的基因也可能是藍果忍冬果實特異性新基因，有待進一步分析研究。

對比兩個測序樣本檢測到的差異表達基因發(fā)現(xiàn)，表達量上調(diào)的基因數(shù)目多于下調(diào)基因數(shù)目，這表明成熟期果實相對于綠熟期果實的樣本有更多的基因表達量。對這些基因進行KEGG 分析后發(fā)現(xiàn)有大量的色素形成相關(guān)途徑被顯著富集，包括與花青素代謝直接相關(guān)的苯丙氨酸、酪氨酸和色氨酸的生物合成途徑、類黃酮生物合成途徑和氨基酸生物合成途徑，表明藍果忍冬果實色澤調(diào)控的復(fù)雜性。GO 富集分析發(fā)現(xiàn)在果實綠熟期到轉(zhuǎn)色期果實成熟過程中，參與細胞過程和代謝的基因較多，同時涉及到的功能最多的基因為結(jié)合和催化功能，從分子水平說明了這一過程中代謝旺盛。藍果忍冬不同發(fā)育時期的差異表達基因在植物激素信號轉(zhuǎn)導(dǎo)、碳代謝、氨基酸生物合成、脂肪代謝等多條代謝途徑中富集，與類黃酮合成緊密相關(guān)，具體受哪些基因控制還需要做進一步信息學(xué)分析。

園藝產(chǎn)品品質(zhì)形成的調(diào)控包括多個方面，在植物激素相關(guān)研究中，表明GA 對果實成熟軟化［31］、有機酸及糖分積累［32］、果實著色［33］、果實大小和形狀［34］等方面都有一定影響。比如，趙榮華等［35］研究表明GA 水平的調(diào)控可同時影響葡萄果實的形態(tài)和品質(zhì)，另有研究發(fā)現(xiàn)生長素和赤霉素的相互作用在辣椒的生長發(fā)育過程中能平衡細胞分裂和細胞擴張［36］。在轉(zhuǎn)錄因子相關(guān)研究中，Medina-Puche等［37］研究表明YB9/FaMYB11、FabHLH3 和FaTTG1 在調(diào)控草莓花青苷合成過程中起關(guān)鍵作用，調(diào)控果皮花青苷的生物合成［38］。本研究發(fā)現(xiàn)藍果忍冬‘中科藍1 號’成熟期果實含有8 種花色苷成分，其中含量最高的花色苷為矢車菊素-3-O-葡萄糖苷，占總花色苷含量的86.28%，而綠熟期果實未檢測出花色苷成分，推測在果實成熟過程中存在調(diào)控花色苷含量的轉(zhuǎn)錄因子，還需進一步研究。

4 結(jié)論

本研究對‘中科藍1 號’藍果忍冬綠熟期和成熟期果實進行了轉(zhuǎn)錄組測序，共有40.32 Gb Clean Data，每個樣品Clean Data 均大于6.07 Gb，Q30 堿基百分比均在93.89%及以上，GC 含量為46.32%-49.38%；共篩選出3 247 個差異表達基因，包含1 642 個上調(diào)基因和1 605 個下調(diào)基因，為研究藍果忍冬果實轉(zhuǎn)色過程的的分子調(diào)控機制提供了參考。