☰
2016年PHM数据集轴承退化多方法分析:特征工程与模型避坑指南
2026/10/22 8:17:39 网站建设 项目流程

简介:面向预测性维护入门者与数据分析学习者,一套围绕2016年PHM数据集的性能特征分析资料包,完整呈现多方法分析过程,可用于理解设备故障预测中的数据预处理、特征划分与模型对比思路。压缩包共759个文件、约17.84MB,以563个CSV数据文件为主体,辅以126个zbak备份、55张PNG可视化图表、5个Python脚本及说明文档,便于对照原始数据、分析代码与结果图进行复盘。目前已有76人学习,适合个人研究或课程实践参考。资源亮点在于涵盖三模态与双模态划分实验,通过不同特征组合探索对故障预测的影响;同时附带README、许可证等文件,明确数据来源与合规使用边界。学习者可借此掌握统计分析与机器学习在PHM数据集上的应用路径,并借鉴作者对设备运行规律的多角度拆解方法,为后续独立开展预测性维护研究提供可操作范本。

1. 多方法分析2016年PHM数据集:为什么单一指标看不透轴承退化

2016年PHM数据挑战赛提供的轴承退化振动数据,一直是故障预测与健康管理领域用来验证剩余使用寿命算法的常用基准。“多方法分析2016年PHM数据集性能特征”这组词背后有两个关键点:一是要回答数据本身表现出哪些退化信号,二是要在时域统计、频段能量、模型预测多个维度上交叉验证,避免单一指标给出片面结论。很多刚接触这个数据集的人喜欢直接套用深度模型,结果发现验证集漂亮、测试集翻车,根本原因是对数据集的性能特征没有建立起整体认知。这篇文章适合做设备健康管理与故障预测的工程师,也适合需要快速上手这个数据集的在校研究者,全文从数据读取讲到残差追溯,每个步骤都能直接复现。

2. 数据读取与预处理:先把2016年PHM数据集的底数摸清

在进入特征工程之前,先把原始文件的格式、量纲和信号完整性检查一遍。这一步看似简单,却是之后所有分析能否成立的基础;数据读取阶段埋下的隐患,会通过模型结果以极其隐蔽的方式暴露出来。这个数据集通常包含多个监测周期的振动信号,每个周期对应设备某一段运行时间,读取时需要特别注意文件名排序与周期编号的对应关系。

2.1 批量读取监测周期,先做完整性与量纲检查

工程上最常见的数据组织方式是每个监测周期一个独立文件。文件格式可能是CSV,也可能是MAT;不管哪种,第一步都是把整批数据读入,用一个字典结构保持周期编号与原始信号的对应关系,不要边读边做数据处理,因为后续的滤波、切窗操作都要反复回到原始数据。

import numpy as np import pandas as pd from pathlib import Path data_dir = Path('./phm2016_run1') cycle_files = sorted(data_dir.glob('*.csv')) raw_series = {} for f in cycle_files: df = pd.read_csv(f) raw_series[f.stem] = df print(f.stem, df.shape, df.columns.tolist())

这里用pathlib的glob获取全部CSV文件,再逐个读取进字典。sorted()对文件名字符串排序可能带来“cycle_10排在cycle_2前面”的问题,所以若文件名包含前导零,要确认补零格式一致;如果格式不统一,建议用正则提取周期编号后按整数排序。

读取完成后需要马上检查三件事:一是列名是否统一,有些文件多一个时间戳列,有些没有,列数不一致会直接影响后面索引;二是量纲是否合理,加速度信号一般在±10g或±50m/s²这个数量级,如果出现10³的异常值,需要额外注意;三是文件长度是否一致,同一批次实验的采样点数应当接近,若某个文件明显短于平均值,可能记录不完整,应标记并决定是否剔除。

2.2 高频随机噪声、基线漂移和脉冲干扰的处理顺序

振动信号预处理有三个目标:去掉直流分量带来的均值偏移,滤掉低频基线漂移,清理传感器受电磁干扰产生的瞬时脉冲。三者顺序不能随意调整,常规顺序是先去除直流偏置,再做高通滤波,最后用滑动中值清理脉冲。

from scipy import signal fs = 25600 # 采样率,以数据集原始记录为基准 sos = signal.butter(10, 5, btype='highpass', fs=fs, output='sos') def preprocess_vib(x, fs=fs): """去直流 -> 高通 -> 中值滤波""" x = x - np.mean(x) x_hp = signal.sosfilt(sos, x) x_med = signal.medfilt(x_hp, kernel_size=5) return x_med

