黑果枸杞果實(shí)發(fā)育過程中轉(zhuǎn)錄組測序分析

孟小偉，牛 赟，2，馬彥軍

（1.甘肅農(nóng)業(yè)大學(xué) 林學(xué)院，甘肅 蘭州 730070；2.淮陰師范學(xué)院，江蘇 淮安 223300）

高通量測序技術(shù)（Illumina測序技術(shù)）具有數(shù)據(jù)量大、準(zhǔn)確性好、快速便捷、運(yùn)行成本低等優(yōu)點(diǎn)[1]。對于缺乏全基因組信息的非模式物種，利用Illumina測序技術(shù)能對研究對象在某一特定狀態(tài)下基因轉(zhuǎn)錄表達(dá)信息進(jìn)行全面地分析[2]，有利于挖掘重要功能基因，而且還可以更加便利地揭示特定生物學(xué)過程的分子機(jī)制[3]，這已成為研究植物優(yōu)良性狀的重要技術(shù)手段[4]。近年來，利用這一技術(shù)已成功地對多種藥用植物在生理活動、生長發(fā)育以及次生代謝產(chǎn)物生物合成的相關(guān)基因進(jìn)行了發(fā)掘和鑒定[5]。

黑果枸杞Lycium ruthenicum為茄科Solanceae枸杞屬Lycium多年生灌木，是典型的鹽生植物，也是一種藥食兩用的優(yōu)良野生植物資源，在我國藏藥和維藥等傳統(tǒng)民族醫(yī)藥中廣泛應(yīng)用[6]，主要分布于中國西北地區(qū)鹽漬化土壤或荒漠環(huán)境區(qū)域，在荒漠區(qū)對維持生態(tài)平衡方面發(fā)揮著重要的作用，是鹽堿地治理的先鋒樹種[7]。近年來黑果枸杞受到人們極大的關(guān)注，主要是由于黑果枸杞成熟果實(shí)中含有大量的花色素，其穩(wěn)定性好，著色性強(qiáng)，是理想的食用天然花色苷[8]，同時黑果枸杞果實(shí)多糖，具有抗疲勞、降血糖等作用[9]，因此黑果枸杞是一種既有生態(tài)價值又有經(jīng)濟(jì)價值的樹種。目前對黑果枸杞果實(shí)研究主要在黑果枸杞花色苷組成[10]、花青素含量及功能[11]、黑果枸杞多糖的功效[12]等方面，對黑果枸杞果實(shí)轉(zhuǎn)錄組及有關(guān)基因方面的研究較少[13]。本研究通過高通量測序技術(shù)，對黑果枸杞果實(shí)發(fā)育過程中青果期、變色期及成熟期進(jìn)行測序，通過對黑果枸杞果實(shí)不同發(fā)育階段測序數(shù)據(jù)進(jìn)行分析，以期從整體水平了解黑果枸杞果實(shí)發(fā)育過程中基因表達(dá)模式和功能，為分析與黑果枸杞多糖和花青素等代謝相關(guān)的結(jié)構(gòu)基因和調(diào)控基因及SSR分子標(biāo)記的引物設(shè)計(jì)提供基礎(chǔ)研究資料。

1 材料與方法

1.1 試驗(yàn)材料

供試材料為黑果枸杞不同發(fā)育階段的果實(shí)，采自甘肅農(nóng)業(yè)大學(xué)校內(nèi)試驗(yàn)地（38°28′N，106°16′E），樹齡為5 a。試驗(yàn)始于2017年6月26日，選取并標(biāo)記長勢優(yōu)良、無病蟲害的黑果枸杞植株（扦插繁育），采摘3個不同發(fā)育時期的果實(shí)，分別為青果期（S1）、變色期（S2）、成熟期（S3）。每3棵枸杞植株的果實(shí)為一個重復(fù)，共重復(fù)3次，果實(shí)采集后用鋁箔紙包裝，置于液氮中速凍，并立即放置到-80 ℃冰箱中保存，備用。

1.2 試驗(yàn)方法

1.2.1 總RNA的提取、文庫構(gòu)建和轉(zhuǎn)錄組測序

用mirVana? miRNA ISOlation Kit (Ambion-1561)試劑盒提取和純化黑果枸杞果實(shí)總RNA。分別用NanoDrop 2000（美國，賽默飛）與Agilent 2100 Bioanalyzer（美國，安捷倫）進(jìn)行總RNA質(zhì)量和純度檢測。將符合質(zhì)量要求的黑果枸杞果實(shí)總RNA送樣至上海歐易生物醫(yī)學(xué)科技有限公司進(jìn)行轉(zhuǎn)錄組測序。

