楊 雪 李文獻

（1武漢中育科技有限公司，湖北武漢 430070；2河南龍鳳山股份有限公司，河南駐馬店463000）

1 引言

1.1 模型設(shè)計背景

我國居民飲食習(xí)慣決定國內(nèi)生豬消費以“鮮肉”為主，生豬產(chǎn)業(yè)鏈從生豬養(yǎng)殖到活豬屠宰最后到豬肉產(chǎn)品送達餐桌，即從供應(yīng)端直接到達需求端，因而影響供需平衡的生產(chǎn)、需求、庫存3個因素里面的庫存因素對供需調(diào)節(jié)的影響很小，由此可以假設(shè)生豬供應(yīng)量約等于消費需求量。隨著我國宏觀經(jīng)濟增長和居民收入水平提高，生豬價格的需求彈性越來越小，目前我國約5 500萬噸的年豬肉消費需求量，已經(jīng)讓供應(yīng)端的收入接近邊際成本，因此需求的剛性增加已經(jīng)不足以引起強有力的價格波動進而推動新一輪豬周期的產(chǎn)生，只有當(dāng)某一因素影響到供應(yīng)端生產(chǎn)效率導(dǎo)致供應(yīng)不足，才有可能引發(fā)一輪豬周期。本輪豬周期就是由于2018年8月以來非洲豬瘟疫情嚴重影響生豬生產(chǎn)效率，從而引發(fā)的生豬價格大幅波動，同時這也是自20世紀90年代以來6輪豬周期中價格波動最強的一次。通過分析引起本輪豬周期的主因非洲豬瘟內(nèi)在邏輯對供應(yīng)影響的規(guī)律，從而達到預(yù)測本輪豬周期價格規(guī)律的目的。

1.2 模型涉及名詞解釋

主因強度因子（Strength，簡稱S）：根據(jù)主因的流行病學(xué)規(guī)律，分析主因發(fā)生后在未來持續(xù)作用的幾個月內(nèi)每個月的作用強度。

主因時間因子（Time，簡稱T）：根據(jù)生豬繁殖周期、MSY等生豬生產(chǎn)技術(shù)指標(biāo)及主因發(fā)生后的復(fù)養(yǎng)時間，分析在不同生產(chǎn)模式下，主因發(fā)生后在未來一段時間內(nèi)的生豬減產(chǎn)復(fù)產(chǎn)周期以及生豬出欄減少量的規(guī)律。

分配函數(shù)（Φ）：根據(jù)主因強度因子（S）和時間因子（T）雙維度分配結(jié)果，總結(jié)主因?qū)е碌纳i出欄減少量分配至未來各月的規(guī)律的函數(shù)。

量價運算框（Operation box）：分析生豬出欄減少量和生豬價格關(guān)系規(guī)律的運算模型，按輸入與輸出變量不同，可分為P-Y量價運算框（輸入實際生豬價格，輸出生豬出欄減少量估計理論值）與量價運算框（輸入預(yù)測生豬處理減少量，輸出預(yù)測生豬價格）。

廣義簡約梯度法運算框（GRG box）：根據(jù)P-Y量價運算框輸出的生豬出欄減少量（Y），結(jié)合分配函數(shù)（Φ）并設(shè)定相關(guān)約束條件，通過廣義簡約梯度法（GRG）輸出最佳估計值：預(yù)測生豬出欄減少量的運算模型。

2 模型結(jié)構(gòu)設(shè)計

2.1 分配函數(shù)的結(jié)構(gòu)及計算

2.1.1 主因強度因子（S）分析與計算方法

強度因子主要是基于主因在一個區(qū)域發(fā)生后，在該區(qū)域的流行病學(xué)作用強度邏輯，即主因發(fā)生后影響持續(xù)M個月，每個月影響強度為βm(m=1,2,3,…,M)，總強度因子（S）如列式（1）：

主因作用產(chǎn)生后，受資金和復(fù)養(yǎng)信心等影響，一部分受主因作用而淘汰的產(chǎn)能在短時間無法復(fù)產(chǎn)，據(jù)此本模型將24個月內(nèi)無法復(fù)產(chǎn)部分的強度因子記為S^sp，如列式（2）：