高通滤波采用10阶巴特沃斯滤波器,截止频率取5Hz。这里选5Hz是考虑到大多数工业轴承的转频在25Hz上下,5Hz截止可以滤掉温度漂移和安装预紧力变化带来的基线,又不会伤到转频本身。如果在实际数据上观察到低频段存在明显的周期性跳动,则需要调高截止频率,但不能超过转频的一半。

滑动中值kernel_size取5比较安全。取值太大会把真正的剥落冲击磨平,导致后续峭度特征失效;取值太小又对连续性脉冲毫无办法。处理完之后,我习惯截取一段信号对比滤波前后的波形和功率谱,确认原始频谱峰没有被滤掉,而不是直接跳到特征计算步骤。

2.3 按固定窗口切分,并把时间维信息保留下来

原始振动信号是连续波形,建模时一般把它切成固定长度的时间窗,每个窗口计算一组特征。窗口长度直接决定频率分辨率和特征序列的长度:窗口太短,FFT分辨不出故障的窄带特征;窗口太长,早期瞬态退化又会被平均掉。实践中常取0.5到2秒。

def slice_windows(x, fs=fs, win_len=1.0, overlap=0.5): win = int(win_len * fs) hop = int(win * (1 - overlap)) slices = [] start = 0 while start + win <= len(x): slices.append(x[start:start + win]) start += hop return np.array(slices)

上面的函数返回一个二维数组,每行是一个时间窗信号。overlap取0.5表示相邻窗口有50%重叠,为的是让特征序列更平滑,同时保留局部变化。需要特别说明的是,有重叠切分出的样本本来就是强相关的,后续建模时不应当做独立样本随机划分,否则会出现非常隐蔽的数据泄漏。正确做法是按周期分组,让同一周期的窗口只在训练或验证集合中出现。

提示:预处理后的窗口信号建议保存为npy文件,并在文件名中保留周期编号与窗口起始时间。后面做任何特征提取或模型验证,都从这份中间结果读取,既省时间也避免重复计算的偏差。

3. 特征提取:从时域到频域再到时频域,三路特征如何暴露退化痕迹

什么样的特征能反映数据集的性能特征?对轴承退化而言,核心是找到跟随损伤演化而变化的统计量,并要求它单调、抗噪、对早期缺陷敏感。时域、频域、时频域三条路线恰好对应不同的敏感性:时域指标擅长捕捉冲击能量变化,频域带能量擅长定位特定频段的共振,时频图辅助诊断异常事件发生的时刻。

3.1 时域指标:RMS与峭度组合,利用形态差异识别早期故障

时域指标的计算非常直接,但“直接”不代表“简单”。RMS是能量指标,适合反映整体磨损程度;峭度和峰值因子对局部冲击敏感,适合捕捉微观裂纹。在同一组特征里,这两类指标对退化阶段的响应曲线并不一致,这正是多方法分析需要的互补性。

def time_domain_features(x): rms = np.sqrt(np.mean(x**2)) peak = np.max(np.abs(x)) std = np.std(x) kurt = np.mean((x - np.mean(x))**4) / (std**4 + 1e-12) skew = np.mean((x - np.mean(x))**3) / (std**3 + 1e-12) crest = peak / (rms + 1e-12) return { 'rms': rms, 'peak': peak, 'kurtosis': kurt, 'skewness': skew, 'crest_factor': crest }

分母加1e-12是防御性写法。轴承在健康阶段振动幅值非常低,标准差接近零,如果不加保护,峭度计算会出现除零或极端值。实际经验是健康段峭度大约在3附近,一旦局部缺陷出现,峭度会迅速跳到10以上,随后在磨损扩大阶段下降到5左右。

RMS与峭度配合使用,可以识别退化阶段。RMS单调上升时,峭度没有明显变化,说明是均匀磨损;RMS平稳但峭度突然抬高,说明出现局部剥落。这些差异如果只用单一指标,很容易被漏掉。此外,峰值因子在早期冲击下迅速升高,中期又会回落,建议与峭度一起看,而不是只取其中一个。

3.2 频域带能量:用功率谱密度统计高低频段的能量迁移

频域特征从另一个视角补充“哪些频率被激发”。轴承局部故障会激励结构共振,产生一个以某个中心频率为主的宽带成分。频域带能量比直接计算FFT峰值更稳定,因为窗口存在非平稳波动,单根谱线的幅值抖动很大,而一个频带内的能量面积受谱线位置偏移的影响较小。

