徐德宏，羅月芳，江靈敏，孟蕾，譚朝陽

(湖南中醫(yī)藥大學(xué) 藥學(xué)院 生物工程實(shí)驗(yàn)室，湖南 長沙 410208)

女貞LigustrumlucidumAit.為木犀科女貞屬植物，其干燥成熟果實(shí)(即女貞子)是一種常用中藥，具有滋補(bǔ)肝腎、明目烏發(fā)的功能，用于肝腎陰虛、眩暈耳鳴、腰膝酸軟、須發(fā)早白、目暗不明、內(nèi)熱消渴、骨蒸潮熱[1]。經(jīng)化學(xué)成分分析和活性檢測發(fā)現(xiàn)，女貞果實(shí)中含有多種具有藥理活性的次級(jí)代謝物[2-3]，如萜類、黃酮類等，因此近年來女貞果實(shí)在醫(yī)藥及保健方面?zhèn)涫苋藗兊年P(guān)注。目前，女貞果實(shí)的研究主要集中在次級(jí)代謝物的成分分析、提取和活性檢測，而關(guān)于其功能基因發(fā)現(xiàn)和注釋、物質(zhì)合成相關(guān)途徑分析、分子標(biāo)記開發(fā)等轉(zhuǎn)錄組學(xué)方面的研究卻未有報(bào)道，故現(xiàn)階段有必要采用適當(dāng)?shù)姆椒▽?duì)以上內(nèi)容展開深入系統(tǒng)的研究。

近年來，各種組學(xué)研究相繼興起，其中第二代高通量測序技術(shù)特別是RNA測序技術(shù)已廣泛應(yīng)用于生物體轉(zhuǎn)錄組基因表達(dá)分析，采用該技術(shù)能全面快速地獲取研究對(duì)象在某一狀態(tài)下基因組轉(zhuǎn)錄信息，從中挖掘重要的功能基因，揭示不同生物學(xué)性狀的分子機(jī)制[4-5]。Illumina/Solexa是目前轉(zhuǎn)錄組測序中最常使用的一種RNA高通量測序技術(shù)，與傳統(tǒng)的EST sanger測序法相比，具有高通量、低成本和高靈敏度等無法比擬的優(yōu)勢，同時(shí)可對(duì)低豐度的基因表達(dá)進(jìn)行檢測[6]。此外，隨著生物信息學(xué)的不斷發(fā)展，特別是序列組裝軟件的開發(fā)，主要蛋白質(zhì)數(shù)據(jù)庫的更新和轉(zhuǎn)錄組分析方法的建立[7-9]，為非模式生物進(jìn)行轉(zhuǎn)錄組高通量測序、從頭拼接組裝和unigene的功能注釋提供了極大的方便，因而Illumina/Solexa第二代高通量測序技術(shù)已廣泛應(yīng)用到了非模式生物的轉(zhuǎn)錄組研究中。

本研究將利用Illumina/Solexa Hiseq2000測序平臺(tái)，對(duì)女貞果實(shí)進(jìn)行轉(zhuǎn)錄組測序，將測序得到的大量數(shù)據(jù)進(jìn)行拼接與組裝，并應(yīng)用生物信息學(xué)方法對(duì)所得序列進(jìn)行功能注釋、功能分類、代謝途徑分析和簡單重復(fù)序列(simple sequence repeat，SSR)位點(diǎn)查找等研究，從而在轉(zhuǎn)錄組水平系統(tǒng)研究女貞果實(shí)，并為進(jìn)一步開展女貞果實(shí)的分子標(biāo)記開發(fā)和重要次級(jí)代謝產(chǎn)物的生物合成奠定基礎(chǔ)。

1 材料與方法

1.1 女貞果實(shí)轉(zhuǎn)錄組文庫構(gòu)建、測序和組裝

