摘 要 为开发木荷分子标记，采用高通量测序技术获得木荷基因组原始数据，经生物信息学软件对木荷基因组序列进行序列拼接、组装和对比，共获得308 418条Contig序列和459 984条Scafford序列。采用MISA软件搜索木荷基因组序列中微卫星（Microsatellite）位点，共得到334 843个SSR序列，总长度5 074 708 bp，占木荷基因组大小的0.98%，木荷基因组SSR序列平均长度为15.2 bp，平均分布频率为644个/Mb。木荷SSR序列中，单核苷酸序列数量最多，共188 217个，占木荷SSR序列总数的56.21%，其次是二核苷酸（23%）>三核苷酸（13%）>四核苷酸（5%）>五核苷酸（2%）>六核苷酸（1%）。木荷全基因组SSR序列中共包括400种重复基元，其中单核苷酸重复基元A和二核苷酸重复基元AT是主要重复基元，分别占总SSR的56%和11%，SSR基元的重复次数分布在4～40次，主要分布在4～25次。本研究丰富了木荷分子标记类型，为进一步群体遗传结构和遗传多样性分析提供了基础数据。
关键词 木荷；SSR标记；生物信息学；基因组；高通量
中图分类号 S792.99 文献标识码 A
Abstract In order to develop the molecular marker of Schima superba， high-throughput sequencing was used to obtain the original data of S. superba. Bioinformatics softwares were used to sequence， assemble and compare the S. superba genomic sequences. A total of 308 418 contig sequences and 459 984 scafford sequences were obtained. In total， 334 843 SSR sequences were obtained by using the bioinformatics analysis software MISA. The total length of the SSR was 5 074 708 bp， accounting for 0.98% of the genome size of S. superba， with an average length of 15.2 bp. The mean distribution frequency was 644/Mb. Among the SSR of S. superba genome， the number of single nucleotides was the largest， with a total number of 188 217， accounting for 56.21% of the total SSR of S. superba， followed by dinucleotide （23%） > trinucleotides （13%） > teranucleotides （5%） > pentanucleotides （2%） > hexanucleotide （1%）. The SSR sequences of S. superba genome was consisted of 400 repeat motifs， in which mononucleotide repeat A and dinucleotide repeat AT were the major repeat motifs， accounting for 56% and 11% of the total SSR respectively. The number of repetitions of SSR motifs was distributed from 4 to 40， mainly between 4 and 25 times. This study would provide basic data for further genetic diversity analysis of S. superba.
Keywords Schima superba； simple sequence repeat； bioinformatics； genomes； high-throughput sequencing
DOI 10.3969/j.issn.1000-2561.2018.09.013
微卫星（Microsatellite），又称简单重复序列（Simple Sequence Repeats，SSRs），在原核生物和真核生物基因組中普遍存在，分布于基因间隔区、UTR区、内含子区等区域。SSR标记为共显性标记，具有数量多、分布广、多态性高和信息量大等优点[1]。SSR序列由1～6个重复并且串联的核苷酸组成，由高度变异的核心序列和高度保守的侧翼序列组成，可根据两端保守的侧翼序列针对性地设计引物，用于多态性的检测[2]。SSR重复类型包括：单核苷酸（Mononucleotide）重复、二核苷酸（Dinucleotide）重复、三核苷酸（Trinucleotide）重复和四核苷酸（Teranucleotide）及以上重复[3]，长度大多为100～200 bp，根据微卫星片段长度的不同又形成简单序列长度多态性（Simple Sequence Length Polymorphism，SSLP）[4]。近年来，SSR分子标记技术作为第二代微卫星分子标记技术，已广泛应用于群体遗传多样性分析[5]、种群内和种群间遗传变异分析[6-7]、连锁遗传图谱构建和群体遗传结构分析等[8-9]。
木荷（Schima superba Gardn.et Champ）属山茶科（Theaceae），木荷属（Schima），常绿阔叶乔木，自然分布于福建、广东、江西、浙江等31N以南，105E以东的地区[10]。木荷为典型的亚热带树种，叶呈厚革质，含水量高达42%，混植于易发生火灾的松树林，可以起到很好的防火阻隔作用[11]。木荷作为南方重要的防火树种，具有很大的研究价值，近几年来已经成为林木遗传育种专家重要的研究树种。目前，从有关木荷SSR位点分析的文献来看，研究者们采用的都是传统的分子生物学方法[12-13]，囿于依赖传统的实验技术方法，具有一定的局限性。
过去，人们对某些疾病现象的揭示局限于单个基因位点的研究，往往不能全面的揭示这些原理，随着分子生物技术的进一步发展，第二代测序技术的出现将研究者们带到了全基因组整体水平的研究层面，打破了传统只对单个基因位点研究的局限，具有通量高、周期短、成本低和信息量全面等优点[14]。然而，面对如此巨大且复杂的数据量，生物信息学的作用不可忽视，生物信息学和第二代测序的出现大大推动了分子生物学的发展。迄今为止，未见有利用高通量测序和生物信息学方法开发木荷基因组SSR的报道。本研究在高通量测序的基础上，利用生物信息学方法搜索木荷基因组SSR位点，分析木荷SSR序列生物信息学特征，为进一步遗传结构和遗传多样性分析提供了基础数据，还有助于利用分子标记辅助育种方法选育我国木荷优良种质资源。
1 材料与方法
1.1 材料
实验材料为单一木荷无性系，来源于广东省英德市林业科学研究所桥头基地木荷优树无性系种质基因库，占地面积75 hm2，相对湿度77%，年降雨量1 906.2 mm，年平均气温21.1 ℃。采集木荷幼嫩叶片5～10片，写上标签并用密封袋装好，带回实验室80 ℃超低温冰箱保存，用于DNA提取。
1.2 方法
1.2.1 木荷DNA的提取 采用常规试剂盒方法提取木荷DNA，经1%琼脂糖凝胶电泳检测，用紫外分光光度计（NanoDrop 2000 spectroph¬oto¬me¬ter）检测DNA浓度和纯度，置于20 ℃低温保存。
1.2.2 木荷基因组的高通量测序 采用高通量测序平台Illumina Hiseq 2500进行测序，测序读长PE250 （Pair-end 250）或PE150（Pair-end 150）。
1.2.3 数据处理 数据处理的过程分为4个步骤：数据质控，筛选高通量测序数据中质量较高的数据（Clean Reads），筛选标准为：去除含N碱基或含有3个以上碱基质量值低于20的测序序列；序列组装，采用Denove[15]序列组装软件，将测序数据组装成scaffold序列；序列标签聚类，采用序列比对方法，为减少比对误差导致的SSR检测错误，要求比对到SSR区域的测序片段为双末端同时比对到参考序列。将相同的序列标签聚类到一起，为后续的多态性检测做准备。最终将得到的数据以fasta格式保存用于后续分析。
1.2.4 数据分析 采用生物信息学分析软件MISA （http：//pgrc.ipk-gatersleben.de/misa） 结合Perl语言，搜索基因组中所有SSR。统计标准设置如下：要求单核苷酸重复不少于10次，二核苷酸重复不少于6次，三核苷酸重复、四核苷酸重复、五核苷酸重复和六核苷酸重复不少于4次。重复序列两边的侧翼序列大小设置为100 bp。
2 结果与分析
2.1 高通量测序数据分析
利用高通量测序技术共获得34 574 976条原始数据（Raw reads），经过序列组装和去除冗余序列得到34 379 120条高质量数据（High-quality reads）（表1），占Raw reads的99.4%，对SSR序列进行拼接分别得到308 418条和459 984条长度大于500 bp的contigs和scafford（图1），碱基最长的contig序列长度为6 108 bp，碱基最长的scaffold序列长度为19 000 bp。
2.2 木荷SSR位点的数量分布
通过SSR生物信息学软件MISA共检测到334 843个SSR序列，木荷基因组SSR平均出现频率为644个/Mb。SSR序列在木荷基因组上种类丰富，单核苷酸至六核苷酸均有分布，但数量上差异较大（表2）。其中，单核苷酸重复位点数量最为丰富（188 217个），占56.21%，二核苷酸重复位点数量高达五分之一（22.36%），三核苷酸重复位点次之（43 631个），占13.03%，四核苷酸重复位点共17 559个，占5.24%，五核苷酸重复位点6 277个，占1.87%，六核苷酸重复位点数量最少，3 027个，仅占0.90%，由此可见，木荷基因组SSR重复基元的数量随着碱基重复次数的增加呈降低趋势（图2，图3）。
2.3 木荷基因组SSR重复基元数量、频率、长度特征
木荷基因组SSR各重复基元的种类和各重复基元的数量分布呈显著差异（表3，图3，图4）， 各重复类别的优势重复基元数量分布情况见图5。从木荷基因组中找到的SSR种类共400种，其中，单核苷酸重复类别共2种，以A重复基元占绝大多数（187 723个），占单核苷酸重复类别的99.7%，占所有SSR总数的两种重复基元的56.1%，出现的频率分别是361.01个/Mb，C重复基元出现频率仅0.59个/Mb。在二核苷酸重复类别中，以AT重复基元数量最多（36 551个），占二核苷酸重复基元的48%，出现频率为70.29个/Mb，AG重复基元数量和频率略低于AT重复基元，分别为34 716个和66.8个/Mb。AC重复基元和CG重复基元出现的次数较低，出现频率分别为9.22个/Mb和0.14个/Mb；在三核苷酸重复类别中，以AAT重复基元出现的次数最多（20 913个），出现频率为40.22个/Mb，然后依次是AAG>ACC>ATC，出现频率分別为15.97个/Mb、8.67个/Mb和5.33个/Mb。
木荷基因组四核苷酸重复类型共32种，AAAT为主要的重复基元（13 468个），出现频率最高，为25.9个/Mb，其次为AAAC重复基元和AAAG重复基元，出现频率分别为1.50个/Mb和1.51个/Mb。五核苷酸重复类型一共93种，数量最多的重复基元是AAAAT，出现频率为5.29个/Mb。六核苷酸重复拷贝类型的种类最多，一共有259种，占所有重复基元数的65%，数量最多的重复类别为AAAAAT，出现频率为1.47个/Mb，其余重复基元不占优势，出现的频率均较低，重复基元种类随重复次数增加。
木荷基因组SSR序列总长度为5 074 708 bp，平均长度为15.12 bp，总长度占木荷总基因组大小的0.98%，单核苷酸重复类型的总长度最长，从单核苷酸至六核苷酸长度依次降低，各重复类别平均长度分布规律与总长度相反，从单核苷酸至六核苷酸呈升高的趋势，依次为12.74、19.62、15.20、17.61、21.32、25.81 bp（表4）。六核苷酸重复基元平均长度最长，其中单核苷酸重复类型共两种，分别为A和C，二核苷酸重复基元类别共四种，三核苷酸、四核苷酸、五核苷酸重复类型以及六核苷酸重复类型种类逐渐递增，依次为10种、32种、93种和259种，六核苷酸重复种类最多，占所有核苷酸种类的65%，由此可见，重复类型种类随着碱基重复次数的增加而增加。
木荷基因组SSR中，单核苷酸重复基元的分布范围在10～40次，主要分布在10～25次，占99%；二核苷酸重复基元分布范围为6～40次，主要分布在6～20次，占97%；三核苷酸重复基元分布在4～21次，主要分布在4～13次，占99%；四核苷酸重复基元分布在4～10次，主要分布在4～7次，占99%；五核苷酸重复基元分布范围为4～9次，主要分布在4～5次，占96%；六核苷酸重复基元范围分布在4～9次，主要分布在4～5次，占95% （图6）。
3 讨论
分析木荷基因组序列共得到334 843个SSR序列，平均每个SSR位点出现的频率为644个/Mb，相比茶树[16]、柑橘[17]、杨树[18]等植物，出现频率更高，说明木荷基因组SSR位点较为丰富。木荷基因组SSR各重复类别中，单核苷酸重复类别出现的频率最高，为361.96个/Mb，其次为二核苷酸重复（146.41个/Mb）>三核苷酸重复（83.91个/Mb）>四核苷酸重复（33.77个/Mb）>五核苷酸重复（12.07个/Mb）>六核苷酸重复（5.82个/Mb）。
目前，SSR标记技术已大量用于农作物和动物基因SSR位点开发[19-20]。木荷基因组SSR重复序列中以单核苷酸重复类型最为丰富，其次是二核苷酸重复类型，这与在模式植物杨树上的研究结果相似[21]。郑燕等[22]在水稻、玉米、高粱这3种禾本科植物的研究中，发现三核苷酸重复出现次数最多，其次为六核苷酸SSR重复序列。孙蛟龙等[23]在浮萍转录组中发现二核苷酸重复类型数量最多，其次为三核苷酸重复类型。由此可见，SSR序列的分布规律在不同的物种上差异较大，这可能是由于物种间基因组大小不同造成优势SSR序列差异。在木荷基因组单核苷酸至六核苷酸重复类别中，数量最多的重复基元分别为A、AT、AAT、AAAT、AAAAT以及AAAAAT，说明这些优势重复基元富含A/T核酸，在木荷基因组中大量存在，这与枣树[24]、杨树[25]、黑麦草[26]等植物基因组SSR序列研究一致。
本研究基于高通量测序技术，首次利用生物信息学分析方法发掘木荷基因组SSR序列，统计分析木荷基因组SSR数量、长度、频率和密度等生物信息学特征，成功开发了大量木荷SSR位点，丰富木荷分子标记类型，为木荷群体遗传结构和遗傳多样性分析提供数据基础，有利于加快建立木荷种质资源评价与保护机制，也有利于进一步木荷种质资源多样性研究和品种选育。
参考文献
[1]张立荣， 徐大庆， 刘大群. SSR和ISSR分子标记及其在植物遗传育种研究中的应用[J]. 河北农业大学学报， 2002（1）： 90-94.
[2]Xu M， Liu X， Wang J W， et al. Transcriptome sequencing and development of novel genic SSR markers for Dendrobium officinale[J]. Molecular Breeding， 2017， 37： 18.
[3]鄢秀芹， 鲁 敏， 安华明. 刺梨转录组SSR信息分析及其分子标记开发[J]. 园艺学报， 2015，42（2）： 341-349.
[4]Xie J K， Chen D Z， Xiao Y Q， et al. Genetic diversity of Dongxiang wild rice （Oryza rufipogon Griff.） detected by SSLP markers[J]. Agricultural Sciences in China， 2003， 2（7）： 12-18.
[5]白 玉. DNA分子标记技术及其应用[J]. 安徽农业科学， 2007， 35（24）： 7 422-7 424.
[6]左力辉， 张 双， 梁海永， 等. 榆树转录组EST-SSR引物开发及遗传多样性分析[J]. 植物遗传资源学报， 2018， 19（1）： 157-166.
[7]杨育峰， 史典义， 王雁楠， 等. 基于转录组测序数据的甘薯SSR标记开发及群体聚类分析[J].分子植物育种， 2018， 16（11）： 3 569-3 579.
[8]洪彦彬， 梁炫强， 陈小平， 等. 花生栽培种SSR遗传图谱的构建[J]. 作物学报， 2009， 35（3）： 395-402.
[9]徐立安， 李新军， 潘惠新， 等. 用SSR研究栲树群体遗传结构[J]. 植物学报， 2001， 43（4）： 409-412.
[10]杨智杰， 陈光水， 谢锦升， 等. 杉木、木荷纯林及其混交林凋落物量和碳归还量[J]. 应用生态学报， 2010， 21（9）： 2 235-2 240.
[11]赵秀军. 木荷的生物学特征及其栽培技术要点[J]. 南方农业， 2017， 11（17）： 50-51.
[12]辛娜娜， 张 蕊， 徐肇友， 等. 木荷1代育种群体遗传多样性分析[J]. 林业科学研究， 2015， 28（3）： 332-338.
[13]金则新， 李钧敏， 李建辉. 木荷种群遗传多样性的ISSR分析[J]. 浙江大学学报（农业与生命科学版）， 2007， 33（3）： 271-276.
[14]周家蓬， 裴智勇， 陈禹保， 等. 基于高通量测序的全基因组关联研究策略[J]. 遗传， 2014， 36（11）： 1 099-1 111.
[15]Li G， Zhao Y， Liu Z， et al. De novo assembly and characterization of the spleen transcriptome of common carp （Cyprinus carpio） using Illumina paired-end sequencing[J]. Fish & Shellfish Immunology， 2015， 44： 420-429.
[16]杨 华， 陈 琪， 韦朝领， 等. 茶树转录组中SSR位点的信息分析[J]. 安徽农业大学学报， 2011， 38 （6）： 882-886.
[17]Chen C， Zhou P， Choi YA， et al. Mining and characterizing microsatellites from citrus ESTs[J]. Theoretical and Applied Genetics， 2006， 112 （7）： 1 248-1 257.
[18]张新叶， 宋丛文， 张亚东， 等. 杨树EST-SSR标记的开发[J]. 林业科学， 2009， 45（9）： 53-59.
[19]吕远大， 李 坦， 石 丽， 等. 基于全基因組重测序信息开发玉米H99自交系特异分子标记[J]. 作物学报， 2014， 40（2）： 191-197.
[20]白俊艳， 王玉琴， 庞有志， 等. 绒山羊EST序列的生物信息学分析[C]//中国畜牧兽医学会信息技术分会. 中国畜牧兽医学会信息技术分会2012年学术研讨会论文集. 北海： 中国畜牧兽医学会. 2012： 5.
[21]Li S， Yin T. Map and analysis of microsatellites in the genome of Populus： The first sequenced perennial plant[J]. Science in China Series C： Life Sciences， 2007， 50（5）： 690-699.
[22]郑 燕， 张 耿， 吴为人. 禾本科植物微卫星序列的特征分析和比较[J]. 基因组学与应用生物学， 2011， 30（5）： 513-520.
[23]孙蛟龙， 方 扬， 靳艳玲， 等. 浮萍转录组数据SSR位点的生物信息学分析[J]. 应用与环境生物学报， 2015， 21（3）： 401-405.
[24]肖 京. 枣基因组SSR位点特征分析及引物开发[D]. 保定： 河北农业大学， 2014.
[25]蒋雪梅， 胡廷章， 向兴胜， 等. 杨树全基因组微卫星序列的统计及其生物信息学分析[J]. 西南农业学报， 2015， 28（2）： 527-533.
[26]Hou S， Sun Z， Linghu B， et al. Genetic diversity of Buckwheat Cultivars （Fagopyrum tartaricum Gaertn.） assessed with SSR markers developed from genome survey sequences[J]. Plant Molecular Biology Reporter， 2016， 34（1）： 233-241.