from scipy import signal def band_energy_features(x, fs=fs, bands=[(500,1000), (1000,2000), (2000,4000), (4000,8000)]): f, psd = signal.welch(x, fs=fs, nperseg=1024, average='median') feats = {} for lo, hi in bands: mask = (f >= lo) & (f <= hi) feats[f'band_{lo}_{hi}'] = np.trapz(psd[mask], f[mask]) feats['center_freq'] = np.sum(f * psd) / np.sum(psd) feats['freq_std'] = np.sqrt(np.sum((f - feats['center_freq'])**2 * psd) / np.sum(psd)) return feats

welch函数自动完成分段平均,nperseg=1024代表每段1024个采样点。average='median'比默认的mean更抗尖峰干扰,这一点在处理含有周期性冲击的信号时差异明显。用np.trapz做梯形积分,得到的是频带内能量的绝对面积,比取峰值更平滑,两个退化窗口之间即使峰值位置发生偏移,带能量仍然有稳定的上升趋势。

关于频段设置,不要照搬这里的数值。工程上的做法是先画全频谱,找到共振峰集中区域,然后把频段边界设在共振峰附近;常见的高频共振带集中在2000到8000Hz,但具体中心频率取决于传感器安装方式和轴承尺寸。

3.3 时频图检查:退化阶段的“表演”用STFT一目了然

时域和频域特征是对整段信号的压缩,丢失了事件发生的时间信息。为了在退化过程中定位异常时刻,最有效的手段是画时频图。短时傅里叶变换在时间轴上滑动窗口,输出时间-频率-幅度三维分布,共振频带何时被激活、冲击能量是否间歇出现,都能清楚地展现在图上。

import matplotlib.pyplot as plt f_stft, t_stft, Zxx = signal.stft(x[:fs*5], fs=fs, nperseg=256, noverlap=128) plt.figure(figsize=(10, 4)) plt.pcolormesh(t_stft, f_stft, np.abs(Zxx), shading='auto') plt.title('退化时段短时傅里叶变换') plt.colorbar(label='幅值')

nperseg=256在25600Hz采样率下对应10ms窗口,频率分辨率大约100Hz。这个分辨率观察宽带共振足够,如果要区分间隔很近的窄带谱峰,则需要增大窗口长度。shading='auto'用来消除网格线的视觉干扰,让色块连续。

对比不同时期的时频图,健康段通常只在转频及其谐波附近出现细线;退化中后期会出现明显的宽带色块,覆盖数kHz范围,并且随冲击发生呈间歇性。若某个时段出现水平细带但幅值突然增高,往往不是故障,而是传感器附近传入的外部振动干扰,应当回到原始信号去定位,而不是把它当作退化证据。

4. 从统计阈值到LSTM:三种建模方法交叉验证数据集性能特征

特征构建完成后,性能特征需要落到模型层面。只用一种模型评估容易产生偶然性,常见的做法是用三个方法做交叉验证——统计阈值法提供一个完全不依赖训练的可解释基线,随机森林验证非线性特征组合是否有信息量,LSTM补充时间上下文建模能力。三层递进关系,能清楚看到性能特征的边界。

4.1 统计阈值基线:用RMS均值加三倍标准差设置报警门限

在建模之前,先用最简单的统计阈值给每个时间窗打标签。这个方法不需要训练样本,只需要用早期健康段的数据估算正常范围。它虽然粗糙,但在评估复杂模型时是必需的参照系——如果后面随机森林、LSTM的效果压不过这个阈值,问题大概率出在特征而不是模型。

rms_values = np.array([time_domain_features(w)['rms'] for w in windows]) threshold = np.mean(rms_values[:20]) + 3 * np.std(rms_values[:20]) alarm_flags = (rms_values > threshold).astype(int) first_alarm = np.argmax(alarm_flags) # 第一次报警的窗口索引

RMS阈值的取值逻辑是取前20个窗口作为健康基线,均值加三倍标准差作为报警门限。这个参数来自经典质量控制理论,在很多轴承退化数据集上不需要修改就能得到合理结果。第一次报警的窗口索引就是预测的“退化开始时间”,可以通过它与真实起始周期的差值,量化这个基线的滞后量。

阈值法的缺点是响应滞后。RMS上升是渐进过程,真正出现微观裂纹时,RMS可能只比基线高10%左右,阈值法要等到RMS明显超出三倍标准差才能报警,这个时间一般比频带能量方法的报警点晚几个周期。所以阈值法适合作为评估下限,而非冠军候选。

4.2 随机森林回归:让非线性特征组合自动生效

随机森林的优势在于能自动捕捉特征之间的非线性交互,且不需要特征独立同分布。把第3章产出的所有时域和频域特征拼成矩阵,标签用归一化的剩余寿命值(第一个窗口为1,最后一个窗口为0),就能训练一个RUL回归模型。

