基于叶绿体和核基因片段序列的甘蒙柽柳谱系地理研究*

温月仙 甘红豪 史胜青 江泽平 吴利禄 褚建民

(1. 中国林业科学研究院林业研究所 国家林业和草原局林木培育重点实验室 北京 100091； 2. 中国林业科学研究院森林生态环境与保护研究所 北京 100091)

物种地理分布区域的变化会影响种群的遗传结构(Ricklefs, 1987； 2010; Herzogetal., 2010)并推动物种的历史演化进程(Petitetal., 1997; Abbottetal., 2000; Hewitt, 2000)。谱系地理学(phylogeography)，是通过重建基因谱系，追溯近缘物种或种群间等位基因的演化关系，科学地解释物种形成现有地理分布格局的成因(Aviseetal., 1987； Avise， 2010)。研究表明，与造山运动相关的地质历史事件对现存植物的分布格局、遗传结构产生了重大影响(于海彬等， 2013； Wenetal., 2014)。例如由于青藏高原的隆升，该地区很多植物类群经历过快速辐射分化(Eatonetal., 2012; 李国栋， 2013； Zhangetal., 2014; Ebersbachetal., 2017)，高大的山体阻碍了居群间的基因交流，遗传分化加剧(张茜， 2008； Xuetal., 2010; Liuetal., 2014)。随着喜马拉雅和青藏高原的隆起，中国的水系也发生了一系列的演变(李吉均，1999； 程捷等， 2001)，而水系的变迁势必会对周边植物的分布格局产生重要影响(孙航等， 2011)。通过对滇榄仁(Terminaliafranchetii)、皱叶醉鱼草(Buddlejacrispa)、冬麻豆(Salweeniawardii)等植物的遗传多样性和谱系地理学研究，发现水系的变迁是导致这些植物类群由间断分布变为连续分布，或由连续分布变为现代的间断分布格局的重要原因(孙航等， 2011)。

黄河是我国第二长河，其形成演化与青藏高原的隆升息息相关。先后经青藏运动C幕、昆黄运动和共和运动，形成了现代黄河(李长安等，1998； 李吉均，1999； 李吉均等， 2001)。关于黄河流域谱系地理的研究多集中于动物方面，例如 Zhou等(2016)通过研究黄河流域4个地区的赤眼鳟(Squaliobarbuscurriculus)，发现赤眼鳟种群在全新世发生了突然的空间扩张，推测是由于气候变暖和黄河流域的变迁引起的； Duan等(2009)研究了黄河鲤科(Cyprinidae)Schizothoracine鱼的谱系地理结构，发现其种群遗传结构与黄河的形成紧密相关。然而，黄河的形成演化对植物遗传结构影响的报道相对较少。

甘蒙柽柳(Tamarixaustromongolica)隶属于柽柳科(Tamaricaceae)柽柳属(Tamarix)，性喜水，耐干旱、盐碱和霜冻，是我国黄河流域特有种，广泛分布于青海、甘肃、宁夏、内蒙古、陕西、山西、河南、山东等省(自治区)(张鹏云等，1990； 刘铭庭， 2012； 李法曾等， 2016)，是研究黄河形成演化影响植物遗传结构的优良材料。目前，对甘蒙柽柳的研究多集中于抗盐性评价及其生理机制等方面，通过对甘蒙柽柳、枸杞(Lyciumchinense)、樟子松(Pinussylvestrisvar.mongolica)、杜松(Juniperusrigida)、白榆(Ulmuspumila)、新疆杨(Populusalbavar.pyramidalis)等树种进行盐胁迫试验，综合评价之后发现甘蒙柽柳的耐盐性最强(刘克东， 2012)。对不同盐浓度下甘蒙柽柳和柽柳(Tamarixchinensis)根系生长及生理特性进行研究，发现高盐条件下，二者均通过增加根系生物量分配比例来缓解盐分对根系的危害(宋香静等， 2017)。但是，目前关于甘蒙柽柳的研究主要集中在宏观生理生态领域，有关甘蒙柽柳进化历史的研究则较少。Liang等(2018)利用2个叶绿体片段(trnL-F,rps16)和1个核基因ITS序列，研究了黄河形成对甘蒙柽柳(青海、甘肃、宁夏、内蒙古、陕西、山西)和柽柳(陕西、山西、河南、山东)扩张的影响，发现甘蒙柽柳居群是从黄河上游迁移至下游的。尽管Liang等(2018)的工作对该种的迁徙有很好的解释，但是该研究的采样地点仅到达黄河上游甘肃省和青海省的交界处(民和县)，而缺乏对黄河溯源加长(黄河支流溯源侵蚀进入共和盆地、兴海盆地等)部分样本的分析，无法完整地揭示甘蒙柽柳居群的迁移过程。

