基于WSET时频能量谱与CNN的滚动轴承故障诊断Matlab实现

做滚动轴承故障诊断做到第三年我踩过最痛的一个坑是把模型复杂度和诊断准确率划了等号。后来把精力放回信号处理端用同步挤压时频能量谱WSET替代普普通通的连续小波变换再把二维时频图直接喂给CNN反而在公开的CWRU滚动轴承数据上稳定跑出高准确率。这篇文章就把整套流程摊开讲数据读取、WSET时频图生成、CNN搭建训练、可视化分析全部在Matlab环境里落地。适合有信号处理基础、但不太熟悉深度学习工程落地的机械专业学生也适合想用公开数据集快速验证新想法的算法工程师。你会看到每个关键步骤我为什么这么选以及实际跑实验时容易翻车的几个地方。1. 为什么用WSET-CNN而不是直接堆网络先说结论滚动轴承故障诊断这个任务特征提取端往往比网络结构更值得花时间。很多论文把网络做得很深但输入还停留在原始波形或者简单频谱。不是不能学问题是原始振动信号里背景噪声、转频谐波和随机冲击混在一起卷积神经网络虽然能自动提取特征却不会自动替你完成“分离调制分量”这种信号处理工作。所以整个思路的核心是把物理信息压缩到一张干净的二维时频能量谱里再让CNN去分类。1.1 时频分析滚动轴承故障信号绕不开的表达方式轴承故障本质上是周期性冲击。内圈、外圈、滚动体每通过损伤位置一次就会激起一次瞬态冲击并在结构固有频率附近激起谐振。这类信号最典型的特征是“非平稳”冲击的位置、幅值随时间变化而且淹没在背景噪声里。时域波形只能看到随机起伏的毛刺不太容易直接区分故障位置。FFT频谱只能给出全局平均根本看不出冲击发生在什么时候。早期工程上常用包络谱通过带通滤波希尔伯特解调把冲击信号提取出来再在频谱里找故障特征频率。这个方法很有效但非常依赖人工经验频带选哪里、边带怎么判不同人做出来的结果可能不一样。时频分析解决了这个问题。它同时保留时间和频率信息把“什么时候激起了什么频带的振动”画在一张二维图上。这一下就直观多了从二维图上能明显看到周期性冲击的分布形态这也是后来深度学习诊断模型普遍选择时频图作为输入的原因。1.2 同步挤压变换把模糊的CWT谱图“对焦”连续小波变换CWT已经能做时频分析但它有一个众所周知的毛病受海森堡不确定性限制能量会沿着频率方向扩散。就是谱图上本来应该是一条细线的地方实际画出来是一小片模糊带。肉眼看着还行但卷积神经网络在提取特征时模糊的带状区域会掩盖不同故障类别之间的细微差异分类边界就不容易学出来。同步挤压变换Synchrosqueezing Transform做的事情可以理解成给CWT谱图“对焦”。它利用小波系数的局部相位信息把扩散出去的能量重新分配搬运到真正的瞬时频率曲线附近。效果就是时频图上的能量团变得更集中、边界更锐利。对一个轴承振动信号来说周期冲击对应的时频脊线会非常清晰转频分量和谐波成分也分得比较开。在Matlab里这一步实现成本极低直接调用wsst函数即可。把输出系数取模平方后就得到了同步挤压时频能量谱我把它记为WSETWavelet-based Synchrosqueezing Energy Time-frequency。后面CNN吃的就是这张能量谱图而不是原始的CWT系数。1.3 为什么特征提取端要选CNN时频图本质上是一张图像CNN天生就是干这个的。卷积核在局部邻域上提取纹理方向、能量团形状、周期性条带结构不需要人为指定“关注哪个频带”一层卷积学冲击条带的方向两层卷积学不同故障模式之间的组合关系分类能力是自动涌现的。有人会问LSTM或者Transformer序列模型行不行时频图里频率轴和时间轴都有明确物理含义序列模型在单张二维图上天然吃亏你总得先手动决定沿哪个轴做序列化信息必然被截断。手工特征均值、峰值因子、峭度、包络谱峰值也很有用但面对不同损伤位置、不同负载工况的组合变化固定特征模板很容易漏掉关键差异。CNN方案的好处是只要WSET图的质量稳定模型可以端到端训练而且Matlab里的Deep Learning Toolbox把预处理、训练、可视化串成了一条完整工具链做验证实验非常顺手。2. CWRU公开轴承数据读取、分段与样本划分2.1 数据集结构与常用工况选择CWRU公开轴承数据集是故障诊断领域最常用的基准数据之一很多人叫它CWRU轴承数据。它提供驱动端和风扇端的加速度计振动信号采样率常见12kHz和48kHz两档健康状态包括正常、内圈故障、外圈故障、滚动体故障故障损伤尺寸从0.007英寸到0.021英寸不等并对应了0、1、2、3 HP四挡负载工况。做四分类实验时常规做法是把所有内圈故障样本统一划为IR类别外圈故障划为OR滚动体故障划为B正常样本划为N。不管0.007英寸还是0.021英寸损伤都当作对应故障类型。如果想做细粒度诊断也可以把损伤尺寸单独变成标签但分类难度会高不少本文按四分类来走。我强烈建议优先用驱动端加速度计、12kHz采样率的数据。每段信号长度足够长样本截取方便而且12kHz数据计算WSET时不会因为采样点太多导致内存紧张。48kHz数据可以作为算法稳定性验证但日常迭代用它明显更慢。下面是常用数据说明速查表参数内容加速度计位置驱动端DE采样率12000 Hz负载档位0 / 1 / 2 / 3 HP故障模式正常(N)、内圈(IR)、外圈(OR)、滚动体(B)损伤尺寸0.007 / 0.014 / 0.021 英寸推荐窗长1024 或 2048 点滑窗步长窗长的50%重叠50%2.2 Matlab读取与滑窗截取数据下载下来是多个.mat文件文件名按工况编号排列。先用load读入变量里一般包含DE、FE和BA等通道。以驱动端数据为例取DE时间序列即可。Matlab读入代码大概是这样的% 以某一段驱动端数据为例 data load(xxx.mat); fs 12000; signal data.DE; % 时间序列 signal signal(:); % 统一成列向量原始信号长度动不动几万点不可能整段喂给网络。需要滑窗截取把长信号切成若干固定长度样本。推荐窗长1024点或2048点原因后面3.1节展开。我通常用一个简单的函数完成截取直接输出Matlab训练需要的四维数组function samples extractWindows(x, winLen, step) x x(:); N floor((length(x) - winLen) / step) 1; samples zeros(winLen, 1, 1, N); for k 1:N idx (k - 1) * step 1; samples(:, 1, 1, k) x(idx : idx winLen - 1); end end输出维度是[winLen, 1, 1, N]分别对应高度、宽度、通道数、样本数这样可以直接喂给trainNetwork不用额外转格式。每次调用函数前把原始信号里跳变的头尾去掉避免引入不自然的高幅值脉冲。2.3 数据划分与标签设计的核心原则这是整个流程里最容易被忽略、也最影响真实结果的一步。很多人直接把同一条长信号的前面部分做训练、后面部分做验证准确率看起来高得离谱但换到新工况数据就露馅。错在信息泄露相邻滑窗之间高度重叠模型相当于已经见过验证集样本的内容。我实测下来最稳妥的做法是“按信号文件划分”。每个工况文件内截取出来的窗口要么全部进训练集要么全部进验证集或测试集绝不让同一条连续信号横跨数据集边界。四分类时每个类别各抽200到300个窗口样本量已经足够撑起一个小型CNN。标签用categorical类型labels categorical(condition, {N, IR, OR, B});然后在划分之前先randperm打乱信号文件顺序保证各故障类型在训练、验证、测试中分布均匀。实际代码里我习惯先按类别抽取固定数量的窗再统一shuffle避免某些类样本多、某些类样本少导致训练不均衡。3. WSET时频能量谱的Matlab实现与可视化分析3.1 调用wsst生成同步挤压系数Matlab小波工具箱提供了现成的wsst函数一行命令就能得到同步挤压系数矩阵[sst, fr] wsst(signal, fs);其中sst是复系数矩阵每个元素对应某个时间位置、某个频率上的小波同步挤压系数fr是该矩阵每一行对应的频率轴单位Hz。矩阵大小默认和时间点、频率分辨率相关具体由信号长度和内部参数决定。对这组系数取模平方就得到能量谱E_abs abs(sst).^2;轴承故障信号的能量动态范围极大。转频附近、高倍谐波区域能量很强冲击引起的宽带谐振区域能量相对弱一些直接用原始能量谱训练模型会被少数强能量像素带偏。我建议做一次对数压缩E_log log1p(E_abs);log1p处理的是log(1 E)既保留小能量区域的相对差异又不会让强能量点过冲。这一步可以明显缩短训练收敛时间。有一点要提前说清楚wsst的计算量不低。对一段十几万点的信号直接算整个时频图内存和时间都吃不消。我的做法是先滑窗截出2048点小片段再对每个片段单独执行wsst。2048点在12kHz采样率下对应约0.17秒信号能覆盖足够多的轴承转频周期。而且片段短计算时间从几十秒级降到秒级训练集迭代速度完全可接受。3.2 从能量谱到CNN输入图尺寸、归一化与避坑同步挤压系数矩阵的形状不是规则正方形直接喂给CNN通常不行。需要先resize到统一尺寸再归一化。经验值用64×64信息保留足够计算量也小最大用到128×128就差不多了再大只会显著拖慢训练并不稳定提升准确率。img imresize(E_log, [64, 64]); img mat2gray(img); % 线性归一化到[0,1]转成四维数组后标签与图像一一对应trainImages zeros(64, 64, 1, numSamples); for i 1:numSamples trainImages(:, :, 1, i) imresize(..., [64 64]); end不少初学者习惯把时频图画成彩色图再保存成jpg然后用三通道彩色图训练。我给的建议是不要这么做。第一jpg是有损压缩会模糊能量谱边缘第二用jet等colormap映射出来的伪彩色会引入额外的颜色边界这些边界和故障类别没有直接物理关系CNN可能学到的是颜色映射产生的伪影。我实测用原始灰度单通道训练分类边界更干净模型泛化也更好。如果确实希望可视化给别人看可以把64×64矩阵丢给imagesc加colormap去展示展示归展示训练归训练两个流程分开。3.3 完整可视化分析步骤时域、频谱、时频对比可视化不只是为了“画张好看的图”它能帮你发现预处理是否出错、特征是否有效。我做完一个实验必看下面这一套图第一张图是原始时域波形。把正常、内圈故障、外圈故障、滚动体故障四条信号放在同一个figure里用subplot排成四行观察幅值包络和冲击节奏是否合理。滚动体故障的冲击往往不如内圈外圈那么规律如果时域图里一点冲击感都看不到应该先怀疑数据读错了而不是急着调模型。第二张图是FFT频谱。主要看边带结构和共振峰分布确认没有明显的直流偏置或传感器噪声异常。这里不是要手工提取频率而是要建立“信号大概长什么样”的直觉。第三张图是关键同一段信号的CWT时频图与WSET时频图并排对比。Matlab里简单写就行subplot(1,2,1); imagesc(t, f_cwt, abs(cwt(signal, fs))); % CWT axis xy; title(CWT); subplot(1,2,2); imagesc(t, fr, E_log); % WSET axis xy; title(WSET);肉眼应该能明显看出WSET的谱线更聚焦能量条带更窄CWT版本则是一团糊状能量带。这一步如果对比不明显多半是wsst参数选择问题可以尝试调整voicePerOctave提高频率细分。第四张图是四类故障的典型WSET图谱拼成网格。你会看到内圈故障的时频图里有密集的、方向倾斜的谱条纹外圈故障相对更均匀滚动体故障则会出现不规则的断续亮带。这种差异就是CNN能学到的核心判别信息。每次训练新模型前我都会认真看一遍这些典型图心里先对“模型靠什么做决策”有个基本预判后面改网络结构才心里有数。4. 基于WSET输入的CNN训练与结果解读4.1 轻量CNN模型结构搭建网络结构不需要太复杂。WSET输入图本身已经包含了很强的时频特征网络太深反而容易过拟合到训练集的噪声细节。我用的是一个非常轻量的CNN两个卷积块、一个全连接层、一个分类层。layers [ imageInputLayer([64 64 1], Normalization, none) convolution2dLayer(3, 16, Padding, same) batchNormalizationLayer reluLayer maxPooling2dLayer(2, Stride, 2) convolution2dLayer(3, 32, Padding, same) batchNormalizationLayer reluLayer maxPooling2dLayer(2, Stride, 2) fullyConnectedLayer(64) reluLayer fullyConnectedLayer(4) softmaxLayer classificationLayer ];这个结构的参数量在十万级别对于几千个训练样本来说很合适。顺便解释几个细节imageInputLayer设Normalization为none因为前面已经做过mat2gray不需要再来一次标准化。卷积核用3×3池化用2×2步长2是图像分类的经典配置。时频图像素之间没有严格的“横平竖直”结构小卷积核更容易提取局部纹理组合。batchNormalizationLayer放在卷积之后、ReLU之前。批归一化能把每层的特征尺度稳住对不同的负载工况数据会稳定很多训练曲线也不容易震荡。4.2 训练选项设置与过程曲线训练配置看起来简单但几个参数直接决定能不能收敛。我常用的配置如下options trainingOptions(sgdm, ... MaxEpochs, 30, ... MiniBatchSize, 64, ... InitialLearnRate, 0.01, ... LearnRateSchedule, piecewise, ... LearnRateDropFactor, 0.1, ... LearnRateDropPeriod, 10, ... Shuffle, every-epoch, ... ValidationData, {valImages, valLabels}, ... ValidationFrequency, 20, ... Plots, training-progress, ... Verbose, false); net trainNetwork(trainImages, trainLabels, layers, options);几个参数的经验是初始学习率0.01对这个小网络刚刚好太大会出现Loss乱跳太小则收敛极慢MiniBatchSize取64显存一般的显卡也扛得住在每个epoch开始前做shuffle能避免同一个连续窗序列反复按同一顺序出现对稳定验证集准确率有帮助。训练过程中打开training-progress图重点看两条曲线的关系。正常情况训练Loss和验证Loss同步下降验证准确率逐步爬升。如果训练Loss降到了很低但验证Loss开始反向上升基本可以判定过拟合优先考虑减小网络尺寸或加Dropout层如果两条曲线一起震荡则先检查学习率是否过大。4.3 混淆矩阵、误分类与改进方向训练完成后用测试集做最终评估YPred classify(net, testImages); figure; confusionchart(testLabels, YPred);混淆矩阵是最直接的诊断结果工具。对角线越亮模型越可靠。我实际跑四分类时最容易出现的误分类组合是滚动体故障被分到内圈或外圈类别。这是因为滚动体故障的冲击位置会在滚动体旋转过程中不断变化冲击特征带有很强的相位调制和内外圈故障的谱结构有时确实很像。遇到这种情况我的处理顺序是先把滚动体故障对应的WSET图再翻出来看一遍确认它和某个其他类别的图谱是否真的接近如果接近考虑把窗长从1024拉到2048让更长的冲击序列把周期特征暴露出来如果还不行再对滚动体故障类别的样本做增强比如轻微叠加高斯噪声、随机裁切不同起始点。实测下来“加长窗长”往往是滚动体召回率恢复的最有效手段。5. 常见问题与故障诊断实操经验5.1 高频问题速查表下面这张表是我在过去几轮实验里反复被问到、或者自己踩过的问题整理成速查表方便对照排查现象可能原因排查方法wsst提示找不到函数未安装Wavelet Toolbox在Matlab命令行用ver确认工具箱状态运行时内存不足信号太长直接算WSET先滑窗截成2048点再逐段计算训练Loss不降未做对数压缩或归一化检查是否用了log1p和mat2gray验证准确率虚高同一条信号两端被分进不同数据集改为按信号文件划分数据验证Loss持续上涨网络过拟合减小网络宽度、降低epoch、增加Dropout彩色图训练反而变差jpg压缩和伪彩色引入伪影改用灰度单通道输入时频图谱线过粗不锐利voicePerOctave设置太低适当提高该参数获得更细频率轴各类样本数量差异大某类别的滑窗数量明显偏少按类别抽取同样数量的样本再训练5.2 几条来自实操的补充经验第一训练集样本不是越多越好而是要覆盖足够多的工况变化。把0HP、1HP、2HP、3HP四档负载全部混进训练集模型确实能学会负载无关的故障特征但如果训练集只放了0HP测试时改成3HP准确率掉下去是正常的这是任务本身的物理特性决定的不是网络结构有问题。第二检查数据标准化时切记不要对整个数据集的均值和方差做全局标准化。因为每条样本来自不同工况全局标准化会把工况差异当作主要变化方向反而掩盖故障特征。我坚持每个样本单独归一化到[0,1]这在工程上最简单也最稳定。第三做WSET可视化时频率轴建议截断到6000Hz以下。12kHz采样的奈奎斯特频率是6000Hz但真正对轴承冲击敏感的结构谐振往往落在低频到中频段画全频段会浪费大量显示区域也会让CNN把注意力分散到无用的服从噪声分布的频谱区域。我个人在实际操作中最感到值回票价的习惯是每次跑完实验先把典型样本的WSET图重新打开看一遍再决定要不要动网络结构。很多“模型效果不佳”的问题根源不在网络而在输入时频图里有明显的坏样本或者标签错位。如果你现在正卡在准确率提不上去的瓶颈我建议你先别急着改卷积层宽度把自己数据里每个类别的WSET图并排画出来看一眼把数据划分规则重新审一遍可能比我在这里推荐任何网络结构都管用。这个流程跑顺之后后面哪怕换成小样本场景或者继续往迁移学习方向走整个底座也是稳的。