王 澳 竹 磊 張海亮 師 睿 劉登科 馬亞賓倪俊卿 郭 剛 蔣桂娥* 王雅春*

(1.中國(guó)農(nóng)業(yè)大學(xué) 動(dòng)物科技學(xué)院/農(nóng)業(yè)農(nóng)村部動(dòng)物遺傳育種與繁殖(家畜)重點(diǎn)實(shí)驗(yàn)室/畜禽育種國(guó)家工程實(shí)驗(yàn)室，北京 100193； 2.河北省畜牧良種工作站，石家莊 050061； 3.北京首農(nóng)畜牧發(fā)展有限公司，北京 100176)

奶牛泌乳曲線是指用于描述奶牛的產(chǎn)奶量隨泌乳天數(shù)變化趨勢(shì)的數(shù)學(xué)模型，對(duì)牧場(chǎng)生產(chǎn)及育種實(shí)踐具有重要的指導(dǎo)意義，例如預(yù)測(cè)總體產(chǎn)奶量及牧場(chǎng)收益、動(dòng)態(tài)調(diào)整奶牛飼養(yǎng)方案、制定牧場(chǎng)生產(chǎn)決策、定制奶牛日常管理系統(tǒng)等。目前，國(guó)內(nèi)外已有大量研究分析了荷斯坦牛日產(chǎn)奶量的變化規(guī)律。研究表明，日產(chǎn)奶量隨泌乳天數(shù)呈規(guī)律性、連續(xù)性變化，具體表現(xiàn)為奶牛自產(chǎn)犢后開始泌乳，日產(chǎn)奶量逐漸上升至高峰后緩慢下降，直至干奶。通?？梢允褂妹谌榍€對(duì)泌乳期內(nèi)日產(chǎn)奶量的變化進(jìn)行描述，其基本特征包括日產(chǎn)奶量上升階段的上升速率，高峰奶，高峰奶出現(xiàn)的泌乳天數(shù)(高峰日)以及日產(chǎn)奶量的下降速率等。Wood模型(Wood’s incomplete gama function model)、AS模型(Ali-Schaeffer model)、Nelder模型(Nelder quadratic polynomial model)、混合線性模型等大量模型均已被成功應(yīng)用于荷斯坦牛泌乳曲線擬合。在各泌乳曲線擬合模型中，應(yīng)用最廣泛且最具代表性是不完全伽馬函數(shù)模型(Wood模型)，該模型由Wood于1967年提出。研究發(fā)現(xiàn)，Wood模型對(duì)荷斯坦牛產(chǎn)奶量的擬合度較高，可以有效解釋荷斯坦牛產(chǎn)奶量93%～98%的變異，故本研究使用Wood模型進(jìn)行泌乳曲線擬合。

受奶牛場(chǎng)生產(chǎn)管理智能化水平的限制，獲取奶牛連續(xù)監(jiān)測(cè)的日產(chǎn)奶量記錄較為困難，國(guó)內(nèi)外關(guān)于奶牛泌乳曲線的研究多利用測(cè)定日DHI(Dairy herd improvement)記錄進(jìn)行，通過(guò)DHI記錄可以建立擬合準(zhǔn)確度較高的模型對(duì)荷斯坦牛日產(chǎn)奶量的變化趨勢(shì)進(jìn)行描述與分析，并驗(yàn)證和比較已有模型在不同群體中的擬合能力。然而DHI記錄在奶牛的每個(gè)泌乳月中只有一個(gè)記錄，在刻畫奶牛的泌乳性能和研究泌乳曲線上存在固有缺陷。隨著我國(guó)規(guī)?；翀?chǎng)智能化擠奶設(shè)備的普及和管理水平的持續(xù)提升，部分牧場(chǎng)已可以實(shí)現(xiàn)全群奶牛每日產(chǎn)奶量的連續(xù)監(jiān)測(cè)，如何有效利用連續(xù)監(jiān)測(cè)日產(chǎn)奶量記錄進(jìn)行泌乳曲線擬合，并用于牧場(chǎng)精細(xì)化管理和奶牛遺傳育種，已成為新的研究熱點(diǎn)與挑戰(zhàn)。

娟姍牛具有良好的乳用性能和對(duì)熱應(yīng)激環(huán)境的高適應(yīng)性，對(duì)于豐富我國(guó)奶業(yè)的品種組成具有重要意義。自20世紀(jì)初娟姍牛被引入我國(guó)以來(lái)，養(yǎng)殖規(guī)模不斷擴(kuò)大，娟荷雜交牛已作為優(yōu)良的新型耐熱群體在我國(guó)南方高溫地區(qū)推廣養(yǎng)殖。國(guó)際上，一些研究已利用測(cè)定日DHI記錄對(duì)娟姍牛的泌乳曲線進(jìn)行了分析。例如，Wasike等研究初步證明了Wood模型對(duì)娟姍牛泌乳曲線的擬合能力；Hansen等、Mohanty等和Khan等的研究均表明了荷斯坦牛和娟姍牛的泌乳曲線在產(chǎn)奶量變化速率及高峰奶上表現(xiàn)出較大差異。然而這一方面的研究在我國(guó)奶牛群體中仍處于空白階段，作為娟姍牛品種特征及與荷斯坦牛種間差異研究的重要一環(huán)，對(duì)娟姍牛泌乳曲線和日產(chǎn)奶量變化規(guī)律進(jìn)行分析并將其與荷斯坦牛的泌乳規(guī)律進(jìn)行比較具有重要意義。

