溪黃草轉(zhuǎn)錄組測序及生物信息學(xué)分析

王繼華，蔡時可，梅瑜，李漢章，楊少海

（1.廣東省農(nóng)作遺傳改良重點實驗室，廣東省農(nóng)業(yè)科學(xué)院作物研究所，廣東廣州 510640；2.廣東建鵬南藥種養(yǎng)有限公司，廣東廣州 510070）

溪黃草Rabdosia serra（Maxim.）Hara 為唇形科香茶菜屬多年生草本植物，喜陰涼濕潤環(huán)境，俗稱溪溝草、黃汁草等，多產(chǎn)于華南地區(qū)的廣東、廣西及臺灣[1-2]。全草均可入藥，味苦、性寒，具有清熱利濕、退黃、涼血散瘀的功效，用于治療濕熱瀉痢、跌打瘀腫、急性黃疸型肝炎、急性膽囊炎、口腔炎、腸炎等病證[3-5]。藥理學(xué)研究表明，溪黃草具有保肝、護肝、抗炎、抗癌作用，是多種中成藥和保健品的主要成分[3，6-7]。溪黃草是廣東省連州市的特產(chǎn)，連州溪黃草在2007 年獲得國家地理標(biāo)志產(chǎn)品（國家質(zhì)量監(jiān)督檢驗檢疫總局，2007年第220號公告）。

溪黃草來源復(fù)雜，其基原植物為唇形科溪黃草、線紋香茶菜、以及線紋香茶菜的變種細花線紋香茶菜和狹基線紋香茶菜[8]。資源鑒定主要根據(jù)根、莖、葉、花、果實等形態(tài)鑒定及味道等進行區(qū)分[9-11]。溪黃草中主要含有萜類、黃酮類、多糖和酚酸類等化學(xué)成分，各種有效成分還具有協(xié)同增效作用[12]。不同的基源植物之間化學(xué)成分相近，但部分藥效成分的含量相差甚遠[13]。溪黃草基因組數(shù)據(jù)仍不完整，不能對有效成分的代謝途徑進行深入的挖掘。目前，高通量測序平臺的RNA-seq技術(shù)在解析藥用成分代謝通路、挖掘關(guān)鍵基因、開發(fā)分子遺傳標(biāo)記等方面得到廣泛應(yīng)用[14-16]。因此，本研究開展溪黃草的轉(zhuǎn)錄組測序及生物信息學(xué)分析，以期為解析其藥用物質(zhì)合成代謝通路、挖掘關(guān)鍵調(diào)控基因以及開發(fā)分子標(biāo)記奠定基礎(chǔ)。現(xiàn)將研究結(jié)果報道如下。

1 材料與方法

1.1實驗材料實驗材料為廣東省農(nóng)業(yè)科學(xué)院作物研究所南藥資源圃栽培的溪黃草。于2018 年12 月份取樣，采集健壯植株的葉片，迅速用錫箔紙包裹并浸入液氮處理，隨后置于超-80 ℃冰箱保存。

1.2溪黃草RNA的提取采用生工生物工程（上海）股份有限公司的總RNA 提取試劑盒（B511311-0025）提取溪黃草總RNA，然后通過1%電泳凝膠檢測提取RNA的完整性。應(yīng)用Invitrogen Qubit?2.0熒光計及試劑盒（Fluorometer Life Tech Invitrogen，Q32886）對總RNA進行定量。

1.3轉(zhuǎn)錄組測序與拼接組裝委托生工生物工程（上海）股份有限公司采用Illumina HiSeq2500 的高通量測序平臺進行轉(zhuǎn)錄組測序。測序得到的原始數(shù)據(jù)，通過FastQC軟件進行質(zhì)量評估和Trimmomatic進行質(zhì)量剪切，過濾掉接頭、低質(zhì)量的序列（reads ＜ 35nt）、帶N 堿基的序列、低質(zhì)堿基（Q值＜20）得到高質(zhì)量的clean data[17]。應(yīng)用Trinily 軟件對clean data 進行de novo拼接組裝，再采用RSeQC軟件去除轉(zhuǎn)錄本中的冗余序列，得到非冗余通用基因（universal gene，unigene）[18]。