由列式（1）和列式（2）可得：

以受主因影響能繁母豬數(shù)（X）為自變量，通過強度因子（S）參數(shù)β可以計算能繁母豬減少數(shù)（Z），

計算模型如下：

其中xi表示自變量在第i個月上的值(i=0,1,2,…,N)，zt表示因變量在第t個月上的值(t=0,1,2,…,N+M-1)，主因強度影響持續(xù)M個月，故因變量可在最新自變量觀測值基礎(chǔ)上預(yù)測至未來M-1個月的值。

該計算公式的矩陣形式可表示為：

2.1.2 主因時間因子（T）分析與計算方法

⑴生豬生產(chǎn)技術(shù)指標(biāo)計算方法

生豬生產(chǎn)具有一定周期性，這也導(dǎo)致主因?qū)ιa(chǎn)產(chǎn)能的減產(chǎn)影響在時間維度上具有周期性規(guī)律，因此生豬繁育周期、MSY等生產(chǎn)技術(shù)指標(biāo)是計算主因時間因子（T）的關(guān)鍵一環(huán)。生豬繁育從母豬繁育到育肥豬出欄具有明顯的階段性，經(jīng)歷階段不同所需的繁育周期常數(shù)也不同，如列式（7） （8）：

MSY指每年每頭母豬出欄育肥數(shù)，本模型基于月度數(shù)據(jù)建立，故定義MSM為每月每頭母豬出欄育肥數(shù)，用來衡量受主因減少的能繁母豬導(dǎo)致生豬出欄減少的效率，計算如列式（9） （10） （11） （12）：

⑵不同生產(chǎn)模式（C）下減產(chǎn)復(fù)產(chǎn)周期(tij)的計算方法

復(fù)養(yǎng)時間和生產(chǎn)模式的不同使得生豬減產(chǎn)復(fù)產(chǎn)時間產(chǎn)生變化，進而影響主因時間因子（T）的計算，記：

B2的情況已包含在強度因子列式（2）的計算中，故時間因子中不再重復(fù)考慮其影響，由列式（8） （13）（14）可得：

⑶主因時間因子（T）的計算方法

時間因子主要分析不同生產(chǎn)模式下受主因影響后生豬減產(chǎn)復(fù)產(chǎn)的時間邏輯，參數(shù)表征主因發(fā)生后，每頭受主因影響減少的能繁母豬在未來一段時間內(nèi)導(dǎo)致的生豬出欄減少量在每個月的分配情況。記：

根據(jù)列式（12） （17）建立Tc的計算模型如下：

其中參數(shù)γij表示第j種生產(chǎn)模式下 (j=1,2,3,4,5)，主因發(fā)生后第i個月 (i=0,1,2,…)的時間因子。當(dāng)減產(chǎn)未開始 (iint(t2j))，時間因子不產(chǎn)生作用；當(dāng)處于減產(chǎn)復(fù)產(chǎn)周期內(nèi)時，時間因子表示每頭受主因影響減少的能繁母豬導(dǎo)致的當(dāng)月生豬出欄減少量（完整月份即為MSM，不足一月的按比例計算）。

由列式（15）和公式2-1可計算總時間因子（T）：

時間因子只考慮主因發(fā)生后24個月內(nèi)可復(fù)產(chǎn)的情況，以可復(fù)產(chǎn)的能繁母豬減少數(shù) (Znsp)為自變量，通過總時間因子（T）參數(shù)γi可以計算標(biāo)豬出欄減少量（Y），計算模型如下：

其中ztnsp表示自變量在第t個月上的值 (t=0,1,2,…,N+M-1)，γi表示時間因子在第i個月上所有生產(chǎn)模式的總值 (i=0,1,2,…,L)，時間因子的作用時長（即i的最大值）取決于不同生產(chǎn)模式中最長減產(chǎn)結(jié)束時間 (L=MAX(int(t2j))，L≤24)。y_n表示因變量在第n個月上的值(n=0,1,2,…,+N+M+L-2)，主因時間影響最長持續(xù)L個月，故因變量可在最新自變量觀測值基礎(chǔ)上預(yù)測至未來L-1個月的值。

公式2-1、2-2、2-3的矩陣形式可表示為：