本文在对我国甘蒙柽柳的野生资源进行系统调查、收集与保存的基础上，根据甘蒙柽柳的枝、叶等形态特征进行了鉴定，比较全面地采集了黄河流域自然分布的甘蒙柽柳样品； 利用叶绿体基因trnQ-rps16片段和核基因ITS序列研究其遗传多样性及谱系地理结构，探讨甘蒙柽柳居群的扩张历史，推测黄河形成对其现有地理分布格局的影响。

1 材料与方法

1.1 样品采集

通过查阅《中国植物志》(张鹏云等，1990)以及各省植物志、文献资料等，对分布于我国黄河流域的甘蒙柽柳进行实地调查，记录其原生生态资料(如经纬度、海拔等信息)，并采集幼嫩叶片。采样植株间距大于50 m，居群之间距离大于15 km，采集的甘蒙柽柳叶片立即干燥保存。目前，已收集到青海、甘肃、宁夏、内蒙古、陕西、山西、河南、山东共8个省(自治区)17个居群266个甘蒙柽柳单株叶片样品(表1)，甘蒙柽柳的形态特征及部分省区的天然分布见图1。以柽柳科水柏枝属(Myricaria)的宽苞水柏枝(M.bracteata)、柽柳属的沙生柽柳(T.taklamakanensis)和刚毛柽柳(T.hispida)为外类群，样品取自甘肃省敦煌市。

表1 甘蒙柽柳采样信息Tab. 1 The sampling information of Tamarix austromongolica

1.2 DNA提取、扩增、测序及引物筛选

植物总DNA提取采用改良的CTAB法(Doyle, 1991)。通过查阅文献，选取3对叶绿体通用引物和核基因ITS序列(表2)，获取其序列信息并合成引物(北京擎科新业生物技术有限公司)，用于甘蒙柽柳多态性叶绿体基因及核基因序列片段的筛选。PCR反应体系为20 μL： 5×Buffer (Mg2+plus) 4 μL，dNTP 1.6 μL(各2.5 mmol·L-1)，Taq DNA聚合酶 0.5 μL (2.5 U·μL-1)(PrimeSTAR® HS DNA Polymerase，TaKaRa)，Primer 0.4 μL (10 μmol·L-1)，模板DNA 1 μL，ddH2O 12.1 μL。PCR反应程序为： 98 ℃预变性1 min； 然后98 ℃变性10 s，58 ℃退火5 s，72 ℃延伸60 s，共35个循环； 最后72 ℃延伸7 min。利用1%的琼脂糖凝胶电泳筛选出条带清晰的PCR扩增产物，并对其进行双向测序(北京擎科新业生物技术有限公司)，测序仪器为ABI 3730XL DNA 分析仪(Applied Biosystems，Foster，CA，USA)。

图1 甘蒙柽柳的形态特征与天然分布Fig.1 Morphological characteristics and natural distribution of T. austromongolicaA. 花； B. 茎叶； C-F: 甘蒙柽柳在山东昌邑市(C)、河南济源市(D)、内蒙古磴口县(E)、青海同德县(F)的分布。A． Flower; B. Cauline and leaf; C-F: The natural distribution of T. austromongolica in Changyi city in Shandong province(C), Jiyuan city in Henan province(D), Dengkou county, Inner Mongolia Autonomous Region(E), and Tongde county in Qinghai province(F).

表2 引物基本信息Tab. 2 Information of primers

1.3 基因片段的单倍型、遗传多样性、遗传结构及历史动态分析

