變壓器運(yùn)行數(shù)據(jù)缺失值修復(fù)研究

秦佳峰，周 超，林 穎，白德盟，鄭文杰

（國網(wǎng)山東省電力公司電力科學(xué)研究院，山東 濟(jì)南 250003）

0 引言

隨著電力系統(tǒng)信息化程度的提高，電網(wǎng)的運(yùn)行和監(jiān)測數(shù)據(jù)均呈指數(shù)級增長［1-3］。巨量的數(shù)據(jù)為智能調(diào)度和決策提供了信息基礎(chǔ)，但也對數(shù)據(jù)處理手段提出了更高的要求［4-6］。變壓器作為電網(wǎng)中的重要一環(huán)，其運(yùn)行數(shù)據(jù)對后續(xù)的設(shè)備狀態(tài)評估至關(guān)重要，但是由于設(shè)備故障、通信故障、操作失誤等原因，實(shí)際取得的運(yùn)行數(shù)據(jù)往往存在缺失的情況，從而降低后續(xù)數(shù)據(jù)挖掘算法的性能。因此，修復(fù)缺失的變壓器運(yùn)行數(shù)據(jù)對提升智能電網(wǎng)的運(yùn)行效率，提高電網(wǎng)系統(tǒng)的可靠性和自愈能力至關(guān)重要。

關(guān)于電力系統(tǒng)缺失數(shù)據(jù)的修復(fù)，目前已有多種方法［7-9］。張曉星等［7］基于聚類的思想計(jì)算電力負(fù)荷曲線的特征函數(shù)，利用該特征函數(shù)修復(fù)缺失數(shù)據(jù)，但這種方法對數(shù)據(jù)的日曲線相似度要求很高。嚴(yán)英杰等［8］將輸變電設(shè)備狀態(tài)數(shù)據(jù)看成時(shí)間序列，用差分整合移動(dòng)平均自回歸模型（Autoregressive Integrated Moving Average，ARIMA）擬合并迭代檢驗(yàn)的方法修復(fù)缺失數(shù)據(jù)，但這種方法利用的信息較少，且不適合連續(xù)缺失的情況。劉沅昆等［9］提出了一種基于Pearson 相關(guān)系數(shù)的歷史數(shù)據(jù)挖掘恢復(fù)方法，但這種方法本質(zhì)上是在做線性回歸，很難刻畫出數(shù)據(jù)的實(shí)際波動(dòng)情況。

通過分析真實(shí)的變壓器運(yùn)行數(shù)據(jù)，希望找到一種缺失數(shù)據(jù)修復(fù)方法，不僅能利用不同日曲線之間的相似性，并且適用不同的缺失點(diǎn)分布情況。函數(shù)型主成分分析（Functional Principal Component Analysis，F(xiàn)PCA）作為一種高維特征提取方法，能夠通過對同一觀測對象的重復(fù)測量，在低維函數(shù)空間上表示出數(shù)據(jù)最主要的波動(dòng)情況［10-11］，能夠很好地?cái)M合出數(shù)據(jù)波動(dòng)的整體趨勢。小波變換作為信號處理的常用方法，能夠通過對函數(shù)的多尺度細(xì)分，聚焦函數(shù)的局部細(xì)節(jié)，提取出其中的有效信息［12］。

因此，結(jié)合FPCA和小波變換處理函數(shù)型數(shù)據(jù)的優(yōu)勢，提出了一種變壓器缺失數(shù)據(jù)修復(fù)方法。該方法能夠自動(dòng)提取數(shù)據(jù)波動(dòng)的特征，并使恢復(fù)所得的數(shù)據(jù)在局部上更符合實(shí)際情形［13］。

1 FPCA原理

1.1 根據(jù)FPCA構(gòu)建模型

具體到變壓器的運(yùn)行數(shù)據(jù)修復(fù)，考慮其中一類運(yùn)行數(shù)據(jù)，可以將第i天內(nèi)變壓器的該類運(yùn)行數(shù)據(jù)視為一個(gè)觀測對象Yi，傳感器在第j個(gè)觀測點(diǎn)對應(yīng)的時(shí)間點(diǎn)tj測得的數(shù)值Yij為該時(shí)刻的觀測值，記第i天的觀測點(diǎn)數(shù)量為ni。由于同一個(gè)變壓器在連續(xù)的幾天內(nèi)運(yùn)行數(shù)據(jù)的波動(dòng)趨勢有一定的相似性，因此可以將N天的數(shù)據(jù)Yij，1 ≤i≤N，1 ≤j≤ni視為一組縱向數(shù)據(jù)進(jìn)行分析。首先，利用FPCA 方法估計(jì)函數(shù)型數(shù)據(jù)的低頻部分Xi(t)。根據(jù)Karhunen-Loève 表示［14］，Xi(t)可表示為

