我國荷斯坦青年公?；蚪M選擇效果分析

蘇丁然，彭 朋，閆青霞，陳紹祜，張勝利，李 姣，劉丑生，孫東曉*

(1. 中國農(nóng)業(yè)大學(xué)動物科技學(xué)院，北京 100193； 2. 中國奶業(yè)協(xié)會，北京 100192； 3. 全國畜牧總站，北京 100125)

基因組選擇是通過覆蓋全基因組的高密度單核苷酸多態(tài)性(single nucleotide polymorphism, SNP)標記，并結(jié)合大量的后裔測定遺傳評估結(jié)果來計算個體的基因組估計育種值(genomic estimated breeding value, GEBV)[1-2]。2001年，Meuwissen等[3]首次提出基因組選擇(genomic selection, GS)的概念。2006年，Schaeffer[4]基于加拿大荷斯坦牛群體的研究結(jié)果表明，通過基因組選擇可以節(jié)約92%的育種成本。畜禽全基因組高密度 SNP 芯片技術(shù)逐漸發(fā)展成熟，2007年首款商業(yè)化牛50K全基因組SNP芯片問世，由美國Illumina公司研發(fā)[5]，由此，基因組選擇技術(shù)在奶牛育種工作中開始得到應(yīng)用，并逐步成為重要的奶牛育種策略[6-8]，其中，美國荷斯坦協(xié)會于2009年1月首次發(fā)布荷斯坦青年公牛的基因組遺傳評定結(jié)果[9]。2012年，中國農(nóng)業(yè)大學(xué)奶牛育種團隊及北京奶牛中心、北京首農(nóng)畜牧、上海奶牛中心、全國畜牧總站和中國奶業(yè)協(xié)會等單位合作，成功建立了中國荷斯坦牛基因組選擇技術(shù)平臺，被農(nóng)業(yè)部指定為我國荷斯坦青年公牛遺傳評定的唯一平臺，并開始在全國推廣應(yīng)用[10]。在基因組選擇提出之前，后裔測定是評定奶牛種用價值最可靠的方法，但后裔測定存在世代間隔長，育種成本高等缺點。基因組選擇相比后裔測定有以下優(yōu)勢：1) 能夠捕獲與目標性狀有關(guān)的所有的遺傳變異，具有較高的選擇準確性(>70%)，高于傳統(tǒng)系譜指數(shù)估計的可靠性[11-18]，且利用生物學(xué)先驗信息可以更好地提高基因組選擇的可靠性[19-21]。2) 可以不依賴表型信息進行選擇[22]。3) 可以大幅度縮短世代間隔，加快遺傳進展，基因組選擇使公牛的世代間隔從7年下降到小于2.5年[23-24]，母牛的世代間隔從4年 左右下降到接近2.5年[25-27]；產(chǎn)奶性狀的年遺傳進展從50%～100%增加到300%，低遺傳力性狀的年遺傳進展從3倍增加到4倍[28-31]，加快了奶牛繁殖性狀、健康性狀等低遺傳力性狀的遺傳進展[32-33]。4) 大大降低了育種成本[34-36]。

我國自2012年開始，在全國所有公牛站實施荷斯坦青年公?；蚪M遺傳評估，截至2019年12月已有3 031頭青年公牛參與基因組遺傳評估，其中1 686頭 公牛目前已獲得了女兒至少一個完整泌乳期的生產(chǎn)性能測定(dairy herd improvement，DHI)和體型鑒定數(shù)據(jù)。基于此，本研究旨在利用我國荷斯坦青年公牛基因組遺傳評估結(jié)果及女兒表型值對基因組選擇的效果進行分析，驗證基因組遺傳評估的準確性，為進一步完善我國荷斯坦青年公牛基因組選擇技術(shù)體系及我國奶牛群體遺傳改良計劃工作提供理論依據(jù)。

1 材料與方法

1.1 數(shù)據(jù)來源與整理

1.1.1 數(shù)據(jù)來源 基于中國農(nóng)業(yè)大學(xué)等單位建立的中國荷斯坦?；蚪M選擇參考群體和技術(shù)平臺，2019年共計對3 031頭荷斯坦公牛進行了基因組遺傳評估，包括基因組性能指數(shù)(genomic China performance index, GCPI)及產(chǎn)奶量、乳蛋白率、乳脂率、體型總分、泌乳系統(tǒng)、肢蹄和體細胞評分等7個 性狀的GEBV；其中，1 686頭公牛目前已經(jīng)獲得后裔測定成績，而包含性能指數(shù)(China perfor-mance index，CPI)及上述7個性狀的估計育種值(estimated breeding value，EBV)。

