許 訓李 鑒孫 罡

（1.中國電子科技集團公司第三十研究所 成都 610000）（2.電子科技大學信息與通信工程學院 成都 612731）

1 引言

故障診斷的目的是根據觀測數據來確定系統(tǒng)是否存在故障節(jié)點，如果存在，需要進一步確認故障節(jié)點的數量和位置。在航空航天、工業(yè)過程控制、汽車系統(tǒng)、計算機系統(tǒng)、無線通信網絡等領域都具有重要應用?？紤]到實際情況，我們能夠得到的觀測數據往往小于系統(tǒng)內的節(jié)點數量，所以故障診斷問題可以抽象成一個欠定方程的求解問題。如果繼續(xù)考慮到節(jié)點發(fā)生故障是一個小概率事件，那么故障狀態(tài)向量則具有稀疏性，可以借助壓縮感知的方法來完成。

近年來，出現了不少利用壓縮感知來實現故障診斷的論文，并且取得了一些研究成果。利用壓縮感知提出了一種針對帶噪聲干擾的情況下電氣系統(tǒng)的故障診斷問題的解決方法［1］；基于壓縮感知框架下壓縮采集的信號，提出了一種能夠在降低平均采樣速率的同時用更少的數據量表現故障特征，實現故障診斷的方法［2］；提出了一種基于壓縮感知和深度學習理論，利用隨機高斯矩陣實現軸承信號的變換域壓縮采集的軸承故障診斷方法［3］，能夠實現對不同故障位置和缺損程度的故障特征自動提取與準確故障診斷；利用壓縮感知實現的建筑電氣系統(tǒng)故障診斷方法，能夠滿足小樣本數據的建筑電氣故障診斷工程應用的需求［4］；基于壓縮感知理論，開展的振動信號稀疏采樣重構方法、故障特征壓縮檢測方法以及多源耦合信號分離與診斷信號的研究［5］。利用壓縮感知實現在旋轉機械遠程故障診斷中的應用［6］。對于系統(tǒng)節(jié)點只具有正常或者故障兩種狀態(tài)的常態(tài)情況時，當前對于二值故障狀態(tài)向量的研究大多只利用單范數考慮到了故障狀態(tài)向量的稀疏性，而忽略了故障狀態(tài)向量的二值性。本文提出一種利用?1，?∞范數混合優(yōu)化的方法，既保證了故障狀態(tài)向量的稀疏性，又同時考慮到了故障狀態(tài)向量的二值性，并且方法也對故障狀態(tài)向量的稀疏度依賴程度也得到了極大的優(yōu)化。

2 理論基礎

壓縮感知理論是由 David L:Donoho［7］，Emmanuel Candes，Justin Romberg 和 Terence Tao［8］等提出的，其基本思想是，原始信號經過某一變換（小波變換、傅立葉變換等）使其在某一正交空間具有稀疏性，再乘以一個投影矩陣（測量矩陣），便得到壓縮后數據，可以用于進行存儲、傳輸等。而解壓過程，則可表示為根據壓縮數據和壓縮矩陣求欠定方程最稀疏解問題。

x代表原始信號的稀疏表示，y代表壓縮數據（觀測信號），A代表投影矩陣，壓縮感知問題即可表示為

這是一個NP難問題，具有指數復雜度，Terence Tao在文獻［9］中指出，當觀測矩陣A滿足限制等距性質（RIP），即存在常數 δs∈(0，1)使

就可將上述非凸問題松弛為

從而可以利用凸優(yōu)化的方法進行求解。

3 數學模型

結合故障診斷的實際問題，本文主要考慮二值故障狀態(tài)向量所對應的數學模型。

假設系統(tǒng)一共有N個節(jié)點，x∈χN表示系統(tǒng)故障狀態(tài)向量，χ為系統(tǒng)分狀態(tài)可能的取值，一般情況下 χ={0，1}，0表示正常工作，1表示存在故障；如xi=1表示系統(tǒng)第i個節(jié)點故障，xi=0表示系統(tǒng)第i個節(jié)點正常工作。同時，節(jié)點發(fā)生故障的概率為 p，且某節(jié)點是否發(fā)生故障不受其他節(jié)點狀態(tài)影響，即xi獨立同分布于Bernuolli(p)。

