西藏牛亚科部分群体线粒体DNA 遗传多样性研究

黄 兴，柴志欣，信金伟，王 会，王吉坤，姬秋梅，钟金城

(1.西南民族大学青藏高原研究院，四川 成都 610041；2.西藏自治区农牧科学院省部共建青稞和牦牛种质资源与遗传改良国家重点实验室，西藏 拉萨 850009)

黄牛作为一种优质畜种，为人类提供肉、奶、毛皮、役力等畜产品，中国黄牛品种资源十分丰富[1]，在中国畜牧业生产及地区农业经济发展中发挥重要作用.针对中国黄牛遗传起源、系统进化，国内外诸多专家学者做了大量研究，《中国牛品种志》中按照地理分布将中国黄牛大体划分为南方黄牛、中原黄牛和北方黄牛[2].黄牛在自然选择下又分为有峰和无峰两种形态，有峰牛叫做瘤牛，大量垂皮，瘤驼，具有抗热特性；无峰牛皮肤多毛，能抗寒，所以因气候差异，中国南方多瘤牛，北方多无峰牛.与国外相比，中国瘤牛瘤驼偏小，而印度瘤牛肩峰高大，差异较大[3].西藏地区海拔高，缺氧，生活在该地区的物种适应性极强，以黄牛、奶牛、牦牛为代表的物种具有很高研究价值，当地政府推进黄牛改良项目，不断加大改良力度，使得养牛产业成为农牧民增收致富的重要方式[4]，伴随着改良进程的加速，出现了杂交牛后代适应性差、产奶量低、繁殖率偏下等问题[5]. 因此，通过分析西藏地区各品种牛遗传多样性及系统进化为优良品种选育，遗传稳定性增强，及西藏地区遗传保种工作提供科学依据.

动物遗传多样性和系统进化分析方法较多，主要包括形态、细胞、生化、分子等遗传标记[6]，而动物线粒体DNA(mtDNA)常被作为研究物种遗传进化的一个分子标记[7-9]，它是一种具有简单双链闭合环型的DNA 分子，除细胞核基因组外唯一核外遗传物质，结构稳定属母系遗传，在不同群体内具有广泛的多态性.作为mtDNA 控制区的D-loop 区，具有变异快和多态性丰富等特点，是研究遗传多样性、亚科内物种间系统发育和进化关系的重要区域. 牛mtDNA 由一段910 bp 左右的D-loop 区以及37 个基因构成，总长度约为16 338 bp[10]. Loftus[11]等在两种独立驯养牛起源和分类证据中根据D-loop 区差异将独立驯养牛划分成瘤牛和普通牛两种类型.Bailey[12]等通过分析牛D-loop 区发现世界牛起源地.Yu[13]等在对12 个南方品种牛mtDNA 研究中，揭示了瘤牛和普通牛是南方黄牛的主要起源. 雷初朝[14]等通过8 个黄牛品种mtDNA 的分析明确了普通牛和瘤牛的起源方向.赖松家[15]等通过四川5 个黄牛品种mtDNA D-loop 区的研究，共发现28 种单倍型，由普通黄牛和瘤牛两种类型组成.张桂香[16]等分析中国地方黄牛品种mtDNA 的遗传变异情况，表明中国黄牛品种遗传多样性丰富.汪琦[17]等通过对三江黄牛mtDNA 的研究，发现三江黄牛为普通牛和瘤牛共同进化而来. 除此之外，钟金城[18]等、李祥龙[19]等、黄勇富[20]等对其他物种如牦牛、山羊、猪等也做了相应mtDNA 限制性内切酶片段长度多态性研究. 西藏地区拥有丰富的地方品种，因技术和交通等不利因素影响，不同牛品种间亲缘关系和遗传多样性尚未系统、全面得到阐释. 本研究旨在对西藏地区部分黄牛品种mtDNA D-loop 区进行研究，为西藏黄牛、奶牛遗传多样性评估和选育提供科学依据，所涉及的黄牛群体个体稀少，部分黄牛群体几乎处于濒危状态，工作人员采样困难，未能采集到大范围的样本，故采集样本量较少.

1 材料与方法

1.1 材料

西藏3 个黄牛群体和1 个奶牛品种，分别是林芝察隅镇巴美牛(ZJ)、林芝察隅镇当地土种牛(TZ)、曲尼帕娟姗牛(JS)、亚东县驼峰牛(TF)，均选取4.5 岁健康雌性个体，采集耳组织，75% 乙醇保存，带回实验室，-80 ℃保存备用.样品采样记录见表1.

表1 牛样品采样记录Table 1 Cattle sample sampling record