式中：μ(t)為均值函數(shù)；φk(t)為協(xié)方差曲面的特征函數(shù)，即主成分；αik為Xi(t)在主成分上的得分。

利用FPCA方法分別估計(jì)均值函數(shù)μ(t)、特征函數(shù)φk(t)、系數(shù)αik。由于通過FPCA方法估計(jì)的Xi(t)忽略了數(shù)據(jù)高頻部分的信息，利用小波刻畫函數(shù)數(shù)據(jù)細(xì)節(jié)的能力，通過對殘差函數(shù)ξi(t)做小波變換，提取出其中高頻部分的信息。因此整體模型可以改寫為

關(guān)于殘差函數(shù)ξi(t)的處理將在下文敘述。同時(shí)，將第i天的第j個(gè)測量時(shí)間點(diǎn)tj記為Tij，將觀測值Yi(tj)記為Yij。

1.2 均值函數(shù)及協(xié)方差曲面估計(jì)

采用加權(quán)最小二乘的方法估計(jì)均值函數(shù)μ(t)。對某一固定的時(shí)間點(diǎn)t，假設(shè)β0(t)和β1(t)為兩個(gè)參數(shù)，其中β0(t)為均值函數(shù)μ(t)的估計(jì)，β1(t)為線性修正函數(shù)的斜率，因此Yij-β1(t-Tij)即為t時(shí)刻Yij對應(yīng)的修正值。取高斯核函數(shù)為權(quán)重值，因此，可通過求解如下極小化問題得到對均值函數(shù)μ(t)的估計(jì)為

式中：β0、β1為關(guān)于時(shí)間t的參數(shù)函數(shù)；hμ是帶寬，為了使算法有更強(qiáng)的自適應(yīng)性，采用廣義交叉檢驗(yàn)的方法選擇帶寬；k(·)是高斯核函數(shù)，其表達(dá)式為

通過求解式（3）中的極小化問題，得到對均值函數(shù)μ(t)的估計(jì)為

對協(xié)方差曲面的估計(jì)同樣可以采用加權(quán)最小二乘的方法。與估計(jì)均值函數(shù)μ(t)一樣，對某一組固定的時(shí)間點(diǎn)(s，t)，β0(s，t)，β11(s，t)和β12(s，t)都是待估計(jì)的參數(shù)，其中β0(s，t)為對協(xié)方差曲面G(s，t)的估計(jì)，β11(s，t)和β12(s，t)分別為兩個(gè)方向上線性修正函數(shù)的斜率。

首先估計(jì)在s≠t時(shí)的協(xié)方差曲面，可以通過求解如下極小化問題得到，即為

式中：hG為帶寬，可通過廣義交叉檢驗(yàn)的方法選擇帶寬；k2(·，·)為一個(gè)二元高斯核函數(shù)。k2(·，·)和Gi(Tij，Til)的表達(dá)式分別為：

通過求解式（6）中的極小化問題，得到在s≠t時(shí)協(xié)方差曲面的估計(jì)為

當(dāng)s=t時(shí)，由于協(xié)方差曲面在垂直對角線方向上的形狀更接近于二次曲線［15］，在垂直對角線方向上用二次函數(shù)修正，而在對角線方向上仍采用線性函數(shù)修正，將坐標(biāo)軸順時(shí)針旋轉(zhuǎn)45°，即為

結(jié)合式（9）和式（11），可得對協(xié)方差曲面的估計(jì)(s，t)。

1.3 主成分的選取及主成分得分的估計(jì)

主成分是協(xié)方差曲面的特征函數(shù)，因此可以通過求解式（13）的積分方程得到對特征函數(shù)的估計(jì)。

具體的，通過將協(xié)方差曲面離散化，得到K組特征值和特征函數(shù)1 ≤k≤K，需要選取其中最重要的k個(gè)主成分。為此，定義方差解釋比（Fraction of Variance Explained，F(xiàn)VE）為

給定一個(gè)閾值θ，選擇前k個(gè)最大的特征值和對應(yīng)的特征函數(shù)，使得AFVE(k)≥θ，完成對主成分的選取。

