基于分形插值的ARIMA大壩預警模型

屠立峰 包騰飛 李月嬌 趙 斌

（1.河海大學 水利水電學院，南京 210098；2.河海大學 水資源高效利用與工程安全國家工程研究中心，南京 210098；3.南京南瑞集團公司 國際公司，南京 210093）

大壩的變形、應力和作用荷載等觀測資料均為典型的時間序列，故可通過過去的觀測值來預估未來的觀測值.對于非平穩的時間序列，差分自回歸移動平均模型（ARIMA模型）具有良好的適應性，常作為經典的時間序列模型［1］.時間序列覆蓋的歷史越長，可靠性也就越高，但也越有可能包含缺失值，所以對該組數據進行分析時，預測和擬合數據與實際監測值會有較大的誤差.對于缺失率較低的情況可用現有的數據進行分析，然而大量信息的丟棄會使數據分布產生偏斜而扭曲數據挖掘和數據分析，從而誤導決策，為此需通過數據填補來解決.通過函數插值可以彌補時間序列存在缺失值而導致預測誤差較大的弊端，然而傳統的插值函數由于兩端的插值區間易出現“龍格現象”而限制了在擬合和數據處理中的應用，基于此學者們提出了分形插值的方法.分形插值是基于分形幾何學提出的迭代函數系統，分形幾何學指出任何一個局部都與整體自相似或統計自相似，因此可將已知數據插值成具有自相似結構的曲線和曲面.研究表明，分形插值與傳統的插值方法相比具有更高的精度和適應性［2］.考慮到大壩監測數據間復雜的動態性和非線性的特點，本文為彌補ARIMA模型對含有時間缺失值的序列預測失準的弊端，對缺失的時間序列進行分形插值，以提高模型的預測精度.

1 分形插值

1.1 分形插值原理

一個數據集合是形如｛（xi，yi）∈R2｜i＝0，1，2，…，N｝的點集，其中

插值函數對應的這組數據是一個連續函數f：［x0，xN］→R，如

點（xi，yi）∈R2稱為插值點，f（x）被稱為插值函數.

傳統的插值函數一般由初等函數的一組基函數線性表觀的，相鄰的插值點只可用直線或光滑的弧線銜接，卻得不到兩點之間的部分情況.分形插值函數是由一類特殊的仿射變換生成的，它可以得到兩個相鄰的信息點之間的局部變化，它為描述了一個不規則、隨機曲線擬合的實驗數據提供了強大的工具［3］.

1.2 分形插值方法

給定數據集｛（xi，yi）｜i＝0，1，2，…，N｝，構造迭代函數系統（IFS）｛R2；Wn，n＝1，2，…，N｝，其中Wn是具有如下型式的仿射變換：

并且

為了確保各小區間的不交迭，令式（3）中的bn＝0.式（4）具體可寫為

式（5）中有4個方程5個參數，所以其中一個是自由參數.在一般情況下，選擇dn為自由參量，稱為變換Wn的垂直比例因子.為保證IFS收斂，令｜dn｜＜1，解方程組（5），并令L＝xN－x0，則

2 基于分形插值的ARIMA模型

2.1 ARIMA模型定義

ARIMA模型是指經過差分將非平穩的時間序列轉化為平穩的時間序列，而后根據轉化后得到的平穩時間序列創建ARIMA模型.該模型認為：系列數據隨著時間的推移，預測形成隨機序列，利用近似數學模型來描述這個隨機序列，識別后的模型就能利用時間序列的過去值對未來進行預測［4］.

ARIMA（p，d，q）模型數學表達式如公式（7）所示：

式中，φm（m＝1，2，…，p）是自回歸模型的系數，θj（j＝1，2，…，q）均滑動模型的系數，p為自回歸階數，q為滑動平均部分的階數，at白噪聲序列.

2.2 基于分形插值的ARIMA模型建模流程

大壩安全監測數據一般為非平穩時間序列，模型建立步驟如下：

1）判斷給定的時間序列是否含有缺失值，若為有缺損值的時間序列需對該組數據進行分形插值，得到數據填補后的序列樣本.反之則跳轉至2）.

2）進行平穩性檢驗得到參數d，即差分階數.平穩性檢驗常用的方法有兩種：一種是圖檢驗法，另一種是構造檢驗統計量進行假設檢驗的方法［5］.圖檢驗法是依據自相關圖的平穩性檢驗方法，因其簡便和運用廣泛的特點，采用該方法將大壩安全監測量轉化為平穩的時間序列.

3）確定模型類型后，使用赤池信息量準則（AIC）或貝葉斯信息準則（BIC）判斷模型的p、q值.由于BIC準則比AIC準則具有更好的收斂性［6］，故本文采用BIC準則計算模型的p、q值.其數學表達式為

式中，k為參數的數量，n為觀察數.

式（8）表明BIC準則主要有兩部分組成，一是參數的數量，隨著階數的增大而增大；二是模型的擬合情況［7］.增加自由參數的數目可提高了擬合的優良性，BIC鼓勵數據擬合的優良性，但是盡量避免出現過度擬合的情況［8］.所以選擇BIC值最小組對應的p、q值.

