者玉琦,武志娟,王吉坤,钟金城,柴志欣,信金伟2
(1. 青藏高原动物遗传资源保护与利用四川省重点实验室,四川 成都 610225;2. 省部共建青稞和牦牛种质资源与遗传改良国家重点实验室,西藏 拉萨 850000;3. 四川省青藏高原草食家畜工程技术中心,四川 成都 610225)
牦牛是青藏高原的特有畜种资源,在我国主要分布于青海、西藏、四川、甘肃、云南、新疆等海拔3500~5000 m的高寒草原草地。西藏自治区作为我国第二大牦牛养殖和主产区,因其特殊的地理特征及生态环境形成了斯布牦牛、江达牦牛、嘉黎牦牛、类乌齐牦牛、桑日牦牛、工布江达牦牛等地方品种或群体,牦牛遗传资源十分丰富。西藏自治区也成为目前我国地方牦牛品种及遗传资源发掘潜力最大的区域,其中西藏高山牦牛、斯布牦牛、帕里牦牛、娘亚牦牛、类乌齐牦牛及查吾拉牦牛已列入国家畜禽品种遗传资源名录,江达牦牛也于2021 年提交申报国家优良牦牛遗传资源。牦牛是经长期自然选择和人工选择形成的能适应高寒、缺氧、强紫外线、牧草期短等极端环境条件的特色牛种,在遗传上是一个极为宝贵的基因库,因此应该在正确分类以及认识的基础上对其合理利用,利用特色优势资源带动畜产经济发展[1]。
线粒体DNA(mitochondrial DNA, mtDNA)作为高等动物唯一的核外遗传物质,因其具有母系遗传、多态性丰富、结构简单等特点被认为是研究物种起源、演化和分类最好的分子遗传标记[2]。动物mtDNA 包含细胞色素C 氧化酶(cytochrome C oxidase,COX)3 个亚基基因(COX1,COX2,COX3)、Cytb、ATP6、ATP8等13 个蛋白基因,它们都是线粒体内膜呼吸链的重要组成因子。细胞色素C 氧化酶是由细胞核基因组和线粒体基因组分别编码的亚基共同组成的复杂复合物[3-4],COX3是线粒体呼吸链电子传递的终末复合物,属亚铁血红素,铜氧化酶的超家族,是线粒体氧化能力的关键调节因子。已有研究表明,包括COX3在内的细胞色素C 氧化酶相对保守,常用来研究物种分子系统演化及分类,被认为是对属内不同种、种内不同亚种或不同地理型之间物种鉴定较为有效的分子标记[5-6]。目前,研究者基于线粒体COX3基因已对四川黑熊[7]、瓢鸡[8]、长爪沙鼠[9]等物种开展了种系发育研究,获得了比较有参考价值的分类结果,但利用COX3基因对牦牛遗传多样性及分子进化方面的研究还鲜有报道。本研究选取西藏7 个特色牦牛群体共140 头个体,测定COX3基因蛋白质编码区(coding sequence,CDS)序列,分析其遗传多样性,探讨不同群体间的系统进化关系,旨在为西藏牦牛遗传资源保护及利用提供参考。
2022 年5 月样品采集于西藏自治区,相关试验于2022 年7-9 月在西南民族大学青藏高原动物遗传资源保护与利用四川省重点实验室完成。
从西藏帕里牦牛、嘉黎牦牛、类乌齐牦牛、工布江达牦牛、斯布牦牛、桑日牦牛和江达牦牛中分别选取成年健康牦牛20 头,共140 头(表1),用耳样采集器采集长宽各1 cm 耳组织,75%乙醇保存带回实验室,于-80 ℃保存备用。
表1 样本信息Table 1 Sample information
剪绿豆大小耳组织样品去除毛发后,利用动物组织基因组DNA 提取试剂盒(TianGen 生物技术公司)提取其基因组DNA,用1%琼脂糖凝胶电泳、超微量核酸蛋白测定仪(NanoDrop One,美国)分别检测DNA 的纯度和浓度,检测符合要求后置于-20 ℃保存。
根据GenBank 中公布的牦牛COX3基因CDS 核苷酸序列(NC_006380.3),利用Primer Premier 5.0 设计特异性引物,引物序列为:F:5′-TAATCGGAGGGGCTACACTT-3′,R:5′-TCGTATGGACTTGTCTTCTCA-3′。PCR 反应体系为25 μL,分别加入上、下游引物各1 μL(10 pmol·L-1),DNA 模板1 μL,2×Long Taq DNA 预混酶12.5 μL,灭菌ddH2O 9.5 μL。PCR 扩增体系:94 ℃预变性4 min,94 ℃变性45 s,54 ℃退火35 s,72 ℃延伸1 min,共32 个循环;72 ℃延伸10 min,4 ℃保存。将PCR 扩增产物经1%琼脂糖凝胶电泳分离,用胶回收试剂盒(OMEGA)回收纯化,纯化后的PCR 产物(每个体3 个重复)送北京擎科生物科技有限公司成都分公司测序。
利用DNAMAN 软件对所测序列进行比对、校正获得西藏牦牛mtDNACOX3基因CDS 序列;采用鼠小弟(https://cloud.keyandaydayup.com/)统计基因序列长度及碱基组成;利用DNASP 5.1 软件进行多态位点、单倍型多样性和核苷酸多样性分析;根据西藏牦牛COX3基因CDS 单倍型序列的变异位点,利用MEGA 7.0 软件计算群体间的Kimura 双参数遗传距离,构建群体间聚类关系,并绘制单倍型网络关系图,采用邻接法(neighbor joining,NJ)构建牦牛系统发育树。
PCR 扩增产物经凝胶成像系统检测可见条带单一,清晰明亮,大小为1200 bp 左右,且无杂带及拖带(图1),符合后续试验需要。
图1 江达牦牛COX3 基因PCR 产物电泳Fig. 1 Electrophoresis of PCR product of COX3 gene of Jiangda yak
2.2.1 核苷酸及氨基酸组成 经测序获得7 个西藏牦牛群体共140 头个体的COX3基因CDS 序列,其长度均为781 bp,T、C、A、G 4种核苷酸的平均比例分别为29.2%(28.9%~29.6%)、29.4%(28.9%~29.6%)、26.2%(26.0%~26.5%)、15.2%(15.0%~15.2%),其中A+T 含 量55.4%,G+C 含 量44.6%,A+T 碱 基含量占比较大。经DNASP 5.1 软件统计分析显示,所有牦牛群体序列共有55 个多态位点,其中单一多态位点5个(位于144、217、757、765、769 位碱基处),约占多态位点总数的9.09%,简约信息位点50 个(位于27、36、51、64、66、69、91、94、97、102、108、119、120、121、123、127、150、177、181、186、198、252、255、273、300、309、327、330、366、375、417、451、453、462、471、510、519、540、549、555、574、576、582、594、621、633、657、667、681、744 位碱基处),约占多态位点总数的90.91%。同义突变45 个,序列变异中包含碱基替换11 个,无碱基缺失、插入(表2)。
表2 西藏牦牛遗传多样性Table 2 Genetic diversity of Tibetan yak
经分析得到7 个西藏牦牛群体COX3基因不同氨基酸的平均含量(表3),可以看出牦牛富含多种氨基酸,其中缬氨酸与亮氨酸平均含量较高,色氨酸平均含量最低(0.40%)。酸性氨基酸、碱性氨基酸的含量分别为4.84%、9.68%;亲水性氨基酸、疏水性氨基酸的含量分别为35.08%、59.35%。
表3 西藏牦牛COX3 基因的氨基酸组成Table 3 Amino acid composition of COX3 gene of Tibetan yak (%)
2.2.2 核苷酸多样性及单倍型组成 利用DNASP 软件统计分析牦牛群体的核苷酸多样性(表4)和单倍型多样性(表5)。西藏牦牛COX3基因CDS 区共有11 种单倍型(表6),平均单倍型多样性(average haplotype diversity,Hd)、平均核苷酸多样性(average nucleotide diversity,Pi)分别为0.665 和0.00480,牦牛群体单倍型多样性变化范围为0.100~0.858;其中类乌齐牦牛最高,为0.858,嘉黎牦牛最低,为0.100;核苷酸多样性为0.00013~0.01802;其中类乌齐牦牛最高,为0.01802,嘉黎牦牛最低,为0.00013。综合核苷酸多样性、单倍型多样性以及单倍型数量来看,西藏牦牛群体遗传多样性丰富,且类乌齐牦牛的遗传多样性最高,而嘉黎牦牛的遗传多样性最低(表4 和表5)。西藏牦牛群体的Tajima’s D 值检测结果显示,工布江达牦牛和桑日牦牛值为正数,其余牦牛群体均为负值,此外除帕里牦牛外,其余牦牛群体差异均不显著(P>0.10)(表4)。
表4 西藏牦牛COX3 基因CDS 区核苷酸多样性Table 4 Nucleotide diversity in CDS region of COX3 of Tibetan yak
表5 西藏牦牛群体的单倍型多样性Table 5 Haplotype diversity of Tibetan yak groups
表6 西藏牦牛群体的单倍型数量分布Table 6 Haplotype distribution of Tibetan yak groups
统计分析各群体单倍型数量及分布情况(表6),Hap_2,Hap_1 及Hap_4 数量较多,分别为72、30 和23 个,几乎分布于本研究的所有牦牛群体,其余单倍型占比较小,仅分布于少数几个牦牛群体。类乌齐牦牛单倍型最多,共享8 种单倍型,其次是帕里牦牛,与其他地方品种共享6 种单倍型。
2.2.3 西藏牦牛群体间的遗传距离及聚类分析 根据西藏牦牛mtDNACOX3基因CDS 区碱基序列,利用MEGA 7.0 基于Kimura-2-parameter 计算7 个牦牛群体间的遗传距离。群体间遗传距离的变异范围为0.001~0.014,类乌齐牦牛和帕里牦牛间遗传距离最大(表7),但本研究7 个西藏牦牛群体间的遗传距离总体上均较小,表明它们之间的亲缘关系较近,且存在基因交流。
表7 西藏牦牛群体间Kimura 双参数遗传距离Table 7 Kimura two-parameter genetic distance among Tibetan yak groups
根据双参数遗传距离,对7 个西藏牦牛群体进行聚类分析,构建系统发育树(图2)。结果表明,斯布牦牛与桑日牦牛、工布江达牦牛、帕里牦牛、嘉黎牦牛及江达牦牛聚为一类,类乌齐牦牛单独聚为一类。
图2 基于Kimura 双参数遗传距离的西藏牦牛类群间NJ 聚类关系Fig. 2 NJ tree based on Kimura 2-parameter distance among groups of Tibetan yak
2.2.4 系统发育树及单倍型网络关系 根据单倍型序列的变异位点,运用MEGA 7.0 软件绘制西藏牦牛COX3基因CDS 区单倍型网络关系图及系统发育树。11 种单倍型可分为2 个聚类簇,说明本研究的7个牦牛群体有两个母系起源(图3)。聚类簇Ⅰ包含9种单倍型,占全部单倍型的81.8%,涵盖本研究中大部分牦牛群体,聚类簇Ⅱ包括2 种单倍型,占单倍型总数的18.2%。对单倍型个体构建系统发育树(图4),发 现Hap_1、Hap_2、Hap_3、Hap_6、Hap_7、Hap_8、Hap_11 聚为一类后才与Hap_4 和Hap_10 聚为一大类,Hap_5 和Hap_9 聚为一类,其进化关系与单倍型聚类及分布状态一致。
图3 西藏牦牛mtDNA COX3 基因单倍型网络关系Fig. 3 Network based on mtDNA COX3gene haplotype of Tibetan yak
图4 西藏牦牛11 种单倍型的系统发育树Fig. 4 Phylogenetic trees based on 11 haplotypes of Tibetan yak
2.2.5 不同牛种间单倍型网络进化关系 将从GenBank 下载的6 条其他牛种(家牦牛、野牦牛、印度野牛、白肢野牛、普通牛、原牛)的COX3基因序列作为对照,利用MEGA 7.0 软件构建了13 个群体COX3基因以及11 个单倍型间的系统发育树。江达牦牛首先与家牦牛聚为一类,接着与本研究的其他牦牛群体聚为一大类,然后又与野牦牛聚为一大类,最后才与印度野牛、白肢野牛聚为一类(图5)。本研究中的西藏牦牛基本上可划分到家牦牛、原牛和普通牛三大单倍型群体中,其中8 个单倍型属于家牦牛支系,2 个属于原牛和普通牛支系,表明家牦牛、原牛和普通牛是西藏牦牛的混合母系起源,但受家牦牛影响较大,与NJ系统发育树聚类结果一致(图6)。
图5 不同牛种COX3 基因的NJ 系统发育树Fig. 5 Phylogenetic trees of COX3 gene based on NJ analysis in different cattle breeds
图6 不同牛种COX3 基因单倍型的NJ 系统发育树Fig. 6 Phylogenetic trees of COX3haplotype based on NJ analysis in different cattle breeds
本研究测定了140 条西藏牦牛mtDNACOX3基因CDS 区序列,长度均为781 bp,序列中A+T 含量较G+C 含量高,说明西藏牦牛mtDNACOX3富含碱基AT,具有一定的碱基偏好性,符合大多数动物mtDNA AT 碱基含量偏高的普遍特征,属于该基因序列中的优势碱基组成,这与杜玉杰等[7]对四川黑熊COXⅢ基因(A+T 含量56.4%),周武等[10]对豫西黑猪COXⅠ基因(A+T 含量59.3%)及谢雯琴等[11]对蛇源裂头蚴COX3基因(A+T 含量65.13%~66.37%)A、T 碱基偏好的研究结果一致,这一结果与一般脊椎动物mtDNA 中 AT 含量在50%~63%碱基组成相吻合[12]。同时,相较于其他3 种碱基,G 碱基的含量最少(15.2%),表现出明显的非G 碱基偏好性,这种非G 偏爱(anti-ghias)的特性是后生动物线粒体基因组L-链碱基组成所独有的特点[13]。本研究西藏牦牛COX3基因富含多种氨基酸,其氨基酸组成与蛋白质的质量密切相关,并受多种因素影响。动物机体中氨基酸平衡程度会直接影响动物的采食量,且动物体内氨基酸供给处于平衡会降低动物尿氮损失从而提高氨基酸的利用效率,氨基酸不平衡会降低动物的生长率[14]。缬氨酸与亮氨酸是影响牦牛肉质及牦牛乳组成特性的主要氨基酸[15-16],且缬氨酸与肉品香味有关[17],本研究中缬氨酸与亮氨酸含量较高,从侧面反映出西藏牦牛肉质风味鲜美且乳成分营养价值较高。此外,通过对西藏牦牛COX3序列的比较分析,共发现11 种单倍型,55 个多态位点,群体间多态位点较丰富,变异较高的原因可能与牦牛所处环境有关,紫外线可能是诱发突变的主要原因,亦或与部分牦牛群体受到的选择压力较小,积累的突变较多有关。通常认为群体中的单倍型多样性(Hd)、核苷酸多样性(Pi)是衡量一个群体线粒体DNA 变异程度的两个关键指标,Hd和Pi值越大,表明群体的遗传多态性越丰富[18],本研究中西藏7 个牦牛群体的Hd和Pi分别为0.665 和0.00480,略低于姬秋梅等[1]对西藏牦牛Cytb基因的研究结果(Hd=0.884),但高于赵上娟等[2]对西藏牦牛CO Ⅲ基因(Pi=0.0018)以及李双等[19]对温岭高峰牛线粒体DNA 全基因组(Hd=0.778,Pi=0.0017)的研究结果,可能与所研究的牦牛群体数量及参照基因不同有关,但类乌齐牦牛的核苷酸多样性和单倍型多样性均高于其他牦牛群体,遗传多样性最丰富。Tajima’s D 值检测结果显示,工布江达牦牛与桑日牦牛值为正数,说明其序列进化方式为平衡选择,存在部分单倍型分化及群体收缩现象;其余牦牛群体均为负值,表明这些牦牛群体在进化过程中存在负向选择或群体扩张现象。嘉黎牦牛、类乌齐牦牛、斯布牦牛以及江达牦牛群体的Tajima’s D 值均为负值,且中性检验结果不显著,表明以上牦牛群体在进化过程中受到自然或者人工选择后,群体发生了扩张现象,同时线粒体DNA序列在进化上遵循中性模型。以上结果均表明,西藏牦牛COX3基因序列的核苷酸替代率高,遗传多样性较为丰富。
本研究根据Kimura 双参数遗传距离对西藏牦牛群体进行聚类分析,根据NJ 聚类关系,可将西藏7 个牦牛群体分为两大类,即类乌齐牦牛单独聚为一类,其余牦牛群体则聚为一类。相关研究表明,西藏牦牛群体的划分与地理位置并不呈现明显的相关性,不能单一地通过区域分布划分,可能与牦牛的“单一起源”有关[20]。类乌齐牦牛在本研究7 个牦牛群体中遗传多样性最丰富,其次是帕里牦牛,此外,类乌齐牦牛也是本研究中单倍型最多的群体,结合西藏发展史和考古学等资料[21-22],类乌齐牦牛作为昌都市具有代表性的牦牛品种,随着古羌人和古鲜卑人的融合,迁移以及现代西藏各民族的迁移,类乌齐牦牛也逐步走向青藏高原各地。本研究两大分类与Chai等[23]利用全基因组测序对全国32 个家牦牛群体的进化和分化研究的分类结果一致,帕里牦牛、嘉黎牦牛、工布江达牦牛、斯布牦牛、桑日牦牛、江达牦牛聚为一类,主要分布于西藏的核心区域,包括西藏的中部、西部和南部;类乌齐牦牛聚为一类,主要分布于青藏高原的边缘地带。与姬秋梅等[1]基于mtDNAcytb基因及张成福等[24]根据mtDNA D-loop 区的分类结果因牦牛群体数量不同而存在一定差异。此外,帕里牦牛与其他地方群体共享6 种单倍型,仅次于类乌齐牦牛,表明帕里牦牛与其他地方品种或群体有着广泛的基因交流,这与帕里牦牛作为西藏传统的三大优势牦牛群体(嘉黎牦牛、帕里牦牛和斯布牦牛)之一,分布于西藏牧业发展的腹心地带,牦牛引进覆盖着周边地区的实际情况一致。单倍型进化关系对比可以看出家牦牛、原牛和普通牛是本研究7 个西藏牦牛群体的混合母系起源,但受家牦牛影响较大,用最大似然法构建的西藏牦牛系统进化树与11 种单倍型的系统发育树结果一致。
此外,研究者对mtDNACOX3基因作为物种系统分类依据的优劣评价有所不同,可能与所研究的物种不同有关[5],也可能是由自然或人工选择压力下,不同DNA 序列的进化速率存在一定差异所致[25-26],有待进一步的深入研究。本研究也为mtDNACOX3基因作为系统分类候选标记基因提供了一定参考。
牦牛作为“全能型”家畜,可为人们提供肉、乳、皮、毛、役、燃料等畜产品,是青藏高原主要的经济畜种。因此对西藏不同地区尤其是优势种的数量、性状指标等进行系统研究分析, 开展特色牦牛群体的种质资源评价,适当开展优良品种选育、杂交优势利用等对牦牛遗传多样性的合理开发和可持续性利用尤为重要。
本研究西藏7 个牦牛群体mtDNACOX3基因CDS 区序列遗传多样性丰富,其中类乌齐牦牛遗传多样性最高;单倍型可聚为2 簇,存在两个母系起源;7 个牦牛群体可分为2 大类,即类乌齐牦牛(LWQ)为一类,其余牦牛为一类。西藏牦牛群体可划分到家牦牛、原牛和普通牛三大单倍型群体中,其中8 个单倍型属于家牦牛支系,2 个属于原牛和普通牛支系。