中国荷斯坦牛是由欧洲荷斯坦牛与中国黄牛杂交培育的奶牛品种，已经在国内饲养60多年^[1]，对其产奶量和乳品质的选择是人们长期的育种目标^[2]。数量性状遗传参数的估计不仅是估计育种值、制订育种计划和了解数量性状遗传机制的必要前提，更是奶牛育种的中心任务^[3]。遗传力作为从数量表型范畴进入遗传范畴的限制因素，对于每一个个体来说至关重要^[4]；完整的系谱记录、足够大的样本量和相对稳定的饲养环境等是准确估计遗传力的必备条件^[5-6]。近年来，国内外有关奶牛遗传参数评估的报道较多，毛杰等^[7]对上海地区荷斯坦奶牛体型性状、产奶性状和体细胞评分进行遗传评估发现，体型性状的遗传力在0.17~0.52之间，产奶性状的遗传力范围是0.13~0.59，体细胞评分遗传力为0.13。任小丽等^[8]研究发现，中国荷斯坦牛头胎测定日产奶量、乳脂率、乳蛋白率和体细胞评分遗传力分别为0.14、0.29、0.37和0.09，其中乳脂率和乳蛋白率估计遗传力与Kramer等^[9]的研究结果相近(瑞士褐牛)，且瑞士褐牛的泌乳速度和乳房深度呈中等遗传力，攻击性为低遗传力性状^[10]。Wethal和Heringstad^[11]对挪威红奶牛进行遗传评估发现，每日盒装时间(BT)和平均挤奶速率(FR)和每分钟挤奶千克数(MEF)为高遗传力性状(0.48、0.27和0.22)，挤奶频率和挤奶间隔的遗传力较低(0.02~0.07)。Burren等^[12]对瑞士奶牛产精性状进行遗传评估发现，射精量、精子浓度和精子活力的估计遗传力分别为0.42、0.25和0.09。随着育种工作的稳步实施及市场需求的逐渐改变，对奶牛育种不仅仅集中在产奶量和乳品质方面，育种者更为关注的是奶牛的性情、寿命和疾病等综合遗传性能的选育，以获得奶牛养殖的最大经济效益^[13]。常瑶等^[14]研究发现，北京地区中国荷斯坦牛性情评分为1分、2分和3分的比例分别为75.00%、19.25%和5.75%，属于低遗传力性状(0.03)，且与重复输精次数为中度遗传相关(r=0.56)，与产犊到首次配种间隔和直肠温度为中度负遗传相关(r =－0.30, r =－0.17)。李想等^[15]通过单性状动物模型和单性状公畜模型估计中国荷斯坦牛生产寿命遗传力分别为0.05和0.04，多性状动物模型估计的遗传力较为稳定，介于0.05~0.07之间。董祎鑫等^[16]收集北京地区27个场中国荷斯坦牛产后0~35 d的25 026条繁殖疾病记录，对子宫疾病进行遗传参数估计发现遗传力仅为0.01。目前，宁夏地区荷斯坦牛遗传评估也已有报道，李欣等^[17]研究发现，宁夏地区荷斯坦牛305 d产奶量和初产日龄(0.31、0.37)为中等遗传力性状；乳脂率、乳蛋白率、体细胞评分、产犊间隔和体型线性评分均为低遗传力性状(0.07~0.19)，与Soufleri等^[18]的研究结果相近(0.04~0.22)；奶牛乳房结构是发挥生产性能的重要基础，有研究报道，奶牛乳房性状与产奶性状的遗传关系密切^[19]，而对于宁夏地区荷斯坦牛乳房性状的遗传参数评估尚未见报道。

宁夏DHI测定中心从2012年以来经过不懈努力，积累了许多宝贵、准确的系谱资料，足够多的样本资料为遗传评估的进行提供了参考，而主成分分析(PCA)用来作为降低数据维数及探索信息损失最小性状之间关系的多元方法。本研究利用多性状动物模型，借助DMU软件及主成分分析方法，在李欣等^[17]的研究基础上对宁夏地区荷斯坦奶牛群体乳房性状进行遗传评估，计算出属于该地区荷斯坦牛的遗传参数；通过分析各个性状的遗传规律，找出各性状间的遗传联系，以期为进一步育种规划作铺垫。

1 材料与方法 1.1 系谱数据来源与整理

