分形插值在風速時間序列中的應用

郭秀婷，朱昶勝，張生財，趙奎鵬

（1.蘭州理工大學計算機與通信學院，蘭州 730050；2.蘭州理工大學理學院，蘭州 730050）

0 引言

可再生能源發電具有不依賴化石燃料、減少溫室氣體排放、保護環境等一系列優點，因此得到了快速發展；其中風能作為一種低成本、分布域廣的清潔能源得到了各國的高度重視［1］。但風的間歇性和隨機性給風力發電的規劃及管理帶來了很多困擾［2-3］，而風速預測和風電功率預測是應對風電大規模并網運行問題的重要手段，因此風速預測和風電功率預測成為全球研究的熱點［4-6］。受天氣因素、地理環境因素、傳感器精度、電磁干擾以及通信故障等因素影響，采集到的風速數據會出現一定的異常值甚至一段時間的缺失值，但在研究過程中這些異常值和缺失值卻至關重要，因此研究適用于風速的異常值檢測及插值算法具有重要的意義［7-8］。

陳偉等［8］在處理風電異常數據時采用三次樣條插值進行數據重構，但自然風速時間序列是一個波動強烈的非線性時間序列，Lagrange、三次樣條和線性插值等傳統插值算法中，任意兩個信息點之間用直線或者光滑的曲線連接，掩蓋了相鄰插值點間的局部變化特征，因此用于波動強烈的風速時間序列效果并不理想；武艷強等［9］提出的多點三次樣條插值算法對于頻率相對穩定的低頻和中頻曲線插值效果比較理想，但對于包含多種頻率混合的非線性風速曲線卻不適合；楊茂等［10］提出用支持向量機重構風電數據，但由于參數的敏感性因素導致誤差較大；陳偉等［11］提出一種粒子群優化的最小二乘支持向量機對風速數據進行重構，但將缺失風速置為零后的時間序列參與模型訓練，如果連續缺失數據較多，則參與模型訓練的數據本身誤差較大，無法保證插值效果，且計算量大；自然風速時間序列具有典型的分形特征［12-13］；文獻［14-15］用分形插值算法及分段分形插值算法處理TurbSim 仿真軟件在特定條件下生成的模擬風速時間序列，這種理想風速時間序列不受瞬時天氣條件的影響，與實測風速有一定的差距。

本文將自適應粒子群算法優化的分形插值算法應用到甘肅省某風電場實測風速時間序列的連續缺失值的重構問題，實驗證明針對波動強烈、連續缺失數據量大的風速時間序列，分形插值算法可以取得比傳統插值算法更高的精度。

1 分形插值原理

分形插值的概念是由美國數學家Barnsley于1986年提出來的。分形插值［16］與多項式插值和曲線插值不同，它是由一類特殊的迭代函數系統（Iterated Function System，IFS）產生的，它是針對集合而非針對點的，能夠反映出相鄰已知相關點之間的局部特征，得到比傳統插值方法更高的精度。

1.1 數據集

假設點集 {(xi，yi)∈R2|i=0，1，…，N}，其中{(xi，yi)∈R2|i=0，1，…，N}，對應于這個數據集合的插值函數是一個連續函數f：[x0，xN]→R，使得：

點(xi，yi)∈R2稱為插值點，函數f內插該數據且f的圖形穿過該插值點。吸引子G={(x，f(x))；x∈[x0，xN]}是內插數據的連續函數的圖像。

1.2 構造IFS

考慮IFS{R2；Wn，n=1，2，…，N}，其中Wn是具有如下形式的仿射變換：

將已知點代入方程（4），不斷迭代，使系統收斂到其構造的吸引子G附近，即完成了插值運算。

1.3 分形維數

分形維數是刻畫分形特征的定量指標，它可以衡量兩個分形的相似程度，是衡量分形物體不規則程度的參數。對于風速時間序列，分形維數表示研究時間段內風速波動的強烈程度，分形維數越大表明風速波動越劇烈［17］。分形維數的計算遵循如下定理：假定數據集{(xi，yi)}不共線且則滿足G的分形維數D是如下方程的唯一常數解：

否則，滿足G的分形維數等于1。

式（5）表明分形維數D不依賴于因變量（數據集），而自由變量垂直比例因子dn較為關鍵，它能影響分形維數D的確定。目前常用盒維數法求解分形維數，即采用邊長為δi(i=1，2，…)的正方形的網格覆蓋分形體，數出覆蓋分形體所需盒子的數目N(δi)，隨著 limδi→0，得到一組數據(lnδi-1，lnN(δi))(i=1，2，…，M)，用最小二乘法擬合，得到的擬合直線的斜率即為分形體的分形維數。

