劉嘉誠，劉 俊，2，趙宏炎，吳澤陽，劉曉明，吳 柳

（1. 陜西省智能電網(wǎng)重點實驗室（西安交通大學電氣工程學院），陜西省西安市 710049；2. 人工智能與數(shù)字經(jīng)濟廣東省實驗室，廣東省廣州市 510320）

0 引言

近年來，隨著節(jié)能減排政策的不斷推進、光伏發(fā)電技術的日益成熟以及發(fā)電成本的不斷下降［1］，中國光伏產(chǎn)業(yè)迅速發(fā)展，裝機容量不斷攀升，截至2019 年底，全國光伏裝機總量已達204.3 GW，且全年光伏發(fā)電量達到224.3 TW ?h［2］。然而，由于光伏發(fā)電本身受天氣等自然環(huán)境因素影響，其出力存在較高的不確定性與波動性，這對系統(tǒng)的安全穩(wěn)定運行、經(jīng)濟調度以及電能質量治理帶來巨大的威脅與挑戰(zhàn)［3-4］，準確預測光伏出力具有重要意義。

光伏出力預測按照時間尺度可分為超短期（未來4 h）、短期（未來72 h）與中長期預測（未來1 月至1 年）［5］；按預測形式又可分為點預測、區(qū)間預測與概率預測。在點預測研究中，時間序列方法、回歸分析方法等統(tǒng)計學方法因其建模簡單且適應性廣得到了廣泛應用［6-7］，但對大量數(shù)據(jù)的依賴性也增加了數(shù)據(jù)處理的難度，可靠有效的特征處理對于挖掘歷史數(shù)據(jù)中的規(guī)律、提升光伏預測精度十分必要。

特征工程從樣本數(shù)據(jù)的特征集出發(fā)，通過對各類特征的構建、提取或算法選擇，挖掘與光伏出力相關性最大的特征子集合，降低無關或弱關聯(lián)氣象特征對于模型的干擾，在簡化模型的同時改善其預測性能。文獻［8-10］構建了風速3 階差分最大值、輻照度最大值比、輻照度趨勢分量等特征；文獻［11］構造包括數(shù)值指標、增量指標、波動指標在內的瞬時氣象模式特征以及時窗趨勢特征；此類特征的構建方法往往依賴大量的專業(yè)知識并且可靠性與有效性難以保證。文獻［12］通過深度神經(jīng)元層將氣象特征進行非線性映射，盡管該方法可以構造出新的特征空間，但是訓練統(tǒng)一適用的深度網(wǎng)絡難度較高，且氣象數(shù)據(jù)與光伏功率之間的關聯(lián)性難以直接解釋。文獻［13-15］采用主成分分析（PCA）算法對特征量進行處理，該算法僅為特征的線性變換與降維，具有一定的局限性。文獻［16-17］分別采用相關系數(shù)和最小絕對值收斂與選擇算子（least absolute shrinkage and selection operator，LASSO）回歸進行特征篩選，然而這兩者只能反映變量之間的線性相關性，與實際情況差距較大。文獻［18］將最大相關最小冗余（max-relevance and min-redundancy，mRMR）算法引入光伏預測的特征篩選環(huán)節(jié)中，以氣象特征與光伏輸出兩兩之間的互信息熵作為特征評價指標，盡管mRMR 算法簡化了多維變量的互信息計算，但仍然具有初始特征固定、特征集需要遍歷尋優(yōu)等缺點。

針對現(xiàn)有光伏出力預測研究中特征相關性分析以及特征篩選方法的不足，本文以互信息為特征評價指標，提出了一種基于擴散核密度估計（diffusion kernel density estimation，DKDE）和改進mRMR 特征選擇算法，并建立了基于XGBoost（eXtreme gradient boosting）算法的高精度光伏發(fā)電預測模型。主要創(chuàng)新性貢獻包括：①針對傳統(tǒng)互信息計算中變量離散化時未考慮變量實際概率分布的問題，提出基于DKDE 的區(qū)間劃分方法，根據(jù)樣本落入各個區(qū)間的概率計算互信息以表征變量間相關性；②提出一種改進的mRMR 特征選擇算法，并行選擇多個初始特征進行特征子集的擴充，并以特征子集的總信息熵作為其權重系數(shù)，以充分利用多特征子集的信息；③在建立多個特征子集的基礎上，基于XGBoost 算法建立短期光伏出力預測模型，與神經(jīng)網(wǎng)絡、支持向量機回歸（SVR）等其他機器學習方法相比，其在超參數(shù)調整、訓練速度等方面都具有顯著優(yōu)勢。最后，選取美國俄勒岡州本德市的一座實際光伏電站數(shù)據(jù)進行模型的訓練與測試，驗證所提出的特征選擇方法與預測模型的有效性與準確性。

