福建内陆山区棘胸蛙遗传多样性分析

王茂元，赖铭勇，黄洪贵，吴妹英，田田，黄柳婷

（福建省淡水水产研究所，福建 福州 350002）

棘胸蛙Rana spinosa（俗称石蛙、石鳞、石鸡、石冻等）主要分布于中国南部以及越南北部丘陵山区[1]。棘胸蛙肉质鲜美、口感爽滑、营养丰富，具有较高的食用和医用价值，是一种经济价值较高的名贵珍品。据《本草纲目》记载：“石蛙主治小儿痨瘦，疳疾最良，产妇尤佳”，市场需求量大[2]。由于人为捕杀和自然环境的恶化，野生棘胸蛙资源遭到严重破坏，已被中国物种红色名录列为易危（Vu）等级[3]。近年来，随着棘胸蛙养殖技术逐步成熟，养殖业者为扩大养殖规模，不断进行频繁引种、近亲繁殖、累代人工养殖等，存在产量下降、品质退化、生物遗传多样性降低等诸多风险，严重危及棘胸蛙种质资源。因此，开展棘胸蛙养殖群体遗传多样性研究，从分子水平分析其种群遗传结构特征和种群分化，对棘胸蛙的种质资源保护与科学开发利用具有十分重要的意义。

测序基因分型技术（Genotyping-by-sequencing，GBS）[4]是近年来常见的单核苷酸多态性标记（Single Nucleotide Polymorphism，SNP）开发最为有效、经济的简化基因组测序方法之一。其原理是采用限制性内切酶加标签的方法，实现多样本高通量平行测序[5]，获取目的物种在基因组层面上的高密度SNP标记。这些标记可用于遗传图谱构建[6,7]、种质鉴定[8,9]、基因多样性[10,11]及群体遗传学[12,13]等多方面的研究。本研究通过利用GBS 简化基因组测序技术分析三明、龙岩和南平三个地区棘胸蛙群体的遗传多样性，旨在从分子水平揭示不同地域棘胸蛙的种群遗传结构，进一步了解福建棘胸蛙的遗传特征和种群分化，为棘胸蛙种质资源保护和遗传改良提供参考。

1 材料与方法

1.1 材料

2017 年10 月，于南平市正绿石蛙养殖有限公司取试验用棘胸蛙23 只、三明市龙源棘胸蛙养殖场20 只以及龙岩市武平鑫旺质胜农业发展有限公司23 只，共计66 个样品，平均体质量为（143.79±26.55）g，均为当地野生蛙群子一代3 龄蛙，分别命名为Rs-N、Rs-S 和Rs-L，采用低温麻醉，使棘胸蛙进入休眠状态，活体解剖取后肢于95%乙醇中保存备用。

1.2 方法

1.2.1 DNA 的提取与纯化

采用SDS 法提取棘胸蛙后肢的肌肉组织全基因组DNA，在1%的琼脂糖凝胶上100 V 电泳40 min，检测样本基因组DNA 的完整性。用核酸浓度测定仪（BioDrop-μLite）测定其浓度和纯度，质检符合要求的高质量DNA 样本进行后期的文库构建和测序。

1.2.2 GBS 文库构建

选择质检合格的基因组DNA 样本，采用限制性内切酶ApekⅠ酶切基因组DNA。酶切后的片断链接带有内切酶序列末端的普通接头和探针接头，再将不同样品的链接产物集中，通过凝胶电泳回收180～480 bp 区间的DNA 片段，再进行PCR 扩增，用磁珠富集纯法，构建测序样本文库。

1.2.3 测序

构建好的文库用Agilent 2100 Bioanalyzer 和QPCR 的方法进行质控，通过质控的文库进行测序。利用第二代测序技术（Next-Generation Sequencing，NGS），在Illumina HiSeq4000 测序平台上进行双末端150 bp 测序。

1.3 数据统计分析

1.3.1 数据处理

由测序所得的原始数据（Raw data）进行质控过滤，包括去除含有接头序列的双末端序列；去除N的含量超过该条序列长度比例的10%时的双末端序列；去除低质量的序列（质量值Q≤5 的碱基占整条序列的50%以上的被定义为低质量序列）；去掉无法依据接头区分样品的序列等。过滤后的数据进行后续分析。

1.3.2 聚类与比对