因此，本研究收集了河北地區(qū)某規(guī)?；翀?chǎng)娟姍牛和荷斯坦牛連續(xù)監(jiān)測(cè)日產(chǎn)奶量記錄，旨在擬合娟姍牛的泌乳曲線，全面分析其305 d產(chǎn)奶量和泌乳曲線相關(guān)指標(biāo)的影響因素，并比較荷斯坦牛和娟姍牛兩個(gè)群體泌乳曲線之間的差異，以期為兩品種奶牛種間差異研究，及有效利用每日連續(xù)監(jiān)測(cè)的日產(chǎn)奶量記錄進(jìn)行生產(chǎn)性能選育提供理論依據(jù)與參考。

1 材料與方法

1.1 數(shù)據(jù)來(lái)源

1

.

1

.

1

日產(chǎn)奶量

本研究使用的數(shù)據(jù)來(lái)源于河北省定州市某規(guī)?；膛?chǎng)，該牧場(chǎng)全群存欄約18 000頭，同時(shí)飼養(yǎng)荷斯坦牛和娟姍牛，采用散欄牛舍進(jìn)行飼養(yǎng)，實(shí)行TMR(Total mixed rations)飼喂、轉(zhuǎn)盤式擠奶臺(tái)集中擠奶，成母牛全群平均單產(chǎn)達(dá)到10 000 kg。該牧場(chǎng)的轉(zhuǎn)盤式擠奶機(jī)可以測(cè)定并連續(xù)記錄全群奶牛的每日產(chǎn)奶量。本研究收集了牧場(chǎng)2017年1月—2021年1月7 814頭泌乳牛日產(chǎn)奶量記錄共1 903 723條，覆蓋了10 524個(gè)泌乳期，包括荷斯坦牛1 056 155條和娟姍牛847 568條，產(chǎn)奶量記錄由奶廳智能監(jiān)測(cè)系統(tǒng)自動(dòng)導(dǎo)出。

1

.

1

.

2

繁殖記錄

試驗(yàn)?zāi)翀?chǎng)的繁殖記錄由牧場(chǎng)管理軟件導(dǎo)出的事件信息整理而來(lái)，提取試驗(yàn)牛的胎次、產(chǎn)犢日期、配種日期、配種次數(shù)、出生日期和初產(chǎn)月齡等用于分析。

1.2 數(shù)據(jù)質(zhì)控

對(duì)原始數(shù)據(jù)進(jìn)行質(zhì)控時(shí)，要求每頭奶牛的泌乳天數(shù)處于1～305 d、日產(chǎn)奶量處于2.5～100 kg、初產(chǎn)日齡處于600～1 800 d；對(duì)每個(gè)泌乳期內(nèi)日產(chǎn)奶量記錄的條數(shù)進(jìn)行質(zhì)控，具體質(zhì)控過(guò)程及每一步驟所對(duì)應(yīng)的個(gè)體數(shù)、泌乳期數(shù)和數(shù)據(jù)量如表1所示。其中，步驟4剔除異常值時(shí)，分別對(duì)兩個(gè)群體每個(gè)泌乳天數(shù)(Days in milk, DIM)內(nèi)的日產(chǎn)奶量記錄進(jìn)行一次質(zhì)控，剔除每個(gè)泌乳天數(shù)內(nèi)產(chǎn)奶量低于平均值減3倍標(biāo)準(zhǔn)差的個(gè)體記錄；步驟5～8保證了日產(chǎn)奶量數(shù)據(jù)在泌乳早、中、后期均勻分布且擁有充足的數(shù)據(jù)量用于奶牛個(gè)體泌乳曲線擬合。質(zhì)控完成后，形成數(shù)據(jù)集1和數(shù)據(jù)集2；數(shù)據(jù)集1包括921頭娟姍牛的435 836條(覆蓋1 517個(gè)泌乳期)日產(chǎn)奶量記錄；數(shù)據(jù)集2包括1 067頭荷斯坦牛的311 413條(覆蓋1 067個(gè)泌乳期)日產(chǎn)奶量記錄。數(shù)據(jù)集1和2中，每個(gè)個(gè)體每胎次內(nèi)產(chǎn)奶量記錄條數(shù)分布如圖1所示，兩群體奶牛在1～305 d內(nèi)大多有300條以上產(chǎn)奶量記錄，然而，僅有1 %的泌乳期包含完整的305條產(chǎn)奶量記錄，大量泌乳期內(nèi)的記錄不完整，需要對(duì)缺失值進(jìn)行填充。

表1 數(shù)據(jù)質(zhì)控過(guò)程
Table 1 Process control of the data quality

步驟Step質(zhì)控項(xiàng)目Item質(zhì)控條件Condition娟姍牛 Jersey荷斯坦牛 Holstein個(gè)體數(shù)/頭Numberof cows泌乳期/個(gè)Number oflactations數(shù)據(jù)量/條Data size個(gè)體數(shù)/頭Numberof cows泌乳期/個(gè)Number oflactations數(shù)據(jù)量/條Data size0原始數(shù)據(jù)Raw data1 7213 850847 5686 0936 6741 056 1551泌乳天數(shù)Days in milk1～305 d1 6623 680735 7405 8036 150876 3722日產(chǎn)奶量Daily milk yield2.5～100 kg1 6573 672721 6695 7816 122863 5493初產(chǎn)日齡Age at first calving600～1800 d1 6573 672721 6695 7816 122863 5494異常值Outliers按DIM對(duì)數(shù)據(jù)進(jìn)行分組,利用3倍標(biāo)準(zhǔn)差去除下側(cè)異常值1 6573 672721 6695 7816 122863 5495305 d內(nèi)記錄條數(shù)Number of records within 305 d >220條1 1352 100574 1661 6341 634448 5286305 d內(nèi)記錄的時(shí)間跨度Time span of records within 305 d >260 d9541 629465 3841 1251 125326 6627首次記錄DIMDIM with the first record within alactation <30 d9501 616461 9321 1031 103320 8578末次記錄DIMDIM with the lastrecord within alactation >270 d9211 517435 8361 0671 067311 413

