多元時間序列缺失數據填補方法

李正欣, 張鳳鳴, 王 瑛, 陶 茜, 李 超

(空軍工程大學裝備管理與無人機工程學院, 陜西 西安 710051)

0 引 言

按時間順序獲取的一系列觀測值xt(j)稱為時間序列,其中,t(t=1,2,…,n)表示第t個時刻,j(j=1,2,…,m)表示第j個變量,xt(j)表示第j個變量在第t個時刻上的記錄值[1]。當m=1時,xt(j)為一元時間序列;當m>1時,xt(j)為多元時間序列。

時間序列廣泛存在于氣象、水文、醫學和經濟等眾多領域[2]。與一元時間序列相比,多元時間序列更具一般性,例如,股票交易記錄包含交易量、開盤價、收盤價、最高價和最低價[3];“黑匣子”記錄的飛行數據包含速度、高度、俯仰角、傾斜角和航向角等。廣義上,任何包含多變量的數據存儲都可以視為多元時間序列[4]。

時間序列中隱含著深層知識,針對時間序列的數據挖掘已成為目前數據挖掘領域中的研究熱點[5]。然而,在實際應用中,數據采集受多種因素的干擾,難免存在缺失數據,造成序列不連續,從而影響后續的分析處理。在數據挖掘領域UCI公開數據集中,超過40%的數據庫都含有缺失數據[6]。因此,研究多元時間序列的缺失數據填補方法是分析處理時間序列的前提、具有重要的應用價值。

1 相關研究

傳統的數據填補方法有均值填補、回歸填補、期望值最大化方法等[7],這些方法計算簡單,但當缺失數據量較多、波動較大時,填補精度大幅度下降。

隨后在傳統方法基礎上,又融入了統計分析、機器學習等方法,在一定程度上提高了填補精度。然而,這些成果主要針對傳統數據類型,多元時間序列在時間維和變量維上,都具有較強的相關性,與傳統數據具有較大的差異,因此,這些方法并不完全適用于多元時間序列。

目前,針對多元時間序列填補方法的研究成果相對較少[7],主要有多元回歸、神經網絡、支持向量機(support vector machine,SVM)和組合填補方法等。

多元自回歸模型[8]是一種經典的時間序列分析方法,它適用于線性模型,但對非線性的時間序列效果不佳。

神經網絡[9]可以對任意非線性函數進行逼近,具有非線性、自學習的特點,可用于時間序列缺失數據填補。但基于經驗風險最小化原則,要求較大的樣本規模,對于小樣本情況,容易出現過擬合、泛化能力不足等問題。

SVM[10]追求結構風險最小化,專門針對有限樣本,得到現有信息下的全局最優解,且計算復雜度與維度無關。然而,大多數基于SVM的填補方法[11-12],僅關注缺失數據所屬的單變量信息,沒有針對多元時間序列的特點、有效利用其他相關的多變量信息。文獻[13]提出了一種基于狀態匹配與SVM的多變量填補方法。首先,判斷存在缺失數據的序列所屬的類;然后,分析該類序列中,同缺失數據所屬變量相關的其他變量;最后,在相關變量中選擇樣本點,使用SVM填補缺失數據。然而,該方法在判斷序列所屬類別及分析相關變量時,使用人工判讀,對經驗依賴性較強。

文獻[14]提出基于最小二乘支持向量機(least squares support vector machine,LSSVM)的組合權重填補方法,同時利用缺失數據所屬的單變量以及與其相關的多變量信息,對單變量填補值與多變量填補值加權求和,作為最終填補結果。然而,該方法用一個權重值處理所有的缺失數據,難以體現各個缺失數據的差異性;此外,權重值在真實數據已知的情況下才能得到最優結果,限制了該方法的使用。

時間序列通常具有非平穩、非線性、非結構化等特征,此外,具有典型模式的序列通常樣本規模較小。針對以上特征,提出一種基于LSSVM的組合閾值填補方法。首先,以存在缺失數據的多元時間序列為對象,搜索與其同類的相似序列;然后,利用LSSVM分別進行多變量填補和單變量填補;最后,根據多變量和單變量填補結果的差異度,進行組合閾值填補。

2 基于LSSVM的組合閾值填補

LSSVM把最小二乘線性系統引入到SVM中,將解二次規劃問題轉化為求解線性方程組,能有效提高求解速度與泛化能力[15]。

