STFT-CNN-BiGRU在旋转机械故障诊断中的完整实践与Matlab实现

我断断续续做了几年旋转机械故障诊断的项目最深的体会是问题往往不在算法本身而在特征到底怎么提。以前靠频谱分析、包络谱凭经验挑频带后来流行手工特征加传统分类器又陷入特征工程的无底洞。直到把短时傅里叶变换STFT和卷积神经网络CNN、双向门控循环单元BiGRU串成一条完整的流水线我才真正感觉到既不用天天手工调特征、又能把时间维度信息充分利用起来的方案是真实存在的。这篇文章就把我基于STFT-CNN-BiGRU做故障诊断的完整过程记录下来包括整套思路为什么这么搭、时频图参数怎么定、卷积和循环网络结构怎么设计、Matlab代码怎么落地以及训练过程中踩过的那些坑。内容适合正在做故障诊断、PHM故障预测与健康管理相关课题的同学也适合想把深度学习这套东西真正用起来的工业信号处理工程师。文章里给的代码框架足够直接改一改就可以迁移到自己的数据上。1. 故障诊断为什么需要STFT-CNN-BiGRU这套组合1.1 传统诊断方法的三个死穴先说我以前的做法。传统的故障诊断流程基本是这样采集振动信号做傅里叶变换看频谱找到某个故障特征频率处的峰值再根据峰值大小判断故障类型。看起来简单直接但真拿到工程数据就知道问题在哪。第一个死穴是不稳定工况下频谱根本不好用。设备转速一变特征频率跟着漂移频谱图上能量到处扩散肉眼根本看不出来。第二个死穴是对早期故障不敏感。微弱故障特征往往被噪声淹没单纯看幅值规律性很差。第三个死穴最要命就是特征提取完全依赖人工经验换了设备、换了工况之前总结的经验规则可能全部作废。后来用过一维CNN直接处理原始波形效果比传统方法好一些但有个问题一直绕不开网络要自己从原始时域信号里学会频率这个概念等于把信号处理的工作全部交给神经网络数据量小了根本学不出来。工业现场能拿到的高质量故障样本通常只有几百到几千条这种纯数据驱动的路子很容易过拟合。1.2 三条链路各自解决了什么STFT-CNN-BiGRU这套组合的思路说白了就是让信号处理、空间特征提取、时序建模三段各司其职。STFT负责把一维振动信号变成二维时频图。它做的事情很简单把整段信号切成一帧一帧的小段每一段做一次傅里叶变换然后把所有帧的频谱按时间排列起来。这样得到的时频图横轴是时间纵轴是频率像素亮度代表该时刻该频率处的能量强度。故障特征在这种表示里就显现出来了。CNN负责在时频图上做空间特征提取。卷积核能在局部区域自动学习能量团块的形态比如某一频带内的周期性冲击会形成一条条水平亮线某类故障会在特定频段出现能量集中。这些模式不需要人工定义CNN能自己找出来。BiGRU负责处理时间轴上的依赖关系。时频图本质上是一个时间序列故障信号在时间轴上往往有一种渐进演化的过程比如早期轻微磨损的冲击能量逐渐增强、频率成分逐渐变化。BiGRU能沿着时间方向捕捉这种变化趋势而且双向结构意味着它既考虑了过去的信息也兼顾了未来的信息诊断决策的依据比单向模型更完整。从整体来看这套组合解决的核心问题就是用STFT的物理先验降低CNN的学习难度再用BiGRU补上CNN不擅长处理的时序全局信息最终得到一个小样本下依然能稳定收敛的诊断模型。2. 方案选型解析为什么是STFT、CNN和BiGRU而不是别的2.1 时频分析方法那么多为什么选STFT这一步我纠结过很久。时频分析工具不止STFT一种常用的还有小波变换、Wigner-Ville分布、经验模态分解EMD。做一个简单对比就清楚了。方法核心思路优势劣势工程实现成本短时傅里叶变换分帧后逐帧FFT参数少、速度快、结果直观窗长固定时频分辨率固定极低连续小波变换变尺度小波基与信号内积低频高频分辨率可变尺度参数难选计算量大中Wigner-Ville分布信号瞬时自相关函数的傅里叶变换时频聚集性最好交叉项干扰严重中EMD将信号分解为本征模态函数自适应性强模态混叠、端点效应高工程上选了STFT核心原因是可控和稳定。窗长确定了时频分辨率就确定了同一个数据集里每张图都用完全相同的参数生成后续模型不会因为特征图尺度不一致出问题。小波变换虽然理论上更优雅但母小波类型和尺度参数的选择本身又是一个新的调参难题而且连续小波变换计算量明显更大训练数据一多预处理时间成本直接上天。2.2 CNN为什么是时频图上最合适的特征提取器时频图本质上是图像CNN本来就是为图像设计的。它相比全连接网络有两个关键优势。第一个是局部感受野。卷积核每次只看一个小邻域比如3x3的区域这恰好对应某个短时间窗口内、某个窄频带内的能量变化。局部模式在图像的不同位置重复出现时卷积核可以共享参数、分别检测。第二个是权值共享。一个卷积核在整张图上滑动参数量大幅减少这对小样本训练非常关键。实际做的时候我发现时频图跟自然图像有一个重要区别它没有平移不变性的明显意义。图像里的猫挪个位置还是猫但时频图里的能量团挪到另一个频率位置就代表完全不同的故障模式了。所以设计时频图处理的CNN时垂直方向频率轴的平移不能当作增强手段来用否则会把模型搞糊涂。这一点我在后文的数据增强部分还会专门说。2.3 BiGRU出现的逻辑特征图仍然是时间序列CNN把时频图提炼成高维特征图之后很多人会直接展平送全连接层分类我当时就是这么干的效果也不错。但有个信息被浪费了就是特征沿时间轴的顺序关系。故障信号的时间序列是有结构的。齿轮裂纹从萌生到扩展冲击能量会越来越强轴承内圈故障在负载变化时特征频率周围会有明显的调制边带这种调制模式随时间展开。CNN在池化过程中会不断压缩空间分辨率时间轴信息逐渐模糊。BiGRU就是在这个位置补位的把CNN输出的特征图按时间方向切成帧序列送入双向循环网络让模型在更高语义层次上理解能量随时间的演化。选择GRU而不是LSTM是综合考虑参数量和效果的结果。GRU比LSTM少一个门控结构参数量更小训练速度更快在中等规模数据上的表现跟LSTM几乎持平。加上双向之后模型能同时利用当前时刻前和后两个方向的上下文信息对瞬态冲击和趋势特征的判别能力都更强。3. 数据预处理STFT时频图参数怎么定直接决定诊断上限3.1 一段振动信号如何变成一张时频图这一步是整个流程的地基。STFT的原理可以这么理解用一个固定长度的窗函数从信号头部开始截取一段数据对这一段做FFT得到频谱然后窗函数向前滑动一个步长再截取、再做FFT如此循环直到覆盖完整段信号。在Matlab里最直接的实现就是用spectrogram函数它帮你把分帧、加窗、FFT、拼接都封装好了。我用的核心代码是这样的% 单段信号的STFT时频图生成 fs 20000; % 采样率根据实际传感器设置 winLen 1024; % 窗长 noverlap floor(0.75 * winLen); % 重叠率75% nfft 1024; % FFT点数 [~, ~, T, P] spectrogram(data, hamming(winLen), noverlap, nfft, fs); % 转成dB刻度加eps防止log0 S 10 * log10(P eps); % 归一化到0-1灰度范围 S mat2gray(S);处理完得到的是一个二维矩阵行对应频率列对应时间每个元素是对应时频点的能量强度。之后把矩阵resize到统一尺寸作为CNN的输入。3.2 窗长、重叠率、窗型到底怎么选STFT有一个核心矛盾窗长越长频率分辨率越高但时间分辨率越低窗长越短时间分辨率越高但频率分辨率越差。机械故障信号通常需要同时兼顾频带的精细结构和冲击事件的时间位置所以平衡是关键。拿上面代码里的参数举例采样率20kHz窗长1024点那么一个窗的时间长度是1024/20000 51.2ms频率分辨率是20000/1024 ≈ 19.5Hz。如果换成4096点窗频率分辨率能到4.9Hz但一个窗覆盖205ms冲击瞬态的时间定位能力就差了。我的经验是先确认你要关注的最高特征频率有多高窗长取一个让频率分辨率小于特征频率间隔的值。比如要区分两个相距30Hz的特征成分频率分辨率至少要达到15Hz左右对应窗长约1400点左右取2的幂就是2048。如果只做模式识别并不需要精确读出峰值频率1024点窗往往已经够用。重叠率影响时频图的时间轴密度。75%重叠意味着窗每次只滑动25%的长度时间帧更密集图像纵向信息更丰富代价是计算量变大和相邻帧冗余度升高。实际项目中我一般固定用75%这是信号处理领域一个比较常规的折中值。窗型上Hamming窗是我的首选。矩形窗旁瓣太高容易引起频谱泄漏Kaiser窗参数beta调起来麻烦Hamming和Hann性能接近但Hamming的旁瓣略低一点处理机械振动信号时背景噪声压制效果更好。3.3 时频图尺寸统一和样本划分的细节不同样本长度不同STFT出来的矩阵列数就不同CNN输入尺寸必须一致。我的做法是固定一个样本长度比如每段信号取4096点统一做STFT后resize到224x224方便后续网络输入。% 批量生成训练样本 numSamples length(signalCell); H 224; W 224; allData zeros(H, W, 1, numSamples); for i 1:numSamples sig signalCell{i}; [~, ~, ~, P] spectrogram(sig, hamming(1024), floor(0.75*1024), 1024, fs); S mat2gray(10 * log10(P eps)); allData(:, :, 1, i) imresize(S, [H, W]); end labels categorical(classIdx);从每段原始信号里取固定长度也有讲究。我当时做的是把长信号切成长度4096点的片段保证每个片段至少覆盖几个到十几个转频周期这样单张时频图里就包含足够多的冲击事件CNN能学到完整的故障模式。切的时候还要注意同一个连续时间段切出来的两个片段不能同时被分到训练集和测试集否则会造成数据泄漏测试指标虚高。4. 网络结构设计CNN部分与BiGRU部分的参数细节4.1 整体网络结构一览我的网络结构可以分成四段输入层、CNN特征提取段、BiGRU时序建模段、分类输出段。具体参数如下表。层级模块参数配置1图像输入224x224x12卷积块1Conv3x3, 32通道Same PaddingBNReLUMaxPool 2x23卷积块2Conv3x3, 64通道Same PaddingBNReLUMaxPool 2x24卷积块3Conv3x3, 128通道Same PaddingBNReLUMaxPool 2x25特征重排展平并按时间轴重组为序列6BiGRUhidden_size128双向输出最后时刻隐藏状态7分类头FC 512ReLUDropout 0.3FC numClassesSoftmax三个卷积块逐步把通道数从32升到128空间尺寸从224x224降到28x28。这么设计的原因是浅层卷积捕捉边缘和局部纹理深层卷积能组合出更抽象的故障模式。通道递增的规律是深度学习里常见的套路信息逐层抽象参数量却不会失控。4.2 卷积模块设计的几个关键决策第一卷积核尺寸选3x3而不是5x5或7x7。3x3是图像处理里验证过的性价比之王两个3x3卷积串联感受野等同5x5但参数量更少、非线性更强。第二每个卷积后面都跟着BatchNorm层。BatchNorm能让每层输入分布稳定大幅加速收敛对故障诊断这种数据分布差异较大的任务尤其重要。第三池化用2x2最大池化而不是平均池化。最大池化保留了区域最强烈的响应对应时频图里能量最强的冲击成分这恰好是故障诊断最关心的信息。4.3 从CNN特征图到BiGRU的维度转换这一环节是实现的难点。CNN输出的特征图是四维的经过三个池化层后变成28x28x128。我不直接展平而是把时间轴单独留下作为序列方向% 将CNN输出的特征图重排为序列数据 % 原始尺寸 [H, W, C, N] 这里 H28, W28, C128 % 把W当作时间方向每一帧特征拼接所有通道 seqInput permute(features, [4, 2, 1, 3]); % 变成 [N, W, H, C] seqSize [W, H * C]; % 序列长度28每帧特征维度3584这一变换的物理含义很直观时频图的每一列对应一个时间点CNN把这一时刻的局部频谱信息编码成一个高维向量BiGRU再按时间顺序读这些向量。Matlab里接BILSTM层时需要把数据整理成对应格式我习惯在预处理阶段就构造好序列数据再用自定义训练循环传入这样逻辑更清晰、不容易出错。4.4 BiGRU内部机制简化理解GRU单元内部有更新门和重置门。更新门决定上一时刻的隐藏状态有多少信息保留到当前时刻重置门决定当前输入和过去状态如何组合。双向GRU就是把两个方向独立的GRU跑一遍每个时刻的最终输出是两个方向的拼接。设隐藏单元为128双向拼接后每个时刻输出256维。在故障诊断里这个机制的价值在于正向GRU能捕捉故障特征从无到有的发育过程反向GRU能感知能量衰退或边带变化的趋势融合后对故障类别判断更稳定。我实测下来在典型故障数据上双向比单向的准确率高3到5个百分点代价是训练时间多花将近一倍这笔账是划算的。5. Matlab代码落地从网络搭建到训练评估的完整流程5.1 网络搭建代码在Matlab里搭网络可以直接用layerGraph逐层拼接代码可读性很好。我用的是以下框架numClasses 4; % 假设4类故障 imageSize [224, 224, 1]; sequenceDim 3584; % 重排后的每帧特征维度 hiddenSize 128; lgraph layerGraph(); % CNN部分 cnnLayers [ imageInputLayer(imageSize, Normalization, none, Name, input) convolution2dLayer([3 3], 32, Padding, same, Name, conv1) batchNormalizationLayer(Name, bn1) reluLayer(Name, relu1) maxPooling2dLayer([2 2], Stride, 2, Name, pool1) convolution2dLayer([3 3], 64, Padding, same, Name, conv2) batchNormalizationLayer(Name, bn2) reluLayer(Name, relu2) maxPooling2dLayer([2 2], Stride, 2, Name, pool2) convolution2dLayer([3 3], 128, Padding, same, Name, conv3) batchNormalizationLayer(Name, bn3) reluLayer(Name, relu3) maxPooling2dLayer([2 2], Stride, 2, Name, pool3) ]; lgraph addLayers(lgraph, cnnLayers); % BiGRU与分类层 gruLayers [ sequenceInputLayer(sequenceDim, Name, seq_input) bilstmLayer(hiddenSize, OutputMode, last, Name, bilstm) dropoutLayer(0.3, Name, dropout) fullyConnectedLayer(512, Name, fc1) reluLayer(Name, relu_fc) fullyConnectedLayer(numClasses, Name, fc_out) softmaxLayer(Name, softmax) classificationLayer(Name, output) ]; lgraph addLayers(lgraph, gruLayers); % 连接CNN输出与BiGRU序列输入 lgraph connectLayers(lgraph, pool3, seq_input); analyzeNetwork(lgraph);这里需要说明的是CNN输出的四维特征图在Matlab中直接接到sequenceInputLayer时需要你自行完成特征重排并喂入数据。你可以把重排放在数据预处理阶段也可以写自定义层我习惯在训练循环里现做重排省去自定义层的调试时间。5.2 训练选项配置训练选项是一个容易被低估的模块。好的训练参数能让模型快速收敛差的参数会让loss反复震荡甚至直接发散。我在项目里用的是options trainingOptions(adam, ... InitialLearnRate, 1e-3, ... MaxEpochs, 50, ... MiniBatchSize, 32, ... ValidationData, {valData, valLabels}, ... ValidationFrequency, 10, ... Shuffle, every-epoch, ... L2Regularization, 1e-4, ... Plots, training-progress);学习率1e-3是Adam优化器在小模型上的常规起点。如果发现训练曲线震荡厉害先把学习率降到3e-4或1e-4再试。MiniBatchSize取32主要受显存限制batch越大梯度估计越稳但对内存的消耗也越大。ValidationFrequency设为10意思是在10次迭代后做一次验证集评估这样能及时看出过拟合趋势。5.3 测试评估与混淆矩阵训练完成后用测试集评估模型效果这一步必须严谨。我一般同时看整体准确率和混淆矩阵因为只在准确率高的时候就看结束可能漏掉某些类别样本较少导致分类偏斜的问题。% 测试集推理 predLabels classify(net, testData); % 准确率 acc mean(predLabels testLabels); fprintf(Test Accuracy: %.2f%%\n, acc * 100); % 混淆矩阵 confMat confusionmat(testLabels, predLabels); figure; confusionchart(confMat, categories(testLabels));混淆矩阵能直接告诉你哪些类别之间容易混淆。我遇到过的情况是滚动体故障和复合故障在时频图上的表现比较接近混淆矩阵里这两个类互相误判的比例明显偏高。这时候就需要检查STFT参数是不是丢失了某些关键频率成分的细节而不是急着改网络结构。6. 训练过程中的常见问题与排查技巧6.1 数据泄漏很多人都会踩的坑关于数据划分我见过太多人栽在这里。故障诊断的数据都是连续采样的同一电机在同一个转速下连续跑的振动信号前后时间段内的数据特征高度相似。如果随机划分训练集和测试集那么模型在训练时已经见过测试集高度相似的数据测试准确率虚高很正常真实场景里根本没法复现。正确做法是把采集时长分成几个独立时间段用前面时间段的数据训练后面的时间段的验证和测试。这既是工程习惯也符合模型评估的原则。同理切片段的时候也要保证同一段长信号切出来的相邻片段不要分别落入训练集和测试集。我通常会先用时间索引做分组再在组内切样本这样就能彻底避开数据泄漏问题。6.2 Loss不降或震荡时先查哪里模型训练过程中我最常遇到的三个问题及排查顺序是这样的。Loss根本不降先看数据预处理。确认时频图做完log和归一化后是不是存在大量全零或接近全零的图这类样本对训练没有意义。再确认标签跟样本是否对齐这个问题低级但极其常见。如果都没问题再怀疑学习率过高把学习率从1e-3降到1e-4试一轮。验证集Loss持续高于训练集Loss且差值越来越大这是过拟合的典型信号。先增加Dropout比例到0.5再做权重衰减最后考虑数据增强。如果增强前准确率已经很高说明模型容量可能过大减少卷积核数量反而更有效。训练Loss先下降后突然变成NaN多数原因是学习率过大导致梯度爆炸。解决办法是降低学习率或者给梯度设置裁剪。Matlab的trainingOptions里没有直接的梯度裁剪参数但可以通过自定义训练循环实现或者干脆把学习率调低并用更小的batch训练。6.3 类别不均衡与数据增强的注意事项故障诊断数据天然存在类别不均衡。正常工况样本最多各种故障样本相对较少。在不均衡严重的情况下直接训练会让模型倾向于把多数类别猜对。处理办法不外乎三种对少数类过采样、对多数类欠采样、改用加权交叉熵损失。工程上我最推荐先做过采样简单直接不会丢失信息。数据增强方面对于时频图我会做两种增强一是添加轻微高斯噪声模拟采集过程中的随机干扰二是沿时间轴做小幅度的随机平移。这两种增强都在逼模型学到更本质的频率特征而不是死记硬背图里的位置信息。但有一种增强我是坚决不做的就是沿频率轴翻转。在自然图像里上下翻转可能不影响语义但时频图不一样高频能量变到低频位置会直接扭曲故障语义。翻转后的图根本不是真实物理过程能产生的信号模型学到的东西也就没有了物理意义。这是时频图增强跟自然图像增强最大的区别。6.4 常见问题速查表现象可能原因处理手段准确率虚高但换数据崩溃数据泄漏划分不独立按时间段分组划分Loss不降预处理错误、标签错位、学习率过高检查时频图调低学习率验证集过拟合模型容量大数据量小加Dropout、权重衰减、数据增强Loss变成NaN梯度爆炸降低学习率减小batch某两类互相混淆时频图特征不够区分调整窗长检查特征频率段增加数据显存不足批量或图像尺寸过大减小batch减小输入尺寸7. 项目扩展与实际操作体会这套组合跑通之后我做了两次扩展效果都不错。第一次是把单传感器改成多传感器把不同测点的振动信号分别做STFT然后当作多通道图像输入CNN。相当于让网络同时看设备不同部位的状态诊断鲁棒性明显提高。第二次是在CNN和BiGRU之间插入一个简单的通道注意力模块自动给故障特征明显的频率通道加权重在诊断准确率上又提升了一点。我个人的体会是这类模型最耗时间的不是网络搭建而是STFT预处理和训练参数的磨合。网络结构直接抄现成的就能跑但窗长、重叠率、样本长度这些参数必须结合自己数据的采样率和特征频率来定没有任何一个通用值能解决所有问题。最后分享一个实用技巧训练之前先打印几张不同类别的时频图用眼睛扫一遍。如果不同故障类别之间的时频图肉眼就几乎分不出差异那后面再怎么调网络都白搭。预处理做到位这个模型的下限就已经很高了预处理偷懒再强的网络都兜不住。做故障诊断永远要把功夫花在信号本质的理解上而不是跟风追新模型。