基于叶绿体基因的文冠果遗传多样性及其遗传结构

刘 莉，赵永军，王 磊，吴 丹，鲁仪增，陆 璐，王 震，庄振杰，解孝满,*

(1.山东师范大学 生命科学学院，山东 济南 250014；2.山东省林草种质资源中心，山东 济南 250102)

文冠果(Xanthocerassorbifolium)又名文官树、崖木瓜、文光果等，为无患子科(Sapindaceae)文冠果属(Xanthoceras)落叶灌木或小乔木，单属单种[1-2]。其不仅适应性强，观赏性好，而且叶、果、树皮具有多种用途；种仁营养价值高，籽油主要脂肪酸为亚油酸和油酸，且富含植物源神经酸，食用历史悠久，是我国北方具有发展潜力的特有木本油料树种[3-4]。随着分子生物学的发展，分子标记早已运用到林木遗传育种研究中[5]。文冠果核基因组于2019年发布共享，2021年初步实现了基于文冠果基因组的SSR标记开发及验证，有效促进了文冠果的遗传多样性与遗传结构研究[6-8]。但是，目前文冠果基于核基因组的研究总体起步较晚，技术尚不十分成熟，参试的样本较少，不足以全面揭示文冠果种质资源遗传多样性特征。另外，由于文冠果交叉引种历史较长，基于核基因组解析当前栽培群体和野生群体共存格局下的文冠果遗传多样性和遗传结构及其迁移路线存在较大的技术瓶颈。

在大多数被子植物中，母系遗传的叶绿体基因组具有高度保守、非编码区进化速度快等特点[9]，而且叶绿体基因组非编码区中的基因间区具有丰富的遗传变异[10]，叶绿体基因单倍型可直接用于植物遗传多样性检测[11]，在系统发育研究方面具有独特优势[12-13]。朱仁斌等[14]利用5条叶绿体基因片段对7个省份的38个文冠果野生群体(399个个体)进行了地理谱系调查，研究表明文冠果具有较高的遗传差异(Fst=0.765)和相对较低的基因流(Nm=0.03)。但该研究的参试样本来自“东南—西北”方向的采样带，缺失内蒙古、黑龙江、北京和山东等区域样本，而且采用的引物数量较少，挖掘的单倍型较少，得出的进化趋势及相关结论仍需要进一步验证。

为进一步深入了解当前文冠果栽培群体和野生群体共存格局下的遗传多样性及其引种迁移格局变化，本研究基于叶绿体全基因组单核苷酸多态性(single nucleotide polymorphism，SNP)及其单倍型特征进行遗传多样性和遗传结构分析，以期为文冠果种质资源的有效保护与高效利用提供理论依据。

1 材料与方法

1.1 材料来源

以生长势强、结果量大为指标，于2020-2021年从文冠果主要集中分布的16个地区(35°43′48″-44°34′48″N，107°39′01″-129°34′48″E)的栽培群体和野生群体进行选优，共计筛选优树172株，其中野生优树81株，栽培优树91株。调查优树生境信息，然后采集叶片，并快速置于硅胶颗粒中干燥备用。为确保文冠果样品具有高度代表性，保持各采集单株距离至少200 m，以确保样本之间无亲缘关系；群体内的样本数量应15份以上，数量不够的群体以实际采集份数为准，各群体材料来源情况见图1。

图1 不同文冠果群体采样基本信息

1.2 样本DNA提取及叶绿体基因组测序组装

采用磁珠法提取文冠果叶片DNA，Nanodrop 2 000分光光度计(IMPLEN,CA,USA)、Qubit 3.0 Flurometer(life technologies,CA,USA)和1%浓度(质量浓度)的琼脂糖凝胶分别进行纯度、浓度和完整性检测。合格样品由安诺优达基因科技(北京)有限公司使用Illumina测序平台测序，读长为PE150，测序深度10×。Bcl-convert软件对原始数据(raw reads)进行质控过滤，使用BWA软件的mem算法从每个样品的全基因质控数据(clean data)中提取出比对文冠果叶绿体参考基因组NC_037448.1的数据[13-14]；使用Get Organelle程序包进行叶绿体基因组的组装[15]；使用PGA软件进行注释[16]。