1.4基因功能注釋采用基于局部比對算法的搜索工具（Basic Local Alignment Search Tool，BLAST）將組裝的unigene 與保守域數(shù)據(jù)庫（Conserved Domain Database，CDD）、真核生物蛋白質(zhì)同源簇數(shù)據(jù)庫/蛋白相鄰的聚類（eu Karyotic ortholog groups/clusters of orthologous groups，KOG/COG）、非冗余（Non-redundant，NR）、核酸序列數(shù)據(jù)庫（Nucleotide Sequence Database，NT）、蛋白結(jié)構(gòu)域預(yù)測（Protein Families Database of Alignments and Hidden Markov Models，PFAM）、Swissprot、TrEMBL 等多個數(shù)據(jù)庫比對得到功能注釋信息。使用Transdecoder進行編碼序列（Coding Sequence，CDS）預(yù)測。根據(jù)unigene 與Swissprot、TrEMBL 的注釋結(jié)果得到基因本體論（Gene Ontology，GO）功能注釋信息，利用京都基因與基因組百科全書（Kyoto Encyclopedia of Genes and Genomes，KEGG）自動注釋服務(wù)器（KEGG Automatic Annotation Server，KAAS）得到KEGG注釋信息。

1.5基因結(jié)構(gòu)分析使用微衛(wèi)星識別工具（Microsatellite Identification Tool，MISA）軟件鑒定溪黃草unigene上存在的SSR 位點，并利用Primer 3軟件（http：//primer3.sourceforge.net/releases.php）設(shè)計SSR引物[19]。

2 結(jié)果與分析

2.1轉(zhuǎn)錄組測序與de novo組裝溪黃草cDNA 文庫的構(gòu)建由生工生物工程（上海）股份有限公司完成。應(yīng)用Illumina Hiseq 2500測序平臺測序，共獲得61 944 850條raw reads，總堿基數(shù)為9 291 727 500 bp。使用Trimmomatic 對原始測序數(shù)據(jù)進行處理，去掉含有帶接頭、低質(zhì)量的序列，共得到60 234 786 clean reads，總堿基數(shù)目為8 704 764 735 bp，GC含量為51.21%，Q30 bases ratio達到95.64%，表明文庫構(gòu)件質(zhì)量良好，測序得到的數(shù)據(jù)準(zhǔn)確可靠。使用Trinity 將clean reads 進行de novo組裝成轉(zhuǎn)錄本，共得到86 204 條轉(zhuǎn)錄本，平均長度為1 251.93 bp，N50為1 894 bp，序列長度大于500 bp的有59 403條，占總序列數(shù)目的68.91%。見表1。對Trinity拼裝得到的轉(zhuǎn)錄本去冗余，共獲得37 418 條unigene，平均長度為1 054.1 bp，N50為1 840 bp，其中20 426條序列長度大于500 bp，占總序列數(shù)目的54.59%，序列長度在1 000 bp以上的有14 226條，占總數(shù)的38.02%見。見表1、圖1。

表1 溪黃草轉(zhuǎn)錄組測序的結(jié)果Table 1 Summary of assembled transcripts and unigenes of the Rabdosia serra（Maxim.）Hara transcriptome

圖1 溪黃草unigene序列長度分布Figure 1 Length distribution of assembled unigenes of Rabdosia serra（Maxim.）Hara

2.2 Unigene功能注釋將溪黃草組裝后的unigene序 列 與 CDD、 KOG、 COG、 NR、 NT、 PFAM、Swissprot、TrEMBL 等多個數(shù)據(jù)庫進行比對，共有23 978 條（64.08%）unigene 在至少一個數(shù)據(jù)庫中獲得功能注釋，2 429 條（6.49%）unigene 在所有數(shù)據(jù)庫種均能獲得注釋，見表2。其中NR 數(shù)據(jù)庫中注釋到的unigene 數(shù)目（23 623 條）最多，占總unigene的63.13%，其次為TrEMBL 數(shù)據(jù)庫，注釋到的unigene 為 23 393 條，占 62.52%。KEGG 數(shù)據(jù)庫中注釋到的unigene 數(shù)目（3 674 條，9.82%）最少。尚有13 440條unigene沒有得到有效注釋，占35.92%。見表2。通過與NR 庫的比對，獲得溪黃草unigene序列與近緣種屬的近似情況并獲得同源序列的功能信息，共有23 616條unigene獲得注釋。匹配較多的物種主要有Sesamum indicum、Erythranthe guttata、Dorcoceras hygrometricum和Salvia miltiorrhiza，分別占58.88%、20.39%、3.16%和1.97%，見圖2。表明溪黃草與唇形目（管狀花目）Sesamum indicum的序列相似度最高。

表2 Unigene 的功能注釋Table 2 Functional annotation of assembled unigenes

