某豬場母豬繁殖性能預測和logistic回歸分析

翟靜靜，劉則彰，蔣思文

（華中農業大學，湖北 武漢 430070）

母豬在養殖生產中扮演著重要的角色，在規模化養豬生產中，母豬的繁殖效率是重中之重，是養殖場的核心，也是決定豬場生產水平和經濟效益的關鍵因素[1]。母豬繁殖性能主要指標包括總產仔數、產活仔數、斷奶仔豬數、斷奶發情間隔（Weaning to estrus interval,WEI）、PSY等[2]，還包括繁殖力（如分娩率和死亡率）。

胎次也是影響母豬繁殖性能的重要因素，一般來說，繁殖性能隨著胎次數的增加而提高，在3～5胎次時達到最高水平[3,4]。在2～4胎時母豬的受胎率會隨胎次升高而增加，4胎后逐漸下降[5]。胎次主要通過影響分娩率和產仔數來影響母豬繁殖性能[6]。朱世平等發現，依據胎次可將母豬繁殖性能劃分為1-2胎、3-6胎和7-8胎這三個階段，并且得出母豬繁殖的最佳胎次為3-6胎，繁殖高峰期是4胎次（P＜0.05）[7]。

WEI是指母豬斷奶到再次發情配種的間隔時間。據國外的各個豬場數據分析，WEI與產仔數呈負相關，WEI時間越長，產仔數越低，而在3～5 d發情的母豬，分娩率高。因此要想提高斷奶母豬7 d內的配種率，應縮短WEI，保持母豬的理想體況，保障產仔舍和配懷舍的飼養管理工作，減少體重損失，做好誘情工作，提高分娩成績[7]。

生產中把繁殖母豬發情期進行配種后沒有懷孕并重新發情的現象稱為返情。常規返情是受胎失敗，非常規返情是妊娠失敗，都會導致低分娩率的出現。母豬一般配種受胎為80%～90%，而配種未分娩的10%～20%就是返情。返情的增加，會導致配種分娩率降低，從而影響養殖場的經濟效益[8,9]。返情后再配種母豬分娩率降低10%左右[10]。實際生產中母豬返情率一般在10%左右[11]。

時間序列法是一種統計分析方法，將研究對象的數量隨著時間的發生而呈現變動的這種特征按照時間順序進行排列，所形成的這個序列便稱為研究對象的時間序列[12]。常用的時間序列法有移動平均法、加權移動平均法和指數平均法，時間序列的平穩性是建模必不可少的環節[13]。時間序列分為兩類：第一類是某時間序列各時期的數值是不變動的，為確定型時間序列。第二類是此時間序列各時期數值變動，充滿不確定性，為隨機型時間序列。一般情況下都是隨機型的時間序列。朱迪國等對全球豬流行性腹瀉疫情進行時間序列分析發現，豬血液制品貿易可能是該病跨國界傳播的潛在途徑[14]。因此，在母豬繁殖生產數據中，可以運用時間序列對性狀的時間發展進行預測。

對于因變量為只具有兩種分類結果的屬性數據，例如母豬生產中的母豬分娩或未分娩，就需要去建立邏輯回歸模型，這是一種經典的分類算法[15]。Logistic模型中的參數是通過極大似然法進行估計的。邏輯回歸分析是通過一組預測變量（自變量）預測一個或多個響應變量（因變量）一種統計分析方法[16]。在使用Logistic模型前，需判斷是否滿足以下幾個條件：(1)因變量是二分類變量；(2)分類變量的分類是全面的且相互獨立；(3) 連續的自變量與因變量的對數轉換值之間存在線性關系；(4)自變量之間不存在共線性。利用建立的回歸模型可以對母豬性能及后期發展進行預測。Logistic模型也可以用來評估預測變量對響應變量的預期效果。Hosmer and Lemeshow運用Logistic模型預測母豬體重發育相關變量與第2胎中的非妊娠和產仔數的關聯，發現是高度相關[17]。Koketsu等人利用Logistic回歸模型分析與分娩率相關的因素時發現，胎次、分娩季節以及泌乳天數與采食量之間的交互作用會對母豬分娩率產生影響[18]。

