青藏高原东北部黄缨菊的谱系地理学

张 阳 马子兰 徐珊珊 苏 旭* 李梅英

（1. 上海科技馆自然史研究中心，上海 200041；2. 上海市静安区科学技术协会，上海 200040；3. 青海师范大学生命科学学院，西宁 810008）

谱系地理学（Phylogeography）主要利用分子标记探讨物种的分布格局和进化历史与地质事件之间的相互关系，进一步推测物种现代分布格局的历史成因。全球气候变化尤其是第四纪冰期和间冰期气候波动，导致反复而剧烈的环境变化，这些变化影响着全球动植物的地理分布和遗传结构。青藏高原（Qinghai-Tibet Plateau）被誉为“世界屋脊”和“地球第三极”，其复杂的地质结构和独特的生境特征，使该地区成为响应全球气候与环境变化最强烈的地区。近年来，诸多学者通过对青藏高原及毗邻地区多种高山植物谱系地理学的研究，揭示了高山植物响应第四纪冰期气候波动的分布格局和进化历史，认为它们主要存在两种“避难”模式：一种模式是某些高山植物冰期或间冰期高原台面上存在多个相互隔离的微型避难所，冰期后向周边扩散，这类植物有西藏沙棘（）、西川红景天（）、露蕊乌头（）等；另一种模式认为青海云杉（）、长花马先蒿（）、条纹狭蕊龙胆（）等高山植物冰期或间冰期从高原台面退缩至边缘避难所，冰期后重新扩散至高原台面。

研究表明，被子植物的叶绿体DNA 为母系遗传，具有基因不重组、引物通用、进化速率较慢等优点，目前已成为植物谱系地理学研究的最佳分子标记，尤其适用于单个物种的群体遗传结构、遗传多样性和进化历史研究。黄缨菊（C.Winkl.）是青藏高原地区一种特有的高寒草甸药用植物，隶属于菊科（Asteraceae）、黄缨菊属（），适生于海拔2 230～4 150 m的草甸、草原及干燥山坡，主要分布在中国青海、西藏、甘肃东南部、四川及云南西北部等地。黄缨菊具有较强的抗寒、抗旱、抗盐碱等特性，对维持青藏高原生态系统的稳定性具有重要作用；同时，该物种全草入药，对吐血、子宫出血、十二指肠溃疡出血、胃溃疡及过敏性紫癜等疗效明显，有较高的药用价值。迄今为止，有关黄缨菊的研究主要集中于化学成分、抗逆生理特性和分子系统学方面。譬如，Wang等和Barres等基于核糖体DNA 内转录间隔区（nrITS）和叶绿体基因（cpDNA）分别对黄缨菊属及其近缘属的系统亲缘关系进行研究，发现黄缨菊属与蝟菊属（）、肋果蓟属（）和山牛蒡属（）等具有较近的亲缘关系，而与飞廉属（）的亲缘关系相对较远。然而，关于黄缨菊谱系地理学的的研究鲜见报道。据此，本研究利用、-和-3 个cpDNA 片段对黄缨菊进行谱系地理学研究，揭示黄缨菊的群体遗传结构、现代地理分布格局以及冰期时的避难所，旨在为今后黄缨菊物种形成机制的研究提供理论依据。

1 材料与方法

1.1材料

在青海东部和甘肃北部，共收集到黄缨菊20个居群、123 个个体（见表1 和图1）。按照群体遗传学采样方法对黄缨菊进行野外调查和采样，即居群之间直线距离至少20 km，每个居群随机采集3～13 个个体，个体之间间隔20 m 以上。野外采集生长良好且健康的黄缨菊新鲜嫩叶，迅速用硅胶干燥带回实验室，并记录各居群的经纬度、海拔、生境、伴生种等信息。凭证标本存放于中国科学院西北高原生物研究所青藏高原生物标本馆（HNWP）。

图1 黄缨菊居群采样图及cpDNA单倍型分布饼状图大小代表每个居群单倍型的频率（审图号：GS（2019）1837号）Fig.1 The sampled location of all populations and distribution of cpDNA haplotypes of each population in X.subacaulisThe size of pie charts represents the frequency of haplotypes in each population（Map number：GS（2019）1837）

