濒危珍稀植物半枫荷的转录组分析

叶兴状 刘 丹 罗佳佳 范辉华 张国防 刘 宝* 陈世品

(1.福建农林大学林学院，福州 350002； 2.福建省林业科学研究院，福州 350012)

转录组(Transcriptome)能够在缺乏基因资源的条件下，运用转录组测序技术获得物种的生长和代谢规律，并揭示其生物学特性与基因内在关联，同时可获得物种绝大多数的转录产物[1]。转录组测序不仅能够发掘得到新的代谢通路及功能基因还是开展EST-SSR分子标记的前提，而EST-SSR分子标记是目前最重要的群体遗传多样性研究方法之一。虽然第三代高通量测序已问世，但第二代高通量测序相较第三代测序具有成本低、操作简便、运行时间短等特点，是一种性价比高的基因序列研究手段。近年来已有多种药用植物完成转录组测序，如郑少华[2]等对淫羊藿(Epimediumbrevicornu)进行转录组分析，找到类黄酮素代谢通路；Roberto A Barrero等[3]对狼毒(Stellerachamaejasme)进行转录组分析，提出了前列腺素合成上游途径的新见解，为批量产生这种化合物奠定理论基础，将可用于HIV研究或治疗患者；李响等[4]对冬虫夏草(Chinese Caterpillar Fungus)进行转录组分析，开发了一种用于合成主要药用化合物(如虫草素)的模型；Xu Guo等[5]对铁皮石斛(Dendrobiumofficinale)茎进行转录组分析，揭示了代表参与生物碱骨架合成的25个基因及69个独特序列。

半枫荷(Semiliquidambarcathayensis)属金缕梅科(Hamamelidaceae)半枫荷属(Semiliquidambar)是我国特有的濒危国家二级保护植物[6]，半枫荷在《中国生物多样性红色名录—高等植物卷》中被评为VU级[7]，具有枫香属(Liquidambar)和蕈树属(Altingia)两属的综合性状特征，异型叶是其重要形态特征，叶子兼具有枫香(Liquidambarformosana)和阿丁枫(Altingiachinensis)特征，半枫荷嫩叶部分紫红，极为美观，是南方为数不多的变色树种，具有良好的园林绿化应用前景。此外其材质优良、树干通直、旋刨性良好、材质优良，可作旋刨制品及家具，不仅是珍贵用材树种[8～9]，而且药用价值极高，其根、枝、叶，都可入药，有祛风、除湿、活血、止痛、消肿等功效，可治疗腰肌劳损、跌打淤积、风湿性关节炎，且外伤常用其叶止血[10～12]。目前半枫荷的研究处于起步阶段，主要涉及种群生态特征和药理研究等领域，鲜有涉及分子生物学领域，胡君等[13]在四川省发现了半枫荷新分布，并调查半枫荷种群生境特征；廖娜等[14]发现半枫荷中的多酚具有抗菌和显著抗氧化能力；周光雄等[15]首次从金缕半枫荷中分离得齐墩果酸、3-羰基齐墩果酸等9种有抗炎药效的化合物，为半枫荷药理研究提供基础数据；此外孙静[16]等还发现半枫荷根的活性成分齐墩果酸对病毒性肝炎的抗原具有抑制活性的作用。由于半枫荷用途广、经济价值高，再加上天然更新困难，对生境要求苛刻，同时因根可入药、树根常被药农挖断，人为破坏严重，目前仅在贵州、湖南、福建、广西、江西、广东等省份零星分布[11～12,17]，已濒临灭绝，亟需保护。

鉴于此，本研究选择福建省沙县、顺昌野生半枫荷做为研究材料，借助高通量测序平台Illumina HiSeq 2500首次对半枫荷叶进行转录组测序，建立半枫荷基因文库，然后利用生物信息学的方法对半枫荷进行功能注释，代谢途径、CDS预测及SNP检测等分析，以期初步了解半枫荷基因表达、分布情况，如参与药效合成的酶的基因分布与表达情况，为半枫荷遗传多样性研究提供基础数据，同时也为今后半枫荷功能基因的研究奠定基础。

1 材料与方法

1.1 实验材料

本实验选取半枫荷叶为实验材料。于2018年3月份到福建省顺昌县、沙县各选取1棵半枫荷，采集一年生叶子放入预先准备的冰盒中带回实验室，将叶片擦干去除污渍，后用锡箔纸将叶片包好做好标记，并用液氮处理10～20分钟，再用超低温零下80℃冰箱保存，用于半枫荷引物开发及遗传多样性分析实验。