女貞果實(shí)于2015年8月采摘于湖南中醫(yī)藥大學(xué)湘杏學(xué)院校區(qū)內(nèi)的女貞樹上，用Trizol試劑分別提取3棵女貞樹成熟果實(shí)的總RNA，每個(gè)樣本取等量混合組成RNA池，然后用磁珠富集含poly(A)的mRNA并將其打斷成短序列，以這些短序列為模板用六堿基隨機(jī)引物依次合成第1條cDNA鏈和第2條cDNA鏈。cDNA鏈經(jīng)試劑盒純化并加EB緩沖液洗脫之后做末端修復(fù)、加poly(A)并連接測序接頭，然后用瓊脂糖凝膠電泳選擇200～500 bp片段進(jìn)行PCR擴(kuò)增，最后對(duì)建好的測序文庫用Illumina/Solexa Hiseq2000測序平臺(tái)進(jìn)行測序。鑒于測序平臺(tái)數(shù)據(jù)錯(cuò)誤率對(duì)結(jié)果的影響，對(duì)原始reads數(shù)據(jù)進(jìn)行質(zhì)量預(yù)處理得到clean reads，即先使用cutadapt軟件去除3′端測序接頭[10]；再利用滑動(dòng)窗口法去除低質(zhì)量片段：質(zhì)量閾值20(錯(cuò)誤率=1%)，窗口大小10 bp，長度閾值35 bp；最后切除reads中含N部分序列：長度閾值35 bp。參照Grabherr 等提出方法[11]，將所得clean reads通過Trinity軟件(版本號(hào)trinityrnaseq_r2013-02-25)進(jìn)行從頭拼接組裝，即通過clean reads之間的重疊信息組裝得到重疊群(Contigs)，然后局部組裝得到轉(zhuǎn)錄本(Transcripts)，最后用 Easycluster軟件對(duì)Transcripts進(jìn)行同源聚類以得到單基因簇(Unigene)[12-13]。

1.2 Unigene功能注釋、GO分類和代謝路徑分析

通過Blastx在e值<0.000 01的條件下，將Unigene序列比對(duì)到蛋白數(shù)據(jù)庫Nr(Non-redundant protein database，非冗余蛋白數(shù)據(jù)庫)、Swiss-Prot(SwissProt protein database，蛋白質(zhì)序列數(shù)據(jù)庫)，并通過Blastn在e值<0.000 01的條件下，將Unigene比對(duì)到核酸數(shù)據(jù)庫Nt(Nucleotide collection，核酸數(shù)據(jù)庫)，得到與給定Unigene具有最高序列相似性且功能已知的蛋白，從而實(shí)現(xiàn)Unigene的功能注釋。根據(jù)Nr數(shù)據(jù)庫的功能注釋結(jié)果，使用Blast2GO軟件得到unigene的GO條目[14]，然后用WEGO軟件對(duì)所有的unigene進(jìn)行GO功能分類統(tǒng)計(jì)[15]。KEGG 通路注釋是在evalue<0.000 01的條件下，利用blastx把unigene比對(duì)到KEGG數(shù)據(jù)庫，然后根據(jù)比對(duì)結(jié)果查詢與次級(jí)代謝產(chǎn)物合成通路相關(guān)的unigene。

1.3 SSR位點(diǎn)搜索和分析

利用MISA軟件在所有unigene 中搜索SSR位點(diǎn)，參數(shù)設(shè)置如下：二核苷酸至少重復(fù)次數(shù)為6，三核苷酸、四核苷酸、五核苷酸和六核苷酸至少重復(fù)次數(shù)均為4[16]。

2 結(jié)果與分析

2.1 測序結(jié)果評(píng)估和讀序列(reads)拼接組裝

采用Illumina Hiseq2000測序平臺(tái)，對(duì)女貞果實(shí)mRNA逆轉(zhuǎn)錄構(gòu)建的cDNA文庫進(jìn)行測序，共獲得32 856 614條clean reads，且本次測序總的cDNA堿基數(shù)約為4 928 492 100 bp(4.93 Gb)，故每條reads的平均長度為150 bp。Q20及GC含量百分比依次為97.94%和43.97%(見表1)。由以上結(jié)果表明，轉(zhuǎn)錄組測序數(shù)據(jù)量和質(zhì)量都較高，為后續(xù)的組裝提供了較好的原始數(shù)據(jù)。