2.1 LSSVM

訓練樣本集D={(xk,yk)}中:xk∈Rn表示輸入、yk∈R表示輸出,原始空間中的函數估計轉化為求解[14]

s.t.yk=ωTφ(xk)+b+ek,k=1,2,…,N

(1)

式中,γ為正則化參數；φ(·):Rn→Rnb是核空間映射函數；權矢量ω∈Rnb；誤差變量ek∈R；b是偏差量。

可以構造拉格朗日方程求解最優化問題為

L(ω,b,e,α)=J(ω,e)-

(2)

拉格朗日乘子αk∈R。令L對ω,b,ek,αk的偏導數分別等于0,消去ek和ω可得

(3)

式中,Ω=K(xk,xl)=φ(xk)Tφ(xl)；y=[y1,y2,…,yN]T;α=[α1,α2,…,αN]T;1v=[1,1,…,1]T。

根據Mercer條件,存在映射函數φ和核函數K(·,·)使得

K(xk,xl)=φ(xk)Tφ(xl)

(4)

則LSSVM的函數估計為

(5)

α、b由式(3)求解得到。

2.2 組合閾值填補方法

利用LSSVM填補缺失數據有兩種基本思路:單變量填補和多變量填補。單變量填補只利用缺失數據所屬的單變量信息,進行支持向量回歸,填補缺失數據。多變量填補利用與缺失數據所屬變量相關的其他多個變量信息,對缺失數據進行填補。

單變量填補依據變量自身的內在規律,多變量填補利用變量間的相互關系。組合閾值填補方法同時利用單變量和其他多個相關變量信息。

不失一般性,設數據集中,序列X1的第i個變量維度上存在缺失數據。組合閾值填補方法的基本流程如圖1所示。

圖1 組合閾值填補方法的基本流程Fig.1 Process of the threshold filling method

(1) 構建訓練集。在多元時間序列數據集中,搜索與X1相似的其他序列,結果記為X2,X3,…,Xn,構成訓練集。

(2) 多變量填補。在X2,X3,…,Xn的每個時刻上,以其他(m-1)個變量維度上的數據作為輸入、第i個變量維度上的數據作為輸出,訓練LSSVM。然后,針對X1的缺失數據段,利用訓練后的LSSVM模型進行填補,得到多變量填補結果。

(3) 單變量填補。在訓練集的每個序列中,針對第i個變量維度,以k個變量值作為輸入,第k+1個變量值作為輸出,訓練LSSVM。然后,針對X1的缺失數據段,利用訓練后的LSSVM模型進行填補,得到單變量填補結果。

(4) 組合閾值填補。以多變量和單變量的填補結果為基礎,按照

(6)

2.3 核函數與參數優化

LSSVM模型中,核函數把樣本投影到高維空間,把非線性問題轉化為線性問題。由于高斯徑向基核函數(radial basis function,RBF)具有良好性能,且只有一個形狀參數,因此,文中LSSVM采用RBF核函數。

LSSVM模型涉及2個參數:正則化參數γ和RBF核函數的形狀參數σ。γ表示擬合誤差權重,γ越大,擬合函數越逼近樣本點。σ影響擬合結果平滑程度,σ越大,擬合結果越平滑,但會出現擬合效果不佳的情況;σ過小,擬合結果會出現起伏。為了達到較好的填補效果,應合理選擇參數,保證擬合曲線比較平滑,同時具有較小的擬合誤差。

3 實驗與分析

實驗環境為:Intel(R) Core(TM) i7-3770 CPU、4G RAM、Windows 7和Matlab R2010a。

用已知分類結果的公開數據集ASL(Australian sign language)[16]作為實驗數據。ASL是一組手語信號,包含22個變量,左、右手動作各用11個變量描述:5個變量表示各手指的彎曲程度;6個變量(6個自由度)表示手的位置。數據集包含95種語意(95個類),每種語意包含27個序列,一共2 565個序列。

不失一般性,分別用前兩種語意alive、answer的第27個序列(記為alive_27、answer_27)作為實驗對象,序列長度分別為57、54。刪除第6個變量維度上的最后10個數據,模擬缺失數據,如圖2所示。比較多變量填補、單變量填補、組合權重填補[14]和組合閾值填補等4種方法的填補效果,以平均相對誤差e作為比較依據。