2.3 KOG功能注釋將組裝得到的溪黃草unigene與KOG數(shù)據(jù)庫比對，共有11 473條unigene獲得注釋，按其功能共分為26 類。其中：信號轉(zhuǎn)導(dǎo)機制注釋到的unigene 數(shù)目最多，共有1 974 條，占15.35%；一般功能預(yù)測、翻譯后修飾，蛋白質(zhì)轉(zhuǎn)運，分子伴侶伴侶、轉(zhuǎn)錄和碳水化合物的轉(zhuǎn)運和代謝注釋到的unigene 占比分別為11.42%、10.26%、5.50%和5.21%；僅有34條和4條unigene注釋對應(yīng)到核結(jié)構(gòu)和細胞運動性，分別占比0.26%和0.03%。見圖3。此外，676 條unigene 注釋到未知功能。由此可見，溪黃草unigene 涉及信號轉(zhuǎn)導(dǎo)功能最多，可為今后溪黃草代謝物質(zhì)的調(diào)控研究提供寶貴資源。

圖2 NR數(shù)據(jù)庫的同源物種分類Figure 2 Species distribution of the top BLAST hits against the NR database for the assembled unigenes

2.4 GO功能注釋GO 數(shù)據(jù)庫是全面描述生物體中基因及其產(chǎn)物屬性的分類系統(tǒng)，主要分為生物過 程（biological process）、 細 胞 組 分（cellular component）及分子功能（molecular funtion）三大類，見圖4。根據(jù)GO 數(shù)據(jù)庫注釋結(jié)果，共有20 222 條溪黃草unigene注釋成功。這些unigene總共被劃分為65 個功能分類，生物過程注釋到的unigene 最多，共有42 872 個，占35.95%，其中注釋較多的功能分別為細胞過程（11 757個，7.88%）、代謝過程（9 935個，6.66%）、刺激響應(yīng)（5 455個，3.66%）和生物調(diào)節(jié)（4 952 個，3.32%）。注釋到分子功能的unigene 數(shù)目為38 425個，占32.22%，其中注釋較多的功能分別為結(jié)合（11 625 個，7.79%）、催化活性（9 633，6.45%）和轉(zhuǎn)運活性（1 276個，0.86%）。37 943 個（31.82%）unigene 注釋到細胞組分，其中注釋較多的功能分別為細胞（14 297 個，9.58%）、細胞組分（14 265 個，9.56%）、細胞器（10 568 個，7.08%）、膜（8 331 個，5.58%）和膜組分（6 481 個，4.34%）。

圖3 溪黃草unigene的KOG功能分類Figure 3 KOG function classification of Rabdosia serra（Maxim.）Hara unigenes

圖4 溪黃草unigene的GO功能分類Figure 4 GO classification of Rabdosia serra（Maxim.）Hara unigenes

2.5 KEGG功能注釋根據(jù)KEGG 數(shù)據(jù)庫注釋結(jié)果，主要分為五大類功能，包括代謝（2 529 個，41.67%）、遺傳信息過程（1 298 個，21.39%）、細胞過程（638 個，10.51%）、環(huán)境信息過程（623 個，10.27%）和有機系統(tǒng)（981 個，16.16%）。根據(jù)unigene 參與的代謝過程，進一步劃分為32個功能分類，涉及289 個代謝通路，具體見圖5。在代謝中，與氨基酸代謝和脂質(zhì)代謝相關(guān)的unigene 分別為305 個（5.03%）和275 個（4.53%）。在有機系統(tǒng)中，涉及環(huán)境適應(yīng)的unigene 有162 個，占2.67%。在遺傳信息過程中注釋到unigene 最多的代謝過程為轉(zhuǎn)錄（212 個，3.49%），而在細胞過程中，與轉(zhuǎn)運和分解代謝的unigene 最多（310個，5.11%）。信號轉(zhuǎn)導(dǎo)在植物發(fā)育過程以及對外界刺激的響應(yīng)過程起重要作用，本研究共檢測到590 條unigene 涉及信號轉(zhuǎn)導(dǎo)，可為下一步鑒定溪黃草的信號因子提供研究基礎(chǔ)。溪黃草主要含萜類、黃酮類、酚類、氨基酸等化學(xué)成分，其中二萜類化合物非常豐富。本研究鑒定到19個unigene涉及黃酮類物質(zhì)的生物合成，12個unigene 涉及倍半萜類化合物和三萜類化合物生物合成，40個unigene 與萜類物質(zhì)骨架生物合成相關(guān)，13個unigene 參與雙萜的生物合成，51個unigene 涉及苯丙素的生物合成。這些與次生代謝相關(guān)unigene 的鑒定結(jié)果為進一步解析溪黃草藥用成分物質(zhì)的生物合成提供了可能。

圖5 溪黃草unigene的KEGG功能分類Figure 5 KEGG classification of Rabdosia serra（Maxim.）Hara unigenes

