趙 征， 南宏鋼， 喬錦濤

(華北電力大學(xué) 控制與計(jì)算機(jī)工程學(xué)院，河北 保定 071003)

0 引 言

截止2019年年底，我國(guó)風(fēng)電累計(jì)并網(wǎng)裝機(jī)容量達(dá)到210.05 GW，占全部發(fā)電裝機(jī)容量的10.4%[1]，風(fēng)電已經(jīng)成為繼火電、水電之后的第三大能源。風(fēng)速作為影響風(fēng)功率的第一大因素，然而風(fēng)速具有的隨機(jī)性、非平穩(wěn)性給電力系統(tǒng)的安全穩(wěn)定運(yùn)行帶來極大的沖擊與挑戰(zhàn)，因此準(zhǔn)確預(yù)測(cè)風(fēng)速對(duì)電網(wǎng)的正常運(yùn)行至關(guān)重要[2]。

根據(jù)風(fēng)速預(yù)測(cè)模型的不同，風(fēng)速預(yù)測(cè)方法主要有物理學(xué)方法、統(tǒng)計(jì)學(xué)方法及組合方法[3]。

物理學(xué)方法需要豐富的氣象數(shù)據(jù)與大量的機(jī)組技術(shù)參數(shù)[4]，建模過程復(fù)雜。而統(tǒng)計(jì)學(xué)方法[5]利用風(fēng)速歷史數(shù)據(jù)，建立系統(tǒng)輸入與輸出之間的映射關(guān)系,從而對(duì)未來風(fēng)速趨勢(shì)做出預(yù)測(cè)，建模過程簡(jiǎn)單可行。統(tǒng)計(jì)學(xué)方法中時(shí)間序列分析法主要有：自回歸滑動(dòng)平均法[6]、自回歸差分滑動(dòng)平均法[7]、廣義自回歸條件異方差法[8]等。組合方法雖然在一定程度上能夠提高模型的預(yù)測(cè)精度，但是加權(quán)系數(shù)的強(qiáng)壯性一直是這類方法的瓶頸[9]。

考慮到時(shí)間序列法采用風(fēng)速歷史數(shù)據(jù)，建模簡(jiǎn)單計(jì)算量少，其本身所具有的時(shí)序性和自相關(guān)性已經(jīng)為建模及預(yù)測(cè)提供了足夠的信息，適用于超短期及短期預(yù)測(cè)。ARIMA模型用于刻畫時(shí)間序列的自相關(guān)性；GARCH模型用于刻畫時(shí)間序列的波動(dòng)性[9]。文獻(xiàn)[10]建立了超短期風(fēng)速預(yù)測(cè)的ARIMA-GARCH模型，并且對(duì)GARCH模型中的條件方差進(jìn)行加權(quán)混合修正，以此來提高預(yù)測(cè)精度，但是考慮到風(fēng)速的波動(dòng)性較強(qiáng)，該模型并未對(duì)原始風(fēng)速數(shù)據(jù)進(jìn)行預(yù)處理。文獻(xiàn)[11]建立了短期風(fēng)速預(yù)測(cè)的CEEMDAN-QGA(quantum genetic algorithm, QGA)-BP(back propagation, BP)模型，但是并未考慮到原始風(fēng)速序列經(jīng)CEEMDAN分解后，高頻分量仍存在較大的隨機(jī)性，對(duì)其進(jìn)行VMD分解可有效地降低隨機(jī)性。二次分解[12]體現(xiàn)了對(duì)原始數(shù)據(jù)處理的“精細(xì)化”思想。

鑒于以上考慮，本文提出了一種基于CEEMDAN-VMD二次分解的ARIMA-GARCH超短期風(fēng)速預(yù)測(cè)模型。首先對(duì)原始風(fēng)速數(shù)據(jù)進(jìn)行CEEMDAN分解，采用樣本熵評(píng)估各分量的復(fù)雜性，對(duì)復(fù)雜性高的分量進(jìn)行VMD分解。二次分解后的各分量作為ARIMA-GARCH預(yù)測(cè)模型的輸入，將所有分量預(yù)測(cè)值線性疊加得到最終預(yù)測(cè)結(jié)果。

1 基本方法原理

1.1 自適應(yīng)噪聲完整集成經(jīng)驗(yàn)?zāi)B(tài)分解