使用DNASTAR Lasergene(V7.1)软件对原始序列进行正反拼接(Swindelletal., 1997)，并参照测序峰图进行人工校正，采用BioEdit(V7.0.5.2)软件比对序列(Hall, 1999)。其中，ITS序列中存在简并碱基，因此首先利用DnaSP(V5.10)软件将比对校正后的ITS序列拆分，再进行后续计算分析(Libradoetal., 2009)。利用DnaSP(V5.10)统计叶绿体序列和ITS序列的单倍型，计算各居群的单倍型多样性指数Hd、核苷酸多态性指数π、基因流指数Nm、Tajima’sD以及Fu’s Fs，并进行失配分布分析。使用NETWORK(V5.0.0.3)软件中的Median-Joining算法构建单倍型间的网络结构图(Bandeltetal., 1999)，并分析各单倍型间的系统进化关系。运用Arlequin(V3.5)的分子方差分析模块计算甘蒙柽柳居群内、居群间的遗传变异(Excoffieretal., 2010)。使用PERMUT程序计算甘蒙柽柳的总遗传多样性HT、居群内平均遗传多样性HS、遗传分化系数GST和NST(Nei, 1973； 1978; Ponsetal., 1996)，检验是否具有明显的谱系地理结构。

2 结果与分析

2.1 引物筛选结果

本研究中trnH-trnK 引物扩增失败，未能获得目标条带；trnK引物扩增成功但无变异位点； 叶绿体基因trnQ-rps16序列以及核基因ITS序列扩增成功且发现变异位点(>3个)，因此采用trnQ-rps16和ITS序列对所有甘蒙柽柳样品进行扩增、测序。

2.2 甘蒙柽柳叶绿体和核基因单倍型变异与遗传多样性

叶绿体trnQ-rps16片段成功测序17个居群207个甘蒙柽柳个体，比对校正后的序列长度为590 bp，其中包含3个碱基突变位点、1个插入缺失位点，得到4个单倍型。单个种群的单倍型多样性(Hd)为0～0.6，核苷酸多样性(π)为0～1.69 ×10-3(表3)，总的单倍型多样性Hd为0.12，总的核苷酸多样性π为0.21 ×10-3。其中，永靖县太极镇(TJ)的Hd(0.509)、π(0.86×10-3)最高。

采用核基因ITS片段对232个个体测序分析，经比对后ITS序列长度为612 bp，包含29个碱基突变位点，无插入缺失，得到32个单倍型。单个种群的单倍型多样性(Hd)为0.108～0.867，核苷酸多样性(π)为0.18 ×10-3～4.45 ×10-3(表3)，总的单倍型多样性Hd为0.76，总的核苷酸多样性π为2.59×10-3。所有单倍型序列均提交至GenBank(https:∥www.ncbi.nlm.nih.gov/WebSub/?tool=genbank)，叶绿体基因trnQ-rps16序列的接收号为MK392589-MK392592，核基因ITS序列的接收号为MK392593-MK392624。

2.3 甘蒙柽柳叶绿体和核基因单倍型分布与种群遗传结构

叶绿体单倍型H1出现频率最高(194/207=93.7%)，分布于所有居群中，且位于网络结构图中心位置，推测为古老单倍型； H2(6/207)、H4(5/207)次之，H2分布于积石山县(JS)、永靖县兰永公路(LY)，H4分布于永靖县太极镇(TJ)、兰永公路(LY)； H3为昌邑市(CY)特有单倍型(图2A，表3)； H5-H8为外类群单倍型。

核基因单倍型R1频率最高(211/464=45.5%)，除永济市(YJ)、滨州市(BZ)外均有分布，R2(24/464)、R3(58/464)、R4(39/464)、R5(44/464)次之(表3)。由核基因ITS单倍型网络结构图可看出，甘蒙柽柳居群可分为R1、R5两个亚群，每个亚群又各自发生分化。单倍型R1分布范围广，频率最高，处于单倍型网络结构图的中心位置，推测R1为古老单倍型(图2B，表3)。单倍型R33-R39为外类群单倍型。