(a)數(shù)據(jù)集1，娟姍牛；(b)數(shù)據(jù)集2，荷斯坦牛(a) Dataset 1 (Jersey cattle); (b) Dataset 2 (Holstein cattle)圖1 數(shù)據(jù)集1和2中每個(gè)泌乳期內(nèi)日產(chǎn)奶量記錄條數(shù)分布Fig.1 Distribution of record number within a lactation in dataset 1 and 2

1.3 數(shù)據(jù)分析

1

.

3

.

1

泌乳曲線擬合本研究選取Wood模型進(jìn)行奶牛個(gè)體泌乳曲線擬合，數(shù)據(jù)集1和2的后續(xù)處理及擬合過(guò)程如圖2所示。首先，步驟a采用簡(jiǎn)單移動(dòng)平均法建立4-d窗口對(duì)每個(gè)泌乳期的日產(chǎn)奶量記錄進(jìn)行預(yù)處理，即對(duì)于奶牛第

n

天的日產(chǎn)奶量，利用第

n

-3天、第

n

-2天、第

n

-1天和第

n

天數(shù)據(jù)的平均值進(jìn)行賦值，第1～3天產(chǎn)奶量記錄保持原值，如式(1)所示。使4-d窗口遍歷整個(gè)泌乳期，完成移動(dòng)窗口平均平滑。

圖2 奶牛個(gè)體泌乳曲線擬合過(guò)程Fig.2 The fitting process of individual lactation curve

(1)

基于步驟a平滑后的數(shù)據(jù)，步驟b采用Wood模型擬合每個(gè)個(gè)體每胎的泌乳曲線，并計(jì)算決定系數(shù)(

R

)和均方根誤差(Root mean squared error, RMSE)作為擬合準(zhǔn)確度的判斷標(biāo)準(zhǔn)。Wood模型如式(2)所示，其求導(dǎo)函數(shù)形式表示了泌乳期內(nèi)日產(chǎn)奶量的變化速率，如式(3)所示。其中，

t

表示泌乳天數(shù)；

Y

表示該泌乳天數(shù)的日產(chǎn)奶量；

a

表示奶牛的泌乳潛力，

b

表示泌乳曲線下降階段的速率，

c

表示曲線達(dá)到高峰的速率；

e

表示自然常數(shù)。

Y

=

at

e

-

(2)

Y

′=

abt

-1

e

--

act

e

-

(3)

為了降低少數(shù)偏離數(shù)據(jù)點(diǎn)對(duì)擬合曲線的影響，在步驟c中計(jì)算了每一條泌乳曲線實(shí)際記錄點(diǎn)與曲線擬合值間的離差，并利用四分位法去除極端低于泌乳曲線的數(shù)據(jù)點(diǎn)(負(fù)偏離點(diǎn))，再次進(jìn)行Wood模型擬合，反復(fù)執(zhí)行曲線擬合和偏離點(diǎn)判定直至不再產(chǎn)生新的負(fù)偏離點(diǎn)，得到個(gè)體最終Wood曲線。

為了保證后續(xù)分析的準(zhǔn)確性，在步驟d中對(duì)曲線的擬合程度進(jìn)行了二次質(zhì)控，質(zhì)控條件為日產(chǎn)奶量擬合值處于0～100 kg、

R

>0.7、RMSE<3。步驟e中利用泌乳曲線擬合值對(duì)相應(yīng)泌乳期內(nèi)日產(chǎn)奶量記錄中的缺失值進(jìn)行填充，最終得到1 607頭奶牛共490 135條日產(chǎn)奶量數(shù)據(jù)，共計(jì)1 607個(gè)泌乳期，用于計(jì)算305 d產(chǎn)奶量。其中，娟姍牛有1 001個(gè)泌乳期，荷斯坦牛有606個(gè)泌乳期。

1

.

3

.

2

泌乳曲線相關(guān)指標(biāo)定義本研究中，共定義了305 d產(chǎn)奶量(MY)、高峰奶(PY)、高峰日(

t

)、日產(chǎn)奶量變化速率(

v

)、日產(chǎn)奶量最高提升速率(

v

)、日產(chǎn)奶量下降階段平均速率(

v

)、各泌乳階段產(chǎn)奶量占比比值(

r

)等7個(gè)泌乳曲線相關(guān)指標(biāo)用于分析娟姍牛的泌乳曲線特征；在比較荷斯坦牛和娟姍牛泌乳曲線差異時(shí)，定義了高峰奶占比(

p

)、日產(chǎn)奶量最高提升速率占比(

p

)和日產(chǎn)奶量下降階段平均速率占比(

p

)3個(gè)指標(biāo)。其中，MY指奶牛產(chǎn)犢后第1 天至第305 天的累計(jì)產(chǎn)奶量；PY 指奶牛產(chǎn)犢后第1 天至第305 天日最高產(chǎn)奶量；

t

指奶牛高峰奶時(shí)的泌乳天數(shù)；

v

指奶牛日產(chǎn)奶量隨泌乳天數(shù)的變化速率；

v

指奶牛日產(chǎn)奶量變化速率的最大值，利用式(3)進(jìn)行計(jì)算，即

t

=1～305 d時(shí)，

Y

′的最大值；

v

指奶牛高峰日后日產(chǎn)奶量變化速率的平均值，利用式(3)進(jìn)行計(jì)算，即

t

=

t

～305 d內(nèi)，

Y

′的平均值；r指奶牛在泌乳早期(1～99 d)、泌乳中期(100～199 d)和泌乳后期(200～305 d)產(chǎn)奶總量占MY的比值；

p

指PY對(duì)MY的比值；

p

指

v

對(duì)MY的比值；

p

指

v

對(duì)MY的比值。

1

.

3

.

3

統(tǒng)計(jì)分析基于SAS 9.4軟件的 GLM過(guò)程，采用一般線性模型進(jìn)行影響因素分析，多重比較采用Bonferroni

