瓠瓜幼叶转录组功能基因测序及分析

许端祥 赵瑞丽 陈中钐 杜文丽 徐同伟 高 山

(福州市蔬菜科学研究所，福建 福州 350111)

瓠瓜[Lagenaria siceraria(Molina) Standl. ]是原产于热带的葫芦科葫芦属的一年生草本植物，又称瓠子、扁蒲、葫芦，其嫩果品质柔嫩，味道鲜美，富含多种营养物质，具有美容、清热、解毒、治疗肺炎等功效，深受人们欢迎[1-2]。 瓠瓜根系发达，耐热性好、抗病性强，常用作西甜瓜等瓜果类蔬菜的砧木[1]。 瓠瓜在中国栽培历史悠久，主要分布在长江及其以南地区，近几年北方也开始引种栽培[3]。

近年来，随着高通量转录组测序技术的快速发展，紫背天葵[4]、党参[5]、青篱柴[6]、黄秋葵[7]、番茄[8]等物种已完成了转录组功能基因组测序和功能注释，并从中发掘出一些重要功能基因加以利用。 戎利勤等[9]完成了小花草玉梅的转录组测序，并从中筛选出12 个与花发育紧密相关的MADS基因，AGL6、SEP3、FUL1、PI2、SEP1 的表达量均为花形态建成的主要指标；林珲等[10]对青梗花椰菜和白梗花椰菜的转录组数据分析，发现有6 个差异基因与类胡萝卜素生物合成有关，类黄酮生物合成途径中得到9 个差异基因以及在叶绿素代谢途径中发现1 个差异基因与花梗颜色相关；叶新如等[11]以冬瓜嫩叶为材料，利用高通量测序技术获得大量冬瓜转录组信息，从中挖掘冬瓜基因数据及SSR 分子标记，有助于从分子水平对冬瓜进行深入研究；吴新义等[12]采用细胞流式技术测定了4 种中国瓠瓜的基因组大小；Wu 等[13]首次对瓠瓜自交系USVL1VR-Ls 基因组序列进行De Novo 测序分析，揭示了葫芦科植物基因组进化史，为葫芦科植物比较基因组学研究和品种改良提供了依据；Wang 等[14]建立瓠瓜基因组信息相关数据库，为瓠瓜基因组学研究提供了帮助。 由此可见，高通量转录组测序已经成为功能基因研究的高效手段。 本研究采用高通量测序技术对福州特有地方品种福州芋瓠瓠瓜叶片进行转录组测序，对获得的数据进行分析、组装、注释，从中发掘瓠瓜功能基因数据和SSR 分子标记，旨在为进一步研究瓠瓜各类功能基因奠定基础，为福州芋瓠品种保护和利用奠定基础。

1 材料与方法

1.1 试验材料

以福州本地瓠瓜品种福州芋瓠为材料，于2017 年11 月上旬播种于32 孔穴盘中，待幼苗长至3 ～4 片真叶时将4 株完整植株剪碎，混合取样，设3 次重复，液氮速冻-80℃保存备用。 样品委托上海凌恩生物科技有限公司进行转录组测序。

1.2 RNA 的提取纯化

采用Trizol 法提取瓠瓜叶片的总RNA，检测RNA的浓度和完整性；采用磁珠法分离mRNA，将得到的mRNA 逆转录成dsDNA，加入dA 碱基，连接接头，PCR富集测序样本。

1.3 转录组测序及组装

采用Illumina Hiseqxten 测序平台进行转录组测序，获得的图像数据经Base Calling 转化为原始序列(reads)，使用Trimmomatic (version 0.36， http:/ /www.usadellab.org/cms/index.php? page ＝trimmomatic)软件(LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:75)参数获得高质量的reads 序列。 再利用Trinity 软件从头组装生成contig，然后拼接组装成Transcript，最后使用Tgicl 和Phrap 软件进行同源聚类和拼接得到单基因簇(Unigene)。

1.4 序列注释和功能分类

通过BLAST 比对工具，将瓠瓜转录组获得的Unigene 与公共数据库进行比对，将E 值设置为≤1E-5，根据基因的相似性进行功能注释，得到与给定Unigene 具有最高序列相似性的蛋白，从而得到该Unigene 的蛋白功能注释信息。 公共数据库包括非冗余蛋白数据库(non-redundant protein database，Nr)、直系同源蛋白质家族数据库 (cluster of orthologous groups，COG)、基因本体论数据库(gene ontology，GO)、东京基因与基金组百科全书(kyoto encyclopedia of genes and genomes，KEGG)。