表3 甘蒙柽柳居群内单倍型分布及多态性①Tab. 3 Distribution and polymorphism of haplotypes in each population of T. austromongolica

①n: 序列数目；H: 所含单倍型数目；Hd: 单倍型多样性；π: 核苷酸多样性。n: Number of sequences sampled;H: Number of haplotypes;Hd: Diversity of haplotypes;π: Nucleotide diversity.

图2 甘蒙柽柳单倍型网络结构Fig.2 The haplotypes network of T. austromongolicaA. 叶绿体基因单倍型网络结构； B. 核基因单倍型网络结构。A. The network of cpDNA haplotypes; B. The network of nDNA haplotypes.

对叶绿体基因片段分析发现，甘蒙柽柳总的遗传多样性HT为0.13，居群内平均遗传多样性HS为0.109； 居群间遗传分化系数GST为0.19、NST为0.15(P>0.05)，基因流Nm=1.38。核基因ITS序列结果显示，甘蒙柽柳的整体遗传多样性HT为0.82，居群内平均遗传多样性HS为0.618； 居群间遗传分化系数GST为0.24、NST为0.22(P>0.05)，基因流Nm=0.89。分子方差分析(AMOVA)结果表明，叶绿体基因76.4%的遗传变异存在于居群内，核基因居群内的遗传变异为75.7%，显著高于居群间的遗传变异(cpDNA:16.53%； ITS:16.92%)(表4)。

表4 甘蒙柽柳种群的分子方差分析Tab. 4 Analysis of molecular variance(AMOVA)in T. austromongolica populations

2.4 种群历史变化动态分析

叶绿体基因的Tajima’sD值为-1.389 49，Fu’s Fs值为-2.813； 核基因的Tajima’sD值为-1.736 39，Fu’s Fs值为-23.317。叶绿体和核基因的Tajima’sD和Fu’s Fs值均为负值，此外，失配分析发现甘蒙柽柳期望扩张群体的分布曲线与实际观测所得分布曲线基本吻合，呈现单峰分布(图3)，表明甘蒙柽柳种群经历过快速扩张。

图3 甘蒙柽柳叶绿体基因(A)和核基因(B)的失配分析Fig.3 Mismatch-distribution analysis of cpDNA (A) and nDNA (B) in T. austromongolica

3 讨论

3.1 甘蒙柽柳的遗传多样性

甘蒙柽柳叶绿体基因trnQ-rps16片段的结果显示，居群总的单倍型多样性Hd为0.12，总的核苷酸多样性π为0.21 ×10-3，总的遗传多样性HT为0.13，低于Petit等(2005)所汇总的170个物种的叶绿体变异的平均值(0.67)，表明甘蒙柽柳叶绿体基因在物种水平上的遗传多样性较低。Liang等(2018)通过研究甘蒙柽柳和柽柳2个物种的谱系地理结构，发现叶绿体基因的遗传多样性较低。张如华(2011)基于3对cpSSR引物对柽柳进行遗传多样性分析，结果显示柽柳的遗传多样性为0.28，同样处于较低水平。这可能是由于叶绿体DNA为母系遗传，靠种子进行传播，而甘蒙柽柳的种子小且被毛，无休眠期，遇水即可萌发，可随风力或水力远距离传播，也可附着于动物身上进行传播(Hamilton, 1962； 王仲礼， 2005)，居群间基因交流较为容易，加之叶绿体基因相对保守，导致甘蒙柽柳叶绿体基因的遗传多样性水平较低。