圖2 存在缺失數據的alive_27Fig.2 Sequence containing missing data (alive_27)

(7)

式中,n為缺失數據個數;si,ri分別表示第i個填補值和真實值。

3.1 填補效果的比較

為了客觀地比較填補效果,4種方法均以相似性搜索[17]的結果作為訓練集,本文重點是數據填補,因此直接用與alive_27、answer_27同類的其他26個序列作為相似性搜索結果、構成訓練集。

3.1.1 針對alive_27的填補

先對alive_27進行填補,多變量填補γ=2,σ=5時,e達到最小值3.6×10-4。針對不同的步長k,進行單變量填補,遍歷γ,σ各種組合,當k=10時,e達到最小值3.3×10-4,對應的γ=38,σ=67。

組合權重、組合閾值填補均以多變量、單變量填補結果為基礎,涉及參數分別為權重ω、閾值t。當ω=0.6時,組合權重填補e=3.16×10-4達到最小值。用式(6)計算每個時刻,多變量與單變量填補的差異度,如圖3所示。圖3(a)為針對alive_27的差異度,在10個填補數據中,時刻50、57上的差異度較大,閾值t=4×10-4時,組合閾值填補用多變量與單變量填補的均值表示兩個時刻的填補數據。4種方法的填補結果如圖4所示,誤差如表1所示。

圖3 多變量與單變量填補的差異度Fig.3 Difference between the results of multivariate and univariate filling methods

圖4 針對alive_27的填補結果Fig.4 Filling results for alive_27

3.1.2 針對answer _27的填補

類似地,對answer_27進行填補,其中,多變量填補γ=38 700、σ=120;單變量填補k=10、γ=130、σ=152;組合權重填補ω=0.6;根據圖3(b)針對answer_27的差異度,組合閾值填補t=1×10-3。填補結果如圖5所示。

圖5 針對answer_27的填補結果Fig.5 Filling results for answer_27

從alive_27、answer_27的填補結果可以看出,單變量填補體現了缺失數據的變化趨勢;多變量填補能夠體現波動情況。與組合權重方法相比,組合閾值方法確定閾值t,更加明確,且填補誤差較小。

3.2 缺失數據段長度對填補的影響

為了進一步比較4種方法的填補效果,增加缺失數據段的長度。分別刪除序列alive_27中第6個變量維度上的最后20、30、40個數據。填補結果如圖6所示,填補誤差如表2所示,對應的參數設置如表3所示。

隨著缺失數據段長度的增加,4種方法的填補效果同上文類似,但多變量填補在部分數據段上的填補偏差變大。原因在于,單變量填補方法在缺失數據所屬的變量維度上,以k個點為初始值進行填補,能夠較好地保持序列趨勢;但每次填補都以上次填補的結果作為輸入、存在累積誤差,難以體現波動細節。多變量填補方法利用其他多個相關變量信息,每次填補都以已知數據作為輸入,不存在累積誤差,能夠體現波動細節;但沒有利用缺失數據所屬變量信息,在某些點上存在較大偏差。

圖6 缺失數據段長度對填補效果的影響Fig.6 Relation between the length of missing data segment and its filling effect

表3 4種方法的參數設置

相比之下,組合權重、組合閾值填補的精度較高,受缺失數據段長度的影響不大。這是因為兩種方法利用單變量和其他相關多變量信息,較好地體現了波動細節,也有效避免了較大偏差。

然而,組合權重填補對單變量和多變量填補結果加權平均,用一個權重值處理所有的缺失數據,難以體現各個缺失數據的差異性;此外,權重值不易確定也限制了其使用。組合閾值填補以多變量填補結果為基礎,當多變量與單變量填補的差異度較大時,以兩者的均值作為填補結果,既體現了填補的差異性,也避免了權重值的選擇。

4 結 論

在小樣本情況下,針對多元時間序列中的缺失數據,利用缺失數據所屬單變量和其他相關多變量信息,提出了一種基于LSSVM的組合閾值填補方法,并進行了實驗驗證。結果表明,它能較好地體現序列的波動細節,且受缺失數據段長度的影響相對較小。

與組合權重填補方法相比,組合閾值填補方法的參數(差異度閾值)較容易確定,且能較好地體現填補的差異性,填補精度較高。

