全基因组重测序解析秦川牛保种群遗传多样性和遗传结构

马钧，樊安平，王武生，张金川，江晓军，马瑞军，贾社强，刘飞，雷初朝，黄永震

研究报告

全基因组重测序解析秦川牛保种群遗传多样性和遗传结构

马钧1，樊安平2，王武生2，张金川2，江晓军2，马瑞军2，贾社强2，刘飞2，雷初朝1，黄永震1

1. 西北农林科技大学动物科技学院，杨凌 712100 2. 陕西省农牧良种场，宝鸡 722203

在畜禽资源保护中，群体遗传多样性和遗传结构是决定保种效果的重要因素。本研究采用全基因组重测序技术检测100头秦川牛(30头公牛、70头母牛)的基因组变异，通过分析群体遗传多样性、连续纯合片段(runs of homozygosity，ROH)分布特征、亲缘关系和家系结构，对秦川牛的保种效果进行了综合评估。结果显示，100头秦川牛共检测到20,968,017个高质量SNPs位点，平均最小等位基因频率为0.191±0.124，平均多态信息含量为0.279±0.131，平均观察杂合度为0.275±0.131，平均期望杂合度为0.279±0.131，表明秦川保种群遗传多样性较为丰富。秦川牛保种群体平均状态同源(identity by state，IBS)遗传距离为0.243±0.020，其中公牛为0.242±0.021，亲缘关系G矩阵结果与IBS距离矩阵结果一致，均显示秦川牛保种群部分个体间亲缘关系较近。100头秦川牛个体共检测到8258个基因组ROH，ROH总长度为9.64 GB，平均ROH长度为1.167±1.203 Mb，69.35%的ROH是长度为0.5～1 Mb的短ROH。个体平均ROH总长度为96.40 Mb。基于ROH的平均近交系数为0.039±0.039，其中30头秦川公牛的平均近交系数为0.044±0.035，表明部分公牛个体存在一定程度的近交积累。进化树结果显示，秦川牛保种群所测个体可分为8个家系，包括7个含公牛家系和1个不含公牛家系。本研究表明，秦川牛保种群的遗传多样性较为丰富，未出现较大程度近交积累，但部分个体间存在近交风险，应强化选配以确保秦川牛资源的可持续发展。

秦川牛；全基因组重测序；遗传多样性；遗传结构；亲缘关系

根据《国家畜禽遗传资源品种名录(2021年版)》，我国目前拥有55个地方普通牛品种，是世界上牛种资源最多的国家之一。然而，近年来引进国外优良牛种杂交改良本土黄牛，虽然在一定程度上提高了国内牛肉生产水平，但也加剧了本土纯种黄牛数量的减少。畜禽种质资源是国家的战略性资源，做好地方品种的保种工作对于落实种业振兴战略至关重要。秦川牛是我国优良地方黄牛品种，具有体躯高大、适应性好、抗病力强、肉质细致、瘦肉率高及大理石纹品质优良等特点，是一个种质特色鲜明的优良地方品种[1]。然而，目前秦川牛保种群体的遗传多样性及其家系结构尚不清楚，严重阻碍了秦川牛种质资源的保护和利用，亟需相关研究为该群体保种及选配工作的科学开展提供理论指导。

遗传多样性是生物适应环境与进化的基石，也是评估种质资源现状的重要参考指标[2]。群体的遗传多样性越高，适应环境变化的能力就越强[3]。衡量群体遗传多样性有助于了解群体的种质资源现状和未来趋势，是畜禽种质资源保护和开发利用的关键[4]。评估基因组亲缘关系及近交系数有助于了解群体内亲缘关系和近交水平，为制定合理的选配计划提供参考。连续纯合片段(runs of homozygosity，ROH)是生物基因组中纯合基因型的连续片段，是子代继承了亲本的相同单倍型而形成的[5]。研究表明，基因重组、近亲交配、自然和人工选择会影响ROH在基因组中的大小和分布[6]。ROH分析可用于评估群体遗传多样性[5,7]、解析近交水平[8]、推测近交历史[9]、人工选择的追踪[10]及优化育种规划等[11]。

单核苷酸多态性(single nucleotide polymer­phism，SNP)是一种广泛存在于基因组中的遗传变异，被广泛应用于种质资源评估、遗传多样性分析和系统发育进化等研究[12～14]。目前，常见的SNP分型方法包括全基因组重测序(whole genome sequ­en­cing，WGS)、简化基因组测序、SNP芯片等。在马身猪保种群中使用50K SNP芯片进行遗传多样性及亲缘关系分析，发现该群体的遗传多样性丰富，但近交程度高[15]。在青峪猪保种群体的遗传多样性和家系结构评估中，发现该群体在闭锁繁育中出现了遗传多样性损失[16]。在合川黑猪保种群中，发现该群体遗传多样性较丰富，但部分个体间存在较大的近交风险[17]。利用GBS(genotyping by sequencing)简化基因组测序技术评估麦洼牦牛3个保种群的保种效果，表明现有保种策略是可行的，研究结果为麦洼牦牛的保护和利用提供了一定遗传学参考[18]。师睿等[19]基于SNP芯片数据分析新疆近交牛的基因组纯合程度，为该种质的育种规划及后续开发利用提供了一定参考依据。李隐侠等[20]基于绵羊SNP 50K v3芯片对柯尔克孜羊进行遗传多样性和遗传结构分析，结果显示柯尔克孜羊遗传多样性较丰富，群体近交程度低。马克岩等[21]基于简化基因组测序技术，发现永登七山羊的遗传多样性较低，与其他3个甘肃地方绵羊品种(兰州大尾羊、滩羊和岷县黑裘皮羊)有明显分化。

