昭通牛mtDNA D-loop 区遗传多样性分析*

朱群媱，丁利民，容 斌，邓琦炜，胡 瑀，范新阳，苗永旺**

(1.云南农业大学 动物科学技术学院，云南 昆明 650201；2.昭通市畜牧兽医技术推广站，云南 昭通 657000；3.昭通市农业科学院 畜禽研究所，云南 昭通 657000)

牛属于脊椎动物亚门(Vertebrata)哺乳纲(Mammalia)偶蹄目(Artiodactyla)反刍亚目(Ruminar)洞角科(Covicornia)牛亚科(Bovedae)牛属(Bos)，可分为普通牛(B.taurus)和瘤牛(B.indicus) 2个类型。近年对世界范围内牛种线粒体DNA 控制区(mitochondrial DNA control region，mtDNA D-loop)序列的研究表明：普通牛和瘤牛各自起源于独立的驯化事件[1-2]。普通牛的驯化中心位于近东地区，欧洲、非洲和东北亚地区也可能存在普通牛的驯化事件；瘤牛则可能起源于印度河流域[3-4]。在普通牛中主要存在5 个mtDNA 世系(T和T1～T4)，且这些世系存在明显的地理区域分布，近东和欧洲地区以T、T2 和 T3 世系为主，非洲地区以T1 世系为主，而T4 世系为东北亚所特有；在瘤牛中存在2 个mtDNA 世系(I1 和I2)，常见于印度次大陆，且以I1 世系为主[1-3]。对中国家牛的mtDNA 世系研究表明：中国家牛普遍存在普通牛和瘤牛2 个血统[2,5]，且普通牛血统构成主要包含T1～T4 mtDNA 世系，瘤牛血统则包含I1 和I2 mtDNA 世系。

昭通牛因产于云南省昭通市而得名，在全市各县(区)均有分布。昭通市位于云南省东北部，地处云、贵、川三省结合部的乌蒙山腹地，群山林立，水资源和饲料资源丰富，为昭通牛遗传资源的形成和养殖提供了良好条件[6-7]。昭通牛属役肉兼用型牛种，体型中等，体质结实，结构匀称，短宽平直，被毛以黄色为主；对当地自然生态条件有良好的适应性，具有耐粗饲、抗逆性和抗病力强等优良种质特性[6]。昭通牛于1980 年全国畜禽品种资源调查时定名，属云南省地方优良品种，于1987 年列入《云南省家畜家禽品种志》，2011 年列入《中国畜禽遗传资源志·牛志》。

畜禽的遗传多样性不仅可以描述群体的遗传背景信息，为其保种和选育提供参考，还对群体的亲缘关系、起源分化和基因渗入等研究具有重要意义[1]。近年来，多个遗传标记被用于遗传多样性的研究中，其中mtDNA D-loop 序列主要用于畜禽母系遗传研究。mtDNA 具有分子结构简单、不重组和高突变等特征，在哺乳动物中遵循严格的母系遗传，是物种鉴定和祖先母系起源研究的最佳遗传标记。mtDNA D-loop 区是线粒体基因组中突变率最高的区域，因此被广泛用于畜禽遗传多样性研究[8]。2000 年前，引入外来商业化牛品种进行的杂交改良对昭通牛的影响较大，导致其血统混杂，纯种牛数量急剧减少[6-7]。虽然已有研究对昭通牛的遗传多样性进行了报道[9-11]，但都是几年甚至10 年前的研究，且当时研究的样品数量有限，不能准确反映当前该牛的群体遗传结构状况。考虑到当前昭通牛保种与选育利用的迫切需要，本研究利用mtDNA D-loop 区序列为遗传标记，采用PCR 产物直接测序技术对当前该牛的群体变异进行检测，进一步结合已发表的云南本地其他牛的mtDNA 数据进行深入的群体遗传学分析，试图从母系遗传的角度阐明昭通牛的群体遗传结构、血统构成和起源分化，为该牛的保种和选育提供依据。

1 材料与方法

1.1 研究材料

