武紅霞　許文天　羅純　姚全勝　王松標(biāo)　馬小衛(wèi)　詹儒林

摘 要 以紅芒6號(hào)為試材，對(duì)生長(zhǎng)發(fā)育各時(shí)期果皮與果肉樣品提取RNA后等量混合進(jìn)行轉(zhuǎn)錄組測(cè)序，共獲得了68 419 722個(gè)reads，包含6 157 744 980個(gè)核苷酸序列信息，平均讀長(zhǎng)90 bp，序列信息全部登陸NCBI數(shù)據(jù)庫(kù)（SRP035450）。對(duì)reads進(jìn)行拼接，共獲得124 002個(gè)Contig片段，進(jìn)一步組裝獲得54 207個(gè)Unigene片段，平均長(zhǎng)度為838 bp，閾值小于10-5的序列中共有42 515（78.43%）個(gè)unigenes被注釋到公共蛋白數(shù)據(jù)庫(kù)。其中，35 198個(gè)被注釋到生物學(xué)過(guò)程、分子功能和細(xì)胞組分GO類別的44個(gè)功能組，14 619個(gè)Unigene匹配到25類COG功能組，23 741個(gè)Unigene被富集到128個(gè)KEGG代謝通路，包括代謝途徑、次生代謝的生物合成、植物-病原物互作、植物激素信號(hào)轉(zhuǎn)導(dǎo)、苯丙氨酸生物合成、淀粉和蔗糖代謝、類黃酮類化合物合成等。本研究將為進(jìn)一步的芒果功能基因克隆、基因表達(dá)分析、分子標(biāo)記開(kāi)發(fā)等研究奠定基礎(chǔ)。

關(guān)鍵詞 芒果；轉(zhuǎn)錄組；功能注釋

中圖分類號(hào) S682.2 文獻(xiàn)標(biāo)識(shí)碼 A

芒果是著名的熱帶水果，素有“熱帶果王”之美譽(yù)。中國(guó)是主要的芒果生產(chǎn)國(guó)之一，2014年栽培面積14萬(wàn)hm2，產(chǎn)量130萬(wàn)t，主要分布在海南、廣西、云南、四川、廣東和福建等?。▍^(qū)）。前期研究主要集中于栽培生理研究，如果實(shí)成熟過(guò)程[1-4]、果實(shí)香氣揮發(fā)物質(zhì)[5-6]、抗氧化活性[7-8]、采后處理和果實(shí)品質(zhì)研究[9-10]，盡管芒果是中國(guó)熱區(qū)農(nóng)業(yè)產(chǎn)業(yè)的經(jīng)濟(jì)支柱，其基因組較小，但基因組信息相當(dāng)缺乏。近年來(lái)，高通量轉(zhuǎn)錄組測(cè)序技術(shù)的不斷發(fā)展和完善為開(kāi)展不同生物功能基因組學(xué)研究提供了全新的思路和方法，轉(zhuǎn)錄組測(cè)序是基于下一代高通量測(cè)序技術(shù)建立的一種高效、快捷分子生物學(xué)研究手段，使科技工作者能在組學(xué)水平上研究基因組序列未知的非模式生物，從整體水平了解植物在特定階段的基因功能及基因結(jié)構(gòu)，更加便利地提示特定生物學(xué)過(guò)程的分子機(jī)制[11-13]。近年來(lái)，已廣泛應(yīng)用于楊梅[14]、蘋果[15]、梨[16]、葡萄[17]、柑橘[18]等果樹(shù)，芒果研究方面，Dautt-Castro等[19]用轉(zhuǎn)錄組測(cè)序技術(shù)研究了與 Kent芒成熟過(guò)程相關(guān)的基因，在青熟和成熟果實(shí)中，2 306個(gè)基因差異表達(dá)。Luria 等[20]用轉(zhuǎn)錄組測(cè)序揭示了熱水處理對(duì)芒果影響的分子機(jī)制，基于基因表達(dá)譜測(cè)序技術(shù)鑒定了與熱處理相關(guān)基因，Azim等[21]報(bào)道了芒果葉綠體基因組，這些研究豐富了芒果生物信息數(shù)據(jù)。有關(guān)芒果果實(shí)生長(zhǎng)發(fā)育過(guò)程的轉(zhuǎn)錄組信息鮮有報(bào)道。本研究利用高通量測(cè)序技術(shù)平臺(tái)illumina Hiseq 2 000對(duì)芒果發(fā)育成熟過(guò)程的果實(shí)進(jìn)行轉(zhuǎn)錄組測(cè)序、拼接組裝，再用生物信息學(xué)的方法對(duì)得到的Unigene進(jìn)行注釋和功能分類，以期為芒果功能基因的發(fā)掘利用、特異miRNA的鑒定及功能分析等奠定基礎(chǔ)。

