基于线粒体COⅠ的南海北部长棘银鲈遗传多样性分析

裴丽梅，宫亚运，黄柏炎，章 群，徐 示，王业磷

（暨南大学生态学系，广州 510632）

长棘银鲈（Gerres filamentosus）隶属于鲈形目（Perciformes）银鲈科（Gerreidae）银鲈属，广泛分布于世界热带和亚热带海域和部分河口区［1-2］，一般生活在浅海沿岸50 m以深的沙质底部，以生活在沙子和泥浆底部的小型生物为食［3-4］。长棘银鲈为小型优质鱼类 ，我国见于南海和东海南部水域，历史上北部湾曾盛产，但近年来由于过度捕捞，野生资源已经严重衰竭，资源密度由1961—1962年的 78.4 kg·km-2下降至 1992—1993年以后的 0.1 kg·km-2［5］。随着长棘银鲈等鱼类资源的衰竭，北部湾海洋生态系统发生了明显变化，如游泳动物和底栖生物食性鱼类种类数日益减少，而浮游生物食性鱼类种类数及其渔获比例明显增加［6］。虽然自1999年以来实行的休渔制度对北部湾部分鱼类资源恢复有明显效果，但是对衰竭型种类效果甚微［5］，长棘银鲈种质资源保护、恢复迫在眉睫。

了解种群的遗传变异是评估资源现状、制定野生群体保护与恢复策略和实施养殖群体种质遗传改良的重要前提［7］。目前国内外对于长棘银鲈的研究大多集中在形态比较、物种分类、生长发育和线粒体全基因测序等方面［7-13］，尚未有长棘银鲈种群遗传相关的研究报道。线粒体DNA（mtDNA）具有进化速度快、母系遗传、几乎不发生重组等特点［14-15］，已被用作脊椎动物的系统发育和群体遗传分析的有效工具［16-17］。其中细胞色素氧化酶Ⅰ（cytochrome oxidase subunitⅠ，COⅠ）序列进化速率适中，是开展动物DNA条形码研究首选分子标记［18］，用于种群遗传学研究时，也可在一定程度上避免因线粒体控制区易变异饱和、有时不能有效检测到种群结构以及学术上对线粒体控制区序列插入和缺失状况处理的共识不足［19-20］。近年来国内外有不少使用COⅠ基因研究鱼类遗传变异的报道，如日本鲭（Scomber japonicus）［21］、 银 鲳 （Pampus argenteus）［22］和乌贼（Sepia officinalis）［23］等。本研究采用线粒体COⅠ基因序列对南海北部7个长棘银鲈群体进行分析，旨在揭示中国长棘银鲈主要分布区的遗传多样性和种群结构，为种质资源保护和开发利用提供科学依据。

1 材料与方法

1.1 材料来源

2014年7月至2019年10月间采自广东碣石（JS）34尾、广东阳江（YJ）34尾、海南琼海（QH）14尾、海南陵水（LS）24尾、海南东方（DF）20尾、海南新盈（XY）26尾、广西东兴（DX）24尾共7个地点的176尾长棘银鲈（表1），用95%的酒精处理并保存。所有样品均为从近海作业的小型渔船上购买的野生鱼类。

1.2 基因组DNA的提取与PCR扩增

DNA提取和PCR扩增条件参照乐小亮等［24］方法，PCR扩增引物参见WARD等［25］，经琼脂糖凝胶电泳检测合格的PCR扩增产物送华大基因公司测序。

1.3 数据处理