1 光伏預測模型概述

本文光伏出力預測的框架如圖1 所示。

圖1 光伏出力預測框架Fig.1 Framework of photovoltaic output forecasting

預測模型建立的主要步驟如下。

1）數(shù)據(jù)預處理

對于給定的光伏電站歷史運行數(shù)據(jù)與氣象監(jiān)測數(shù)據(jù)，根據(jù)數(shù)據(jù)標簽對異常日數(shù)據(jù)進行剔除，對于缺失值采用前后采樣點的均值替代，最后對所有數(shù)據(jù)進行最大最小歸一化處理，避免數(shù)量級差異對預測造成的不利影響。

2）特征集篩選

基于DKDE 對互信息的計算進行改進，并以此作為氣象特征與光伏出力間相關性的量度；通過改進的mRMR 特征選擇算法并行構建多個特征子集，同時計算出各個子集的總互信息熵作為其得分，根據(jù)得分保留部分最優(yōu)子集。

3）模型構建與光伏預測

分別對各個特征集構建基于XGBoost 算法的光伏預測模型。然后，輸入待預測日各時刻氣象信息得到各模型的預測輸出。最后，將各輸出值按特征集得分加權求和，即為待預測日的光伏發(fā)電功率。

2 基于DKDE 的互信息計算

作為變量之間相關性的有效量度，互信息能夠表征某一特征對另一變量概率分布的可解釋性，而不僅僅局限于線性或某種特定的非線性關系。因此，特別適用于描述氣象參數(shù)與光伏出力之間復雜關聯(lián)性的強弱程度。

給定變量X和Y的聯(lián)合概率分布為p(x，y)，其中，x和y分別表示X和Y中的任意變量，則變量X和Y之間的互信息熵可表示為：

式中：p(x)和p(y)分別為x和y的概率密度函數(shù)（probability density function，PDF）。

在絕大多數(shù)情況下，p(x)、p(y)與p(x，y)都難以精確獲取，故需要對變量取值進行區(qū)間劃分使其離散化，即使連續(xù)變量的取值數(shù)量為有限個，則互信息熵計算式可轉換為：

對于離散變量，互信息可通過式（2）直接計算，而當變量X和Y中有一個取值連續(xù)時，首先需要獲得該變量的概率密度以及聯(lián)合PDF。對此問題通常有2 種處理方法：①將連續(xù)變量分箱處理，將變量取值劃分為有限個區(qū)間，分別統(tǒng)計變量落入各個區(qū)間的概率，該方法結構簡單、計算量小，但是區(qū)間劃分難以界定，大多數(shù)采用等距劃分；②通過概率密度估計或者窗函數(shù)［19-21］估算變量之間的聯(lián)合PDF，但是這一過程中核函數(shù)帶寬通常難以選取，從而導致結果往往與實際有較大偏差。對此，本文提出一種基于DKDE 的計算方法，首先通過DKDE 分別估計出單個變量的概率密度分布，然后從PDF 中提取波峰作為變量的典型模式，以波谷作為取值區(qū)間邊界將變量離散化，從而通過式（2）計算氣象參數(shù)與光伏出力的互信息。

給定一維隨機變量X的N個采樣樣本x1，x2，…，xi，…，xN，X的高斯核密度估計定義為：

式中：φ(?)為高斯核函數(shù)，其定義如式（4）所示；t為核函數(shù)帶寬的平方。

文獻［22］指出，核函數(shù)也可以表示為如式（5）所示的線性擴散偏微分方程的解：