1.2 方法

1.2.1 半枫荷总RNA的提取

将采集的2棵半枫荷叶样品混合在一起，然后用生工生物(上海)有限公司的EZ-10DNA away RNA Mini-prep Kit试剂盒提取总RNA，使用杭州奥盛仪器有限公司的Nano-200微量分光光度计测定RNA的浓度和质量，如果OD260NM/OD280NM的范围为1.8～2.2，28S/18S>2，RNA完整值(RIN)>8.5，OD260NM/OD230NM>1.8，则认为提取的RNA符合建库要求。

1.2.2 半枫荷测序与拼接组装

采用Illumina Hi Seq2500的高通量测序技术对半枫荷的转录组进行测序，测序得到的原始数据，里面含有带接头的、低质量的序列。为了保证信息分析质量，对原始数据过滤，得到Clean数据，参数如下：

①去除带N碱基的序列；

②去除reads中的接头序列；

③从reads 3′到5′方向开始去除低质碱基(Q值<20)；

④使用滑窗法去除reads尾部质量值在20以下的碱基(窗口大小为5 bp)；

⑤去除reads长度小于35nt的reads本身及其配对reads。

⑥从reads 5′到3′方向开始去除低质碱基(Q值<20)；

将reads回帖到全部转录本，并归一化后得到reads在转录本上的平均分布，发现有少许偏移，但未发生明显偏移，符合建库要求。然后将Clean数据利用Trinity[18]软件(参数min_kmer_cov 2，其余默认)进行拼接组装。建库过程中的PCR扩增会产生冗余序列(duplicate reads)，指的是碱基排列完全相同的序列，如果建库扩增时不均匀，可能会导致部分基因冗余序列过高，造成冗余序列分布的异常，使用RSeQC软件去除冗余序列，得到Unigenes。

1.2.3 基因功能注释

运用NCBI Blast+将转录本与CDD、KOG、NR、NT、GO、Swissprot、PFAM、KEGG、TrEMBL等9个数据库进行比对，得到其功能注释信息。根据转录本与Swissprot、TrEMBL的注释结果得到GO功能注释信息。利用KAAS[19]与基因和蛋白质数据库进行比对，得到转录本KEGG注释信息。利用以上数据库的分析结果，筛选Unigenes最优比对片段作为该Unigenes的CDS,使用getorf[20]检测Unigenes的开放阅读框(Open Reading Frame，ORF)，采用hmmsearch[21]将ORF比对到转录因子蛋白结构域(数据来源于PlantfTFDB)，然后根据PlantTFDB描述的转录因子家族特征对Unigenes进行能力鉴定。最后使用BCFtools[22]根据Mapping结果找出可能的单核苷酸多态性(Single Nucleotide Polymorphsims，SNP)位点，进行SNP分析。

2 结果与分析

2.1 半枫荷转录本组装与CDS预测分析

通过使用软件将clean数据Denovo组装成转录本共得到127 249个Transcripts。对半枫荷转录本长度进行分析(图1A,表1)，经组装、去冗长等处理得到77 629个Unigenes，200～500 bp范围内的Unigenes最多，为36 630条，其次是300～400 bp、400～500 bp和大于2 000 bp。最短长度为201 bp，最长10 975 bp，N50长度1 105 bp，总长度48 512 368 bp，Unigene平均长度为624.93 bp，其中长度分布在200～500 bp的Unigenes占比67.88%，分布在500～1 000 bp的有11 296个，占比14.55%，长度大于1 000 bp的有13 635个，占比17.57%。同时，为了更好的发掘半枫荷功能基因和便于遗传指纹图谱的构建，本研究对半枫荷编码序列(Coding Sequence，CDS)进行预测，获取57 671个Unigenes作为CDS。CDS序列1 000 nt以上的有8 050个，仅占13.96%，100～500 nt的序列有40 608个，占70.41%，其中200～300 nt的序列有20 615个，占总CDS序列的37.75%，500～1 000 nt的序列仅占15.63%。将图1B与图1A对比，序列分布更加集中于100～600 nt。

图1 半枫荷Unigene的长度和CDS长度的分布图Fig.1 Distribution of unigene lengths and CDS length for S.cathayensis