1 686頭荷斯坦公牛的女兒DHI數(shù)據(jù)來自全國27個省市自治區(qū)2012—2019年期間2 018個牛場的416 086頭女兒的720 393條DHI以及1 216頭荷斯坦公牛的56 902頭女兒的113 710條體型鑒定記錄。DHI數(shù)據(jù)包括305 d產(chǎn)奶量、305 d乳脂率、305 d乳蛋白率與305 d體細胞評分；體型鑒定數(shù)據(jù)包括體型總分、肢蹄評分與泌乳系統(tǒng)評分。其中3個 305 d產(chǎn)奶性狀和305 d體細胞評分的表型值是基于DHI測定日記錄，根據(jù)泌乳曲線擬合的數(shù)學(xué)模型計算得到的；體型總分、肢蹄評分和泌乳系統(tǒng)評分為各部位中所包含的體型性狀按各自權(quán)重計算得到。全部表型數(shù)據(jù)由中國奶牛數(shù)據(jù)中心提供。

1.1.2 數(shù)據(jù)整理 通過對荷斯坦青年公牛女兒生產(chǎn)性能記錄原始數(shù)據(jù)的分布進行分析，結(jié)合奶牛的生理規(guī)律，并參照行業(yè)標準《中國荷斯坦牛生產(chǎn)性能測定技術(shù)規(guī)范》，對生產(chǎn)性能數(shù)據(jù)按以下標準進行篩選：1) 牛號編碼正確；2) 第1、2、3泌乳期產(chǎn)奶性狀記錄完整；3) 無重復(fù)記錄；4) 305 d產(chǎn)奶量為1 500～18 000 kg；5) 305 d乳脂率為2.0%～6.0%；6) 305 d乳蛋白率為2.0%～5.0%；7) 305 d體細胞評分為0～9；8) 系譜完整。

按照國家標準《中國荷斯坦牛體型鑒定規(guī)程》，對體型鑒定數(shù)據(jù)按以下質(zhì)控標準進行篩選：1) 母牛的牛號正確；2) 第1泌乳期體型測定記錄完整；3) 牛 只體型測定記錄完整；4) 各部位評分和體型總分介于0～100分之間；5) 無重復(fù)記錄；6) 系譜完整。

1.2 統(tǒng)計分析方法

1.2.1 荷斯坦公牛的基因組遺傳評估結(jié)果與后裔測定結(jié)果相關(guān)性分析 利用R軟件，對基因組遺傳評估結(jié)果與后裔測定結(jié)果進行相關(guān)性分析，分別計算以下指標間的Spearman相關(guān)系數(shù)：1)荷斯坦青年公牛GCPI與后裔CPI；2)305 d產(chǎn)奶量、乳蛋白率、乳脂率、體細胞評分和體型總分、泌乳系統(tǒng)評分、肢蹄評分等7個性狀各自的GEBV與EBV；3)對 荷斯坦青年公牛的GCPI、CPI分別由高到低進行排名，計算GCPI、CPI排名均在前 10%、前 20%、前 30%、前 40%、前 50%的公牛以上7個性狀GEBV 與 EBV。

1.2.2 基于女兒性能數(shù)據(jù)分析基因組選擇效果 根據(jù)305 d產(chǎn)奶量、乳蛋白率、乳脂率、體細胞評分和體型總分、泌乳系統(tǒng)、肢蹄評分等7個性狀的GEBV分別對公牛由高到低進行排名并劃分為2種分組，分組1劃分標準見表1，分組2劃分標準見表2。

表1 荷斯坦青年公牛GEBV劃分標準

表2 荷斯坦青年公牛高低組GEBV劃分標準

通過R軟件對各組公牛女兒表型值進行方差分析，并利用Bonferroni法進行多重比較。

產(chǎn)奶性狀分析模型：y=μ+F+P+G+e，體型性狀分析模型：y=μ+F+G+e，式中，y為女兒表型值，μ為群體均值，F(xiàn)為場-固定效應(yīng)，P為胎次效應(yīng)-固定效應(yīng)，G為公牛GEBV分組-固定效應(yīng)，e為隨機殘差。