利用Trinity軟件按paired-end方法對(duì)reads進(jìn)行拼接組裝，獲得229 286個(gè)重疊群(contig)，平均長為322 bp，總長為73 933 206 bp，N50為475 bp(見表1)，其中小于200 bp 的重疊群數(shù)量占59.32%，200～3000 bp以及3000 bp以上各占40.16%和0.52%，具體長度分布見圖1。所得重疊群再次組裝共獲得163 092條unigene，平均長為665.68 bp，總長為108 567 136 bp(108.57 Mb)，N50為1345 bp(見表1)，其中100～500、500～2000、2000 bp以上的各占64.01%、28.67%、7.32%，具體長度分布見圖2。

表1 女貞果實(shí)轉(zhuǎn)錄組測序和組裝結(jié)果

圖1 Contig長度分布圖

圖2 Unigene長度分布圖

2.2 序列比對(duì)和功能注釋

將所獲得的unigene與公共數(shù)據(jù)庫Nr和Swiss-Prot進(jìn)行Blastx比對(duì)，并利用Blastn比對(duì)Nt核酸數(shù)據(jù)庫，結(jié)果顯示共有83 903條unigene獲得了基因注釋，占unigene總數(shù)的51.45%。

2.3 功能分類和代謝途徑分析

2.3.1 GO和COG功能分類 為進(jìn)一步揭示女貞果實(shí)中unigene的功能分類，組裝產(chǎn)生的163 092條unigene在GO(gene ontology)和COG(cluster of orthologous groups)數(shù)據(jù)庫中分別進(jìn)行比對(duì)分析。在GO分類中，共有16 876 條unigene被注釋，其中12 203條unigene歸入生物學(xué)過程類別，3015條unigene歸入分子功能類別以及1658條unigene歸入細(xì)胞組分類別(見圖3)。在這3個(gè)類別中，生物學(xué)過程類別包括22個(gè)功能組，分子功能類別包括15個(gè)功能組，細(xì)胞組分類別則包括16個(gè)功能組(見圖3)。代謝過程是生物學(xué)過程類別中最大功能組(9468條)；細(xì)胞和細(xì)胞成分是細(xì)胞組分類別中最大功能組(都為5690條)；催化活性是分子功能類別中最大功能組(9519條)。通過以上GO分析，使我們?cè)谵D(zhuǎn)錄水平對(duì)女貞果實(shí)功能基因的分類情況有了一個(gè)系統(tǒng)的了解。在COG分類中可將unigene分為25個(gè)類別(見圖4)，其中unigene條數(shù)最多的10個(gè)類別依次是一般功能預(yù)測類(6875條)；復(fù)制、重組與修復(fù)類(3969條)；轉(zhuǎn)錄類(3733條)；信號(hào)轉(zhuǎn)導(dǎo)機(jī)制類(2926條)；蛋白質(zhì)翻譯后修飾、折疊和分子伴侶類(2632條)；碳水化合物運(yùn)輸與代謝類(2198條)；氨基酸運(yùn)輸與代謝類(1667條)；功能未知類(1371條)；細(xì)胞周期調(diào)控與分裂、染色體重排類(1201條)和無機(jī)離子運(yùn)輸與代謝類(1197條)，而細(xì)胞外結(jié)構(gòu)則是最小的類，僅包含4條unigene。通過以上COG分析，我們可以看到在女貞果實(shí)中，unigene涉及的功能類別主要是基因表達(dá)、蛋白質(zhì)合成和物質(zhì)代謝相關(guān)方面，這就正好符合果實(shí)作為一種營養(yǎng)器官，利用以上功能活動(dòng)為自身積累營養(yǎng)物質(zhì)這一基本功能。

2.3.2 KEGG代謝途徑分析 KEGG(Kyoto encyclopedia of genes and genomes)可對(duì)細(xì)胞內(nèi)unigene參與的代謝途徑進(jìn)行分析，本研究將獲得的unigene比對(duì)到此數(shù)據(jù)庫中，結(jié)果顯示女貞果實(shí)的unigene主要參與了21類代謝途徑(見圖5)，每條途徑中又包括下一級(jí)具體的代謝過程，其中涉及主要的次級(jí)代謝產(chǎn)物黃酮類合成途徑(Ko00941)含有258條unigene，萜類合成途徑(Ko00900，Ko00902，Ko00904，Ko00909)含有929條unigene。以上這些unigene首次在女貞代謝途徑中被發(fā)現(xiàn)，這為日后進(jìn)一步驗(yàn)證它們的功能、并利用它們?cè)隗w外生物合成相關(guān)的藥用活性物質(zhì)奠定了基礎(chǔ)。