项目Item数量No.≥500bp≥1000bpN50N90最长值Max len最小值Min Len总长度Total len平均长度Average lenTranscript127249572803214413052851079520197669341767.55Unigene77629249311363511052451079520148512368624.93

2.2 功能注释及COG分类

将denovo组装得到的77 629个基因序列分别与CDD、KOG、NR、NT、GO、Swissprot、PFAM、KEGG、TrEMBL这9个数据库进行比对，成功注释到9个数据库中(表2)；在9个数据库中至少1个数据库注释成功的Unigenes数为45 293，占总Unigenes数的58.35%；共有2 692个Unigenes与9个数据库都能匹配成功，占总Unigenes数的3.47%(表2)；而只能注释到NR、KEGG、COG、Swissprot蛋白库的Unigenes分别只有6 971、0、120、3 362个(图2)；此外，还有32 336(41.65%)个Unigenes序列未能获得功能注释信息，这有可能是测序技术不同引起的。

表2 Unigene注释结果统计表

图2 半枫荷Unigene的功能注释维恩图Fig.2 Vean diagram functional annotation distribution of unigenes of transcriptome for S.cathayensis

2.3 KOG功能注释

将组装得到的Unigenes比对到KOG蛋白质库中，半枫荷KOG功能注释获信息结果显示细分为25个子类、共有22 802个Unigenes被注释(占Unigenes总数的29.37%)，其各个子类的基因表达丰度差异较大，共包括25 253个功能注释信息，基本涵盖了半枫荷大部分的生命活动。其中，信号传递机制，翻译后修饰、蛋白质折叠和分子伴侣，一般功能基因最多，分别有2 858(11.32%)、2 699(10.69%)和2 552(10.11%)。而核结构和细胞活性较少，分别仅有55(0.22%)和6(0.02%)个。另有971(占3.85%)个功能注释信息未明确其生物学功能(图3)。

图3 半枫荷转录组Unigenes的KOG功能分布图 A.RNA加工和修饰；B.染色体结构和动力学；C.能源产生与转化；D.细胞周期调控，细胞分裂，染色体分离；E.氨基酸转运和代谢；F.核酸转运和代谢；G.碳水化合物转运和代谢；H.辅酶转运和代谢；I.脂类转运和代谢；J.翻译，核糖体结构和生物发生；K.转录；L.复制，重组和修饰；M.细胞壁/细胞膜生物发生；N.细胞活性；O.翻译后修饰,蛋白翻转,伴侣；P.无机离子转运和代谢；Q.次生代谢物合成,转运和代谢；R.只有一般功能预测；S.未知功能；T.信号传递机制；U.细胞间运输,分泌物和囊泡运动；V.防御机制；W.细胞外结构；Y.核结构；Z.细胞骨架Fig.3 KOG functional annotation distribution of unigenes of transcriptome for S.cathayensis A. RNA processing and modification; B. Chromatin structure and dynamics; C. Energy production and conversion; D. Cell cycle control,cell division,chromosome partitioning; E. Amino acid transport and metabolism; F. Nucleotide transport and metabolism; G. Carbohydrate transport and metabolism; H. Coenzyme transport and metabolism; I. Lipid transport and metabolism; J. Translation,ribosomal structure and biogenesis; K. Transcription; L. Replication,recombination and repair; M. Cell wall/membrane/envelope biogenesis; N. Cell motility; O. Posttranslational modification,protein turnover,chaperones; P. Inorganic ion transport and metabolism; Q. Secondary metabolites biosynthesis,transport and catabolism; R. General function prediction only; S. Function unknown; T. Signal transduction mechanisms; U:Intracellular trafficking,secretion,and vesicular transport; V. Defense mechanisms; W. Extracellular structures; Y. Nuclear structure; Z. Cytoskeleton

图4 半枫荷转录组Unigenes的NR注释物种分布图Fig.4 NR annotated species distribution of Unigenes of transcriptome for S.cathayensis

2.4 NR功能注释

将获得的半枫荷Unigenes比对到NR库中(E≤1e-5)，NR功能注释获信息结果显示有39 076个Unigenes被注释(占50.34%)，且与其它物种已知基因序列具有不同程度的匹配性。如注释物种分布图所示(图4)：注释为葡萄(Vitisvinifera)相关基因的序列最多，达10 111条，占NR库中被注释Unigenes的25.89%，其次，序列同源性大于2%的近缘物种还有大麦亚种(Hordeumvulgaresubsp. Vulgare)、荷花(Nelumbonucifera)、可可(Theobromacacao)、枣(Ziziphusjujuba)、麻风树(Jatrophacurcas)、蓖麻(Ricinuscommunis)、甜橙(Citrussinensis)、毛果杨(Populustrichocarpa)和梅(Prunusmume)，其余46.07%的注释Unigenes分布于其他648个物种中。由于缺乏半枫荷基因组信息和技术局限性，仍有29(0.074%)个Unigenes在NR库中未能获得注释。