本研究以陕西省农牧良种场—国家秦川牛保种场的100头秦川牛保种个体为研究对象，基于全基因组重测序数据评估了该群体的遗传多样性、近交系数、亲缘关系及家系结构，以解析秦川牛保种群的种质资源现状，旨在为该群体优质保种、科学选配及开发利用提供遗传学依据。

1 材料与方法

1.1 实验材料

本研究选取来自陕西省农牧良种场(国家秦川牛保种场)100头年龄在2～4岁之间的种牛为研究对象，其中公牛30头，母牛70头(表1)。为了方便后续数据分析和对比研究，采用“真实牛号＋性别”的方式标识牛只个体号。所有牛只使用一次性含EDTA抗凝剂的真空采血管颈静脉采血，每头牛采集血液6 mL，与管内EDTA抗凝剂充分混匀后，低温运输至实验室，于-20℃冰箱冷冻保存。血液基因组DNA使用天根生化科技(北京)有限公司的血液基因组DNA提取试剂盒(DP348-03)进行提取，将提取的DNA样品用琼脂糖凝胶电泳分析DNA的纯度和完整性，Nanodrop检测DNA的纯度(260/280比值)，并使用Qubit 2.0(美国Invitrogen公司)对DNA浓度进行精确定量。根据美国Illumina公司推荐的基因组重测序文库构建方法进行文库制备，使用Illumina NovaSeq平台生成150 bp双末端reads。

1.2 基因组SNPs检测