1 材料與方法

1.1 材料

試驗(yàn)于2012年3～7月在廣東省雷州市覃斗鎮(zhèn)芒果基地進(jìn)行，試驗(yàn)材料為15年生的紅芒6號(hào)芒果樹(shù)，試驗(yàn)地土壤屬磚紅壤，株行距4 m×4 m，樹(shù)體生長(zhǎng)良好。分別于幼果期（花后50 d）、膨大期（花后80 d）、采收期（青熟）與成熟期（常溫放置后成熟）采集芒果果實(shí)，每個(gè)時(shí)期采10～15個(gè)果，采集后迅速帶回實(shí)驗(yàn)室，將果皮與果肉分開(kāi)并用液氮速凍后保存于-80 ℃冰箱中備用。

1.2 方法

1.2.1 芒果果皮與果肉RNA提取 RNA提取參照Shan等[22]建立的方法。

1.2.2 芒果果實(shí)轉(zhuǎn)錄組測(cè)序 基于轉(zhuǎn)錄組學(xué)研究的優(yōu)勢(shì)，以紅芒6號(hào)芒果為材料，取不同發(fā)育階段的果皮與果肉樣本，分別提取總RNA后等量混合用于測(cè)序，轉(zhuǎn)錄組測(cè)序委托深圳華大基因公司完成。提取的總RNA使用DNaseI消化DNA后，用帶有Oligo（dT）的磁珠富集真核生物mRNA；加入打斷試劑在Thermomixer中適溫將mRNA打斷成短片段，以打斷后的mRNA為模板合成一鏈cDNA，然后配制二鏈合成反應(yīng)體系合成二鏈cDNA，并使用試劑盒純化回收、粘性末端修復(fù)、cDNA的3′末端加上堿基“A”并連接接頭，然后進(jìn)行片段大小選擇，最后進(jìn)行PCR擴(kuò)增；構(gòu)建好的文庫(kù)用Agilent 2 100 Bioanalyzer和ABI StepOnePlus Real-Time PCR System質(zhì)檢合格后，使用Illumina HiSeqTM 2000進(jìn)行測(cè)序。

1.2.3 數(shù)據(jù)組裝及基因功能注釋 將測(cè)序獲得的原始數(shù)據(jù)去除接頭序列、或N含量過(guò)高或低質(zhì)量的序列，得到凈讀長(zhǎng)，用Trinity軟件按Min contig length 100，Group pairs distance 250，path reinforcement distance 85，Min kmer cov 2參數(shù)進(jìn)行數(shù)據(jù)組裝。將具有重疊區(qū)域的reads連成更長(zhǎng)的Contig，采用Overlap的方法進(jìn)一步拼接成Unigene。

通過(guò)Blastx程序?qū)nigene序列比對(duì)到蛋白數(shù)據(jù)庫(kù)NR、Swiss-Prot、KEGG和COG（E值<0.00 001），并通過(guò)Blastn程序?qū)nigene與核酸數(shù)據(jù)庫(kù)NT（E值<0.00 001）進(jìn)行比對(duì)，獲得Unigene最高序列相似性的蛋白，便得到該Unigene的蛋白功能注釋信息。將Unigene與COG數(shù)據(jù)庫(kù)進(jìn)行比對(duì)，獲得Unigene可能的功能注釋及功能分類。根據(jù)nr注釋信息，利用Blast2GO軟件得到Unigene的GO注釋信息，并對(duì)所有Unigene用WEGO軟件進(jìn)行GO功能分類統(tǒng)計(jì)，根據(jù)KEGG注釋信息進(jìn)一步分析得到Unigene的Pathway途徑注釋。

預(yù)測(cè)編碼蛋白框：按NR、Swiss-Prot、KEGG和COG順序?qū)nigene序列與以上數(shù)據(jù)庫(kù)Blastx比對(duì)（E值<10-5），取Blast比對(duì)結(jié)果中排列最高的蛋白確定該Unigene的編碼區(qū)序列，將編碼區(qū)序列翻譯成氨基酸序列，得到該基因編碼區(qū)的核苷酸序列（序列方向5′→3′）和氨基酸序列。若Unigene比對(duì)不上以上數(shù)據(jù)庫(kù)，用軟件ESTScan[23]預(yù)測(cè)編碼區(qū)，得到其編碼區(qū)的核苷酸序列和氨基酸序列。