最后用條件期望來估計(jì)主成分得分αik，記

其中δjl在j=l時(shí)為1，j≠l時(shí)為0，Γ為t的定義域?yàn)?/p>

該函數(shù)能夠較好地近似函數(shù)型數(shù)據(jù)的整體波動(dòng)趨勢，但是在局部細(xì)節(jié)方面存在過于平滑的現(xiàn)象，因此需要分析真實(shí)值與估計(jì)函數(shù)做差得到的殘差函數(shù)ξi(t)，通過用小波變換提取其中的高頻有效信息，得到更為精確的估計(jì)。

2 小波變換原理

2.1 根據(jù)小波變換構(gòu)建模型

假設(shè)殘差函數(shù)ξi(t)是一個(gè)包含噪聲的一維信號，即ξi(t)可以表示為

式中：?i(t)為真實(shí)信號；ei(t)為噪聲。將ξi(t)看成是一個(gè)閉區(qū)間上的有界連續(xù)函數(shù)，因此可進(jìn)一步假設(shè)ξi(t)是平方可積的，根據(jù)離散小波變換的理論［16］，ξi(t)展開為

式中：Φ(t)和ψ(t)分別為尺度函數(shù)和小波基函數(shù)；cj0(k)和dj(k)分別為對應(yīng)尺度下的系數(shù)。

由于實(shí)際應(yīng)用中，得到的ξi(t)是一列離散的采樣點(diǎn)，根據(jù)采樣定理，尺度無法無限細(xì)分下去，因此具體計(jì)算中將ξi(t)展開為

式中：j0=0；j1為分解層數(shù)。

2.2 去噪方案的選取

由于小波基、分解層數(shù)和閾值方案的選擇會對去噪效果產(chǎn)生影響，對不同類型的數(shù)據(jù)可能需要選擇不同的參數(shù)進(jìn)行去噪，隨機(jī)選取真實(shí)數(shù)據(jù)點(diǎn)中的10%作為測試集，根據(jù)不同參數(shù)組合在測試集上的表現(xiàn)選擇最優(yōu)的去噪方案。

可供 選擇的 小波基［17-19］包 括db1-8、sym2-6、coif1-4，分解層數(shù)為1～8，閾值方案包括閾值選取規(guī)則和小于閾值小波系數(shù)的處理方法兩個(gè)部分。

2.2.1 閾值選取規(guī)則

閾值選取規(guī)則主要包括四種方式：

1）無偏風(fēng)險(xiǎn)估計(jì)閾值（rigrsure），算法為：

b）若Aeta＜Acrit，則λheur=λsqt；若Aeta≥Acrit，則λheur=min{λrigr，λsqt}。

4）極大極小閾值（minimax）為

2.2.2 小于閾值小波系數(shù)的處理方法

小于閾值小波系數(shù)的處理方法主要包括兩種：

1）硬閾值處理為

2）軟閾值處理為

3 變壓器運(yùn)行數(shù)據(jù)缺失點(diǎn)修復(fù)方案

變壓器運(yùn)行數(shù)據(jù)缺失點(diǎn)修復(fù)方案主要包括三個(gè)部分，首先用FPCA方法得到觀測對象的初次估計(jì)函數(shù)；然后估計(jì)殘差函數(shù)，并用交叉驗(yàn)證的方法找到最優(yōu)的小波去噪方案，得到對殘差函數(shù)的估計(jì)；最后將殘差函數(shù)的估計(jì)作為修正函數(shù)，與FPCA初次估計(jì)函數(shù)結(jié)合，得到觀測對象在整個(gè)定義域中的估計(jì)值，用該估計(jì)值完成缺失點(diǎn)的修復(fù)。

以變壓器第i天內(nèi)某一類運(yùn)行數(shù)據(jù)作為一個(gè)觀測對象Yi，考慮連續(xù)n天內(nèi)同一變電站的同一類數(shù)據(jù)，其中1 ≤d≤n。將第d天的觀測對象Yd的真實(shí)觀測值集合記為Sd={Yd1，Yd2，···，Ydnd}，nd為連續(xù)觀測n天中的第d天，計(jì)算第d天的估計(jì)函數(shù)，包括9個(gè)步驟。