2.5 GO功能注释及分类

根据转录本与Swissprot、TrEMBL的注释结果得到GO功能注释信息。综合描述半枫荷叶片中相关基因的生物学特征，从宏观上了解其表达基因的功能分布情况，对获得相应GO条目的Unigenes进行统计分析(图5)，共有39 293个Unigenes获得317 475个GO注释，平均每条8.08个。把获得GO注释的Unigenes划分为细胞组分、分子功能和生物学过程三个类型，分别各获得133 668(42.10%)、53 236(16.77%)和130 571(41.13%)个注释，可细分为65个子类，其中，细胞、细胞组分、结合、细胞器、催化活性、细胞过程、代谢过程和单一有机体过程等获得注释较多。

2.6 KEGG代谢通路分析

将半枫荷的Unigenes序列映射到KEGG代谢库(E≤1e-5)，根据注释信息对其可能参与或涉及的代谢途径进行统计分析。结果表明，4 029条Unigenes获得注释(5.19%)；其参与的代谢通路可归为5大类别、33个子类。由图6可知，5种代谢通路大类中，代谢相关的通路获得2 786个注释，占41.99%，其次是遗传信息处理和生物系统分别获得1 562(23.55%)、1 036(15.62%)，再次是细胞过程与环境信息处理分别获得650(9.80%)、600(9.04%)。进一步细分为33子类代谢途径，其中，翻译获得注释最多，为722(17.92%)，其次为信号传导、折叠排序和退化、碳水化合物代谢。共有8 923条次Unigenes归入286条代谢途径，按数量排列，靠前的有15个代谢通路(表3)。核糖体代谢途径注释到的Unigenes数量最多，有389(4.36%)条Unigenes，其次为碳代谢(221，2.48%)、氨基酸的生物合成(193，2.16%)和内质网蛋白处理(171，1.92%)等途径。

在半枫荷的转录组数据中，有64条Unigenes映射到苯丙素生物合成(ko00940)；有65条Unigenes映射到萜类物质生物合成(ko00904，ko00902，ko00909，ko00900)；有27条Unigenes映射到黄酮类物质生物合成(ko00941，ko00942，ko00944，ko00965)；有2条Unigenes映射到硫甙生物合成(ko00966)，可合成糖苷类物质；有18条Unigenes映射到了合成生物碱类(ko00232、ko00950)通路；有1条unigene映射到了黄酮和黄酮醇的生物合成(ko00944)通路，有1条Unigene映射到了异黄酮的生物合成(ko00943)通路(表4)。半枫荷中多条次生代谢途径及相关Unigenes的发现，说明其次生代谢生物合成途径的复杂性，通过对半枫荷化学成分的分析研究，为后续从半枫荷中分离新化合物提供了线索，也为阐述半枫荷药效物质基础提供了理论依据。