1.2.2 De novo拼接

De novo拼接是指在無參考基因組的情況下，將有overlap的reads連接成一個更長的序列，通過不斷的延伸，最后拼接出transcript。本研究采用Trinity（version:trinityrnaseq_r20131110）[14]軟件paired-end的拼接方法得到Transcript序列，依據(jù)序列的相似性及其長度，篩選出最長的一條Transcript作為Unigene，以此作為后續(xù)分析的參考序列[15]。

1.2.3 Unigene功能注釋

通過blastx將獲得的Unigene序列分別與非冗余數(shù)據(jù)庫（non-redundant，NR）、SwissProt和真核生物蛋白相鄰類的聚簇（Clusters of orthologous groups for eukaryotic complete genomes，KOG）庫進(jìn)行比對，取e＜1e-5的注釋基于BLAST（Basic Local Alignment Search Tool）尋找蛋白質(zhì)或核酸的相似序列。有兩個值，一個是S值，一個是E值。S值表示兩序列的同源性，分值越高表明它們之間相似的程度越大。E值就是S值可靠性的評價。它表明在隨機(jī)的情況下，其它序列與目標(biāo)序列相似度要大于S值的可能性，所以它的分值越低越好[16]。得到跟給定Unigene具有最高序列相似性的蛋白，通過利用KAAS（http://www.genome.jp/kaas-bin/kaas_main）得到Unigene的京都基因與基因組百科全書（KEGG）注釋信息，基于SwissProt的注釋結(jié)果，根據(jù)Uniprot ID映射GO term，從而得到該Unigene的基因本體（gene ontology，GO）蛋白功能注釋信息。

2 結(jié)果與分析

2.1 總RNA質(zhì)量檢測和組裝結(jié)果分析

由表1可看出，黑果枸杞果實(shí)在3個發(fā)育時期9個樣本，通過轉(zhuǎn)錄組測序共獲得49019994～50537322原始讀序（raw reads），各樣本的clean data均達(dá)到6 Gb以上，Q30（測序錯誤率小于0.1%）均在89.85%以上，GC含量均在40%以上。利用Trinity軟件對測序所得數(shù)據(jù)進(jìn)行合并組裝（表2），共獲得43 573個Unigene，其中長度在1 kb上的Unigene有19 453條，這些Unigene可作為后續(xù)實(shí)驗(yàn)研究的重點(diǎn)對象；Unigene平均長度為1 262.65 bp。一般N50長度超過800 bp就可認(rèn)為組裝序列完整性較好，本實(shí)驗(yàn)所組裝得到的序列N50為1 743 bp，這說明本次轉(zhuǎn)錄組測序數(shù)據(jù)量和質(zhì)量都較高，可以用于后續(xù)分析。

表1 有效數(shù)據(jù)評估統(tǒng)計(jì)Table 1 Valid data evaluation statistics

表2 組裝結(jié)果統(tǒng)計(jì)分析Table 2 Statistical analysis of assembly results

2.2 基因功能注釋

由表3可知，黑果枸杞果實(shí)發(fā)育過程中43 573個Unigene中有23 723個Unigene在不同數(shù)據(jù)庫中得到了注釋，占總Unigene的54.44%，還有19 850個Unigene在這些數(shù)據(jù)庫中沒有得到注釋，其中注釋到NR數(shù)據(jù)庫的Unigene最多，達(dá)到23 559個，占54.07%。共有3 726個Unigene在所有數(shù)據(jù)庫中都得以注釋（圖1）。

表3 單基因序列注釋統(tǒng)計(jì)Table 3 Statistics of single gene sequence annotation

2.2.1 Unigene的Nr及SwissProt數(shù)據(jù)庫比對分析

由表3～4可看出，黑果枸杞果實(shí)轉(zhuǎn)錄組所得Unigene有23 559個Unigene在NR數(shù)據(jù)庫中得到注釋。在NR數(shù)據(jù)庫注釋中與馬鈴薯Solanum tuberosum同源序列最多，為4 586個，占注釋Unigene的19.47%；與擬南芥Arabidopsis thaliana同源序列最少，為248個，占注釋Unigened的1.05 %；與其他物種的為2 154個，占注釋Unigene的9.14%。E值等于0的Unigene有2 993個，占注釋到NR數(shù)據(jù)庫Unigene的12.70%，E值在0～1e-100的Unigene最多，為9 024，占注釋到NR數(shù)據(jù)庫Unigene的38.30%。