研究根據豬場的生產記錄，時間序列和Logistic回歸模型這兩種方法，把所得的豬場數據進行匯總處理，運用時間序列去進行短期預測母豬繁殖性能高低的時間段和母豬狀態，來判斷豬場中繁殖母豬在某一時間段可能會出現的情況；通過Logistic回歸模型進一步剖析影響因素，對影響該豬場母豬繁殖性能的因素進行調查研究，為該豬場采取相應措施提高母豬繁殖性能提供科學依據，有效提高豬場經濟效益。

1 材料與方法

1.1 數據來源

此研究選擇某公司的2016-2019年7個分場的母豬繁殖性能記錄數據，包括產仔性能、配種分娩明細及淘汰死亡等。以母豬繁殖性能為主線，對規模豬場不同生產指標的各種數據進行處理分析，根據對生產數據的分析，找到影響繁殖性能的因素，制定行之有效的解決方案，以此展開相關研究。

1.2 數據統計分析

對數據采用Excel軟件進行統計，數據以平均值±標準差表示；運用GraphPad Prism 7.0作圖，使用SPSS 19.0建立邏輯回歸模型和對數據進行顯著性檢驗，P＜0.05表示差異顯著，P＜0.01表示差異極顯著；運用SAS 9.4軟件構建時間序列模型ARIMA進行預測。

1.3 Logistic模型的構建

此研究采用二元Logistic回歸分析方法對母豬胎次和年份進行了分析。母豬繁殖性能相關因素分析研究的變量與賦值見表1。根據收集到的數據，將數據分為兩類情況，如受胎：1和未受胎：0；首先采用多重共線性模型篩選變量（|r|＞0.7）。然后，對每個解釋變量進行單因素分析，將上述變量P＜0.1的解釋變量納入logistic回歸分析模型中，運用基于最大似然估計的向前逐步回歸法（Forward: LR法）建 立logistic回歸模型，確定與母豬繁殖性能相關風險因素的依據是似然比統計量的P值＜0.05。回歸系數用95%置信區間（Conf idence intervals，CI）的似然比（Odds ratios，OR）表示。OR值即是相對危險度的精確估計值，如果OR大于1，變量施加的風險增加。Logistic回歸模型的公式如下：

其 中，β0為 常 量，A、B分別表示因素斷奶或配種年份、母豬斷奶或配種胎次。把胎次、年份等影響7 d斷配率、返情率、受胎率、分娩率這些指標的因素進行定義變量和賦值（見表1）。

表1 母豬繁殖性能相關因素分析研究的變量與賦值

2 結果與分析

2.1 時間序列分析預測

2.1.1 時序圖檢驗

此研究收集了2016-2019.11、2016-2019.08、2016-2019.07、2016-2019.07的豬場健仔數、斷奶頭數、7 d斷配率、分娩率的數據。選取了其中80%的數據用于時間序列分析，留下20%的數據作為實際值與預測值對比。首先對原始序列進行平穩性檢驗，分別描繪平均每頭母豬健仔數（圖1a）、斷奶頭數（圖1b）、7 d斷配率（圖1c）、分娩率（圖1d）的時序圖，從時序圖中可以看出該時間序列具有明顯的遞增趨勢性和下降趨勢性，為典型的非平穩性序列。

2.1.2 一階差分時序圖

時間序列預處理，可知這四個指標均為不平穩序列，ARIMA模型基于平穩時間序列的基礎上，因此根據數據特征對數據進行差分處理，一階差分后的時序圖（圖2）顯示，各觀測值圍繞著一個常數附近上下波動，數據逐漸趨于平穩，序列表現平穩的性質。