1.5 转录组SSR 位点分析

采用MISA 软件对获得的Unigene 进行SSR 位点搜索，以一、二、三、四、五和六核苷酸的重复次数分别不小于11、6、5、5、5 和5 次为搜索标准，运用MISA 软件识别SSR 及进行引物设计。

2 结果与分析

2.1 转录组数据组装分析

福州芋瓠是福州特有地方保护品种，其种植历史悠久，为更好地研究保护利用福州芋瓠，首先以Wu等[13]拼接的USVL1VR-Ls 基因组序列为参考基因组，进行有参转录组分析，共获得664 252 268 个原始reads片段，经过滤后得到高质量的测序数据525 959 048 个reads 片段，共计77 337 803 108 bp(77.34 Gb)个核苷酸，碱基Q20 大于98%，Q30 大于94%，说明测序质量较高，该测序结果可用于后续分析，但mapping 到参考基因组上的比对率很低，绝大部分的样品比对率在50%～60%，均低于70%，猜测导致这种结果的原因可能是本研究测序的福州芋瓠与Wu 等[13]所测序的瓠瓜基因组的品种USVL1VR-Ls 遗传差异较大，为此转而实施无参转录组分析。

将获得的高质量的数据经组装后共获得107 600个Transcript，总长度达127 115 912 bp，Transcript 序列再次聚类和组装后共获得87 518 个Unigene，总长度高达91 405 320 bp，平均长度为1 044 bp，N50 为1 386 bp， 其中最长的Unigene 为60 483 bp。 87 518个Unigene，长度主要集中在1 ～1 000 bp 的有66 538个，占总Unigene 的76.03%，1 000 ～2 000、2 000 ～4 000 bp、＞4 000 bp 的Unigene 分别占15.77%、5.01%、3.20%(图1)。

图1 瓠瓜Unigene 的长度分布Fig.1 Length distribution of bottle gourd Unigenes

2.2 瓠瓜Unigene 的功能注释

2.2.1 Nr 功能注释 将组装后得到的87 518 个Unigene 在Nr 数据库进行注释，共注释到55 725 个Unigene。 在Nr 注释中，E 值介于1E-50～1E-10 之间的Unigene 数量最多，占总Unigene 的41.71%(23 242个)，其次是1E-100 ～1E-50，为25.38%(14 144 个)(表1)。 Nr 注释中相似度相对较高，高于90%的Unigene 有23 412 个(42.01%)。 Nr 功能注释到葫芦科的Unigene 数量较多，有26 354 个，占总Unigene 的47.29%，其中注释到甜瓜上的Unigene 最多，有10 505个(18.85%)，注释到黄瓜上的Unigene 有9 166 个(16.45%)，注释到南瓜上的Unigene 有4 888 个(8.77%)，注释到苦瓜上的Unigene 有1 795 个(3.22%)。 此外，Unigene 在其他植物中也有注释到，如注释到蓖麻上的Unigene 有4 785 个(8.59%)，注释到水桔梗上的Unigene 有2 818 个(5.06%)，注释到牛头蛎球菌上的Unigene 有1 621 个(2.91%)，注释到枣上的Unigene 有1 276 个(2.29%)，注释到水稻上的Unigene 有1 227 个(2.20%)(图2)。

表1 瓠瓜Unigene 在Nr 数据库中的E 值分布和相似度分布Table 1 E-value and identity distribution of bottle gourd in Nr database

图2 瓠瓜Unigene 的Nr 库物种分布Fig.2 Species distribution in Nr library of bottle gourd Unigene

2.2.2 GO 功能分析 由表2 可知，瓠瓜Unigene 在GO 数据库中共注释到18 278 个Unigene，占总Unigene 的20.93%。 注释到的Unigene 根据其功能可以分为分子功能、细胞组分和生物学过程三大类和55个亚类，其中分子功能包含14 个亚类，结合活性、催化活性的Unigene 较多，分别有8 059(44.09%)、9 517个(52.07%)，而归类到蛋白标签、金属蛋白活性、营养储存活性、抗氧化活性、分子转导活性、转录因子活性-蛋白结合的Unigene 较少。 在细胞组分中，归类到细胞器、膜部分、膜结构、细胞、细胞部分的Unigene 较多，分别为3 075(16.82%)、4 269(23.36%)、4 974(27.21%)、4 723(25.84%)、4 598 个(25.16%)，归类到其他亚类的Unigene 较少。 在生物学过程功能中，归类到单一有机体进程、细胞过程、代谢进程的Unigene较多，分别为5 666 (31.00%)、 8 307(45.45%)、9 118 个(49.89%)，归类到其他类别的Unigene 均较少。