变异分析软件GATK检测叶绿体基因组中所有的潜在的SNP位点[17]，使用Variscan软件计算文冠果各群体的Nei多样性指数(Nei)、香浓指数(H)、多态性信息含量(PIC)、核苷酸多样性(Pi)以及群体间的遗传分化指数(Fst)，Tajima's D值。其中：根据多态性信息含量的衡量标准：PIC≥0.5 时为高度多态性；0.25≤PIC<0.5 时为中度多态性；PIC<0.25 时为低度多态性[18]；种群分化程度用Fst来衡量，其取值从0到1，0表示2个种群间是随机交配的，基因型完全相似，1则表示2个种群完全隔离[19]。

表1 SNP突变模式分布统计

172个样本共有基因序列排序连接后，使用MAFFT(v7.205)比对后得到一致性比对文件[20]；使用DnaSP(v5.10)软件计算单倍型参数、PopART(v1.7)软件构建单倍型网络图(TCS network)、Admixture软件构建遗传结构(structure)图、GVTA软件构建亲缘关系(Kinship)图[21-22]。

2 结果与分析

2.1 文冠果叶绿体基因组特征及其差异分析

2.1.1 基因组特征分析 由图2可以看出，172个文冠果优树样本的叶绿体基因组均为共价闭环的双链环状分子，典型的4部分结构，包含有大单拷贝区(LSC，large single-copy)、小单拷贝区(SSC，small single-copy)和2个反向重复区(IRA和IRB，inverted repeats)。其序列全长位于159 048～160 069 bp，均小于参考叶绿体基因组(161 231 bp)；LSC 区段长度为86 190～86 946 bp，均大于参考叶绿体基因组的LSC长度(85 299 bp)；SSC区段长度为18 097～18 988 bp，参考叶绿体基因组的SSC长度位于样本SSC区段长度之间；IRs区段长度为27 243～27 257 bp，均小于叶绿体参考基因组的IRs区段长度(28 620 bp)；且其GC含量均为37.7%，与参考叶绿体基因组的GC含量一致[23]。

图2 以QS101为例文冠果叶绿体基因组环状图

注释结果中发现172个文冠果优树资源的叶绿体基因组中的基因种类和数量没有明显差异，其中1个蛋白编码基因(rps12)含有3个拷贝；21个基因含有2 个拷贝，包括 10个蛋白编码基因(ndhB、rpl2、rp123、rps7、trnA、trnA-TGC、trnI-GAT、trnI-CAU、ycf15、ycf2)、6个tRNA基因(trnI-GAU、trnL-CAA、trnV-GAC、trnA-UGC、trnR-ACG和trnN-GUU)和4个rRNA基因(16SrRNA、23SrRNA、4.5SrRNA和5SrRNA)。

2.1.2 叶绿体基因组SNP差异分析 与参考叶绿体基因组相比，共计从172个参试叶绿体基因组中检测出6类突变，共有52 460个SNP位点。其中，C-T(G-A)突变类型数量最多,其次为T-C(A-G)变异数量为12 212个，C-G(G-C)突变类型最少为170个。由图3可以看出，非编码区的突变位点较多，其中atpI、petL、petN、rpoC2、petA、psbE、psbM、psbJ、trnE-UUC和trnT-GGU突变位点相较于其他基因多(Pi>2%)。编码区ycf1基因的突变相较于其他基因也较高。总体来看，非编码区具有高突变性，文冠果高突变位点主要在控制花器官发育和光合作用等有关基因上，为探究文冠果多样性提供具体研究价值。

图3 SNP位点的(Pi%)的百分比

2.2 不同群体文冠果遗传多样性和遗传结构分析

2.2.1 遗传多样性 遗传多样性分析显示(表2)，11个文冠果群体的叶绿体基因组SNP位点多态性信息含量(PIC)在0.12(RG8)～0.35(RG4)，平均为0.27，看出多态性处于中低度水平。香浓指数(H)为0.25(RG8)～0.69(QS1)，平均为0.50。Nei多样性指数(Nei)在0.14(RG8)～0.48(RG4)，平均为0.3，均<0.5，说明文冠果基因多样性水平偏低。核苷酸多样性指标(Pi)平均值为4.00×10-3，其中最高的是QS7群体，为6.12×10-3，核苷酸多样性指标(Pi)总体呈现比较低的状态，其中QS6群体和RG8群体的Pi值最低。给每项指标赋值，每项排名相加总分按照多样性综合指标从大到小排序，群体遗传多样性水平为QS9>QS5>QS10>RG4>QS7>QS1>RG3>QS2>RG1>RG6>RG2>QS4>RG5>QS8>QS6>RG8。大部分野生群体的遗传多样性水平较高，而人工引种种植群体的多样性较低，存在引种单一问题。