t

檢驗(yàn)，顯著水平設(shè)為0.05?；跀?shù)據(jù)集1利用模型(1)分析娟姍牛305 d產(chǎn)奶量及泌乳曲線相關(guān)指標(biāo)的影響因素，基于數(shù)據(jù)集1和數(shù)據(jù)集2的并集利用模型(2)分析娟姍牛和荷斯坦牛泌乳曲線之間的差異。此外，采用SAS 9.4軟件的CORR過(guò)程分析娟姍牛MY、PY、

t

、

v

和

v

等指標(biāo)之間的相關(guān)。

Y

=

μ

+

P

+

S

+

N

+

A

+

e

(1)

式中：

Y

表示泌乳曲線

R

、MY、PY、

t

，

μ

表示群體平均值，

P

表示胎次(

i

=1～4)，

S

表示產(chǎn)犢季節(jié)(

j

=1～4)，

N

表示配種次數(shù)(

k

=1、2)，

A

表示初產(chǎn)月齡(

m

=1～4)，

e

表示隨機(jī)殘差。

Y

=

μ

+

P

+

S

+

N

+

A

+

B

+

e

(2)

式中：

Y

表示曲線模型

R

及泌乳曲線相關(guān)指標(biāo)，

B

表示品種效應(yīng)(

n

=1、2)，其他因子的含義同模型1。

本研究中，將胎次劃分為4個(gè)水平，包括1胎、2胎、3胎和4胎及以上；產(chǎn)犢季節(jié)劃分為4個(gè)水平，包括春季(3～5月)、夏季(6～8月)、秋季(9～11月)和冬季(12～次年2月)；配種次數(shù)劃分為3個(gè)水平，包括1次、2次和3次及以上；初產(chǎn)月齡劃分為4個(gè)水平，包括<23、23～24、25～26和27月齡及以上。

2 結(jié)果與分析

2.1 娟姍牛各泌乳曲線相關(guān)指標(biāo)分布

本研究采用Wood模型對(duì)娟姍牛個(gè)體泌乳曲線進(jìn)行擬合，平均

R

為0.86±0.07，模型擬合度可以滿足后續(xù)分析。娟姍牛個(gè)體泌乳曲線相關(guān)指標(biāo)的分布如圖3所示，包括305 d產(chǎn)奶量(圖3(a))、高峰奶(圖3(b))、高峰日(圖3(c))、日產(chǎn)奶量最高提升速率(圖3(d))、日產(chǎn)奶量下降階段平均速率(圖3(e))、

R

(圖3(f))。MY、PY和

v

呈正態(tài)分布，

t

和

v

分別呈右偏態(tài)和左偏態(tài)分布，而模型

R

在0.7～1.0間分布較為均勻。

(a)305 d產(chǎn)奶量;(b)高峰奶;(c)高峰日;(d)日產(chǎn)奶量最高提升速率;(e)日產(chǎn)奶量下降階段平均速率;(f)R2(a) 305 d milk yield; (b) Peak yield; (c) Peak day; (d) Maximum rate of increase of daily milk yield; (e) Average rate of decline of daily milk yield; (f) R2圖3 娟姍牛泌乳曲線相關(guān)指標(biāo)和模型R2分布Fig.3 Distribution of lactation curve indicators and R2 in Jersey cattle

娟姍牛泌乳曲線相關(guān)指標(biāo)的描述性統(tǒng)計(jì)如表2所示。娟姍牛MY的平均值為6 387.60 kg，變異系數(shù)為26.52%，最大值為11 622.89 kg，最小值為1 984.76 kg。PY為31.25 kg，

t

為59.56 d，在泌乳前、中、后期的產(chǎn)奶總量比為1∶0.873∶0.628。娟姍牛

v

平均為2.78 kg/d，99.0 %的

v

出現(xiàn)于泌乳第1天至第2天內(nèi)，即泌乳第2天的產(chǎn)奶量較第1天的高2.78 kg。在日產(chǎn)奶量下降階段，娟姍牛

v

為-0.06 kg/d，即從泌乳高峰日(泌乳第59天)至第305天，日產(chǎn)奶量下降約14.76 kg。

表2 娟姍牛泌乳曲線相關(guān)指標(biāo)描述性統(tǒng)計(jì)量
Table 2 Describe statistics of lactation curve indicators in Jersey cattle

項(xiàng)目Item平均值Mean標(biāo)準(zhǔn)差SD最大值Max最小值MinR20.860.070.990.70MY305/kg6 387.601 694.1711 622.891 984.76PY/kg31.256.9957.2512.55tPY/d59.5628.002919泌乳早期(1～99 d)產(chǎn)奶量占比/%Proportion of milk yield from DIM 1 to DIM 9939.975.0158.9527.17泌乳中期(100～199 d)產(chǎn)奶量占比/%Proportion of milk yield from DIM 100 to DIM 19934.912.3444.5726.13泌乳后期(200～305 d)產(chǎn)奶量占比/%Proportion of milk yield from DIM 200 to DIM 30525.114.8636.418.32vmax/(kg/d)2.780.925.51-0.04vdown/(kg/d)-0.060.030.00-0.17

注：MY,305 d產(chǎn)奶量;PY,高峰奶;,高峰日;,日產(chǎn)奶量最高提升速率;,日產(chǎn)奶量下降階段平均速率。下同。
Note: MY, 305 d milk yield; PY, Peak yield; , Peak day; , Maximum rate of increase of daily milk yield; , Average rate of decline of daily milk yield. The same below.

2.2 娟姍牛各泌乳曲線指標(biāo)間相關(guān)性

本研究對(duì)娟姍牛泌乳曲線各指標(biāo)進(jìn)行相關(guān)性分析，結(jié)果如表3所示。MY與PY、

t

、

v