表1 黄缨菊cpDNA片段6种单倍型序列的变异位点Table 1 Variable sites of six haplotype sequences base on cpDNA fragments in X.subacaulis

1.2 DNA 提取和PCR 扩增

用植物基因组DNA 提取试剂盒（TIANGEN）提取硅胶干燥的叶片总DNA。利用-、和-3 个cpDNA 片段 对黄 缨菊 所有个体进行PCR 扩增。PCR 扩增体系为25.00 μL，包含17.25 μL 的HO、2.50 μL 的10×PCR buf‑fer（1.5 mmol·LMg）、2.00 μL 的dNTPs 混 合 液（10 mmol·L）、正反引物各1.00 μL（5 μmol·L）、0.25 μL 的DNA 聚合酶（5 U·μL）、DNA 模板1.00 μL（约100 ng）；反应程序如下：94 ℃预变性5 min；94 ℃变性50 s，55～59 ℃退火1 min，72 ℃延伸1 min，循环36 次；72 ℃延伸10 min。PCR 扩增产物经1%琼脂糖凝胶电泳检测后送至北京六合华大基因科技有限公司进行测序。

1.3数据分析

利用MEGA 6.0 软件对黄缨菊所有测序个体的序列进行比对和手工校正，然后将cpDNA 的3 个片段（-、和-）串联。采用DnaSP 5.0 软件统计黄缨菊单倍型数目和变异位点；应用PERMUT 2.0软件计算总的遗传多样性（）、居群内平均遗传多样性（），以及居群遗传分化系数和，对和值进行1 000 次重复置换检验；运用ARLEQUIN 3.5 软件计算每个居群的遗传多样性（）、核苷酸多样性（）及单倍型组成和频率，并利用分子变异（AMOVA）分析方法检测黄缨菊居群内和居群间的遗传变异水平及分化系数（），同时依据公式=（1-）/2计算物种水平上居群间平均基因流；利用Ge‑nAlEx 6.5中的Mantel 统计学检验比较黄缨菊居群平均遗传距离与地理距离之间的相关性。

以葵花大蓟（）为外类群，利用NETWORK 10.2软件中的Median-joining 算法构建单倍型网络图；采用PhyloSuite v1.2.2软件，选用最佳模型F81，利用最大似然法（Maximum likeli‑hood，ML）和贝叶斯法（Bayesian inference，BI）构建黄缨菊单倍型系统发育关系。ML 树利用启发式搜索法运行5 000次自展重复（Bootstrap）检验各分支的置信度，BI 树利用马尔科夫链蒙特卡洛方法（Markov chain Monte Carlo process，MCMC）运 行200 万代获得各个分支的后验概率（Posterior prob‑ability，PP），其余为默认参数；使用DnaSP 5.0软件计算中性检验值Tajima’s和Fu and Li’s，并进行错配分布分析。

2 结果与分析

2.1序列变异及遗传多样性

通过对黄缨菊所有个体cpDNA 片段、和-分析，获得序列长度为2 074 bp，（G+C）和（A+T）含量分别为37.6%、62.4%；cpDNA 片段含有26个变异位点，其中有16个碱基插入和缺失、8 个碱基颠换和2 个碱基转换（见表1）。基于叶绿体DNA 片段的单倍型分析显示，黄缨菊居群共检测到6个单倍型（H1～H6），其中单倍型H1的出现频率最高，单倍型H3、H5和H6为特有单倍型（见图1）；将-、和-片段的单倍型序列分别提交至GenBank（MZ274274-274277、MZ312254-312255 和MZ312260-312263）。黄缨菊居群总的遗传多样性（）为3.101，范围为0～0.533；总的核苷酸多样性（）为0.008 903，范围为0～0.003 525；参试居群中P1、P8～P9、P11～P13和P17 的遗传多样性水平高于其他居群（见表2）。每个居群的、、单倍型组成数量及频率详见表2。