将测序后得到的序列在 Bioedit 7.0.9.0［26］软件上进行人工校对和多重比较。运用MEGA 7.0［27］软件计算核苷酸组成、碱基变异情况、转换和颠换值。通过 DnaSP 6.12.03［28］软件计算单倍型数、单倍型多样性（haplotype diversity，Hd）、核苷酸多样性（nucleotide diversity，π）。使用Network 5.0［29］构 建 单 倍 型 网 络 图。通 过Arlequin 3.5.2.2［30］软件计算遗传分化系数（Fst）和群体间基因流值（Nm）；进行AMOVA分子方差分析遗传变异的分布情况，并以1 000次置换分析统计学显著性；计算Tajima’sD和Fu’sFs值并进行中性检验，得到SSD值、Rg值、τ值参数。根据公式T＝（τ／2μk）×代时［31］计算种群扩张时间。其中，τ是种群扩张时间参数；μ是序列的变异速率 ［采用 1% ～3%／MY［32］（million year，MY）］；k是序列长度；T是种群扩张时间。

2 结果与分析

2.1 COⅠ基因片段序列特征及遗传多样性

在长棘银鲈COⅠ基因5′端652 bp序列中，T、C、A和 G含量分别为 29.0%、28.7%、23.3%和19.1%。A+T的含量（52.3%）略高于C+G的含量（47.8%），与大多数硬骨鱼类相似［33］。共发现单倍型25个，变异位点155个，其中简约信息位点87个。转换和颠换比为2.41，表明序列未达到饱和，适合进行种群遗传发育分析［34］。7个群体整体呈现较高的单倍型多样性（Hd＝0.767±0.018）与较低的核苷酸多样性（π＝0.0550±0.003 7）。其中遗传多样性指数最高的是碣石群体（Hd＝0.606±0.092，π＝0.059 2±0.006 6），其次为阳江群体（Hd＝0.323±0.102，π＝0.027 0±0.010 7）；琼海群体遗传多样性指数最低（Hd＝0.560±0.125，π＝0.001 0±0.000 3），陵水群体次低（Hd＝0.685±0.077，π＝0.0015±0.000 3）。

2.2 长棘银鲈遗传分化

176条COⅠ序列中共检测到25个单倍型，单倍型频率最高的为 Hap1（60条），其次是Hap15（50条）。共享单倍型4个，其中 Hap1、Hap2和Hap3为所有群体共享；Hap15为碣石、阳江、新盈3个群体共享；其它单倍型为各群体特有。

用 Arlequin 3.5.2.2软件以 Kimura双参数模型构建的单倍型网络中呈现非典型的星状结构，没有较为明显的地理聚群和谱系结构（图1）。7个群体间的遗传分化情况（Fst＝-0.025 5～0.851 4，Nm＝-68.993～56.971）详见表 2。其中，碣石和阳江群体分别与琼海、陵水、东方、新盈、东兴5个群体间出现极显著性遗传分化（Fst＞0.25，P＜0.001）［35］，且碣石和阳江群体分别与琼海、陵水、东方、新盈、东兴5个群体之间基因交流匮乏（0＜Nm＜1）［36-37］。碣石和阳江群体间虽然为随机交配（Nm＞4），但也出现了低等程度分化（0.05＜Fst＝0.111 1＜0.15，P＜0.05），其余组别之间无明显分化。

将7个群体通过4种不同方式进行AMOVA分析：①7个地理群体；②大陆群组：碣石、阳江、东兴；海岛群组：琼海、陵水、东方、新盈；③碣石和阳江为一组；琼海、陵水、东方、新盈、东兴为另一组；④碣石群体；阳江群体；其他群体（琼海、陵水、东方、新盈）。AMOVA分析（表3）表明4种分组方式的主要变异来源均发生在群体内，变异比例分别为 74.13%、68.87%、61.69% 和64.76%，同时4种分组方法下4种群组间存在显著性分化（Fst＞0.25，P＜0.001），原因可能与长棘银鲈只进行短距离洄游有关［5］。

表1 南海长棘银鲈群体的遗传多样性Tab.1 Genetic diversity of G.filamentosus in coastal waters of the South China Sea

图1 长棘银鲈单倍型网络图Fig.1 Parsimony network of G.filamentosus hap lotypes

2.3 长棘银鲈种群历史动态