其中，a(?)和p(?)均為在定義域內二次可導的正值函數(shù)，且g(x；t)需滿足如下初始條件。

式中：δ(?)為Dirac 函數(shù)。

當a(?)和p(?)均取值為常數(shù)1 時，即得到擴散偏微分方程的一個特例，即

在偏微分方程框架下，式（3）即為式（7）在初始條件（式（6））下的唯一解。若考慮x的取值范圍為有界，假定其定義域限制為[0，1]，則對式（7）的求解還需要添加Neumann 邊界條件，即

求解后得到的擴散高斯核密度估計仍為核函數(shù)累加的形式，即

擴散高斯核密度估計過程如下。

DKDE 能夠自適應地根據(jù)采樣樣本自動計算最優(yōu)帶寬t*，避免了人為隨機選取的不確定性，且文獻［22-24］的研究表明，DKDE 能夠準確擬合變量的實際概率分布。圖2 以溫度這一氣象參數(shù)為例，展示了其概率分布的核密度估計與區(qū)間劃分情況。通過DKDE 計算得到的概率分布為近似正態(tài)的多峰分布。當帶寬t*選擇過大時，容易忽略局部信息而近似為正態(tài)分布；而當t*選擇過小時，則會導致采樣的偶然性與隨機性被放大。

圖2 2019 年1 月美國本德市溫度概率分布與區(qū)間劃分情況Fig.2 Probability distribution and interval division oftemperature in Bend,US,in January 2019

離散互信息用于反映變量處于不同模式的關聯(lián)性，其值取決于變量離散化的方式。考慮到典型模式下的變量取值通常具有統(tǒng)一集中的分布規(guī)律，反映在概率密度分布圖像上便是連續(xù)區(qū)間內呈現(xiàn)出概率峰，而在模式邊緣處取值的概率相對較低。因此，本文的區(qū)間劃分方法上充分考慮了特征的概率分布情況，即以概率峰作為變量的典型模式，以波谷作為區(qū)間邊界。相比之下，傳統(tǒng)的等距劃分只是為了方便互信息計算所采取的簡單離散化，不能體現(xiàn)這一規(guī)律性。

3 基于改進mRMR 算法的特征篩選

特征選擇是光伏發(fā)電預測數(shù)據(jù)處理與模型構建的重要步驟，維度過高、冗余過多的特征空間容易導致氣象信息與輸出功率的映射規(guī)律難以挖掘，增加模型的訓練難度與計算復雜度。從信息角度出發(fā)，特征選擇的目的在于選擇包含輸出變量盡可能多的信息且維數(shù)最少的特征子集（X1，X2，???，Xn），即同時滿足maxI(X1，X2，???，Xn；Y)與minn，其中多變量的互信息可通過式（12）計算［19］。

式中：Sn為n個特征組成的特征集合。

由于采樣數(shù)量限制以及高維矩陣求逆計算過于復雜，聯(lián)合PDF（p(x1，x2，???，xn)）無法精確獲取，因此直接計算多變量間的互信息難以實現(xiàn)，對此mRMR 算法將特征篩選近似等效為如下優(yōu)化問題。

式（13）中，第1 項代表所選氣象特征集Sn與光伏輸出Y的相關性，第2 項代表氣象特征之間的信息冗余度。通常mRMR 采用增量搜索算法［19］，每次挑選當前的最優(yōu)特征加入候選特征子集，即新增特征Xq滿足：

式中：SSh?1為全體特征集S中剔除已選特征集Sh?1后剩余的待篩選特征集合。

mRMR 特征選擇過程如下。

1）給定一個足夠大的整數(shù)n，根據(jù)式（14）逐次選擇特征生成n個序列特征集S1?S2?…?Sn?1?Sn。

2）根據(jù)n個特征集分別訓練模型進行預測，得到預測誤差序列E={e1，e2，???，en}，選擇相鄰l個均值、方差均較小的誤差組成誤差集Ω。

3）在Ω中搜索最小值，其對應的特征集即為最優(yōu)特征集。