圖1 Unigene中NR,SwissProt,KOG,GO,KEGG注釋交并集Fig.1 Unigene NR, SwissProt, KOG, GO, KEGG annotation intersection

由表3～5可看出，黑果枸杞果實(shí)轉(zhuǎn)錄組所得Unigene在SwissProt數(shù)據(jù)庫中得以注釋的Unigene有17 212個，其中E值為0注釋的Unigene有3 393個，占19.71%；E值 介 于0到1e-100的Unigene有2 970個，占注釋Unigene的17.26%。黑果枸杞果實(shí)轉(zhuǎn)錄組所得Unigene在SwissProt數(shù)據(jù)庫注釋中與擬南芥同源序列最多，為11 159個，占注釋Unigene的64.83%；與爪蟾Xenopus laevis同源序列最少，為195，占注釋Unigene的1.13%。

2.2.2 GO注釋和分類

通過GO數(shù)據(jù)庫對比分析，共有15 064個Unigene在GO數(shù)據(jù)庫中得到了注釋，獲得的GO數(shù)據(jù)庫注釋的Unigene 可分為分子功能（Molecular function）、細(xì) 胞 組 分（Cellular component）和生物過程（Biological process）3大類別，進(jìn)一步可細(xì)分為58個功能組。分子功能涉及20 197個Unigene，分為16個功能組，其中結(jié)合（Binding，8 959個）和催化活性（Catalytic activity，7 557個）含Unigene較多；細(xì)胞組分涉及51 077個Unigene，分為19個功能組，其中細(xì)胞（Cell，11 571個）、細(xì)胞部分（Cell part，11 561個）、細(xì)胞器（Organelle，8 998）涉及的Unigene較多；生物過程涉及56 167個Unigene，分為23個功能組，其中細(xì)胞進(jìn)程（cellular process，10 157）和代謝進(jìn)程（metabolic process，8 695）涉及的Unigene較多。

表4 黑果枸杞Unigene的NR比對分析Table 4 Nr mapping analysis of Lycium ruthenicum Unigene

表5 黑果枸杞Unigene的SwissProt比對分析Table 5 SwissProt mapping analysis of Lycium ruthenicum Unigene

表6 黑果枸杞的Unigene GO功能分類Table 6 GO functional categories of Lycium ruthenicum Unigene

2.2.3 KOG注釋和分類

為了進(jìn)一步分析黑果枸杞果實(shí)發(fā)育時期轉(zhuǎn)錄組Unigene的功能，進(jìn)行了KOG功能分類分析，共有13 128個Unigene獲得14 487個KOG注釋，平均每條1.1個。分類結(jié)果如表7所示，共獲得25個不同的功能分類。一般功能預(yù)測（General function prediction only）的Unigene為4 180個，是最大的功能類群；其次是翻譯后修飾、蛋白質(zhì)翻轉(zhuǎn)和分子伴侶（Posttranslation almodification,protein turnover, chaperones）有1 391個，信號轉(zhuǎn)導(dǎo)類機(jī)制（Signal transduction mechanisms）次之，有1 224條，最少的是細(xì)胞運(yùn)動（Cell motility）功能類別，僅7個。黃酮類代謝途徑所屬的Q類（次生代謝產(chǎn)物生物合成、運(yùn)輸和代謝）共獲得了520個Unigene注釋。

2.2.4 KEGG代謝通路分析

KEGG是系統(tǒng)分析基因產(chǎn)物在細(xì)胞中的代謝途徑以及基因產(chǎn)物功能的數(shù)據(jù)庫。根據(jù)KEGG數(shù)據(jù)庫的注釋信息能進(jìn)一步得到Unigene的pathway注釋。結(jié)合KEGG數(shù)據(jù)庫，黑果枸杞果實(shí)轉(zhuǎn)錄組注釋到KEGG的4 951個Unigenes，獲得了9 754個KEGG注釋。注釋的Unigene參與的代謝通路可分為4大類別23個子類。由表8可看出，4大類別代謝通路中，與代謝（Metabolism）相關(guān)的通路獲得4 351個Unigene注釋，遺傳信息處理（Genetic information processing）相關(guān)的通路獲得1 841個Unigene注釋，細(xì)胞過程（Cellular processes）相關(guān)的通路獲得1 694個Unigene注釋，環(huán)境信息處理（Environmental information processing）相關(guān)的通路獲得1 868個Unigene注釋。進(jìn)一步細(xì)分為23個子類代謝通路，其中信號傳導(dǎo)（Signal transduction）獲得Unigene注釋最多，為1 833個，其次為碳水化合物代謝（Carbohydrate metabolism），為1 068個。次生代謝物生物合成（Biosynthesis of other secondary metabolites） 有321個Unigene。