1.2 基因组DNA 的提取及检测

利用动物组织基因组DNA 提取试剂盒(Tian Gen生物技术公司)提取DNA，紫外分光光度计检测其浓度，1.5%的琼脂糖检测，并置于-20 ℃保存备用.

1.3 引物设计和PCR 扩增

根据雷初朝[14]等公布的特异性引物，对4 个群体牛mtDNA D-loop 区进行扩增，引物由英潍捷基(上海)生物技术有限公司合成. 上游引物为:5’-CTGCAGTCTCACCATCAACC -3’， 下 游 引 物 为: 5’ -GGGGTGTAGATGCTTGC-3’.PCR 扩增体系为25 μL:ddH2O 9.5 μL、上下游引物各1 μL、DNA 模板1 μL、Dream Taq Green PCR Master Mix(2 ×) 12.5 μL.PCR反应程序为:95 ℃预变性4 min；95 ℃变性30 s，57.4℃退火45 s，72 ℃延伸1 min，循环35 次；72 ℃延伸10 min，4 ℃保存.

经1.5%琼脂糖凝胶电泳检测PCR 产物后用凝胶成像系统观察扩增情况，并保存结果. 扩增成功个体送上海生物工程技术有限公司测序.

1.4 数据处理

测序结果应用DNAMAN 软件与原始序列进行比对校正，Clustal X 1. 83、MegAlign、BioEdit 等软件进行序列比对.利用MEGA 5. 0、DnaSP 4. 0 等软件构建系统进化树，以及计算碱基组成、核苷酸多态性、单倍型多态性等遗传指标. 利用MEGA6 软件进行起源进化关系分析，并添加2 头欧洲普通牛(GenBank 登录号:L27712，GenBank 登录号:L27716)，2 头印度瘤牛(GenBank 登录号:L27722，L27732)，1 头非洲瘤牛(GenBank 登录号:L27728). 采用邻接法(NJ 系统发育树).

2 结果

2.1 mtDNA D-loop 区核苷酸序列及其遗传多样性

将测序结果进行人工比对、剪切，得到4 个群体对应的mtDNA D-loop 区，18 个个体D-loop 区长度分布在892 bp ～912 bp，品种间存在长度差异，4 个群体A、G、T、C、A ＋T、C ＋G 的平均含量约为:32.70%、13.80%、28.80%、24.70%、61.50%、38.50%，表明西藏黄牛和奶牛D-loop 序列富含A ＋T 碱基，此外Dloop 区末端发现一段连续的碱基C.

西藏3 个黄牛群体和一个奶牛品种中，共检测到颠换、转换、插入、缺失以及颠换和转换共存同一位点等5 种突变类型(图1).利用生物软件DNASPv5 比对，因序列间对齐间隙或缺少数据等情况存在，共缺失30 bp，总多态位点为125 个，单倍型数量为17，突变总数达到132，平均单倍型多样性约为0.992，平均核苷酸多样性约为0.039，平均核苷酸差异数约为34.767.其中单一多态位点25 个，包含24 个单碱基突变和一个三碱基突变，位置分别位于:25 110 117 135 152 198 259 260 266 286 303 306 336 374 409 410 519 588 712 739 750 831 875 911 和204.简约信息位点为100 个，包含94 个两碱基突变和6 个三碱基突变的简约信息位点，分别位于:1 4 27 28 38 40 55 56 63 80 82 87 90 124 125 130 131 134 148 153 157 158 160 163 165 166 168 172 174 183 203 231 256 258 267 268 277 283 292 294 295 304 312 318 320 322 323 325 326 327 329 331 332 335 347 348 351 357 406 415 416 439 440 443 458 459 466 471 475 505 511 512 513 557 655 715 718 721 722 728 755 773 782 798 799 805 818 829 842 846 850 865 877 879 和162 170 284 319 912 913.

将4 个群体分离，单独比对发现，巴美牛单倍型数为4，单倍型多样性(Hd)为1.000；当地土种牛单倍型数为3，单倍型多样性(Hd)为1.000；娟姗牛单倍型数为5，单倍型多样性(Hd)为1.000；驼峰牛单倍型数为5，TF15 与TF16 共享一种单倍型，单倍型多样性(Hd)为0.933，巴美牛，娟姗牛，当地土种牛单倍多样性最高，驼峰牛单倍型多样性最低，具体见表2.

图1 西藏牛亚科部分群体mtDNA D-loop 序列比对Fig.1 The mtDNA D-loop sequence alignment of some subfamilies of cattle in Tibet

表2 西藏牛亚科部分群体的单倍型多样性Table 2 Haplotype diversity of some subfamilies of cattle in Tibet