y∈RM表示觀測數據，M＜N；ω∈RM表示觀測噪聲，通常為高斯噪聲，ωi相互獨立且和x也相互獨立，可以假設ωi～Normal(0，σ2ω)；A∈RM×N表示故障觀測矩陣，為了使矩陣滿足RIP性質，我們選擇Gaussian隨機矩陣作為觀測矩陣，即A所有元素aij獨立同分布于Normal(0，1/M)。故障診斷的線性測量模型可以表示為

觀測數據 y和觀測數據A已知，需要求解狀態(tài)故障向量x。

考慮到系統(tǒng)節(jié)點發(fā)生故障是一個小概率事件(p＜0.1)，故障狀態(tài)向量x具有稀疏性，進一步考慮到xi智能取0或1，可以設計如下目標函數:

其中，||x||1保證了 x的稀疏性，而||x-0.5||∞保證了x的二值性［10］，顯然目標函數是一個凸函數，可以利用凸優(yōu)化的方法解決該問題，具體實現步驟將在第四節(jié)中詳細論述。

4 算法實現

考慮用內點法［11］來處理式（5），將不等式約束的優(yōu)化問題轉換為形式如下的無約束優(yōu)化問題:

從式（6）中可以看出，當x在式（5）的約束范圍以外時，目標函數的值會趨于無窮。γ是懲罰因子，在迭代過程中會逐漸變小，當γ取很小時，式（6）和式（5）的解十分接近。進一步考慮用?p范數逼近的方法來處理 λ||x-0.5||∞［12］，式（6）又可以進一步寫成如下形式:

隨著 p的逐漸增大，式（6）的解將會逐漸逼近式（5）的解。而式（6）則可以通過牛頓迭代來求解，需要求取優(yōu)化目標的梯度以及Hessian矩陣。令

優(yōu)化目標可以寫成

梯度如下:

具體實現步驟如下:

1）數據初始化:迭代初值xiter=AT(AAT)-1y，p的初始值、步長、上限值分別為 pmin、Δp、pmax，懲罰因子 γ的初始值為 γ0，衰減系數 γ-∈(0，1)；

2）進行n次牛頓迭代（無論是否收斂），具體步驟如下:

3）更新:xiter=x*，p=p+Δp ；

4）判斷條件:p＜pmax，如果成立則返回2），否則繼續(xù)進行下一步；

5）牛頓迭代，具體步驟如下:

上面方法步驟2）中，牛頓迭代沒有必要一定要出目標函數的最優(yōu)解，只要進行固定次迭代，能夠保證算法正常進行即可。參數λ表示?∞范數項的權重，λ越小，表示目標函數中?∞范數項所占比重越小，?1范數項所占比重越大，得到的優(yōu)化結果就越偏重于稀疏性約束。λ越大，表示目標函數中范數項?∞所占比重越大，?1范數項所占比重越小，得到的優(yōu)化結果就越偏重于有限字符集約束。

本文希望得到一個同時受稀疏性和有限字符集約束的優(yōu)化結果，并且兩個約束程度應該比較接近。將原始信號x帶入目標函數，可以得到 ?1范數項的值約為幾十（跟故障概率有關），?∞范數項的值為0.5λ，希望二者是一個數量級，故本次試驗考慮取λ=100。

經過上面一系列計算，可以得到一個初步的故障狀態(tài)向量xiter，這并不是一個二值向量，其中幾乎都是小數，嚴格為0或者嚴格為1的部分很少，所以需要想辦法將xiter進行二值化，從而得到最終結果。

本文首先考慮了兩種二值化的辦法。最簡單的一種是，直接設定一個硬閥值如0.5，將xiter大于閥值的部分置1，小于的部分置0。這種做法很簡單，在噪聲較小如SNR=1000的時候，能取得不錯的效果。但噪聲增大以后，效果很差，并且硬閥值具體取多少很難確定。

另一種是，簡單的認為xiter中值越大的位置，為1概率越大，如果原始信號的稀疏度K已知，那么可以將xiter中最大的K位置1，其余的全部置0。這種方法比第一種效果要好很多，在仿真試驗中，采用這種方法進行二值恢復，大概能夠有90%的概率完全恢復出原始信號。（在故障概率 p=0.05，SNR=10的情況下，重復仿真100次，能夠完全恢復x90次）但是這種方法需要對先對稀疏度即故障個數進行估計，雖然稀疏度K的估計并不難求，且在能夠獲得多次觀測數據的情況下估計的準確率很高。