自適應(yīng)噪聲完整集成經(jīng)驗(yàn)?zāi)B(tài)分解[13]在原始信號(hào)x[n]中添加滿足標(biāo)準(zhǔn)正態(tài)分布的高斯白噪聲ω[n]，則第i次的信號(hào)可以表示為xi[n]=x[n]+ωi[n](i=1,…,I)，其中，I為實(shí)驗(yàn)次數(shù)，CEEMDAN算法步驟如下：

(2)在r1[n]中添加白噪聲，即分解r1[n]+ω1[n],(i=1,…,I)，得到第2個(gè)模態(tài)分量IMF2及余量r2[n]。

(3)重復(fù)步驟(1)、(2)。當(dāng)分解出的剩余分量不適合被分解時(shí)，停止分解。最終原始信號(hào)可表示為

(1)

式中：K為分解后本征模函數(shù)的數(shù)量。

1.2 變分模態(tài)分解

變分模態(tài)分解[14]的核心是變分問題的構(gòu)造與求解。變分問題可表述為：經(jīng)VMD分解后得到若干個(gè)模態(tài)的估計(jì)帶寬最小，并且所有模態(tài)可以重構(gòu)為原信號(hào)。

VMD算法步驟如下：

(1)計(jì)算調(diào)制信號(hào)梯度的平方L2范數(shù)，對(duì)各模態(tài)信號(hào)帶寬進(jìn)行估計(jì)，即

(2)

式中：{uk}表示各模式的集合；{ωk}表示各模式的中心頻率；δ(t)表示沖激函數(shù)；f表示原始信號(hào)。

(2)加入二次懲罰因子α和拉格朗日算子λ(t)將原問題轉(zhuǎn)換為非約束變分問題。即

(3)

式中：二次懲罰項(xiàng)提高了懲罰項(xiàng)有限權(quán)重的收斂性，拉格朗日乘子用于執(zhí)行約束。

(4)

(5)

(6)

收斂準(zhǔn)則：

(7)

2 CEEMDAN-VMD-ARIMA-GARCH預(yù)測(cè)模型

2.1 自回歸差分滑動(dòng)平均模型

自回歸模型AR(p)，是利用時(shí)間序列當(dāng)前時(shí)刻值及其滯后p階的歷史信息建立多項(xiàng)式來刻畫時(shí)間序列的特征，具體表達(dá)式如下：

ut=c+φ1ut-1+φ2ut-2+…+φput-p+εtt=1,2,…,T

(8)

式中：c為常數(shù)項(xiàng)；φ1,φ2,…,φp為多項(xiàng)式系數(shù)；εt為服從εt～N(0,σ2)的白噪聲。

滑動(dòng)平均模型MA(q)，是利用時(shí)間序列當(dāng)前時(shí)刻值及其滯后q階的歷史信息建立加權(quán)多項(xiàng)式來刻畫時(shí)間序列的特征，具體表達(dá)式如下

ut=μ+θ1εt-1+θ2εt-2+…+θqεt-q+εtt=1,2,…,T

(9)

式中：μ為常數(shù)項(xiàng)；θ1,θ2,…,θq為多項(xiàng)式系數(shù)；εt為服從εt～N(0,σ2)的白噪聲。

自回歸滑動(dòng)平均模型ARMA(p,q)是AR(p)模型與MA(q)模型的組合，具體表達(dá)式如下：

(10)

ARMA模型適用于解決平穩(wěn)過程的時(shí)間序列，對(duì)于非平穩(wěn)時(shí)間序列需要進(jìn)行差分運(yùn)算，該模型變?yōu)樽曰貧w差分滑動(dòng)平均模型。

dYt=Wt

(11)

式中：Wt為非平穩(wěn)時(shí)間序列Yt經(jīng)過d階差分后得到的平穩(wěn)時(shí)間序列；d為d階差分算子。

2.2 模型定階

模型定階用于確定平穩(wěn)時(shí)間序列模型(ARMA)的階數(shù)。本文采用貝葉斯信息準(zhǔn)則[15](bayesian information criterion，BIC)確定模型階數(shù)。BIC信息準(zhǔn)則相較于AIC(akaike information criterion, AIC)信息準(zhǔn)則具有更低的階數(shù)。BIC信息準(zhǔn)則為

(12)

2.3 參數(shù)估計(jì)