等均呈極顯著正相關(guān)(

P

<0.01)，而與

v

呈極顯著負(fù)相關(guān)(

P

<0.01)；其中，MY分別與PY、

v

呈強(qiáng)正相關(guān)關(guān)系，即娟姍牛的高峰奶越高，日產(chǎn)奶量初期提升越快，全泌乳期產(chǎn)奶量越高。

v

分別與PY、

v

呈極顯著中等負(fù)相關(guān)關(guān)系，即娟姍牛高峰日前產(chǎn)奶量越高，下降階段日產(chǎn)奶量下降速率越大。

表3 娟姍牛泌乳曲線各指標(biāo)間相關(guān)系數(shù)
Table 3 Correlation coefficient of lactation curve indicators in Jersey cattle

項(xiàng)目ItemMY305/kgPY/kgtPY/dvmax/(kg/d)vdown/(kg/d)MY305/kg0.875**0.154**0.602**-0.164**PY/kg<0.01-0.0430.700**-0.505**tPY/d<0.010.173-0.101**0.024vmax/(kg/d)<0.01<0.01<0.01-0.420**vdown/(kg/d)<0.01<0.010.451<0.01

注：上三角為相關(guān)系數(shù)，下三角為值；右上角無(wú)上標(biāo)代表相關(guān)關(guān)系不顯著，上標(biāo)*代表相關(guān)關(guān)系顯著(<0.05)，上標(biāo)**代表相關(guān)關(guān)系極顯著(<0.01)。
Note: The upper triangle is the correlation coefficient and the lower triangle is the -value. No superscript in the upper right corner means the correlation is not significant, * means the correlation is significant (<0.05), ** means the correlation is highly significant (<0.01).

2.3 娟姍牛各指標(biāo)影響因素

利用模型1分析胎次、產(chǎn)犢季節(jié)、初產(chǎn)月齡和配種次數(shù)對(duì)

R

和泌乳曲線相關(guān)指標(biāo)的影響，各因素不同水平的最小二乘均值及其標(biāo)準(zhǔn)誤如表4所示。

2

.

3

.

1

模型擬合準(zhǔn)確度(

R

)本研究中，胎次、產(chǎn)犢季節(jié)和初產(chǎn)月齡對(duì)

R

均有顯著影響(

P

<0.05)。隨胎次的升高，

R

呈現(xiàn)上升的趨勢(shì)，1胎娟姍牛的

R

最低，為0.836，4胎及以上娟姍牛的

R

最高，為0.870。不同產(chǎn)犢季節(jié)中，夏季產(chǎn)犢奶牛的擬合效果顯著(

P

<0.05)高于冬季產(chǎn)犢的奶牛。不同初產(chǎn)月齡中，23～24月齡首次產(chǎn)犢的娟姍牛的擬合效果顯著(

P

<0.05)高于其他初產(chǎn)月齡的牛；過(guò)早(小于23月齡)或過(guò)晚(大于27月齡)首次產(chǎn)犢的奶牛，其泌乳曲線的擬合效果均較差。

2

.

3

.

2

305 d產(chǎn)奶量、高峰奶、高峰日由表4可知，胎次、產(chǎn)犢季節(jié)和初產(chǎn)月齡對(duì)娟姍牛MY和PY均有顯著影響(

P

<0.05)，胎次、產(chǎn)犢季節(jié)和配種次數(shù)對(duì)娟姍牛

t

影響顯著(

P

<0.05)。隨著胎次的升高，MY呈現(xiàn)逐漸增加的趨勢(shì)，在第1胎時(shí)最低(5 190.64 kg)，4胎及以上時(shí)達(dá)到最高(6 977.06 kg)，PY隨胎次的變化趨勢(shì)與MY相同。在各產(chǎn)犢季節(jié)中，MY由高到低為春季(6 727.30 kg)、冬季(6 724.76 kg)、秋季(6 001.00 kg)和夏季(5 893.55 kg)。不同產(chǎn)犢月齡分組中，23～24月齡及27月齡以上初次產(chǎn)犢的娟姍牛，其MY顯著高于其他月齡(

P

<0.05)。

表4 各因素對(duì)娟姍牛泌乳曲線擬合準(zhǔn)確度及相關(guān)指標(biāo)的影響
Table 4 Effects of various factors on and lactation curve indicators in Jersey cattle

項(xiàng)目Item頭數(shù)NumberR2MY305/kgPY/kgtPY/d胎次 Parity 1胎 Lactation 12150.84 b5 190.64±109.04 c24.40±0.42 c67.26±1.95 a 2胎 Lactation 22710.87 a6 602.14±94.10 ab32.21±0.36 b54.87±1.68 b 3胎 Lactation 32950.86 a6 576.78±91.02 b32.46±0.35 b58.43±1.63 b 4胎及以上 Lactation 4+2200.87 a6 977.06±106.74 a34.25±0.41 a60.80±1.91 ab產(chǎn)犢季節(jié) Calving season 春季(3～5月) Spring2550.87 b6 727.30±99.24 a33.51±0.38 a60.65±1.77 b 夏季(6～8月) Summer3130.89 a5 893.55±87.84 b30.07±0.34 b56.11±1.57 bc 秋季(9～11月) Fall1890.85 c6 001.00±109.88 b29.39±0.42 b52.32±1.96 c 冬季(12～次年2月) Winter2440.83 d6 724.76±99.53 a30.35±0.38 b72.28±1.78 a初產(chǎn)月齡 Age at first calving <231770.86 ab6 295.79±121.79 ab30.34±0.47 ab64.50±2.18 23～243330.87 a6 482.99±85.71 a31.43±0.33 a60.55±1.53 25～262380.86 ab6 048.78±102.23 b29.74±0.39 b57.34±1.83 27月齡及以上2530.85 b6 519.07±98.76 a31.81±0.38 a58.97±1.76配種次數(shù) Number of insemination 1次4470.866 301.68±75.3630.94±0.2957.75±1.35 b 2次2140.876 287.91±104.4030.62±0.4059.92±1.86 ab 3次及以上3400.866 420.37±86.3330.94±0.3363.34±1.54 a