2.2.3 COG 功能分类 由表3 可知，瓠瓜叶片转录组测序中有41 635 个Unigene 在COG 数据库中获得注释，占总Unigene 的47.67%，可分为信息存储与处理、细胞过程和信号传递、新陈代谢、无特征基因四大类，其中新陈代谢、细胞过程和信号传递类的Unigene较多，分别为19 438(46.69%)、8 679 个(20.85%)，信息存储与处理类有7 620 个(18.30%)，无特征基因类5 898 个(14.17%)。 四大类中的Unigene 可根据其功能分为25 个功能区域，其中参与氨基酸运输和代谢的Unigene 数量最多，占10.82%，翻译、核糖体结构和生物合成(6.03%)、重复/重组和修饰(7.37%)、能量生成和转换(8.26%)、 碳水化合物运输和代谢(7.18%)、无机离子运输和代谢(7.39%)功能区域的Unigene 也较为丰富，而其他功能区域的Unigene 数量相对较少。 此外，还有2 584 个Unigene 的功能未知。

2.2.4 KEGG 代谢通路分析 将获得的瓠瓜Unigene与KEGG 数据库进行比对， 由表4 可知，有24 770 个Unigene 获得注释，占总Unigene 的28.4%。 统计分析显示，注释到的Unigene 共涉及到220 个KEGG 代谢途径，其中参与代谢途径的Unigene 最多，为9 456 个，占比38.18%；参与次生代谢产物的生物合成的Unigene 有4 841 个，占比19.54%；生态环境下的微生物代谢作用的Unigene 有2 528 个，占比10.21%；氨基酸生物合成的Unigene 有1 990 个，占比8.03%；碳代谢的Unigene 有1 828 个，占比7.43%；参与其他代谢途径的Unigene 相对较少。 上述结果说明，瓠瓜的代谢活动较为活跃。

表2 瓠瓜Unigene 的GO 功能分类Table 2 GO functional categories of bottle gourd Unigene

表2(续)

表3 瓠瓜Unigene 的COG 功能分类Table 3 COG functional categories of bottle gourd Unigene

表3(续)

表4 瓠瓜Unigene 的KEGG 代谢通路分析Table 4 KEGG pathway analysis of bottle gourd Unigene

表4(续)

表4(续)

表4(续)

表4(续)

表4(续)

2.3 瓠瓜转录组的SSR 位点分析

MISA 软件检索结果显示，瓠瓜87 518 条Unigene中有8 617 条Unigene 中存在SSR 序列，共检索到11 029 个SSR 位点，发生频率为9.846%，平均每8.16 kb 出现1 个SSR。 其中1 858 条Unigene 含有2 个或2 个以上的SSR 位点，920 条Unigene 属于复合SSR。在6 种 SSR 重复类型中， 单核苷酸比例最高(55.51%)，其次是三核苷酸(25.41%)，二核苷酸占17.07%(表5)。 瓠瓜6 种SSR 重复类型中包含有87种重复基序，单、二、三、四、五、六核苷酸6 种重复类型分别包含2、4、10、20、18、33 种重复基序类型。 其中，单核苷酸以A/T 为主，发生频率为98.22%，总体占比54.52%；二核苷酸以AG/CT，其次为AT/AT；三核苷酸以AAG/CTT 为主，四、五、六核苷酸重复基序类型多，数量少(表6)。

表5 瓠瓜转录组的SSR 类型、数量及分布频率Table 5 Type， number and frequency of SSRs in bottle gourd

3 讨论

随着瓠瓜商品化，其种植面积逐年扩大，瓠瓜种质繁育和病虫害方面的研究已受到广泛关注[15-16]，关于分子生物技术的研究主要集中在种质资源的遗传多样性方面[17-19]，对瓠瓜功能表达基因克隆方面的研究也有部分报道[3]。 2016 年吴新义等[12]采用细胞流式技术测定了4 种中国瓠瓜的基因组大小为29.11 ～344.56 Mb；Wu 等[13]首次对瓠瓜自交系USVL1VR-Ls进行基因组De Novo 测序分析，揭示了葫芦科植物基因组进化史；Wang 等[14]建立了瓠瓜基因组信息相关数据库，为瓠瓜基因组学研究提供帮助。

