李正欣, 張鳳鳴, 王 瑛, 陶 茜, 李 超
(空軍工程大學裝備管理與無人機工程學院, 陜西 西安 710051)
按時間順序獲取的一系列觀測值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)為多元時間序列。
時間序列廣泛存在于氣象、水文、醫(yī)學和經(jīng)濟等眾多領(lǐng)域[2]。與一元時間序列相比,多元時間序列更具一般性,例如,股票交易記錄包含交易量、開盤價、收盤價、最高價和最低價[3];“黑匣子”記錄的飛行數(shù)據(jù)包含速度、高度、俯仰角、傾斜角和航向角等。廣義上,任何包含多變量的數(shù)據(jù)存儲都可以視為多元時間序列[4]。
時間序列中隱含著深層知識,針對時間序列的數(shù)據(jù)挖掘已成為目前數(shù)據(jù)挖掘領(lǐng)域中的研究熱點[5]。然而,在實際應用中,數(shù)據(jù)采集受多種因素的干擾,難免存在缺失數(shù)據(jù),造成序列不連續(xù),從而影響后續(xù)的分析處理。在數(shù)據(jù)挖掘領(lǐng)域UCI公開數(shù)據(jù)集中,超過40%的數(shù)據(jù)庫都含有缺失數(shù)據(jù)[6]。因此,研究多元時間序列的缺失數(shù)據(jù)填補方法是分析處理時間序列的前提、具有重要的應用價值。
傳統(tǒng)的數(shù)據(jù)填補方法有均值填補、回歸填補、期望值最大化方法等[7],這些方法計算簡單,但當缺失數(shù)據(jù)量較多、波動較大時,填補精度大幅度下降。
隨后在傳統(tǒng)方法基礎(chǔ)上,又融入了統(tǒng)計分析、機器學習等方法,在一定程度上提高了填補精度。然而,這些成果主要針對傳統(tǒng)數(shù)據(jù)類型,多元時間序列在時間維和變量維上,都具有較強的相關(guān)性,與傳統(tǒng)數(shù)據(jù)具有較大的差異,因此,這些方法并不完全適用于多元時間序列。
目前,針對多元時間序列填補方法的研究成果相對較少[7],主要有多元回歸、神經(jīng)網(wǎng)絡(luò)、支持向量機(support vector machine,SVM)和組合填補方法等。
多元自回歸模型[8]是一種經(jīng)典的時間序列分析方法,它適用于線性模型,但對非線性的時間序列效果不佳。
神經(jīng)網(wǎng)絡(luò)[9]可以對任意非線性函數(shù)進行逼近,具有非線性、自學習的特點,可用于時間序列缺失數(shù)據(jù)填補。但基于經(jīng)驗風險最小化原則,要求較大的樣本規(guī)模,對于小樣本情況,容易出現(xiàn)過擬合、泛化能力不足等問題。
SVM[10]追求結(jié)構(gòu)風險最小化,專門針對有限樣本,得到現(xiàn)有信息下的全局最優(yōu)解,且計算復雜度與維度無關(guān)。然而,大多數(shù)基于SVM的填補方法[11-12],僅關(guān)注缺失數(shù)據(jù)所屬的單變量信息,沒有針對多元時間序列的特點、有效利用其他相關(guān)的多變量信息。文獻[13]提出了一種基于狀態(tài)匹配與SVM的多變量填補方法。首先,判斷存在缺失數(shù)據(jù)的序列所屬的類;然后,分析該類序列中,同缺失數(shù)據(jù)所屬變量相關(guān)的其他變量;最后,在相關(guān)變量中選擇樣本點,使用SVM填補缺失數(shù)據(jù)。然而,該方法在判斷序列所屬類別及分析相關(guān)變量時,使用人工判讀,對經(jīng)驗依賴性較強。
文獻[14]提出基于最小二乘支持向量機(least squares support vector machine,LSSVM)的組合權(quán)重填補方法,同時利用缺失數(shù)據(jù)所屬的單變量以及與其相關(guān)的多變量信息,對單變量填補值與多變量填補值加權(quán)求和,作為最終填補結(jié)果。然而,該方法用一個權(quán)重值處理所有的缺失數(shù)據(jù),難以體現(xiàn)各個缺失數(shù)據(jù)的差異性;此外,權(quán)重值在真實數(shù)據(jù)已知的情況下才能得到最優(yōu)結(jié)果,限制了該方法的使用。
時間序列通常具有非平穩(wěn)、非線性、非結(jié)構(gòu)化等特征,此外,具有典型模式的序列通常樣本規(guī)模較小。針對以上特征,提出一種基于LSSVM的組合閾值填補方法。首先,以存在缺失數(shù)據(jù)的多元時間序列為對象,搜索與其同類的相似序列;然后,利用LSSVM分別進行多變量填補和單變量填補;最后,根據(jù)多變量和單變量填補結(jié)果的差異度,進行組合閾值填補。
LSSVM把最小二乘線性系統(tǒng)引入到SVM中,將解二次規(guī)劃問題轉(zhuǎn)化為求解線性方程組,能有效提高求解速度與泛化能力[15]。
訓練樣本集D={(xk,yk)}中:xk∈Rn表示輸入、yk∈R表示輸出,原始空間中的函數(shù)估計轉(zhuǎn)化為求解[14]
s.t.yk=ωTφ(xk)+b+ek,k=1,2,…,N
(1)
式中,γ為正則化參數(shù);φ(·):Rn→Rnb是核空間映射函數(shù);權(quán)矢量ω∈Rnb;誤差變量ek∈R;b是偏差量。
可以構(gòu)造拉格朗日方程求解最優(yōu)化問題為
L(ω,b,e,α)=J(ω,e)-
(2)
拉格朗日乘子αk∈R。令L對ω,b,ek,αk的偏導數(shù)分別等于0,消去ek和ω可得
(3)
式中,Ω=K(xk,xl)=φ(xk)Tφ(xl);y=[y1,y2,…,yN]T;α=[α1,α2,…,αN]T;1v=[1,1,…,1]T。
根據(jù)Mercer條件,存在映射函數(shù)φ和核函數(shù)K(·,·)使得
K(xk,xl)=φ(xk)Tφ(xl)
(4)
則LSSVM的函數(shù)估計為
(5)
α、b由式(3)求解得到。
利用LSSVM填補缺失數(shù)據(jù)有兩種基本思路:單變量填補和多變量填補。單變量填補只利用缺失數(shù)據(jù)所屬的單變量信息,進行支持向量回歸,填補缺失數(shù)據(jù)。多變量填補利用與缺失數(shù)據(jù)所屬變量相關(guān)的其他多個變量信息,對缺失數(shù)據(jù)進行填補。
單變量填補依據(jù)變量自身的內(nèi)在規(guī)律,多變量填補利用變量間的相互關(guān)系。組合閾值填補方法同時利用單變量和其他多個相關(guān)變量信息。
不失一般性,設(shè)數(shù)據(jù)集中,序列X1的第i個變量維度上存在缺失數(shù)據(jù)。組合閾值填補方法的基本流程如圖1所示。
圖1 組合閾值填補方法的基本流程Fig.1 Process of the threshold filling method
(1) 構(gòu)建訓練集。在多元時間序列數(shù)據(jù)集中,搜索與X1相似的其他序列,結(jié)果記為X2,X3,…,Xn,構(gòu)成訓練集。
(2) 多變量填補。在X2,X3,…,Xn的每個時刻上,以其他(m-1)個變量維度上的數(shù)據(jù)作為輸入、第i個變量維度上的數(shù)據(jù)作為輸出,訓練LSSVM。然后,針對X1的缺失數(shù)據(jù)段,利用訓練后的LSSVM模型進行填補,得到多變量填補結(jié)果。
(3) 單變量填補。在訓練集的每個序列中,針對第i個變量維度,以k個變量值作為輸入,第k+1個變量值作為輸出,訓練LSSVM。然后,針對X1的缺失數(shù)據(jù)段,利用訓練后的LSSVM模型進行填補,得到單變量填補結(jié)果。
(4) 組合閾值填補。以多變量和單變量的填補結(jié)果為基礎(chǔ),按照
(6)
LSSVM模型中,核函數(shù)把樣本投影到高維空間,把非線性問題轉(zhuǎn)化為線性問題。由于高斯徑向基核函數(shù)(radial basis function,RBF)具有良好性能,且只有一個形狀參數(shù),因此,文中LSSVM采用RBF核函數(shù)。
LSSVM模型涉及2個參數(shù):正則化參數(shù)γ和RBF核函數(shù)的形狀參數(shù)σ。γ表示擬合誤差權(quán)重,γ越大,擬合函數(shù)越逼近樣本點。σ影響擬合結(jié)果平滑程度,σ越大,擬合結(jié)果越平滑,但會出現(xiàn)擬合效果不佳的情況;σ過小,擬合結(jié)果會出現(xiàn)起伏。為了達到較好的填補效果,應合理選擇參數(shù),保證擬合曲線比較平滑,同時具有較小的擬合誤差。
實驗環(huán)境為:Intel(R) Core(TM) i7-3770 CPU、4G RAM、Windows 7和Matlab R2010a。
用已知分類結(jié)果的公開數(shù)據(jù)集ASL(Australian sign language)[16]作為實驗數(shù)據(jù)。ASL是一組手語信號,包含22個變量,左、右手動作各用11個變量描述:5個變量表示各手指的彎曲程度;6個變量(6個自由度)表示手的位置。數(shù)據(jù)集包含95種語意(95個類),每種語意包含27個序列,一共2 565個序列。
不失一般性,分別用前兩種語意alive、answer的第27個序列(記為alive_27、answer_27)作為實驗對象,序列長度分別為57、54。刪除第6個變量維度上的最后10個數(shù)據(jù),模擬缺失數(shù)據(jù),如圖2所示。比較多變量填補、單變量填補、組合權(quán)重填補[14]和組合閾值填補等4種方法的填補效果,以平均相對誤差e作為比較依據(jù)。
圖2 存在缺失數(shù)據(jù)的alive_27Fig.2 Sequence containing missing data (alive_27)
(7)
式中,n為缺失數(shù)據(jù)個數(shù);si,ri分別表示第i個填補值和真實值。
為了客觀地比較填補效果,4種方法均以相似性搜索[17]的結(jié)果作為訓練集,本文重點是數(shù)據(jù)填補,因此直接用與alive_27、answer_27同類的其他26個序列作為相似性搜索結(jié)果、構(gòu)成訓練集。
3.1.1 針對alive_27的填補
先對alive_27進行填補,多變量填補γ=2,σ=5時,e達到最小值3.6×10-4。針對不同的步長k,進行單變量填補,遍歷γ,σ各種組合,當k=10時,e達到最小值3.3×10-4,對應的γ=38,σ=67。
組合權(quán)重、組合閾值填補均以多變量、單變量填補結(jié)果為基礎(chǔ),涉及參數(shù)分別為權(quán)重ω、閾值t。當ω=0.6時,組合權(quán)重填補e=3.16×10-4達到最小值。用式(6)計算每個時刻,多變量與單變量填補的差異度,如圖3所示。圖3(a)為針對alive_27的差異度,在10個填補數(shù)據(jù)中,時刻50、57上的差異度較大,閾值t=4×10-4時,組合閾值填補用多變量與單變量填補的均值表示兩個時刻的填補數(shù)據(jù)。4種方法的填補結(jié)果如圖4所示,誤差如表1所示。
圖3 多變量與單變量填補的差異度Fig.3 Difference between the results of multivariate and univariate filling methods
圖4 針對alive_27的填補結(jié)果Fig.4 Filling results for alive_27
3.1.2 針對answer _27的填補
類似地,對answer_27進行填補,其中,多變量填補γ=38 700、σ=120;單變量填補k=10、γ=130、σ=152;組合權(quán)重填補ω=0.6;根據(jù)圖3(b)針對answer_27的差異度,組合閾值填補t=1×10-3。填補結(jié)果如圖5所示。
圖5 針對answer_27的填補結(jié)果Fig.5 Filling results for answer_27
從alive_27、answer_27的填補結(jié)果可以看出,單變量填補體現(xiàn)了缺失數(shù)據(jù)的變化趨勢;多變量填補能夠體現(xiàn)波動情況。與組合權(quán)重方法相比,組合閾值方法確定閾值t,更加明確,且填補誤差較小。
為了進一步比較4種方法的填補效果,增加缺失數(shù)據(jù)段的長度。分別刪除序列alive_27中第6個變量維度上的最后20、30、40個數(shù)據(jù)。填補結(jié)果如圖6所示,填補誤差如表2所示,對應的參數(shù)設(shè)置如表3所示。
隨著缺失數(shù)據(jù)段長度的增加,4種方法的填補效果同上文類似,但多變量填補在部分數(shù)據(jù)段上的填補偏差變大。原因在于,單變量填補方法在缺失數(shù)據(jù)所屬的變量維度上,以k個點為初始值進行填補,能夠較好地保持序列趨勢;但每次填補都以上次填補的結(jié)果作為輸入、存在累積誤差,難以體現(xiàn)波動細節(jié)。多變量填補方法利用其他多個相關(guān)變量信息,每次填補都以已知數(shù)據(jù)作為輸入,不存在累積誤差,能夠體現(xiàn)波動細節(jié);但沒有利用缺失數(shù)據(jù)所屬變量信息,在某些點上存在較大偏差。
圖6 缺失數(shù)據(jù)段長度對填補效果的影響Fig.6 Relation between the length of missing data segment and its filling effect
表3 4種方法的參數(shù)設(shè)置
相比之下,組合權(quán)重、組合閾值填補的精度較高,受缺失數(shù)據(jù)段長度的影響不大。這是因為兩種方法利用單變量和其他相關(guān)多變量信息,較好地體現(xiàn)了波動細節(jié),也有效避免了較大偏差。
然而,組合權(quán)重填補對單變量和多變量填補結(jié)果加權(quán)平均,用一個權(quán)重值處理所有的缺失數(shù)據(jù),難以體現(xiàn)各個缺失數(shù)據(jù)的差異性;此外,權(quán)重值不易確定也限制了其使用。組合閾值填補以多變量填補結(jié)果為基礎(chǔ),當多變量與單變量填補的差異度較大時,以兩者的均值作為填補結(jié)果,既體現(xiàn)了填補的差異性,也避免了權(quán)重值的選擇。
在小樣本情況下,針對多元時間序列中的缺失數(shù)據(jù),利用缺失數(shù)據(jù)所屬單變量和其他相關(guān)多變量信息,提出了一種基于LSSVM的組合閾值填補方法,并進行了實驗驗證。結(jié)果表明,它能較好地體現(xiàn)序列的波動細節(jié),且受缺失數(shù)據(jù)段長度的影響相對較小。
與組合權(quán)重填補方法相比,組合閾值填補方法的參數(shù)(差異度閾值)較容易確定,且能較好地體現(xiàn)填補的差異性,填補精度較高。
[2] 韓敏, 許美玲, 任偉杰. 多元混沌時間序列的相關(guān)狀態(tài)機預測模型研究[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] 馬茜,谷峪,李芳芳,等.順序敏感的多源感知數(shù)據(jù)填補技術(shù)[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] 呂政, 趙珺, 劉穎, 等. 基于最大方差權(quán)信息系數(shù)的煤氣數(shù)據(jù)填補[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] 雷春麗, 芮執(zhí)元. 基于多元自回歸模型的電主軸熱誤差建模與預測[J]. 機械科學與技術(shù), 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] 許磊, 張鳳鳴. 缺失飛參數(shù)據(jù)填補的組合方法研究[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] 謝川,倪世宏,張宗麟.一種缺失飛行參數(shù)預處理的新方法[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的缺失飛行數(shù)據(jù)組合填補方法[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.