表2 供试119个文冠果优树多态性信息

对16个群体进行Tajima's D计算，结果在-0.115～0.258，Tajima'sD值为负时，表明种群经历过快速扩张，其中QS6、QS8、RG1和RG8群体的Tajima's D值为负值，说明这3个群体经历了快速扩张，QS9和QS10群体的Tajima's D值为0，符合中性进化假说。

2.2.2 遗传分化分析 遗传分化系数(Fst)分析显示(表3)，文冠果不同野生群体间的遗传分化系0.012～0.933，其中最小的是QS9群体和QS6群体，说明这2个群体存在较大的基因交流，遗传分化不明显，QS1、QS2、QS4和QS5与QS7、QS8、QS9和QS10群体两两之间Fst值较大，说明黄河为界东西部文冠果资源交流较少。在人工种植群体中，几个种植群体与QS10和QS8群体的Fst值都>0.5，几个种植群体与QS7和QS9群体的Fst值也在0.5左右，表示群体间存在着很大的遗传分化。

2.2.3 群体结构分析 设置假定祖先数K为1～8，并分别计算其CV error值，并制作百分比堆积图，当K等于n时CV error值最小，则认为群体中包含n个亚群较为合理。由图4A可知，当K=4时为最优结果，因此172份文冠果资源分为4个祖先成分最为接近于群体的真实情况；如图4B所示，本次分析的4个祖先成分主色段分别为蓝色、黄色、橘色以及紫色。第Ⅰ群组主色段为黄色，包括7个野生群体(QS2、QS4、QS6、QS7、QS8、QS9、QS10)，共50份资源；第Ⅱ群组主色段为蓝色，包括2个野生群体(QS1和QS2)，7个种植群体(RG1、RG2、RG3、RG4、RG5、RG6和RG8)其叶绿体遗传背景较纯；第Ⅲ群组主色段为橘色，4个野生群体(QS1、QS2、QS4和QS5)，7个种植群体(RG1、RG2、RG3、RG4、RG5、RG6和RG8)；第Ⅳ组为紫色，为QS2群体的4份文冠果资源。图4C、图4D与Structure图结果吻合。

注：文冠果不同群体数量CV error曲线(A)；K=4时，文冠果Structure图(B)；文冠果KinShip图(C)；文冠果主成分(PCA)分析图(D)。

2.3 文冠果叶绿体DNA单倍型多态性及其谱系地理结构分析

2.3.1 单倍型变异、分布及其多态性特征 由表4可以看出172份文冠果资源共有65种单倍型。QS6群体有14种单倍型，单倍型数量最多，QS2和QS4有8种单倍型，QS8有7种单倍型，RG4、RG5、RG6群体都有5种单倍型。QS9群体单倍型数量只有1个，是含有单倍型最少的群体。16个群体中有8个群体含有Hap_18，其中群体内有37份资源含有该单倍型；7个群体含有Hap_11单倍型，其中有24份资源含有该单倍型。RG8群体(鄂尔多斯)中含有6种单倍型，在种植群体中含有单倍型种类数量最多；RG3群体(承德市)含有3种单倍型，在种植群体中含有单倍型种类数量最少。由表4可以看出，Hap_11、Hap_18单倍型在各个群体中分布广泛，这2个单倍型是文冠果亲缘关系远近的重要的单倍型。

表4 各群体所含单倍型类型分布

由表5可以看出，QS1、QS2、QS5、QS7、QS8和QS9群体单倍型多样性(Hd)达到1，其次是QS6和QS4群体的单倍型多样性(Hd)达到0.98；最低是QS10群体的0.6。可见，人工种植群体的多样性都普遍低，人为影响较大。QS5群体的单倍型多样性方差(Vh)与单倍型多样性标准差(Sh)均为最高，分别为0.009 26和0.096，QS6群体的单倍型多样性方差(Vh)与单倍型多样性标准差(Sh)最低，分别为0.000 44和0.021。

表5 供试172个文冠果优树单倍型多样性信息