注：表中的標(biāo)準(zhǔn)誤均小于0.01，故未在表中列出；同列數(shù)據(jù)后所標(biāo)字母相異表示差異顯著(<0.05)，所標(biāo)字母相同或無(wú)字母表示差異不顯著(>0.05)。
Note: The standard errors of are less than 0.01 which are not listed in table 4. In the same column, values with different right letter superscripts mean significant difference (<0.05), while the same or no letter superscripts mean no significant difference (>0.05).

對(duì)于高峰日，1胎娟姍牛

t

(67.26 d)顯著(

P

<0.05)遲于其他胎次，

t

在2胎時(shí)降至最低，而后逐漸上升。在各產(chǎn)犢季節(jié)中，

t

由高到低依次為冬季(72.28 d)、春季(60.65 d)、夏季(56.11 d)和秋季(52.32 d)。初產(chǎn)月齡對(duì)

t

無(wú)顯著影響，但呈現(xiàn)出隨初產(chǎn)月齡的增加，

t

逐漸降低的趨勢(shì)。隨著配種次數(shù)的增加，

t

呈上升的趨勢(shì)，在3次及以上時(shí)達(dá)到最高(63.34 d)。

2.4 娟姍牛與荷斯坦牛泌乳曲線比較

利用模型2分析了娟姍牛和荷斯坦牛泌乳曲線間的差異，結(jié)果如表5和圖4所示。娟姍牛泌乳曲線的擬合程度顯著高于荷斯坦牛(

P

<0.05)，娟姍牛的

t

較荷斯坦牛提前了39.84 d，

p

較荷斯坦牛高0.036%。在泌乳初期，娟姍牛的初始產(chǎn)奶量(泌乳第1天時(shí)的平均產(chǎn)奶量)較荷斯坦牛更高，且

p

顯著高于荷斯坦牛；在日產(chǎn)奶量下降階段，娟姍牛的

v

約為0.063 kg/d，約占305 d產(chǎn)奶量的0.01‰，顯著(

P

<0.05)高于荷斯坦牛。娟姍牛在泌乳前、中、后期產(chǎn)奶量在305 d產(chǎn)奶量占比的比值為1∶0.873∶0.628，娟姍牛泌乳主要集中于泌乳前期，在泌乳后期總奶量顯著(

P

<0.05)降低，而荷斯坦牛為1∶1.106∶0.956，在泌乳前、中、后期產(chǎn)奶總量較為穩(wěn)定。

表5 娟姍牛和荷斯坦牛泌乳曲線相關(guān)指標(biāo)之間的差異
Table 5 Differences of lactation curve indicators between Jersey and Holstein cattle

品種BreedR2MY305/kgPY/kgpPY/%tPY/d娟姍牛Jersey0.86 a6 322.42±55.06 b31.42±0.13 b0.497 a60.03±1.27 b荷斯坦牛Holstein0.84 b9 517.07±112.83 a43.87±0.38 a0.461 b99.87±2.60 a品種Breedvmax/(kg/d)pmax/%vdown/(kg/d)pdown/‰rlac娟姍牛Jersey2.72±0.03 b0.043 a-0.063 b-0.010 a1∶0.873∶0.628荷斯坦牛Holstein3.35±0.06 a0.035 b-0.065 a-0.007 b1∶1.106∶0.956

注：表中、、、和的標(biāo)準(zhǔn)誤均小于0.01，故未在表中列出；同列數(shù)據(jù)后所標(biāo)字母相異表示差異顯著(<0.05)，所標(biāo)字母相同表示差異不顯著(>0.05)。MY,305 d產(chǎn)奶量;PY,高峰奶;,PY對(duì)MY的比值;,高峰日;,日產(chǎn)奶量最高提升速率;,對(duì)MY的比值;,日產(chǎn)奶量下降階段平均速率;,對(duì)MY的比值;,奶牛在泌乳早、中、后期產(chǎn)奶總量占MY的比值。
Note: The standard errors of , 、、 and are less than 0.01 which are not listed in table 5. In the same column, values with different right letter superscripts mean significant difference (<0.05), while the same or no letter superscripts mean no significant difference (>0.05). MY, 305 d milk yield; PY, peak yield; , the ratio of PY to MY; , peak day; , maximum rate of increase of daily milk yield; , the ratio of to MY; , average rate of decline of daily milk yield; , the ratio of to MY; , the ratio of total milk yield in early, middle and late lactation stage to MY.

(a)娟姍牛和荷斯坦牛每個(gè)泌乳天數(shù)下日產(chǎn)奶量占305 d產(chǎn)奶量的比值；(b)娟姍牛和荷斯坦牛每個(gè)泌乳天數(shù)下的平均日產(chǎn)奶量(a) Proportion of daily milk yield per lactation for Jersey and Holstein; (b) Daily milk yield per lactation for Jersey and Holstein圖4 娟姍牛和荷斯坦牛泌乳曲線之間的差異Fig.4 Differences of lactation curve in Jersey and Holstein cattle

3 討 論

3.1 連續(xù)監(jiān)測(cè)數(shù)據(jù)的處理方法