基于核基因ITS序列的研究发现甘蒙柽柳物种水平上的遗传多样性较高，居群总的单倍型多样性Hd为0.76，总的核苷酸多样性π为2.59 ×10-3，总的遗传多样性HT为0.82，高于Petit等(2005)基于77个物种的核基因标记得出的被子植物核基因变异的平均值(0.137)。张如华(2011)采用SSR分子标记研究了10个柽柳居群的遗传多样性，结果显示柽柳物种水平的基因多样度为0.590 6； 李锐(2007)对6个柽柳居群研究发现其Nei多样性指数(h)为0.456； 赵景奎等(2008)通过RAPD技术，研究的黄河三角洲3个自然柽柳居群的种水平基因多样度为0.406 1。本研究结果与上述研究结果一致，均表明柽柳属在物种水平上具较高的遗传多样性。物种的繁育系统是影响其遗传多样性的重要因素，广泛分布的、异交长寿木本植物具有较高的遗传多样性(Hamricketal., 1989; Nybometal., 2000)。甘蒙柽柳为多年生木本植物，广泛分布于黄河流域，适应性强，性喜水也耐盐碱，一般为优势种和建群种，群体数量较大，保留了较丰富的遗传多样性。另外，核基因为双亲遗传，通过连锁互换、自由组合等形式发生基因重组，增加了居群的遗传变异。甘蒙柽柳的繁育系统为兼性异交，需要靠蝇类和蜜蜂等多种昆虫为其传粉，大大增加了异交的机会； 此外，甘蒙柽柳一年三季开花，花期5—9月，花朵和花序数量众多，且每朵花中含10～20个胚珠，为产生大量种子提供了基础，果实成熟期和种子萌发期较短(王仲礼， 2005)，庞大的种子库为其遗传多样性提供了物质基础。

3.2 甘蒙柽柳的遗传结构

居群遗传分化系数GST是基于单倍型频率的遗传指数(Nei, 1973； 1978)，而NST是基于单倍型的频率和遗传距离的遗传指数，因此，当NST显著大于GST时，表明居群间存在明显的谱系地理结构，当NST等于或小于GST时，表明居群间遗传分化不明显(Ponsetal., 1996)。本研究中，通过对叶绿体基因片段研究发现，居群遗传分化系数GST为0.19，NST为0.15(P>0.05)，小于Petit等(2005)总结的124种被子植物母系遗传的叶绿体遗传分化系数GST(0.637)，表明甘蒙柽柳居群间遗传分化较小，且NST、GST均不显著(P>0.05)，NST小于GST，表明居群间无明显的谱系地理结构。居群总的遗传多样性HT(0.13)与居群内平均遗传多样性HS(0.109)数值相近，说明居群间基因交流较为频繁，本研究中甘蒙柽柳居群间具有较高叶绿体基因流(Nm=1.38)也证实了这一点。与其他使用叶绿体基因序列对广布种天然居群的研究相比，如朱仁斌(2016)利用5条叶绿体DNA对自然分布区的文冠果(Xanthocerassorbifolium)研究发现，其遗传分化系数GST为0.783，基因流为0.03； 孙荣喜(2017)对枫香树(Liquidambarformosana)的研究结果显示GST为0.645，Nm为0.18； 本研究中的遗传分化系数GST较低，基因交流相对频繁，这可能与甘蒙柽柳特殊的生境有关。甘蒙柽柳为黄河流域特有种，大多生长在黄河周边，杨维康等(2002)调查研究了新疆境内柽柳属植物的分布与群落类型，表明水在柽柳属植物现代地理分布形成中起决定性作用，甘蒙柽柳的种子可随风力或水力远距离传播，也可附着于动物身上进行传播(Hamilton, 1962； 王仲礼， 2005)，居群间较为频繁的基因交流，可能是导致遗传分化较小的主要原因。

核基因ITS序列的结果显示，居群间遗传分化系数GST为0.24，NST为0.22(P>0.05)，虽高于Petit等(2005)总结的基于双亲遗传标记的22种植物的平均遗传分化程度(GST=0.15)，但GST和NST均不显著，且NST小于GST，说明居群间无明显的谱系地理结构。AMOVA结果亦显示甘蒙柽柳的遗传变异主要发生在居群内(75.68%)，明显高于居群间(16.92%)的遗传变异。张如华(2011)通过SSR分子标记对10个柽柳居群的遗传变异进行研究，SSR结果显示柽柳居群的遗传分化系数GST为0.050 3，基因流为4.736 5； 赵景奎等(2008)基于RAPD标记对黄河三角洲的3个柽柳居群进行研究，发现GST为0.050 7，基因流为9.356 4。本研究中的GST(0.24)值大于上述研究结果，基因流相对较小，一方面可能是由于分子标记的不同所导致的差异，另一方面因为本研究中所采集的样本范围广泛，西起青海省，东至山东省，地理距离在一定程度上限制了居群间的基因交流。但甘蒙柽柳的结果亦显示居群间遗传分化较小，遗传变异主要存在于居群内，这可能与甘蒙柽柳的生境以及种子易于传播有关，河流促进了各居群间的基因交流，致使其遗传分化系数降低。