from sklearn.ensemble import RandomForestRegressor from sklearn.preprocessing import StandardScaler feature_names = list(feature_matrix.columns) X = feature_matrix.values y = rul_values scaler = StandardScaler() X_scaled = scaler.fit_transform(X) rf_model = RandomForestRegressor(n_estimators=300, max_depth=12, min_samples_leaf=2, random_state=42) rf_model.fit(X_scaled, y) mae_rf = mean_absolute_error(y, rf_model.predict(X_scaled))

n_estimators=300、max_depth=12、min_samples_leaf=2这三组参数在滚动轴承RUL任务上是一个比较稳定的起点。树的深度限制很重要,因为样本数不多,不限制深度会完全记住训练集,测试集表现很差;min_samples_leaf=2让叶子节点至少有2个样本,减少单点噪声的影响。

随机森林在这个数据集上一般能拿到比阈值法更低的MAE,但有一个需要注意的现象:它对早期窗口的预测误差小而后期窗口误差大。因为在退化后期RMS曲线的斜率越来越陡,树模型只能记住训练时的变化模式,一旦测试序列的退化速度与训练不同,就会发生系统性偏差。这个残差规律提示我们,后续不能只看MAE,要看残差曲线。

4.3 LSTM时间序列模型:用相邻窗口序列刻入退化轨迹

随机森林的缺憾是无法感知时间上的形态。LSTM通过循环结构把前面20个窗口的信息传递给当前时刻,理论上能学到“这个退化过程是在加速还是减速”,从而提供更稳健的中期预测。输入形状为(样本数,时间步长,特征数),时间步长是本次引入跟随机森林对比的关键超参数。

from tensorflow.keras.models import Sequential from tensorflow.keras.layers import LSTM, Dense def build_sequences(X, y, seq_len=20): Xs, ys = [], [] for i in range(len(X) - seq_len + 1): Xs.append(X[i:i+seq_len]) ys.append(y[i+seq_len-1]) return np.array(Xs), np.array(ys) X_seq, y_seq = build_sequences(X_scaled, y, seq_len=20) model_lstm = Sequential([ LSTM(32, activation='relu', return_sequences=True, input_shape=(20, X.shape[1])), LSTM(16, activation='relu'), Dense(8, activation='relu'), Dense(1, activation='linear') ]) model_lstm.compile(optimizer='adam', loss='mse') model_lstm.fit(X_seq, y_seq, epochs=30, batch_size=64, validation_split=0.2, verbose=1)

LSTM训练时,validation_split=0.2从尾部按顺序取最后20%样本做验证,不是随机打乱。这样能模拟真实的测试场景,因为预测时只能利用当前时刻以前的数据。seq_len=20在窗口长度为1秒时相当于20秒的观测历史,一般够用;历史过长反而引入早期无信息段,造成信息稀释。

LSTM的坑在于结果不稳定。首次训练往往有较好的验证loss,第二次训练因为初始化不同,loss可能就明显变差。这是小样本时序模型的正常现象,不一定是代码问题。处理的方法是固定随机种子,同时用3次独立训练的平均结果评估,而不是取单次结果。

4.4 多方法对比结果:复杂模型并不完全占优

将三种方法在2016年PHM数据集上的关键表现汇总成对比,便于直接看出性能特征的边界条件:

方法是否监督常用特征主要优点典型局限
统计阈值否RMS完全可解释报警滞后
随机森林是全部特征非线性交互后期残差有偏
LSTM是特征+时间序列捕捉趋势末期突变难匹配

实际运行中常见的现象是:随机森林的整体MAE优于LSTM,但在“退化中期”的20个窗口内,LSTM的MAE更低;而到最后10个窗口,LSTM的误差又会反弹。原因是退化末期信号包含大量冲击瞬变,LSTM学到的平滑趋势反而无法匹配这种突变。

这说明“模型更复杂”在这个数据集上并不必然带来更好的结果。性能特征分析的目的不是选一个冠军模型,而是搞清楚每种方法在哪些退化阶段有效、误差从哪里来。这个结论会直接引导下一章的排查方向。

5. 避坑笔记:2016年PHM数据集性能特征分析中的翻车现场

这一章按“现象、原因、解决”的顺序逐一记录,覆盖数据处理和建模中最常见的五个隐患。这些问题在初版分析里几乎都会遇到,提前了解可以省下大量连续调参的时间。

5.1 训练集和验证集随机打乱,把时序泄漏喂给了模型

现象:随机森林在验证集上MAE非常低;切换到测试集后误差暴增。LSTM的validation曲线同样漂亮,但换数据后失效。