本试验收集来自宁夏畜牧工作站2012年12月至2020年11月24个牧场的13 031头荷斯坦奶牛的乳房性状数据，各个牧场的数据均是请专业的体型线性评分人员进行测量、评分和记录。经过源数据筛选，剔除父母号缺失、空值、错误值和极端值后，得到系谱完整、有准确记录的牛只12 415头，共有111 735个乳房性状观测值，表 1为源数据的筛选结果。

1.2 乳房性状影响因素

场(field)：24个场，每个场为1个处理水平；出生年度(birth year)：2007—2018年，共12个水平；出生季节(birth season)：本试验将1年划分为春(3~5月)、夏(6~8月)、秋(9~11月)、冬(12月至次年2月)4个季节，每个季节为1个水平；鉴定月龄(month of age)：每个月龄段为1个水平；胎次(parity)：1~4胎，每胎次为1个处理水平。

1.3 乳房性状影响因素的显著性检验

影响因素的显著性检验使用SAS 9.2的线性混合模型(linear mixed model)：

$ \begin{array}{*{20}{l}} {\;\;\;\;\;\;{Y_{ijklmno}} = \mu + {\rm{ }}\mathit{Far}{\mathit{m}_i} + {\rm{ }}\mathit{Yea}{\mathit{r}_j} + {\rm{ }}\mathit{Seaso}{\mathit{n}_k} + {\rm{ }}\mathit{Mont}{\mathit{h}_l} + }\\ {\mathit{ Parit}{\mathit{y}_m} + {a_{ijklmn}} + {e_{ijklmno}}} \end{array} $

式中，Y_ijklmno为奶牛乳房深度评分/乳房质地评分/悬韧带评分/前乳房附着评分/前乳头位置评分/前乳头长度评分/后乳房高评分/后乳房宽评分/后乳头位置评分；μ为总平均数；Farm_i为出生场次固定效应；Year_j为出生年份固定效应；Season_k为出生季节固定效应；Month_l为生长月龄固定效应；Parity_m是胎次固定效应；a_ijklmn为个体的随机遗传效应；e_ijklmno为随机残差效应。

1.4 乳房性状遗传参数估计模型

使用Microsoft Excel 2016软件进行源数据整理，采用动物模型DMU(v6.0)软件DMUAI模块的AI-REML结合EM算法计算方差组分估计值。规定模型中个体动物效应和随机残差效应为随机效应。对于遗传参数的估计，根据数据的具体情况，配合下列动物模型进行最小二乘方差分析：

$ y=X \beta+Z \mu+e $

式中，y为观察值向量；X为固定效应结构矩阵；β为固定效应向量；Z为随机效应结构矩阵；μ为随机效应向量；e为随机残差向量。

模型的假设前提条件：E(y) = xβ, E(h) = 0, E(α) = 0

Var(h) = Iσ_h², Var(α) = Aσ_α², Var(e) = Iσ_e²

得到的方差组分带入下列公式计算遗传参数：

遗传力$h^{2}=\frac{V_{A}}{V_{A}+V_{E}}$

遗传相关${r_A} = \frac{{\mathit{Cov}\left( {{a_1}, {a_2}} \right)}}{{\sqrt {{\sigma ^2}{a_1}{\sigma ^2}{a_2}} }}$

表型相关${r_p} = \frac{{\mathit{Cov}\left( {{P_1}, {P_2}} \right)}}{{\sqrt {{\sigma ^2}{P_1}{\sigma ^2}{P_2}} }}$

式中，h²为遗传力；V_A为性状的加性遗传方差(育种值方差)；V_E为性状的残差方差；r_A为性状间的遗传相关；r_P为性状间的表型相关；Cov(a₁, a₂)为性状a₁与a₂的加性效应协方差；σ²a₁为性状α₁的加性方差；σ²a₂为性状a₂的加性方差；Cov(P₁, P₂)为性状P₁与P₂的表型协方差；σ²P₁为性状P₁的表型方差；σ²P₂为性状P₂的表型方差。

育种值是指从基因型效应中剔除了显性效应和上位效应后的加性效应部分，在世代传递中能够稳定遗传。本研究采用动物模型BLUP，将动物个体本身的加性遗传效应(即育种值)作为随机效应放在模型中计算。假设个体所考察性状不考虑显性效应和上位效应，则观察值y 通常可以用如下模型来描述：