注：1.生物黏附；2.生物調(diào)節(jié)；3.細(xì)胞成分組織或生物發(fā)生；4.細(xì)胞進(jìn)程；5.發(fā)育進(jìn)程；6.定位的建立；7.生長；8.免疫系統(tǒng)進(jìn)程；9.定位；10.運(yùn)動(dòng)；11.代謝進(jìn)程；12.多有機(jī)體進(jìn)程；13.多細(xì)胞有機(jī)體進(jìn)程；14.生物進(jìn)程的負(fù)調(diào)控；15.生物進(jìn)程的正調(diào)控；16.生物進(jìn)程調(diào)控；17.復(fù)制；18.生殖進(jìn)程；19.刺激應(yīng)答；20.節(jié)律進(jìn)程；21.信號(hào)；22.單有機(jī)體進(jìn)程；23.細(xì)胞；24.細(xì)胞連接；25.細(xì)胞成分；26.細(xì)胞外基質(zhì)；27.細(xì)胞外基質(zhì)成分；28.胞外區(qū)域；29.大分子復(fù)合物；30.膜；31.膜成分；32.膜關(guān)閉內(nèi)腔；33.核區(qū)；34.細(xì)胞器；35.細(xì)胞器成分；36.共質(zhì)體；37.病毒體；38.病毒體成分；39.抗氧化活性；40.結(jié)合；41.催化活性；42.電子載體活性；43.酶調(diào)節(jié)活性；44.鳥苷酸交換因子活性；45.金屬伴侶活性；46.分子轉(zhuǎn)導(dǎo)活性；47.核酸結(jié)合的轉(zhuǎn)錄因子活性；48.營養(yǎng)庫活性；49.蛋白結(jié)合的轉(zhuǎn)錄因子活性；50.蛋白尾；51.受體活性；52.結(jié)構(gòu)分子活性；53.轉(zhuǎn)運(yùn)活性。圖3 女貞果實(shí)unigene的GO分類

注：1.翻譯，核糖體結(jié)構(gòu)與生物發(fā)生；2.轉(zhuǎn)錄；3.信號(hào)傳導(dǎo)機(jī)制；4.次級(jí)產(chǎn)物合成，運(yùn)輸及代謝；5.加工與修飾；6.復(fù)制、重組與修復(fù)；7.蛋白質(zhì)翻譯后修飾與轉(zhuǎn)運(yùn)，分子伴侶；8.核苷酸運(yùn)輸與代謝；9.核結(jié)構(gòu)；10.脂類轉(zhuǎn)運(yùn)與代謝；11.細(xì)胞內(nèi)轉(zhuǎn)運(yùn)，分泌和小泡運(yùn)輸；12.無機(jī)離子轉(zhuǎn)運(yùn)和代謝；13.一般功能預(yù)測；14.功能未知；15.胞外結(jié)構(gòu)；16.能量生成和轉(zhuǎn)換；17.防御機(jī)制；18.細(xì)胞骨架；19.輔酶運(yùn)輸和代謝；20.染色質(zhì)結(jié)構(gòu)與變化；21.細(xì)胞壁/膜生物發(fā)生；22.細(xì)胞運(yùn)動(dòng)；23.細(xì)胞周期控制，細(xì)胞分裂，染色體區(qū)分；24.碳水化合物運(yùn)輸和代謝；25.氨基酸運(yùn)輸與代謝。圖4 女貞果實(shí)unigene的COG分類