步驟1）隨機(jī)選取其中10%的點(diǎn)作為測試集，記Ud={Ydt1，Ydt2，···，Ydtp}，p=[0.1×nd]，其 余90% 的點(diǎn)作為訓(xùn)練集，記為Vd={Ydtp+1，Ydtp+2，···，Ydtnd}。以第d天的訓(xùn)練集Vd，以及除第d天外所有的真實(shí)測量值作為訓(xùn)練集，用FPCA 方法建模，其中t的單位為分鐘，考慮到實(shí)際情況，所有可能的觀測點(diǎn)都是整數(shù)分鐘。

步驟3）由式（6），得到s≠t時(shí)的協(xié)方差曲面的估計(jì)函數(shù)，再由式（11），得到s=t時(shí)協(xié)方差曲面的估計(jì)函數(shù)，從而得到整個(gè)協(xié)方差曲面的估計(jì)函數(shù)(s，t)。

步驟4）根據(jù)式（13），通過將協(xié)方差曲面離散化的方法，得到K組特征函數(shù)和對應(yīng)的特征值的估計(jì)選取前k個(gè)最大的特征值及其對應(yīng)的特征函數(shù)，使得式（14）定義的AFVE(k)大于給定的閾值θ。得到k個(gè)主成分的估計(jì)

步驟5）根據(jù)式（16），估計(jì)出觀測對象Yd在主成分(t)上的得分。綜合步驟1）—步驟5），得到對觀測對象Yd的初次估計(jì)函數(shù)，即低頻部分的估計(jì)為

步驟6）估計(jì)殘差函數(shù)(t)。為方便敘述，將前述Sd中的觀測時(shí)間點(diǎn)記為，同理，記，tdnd}。當(dāng)t∈時(shí)，則殘差函數(shù)為當(dāng)t?時(shí)，補(bǔ)充定義td0=1，td(nd+1)=1440，則必存在p，使得t∈(tdp，td(p+1))。定義則當(dāng)Agap較小時(shí)，可用移動(dòng)平均來估計(jì)該點(diǎn)的殘差函數(shù)的值，當(dāng)Agap較大時(shí)，對殘差函數(shù)的估計(jì)意義不大，反而可能造成最終估計(jì)的錯(cuò)誤，因此直接置為零。綜上，對殘差函數(shù)(t)的估計(jì)為

步驟7）利用MATLAB軟件中的小波工具箱對殘差函數(shù)進(jìn)行小波變換，并選擇最優(yōu)的去噪方案。

首先選定一個(gè)去噪方案，用wden 函數(shù)進(jìn)行小波去噪。該函數(shù)的四個(gè)參數(shù)“wname”、“n”、“tptr”、“sorh”分別對應(yīng)于小波基、分解層數(shù)、閾值選取規(guī)則以及軟硬閾值。最后一個(gè)參數(shù)“scal”表示定義的閾值是否需要重新調(diào)整，默認(rèn)為“one”，即不用重新調(diào)整。輸入去噪方案對應(yīng)的參數(shù)以及信號(t)之后，就可得到去噪后的殘差函數(shù)(t)。

選取不同的去噪方案，重復(fù)步驟7），記錄最小的RRMSEi及對應(yīng)的去噪方案。

步驟8）將訓(xùn)練集擴(kuò)充為整個(gè)觀測集，重復(fù)步驟1）—步驟5），更新初次估計(jì)函數(shù)(t)，再更新對殘差函數(shù)(t)的估計(jì)。即用重復(fù)步驟6），得到更新后的殘差函數(shù)為

這樣就完成了對觀測對象Yd的修復(fù)。

4 算例分析

以某220 kV 主變壓器某月的高壓側(cè)-I 數(shù)據(jù)為例，對本文所提缺失點(diǎn)修復(fù)方法進(jìn)行測試?？紤]到實(shí)際應(yīng)用中既存在零散的隨機(jī)缺失點(diǎn)，又存在連續(xù)缺失的情況，因此分別對這兩種情況做了測試。

首先考慮缺失點(diǎn)零散分布的情況。以該月前十天的數(shù)據(jù)為例進(jìn)行分析，1 日實(shí)際獲取了997 個(gè)觀測點(diǎn)，隨機(jī)選取其中100 個(gè)點(diǎn)作為測試集，假設(shè)這些觀測點(diǎn)缺失，以此為參照，評價(jià)修復(fù)效果。

經(jīng)過步驟1）—步驟5）后，得到的初次估計(jì)函數(shù)如圖1 所示。真實(shí)測量值及其變化趨勢用黑點(diǎn)及黑色虛線連接，假設(shè)缺失的100 個(gè)點(diǎn)用藍(lán)色星狀點(diǎn)表示，初次估計(jì)函數(shù)用紅色實(shí)線表示。