參數(shù)估計(jì)用于確定ARMA模型中的未知系數(shù)φ1,φ2,…,φp和θ1,θ2,…,θq。常用的參數(shù)估計(jì)方法有矩估計(jì)法、最大似然估計(jì)法及最小二乘估計(jì)法。矩估計(jì)法運(yùn)算繁瑣且無法得到全局最優(yōu)解。最大似然估計(jì)法難以將似然函數(shù)轉(zhuǎn)化為參數(shù)的顯示函數(shù)。本文采用最小二乘估計(jì)法估計(jì)模型未知參數(shù)。

2.4 異方差性檢驗(yàn)

若時(shí)間序列模型方差不是常數(shù)，其大小隨時(shí)間變化且與過去的誤差有關(guān)，表明該模型的殘差具有異方差性。自回歸條件異方差模型(autoregressive conditional heteroskedastictiy，ARCH)通過描述方差的變化規(guī)律，對(duì)原模型的殘差序列進(jìn)行優(yōu)化。在此基礎(chǔ)上波勒斯萊文進(jìn)行改進(jìn)，建立廣義自回歸條件異方差模型。標(biāo)準(zhǔn)的GARCH(1,1)模型為

(13)

ARCH LM檢驗(yàn)[16]經(jīng)常被用來檢驗(yàn)?zāi)Ｐ蜌埐钪惺欠窬哂挟惙讲钐匦裕僭O(shè)為

H0:α1=α2=…=αq=0

(14)

LM檢驗(yàn)統(tǒng)計(jì)量TR2為樣本長(zhǎng)度乘以回歸檢驗(yàn)R2，漸進(jìn)服從χ2(p)分布，其中p為χ2分布的自由度。給定置信區(qū)間，當(dāng)TR2>χ2(p)對(duì)應(yīng)的臨界值時(shí)，則拒絕原假設(shè)，證明存在異方差效應(yīng)。

2.5 CEEMDAN-VMD-ARIMA-GARCH預(yù)測(cè)模型

針對(duì)風(fēng)速時(shí)間序列的隨機(jī)性、非平穩(wěn)性，本文提出了CEEMDAN-VMD-ARIMA-GARCH預(yù)測(cè)模型。圖1為預(yù)測(cè)模型流程圖。模型預(yù)測(cè)步驟如下：

(1)原始風(fēng)速序列經(jīng)過CEEMDAN分解，獲得含剩余分量在內(nèi)的n個(gè)分量。

(2)采用樣本熵[17]評(píng)估各分量的復(fù)雜性。針對(duì)復(fù)雜度相對(duì)較高分量進(jìn)行二次分解-VMD分解。判斷經(jīng)二次分解后各分量的平穩(wěn)性及模型的殘差特性，分別建立對(duì)應(yīng)的ARIMA-GARCH預(yù)測(cè)模型。

(3)對(duì)復(fù)雜度相對(duì)較低各分量分別判斷平穩(wěn)性及殘差特性，建立ARIMA-GARCH預(yù)測(cè)模型。

(4)將所有分量的預(yù)測(cè)值疊加，得到最終風(fēng)速預(yù)測(cè)值。

圖1 CEEMDAN-VMD-ARIMA-GARCH預(yù)測(cè)模型流程圖Fig.1 CEEMDAN-VMD-ARIMA-GARCH forecast model flow chart

3 算例實(shí)現(xiàn)

本文選用西北某風(fēng)電場(chǎng)1月1日至3日連續(xù)3天的風(fēng)速數(shù)據(jù)作為實(shí)驗(yàn)樣本，采樣時(shí)間間隔為5 min，共有864個(gè)數(shù)據(jù)點(diǎn)，預(yù)測(cè)未來3 h風(fēng)速。

3.1 CEEMDAN分解風(fēng)速序列

CEEMDAN能夠較好地解決EMD(empirical mode decomposition, EMD)、EEMD(ensemble empirical mode decomposition, EEMD)中出現(xiàn)的模態(tài)混疊現(xiàn)象。CEEMDAN分解樣本數(shù)據(jù)，原始風(fēng)速分解結(jié)果如圖2所示。

圖2 CEEMDAN分解結(jié)果圖Fig.2 CEEMDAN decomposition result chart