表2 黄缨菊采样信息及各居群遗传多样性Table 2 Sampling information and genetic diversity of X.subacaulis

2.2居群遗传结构

研究结果表明，黄缨菊居群内平均遗传多样性=0.155，总的遗传多样性=0.451，居群间遗传分化系数（0.727）大于（0.656），但不显著（>0.05），说明黄缨菊在现有分布区域内没有明显的谱系地理结构。AMOVA 分析结果显示，黄缨菊居群间遗传变异为68.98%，而居群内遗传变异为31.02%，=0.689 85（<0.01），表明黄缨菊的遗传变异主要存在于居群间，且遗传分化水平较高（见表3）。同时，黄缨菊居群间的基因流=0.025，显示居群间基因流较小。此外，基于Mantel检验表明黄缨菊居群间平均遗传距离与地理距离显著正相关（=0.214，=0.01），即居群间地理距离越大其遗传分化明显（见图2）。

表3 黄缨菊居群单倍型的AMOVA分析Table 3 Analyses of molecular variance（AMOVA）of haplotypes from populations in X.subacaulis

图2 黄缨菊居群地理距离与遗传距离的Mantel检验Fig.2 Mantel test between geographic and genetic distance of populations in X.subacaulis

2.3居群历史动态

黄缨菊错配分布分析显示，其错配分布曲线为单峰曲线（见图3）；同时，Sum of squared devia‑tion（SSD）值 为0.038（=0.109）和Harpending’s raggedness index值为0.167（=0.247），均不显著（>0.05），符合种群扩张模型假设，说明黄缨菊在现有分布区域内经历过扩张。然而，中性检验结果发现，Tajima’s=-1.070 35、Fu and Li’s=0.634 08，均不显著（>0.05），表明黄缨菊居群没有经历过扩张现象，这与错配分布分析的结果不同。此外，AMOVA 分析显示，黄缨菊种群暴发前的种群大小=0，群暴发后的种群大小=99 999，有效种群大小（effective population size）的变化范围（～=0～999 99）非常高，认为黄缨菊居群近期经历过扩张。

图3 黄缨菊居群cpDNA序列的错配分布分析Fig.3 Mismatch distribution of cpDNA sequences from X.subacaulis

2.4单倍型系统发育分析

以葵花大蓟（MZ312258、MZ312256、MZ312264）和 飞 廉（，MZ312259、MZ312257、MZ312265）为外类群，利用PhyloSuite 软件基于最大似然法（ML）和贝叶斯法（BI）构建黄缨菊单倍型系统发育树，结果显示黄缨菊6 个单倍型聚为3支（Clade Ⅰ、Ⅱ、Ⅲ），支持率达到80%以上。其中，单倍型H1、H2、H6 聚为一支（Clade Ⅰ），位于系统发育树顶端，亲缘关系较近；特有单倍型H3、H5 聚为一支（Clade Ⅱ），位于系统发育树中部，具有较近的系统发育关系；单倍型H4 位于发育树的基部（Clade Ⅲ），独自构成一支（见图4A）。同时，基于NETWORK 软件构建了黄缨菊的单倍型网络图，发现单倍型H1 位于网络图中央，占比最大（67%）。其中，单倍型H1 与H3、H4、H5 间以较多的突变步和“隐单倍型”（目前未检测到的单倍型或已灭绝的古老单倍型）相连，表明本研究检测到的单倍型应该具有的独立发育历史（如单倍型H3、H4、H5），并且单倍型之间的网络关系与系统发育树显示的亲缘关系一致（见图4B）。

图4 黄缨菊单倍型的系统发育树（A）和中央网络连接图（B）A.分支上的数字表示BI和ML的支持率，支持率小于80%未标注；B.圆圈大小代表单倍型频率，实心方块代表“隐单倍型”Fig.4 Phylogenetic tree（A）and median-joining network（B）of haplotypes in X.subacaulisA.Numbers above branches are bootstrap values of BI and ML，and the support rate is less than 80% without marking；B.The size of circles repre‑sents the haplotype frequency，and parallelogram represents missing haplotypes

3 讨论

