分解原理與MATLAB實(shí)現(xiàn):參數(shù)設(shè)置、避坑指南及工程應(yīng)用)
簡(jiǎn)介VMD變分模態(tài)分解MATLAB程序包面向信號(hào)處理研究者和工程師用于對(duì)非線性、非平穩(wěn)信號(hào)進(jìn)行自適應(yīng)分解提取不同頻率成分可應(yīng)用于噪聲去除、特征提取、模式識(shí)別、金融時(shí)間序列分析等場(chǎng)景也可延伸至機(jī)械故障診斷與生物醫(yī)學(xué)信號(hào)處理方向。壓縮包體積僅4KB包含2個(gè)m源文件與1個(gè)txt文本文件兼顧核心算法實(shí)現(xiàn)、示例測(cè)試與許可說(shuō)明m文件分別提供主函數(shù)和演示腳本txt文件為許可協(xié)議說(shuō)明便于用戶(hù)依法使用和二次修改。已有3433人學(xué)習(xí)下載。借助該程序包讀者可直接在MATLAB中調(diào)用VMD算法完成信號(hào)分解依托測(cè)試腳本快速掌握調(diào)用流程與參數(shù)設(shè)置并結(jié)合自身數(shù)據(jù)調(diào)整模態(tài)數(shù)量、懲罰系數(shù)等關(guān)鍵參數(shù)從而深入理解信號(hào)內(nèi)在結(jié)構(gòu)開(kāi)發(fā)者也可參考許可條款合規(guī)地將其集成到自己的研究或工程項(xiàng)目中。1. VMD能用在哪從模態(tài)混疊說(shuō)起的一個(gè)信號(hào)分解方案做滾動(dòng)軸承故障診斷或者結(jié)構(gòu)振動(dòng)分析的人多半經(jīng)歷過(guò)這樣的場(chǎng)景傳感器采回來(lái)一段振動(dòng)信號(hào)想從時(shí)域波形里分辨出故障特征頻率結(jié)果發(fā)現(xiàn)幾個(gè)頻率成分在時(shí)域上纏成一團(tuán)頻譜圖上也是糊成一片。這時(shí)候常見(jiàn)的做法是拿經(jīng)驗(yàn)?zāi)B(tài)分解去拆信號(hào)但拆出來(lái)的固有模態(tài)函數(shù)經(jīng)常出現(xiàn)模態(tài)混疊同一頻率碎在好幾個(gè)分量里或者兩個(gè)不同的頻率被并進(jìn)一個(gè)分量分析起來(lái)讓人頭疼。VMD變分模態(tài)分解就是在這種背景下被提出來(lái)的一套信號(hào)分解方法它把信號(hào)分解問(wèn)題從「遞歸篩分」改寫(xiě)成了「變分約束求解」把模態(tài)看成有限帶寬的子信號(hào)通過(guò)迭代更新中心頻率和帶寬把原始信號(hào)拆成K個(gè)分量。這篇文章講的不是理論推導(dǎo)本身而是圍繞一套MATLAB程序包把VMD從原理到程序?qū)崿F(xiàn)、參數(shù)設(shè)置再到常見(jiàn)的坑完整梳理一遍適合已經(jīng)有基本信號(hào)處理基礎(chǔ)、想在手頭數(shù)據(jù)上用VMD做特征提取或去噪的工程師和數(shù)據(jù)科學(xué)方向的研究者。2. VMD的數(shù)學(xué)內(nèi)核與MATLAB程序包的實(shí)現(xiàn)骨架2.1 三個(gè)約束條件為什么VMD能把信號(hào)拆得比EMD更規(guī)整理解VMD的變分模型關(guān)鍵是抓住三個(gè)約束條件。第一個(gè)約束是每個(gè)模態(tài)圍繞一個(gè)中心頻率展開(kāi)且在頻域上具有有限的帶寬這是VMD區(qū)別于經(jīng)驗(yàn)?zāi)B(tài)分解的根本設(shè)計(jì)——EMD是靠極值點(diǎn)包絡(luò)迭代篩出本征模態(tài)函數(shù)而VMD把每個(gè)模態(tài)建模成調(diào)幅調(diào)頻信號(hào)它的帶寬由解調(diào)譜的梯度范數(shù)來(lái)衡量。第二個(gè)約束是所有模態(tài)之和要能重構(gòu)原始信號(hào)這是信號(hào)分解問(wèn)題的保真性條件保證分解不會(huì)丟信息。第三個(gè)約束是模態(tài)之間在頻域的分離性VMD通過(guò)交替方向乘子法迭代求解構(gòu)造出的變分問(wèn)題在每一步更新模態(tài)、中心頻率和拉格朗日乘子讓每個(gè)模態(tài)在迭代中逐步收斂到它應(yīng)當(dāng)占據(jù)的頻帶中心。這套模型落到MATLAB程序包里最常見(jiàn)的實(shí)現(xiàn)方式是主函數(shù)加若干嵌套子函數(shù)的結(jié)構(gòu)。主函數(shù)接收信號(hào)、分解層數(shù)K、懲罰因子alpha等參數(shù)內(nèi)部初始化一組模態(tài)分量和中心頻率然后進(jìn)入交替方向乘子法迭代循環(huán)。循環(huán)體里幾個(gè)關(guān)鍵更新式對(duì)應(yīng)的是頻域運(yùn)算所以VMD在MATLAB里天然要做FFT和IFFT的往返。一個(gè)典型的迭代步長(zhǎng)內(nèi)要做的事是對(duì)每個(gè)模態(tài)在當(dāng)前頻域估計(jì)上應(yīng)用維納濾波更新模態(tài)再更新模態(tài)的中心頻率為頻率響應(yīng)的重心位置最后更新拉格朗日乘子以強(qiáng)制重構(gòu)約束。從程序結(jié)構(gòu)上說(shuō)這套迭代不難讀難的是把邊界條件和收斂判據(jù)處理好。2.2 讀懂程序包里的主循環(huán)各步驟對(duì)應(yīng)的數(shù)學(xué)含義我一般拿到一個(gè)VMD程序包先不看那些參數(shù)解析函數(shù)直奔主循環(huán)。一個(gè)標(biāo)準(zhǔn)的主循環(huán)大致長(zhǎng)這樣% 主循環(huán)ADMM迭代求K個(gè)模態(tài) u_hat zeros(N, K); % 頻域模態(tài)存儲(chǔ) omega zeros(K, maxiter); % 中心頻率迭代記錄 tau 0; % 噪聲容忍度0為強(qiáng)制重構(gòu) for n 1:maxiter for k 1:K % 更新第k個(gè)模態(tài)維納濾波步 sum_u_hat sum(u_hat, 2) - u_hat(:, k); u_hat(:, k) (f_hat - sum_u_hat - lambda_hat/2) ... ./ (1 2*alpha*(freq - omega(k, n)).^2); % 更新中心頻率頻域能量重心 omega(k, n1) (freq * (abs(u_hat(:, k)).^2)) ... / sum(abs(u_hat(:, k)).^2); end % 更新拉格朗日乘子 lambda_hat lambda_hat tau*(f_hat - sum(u_hat, 2)); % 收斂判定 e sum(abs(sum(u_hat, 2) - f_hat).^2)/sum(abs(f_hat).^2); if e tol break; end end這套代碼邏輯說(shuō)穿了就是三步循環(huán)更新模態(tài)、更新中心頻率、更新拉格朗日乘子。注意第8行的維納濾波更新里分母上的(12alpha(freq - omega).^2)就是頻域里的帶寬約束項(xiàng)alpha越大模態(tài)的頻帶被壓得越窄。程序包里這段代碼的循環(huán)上限maxiter通常設(shè)500但絕大多數(shù)信號(hào)在幾十步內(nèi)就滿足收斂判據(jù)了。如果看到循環(huán)跑滿次數(shù)還沒(méi)收斂先檢查信號(hào)里是否包含明顯高于奈奎斯特頻率附近的分量或者K設(shè)得和實(shí)際模態(tài)數(shù)不匹配。2.3 主函數(shù)簽名的參數(shù)逐一拆解alpha、tau、DC、init的玩法MATLAB程序包里常見(jiàn)的函數(shù)簽名是[u, u_hat, omega] VMD(f, alpha, tau, K, DC, init, tol)不同程序包的參數(shù)順序可能不一樣但核心參數(shù)就那么幾個(gè)。alpha是懲罰因子控制模態(tài)帶寬的緊致程度取值在200到3000之間比較常用alpha越小每個(gè)模態(tài)的帶寬越寬越能容納頻率漂移alpha越大模態(tài)越窄但容易在某些信號(hào)里把真實(shí)分量劈成好幾段。tau是噪聲容忍參數(shù)為0時(shí)要求嚴(yán)格重構(gòu)信號(hào)里的每個(gè)細(xì)節(jié)都要保留設(shè)為0.5這種中間值會(huì)讓程序允許一定重構(gòu)誤差等價(jià)于給分解引入平滑效果適用于強(qiáng)噪聲場(chǎng)景。DC參數(shù)設(shè)為1時(shí)把第一個(gè)模態(tài)的中心頻率固定在0 Hz這等效于強(qiáng)制把信號(hào)的直流分量或者極低頻趨勢(shì)項(xiàng)分離成單獨(dú)一個(gè)模態(tài)。init決定中心頻率的初始化方式init1用均勻分布的頻率作為初始值適合完全不清楚信號(hào)頻率成分的情況init0則全部從0開(kāi)始迭代收斂慢但對(duì)初始猜測(cè)不敏感。這四個(gè)參數(shù)是程序的四個(gè)旋鈕調(diào)參時(shí)不要一次性全動(dòng)先固定alpha和tau掃K值找到模態(tài)數(shù)的合理范圍再回頭細(xì)調(diào)alpha。程序包里還經(jīng)常會(huì)見(jiàn)到一個(gè)輔助函數(shù)用于頻域坐標(biāo)生成和信號(hào)長(zhǎng)度補(bǔ)齊主要是為了FFT計(jì)算方便在信號(hào)長(zhǎng)度不是2的冪時(shí)也不會(huì)報(bào)錯(cuò)性能只會(huì)有微小差異。3. 跑通你的第一組VMD分解準(zhǔn)備數(shù)據(jù)、調(diào)用程序包3.1 用仿真信號(hào)驗(yàn)證程序包是否裝好三個(gè)正弦疊加的指令序列拿到程序包先別急著處理真實(shí)數(shù)據(jù)用一個(gè)構(gòu)造的仿真信號(hào)把程序跑通同時(shí)驗(yàn)證程序包里各功能模塊是否正常。常見(jiàn)的做法是構(gòu)造三段正弦疊加信號(hào)頻率選在比如50 Hz、120 Hz、200 Hz采樣率設(shè)為1000 Hz采樣點(diǎn)數(shù)1024。代碼塊如下fs 1000; t (0:1023)/fs; f1 50; f2 120; f3 200; sig sin(2*pi*f1*t) 0.8*sin(2*pi*f2*t) 0.5*sin(2*pi*f3*t); % 調(diào)用VMD分解 [K, alpha] deal(3, 2000); [u, ~, omega] VMD(sig, alpha, 0, K, 0, 1, 1e-7);運(yùn)行完以后omega矩陣?yán)镒詈笫諗康娜兄行念l率應(yīng)該基本收斂到50、120、200附近。我把這組仿真當(dāng)作程序包的自檢用例如果中心頻率偏離超過(guò)正負(fù)2 Hz先檢查信號(hào)是否加了窗因?yàn)榧哟皶?huì)在頻域展寬主瓣VMD對(duì)這類(lèi)效應(yīng)敏感再檢查alpha是否過(guò)小導(dǎo)致模態(tài)帶寬太大把相鄰兩個(gè)頻率的成分吸到了一起。仿真信號(hào)的好處是模態(tài)數(shù)是已知的可以用它來(lái)反向驗(yàn)證程序包的收斂精度實(shí)測(cè)經(jīng)驗(yàn)中tol從1e-7放寬到1e-6對(duì)分解結(jié)果影響不大但迭代步數(shù)能省掉接近一半。3.2 真實(shí)數(shù)據(jù)場(chǎng)景長(zhǎng)度選擇、去趨勢(shì)、歸一化的標(biāo)準(zhǔn)流程真實(shí)數(shù)據(jù)比仿真信號(hào)復(fù)雜的地方在于噪聲、趨勢(shì)項(xiàng)和幅值動(dòng)態(tài)范圍。我在處理實(shí)測(cè)振動(dòng)數(shù)據(jù)時(shí)習(xí)慣先做一次去均值去趨勢(shì)用detrend函數(shù)把信號(hào)里的線性趨勢(shì)消掉再做幅值歸一化。這樣做的原因很實(shí)際VMD內(nèi)部的重構(gòu)約束直接和信號(hào)幅值掛鉤如果信號(hào)幅值跨度太大比如前一段接近0.2而后一段超過(guò)5重建誤差的判據(jù)會(huì)被大量級(jí)片段主導(dǎo)低頻模態(tài)的收斂精度很受影響。采樣點(diǎn)數(shù)方面VMD對(duì)點(diǎn)數(shù)不敏感從幾百點(diǎn)到上百萬(wàn)點(diǎn)都能算但點(diǎn)數(shù)越多迭代每步的FFT開(kāi)銷(xiāo)越大。常見(jiàn)的做法是截取2的整數(shù)次冪長(zhǎng)度附近的一段數(shù)據(jù)比如4096或8192點(diǎn)既保證頻域分辨率足夠又避免零填充帶來(lái)的頻譜泄漏。開(kāi)始正式分解之前還有個(gè)工程習(xí)慣先對(duì)信號(hào)做一次粗略的功率譜估計(jì)標(biāo)出明顯的譜峰位置把譜峰數(shù)量當(dāng)作K值的初始參考。一個(gè)反直覺(jué)的結(jié)論是K設(shè)得比譜峰數(shù)量略大一到兩個(gè)往往效果更好因?yàn)檎鎸?shí)信號(hào)里除了明顯譜峰之外還有邊頻帶和低頻漂移成分它們不是噪聲但它們需要獨(dú)立的模態(tài)來(lái)承載。這套預(yù)處理流程通常能把VMD的參數(shù)搜索范圍縮小一大半。3.3 從程序包輸出里提取模態(tài)時(shí)域波形和瞬時(shí)頻率的重組邏輯VMD程序包輸出最常見(jiàn)的是一組u和一個(gè)omega矩陣。u矩陣的列就是分解出來(lái)的模態(tài)分量時(shí)域波形omega矩陣記錄的是每次迭代的中心頻率軌跡。分析階段該關(guān)注的不是模態(tài)的幅值本身而是模態(tài)的時(shí)域波形和功率譜的耦合關(guān)系。具體做法是先看u矩陣?yán)锼牡轿鍌€(gè)模態(tài)的功率譜是否各自集中在一個(gè)窄帶內(nèi)如果是說(shuō)明分解質(zhì)量高如果某個(gè)模態(tài)的譜峰拖了很長(zhǎng)的尾巴說(shuō)明alpha太大或者K太多一部分頻率成分被劈裂了。瞬時(shí)頻率的提取通常是在分解完成的模態(tài)上再走一遍希爾伯特變換得到解析信號(hào)后取相位導(dǎo)數(shù)。這樣做的目的是觀察非平穩(wěn)工況下中心頻率的漂移情況比如變轉(zhuǎn)速設(shè)備在升速過(guò)程中某個(gè)模態(tài)的中心頻率應(yīng)該是平滑上升的。如果提取出的瞬時(shí)頻率曲線出現(xiàn)明顯的鋸齒狀跳躍不用急著怪VMD程序包多半是模態(tài)之間發(fā)生了邊界混疊K取小了或者alpha取大了先回去調(diào)參比在瞬時(shí)頻率上做平滑更治本。4. 把參數(shù)調(diào)到恰好K值選多少、alpha怎么設(shè)、收斂怎么判4.1 K值三準(zhǔn)則周期圖峰值計(jì)數(shù)、模態(tài)能量準(zhǔn)則和過(guò)分解檢查K值是VMD里最容易翻車(chē)的一個(gè)旋鈕設(shè)少了模態(tài)被合并設(shè)多了模態(tài)被劈碎。我的經(jīng)驗(yàn)是按三層準(zhǔn)則來(lái)做。第一層是周期圖上的峰值數(shù)把功率譜里明顯高于噪聲底的獨(dú)立譜峰數(shù)量出來(lái)作為K值的最低下限。第二層是檢查分解結(jié)果里是否出現(xiàn)「中心頻率接近且波形高度相似的相鄰模態(tài)」如果有說(shuō)明K偏大因?yàn)閂MD不會(huì)自動(dòng)合并兩個(gè)頻率接近的模態(tài)只會(huì)把同一個(gè)物理分量拆成兩個(gè)帶寬更窄的分量。第三層是模態(tài)的歸一化能量占比如果一個(gè)模態(tài)的能量占總信號(hào)能量不到1%且頻譜形態(tài)接近于噪聲包絡(luò)這個(gè)模態(tài)多半是過(guò)分解的產(chǎn)物可以認(rèn)為K多設(shè)了。一個(gè)更容易操作的判定是觀察omega矩陣最后一列的收斂情況如果K取6時(shí)的某幾個(gè)中心頻率分別在迭代結(jié)束時(shí)停在彼此相差不到5 Hz的范圍內(nèi)說(shuō)明這組模態(tài)在頻域上幾乎重疊K值可以降一檔。反過(guò)來(lái)如果K取3時(shí)出現(xiàn)一個(gè)模態(tài)的頻帶同時(shí)覆蓋了兩個(gè)譜峰K值應(yīng)該升到4或5。在實(shí)際項(xiàng)目里我遇到過(guò)振動(dòng)信號(hào)實(shí)際只有兩個(gè)轉(zhuǎn)動(dòng)頻率分量但K取4效果反而更好多出來(lái)的兩個(gè)模態(tài)分別承載了噪聲和沖擊分量讓前兩個(gè)模態(tài)更干凈。4.2 alpha與帶寬的綁定關(guān)系什么時(shí)候用200什么時(shí)候用5000以上alpha參數(shù)的物理意義是懲罰因子它直接控制了信號(hào)的頻域帶寬約束強(qiáng)度。程序包里傳給VMD的alpha如果太小模態(tài)會(huì)為了擬合信號(hào)細(xì)節(jié)而把帶寬拉得很大結(jié)果就是不同模態(tài)在頻域上大量重疊重構(gòu)誤差小但模態(tài)可解釋性差。alpha如果太大模態(tài)被壓成很窄的一條譜線一個(gè)輕微頻率漂移都會(huì)導(dǎo)致模態(tài)分裂。我在運(yùn)行中總結(jié)出一張粗略的參照表采樣率1000 Hz的軸承信號(hào)alpha在500到1500之間效果較好采樣率5000 Hz或更高的聲發(fā)射信號(hào)alpha要升到3000到8000才壓得住邊界。這張表的價(jià)值在于給調(diào)參一個(gè)出發(fā)點(diǎn)而不是死記數(shù)值。實(shí)際調(diào)試時(shí)我有個(gè)習(xí)慣固定K值后把a(bǔ)lpha按對(duì)數(shù)刻度掃三到四次觀察omega矩陣迭代軌跡的穩(wěn)定性。如果某個(gè)alpha區(qū)間內(nèi)omega的收斂軌跡平滑且中心頻率彼此分離清楚這個(gè)區(qū)間就是可用的工作區(qū)間。特別提醒一點(diǎn)程序包如果寫(xiě)成內(nèi)嵌維納濾波形式alpha直接寫(xiě)進(jìn)頻域更新式里這意味著alpha對(duì)結(jié)果的影響可以用Parseval定理?yè)Q算到時(shí)域去理解。調(diào)試順序上永遠(yuǎn)先粗調(diào)alpha再微調(diào)K值因?yàn)閍lpha決定單模態(tài)的質(zhì)量K值決定模態(tài)的條數(shù)前者出錯(cuò)會(huì)讓后者怎么調(diào)都不對(duì)。4.3 收斂容差和迭代次數(shù)tol放寬多少不影響結(jié)果收斂容忍度tol是VMD程序包里容易被忽略的參數(shù)。它的含義是分解模態(tài)重構(gòu)信號(hào)與原始信號(hào)之間的歸一化均方誤差閾值默認(rèn)給1e-7是非常嚴(yán)格的設(shè)置迭代步數(shù)會(huì)明顯偏多。考慮實(shí)際工程信號(hào)本身含有測(cè)量噪聲重構(gòu)誤差收斂到1e-6附近的曲線變化已經(jīng)極其平緩此時(shí)模態(tài)分解結(jié)果幾乎不再變化。我自己在實(shí)測(cè)數(shù)據(jù)里把tol從1e-7放寬到1e-5做過(guò)對(duì)比逐條模態(tài)波形的最大幅值差異不到3%而計(jì)算時(shí)間節(jié)約了40%以上。如果你的場(chǎng)景是大量批量數(shù)據(jù)的離線的分析任務(wù)建議從1e-6起步如果是需要實(shí)時(shí)性或者嵌入式轉(zhuǎn)譯的場(chǎng)景1e-4也勉強(qiáng)可用但前提是模態(tài)數(shù)量準(zhǔn)確和收斂軌跡不發(fā)生漂移。迭代上限maxiter默認(rèn)500一般夠用但我見(jiàn)過(guò)一個(gè)案例某信號(hào)的主頻成分在迭代初期頻繁跳變中心頻率導(dǎo)致前200步都沒(méi)進(jìn)入單調(diào)收斂軌道最終在接近450步才達(dá)到tol。這種情況不用急著加maxiter先觀察omega矩陣前幾十步的軌跡看跳變來(lái)自初始頻率拉拽還是真實(shí)頻率漂移。康德的經(jīng)驗(yàn)是與其改maxiter不如改init1使中心頻率從均勻分布開(kāi)始迭代通常能讓收斂時(shí)間縮短一半。注意tol和maxiter一起改的場(chǎng)景很少在VMD這個(gè)程序包里收斂判據(jù)是硬性誤差閾值不像其他優(yōu)化問(wèn)題還常需要限制最大步長(zhǎng)。4.4 逐一核驗(yàn)?zāi)B(tài)正交性的三個(gè)觀察指標(biāo)分解做完后要花幾分鐘確認(rèn)模態(tài)質(zhì)量有三個(gè)快速觀察指標(biāo)。第一是相關(guān)系數(shù)矩陣計(jì)算各模態(tài)之間的互相關(guān)系數(shù)理想情況下非對(duì)角元絕對(duì)值應(yīng)小于0.1如果某兩個(gè)模態(tài)的相關(guān)系數(shù)超過(guò)0.3幾乎可以斷定發(fā)生了模態(tài)混疊K值或alpha需要重設(shè)。第二是模態(tài)功率譜的3 dB帶寬與中心頻率之比比值超過(guò)0.25甚至0.5說(shuō)明模態(tài)帶寬已經(jīng)寬得失去窄帶信號(hào)的意義。第三是重構(gòu)殘差的頻譜形態(tài)殘差主要由噪聲和未建模成分組成如果殘差能量集中在某個(gè)窄帶以?xún)?nèi)說(shuō)明分解不充分那個(gè)窄帶里還有分量沒(méi)被注意到。這三個(gè)指標(biāo)用起來(lái)也快寫(xiě)個(gè)幾行MATLAB腳本就能在分解后自動(dòng)輸出。相關(guān)的經(jīng)驗(yàn)是相關(guān)系數(shù)矩陣對(duì)alpha的敏感度遠(yuǎn)高于對(duì)K的敏感度如果你遇到模態(tài)混疊但不想動(dòng)K先試試加增大alpha反過(guò)來(lái)如果你覺(jué)得模態(tài)條數(shù)不對(duì)且alpha已經(jīng)較高這時(shí)候減K會(huì)比動(dòng)alpha更能解決問(wèn)題。不同的出發(fā)病根對(duì)應(yīng)不同的旋鈕這個(gè)判斷本身比調(diào)參操作更值錢(qián)。5. VMD程序包避坑記錄模態(tài)分裂、初始化陷阱和邊界效應(yīng)5.1 現(xiàn)象alpha設(shè)到飽和區(qū)間后端某兩個(gè)模態(tài)的譜峰幾乎重疊一個(gè)相對(duì)常見(jiàn)的翻車(chē)現(xiàn)場(chǎng)是模態(tài)譜峰重疊?,F(xiàn)象很直觀用程序包自帶的繪圖函數(shù)畫(huà)出模態(tài)頻譜K4的結(jié)果中有兩個(gè)模態(tài)的中心頻率只相差不到2 Hz譜峰疊加成雙峰樣子時(shí)域波形上這兩個(gè)模態(tài)看起來(lái)像同一段振蕩信號(hào)被鏡像了一下。原因分析上這種情況最常出在alpha取值在飽和區(qū)間后端比如采樣率1000 Hz的信號(hào)把a(bǔ)lpha設(shè)到8000以上帶寬約束過(guò)強(qiáng)導(dǎo)致VMD在迭代中把單個(gè)物理分量在頻域上強(qiáng)制劃分成多個(gè)窄帶分量以更精確地?cái)M合分量?jī)?nèi)部的幅值調(diào)制細(xì)節(jié)。解決路徑是把a(bǔ)lpha往回調(diào)到3000以?xún)?nèi)或者把K值同時(shí)降一檔讓出帶寬給真實(shí)分量。注意單純?cè)黾覭值并不能解決這類(lèi)問(wèn)題反而會(huì)加劇模態(tài)分裂。另一種導(dǎo)致譜峰重疊的原因是信號(hào)的幅值調(diào)制深度太高調(diào)幅效應(yīng)會(huì)產(chǎn)生較寬的邊帶VMD會(huì)傾向于把它們拆成多個(gè)窄帶模態(tài)看起來(lái)像同一個(gè)物理來(lái)源被劈成了兩條譜。這種情況最有效的處理路徑是在VMD分解之前對(duì)信號(hào)做包絡(luò)解調(diào)讓分解算法的輸入從帶邊帶的載波信號(hào)變成低頻包絡(luò)信號(hào)分解復(fù)雜度直接從四五個(gè)模態(tài)降成兩個(gè)。5.2 現(xiàn)象信號(hào)長(zhǎng)度只有幾百個(gè)采樣點(diǎn)分解結(jié)果在數(shù)據(jù)兩端出現(xiàn)劇烈振蕩短數(shù)據(jù)的邊界效應(yīng)殺傷力很大。我在處理一段只有三百個(gè)點(diǎn)的瞬態(tài)振動(dòng)數(shù)據(jù)時(shí)VMD分解出的模態(tài)在兩端出現(xiàn)明顯比中間幅值大幾倍的尖峰看起來(lái)非常可疑但中心頻率又合理。原因出在VMD內(nèi)部的維納濾波更新涉及頻域處理信號(hào)在時(shí)間域兩端經(jīng)FFT隱含的周期延拓造成了不連續(xù)的跳變這種邊界效應(yīng)在數(shù)據(jù)太短時(shí)會(huì)被放大到誤導(dǎo)模態(tài)形態(tài)。解決這個(gè)問(wèn)題的標(biāo)準(zhǔn)做法是信號(hào)延拓常見(jiàn)的做法是在原始信號(hào)兩端分別拼接一段鏡像或者線性預(yù)測(cè)延拓的數(shù)據(jù)分解完成后再把延拓段裁掉。在MATLAB程序包里我一般會(huì)構(gòu)造一個(gè)帶welch窗的平滑截?cái)?。具體而言先取出原始數(shù)據(jù)做均值歸一化然后在其前后各拼接20%長(zhǎng)度的鏡像數(shù)據(jù)再調(diào)VMD之后把輸出的u矩陣裁回原始長(zhǎng)度。另一種可接受的思路是假設(shè)信號(hào)本身是一個(gè)周期過(guò)程的片段直接給VMD輸入整周期的數(shù)據(jù)做整周期采樣來(lái)保證FFT假設(shè)成立。實(shí)測(cè)經(jīng)驗(yàn)里信號(hào)長(zhǎng)度不低于300個(gè)點(diǎn)且端點(diǎn)幅值接近零時(shí)邊界效應(yīng)基本可以接受不用處理但如果現(xiàn)有數(shù)據(jù)長(zhǎng)度到不了這個(gè)量級(jí)就必須做延拓。5.3 現(xiàn)象程序包里某輔助函數(shù)報(bào)錯(cuò)說(shuō)矩陣維度不一致這個(gè)問(wèn)題通常不是VMD分解算法的邏輯錯(cuò)誤而是程序包接口和用戶(hù)數(shù)據(jù)的維度約定不一致。最常見(jiàn)場(chǎng)景是把一個(gè)行向量信號(hào)傳給VMD程序內(nèi)部卻在頻域初始化里按列向量處理導(dǎo)致維度廣播失敗。原因在于很多從公開(kāi)渠道流傳的VMD程序包默認(rèn)輸入是單精度列向量你用load命令讀進(jìn)來(lái)的結(jié)構(gòu)體字段可能是個(gè)帶時(shí)間的多列矩陣。解決這個(gè)問(wèn)題很簡(jiǎn)單進(jìn)入調(diào)用處把信號(hào)強(qiáng)制轉(zhuǎn)為double類(lèi)型的列向量即可。我通常會(huì)在主腳本開(kāi)頭加一行signal double(signal(:))這行操作雖然不起眼卻能消滅大部分接口報(bào)錯(cuò)。另一類(lèi)和矩陣維度有關(guān)的報(bào)錯(cuò)出現(xiàn)在u_hat的初始化上。如果程序包按信號(hào)長(zhǎng)度假設(shè)FFT點(diǎn)數(shù)為偶數(shù)過(guò)短的奇數(shù)長(zhǎng)度信號(hào)會(huì)導(dǎo)致頻域數(shù)組賦值尺寸不匹配。遇到這種報(bào)錯(cuò)時(shí)不要試圖去改程序包的內(nèi)部數(shù)組分配代碼而是先把信號(hào)補(bǔ)零到下一個(gè)2的冪長(zhǎng)度。VMD對(duì)補(bǔ)零的容忍度比一般頻譜分析高因?yàn)檠a(bǔ)零的主要作用是調(diào)整FFT分辨率模態(tài)中心頻率收斂位置基本不受影響。5.4 現(xiàn)象分解出的第一個(gè)模態(tài)總是跟著直流分量走但DC參數(shù)明明設(shè)了0另一種容易遇到的坑是DC參數(shù)設(shè)了0但第一個(gè)模態(tài)仍然明顯包含低頻趨勢(shì)。原因比較隱蔽程序包里的DC參數(shù)控制的是第一個(gè)模態(tài)的中心頻率是否固定為0 Hz但如果信號(hào)本身含有較強(qiáng)的低頻漂移比如采樣前放大器沒(méi)有完全隔直即使DC設(shè)0VMD也會(huì)在迭代中讓某個(gè)模態(tài)自動(dòng)漂向低頻區(qū)把這部分能量吸收進(jìn)去。這時(shí)候應(yīng)該做的是預(yù)處理時(shí)徹底去掉趨勢(shì)項(xiàng)而不是去調(diào)整DC參數(shù)。我的判斷路徑是先對(duì)信號(hào)做一次平滑濾波求出趨勢(shì)基線減去后再做分解看第一個(gè)模態(tài)的頻譜是否回到目標(biāo)頻率附近。還有一種少見(jiàn)但真實(shí)的情況是程序包里DC參數(shù)寫(xiě)法是判等字符而非邏輯值你傳入0表示「不固定直流」但某個(gè)版本寫(xiě)法用if DC 0來(lái)判斷導(dǎo)致行為被顛倒。遇到這種問(wèn)題時(shí)最快的驗(yàn)證方式是跑一段帶直流偏移的方法看第一個(gè)模態(tài)的頻譜集中在0附近還是集中在目標(biāo)頻率附近就知道程序包的DC語(yǔ)義是否如你預(yù)期。5.5 現(xiàn)象程序運(yùn)行時(shí)間幾十秒起步比想象中慢得多一個(gè)讓人抓狂的問(wèn)題是明明信號(hào)只有幾千個(gè)點(diǎn)分解K3卻跑了快一分鐘。這既不是死循環(huán)也不是算法被寫(xiě)壞而是某些VMD程序包內(nèi)部對(duì)每個(gè)模態(tài)都做了逐點(diǎn)循環(huán)的頻域更新。比如在幾百個(gè)頻域點(diǎn)里套了reshape和循環(huán)代碼寫(xiě)法笨重但結(jié)果算法功能沒(méi)問(wèn)題。我處理這種慢速程序包的經(jīng)驗(yàn)很簡(jiǎn)單先跑一次profiler或者用tic/toc看時(shí)間都消耗在哪如果瓶頸在主循環(huán)里每個(gè)模態(tài)的頻域更新可以把頻域網(wǎng)格計(jì)算提到循環(huán)外面做一次省掉重復(fù)計(jì)算。另一個(gè)顯著的加速點(diǎn)是矩陣化替代循環(huán)把頻域更新式的分母計(jì)算向量化K個(gè)模態(tài)一次算完。改完之后這段程序包也能用實(shí)測(cè)速度通常能提高五到十倍。5.6 現(xiàn)象同一段信號(hào)跑兩次得到不完全相同的模態(tài)結(jié)果這不算玄學(xué)但要說(shuō)明原因。某些VMD程序包內(nèi)部有隨機(jī)數(shù)參與初始化實(shí)際就是init0時(shí)用rand初始化中心頻率初值或算法里含有一個(gè)隨機(jī)的擾動(dòng)項(xiàng)。如果兩次分解只在微小數(shù)值上有差異而譜峰位置一致可以認(rèn)為結(jié)果是穩(wěn)定的如果模態(tài)分裂、合并這些宏觀特征都變了那就是程序包的隨機(jī)因素在起作用。解決方式是把中心頻率初始化改成確定性方案比如按等分頻帶初始化。保留一個(gè)可復(fù)現(xiàn)結(jié)果的隨機(jī)種子也是工程上常用的做法。實(shí)際中我見(jiàn)到有人因?yàn)橄嗤盘?hào)兩次分解結(jié)果不一致而放棄VMD轉(zhuǎn)向別的算法這是應(yīng)急的過(guò)度反應(yīng)——絕大多數(shù)情況下這只是初始化策略問(wèn)題與VMD方法本身無(wú)關(guān)。6. VMD的進(jìn)階玩法把分解結(jié)果用到特征工程里去VMD程序包用順手之后真正的價(jià)值在于把分解出的模態(tài)轉(zhuǎn)化成工程上可用的特征。一個(gè)簡(jiǎn)單有效的做法是把每個(gè)模態(tài)的瞬時(shí)幅值包絡(luò)做頻譜分析得到包絡(luò)譜再把幾個(gè)主要包絡(luò)譜峰頻率和幅值拼成一個(gè)特征向量交給分類(lèi)器做故障模式識(shí)別。這在程序包里實(shí)現(xiàn)起來(lái)不難對(duì)每個(gè)模態(tài)先經(jīng)過(guò)希爾伯特變換取包絡(luò)再做FFT最后提取顯著譜峰整個(gè)流程可以寫(xiě)成一個(gè)循環(huán)批量處理大批量歷史數(shù)據(jù)。我做過(guò)的模擬項(xiàng)目里這種特征比直接對(duì)原始信號(hào)提取時(shí)域統(tǒng)計(jì)量在分類(lèi)準(zhǔn)確率上穩(wěn)定提升百分之15到20尤其是在故障特征微弱、被強(qiáng)背景噪聲掩蓋的場(chǎng)景下。另一個(gè)值得嘗試的方向是模態(tài)篩選重建。分解得到的K個(gè)模態(tài)中有些是故障信息載體有些是噪聲選哪些分量用來(lái)重建取決于應(yīng)用需求。經(jīng)驗(yàn)性的做法是計(jì)算每個(gè)模態(tài)與原始信號(hào)的相關(guān)系數(shù)保留相關(guān)系數(shù)較高的模態(tài)或者更精細(xì)地用每個(gè)模態(tài)的時(shí)域包絡(luò)和某個(gè)已知故障特征頻率做相關(guān)分析來(lái)確定去留。這樣重建的信號(hào)去噪效果比硬閾值或小波閾值更干凈因?yàn)槟闶窃诒A粑锢砜山忉尩姆至慷皇前捶登幸坏?。附帶的價(jià)值是你在做特征篩選時(shí)能解釋每個(gè)特征來(lái)自哪個(gè)頻帶、哪個(gè)模態(tài)而不只是依賴(lài)黑盒模型。如果處理的是非平穩(wěn)變轉(zhuǎn)速信號(hào)建議把VMD和階次跟蹤配合使用先通過(guò)轉(zhuǎn)速計(jì)脈沖算出參考軸轉(zhuǎn)頻再對(duì)每個(gè)模態(tài)做角域重采樣。實(shí)測(cè)效果比直接對(duì)原始信號(hào)做階次跟蹤更好因?yàn)槟B(tài)已經(jīng)剔除了噪聲和干擾頻率重采樣得到的階次譜更干凈。這類(lèi)流程看起來(lái)步驟多但每一步都是現(xiàn)成函數(shù)組合起來(lái)不復(fù)雜。這里也有個(gè)反面經(jīng)驗(yàn)有個(gè)同事直接對(duì)模態(tài)波形做常規(guī)FFT忽略了變轉(zhuǎn)速帶來(lái)的頻率模糊特征頻率被抹成了一個(gè)寬峰白費(fèi)了VMD分離的效果后來(lái)重寫(xiě)成角域重采樣才把故障特征拉出來(lái)。收個(gè)尾。我在VMD上踩過(guò)最深的坑就是把程序包當(dāng)黑匣子參數(shù)只抄默認(rèn)值分解出來(lái)一堆模態(tài)卻搞不清每條成分對(duì)應(yīng)什么物理量最后被領(lǐng)導(dǎo)一問(wèn)就露怯了?,F(xiàn)在每次跑VMD都先花五分鐘看頻譜、估計(jì)模態(tài)數(shù)、設(shè)一組初始參數(shù)再按本文的驗(yàn)證流程逐個(gè)檢核心里才踏實(shí)。希望這篇筆記能幫你在自己的信號(hào)數(shù)據(jù)里少走幾步彎路用最少的時(shí)間把變分模態(tài)分解真正用起來(lái)。本文還有配套的精品資源點(diǎn)擊獲取