2 自適應變異粒子群優化算法

2.1 粒子群優化算法

粒子群優化（Particle Swarm Optimization，PSO）算法是一種群體智能的優化算法，源于對鳥類捕食行為的研究。這種算法最早由Kennedy 和Eberhart 在1995 年提出，是一種非常成熟的求解優化問題的方法。該算法及其改進算法已在各種行業得到廣泛應用［18-20］，并取得了很好的效果。粒子群算法收斂速度快，具有很強的通用性。算法中每個粒子都代表問題的一個潛在最優解，用位置、速度和適應度值三項指標表示粒子的特征。位置用來存儲用來優化的變量，在本文中表示分形插值的垂直比例因子參數值，適應度值表示粒子的優劣，粒子的速度決定粒子移動的方向和距離，速度隨著自身及其他粒子的移動經驗進行動態調整。最初隨機初始化粒子，計算每個粒子的適應度值，根據適應度值更新個體極值和群體極值，然后粒子在解空間中運動，通過跟蹤個體極值和群體極值更新個體的位置和速度，粒子每更新一次位置，就重新計算一次適應度值，通過比較新粒子的適應度值和個體極值、群體極值的適應度值更新個體極值和群體極值，最終的群體極值即最優解。

每次迭代過程中，粒子通過個體極值和群體極值更新自身的位置和速度，公式如下：

其中：Vi為粒子的速度，Xi為粒子的位置，Pi為個體極值，Pg為群體極值，ω為慣性權重，c1和c2為加速度因子，是非負的常數，k為當前迭代次數，r1和r2為分布在［0，1］區間的隨機數，γ為約束因子。

其中ω的計算公式如下：

其中：ωstart和ωend分別為初始慣性權重和迭代至最大次數時的慣性權重，分別取值ωstart=0.9、ωend=0.4，Tmax為最大迭代次數。

本文中參數取值為c1=c2=1.49，γ=0.5，最大迭代次數Tmax=100，種群大小p=15。

2.2 自適應變異粒子群算法

粒子群算法雖然具有收斂快、通用性強的優點，但同時存在易早熟收斂、搜索精度低等缺點，因此，在本文中提出將遺傳算法中的變異思想引入PSO算法，即在粒子每次更新之后，引入變異算子，以一定的概率重新初始化粒子。變異操作拓展了在迭代中不斷縮小的搜索空間，使粒子能夠跳出當前搜索到的最佳位置，在更大的空間中開展搜索，保持了粒子的多樣性，同時提高算法搜索到更優值的可能性。本文中采用的變異算子公式如下：

其中r為在［0，1］區間的隨機數。

3 衡量分形插值結果的指標

3.1 插值誤差

插值曲線與實際曲線之間的誤差稱為插值誤差。本次實驗采用均方根誤差（Root Mean Square Error，RMSE）對插值結果進行評價：

3.2 納什效率系數

納什效率系數（Nash-Sutcliffe Efficiency coefficient，NSE），用于表示模擬結果與實際結果的吻合程度，本文用于驗證插值曲線與實際曲線的吻合程度，吻合程度越高表示插值效果越好，公式如下：

4 實驗結果與討論

所有實驗運行環境均為：1.60 GHz CPU（8 CPUs），8.00 GB RAM，Windows 10，Matlab R2019a。

4.1 實驗數據統計分析

采用甘肅省酒泉風電場2017 年1 月9 日09：55 至22：20的間隔5 min 的150條數據集和2017 年1 月10 日00：00 至14：30 間隔5 min 的175 條數據集進行實驗，簡稱Dataset A 和Dataset B，這兩組數據統計指標即：均值（Mean），最小值（Min），最大值（Max）和標準方差（Standard deviation，Std.），如表1和圖1所示，由此可知該組數據表現出非線性和較強的波動性。在兩組數據集中分別抽取30 個數據點和40 個數據點作為插值初始點。

表1 實驗數據集的統計指標Tab.1 Statistical indexes of experimental datasets

圖1 實驗數據集的曲線Fig.1 Curves of experimental datasets

4.2 參數的設定

4.2.1 垂直比例因子

由式（4）可知，有一個自由變量dn(|dn|＜1)，即垂直比例因子，它是影響分形插值精度的關鍵因素，當dn=0 時，分形插值函數就變成分段線性插值函數。圖2 為選取不同的垂直比例因子時圖形呈現的不同形態，由圖2可以看出隨著|dn|的增大，分形插值擬合曲線越來越粗糙，波動幅度越來越大，曲線的局部特征也更加明顯。分形維數是刻畫風速序列波動的強烈程度的指標，圖3 驗證了|dn|越大，分形維數D越大，即風速曲線波動越強烈。