2.1.3 模型識別

為了驗證一階差分后此序列的平穩性，依次做出該序列的自相關圖和偏相關圖（圖3到圖6）。觀察ACF（自相關）圖和拖尾效應確定ARIMA模型的P值，觀察PACF（偏相關）和截尾效應確定ARIMA模型的q值。其中差分的次數就是模型ARIMA（p， d，q）中的階數d，因此通過最佳準則函數定價方法選擇ARIMA（0，1，1）、ARIMA（5，1，5）、ARIMA（0，1，0）、ARIMA（0，1，0）分別為健仔數、斷奶頭數、7 d斷配率、分娩率的相對最優模型。

2.1.4 模型的預測

根據已經建立的ARIMA模型對母豬健仔數、斷奶頭數、7 d斷配率及分娩率的數據進行預測。預測時間為9個月，健仔數的預測時間從2019年3月到11月，斷奶頭數預測時間為2018年12月到2019年8月，7 d斷配率和分娩率預測時間都是從2018年11月到2019年7月。為評估模型的預測效果，引入均方誤差（Mean square error，MSE）和平均絕對百分比誤差（Mean absolute percentage error，MAPE）兩個指標。預測的偏差較小，與實際觀察值基本吻合（表2到表5）。預測值的變化趨勢與實際值大致保持一致，實際值基本落在預測曲線95%的可信區間內，該模型擬合預測成功，具有一定的實際價值。由表6誤差評估表可見，MSE值約小，說明預測模型具有更好的精確度，MAPE＜20%，可以從中判斷該模型擬合預測很成功。

2.2 因素分析

2.2.1 受胎率的logistic回歸分析

原始數據為192 874條配種記錄，剔除掉胎次為0的數據，剩下143 464條記錄用于后續的因素分析。由表7可知，配種年份和配種胎次會影響母豬受胎率。高胎次（≥7）（OR值：0.889，95%CI：0.810～0.975）的母豬跟低胎次（1～2）的母豬相比受胎率小于0.889倍；在2019年（OR值：1.247，95%CI：1.157～1.344）配種的母豬比2016年配種的母豬成功受胎的概率高1.247倍。

2.2.2 返情率的logistic回歸分析

收集了7個分區豬場192 874條配種記錄，剔除掉胎次為0的數據，剩下143 464條記錄用于后續的因素分析。由表8可知，配種年份和配種胎次會影響母豬返情率。3～6胎 次（OR值：1.272，95%CI：1.208～1.339）的母豬比1～2胎次的母豬出現返情現象的概率高1.272倍；在2019年（OR值：0.802，95%CI：0.744～0.864）配種的母豬比2016年配種的母豬發生返情的概率小。

2.2.3 分娩率的logistic回歸分析

表2 ARIMA模型預測母豬健仔數效果表 頭

表3 ARIMA模型預測母豬斷奶頭數效果表 頭

表4 ARIMA模型預測母豬7 d斷配率效果表 %

原始數據為192 874條配種記錄，剔除掉胎次為0的數據，剩下143 464條記錄用于后續的因素分析。由表9可知，配種年份和配種胎次會影響母豬分娩率。2017年（OR值：0.782，95%CI：0.747～0.817）的母豬比2016年的母豬分娩率低；高胎次（≥7）（OR值：0.851，95%CI：0.808～0.897）的母豬成功分娩的概率比低胎次（1～2）的母豬小85%。

2.2.4 7 d斷配率的logistic回歸分析

收集了7個分區豬場2016-2019年間142 478條斷奶配種記錄，剔除掉WEI＞100的數據，剩下138 114條記錄用于后續的因素分析。由表 10可知，斷奶年份和斷奶胎次會影響母豬7 d斷配率。2019年（OR值：1.610，95%CI：1.526～1.698）的母豬斷奶發情間隔大于7 d的發生概率是2016年的1.610倍；3～6胎次（OR值：0.939，95%CI：0.904～0.975）的母豬斷奶發情間隔大于7 d的發生概率比低胎次(1～2)的母豬低。