2.1.3 分配函數(shù)（Φ）計算方法

主因在某一區(qū)域發(fā)生后，在該區(qū)域流行持續(xù)M個月，在這M個月內(nèi)，每個月不同強度的主因同時又作為獨立的1起疫情以該月為起點對未來i個月的生豬出欄減少量造成影響。主因以上述方式構(gòu)成了強度-時間雙維度的作用矩陣 (φim)，由列式（1） （2） （3） （12）和公式2-2可得：

不可復(fù)產(chǎn)部分最長持續(xù)時間為B2≥24，可復(fù)產(chǎn)部分最長持續(xù)時間為L≤24，故φim的i值的最大值取B2。分配函數(shù)（Φ）為同一月份的合計影響值，需將矩陣（φim）轉(zhuǎn)換為數(shù)據(jù)序列φt，表示自第1起主因發(fā)生后第t個月的分配函數(shù)值，計算如下：

以受主因影響能繁母豬數(shù)（X）為自變量，通過分配函數(shù)（Φ）參數(shù)φ_t可以計算標(biāo)豬出欄減少量（Y），計算模型如下：

其中xi表示自變量在第i個月上的值 (i=0,1,2,…,N)，φt表示分配函數(shù)在第t個月上的值 (t=0,1,2,…,B2+M-1)。yn表示因變量在第n個月上的值 (n=0,1,2,…,N+B2+M-2)，分配函數(shù) (Φ)最長分配時間為B2+M-1個月，故因變量可在最新自變量觀測值基礎(chǔ)上預(yù)測至未來B2+M-2個月的值。

公式3-1、3-2、3-3的矩陣形式可表示為：

2.2 量價運算框的結(jié)構(gòu)及計算

在2.1章節(jié)分配函數(shù)的設(shè)計中，已建立主因發(fā)生后標(biāo)豬出欄減少量（Y）在之后一段時間內(nèi)每個月的分配規(guī)律模型。本章節(jié)的量價運算框則構(gòu)建了標(biāo)豬出欄減少量（Y）與生豬價格（P）之間的運算模型，按輸入與輸出變量不同，可分為P-Y量價運算框（輸入實際生豬價格，輸出生豬出欄減少量估計理論值）與量價運算框（輸入預(yù)測生豬處理減少量，輸出預(yù)測生豬價格）。

生豬價格需求函數(shù)Q=f(P)以及豬周期起點標(biāo)準(zhǔn)豬肉年供應(yīng)量的變化函數(shù)Qt=f(i)是實現(xiàn)量價運算的重要函數(shù)：

2.3 廣義簡約梯度法（GRG）運算框結(jié)構(gòu)及計算

根據(jù)分配函數(shù)Φ（公式3-3）構(gòu)建規(guī)劃求解的Y-X函數(shù)方程，采用廣義簡約梯度法（GRG）求自變量X變化情況下因變量Y的最優(yōu)解，即得出預(yù)測標(biāo)豬出欄減少量，使得P≥Pthr（Pthr為設(shè)定有效性篩選閾限值）條件下，GRG法計算出的最優(yōu)解與2.2中由量價運算框計算出的理論值Y的相關(guān)系數(shù)最大，運算邏輯如圖2：

2.4 價格預(yù)測技術(shù)路線及結(jié)果評估方法

量價運算框與GRG運算框是本模型價格預(yù)測技術(shù)路線的核心成分，以實際價格序列數(shù)據(jù)作為輸入值，經(jīng)過3步運算，可得到未來預(yù)測價格的序列數(shù)據(jù)。將輸入的序列數(shù)據(jù)分為3個不相交的集合：

訓(xùn)練集：輸入模型進行運算，獲得訓(xùn)練誤差；

驗證集：與訓(xùn)練集相對獨立，獲得模型在該集上的預(yù)測誤差；

測試集：與訓(xùn)練集和驗證集相對獨立，獲得模型在該集上的測試誤差。

本文采用3項擬合指標(biāo)對模型預(yù)測結(jié)果進行評估，具體評估方案如表1：

2.5 預(yù)測結(jié)果迭代計算