$ y=\sum\limits_{j=1}^{r} b_{j}+a+p+e $

式中，b_j是第j个系统环境效应，一般为固定效应；a为该个体的加性遗传效应(育种值)，即随机效应；p为随机永久性环境效应；e为残差效应。

1.5 乳房性状主成分分析

主成分分析采用SAS 9.2的Princomp程序完成，前m个主成分的累计贡献率为：($\sum_{i=1}^{m} \lambda_{i} / p$)× 100%。本研究只考虑特征值≥1的主成分。第一个主成分解释了数据集中总方差的最大百分比。第二个主成分解释占比第二的百分比，以此类推，直到解释数据集中的所有差异。

2 结果 2.1 乳房性状基本表型参数

荷斯坦牛乳房性状表型参数统计结果见表 2，线性评定的9个乳房性状平均值在4.99~6.24之间。体型性状评分并非是分值越高越好，不能只单纯考虑其测定值，要与最优分作比对。当理想得分在中间值时，可能会出现理想分两端求平均值得到理想得分的情况，无疑会造成错误，所以本研究是用每个得分与理想得分求出差值的绝对值，然后再求其均值。其中，后乳房高的最优分相差最大，为3.90；悬韧带、乳房质地、前乳房附着和后乳房宽与理想体型得分的差距较大，在2.76~3.71之间，说明宁夏地区荷斯坦牛以上5种乳房性状还不是很理想，之后应重点关注。而前乳头位置(0.78)、前乳头长度(0.64)、后乳头位置(0.75)和乳房深度(0.88)与理论最优分相差较小，说明这4个乳房性状较为理想。

2.2 乳房性状影响因素的显著性检验

由表 3可知，场和出生年度对9个乳房性状影响均极显著(P < 0.01)；出生季节对乳房深度、乳房质地、前乳房附着、前乳头位置、后乳房高、后乳房宽有极显著影响(P < 0.01)，对悬韧带、前乳头长度和后乳头位置影响显著(P < 0.05)；胎次对前乳头长度和后乳房高影响不显著，对其他性状的影响均达到极显著水平(P < 0.01)。鉴定月龄对乳房深度、悬韧带、前乳房附着、前乳头位置、前乳头长度、后乳头位置有显著或极显著影响(P < 0.05或P < 0.01)。因此，运用DMU模型进行每个性状遗传参数计算时需剔除对乳房性状影响不显著的因素。

2.3 乳房性状的估计遗传力

表 4为各性状的加性方差、残差方差及遗传力的估计值。可以看出，宁夏地区荷斯坦牛乳房性状的遗传力在0.080~0.378之间。其中，后乳房宽的遗传力最高(0.378)，其次为后乳房高(0.272)和乳房深度(0.242)；乳房质地的遗传力最低，为0.080；其余乳房性状遗传力均在0.169~0.231之间，总体来说，宁夏地区荷斯坦牛乳房性状具有中等偏低的遗传力。

2.4 各乳房性状间的遗传相关和表型相关

宁夏地区荷斯坦牛乳房性状遗传相关整体在－0.277~0.706之间，表型相关在－0.084 ~ 0.316之间，整体上较遗传相关低(表 5)，但是，前乳房附着和乳房深度、前乳房附着和乳房质地、前乳房附着和悬韧带、乳房质地和悬韧带、后乳头位置和前乳头位置、后乳头位置和悬韧带遗传相关和表型相关的趋势一致，并且均为弱相关或中度相关，相关系数在0.210~0.706之间。其次，乳房深度和悬韧带、乳房质地和前乳头长度、前乳房附着和后乳房宽、前乳房附着和后乳头位置均为中等遗传相关，相关系数分别为0.306、0.443、0.437、0.477。其余性状遗传和表型均为弱相关或负相关。

2.5 各乳房性状的估计育种值

各性状育种值排名前50的奶牛见表 6。乳房深度、悬韧带、前乳房附着、前乳头长度、后乳房宽和后乳头位置育种值排名前50数量最多的牧场均是640190，分别为23、31、30、19、16和37头；乳房质地和后乳房高育种值排名前50数量最多的牧场均是640004，分别为13和49头；在640251牧场中前乳头位置育种值排名前50的数量最多，为35头。对于9个乳房性状整体而言，育种值排名前50的奶牛个体数最多的牧场为640190(163头)，其次为640004(75头)和640251(69头)；640253、641009、640166、640173和640176牧场中没有育种值排名前50的牛。