本文進一步考慮，提出一種更好的且不用顧忌系數度K的二值化辦法。可以借助正交匹配追蹤算法［13～15］將現行觀測方程 (4)寫成如下形式:

其中Γ表示x的支撐集，ai為觀測矩陣A的列向量，也稱原子。

上面說過，可以簡單地認為xiter中值越大的位置，為1的概率越大，可以先構造一個所有元素均為0的N維列向量x*。

將 x*的 k1位置1，假設ki表示 xiter中第i大元素所在位置，k1即表示xiter中最大元素所在位置。如果k1∈Γ，

由于aij獨立同分布于 Normal(0，1/M)，wi獨立同分布于Normal(0，σ2w)，根據式(20)可以推出此時r中每一個元素ri獨立同分布于Normal{0，(k-1)/M+a2w}。

如果k1?Γ

同理可知，r中每一個元素ri獨立同分布于

可以發(fā)現，如果將支撐集中的元素置1，r的方差會變小，由于r中每一個元素ri獨立同分布于某一均值為0的正態(tài)分布，可以將r看作該正態(tài)分布的觀測樣本，r的方差看作樣本方差，為分布方差的無偏估計量，可以依此設計如下二值化算法。

1）數據初始化:x*=0，ki，i=1，其中i只取到n，n表示xiter大于0的元素個數；

2）計算置1前的方差sbefore=sd(y-Ax*)，sd表示計算方差函數；3）置1后并計算置1后的方差:

4）判斷條件:safter＜sbefore，若滿足則繼續(xù)下一步，若不滿足則令=0后繼續(xù)下一步；

5）如果 i=＜n，i=i+1，返回步驟2），否則停止跌單。

該算法經測試十分有效，完全恢復率在97%左右，（在故障概率 p=0.05，SNR=10的情況下，重復仿真100次，能夠完全恢復x97次），比稀疏度已知的時候還要高。

將?p逼近的牛頓迭代算法和上述二值化辦法結合，即可通過觀測值求出故障狀態(tài)向量x。

5 仿真與結果分析

本節(jié)將利用Matlab會對前文所設計的算法HybirdOptimize，進行仿真測試并且將與利用OMP和MAP-OMP算法進行壓縮感知的支撐恢復實現故障診斷方案進行對比。

圖1表示不同故障概率下，算法對故障狀態(tài)向量的恢復效果（信噪比SNR=10）。其中a的橫坐標表示發(fā)生故障的概率，縱坐標表示完美診斷的概率，完美診斷代表通過算法推斷出的狀態(tài)故障向量和真實的狀態(tài)故障向量完全一樣。b的縱坐標表示平均每次試驗中，發(fā)生錯誤診斷的節(jié)點個數，可以用下式表示。

n表示試驗次數，實驗共進行了1000次仿真，n=1000。通常情況下，信號的恢復效果，還需要用均方誤差（MSE）來衡量，但是考慮到x的二值性，均方誤差（MSE）和錯誤診斷數量（QWD）只相差一個N的倍數而已。

圖1 故障概率對故障診斷效果的影響

由圖1可以看出?1，?∞范數混合優(yōu)化的效果最好，尤其是在稀疏度較高的時候，其對故障狀態(tài)向量的恢復效果要優(yōu)于用OMP和MAP-OMP來進行支撐恢復的方法。

圖2表示不同信噪比下，算法對故障狀態(tài)向量的恢復效果（故障概率p-fault=0.05）。其中a，b縱坐標含義與圖1相同，橫坐標則表示信噪比。

圖2 信噪比對故障診斷效果的影響

由仿真結果可以看出，?1，?∞范數混合優(yōu)化在噪聲較大的情況下，效果并不理想。其原因是在設計目標函數的時候，并沒有考慮對噪聲的抑制，只有在二值化的過程中稍微考慮的噪聲的影響，所以其抗噪能力較弱，這也是需要繼續(xù)研究提高的地方。

6 結語

對于二值故障狀態(tài)向量來說，?1，?∞范數混合優(yōu)化能夠取得非常好的效果。并且，?1，?∞范數混合優(yōu)化算法有一個很大的優(yōu)勢是，它不需要事先知道故障狀態(tài)向量x的稀疏度K就能取得很好的效果，這是MAP-OMP和OMP很難做到的。從空間復雜度和時間復雜度來說，本文所涉及的算法比起支撐恢復算法也要稍顯復雜一點，因此算法很多部分還值得繼續(xù)探究。