本模型預(yù)測結(jié)果并不是靜態(tài)值，而是進行迭代計算，即模型運算結(jié)果隨輸入的實際價格數(shù)據(jù)序列的更新而更新。這種迭代設(shè)計使本模型可以盡可能接受最新的價格波動信息，對歷史預(yù)測進行糾偏。記：

3 模型測試及預(yù)測結(jié)果

3.1 樣本數(shù)據(jù)與抽樣說明

本輪豬周期受2018年8月國內(nèi)暴發(fā)的非洲豬瘟疫情影響，故為驗證本模型預(yù)測效果，本文選取2018年8月1日至2020年7月31日（不含節(jié)假日）的全國生豬均價日度數(shù)據(jù)作為原始數(shù)據(jù)總體（F），數(shù)據(jù)取自卓創(chuàng)資訊，見圖4。從總體F中抽取每月的最高價形成月度數(shù)據(jù)樣本（P，即當(dāng)月生豬最高日均價）作為研究樣本，記：

價格越高，則受主因影響導(dǎo)致的減產(chǎn)越嚴重，主因?qū)r格波動的解釋力就越強，而均價則會在一定程度上削弱主因的解釋力，故本文研究樣本的選擇使用最高價而不是均價。

3.2 參數(shù)α與β的確定及分配函數(shù)計算結(jié)果

對多組α與β取值進行評估篩選，最終確定值為

表1 預(yù)測評估指標(biāo)及方法

生豬繁育周期常數(shù)及生產(chǎn)技術(shù)指標(biāo)如下：

表2 生產(chǎn)技術(shù)指標(biāo)

將上述參數(shù)代入公式3-3得出分配函數(shù)，分配函數(shù)曲線近似正偏態(tài)分布曲線，上升段短而快，下降段長而緩，如圖5。

3.3 預(yù)測價格計算

依據(jù)生豬行業(yè)的歷史數(shù)據(jù)及國家統(tǒng)計公布的人口數(shù)據(jù)，可得量價運算框涉及常數(shù)q＝75%，w＝110 kg，Pop＝14億，同時對參數(shù)θ多組取值進行評估篩選，最終確定最優(yōu)值θ＝2.5。公式4-1、4-2和4-3引用本團隊另一研究結(jié)果：

將生豬價格樣本數(shù)據(jù)（P）與上述參數(shù)值輸入P-Y量價運算框，可計算出標(biāo)豬出欄減少量（Y）。但是本模型是建立在主因?qū)е庐a(chǎn)能減少進而影響價格波動的理論基礎(chǔ)上的，忽略了影響價格波動的其他因素，如季節(jié)性需求變化和突發(fā)事件等。因此生豬價格（P）通過P-Y量價運算框計算出的每月標(biāo)豬出欄減少量（Y），與本模型涉及的受主因影響的每月標(biāo)豬出欄減少量之間存在偏差。當(dāng)主因的解釋力越強時，這種偏差就越小，本模型的預(yù)測性就越強，因此需要對Y進行有效性篩選，進而得出受主因影響的標(biāo)豬減產(chǎn)量的最佳估計。

有效性篩選以生豬價格閾限為標(biāo)準(zhǔn)，假設(shè)P≥Pthr時，主因?qū)r格波動的解釋力最高，其他因素的影響幾乎可以忽略不計。分析2006年6月至2018年7月3輪豬周期的價格波動，生豬價格的最高點都在20元/kg左右，而本輪豬周期生豬價格超過歷史高位20元/kg，則假設(shè)Pthr=20。

將Y與有效性篩選閾值Pthr=20輸入GRG運算框，參考國內(nèi)能繁母豬存欄情況限定自變量X的取值范圍為10≤X≤800（單位：萬頭），計算出Y的最佳估計。

為比較多次迭代結(jié)果變化，將不同時間更新的樣本數(shù)據(jù)(Pi)輸入模型，得出2019年12月（i＝16）至2020年7月（i＝23）共8次迭代計算的結(jié)果，如圖7。

3.4 預(yù)測結(jié)果檢驗與評估

多次迭代計算結(jié)果在MAE、MAPE、MSE和RMSE指標(biāo)上均達到檢驗標(biāo)準(zhǔn)，在訓(xùn)練集上計算的r和R2也均大于0.95，驗證集由于受樣本量限制（n＜8）沒有表現(xiàn)出線性相關(guān)性，計算的r和R2參考性不大。見表3。