表7 黑果枸杞Unigene的KOG功能分類Table 7 KOG functional categories of Lycium ruthenicum Unigene

以KEGG pathway數(shù)據(jù)庫作為參考，可將注釋到KEGG數(shù)據(jù)庫中的9 754個Unigene定位到215個具體的代謝途徑分支。表9列出注釋到KEEG前10的代謝途徑以及類黃酮生物合成、苯丙醇生物合成和花青素生物合成途徑。由表9可看出，核糖體（Ribosome）代謝途徑注釋到的Unigenes數(shù)量最多，有263個Unigenes；其次為內(nèi)質(zhì)網(wǎng)蛋白加工途徑（Protein processing in endoplasmic reticulum），為254個Unigene。注釋到類黃酮生物合成途徑中的Unigene有45個，苯丙醇生物合成途徑的Unigene有133個，花青素生物合成途徑的Unigene有1個，黃酮和黃酮醇的生物合成途徑（Flavone and flavonol biosynthesis）Unigene有2個。

表8 黑果枸杞Unigene的KEGG功能分類Table 8 KEGG functional categories of Lycium ruthenicum Unigene

表9 黑果枸杞Unigene的KEGG代謝途徑分析Table 9 KEGG analysis of metabolic pathway of Lycium ruthenicum Unigene

2.3 SSR信息分析

利用軟件MISA對黑果枸杞果實(shí)發(fā)育過程中轉(zhuǎn)錄組測序所獲得的Unigene進(jìn)行SSR預(yù)測，結(jié)果見表10。由表10可知，共有16 815個SSR位點(diǎn)：單核苷酸SSR最多，為12 262個，占72.92%；6核苷酸SSR最少，為0.04%。重復(fù)單元重復(fù)出現(xiàn)的次數(shù)大于11次以上最多，為5 297，占31.50%；重復(fù)單元重復(fù)出現(xiàn)9次的最少，為287，占1.71%。

表10 SSR不同重復(fù)基序分布Table 10 Distribution of different repeat motifs in SSR

3 討 論

第二代高通量測序技術(shù)因測序時間短、成本低和所獲得數(shù)據(jù)量大等優(yōu)點(diǎn)，被廣泛應(yīng)用于非模式生物分子生物學(xué)研究中。例如對文冠果[17]、杜仲[18]、藍(lán)靛果忍冬[19]等非模式植物果實(shí)發(fā)育過程中轉(zhuǎn)錄組測序分析，全面了解了文冠果、杜仲、藍(lán)靛果忍冬等果實(shí)基因表達(dá)情況。本研究利用第二代高通量測序技術(shù)對黑果枸杞果實(shí)發(fā)育不同階段進(jìn)行轉(zhuǎn)錄組測序，建立了黑果枸杞的轉(zhuǎn)錄組數(shù)據(jù)庫，獲得了大量基礎(chǔ)數(shù)據(jù)。對這些數(shù)據(jù)進(jìn)行分析、序列組裝，所得序列在不同數(shù)據(jù)庫的功能注釋及分類、代謝途徑等分析，揭示黑果枸杞果實(shí)不同發(fā)育階段整體基因表達(dá)特征，為黑果枸杞果實(shí)發(fā)育過程中的分子生物學(xué)研究提供了基礎(chǔ)資料。

本研究所得到的黑果枸杞果實(shí)發(fā)育過程中的43 573個Unigene有23 723個Unigene在NR、SwissProt、GO、KOG、KEGG等不同數(shù)據(jù)庫中得到了注釋，占總Unigened的54.44%，還有19 850個Unigene在這些數(shù)據(jù)庫中沒有得到注釋。這一結(jié)果在許多非模式生物轉(zhuǎn)錄組測序中都存在[20]，這主要是由于非模式生物缺乏基因組方面研究的基礎(chǔ)資料，使得部分Unigene在NR、SwissProt、GO、KOG、KEGG等數(shù)據(jù)庫中無法得以注釋，而這些未得到注釋的Unigene有可能是非模式生物特有的基因，因此需要對這些未得到注釋的Unigene進(jìn)一步從結(jié)構(gòu)、功能等方面來深入研究，從而揭示他們在該生物生長發(fā)育過程中所發(fā)揮的作用。由于黑果枸杞現(xiàn)有遺傳信息量少的原因，測序結(jié)果中所得到的Unigene即使能在NR、SwissProt、GO、KOG、KEGG等數(shù)據(jù)庫中得以注釋，也有一部分Uuigene是注釋到動物和人類基因數(shù)據(jù)庫中。本研究中SwissProt注釋時就有2.04%的Unigene與智人同源，1.55%的Unigene與小鼠同源，還有部分Unigene同黑腹果蠅和爪蟾同源。在KEGG代謝通路pathway分析中，有部分Unigenere與人類疾病路徑相關(guān)。這一結(jié)果在許多其他生物測序中都存在[21]，這主要是因?yàn)檫@些非模式生物現(xiàn)有EST數(shù)據(jù)很少，而為了確定這些非模式生物在測序過程中獲得大量Unigene的功能，只能與其他植物甚至動物和人類的基因數(shù)據(jù)庫進(jìn)行比對，因此產(chǎn)生了這樣的結(jié)果。

