本輪豬周期主因非洲豬瘟量化模型設(shè)計(jì)及價(jià)格預(yù)測效果評估

楊 雪 李文獻(xiàn)

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

1 引言

1.1 模型設(shè)計(jì)背景

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

1.2 模型涉及名詞解釋

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

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

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

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

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

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

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

2.1.1 主因強(qiáng)度因子（S）分析與計(jì)算方法

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

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

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

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

計(jì)算模型如下：

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

該計(jì)算公式的矩陣形式可表示為：

2.1.2 主因時(shí)間因子（T）分析與計(jì)算方法

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

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

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

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

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

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

⑶主因時(shí)間因子（T）的計(jì)算方法

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

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

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

由列式（15）和公式2-1可計(jì)算總時(shí)間因子（T）：

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

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

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

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

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

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

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

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

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

2.2 量價(jià)運(yùn)算框的結(jié)構(gòu)及計(jì)算

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

表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ù)測價(jià)格計(jì)算

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

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

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

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

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

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

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

4 討論

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

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

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

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

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

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

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

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