3 小結與討論

數據的時間序列大多是不平穩的，數據本身具有季節趨勢，若采用簡單的模型進行分析，往往不能得到準確的預測結果。隨著對時間研究的發展，時間序列研究的方法更加多樣，其中ARIMA模型是常用的非平穩的模型[13]。時間序列預測法對于中短期預測的效果要優于長期預測的效果，因為時間過長發生外界因素變化的可能性會加大。世界上沒有一個模型可以做到100%的精確，只能隨著時間發展進行改善，不斷添加新的數值和因素。此研究中以健仔數、斷奶頭數、7 d斷配率和分娩率的實際數據為基礎，擬合ARIMA模型，同時對未來9個月的數據進行了預測和比較。研究結果表明健仔數在未來幾個月呈上升趨勢，與實際值相比最大偏差為0.8頭；斷奶頭數預測值具有一定的波動性，與實際值最大偏差為1.67；7 d斷配率預測值持續上升，實際值時而下降時而上升，最大偏差為11.55個百分點；分娩率預測結果呈下降趨勢，實際值具有波動，最大偏差為4.26個百分點。這幾個指標雖有一定偏差，但實際值基本落在預測曲線95%的置信區間內，說明模型對原始序列信息的提取比較充分。通過預測效果表和MAPE值（值越小，代表預測的準確度越強）來看，利用ARIMA模型對這些數據進行擬合和預測是合理可行的，但也會有一定的局限性，生產過程中可能會因環境變化或疾病、飼養管理方案改變引起差異。在短期預測當中，豬場未來一段時間的母豬生產情況將會高于前段時間的均值，或是有走低的趨勢，這值得生產者注意調整適當的養殖方法。時間序列分析方法可以用在養豬生產中，可以預測未來一段時間的生產變化發展。

表5 ARIMA模型預測母豬分娩率效果表 %

表6 ARIMA模型誤差評估

表7 與母豬受胎率相關因素的logistic回歸結果

表8 與母豬返情率相關因素的logistic回歸結果

對受胎率、返情率、分娩率和7 d斷配率的二元logistic回歸分析，表明受胎率、返情率、分娩率和7 d斷配率均與年份和胎次有回歸關系。此研究發現在2017年受胎率高于2016年，返情率、分娩率和7 d斷配率低于2016年，說明2017年該場母豬配種要比2016年更好，而分娩率低可能是受胎后其他因素導致的。2018和2019年受胎率和7 d斷配率高于2016年，返情率和分娩率低于2016年，分娩率進一步降低可能和非洲豬瘟的暴發有關。而隨著胎次的增加受胎率、分娩率和7 d斷配率下降，返情率上升。此結果與朱世平等得出母豬繁殖的最佳胎次為3～6胎相反[7]，說明該公司繁殖母豬的繁殖性能仍有待提高。該公司豬場年份之間具有差異，特別是2017年性能比較好，后來隨著胎次增加繁殖性能不增反減，可能是有新建場區，母豬養殖水平和設施等還未跟上，也可能是豬場管理制度不嚴，需要從飼養管理方案到技術設備、母豬淘汰制度、發情鑒定和足夠的后備母豬儲備、專業的人員操作進行完善[19]，細化采精、配種和分娩等記錄，提高豬場效率。

表9 與母豬分娩率相關因素的logistic回歸結果

表10 與母豬7 d斷配率相關因素的logistic回歸結果

時間序列分析方法對于母豬繁殖性能方面有一定的短期預測效果，可進一步探索用于其他生產指標。相較于之前的研究，此次研究根據數據更全面分析了母豬的繁殖性能。并應用logistics模型和時間序列ARIMA模型對母豬繁殖性能進行因素分析及短期預測，試驗數據量大，為生產提供了更具參考價值的依據。