原因:轴承退化是连续演变过程,相邻时间窗的特征高度相似。随机打乱之后,同一段退化前后期的窗口分别进入训练集和验证集,模型相当于从“未来数据”里学会了预知当前状态,评价分自然失真。

解决:按周期编号分组划分数据,用GroupShuffleSplit或手动分组,确保训练集出现的周期不会出现在验证集;同时验证集必须是训练集之后的时间段,避免窗口重叠造成泄漏。

5.2 未清理传感器脉冲,峰度特征被单个尖峰拉爆

现象:峰度和峰值因子曲线在健康段出现孤立的高峰,数值甚至超过退化段;建模后这两个特征的重要性排名虚高,测试效果却差。

原因:现场信号混有电磁脉冲干扰,脉冲幅值可能是正常振动的5倍以上。峭度和峰值因子对极端样本极其敏感,一个脉冲就能让窗口特征完全偏离真实状态。

解决:在预处理阶段用滑动中值滤波把尖峰削掉,或对峭度、峰值因子做分位数截尾,把超过99.7%分位的数值替换为边界值。同时要注意不能只清理特征而让原始信号保留尖峰,因为后续还要做STFT诊断。

5.3 每个时间窗独立做标准化,把RMS的增长趋势抹掉了

现象:特征输入LSTM前做过标准化,训练时loss下降正常,但预测输出几乎变成常数;随机森林的预测也不如阈值法。

原因:如果对每个时间窗单独做z-score标准化,RMS从0.1升到1.0的过程会被压缩成均值为0、方差为1的序列,退化最重要的趋势信息就此消失。RUL预测的核心是振幅的增长曲线,单窗归一化等于主动删除了目标信息。

解决:全局标准化只允许用一份scaler在训练集上拟合统计量,再用这份scaler转换验证集和测试集。测试集不能重新计算mean和std,这一点常被忽略。

5.4 LSTM时间步长照搬默认值,训练变慢且效果变差

现象:时间步长设到50甚至100时,训练时间成倍增加,loss下降变慢,最后的MAE比seq_len=20更差;预测曲线整体滞后。

原因:时间步长过长时,LSTM需要记忆的序列跨度太大,前期健康窗口的特征占比太高,把后期退化信息稀释了;同时退化较短的测试序列上,过长的时间步长会截掉尾部样本。

解决:在10、20、30三个值之间做网格搜索。在50%重叠切窗策略下,20通常是记忆效果与计算效率的较好折中;序列样本数不足时,缩到10更稳妥。

5.5 健康段窗口加入回归训练,让模型输出“负寿命”

现象:训练数据中混入未知阶段的健康窗口,模型输出负的RUL;回归没有约束,越界了。

原因:健康段振动幅值接近零,回归模型为了让退化末端的RMSE最小化,会把健康段映射为负值,这是一种典型的边界外推过拟合。

解决:对输出层加ReLU或自定义激活把预测限制在[0,1]区间;同时把纯健康窗口单独挑出来作为辅助分类任务,在回归损失之外加一项健康/非健康分类损失,防止模型把健康段强行映射到0以下。

6. 进阶:用特征重要度和残差曲线给多方法分析收尾

走到这里,数据已经完成了预处理、特征提取、三路建模和坑排查,怎样判断分析是否真正做完?我的习惯是把最后一步放在特征重要度排序和残差曲线追溯上。这个收尾动作保证了整套多方法分析可信、可解释、可复用。

先看随机森林的feature_importances_,把排名前5的特征列出来,检查它们是否来自不同的特征族。常见情况是排名靠前的有RMS、某个高频频带能量、峭度和重心频率,这说明时域和频域都在贡献信息,性能特征不是靠某一个指标独自撑起来的。

importance = pd.Series(rf_model.feature_importances_, index=feature_names) top5 = importance.sort_values(ascending=False).head(5) print(top5)

接下来画残差曲线,即预测RUL与真实RUL的差值随时间的变化。若残差在早期为正、中后期为负,说明模型系统性高估早期寿命、低估后期寿命,根源往往是训练集中的退化速度分布不匀;若残差呈现锯齿状波动,则要回头检查异常脉冲是否清理干净。对残差偏差很大的样本,我会逐一回到STFT图上确认原因,这个习惯帮我在不止一个项目里发现过传感器松动和外部干扰,而不是把问题简单归咎于模型。

多方法分析最终要不要在某个方向上继续深入,主要看残差追溯是否有新解释;如果特征重要度稳定、残差没有系统性偏置,这套流程就可以作为后续算法选型的可信边界。希望这篇经验能帮到你,让你在动手拆这个数据集时少走弯路。

本文还有配套的精品资源,点击获取

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询