图5 半枫荷转录组Unigenes的GO功能分类统计图 1.抗氧化活性；2.结合；3.催化活性；4.电子载体活动；5.metallochaperone activity；6.金属伴活动分子功能调节器；7.分子转导活性；8.核酸结合转录因子活性；9.营养库活性类；10.蛋白质标记；11.信号转导活性；12.结构分子活性；13.转录因子活动与蛋白质结合；14.传译调治活性；15.转运蛋白活性；16.受体活性；17.酶调节活性；18.受体调节活性；19.运转状态；20.生物粘附；21.生物调节；22.细胞聚集；23.细胞杀伤；24.细胞成分组织或生物合成；25.细胞过程；26.解毒；27.发育过程；28.生长；29.免疫系统过程；30.定位；31.运转；32.代谢过程；33.多组织进程；34.多细胞生物过程；35.生物过程负调控；36.生物过程的积极调节；37.生物过程的调节；38.繁殖；39.reproductive process；40.生殖过程对刺激的反应；41.节律过程；42.信号传导；43.单一有机体过程；44.建立斑点定位；45.细胞；46.细胞连接；47.细胞部分；48.细胞外基质；49.胞外基质成分；50.胞外区域；51.胞外区域部分；52.大分子复合体；53.膜；54.膜部分；55.膜封闭腔；56.拟核；57.细胞器；58.细胞器部分；59.其他有机体；60.其他有机体部分；61.超分子纤维；62.共质体；63.突触；64.突触部分；65.病毒；66.病毒部分Fig.5 GO functional classification of unigenes of transcriptome for S.cathayensis 1.antioxidant activity； 2.binding； 3.catalytic activity； 4.electron carrier activity； 5.metallochaperone activity； 6.molecular function regulator； 7.molecular transducer activity； 8.nucleic acid binding transcription factor activity； 9.nutrient reservoir activity； 10.protein tag； 11.signal transducer activity； 12.structural molecule activity； 13.transcription factor activity,protein binding； 14.translation regulator activity； 15.transporter activity； 16.receptor activity； 17.enzyme regulator activity； 18.receptor regulator activity； 19.behavior； 20.biological adhesion； 21.biological regulation； 22.cell aggregation； 23.cell killing； 24.cellular component organization or biogenesis； 25.cellular process； 26.detoxification； 27.developmental process； 28.growth； 29.immune system process； 30.localization； 31.locomotion； 32.metabolic process； 33.multi-organism process； 34.multicellular organismal process； 35.negative regulation of biological process； 36.positive regulation of biological process； 37.regulation of biological process； 38.reproduction； 39.reproductive process； 40.response to stimulus； 41.rhythmic process； 42.signaling； 43.single-organism process； 44.establishment of localization； 45.cell； 46.cell junction； 47.cell part； 48.extracellular matrix； 49.extracellular matrix component； 50.extracellular region； 51.extracellular region part； 52.macromolecular complex； 53.membrane； 54.membrane part； 55.membrane-enclosed lumen； 56.nucleoid； 57.organelle； 58.organelle part； 59.other organism； 60.other organism part； 61.supramolecular fiber； 62.symplast； 63.synapse； 64.synapse part； 65.virion； 66.virion part

图6 半枫荷转录组Unigenes的KEGG功能注释分布统计图 1.细胞生长和死亡；2.细胞运动；3.细胞群体；4.运输和分解代谢；5.膜运输；6.信号转导；7.信号分子和相互作用；8.折叠排序与退化；9.复制和修复；10.转录；11.翻译；12.氨基酸代谢；13.其他次级代谢产物的生物合成；14.碳水化合物代谢；15.能量代谢；16.多糖合成与代谢；17.类脂(化合)物代谢作用；18.辅助因子和维生素的代谢；19.其他氨基酸的代谢；20.萜类化合物和聚酮化合物的代谢；21.核苷酸代谢；22.概观；23.异生素生物降解和新陈代谢； 24.老化；25.循环系统；26.进化；27.消化系统；28.内分泌系统；29.环境适应；30.排泄系统；31.免疫系统；32.神经系统；33.感官系统Fig.6 KEGG functional annotation distribution of unigenes of transcriptome for S.cathayensis 1.Cell growth and death； 2.Cell motility； 3.Cellular community； 4.Transport and catabolism； 5.Membrane transport； 6.Signal transduction； 7.Signaling molecules and interaction； 8.Folding sorting and degradation； 9.Replication and repair； 10.Transcription； 11.Translation； 12.Amino acid metabolism； 13.Biosynthesis of other secondary metabolites； 14.Carbohydrate metabolism； 15.Energy metabolism； 16.Glycan biosynthesis and metabolism； 17.Lipid metabolism； 18.Metabolism of cofactors and vitamins； 19.Metabolism of other amino acids； 20.Metabolism of terpenoids and polyketides； 21.Nucleotide metabolism； 22.Overview； 23.Xenobiotics biodegradation and metabolism； 24.Aging； 25.Circulatory system； 26.Development； 27.Digestive system； 28.Endocrine system； 29.Environmental adaptation； 30.Excretory system； 31.Immune system； 32.Nervous system； 33.Sensory system

Table3TopsixteenmetabolicpathwaysinvolvedinS.cathayensisunigenes