本研究對黑果枸杞果實(shí)發(fā)育過程中3個階段的43 573條Unigene進(jìn)行SSR搜索，得到16 815個SSR位點(diǎn)，出現(xiàn)頻率為38.59%，高于野三七（頻率為16.86%）[22]、臘梅（頻率為12.35%）[23]、紅豆杉（頻率為2.07%）[24]等。這表明不同植物轉(zhuǎn)錄組測序結(jié)果中SSR的發(fā)生頻率差異較大，產(chǎn)生這種差異可能與所研究物種、轉(zhuǎn)錄組測序所得數(shù)據(jù)量以及檢索標(biāo)準(zhǔn)等因素有關(guān)[25]。在黑果枸杞果實(shí)轉(zhuǎn)錄組中，SSR以單核苷酸重復(fù)基序（72.92%）為主，這與‘芙蓉李’研究中SSR主要類型為單核苷酸重復(fù)基序（42.19%）為主一致[26]。目前大量研究結(jié)果表明利用轉(zhuǎn)錄組數(shù)據(jù)開發(fā)SSR標(biāo)記是可行的[27]。本研究基于黑果枸杞果實(shí)發(fā)育過程中轉(zhuǎn)錄組的高通量測序數(shù)據(jù)，從RNA水平上有針對性地進(jìn)行了黑果枸杞特異性SSR位點(diǎn)的檢索和評價，為進(jìn)一步開發(fā)新的黑果枸杞功能基因及調(diào)控基因的SSR分子標(biāo)記、黑果枸杞遺傳多樣性研究、黑果枸杞遺傳圖譜的構(gòu)建、分子標(biāo)記輔助育種等研究提供了基礎(chǔ)資料。

4 結(jié) 論

本研究利用Illumina HiSeqTM2500測序儀對黑果枸杞果實(shí)發(fā)育過程中轉(zhuǎn)錄組進(jìn)行測序，拼接出Unigenes 43 573條。注釋到GO、KEGG、KOG、NR和Swiss-Prot數(shù)據(jù)庫上的Unigenes總數(shù)為23 723條，占總Unigened的54.44%，共有3 726個Unigene在所有數(shù)據(jù)庫中都得以注釋，還有19 850個Unigene在這些數(shù)據(jù)庫中沒有得到注釋。獲得的GO數(shù)據(jù)庫注釋的Unigene可分為細(xì)胞組分、分子功能和生物過程3大類別，進(jìn)一步可細(xì)分為62個功能組。KOG功能分類中共獲得25個不同的功能分類，黃酮類代謝途徑所屬的Q類（次生代謝產(chǎn)物生物合成、運(yùn)輸和代謝）共獲得了520個Unigene注釋。黑果枸杞果實(shí)轉(zhuǎn)錄組注釋到KEGG的4 951個Unigenes獲得了9 754個KEGG注釋，參與的代謝通路可歸為4大類別23個子類。將注釋到KEGG數(shù)據(jù)庫中的9 754個Unigene定位到215個具體的代謝途徑分支，其中注釋到類黃酮生物合成途徑中的Unigene有45個，苯丙醇生物合成途徑的Unigene有133個，花青素生物合成途徑的Unigene有1個，黃酮醇的生物合成途徑的Unigene有2個。在黑果枸杞果實(shí)發(fā)育過程的轉(zhuǎn)錄組中發(fā)現(xiàn)16 815個SSR位點(diǎn)，最多的為單核苷酸SSR，占72.92%，重復(fù)單元重復(fù)出現(xiàn)的次數(shù)大于11次以上最多，占31.50%。