在云南省昭通市昭阳区旧圃镇和彝良县洛泽河镇采集了98 份昭通牛耳组织样品(45 头母牛和53 头公牛)，置于装有500 µL 95%酒精的1.5 mL 离心管中带回实验室，—20 ℃保存，用于基因组DNA 的提取。样品的采集由当地畜牧部门的专家负责纯种牛的鉴定，要求具备典型的品种特征、无直接血缘关系且在分布上具有代表性。为了更好地评价昭通牛的遗传多样性特征并弄清其与云南其他地方牛种的遗传差异，从NCBI数据库(https://www.ncbi.nlm.nih.gov/)下载了189条云南地方牛的mtDNA D-loop 区序列(GenBank登录号：MF410466～MF410654)，其中包括53 头文山牛、34 头德宏高峰牛、49 头迪庆牛和53 头滇中牛，用于与本研究数据进行比较分析。

1.2 基因组DNA 的提取

采用酚—氯仿法提取基因组DNA，并用TE 液溶解，再用1.5%琼脂糖凝胶电泳检测DNA原液的完整性，用Nano Drop 2000 紫外可见分光光度计(Thermo Fisher Scientific，Waltham，MA，USA)检测其质量浓度和OD 值，并稀释成质量浓度为100 ng/µL 的DNA 工作液备用。

1.3 目的片段扩增及序列测定

所用引物为LOFTUS 等[12]公布的特异性引物扩增D-loop 全序列片段(Cattle-F_L15757：5 ′ -CTGCAGTCTCACCATCAACC-3 ′，Cattle-R_H496：5 ′-GGGGTGTAGATGCTTGC-3′)，目的片段扩增及测序所需引物相同。PCR 反应体系为20.0 μL，包括10×buffer (Mg2+15 mmol/mL) 2.0 μL，dNTP 1.6 μL，10 μmol/μL 上、下游引物各0.4 μL，5 U/μL rTaqDNA 聚合酶0.1 μL 和100 ng/µL 的 DNA 工作液模板2.0 μL，加ddH2O 补足体系。反应程序为：95 ℃预变性3 min→34 个循环(94 ℃变性30 s→60 ℃退火40 s→72 ℃延伸40 s)→72 ℃后延伸5 min。PCR 产物由生工生物工程(上海)股份有限公司进行双向直接测序。

1.4 数据统计与分析

采用Lasergene 软件包(DNA Star Inc.，USA)[13]中的SeqMan 程序对昭通牛mtDNA D-loop 序列进行人工校正和输出；采用ClustalX 软件[14]进行序列比对，采用BioEdit 软件[15]调整和编辑比对好的序列；采用MEGA6 软件[16]对序列进行单倍型输出和碱基组成分析，计算群体内遗传距离和群体内遗传多样性指数等；群体间遗传距离采用最大似然法基于T92+G 模型进行估算，采用Bootstrap 法(5 000 次)估算稳健的标准误差，然后采用邻接法根据群体间遗传距离构建群体间遗传关系树，并利用在线软件iTOL (https://itol.embl.de/#)对其进行美化；利用DNASP 6 软件[17]计算群体遗传多样性参数(单倍型数、单倍型多样度和核苷酸多样度)，估算群体间的遗传分化指数(genetic differentiation index，Fst值)，根据公式Nm=(1/Fst—1)/4[18]计算基因流(gene flow，Nm值)；采用Network 10 软件(https://www.fluxus～engineering.com/)的中接法构建昭通牛群体各单倍型遗传关系的中介网络图，其中赋予转换和插入/缺失的权重分别为10 和15；采用MEGA6 软件[16]构建昭通牛群体单倍型的系统发育树。

2 结果与分析

2.1 序列变异

昭通牛mtDNA D-loop 区的突变位点及定义的单倍型如图1 所示。98 头昭通牛的mtDNA Dloop 区序列扩增片段长915 bp，去掉两端测序质量较差的序列后，用于分析的序列长度为898 bp。序列中A、T、C 和G 碱基的平均含量分别为33.30%、29.00%、23.90%和13.80%，序列A+T 平均含量为62.30%，共检测到73 个核苷酸替换位点和3 个插入/缺失，核苷酸替换位点约占检测核苷酸位点总数的8.13%，核苷酸替换主要为转换。在发现的73 个突变位点中，包含18 个单一信息位点和55 个简约信息位点。基于发现的核苷酸替换位点及插入/缺失信息，在昭通牛群体中共确定36 种单倍型，单倍型H02 和H22 在群体中所占比例较高，其频率分别为29.59%和11.22%。昭通牛群体的单倍型多样度为0.880±0.027，核苷酸多样度为0.024 43±0.000 94，群体内遗传距离为0.026±0.004，群体内遗传多样性指数为0.027±0.004。

图1 昭通牛mtDNA D-loop 区的突变位点及定义的单倍型Fig.1 Mutation sites and defined haplotypes based in mtDNA D-loop sequences of Zhaotong cattle

2.2 群体遗传多样性

云南5 个本地牛群体共287 条mtDNA D-loop区序列中A、T、C 和G 碱基的平均含量分别为33.40%、28.80%、24.00%和13.80%，A+T 和C+G的平均含量为62.20%和37.80%，共发现93 个变异位点，其中单一信息位点24 个，简约信息位点69 个，共确定86 个单倍型；单倍型多样度为0.882±0.017，核苷酸多样度为0.024 82±0.000 30，群体内遗传距离为0.027±0.004，群体内遗传多样性指数为0.027±0.005 (表1)。

2.3 群体间遗传距离和遗传分化指数

由表2 可知：昭通牛与迪庆牛间的遗传距离最小(0.024)，与德宏高峰牛间的遗传距离最大(0.031)。根据群体间遗传距离数据构建的遗传关系树(图2)可知：昭通牛与迪庆牛聚为一大支，德宏高峰牛与文山牛、滇中牛聚为另一支。

表2 云南地方牛群体间遗传距离Tab.2 Genetic distances between local cattle in Yunnan

图2 基于群体间遗传距离构建的聚类树Fig.2 Cluster tree based on the genetic distances between populations

由表3 可知：各群体间的遗传分化指数(Fst值)在—0.015 34～0.671 32 之间，昭通牛与其他地方牛群体间的Fst值在0.023 32～0.552 10 之间，其中，与迪庆牛的Fst值最小，与德宏高峰牛的Fst值最大；昭通牛与其他地方牛群体间的基因流(Nm值)在0.20～10.47 之间。

表3 云南地方牛群体间遗传分化指数和基因流Tab.3 Genetic differentiation index and gene flow among of local cattle in Yunnan

2.4 世系组成

按照近年对普通牛和瘤牛mtDNA 世系划分的标准[2-3]，云南5 个本地牛的mtDNA 世系划分结果(表4)显示：昭通牛存在普通牛和瘤牛2 个mtDNA 世系。就瘤牛世系而言，昭通牛在目前已发现的2 个瘤牛mtDNA 世系中都有分布，其中I1 世系占37.76%，I2 世系占1.02%；就普通牛世系而言，昭通牛群体中只含有T2、T3 和T4世系，分别占5.10%、43.88%和12.24%。

表4 云南本地牛群体的mtDNA 世系及占比Tab.4 mtDNA lineages and proportions of the Yunnan native cattle

2.5 单倍型遗传关系的系统发育树与中介网络图

由图3a 可知：36 个单倍型中，H01～H09 属于瘤牛血统，其中H01～H07 和H09 属于瘤牛血统的I1 世系，H08 属于瘤牛血统的I2 世系；H10～H36 属于普通牛血统，其中H11 和H12 属于普通牛血统的T2 世系，H10、H13～ H24、H27～H30 和H34～H36 属于普通牛血统的T3 世系，H25、H26和H31～H33 属于普通牛的T4 世系。昭通牛的群体遗传构成是普通牛和瘤牛mtDNA 世系的混合组成，普通牛和瘤牛世系的占比分别为61.22%和38.78%。中介网络图(图3b)分析表明：普通牛mtDNA 单倍型世系与瘤牛mtDNA 单倍型世系间差异较大，两者间存在至少25 个突变位点。

图3 昭通牛群体单倍型遗传关系的系统发育树(a)和中介网络图(b)Fig.3 Phylogenetic tree (a) and mediation network diagram (b) of genetic relationships among haplotypes in Zhaotong cattle

3 讨论

中国多个地方黄牛群体的mtDNA D-loop 区遗传多样性研究表明：中国地方黄牛群体具有丰富的遗传多样性[2,11,19]。本研究对98 头昭通牛mtDNA D-loop 区序列进行测定，共发现73 个突变位点，定义了36 种单倍型，单倍型多样度为0.880±0.027，核苷酸多样度为0.024 43±0.000 94，群体内遗传距离为0.026±0.004，群体内遗传多样性指数为0.027±0.004。其中，单倍型多样度和核苷酸多样度是衡量群体遗传变异的重要指标[19]，它们的值越大，说明群体的遗传多样性越丰富。核苷酸多样度不受样本量和序列长度的影响，因此是衡量群体遗传多样性水平的首要指标。对中国57 个地方牛群体的mtDNA D-loop 序列进行变异分析发现：所有序列的单倍型多样度为0.904±0.008，核苷酸多样度为0.025 7±0.000 1；其中，26个牛种的单倍型多样度低于昭通牛(0.880)，36个牛种的核苷酸多样度低于昭通牛[2]。LEI 等[5]对33 个中国地方黄牛群体的研究表明：中国黄牛的mtDNA D-loop区序列的单倍型多样度为0.914 8±0.013 1，核苷酸多样度为0.025 8±0.012 6；其中19个黄牛群体的核苷酸多样度小于昭通牛。XIA 等[20]对广西3 个本地牛品种的mtDNA D-loop 序列进行测定及分析发现：隆林牛、南丹牛和涠洲牛的单倍型多样度分别为0.957±0.030、0.876±0.063和0.949±0.037，核苷酸多样性分别为0.004 09±0.001 63、0.006 40±0.001 43 和0.000 24±0.000 03，3 个牛种的核苷酸多样度都低于昭通牛。昭通牛与本研究对比的云南其他4 个地方牛群体相比，核苷酸多样度最高。以上结果表明昭通牛群体的群体遗传多样性在中国地方牛中处于相对较高的水平。

遗传距离可以衡量群体间的遗传关系和遗传差异[21]。本研究表明：昭通牛与迪庆牛间的遗传距离最小(0.024)，与德宏高峰牛间的遗传距离最大(0.031)，表明昭通牛与迪庆牛间的遗传关系最近，而与德宏高峰牛间的遗传关系最远，遗传差异较大。遗传分化指数(Fst值)是估测种群间和种群内遗传相似性的指数，以亚种群间的遗传分化占总的遗传多样性的比例来表示[22]，是群体间遗传分化的度量值。基因流(Nm值)是群体间的基因流动，是影响群体内和群体间遗传变异程度的重要因素[18]。本研究显示：昭通牛与迪庆牛间遗传分化最小，Fst值为0.023 32，Nm值为10.47；与德宏高峰牛间遗传分化最大，Fst值为0.552 10，Nm值为0.20。以上结果表明昭通牛与迪庆牛间的遗传分化最小，而与德宏高峰牛间的遗传分化最大。

本研究表明：昭通牛存在普通牛和瘤牛2 个母系世系，普通牛血统占61.22%，瘤牛血统占38.78%。在遗传组成上，昭通牛与云南其他4 个牛种存在一定的差异。迪庆牛、昭通牛、文山牛和滇中牛的普通牛mtDNA 世系占比分别为73.47%、61.22%、35.85%和35.85%，其瘤牛mtDNA 世系占比分别为26.53%、38.78%、64.15%和64.15%；而德宏高峰牛的血统构成有别于这4 个牛种，仅存在瘤牛血统，占比为100.00%。形成这种遗传组成地理格局的原因可能与普通牛和瘤牛2 种血缘世系在中国的扩散路线及其相遇混合的程度有较大关系。中国云南省的西部和南部与南亚瘤牛的主要扩散地缅甸、老挝和越南等东南亚国家接壤，这一地区的牛种受瘤牛血统的影响较大。就本研究中的5 个牛种而言，德宏高峰牛、文山牛和滇中牛的瘤牛mtDNA 世系占比较高，其分布地与东南亚地区相邻或较近，特别是德宏高峰牛，其分布地与缅甸相邻，离印度(瘤牛分布地)较近，其瘤牛血统占比最高(100.00%)；而分布于云南省北部的迪庆牛和东北部的昭通牛，其普通牛mtDNA 世系占比较高，它们的分布地距离东南亚地区相对较远，与四川和贵州相邻，推测它们受中国内地牛种普通牛血统的影响更大。本研究发现：昭通牛的普通牛T2 世系占5.10%，T3 世系占43.88%，T4 世系占12.24%；瘤牛I1 世系占37.76%，I2 占1.02%。推测昭通牛中的普通牛T2、T3 和T4 世系可能是在中原文化向云南的传播过程中由中国北方牛种混入，这些普通牛mtDNA 世系最初起源于近东、欧洲和东北亚地区，而昭通牛中的瘤牛I1 和 I2 世系则来自由东南亚地区传入的南亚瘤牛血统。在一系列的文化传播与交流中，来自2 个方向的普通牛血统和瘤牛血统在当地汇合，经过长期的人工选育和风土驯化，逐渐形成了现在的昭通牛。

周艳等[10]和JIA 等[1]分别对7 头和21 头昭通牛mtDNA D-loop 区序列进行了测定，结果显示：序列的单倍型多样度为0.952±0.096，核苷酸多样度为0.025 04；普通牛和瘤牛血统分别占71.43%和28.57%。LI 等[11]对云南地方牛mtDNA 多样性及系统地理结构的研究结果显示：昭通牛群体的单倍型多样度为0.974，核苷酸多样度为0.025 47；普通牛血统和瘤牛血统分别占60.9%和39.1%。与以上结果相比，本研究揭示目前的昭通牛群体遗传多样性水平有所下降，群体中普通牛和瘤牛mtDNA 世系占比有一定变化，普通牛血统比例下降，瘤牛血统比例上升，说明该牛种近年可能存在更多瘤牛血统的基因渗入。昭通牛与瘤牛血统占比较高的云南文山牛、滇中牛群体间基因流较大，该牛近年是否与这2 个牛种或云南其他牛种有过基因交流，导致其瘤牛血统比例上升，有必要进一步研究。

4 结论

昭通牛群体遗传多样性丰富，存在普通牛和瘤牛2 个母系起源，其血统占比分别为61.22%和38.78%，血统构成有别于云南其他地方牛种。本研究推测昭通牛群体中的普通牛T2、T3 和T4 世系可能是在中原文化的传播过程中由中国北方牛种中的普通牛世系混入所致；而瘤牛I1 和I2 世系则来自东南亚地区传入的瘤牛血统。现今的昭通牛与10 年前的昭通牛在血统构成上有较大差异，推测近年可能存在一定瘤牛血统的基因渗入。