第4步：在上述模型識別的基礎上，采用極大似然估計或最小二乘估計法進行參數估計.

第5步：進行假設檢驗，診斷at是否為白噪聲序列.

第6步：根據建立好的ARIMA模型對大壩變形參數進行擬合和預測.

3 工程應用實例

小灣水電站為混凝土雙曲拱壩，壩頂高程1 245 m，最大壩高293.5m，總裝機容量420萬kW.以小灣壩頂20120701～20120919觀測得到的數據為例，其中 20120706～20120710，20120726～20120730，20120816～20120825觀測數據缺失.首先利用Matlab數學軟件對缺失的數據進行分形插值，再利用SPSS軟件對得到的連續時間序列建立ARIMA模型進行擬合預測，最后與原始的ARIMA模型進行對比.本文采用小灣拱壩壩頂位移的原始觀測值前70組進行擬合，后10組數據進行預測.

3.1 分形插值

通過Matlab數學軟件對觀測值中的缺失值進行分形插值，啟動Matlab軟件，在命令窗口中做如下操作：

通過計算可以得到20120706～20120710，20120726～20120730，20120816～20120825值插值后的位移值.表1為經過分形插值后的數據.

表1分形插值數據表

3.2 基于分形插值的ARIMA模型

通過差分將壩頂觀測位移平穩化，下面利用圖檢驗法分析序列的平穩性以確定模型的階數.一階差分后的自相關系數（ACF）和偏自相關系數（PACF）如圖1所示.

圖1 ACF（左）和PACF（右）系數圖

由圖1可看出，當滯后步數大于等于2時，自相關系數和偏自相關系數均落在置信區間內，因此可以認為經一階差分后時間序列是平穩的.表2為從低階到高階不同的BIC值.

表2 ARIMA模型BIC準則值

由表2可知，ARIMA（1，1，3）模型對應的BIC值最小，故其相應的精度就越高.運用相同的方法對未經分形插值的時間序列建立ARIMA（0，1，1）模型，分別對原始監測數據進行擬合，并計算相應的誤差，擬合結果見表3.

表3 兩種模型的擬合值與相對誤差

續表3 兩種模型的擬合值與相對誤差

由表3可看出，除個別觀測日期外，基于分形插值的ARIMA模型的擬合值更接近于實測值，利用基于分形插值的ARIMA模型的相對誤差小于原始的ARIMA模型.分別利用ARIMA模型和基于分形插值ARIMA模型對實測值展開預測，預測結果見表4.

表4 兩種模型的預測值與相對誤差

由表4可看出，基于分形插值的ARIMA模型的預測值更接近于實測值，利用分形插值的ARIMA模型的相對誤差都小于原始ARIMA模型的相對誤差.

將表3～4計算所得兩種模型的擬合數據和預測值與實測值進行比較，比較結果如圖2所示.

圖2 兩種模型的擬合值與實測值的比較

通過表3～4以及圖2可以得到，對于含有缺失值的時間序列，基于分形插值的ARIMA模型的預測值的平均相對誤差由原先的3.56%降低至1.34%，表明該模型的預測精度顯著提高.

4 結 論

1）在實際監測過程中，數據缺失是常見的現象，本文以小灣水電站壩頂位移預測為例，建立了基于分形插值的ARIMA模型，利用分形插值法對含有數據缺失的時間序列進行插值，從而建立ARIMA時間序列模型，得到壩頂位移預的擬合預測.計算結果表明，改進后的ARIMA模型的預測精度明顯提高.

2）通過預測值和實際值比較發現，由于ARIMA模型對于單調數據處理的不足，易出現預測值單調迅速增長的現象，導致預測值均大于實測值.針對ARIMA模型的不足，近年來學者們已經陸續提出了基于ANN－ARIMA、GM－ARIMA等應用模型.基于此，對于含有一定數據缺失的時間序列可將分形插值與優化后的ARIMA模型相結合，進一步的提高預測精度，不失為一種新的研究方向.

［1］ 岳莉莉.基于時間序列分析的風速短期預測方法研究［D］.北京：華北電力大學，2012.

［2］ 李信富，李小凡.分形插值與拉格朗日插值的比較研究［J］.黑龍江大學自然科學學報，2008（3）：323－326，331.

［3］ 孫洪泉.分形幾何與分形插［M］.北京：科學出版社，2011.

［4］ 王正宇，王紅玲.基于ARIMA模型的我國GDP分析預測［J］.對外經貿，2011，12：107－108.

［5］ 王 燕.應用時間序列分析［M］.北京：中國人民大學出版社，2005.

［6］ Olivier C，Courtellemont P，Colot O.Comparison of Histograms：a Tool for Detection［J］.European Journal of Diagnosis and Safety in Automation，1994，4（3）：335－355.

［7］ 胡效雷，何祖威.基于GM－ARMA組合模型的年電力需求預測［J］.廣東電力，2007（2）：10－13.

［8］ 馮龍龍，李 星，李曉晨，等.GM－ARIMA模型在大壩安全監測中的應用［J］.三峽大學學報：自然科學版，2013，35（5）：7－10.