2.6 乳房性状主成分分析

表 7给出了宁夏地区荷斯坦牛乳房性状的特征值、相邻特征值之差、解释的方差比例及累计贡献率。特征值≥1的4个主成分占总方差的64.19%。在9个乳房性状中，PRIN1占总变异的比例最大(23.39%)，表明第一主成分反映的综合信息含量最多，提取主成分的荷载系数如表 8所示，其信息含量的表达式为：PRIN1 = 0.356(乳房深度)+ 0.305(乳房质地)+0.424(悬韧带)+ 0.467(前乳房附着)+0.297(前乳头位置)+ 0.260(前乳头长度)+ 0.157(后乳房高)+ 0.291(后乳房宽)+0.342(后乳头位置)，RIN1在悬韧带和前乳房附着上具有更高的负荷，可能代表悬韧带和前乳房附着的特征。PRIN2占17.47%的变异性，PRIN2 = －0.223(乳房深度)+ 0.463(乳房质地)－0.180(悬韧带)+ 0.123(前乳房附着)－0.477(前乳头位置)+ 0.143(前乳头长度)+ 0.367(后乳房高)+ 0.428(后乳房宽)－0.350(后乳头位置)，在乳房质地、前乳头位置、后乳房宽方面负荷较大；而PRIN3占12.67%的比例变异性，PRIN3 = －0.407(乳房深度)－0.069(乳房质地)+ 0.230(悬韧带)－0.334(前乳房附着)+ 0.123(前乳头位置)－0.417(前乳头长度)+ 0.480(后乳房高)+ 0.264(后乳房宽)+ 0.421(后乳头位置)，在乳房深度、前乳头长度、后乳房高、后乳头位置方面负荷较高。PRIN4占10.65%的比例变异性，PRIN4 = －0.495(乳房深度)－0.009(乳房质地)+ 0.089(悬韧带)－0.286 (前乳房附着)+ 0.201(前乳头位置)+ 0.651(前乳头长度)－0.360(后乳房高)+ 0.254(后乳房宽)+ 0.082(后乳头位置)，在乳房深度、前乳头长度方面负荷较大。由于特征值小于1的主成分解释的方差比例较小，无法清楚地描述一组性状，因此不考虑进一步分析。

3 讨论 3.1 乳房性状的遗传力分析

遗传力主要反映基因加性效应的大小，是性状、群体和环境的综合体现，因此性状估计遗传力一般均存在一定的偏差，但仍有相对的稳定性，在育种规划中提高遗传力估计值准确性是非常有必要的^[20]。宁夏地区荷斯坦牛乳房性状的遗传力在0.080~0.378之间，整体上属于中等偏低遗传力性状。其中，后乳房宽的遗传力最高(0.378)，为中等遗传力，与瑞士褐牛(0.42)^[10]的遗传力相近，高于上海^[21]和武汉地区^[22]的荷斯坦牛(0.11)，可能与出生年度和环境因素等有关；乳房深度(0.242)、悬韧带(0.183)、前乳房附着(0.230)、前乳头位置(0.231)、前乳头长度(0.173)、后乳房高(0.272)和后乳头位置(0.169)为低遗传力性状，与Campos等^[23]、Setati等^[24]、Haile-mariam和Pryce^[25]、Boettcher等^[26]的研究结果基本相一致，但高于Olasege等^[1]的研究结果(0.05~0.16)，可能是由于宁夏地处西北沙漠性干旱地区，永久性的环境因素对奶牛的影响较大，另外不同的体型线性评分人员之间也会存在一定偏差。有趣的是，本研究中乳房质地的遗传力仅为0.080，低于其他人的研究结果，可能在估计这些性状时使用了不同的评分系统、种群大小、分类系统和模型^[27]，其次，由于不同牧场牛舍卫生等级不同，牧场管理水平对乳房质地的影响也较大。遗传力估计是群体特异性的^[28]，宁夏地区荷斯坦牛乳房性状的中等偏低遗传力进一步揭示了非加性遗传和环境效应在这些性状总变异中的重要性。后期课题组会收集近些年更多的生产数据，对宁夏地区荷斯坦奶牛系谱进行详细的记录和追踪，完善乳房性状的遗传评估，为奶牛的选种选配和育种规划提供重要且有意义的资料。