注：1.運(yùn)輸與代謝；2.膜轉(zhuǎn)運(yùn)；3.信號(hào)轉(zhuǎn)導(dǎo)；4.折疊，分類和降解；5.復(fù)制和修復(fù)；6.轉(zhuǎn)錄；7.翻譯；8.藥物抵抗；9.內(nèi)分泌代謝?。?0.氨基酸代謝；11.其他次生物質(zhì)代謝；12.碳水化合物代謝；13.能量代謝；14.總覽圖；15.糖生物合成和代謝；16.脂質(zhì)代謝；17.輔助因子和維生素代謝；18.其他的氨基酸代謝；19.萜類和聚酮類化合物的代謝；20.核苷酸代謝；21.環(huán)境適應(yīng)。圖5 女貞果實(shí)unigene的KEGG分類

2.4 SSR分析

利用SSR分析軟件，從14 270條unigene中共搜索到15 925個(gè)SSR位點(diǎn)。搜索到的SSR位點(diǎn)類型豐富，單核苷酸至六核苷酸類型均有(見表2)。其中，二核苷酸重復(fù)所占比例最高，達(dá)到了43.95%；比例最低的是四核苷酸重復(fù)，僅為2.37%。在檢測到的SSR中，出現(xiàn)頻率最高的10類基序?yàn)椋篈G/GT(3856個(gè))、AT/AT(183個(gè))、AC/GT(125個(gè))、AAT/ATT(940個(gè))、AAG/CTT(736個(gè))、ATC/ATG(501個(gè))、ACC/GGT(491個(gè))、AGC/CTG(316個(gè))、AGG/CCT(312個(gè))、CCG/CGG(133個(gè))。

表2 女貞果實(shí)unigene中的SSR位點(diǎn)數(shù)量與分布

3 討論

中藥材是人類醫(yī)藥衛(wèi)生領(lǐng)域的一個(gè)瑰寶，如何開發(fā)利用它使其不斷的造福當(dāng)代世人，一直是現(xiàn)今許多中醫(yī)藥科研工作者不斷追求的目標(biāo)。隨著生物醫(yī)藥科研技術(shù)的不斷發(fā)展，許多新技術(shù)新方法應(yīng)用于中藥材研究，使得以上目標(biāo)逐步被實(shí)現(xiàn)，其中基于高通量測序技術(shù)的轉(zhuǎn)錄組學(xué)研究便是一個(gè)典型例子[17-18]。目前，轉(zhuǎn)錄組學(xué)在中藥材研究領(lǐng)域主要集中于藥用植物功能基因組的研究，即通過轉(zhuǎn)錄組學(xué)發(fā)現(xiàn)藥用植物次生代謝產(chǎn)物生物合成關(guān)鍵酶基因、闡明次生代謝途徑及其調(diào)控機(jī)制等，這將為利用功能基因和代謝途徑體外生物合成中藥材中藥用活性成分，甚至是為創(chuàng)造性的研發(fā)新型藥物奠定堅(jiān)實(shí)的基礎(chǔ)[19-20]。本研究正是基于以上描述，采用Illumina/Solexa高通量測序技術(shù)對(duì)女貞果實(shí)這一中藥材進(jìn)行了轉(zhuǎn)錄組測序，所得結(jié)果顯示共有4.93 Gb有效數(shù)據(jù)產(chǎn)生，可組裝拼接得到163 092個(gè)unigene，其中有79 189條unigene無法獲得功能注釋，究其原因：一是因?yàn)榕懭狈蚪M、EST和蛋白質(zhì)序列信息，不能進(jìn)行基因組條件下的比對(duì)分析，故導(dǎo)致較多的unigene無法匹配到合適的同源序列上；二是由于不能注釋的unigene多為小于1000 bp的小片段，因而它們很難與公共數(shù)據(jù)庫中已知的序列獲得良好的比對(duì)；三是小片段的unigene有些可能是短的非編碼序列或者新基因的某些不完整片段，從而導(dǎo)致它們無法被已有數(shù)據(jù)庫所注釋。另外，由測序數(shù)據(jù)組裝拼接成unigene的過程中，N50值(指從組裝最長的unigene依次向下求長度的總加和，當(dāng)累加長度達(dá)到組裝長度的一半時(shí)對(duì)應(yīng)的unigene長度)是一個(gè)判斷組裝程度好壞的重要指標(biāo)，其值越大說明組裝得到的長片段越多，組裝效果越好。本次研究unigene的N50值為1345 bp，大于1000 bp，說明本次序列組裝的質(zhì)量和長度可以滿足轉(zhuǎn)錄組分析的基本要求，這為研究中進(jìn)行以序列為基礎(chǔ)的相關(guān)轉(zhuǎn)錄組分析提供了可靠保障。

