BP神经网络做齿轮箱故障诊断：从振动信号特征提取到Python实战

简介基于反向传播BP神经网络的齿轮箱故障诊断MATLAB源码面向设备故障诊断与机器学习交叉应用场景适合机械工程、自动化等专业学生及从事状态监测的工程技术人员。资源包内仅含1个.m脚本压缩后仅约1KB以极简代码覆盖BP网络从输入层设计、隐藏层计算、误差反向传播到故障类型输出的核心流程便于直接阅读、修改和在MATLAB中快速运行验证。该脚本可作为齿轮磨损、断齿、轴承损伤等典型故障识别的基础分类模板也可与振动信号采集、傅里叶变换等特征提取环节衔接帮助理解非线性映射与自学习在工程诊断中的落地方式通过调整隐藏层节点数、学习率与迭代次数还能进一步适配不同工况数据在少量样本条件下验证网络泛化能力。目前已有196人学习下载适用于课程设计、故障诊断入门或算法对比实验尤其适合需要在MATLAB环境中快速验证BP诊断流程的读者。1. BP神经网络做齿轮箱故障诊断先回答“值不值得学”BP_net.rar 这个资源包在齿轮箱故障诊断圈子里流传了很久核心就是拿 bp神经网络 对齿轮箱振动信号做分类。很多工程师下载后第一件事是跑通 demo换成自己现场采的数据却立刻翻车原因多半不在网络而在特征和参数。这篇笔记先把原理讲明白再把一套能复现的 Python 流程给你从振动信号提特征、构造网络、训练诊断到部署时容易踩的坑。适合刚接手设备点检、需要给齿轮箱做状态识别或者正在复现文献模型的人。看完你会知道BP 不是最先进的算法但作为故障诊断基线它依然值得投入。2. 从振动信号到BP输入齿轮箱故障诊断特征提取的三种落地做法BP 网络不能直接吃原始波形至少不能直接吃上万个采样点。齿轮箱振动信号里混杂着啮合频率、边频带、随机冲击和环境噪声把原始时域波形直接灌进网络隐层节点数会被迫膨胀训练数据需求量也成倍增加。所以第一件事是把一段振动信号压缩成一组有物理意义的特征。2.1 为什么用振动信号而不是温度或声音齿轮箱故障诊断的首选信号源是加速度传感器采集的壳体振动。温度信号变化太慢齿轮点蚀或断齿发生到温度可观测往往已经进入严重损坏阶段声音信号受现场环境干扰极大同一条产线上不同机台的噪声混在一起很难分离出目标齿轮箱的故障成分。振动信号直接反映齿轮啮合状态冲击性的故障会在时域波形上产生周期性脉冲在频域上产生以啮合频率为中心的边频带这些特征能被 BP 网络有效学习。需要注意的是这套信号选择思路同样适用于轴承故障诊断只是轴承故障特征频率的计算方式不同。2.2 时域特征均值、峰值、峭度、波形指标怎么算时域统计特征计算成本最低适合嵌入式平台和在线监测。常用的有均值、均方根值、峰值、峭度、波形因子、峰值因子、脉冲因子和裕度因子。均值反映直流分量均方根值代表振动能量的总体水平峰值对突发冲击敏感峭度对早期故障非常敏感正常滚动表面振动近似高斯分布峭度接近 3一旦出现点蚀或剥落峭度会明显上升。下面这份 Python 代码可以一次性算出这些量import numpy as np def time_features(x): # x 是一段时域振动信号长度为 n mean_val np.mean(x) # 均值 rms np.sqrt(np.mean(x**2)) # 均方根值 peak np.max(np.abs(x)) # 峰值 std np.std(x) # 峭度四阶中心矩除以方差的平方 kurt np.mean((x - mean_val)**4) / (std**4 1e-12) # 波形因子均方根值与绝对均值之比 waveform rms / (np.mean(np.abs(x)) 1e-12) # 峰值因子峰值与均方根值之比 crest peak / (rms 1e-12) # 脉冲因子峰值与绝对均值之比 impulse peak / (np.mean(np.abs(x)) 1e-12) # 裕度因子峰值与方根幅值之比 margin peak / (np.mean(np.sqrt(np.abs(x)))**2 1e-12) return np.array([mean_val, rms, peak, kurt, waveform, crest, impulse, margin])注意每个除法分母都加了极小量 1e-12防止信号全为静默时除零。这八个特征量级差异很大后续必须做归一化否则 BP 网络会把大数值特征当成主变量小数值特征直接被淹没。2.3 频域特征FFT谱与边频带不学这个BP白搭齿轮箱故障诊断里频域特征比时域特征更有区分度。正常齿轮啮合时频谱上主要成分是啮合频率及其谐波当齿轮出现局部缺陷啮合频率附近会出现边频带边频带间隔等于故障轴的转频。只看时域峭度很难区分“断齿”和“点蚀”但频域的边频带分布能提供类别线索。常用的频域特征包括重心频率、均方频率、频率方差以及前几个主导频率幅值占比。def freq_features(x, fs): # x: 时域信号, fs: 采样率单位 Hz n len(x) # rfft 节省一半计算量只取正频率 spectrum np.fft.rfft(x) amp np.abs(spectrum) / n # 幅值谱 freqs np.fft.rfftfreq(n, 1/fs) # 对应的频率轴 total_amp np.sum(amp) 1e-12 # 重心频率幅值加权平均频率 fc np.sum(freqs * amp) / total_amp # 均方频率 msf np.sum((freqs**2) * amp) / total_amp # 频率方差频谱围绕重心频率的分散程度 vf np.sum(((freqs - fc)**2) * amp) / total_amp # 前 5 个主峰幅值占总幅值比例 top_idx np.argsort(amp)[::-1][:5] top_ratio np.sum(amp[top_idx]) / total_amp return np.array([fc, msf, vf, top_ratio])重心频率如果明显偏离正常状态说明啮合能量发生了迁移频率方差增大说明频谱结构变得分散。如果做的是行星齿轮箱啮合频率的边频带间隔通常等于行星架转频特征提取逻辑不变但故障频率公式要换。2.4 特征组合用Python写一个特征提取脚本实际工程中我会把时域和频域特征拼成一个复合向量作为 BP 网络的输入。下面这个函数对一段信号提取 12 维特征8 个时域 4 个频域每个训练样本就是这样一个向量。def extract_feature_vector(x, fs): tf time_features(x) ff freq_features(x, fs) return np.concatenate([tf, ff]) # 示例用法假设每段信号 2048 点采样率 5120 Hz # sample np.loadtxt(gear_normal_001.csv) # feature extract_feature_vector(sample, 5120) # print(feature.shape) # (12,)这里没有选原始时间序列作为输入而是用统计特征。优点是训练快、对振幅偏移不敏感缺点是会丢掉相位信息。如果你需要诊断的故障类型非常接近比如“齿面轻微磨损”和“齿面中度磨损”12 维统计特征可能不够后续可以继续加小波包能量特征或包络谱特征。但作为基线这套组合已经能区分大多数齿轮箱故障类别。3. 搭一个能跑的BP神经网络故障诊断模型结构、代码与参数很多故障诊断代码包里的 BP 网络都是三层结构网上流传的 BP_net.rar 里大概率也跑的是这类网络。三层 BP 的原理解释起来很直观输入层接收特征向量隐层做非线性变换输出层输出每个类别的概率。下面这份 bp神经网络python代码我手写了一个可运行的版本不依赖深度学习框架方便你看清每一步在做什么。3.1 三层BP网络为什么够用输入层、隐层、输出层的尺寸选择BP 神经网络的万能逼近定理告诉我们一个包含足够多隐层节点的单隐层网络可以逼近任意连续函数。齿轮箱故障诊断本质是个多分类问题12 维特征映射到 4 类或 6 类故障三层结构完全够用。输入层节点数等于特征维度比如上面的 12输出层节点数等于故障类别数例如“正常、磨损、断齿、点蚀”就是 4 个。隐层节点数没有一个解析解常见经验公式是取输入层和输出层节点数的均值附近再向上尝试。比如输入 12、输出 4隐层可以先取 8然后分别试 6、10、15用验证集精度决策。bp神经网络结构图大致是输入层 12 个圆中间隐层 8 个圆输出层 4 个圆全连接。3.2 训练代码反向传播、学习率、动量项下面是一个最小可用的 BP 网络实现包含了前向传播、反向传播、动量更新和 softmax 输出。注意这里用的是批量梯度下降每轮迭代用全部训练样本计算梯度梯度比随机梯度下降更平稳适合小规模故障诊断样本集。class BPNN: def __init__(self, n_in, n_hidden, n_out, lr0.01, momentum0.9): self.lr lr self.momentum momentum # 权重初始化不能全零否则隐层节点对称 self.w1 np.random.randn(n_in, n_hidden) * 0.1 self.b1 np.zeros(n_hidden) self.w2 np.random.randn(n_hidden, n_out) * 0.1 self.b2 np.zeros(n_out) # 动量缓存 self.vw1 np.zeros_like(self.w1) self.vb1 np.zeros_like(self.b1) self.vw2 np.zeros_like(self.w2) self.vb2 np.zeros_like(self.b2) def sigmoid(self, z): return 1 / (1 np.exp(-z)) def softmax(self, z): # 减去最大值防溢出 expz np.exp(z - np.max(z, axis1, keepdimsTrue)) return expz / np.sum(expz, axis1, keepdimsTrue) def forward(self, X): self.z1 X self.w1 self.b1 self.a1 self.sigmoid(self.z1) self.z2 self.a1 self.w2 self.b2 self.a2 self.softmax(self.z2) return self.a2 def train_one_step(self, X, y): # y 是 one-hot 矩阵形状 (样本数, n_out) out self.forward(X) m len(X) # 交叉熵损失对 z2 的梯度softmax cross-entropy 简化结果 delta2 (out - y) / m # 隐层误差乘权重转置再乘 sigmoid 导数 delta1 (delta2 self.w2.T) * self.a1 * (1 - self.a1) # 动量更新 self.vw2 self.momentum * self.vw2 - self.lr * (self.a1.T delta2) self.vb2 self.momentum * self.vb2 - self.lr * np.sum(delta2, axis0) self.vw1 self.momentum * self.vw1 - self.lr * (X.T delta1) self.vb1 self.momentum * self.vb1 - self.lr * np.sum(delta1, axis0) self.w1 self.vw1 self.b1 self.vb1 self.w2 self.vw2 self.b2 self.vb2学习率是最重要的超参数。经验上取 0.01 起手如果损失下降太慢再调到 0.05如果损失振荡就调回 0.005。动量项取 0.9 能让梯度方向平滑避免陷入局部极小点。权重初始化乘以 0.1 是控制初始输出不要太大否则 sigmoid 容易饱和。这份代码里的隐层激活函数是 sigmoid输出层是 softmax千万不要把输出层也写成 sigmoid——多分类问题的标签是 one-hot 编码softmax 才能让所有类别的概率之和为 1。3.3 数据归一化mapminmax的替代实现用过 MATLAB 的人都知道 mapminmax 函数把特征映射到 [-1, 1]。Python 里需要自己实现这个逻辑而且要特别注意训练集和测试集必须使用同一组 min、max 参数否则测试数据被各自归一化后分布完全对不上。归一化后不同量纲的特征才能被 BP 网络平等对待。def normalize_with_params(X, min_valNone, max_valNone): # X: 二维数组行是样本列是特征 X np.asarray(X, dtypefloat) if min_val is None: min_val np.min(X, axis0) if max_val is None: max_val np.max(X, axis0) # 避免出现零跨度特征 denom max_val - min_val 1e-12 X_norm (X - min_val) / denom * 2 - 1 return X_norm, min_val, max_val # 训练阶段 X_train_norm, mn, mx normalize_with_params(X_train) # 测试阶段必须传训练阶段算出的 mn, mx X_test_norm, _, _ normalize_with_params(X_test, mn, mx)特别提醒如果某个特征在所有训练样本里都是常数max_val - min_val 等于 0加 1e-12 只是兜底实际应该检查特征是否有区分度。我在实际项目里见过有人把归一化函数写在数据读取函数里导致训练和测试各自调用结果验证集准确率 0.9换了新数据直接掉到 0.5。这是血泪经验。3.4 训练集/验证集划分与早停不要把所有数据都拿去训练否则模型会死记硬背。我会把每个类别的样本按比例随机拆成 70% 训练、15% 验证、15% 测试而且要保证每类样本在三个集合中比例一致。早停是防止过拟合最实用的一招每轮训练后计算验证集损失如果连续若干轮没有下降就停止。# 假设已经得到 X_norm 和 y_onehot from sklearn.model_selection import train_test_split X_tr, X_tmp, y_tr, y_tmp train_test_split( X_norm, y_onehot, test_size0.3, stratifyy_label, random_state42) X_va, X_te, y_va, y_te train_test_split( X_tmp, y_tmp, test_size0.5, stratifyNone, random_state42) model BPNN(n_in12, n_hidden8, n_out4, lr0.01, momentum0.9) best_loss np.inf wait 0 early_stop_rounds 15 for epoch in range(500): model.train_one_step(X_tr, y_tr) val_out model.forward(X_va) val_loss -np.mean(np.sum(y_va * np.log(val_out 1e-12), axis1)) if val_loss best_loss: best_loss val_loss wait 0 else: wait 1 if wait early_stop_rounds: print(f早停在第 {epoch1} 轮) break早停阈值设置 15 轮比较稳妥。如果是小样本阈值可以设小点比如 10如果训练过程噪声大阈值设 20 更合适。stratify 参数用于保证划分后类别比例一致但第二次划分时 y_tmp 已经打乱想严格分层可以自定义函数这里为了简洁用随机划分实际数据量够大时影响不大。4. 模型评估与阈值设定诊断结果怎么看才靠谱训练完模型只是第一步真正难的是判断模型在未知数据上能不能用。很多网上流传的故障诊断代码只打印一个“准确率 98%”但那个数字是怎么算出来的、哪些类容易混淆、置信度低的时候该不该报警全都语焉不详。4.1 混淆矩阵、准确率与召回率准确率容易骗人尤其是故障样本占比很低的场景。齿轮箱故障诊断里正常状态样本可能占 80%故障样本只有 20%模型全部预测正常也能拿到 80% 准确率。所以必须看混淆矩阵和每一类的召回率。下面这个函数可以手工打印混淆矩阵def confusion_matrix(y_true, y_pred, n_classes): cm np.zeros((n_classes, n_classes), dtypeint) for t, p in zip(y_true, y_pred): cm[t, p] 1 return cm y_pred_class np.argmax(model.forward(X_te), axis1) cm confusion_matrix(y_label_te, y_pred_class, n_classes4) print(混淆矩阵行真实列预测) print(cm)从混淆矩阵里能看到哪些类互相混淆。比如“断齿”和“点蚀”如果经常互相误判说明特征里缺少能区分这两种故障的分量需要补特征而不是调网络复杂度和迭代次数。4.2 多故障类别下的softmax与one-hotBP 网络做多分类时输出层必须用 softmax标签必须用 one-hot 编码。比如“正常、磨损、断齿、点蚀”四类某个样本属于“断齿”标签向量就是 [0, 0, 1, 0]。代码里 forward 函数已经包含 softmax训练时的损失函数用交叉熵梯度计算在 train_one_step 中做了简化这是标准做法。如果你从别处拷来的代码里输出层是 sigmoid、损失是均方误差跑多分类只能得到一组 0 到 1 的值它们并不保证加和为 1归属判断会出问题。4.3 诊断决策概率输出如何转成报警实际部署时我不建议直接把最大概率对应的类别作为最终诊断。如果两个类别的概率分别是 0.48 和 0.47最大概率类的胜出纯属偶然这时候给维修人员报“断齿”会误导人。更稳妥的做法是设置置信度阈值低于阈值就输出“未知/可疑”要求重新采样或人工复核。prob model.forward(X_te_norm) pred np.argmax(prob, axis1) conf np.max(prob, axis1) # 置信度低于 0.6 的样本标记为 -1未知类 for i in range(len(pred)): if conf[i] 0.6: pred[i] -1这个 0.6 不是拍脑袋定的。我会拿验证集画出置信度分布观察正常样本和故障样本对应概率的差距。如果大部分样本置信度都在 0.9 以上阈值设 0.7 就能筛掉很多模糊样本如果置信度普遍偏低阈值降 0.5避免漏报。阈值设置的原则是宁可多报“需人工确认”也不要把低置信度故障错报成正常。5. BP_net调试避坑5个让故障诊断翻车的常见问题网上资源包里能跑通拿到自己现场数据就失效问题大多不在 BP 本身而在数据组织和特征处理。这一章我列五条高频坑每一条都是我在项目里真实遇到过的。5.1 样本不均衡正常样本远多于故障模型全预测正常现象训练后准确率高达 90%但查看混淆矩阵发现所有故障类全被预测成“正常”故障检出率为 0。原因正常样本占 90%模型只要把所有样本都判成正常损失已经很低梯度对少数类的修正被淹没。解决先统计每类样本数对故障类做过采样重复采样或滑动窗口重叠切分或者对正常类降采样。另一个有效手段是给损失函数加类别权重让少数类样本的误差在反向传播时乘以更大的系数。最简单做法是先用 class_weight 思想在 train_one_step 里把 delta2 按类别权重缩放。# 示例每类样本数 counts np.bincount(y_label) weights 1.0 / (counts.astype(float) 1e-12) # 训练时把每个样本的梯度乘上其类别权重 sample_weight weights[y_label] # 在 train_one_step 的 delta2 计算中乘上 sample_weight 的列形状5.2 滑动窗口切分不当训练集和测试集重叠导致信息泄漏现象同一段连续振动信号里先用前 2048 点做训练样本再用紧挨着的后 2048 点做测试样本测试准确率高得离谱但现场实测不灵。原因相邻窗口来自同一段运行状态高度相关测试样本几乎就是从训练样本“平移”出来的模型实际上记住了噪声模式而不是故障特征。解决滑动窗口切分前先按时间把信号分成若干段比如一段信号 10 秒前 7 秒用于训练、后 3 秒用于测试或者在不同运行转速/不同载荷下分别采集训练和测试数据。窗口之间可以重叠但训练集和测试集必须来自不同时间段。5.3 特征量级差几个数量级网络训练不收敛现象训练损失忽大忽小甚至变成 NaN准确率永远停在某个固定值。原因峭度可能高达几十重心频率也可能上千直接喂给网络梯度计算时大数值特征对应的权重更新幅度巨大反复振荡。解决严格按 3.3 节做最小最大归一化并且用训练集算参数、测试集沿用同一参数。我见过有人对每列特征单独归一化时把一列常数值数据加 1e-12导致其归一化后约等于 -1这列特征基本失效所以还要检查特征是否有方差。5.4 隐层节点数抄别人的自己的数据精度上不去现象照搬网络达人的结构隐层节点取 30自己只有 200 个样本训练后验证损失居高不下网络还明显过拟合。原因隐层节点越多参数量越大需要的训练样本呈指数级增长。齿轮箱故障诊断本来就是小样本场景用 8 到 15 个隐层节点往往比 30 个效果好。解决从 6 开始每次增加 2画验证集准确率曲线。选择曲线峰值对应的节点数再往大走一格看看是否继续上升。如果上升明显说明特征复杂可以加如果持平或下降说明已经过拟合。5.5 早停阈值设得太松模型还没学到就停了现象训练日志显示第 20 轮就早停验证损失还在下降但连续 15 轮没创新低最终准确率只有 70%。原因早停阈值设置成“验证损失不下降立即停止”但训练初期损失下降慢很容易触发连续不下降的误判。算法学习本来就有平台期平台的持续时间可能超过阈值。解决把 early_stop_rounds 设成 15 到 30并且只在训练轮数超过 100 后才启用早停。另外记录最佳验证损失对应的权重而不是最后一轮权重。我在项目里习惯额外存一份 best_model 的权重副本训练结束后还原它。# 保存最佳验证损失时的权重 if val_loss best_loss: best_loss val_loss best_w1 model.w1.copy() best_b1 model.b1.copy() best_w2 model.w2.copy() best_b2 model.b2.copy() # 训练结束后恢复 model.w1, model.b1 best_w1, best_b1 model.w2, model.b2 best_w2, best_b26. 进阶用滚动验证和特征扰动确认模型的真实诊断能力模型部署前我一般还要做两个验证。第一是滚动验证把采集到的连续振动信号按时间顺序切成 5 段每次拿其中 4 段训练、1 段测试轮流五次。这样做比随机划分更贴近真实工况因为设备状态随时间漂移同一时段的数据相关性高滚动验证能暴露跨时段泛化能力。第二是特征扰动验证对每个特征列逐个加上随机噪声观察模型输出概率的变化幅度。变化越大的特征对诊断结果影响越强。# 特征扰动对第 j 列特征加正态噪声比较预测类别变化比例 import numpy as np base_pred np.argmax(model.forward(X_val_norm), axis1) stability np.zeros(n_features) for j in range(n_features): X_perturbed X_val_norm.copy() X_perturbed[:, j] np.random.normal(0, 0.05, sizelen(X_val_norm)) pert_pred np.argmax(model.forward(X_perturbed), axis1) stability[j] np.mean(pert_pred base_pred) print(各特征保持分类不变的比例:, stability)如果某个特征加噪声后正确率大幅下降说明网络严重依赖它如果所有特征加噪声后预测都不变说明模型可能学成了“输出恒为某一类”这时候要回到混淆矩阵检查。我现在的习惯是在部署前先把每个故障模式的召回率单独打印出来而不是只看总准确率。这个习惯救过我很多次因为总有某个故障被平均分数掩盖。BP 网络做齿轮箱故障诊断做到这一步已经具备上线条件剩下的就是持续用现场数据重训把新工况吸进模型里。希望帮到你。本文还有配套的精品资源点击获取