2.6 SSR分析見圖6。采用MISA 對組裝的unigene 進行SSR 檢測，并對SSR 的類型和密度進行統(tǒng)計。結(jié)果表明，在7 809條unigene中共鑒定到9 489 個 SSR 位點。其中，921 條 unigene 中檢測到944 個（9.95%）復(fù)雜重復(fù)類型的SSR 位點。最豐富的重復(fù)類型是雙堿基重復(fù)，共檢測到4 208 個位點，占44.35%；其次為單堿基重復(fù)（2 524 個，26.60%）、三堿基重復(fù)（1 685 個，17.76%）、六堿基重復(fù)（70 個，0.74%）和四堿基重復(fù)（42 個，0.44%）；最少的為五堿基重復(fù)，僅檢測到16 個位點，占0.17%。見圖6-A。在SSR 分布的密度上，雙堿基重復(fù)最高，達到142.2個/Mbp，五堿基重復(fù)最低，僅為0.48 個/Mbp。見圖6-B。基于9 489 個SSR 位點，使用Primer 3.0 設(shè)計引物，為進一步開發(fā)溪黃草的遺傳標(biāo)記和近緣種屬的遺傳圖譜提供了基礎(chǔ)數(shù)據(jù)。

圖6 溪黃草unigene的SSR位點分析Figure 6 The analysis of SSR sites of Rabdosia serra unigenes

3 討論

近年來，隨著高通量測序技術(shù)的發(fā)展，憑借其產(chǎn)出數(shù)據(jù)量大、成本低、不需要參考基因組等優(yōu)勢，在新基因發(fā)掘、功能基因鑒定和分子標(biāo)記開發(fā)上的應(yīng)用越來越廣泛。本研究基于Illumina HiSeq2500的高通量測序平臺對溪黃草葉片進行轉(zhuǎn)錄組測序，共獲得60 234 786 條clean reads，總堿基數(shù)目為8 704 764 735 bp，GC 含量為51.21%，Q30 bases ratio 達到95.64%，表明文庫構(gòu)件質(zhì)量良好。N50 是評價組裝序列完整性的重要指標(biāo)。通過de nove組裝，本研究共獲得37 418 條unigene，平均長度為1 054.1 bp，N50為1 840 bp，與已經(jīng)構(gòu)建云南松（1 818 bp）、香榧（1 702 bp）的N50 長度相接近，但較中藥黃芩（797.64 bp）的N50 長，表明溪黃草轉(zhuǎn)錄組序列組裝質(zhì)量較高[20-22]。

溪黃草組裝后的unigene 序列與多個公共數(shù)據(jù)庫比對，共有 23 978 條（64.08%）unigene 在至少1 個數(shù)據(jù)庫中獲得功能注釋，2 429 條（6.49%）unigene在所有數(shù)據(jù)庫中均能獲得注釋，但仍有約35.92%的unigene 沒有獲得注釋信息，可能是由于部分序列組裝長度過短，缺少保守的核心序列以及溪黃草基因組信息匱乏暫時無法獲得準(zhǔn)確的功能注釋。在KOG 數(shù)據(jù)庫和KEGG 數(shù)據(jù)庫中分別鑒定到1 974條和590條unigene涉及信號轉(zhuǎn)導(dǎo)，可為下一步解析溪黃草響應(yīng)外界刺激以及代謝物質(zhì)的調(diào)控提供靶標(biāo)基因。同時，還鑒定到19個unigene涉及黃酮類物質(zhì)的合成，51個unigene 涉及苯丙素的生物合成，以及大量與萜類物質(zhì)合成相關(guān)的unigene，這為明晰溪黃草黃酮類和萜類物質(zhì)的合成途徑和代謝網(wǎng)絡(luò)提供了數(shù)據(jù)支撐。

SSR標(biāo)記操作簡單、重復(fù)性好，與傳統(tǒng)方法相比，高通量能夠挖掘出大量的SSR 位點。甘草、鐵皮石斛、枸杞等中藥材均利用轉(zhuǎn)錄組測序開發(fā)了SSR 標(biāo)記用于分子輔助育種和遺傳圖譜的構(gòu)建[23-25]。本研究應(yīng)用MISA 軟件共檢測到9 489 個SSR 位點，最豐富的重復(fù)類型是雙堿基重復(fù)（4 208 個位點），占44.35%，達到142.2 個/Mbp，該研究結(jié)果與中藥黃芩類似[20]。本研究對溪黃草轉(zhuǎn)錄組進行了初步的探究，彌補了溪黃草基因組信息的不足，為解析次級代謝物質(zhì)合成通路及分子生物學(xué)方面的研究打下了基礎(chǔ)。溪黃草SSR 位點的發(fā)掘，可為溪黃草分子標(biāo)記的開發(fā)、遺傳多樣性分析、遺傳圖譜構(gòu)建等奠定理論基礎(chǔ)，為利用分子手段鑒定和區(qū)分溪黃草及其基源植物提供了依據(jù)。