3.2 乳房性状间的遗传相关和表型相关分析

由于生物体各性状间存在着一定联系，当对某一个性状进行选择时，也会对其他性状间接地产生选择效果。宁夏地区荷斯坦牛乳房性状间的遗传相关(－0.277 ~ 0.706)整体高于表型相关(－0.084 ~ 0.316)，大部分性状的表型相关与遗传相关趋势相似^[29]。其中，表现出极显著中度遗传相关的性状为乳房深度和前乳房附着、乳房质地和悬韧带、前乳房附着和悬韧带、后乳头位置和悬韧带，其余性状间多为弱相关；表型上多数性状表现出极显著弱相关，本研究结果与Boettcher等^[26]和Degroot等^[30]的研究结果相一致，表明当选择是基于这些性状中的某一个时，将对其他的乳房性状有负面的连锁效应。然而，在实践中多个乳房性状在评分时呈现同一趋势，因此各性状间强的正相关也有利于保持奶牛在乳房形态上的一致性。

3.3 乳房性状的育种值分析

动物个体育种值是可以稳定遗传的，因此根据它进行种畜选择就可以获得最大的选择进展^[31]。本研究对于9个乳房性状整体而言，乳房性状育种值排名前50的奶牛个体数最多的牧场为640190(163头)，其次为640004(75头)和640251(69头)，说明以上3个牧场乳房性状整体评分处于中上水平，从侧面也可以反映出牧场比较好的管理水平、齐备的硬件设施及整个牛群的生产性能；在今后的选育过程中，可以将这些乳房性状育种值整体靠前的牛只作为核心育种群，培育优秀奶牛，结合胚胎移植，扩大核心群规模；也可以在有单个性状优势的牧场有目的地选育，对某一个单一性状进行提升。

3.4 乳房性状的主成分分析

主成分分析是一种多元统计方法，它采用降维的思想，在丢失少量信息的情况下，通过线性变换将反映事物内部特征的多个指标转化为多个综合指标^[32]。通常，变换后的几个综合指标称为主成分。每个主成分都是原始变量的线性组合，彼此线性无关。在解决较为复杂的问题时，我们可以通过考虑几个主要成分来反映数据的内部特征，从而简化问题，提高分析效率，减少信息损失^[33]。因此，主成分分析在畜禽重要经济性状的选择中广泛应用。Mello等^[34]对巴西红信德乳牛的5个生产性状(初产犊龄、产犊间隔期、繁殖效率、总产奶量、泌乳期)进行主成分分析，筛选出的3个主成分占总方差的90.79%；商涛和程郁昕^[35]报告的奶牛乳用特征性状4个主成分占总方差的87.39%，程郁昕等^[36]报告的占总方差为93.46%。本研究特征值≥1的4个主成分占总方差的64.19%，低于上述研究中的总变异性，然而，高于Chu和Shi^[37]和陆呈委等^[38]的报道，尽管有约35%的信息丢失了，但这4个主成分综合地反映了各个性状的信息。在9个乳房性状中，PRIN1占总变异的比例最大(23.39%)，反映的综合信息含量也最多，主要反映奶牛的悬韧带和前乳房附着；乳房质地、前乳头位置和后乳房宽主要在PRIN2中体现，PRIN3在剩余部分乳房性状中负荷均较高；De Haas等^[39]和Kern等^[40]研究发现，遗传上乳房较深、较高、较宽、质地柔软对奶牛的产奶量有一定的提高，因此实际生产中在提高产奶量方面可以倾向于第二和第三主成分。

4 结论

本研究发现，宁夏地区荷斯坦牛乳房性状处于中等偏低遗传力范围(0.080~0.378)；其中后乳房宽和后乳房高的遗传力在0.25以上，同样是排名前50头奶牛育种值较高的性状(1.186以上)，其余性状育种值在0.382~1.129之间。乳房深度和前乳房附着、乳房质地和悬韧带、前乳房附着和悬韧带、后乳头位置和悬韧带分别呈极显著中度遗传相关，其余性状之间多为弱相关(遗传或表型相关)。主成分分析共得到4个新的成分，占总方差的64.19%，尽管有约35%的信息丢失了，但这4个主成分综合地反映了各个性状的信息，从而避免了分析大量相关性状，提高了分析大型复杂数据的精度，减轻了计算负担。