2 結(jié) 果

2.1 產(chǎn)奶與體型性狀的基本統(tǒng)計量

由表3可知，產(chǎn)奶量、乳脂率和乳蛋白率的均值分別為8 382.99 kg、3.85%和3.31%，變異系數(shù)均小于0.30；體細胞評分的均值為4.02，變異系數(shù)為0.36；體型總分、肢蹄和泌乳系統(tǒng)評分的均值分別為85.68、85.46和85.91，變異系數(shù)均小于0.10。

表3 產(chǎn)奶與體型性狀的基本統(tǒng)計量

圖1為產(chǎn)奶性狀和體型性狀的表型值分布直方圖，可以看出產(chǎn)奶量、乳脂率、乳蛋白率和SCS基本上都服從正態(tài)分布，而體型總分、肢蹄評分和泌乳系統(tǒng)評分均偏離正態(tài)分布，正態(tài)分布檢驗也顯示，產(chǎn)奶性狀的4個性狀服從正態(tài)分布，而體型性狀的3個性狀偏離正態(tài)分布。

圖1 產(chǎn)奶性狀和體型性狀表型值分布直方圖

2.2 公牛的基因組遺傳評估與后裔測定成績相關(guān)性分析

荷斯坦公牛GCPI與后裔CPI的Spearman相關(guān)系數(shù)為0.47。各性狀的GEBV與EBV的相關(guān)性分析結(jié)果如表4所示，產(chǎn)奶量和體細胞評分具有較強的正相關(guān)(0.4

表4 荷斯坦青年公牛各性狀GEBV與EBV的相關(guān)性

本研究進一步比較分析了GCPI、CPI排名均在前10%、前20%、前30%、前40%、前50%的公牛各性狀GEBV和EBV的相關(guān)性，結(jié)果如表5所示，排名前10%的公牛各性狀育種值Spearman相關(guān)系數(shù)除乳蛋白率和乳脂率外都較高(rs>0.7)，排名前20%、前30%、前40%、前50%的公牛及總樣本的各性狀育種值Spearman相關(guān)系數(shù)均低于前10%。

表5 荷斯坦青年公牛各性狀GEBV與EBV的相關(guān)系數(shù)

2.3 基于女兒生產(chǎn)性能測定(DHI)數(shù)據(jù)分析公?；蚪M選擇效果

2.3.1 公牛GEBV連續(xù)分組1 公牛女兒表型值和公牛GEBV的變化趨勢如圖2所示。利用全國的生產(chǎn)性能數(shù)據(jù)分析基因組選擇效果，發(fā)現(xiàn)公牛產(chǎn)奶量、乳蛋白率GEBV與公牛女兒表型值的變化趨勢一致，隨著表型值下降而下降；乳脂率中除低組外，GEBV與表型值變化趨勢基本一致；體細胞評分高組、中高組與中等組的GEBV與表型值變化趨勢一致，中低組與低組的表型值呈上升趨勢，與GEBV變化趨勢不一致。

圖2 各公牛組GEBV及其女兒表型值變化趨勢

通過對北京和上海地區(qū)及其他省市(地區(qū))生產(chǎn)性狀的基因組選擇效果進行比較分析發(fā)現(xiàn)，在產(chǎn)奶量中，北京及上海地區(qū)公牛GEBV與公牛女兒表型值的變化幅度較其他地區(qū)的變化幅度更為一致；在乳蛋白率中，其他地區(qū)公牛女兒表型值呈先上升后下降的趨勢，與公牛GEBV的下降趨勢不一致；乳脂率中除高組外，北京及上海地區(qū)公牛GEBV與公牛女兒表型值變化趨勢較為一致均呈下降趨勢，其他地區(qū)各組公牛GEBV呈下降趨勢，而與公牛女兒表型值先下降后上升的趨勢不一致；體細胞評分中，其他地區(qū)公牛女兒表型值在低組呈上升趨勢，與公牛GEBV的下降趨勢不一致。

生產(chǎn)性狀各GEBV分組的最小二乘均值與多重比較結(jié)果如表6所示，利用全國生產(chǎn)性能數(shù)據(jù)對產(chǎn)奶性狀基因組選擇效果進行分析，結(jié)果顯示，公牛各生產(chǎn)性狀GEBV從高到低各組間女兒相應(yīng)性狀的表型值大部分存在極顯著差異(P<0.01)。其中，公牛乳蛋白率GEBV低組與中低組組間差異不顯著(P>0.05)；乳脂率GEBV排名后5%的公牛女兒的乳脂率異常偏高，且和中等組、中高組相比，差異不顯著(P>0.05)；體細胞評分中等組、中低組、低組的組間差異不顯著(P>0.05)。

表6 產(chǎn)奶性狀表型值及各GEBV組的比較(最小二乘均值±標準誤)

利用北京及上海地區(qū)生產(chǎn)性能數(shù)據(jù)對產(chǎn)奶性狀基因組選擇效果進行分析，結(jié)果顯示，產(chǎn)奶量和乳蛋白率各組之間均存在極顯著差異(P<0.01)；在乳脂率中，高組低于中高組，低組略高于中低組，差異均不顯著(P>0.05)；在體細胞評分中，除中組與中低組之間差異不顯著(P>0.05)，其他各組均差異顯著且由高組到低組呈下降趨勢。

利用其他地區(qū)生產(chǎn)性能數(shù)據(jù)對產(chǎn)奶性狀基因組選擇效果進行分析，結(jié)果顯示，產(chǎn)奶量中各組由高到低呈下降趨勢，且各組間差異極顯著(P<0.01)；在乳蛋白率中，高組乳蛋白率水平異常偏低，低于中高組和中等組，且差異極顯著(P<0.01)；在乳脂率中，高組低于中高組，但差異不顯著(P>0.05)，而低組高于中低組，且差異極顯著(P<0.01)。

2.3.2 公牛GEBV兩尾分組2 對各性狀GEBV排名高低組公牛在女兒表型值上的差異進行對比分析，結(jié)果如圖3所示。利用全國的生產(chǎn)性能數(shù)據(jù)高低組進行分析，結(jié)果顯示，產(chǎn)奶量、乳脂率和體細胞評分的前5%和后5%差值分別為999 kg、0.11%和0.34分，前25%與后25%差值分別為452 kg、0.02%和0.21分，前5%與后5%的差值均大于前25%與后25%的差值；前5%和后5%公牛的女兒在乳蛋白率上僅相差0.04%，與前25%和后25%的差值(0.04%)相同。

圖3 公牛女兒高組與低組表型值對比圖

北京及上海地區(qū)產(chǎn)奶量的前5%與后5%、前25%與后25%的差值分別為566 kg和419 kg，其他省市(地區(qū))的差值分別為1 008 kg和473 kg，其他省市(地區(qū))的高低組差值要高于北京及上海地區(qū)；在乳蛋白率中，北京及上海地區(qū)的高組與低組差值分別為0.10%和0.04%，均高于其他省市(地區(qū))的0.01%和-0.01%；北京及上海地區(qū)的乳脂率與體細胞評分的前5%與后5%差值分別為0.05%和0.48，均高于其他省市(地區(qū))的0.02%和0.28，而前25%與后25%的差值要低于其他省市(地區(qū))。

2.4 基于女兒體型鑒定數(shù)據(jù)分析公牛因組選擇效果

2.4.1 公牛GEBV連續(xù)分組1 公牛GEBV及其女兒表型值變化趨勢如圖4所示。利用全國體型鑒定數(shù)據(jù)分析基因組選擇效果，結(jié)果發(fā)現(xiàn)，公牛女兒的體型總分、泌乳系統(tǒng)表型值與公牛GEBV變化趨勢相反，肢蹄評分表型值與公牛GEBV趨勢相同，呈下降趨勢。

圖4 公牛GEBV及其女兒表型值變化趨勢

進一步比較分析北京和上海地區(qū)及其他省市(地區(qū))的體型性狀基因組選擇效果發(fā)現(xiàn)，在北京和上海地區(qū)，3種體型性狀公牛GEBV與公牛女兒表型值變化趨勢基本一致，均呈下降趨勢；在其他地區(qū)中，公牛女兒體型總分與泌乳系統(tǒng)評分呈上升趨勢，與公牛GEBV的下降趨勢相反，公牛女兒的肢蹄評分為先上升后下降，也與公牛GEBV的下降趨勢不一致。

體型性狀各GEBV分組的最小二乘均值與多重比較結(jié)果如表7所示。利用全國體型鑒定數(shù)據(jù)分析基因組選擇效果，結(jié)果表明，在體型總分和泌乳系統(tǒng)表型值中，高組比低組分別低1.04分和1.09分，且差異均極顯著(P<0.01)；在肢蹄評分表型值中，高組、中高組和中等組公牛間女兒肢蹄評分差異不顯著(P>0.05)，而高組比低組高0.81分，差異極顯著(P<0.01)。

表7 體型性狀表型值及各GEBV組的比較(最小二乘均值±標準誤)

利用北京及上海地區(qū)的體型鑒定數(shù)據(jù)分析基因組選擇效果，結(jié)果表明，在體型總分中，除中組外其他各組均呈下降趨勢，高組與低組相差0.70分，差異極顯著(P<0.01)；泌乳系統(tǒng)各組均呈下降趨勢，高組與低組相差1.80分，差異極顯著(P<0.01)；肢蹄評分中，各組呈下降趨勢，高組與低組相差1.20分，差異極顯著(P<0.01)。

利用其他省市(地區(qū))的體型鑒定數(shù)據(jù)分析基因組選擇效果，結(jié)果表明，體型總分與泌乳系統(tǒng)中各組呈上升趨勢，高組比低組分別低0.55分與1.10分，且差異極顯著(P<0.01)；肢蹄評分各組也呈上升趨勢，但高組與低組之間差異不顯著(P>0.05)。

2.4.2 公牛GEBV雙尾分組2 對各體型性狀GEBV排名高低組的公牛在女兒表型值上的差異進行對比分析，結(jié)果如圖5所示。在體型總分與泌乳系統(tǒng)中，前5%均低于后5%且分別相差-1.04和-1.09，而前25%與后25%分別相差0.30和0.46；在肢蹄評分中，高低組分別相差0.81和0.50，前5%與后5%的差值要大于前25%與后25%的差值。

圖5 公牛女兒高組與低組表型值對比圖

北京及上海地區(qū)的體型總分、肢蹄評分和泌乳系統(tǒng)評分前5%與后5%的差值分別為0.70、1.20和0.80，前25%與后25%的差值分別為0.30、0.82和0.46。其他省市(地區(qū))的前5%與后5%的差值分別為-0.55、0.13和-1.10，前25%與后25%的差值分別為-2.42、0.02和-4.44，體型總分與泌乳系統(tǒng)評分的高組均低于低組。

3 討 論

3.1 產(chǎn)奶與體型性狀表型分析

本研究中，產(chǎn)奶性狀的表型數(shù)據(jù)均服從正態(tài)分布，而體型性狀的3個性狀有些偏離正態(tài)分布。體型性狀是由體型鑒定員按照9分制(1～9的整數(shù))對頭胎產(chǎn)犢后30～180 d內(nèi)荷斯坦母牛20個性狀進行線性評分，再按照各個性狀所占權(quán)重計算得到部位評分，然后各部位評分加權(quán)得到體型總分。因此，體型性狀表型值不同于數(shù)量性狀的連續(xù)分布，使其偏離正態(tài)分布。

3.2 基因組遺傳評估結(jié)果與后裔測定成績的相關(guān)性

整體上看，我國荷斯坦青年公?；蚪M遺傳評估結(jié)果與后裔測定結(jié)果一致性較好。公牛GCPI與CPI的相關(guān)系數(shù)為0.47，呈正相關(guān)，與國內(nèi)其他研究一致[37-39]。本研究中，產(chǎn)奶量GEBV與EBV的相關(guān)系數(shù)為0.69，要低于Wiggans等[12]報道的產(chǎn)奶量育種值的可靠性(0.70～0.80)，而體細胞評分的相關(guān)系數(shù)為0.46，也低于新西蘭家畜遺傳改良公司(Livestock Improvement Corporation Ltd.，LIC)公布的體細胞評分GEBV估計可靠性(0.50～0.67)[39]，主要原因可能為我國參考群體規(guī)模較小，且奶牛的表型數(shù)據(jù)準確性也低于發(fā)達國家[40-41]；而乳蛋白率、乳脂率的相關(guān)系數(shù)分別為0.14、0.18，要低于產(chǎn)奶量和體細胞評分的相關(guān)性，原因可能是乳脂率與乳蛋白率的測定易受DHI奶樣采集和樣品測定環(huán)節(jié)中操作不規(guī)范的影響。體型性狀(體型總分、泌乳系統(tǒng)和肢蹄評分)GEBV與EBV的相關(guān)性均較低，分別為0.19、0.24和0.26，體型性狀評分是由體型鑒定員通過觀察來對每個部位進行打分，在鑒定過程中易受鑒定員的主觀因素和鑒定環(huán)境的影響，而導(dǎo)致體型鑒定數(shù)據(jù)的準確性降低；相較于產(chǎn)奶性狀，體型性狀的數(shù)據(jù)量較少，也可能是導(dǎo)致體型性狀準確率較產(chǎn)奶性狀低的原因。有研究表明，通過擴大奶牛參考群體規(guī)模，可以提高公牛GEBV的準確性[42]，而我國的奶牛參考群體相較于奶業(yè)發(fā)達國家數(shù)量較少，這可能是導(dǎo)致我國基因組選擇準確性低的原因。

3.3 基于女兒表型值分析荷斯坦青年公?；蚪M評估效果

產(chǎn)奶性狀及體細胞評分的GEBV與表型值的變化趨勢基本一致，這與德國育種數(shù)據(jù)處理中心(Vereinigte Informations systeme Tierhaltung，VIT)公布的結(jié)果相一致，平均基因組育種值隨著表型值升高而升高，且生產(chǎn)性狀的GEBV可靠性達75%；加拿大奶業(yè)信息網(wǎng)(Canada Dairy Network，CDN)公布的基因組育種值考慮系譜信息與直接基因組評估值(direct genomic value，DGV)，平均可靠性達60%?；谌珖鴶?shù)據(jù)的乳脂率GEBV排名后5%的公牛，其女兒乳脂率反而偏高，除北京和上海以外的其他省市(地區(qū))排名后25%的公牛其女兒的乳脂率水平也異常偏高，其原因可能是有些省市(地區(qū))的部分牛場在DHI奶樣采集和測定環(huán)節(jié)的規(guī)范性存在不足，導(dǎo)致乳脂率、乳蛋白率數(shù)據(jù)的可靠性偏低，不能很好地反映出基因組選擇的效果。另外體型總分、泌乳系統(tǒng)評分女兒表型值變化趨勢同公牛GEBV變化趨勢相反。荷斯坦奶牛的體型鑒定是將20個線性性狀，按生物學(xué)特性的變異范圍，并由各省市(地區(qū))的鑒定員按照1～9的整數(shù)來進行線性評定。因此，體型數(shù)據(jù)易受鑒定員因素影響，使得不同省市(地區(qū))的體型鑒定數(shù)據(jù)存在一定的誤差，由此，全國的體型鑒定數(shù)據(jù)不能很好地反映公?；蚪M育種值的高低。