4 討論

本模型試驗測試過程使用的數(shù)據(jù)來源于卓創(chuàng)資訊，目前國內(nèi)生豬價格統(tǒng)計數(shù)據(jù)來源眾多，不同數(shù)據(jù)來源可以得到不同的預(yù)測結(jié)果，但只要原始數(shù)據(jù)統(tǒng)計標(biāo)準(zhǔn)一致就不影響模型趨勢的預(yù)測。由于當(dāng)前數(shù)據(jù)樣本有限，模型測試時使用的訓(xùn)練集為全部樣本數(shù)據(jù)，暫無其他數(shù)據(jù)集作為驗證集與測試集進一步檢驗?zāi)Ｐ托Ч?，因此由?xùn)練集計算出的評估指標(biāo)值存在“過度擬合”的可能（即訓(xùn)練誤差低估實際預(yù)測誤差），無法完全表征本模型對未來價格的預(yù)測效果好壞。但是后期的模型迭代計算設(shè)計將彌補這一問題，使預(yù)測值隨迭代更新。

表3 預(yù)測結(jié)果檢驗與評估

本模型使用的參數(shù)α、β和θ取值的篩選過程，限于本文篇幅，不在此詳細闡述，后續(xù)再對參數(shù)估計過程進行解釋。

本模型量價運算框中應(yīng)用需求函數(shù)及其反函數(shù)（公式4-1、4-2）的價格范圍是12≤P≤40，當(dāng)生豬價格超過這一數(shù)據(jù)范圍時，會產(chǎn)生模型鈍化。本輪豬周期實際預(yù)測過程中，超出定義范圍的價格出現(xiàn)次數(shù)很少，因此不影響模型整體的預(yù)測效果。

本模型量價運算框在推導(dǎo)標(biāo)豬數(shù)量對應(yīng)豬肉產(chǎn)量時，標(biāo)豬的體重（w）按110 kg計算，模型雖然最后通過2次量價運算框的計算（P-Y量價運算框與量價運算框）將這一因素抵消，從而不影響預(yù)測價格的趨勢，但w會使預(yù)測價格的絕對值偏大或偏小。我們將繼續(xù)研究價格與銷售體重的最佳銷售經(jīng)濟模型，以進一步消除預(yù)測誤差。

本模型GRG運算框采用的有效性篩選閾值Pthr=20，該閾值的選擇是存在決策誤差的，但后期我們將根據(jù)數(shù)據(jù)量的增加，不斷地校正這一因素，以減少決策誤差。

本模型除主因影響外，只考慮宏觀剛性需求的變化（公式4-3），未考慮季節(jié)性因素和突發(fā)性因素（如2020年新冠肺炎疫情）對豬肉消費的影響，因此預(yù)測價格對這2個因素產(chǎn)生的波動不敏感，當(dāng)相關(guān)事件發(fā)生時預(yù)測值的絕對誤差會偏大，但當(dāng)這一事件結(jié)束時，模型將自行校正，因此不影響長期趨勢的預(yù)測。如圖7所示的多次迭代計算結(jié)果，受新冠肺炎疫情影響，2019年底至2020年初期間，雖然主因?qū)е碌臏p產(chǎn)因素仍然存在，但這一時期需求同時被抑制，導(dǎo)致生豬價格未能表現(xiàn)出主因的真實影響，主因影響被低估，于是預(yù)測價格出現(xiàn)整體偏低。隨著新冠肺疫情得到控制，需求恢復(fù)，模型自動對前期預(yù)測進行校正，使得2020年7月迭代計算結(jié)果高于前期預(yù)測，因此本模型應(yīng)用的預(yù)測價格通常采用最新迭代計算結(jié)果。我們后期將繼續(xù)研究季節(jié)因素對需求的影響，以減少該因素產(chǎn)生的預(yù)測誤差。

本模型只適用于分析本輪豬周期起因（非洲豬瘟）引起的減產(chǎn)因素對價格波動的影響，當(dāng)非洲豬瘟因素結(jié)束或者生豬價格低于20元/kg時，本模型就不具備預(yù)測性，則需應(yīng)用其他模型進行價格預(yù)測。