受奶牛場(chǎng)日產(chǎn)奶量監(jiān)測(cè)設(shè)備和技術(shù)條件限制，國(guó)內(nèi)外大多利用DHI記錄對(duì)奶牛泌乳曲線進(jìn)行研究，然而這種方式存在數(shù)據(jù)量小、精細(xì)化程度低的缺陷，無(wú)法對(duì)奶牛日產(chǎn)奶量變化規(guī)律進(jìn)行深入描述與刻畫。隨著智能化擠奶設(shè)備在我國(guó)牧場(chǎng)的推廣及應(yīng)用，如何對(duì)大量日產(chǎn)奶量連續(xù)監(jiān)測(cè)數(shù)據(jù)進(jìn)行研究與分析，并將其有效應(yīng)用于奶牛遺傳育種及牧場(chǎng)實(shí)踐中，必將在今后成為領(lǐng)域內(nèi)新的問(wèn)題與挑戰(zhàn)。本研究針對(duì)連續(xù)監(jiān)測(cè)記錄中數(shù)據(jù)波動(dòng)大、噪聲數(shù)據(jù)多的特點(diǎn)，創(chuàng)新了數(shù)據(jù)處理方式，取得了良好的效果。首先，本研究利用簡(jiǎn)單移動(dòng)平均法建立4 d窗口對(duì)數(shù)據(jù)集進(jìn)行預(yù)處理。簡(jiǎn)單移動(dòng)平均法能夠?qū)r(shí)間序列數(shù)據(jù)進(jìn)行逐項(xiàng)推移，依次計(jì)算出給定項(xiàng)數(shù)的均值，是一種簡(jiǎn)單高效的預(yù)測(cè)方法，能夠有效消除數(shù)據(jù)隨機(jī)波動(dòng)的影響，突顯出數(shù)據(jù)的變化趨勢(shì)。Poppe等研究證明了利用滑動(dòng)窗口消除連續(xù)日產(chǎn)奶量記錄中隨機(jī)波動(dòng)的有效性。在擬合泌乳曲線時(shí)，對(duì)每一個(gè)泌乳期，采用邊質(zhì)控邊擬合的方式，反復(fù)利用四分位法選擇并去除噪聲數(shù)據(jù)，提升了曲線擬合質(zhì)量，保證曲線能夠真實(shí)還原奶牛日產(chǎn)奶量的變化趨勢(shì)。最后，為了保證后續(xù)分析的準(zhǔn)確性，本研究對(duì)泌乳曲線進(jìn)行了二次質(zhì)控。綜上，本研究數(shù)據(jù)處理方法為今后奶牛日產(chǎn)奶量及其他連續(xù)監(jiān)測(cè)數(shù)據(jù)的處理及分析提供了參考。

3.2 娟姍牛305 d產(chǎn)奶量及其影響因素

娟姍牛產(chǎn)奶量的遺傳力為0.08～0.28，與荷斯坦牛相近(0.30左右)。因此，娟姍牛的產(chǎn)奶量同樣受遺傳和環(huán)境因素共同影響。本研究發(fā)現(xiàn)，胎次對(duì)娟姍牛305 d產(chǎn)奶量有顯著影響，頭胎產(chǎn)奶量顯著低于其他胎次，與Mohanty等和Wasike等研究結(jié)果一致。有研究指出，荷斯坦牛初產(chǎn)牛由于仍處于生長(zhǎng)階段，采食量較低、乳腺分泌組織發(fā)育不完善，305 d產(chǎn)奶量較經(jīng)產(chǎn)牛更低；娟姍牛初產(chǎn)牛的生長(zhǎng)趨勢(shì)與荷斯坦牛類似，故胎次對(duì)產(chǎn)奶量表現(xiàn)出相同影響。在娟姍牛的日常管理中，應(yīng)注意分群飼養(yǎng)，注重低胎次奶牛的生長(zhǎng)需求和日糧的營(yíng)養(yǎng)水平，以提高牧場(chǎng)的生產(chǎn)效益。

產(chǎn)犢季節(jié)對(duì)娟姍牛305 d產(chǎn)奶量影響顯著，夏季產(chǎn)犢的奶牛305 d產(chǎn)奶量最低，春季最高。研究表明，溫濕度指數(shù)的變化對(duì)奶牛日產(chǎn)奶量存在顯著影響，夏季產(chǎn)犢的奶牛在泌乳高峰期將經(jīng)歷嚴(yán)重的熱應(yīng)激，影響到泌乳曲線上升階段日產(chǎn)奶量的提升，從而大幅降低了其305 d產(chǎn)奶量。李朝明等研究表明，娟姍牛在受熱應(yīng)激影響后，產(chǎn)奶量約下降了14.08%，同時(shí)，有研究表明，奶牛不同季節(jié)的采食量和營(yíng)養(yǎng)水平也存在較大差異，這可能也是影響娟姍牛不同產(chǎn)犢季節(jié)下305 d產(chǎn)奶量的重要因素之一。

奶牛的配種次數(shù)相對(duì)反映了其妊娠能力，是一種重要的功能性狀；初產(chǎn)月齡反映了奶牛的生長(zhǎng)水平及受孕能力，是一種重要的生產(chǎn)指標(biāo)。對(duì)于荷斯坦牛而言，許多研究表明，初產(chǎn)月齡在24月齡時(shí)可獲得最佳的生產(chǎn)效益，而配種次數(shù)與產(chǎn)奶量之間的關(guān)系并不明確。然而，對(duì)于娟姍牛的初產(chǎn)月齡、配種次數(shù)等繁殖性狀與泌乳性能之間的關(guān)系，國(guó)內(nèi)外均鮮有研究。本研究結(jié)果表明，隨著初產(chǎn)月齡和配種次數(shù)的增加，娟姍牛的305 d產(chǎn)奶量呈現(xiàn)逐漸上升的趨勢(shì)。大量研究報(bào)道了奶牛泌乳性能與繁殖性能之間存在負(fù)相關(guān)關(guān)系，高產(chǎn)奶牛往往表現(xiàn)出產(chǎn)犢間隔更長(zhǎng)、初產(chǎn)月齡更高及配種次數(shù)增加的趨勢(shì)，與本研究結(jié)果相符。因此，在牧場(chǎng)生產(chǎn)管理中，要注重處理娟姍牛泌乳性能和繁殖性能之間的關(guān)系，保障牧場(chǎng)整體收益。