2.3.2 单倍型系统发育及其遗传结构 根据组装好的国内172个文冠果叶绿体基因组进行单倍型分型并绘制单倍型网络图。图5为野生群体文冠果资源单倍型网络图，每个彩色圆代表着一种单倍型，圆的大小表示该单倍型所含有的样本总量，2个圆之间的连接代表2个单倍型彼此相关，而连线上的短线代表由一种单倍型变成另一种单倍型所需要经过的碱基替换个数。另外，每个黒色的圆点表示推断缺失样品可能存在的单倍型。该结果与叶绿体基因组遗传结构分组相对应。第1支是由Hap_55变异形成的分支，主要包括28个单倍型，Hap_55是野生群体QS9的单倍型，属于宁夏固原群体，由此单倍型变异而形成的野生资源单倍型主要包括QS4、QS6、QS7和QS8群体。除Hap_55到Hap_30和Hap_39的变异较大些外，其他单倍型之间的突变位点小，差异较小，这部分单倍型地理位置主要是在黄河几字形西南侧，子午岭、黄龙山南部。这部分区域地形有河流与山脉，例如，QS9和QS8群体之间的六盘山、泾河、渭河，可能山的运动促使QS9群体的Hap_55的资源向四周扩散。对应的第2部分是缺失型单倍型形成的分支，包括QS1群体的资源和QS2群体部分资源，这部分资源位于陕晋段黄河东侧、太行山南侧区域；第3支是由Hap_11形成的分支，是由太行山西北侧的QS1群体中的单倍型变异而来，形成QS4群体和QS5群体的单倍型，这3个群体都分布在太行山、吕梁山、黄龙山和子午岭形成的山脉的北侧，因此推测该单倍型的分布与山脉有关；第4支是由缺失型单倍型变异而形成的4个单倍型，分别为Hap_14、Hap_15、Hap_16和Hap_17，这几个单倍型对应的是QS2群体中的QS201、QS202、QS203和QS204，为山西临汾大宁县群体资源，这几部分为该区域特有单倍型，推测该区域有文冠果避难所。

图5 野生文冠果资源单倍型网络

文冠果引种时间较短，因此人为种植的文冠果叶绿体基因组变化不大。将野生群体与人为引种群体，共同建立单倍型网络图(图6)，通过单倍型变化探究人为引种的资源从哪些地区引种。在Hap_11形成的第3分支上，RG1、RG3、RG4和RG5群体可以确定从QS1群体引种，有共同的单倍型Hap_11。RG8资源与QS5群体共同单倍型为Hap_27。Hap_29的资源有QS5、RG4、RG5和RG8群体，且RG8群体资源的单倍型由Hap_29变异而来，基本可以确定RG8群体的文冠果资源是由QS5引种而来。Hap_29的资源有RG4、RG5和RG6，而这几份资源是由QS10群体的资源引种而来。

图6 供试172个文冠果优树单倍型网络

3 结论与讨论

遗传多样性作为生物多样性的重要基础，是物种长期生存、进化和适应的结果[24]。文冠果长期处于野生状态，交叉引种繁杂，采用母系遗传的叶绿体基因组SNP变异位点及单倍型可以较好地解决文冠果种质资源的遗传多样性及其遗传结构。保护物种，换个层面说就是保护物种的进化潜力或遗传多样性[25]。若物种内的遗传多样性水平越大，则说明对环境变化的适应性就越大，其进化的潜力也就越大[26]。但是，保护物种应该有对物种多样性水平和种群遗传结构全面的了解的基础，因此了解物种的遗传变异的空间分布模式对于科学有效地保护物种起着非常重要的作用[27]。目前文冠果野生资源，由于人为砍伐，已经逐渐减少，因此需要提高文冠果资源保护意识，以便充分发挥文冠果的经济、生态、美学价值。总体来说，文冠果分布范围窄，适生区内地形变化小，致使文冠果遗传变化较小。本研究发现，野生区域的文冠果优树资源遗传多样性较高，而栽培区域文冠果遗传多样性普遍较低，因此在当前野生群体与栽培共同存在的格局下，优先保护野生群体，重点保护部分优株群体，可促进文冠果种质资源的有效保护与高效利用。

3.1 基于SNP的多态性分析与群体结构