巴美牛mtDNA D-loop 区转换45 个，颠换3 个，插入及缺失为3 个，分别占4.9%、0.3%、0.3%，核苷酸多样性(Pi)为0.03528；当地土种牛mtDNA D-loop区转换89 个，颠换15 个，插入及缺失为11 个，并存在一处转换和颠换共存位点，分别占9.9%、1.7%、1.2%、0.1%，核苷酸多样性(Pi)为0.07929；娟姗牛mtDNA D-loop 区转换48 个，颠换4 个，插入及缺失为7 个，并存在1 处转换和颠换共存位点，分别占5.3%、0.4%、0.7%、0.1%，核苷酸多样性(Pi)为0.03665，Tajima's D 中性检验D 值为1.31274(P ＞0.1)，不显著；驼峰牛mtDNA D-loop 区转换3 个，颠换3 个，插入及缺失为4 个，并存在1 处转换和颠换共存位点，分别占0. 3%、0. 3%、0. 4%、0. 1%，核苷酸多样性(Pi)为0.00353，Tajima's D 中性检验D 值为-0.51093(P ＞0.1)，不显著，表明娟姗牛和驼峰牛可能未出现群体扩张模式.4 个群体，核苷酸多样性变化范围在0.00353 ～0.07929，当地土种牛的核苷酸多样性最高，巴美牛和娟姗牛核苷酸多样性差异不大，驼峰牛的核苷酸多样性最低.其结果如表3、表4.

表3 西藏牛亚科部分群体mtDNA D-loop 核苷酸多态性Table 3 Nucleotide polymorphism of mtDNA D-loop of some subfamilies of cattle in Tibet

表4 西藏牛亚科部分群体mtDNA D-loop 区平均核苷酸变异率(%)Table 4 Average nucleotide variation rate of mtDNA D-loop region of some subfamilies of cattle in Tibet(%)

2.2 西藏牛亚科部分群体系统进化关系分析

利用MegAlign 软件对4 个群体mtDNA D-loop 区两两比对发现(见图2)，4 个群体牛歧异度在0.1%～11.9%之间变化，同源性在88.7% ～99.9%；驼峰牛和巴美牛的歧异度最小为0.1%，同源性为99.9%，推测巴美牛与驼峰牛有相似血统；娟姗牛和当地土种牛的歧异度最大，达到11.9%，同源性为88.7%，表明当地土种牛作为当地固有的群体和外来的娟姗牛品种有着明显的遗传分化.

图2 西藏牛亚科部分群体mtDNA D-loop 区歧异度及差异百分比Fig.2 Divergence and percent identity of mtDNA D-loop region of some subfamilies of cattle in Tibet

利用DNAsp 软件进行比对发现，4 个西藏牛亚科群体mtDNA D-loop 区序列由17 种单倍型组成，其中驼峰牛3 号和驼峰4 号共享一种单倍型；其余每一头对应一种单倍型，表明西藏牛亚科群体mtDNA 遗传多样性丰富.由图3 可知，17 种单倍型大致可以分为三大类，第一类包括2 号当地土种牛(TZ2)和3 号当地土种牛(TZ3)；第二类包括3、4 号巴美牛，2、4、5 号娟姗牛和两种印度瘤牛，且发现Sahiwal 牛(L27732)与3、4 号巴美牛单倍型更相似，Hariana 牛(L27722)与2 号娟姗牛单倍型较相似；第三类将剩余单倍型个体和2 头欧洲普通牛和1 头非洲瘤牛划归为一大类，其中非洲瘤牛和娟姗牛3 号单倍型更一致.4 个群体牛系统发育进化树与同源性比对情况基本一致，除巴美牛外，其余群体自成一类. 娟姗牛与当地土种牛以及驼峰牛亲缘关系较远，但巴美牛与娟姗牛、驼峰牛、当地土种牛分歧时间较近.

图3 西藏牛亚科部分群体mtDNA D-loop 全序列NJ 系统发育树(普通牛与瘤牛为对照)Fig.3 NJ tree of mtDNA D-loop region sequences of some subfamilies of cattle in Tibet(Bos taurus and Bos indicus as control)

3 讨论

3.1 mtDNA D-loop 区遗传多样性分析