圖2中，原始風(fēng)速序列被分解為8個(gè)波動(dòng)較小的IMF分量和1個(gè)剩余分量RES。結(jié)合時(shí)域圖與頻譜圖，不同分量的主頻率成分不同，對(duì)應(yīng)時(shí)域的波動(dòng)性不同。前3個(gè)分量的主頻率成分不單一，說明它們?nèi)匀痪哂休^高的復(fù)雜度。為了提高模型預(yù)測(cè)精度，本文采用樣本熵評(píng)估分量的復(fù)雜性，具體結(jié)果如圖3所示。

圖3 各子分量樣本熵值Fig.3 Entropy of each subcomponent sample

從圖3各子分量樣本熵值可以看出，前3個(gè)分量的樣本熵值較高，熵值均高于1。從第4個(gè)分量開始熵值急劇減小，其中第4個(gè)分量熵值約為0.6，遠(yuǎn)小于第三個(gè)分量的熵值。綜上所述，認(rèn)為前三個(gè)分量是復(fù)雜度相對(duì)較高分量。

3.2 VMD二次分解

本文對(duì)IMF1、IMF2、IMF3三個(gè)樣本熵值較高的分量進(jìn)行VMD二次分解。

首先確定IMF1分量VMD分解的層數(shù)K。由于VMD分解層數(shù)需要人為確定，文獻(xiàn)[19]提出了一種較為便捷的模態(tài)數(shù)K的選取方法，其基本思想是考慮VMD分解之后各模態(tài)的中心頻率互不相同，如果出現(xiàn)相鄰兩個(gè)模態(tài)之間中心頻率相近，則認(rèn)為出現(xiàn)過分解。文獻(xiàn)[19]中認(rèn)為當(dāng)前分解層數(shù)減1作為最佳分解層數(shù)。本文選用該方法來選取VMD的分解層數(shù)。此外，VMD分解的初始其他參數(shù)中，α取值1 000，τ取0。IMF1分量VMD分解不同K值下的中心頻率如表1所示。

從表1可以看出，當(dāng)K=4時(shí)，相鄰兩個(gè)模態(tài)分量中心頻率的距離最小值為198×10^(-6)。與模態(tài)數(shù)K=2，K=3的中心頻率距離值有明顯的區(qū)別?？梢哉J(rèn)為此時(shí)屬于模態(tài)“接近”情況。此外，考慮到的運(yùn)算效率，VMD的分解層數(shù)也不應(yīng)該過大。因此，我們認(rèn)為模態(tài)數(shù)K=4時(shí)，VMD分解屬于過分解。因此，模態(tài)數(shù)K選擇3。K選擇3時(shí)，IMF1分量的VMD分解圖如圖4所示。

表1 IMF1分量VMD分解不同K值下的中心頻率

圖4 IMF1模態(tài)數(shù)K為3時(shí)VMD分解圖Fig.4 VMD decomposition diagram when IMF1 modal number K is 3

同理，對(duì)IMF2、IMF3分別進(jìn)行VMD分解，經(jīng)過分析不同分解層數(shù)K對(duì)應(yīng)的中心頻率，IMF2分量VMD分解層數(shù)K選擇2較為合適，IMF3分量VMD分解層數(shù)K選擇2較為合理。

至此，原始風(fēng)速數(shù)據(jù)經(jīng)過CEEMDAN-VMD二次分解后共有13個(gè)分量，其對(duì)應(yīng)的樣本熵值如圖5所示。

圖5 CEEMDAN-VMD二次分解各分量樣本熵值Fig.5 CEEMDAN-VMD quadratic decomposition of each component sample entropy

3.3 模型定階

建立預(yù)測(cè)模型之前，每個(gè)時(shí)間序列需要進(jìn)行平穩(wěn)性檢驗(yàn)。經(jīng)ADF(augmented dickey-fuller test, ADF)檢驗(yàn)，13個(gè)分量均為平穩(wěn)時(shí)間序列，分別建立ARMA模型。

表2 根據(jù)BIC準(zhǔn)則確定ARMA階數(shù)表

表3 CEEMDAN-VMD二次分解各分量模型階數(shù)表

3.4 異方差性檢驗(yàn)

表4 ARCH LM檢驗(yàn)表

表4中，統(tǒng)計(jì)量TR2服從卡方分布，對(duì)應(yīng)的值為15.883 9，小于置信水平5%對(duì)應(yīng)卡方分布的臨界值，不拒絕原假設(shè)，表明序列具有異方差性。