本研究采用Illumina Hiseqxten 测序平台对瓠瓜叶片进行转录组测序分析，获得的Q20 和Q30 分别为98%和94%，N50 为1 386 bp，说明本次测序质量较高，满足后续分析的质量要求，但以Wu 等[13]拼接的USVL1VR-Ls 基因组序列为参考基因组，进行有参转录组分析时的比对率很低，绝大部分的样品比对率为50%～60%，均低于70%，猜测导致这种结果的原因可能是本研究所测序的福州芋瓠与Wu 等[13]所测序的瓠瓜品种USVL1VR-Ls 遗传差异较大，因此，本研究进行无参转录组分析来分析福州芋瓠的特异性。 将获得的序列进行拼接组装后共获得87 518 个Unigene，平均长度为1 044 bp，序列相对较长，其中≥1 kb 的Unigene 有20 980 个，远高于紫色黄秋葵[20]、黑油椿香椿[21]和紫背天葵[4]，其原因可能是不同物种中大片段Unigene 的数量不同。

将组装后获得的全部Unigene 与Nr 公共数据库比对，有55 725 个Unigene 获得功能注释，占63.81%，36.19%的Unigene 未获得相关注释信息，此现象在其他植物中也存在[22-23]，其原因可能是Unigene 片段太小未匹配到单数据序列，也可能是瓠瓜中存在新的功能基因。 在GO 数据库中有18 278 个Unigene 获得注释，占总Unigene 的20.93%，相对注释率较低，这与GO 数据库信息不够完善有关[4]。 COG 数据库注释到的Unigene 较多，占总Unigene 的72.24%，功能涉及到信息存储与处理、细胞过程和信号传递、新陈代谢等各类生命活动，其中与新陈代谢有关的Unigene 占46.69%，说明瓠瓜自身具有较强的代谢能力。 KEGG全局通路注释到24 770 个Unigene 涉及到220 个KEGG 代谢途径，参与次生代谢产物的生物合成、微生物代谢、氨基酸生物合成、碳代谢等代谢途径的Unigene 最多，这与COG 分析结果相一致，此结果为研究瓠瓜代谢途径相关功能基因提供了依据。

根据葫芦科植物中部分SSR 引物具有通用性的特点，前人已成功筛选出一些适用于瓠瓜种子纯度鉴定的SSR 引物[24]，但仍不能满足研究和应用的需求，许多研究表明，可以从转录组序列数据中发掘、开发出丰富的SSR 标记[25-27]。 本研究从瓠瓜转录组87 518条Unigenes 中筛选出11 029 个SSR 位点，发生频率为9.846%，该结果高于南瓜[25]和杏鲍菇[26]，但低于丝瓜[27]，这可能与物种基因组差异、测序方法、SSR 筛选条件等因素有关。 本试验中，瓠瓜6 种SSR 重复类型中单、二、三核苷酸出现频率较高，四、五、六核苷酸出现频率相对较低，这与李荣华等[28]、李海波等[29]的研究结果相似。 本研究中，瓠瓜包含87 种SSR 重复基序类型，以A/T、AG/CT 和AAG/CTT 含量最为丰富，这与朱海生等[27]、李海波等[29]的研究结果相似。

表6 瓠瓜转录组的SSR 重复基序类型及频率Table 6 Motif type and frequency of SSRs in bottle gourd

4 结论

本研究利用Illumina Hiseq 高通量测序平台对福州芋瓠瓠瓜叶片进行无参转录组分析，共获得87 518个Unigene， 55 725 个Unigene 在Nr 数据库获得注释，其中26 354 个Unigene 注释到葫芦科作物中；GO 数据库中注释到18 278 个Unigene 可以分为分子功能、细胞组分和生物学过程功能三大类和55 个亚类；COG数据库中有41 635 个Unigene 获得注释，分布在信息存储与处理、细胞过程和信号传递、新陈代谢、无特征基因四大类的25 个功能区域；KEGG 数据库中有注释到的24 770 个Unigene 涉及到220 个KEGG 代谢途径，可以全面了解瓠瓜叶片的代谢途径信息；在瓠瓜叶片转录组中存在11 029 个SSR 位点，发生频率为9.846%。 本研究得到大量福州芋瓠瓠瓜叶片基因序列，了解了瓠瓜叶片基因的表达情况，为后续深入开展瓠瓜功能基因研究及SSR 分子标记开发奠定了基础。