原始数据经过过滤后获得的数据为有效数据（Clean data）。将每个样本的有效数据分别与拟参考基因组非洲爪蟾Xenopus laevis 进行BWA 0.5.7 比对分析[14,15]，bwa mem 使用参数：-M-R-R"@RG；ID：wHAIPI016410-33；PL：Illumina；PU：150211_I19 1_FCC6L7EANXX_L5_wHAIPI016410-33；LB：wHAI PI016410-33；SM：TEST1；CN：BGI"all.chrs.con.faread1.fq.gz read2.fq.gz ＞aligned_reads.sam。统计比对完成后，应用samtools eference{samtools1}软件[16,17]根据基因组上的坐标对SAM 文件进行排序，并将其转化成BAM 文件。

1.3.3 SNP 检测

应用GATK 检测SNPs，应用BGI 自主开发软件对其进行注释和统计[18,19]，SNP 过滤参数为："QD＜2.0||FS＞60.0||MQ＜40.0||MQRankSum＜-12.5||Read-PosRankSum＜-8.0"。SNP 检测的过程如下：（1）将所有样本的比对文件bam 放到一个list 文件里面；（2）使用GATK 得到所有样本的SNP；（3）使用GATK过滤得到的变异结果，选取可靠的变异结果。

1.3.4 群体遗传学分析测定

基于各样品的SNP 使用p-distance 方法计算各样品的遗传距离矩阵[20]，用NJ 法构建群体间的系统发育树。应用Admixture 软件基于最大似然法对66 份棘胸蛙样品进行主成分分析。

2 结果与分析

2.1 总DNA 的检测

用1%琼脂糖凝胶电泳检测棘胸蛙样本的总DNA（图1）。结果显示：基因组DNA 条带单一清晰可见，无拖尾，表明基因组DNA 提取质量良好，无明显降解、无RNA 和蛋白质污染，能够满足后期建库测序要求。

2.2 数据分析

通过Illumina HiSeq 测序平台，下机原始数据（Raw data）按照探针拆分、过滤后获得有效数据（Clean data）共66.79 G。三区市66 个棘胸蛙样本测序数据统计见表1。由表1 计算得出，66 个样品的平均原始数据约为6.98 M，平均有效数据约为3.69 M，且Q20 均在97%以上，Q30 均高于93%，GC 含量平均为46.70%，表明测序数据质量值高，能满足后续信息分析。

2.3 比对与SNP 检测

以非洲爪蟾为参考基因组，应用BWA 软件将所有的测序数据比对到参考基因组上。其基因组大小为2 734 636 505 bp，有效基因组大小为2 408 724 787 bp（参考序列中不含N），参考基因组GC 含量为34.34%。所有样本的比对率介于58.40%和60.76%之间。

对3 个地区棘胸蛙所获得的SNP 进行统计汇总（表2），为确保鉴定得到的SNP 的质量，对检测到的SNP 进行过滤：（1）异置信度与质量深度的比值≥2；（2）使用Fisher's 精确检验来检测链偏倚现象而得到的Fhred 格式的p 值，本试验设置为≤60；（3）均方根值比对质量值≥40；（4）变异非参检验、秩和检验的零-均值（Zscore）和比对序列质量值的比值≥-12.5；（5）变异非参检验、秩和检验的零-均值（Zscore）和比对序列位置偏差的比值≥-8.0。最终经过过滤，3 个地区的棘胸蛙样本共获得396 969个SNP 位点，其中龙岩群体获得SNP 位点数最多，为139 169 个；三明群体最少，为123 437 个SNP 位点。经过再次筛选检测，3 个地区的棘胸蛙样品共获得8 615 个共有SNP 位点用于高级分析。

表3 3 个地区的棘胸蛙样本基因分型结果统计Tab.3 Test result of genotypes of all samples spine frog Rana spinosa collected from three areas

由表3 可知，南平棘胸蛙群体（Rs-N）纯合基因型的占比平均为99.77%，杂合基因型的占比平均为0.23%；三明棘胸蛙群体（Rs-S）纯合基因型的占比平均为99.66%，杂合基因型的占比平均为0.34%；龙岩棘胸蛙群体（Rs-L）纯合基因型的占比平均为99.65%，杂合基因型的占比平均为0.35%。3 个群体的纯合基因型的占比都在99%以上，表明3 个群体种质都比较纯，遗传多样性不够丰富；南平群体的纯合性占比最高，杂合性最低，表明南平群体种质更加纯合，受到外源种质的污染更低，种质保存更好；龙岩群体杂合基因型占比较高，说明其群体的种质不够纯，但遗传多样性相对丰富。