圖1 初次估計(jì)函數(shù)

可以看到，初次估計(jì)函數(shù)已經(jīng)能夠較好地?cái)M合真實(shí)觀測值，但是在局部細(xì)節(jié)上還有欠缺。在步驟1）—步驟5）的基礎(chǔ)上，繼續(xù)實(shí)施步驟6）—9）后，得到的最終估計(jì)函數(shù)如圖2 所示。黑點(diǎn)和藍(lán)色星狀點(diǎn)的含義同圖1，紅色實(shí)線表示最終估計(jì)函數(shù)。

圖2 最終估計(jì)函數(shù)

結(jié)果表明，在對初次估計(jì)函數(shù)加上修正函數(shù)之后，最終估計(jì)函數(shù)能夠很好地?cái)M合真實(shí)測量值。以測試集上的RRMSEd作為指標(biāo)，將這種方法同三次樣條插值方法進(jìn)行對比。重復(fù)進(jìn)行了10 次實(shí)驗(yàn)，每次的測試集都是隨機(jī)選取，對比兩種方法的均方根誤差（Root Mean Squard Error，RMSE）如表1所示。

表1 高壓側(cè)-I上的RMSE對比

結(jié)果顯示，經(jīng)本文方法得出的估計(jì)函數(shù)在測試集上的RMSE 要明顯小于三次樣條插值法的RMSE，10次實(shí)驗(yàn)平均能夠使RMSE下降21.3%，本文的缺失點(diǎn)修復(fù)方法更加穩(wěn)定可靠。

當(dāng)缺失點(diǎn)為大塊的連續(xù)缺失時(shí)，普通的插值方法已經(jīng)沒法給出一個(gè)合理的估計(jì)，反而會添加錯(cuò)誤數(shù)據(jù)。因此仍然以該月1 日至10 日的數(shù)據(jù)作為一組進(jìn)行分析，假設(shè)第700～799 個(gè)觀測值缺失，利用所述的方法，對該段的缺失值進(jìn)行恢復(fù)?；謴?fù)的結(jié)果如圖3所示。

圖3 連續(xù)缺失時(shí)的恢復(fù)曲線

為了檢驗(yàn)在所考慮數(shù)據(jù)集上的整體效果，測試該方法應(yīng)用于其他11 類數(shù)據(jù)上的效果，同高壓側(cè)-I的做法一樣，隨機(jī)選取100 個(gè)真實(shí)觀測值作為測試集，重復(fù)實(shí)驗(yàn)10 次，分別計(jì)算本文方法和三次樣條方法在測試集上的平均RMSE，結(jié)果如表2 所示。測試結(jié)果顯示，本文方法在不同類型的變壓器運(yùn)行數(shù)據(jù)恢復(fù)上均有較好的表現(xiàn)。對三次樣條差值方法有很大程度的提高。除高壓側(cè)-Q與中壓側(cè)-Q外，RMSE下降的百分比均超過15%。

表2 12類數(shù)據(jù)上的RMSE對比

5 結(jié)語

利用FPCA 對數(shù)字信號整體特征的提取能力和小波變換對局部信息的提取能力，對變壓器運(yùn)行數(shù)據(jù)的缺失數(shù)據(jù)修復(fù)這一實(shí)際問題，提出了基于FPCA和小波變換的變壓器運(yùn)行數(shù)據(jù)缺失數(shù)據(jù)修復(fù)方法。

使用FPCA 方法能夠?qū)崿F(xiàn)對需修復(fù)數(shù)據(jù)集的初次估計(jì)函數(shù)；并利用交叉驗(yàn)證方法找到最優(yōu)的小波去噪方案，得到對殘差函數(shù)的估計(jì)；通過將殘差函數(shù)的估計(jì)作為修正函數(shù)，與FPCA 初次估計(jì)函數(shù)結(jié)合，完成缺失點(diǎn)的修復(fù)。

該方法能夠根據(jù)少量歷史數(shù)據(jù)自動(dòng)尋找數(shù)據(jù)的波動(dòng)特征，并在此基礎(chǔ)上對局部細(xì)節(jié)進(jìn)行優(yōu)化，從而完成對缺失數(shù)據(jù)的恢復(fù)。該方法對缺失點(diǎn)是連續(xù)分布還是離散分布沒有要求，具有適用性廣、恢復(fù)精度高的特點(diǎn)。