圖2 垂直比例因子對曲線的影響Fig.2 Influence of vertical scaling factor on curve

圖3 分形維數隨垂直比例因子的變化Fig.3 Change of fractal dimension with vertical scaling factor

但|dn|并不是越大越好，即不代表插值曲線局部特征越復雜越好，圖4 為RMSE 隨dn取值的變化，由此可知，dn嚴重影響分形插值的誤差。

圖4 RMSE隨垂直比例因子的變化Fig.4 Change of RMSE with vertical scaling factor

本文采用自適應變異粒子群算法對垂直比例因子dn進行尋優，粒子群中每個粒子包含垂直比例因子，且以式（9）RMSE為適應度值，即最優垂直比例因子取值滿足以下公式：

4.2.2 迭代次數

表2 給出了迭代次數對分形維數、插值誤差、運行時間及插值點數的影響。由表2 可知，隨著迭代次數的增多，分形維數基本不變，即曲線波動的強烈程度不隨迭代次數的改變而改變，插值誤差也基本不受迭代次數的影響，但迭代次數越多，插值得到的數據點數越大，曲線細節會更豐富，但不會變光滑。圖5 給出了迭代次數分別為2 和4 時的圖形，由局部放大圖可以明顯看到迭代2 次和迭代4 次插值得到的數據點的密集程度的差別。迭代次數n與生成的插值點數m之間的關系［17］為：

其中q為初始插值點的個數。雖然迭代次數越多，插值得到的曲線細節越豐富，但該參數不宜取得過大，符合實際插值要求即可，迭代次數過大會增加計算時間，降低計算效率，本次實驗迭代次數選擇2次即可滿足要求。

圖5 迭代次數對插值點數的影響Fig.5 Influence of iteration number on interpolation point number

表2 迭代次數對各指標的影響Tab.2 Influence of iteration number on different indexes

4.3 分形插值與傳統插值的插值結果的對比研究

分別采用三次樣條插值、Lagrange 插值等傳統插值方法和分形插值算法進行比較分析，Dataset A 數據集的插值結果如圖6 所示。結果表明：分形插值較好地保持了圖形的形狀以及粗糙程度，能更好地凸顯圖形的局部細節特征；而Lagrange 插值和三次樣條插值只是保持了圖形的大致形狀，相鄰點間用光滑曲線連接，忽略了圖形局部豐富的信息。

表3是對Dataset A數據集插值結果的各種量化指標，三種插值結果的均值與原始數據相差均不大，對于表征曲線波動強烈程度的分形維數，Lagrange插值與分形插值相當，三次樣條插值最差；分形插值誤差最小為0.501 8，與Lagrange插值算法比較，減小了66.52%，與三次樣條插值算法比較，減小了58.57%；分形插值序列與原始序列的吻合程度同樣最佳，達到0.812 0，遠遠超過Lagrange插值和三次樣條插值兩種傳統插值方法。

圖6 Dataset A插值效果比較Fig.6 Interpolation effect comparison of Dataset A

由以上分析可知，對于波動強烈并且連續缺失數據較多的風速時間序列，分形插值不僅可以保持圖形的整體形狀和波動特性，而且插值精度更高，與原始曲線吻合更好，與傳統的插值算法相比，分形插值具有明顯的優勢。為了進一步證實上述結論的可靠性，對Dataset B數據集插值計算，得到插值結果如圖7所示。

表3 Dataset A 插值性能對比Tab.3 Interpolation performance comparison of Dataset A

圖7 Dataset B插值效果比較Fig.7 Interpolation effect comparison of Dataset B

表4 給出了Dataset B 數據集的各插值結果。Dataset B 數據集的各插值誤差RMSE 和納什效率系數NSE 的規律與Dataset A 的規律相同，進一步說明了分形插值比傳統的插值方法更適合連續缺失數據較多且波動強烈的風速序列。

表4 Dataset B 插值性能對比Tab.4 Interpolation performance comparison of Dataset B

5 結語

自適應變異粒子群優化的分形插值算法為風速序列的連續缺失數據插值提供了一種新途徑，并對分形插值算法中的兩個重要參數進行討論分析，利用自適應變異粒子群算法對垂直比例因子尋優，最后對分形插值算法與傳統插值算法的進行了對比研究。通過計算和驗證得出如下結論：

1）垂直比例因子是影響分形插值精度的重要因素，采用自適應變異粒子群優化垂直比例因子，分形插值取得了很好的插值效果，證明了其有效性。

2）針對連續缺失的波動強烈的風速序列，分形插值比傳統插值方法更有優勢，分形插值不僅可以保持風速曲線的整體變化趨勢，而且可以反映數據的局部波動特征，并且具有更高的插值精度和吻合度。