轉(zhuǎn)錄組研究的一個(gè)重要目標(biāo)是發(fā)現(xiàn)功能基因，并闡明功能基因編碼產(chǎn)物所處的代謝途徑。要實(shí)現(xiàn)這一目標(biāo)，研究者需利用GO數(shù)據(jù)庫、COG數(shù)據(jù)庫和KEGG數(shù)據(jù)庫對(duì)組裝得到的unigene進(jìn)行分析。本研究中筆者通過GO數(shù)據(jù)庫的分析，首次系統(tǒng)地闡明了女貞果實(shí)功能基因的分類情況，并利用COG數(shù)據(jù)庫對(duì)女貞果實(shí)unigene進(jìn)行了功能注釋，使其從分子水平找尋自己的直系同源體，從而預(yù)測它們的生物學(xué)功能，最后經(jīng)KEGG數(shù)據(jù)庫把一些預(yù)測的功能基因定位在了涉及主要次級(jí)代謝產(chǎn)物黃酮類和萜類合成的途徑中。在此基礎(chǔ)上，為進(jìn)一步驗(yàn)證這些預(yù)測基因的真正功能，筆者將在后續(xù)的研究中對(duì)它們進(jìn)行體外克隆和異源表達(dá)實(shí)驗(yàn)，并在體外和體內(nèi)兩個(gè)方面對(duì)表達(dá)產(chǎn)物進(jìn)行功能活性檢測，從而用實(shí)驗(yàn)結(jié)果檢驗(yàn)轉(zhuǎn)錄組分析結(jié)果的正確與否。

SSR是真核生物基因組非編碼區(qū)中存在的一些簡單重復(fù)序列，在物種進(jìn)化過程中呈現(xiàn)較為保守的特性，已廣泛應(yīng)用于中藥材物種、種質(zhì)資源以及農(nóng)家種鑒定[21]。目前，女貞缺乏可使用的分子標(biāo)記，此次果實(shí)轉(zhuǎn)錄組研究所獲大量數(shù)據(jù)，可為開發(fā)女貞SSR標(biāo)記提供豐富的篩選資源。在今后的研究中，可利用已獲得的15 925個(gè)SSR位點(diǎn)，通過引物設(shè)計(jì)軟件設(shè)計(jì)合理的SSR引物，然后利用SSR引物進(jìn)行擴(kuò)增檢測，從中篩選出擴(kuò)增穩(wěn)定、條帶清晰、多態(tài)性好的引物，為進(jìn)一步開發(fā)女貞的SSR標(biāo)記奠定基礎(chǔ)[22]。

本研究利用Illumina/Solexa Hiseq2000測序平臺(tái)首次對(duì)女貞果實(shí)進(jìn)行了轉(zhuǎn)錄組測序研究，獲得了大量的轉(zhuǎn)錄組數(shù)據(jù)。對(duì)這些數(shù)據(jù)從拼接組裝、功能注釋、代謝途徑和SSR位點(diǎn)查詢4個(gè)方面進(jìn)行分析研究，使我們?cè)谵D(zhuǎn)錄水平對(duì)女貞果實(shí)有了比較詳實(shí)的認(rèn)識(shí)，同時(shí)也為今后開發(fā)利用女貞果實(shí)的藥用次生代謝產(chǎn)物和分子標(biāo)記物，以及為女貞基因組的測序組裝提供了具有參考價(jià)值的基礎(chǔ)數(shù)據(jù)。

[1] 國家藥典委員會(huì).中華人民共和國藥典：一部[S].北京：中國醫(yī)藥科技出版社，2015：45-46.

[2] 尹輝.中藥女貞子化學(xué)成分研究綜述[J].九江學(xué)院學(xué)報(bào)(自然科學(xué)版)，2015，30(1)：74-75.