盡管mRMR 算法已將互信息計算大大簡化，但在選擇特征時需要進行多次嘗試且不斷訓練模型，增加了其計算復雜程度。而且，由于增量搜索過程中始終固定某一氣象參數(shù)作為初始特征，可能潛在地導致后續(xù)特征選擇過于局限，因此本文對mRMR算法作出以下改進。

1）對于初始的特征空間S，分別以每一個特征作為初始特征進行增量搜索，構建m個特征子集S1，S2，???，Sm。

2）增量搜索過程在新特征無法提供新的信息時停止，即式（14）取值非正時停止搜索。

3）以特征集構建過程中的總互信息增量作為特征集評價標準，歸一化后經(jīng)softmax 函數(shù)處理作為最終得分，最后保留得分高于平均值的v個特征集并重新標記為S1，S2，…，Sv，對應的得分為V1，V2，…，Vv。

4）分別對v個特征集進行模型訓練，得到v個預測器P1，P2，???，Pv，對于任意給定的輸入樣本x=[x1，x2，???，xn]，選擇特征集對應的特征輸入各個預測器中，得到預測輸出P1(x)，P2(x)，…，Pv(x)，模型的最終預測結果由預測器得分加權得到，即

4 基于XGBoost 算法的光伏預測模型

XGBoost 算法［25-26］是一種Boosting 類的集成學習算法，通過多輪迭代訓練多個弱分類器并將其組合得到最終的強分類器。XGBoost 算法一般采用分類與回歸樹（classification and regression tree，CART）作為基本模型，對于一個回歸問題，可構建如式（16）所示的集成映射模型。

式中：x∈Rm為m維的樣本輸入；fr為第r棵回歸樹，其輸出fr(x)等于樣本x歸入的葉節(jié)點得分。R棵級聯(lián)樹得分累加即得到最終的預測輸出y?R。

式中：xi和yi分別為第i個樣本輸入與輸出；loss(?)為單個樣本預測輸出與實際輸出的損失函數(shù)，通常選擇為平方誤差；TR為第R棵回歸樹所包含的葉節(jié)點數(shù)；ωs為第s個葉節(jié)點的得分；γ和λ為自定義權系數(shù)，一般在［0，1］區(qū)間內取值。

式（17）的第1 項表示預測與實際輸出的誤差，第2 項則用于限制樹的復雜度，避免過擬合。將損失函數(shù)LR在函數(shù)空間泰勒展開并省略高次項可得：

式（18）將損失函數(shù)轉換為葉節(jié)點得分的凸函數(shù)，為使回歸樹的損失達到最小，第s個葉節(jié)點的得分ωs應取為：

考慮到回歸樹的生成過程中，無法對其所有結構進行遍歷，只能逐層進行根節(jié)點的二分裂，則分裂特征與分裂閾值應使得式（20）所示的目標函數(shù)達到最大。

式中：SL為二分裂后歸屬于左半枝的樣本集合；SM為歸屬于右半枝的樣本集；Sroot=SL∪SM。值得注意的是，特征與閾值可通過遍歷的方式進行選取，且分裂后的葉節(jié)點得分仍可通過式（19）計算，分裂在回歸樹達到預設的最大深度（一般不超過12）時停止。XGBoost 算法模型示例見附錄A 圖A1。

5 算例分析

本文選取2018—2019 年美國俄勒岡州本德市的一座光伏發(fā)電站歷史監(jiān)測數(shù)據(jù)［27］為研究對象，時間間隔為1 min?？紤]到不同季節(jié)的復雜氣象情況，分別對各月進行單獨建模。在模型構建過程中，采用2018 年全月以及2019 年對應月前3 個星期的日間數(shù)據(jù)進行特征提取與XGBoost 模型訓練，并以2019 年當月最后1 個星期的日間實測數(shù)據(jù)作為測試集，根據(jù)數(shù)值天氣預報（numerical weather prediction，NWP）數(shù)據(jù)預測光伏功率并與實測出力對比分析?？紤]到光伏預測時氣象數(shù)據(jù)需要由氣象站或電站的氣象工作人員預測給出，必然存在一定誤差，故本文對當月最后1 個星期測試集的輸入氣象參數(shù)均疊加?10%～10%的隨機擾動。為驗證所提出的特征提取方法與基于XGBoost 算法的預測模型的有效性，對比研究了多種不同的特征處理方法與不同預測算法對預測精度的影響。