SNP标记往往是二等位多态性的，从理论上看，每一个SNP位点都可以有4 种不同的变异形式，但实际上发生的只有2种，即转换和颠换。由于这种形式只有2种等位碱基，可以直接估计等位基因频率，因此 SNP具有二态型[28]。叶绿体基因组SNP多态性分析，都表现出了野生群体的多态性大于人为栽培群体。研究表明，宁夏固原、陕西咸阳、宝鸡、山西大同、忻州、朔州群体的文冠果遗传多样性处于较高水平，辽宁省朝阳市和河北省承德市的种植群体文冠果遗传多样性水平较高，指导我们在保护搜集文冠果野生材料时，尽量加大这几个地区的取样量和优先对这些地区进行保护。文冠果虽然在内蒙古地区栽培面积大，但发现该区域遗传多样性水平并没有很高，因此，内蒙古引种文冠果问题上应注意多区域大范围引种，可以引种亲缘关系较远的文冠果资源作亲本，其后代就能获得更大的遗传变异，为文冠果育种选择提供更多的材料。文冠果野生资源人为破坏严重，应当对野外优树进行保护，以此保护文冠果的多样性。多态性信息含量(PIC)、Nei多样性指数(Nei)和核苷酸多样性指标(Pi)显示文冠果群体内为低度多态性水平，这比基于SSR标记的结果偏低，主要原因为叶绿体基因相较于核基因变异较小[29-31]。

根据叶绿体基因组结果将文冠果分为4个组，可以划为西部、中部、东北部和东南部。根据地形位置进行分析，影响分组的因素为太行山、吕梁山、黄龙山和子午岭形成的山系，陕晋段黄河和吕梁山围成的区域，及陕晋段黄河，各个区域的亲缘关系与Fst结果相同。

3.2 单倍型多态性及其分布特征

文冠果是古老的树种，在漫长进化过程中逐步积累了丰富的遗传变异，因此形成了多样的单倍型。本研究中，对119个文冠果的全叶绿体基因组进行测序，从整个基因组层面充分挖掘多态性位点，结果显示其单倍型变异丰富，多样性高达0.9，比文献报道的叶绿体单倍型多样性高(0.76)[32]。进而说明，采用叶绿体基因片段，虽然基于检测到的单倍型能够适度解析文冠果的谱系地理信息，但是标记数量有限，获得标记位点较少，不能准确地了解文冠果单倍型情况。

3.3 系统发育及引种迁移趋势

遗传变异在空间分布上产生遗传分化，造就了一个物种种群的遗传结构[33]。群体的遗传结构改变是由于随机的遗传漂变及基因流的隔离[34]。叶绿体基因比较稳定，不易变化，因此可以看出各群体间亲缘关系的远近。甘肃合水和鄂尔多斯群体间的分化系数大，说明这2个群体间的遗传分化程度最高，即其亲缘关系最远。承德和临汾群体遗传分化系数低，说明这2个群体亲缘关系比较近。大多数叶绿体基因组是母系遗传，文冠果通过种子传播来促进群体间的基因流动。从单倍型网络图可以看出，合水、榆林和临汾是野生资源主要分布区，这3个群体单倍型差异较大，结合野外观察，文冠果种子只能散落在母树四周，果实靠啮齿动物搬运进行传播，且3个群体被子午岭、吕梁山及黄河阻碍了种子传播，对群体间的基因交流产生了重要影响。但是由于野生资源样本量少，很多未知单倍型尚未挖掘，因此今后将进一步增加野生型资源进系统发育研究，进而探明文冠果引种驯化及其迁移历史，为文冠果种质资源的保护利用提供支持。

不论是从单倍型进行聚类还是SNP聚类，都表明文冠果分为4个亚群。从文冠果65种单倍型Network网络图中可以得出，文冠果有2次大规模扩散，分别为Hap_27、Hap_11和Hap_18。其中Hap_27是重要的单倍型，分布在陕西榆林和延安一带。Hap_11和Hap_18的扩张与人为因素有关。Hap_11和Hap_55位于网状图的中心位置为古老型单倍型，位于固原市的 Hap_55单倍型扩散的原因可能是由于六盘山、泾河、渭河的影响；而位于山西大同的Hap_11单倍型扩散原因不清楚，需要增加样本量进一步探讨。

在对人工栽培的文冠果与野生型文冠果的单倍型分析时，明确了现在人工栽培的文冠果来源地区。根据各个群体Fst值、共有单倍型及单倍型变异趋势，看出本研究中赤峰市翁牛特旗和巴林左旗、北京东城区、淄博临淄区、牡丹江、辽宁省朝阳市和河北承德的资源是从山西大同引种而来；内蒙古鄂尔多斯的文冠果是由吕梁市中阳县和榆林靖边县引种而来。