代谢通路PathwayUnigenes数量Number of unigenesUnigene比例Percent of unigenes(%)代谢通路编号Pathway ID核糖体Ribosome3894.36ko03010碳代谢Carbon metabolism2212.48ko01200氨基酸合成Biosynthesis of amino acids1932.16ko01230内质网蛋白加工Protein processing in endoplasmic reticulum1711.92ko04141剪接Spliceosome1581.77ko03040RNA转运 RNA transport1541.73ko03013氧化磷酸化Oxidative phosphorylation1501.68ko00190嘌呤代谢Purine metabolism1361.52ko00230植物激素信号转导Plant hormone signal transduction1281.43ko04075胞吞作用Endocytosis1101.23ko04144泛素介导的蛋白水解Ubiquitin mediated proteolysis1081.21ko04120植物病原体相互作用Plant-pathogen interaction1081.21ko04626淀粉和蔗糖代谢Starch and sucrose metabolism1001.12ko00500糖酵解/葡糖糖生成Glycolysis/Gluconeogenesis991.11ko00010mRNA监测通路mRNA surveillance pathway961.08ko03015嘧啶代谢Pyrimidine metabolism941.05ko00240

2.7 TF编码能力及SNP检测

本研究共预测出88个基因家族共1 547个编码转录因子(Transcription Factor，TF)的Unigenes(图7)，另有57个转录因子未知具体基因家族；其中锌指蛋白C2H2类TF家族占比最多，为151个(9.76%)，其次为AP2/ERF-ERF、MYB-related和bHLH家族TF，分别占75(4.85%)、69(4.46%)和68(4.40%)个TF，而HRT、SOH1、STAT、ULT所占比例最少，均为1(0.065%)个。同时预测出同源异型盒(Homeobox，HB)家族TF、包括HB-HD-ZIP、HB-PHD、HB-KNOX、HB-WOX、HB-BELL、HB-other；可能参与调控药效成分合成有关的TF有AP2/ERF-ERF、WRKY、bHLH、bZIP、MYB等；还发现可能参与调节种子发育、萌发有关的TF有bZIP、MYB、MYB-related、C2C2-Dof、NAC、NF-YB、NF-YA、NF-YC、MADS-M-type、MADS-MIKC；与逆境抗性相关TF主要有：bZIP、AP2/ERF-ERF、AP2/ERFWRKY、HB类、MYB、MYB-related和HSF。这些TF功能的发现表明半枫荷转录调控的复杂性、交叉性和多样性，也为后续研究半枫荷濒危原因、药效药理研究、遗传多样性、指纹图谱、连锁图谱构建奠定基础。

表4 半枫荷药用活性成分相关的次生代谢产物生物合成途径

图7 半枫荷转录组Unigene主要的转录因子家族分类图Fig.7 The main transcription factor family classification of unigenes of transcriptome for S.cathayensis

图8 半枫荷转录组Unigene的SNP变异类型分布图Fig.8 SNP variant type distribution of unigenes of transcriptome for S.cathayensis

此外，本研究检测出125 917个单核苷酸多态性(Single Nucleotide Polymorphsims，SNP)位点(图8)，分析发现半枫荷Unigenes序列上SNP分布不均匀，转换突变类型(蓝色)的SNP数量比颠换突变(红色)的2倍略多。其中，转换突变类型(A→G、C→T、G→A、T→C)有84 445个，占67.06%，颠换突变类型(A→C、A→T、C→A、C→G、G→C、G→T、T→A、T→G)有41 472个，占32.94%。半枫荷SNP分析，结合Unigenes的表达量及表型，可在mRNA层面的基因型差异区别于其它物种，通过构建半枫荷SNP图谱，对了解其亲属关系和物种起源进化关系具有重要的生物学意义。

3 讨论

本研究采用RNA-Seq技术首次对半枫荷叶片进行了高通量转录组测序分析，拼接组装后获得77 629条Unigene，N50长度为1 105，平均长度为624.93 bp，均比川芎(Ligusticumchuanxiong)[23]、鹰嘴豆(Cicerarietinum)[24]、葫芦巴(Trigonellafoenum-graecum)[25]更长，说明本研究所得序列中长片段较多，拼接效果较好。半枫荷SSR出现频率为13.74%、比同是中药材的野三七(Panaxstipuleanatus)[26](16.82%)、半夏(Pinelliaternata)[27](16.24%)更低，但比党参(Codonopsispilosula)[28](12.22%)、莲雾[29](Syzygiumsamarangense)(9.27%)、川芎(Ligusticumchuanxiong)[23](8.11%)更高，其中有13 635条Unigenes的长度在1 kb以上，进一步说明半枫荷SSR位点分布较密集，可开发位点多。