5.1 特征集篩選結果

光伏電站實測特征數(shù)據(jù)包括太陽天頂角（solar zenith angle，SZA）、太陽方位角（solar azimuthal angle，AZM）、地外輻射（extraterrestrial irradiance，ETR）、水平總輻射（global horizontal irradiance，GHI）、傾斜面總輻照度（global tilted irradiance，GTI）、風速、溫度以及電池板溫度共計8 維特征。由于光伏出力可能與時間存在潛在關聯(lián)性，額外增加年、日、時、分4 維特征數(shù)據(jù)，構成總計12 維的特征空間。

以2018 年1 月為例，各特征變量與光伏功率分布的互信息計算結果如附錄B 圖B1 所示，主對角線為變量本身所包含的信息熵，其值越大表明變量分布的不確定性越大；其余元素為變量之間的相關性量度，對應網(wǎng)格顏色越深表示變量之間的可解釋性越強。

在不同月份，由于氣候差異，不同氣象特征對光伏輸出功率的影響可能不同，且改進的mRMR 特征選擇過程需要并行構建多個特征子集，因此建模過程中各個特征的重要性存在較大差異。附錄B 圖B2 統(tǒng)計了各個月份、各氣象參數(shù)的被選擇頻次數(shù)，其中日標簽、ETR、GHI、GTI 以及電池板溫度與光伏輸出功率密切相關，廣泛應用于各個模型中；時、分與SZA 在夏季與秋季應用廣泛，對于預測具有重要作用，而在春、冬兩季較少被選入特征子集中；月標簽、AZM、風速以及溫度與光伏功率的相關性較弱，只在個別月份有利于預測。

5.2 不同方法的預測精度對比分析

為驗證所提出模型構建方法的有效性，本文對比分析了采用不同特征提取方法與預測算法時的光伏預測結果，特征方法包括本文提出的改進mRMR與DKDE 結合、傳統(tǒng)mRMR、無特征處理以及Lasso特征選擇方法，預測算法則包括XGBoost、深度神經(jīng)網(wǎng)絡（DNN）、SVR 以及最近鄰回歸（kNN）。光伏預測模型的預測精度可通過平均絕對百分比誤差（MAPE）、平均絕對誤差（MAE）與均方根誤差（RMSE）3 個指標（eMAPE、eMAE、eRMSE）進行評價，各指標計算如式（21）—式（23）所示。

圖3 展示了采用不同特征方法與預測算法組合時各月最后1 個星期預測結果的eMAPE分布。在采用XGBoost 算法進行預測時，本文所提出的改進mRMR 與DKDE 結合的特征方法在各月均有較好表現(xiàn)。在夏、秋兩季，eMAPE一般在6%左右，相比于未經(jīng)特征處理時下降約1.5%；而在春季與冬季等氣象情況較為復雜的月份，eMAPE不超過16%；且本文所提出特征方法優(yōu)勢明顯。在預測精度最低的2 月，改進mRMR 與DKDE 結合可使eMAPE從19.3%（未經(jīng)特征處理）降低至15.2%，降幅達4.1%。圖3同時也對比了改進mRMR 與DKDE 單獨應用的場景，結果表明改進mRMR 與DKDE 的應用均對模型精度的提升有益，如在1 月僅應用改進mRMR 可使eMAPE從無特征處理時的13.5%降低至12.1%，在此基礎上再引入DKDE 可使其進一步降低至9.8%。對比DNN、SVR 與kNN 等機器學習預測算法，基于XGBoost 算法的模型預測誤差在大多數(shù)月份都顯著低于其他對比模型的預測誤差。除預測效果外，附錄B 表B1 還展示了不同預測算法的超參數(shù)數(shù)量以及模型的平均訓練時間。在相同的軟硬件配置情況下，XGBoost 算法的訓練時間相較于SVR 與DNN可縮短一個數(shù)量級。盡管kNN 無須進行參數(shù)訓練，但在每次預測時卻需要計算待預測日氣象參數(shù)與訓練集中所有樣本的距離并排序，在訓練集樣本數(shù)量較多以及預測較為頻繁時計算量十分巨大。故綜上所述，XGBoost 算法在預測精度、超參數(shù)調整及訓練速度方面均具有較為顯著的優(yōu)勢。

