棉花是一種重要的紡織原材料,也是一種戰(zhàn)略物資,從事棉花相關(guān)行業(yè)的人口有2000多萬。同其他農(nóng)產(chǎn)品如豬肉、大豆、大蒜等產(chǎn)品類似,棉花產(chǎn)業(yè)也有一定的波動性,這種波動性既影響農(nóng)民的利益,也影響消費(fèi)者的利益,更影響決策者的政策制定,因此提前預(yù)測棉花的產(chǎn)量、種植面積等信息與多方人員關(guān)系密切,只有及時(shí)、準(zhǔn)確地對棉花產(chǎn)業(yè)進(jìn)行預(yù)測,制定合理正確的決策,才能降低其波動性,保護(hù)產(chǎn)業(yè)鏈上諸多部門和人員的利益,實(shí)現(xiàn)棉花產(chǎn)業(yè)的健康、持續(xù)發(fā)展。
在眾多預(yù)測算法中,灰色系統(tǒng)理論由于對統(tǒng)計(jì)數(shù)據(jù)少、信息不完全數(shù)據(jù)序列的擬合與建模獨(dú)具特色,已在經(jīng)濟(jì)、科學(xué)等諸多領(lǐng)域的預(yù)測、決策及評估中得到了廣泛應(yīng)用。陳啟軍等建立改進(jìn)型GM(1,1)灰色預(yù)測模型,對溫州市2018—2022年農(nóng)業(yè)灌溉水利用系數(shù)進(jìn)行預(yù)測,結(jié)果顯示預(yù)測結(jié)果基本符合溫州市實(shí)際情況[1]。劉超等建立GM(1,1)警度預(yù)測模型,利用2009—2016年的數(shù)據(jù)進(jìn)行檢驗(yàn),結(jié)果顯示模型的可靠性良好[2]。然而,GM(1,1)模型中只有兩個(gè)可調(diào)參數(shù),這會對模型的精度產(chǎn)生一定的影響。提高模型精度的一種比較好的辦法是增加可調(diào)參數(shù),本文基于GM(1,1)模型設(shè)計(jì)了6參數(shù)模型和9參數(shù)模型,然而,參數(shù)過多會給求解帶來困難,甚至無法求解。近年來出現(xiàn)的仿生計(jì)算方法如遺傳算法、蟻群算法及粒子群算法[3-4]在數(shù)值優(yōu)化及參數(shù)尋優(yōu)方面有巨大優(yōu)勢,尤其是粒子群算法在數(shù)值尋優(yōu)方面算法相對簡單且易于實(shí)現(xiàn),因而本文采用粒子群算法求解模型中的未知參數(shù)。
本研究以2000—2019年全國棉花產(chǎn)量為基礎(chǔ)數(shù)據(jù)序列,基于GM(1,1)模型[GM(1,1)中只含有兩個(gè)變量a、b,為便于描述將其稱為“2變量ab擬合”模型]設(shè)計(jì)了“6變量擬合”模型及“9變量擬合”模型,并使用粒子群算法求解這些變量,找出最優(yōu)模型,再使用最優(yōu)模型預(yù)測2020年的棉花產(chǎn)量,以期為棉花產(chǎn)業(yè)的生產(chǎn)、交易、加工及決策提供依據(jù)。
灰色預(yù)測模型(Gray Forecast Model)需要的建模信息少,特別適合對小樣本進(jìn)行建模并預(yù)測未來發(fā)展趨勢,且精度較高,是處理小樣本預(yù)測問題的一種有效算法。
研究基于2000—2019年我國棉花產(chǎn)量(數(shù)據(jù)來源:國家統(tǒng)計(jì)局網(wǎng)站),各年產(chǎn)量如圖1所示。本文在基礎(chǔ)GM(1,1)模型[為與另外兩種模型相區(qū)別,稱基本GM(1,1)模型為“2變量ab擬合”]基礎(chǔ)上設(shè)計(jì)了“6變量擬合”模型及“9變量擬合”模型,需要先測試這3個(gè)模型中哪一個(gè)是最優(yōu)的,再使用此最優(yōu)的模型擬合歷史數(shù)據(jù)并預(yù)測2020年的棉花產(chǎn)量。
圖1 2000—2019年歷年全國棉花產(chǎn)量
灰色模型的本質(zhì)是利用將原始數(shù)據(jù)序列累加后一般會呈現(xiàn)指數(shù)規(guī)律的特性,先對原始數(shù)據(jù)序列求累加序列,再擬合此累加序列的方程,根據(jù)此方程可以預(yù)測未來的數(shù)據(jù)(仍然是累加序列),經(jīng)過反變換后可以得到原始數(shù)據(jù)序列的預(yù)測值。
設(shè)一個(gè)數(shù)列X(0)=[x(0)(1),x(0)(2),…,x(0)(n)],有n個(gè)元素,其中x(0)(k)≥0,k=1,2,…,n;X(1)=[x(1)(1),x(1)(2),…,x(1)(n)]稱為X(0)的一階累加序列(1-AGO),其中x(1)(k) =Σki=1x(0)(i)=x(1)(k-1)+x(0)(k),k=1,2,…,n;圖1中展示的最近20年棉花產(chǎn)量的一階累加和結(jié)果如圖2所示。
灰微分方程為式(1)
式(1)中a稱為發(fā)展系數(shù),b稱為灰作用量,參數(shù)向量a=[a,b]T可以運(yùn)用最小二乘法估計(jì)
其中Y,B分別為式(3)所示:
當(dāng)k≥n時(shí),所求得的x(1)(k+1)就是X(1)序列的預(yù)測值。
因此求得k≥n時(shí)的x(1)(k+1)序列的值,就可以利用式(5)求出X(0)序列的預(yù)測值。式(4)中只有兩個(gè)變量a和b,可由式(2)求出具體的值,為和另外兩個(gè)模型相區(qū)別,將此模型稱為“2變量ab擬合”。
圖2 2000—2019年20年棉花產(chǎn)量的一階累加和
從式(4)中可以看出其中包含一個(gè)自變量k和兩個(gè)待定變量a、b;它是式(6)的一種特殊情況,當(dāng)式(6)中 的a1、a4及a5等于0時(shí),式(6)就退化為式(4),由于式(6)包含的可能性更多,因而將式(6)稱作增廣灰色模型。為與其他模型相區(qū)別,考慮到其中包含6個(gè)待定變量,將此模型命名為“6變量擬合”模型。
同理,可以定義包含范圍更廣的模型,如式(7)所示,其中包含9個(gè)待定變量,故將此模型命名為“9變量擬合”模型。當(dāng)a1、a4、a7等于0時(shí)退化為式(6)。
如前所述,式(4)中的待定參數(shù)a和b可由式(2)計(jì)算出來,然而式(6)的待定參數(shù)a1至a6,以及式(7)中的待定參數(shù)卻難以通過少數(shù)矩陣運(yùn)算求解出來;考慮到仿生計(jì)算方法如遺傳算法、蟻群算法及粒子群算法在數(shù)值優(yōu)化及參數(shù)尋優(yōu)方面有巨大優(yōu)勢,尤其是粒子群算法在數(shù)值尋優(yōu)方面算法相對簡單且易于實(shí)現(xiàn),因而本文采用粒子群算法求解式(6)和式(7)中的參數(shù)。
粒子群算法(Particle Swarm Optimization, PSO)由Eberhart和Kennedy基于鳥群覓食行為的模擬于1995年提出。用一個(gè)粒子(實(shí)際是一個(gè)向量)來模擬一個(gè)鳥類個(gè)體,粒子的當(dāng)前值(鳥的當(dāng)前位置)即為所求問題的一個(gè)候選解,粒子不停地調(diào)整過程就是鳥的飛行過程。粒子具有速度和位置兩個(gè)屬性。粒子的飛行速度根據(jù)自己的歷史最優(yōu)位置和種群中歷史最優(yōu)粒子的位置按一定比例關(guān)系進(jìn)行修正,這種過程不斷循環(huán),不斷修正速度和位置,最終得到問題的最優(yōu)解。
算法流程如下:
(1)初始化。設(shè)置最大循環(huán)次數(shù)或終止條件,粒子個(gè)數(shù)M(鳥群中鳥的數(shù)量),粒子的最大飛行速度(速度區(qū)間),粒子中每個(gè)標(biāo)量的取值范圍(搜索區(qū)間),在速度區(qū)間中隨機(jī)初始化每個(gè)粒子的速度值,在搜索區(qū)間中隨機(jī)初始化每個(gè)粒子向量中各個(gè)標(biāo)量的值。
(2)個(gè)體最優(yōu)解與全局最優(yōu)解。定義適應(yīng)度函數(shù),個(gè)體最優(yōu)解(Pi)為每個(gè)粒子個(gè)體在飛行過程中搜索到的歷史最優(yōu)解,從所有粒子中找到的最優(yōu)解叫全局最優(yōu)解(Pg),每個(gè)粒子要朝著自己的個(gè)體最優(yōu)解和全局最優(yōu)解飛行,不斷更新。
(3)更新速度和位置。更新粒子速度的公式為式(8):
其中,ω稱為慣性因子,只能為正值或者0,ω值大時(shí)全局尋優(yōu)能力強(qiáng),反之,全局尋優(yōu)能力弱;c1、c2稱為加速常數(shù),其中c1表示個(gè)體學(xué)習(xí)因子,c2表示全局學(xué)習(xí)因子。本文中設(shè)置ω=0.8、c1=0.2、c2=0.2。random(0,1)表示生成區(qū)間[0,1]上隨機(jī)數(shù)的函數(shù);Vid表示第i個(gè)粒子速度的第d個(gè)標(biāo)量,Pid表示第i個(gè)粒子個(gè)體最優(yōu)解中的第d個(gè)標(biāo)量,Pgd表示全局最優(yōu)解中的第d個(gè)標(biāo)量。
更新粒子位置的公式為式(9):
其中,Xid表示第i個(gè)粒子的第d個(gè)標(biāo)量。
本文試驗(yàn)主要分三部分,分別對比“2變量ab擬合”模型、“6變量擬合”模型及“9變量擬合”模型的擬合精度,并從中選出一個(gè)最好的模型預(yù)測2020年的棉花產(chǎn)量?;疑P偷谋举|(zhì)是建立指數(shù)模型,在建模前先要對歷年棉花產(chǎn)量X(0)序列求一階累加和,得到X(1)序列,即對X(1)中每個(gè)元素,有x(1)(k)= Σki=1x(0)(i),k=1,2,…,n;結(jié)果如表1中的第4列所示。
3.1 2 變量ab擬合
根據(jù)X(0),X(1)構(gòu)造公式(3)中的矩陣Y、B,根據(jù)最小二乘法解得參數(shù)向量a=[a,b]T=[-0.0117,516.7069]T,從而可得模型
由式(10)可以求出X(1)的模型擬合值,具體數(shù)值如表3中第5列所示。圖3直觀對比了X(1)及模型的擬合值,由圖中可見兩條線前面8個(gè)數(shù)據(jù)差異較大,后面12個(gè)數(shù)據(jù)差異較小。
表1 近20年棉花產(chǎn)量、一階累加和及三種模型對一階累加和擬合的結(jié)果
圖3 一階累加及其“2變量ab擬合”結(jié)果對比
式(6)中含有6個(gè)待求變量,用普通的最小二乘法難以求出結(jié)果,改用粒子群算法根據(jù)X(1)求得參數(shù)a1=54.3377,a2=-46.7662,a3=0.0039,a4=-2.4487,a5=43.7599,a6=4.5330;因而公式(6)的擬合表達(dá)式為:
根據(jù)式(11)求得X(1)的擬合值,列于表3中第5列;圖4是X(1)及其擬合值的對比圖,從圖4中可以看出公式(11)的擬合結(jié)果,不論在序列的首部還是尾部都吻合得很好,效果優(yōu)于“2變量ab擬合”方法。
圖4 一階累加及其“6變量擬合”結(jié)果對比
與“6變量擬合”類似,同樣可用粒子群算法根據(jù)X(1)求得式(7)中的9個(gè)參數(shù),a1=44.3481,a2=37.4067,a3=217.1029,a4=0.0057,a5=0.0102,a6=-0.2057,a7=22.2540,a8=60.4856,a9=74.9205,因而式(7)的擬合表達(dá)式為:
根據(jù)式(12)求得X(1)的擬合值,列于表3中第6列;圖5是X(1)及其擬合值的對比圖,從圖5中可以看出式(12)的擬合結(jié)果在序列的首部和尾部都吻合得很好,效果與“6變量擬合”方法類似。
圖5 一階累加及其“9變量擬合”結(jié)果對比
表2列出了3種方法誤差對比??梢钥闯觥?變量擬合”方案對X(1)序列擬合的絕對誤差、相對誤差及平均誤差都是最小的。但是,對X(1)擬合效果好,不等于對X(0)序列擬合效果也好。表3是由對X(1)序列的3種擬合值反推出的2018、2019年棉花年產(chǎn)量,可以看出“2變量ab擬合”方法對2018、2019年產(chǎn)量估計(jì)值的相對誤差均超過了5%,“9變量擬合”方法對2018年產(chǎn)量擬合的誤差小于5%,而對2019年產(chǎn)量的擬合誤差大于5%;只有“6變量擬合”方法對這兩年的擬合誤差均小于5%;這說明擬合所用多項(xiàng)式的次數(shù)以及參數(shù)的個(gè)數(shù)并不是越多越好,因而比較合理的方案是選擇“6變量擬合”,即公式(11)方法預(yù)測2020年的產(chǎn)量。
表2 3種擬合方法的誤差對比
由式(11)可以算出X(1)序列的第21個(gè)擬合值,X(1)(21)=12587.7,2020年的產(chǎn)量估計(jì)值X(0)(21)=X(1)(21)-X(1)(20)=12587.7-12008=579.7(萬噸)。
表3 3種擬合方法對2018、2019年產(chǎn)量的擬合結(jié)果及誤差對比
在眾多預(yù)測算法中,灰色系統(tǒng)理論由于對統(tǒng)計(jì)數(shù)據(jù)少、信息不完全數(shù)據(jù)序列的擬合與建模獨(dú)具特色,已在經(jīng)濟(jì)、科學(xué)等諸多領(lǐng)域的預(yù)測、決策及評估中得到了廣泛應(yīng)用。本文以2000—2019年全國棉花產(chǎn)量為基礎(chǔ)數(shù)據(jù)序列,基于GM(1,1)模型設(shè)計(jì)了“6變量擬合”模型及“9變量擬合”模型,并使用粒子群算法求解模型中的變量,發(fā)現(xiàn)“6變量擬合”模型精度最高,是最優(yōu)模型,再使用“6變量擬合”模型預(yù)測2020年的棉花產(chǎn)量。