所有的Unigene与CDD、KOG、NR、NT、GO、Swissprot、PFAM、KEGG、TrEMBL等9个数据库进行比对，共有45 293条(58.35%)成功获得功能注释，但仍有32 337条(41.65%)序列定位不清楚，得到功能注释基因相对较少，这与香樟(Cinnamomumcamphora)[30]，紫背天葵(Begoniafimbristipula)[31]，柳树(Salixbabylonica)[32]均出现大量未获得注释Unigenes的结果相似，可能是由于这些Unigenes本身为为非编码序列，或是技术上的限制及已测序的半枫荷近缘种较少有关[33]。成功注释到NR数据库中的Unigenes中，注释为葡萄同源种的有25.89%，远高于其他648个被注释物种，出现这种情况，可能是由于葡萄与半枫荷的进化史和生活史较为接近的缘故，亦可能因为葡萄具有参考基因组。而与半枫荷能比对上仅有6个同科植物，共有43条Unigenes，其中枫香(Liquidambarformosana)、北美枫香(Liquidambarstyraciflua)、缺萼枫香(Liquidambaracalycina)、马蹄荷(Exbucklandiapopulnea)、北美金缕梅(Hamamelisvirginiana)、红花檵木(Loropetalumchinense)，分别有31、8、1、1、1和1，之所以会出现这种情况，可能是由于金缕梅科植物转录组、基因组信息严重缺乏、也可能是因为本次测序得到的Unigenes中半枫荷特有的新功能基因较多导致未能获得相应匹配。

对拼接组装获得的45 293条Unigenes进行了代谢途径分析和功能分类，并根据KEGG数据库分析其可能参与或涉及的代谢通路，共获得4 029(5.19%)个Unigenes注释，注释类别可细分为5个代谢通路大类、33类代谢途径，其中定位到代谢相关通路的基因数最多，约占总注释量的42%，进一步证实半枫荷具有较强的代谢活动能力。值得一提的是，本研究共发现286个具体KEGG代谢途径分支参与到半枫荷碳水化合物代谢、生命活动代谢、次生物质代谢和遗传信息处理过程中，其中，232个Unigenes参与到合成药效次生代谢产物中、包括生物碱类、萜类、苯丙素类、黄酮类、醌类和酚类等，基于以上半枫荷药用代谢产物合成途径，初步构建了半枫荷活性成分基因调控网络。我们知道评价中药品质的优劣，不止是看某一成分含量的高低，而是成分与成分间相互影响的整体特征，因此本研究的结果及后续将要开展的EST-SSR分子标记为半枫荷药效质量评价标准提供基础数据。另外，本研究用半枫荷叶片做转录组分析，而半枫荷根药效比叶更大，可以预测半枫荷根相关药效次级代谢产物代谢强度更大，但需进一步验证。本研究根据组装结果预测出88个基因家族共1 547个编码转录因子的Unigenes，这些家族在半枫荷生长发育过程中发挥着不同的作用，其中C2H2锌指蛋白基因家族数量最多，可能参与激素的信号传导、生长发育、抗病信号途径、细胞分化、胁迫应答等相关活动[34]。此外，由于条件限制，本研究仅仅是发现了转录因子与药效成分合成有关，但具体调控的靶基因和作用机制还不清楚，需要进一步研究。本研究通过转录组测序分析，共检测出125 917个SNP多态位点，包括84 445个转换突变类型，41 472个颠换突变类型，即转换突变类型约占SNP的三分之二，与杜玮南[35]等得出的结论一致。通过SNP研究将为今后开展半枫荷遗传图谱构建、物种起源进化、遗传多样性分析、品种鉴定以及分子标记辅助选择育种等研究提供重要参考。

通过生物信息学方法对半枫荷进行CDS预测、功能注释、代谢相关通路、SNP检测和TF编码能力预测，综合KOG、GO和KEGG三大功能数据库的注释结果，表明半枫荷转录、复制和翻译等基因表达丰富度较高，具有较强的代谢活动和遗传信息处理能力，同时找出半枫荷药效物质可能相关的次生代谢产物生物合成途径及转录因子。所得转录组信息可作为金缕梅科植物基因组序列的重要组成部分，进一步丰富了该科植物的基因信息库，为进一步研究药用植物半枫荷的代谢途径，以及该物种的分子遗传育种和基因工程提供基础数据，也将为保护和利用这一濒危药用植物资源奠定基础。