北京及上海地區(qū)的乳蛋白率、乳脂率和體細胞評分的高組與低組差值均大于其他省市(地區(qū))，體型性狀表型值高組均顯著高于低組，而其他省市(地區(qū))的高組卻低于低組。北京及上海地區(qū)相對于其他省市(地區(qū))更能反映公牛的基因組選擇效果。主要原因可能是北京及上海地區(qū)的奶牛育種歷史悠久，育種水平居全國之首，牛場DHI和體型鑒定較規(guī)范，產(chǎn)奶、體細胞評分和體型性狀的表型數(shù)據(jù)準確性較高；另外，我國的奶牛基因組選擇參考群牛只主要來自北京及上海地區(qū)的規(guī)范牛場，這兩個地區(qū)的荷斯坦青年公牛和參考群的遺傳聯(lián)系更高，因此，北京和上海地區(qū)的數(shù)據(jù)可以更好地反映基因組選擇效果。

4 結(jié) 論

基于1 686頭荷斯坦公牛基因組選擇及其女兒表型數(shù)據(jù)的分析結(jié)果表明，公牛GCPI與CPI呈正相關(guān)(rs>0.3)，產(chǎn)奶量和體細胞評分的基因組育種值(GEBV)與估計育種值(EBV)呈較強的正相關(guān)(0.4