3.1居群遗传结构分析

基于叶绿体DNA 群体遗传学数据，先前学者通过对青藏高原地区多种植物遗传结构的研究，表明该地区物种具有高的和低的，如西藏红豆杉（；0.884，0.271）、西藏沙棘（0.956，0.372）、长花马先蒿（0.770，0.332）、露蕊乌头（0.739，0.207）、菊叶红景天（；0.970，0.424）、偏花报春（；0.966，0.178）等。本研究结果显示，黄缨菊居群也具有高的和低的（=0.451，=0.155），表明黄缨菊居群间叶绿体基因组分化明显，进一步验证了先前研究结果的正确性和合理性。

同时，AMOVA 分析表明，黄缨菊68.98%的遗传变异来源于居群间，而居群内遗传变异仅为31.02%，=0.689 85，=0.225，说明黄缨菊居群间遗传分化程度较高。笔者认为这是由于黄缨菊cpDNA 为母系遗传，种子是其繁殖的主要方式，而该物种为垫状草本植物，叶边具刺，种子无利于传播的附属结构，加之黄缨菊主要生长于青藏高原地区地形复杂的草甸、草原及干燥山坡，使居群间基因交流受阻，主要发生于居群内；并且Mantel 检验的结果也进一步支持这一结论。然而，通过分析黄缨菊群体遗传分化系数和，发现>（>0.05），该物种在分布范围内不存在明显的谱系地理结构，这与其他许多高山植物如菊叶红景天（=0.655，=0.652）、东俄洛橐吾（；=0.614，=0.595）、银露梅（；=0.864，=0.720）等的谱系地理结构一致。

3.2冰期避难所及居群进化历史推测

第四纪冰期气候动荡对植物地理分布产生重大影响，尤其对青藏高原植物类群的群体遗传结构和进化历史影响更为剧烈。于海彬等通过对青藏高原地区36种高山植物谱系地理分布格局的综合分析，提出该地区植物主要具有“退缩—回迁”和“微型避难所”两种模式。他们认为，符合“退缩—回迁”模式的植物第四纪冰期时常退缩到高原边缘避难所，冰期后回迁到高原台面，即遗传多样性较高或古老单倍型较多的居群通常聚集于高原边缘避难所，高原边缘居群的遗传多样性较高，而高原台面居群的遗传多样性较低；而符合“微型避难所”模式的植物第四纪冰期时，并非所有居群都向高原边缘迁移，有部分居群可在本地适宜的微环境中生存，从而在高原台面形成一个或多个不连续的微型避难所，冰期后从高原台面原生避难所向外扩张，即高原台面避难所居群的遗传多样性较高。

溯祖理论认为，位于单倍型网络图中心的单倍型通常为古老单倍型，因而本研究中的单倍型H1 应为黄缨菊居群的古老单倍型，因此拥有单倍型H1 的居群所在的区域可以推断为冰期避难所；同时，本研究结果显示分布于青海湖流域的黄缨菊居群遗传多样性和核苷酸多样性偏高，而其他区域的遗传多样性和核苷酸多样性较低。据此，笔者认为青海湖流域很可能是黄缨菊的冰期避难所。一方面，青海湖流域四面环山的有利地形为其成为黄缨菊第四纪冰期的避难所提供了地理条件；另一方面，青海湖流域充足的光照和水分为黄缨菊生存提供必备的环境条件。同时，本研究还发现甘肃南部的黄缨菊居群P18 也具有丰富的遗传多样性和单倍型多样性，认为该区域可能是黄缨菊第四纪冰期的另一个避难所。这与先前诸多高山植物如露蕊乌头、西川红景天等的第四纪冰期避难所的模式相同。此外，黄缨菊错配分布曲线为单峰曲线，偏离总方差（SSD）和粗糙指数（）不显著，有效种群大小为0～999 99，均支持黄缨菊居群经历过近期扩张。

总之，本研究认为冰期气候动荡和青藏高原隆升共同塑造了黄缨菊的现代地理分布格局，这为今后系统研究黄缨菊的进化历史和物种形成机制奠定了基础。