2 結(jié)果與分析

2.1 轉(zhuǎn)錄組組裝結(jié)果

2.1.1 果實(shí)轉(zhuǎn)錄組組裝概要 轉(zhuǎn)錄組測(cè)序共獲得82 817 950讀段，經(jīng)過(guò)去除雜質(zhì)和冗余處理，得到68 419 722條凈讀段，6.1Gb（6 157 744 980 bp）的有效數(shù)據(jù)，13×的測(cè)序深度，平均讀長(zhǎng)90 bp。其中讀長(zhǎng)大于20個(gè)堿基比例為95.27%，GC%值為45.52%，可以看出此次轉(zhuǎn)錄組測(cè)序結(jié)果較好，可為后續(xù)的數(shù)據(jù)組裝提供很好的原始數(shù)據(jù)。轉(zhuǎn)錄組測(cè)序組裝得到124 002個(gè)Contig，平均長(zhǎng)度為338 bp，N50值為623 bp。進(jìn)一步拼接得到54 207個(gè)Unigene，平均長(zhǎng)度838 bp，N50值為1 328。按序列相似度對(duì)Unigene進(jìn)行基因家族聚類，相似度大于70%的命為Clusters，共26 413個(gè)，Singleton序列有27 794個(gè)，其中E值小于10-5的序列共43 751個(gè)（表1）。54 207個(gè)Unigene的注釋信息見(jiàn)http：//www.sciencedirect.com/science/article/pii/S1874391914001523。

2.1.2 Contig和Unigene序列長(zhǎng)度分布 組裝的124 002個(gè)Contig中，核苷酸長(zhǎng)度在100～200 bp的序列數(shù)有74 842，比例達(dá)到了60.36%；200～300 bp有17 520，占14.13%；300～400 bp有8 714，占7.03%；400～500 bp 的有4 538條，占3.66%；而≥ 500 nt的共有18 388條，比例為14.83%。由此可見(jiàn)，Contig片斷長(zhǎng)度以100～200 bp為主（表2）。在Contig數(shù)據(jù)的基礎(chǔ)上，采用over-lap進(jìn)一步拼接，共獲得了54 207個(gè)unigene片段，其中核苷酸長(zhǎng)度在100～500 bp的序列數(shù)有25 509，比例高達(dá)47.06%，500～1 000 bp有12 727，比例為23.48%；1 000～1 500 bp有7 312，占13.49%；1 500～2 000 bp有 4 334條，占8.01%；≥2 000 bp有4 315，占7.96%（表3）。

2.1.3 Unigene的序列同源性分析 經(jīng)Blast比對(duì)無(wú)冗余核酸數(shù)據(jù)庫(kù)NR后，對(duì)Unigene進(jìn)行分析，從比對(duì)的E值來(lái)看（圖1-A），發(fā)現(xiàn)29.68%的序列具有一定的同源性（10-30

序列的相似性比較結(jié)果表明，相似性高于80%的序列高達(dá)27.8%，相似性在60%～80%的序列高達(dá)45.7%，而26.5%的序列相似性在18%～60%間（圖1-B）。對(duì)獲得的 Unigene進(jìn)行Blastn分析（圖1-C），大部分Unigene序列和葡萄、蓖麻、毛果楊、大豆、苜蓿、擬南芥的同源序列比對(duì)上。統(tǒng)計(jì)顯示29.9%的Unigene可比對(duì)到葡萄；27.6%的Unigene與蓖麻同源；依次為毛果楊（22.9%），大豆（6.3%），擬南芥（1.1%）等，這說(shuō)明芒果與葡萄的進(jìn)化關(guān)系較近。

2.2 Unigene功能注釋

通過(guò)Blastx程序在不同的數(shù)據(jù)庫(kù)進(jìn)行搜索比對(duì)（表4），42 515條unigene比對(duì)到NR數(shù)據(jù)庫(kù)，26 380條Unigene比對(duì)到Swissprot數(shù)據(jù)庫(kù)，23 741條Unigene用于KEGG途徑分析，14 619條進(jìn)行COG分析，35 198條被用于GO分類。