[3] 劉亭亭，王萌.女貞子化學(xué)成分與藥理作用研究進(jìn)展[J].中國實(shí)驗(yàn)方劑學(xué)雜志，2014，20(14)：228-234.

[4] Nagalakshmi U，Waern K，Snyder M.RNA-Seq：A Method for Comprehensive Transcriptome Analysis[J].Curr Protoc Mol Biol，2010(suppl 89):Unit 4.11.1-13.

[5] Hoeijmakers W A，Bártfai R，Stunnenberg H G.Transcriptome analysis using RNA-Seq[J].Methods Mol Biol，2013，923：221-239.

[6] 田英芳，張曉政，周錦龍.轉(zhuǎn)錄組學(xué)研究進(jìn)展及應(yīng)用[J].中學(xué)生物教學(xué)，2013(12)：29-31.

[7] Martin J A，Wang Z.Next-generation transcriptome assembly[J].Nature Reviews Genetics，2011，12(10)：671-682.

[8] Chen C，Huang H，Wu C H.Protein Bioinformatics Databases and Resources[J].Methods Mol Biol，2011，694:3-39.

[9] Collins L J，Biggs P J，Voelckel C，et al.An approach to transcriptome analysis of non-model organisms using short-read sequences[J].Genome Informatics，2008，21：3-14.

[10] Martin M.Cutadapt removes adapter sequences from high-throughput sequencing reads[J].Embnet Journal，2011，17(1)：10-12.

[11] Grabherr M G，Haas B J，Yassour M，et al.Full-length transcriptome assembly from RNA-Seq data without a reference genome[J].Nat Biotechnol，2011，29：644-652.

[12] Picardi E，Mignone F，Pesole G.EasyCluster：a fast and efficient gene-oriented clustering tool for large-scale transcriptome data[J].BMC Bioinformatics，2009，10(6)：1-12.

[13] Bevilacqua V，Pietroleonardo N，Giannino E，et al.EasyCluster2：an improved tool for clustering and assembling long transcriptome reads[J].BMC Bioinformatics，2014，15(15)：1-10.

[14] Conesa A，G?tz S，García-Gómez J M，et al.Blast2GO：A universal tool for annotation，visualization and analysis in functional genomics research[J].Bioinformatics，2005，21：3674-3676.

[15] Ye J，F(xiàn)ang L，Zheng H，et al.WEGO：A web tool for plotting GO annotations[J].Nucleic Acids Res，2006，34：293-297.

[16] Wang H，Jiang J，Chen S，et al.Next-Generation Sequencing of the Chrysanthemum nankingense(Asteraceae)Transcriptome Permits Large-Scale Unigene Assembly and SSR Marker Discovery[J].PloS One，2013，8(4)：e62293.

[17] 張曉萌，李健春，王瓊，等.轉(zhuǎn)錄組測序技術(shù)在中醫(yī)藥領(lǐng)域的應(yīng)用[J].中國現(xiàn)代中藥，2016，18(8)：1084-1087.

[18] 陳惠，辛麗麗，龔婕寧.基于全轉(zhuǎn)錄組測序技術(shù)的轉(zhuǎn)錄組學(xué)在中醫(yī)藥領(lǐng)域的應(yīng)用前景分析[J].環(huán)球中醫(yī)藥，2013，6(10)：759-763.

[19] 王堯龍，黃璐琦，袁媛，等.藥用植物轉(zhuǎn)錄組研究進(jìn)展[J].中國中藥雜志，2015，40(11)：2055-2061.

[20] 吳瓊，孫超，陳士林，等.轉(zhuǎn)錄組學(xué)在藥用植物研究中的應(yīng)用[J].世界科學(xué)技術(shù)—中醫(yī)藥現(xiàn)代化，2010，12(3)：457-462.

[21] 周駿輝，袁媛，黃璐琦.SSR標(biāo)記在中藥材分子身份證體系構(gòu)建中的應(yīng)用[J].中國現(xiàn)代中藥，2016，18(10)：1233-1236.

[22] 劉峰，王運(yùn)生，田雪亮，等.辣椒轉(zhuǎn)錄組SSR挖掘及其多態(tài)性分析[J].園藝學(xué)報(bào)，2012，39(1)：168-174.