中性检验结果见表4。除陵水 Tajima’sD（-1.518 6，0.001＜P＜0.05）和 Fu’sFs（-3.293 6，P＜0.05）、东 方 Tajima’sD（-2.590 6，P＜0.001）、新 盈 Tajima’sD（-2.650 4，P＜0.001）为显著性负值外，其余群体Tajima’sD和Fu’sFs均为不显著性负值或正值。结合陵水、东方、新盈群体SSD值及Rg值数值均较小且不具有统计显著性，表明上述观测值符合扩张模型。因Fu’sFs检验比Tajima’sD对近期扩张更敏感［38］，再结合上述3个群体的核苷酸错配图（图2-a～c）为明显单峰以及单倍型网络图为非典型星型结构，可认为陵水群体经历了近期的快速扩张，而东方和新盈群体经历过较为早期的扩张［39-40］。7个群体整体的 Tajima’sD＝-0.765 9（P＞0.1）、Fu’sFs＝6.849 1（P＞0.5），整体上未呈现扩张现象，结合核苷酸错配图（图2-d）呈现明显多峰，说明可能是种群经历多次扩张叠加的结果［40］。陵水、东方、新盈3个群体扩张参数 τ分别为 1.068 4、0.949 2和1.175 8，以线粒体COⅠ基因进化速率为（1%～3%）／百万年［32］，性成熟期 3年，估算出陵水、东方、新盈3个群体扩张大约分别发生在0.245 8～0.081 9、0.218 4～0.072 8、0.270 5～0.090 2百万年之前，均处于晚更新世时期。

表2 长棘银鲈7群体的F st及N m值Tab.2 F st and N m among 7 G.filamentosus populations

表3 长棘银鲈种群结构分子方差分析Tab.3 AMOVA analysis of G.filamentosus

表4 长棘银鲈中性检验和不对称分析Tab.4 Neutrality tests and m ismatch analysis of G.filamentosus

图2 长棘银鲈核苷酸不对称分布Fig.2 Nucleotidem ismatch distribution of G.filamentosus

3 讨论

3.1 长棘银鲈的遗传多样性与种群历史动态

种群遗传多样性的高低很大程度上与物种的适应能力、分化能力以及对环境影响的耐受力相关［41］，而单倍型多样性（Hd）和核酸多样性（π）是衡量某一物种遗传多样性的2个重要指标［42］。长棘银鲈COⅠ序列的遗传多样性（Hd＝0.767±0.018，π＝0.055 0±0.003 7）呈现的高Hd值低π值的特点，与同海域基于COⅠ序列的银鲳 （Pampus argenteus）［22］和 湖 栖 虾 虎 鱼（Gobiopterus lacustris）［43］等相似。南海北部［碣石（0.059 2±0.006 6）和阳江（0.027 0±0.010 7）］的遗传多样性明显高于南海南部［琼海（0.001 0±0.000 3）和陵水（0.001 5±0.000 3）］，推测可能与碣石和阳江所处的海域有关。在冰期大陆架露出时，南海变成边缘海，阳江和碣石群体濒临该边缘海，种群相对稳定，受瓶颈效果影响小；另一方面，碣石与阳江所处海域会同时受到南海暖流、黑潮及沿岸流的影响，使得碣石与阳江群体具有更强的分化能力和环境耐受力。而琼海、陵水、东方、新盈的潮汐汊道港湾特征使内湾水域处于比较稳定的半封闭海洋环境［43］；东兴则位于半封闭的北部湾，海洋环境也相对稳定，遗传多样性偏低。此外，冰期时海南岛及北部湾所在的大陆架露出，间冰期时海平面上升［44-45］，海南岛和北部湾的长棘银鲈群体发生扩张，也致使琼海、陵水、东方、新盈和东兴群体遗传多样性偏低。