2.2.1 Unigene的GO功能分類分析 12 923條Unigene與GO數(shù)據(jù)庫(kù)中的基因具有相似性，且較多的單條Unigene與多種基因相對(duì)應(yīng)，建立了302 479條對(duì)應(yīng)關(guān)系，從而得到盡可能多的注釋和分類。芒果果實(shí)轉(zhuǎn)錄組中的Unigene根據(jù)GO功能可分為生物學(xué)過(guò)程、細(xì)胞組分和分子功能3大類58分支（表5）。生物學(xué)過(guò)程類別中，細(xì)胞代謝過(guò)程和代謝過(guò)程涉及的基因最多，分別有23 131和22 003 條，細(xì)胞組分類別中，組成細(xì)胞和細(xì)胞部分涉及的基因最多，均為28 286條，分子功能類別中，結(jié)合功能涉及的基因最多，有17 458條，其次是具有催化活性功能的基因，達(dá)到了16 842條，其他種類基因的表達(dá)豐度不盡相同。

2.2.2 Unigene的COG功能分類 將Unigene和COG數(shù)據(jù)庫(kù)進(jìn)行比對(duì)，27 869條Unigene與數(shù)據(jù)庫(kù)中的基因具有相似性，且較多的單條Unigene能夠與多種基因相對(duì)應(yīng)，建立了65 536條對(duì)應(yīng)關(guān)系（表6）。Unigene根據(jù)功能大致可分為25類，Unigene的COG功能種類比較全面，涉及了大多數(shù)的生命活動(dòng)，其中整體功能類的基因數(shù)量最多，有4 774 條，其次是轉(zhuǎn)錄功能（2 596），復(fù)制、重組和修復(fù)功能（2 268），翻譯后加工、蛋白質(zhì)折疊的伴侶蛋白（2 098），碳水化合物運(yùn)輸和代謝（1 627），核結(jié)構(gòu)相關(guān)基因類的數(shù)量最少，只有4條，而未知功能的序列多達(dá)1 361條。COG分類結(jié)果表明這些Unigene涉及許多的生物學(xué)功能，這些未知基因可能參與芒果果實(shí)生長(zhǎng)發(fā)育和成熟的調(diào)控。

2.2.3 Unigene的KEGG代謝途徑分類 芒果果實(shí)轉(zhuǎn)錄組數(shù)據(jù)有23 741個(gè)Unigene能注釋到128類代謝途徑中（表7），包括代謝途徑、次生代謝的生物合成通路、植物-病原互作、植物激素信號(hào)轉(zhuǎn)導(dǎo)、DNA剪切、RNA運(yùn)輸、蛋白質(zhì)在內(nèi)質(zhì)網(wǎng)的加工、淀粉和蔗糖代謝、苯丙氨酸的生物合成、萜類化合物合成，脂類代謝，RNA降解、黃酮和黃酮醇的生物合成等，其中代謝途徑和次生代謝的生物合成所注釋的Unigene數(shù)最多，分別為4 969和2 450。

2.3 編碼蛋白框（CDS）的核苷酸和氨基酸序列分析

[5] Pandit， S S， Chidley， H G， Kulkarni， R S， et al. Cultivar relationships in mango based on fruit volatile profiles[J]. Food Chem， 2009， 114（1）： 363-372.

[6] Pandit S S， Kulkarni R S， Chidley H G， et al. Changes in volatile composition during fruit development and ripening of ‘Alphonso mango[J]. J Sci Food Agr， 2009， 89（12）： 2 071-2 081.

[7] Ma X W， Wu H X， Liu L Q， et al. Polyphenolic compounds and antioxidant properties in mango fruits[J]. SciHortic， 2011， 129（1）： 102-107.

[8] MantheyJ A， Perkins-Veazie P. Influences of harvest date and location on the levels of beta-carotene， ascorbic acid， total phenols， the in vitro antioxidant capacity， and phenolic profiles of five commercial varieties of mango（Mangifera indica L.）[J]. J Agric Food Chem， 2009， 57（22）： 10 825-10 830.

[9] Ornelas-Paz J D， Yahia E M， Gardea A A. Changes in external and internal color during postharvest ripening of ‘Manila and ‘Ataulfo mango fruit and relationship with carotenoid content determined by liquid chromatography-APcI（+）-time-of-flight mass spectrometry[J]. Postharvest Biol Tec， 2008， 50（2-3）： 145-152.

[10] Jiang Y， Joyce D C. Effects of 1-methylcyclopropene alone and in combination with polyethylene bags on the postharvest life of mango fruit[J]. Ann Appl Biol， 2000， 137（3）： 321-327.

[11] Schuster S C. Next-generation sequencing transforms today's biology[J]. Nat Methods， 2008， 5（1）： 16-18.

[12] Shendure J， Ji H. Next-generation DNA sequencing[J]. Nat Biotechnol， 2008， 26（10）： 1 135-1 145.