[2] 韓敏, 許美玲, 任偉杰. 多元混沌時間序列的相關狀態機預測模型研究[J]. 自動化學報, 2014, 40(5): 822-829.

HAN M, XU M L, REN W J. Research on multivariate chaotic time series prediction using mRSM model[J]. Acta Automatica Sinica, 2014, 40(5): 822-829.

[3] 吳虎勝, 張鳳鳴, 張超, 等. 多元時間序列的相似性匹配[J]. 應用科學學報, 2013, 31(6): 643-649.

WU H S, ZHANG F M, ZHANG C, et al. Matching similar patterns for multivariate time series[J]. Journal of applied sciences, 2013, 31(6): 643-649.

[4] PREE H, HERWIG B, GRUBER T, et al. On general purpose time series similarity measures and their use as kernel functions in support vector machines[J]. Information Sciences, 2014, 281(4): 478-495.

[5] LI H L,GUO C H.Piecewise cloud approximation for time series mining[J]. Knowledge-Based Systems, 2011, 24(4): 492-500.

[6] 馬茜,谷峪,李芳芳,等.順序敏感的多源感知數據填補技術[J].軟件學報,2016,27(9): 2332-2347.

MA Q, GU Y, LI F F, et al. Order-sensitive missing value imputation technology for multi-source sensory data[J]. Journal of Software, 2016, 27(9): 2332-2347.

[7] 呂政, 趙珺, 劉穎, 等. 基于最大方差權信息系數的煤氣數據填補[J]. 控制理論與應用, 2015, 32(5): 646-654.

Lü Z, ZHAO J, LIU Y, et al. Missing data imputation based on maximal variance weight information coefficient for gas flow in steel industry[J].Control Theory & Applications,2015,32(5):646-654.

[8] 雷春麗, 芮執元. 基于多元自回歸模型的電主軸熱誤差建模與預測[J]. 機械科學與技術, 2012, 31(9): 1526-1529.

LEI C L, RUI Z Y. Thermal error modeling and forecasting based on multivariate autoregressive model for motorized spindle[J].Mechanical Science and Technology for Aerospace Engineering, 2012, 31(9): 1526-1529.

[9] ZOUCAS F A M, BELFIORE P. An empirical analysis of a neural network model for the time series forecasting of different industrial segments[J]. International Journal of Applied Decision Sciences, 2015, 8(3): 261-283.

[10] VAPNIK V. The nature of statistical learning theory[M]. NewYork: Springer-Verlag, 1995.

[11] 許磊, 張鳳鳴. 缺失飛參數據填補的組合方法研究[J]. 計算機工程與應用, 2010, 46(21): 210-212.

XU L, ZHANG F M. Research on combined method in missing flight data complement[J]. Computer Engineering and Applications, 2010, 46(21): 210-212.

[12] CHEN T T, LEE S J. A weighted LS-SVM based learning system for time series forecasting[J]. Information Sciences, 2015, 299(C): 99-116.

[13] 謝川,倪世宏,張宗麟.一種缺失飛行參數預處理的新方法[J].計算機仿真,2005,22(4): 27-30.

XIE C, NI S H, ZHANG Z L. A new method for preprocessing in absent flight parameter[J].Computer Simulation,2005,22(4): 27-30.

[14] 楊軻, 張曉豐, 趙錄峰, 等. 基于LSSVM的缺失飛行數據組合填補方法[J]. 火力與指揮控制, 2013, 38(1): 84-90.

YANG K, ZHANG X F, ZHAO L F, et al. Combined method of complementing missing flight data based on LSSVM[J]. Fire Control & Command Control, 2013, 38(1): 84-90.

[15] LANGONE R, ALZATE C, KETELAERE B D, et al. LS-SVM based spectral clustering and regression for predicting maintenance of industrial machines[J]. Engineering Applications of Artificial Intelligence, 2015, 37(37): 268-278.

[16] MOHAMMED W K. High-quality recordings of australian sign language signs[EB/OL].[2017-01-20].http :∥kdd.ics.uci.edu/databases/ High-quality Australian Sign Language/ High-quality Australian Sign Language. html.

[17] KEOGH E, RATANAMAHATANAC. Exact indexing of dynamic time warping[J]. Knowledge and Information Systems, 2005, 7(3): 358-386.