Fst值显示阳江和碣石群体分别与琼海、陵水、东方、新盈、东兴5个群体间发生极显著性遗传分化（Fst＞0.25，P＜0.001），且 AMOVA分析结果显示长棘银鲈4种划分方式的群组间都有显著性分化。长棘银鲈在河口生活，繁殖期迁移至沿海地区达到成熟的最后阶段，然后进行产卵，不进行远距离洄游［3］。琼州海峡是在全新世中期全球性海侵淹没峡区原来的常态低地而形成，主要形成时段为（10 570±560）～（7 125±96）aBP［46-47］。琼州海峡形成前，黑潮从巴士海峡进入南海，同沿岸流一样在沿广东海岸流向西南过程中受阻，导致碣石和阳江群体与海南各群体及东兴群体间发生遗传分化；琼州海峡形成后，黑潮从巴士海峡进入南海，同沿岸流一样在沿广东海岸流向西南过程中受琼州海峡和南海大陆架特殊地况的影响，也在一定程度上阻碍了碣石和阳江群体的远距离传播，导致碣石和阳江群体与海南各群体及东兴群体间有明显的遗传分化。

大部分中国近海鱼类在近期历史上因冰期海平面变化而出现种群快速扩张［45］。长棘银鲈Tajima’sD和Fu’sFs中性检验和核苷酸不对称分析图呈现单峰，显示陵水、东方和新盈群体都经历过种群扩张，推算扩张时间均处于晚更新世时期。海南岛为南海大陆架岛屿，在晚更新世构造运动频繁，火山活动活跃，使不少海湾底部抬升为阶地；全新世气候变暖，海平面上升，海水侵入河口、海湾、洼地等，使长棘银鲈群体随着海侵发生扩张［44］。

3.2 长棘银鲈种质资源保护

长棘银鲈生活在浅海沿岸水域，是包括珊瑚礁、红树林等生境生态系统的重要组分［2-3］。近年来捕捞技术的发展与革新使得长棘银鲈资源快速衰退［5-6］，长棘银鲈野生资源的保护迫在眉睫。碣石和阳江群体分别与琼海、陵水、东方、新盈群体间有极显著性分化，可以将碣石和阳江以及其他群体作为2个独立的管理保护单位。其中碣石群体遗传多样性最高（Hd＝0.606±0.092，π＝0.059 2±0.006 6），应作为长棘银鲈种质资源保护和开发利用的核心群体和优先对象；琼海（Hd＝0.560±0.125，π＝0.001 0±0.000 3）和陵水（Hd＝0.685±0.077，π＝0.001 5±0.000 3）群体遗传多样性较低，更容易受到环境变化的影响，建议及时采取措施，避免捕捞过度使其遗传多样性难以恢复，重蹈大黄鱼、带鱼等种质资源衰退的覆辙。

随着分子生物学技术的进步，特别是高通量测序技术的进步和测序成本的下降，种群遗传学研究开始进入基因组时代。线粒体DNA有效群体大小不仅为双亲遗传核基因组的1／4，可提供核基因组不能提供的母系遗传信息，开展种群遗传研究的成本远低于种群全基因组分析，仍是种群遗传的首选分子标记［48］。基于线粒体COⅠ的长棘银鲈种群遗传研究发现，南海北部琼州海峡两侧的群体间存在极显著的遗传分化，可作为2个不同的管理保护单位；广东沿海群体的遗传多样性明显高于海南和广西群体，其中碣石群体遗传多样性最高，可优先保护。研究结果为长棘银鲈种群的保护和开发利用措施的制订与实施提供了科学依据。但受采样范围的限制，研究结果不能排除未采样地区存在种群分化的可能，今后需要进行物种全分布区的采样调查。另外，相对于其他线粒体基因，COⅠ基因的进化速度较慢，较适合发现历史相对久远的种群扩张，对近期的种群扩张不太敏感，在后续有关长棘银鲈种群遗传的研究中，可考虑使用高多态性的微卫星指纹图谱分析技术以及基于二代测序的简化基因组乃至全基因组分析技术，以揭示双亲遗传的核基因遗传信息，反映更为近期的全球气候变化与人类活动的影响，从而制订更为精准的保护与开发利用策略。