下机的原始数据利用Trimmomatic软件去除接头和低质量reads[22]，然后用BWA-MEM软件比对至牛参考基因组(ARS-UCD 1.2)[23]。随后，利用 SAMtools软件对比对后的bam文件进行排序[24]，排序后使用Picard (https://broadinstitute.github.io/picard/)的MarkDuplicates模块去除重复的reads。利用GATK v4.3.0.0的BaseRecalibrator 模块进行碱基质量重校正[25]。对于校正后的bam文件，利用 Qualimap 2软件统计样本的测序深度[26]。使用GATK的HaplotypeCaller模块分别对每条常染色体进行变异检测，并用GenotypeGVCFs模块合并得到总VCF文件，最后使用SelectVariants模块处理得到SNP数据集。

表1 秦川牛样品采集信息

个体号中的M和F分别表示公牛和母牛。

1.3 基因型数据质控

为得到高质量基因组变异数据，首先利用GATK v4.3.0.0的VariantFiltration模块对SNP进行过滤，过滤参数采用GATK官方推荐的标准：①Variant Confidence/Quality by Depth(QD)<2.0；②RMS Mapping Quality(MQ)< 40.0；③Phred-scaled-value using Fisher’s exact test to detect strand bias(FS)> 60；④Z-score according to the Wilcoxon rank sum test of Alt. Ref read mapping qualities(MQRankSum)< –12.5；⑤Z-score according to the Wilcoxon rank sum test of Alt. Ref read position bias(ReadPos­Rank­Sum)<–8；⑥Symmetric Odds Ratio of 2×2 contin­gency table to detect strand bias(SOR)>3.0[25]。然后，使用BCFtools v1.8过滤得到双等位SNP集合[27]。最后，使用PLINK v1.9软件去除不符合哈代温伯格平衡、MAF值小于0.05和性染色体上的SNP位点，以及基因型缺失大于10%的个体[28]，获得20,968,017个常染色体SNPs位点用于后续分析。最后使用PLINK基于位点间的连锁不平衡进行过滤，参数为“--indep-pairwise 50 5 0.2”，通过过滤后的1,052,677个SNPs位点用于亲缘关系和家系结构分析。

1.4 遗传多样性分析

为了解现有秦川牛保种群体的基因组遗传多样性现状。本研究利用PLINK v1.9软件通过计算期望杂合度(expected heterozygosity，He)、观察杂合度(observed heterozygosity，Ho)、最小等位基因频率(minor allele frequency，MAF)和多态信息含量(polymorphism information content，PIC)，对秦川牛保种群体的遗传多样性进行评估[28]。

1.5 全基因组ROH统计分析

本研究利用PLINK v1.9软件统计每个样本的基因组纯合性片段，参数设置为：滑窗大小设为50个SNP(--homozyg-window-snp 50)，一段ROH中每50 kb必须有1个SNP(--homozyg-density 50)，滑动窗口允许3个SNP位点杂合(--homozyg-window- het 3)，窗口允许缺少数设为5(--homozyg-window- missing 5)[28]。对ROH在秦川牛群体中的分布、长度和数目情况等进行统计分析。为更加了解群体历史，进一步将其根据物理长度分为短片段(0.5～1 Mb)、中片段(1～5 Mb)和长片段(>5 Mb)三类进行分类统计。基于ROH的基因组近交系数(F)由每个个体检测到的常染色体ROH片段总长度占基因组常染色体总长度(2.5 Gb)的比值进行计算[29]，公式为：

1.6 秦川牛保种群体亲缘关系分析

为了评估个体间的亲缘关系，采用IBS遗传距离和G矩阵两种方法进行分析[30]。使用PLINK v1.9软件计算个体间的遗传距离，构建IBS距离矩阵。使用GCTA v1.94.1软件构建个体间的亲缘关系G矩阵。最后，使用R语言绘制热图可视化个体之间的IBS遗传距离和G矩阵结果，并分析保种群个体间的亲缘关系。

1.7 秦川牛群体家系构建

为了解析试验群体的家系结构，采用PHYLIP v3.69软件基于个体之间IBS距离矩阵使用邻接法(neighbor-joining，NJ)构建系统发育树[31]，并使用FigTree v1.3.1软件进行可视化[32]。然后，结合亲缘关系分析结果，以个体间亲缘关系系数大于等于0.1为标准划分家系，分析秦川牛保种群体的家系结构[17,20,30]。

2 结果与分析

2.1 秦川牛基因组遗传多样性分析结果

本研究对100头秦川牛进行基因组重测序，产生22,024,263,265个高质量reads，其中99.84%的reads可比对到参考基因组上，所有样本的平均测序深度达到11.68×。通过GATK v4.3.0.0、BCFtools v1.8和PLINK v1.9等一系列软件对实验样本的基因组变异进行质控，最终获得一个包含20,968,017个SNPs位点的秦川牛高质量基因组SNP变异集合。在常染色体中，1号染色体上鉴定到的SNPs数目最多，为1,316,538个；25号染色体的SNPs数目最少，为373,259个(图1)。功能注释结果如图2所示，这些SNPs主要分布在基因间区(12,424,579，59.25%)或内含子区(7,924,936，37.80%)；外显子仅占总SNPs数量的0.74%(156,185)，包括51,581个非同义SNPs和100,071个同义SNPs。

所有个体的基因型平均检出率为0.997。基于获得的高质量SNPs评估秦川牛保种群体的遗传多样性水平，结果显示，秦川牛群体基因组SNP变异的MAF介于0.05～0.1的比例最高，为31.54%，介于0.4～0.5的占比最低，为9.2%，平均MAF为0.191±0.124(图3A，表2)。基因组SNPs的多态信息含量范围介于0.098～0.500之间，平均多态信息含量为0.279±0.131，具体分布如图3B所示。其中，多态信息含量在0～0.25之间的SNPs有10,145,311个，占48.38%；多态信息含量在0.25～0.5之间的SNPs有10,822,706个，占比51.62%。秦川牛保种群平均观察杂合度为0.275±0.131，平均期望杂合度为0.279±0.131(表2)，平均观察杂合度略小于平均期望杂合度，但二者相差不大，推测其可能受到了选择或近交。

图1 秦川牛基因组SNPs在每条染色体上的分布情况

图2 秦川牛基因组SNPs功能注释

图3 秦川牛保种群体遗传多样性分析结果

A：最小等位基因频率占比图；B：多态信息含量占比图。

表2 秦川牛群体遗传多样性

2.2 秦川牛保种群体IBS距离矩阵分析结果

利用Plink v1.9软件计算秦川牛个体间IBS遗传距离，结果显示100头秦川牛个体间IBS距离值介于0.128～0.278之间，平均IBS值为0.243±0.020，表明个体间的平均遗传距离较远，且差异较大；30头种公牛之间IBS遗传距离在0.128～0.269之间，平均遗传距离为0.242±0.021，表明种公牛之间的遗传距离均较远。秦川牛保种群体IBS距离矩阵可视化结果如图4所示。大部分秦川牛个体间IBS遗传距离较远，呈中等程度的亲缘关系(图4颜色较深的方格)；部分个体间的IBS遗传距离较近，存在较高的亲缘关系(图4颜色较浅的方格)，表明该保种群存在较大的潜在近交风险。

2.3 秦川牛保种群体G矩阵分析结果

基于全基因组SNPs位点构建的G矩阵能真实反映个体间的亲缘关系。利用GCTA软件构建基因组关系G矩阵，进一步分析秦川牛保种群体的亲缘关系。结果如图5所示，颜色越接近红色表示亲缘关系越近。同IBS距离矩阵结果相似，大部分个体间的亲缘关系呈中等程度(图5中颜色较浅的方格)，但部分牛之间的亲缘关系较近(图5中颜色较深的方格)，均表明秦川牛保种群体存在较大的近交风险，必须强化选配措施。

横、纵坐标为本研究所使用的100头秦川牛个体，每一个小方格代表两个个体间的遗传距离值，颜色越接近红色表明两个个体间的遗传距离越大，反之亦然。

2.4 秦川牛保种群体基因组ROH基本统计及近交系数结果

本研究在100头秦川牛保种群体中共检测到8258个ROH片段，ROH总长度为9.64 GB；平均每头秦川牛有82.58个ROH片段，个体平均ROH总长度为96.40 Mb，单个ROH片段平均长度为1.167±1.203 Mb。图6A展示了秦川牛群体单个ROH长度分布情况，发现大部分ROH长度都在4 Mb以下，长度为0.5～1 Mb的ROH最多，占比69.35%，其次是长度为1～2 Mb的ROH，比例为17.89%，大于4 Mb的占比很小。在整个研究群体中，不同牛只检测出的ROH数目变化较大，从最少的25个ROH片段到最多的328个ROH片段。个体基因组ROH总长度范围在21.50～880.90 Mb之间，在100头秦川牛中，79头牛的ROH长度在50～100 MB之间(图6B)。ROH在染色体上的数量分布如图6C所示，在秦川牛的29条常染色体上，1号染色上分布数目最多，为637个，28和29号染色体上分布数目最少，均为85个。从染色体上ROH分布长度看，1号染色体ROH总长度最长，为775.71 Mb，而29号染色体分布的ROH长度最短，为110.07 Mb(图6D)。为深入了解群体历史，进一步将ROH依据物理长度将其分成长、中、短三类，结果如图6E所示。0.5～1 Mb之间的短ROH片段数目和长度分别是5727个和3764.76 Mb，1～5 Mb之间的中ROH片段数目和长度分别为2344个和4594.06 Mb，大于5 Mb的长ROH片段数目和长度分别为187个和1281.38 Mb。短ROH数目在总ROH数目中比例最大，但其长度在总长度中占比并非最大，在ROH总长度中占比最多的是1～5 Mb的ROH。基于ROH的近交系数F分析结果显示(图6F)，现有保种群的个体近交系数范围在0.009～0.352之间，平均近交系数为0.039± 0.039，其中母牛个体号为1294F的个体近交系数最大为0.352。30头秦川公牛F为0.044±0.035，略高于保种群整体平均水平，表明部分公牛个体已出现一定程度的近交积累。

图5 秦川牛保种群G矩阵分析可视化结果

横、纵坐标为本研究所使用的100头秦川牛个体，每一个小方格代表两个个体间的亲缘关系，颜色越接近红色表明两个个体间的亲缘关系越近，反之亦然。

图6 秦川牛基因组ROH统计分析

A：秦川牛基因组单个ROH占比统计；B：秦川牛个体ROH长度的样本数分布；C：ROH在不同染色体上的数量分布情况；D：ROH在不同染色体上的长度分布情况；E：短、中和长三类ROH的数量和长度统计图；F：基于ROH的近交系数可视化。

2.5 秦川牛保种群体家系结构分析

由于种公牛在保种群中的重要性，本研究使用NJ法对30头公牛进行了聚类分析。结合基因组亲缘关系分析结果，以公牛间基因组亲缘系数大于等于0.1为标准进行聚类分析[17,20,30]。结果如图7所示，现有种公牛样本可分为7个家系，分别命名为家系A、家系B、家系C、家系D、家系E、家系F和家系G。进化树中用同一颜色标注的样本被评估为同一家系(表3)。其中，存在部分交叉个体。例如，样本1329M既与家系E的个体亲缘关系较近，又与家系D的部分个体(1278M)亲缘关系较近，所以同时被归并到这两个家系中。最终，这7个家系分别含有10、7、2、3、3、2和4个公牛。根据所有个体间的G矩阵亲缘关系分析结果，将与公牛的亲缘关系值大于0.1的母牛与公牛划分为同一家系，最终将保种群所测个体分为7个含有公牛的家系。在这7个家系外，还有13头母牛与所检测的公牛亲缘关系都比较远，将其划分至1个不包含公牛血缘的家系H(表3)。此外，部分家系中的母牛存在交叉现象，同时被归并到不同的家系中。其中母牛1109F同时被划分到B和E家系中，母牛1409F同时被划分到家系D和E，母牛1158F同时被划分到家系E和F，母牛1185F同时被划分到家系A和G，母牛1144F同时被划分到家系A、B、C和F中。

图7 秦川牛保种群体公牛进化树结果

3 讨论

畜禽种质资源是保障国家畜产品供给的战略性资源，是种业振兴的物质基础。保种工作旨在维护种群内的基因多样性，减少遗传漂变和基因缺失等的不良影响，保障群体的健康发展和适应环境能力，促进畜禽品种改良和创新[33]。科学评估保种效果可以确定现有群体是否存在遗传风险以及采取何种措施进行改善。秦川牛是我国优良的地方黄牛品种，在我国畜牧业发展中发挥着极其重要的作用，但现有保种群的全基因组遗传多样性和群体结构尚未得到系统性研究。本研究采用全基因组重测序技术检测秦川牛保种场的100头秦川牛，鉴定到20,968,017个SNPs位点，并基于这些高质量SNPs从遗传多样性、亲缘关系、近交水平和家系结构等方面对保种群的保种效果进行评估，为秦川牛保种群制定保种计划提供科学依据。

遗传多样性反映了物种或种群在进化过程中的适应性和生存能力，种群内的基因多样性越高，意味着群体更具适应性和生存能力[21,34]。杂合度水平对生物体的机能表现有深远的影响，有研究发现生物体的杂合度与其适应性(例如生长速率、存活率或繁殖力)之间存在正相关[35]，而低杂合度与近亲繁殖抑郁症和适应性降低有关[36]。Gao等[37]基于BovineSNP50芯片发现包括32头秦川牛在内45个世界范围牛品种的Ho为0.142～0.327，其中32头秦川牛的Ho为0.297。本研究中，秦川牛基因组的Ho和He分别为0.275±0.131和0.279±0.131，属于中等杂合度水平，但杂合度有所降低。同时，与Zhang等[38]基于7003个SNPs在30头秦川牛的研究结果(Ho：0.370；He：0.364)也相差较大，这可能与研究所用位点数目、位点选择、采样地区以及样本数目有关。本研究采用全基因组重测序进行分析，相对于商业芯片，获得的基因组SNP变异更加全面。然而，试验所用保种群个体经长期闭群繁育，可能存在遗传多样性有所降低的问题。此外，秦川牛基因组的平均观察杂合度略低于平均期望杂合度，但二者相差不大，表明该保种群具有较低的杂合子丧失，推测其可能受到了一定程度选择压力，但未出现显著的遗传漂移或近交等。最小等位基因频率是对每个群体中不常见等位基因频率的估计[39]。秦川牛群体MAF分布结果显示，该群体低等频率的等位基因数目较多，而较高MAF值的SNPs相对较少，表明该群体的遗传多样性程度相对较高。基因组多态信息含量反映某一位点在群体中的多态程度，常用来评估群体遗传多样性。根据Botstein等[40]报道，对于微卫星位点，PIC值大于0.5时，可被视为高多态性位点；反之，当PIC值小于0.25时，归类为低多态性位点；而当PIC值处于0.25～0.5之间时，则属于中度多态性位点。本研究发现秦川牛保种群基因组SNPs位点的PIC值分布相对均匀，平均PIC值为0.279±0.131，大部分位点PIC值都分布在0.1～0.5之间，其中51.62%的SNPs位点多态信息含量在0.25～0.5之间，说明该群体的遗传多样性较为丰度，具有中度遗传变异，尤其在较高的PIC值区间内。以上结果表明，近年来秦川牛保种群体所采用的保种措施有效，群体遗传多样性并未出现明显丢失。

表3 秦川牛保种群体家系划分结果

个体号中的M和F分别表示公牛和母牛。

ROH在基因组上的分布情况反映了群体遗传多样性、近交历史及近交水平等[5,7～9]，并且研究发现基于ROH计算的近交系数(F)与基于系谱的近交系数(F)具有较高相关性，且比F、F等近交系数计算方法更能反映群体真实近交水平[41]。短ROH表示发生在较远时期的亲缘关系，长ROH反映近期发生的近交[42]。长ROH数目越多，家系内存在近交的可能性越高。本研究在100头秦川牛保种个体中共检测到8258个ROH片段，其中长度为0.5～1 Mb的短ROH数目最多，长度为5 Mb以上的ROH数目最少，揭示秦川牛保种群在较近时期并未出现较大程度近交。基于ROH的近交系数分析显示，秦川牛保种群的近交系数为0.039±0.039，处于相对较低的近交水平，但高于先前研究报道的麦洼牦牛保种群的平均近交系数(0.0168)[18]。30头公牛的平均近交水平略高于整体平均水平，表明部分公牛个体已出现一定程度的近交积累，在后续种公牛选留中建议保留近交水平较低的个体，以避免保种群遗传多样性丧失。此外，研究还发现部分个体具有较高的近交水平，其中母牛1294F个体的近交系数达到0.352，出现较高程度的近交。在后续保种工作中，对于这些近交系数相对较高的个体要适当引入外血，或使用与其关系较远家系的公牛配对以降低其后代的近交系数。

目前，保种场普遍采用闭群繁育方式来保护畜禽遗传资源，种群遗传多样性和群体结构极易受到人工选择和杂交等因素的影响[43]。因此，理清保种群体个体间的亲缘关系和群体遗传结构对于保存群体遗传多样性及保障保种群体的可持续发展至关重要。基于IBS遗传距离和G矩阵的亲缘关系分析结果显示，部分秦川牛个体间的亲缘关系相对较近，这些个体间具有较大的近交风险，在后续配种工作中，要避免这些亲缘关系较近的公母牛进行配对。此外，本研究中秦川牛保种群的平均IBS遗传距离为0.243±0.020，公牛平均IBS遗传距离为0.242±0.021，公牛的平均IBS遗传距离略低于整个秦川牛保种群，这与先前蔡春波等[15]和刘彬等[16]分别在马身猪和青峪猪保种群的研究结果有所不同，这表明部分秦川公牛之间亲缘关系相对较近，因此在后期选留种公牛时需控制好亲缘关系系数。由于家系构成在保存畜禽遗传多样性中的重要作用，基于30头种公牛的IBS遗传距离构建了系统进化树，并以个体间的亲缘关系远近程度，将本研究所测秦川牛划分了7个拥有种公牛的家系和1个不含公共公牛的家系，并且部分家系牛只存在交叉现象。

为了验证基因组分析在保种实践中的可靠性，本研究查阅了保种场的系谱记录。结果发现，母牛1294F的父亲和爷爷是同一头公牛，证实了该牛只基因组近交系数过大是配种失误产生的近交造成的。本研究进一步将家系划分结果与保种场系谱记录进行比较，发现在同一家系中，亲缘关系热图上位置相近且颜色相似的个体大部分是同一头种公牛或该家系其他公牛的后代。此外，一些个体与两个或多个家系的部分牛只都具有较高的亲缘关系值，是由于其拥有相同的母亲或祖代。因此，基于基因组的家系研究结果可有效填补个体系谱记录不全的缺陷，优化保种群的家系构成，并发现畜禽保种中的潜在问题。在后续保种选配过程中，可以结合家系结构和保种场的系谱记录，制定合理的配种方案，减少同一家系公母个体之间的配种，从而科学地开展保种工作，避免种群出现近交和品种内分化。

本研究基于高通量全基因组重测序数据从遗传多样性、ROH分布特征、亲缘关系和家系结构等方面对秦川牛的保种效果进行综合评估，发现秦川牛保种群的遗传多样性较为丰富；保种群体未出现大面积近交，但部分个体间存在较大的近交风险。本研究结果可为秦川牛的保种选配工作提供一定的科学依据。

[1] Wang B, Ma GJ, Zhang Y, Zhou FL, Cheng ZG, Hao P. External factors and improvement model of Qinchuan cattle hybrid improvement., 2020, (21): 39–40.王勃, 马国际, 张岩, 周福来, 程助国, 郝鹏. 秦川牛杂交改良外部因素及改良模式. 畜牧兽医科学(电子版), 2020, (21): 39–40.

[2] Yan BY. Population genomics research on Chaidamu yellow cattle[Dissertation]. Lanzhou University, 2019.闫碧瑶. 柴达木黄牛的群体基因组学研究[学位论文]. 兰州大学, 2019.

[3] Bravo S, Larama G, Quiñones J, Paz E, Rodero E, Sepúlveda N. Genetic diversity and phylogenetic relationship among araucana creole sheep and Spanish sheep breeds., 2019, 172: 23–30.

[4] Jia XB, Ding P, Chen SY, Zhao SK, Wang J, Lai SJ. Analysis of MC1R, MITF, TYR, TYRP1, and MLPH genes polymorphism in four rabbit breeds with different coat colors., 2021, 11(1): 81.

[5] Liu JX, Wei X, Deng TY, Xie R, Han JL, Du LX, Zhao FP, Wang LX. Genome-wide scan for run of homozygosity and identification of corresponding candidate genes in sheep populations., 2019, 50(8): 1554–1566.刘家鑫, 魏霞, 邓天宇, 谢锐, 韩建林, 杜立新, 赵福平, 王立贤. 绵羊全基因组ROH检测及候选基因鉴定. 畜牧兽医学报, 2019, 50(8): 1554–1566.

[6] Broman KW, Weber JL. Long homozygous chromosomal segments in reference families from the centre d'Etude du polymorphisme humain., 1999, 65(6): 1493–1500.

[7] Al-Mamun HA, Clark SA, Kwan P, Gondro C. Genome- wide linkage disequilibrium and genetic diversity in five populations of Australian domestic sheep., 2015, 47: 90.

[8] Purfield DC, Berry DP, McParland S, Bradley DG. Runs of homozygosity and population history in cattle., 2012, 13: 70.

[9] Browning SR, Browning BL. Identity by descent between distant relatives: detection and applications., 2012, 46: 617–633.

[10] Liu G, Sun FZ, Zhu FX, Feng HY, Han X. Runs of homozygosity and its application on livestock genome study., 2019, 41(4): 304–317.刘刚, 孙飞舟, 朱芳贤, 冯海永, 韩旭. 连续性纯合片段在畜禽基因组研究中的应用. 遗传, 2019, 41(4): 304–317.

[11] Marras G, Gaspa G, Sorbolini S, Dimauro C, Ajmone- Marsan P, Valentini A, Williams JL, Macciotta NPP. Analysis of runs of homozygosity and their relationship with inbreeding in five cattle breeds farmed in Italy., 2015, 46(2): 110–121.

[12] Ma J, Gao X, Li JY, Gao HJ, Wang ZZ, Zhang LP, Xu LY, Gao H, Li HW, Wang YH, Zhu B, Cai WT, Wang CY, Chen Y. Assessing the genetic background and selection signatures of Huaxi cattle using high-density SNP array., 2021, 11(12): 3469.

[13] Van der Nest MA, Hlongwane N, Hadebe K, Chan WY, van der Merwe NA, De Vos L, Greyling B, Kooverjee BB, Soma P, Dzomba EF, Bradfield M, Muchadeyi FC. Breed ancestry, divergence, admixture, and selection patterns of the Simbra crossbreed., 2020, 11: 608650.

[14] Zhang SJ, Yao Z, Li XM, Zhang ZJ, Liu X, Yang P, Chen NB, Xia XT, Lyu SJ, Shi QT, Wang EY, Ru BR, Jiang Y, Lei CZ, Chen H, Huang YZ. Assessing genomic diversity and signatures of selection in Pinan cattle using whole-genome sequencing data., 2022, 23(1): 460.

[15] Cai CB, Zhang XL, Zhang WF, Yang Y, Gao PF, Guo XH, Li BG, Cao GQ. Evaluation of genetic structure in Mashen pigs conserved population based on SNP chip., 2021, 52(4): 920–931.蔡春波, 张雪莲, 张万峰, 杨阳, 高鹏飞, 郭晓红, 李步高, 曹果清. 运用SNP芯片评估马身猪保种群体的遗传结构. 畜牧兽医学报, 2021, 52(4): 920–931.

[16] Liu B, Shen LY, Chen Y, Li Q, Liao K, Guo ZY, Zhang SH, Zhu L. Analysis of genetic structure of conservation population in Qingyu pig based on SNP chip., 2020, 51(2): 260–269.刘彬, 沈林園, 陈映, 李强, 廖坤, 郭芝忺, 张顺华, 朱砺. 基于SNP芯片分析青峪猪保种群体的遗传结构. 畜牧兽医学报, 2020, 51(2): 260–269.

[17] Long X, Chen L, Wu PX, Zhang TH, Pan HM, Zhang L, Wang JY, Guo ZY, Chai J. Evaluation of the genetic structure and selection signatures in Hechuan Black pigs conserved population., 2023, 54(5): 1854–1867.龙熙, 陈力, 吴平先, 张廷焕, 潘红梅, 张亮, 王金勇, 郭宗义, 柴捷. 合川黑猪保种群遗传结构及选择信号分析. 畜牧兽医学报, 2023, 54(5): 1854–1867.

[18] Ma SL, Li XW, Li X, Xie SQ, Liu YL, Tang J, Jiang MF. Assessment of genetic structure of 3 Maiwa yak preserved populations based on genotyping-by-sequencing techno­logy., 2022, 31(9): 183–194.马士龙, 李小伟, 李响, 谢书琼, 刘益丽, 唐娇, 江明锋. 基于GBS简化基因组测序评估3个麦洼牦牛保种群的遗传结构研究. 草业学报, 2022, 31(9): 183–194.

[19] Shi R, Zhang Y, Wang YC, Huang T, Lu GC, Yue T, Lu ZX, Huang XX, Wei XP, Feng ST, Chen J, Kagedeer Wulan, Abulizi Ruxianguli, Abulizi Ruxianguli. The evaluation of genomic homozygosity for Xinjiang inbred population by SNP panels., 2020, 42(5): 493–506.师睿, 张毅, 王雅春, 黄涛, 卢国昌, 岳涛, 卢振西, 黄锡霞, 卫新璞, 冯书堂, 陈军, 乌兰·卡格德尔, 茹先古丽·阿不力孜, 努尔胡马尔·木合塔尔. 利用SNP芯片信息评估新疆近交牛基因组纯合度. 遗传, 2020, 42(5): 493–506.

[20] Li YX, Yashengjiang Nasier, Sai Like Duman, Qian Y, Cao SX, Wang WL, Meng CH, Zhang J, Zhang JL. Evaluation of genetic diversity and genetic structure in Kirgiz sheep population cased on SNPs chip., 2023, 54(2): 572–583.李隐侠, 牙生江·那斯尔, 赛里克·都曼, 钱勇, 曹少先, 王伟列, 孟春花, 张俊, 张建丽. SNP芯片评估柯尔克孜羊群体遗传多样性和遗传结构. 畜牧兽医学报, 2023, 54(2): 572–583.

[21] Ma KY, Han JT, Bai YQ, Li TT, Ma YJ. Genetic diversity analysis of Yongdeng seven goats based on specific-locus amplified fragment sequencing., 2023, 54(5): 1939–1950.马克岩, 韩金涛, 白雅琴, 李讨讨, 马友记. 基于简化基因组测序的永登七山羊遗传多样性分析. 畜牧兽医学报, 2023, 54(5): 1939–1950.

[22] Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data., 2014, 30(15): 2114–2120.

[23] Li H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM., 2013.

[24] Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G, Durbin R, 1000 Genome Project Data Processing Subgroup. The Sequence Alignment/Map format and SAMtools., 2009, 25(16): 2078–2079.

[25] McKenna A, Hanna M, Banks E, Sivachenko A, Cibulskis K, Kernytsky A, Garimella K, Altshuler D, Gabriel S, Daly M, DePristo MA. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data., 2010, 20(9): 1297– 1303.

[26] Okonechnikov K, Conesa A, García-Alcalde F. Qualimap 2: advanced multi-sample quality control for high- throughput sequencing data., 2016, 32(2): 292–294.

[27] Li H. A statistical framework for SNP calling, mutation discovery, association mapping and population genetical parameter estimation from sequencing data., 2011, 27(21): 2987–2993.

[28] Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MAR, Bender D, Maller J, Sklar P, de Bakker PIW, Daly MJ, Sham PC. PLINK: a tool set for whole-genome association and population-based linkage analyses., 2007, 81(3): 559–575.

[29] Zhao GY. Association and prediction of traits based on genomic homozygous segments in beef cattle[Disserta­tion]. Chinese Academy of Agricultural Sciences, 2021.赵国耀. 基于肉牛基因组纯合片段的性状关联与预测 [学位论文]. 中国农业科学院, 2021.

[30] Liu CL, Lu D, Zhou QY, Wan MC, Hao XD, Zhang XH, Zheng RQ, Ji HY. Analysis of population genetic structure of Hang pigs by high density SNP chips., 2022, 53(8): 2502–2513.刘晨龙, 卢丹, 周泉勇, 万明春, 郝晓东, 张献贺, 郑瑞强, 季华员. 利用高密度SNP芯片分析杭猪的群体遗传结构. 畜牧兽医学报, 2022, 53(8): 2502–2513.

[31] Baum BR. PHYLIP: phylogeny inference package. Version 3.2. Joel felsenstein., 1989, 64(4): 539–541.

[32] Xu LY, Bickhart DM, Cole JB, Schroeder SG, Song JZ, Tassell CPV, Sonstegard TS, Liu GE. Genomic signatures reveal new evidences for selection of important traits in domestic cattle., 2015, 32(3): 711–725.

[33] Edwards CE, Tessier BC, Swift JF, Bassüner B, Linan AG, Albrecht MA, Yatskievych GA. Conservation genetics of the threatened plant species(Missouri bladderpod) reveals strong genetic structure and a possible cryptic species., 2021, 16(3): e0247586.

[34] Liu F, Qu YK, Geng C, Wang AM, Zhang JH, Li JF, Chen KJ, Liu B, Tian HY, Yang WP, Yu YB. Analysis of the population structure and genetic diversity of the red swamp crayfish () in China using SSR markers., 2020, 47: 59–71.

[35] Reed DH, Frankham R. Correlation between fitness and genetic diversity., 2003, 17(1): 230–237.

[36] Eldridge MDB, King JM, Loupis AK, Spencer PBS, Taylor AC, Pope LC, Hall GP. Unprecedented low levels of genetic variation and inbreeding depression in an island population of the Black-Footed Rock-Wallaby., 1999, 13(3): 531–541.

[37] Gao YH, Gautier M, Ding XD, Zhang H, Wang YC, Wang X, Faruque MO, Li JY, Ye SH, Gou X, Han JL, Lenstra JA, Zhang Y. Species composition and environmental adap­tation of indigenous Chinese cattle., 2017, 7(1): 16196.

[38] Zhang WG, Gao X, Zhang Y, Zhao YM, Zhang JB, Jia YT, Zhu B, Xu LY, Zhang LP, Gao HJ, Li JY, Chen Y. Genome-wide assessment of genetic diversity and popu­lation structure insights into admixture and introgression in Chinese indigenous cattle., 2018, 19(1): 114.

[39] Yuan J, Xu GQ, Zhou X, Xu SP, Li S, Li WM, Liu B. SNP chip-based monitoring of population conservation effect of Tongcheng pigs., 2022, 53(08): 2514–2523.袁娇, 徐国强, 周翔, 徐三平, 李胜, 黎望明, 刘榜. 基于SNP芯片监测通城猪的保种效果. 畜牧兽医学报, 2022, 53(08): 2514–2523.

[40] Botstein D, White RL, Skolnick M, Davis RW. Constru­ction of a genetic linkage map in man using restriction fragment length polymorphisms., 1980, 32(3): 314–331.

[41] Forutan M, Ansari Mahyari S, Baes C, Melzer N, Schenkel FS, Sargolzaei M. Inbreeding and runs of homozygo­sity before and after genomic selection in North American Holstein cattle., 2018, 19(1): 98.

[42] Keller MC, Visscher PM, Goddard M. Quantification of inbreeding due to distant ancestors and its detection using dense single nucleotide polymorphism data., 2011, 189: 237–249.

[43] Wang Y, Dong RL, Li X, Cui C, Yu GH. Analysis of the genetic diversity and family structure of the Licha Black pig population on Jiaodong Peninsula, Shandong Province, China., 2022, 12(8): 1045.

Analysis of genetic diversity and genetic structure of Qinchuan cattle conservation population using whole-genome resequencing

Jun Ma1, Anping Fan2, Wusheng Wang2, Jinchuan Zhang2, Xiaojun Jiang2, Ruijun Ma2, Sheqiang Jia2, Fei Liu2, Chuchao Lei1, Yongzhen Huang1

In the conservation of livestock and poultry resources, population genetic diversity and genetic structure of the conservation population are important factors affecting the effectiveness of conservation. In this study, whole-genome resequencing technology was used to detect genomic variation in 100 Qinchuan cattle (30 bulls and 70 cows). By analyzing population genetic diversity, runs of homozygosity (ROH) distribution features, kinship relationships, and family structure, the conservation effectiveness of Qinchuan cattle was comprehensively evaluated. The results showed that a total of 20,968,017 high-quality SNPs were detected in 100 Qinchuan cattle, the average minimum allele frequency was 0.191±0.124, the average polymorphic information content was 0.279±0.131, and the average observed heterozygosity was 0.275±0.131, the average expected heterozygosity is 0.279±0.131, indicating that the genetic diversity of the Qinchuan cattle conservation population is relatively rich. The average identity by state (IBS) distance of the Qinchuan conservation population was 0.243±0.020, with a value of 0.242±0.021 for the bulls. The results of the kinship G-matrix were consistent with the results of the IBS distance matrix, both showing that some individuals in the conservation population had close kinship. A total of 8258 genomic ROH were detected in 100 Qinchuan cattle, with a total length of 9.64 GB. The average length of ROH fragments was 1.167±1.203 Mb, 69.35% of the ROH were short ROH with a length of 0.5～1 Mb, and the average total length of ROH per individual was 96.40 Mb. The average inbreeding coefficient based on ROH was 0.039±0.039, with a value of 0.044±0.035 for the bulls, indicating that some bulls had a certain degree of inbreeding accumulation. The results of the phylogenetic tree combined with kinship analysis showed that the individuals in the Qinchuan cattle conservation population could be divided into eight families, including seven families with bulls and one family without bulls. This study demonstrated that the genetic diversity of the Qinchuan conservation population is relatively rich, with no significant inbreeding accumulation, but there is a risk of inbreeding among some individuals. Therefore, it is necessary to strengthen selection and mating to ensure the sustainable development of Qinchuan cattle resources.

Qinchuan cattle; whole-genome resequencing; genetic diversity; genetic structure; kinship

2023-04-25;

2023-06-11;

2023-06-30

陕西省重点研发计划项目(编号：2023-YBNY-141)和财政部与农业农村部-国家现代农业产业技术体系(编号：CARS-37)资助[Sup­portedby the Key R&D project of Shaanxi Province (2023-YBNY-141) and Ministry of Finance and Ministry of Agriculture and Rural Affairs-National Modern Agricultural Industry Technology System (CARS-37)]

马钧，在读博士研究生，专业方向：动物遗传育种与繁殖。E-mail: Junma96@163.com

黄永震，博士，副教授，研究方向：动物遗传育种与繁殖。E-mail: hyzsci@nwafu.edu.cn

10.16288/j.yczz.23-115

(责任编委: 赵要风)