現(xiàn)多輸入時(shí)序預(yù)測全流程)
1. 項(xiàng)目背景與整體思路這兩年做時(shí)序預(yù)測的人越來越多但大多數(shù)人一開始都是從單一模型入手的——拿個(gè)歷史序列直接扔進(jìn)LSTM、BP或者GRU里出一組誤差指標(biāo)就算交差了。等真正接觸到工業(yè)現(xiàn)場的實(shí)測數(shù)據(jù)才發(fā)現(xiàn)理想很豐滿現(xiàn)實(shí)很骨感現(xiàn)場數(shù)據(jù)幾乎都是非平穩(wěn)的噪聲多、模態(tài)混疊、特征維度高、變量之間還有復(fù)雜的非線性耦合。這時(shí)候再去套單一神經(jīng)網(wǎng)絡(luò)模型不是不收斂就是收斂到一半訓(xùn)練崩了精度更是慘不忍睹。我這次做的這個(gè)項(xiàng)目就是把三條路擰成一股繩來解決上述問題VMD變分模態(tài)分解 KPCA核主成分分析 PINN物理信息神經(jīng)網(wǎng)絡(luò)全部用MATLAB實(shí)現(xiàn)目標(biāo)是很常見的多輸入單輸出時(shí)序預(yù)測任務(wù)。所謂多輸入單輸出就是模型接收多個(gè)相關(guān)變量的歷史觀測值只預(yù)測未來一個(gè)關(guān)鍵變量的值。這種結(jié)構(gòu)在企業(yè)里的需求量極大比如用環(huán)境溫度、濕度、風(fēng)速等多個(gè)氣象變量去預(yù)測光伏功率或者用轉(zhuǎn)速、負(fù)載電流、軸承溫度去預(yù)測設(shè)備剩余壽命本質(zhì)上都是同一個(gè)問題。整套方案的核心思路不復(fù)雜VMD先把原始信號拆解成若干個(gè)不同頻段的子序列相當(dāng)于把一個(gè)糾纏不清的混合信號分層歸位KPCA再把這些子序列連同其他外部變量一起做非線性特征壓縮去掉冗余、留下最有解釋力的成分最后交給PINN做回歸預(yù)測。PINN和普通的神經(jīng)網(wǎng)絡(luò)不太一樣它在損失函數(shù)里嵌入了物理方程約束預(yù)測曲線不容易跑偏特別適合樣本量小但規(guī)律性強(qiáng)的工業(yè)場景。這套組合拳最適合三類人一是做風(fēng)電、光伏、電力負(fù)荷預(yù)測的工程師天天被非平穩(wěn)數(shù)據(jù)折磨二是做故障診斷與剩余壽命預(yù)測的研究生手里數(shù)據(jù)源多但不知道怎么做特征工程三是對MATLAB熟悉、想在圖神經(jīng)網(wǎng)絡(luò)之外找一種更穩(wěn)妥預(yù)測路線的從業(yè)者。看完這篇文章你至少能搞清楚三件事這三個(gè)算法為什么要按這個(gè)順序串聯(lián)、每一步在MATLAB里具體怎么落地、以及實(shí)際調(diào)試中那些文檔里不會(huì)寫但能救命的坑。2. 為什么一定要VMD-KPCA-PINN這套組合2.1 單模型處理不了非平穩(wěn)信號先把最基礎(chǔ)的問題聊透為什么不用純神經(jīng)網(wǎng)絡(luò)直接干你看看典型的工業(yè)實(shí)測數(shù)據(jù)就明白了。一條風(fēng)機(jī)功率曲線可能同時(shí)包含慢變的趨勢分量、因?yàn)轱L(fēng)速突變產(chǎn)生的間歇毛刺、還有齒輪箱周期性振動(dòng)帶來的高頻噪聲。神經(jīng)網(wǎng)絡(luò)本身其實(shí)不具備自動(dòng)分離頻段的能力它更像一個(gè)把所有輸入特征亂燉的廚師——如果輸入里有不同時(shí)間尺度的模式混雜在一起模型會(huì)試圖用同一組參數(shù)去擬合所有模式結(jié)果就是什么都擬合不好。這和你在MATLAB里用單一LSTM做風(fēng)電功率預(yù)測最終誤差指標(biāo)總在某個(gè)閾值下不去的現(xiàn)象是吻合的。VMD就是來解決這個(gè)問題的。它能把原始信號f(t)分解成K個(gè)有限帶寬的固有模態(tài)函數(shù)u_k(t)每個(gè)模態(tài)都有一個(gè)自己的中心頻率ω_k并且這些模態(tài)的帶寬之和盡可能小。數(shù)學(xué)上說它是在求解一個(gè)變分約束問題[ \min_{{u_k},{\omega_k}} \sum_k \left| \partial_t \left[ (\delta(t)\frac{j}{\pi t}) * u_k(t) \right] e^{-j\omega_k t} \right|_2^2 ]你不需要逐字推導(dǎo)這個(gè)公式只要抓住關(guān)鍵就行VMD輸出的K條子序列在頻域上是分開的每條子序列保持了自己的時(shí)間域特征但又沒有了混疊帶來的干擾。分解完之后高頻噪聲進(jìn)了某個(gè)IMF慢變趨勢進(jìn)了另一個(gè)IMF模型就能真正看到干凈的結(jié)構(gòu)。2.2 KPCA比PCA更合適高維非線性特征VMD把序列分解之后輸入維度不降反升。假設(shè)原始序列有5個(gè)變量每個(gè)變量分解出4個(gè)模態(tài)再加上溫度和風(fēng)速之類的外部變量輸入空間立刻膨脹到20多維度。這些維度之間還有很強(qiáng)的相關(guān)性比如相鄰模態(tài)之間可能存在信息冗余。直接全塞給神經(jīng)網(wǎng)絡(luò)倒不是不行但訓(xùn)練效率會(huì)肉眼可見地下降表現(xiàn)在訓(xùn)練曲線震蕩劇烈甚至出現(xiàn)過擬合。這時(shí)候傳統(tǒng)PCA能幫上忙嗎能但不夠。PCA本質(zhì)上是在原始空間里尋找最大方差方向的線性投影對付線性相關(guān)的特征很有效一旦特征之間有非線性耦合關(guān)系就抓瞎了。比如風(fēng)速和光伏功率在低風(fēng)速區(qū)接近線性、在高風(fēng)速區(qū)卻進(jìn)入飽和區(qū)這種非線性關(guān)系PCA是看不出來的。KPCA的做法是先通過核函數(shù)將原始數(shù)據(jù)映射到高維特征空間再在這個(gè)高維空間里做PCA。簡單理解就是把非線性問題通過一個(gè)核技巧拉直變成線性問題再降維。這樣做的好處是既保留了非線性結(jié)構(gòu)又拿到了降維后的低維表示。這個(gè)先拉直、再降維的思路和你在坐標(biāo)系里處理彎曲數(shù)據(jù)時(shí)先做變換再用直線擬合是同一個(gè)道理。2.3 PINN讓預(yù)測結(jié)果不脫離物理規(guī)律代碼跑到最后一步模型如果只是個(gè)普通全連接網(wǎng)絡(luò)或者LSTM其實(shí)前面的VMD和KPCA已經(jīng)把特征工程做得相當(dāng)好了精度也不會(huì)差。但PINN的價(jià)值在于另一個(gè)層面——它能把你對物理系統(tǒng)的先驗(yàn)知識直接寫進(jìn)網(wǎng)絡(luò)訓(xùn)練過程。舉個(gè)具體的例子。做儲(chǔ)能電池SOC預(yù)測時(shí)你很清楚SOC的變化率不能無限大它受最大充放電倍率約束SOC的取值一定落在0到1之間不可能出現(xiàn)1.3這種荒謬結(jié)果。這些規(guī)律是物理硬約束普通神經(jīng)網(wǎng)絡(luò)是不知道的它只會(huì)根據(jù)訓(xùn)練集的統(tǒng)計(jì)慣性來輸出。如果測試集的工況發(fā)生了偏移純數(shù)據(jù)模型的預(yù)測就會(huì)飛出物理邊界這在工程上是不可接受的。PINN的解決方式非常直接把物理方程的殘差作為一個(gè)懲罰項(xiàng)加進(jìn)損失函數(shù)。設(shè)網(wǎng)絡(luò)預(yù)測值為(\hat{y})真實(shí)值為(y)物理約束項(xiàng)的殘差為(R(\hat{y}))那么總損失就是[ Loss \frac{1}{N}\sum_{i1}^{N}(\hat{y}i - y_i)^2 \lambda \cdot \frac{1}{N}\sum{i1}^{N}R^2(\hat{y}_i) ](\lambda)是物理約束的強(qiáng)度系數(shù)(R(\hat{y}))可以是一階導(dǎo)數(shù)約束(d\hat{y}/dt - f(\hat{u}))也可以是上下邊界約束(\max(0, \hat{y}-1) \max(0, -\hat{y}))。數(shù)學(xué)形式不固定完全依賴你對系統(tǒng)規(guī)律的理解深度。訓(xùn)練過程中網(wǎng)絡(luò)不僅要讓預(yù)測值逼近真實(shí)值還要讓自己輸出的解盡可能滿足物理方程這相當(dāng)于給模型上了雙保險(xiǎn)。3. 環(huán)境準(zhǔn)備與數(shù)據(jù)鏈路設(shè)計(jì)3.1 工具箱選擇和版本適配問題實(shí)話說這套流程對環(huán)境的要求不算苛刻。VMD在MATLAB里沒有官方內(nèi)置函數(shù)但File Exchange上有原作者的VMDFunction實(shí)現(xiàn)下載后放到搜索路徑即可調(diào)用KPCA如果沒有統(tǒng)計(jì)機(jī)器學(xué)習(xí)的工具箱許可用File Exchange上的kPCA代碼也能跑PINN需要Deep Learning Toolbox建議2021a以上的版本因?yàn)閺倪@開始自定義訓(xùn)練循環(huán)的語法才比較全。我在實(shí)際部署時(shí)用的是MATLAB R2023bDLL和C編譯器都適配得很好沒有遇到舊版本里dlnetwork與自定義損失函數(shù)之間接口不一致的問題。如果你是2020a或者更早的版本建議先把工具箱升級到2021a之后再動(dòng)手否則在自定義訓(xùn)練循環(huán)那一步會(huì)卡很久。新版MATLAB對GPU訓(xùn)練的支持也更友好在訓(xùn)練PINN的時(shí)候只有一個(gè)GPU和沒有GPU完全是兩個(gè)訓(xùn)練速度強(qiáng)烈建議有條件就把支持CUDA的顯卡開起來。這里要特別提醒一句不要讓所有算法模塊在同一個(gè)腳本里一鍋燉。VMD分解和KPCA降維這兩步建議拆成獨(dú)立的函數(shù)文件單獨(dú)運(yùn)行、保存中間結(jié)果到.mat文件。因?yàn)閂MD求解是個(gè)迭代過程一旦K值選得不合適重新分解很費(fèi)時(shí)間把中間結(jié)果緩存下來后續(xù)調(diào)PINN模型參數(shù)時(shí)不需要重新跑前兩步能省下一晚上調(diào)試時(shí)間。3.2 多輸入單輸出問題的數(shù)據(jù)組織方式多輸入單輸出建模的第一步是把原始的三維世界變成二維矩陣。以我的實(shí)驗(yàn)數(shù)據(jù)為例原始數(shù)據(jù)是一個(gè)3065×7的矩陣7列分別是環(huán)境溫度、風(fēng)速、濕度、氣壓、光照強(qiáng)度、上一時(shí)刻功率和當(dāng)前時(shí)刻功率。目標(biāo)是預(yù)測當(dāng)前時(shí)刻功率所以真正參與建模的輸入是前6列輸出是第7列。這里有個(gè)時(shí)序建模的基本操作構(gòu)造滑窗樣本。假設(shè)滑窗長度lookback5步長step1那么第1個(gè)樣本的輸入是第1到5行的6個(gè)變量共30個(gè)特征輸出是第6行的功率值第2個(gè)樣本輸入是第2到6行的6個(gè)變量輸出是第7行的功率值。最后得到一個(gè)(n_samples, 30)的輸入矩陣和一個(gè)(n_samples, 1)的輸出矩陣。為什么滑窗長度取5而不是直接取全部歷史這背后是自相關(guān)分析的邏輯我先把功率序列做了自相關(guān)函數(shù)圖發(fā)現(xiàn)滯后5步之后自相關(guān)系數(shù)掉到0.3以下說明超過5個(gè)時(shí)刻的歷史信息對當(dāng)前時(shí)刻的貢獻(xiàn)已經(jīng)弱了。取長了模型參數(shù)增多訓(xùn)練時(shí)間變長精度提升極少取短了信息不夠誤差明顯上抬。你可以用MATLAB的autocorr函數(shù)跑一遍自己的數(shù)據(jù)按0.3這個(gè)閾值確定滑窗長度。3.3 訓(xùn)練集測試集劃分不能亂切時(shí)序預(yù)測里的訓(xùn)練集測試集劃分和常規(guī)機(jī)器學(xué)習(xí)里的隨機(jī)劃分完全是兩回事。隨機(jī)劃分會(huì)直接導(dǎo)致數(shù)據(jù)泄漏訓(xùn)練集里包含測試集時(shí)刻之后的未來數(shù)據(jù)模型在訓(xùn)練時(shí)偷看到了未來驗(yàn)證誤差會(huì)異常好看一到真正的未來時(shí)刻就原形畢露。正確的做法是按時(shí)間順序切分。我的實(shí)驗(yàn)里取了前80%的樣本作為訓(xùn)練集后20%作為測試集中間不留交叉。歸一化的參數(shù)、KPCA的投影矩陣、VMD的分解邊界全部只從訓(xùn)練集上估計(jì)測試集的數(shù)據(jù)要用訓(xùn)練集算好的參數(shù)轉(zhuǎn)換千萬不能在整個(gè)數(shù)據(jù)集上先歸一化再去劃分否則訓(xùn)練集和測試集的信息就互相污染了。4. VMD分解的MATLAB實(shí)現(xiàn)與K值確定4.1 VMD核心調(diào)用參數(shù)詳解在MATLAB里調(diào)用VMD的核心代碼很簡單但每個(gè)參數(shù)背后的含義值得花點(diǎn)時(shí)間說清楚。常用的調(diào)用格式是[u, u_hat, omega] VMD(f, alpha, tau, K, DC, init, tol);其中f是被分解的一維信號alpha是懲罰因子tau是噪聲容忍度K是模態(tài)數(shù)量DC表示是否保留直流分量init是初始化方式tol是收斂容差。輸出u是K行N列的矩陣每一行代表一個(gè)IMF子序列omega是各模態(tài)的中心頻率。alpha這個(gè)參數(shù)最容易讓人困惑。它控制的是模態(tài)帶寬的懲罰強(qiáng)度alpha越大每個(gè)模態(tài)的帶寬越窄各模態(tài)的中心頻率分得越開分解結(jié)果越純凈但也容易出現(xiàn)模態(tài)丟失alpha越小模態(tài)帶寬越寬模態(tài)之間可能發(fā)生重疊。經(jīng)驗(yàn)值在2000附近但這不是死規(guī)矩。我的做法是固定K不變把a(bǔ)lpha從500拉到5000觀察各模態(tài)中心頻率的分布找到一個(gè)既不重疊又不丟失分量的值。tau在信號含噪比較嚴(yán)重時(shí)推薦取0此時(shí)VMD等效于強(qiáng)制噪聲抑制如果信號比較干凈取默認(rèn)值0.25就行。DC通常在信號有直流偏置時(shí)設(shè)為1否則設(shè)為0。init推薦用1表示全部初始化為均勻頻率分布這對大多數(shù)信號來說是最穩(wěn)妥的開始方式。4.2 K值選大了還是選小了怎么看K值模態(tài)數(shù)是整個(gè)VMD環(huán)節(jié)里最需要經(jīng)驗(yàn)的地方。K選小了分解不充分兩個(gè)不同頻率的成分?jǐn)D在一個(gè)模態(tài)里模態(tài)混疊的問題依然存在K選大了出現(xiàn)過度分解出現(xiàn)虛假模態(tài)——也就是某個(gè)模態(tài)和另一個(gè)模態(tài)在頻域上幾乎重合物理意義不清晰。先說我試過且有效的兩種判定方法。第一種是用VMD輸出的中心頻率omega直接判斷在同一alpha下如果中心頻率從小到大排列之后末尾兩個(gè)頻率的差值小于前面相鄰差值的一半十有八九是K取大了。第二種是觀察各模態(tài)與原始信號的相關(guān)系數(shù)正常情況下各模態(tài)與原始信號的相關(guān)性應(yīng)該從高頻到低頻遞減如果突然跳出來一個(gè)模態(tài)相關(guān)性異常低或者和某個(gè)鄰近模態(tài)相關(guān)性極高說明過分解了。實(shí)際操作起來我建議直接跑K3到K8的多組分解每組用同一alpha畫出各模態(tài)的中心頻率折線圖。要的就是中心頻率從低到高分布均勻、末尾不扎堆、無空洞的那組K。我在實(shí)測中分解風(fēng)機(jī)功率數(shù)據(jù)K5時(shí)有明顯的中心頻率扎堆現(xiàn)象K4就非常均勻最后選K4預(yù)測效果提升了大約10%的RMSE。4.3 VMD分解后的數(shù)據(jù)重組VMD分解完千萬別直接把每條IMF單獨(dú)丟給模型。正確的做法是把分解后的子序列作為新的特征列和原始的其他變量拼在一起構(gòu)成一個(gè)擴(kuò)維矩陣。我這里的做法是把原始功率信號分解成4條IMF將這4條信號按列拼到7個(gè)原始變量之后形成11列的中間矩陣。然后進(jìn)入KPCA降維。這個(gè)重組過程在MATLAB里就是一條矩陣拼接命令的事IMF_all u; % u是VMD輸出的K行N列矩陣 features_combined [X_original, IMF_all]; % 拼接成新特征矩陣需要注意的是每條IMF在送入KPCA之前要不要再做一次差分、取對數(shù)之類的變換取決于數(shù)據(jù)的實(shí)際特性。如果信號本身波動(dòng)劇烈可以先做對數(shù)變化讓幅度更平穩(wěn)但別做過頭否則KPCA提取的非線性主成分可能只反映了變換噪聲。5. KPCA降維的實(shí)戰(zhàn)細(xì)節(jié)5.1 核函數(shù)選擇與參數(shù)標(biāo)定KPCA和PCA在原理上最大的區(qū)別在于核函數(shù)的選擇。我試過線性核、多項(xiàng)式核和高斯核RBF核從預(yù)測結(jié)果來看高斯核的效果最穩(wěn)。原因也好理解高斯核對特征空間的非線性映射表達(dá)力最強(qiáng)能夠捕捉風(fēng)力功率曲線中的飽和特性而線性核和多項(xiàng)式核在這個(gè)場景里映射能力有限。高斯核的寬度參數(shù)sigma是唯一的自由參數(shù)。sigma太小核矩陣過于尖銳自動(dòng)退化成近鄰匹配降維后的特征基本就是原始特征的翻版sigma太大核矩陣趨近于均勻拉不開點(diǎn)與點(diǎn)之間的差異降維效果等同于PCA。標(biāo)定sigma的實(shí)用辦法是用交叉驗(yàn)證在訓(xùn)練集內(nèi)部再切一個(gè)驗(yàn)證子集對sigma做對數(shù)網(wǎng)格搜索取驗(yàn)證集上重構(gòu)誤差最小的sigma。MATLAB里用File Exchange上的kPCA函數(shù)時(shí)需要先自己構(gòu)造核矩陣并中心化。核矩陣中心化這一步特別容易漏但漏了的話結(jié)果完全是錯(cuò)的。原理是PCA的特征分解必須作用于數(shù)據(jù)中心化后的協(xié)方差矩陣KPCA就必須在特征空間里做同樣的中心化而這個(gè)中心化操作反映在核矩陣上是K_c K - mean(K, 2) - mean(K, 1) mean(K, :);5.2 累計(jì)貢獻(xiàn)率怎么定才合理處理核矩陣特征值分解后每個(gè)特征值對應(yīng)的特征向量就是一個(gè)主方向。主方向該保留多少個(gè)取決于累計(jì)貢獻(xiàn)率也就是前r個(gè)特征值之和占總特征值之和的比例。通用經(jīng)驗(yàn)是累計(jì)貢獻(xiàn)率達(dá)到95%但這是常規(guī)PCA的經(jīng)驗(yàn)KPCA加核變換后特征值衰減通常更快我在實(shí)際運(yùn)行中發(fā)現(xiàn)累計(jì)貢獻(xiàn)率到90%左右預(yù)測精度就已經(jīng)穩(wěn)定了再往上加維度只會(huì)增加無效計(jì)算。還有一點(diǎn)值得單獨(dú)說KPCA降維描述的是整個(gè)特征集的綜合投影但降維后的特征已經(jīng)不具備物理含義了你沒法說第3個(gè)主成分代表風(fēng)速。所以在工程解釋性要求高的場景里建議降維后的維度控制在5到8之間既保留了統(tǒng)計(jì)信息又不至于讓后續(xù)分析完全喪失可解釋線索。5.3 投影系數(shù)一致性訓(xùn)練集和測試集必須共用一個(gè)投影矩陣這一條是KPCA實(shí)操里最容易踩的坑也是很多網(wǎng)上代碼里出錯(cuò)的地方。KPCA的投影方向和主成分系數(shù)只應(yīng)該從訓(xùn)練集上計(jì)算測試集數(shù)據(jù)要做的是套用訓(xùn)練集已經(jīng)算好的投影系數(shù)來獲得低維表示而不是在測試集上重新做一次KPCA。我在代碼里是這樣組織的先對訓(xùn)練集構(gòu)建核矩陣、中心化、特征分解存下需要的投影矩陣和主成分方向測試集數(shù)據(jù)進(jìn)入時(shí)先計(jì)算測試集與訓(xùn)練集樣本之間的核矩陣再做同樣的中心化最終乘以訓(xùn)練集的特征向量矩陣得到測試集的低維表示。這么做才能保證訓(xùn)練集和測試集處于同一特征空間。6. PINN網(wǎng)絡(luò)的搭建與訓(xùn)練細(xì)節(jié)6.1 網(wǎng)絡(luò)結(jié)構(gòu)設(shè)計(jì)不是越深越好PINN的骨干網(wǎng)絡(luò)我采用的是典型的三層全連接結(jié)構(gòu)輸入層接收KPCA降維后的特征序列比如輸出維度設(shè)置為8的窗口特征中間層用兩個(gè)全連接層每層32個(gè)神經(jīng)元激活函數(shù)用tanh輸出層一個(gè)神經(jīng)元對應(yīng)預(yù)測值。整體參數(shù)數(shù)量大約在1300個(gè)左右。為什么不用LSTM或GRU這種循環(huán)結(jié)構(gòu)不是不好而是這個(gè)場景里輸入序列已經(jīng)在VMD和KPCA階段被做了充分的時(shí)頻展開和信息壓縮剩余的時(shí)間動(dòng)態(tài)并不復(fù)雜用循環(huán)網(wǎng)絡(luò)的訓(xùn)練成本高調(diào)參門檻也高。全連接網(wǎng)絡(luò)在這一步反而簡單、穩(wěn)定、收斂快。當(dāng)然如果你把VMD-KPCA的輸出直接接一個(gè)LSTM我也試過最終結(jié)論是精度提升不到2%訓(xùn)練時(shí)間卻翻了4倍性價(jià)比太低。tanh激活函數(shù)在這里優(yōu)于ReLU。因?yàn)槲锢硐到y(tǒng)的輸出范圍通常是有界的比如功率值不可能為負(fù)tanh自帶對稱性且導(dǎo)數(shù)連續(xù)對包含物理約束的損失函數(shù)優(yōu)化起來更友好。ReLU在0點(diǎn)不可導(dǎo)放到PINN的梯度計(jì)算中容易在某些迭代步突然爆炸。6.2 自定義訓(xùn)練循環(huán)PINN的核心代碼MATLAB R2021a之后可以用dlnetwork和自定義訓(xùn)練循環(huán)實(shí)現(xiàn)任意損失函數(shù)。PINN的獨(dú)特之處在于損失函數(shù)里要同時(shí)包含數(shù)據(jù)項(xiàng)和物理項(xiàng)這需要我們在每次迭代時(shí)手動(dòng)計(jì)算梯度。核心代碼框架如下% 網(wǎng)絡(luò)輸入和輸出使用dlarray X_dl dlarray(X_train, CB); Y_dl dlarray(Y_train, CB); for epoch 1:maxEpochs % 前向傳播 Y_pred forward(net, X_dl); % 數(shù)據(jù)損失 lossData mse(Y_pred, Y_dl); % 物理損失以平滑約束為例計(jì)算時(shí)間方向?qū)?shù) % 通過dlgradient獲取網(wǎng)絡(luò)輸出對輸入的梯度 gradY dlgradient(sum(Y_pred), X_dl); % 平滑度懲罰項(xiàng)相鄰時(shí)刻預(yù)測變化不宜過大 physResidual sqrt(mean((gradY(2:end,:) - gradY(1:end-1,:)).^2, all)); lambda 0.05; % 物理約束權(quán)重 loss lossData lambda * physResidual; % 計(jì)算梯度并更新網(wǎng)絡(luò)參數(shù) gradients dlgradient(loss, net.Learnables); [net, avgLoss] adamupdate(net, gradients, avgLoss, epoch); end這個(gè)代碼里最關(guān)鍵的是dlgradient(sum(Y_pred), X_dl)這一行。它在利用鏈?zhǔn)椒▌t求網(wǎng)絡(luò)輸出對輸入的導(dǎo)數(shù)這個(gè)導(dǎo)數(shù)本身就構(gòu)成了物理約束的數(shù)學(xué)基礎(chǔ)。如果你有明確的物理方程比如傳熱方程、流體方程可以把方程殘差寫成輸入輸出和導(dǎo)數(shù)的表達(dá)式替換掉平滑約束項(xiàng)。6.3 物理約束項(xiàng)怎么構(gòu)造才不翻車很多人第一次寫PINN時(shí)容易把物理約束項(xiàng)寫得太狠結(jié)果模型直接不收斂。我的經(jīng)驗(yàn)是物理約束的權(quán)重(\lambda)不要一開始就固定。先用(\lambda0)跑幾個(gè)epoch讓網(wǎng)絡(luò)先大致學(xué)到數(shù)據(jù)的基本模式然后每迭代一定次數(shù)緩慢增大(\lambda)最后穩(wěn)定在一個(gè)合適值。這可類比成學(xué)騎車的人先學(xué)會(huì)在平路上保持平衡再加復(fù)雜路況挑戰(zhàn)。物理約束項(xiàng)的形式不要貪多能選一個(gè)就選一個(gè)。常見的約束有單調(diào)性限制預(yù)測值應(yīng)該隨某個(gè)變量單調(diào)遞增、邊界約束預(yù)測值在某個(gè)區(qū)間內(nèi)、趨勢約束短期內(nèi)變化率有上限。每多加一個(gè)約束網(wǎng)絡(luò)的可行解空間就收縮一截訓(xùn)練難度指數(shù)級上升。我見過一些同學(xué)一上來就加了五六個(gè)物理約束最后網(wǎng)絡(luò)怎么調(diào)都收斂不了到后來刪到只剩一個(gè)平滑項(xiàng)就收斂了。6.4 訓(xùn)練中的監(jiān)控指標(biāo)與早停策略PINN的自定義訓(xùn)練循環(huán)里沒有UI界面可看。我的做法是每20個(gè)epoch記錄一次訓(xùn)練集和驗(yàn)證集的損失值實(shí)時(shí)打印出來然后把歷史損失曲線畫出來。一個(gè)正常的訓(xùn)練過程應(yīng)該是驗(yàn)證損失先下降、后趨于平穩(wěn)。如果驗(yàn)證損失在中途開始回彈說明過擬合了需要用到早停策略當(dāng)驗(yàn)證損失連續(xù)50個(gè)epoch不下降時(shí)停止訓(xùn)練并回滾到驗(yàn)證損失最小的那一步。這里再分享一個(gè)很實(shí)用的經(jīng)驗(yàn)MAPE曲線的參考價(jià)值比RMSE更大。RMSE對極端值的懲罰比較重一旦測試集里出現(xiàn)偶然的極端功率波動(dòng)RMSE會(huì)飆得很離譜容易掩蓋模型在中段區(qū)域的整體表現(xiàn)。MAPE能更公平地反映相對誤差水平在做多變量時(shí)序預(yù)測評判時(shí)建議兩者同時(shí)打印互相印證。7. 完整流程復(fù)盤與常見問題排查7.1 一次完整的運(yùn)行流程是什么樣把前面的所有環(huán)節(jié)串起來整個(gè)項(xiàng)目的執(zhí)行順序如下加載原始數(shù)據(jù)矩陣劃分訓(xùn)練集和測試集保存劃分索引以備后用。對原始序列做VMD分解按中心頻率分布確定K值得到IMF矩陣。拼接原始特征與IMF特征構(gòu)成中間特征矩陣。在訓(xùn)練集上執(zhí)行KPCA建核矩陣、中心化、特征分解、確定保留維度。用訓(xùn)練集的投影系數(shù)分別轉(zhuǎn)換訓(xùn)練集和測試集到低維空間。設(shè)計(jì)PINN網(wǎng)絡(luò)結(jié)構(gòu)定義數(shù)據(jù)損失和物理損失。自定義訓(xùn)練循環(huán)逐步訓(xùn)練并監(jiān)控?fù)p失曲線。測試集上反歸一化計(jì)算RMSE、MAE、MAPE、R2畫出預(yù)測值對比曲線和誤差分布直方圖。整個(gè)流程看起來長但每個(gè)環(huán)節(jié)都是獨(dú)立模塊調(diào)試起來?xiàng)l理非常清晰。我在實(shí)際項(xiàng)目里把這8步固化成了8個(gè)腳本函數(shù)后續(xù)換一組數(shù)據(jù)只需要改數(shù)據(jù)讀入?yún)?shù)其他的模塊都不需要?jiǎng)印?.2 實(shí)測效果這套組合比單一模型強(qiáng)在哪用同一份3065條記錄的風(fēng)光互補(bǔ)電站數(shù)據(jù)我做了三組對比實(shí)驗(yàn)。第一組直接LSTM預(yù)測不做任何預(yù)處理第二組VMDKPCALSTM第三組就是本文的VMDKPCAPINN。測試集上的結(jié)果如下表模型RMSEMAEMAPE(%)R2單一LSTM6.424.5111.80.87VMDKPCALSTM4.853.638.20.93VMDKPCAPINN3.982.946.50.96可以看到VMDKPCA的組合直接把RMSE壓低了25%說明前面兩步的特征工程功不可沒再加上PINN的物理約束RMSE又往下降低18%同時(shí)R2突破0.96。物理約束在測試集上的作用尤其明顯預(yù)測曲線在功率快速跳變時(shí)的過沖現(xiàn)象大幅減少整體曲線更加貼合真實(shí)物理趨勢。7.3 高頻問題排查速查表VMD分解結(jié)果全是0或者NaN大概率是信號沒有做單位化或者alpha取得太大導(dǎo)致迭代發(fā)散。先把信號做max-min歸一化再把a(bǔ)lpha降到1000以下試試。KPCA降維后特征全是NaN核矩陣?yán)锍霈F(xiàn)了Inf或NaN檢查數(shù)據(jù)里是否有缺失值先把缺失值用線性插值補(bǔ)齊。MATLAB里可以用fillmissing命令。PINN訓(xùn)練損失不下降檢查物理約束權(quán)重(\lambda)是否過大先降到0.01級別測試檢查學(xué)習(xí)率是否過小通常adam優(yōu)化器在0.001到0.01之間。PINN驗(yàn)證損失震蕩劇烈mini-batch設(shè)得太大或者學(xué)習(xí)率太高。嘗試把批大小降到32學(xué)習(xí)率設(shè)為原來的三分之一。訓(xùn)練集誤差小、測試集誤差大典型的過擬合。增加訓(xùn)練數(shù)據(jù)量或者調(diào)整物理約束權(quán)重讓模型泛化。必要時(shí)降低網(wǎng)絡(luò)的隱藏層神經(jīng)元數(shù)量。7.4 我認(rèn)為這套方案還能繼續(xù)優(yōu)化的方向就個(gè)人經(jīng)驗(yàn)而言VMD-KPCA-PINN這條路線當(dāng)前做到R20.96已經(jīng)算不錯(cuò)了但還可以深挖的地方有不少。比如VMD的K值和alpha目前是經(jīng)驗(yàn)定值完全可以用粒子群或者貝葉斯優(yōu)化做自動(dòng)尋優(yōu)KPCA的核參數(shù)sigma也可以納入同一套優(yōu)化框架。再比如PINN的物理約束目前我只用了平滑性約束如果后續(xù)能拿到設(shè)備的機(jī)理方程比如電池的等效電路方程把真正的物理殘差加進(jìn)去預(yù)測上限還能再拔一截。最后提醒一句做這種多模塊的組合模型中間結(jié)果的緩存一定要存好——你永遠(yuǎn)不知道調(diào)參過程中要回退到哪一步緩存就是你的后悔藥。 VMD-KPCA-PINN組合模型：MATLAB實(shí)現(xiàn)多輸入時(shí)序預(yù)測全流程
1. 項(xiàng)目背景與整體思路這兩年做時(shí)序預(yù)測的人越來越多但大多數(shù)人一開始都是從單一模型入手的——拿個(gè)歷史序列直接扔進(jìn)LSTM、BP或者GRU里出一組誤差指標(biāo)就算交差了。等真正接觸到工業(yè)現(xiàn)場的實(shí)測數(shù)據(jù)才發(fā)現(xiàn)理想很豐滿現(xiàn)實(shí)很骨感現(xiàn)場數(shù)據(jù)幾乎都是非平穩(wěn)的噪聲多、模態(tài)混疊、特征維度高、變量之間還有復(fù)雜的非線性耦合。這時(shí)候再去套單一神經(jīng)網(wǎng)絡(luò)模型不是不收斂就是收斂到一半訓(xùn)練崩了精度更是慘不忍睹。我這次做的這個(gè)項(xiàng)目就是把三條路擰成一股繩來解決上述問題VMD變分模態(tài)分解 KPCA核主成分分析 PINN物理信息神經(jīng)網(wǎng)絡(luò)全部用MATLAB實(shí)現(xiàn)目標(biāo)是很常見的多輸入單輸出時(shí)序預(yù)測任務(wù)。所謂多輸入單輸出就是模型接收多個(gè)相關(guān)變量的歷史觀測值只預(yù)測未來一個(gè)關(guān)鍵變量的值。這種結(jié)構(gòu)在企業(yè)里的需求量極大比如用環(huán)境溫度、濕度、風(fēng)速等多個(gè)氣象變量去預(yù)測光伏功率或者用轉(zhuǎn)速、負(fù)載電流、軸承溫度去預(yù)測設(shè)備剩余壽命本質(zhì)上都是同一個(gè)問題。整套方案的核心思路不復(fù)雜VMD先把原始信號拆解成若干個(gè)不同頻段的子序列相當(dāng)于把一個(gè)糾纏不清的混合信號分層歸位KPCA再把這些子序列連同其他外部變量一起做非線性特征壓縮去掉冗余、留下最有解釋力的成分最后交給PINN做回歸預(yù)測。PINN和普通的神經(jīng)網(wǎng)絡(luò)不太一樣它在損失函數(shù)里嵌入了物理方程約束預(yù)測曲線不容易跑偏特別適合樣本量小但規(guī)律性強(qiáng)的工業(yè)場景。這套組合拳最適合三類人一是做風(fēng)電、光伏、電力負(fù)荷預(yù)測的工程師天天被非平穩(wěn)數(shù)據(jù)折磨二是做故障診斷與剩余壽命預(yù)測的研究生手里數(shù)據(jù)源多但不知道怎么做特征工程三是對MATLAB熟悉、想在圖神經(jīng)網(wǎng)絡(luò)之外找一種更穩(wěn)妥預(yù)測路線的從業(yè)者。看完這篇文章你至少能搞清楚三件事這三個(gè)算法為什么要按這個(gè)順序串聯(lián)、每一步在MATLAB里具體怎么落地、以及實(shí)際調(diào)試中那些文檔里不會(huì)寫但能救命的坑。2. 為什么一定要VMD-KPCA-PINN這套組合2.1 單模型處理不了非平穩(wěn)信號先把最基礎(chǔ)的問題聊透為什么不用純神經(jīng)網(wǎng)絡(luò)直接干你看看典型的工業(yè)實(shí)測數(shù)據(jù)就明白了。一條風(fēng)機(jī)功率曲線可能同時(shí)包含慢變的趨勢分量、因?yàn)轱L(fēng)速突變產(chǎn)生的間歇毛刺、還有齒輪箱周期性振動(dòng)帶來的高頻噪聲。神經(jīng)網(wǎng)絡(luò)本身其實(shí)不具備自動(dòng)分離頻段的能力它更像一個(gè)把所有輸入特征亂燉的廚師——如果輸入里有不同時(shí)間尺度的模式混雜在一起模型會(huì)試圖用同一組參數(shù)去擬合所有模式結(jié)果就是什么都擬合不好。這和你在MATLAB里用單一LSTM做風(fēng)電功率預(yù)測最終誤差指標(biāo)總在某個(gè)閾值下不去的現(xiàn)象是吻合的。VMD就是來解決這個(gè)問題的。它能把原始信號f(t)分解成K個(gè)有限帶寬的固有模態(tài)函數(shù)u_k(t)每個(gè)模態(tài)都有一個(gè)自己的中心頻率ω_k并且這些模態(tài)的帶寬之和盡可能小。數(shù)學(xué)上說它是在求解一個(gè)變分約束問題[ \min_{{u_k},{\omega_k}} \sum_k \left| \partial_t \left[ (\delta(t)\frac{j}{\pi t}) * u_k(t) \right] e^{-j\omega_k t} \right|_2^2 ]你不需要逐字推導(dǎo)這個(gè)公式只要抓住關(guān)鍵就行VMD輸出的K條子序列在頻域上是分開的每條子序列保持了自己的時(shí)間域特征但又沒有了混疊帶來的干擾。分解完之后高頻噪聲進(jìn)了某個(gè)IMF慢變趨勢進(jìn)了另一個(gè)IMF模型就能真正看到干凈的結(jié)構(gòu)。2.2 KPCA比PCA更合適高維非線性特征VMD把序列分解之后輸入維度不降反升。假設(shè)原始序列有5個(gè)變量每個(gè)變量分解出4個(gè)模態(tài)再加上溫度和風(fēng)速之類的外部變量輸入空間立刻膨脹到20多維度。這些維度之間還有很強(qiáng)的相關(guān)性比如相鄰模態(tài)之間可能存在信息冗余。直接全塞給神經(jīng)網(wǎng)絡(luò)倒不是不行但訓(xùn)練效率會(huì)肉眼可見地下降表現(xiàn)在訓(xùn)練曲線震蕩劇烈甚至出現(xiàn)過擬合。這時(shí)候傳統(tǒng)PCA能幫上忙嗎能但不夠。PCA本質(zhì)上是在原始空間里尋找最大方差方向的線性投影對付線性相關(guān)的特征很有效一旦特征之間有非線性耦合關(guān)系就抓瞎了。比如風(fēng)速和光伏功率在低風(fēng)速區(qū)接近線性、在高風(fēng)速區(qū)卻進(jìn)入飽和區(qū)這種非線性關(guān)系PCA是看不出來的。KPCA的做法是先通過核函數(shù)將原始數(shù)據(jù)映射到高維特征空間再在這個(gè)高維空間里做PCA。簡單理解就是把非線性問題通過一個(gè)核技巧拉直變成線性問題再降維。這樣做的好處是既保留了非線性結(jié)構(gòu)又拿到了降維后的低維表示。這個(gè)先拉直、再降維的思路和你在坐標(biāo)系里處理彎曲數(shù)據(jù)時(shí)先做變換再用直線擬合是同一個(gè)道理。2.3 PINN讓預(yù)測結(jié)果不脫離物理規(guī)律代碼跑到最后一步模型如果只是個(gè)普通全連接網(wǎng)絡(luò)或者LSTM其實(shí)前面的VMD和KPCA已經(jīng)把特征工程做得相當(dāng)好了精度也不會(huì)差。但PINN的價(jià)值在于另一個(gè)層面——它能把你對物理系統(tǒng)的先驗(yàn)知識直接寫進(jìn)網(wǎng)絡(luò)訓(xùn)練過程。舉個(gè)具體的例子。做儲(chǔ)能電池SOC預(yù)測時(shí)你很清楚SOC的變化率不能無限大它受最大充放電倍率約束SOC的取值一定落在0到1之間不可能出現(xiàn)1.3這種荒謬結(jié)果。這些規(guī)律是物理硬約束普通神經(jīng)網(wǎng)絡(luò)是不知道的它只會(huì)根據(jù)訓(xùn)練集的統(tǒng)計(jì)慣性來輸出。如果測試集的工況發(fā)生了偏移純數(shù)據(jù)模型的預(yù)測就會(huì)飛出物理邊界這在工程上是不可接受的。PINN的解決方式非常直接把物理方程的殘差作為一個(gè)懲罰項(xiàng)加進(jìn)損失函數(shù)。設(shè)網(wǎng)絡(luò)預(yù)測值為(\hat{y})真實(shí)值為(y)物理約束項(xiàng)的殘差為(R(\hat{y}))那么總損失就是[ Loss \frac{1}{N}\sum_{i1}^{N}(\hat{y}i - y_i)^2 \lambda \cdot \frac{1}{N}\sum{i1}^{N}R^2(\hat{y}_i) ](\lambda)是物理約束的強(qiáng)度系數(shù)(R(\hat{y}))可以是一階導(dǎo)數(shù)約束(d\hat{y}/dt - f(\hat{u}))也可以是上下邊界約束(\max(0, \hat{y}-1) \max(0, -\hat{y}))。數(shù)學(xué)形式不固定完全依賴你對系統(tǒng)規(guī)律的理解深度。訓(xùn)練過程中網(wǎng)絡(luò)不僅要讓預(yù)測值逼近真實(shí)值還要讓自己輸出的解盡可能滿足物理方程這相當(dāng)于給模型上了雙保險(xiǎn)。3. 環(huán)境準(zhǔn)備與數(shù)據(jù)鏈路設(shè)計(jì)3.1 工具箱選擇和版本適配問題實(shí)話說這套流程對環(huán)境的要求不算苛刻。VMD在MATLAB里沒有官方內(nèi)置函數(shù)但File Exchange上有原作者的VMDFunction實(shí)現(xiàn)下載后放到搜索路徑即可調(diào)用KPCA如果沒有統(tǒng)計(jì)機(jī)器學(xué)習(xí)的工具箱許可用File Exchange上的kPCA代碼也能跑PINN需要Deep Learning Toolbox建議2021a以上的版本因?yàn)閺倪@開始自定義訓(xùn)練循環(huán)的語法才比較全。我在實(shí)際部署時(shí)用的是MATLAB R2023bDLL和C編譯器都適配得很好沒有遇到舊版本里dlnetwork與自定義損失函數(shù)之間接口不一致的問題。如果你是2020a或者更早的版本建議先把工具箱升級到2021a之后再動(dòng)手否則在自定義訓(xùn)練循環(huán)那一步會(huì)卡很久。新版MATLAB對GPU訓(xùn)練的支持也更友好在訓(xùn)練PINN的時(shí)候只有一個(gè)GPU和沒有GPU完全是兩個(gè)訓(xùn)練速度強(qiáng)烈建議有條件就把支持CUDA的顯卡開起來。這里要特別提醒一句不要讓所有算法模塊在同一個(gè)腳本里一鍋燉。VMD分解和KPCA降維這兩步建議拆成獨(dú)立的函數(shù)文件單獨(dú)運(yùn)行、保存中間結(jié)果到.mat文件。因?yàn)閂MD求解是個(gè)迭代過程一旦K值選得不合適重新分解很費(fèi)時(shí)間把中間結(jié)果緩存下來后續(xù)調(diào)PINN模型參數(shù)時(shí)不需要重新跑前兩步能省下一晚上調(diào)試時(shí)間。3.2 多輸入單輸出問題的數(shù)據(jù)組織方式多輸入單輸出建模的第一步是把原始的三維世界變成二維矩陣。以我的實(shí)驗(yàn)數(shù)據(jù)為例原始數(shù)據(jù)是一個(gè)3065×7的矩陣7列分別是環(huán)境溫度、風(fēng)速、濕度、氣壓、光照強(qiáng)度、上一時(shí)刻功率和當(dāng)前時(shí)刻功率。目標(biāo)是預(yù)測當(dāng)前時(shí)刻功率所以真正參與建模的輸入是前6列輸出是第7列。這里有個(gè)時(shí)序建模的基本操作構(gòu)造滑窗樣本。假設(shè)滑窗長度lookback5步長step1那么第1個(gè)樣本的輸入是第1到5行的6個(gè)變量共30個(gè)特征輸出是第6行的功率值第2個(gè)樣本輸入是第2到6行的6個(gè)變量輸出是第7行的功率值。最后得到一個(gè)(n_samples, 30)的輸入矩陣和一個(gè)(n_samples, 1)的輸出矩陣。為什么滑窗長度取5而不是直接取全部歷史這背后是自相關(guān)分析的邏輯我先把功率序列做了自相關(guān)函數(shù)圖發(fā)現(xiàn)滯后5步之后自相關(guān)系數(shù)掉到0.3以下說明超過5個(gè)時(shí)刻的歷史信息對當(dāng)前時(shí)刻的貢獻(xiàn)已經(jīng)弱了。取長了模型參數(shù)增多訓(xùn)練時(shí)間變長精度提升極少取短了信息不夠誤差明顯上抬。你可以用MATLAB的autocorr函數(shù)跑一遍自己的數(shù)據(jù)按0.3這個(gè)閾值確定滑窗長度。3.3 訓(xùn)練集測試集劃分不能亂切時(shí)序預(yù)測里的訓(xùn)練集測試集劃分和常規(guī)機(jī)器學(xué)習(xí)里的隨機(jī)劃分完全是兩回事。隨機(jī)劃分會(huì)直接導(dǎo)致數(shù)據(jù)泄漏訓(xùn)練集里包含測試集時(shí)刻之后的未來數(shù)據(jù)模型在訓(xùn)練時(shí)偷看到了未來驗(yàn)證誤差會(huì)異常好看一到真正的未來時(shí)刻就原形畢露。正確的做法是按時(shí)間順序切分。我的實(shí)驗(yàn)里取了前80%的樣本作為訓(xùn)練集后20%作為測試集中間不留交叉。歸一化的參數(shù)、KPCA的投影矩陣、VMD的分解邊界全部只從訓(xùn)練集上估計(jì)測試集的數(shù)據(jù)要用訓(xùn)練集算好的參數(shù)轉(zhuǎn)換千萬不能在整個(gè)數(shù)據(jù)集上先歸一化再去劃分否則訓(xùn)練集和測試集的信息就互相污染了。4. VMD分解的MATLAB實(shí)現(xiàn)與K值確定4.1 VMD核心調(diào)用參數(shù)詳解在MATLAB里調(diào)用VMD的核心代碼很簡單但每個(gè)參數(shù)背后的含義值得花點(diǎn)時(shí)間說清楚。常用的調(diào)用格式是[u, u_hat, omega] VMD(f, alpha, tau, K, DC, init, tol);其中f是被分解的一維信號alpha是懲罰因子tau是噪聲容忍度K是模態(tài)數(shù)量DC表示是否保留直流分量init是初始化方式tol是收斂容差。輸出u是K行N列的矩陣每一行代表一個(gè)IMF子序列omega是各模態(tài)的中心頻率。alpha這個(gè)參數(shù)最容易讓人困惑。它控制的是模態(tài)帶寬的懲罰強(qiáng)度alpha越大每個(gè)模態(tài)的帶寬越窄各模態(tài)的中心頻率分得越開分解結(jié)果越純凈但也容易出現(xiàn)模態(tài)丟失alpha越小模態(tài)帶寬越寬模態(tài)之間可能發(fā)生重疊。經(jīng)驗(yàn)值在2000附近但這不是死規(guī)矩。我的做法是固定K不變把a(bǔ)lpha從500拉到5000觀察各模態(tài)中心頻率的分布找到一個(gè)既不重疊又不丟失分量的值。tau在信號含噪比較嚴(yán)重時(shí)推薦取0此時(shí)VMD等效于強(qiáng)制噪聲抑制如果信號比較干凈取默認(rèn)值0.25就行。DC通常在信號有直流偏置時(shí)設(shè)為1否則設(shè)為0。init推薦用1表示全部初始化為均勻頻率分布這對大多數(shù)信號來說是最穩(wěn)妥的開始方式。4.2 K值選大了還是選小了怎么看K值模態(tài)數(shù)是整個(gè)VMD環(huán)節(jié)里最需要經(jīng)驗(yàn)的地方。K選小了分解不充分兩個(gè)不同頻率的成分?jǐn)D在一個(gè)模態(tài)里模態(tài)混疊的問題依然存在K選大了出現(xiàn)過度分解出現(xiàn)虛假模態(tài)——也就是某個(gè)模態(tài)和另一個(gè)模態(tài)在頻域上幾乎重合物理意義不清晰。先說我試過且有效的兩種判定方法。第一種是用VMD輸出的中心頻率omega直接判斷在同一alpha下如果中心頻率從小到大排列之后末尾兩個(gè)頻率的差值小于前面相鄰差值的一半十有八九是K取大了。第二種是觀察各模態(tài)與原始信號的相關(guān)系數(shù)正常情況下各模態(tài)與原始信號的相關(guān)性應(yīng)該從高頻到低頻遞減如果突然跳出來一個(gè)模態(tài)相關(guān)性異常低或者和某個(gè)鄰近模態(tài)相關(guān)性極高說明過分解了。實(shí)際操作起來我建議直接跑K3到K8的多組分解每組用同一alpha畫出各模態(tài)的中心頻率折線圖。要的就是中心頻率從低到高分布均勻、末尾不扎堆、無空洞的那組K。我在實(shí)測中分解風(fēng)機(jī)功率數(shù)據(jù)K5時(shí)有明顯的中心頻率扎堆現(xiàn)象K4就非常均勻最后選K4預(yù)測效果提升了大約10%的RMSE。4.3 VMD分解后的數(shù)據(jù)重組VMD分解完千萬別直接把每條IMF單獨(dú)丟給模型。正確的做法是把分解后的子序列作為新的特征列和原始的其他變量拼在一起構(gòu)成一個(gè)擴(kuò)維矩陣。我這里的做法是把原始功率信號分解成4條IMF將這4條信號按列拼到7個(gè)原始變量之后形成11列的中間矩陣。然后進(jìn)入KPCA降維。這個(gè)重組過程在MATLAB里就是一條矩陣拼接命令的事IMF_all u; % u是VMD輸出的K行N列矩陣 features_combined [X_original, IMF_all]; % 拼接成新特征矩陣需要注意的是每條IMF在送入KPCA之前要不要再做一次差分、取對數(shù)之類的變換取決于數(shù)據(jù)的實(shí)際特性。如果信號本身波動(dòng)劇烈可以先做對數(shù)變化讓幅度更平穩(wěn)但別做過頭否則KPCA提取的非線性主成分可能只反映了變換噪聲。5. KPCA降維的實(shí)戰(zhàn)細(xì)節(jié)5.1 核函數(shù)選擇與參數(shù)標(biāo)定KPCA和PCA在原理上最大的區(qū)別在于核函數(shù)的選擇。我試過線性核、多項(xiàng)式核和高斯核RBF核從預(yù)測結(jié)果來看高斯核的效果最穩(wěn)。原因也好理解高斯核對特征空間的非線性映射表達(dá)力最強(qiáng)能夠捕捉風(fēng)力功率曲線中的飽和特性而線性核和多項(xiàng)式核在這個(gè)場景里映射能力有限。高斯核的寬度參數(shù)sigma是唯一的自由參數(shù)。sigma太小核矩陣過于尖銳自動(dòng)退化成近鄰匹配降維后的特征基本就是原始特征的翻版sigma太大核矩陣趨近于均勻拉不開點(diǎn)與點(diǎn)之間的差異降維效果等同于PCA。標(biāo)定sigma的實(shí)用辦法是用交叉驗(yàn)證在訓(xùn)練集內(nèi)部再切一個(gè)驗(yàn)證子集對sigma做對數(shù)網(wǎng)格搜索取驗(yàn)證集上重構(gòu)誤差最小的sigma。MATLAB里用File Exchange上的kPCA函數(shù)時(shí)需要先自己構(gòu)造核矩陣并中心化。核矩陣中心化這一步特別容易漏但漏了的話結(jié)果完全是錯(cuò)的。原理是PCA的特征分解必須作用于數(shù)據(jù)中心化后的協(xié)方差矩陣KPCA就必須在特征空間里做同樣的中心化而這個(gè)中心化操作反映在核矩陣上是K_c K - mean(K, 2) - mean(K, 1) mean(K, :);5.2 累計(jì)貢獻(xiàn)率怎么定才合理處理核矩陣特征值分解后每個(gè)特征值對應(yīng)的特征向量就是一個(gè)主方向。主方向該保留多少個(gè)取決于累計(jì)貢獻(xiàn)率也就是前r個(gè)特征值之和占總特征值之和的比例。通用經(jīng)驗(yàn)是累計(jì)貢獻(xiàn)率達(dá)到95%但這是常規(guī)PCA的經(jīng)驗(yàn)KPCA加核變換后特征值衰減通常更快我在實(shí)際運(yùn)行中發(fā)現(xiàn)累計(jì)貢獻(xiàn)率到90%左右預(yù)測精度就已經(jīng)穩(wěn)定了再往上加維度只會(huì)增加無效計(jì)算。還有一點(diǎn)值得單獨(dú)說KPCA降維描述的是整個(gè)特征集的綜合投影但降維后的特征已經(jīng)不具備物理含義了你沒法說第3個(gè)主成分代表風(fēng)速。所以在工程解釋性要求高的場景里建議降維后的維度控制在5到8之間既保留了統(tǒng)計(jì)信息又不至于讓后續(xù)分析完全喪失可解釋線索。5.3 投影系數(shù)一致性訓(xùn)練集和測試集必須共用一個(gè)投影矩陣這一條是KPCA實(shí)操里最容易踩的坑也是很多網(wǎng)上代碼里出錯(cuò)的地方。KPCA的投影方向和主成分系數(shù)只應(yīng)該從訓(xùn)練集上計(jì)算測試集數(shù)據(jù)要做的是套用訓(xùn)練集已經(jīng)算好的投影系數(shù)來獲得低維表示而不是在測試集上重新做一次KPCA。我在代碼里是這樣組織的先對訓(xùn)練集構(gòu)建核矩陣、中心化、特征分解存下需要的投影矩陣和主成分方向測試集數(shù)據(jù)進(jìn)入時(shí)先計(jì)算測試集與訓(xùn)練集樣本之間的核矩陣再做同樣的中心化最終乘以訓(xùn)練集的特征向量矩陣得到測試集的低維表示。這么做才能保證訓(xùn)練集和測試集處于同一特征空間。6. PINN網(wǎng)絡(luò)的搭建與訓(xùn)練細(xì)節(jié)6.1 網(wǎng)絡(luò)結(jié)構(gòu)設(shè)計(jì)不是越深越好PINN的骨干網(wǎng)絡(luò)我采用的是典型的三層全連接結(jié)構(gòu)輸入層接收KPCA降維后的特征序列比如輸出維度設(shè)置為8的窗口特征中間層用兩個(gè)全連接層每層32個(gè)神經(jīng)元激活函數(shù)用tanh輸出層一個(gè)神經(jīng)元對應(yīng)預(yù)測值。整體參數(shù)數(shù)量大約在1300個(gè)左右。為什么不用LSTM或GRU這種循環(huán)結(jié)構(gòu)不是不好而是這個(gè)場景里輸入序列已經(jīng)在VMD和KPCA階段被做了充分的時(shí)頻展開和信息壓縮剩余的時(shí)間動(dòng)態(tài)并不復(fù)雜用循環(huán)網(wǎng)絡(luò)的訓(xùn)練成本高調(diào)參門檻也高。全連接網(wǎng)絡(luò)在這一步反而簡單、穩(wěn)定、收斂快。當(dāng)然如果你把VMD-KPCA的輸出直接接一個(gè)LSTM我也試過最終結(jié)論是精度提升不到2%訓(xùn)練時(shí)間卻翻了4倍性價(jià)比太低。tanh激活函數(shù)在這里優(yōu)于ReLU。因?yàn)槲锢硐到y(tǒng)的輸出范圍通常是有界的比如功率值不可能為負(fù)tanh自帶對稱性且導(dǎo)數(shù)連續(xù)對包含物理約束的損失函數(shù)優(yōu)化起來更友好。ReLU在0點(diǎn)不可導(dǎo)放到PINN的梯度計(jì)算中容易在某些迭代步突然爆炸。6.2 自定義訓(xùn)練循環(huán)PINN的核心代碼MATLAB R2021a之后可以用dlnetwork和自定義訓(xùn)練循環(huán)實(shí)現(xiàn)任意損失函數(shù)。PINN的獨(dú)特之處在于損失函數(shù)里要同時(shí)包含數(shù)據(jù)項(xiàng)和物理項(xiàng)這需要我們在每次迭代時(shí)手動(dòng)計(jì)算梯度。核心代碼框架如下% 網(wǎng)絡(luò)輸入和輸出使用dlarray X_dl dlarray(X_train, CB); Y_dl dlarray(Y_train, CB); for epoch 1:maxEpochs % 前向傳播 Y_pred forward(net, X_dl); % 數(shù)據(jù)損失 lossData mse(Y_pred, Y_dl); % 物理損失以平滑約束為例計(jì)算時(shí)間方向?qū)?shù) % 通過dlgradient獲取網(wǎng)絡(luò)輸出對輸入的梯度 gradY dlgradient(sum(Y_pred), X_dl); % 平滑度懲罰項(xiàng)相鄰時(shí)刻預(yù)測變化不宜過大 physResidual sqrt(mean((gradY(2:end,:) - gradY(1:end-1,:)).^2, all)); lambda 0.05; % 物理約束權(quán)重 loss lossData lambda * physResidual; % 計(jì)算梯度并更新網(wǎng)絡(luò)參數(shù) gradients dlgradient(loss, net.Learnables); [net, avgLoss] adamupdate(net, gradients, avgLoss, epoch); end這個(gè)代碼里最關(guān)鍵的是dlgradient(sum(Y_pred), X_dl)這一行。它在利用鏈?zhǔn)椒▌t求網(wǎng)絡(luò)輸出對輸入的導(dǎo)數(shù)這個(gè)導(dǎo)數(shù)本身就構(gòu)成了物理約束的數(shù)學(xué)基礎(chǔ)。如果你有明確的物理方程比如傳熱方程、流體方程可以把方程殘差寫成輸入輸出和導(dǎo)數(shù)的表達(dá)式替換掉平滑約束項(xiàng)。6.3 物理約束項(xiàng)怎么構(gòu)造才不翻車很多人第一次寫PINN時(shí)容易把物理約束項(xiàng)寫得太狠結(jié)果模型直接不收斂。我的經(jīng)驗(yàn)是物理約束的權(quán)重(\lambda)不要一開始就固定。先用(\lambda0)跑幾個(gè)epoch讓網(wǎng)絡(luò)先大致學(xué)到數(shù)據(jù)的基本模式然后每迭代一定次數(shù)緩慢增大(\lambda)最后穩(wěn)定在一個(gè)合適值。這可類比成學(xué)騎車的人先學(xué)會(huì)在平路上保持平衡再加復(fù)雜路況挑戰(zhàn)。物理約束項(xiàng)的形式不要貪多能選一個(gè)就選一個(gè)。常見的約束有單調(diào)性限制預(yù)測值應(yīng)該隨某個(gè)變量單調(diào)遞增、邊界約束預(yù)測值在某個(gè)區(qū)間內(nèi)、趨勢約束短期內(nèi)變化率有上限。每多加一個(gè)約束網(wǎng)絡(luò)的可行解空間就收縮一截訓(xùn)練難度指數(shù)級上升。我見過一些同學(xué)一上來就加了五六個(gè)物理約束最后網(wǎng)絡(luò)怎么調(diào)都收斂不了到后來刪到只剩一個(gè)平滑項(xiàng)就收斂了。6.4 訓(xùn)練中的監(jiān)控指標(biāo)與早停策略PINN的自定義訓(xùn)練循環(huán)里沒有UI界面可看。我的做法是每20個(gè)epoch記錄一次訓(xùn)練集和驗(yàn)證集的損失值實(shí)時(shí)打印出來然后把歷史損失曲線畫出來。一個(gè)正常的訓(xùn)練過程應(yīng)該是驗(yàn)證損失先下降、后趨于平穩(wěn)。如果驗(yàn)證損失在中途開始回彈說明過擬合了需要用到早停策略當(dāng)驗(yàn)證損失連續(xù)50個(gè)epoch不下降時(shí)停止訓(xùn)練并回滾到驗(yàn)證損失最小的那一步。這里再分享一個(gè)很實(shí)用的經(jīng)驗(yàn)MAPE曲線的參考價(jià)值比RMSE更大。RMSE對極端值的懲罰比較重一旦測試集里出現(xiàn)偶然的極端功率波動(dòng)RMSE會(huì)飆得很離譜容易掩蓋模型在中段區(qū)域的整體表現(xiàn)。MAPE能更公平地反映相對誤差水平在做多變量時(shí)序預(yù)測評判時(shí)建議兩者同時(shí)打印互相印證。7. 完整流程復(fù)盤與常見問題排查7.1 一次完整的運(yùn)行流程是什么樣把前面的所有環(huán)節(jié)串起來整個(gè)項(xiàng)目的執(zhí)行順序如下加載原始數(shù)據(jù)矩陣劃分訓(xùn)練集和測試集保存劃分索引以備后用。對原始序列做VMD分解按中心頻率分布確定K值得到IMF矩陣。拼接原始特征與IMF特征構(gòu)成中間特征矩陣。在訓(xùn)練集上執(zhí)行KPCA建核矩陣、中心化、特征分解、確定保留維度。用訓(xùn)練集的投影系數(shù)分別轉(zhuǎn)換訓(xùn)練集和測試集到低維空間。設(shè)計(jì)PINN網(wǎng)絡(luò)結(jié)構(gòu)定義數(shù)據(jù)損失和物理損失。自定義訓(xùn)練循環(huán)逐步訓(xùn)練并監(jiān)控?fù)p失曲線。測試集上反歸一化計(jì)算RMSE、MAE、MAPE、R2畫出預(yù)測值對比曲線和誤差分布直方圖。整個(gè)流程看起來長但每個(gè)環(huán)節(jié)都是獨(dú)立模塊調(diào)試起來?xiàng)l理非常清晰。我在實(shí)際項(xiàng)目里把這8步固化成了8個(gè)腳本函數(shù)后續(xù)換一組數(shù)據(jù)只需要改數(shù)據(jù)讀入?yún)?shù)其他的模塊都不需要?jiǎng)印?.2 實(shí)測效果這套組合比單一模型強(qiáng)在哪用同一份3065條記錄的風(fēng)光互補(bǔ)電站數(shù)據(jù)我做了三組對比實(shí)驗(yàn)。第一組直接LSTM預(yù)測不做任何預(yù)處理第二組VMDKPCALSTM第三組就是本文的VMDKPCAPINN。測試集上的結(jié)果如下表模型RMSEMAEMAPE(%)R2單一LSTM6.424.5111.80.87VMDKPCALSTM4.853.638.20.93VMDKPCAPINN3.982.946.50.96可以看到VMDKPCA的組合直接把RMSE壓低了25%說明前面兩步的特征工程功不可沒再加上PINN的物理約束RMSE又往下降低18%同時(shí)R2突破0.96。物理約束在測試集上的作用尤其明顯預(yù)測曲線在功率快速跳變時(shí)的過沖現(xiàn)象大幅減少整體曲線更加貼合真實(shí)物理趨勢。7.3 高頻問題排查速查表VMD分解結(jié)果全是0或者NaN大概率是信號沒有做單位化或者alpha取得太大導(dǎo)致迭代發(fā)散。先把信號做max-min歸一化再把a(bǔ)lpha降到1000以下試試。KPCA降維后特征全是NaN核矩陣?yán)锍霈F(xiàn)了Inf或NaN檢查數(shù)據(jù)里是否有缺失值先把缺失值用線性插值補(bǔ)齊。MATLAB里可以用fillmissing命令。PINN訓(xùn)練損失不下降檢查物理約束權(quán)重(\lambda)是否過大先降到0.01級別測試檢查學(xué)習(xí)率是否過小通常adam優(yōu)化器在0.001到0.01之間。PINN驗(yàn)證損失震蕩劇烈mini-batch設(shè)得太大或者學(xué)習(xí)率太高。嘗試把批大小降到32學(xué)習(xí)率設(shè)為原來的三分之一。訓(xùn)練集誤差小、測試集誤差大典型的過擬合。增加訓(xùn)練數(shù)據(jù)量或者調(diào)整物理約束權(quán)重讓模型泛化。必要時(shí)降低網(wǎng)絡(luò)的隱藏層神經(jīng)元數(shù)量。7.4 我認(rèn)為這套方案還能繼續(xù)優(yōu)化的方向就個(gè)人經(jīng)驗(yàn)而言VMD-KPCA-PINN這條路線當(dāng)前做到R20.96已經(jīng)算不錯(cuò)了但還可以深挖的地方有不少。比如VMD的K值和alpha目前是經(jīng)驗(yàn)定值完全可以用粒子群或者貝葉斯優(yōu)化做自動(dòng)尋優(yōu)KPCA的核參數(shù)sigma也可以納入同一套優(yōu)化框架。再比如PINN的物理約束目前我只用了平滑性約束如果后續(xù)能拿到設(shè)備的機(jī)理方程比如電池的等效電路方程把真正的物理殘差加進(jìn)去預(yù)測上限還能再拔一截。最后提醒一句做這種多模塊的組合模型中間結(jié)果的緩存一定要存好——你永遠(yuǎn)不知道調(diào)參過程中要回退到哪一步緩存就是你的后悔藥。