具有世界屋脊之称的青藏高原生物遗传资源丰富，保护和有效利用青藏高原遗传资源意义重大，而遗传多样性研究有助于了解对不同品种或类群进化历史以及生物地位. 本研究对林芝当地土种牛、巴美牛以及来自曲尼帕的娟姗牛，亚东县驼峰牛4 个群体mtDNA D-loop 区研究发现，所有个体D-loop 区均在910 bp 左右，与Loftus[11]、Bailey[12]、雷初朝[14]、赖松家[15]张桂香[16]等测得黄牛mtDNA D-loop 区长度接近.在D-loop 区末端发现一段连续的碱基C，与Rand[21]、Anderson[10]、Hauswirth[22]等发现脊椎动物mtDNA D-loop 区包含一段连续的碱基为中止序列的结果一致.统计发现，西藏牛亚科部分群体mtDNA Dloop 区富含A ＋T 碱基，这与张桂香[16]、林瑞意[23]分别对我国部分黄牛品种线粒体D-loop 区遗传多样性与起源分化研究和贵州务川黑牛线粒体DNA D-loop全序列分析中发现D-loop 区富含A、T 碱基的结果一致，符合线粒体DNA D-loop 区多A ＋T 碱基的特性.

研究发现，18 个个体共发现17 种单倍型，单倍型数占总检个体数的94.4%，表明其遗传多样性丰富.4个群体突变类型主要以转换为主，大量研究表明线粒体DNA 在进化过程中，转换频率高于颠换.在突变类型中，发现了同一位点存在转换和颠换共存现象，这在雷初朝等对中国部分黄牛品种mtDNA 遗传多态型研究中已报道[14]，说明不同黄牛品种mtDNA 存在进化保守性及突变相似性. 4 个群体发现大量突变位点，推测与西藏海拔、光照强度等自然环境具有一定关系，环境诱变，紫外辐射均可诱发突变，西藏林芝平均海拔3 100 m，西藏亚东县平均海拔在3 400 m，西藏曲尼帕位于拉萨，其平均海拔更是在3 600 m 左右，以上地区气压低，紫外线辐射强，长期生活在此环境，易发生突变积累. 当地土种牛和巴美牛因样品数量较少，未能进行Tajima'D 中性检验，娟姗牛Tajima's D 值为1.31274，不显著(P ＞0.1)，驼峰牛Tajima's D值为-0.51093，不显著(P ＞0.1)，均符合中性进化，未出现过群体扩张或瓶颈效应现象.种群增长曲线表明:食物和空间充足及无敌害等理想条件下，数量会连续增加，根据相应的报道显示雪灾是西藏最常见的自然灾害，常造成牲畜死亡[24-25]，作为家养牛，排除敌害因素，漫长的冬季，难以持续有效补饲可能制约了种群发展.平均核苷酸多样性为0.039，比张桂香[16]研究我国部分黄牛(Pi:0.0227)和汪琦[17]报道三江黄牛(Pi:0.00625)的平均核苷酸多样性均要高，这表明西藏黄牛具有丰富的遗传多样性.

3.2 西藏牛亚科部分群体的系统进化关系

西藏黄牛作为农牧区的主要畜种，体型偏小，生长发育缓慢，生产性能偏低，但抗逆性强，随着选育技术日益成熟，陆续出现一些改良品种，其遗传多样性研究受到人们密切关注.研究调查分析不同品种之间的亲缘关系、起源进化有助于合理利用和保护遗传资源.将所有个体mtDNA D-loop 区进行比对发现，当地土种牛和娟姗牛歧异度最大，同源性仅为88.7%. 娟姗牛作为外来引进品种，而当地土种牛作为当地固有品种，本身亲缘关系较其他远. 以2 头欧洲普通牛，1头非洲瘤牛，2 头印度瘤牛作为对照，采用邻接法绘制系统发育树，4 个品种中娟姗牛和部分巴美牛均和两种瘤牛类型归为一类，西藏其他牛种不含有任一瘤牛单倍型，均为普通牛，这与张桂香[16]、陈智华[26]、罗勇发[27]等证明西藏牛为普通牛类型观点一致. 土种牛2 号和土种牛3 号单独归为一类，表明来自林芝下察隅镇的当地土种牛受到其他牛种影响较小，结合核苷酸多样性和单倍型数分析发现，该品种牛遗传多样性丰富.而巴美牛与其他牛种有一定的亲缘关系，人工选择对其影响较大.

4 结论

本研究从mtDNA 基因组出发，对4 个群体mtDNA D-loop 区进行克隆测序，并初步分析了序列结构，探究了西藏牛亚科部分品种的遗传多样性和进化关系，西藏黄牛具有丰富的遗传多样性，不同品种间遗传关系复杂.要彻底弄清西藏黄牛的起源演化，需大量收集更多地区的黄牛品种，扩大研究范围，进一步深入研究也许才能解决这一问题，使得种质资源得到有效保护和利用.