彭沛,趙永平,*,王雨瑋
1.南京航空航天大學 能源與動力學院,南京 210016
2.北京航空工程技術研究中心,北京 100076
在大數據時代,航空發(fā)動機采集的傳感器信息的應用由最初的事故原因調查擴展到發(fā)動機設計性能評估、飛行訓練考核輔助評價、協(xié)助地勤人員進行狀態(tài)監(jiān)控及視情維護等方面。伴隨著飛機整體性能的提升和發(fā)動機結構的復雜化,人們對飛行安全及發(fā)動機可靠性方面將會提出更高的要求,因而充分挖掘并利用飛行中發(fā)動機采集到的數據信息將會成為后續(xù)研究的必然趨勢。傳感器的測量過程由于發(fā)動機工作的復雜環(huán)境而易受隨機噪聲的影響,這就導致地面研究人員對原始獲取的數據難以直觀地分析研究。因此,如何將外場獲取的數據最大程度還原出發(fā)動機在飛機飛行時真實的工作狀態(tài),這對后續(xù)基線模型的建立、實施狀態(tài)監(jiān)控及維護決策有重要的意義[1-4]。但研究人員處理發(fā)動機獲取的傳感器數據目前仍采用人工判讀的方法,尤其在批量分析時會因判讀數據量大而耗時耗力。另外,發(fā)動機環(huán)境參數的不同和個體制造的差異會導致判讀規(guī)律的變化,這也增加了研究人員對工作模式判讀的難度。
國內外為解決上述問題做出很多嘗試。王奕惟等[5]將快速存儲記錄器(Quick Access Re‐coder,QAR)的數據圖像化并采用卷積神經網絡來判別狀態(tài)。該方法本質上是對時間序列人為采樣的單點識別而未考慮到傳感器數據的動力學性質,這割裂了時序上的連續(xù)性,并且在識別過程中該方法也未考慮過渡態(tài)。CWT-CNN 方法[6]結合了小波變換和卷積神經網絡,但其識別的單維時序樣本并非含多種模式的信號,而是人為分割后單一模式的信號。此外,Cynthia 等[7]提出了將剪枝精確線性時間法[8]和自組織神經網絡結合并給出過渡態(tài)的判別方法。該模型成功地識別了飛機的多個飛行狀態(tài),但僅涉及單維時序信號如油門桿角度,這樣就忽視了真實數據中的時滯效應,即按油門桿角度切換模式時,其他監(jiān)視參數如轉速、瞬時耗量等不會同步地發(fā)生改變。
多維時序數據的挖掘在機械零件的傳感器分析、金融財政等早有廣泛的應用。此類問題可細分為時序分割、模式發(fā)現、模式匹配等任務。例如基于動態(tài)規(guī)劃的算法[9]、剪枝精確線性時間法、基于模式的隱馬爾可夫模型[10]、基于高斯假設的貪婪分割模型[11]。這幾類方法僅能由不同模式來分割序列,并不能直接識別出模式,而且這些方法對參數調節(jié)很敏感,例如片段數和相關閾值等。在航空發(fā)動機收集到的大批量多維的時序數據下,不同航段采集的時序數據包含的工作模式數量和類別均不同,不同時序數據的要調節(jié)的參數有很大差別且隨著數據量的增加調節(jié)參數往往會耗費大量的時間。在這種場景下理想的方法應該能尋找任意的模式并不受參數的影響,即無參數化。AutoPlait[12]是基于多級鏈隱馬爾可夫模型,它利用無損壓縮原則實現模型無參數化。Auto-Cyclone[13]也采用相似的思路實現模式季節(jié)性挖掘。但上述2 種方法并不適合于發(fā)動機的工作狀態(tài)的挖掘,它們更傾向于挖掘有季節(jié)性的模式而模糊了發(fā)動機的穩(wěn)態(tài)模式和過渡態(tài)模式的界限,并且復雜的過渡態(tài)模式也無法給出很好的處理方法。事實上,上述2 種無參數方法更適合于飛行動作識別,但模型訓練也需要豐富的飛行姿態(tài)數據。
考慮到現有方法在實際工程中遇到的諸多問題,本文針對航空發(fā)動機工作時產生的時序數據,設計了一個快速自動挖掘工作模式的方法,即AutoMiner。該方法自然地集成了時序分割、模式發(fā)現和模式識別任務,它以最原始的全航段多維時序數據為研究對象,采用最小編碼代價的原則進行模型無參數化,以應對多變的時序數據。并且該方法還借鑒了標記傳播的思想[14]來解決過渡態(tài)的識別問題。在模型訓練方面,該方法支持并行化并利用儲存單元省去多次重復計算,這大大縮短了訓練時間。本文模型的設計初衷為航空發(fā)動機工作模式的挖掘,但由于模型的可遷移性,還能被應用于更高級模式的識別,如飛行動作。最后,在準確度、迭代效率、可遷移性等方面,本文開展一系列實驗,通過與多種方法對比可視化地展示了AutoMiner 方法的優(yōu)越性。
該節(jié)介紹工作模式的線性模型,由最小化編碼損失的方式提出優(yōu)化目標,并給出利用標記傳播法識別過渡態(tài)在內的發(fā)動機工作模式的方法。
從航空發(fā)動機上獲取指定的傳感器數據,可組成多維時間序列X=[x1,x2,…,xn],xi為第i時刻的d維列向量。假設模型根據X的本征模式能將其分割成K段,即S={s1,s2,…,sK},第i段si可視為X中由斷點τi和τi+1分割的部分,則斷點序號集組成τ={τ1,τ2,…,τK+1},其中τ1=1,τK+1=n。S中每段按照工作模式可以賦予相應標記并組成F={f1,f2,…,fK},fi∈{1,2,…,c}表示第i段si屬于工作模式fi,c為工作模式總數。
當航空發(fā)動機處于相同工作模式時,內在物理規(guī)律應是相似的,工作模式的轉變會引起這種內在物理規(guī)律的變化,而這種規(guī)律體現到可測傳感器數據上就是多維時序數據的變化。由此思路,對S中每個片段建立數學模型,自然地認為相同工作模式的數學模型是相似的,而不同工作模式的模型會有顯著區(qū)別。對第i段si可建模為
在計算過程中,編碼描述數據X及儲存算法模型需要付出一定的代價,而最小編碼代價對應的模型是最理想的,這個思想與經典的奧卡姆剃刀定律相似。于是,可建立優(yōu)化目標為
式中:CT表示編碼的總代價,即本文所需最小化的優(yōu)化目標,分別由描述算法模型的編碼代價和算法模型應用到數據X后產生的編碼代價組成。
假設本文討論的編碼過程均是無損的,定義符號log*(a) 為對整數a的編碼長度[15],即log*(a)=log2(a)+log2(log2(a))+…,按此規(guī)律求和且只有正數項包含在求和項中。AutoMiner算法模型簡式為:{Θ1,Θ2,…,ΘK,S,F,K,r},r代表S中片段的總類別數,區(qū)別于1.1 節(jié)中提及的工作模式數c。例如本文涉及到8 個穩(wěn)態(tài)模式,則包含過渡態(tài)后總的模式數c為9,在優(yōu)化目標中總類別數r的計算方法為
式中:IF(?)為判斷函數,例如對IF(l∈F)而言,若語句l∈F為真,則取值1,反之為0。式(3)中r將前8 個穩(wěn)態(tài)按是否出現計算,第9 個表示的過渡態(tài)按出現次數計算。
式中:cF為浮點系數[16],其默認取值為4×8,即32 bit;q為線性模型的階數。
由①~⑤便可獲得算法模型編碼所需的代價,即式(2)第1 項,而式(2)求和中的第2 項,體現模型在給定數據X完成任務的能力,即序列分割及模式挖掘的準確度。本文采用哈夫曼編碼,將編碼代價轉化為線性模型的似然概率,得到
式中:Θ={Θ1,Θ2,…,ΘK};P(X|Θ) 為似然函數值。
考慮到同一工作模式中線性模型的相關參數應盡可能相似,利用式(1)中偏差項服從正態(tài)分布的性質,構造對數似然函數為
式中:det 表示對矩陣取行列式;l(i)表示由第i個片段計算得到的對數似然函數。
在序列分割過程中,當片段長度|si|小于維數d時會導致矩陣求逆出現病態(tài),本文引入正則化的方法將式(6)進行改造,如式(7)所示,這種正則化方法常用于解決高維時間序列的問題[17]。
式中:λ為正則化系數,常取較小值,本文取10?5;tr 表示矩陣的求跡運算。
通過極大化對數似然函數能求解出參數Θ,令式(7)對Θ的偏導數分別為0,得
式中:Ιd表示d階單位矩陣。當不考慮正則化時即λ=0,Σi與相等。式(7)求偏導數涉及到諸多矩陣運算,對式(8)更詳細的推導可見附錄A。
于是,將式(8)結果重新代入式(7),可得到最大似然值的結果,并結合式(5)和式(6)可得
綜合式(5)、式(9)可得式(2)第2 項的表達式為
式中:為便于后續(xù)對編碼代價迭代求解斷點,將結果化為?(τi,τi+1)的求和,且該項只與斷點有關。
最后可以通過①~④及式(4)、式(10)獲得編碼的總代價CT,即模型的優(yōu)化目標
式中:優(yōu)化目標僅為斷點集τ的函數,隨后便可尋優(yōu)出合適的斷點組合,并使總編碼代價最小。
在計算式(11)中編碼的總代價CT前,還需在尋優(yōu)迭代過程中對片段集S賦予合適的標記集F,即找到每段隸屬的工作模式。在計算出F后便可由式(3)計算出總類別數r。
根據本文含過渡態(tài)模式的工程場景可做出如下描述:在F標記集中,第i段si可以劃分為模式fi∈{1,2,…,9},其中前8 個模式為穩(wěn)定態(tài),第9 個模式為過渡態(tài)。假設已獲得待識別片段序列為s0,對該段序列建立線性模型可獲得參數Θ0并組成樣本點。參照行業(yè)工程手冊以及專家建議可按同樣方式選取特征生成帶標記的參照樣本集表示采集第i個參照片段數據后維的特征向量,實驗中為Θ0包含的參數總數。記參照樣本集的對應標記矩陣為YL=[y1,y2,…,ym],9 維的標記向量yi中除第l個元素yil為1,其他項均為0,即表示參照數據屬于第l類模式。由于可獲取的參照樣本僅為前8類的穩(wěn)態(tài)模式,則參考樣本的標記中yi9均為0。
式中:迭代過程中的初始標記值Y(0)被賦值為先驗標記矩陣Y,其值被預先人為設置。式(13)中第2 項表示先驗信息在傳播時的參與程度。為讓傳播過程中僅改變標記y0,比例矩陣Dα設置為分塊矩陣,α越大即待標記數據會接收更多其他點的信息,本文默認取值為0.99[18]。
不難證明式(13)的迭代必收斂于式(14)。
由式(14)獲得軟標記y0,其每個元素表示樣本隸屬對應模式的概率,選取最大概率對應的序號l即為所屬的模式。本質上,參照樣本集雖只含穩(wěn)態(tài)模式,但該方法利用標記傳播更新生成了軟標記,使得過渡態(tài)模式對應的軟標記值在這里充當著剩余概率的角色,即樣本不屬于前8 種穩(wěn)態(tài)模式的概率。
時序模式識別步驟如算法1 所示。為盡可能降低計算復雜度,在采用式(13)計算W時對已知的距離矩陣無需重復計算。于是由算法1 步驟2 可將算法1 的計算復雜度由O(m2)降至O(m)。
算法1 工作模式時序挖掘算法Algorithm 1 Working mode time-series mining algorithm
式(12)得到總編碼代價僅與斷點τ有關,問題轉化為尋找最優(yōu)的斷點組合使總編碼代價最小,在該斷點組合下對應的片段集S以及標記集F即是AutoMiner 方法挖掘時序數據X獲得的最優(yōu)方案。在AutoMiner 方法中提出的斷點優(yōu)化算法是種啟發(fā)式算法,它無需全局遍歷所有的解而是在迭代范圍內找尋局部最優(yōu)的組合。每一次迭代時添加一個斷點,并局部調整斷點組合以保證總編碼代價總是維持著最小值的狀態(tài),直到找到總編碼代價不再減小的斷點組合。這個思路類似標準的“自上而下”方法[19],但該方法容易陷入局部最優(yōu)并且承擔著較高的計算成本。本文將提出的啟發(fā)式算法細分為圖1 所示的3 個步驟來分別介紹。
圖1 斷點尋優(yōu)示意圖Fig.1 Schematic diagram of breakpoint optimization
啟發(fā)式方法如圖1(a)所示。從式(12)中可以看出,每次迭代添加的斷點僅對斷點所在的片段產生影響。于是做出下述約定:對于片段si而言,在其內添加的一個斷點t并將片段si分別分割為左右2 個子片段,其總編碼代價的變化量為
式中:g是K0、r0、K1和r1的函數,它由式(4)和式(11)計算獲得。式(10)中展示了與斷點有關的函數?的表達式。K0和r0為添加斷點t前的總斷點數和總類別數,同理添加斷點后可表示為K1和r1。易知K1=K0+1,而r0和r1可由式(3)獲得。
在si內尋找最優(yōu)斷點tc的思路如算法2 所示,其中算法中總斷點數K、添加斷點前的標記集F等均默認已知。算法2 步驟1 表示從τi+T到τi+1?T以尋優(yōu)步長T來迭代選取t,設置較低步長有利于精確的搜索但會增加計算時間。
此外,算法2 還采取了其他減少計算開銷的方法。由于對于算法2 步驟1~步驟7 的循環(huán)之間并未存在聯(lián)系,即可支持并行化計算來節(jié)省大量計算時間。引入儲存單元Pi是為了將此次尋優(yōu)結果保存下來,防止后續(xù)迭代過程循環(huán)調用算法2。而選擇總代價的變化量作為判斷指標,也從一定程度上,避免了直接比較總編碼代價帶來的計算高額量。若考慮并行化,算法2 對于si的計算復雜度可由降為O(md2)。
算法2 時序片段內部最優(yōu)斷點選擇策略Algorithm 2 Internal optimal breakpoint selection strategy for time-series segments
由于尋優(yōu)算法的策略是啟發(fā)式的,為了最大程度避免陷入局部最優(yōu)解,在每次迭代添加斷點后需要動態(tài)回顧全局,再次更新斷點組合以保證求解的斷點組合的總編碼代價保持最低值。參考圖1(b),更新策略如算法3 所示。當K不大時更新次數將是適度的,而隨著K的增大,每次更新只影響局部斷點且多次迭代后斷點集趨于最優(yōu)組合,所以在容錯閾值范圍ξ內更新一定會在適度次數后收斂,實驗中ξ=5。在算法3 步驟2中,棧儲存添加點或變化點,隨后每次僅對棧彈出點的鄰域最近兩斷點進行更新。在斷點改變后也會同步更新原來的片段,這時與片段對應儲存單元中的內容將會清除并等待下次調用算法2 時再次賦值。
算法3 添加斷點后的動態(tài)更新策略Algorithm 3 Dynamic update strategy after addition of breakpoints
依照圖1(c),算法4 中展示了尋優(yōu)的總體思路。算法4 中步驟3~步驟12 的循環(huán)是獨立運行的并支持并行化運算。在算法4 步驟2~步驟27中,若片段si未因斷點更新改變且該時序片段的最優(yōu)斷點未被此次迭代選中,下一次迭代開始時便無需對si重復尋優(yōu),這也是設置存儲單元Pi的目的,因此每個存儲單元總是唯一地指向片段si。若斷點平均更新次數為L,算法3 的計算復雜度為,K為求解獲得的最優(yōu)斷點數??紤]并行化后,計算復雜度為。
算法4 斷點尋優(yōu)的啟發(fā)式算法Algorithm 4 Breakpoint-seeking heuristic algorithm
在某航段上截取部分簡單的三維時間序列并可視化地展示AutoMiner 方法斷點尋優(yōu)的工作流程圖,如圖2 所示,其中M2 為過渡態(tài)模式,而M1、M3、M4對應不同的穩(wěn)態(tài)模式。算法初始目標是添加單個合適的斷點使得總代價下降,這僅需對時序數據掃描一次便可得圖2(a)。隨后,圖2(b)~圖2(d)為逐步添加斷點的過程。由于算法4步驟22 加入了終止條件,算法將在第6 次迭代時終止并將第5 次的迭代結果作為最優(yōu)斷點組合輸出。如圖2(e)所示,若不加以約束添加斷點總代價反而會隨之增加。
圖2 AutoMiner 方法尋優(yōu)流程Fig.2 Optimization process of AutoMiner method
圖3 是將AutoMiner 方法運用于實際工程場景中的示意圖??梢钥闯鲈摲椒▽耐鈭霾杉l(fā)動機上原始的傳感器數據進行自動處理以獲得大批量數據能進行后續(xù)的分析。參考樣本集按研究的時序對象而變化,??赏ㄟ^專家研究設計手冊或借鑒歷史已標注工作模式的數據來獲得。此處歷史數據中包含了不同退化狀態(tài)下的工作模式數據信息,但實際中退化量與出廠的理論值通常偏差不大,選擇手冊的信息作為參照也是合理的。
圖3 實際工程場景應用的示意圖Fig.3 Application diagram of engineering scenario
本文以某型機載航空發(fā)動機在飛行過程中采集的多組傳感器數據為研究對象,實驗中涉及3 組飛機全航段工作數據,此處一個航段記錄為飛機的一次長時間飛行任務。這些數據囊括了實驗考慮的全部工作模式,在飛機飛行時由4 臺同型號發(fā)動機的傳感器進行同步測量。
對于實驗參數的篩選需要參考實驗手冊,選取隨工作模式的變化量最大的多個物理指標,比如絕對氣壓高度體現了飛機起飛降落及巡航狀態(tài),從而也反映發(fā)動機工作模式的變化。在工作模式需要改變時,飛行員會通過外部干預調節(jié)油門桿角度,這會引起發(fā)動機轉子轉速、瞬時耗量等物理量的變化。為解決參數的時滯效應,本文選擇上述4 類物理指標作為不同工作模式的分辨依據,其中油門桿角度、轉速和瞬時耗量均由2 組以上傳感器參與測量,收集以上數據可共同組成描述發(fā)動機工作狀態(tài)的多維時間序列。后續(xù)實驗將證明,考慮發(fā)動機傳感器的多維時序數據比常規(guī)只考慮單組油門桿角度作為辨別依據的方案更能適應復雜的數據情況。本文涉及的發(fā)動機工作模式的人為編號如表1所示。
表1 某型航空發(fā)動機工作模式說明Table 1 Description of operating mode of aero-engine
飛機常需要應對復雜的飛行條件甚至完成某些高難度的動作,這使得發(fā)動機上測量的參數不可避免地會產生波動,并且熱力及電磁環(huán)境干擾和不穩(wěn)定的傳感器測量狀態(tài)也會引入噪聲。為了減少噪聲以提高挖掘時序數據信息的準確度,便需要對采集到的數據進行平滑降噪。選取某航段上一臺發(fā)動機的油門桿角度進行比較,實驗方法包括滑動中值濾波[20]、滑動均值濾波以及小波軟閾值濾波[21]。對比實驗結果如圖4 所示。數據點之間采樣間隔為1 s,考慮到濾波時盡可能保留過渡態(tài)模式信息,選擇滑動窗口的長度為60 s,同時采用8 層小波分解來計算近似系數。
圖4 3 種濾波方法的對比結果Fig.4 Comparative results of three filtering methods
由于中位數面對異常數據以及高頻噪聲更加魯棒,從圖4 中可以看出相比均值濾波而言,中值濾波能得到更適合識別工作模式的時序數據。在小波軟閾值濾波中,不同分解的層數及合適的小波基會隨著時序信號的改變而變化,因此不適用于大批量自動處理數據的場景。
為展示AutoMiner 方法在航空發(fā)動機工作模式挖掘場景下的優(yōu)越性,本節(jié)提取部分航段中傳感器采集的時序數據與多種時序分割算法對比,包括基于高斯核的KCP 方法[22]和基于高斯模型的GSS 方法[11]。為在模式挖掘和無參數化算法方面進一步比較,對比算法中還包括成功運用于飛行模式挖掘的PELT-SOM 算法[7],以及基于多級鏈隱馬爾可夫模型的AutoPlait無參數化算法[9]。
3.3.1 時序分割能力的對比評估
在時序分割方面,合適的評價指標至今仍是廣泛討論的內容。例如,專家評判分割斷點在第1 000 s,模型求解出斷點為第1 005 s 是否有效。若容錯閾值為10,模型得到的斷點為第1 011 s是否應該被懲罰。嚴格的懲罰制度得到準確度、混淆矩陣等指標并不適合時序分割問題[23]。上述問題歸結于對評判邊界并未做出合適的處理,本文借鑒文獻[24]的思路選擇評價指標為R分數,其表達式為
式中:τ(1)為模型求解的斷點集。數據真實的斷點集τ(2)中距離斷點最近的斷點表示為,此處討論不包含時序數據的起點和終點。指標R的取值范圍為0~1,其值越大說明該模型分割效果越好。
采集部分飛行航段上發(fā)動機傳感器數據進行了2 組實驗,實驗中各物理量僅展示1 個傳感器的測量結果并通過一定比例地無量綱化操作繪制出如圖5、圖6 所示的曲線,其中曲線圖中的分割線為人工結合工程經驗給出,即可視為真實的分割斷點位置。另外,圖5(b)、圖5(c)、圖5(f)和圖6(b)、圖6(c)、圖6(f)涉及算法均可支持工作模式的識別,而圖5(d)、圖5(e)和圖6(d)、圖6(e)中僅支持對時序數據的分割,其中對比算法中參數β為懲罰系數,λ為正則化系數,M 表示某些對比算法對該工程場景失效而僅識別出的未知模式。
圖5 在航段1 中工作模式挖掘的可視化圖Fig.5 Visualization of pattern mining in Flight segment 1
圖6 在航段2 中工作模式挖掘的可視化圖Fig.6 Visualization of pattern mining in Flight segment 2
將時序分割結果代入式(16)中可得到表2。理想的算法應該具有較高的R分數且生成斷點數接近于真實斷點數。由于圖5、圖6 中可以看出航段2 相比航段1 包含更豐富的工作模式,在復雜的數據組成表2 中航段2 的R分數也會相較航段1 偏低。對于PELT-SOM 和KCP 算法而言,懲罰系數越小,意味著允許生成更多的斷點,一方面會使R分數在一定程度上增大,但過多的斷點會增加計算時長,同時使模型結果過于復雜而不利于后續(xù)分析。從圖5(c)和圖6(c)中可以看出,若懲罰系數過小,算法會忽視某些必要的過渡模式。就GSS 算法而言,相比上述2 種算法表現相對更好,但其不支持模式挖掘功能并且需要提前調節(jié)超參數,它并不是此場景下最好的選擇。由于AutoPlait 算法設計初衷為捕捉更高級復雜的模式,且并不能適應過渡態(tài)和數據的無規(guī)律變化,這導致該算法在2 個航段上的實驗并不理想。綜合表2 和圖5、圖6 不難看出,AutoMiner 方法在航段1 和航段2 上較其他幾類方法表現更優(yōu)越。
表2 在航段1 和航段2 上分割評價指標R 分數的對比Table 2 Comparison of split evaluation indicators R score on Flight segments 1 and 2
3.3.2 算法可遷移性的實驗解釋
為研究AutoMiner 算法在更高級模式的挖掘能力,隨機選擇有規(guī)律性的時序數據,改變模型階數可獲得圖7 中的結果。由航段3 的實驗可以看出,復雜的模型可以挖掘出由多個工作模式組成的更高級的模式。圖7(c)、圖7(d)中斜線塊所示為高級的模式,由于未提供相應高級模式的參照數據集,其自然地被視為過渡態(tài)模式。模型階次越高意味著算法采用了更復雜的模型,這可以挖掘出更復雜的模型狀態(tài)。如圖7(e)的算法中選擇多級鏈馬爾可夫模型便可以挖掘出圖中M1~M5 所示的更高級的工作模式。
圖7 在航段3 中算法遷移性實驗Fig.7 Algorithm transportability experiments in Flight segment 3
3.3.3 實驗迭代過程的討論
在未采用并行化措施時,AutoMiner 算法計算復雜度為。從航段2 的實驗可以看出,隨著迭代次數的增加,每步迭代需要的時間將快速減少,這是因為隨著斷點的增多對應時序片段長度較少,加上引入存儲單元也大大縮短了后續(xù)迭代所需的計算時間。圖8 中T為采樣步長,圖例中括號為與之對應的評價指標R分數。從圖8 中可以看出適當增大采樣步長會使算法總耗時急劇減少,但在一定程度內分割指標的減少量并不大,所以可以根據實際情況適當增加采樣步長犧牲較少的精度以換取較低計算成本。
圖8 航段2 中算法總耗時隨迭代次數的變化Fig.8 Variation of total time consumed by algorithm with number of iterations in Flight segment 2
圖9 為AutoMiner 算法在航段1~航段3 上求解斷點組合時總代價隨迭代次數的變化圖,其中選用的模型默認為一階線性模式。為獲得完整的總代價趨勢變化圖,將算法4 步驟22 的中止功能暫時關閉。在圖9 中變化曲線存在明顯的拐點,在迭代達到曲線極小值點后,進一步增加斷點數量反而會使總代價增大。另外,通過式(10)的分析不難發(fā)現,當模型和時序數據間的模型誤差逐漸減小,即足夠小時,編碼代價會在數值上變?yōu)樨摂?,但這與總代價最小化的尋優(yōu)策略并不矛盾。
圖9 AutoMiner 算法在迭代過程中總代價的變化Fig.9 Variation of total cost of AutoMiner algorithm during iteration process
3.3.4 模式挖掘效果的對比評估
現今模式挖掘效果的評估方法常采用準確度、F1分數,其中F1分數是由混淆矩陣得到的精確率和召回率的調和平均數。在航段1~航段3來源更完整的3 組航段上繼續(xù)開展實驗,將Au‐toMiner 算法挖掘的工作模式與PELT-SOM 算法的結果進行比較可以得到表3 所示結果。從中可以看出過渡態(tài)模式上2 種方法的各項指標均弱于穩(wěn)態(tài)模式,而相較之下AutoMiner 方法會獲得更好的結果。若時序數據位于2 種工作模式之間,例如圖5 中前1 000 s 對應發(fā)動機停轉,PELT-SOM 會將其與最相似的穩(wěn)態(tài)模式即慢車模式優(yōu)先匹配從而導致判斷失誤。AutoMiner 方法通過標記傳播獲得軟標記,比較屬于各模式的概率并把不屬于穩(wěn)態(tài)模式的概率視為過渡態(tài)從而解決了此類問題。
表3 模式挖掘評價指標的對比Table 3 Comparison of evaluation metrics for pattern mining
為了對模型獲取時序片段的工作模式開展可視化分析,選用t-SNE 方法[25]將從原始時間序列中提取的多維特征降至二維的特征圖上并繪制如圖10 所示散點圖,圖中按模式總量排布并且括號中表示通過不同模型挖掘出的對應工作模式的時序片段數量。圖10(b)中軟標記最大值代表對應工作模式的概率值,它在圖10(b)中表示為散點的形狀大小。由圖10 可以看出AutoMiner 模型獲得的工作模式更加接近于真實的工作模式。
圖10 所有時序片段的工作模式的可視化分布Fig.10 Visualization of pattern distribution of whole time-series segmentation
1)AutoMiner 方法通過對模型進行編碼,利用編碼代價最小的思想來獲取最優(yōu)的時序片段組合,并在分段時利用標記傳播來解決過渡態(tài)模式識別的問題。相較單變量的方法,本文采用的多變量方法能更好地解決狀態(tài)參數的時滯問題。
2)通過在某型發(fā)動機全航段數據上的實驗表明,在與諸多算法對比后,AutoMiner 方法在分割指標以及模式挖掘指標上均表現優(yōu)異。
3)AutoMiner 方法支持并行化計算并采取諸多措施加速迭代過程,實驗還證明了該方法在更復雜工作模式方面具有很好的可遷移性。
附錄A:
1)式(7)對θi求偏導數并令之為0,求解θi。
求解對θi的偏導數時,將Σi視為已知常數矩陣,并對式(7)的跡取微分可得:
依照微分和偏導數轉化公式并結合跡的性質可將式(A1)化為
對式(A2)化簡后得到
令式(A2)為0,即可化簡出式(8)中θi的結果。
2)式(7)對Σi求偏導數并令之為0,求解Σi。
同樣,求解對Σi的偏導數時,將θi視為已知常數矩陣,并對式(7)的跡取微分可得:
式(A4)可化為含A矩陣的形式,可以直觀看出A矩陣為對稱矩陣,便得到所求偏導數如式(A5)所示,A矩陣在此情形下可以視為雅可比矩陣。
式中:Id為d階單位矩陣;符號“°”表示矩陣的哈達瑪積,即矩陣點乘。若矩陣A由元素aij組成,令式(A5)結果為0,可得
式(A6)反映出A=0,化簡可得式(8)中Σi的結果。