考慮到高階GARCH模型的參數(shù)估計(jì)過程較為復(fù)雜，因此對(duì)各分量ARMA模型擬合殘差建立GARCH(1,1)模型進(jìn)行修正。其余12個(gè)分量采用相同的方法進(jìn)行ARCH LM檢驗(yàn)。

3.5 誤差評(píng)價(jià)指標(biāo)

(15)

(16)

(17)

(18)

3.6 不同模型預(yù)測(cè)結(jié)果比較

本文模型CEEMDAN-VMD-ARIMA-GARCH處理原始風(fēng)速時(shí)間序列，得到3 h風(fēng)速預(yù)測(cè)值。為了對(duì)比驗(yàn)證所提方法的有效性，本文建立了對(duì)比模型一：CEEMDAN-ARIMA-GARCH預(yù)測(cè)模型；對(duì)比模型二：原始風(fēng)速-ARIMA-GARCH模型。原始風(fēng)速作為ARIMA-GARCH預(yù)測(cè)模型的輸入，得到風(fēng)速的預(yù)測(cè)結(jié)果。三種模型預(yù)測(cè)結(jié)果如圖6所示。

圖6 三種模型預(yù)測(cè)結(jié)果對(duì)比Fig.6 Comparison of prediction results of three models

圖6中，采用三種不同的數(shù)據(jù)分解策略，構(gòu)建相同預(yù)測(cè)模型得到不同的預(yù)測(cè)結(jié)果。觀察圖中曲線，采用原始風(fēng)速序列直接建立預(yù)測(cè)模型，雖然預(yù)測(cè)值能夠整體反映真實(shí)值的變化趨勢(shì)，但其結(jié)果始終存在較大的滯后；采用數(shù)據(jù)分解時(shí)，CEEMDAN分解、CEEMDAN-VMD二次分解的預(yù)測(cè)曲線都能夠有效地克服滯后，預(yù)測(cè)值始終能夠及時(shí)地反映出真實(shí)值的變化。但當(dāng)真實(shí)風(fēng)速值發(fā)生較大的突變時(shí)，采用CEEMDAN單一數(shù)據(jù)分解方法所得預(yù)測(cè)結(jié)果與真實(shí)風(fēng)速之間存在較大的偏差，在時(shí)刻點(diǎn)3、15處表現(xiàn)尤為明顯。比較前兩種方法，利用CEEMDAN-VMD二次分解所得預(yù)測(cè)結(jié)果既能克服預(yù)測(cè)存在的滯后現(xiàn)象，又能在風(fēng)速存在大幅變化的時(shí)刻有較為良好的預(yù)測(cè)效果。表5列出了本文所提模型的誤差指標(biāo)。為了更直觀地體現(xiàn)本文所提模型的有效性，在表5中用預(yù)測(cè)模型代替ARIMA-GARCH模型。

表5 模型預(yù)測(cè)誤差

從表5可以看出，直接采用原始風(fēng)速序列建立預(yù)測(cè)模型，其效果最差；CEEMDAN-ARIMA-GARCH模型預(yù)測(cè)效果次之；本文所提模型CEEMDAN-VMD-ARIMA-GARCH模型預(yù)測(cè)效果最好，能夠有效地提高預(yù)測(cè)精度。

4 結(jié) 論

本文針對(duì)風(fēng)速具有隨機(jī)性、非平穩(wěn)性的特點(diǎn)，提出了CEEMDAN-VMD-ARIMA-GARCH預(yù)測(cè)模型，有效地提高了風(fēng)速時(shí)間序列預(yù)測(cè)的準(zhǔn)確性。通過仿真對(duì)比實(shí)驗(yàn)，得到以下結(jié)論：

(1) CEEMDAN分解能夠消除模態(tài)混疊現(xiàn)象，具有優(yōu)良的分解效果；

(2) 以樣本熵判斷各分量的復(fù)雜程度，將復(fù)雜性較高的分量進(jìn)行VMD二次分解，可以進(jìn)一步降低信號(hào)的復(fù)雜性,與CEEMDAN分解和ARIMA-GARCH相結(jié)合的模型比較，平均絕對(duì)誤差減少了7.1%；

(3) 本文考慮了對(duì)原始風(fēng)速序列二次分解處理，引入GARCH模型用于修正時(shí)間序列模型殘差序列的異方差性，但是并沒有考慮GARCH模型的修正。因此，今后重點(diǎn)考慮對(duì)GARCH模型的修正進(jìn)而提高預(yù)測(cè)精度。