圖3 不同特征方法與預測模型對應的各月MAPEFig.3 Monthly MAPE corresponding to different feature methods and forecasting models

各個季度首月預測結果的MAE 與RMSE 值如附錄B 表B2 所示，均為光伏出力歸一化至［0，1］范圍后的計算誤差。與eMAPE分布情況類似，采用本文提出的DKDE 與改進mRMR 結合的特征提取方法以及XGBoost 算法所構建的光伏預測模型在所有對比模型中的eMAE與eRMSE均為最低。

5.3 DKDE 帶寬選擇對預測精度的影響分析

考慮到核函數(shù)帶寬的選擇直接影響到變量概率分布的估計結果，從而導致區(qū)間劃分、互信息計算以及特征篩選結果不同，并最終影響光伏預測的精度。本節(jié)對不同帶寬下的光伏預測結果進行驗證分析，圖4 以DKDE 得到的帶寬t*作為基準，對比了帶寬分別為100t*、0.01t*以及等距劃分區(qū)間條件下的各月預測eMAPE；其中等距劃分可視作核函數(shù)帶寬取值為無窮大的情形，此時任意變量均為取值范圍內的均勻分布。由圖4 可以看出，DKDE 的帶寬選擇在大多數(shù)月份均優(yōu)于其他情況，且在春、冬季等預測誤差較大的月份優(yōu)勢更為顯著。

圖4 不同帶寬估計對應的光伏預測精度Fig.4 Photovoltaic forecasting accuracy corresponding to different bandwidth estimation

5.4 單日光伏出力預測分析

在春冬等天氣復雜多變的季節(jié)，光伏預測精度通常較低，而在夏秋等天氣狀況良好的季節(jié)，模型通常表現(xiàn)較好。為驗證本文光伏預測模型在不同天氣狀況下的可靠性，在不同月份選取了2 類典型的天氣狀況進行測試。圖5（a）展示了天氣晴朗時（8 月24 日）的光伏預測情況，采用改進mRMR 與改進互信息計算構建的預測模型全日eMAPE為4.9%，相比于傳統(tǒng)mRMR 與不作特征處理分別有0.3% 與0.5%的改善。對于多云天氣，圖5（b）展示了預測精度最低的2 月中的一天（2 月23 日）各模型的預測結果，采用全部特征所構建的模型單日預測eMAPE為22.2%，而通過本文的特征篩選方法，eMAPE降低至14.8%，本文所提出的方法在復雜多云天氣的優(yōu)勢更為顯著。

圖5 典型日光伏出力預測曲線Fig.5 Photovoltaic power output forecasting curves on typical days

6 結語

本文重點研究了光伏出力預測中的特征選擇問題，并在此基礎上建立了適用于不同季節(jié)/月份與各類天氣的光伏預測模型。針對特征相關性量度——互信息難以直接計算的問題，采用DKDE 計算各個特征的PDF，根據(jù)其在定義域內的實際分布情況劃分取值區(qū)間并將連續(xù)互信息計算離散化。針對傳統(tǒng)mRMR 特征算法需要大量嘗試且初始特征選擇唯一的缺陷，本文對增量搜索過程中的終止條件進行了改進，提出一種可篩選多特征集的改進mRMR 算法。最后，采用XGBoost 集成學習算法對短期光伏出力進行預測。光伏電站實際運行數(shù)據(jù)測試實驗表明，在不同月份，本文所提出的特征篩選方法相比于其他方法均具有較大優(yōu)勢，預測精度得到顯著提升。

在后續(xù)研究工作中，可進一步考慮對天氣狀況的精細化聚類劃分，針對各種復雜天氣類型分別建模，進一步提升各類天氣條件下的光伏預測精度。

附錄見本刊網(wǎng)絡版（http：//www.aeps-info.com/aeps/ch/index.aspx），掃英文摘要后二維碼可以閱讀網(wǎng)絡全文。