2.4 群体分化分析

为探究3 个棘胸蛙群体之间的遗传分化程度，基于各样品的SNP 使用p-distance 方法计算各样品的遗传距离矩阵，使用TreeBest 软件中的邻接算法，经过1 000 次迭代绘制系统NJ 发育树（图2）。

由图2 可知：南平（Rs-N）23 个棘胸蛙群体中绝大部分个体聚成1 支，可以很好地与三明群体（Rs-S）和龙岩群体（Rs-L）区分开来，说明与这两个群体的进化差异显著，群体遗传距离近，种质比较纯，更适合作为种质选育材料；三明群体和龙岩群体个体间都不能完全聚成1 支，表现为群体中的个体相互掺杂聚成多个分支，呈现无规律分散，说明这两个地区的棘胸蛙群体各自进化差异不显著，群体种质不如南平地区的种纯。

2.5 群体主成分分析

根据个体水平的遗传距离对3 个棘胸蛙群体进行主成分分析（图3）。由图3 可知，3 个地区的棘胸蛙群体区分度不明显，群体遗传差异不显著。

3 讨论

群体遗传学及其分子系统进化等研究的基础是遗传标记技术。从20 世纪80 年代至今，DNA分子标记技术快速发展，由最初的限制性片段长度多态性（Restriction fragment length polymorphism，RFLP）[21]到随机性扩增DNA（Random amplified polymorphic DNA，RAPD）[22,23]，到现在应用十分广泛的微卫星DNA（Microsatellite DNA）[24]以及单核苷酸多态性（Single nucleotide polymorphism，SNP）[25]。这些传统的分子标记方法工作量大，比较繁琐和耗时。近些年，随着高通量测序技术（Next generation sequencing，NGS）的飞速发展，简化基因组测序（Reduced-representation sequencing）技术避免了上述方法的不足，开始逐渐被应用于非模式生物的研究中[26-28]，成为一种十分有前景的分子标记技术。

棘胸蛙作为一种经济价值较高的易危蛙类，逐渐被人们所认识。国内外关于棘胸蛙的研究主要集中在分类与分布、形态学、繁殖与发育生物学、细胞遗传学以及人工养殖技术等方面，而从分子生物学角度探讨棘胸蛙种群遗传结构和系统发生的研究还不多，并且多停留在线粒体DNA 标记技术[29,30]以及传统的分子遗传标记技术[31-33]。本研究突破传统分子标记技术，首次将GBS 技术应用于棘胸蛙遗传多样性分析，在3 个地区的棘胸蛙群体中开发获得8 615 个高质量SNP，其中以龙岩群体的SNP 最多，为3 589 个，南平群体的SNP 最少，为2 153 个，说明不同地区棘胸蛙群体SNP 数存在较大差异，这可能与其栖息的生态环境或人工干预不同有关；在同一群体不同个体间SNP 也不相同，这为今后开展品种选育提供更多的遗传信息。

本研究中的群体进化树表明，南平棘胸蛙群体主要聚在一起，与其地理分布一致，不存在与三明和龙岩个体聚类时混合交叉的现象，能很好地与三明和龙岩群体区分开。三明和龙岩群体在聚类时，个体混合交叉聚类现象明显，不能很好地相互区分开。PCA 主成分分析表明，3 个群体的区分度均不明显，虽然在聚类时呈现一定程度的分化，但未形成一定的地理隔离现象。

总之，本研究首次采用GBS 简化基因组测序技术，揭示了福建省南平市、三明市和龙岩市3 个地区棘胸蛙群体的遗传多样性与遗传结构关系。南平种群的遗传多样性已表现出杂合度低、遗传多样性下降等趋势。系统发生树以及PCA 主成分分析发现，南平种群与其他2 个种群的种群分化表现出一定程度的遗传分化，但还未形成地理隔离。此次研究结果为制定有效的棘胸蛙种质保护措施，更好地保护与利用棘胸蛙种质资源提供了参考。但由于本次试验受采样困难、覆盖度不高等客观原因限制，研究结果还不够深入，下一步将提高采样密度，结合核基因标记技术，进行更深入、精细的研究分析，以期获得更为详细的棘胸蛙群体遗传信息资料。