3.3 娟姍牛和荷斯坦牛泌乳曲線比較

本研究選取Wood模型作為個(gè)體泌乳曲線的擬合模型，以比較娟姍牛和荷斯坦牛泌乳曲線間的差異。許多研究報(bào)道了不同品種奶牛泌乳曲線形狀的差異，包括荷斯坦牛、娟姍牛、娟荷雜交牛、丹麥紅牛等，其差異主要表現(xiàn)在高峰奶、高峰日、泌乳初期產(chǎn)奶量以及泌乳持續(xù)力等指標(biāo)。泌乳持續(xù)力是指奶牛在日產(chǎn)奶量達(dá)到峰值后維持產(chǎn)量的能力，在DHI報(bào)告中通常使用測(cè)定日記錄進(jìn)行計(jì)算。在本研究中，高峰日后日產(chǎn)奶量下降階段的平均速率(

v

)能夠間接表現(xiàn)泌乳持續(xù)力；

v

越小，泌乳持續(xù)力越低。Mohanty等研究指出，奶牛的高泌乳性能與其泌乳初期高初始產(chǎn)量、高持續(xù)力顯著相關(guān)。本研究中，娟姍牛泌乳初期的產(chǎn)奶量及提升速率顯著高于荷斯坦牛，在兩群體日產(chǎn)奶量達(dá)到高峰時(shí)，娟姍牛的高峰奶占比更高，高峰日較荷斯坦牛顯著提前，且日產(chǎn)奶量的下降速率顯著快于荷斯坦牛，故娟姍牛在泌乳中后期的總產(chǎn)奶量顯著降低，從而使娟姍?？偯谌樾阅艿陀诤伤固古?，各泌乳階段產(chǎn)奶量占比間的差異也更大。綜上，在生產(chǎn)管理和奶牛育種中，應(yīng)注重奶牛高峰奶和持續(xù)力的綜合提升，以發(fā)揮奶牛泌乳潛力，提高泌乳期總產(chǎn)奶水平。

3.4 Wood模型對(duì)娟姍牛泌乳曲線的擬合能力

Wood模型強(qiáng)調(diào)了泌乳曲線的基本特征，對(duì)奶牛群體產(chǎn)奶量的擬合度較高，因此成為研究和應(yīng)用中最廣泛使用的模型；同時(shí)，Wood模型也存在一些不足。Dematawewa等和Dijkstra等指出，Wood模型所估計(jì)的泌乳曲線易高估泌乳前期和后期的產(chǎn)奶量，而對(duì)泌乳中期的產(chǎn)奶量有所低估；此外，Wood模型的3個(gè)參數(shù)之間具有高度相關(guān)性(0.70～0.90)，這種參數(shù)特征不可避免地會(huì)造成模型對(duì)原數(shù)據(jù)分布的高敏感性。本研究表明，Wood模型對(duì)娟姍牛泌乳曲線擬合效果較優(yōu)，模型

R

可達(dá)0.86，顯著高于荷斯坦牛。Wasike等和Khan等表明，相較于其他模型，Wood模型對(duì)娟姍牛日產(chǎn)奶量的適應(yīng)性更強(qiáng)，表現(xiàn)出更優(yōu)的擬合精度。不同模型在不同品種奶牛泌乳曲線擬合的效果存在差異，這種差異可能與品種效應(yīng)、模型適用范圍等有關(guān)。本研究發(fā)現(xiàn)，荷斯坦牛的305 d產(chǎn)奶量顯著高于娟姍牛；若將本研究中所有泌乳曲線相關(guān)指標(biāo)由比例形式轉(zhuǎn)換為原形式(圖4(b))，則娟姍牛的高峰奶和

v

均顯著低于荷斯坦牛，即娟姍牛在全泌乳期內(nèi)日產(chǎn)奶量的變化幅度較荷斯坦牛更小，這可能也是導(dǎo)致Wood模型對(duì)娟姍牛擬合度更高的原因之一，即Wood模型可能更適于擬合泌乳期內(nèi)日產(chǎn)奶量波動(dòng)較小的奶牛。Melzer等研究指出，與AS模型及非參數(shù)模型等參數(shù)計(jì)算較為復(fù)雜的模型相比，Wood模型對(duì)變化速率較快或波動(dòng)較大的日產(chǎn)奶量數(shù)據(jù)擬合效果相對(duì)較差，與本研究推測(cè)相符。

4 結(jié) 論

本研究利用娟姍牛和荷斯坦牛日產(chǎn)奶量的每日連續(xù)監(jiān)測(cè)記錄進(jìn)行了泌乳曲線擬合，分析了娟姍牛泌乳曲線相關(guān)指標(biāo)的影響因素，并比較了娟姍牛和荷斯坦牛泌乳曲線間的差異。研究發(fā)現(xiàn)，胎次、產(chǎn)犢季節(jié)和初產(chǎn)月齡等因素對(duì)娟姍牛305 d產(chǎn)奶量有顯著影響，娟姍牛與荷斯坦牛在高峰奶、高峰日和日產(chǎn)奶量變化速率等方面均表現(xiàn)出顯著差異，娟姍牛泌乳集中于泌乳前期(1～99 d)且高峰日較荷斯坦牛提前39.84 d。綜上，本研究為揭示我國(guó)飼養(yǎng)條件下娟姍牛的泌乳特征，比較相同管理?xiàng)l件下荷斯坦牛與娟姍牛泌乳曲線的種間差異提供了理論依據(jù)，并為利用每日連續(xù)監(jiān)測(cè)的日產(chǎn)奶量記錄進(jìn)行牧場(chǎng)精準(zhǔn)管理和產(chǎn)奶性能選育提供了參考。