3.3 甘蒙柽柳种群的历史动态分析

河流可作为非生物物理屏障，也可为有些物种分布和基因流提供地理便利或通道(刘慎谔，1985； Gengetal., 2015)。刘慎谔(1985)在河流流向对植物分布的影响中提到，缫丝花(Rosaroxburghii)多分布于四川、贵州，而云南仅大理有分布，推测该物种是沿着金沙江古道扩散至大理的。Zhang等(2011)探讨了河流被袭夺或改道对滇榄仁分布格局和遗传多样性变化的影响，发现滇榄仁现存的间断分布模式与红河改道造成河流重建的地质事件有关。本研究中发现甘蒙柽柳居群经历过快速扩张(图3)，这可能与黄河的贯通有关。

单倍型网络结构图可以反映单倍型的进化关系，研究表明，一般位于网络结构图的中间位置、频率较高、地理分布广泛的单倍型更为古老(Freelandetal., 2012)，本研究中叶绿体单倍型H1、核基因单倍型R1的频率最高，位于网络结构图的中心位置，且地理分布最为广泛(图2)，可能为古老单倍型。基于奠基者效应，原始居群的遗传多样性要比迁移扩散形成的居群高，具有更多的特有单倍型(Hewitt, 2000)。本研究中，叶绿体基因结果显示，永靖县太极镇(TJ)的Hd、π最高，积石山县(JS)次之，具特有单倍型H2、H4； 核基因ITS的结果显示，永靖县兰永公路(LY)的π最高，Hd较高，具特有单倍型R31、R32(表3)，而永靖县与积石山县以黄河为界，隔河相望，因此，认为该地区附近可能为甘蒙柽柳居群的起源中心。种群遗传结构与物种地理分布区的变化息息相关(Ricklefs, 1987; Abbottetal., 2000; Hewitt, 2000)。青藏运动C幕时期，湟水谷底和黄河兰州段首次形成了连续分布的河流冲积层，在兰州发现了黄河最高级阶地，此处的黄河最为古老(李吉均等，1996； 李长安等，1998)。永靖县临近于湟水与黄河的交汇处，东北与兰州市接壤，这也支持了甘蒙柽柳起源中心在永靖县与积石山县附近的推测。

叶绿体基因为母性遗传，无重组，具有独立的进化历史，相对于核基因更能反映居群的历史动态(Kingetal., 1998; Schaaletal., 1998; Petitetal., 2005)，本研究中的叶绿体结果显示，TJ、LY、JS的单倍型种类丰富，多态性高，而其他居群大多只有1种单倍型(H1)，推测甘蒙柽柳居群是从甘肃省分别向上游的青海省和中下游省份迁移的，奠基者效应造成新建成的居群遗传多样性较低。黄河经昆黄运动、共和运动，古黄河一支流——刘家峡水系向高原溯源侵蚀，先后进入青海省的循化盆地、共和盆地、兴海盆地等，并贯通了中下游段的一系列古湖泊、峡谷，最终形成现代统一的黄河(李长安等，1998； 李吉均，1999； 李吉均等， 2001)。因此，推测黄河的形成贯通是导致甘蒙柽柳居群迁移扩散的重要因素。

4 结论

本研究探讨了我国特有种甘蒙柽柳的谱系地理结构及居群迁移演化进程。甘蒙柽柳叶绿体基因的遗传多样性水平较低，核基因的遗传多样性较高。居群间的遗传分化小，遗传变异主要发生在居群内，无明显的谱系地理结构。该居群经历过快速扩张，这可能与黄河的贯通有关，单倍型H1、R1为古老单倍型，甘蒙柽柳的起源中心可能为甘肃省积石山县、永靖县附近，随着黄河的形成从甘肃省分别向黄河上游(青海省)和中下游迁移。