[13] Wang Z， Gerstein M， Snyder M. RNA-Seq： a revolutionary tool for transcriptomics[J]. Nat Rev Genet， 2009， 10（1）： 57-63.

[14] Feng C， Chen M， Xu C J， et al. Transcriptomic analysis of Chinese bayberry （Myrica rubra） fruit development and ripening using RNA-Seq[J]. BMC Genomics， 2012， 13： 19.

[15] Vimolmangkang S， Zheng D， Han Y， et al. Transcriptome analysis of the exocarp of apple fruit identifies light-induced genes involved in red color pigmentation[J]. Gene， 2014， 534（1）： 78-87.

[16] Liu G Q， Li W S， Zheng P H， et al. Transcriptomic analysis of ‘Suli pear （Pyrus pyrifolia white pear group） buds during the dormancy by RNA-Seq[J]. BMC Genomics， 2012， 13（1）： 1-18.

[17] Sweetman C， Wong D C， Ford C M， et al. Transcriptome analysis at four developmental stages of grape berry （Vitis vinifera cv. Shiraz） provides insights into regulated and coordinated gene expression[J]. BMC Genomics， 2012， 13（1）： 105.

[18] Yu K Q， Xu Q， Da X L， et al. Transcriptome changes during fruit development and ripening of sweet orange（Citrus sinensis）[J]. BMC Genomics， 2012， 13（1）： 10.

[19] Dautt-Castro M， Ochoa-Leyva A， Contreras-Vergara C A， et al. Mango （Mangifera indica L.） cv. Kent fruit mesocarp de novo transcriptome assembly identifies gene families important for ripening[J]. Front plant sci， 2015， 6（10）： 62.

[20] Luria N， Sela N， Yaari M， et al. De-novo assembly of mango fruit peel transcriptome reveals mechanisms of mango response to hot water treatment[J]. BMC Genomics， 2014， 15（1）： 957.

[21] Azim M K， Khan I A， Zhang Y. Characterization of mango （Mangifera indica L.） transcriptome and chloroplast genome[J]. Plant mol biol， 2014， 85（1-2）： 193-208.

[22] Shan L L， Li X， Wang P， et al. Characterization of cDNAs associated with lignification and their expression profiles in loquat fruit with different lignin accumulation[J]. Planta， 2008， 227（6）： 1 243-1 254.

[23] Iseli C， Jongeneel C V， Bucher P. ESTScan： a program for detecting， evaluating， and reconstructing potential coding regions in EST sequences.[C]// Proc Int Conf Intell Syst Mol Biol， 1999： 138-148.

[24] Li X Y， Sun H Y， Pei J B， et al. De novo sequencing and comparative analysis of the blueberry transcriptome to discover putative genes related to antioxidants[J]. Gene， 2012， 511（1）： 54-61.

[25] Ong W D， Voo L Y， Kumar V S. De novo assembly， characterization and functional annotation of pineapple fruit transcriptome through massively parallel sequencing[J]. PLOS One， 2012， 7（10）： e46937.

[26] Pandit S S， Kulkarni R S， Giri A P， et al. Expression profiling of various genes during the fruit development and ripening of mango[J]. Plant Physiol Biochem， 2010， 48（6）： 426-433.

[27] Sherman A， Rubinstein M， Eshed R， et al. Mango （Mangifera indica L.） germplasm diversity based on single nucleotide polymorphisms derived from the transcriptome[J]. BMC Plant Biol， 2015， 15： 277.

[28] 羅 純， 武紅霞， 姚全勝， 等. 芒果轉(zhuǎn)錄組中SSR位點(diǎn)信息分析與引物篩選[J]. 熱帶作物學(xué)報(bào)， 2015， 36（7）： 1 261-1 266.

[29] 王 明， 應(yīng)東山， 王琴飛， 等. 基于454 GSFLX高通量測(cè)序的芒果微衛(wèi)星標(biāo)記特征分析[J]. 南方農(nóng)業(yè)學(xué)報(bào)， 2015（10）： 1 745-1 749.

[30] AltermannE， Klaenhammer T R. PathwayVoyager： pathway mapping using the Kyoto Encyclopedia of Genes and Genomes（KEGG） database[J]. BMC Genomics， 2005， 6： 60.

[31] Blencowe B J， Ahmad S， Lee L J. Current-generation high-throughput sequencing： deepening insights into mammalian transcriptomes[J]. Genes Dev， 2009， 23（12）： 1 379-1 386.