格拉姆角场与轴承故障诊断:从时序信号到图像识别的数据预处理实战

发布时间:2026/8/2 10:53:43
格拉姆角场与轴承故障诊断:从时序信号到图像识别的数据预处理实战 1. 项目概述从代码到数据理解故障诊断的基石拿到一份名为“格拉姆角场东南大学轴承故障诊断代码解读——数据集解读”的代码很多朋友可能会直接一头扎进模型构建和训练的部分急切地想看到诊断准确率。但根据我多年的工业数据分析经验这恰恰是新手最容易踩坑的地方。一个故障诊断项目的成败在模型跑起来之前就已经被数据决定了七八成。这份代码的标题将“数据集解读”放在后半部分而我认为它应该是我们打开任何类似项目时第一个、也是最需要花时间吃透的环节。格拉姆角场Gramian Angular Field, GAF是一种将一维时间序列转换为二维图像矩阵的编码方法它通过保留时间序列的绝对时序关系和数值信息为后续使用成熟的图像分类模型如CNN处理振动信号铺平了道路。而“东南大学轴承数据集”则是在机械故障诊断领域一个非常经典且公开的基准数据集。这个项目本质上就是利用GAF技术将东南大学轴承的振动时序信号“翻译”成图像然后用计算机视觉的方法来识别轴承的健康状态和故障类型。所以在动手修改任何一行模型代码之前我们必须彻底搞清楚我们喂给模型的是什么“粮食”这些“粮食”是怎么从原始的振动信号加工而来的数据里有没有“杂质”或“偏见”只有把数据集这第一道关把好了后续的模型调优、结果分析才有意义。否则很可能出现模型在训练集上表现完美一到实际场景就“翻车”的情况。接下来我就带大家深入这个项目的“后厨”看看数据是如何被准备和处理的。2. 核心数据集东南大学轴承数据深度解析在解读任何相关代码前我们必须先独立于代码理解数据本身的来源、结构和物理意义。这是避免被代码实现带偏、形成自己判断力的关键。2.1 数据来源与采集背景东南大学SEU的轴承数据集是在实验室环境下通过转子实验台采集的。实验台通常包含电机、转轴、支撑轴承、加载装置等部分通过在健康轴承和预设故障的轴承上安装振动加速度传感器来收集数据。故障类型通常包括内圈故障、外圈故障、滚动体故障并且每种故障会有不同尺寸如0.007英寸0.014英寸0.021英寸的模拟损伤。这个数据集之所以经典是因为它工况相对可控负载、转速通常是固定的这减少了变量便于初学者聚焦于故障特征本身。故障模式典型涵盖了旋转机械中最常见的几种轴承故障类型。数据格式规整通常以.matMATLAB数据文件或文本文件形式提供采样频率、数据长度标注清晰。注意实验室数据与现场数据存在“鸿沟”。实验室数据信噪比高故障特征明显而现场数据受背景噪声、工况波动、多源耦合振动影响巨大。因此在实验室数据集上表现优异的模型直接部署到工厂可能需要大幅调整。2.2 数据结构与文件组织通常下载到的SEU数据集文件夹结构如下SEU_Bearing_Dataset/ ├── 正常/ │ ├── normal_1.mat │ ├── normal_2.mat │ └── ... ├── 内圈故障/ │ ├── IR007_1.mat │ ├── IR007_2.mat │ ├── IR014_1.mat │ └── ... ├── 外圈故障/ │ └── ... └── 滚动体故障/ └── ...每个.mat文件里通常存储着一个或多个通道的振动加速度时序数据。关键参数需要从数据说明或代码中提取采样频率Fs如12kHz 24kHz。这决定了信号能捕获的最高频率奈奎斯特频率为Fs/2。数据长度每个文件可能包含几十万甚至上百万个数据点代表一段连续采样的振动信号。转速与负载这些信息对于理解故障特征频率至关重要但数据集有时未必直接附带需要从实验描述中查找。2.3 故障的物理特征与在信号中的体现轴承的故障会在振动信号中产生周期性冲击其频率由故障类型和几何参数决定称为故障特征频率。内圈故障频率BPFI与转频相关频率较高且由于载荷方向变化振幅会有调制现象。外圈故障频率BPFO频率相对固定振幅稳定。滚动体故障频率BSF频率通常低于内外圈故障。保持架故障频率FTF频率最低。在时域波形上健康信号相对平稳而故障信号会出现明显的、周期性的冲击脉冲。在频域通过傅里叶变换我们可以在相应的故障特征频率及其倍频处看到突出的谱线。理解这些物理背景至关重要因为GAF将一维信号转为图像后这些时域和频域的特征会以某种纹理、形状或亮度的模式体现在图像中。我们的目标就是让CNN模型学会识别这些与特定故障对应的图像模式。3. 格拉姆角场原理与代码实现拆解理解了“原材料”原始振动数据后我们来看“烹饪方法”——格拉姆角场。代码中实现GAF的部分是核心我们需要明白每一步的数学意义和工程考量。3.1 GAF转换的核心步骤GAF主要分为两种格拉姆角和场GASF和格拉姆角差场GADF。项目代码中通常使用其中一种或两者结合。其转换流程可分解为以下四步我结合代码中可能出现的函数进行解释第一步数据归一化将原始振动信号X [x1, x2, ..., xn] 缩放到区间[-1, 1]或[0, 1]。# 常见代码片段示例 from sklearn.preprocessing import MinMaxScaler scaler MinMaxScaler(feature_range(-1, 1)) X_normalized scaler.fit_transform(X.reshape(-1, 1)).flatten()为什么必须归一化因为GAF基于角度计算而归一化到[-1,1]区间后数据点可以映射到单位圆上的余弦值。如果数据量纲不统一例如不同通道、不同实验的数据绝对值大的信号会主导角度计算导致信息失真。第二步将数值转换为角度通过反余弦函数将归一化后的值映射为角度弧度制。import numpy as np phi np.arccos(X_normalized) # 此时 phi 在 [0, pi] 区间内这一步是GAF的精华。每个数据点不再是一个孤立的振幅值而是单位圆上的一个点由其与横轴的夹角φ来表征。时序信息被巧妙地编码进了这个角度序列中。第三步计算格拉姆矩阵核心这是生成二维图像的关键。格拉姆矩阵的元素由每两个点之间的三角和差关系构成。格拉姆角和场GASF计算角度之和的余弦。它更侧重于捕捉信号之间的“和”关系反映的是整体相关性。# 伪代码逻辑 GASF np.cos(phi_i phi_j) # 其中 i, j 遍历所有数据点 # 实际代码利用向量化操作避免低效循环格拉姆角差场GADF计算角度之差的正弦。它更侧重于捕捉信号之间的“差”或相对变化关系。GADF np.sin(phi_i - phi_j)生成的GASF或GADF矩阵是一个n x n的对称矩阵对于GASF或反对称矩阵对于GADF其中n是输入时序片段的长度。这个矩阵就是我们要的“图像”。第四步图像化与裁剪生成的n x n矩阵可能很大例如1000x1000。直接作为CNN输入可能计算量过大。因此代码中通常会有以下操作降采样在计算GAF前先对长时序信号进行切片或降采样使n控制在一个合理大小如224适配ImageNet预训练模型。图像缩放生成GAF矩阵后使用cv2.resize或PIL.Image.resize将其缩放到统一尺寸如224x224。伪彩色映射GAF矩阵是单通道的每个像素一个值。为了适配通常输入为3通道的CNN如ResNet需要将其转换为“伪彩色”图像。常见方法是使用cv2.applyColorMap如COLORMAP_JET或简单地将同一矩阵复制到三个通道。3.2 代码中的关键参数与选择在解读data_preprocessing.py或generate_gaf_images.py这类文件时要重点关注以下参数切片长度segment_length如1024个点。这决定了每张“图像”代表多长时间的振动信号。太短可能包含不了一个完整的故障冲击周期太长则图像分辨率过高且可能混合多种状态。经验上这个长度应能覆盖至少2-3个故障特征周期。例如转速为1800 rpm30 HzBPFI约为5.4倍转频即162 Hz周期约6.2毫秒。若采样频率为12 kHz则一个周期约74个点。因此切片长度取512或1024点是合理的。重叠率overlap_rate如0.550%。滑动窗口切分信号时相邻切片之间的重叠比例。提高重叠率可以增加生成图像的数量缓解数据量不足的问题但也会引入更强的样本相关性可能影响模型泛化能力评估的准确性。GAF类型选择只用GASF还是GADF或者将两者合并为双通道图像不同的故障特征在不同场中的表现可能不同。我个人的经验是对于轴承的周期性冲击故障GASF往往能更好地保留冲击的时序相关性效果更稳定。可以尝试融合但会增加模型输入通道和计算量。缩放方法与插值算法将GAF矩阵缩放到目标尺寸如224x224时cv2.INTER_LINEAR双线性插值是常用选择。应避免使用INTER_NEAREST最近邻插值因为它可能在图像中引入块状伪影破坏连续的特征模式。4. 数据预处理流程全链路实操现在我们把数据集和GAF原理串联起来看一个完整的、可复现的数据预处理流水线应该如何构建。这是项目代码的核心骨架。4.1 步骤一原始数据加载与探查在写任何处理代码之前先用Jupyter Notebook或脚本进行数据探查。import scipy.io as sio import numpy as np import matplotlib.pyplot as plt # 1. 加载一个.mat文件示例 data_dict sio.loadmat(‘path/to/IR007_1.mat‘) # 打印所有键查看数据结构 print(data_dict.keys()) # 通常振动数据在 ‘data‘, ‘vibration‘, ‘X‘ 等键下 vibration_signal data_dict[‘X‘].flatten() # 假设键名为‘X‘并转换为一维数组 # 2. 绘制时域波形 plt.figure(figsize(12, 4)) plt.plot(vibration_signal[:5000]) # 只看前5000个点 plt.title(‘Raw Vibration Signal (Time Domain)‘) plt.xlabel(‘Sample Points‘) plt.ylabel(‘Amplitude‘) plt.grid(True) plt.show() # 3. 计算并绘制频谱快速傅里叶变换 from scipy.fft import fft, fftfreq Fs 12000 # 假设采样频率为12kHz N len(vibration_signal) yf fft(vibration_signal) xf fftfreq(N, 1/Fs)[:N//2] # 取正频率部分 plt.figure(figsize(12,4)) plt.plot(xf, 2.0/N * np.abs(yf[0:N//2])) plt.title(‘Frequency Spectrum‘) plt.xlabel(‘Frequency (Hz)‘) plt.ylabel(‘Magnitude‘) plt.grid(True) plt.xlim([0, Fs/2]) # 显示到奈奎斯特频率 plt.show()这个探查步骤能帮你确认信号质量观察是否有明显的故障冲击并验证采样频率。4.2 步骤二数据切片与标签生成这是为后续GAF转换准备输入片段和对应标签。def segment_signal(signal, label, segment_length, overlap_ratio): 将一维信号切分成固定长度的片段并分配标签。 参数: signal: 一维振动信号数组。 label: 该信号对应的整数型标签如0:正常1:内圈故障...。 segment_length: 每个片段的长度。 overlap_ratio: 重叠率0-1之间。 返回: segments: 片段列表形状为 (num_segments, segment_length)。 labels: 标签列表形状为 (num_segments,)。 segments [] labels [] step int(segment_length * (1 - overlap_ratio)) if step 0: step 1 # 避免死循环 num_segments (len(signal) - segment_length) // step 1 for i in range(num_segments): start i * step end start segment_length segment signal[start:end] # 可选这里可以添加片段能量检查过滤掉能量过低的无效片段 segments.append(segment) labels.append(label) return np.array(segments), np.array(labels) # 遍历所有数据文件夹收集所有片段和标签 all_segments [] all_labels [] class_folders {‘normal‘: 0, ‘IR‘: 1, ‘OR‘: 2, ‘Ball‘: 3} # 示例映射 for class_name, label_id in class_folders.items(): folder_path os.path.join(‘dataset‘, class_name) for file_name in os.listdir(folder_path): if file_name.endswith(‘.mat‘): file_path os.path.join(folder_path, file_name) signal load_signal_from_mat(file_path) # 自定义加载函数 segments, labels segment_signal(signal, label_id, segment_length1024, overlap_ratio0.5) all_segments.extend(segments) all_labels.extend(labels) # 转换为NumPy数组 all_segments np.array(all_segments) all_labels np.array(all_labels) print(f“总片段数 {all_segments.shape[0]} 片段长度 {all_segments.shape[1]} 类别数 {len(np.unique(all_labels))}“)4.3 步骤三GAF图像批量生成将上一步得到的所有信号片段批量转换为GAF图像。from pyts.image import GramianAngularField # 可以使用pyts库也可以自己实现 import cv2 # 方法1使用pyts库推荐稳定且高效 def generate_gaf_images_pyts(segments, image_size224, method‘summation‘): 使用pyts库生成GAF图像。 method: ‘summation‘ (GASF) 或 ‘difference‘ (GADF) gasf GramianAngularField(image_sizeimage_size, methodmethod) # pyts要求输入形状为 (n_samples, n_timestamps) images_gasf gasf.fit_transform(segments) # 输出形状 (n_samples, image_size, image_size) # 将值域从[-1,1]或[0,1]映射到[0, 255]的uint8并应用伪彩色 images_uint8 ((images_gasf 1) * 127.5).astype(np.uint8) # 假设值域为[-1,1] colored_images [] for img in images_uint8: colored cv2.applyColorMap(img, cv2.COLORMAP_JET) colored_images.append(colored) return np.array(colored_images) # 形状 (n_samples, image_size, image_size, 3) # 方法2手动实现更灵活便于理解原理 def gramian_angular_field(series, method‘summation‘): 手动计算单一样本的GAF矩阵 # 归一化 min_val, max_val series.min(), series.max() scaled_series (2 * (series - min_val) / (max_val - min_val)) - 1 # 归一化到[-1,1] scaled_series np.clip(scaled_series, -1, 1) # 防止反余弦计算溢出 # 计算角度 phi np.arccos(scaled_series) # 计算格拉姆矩阵 if method ‘summation‘: # GASF cos(φ_i φ_j) cos_sum np.cos(np.add.outer(phi, phi)) return cos_sum elif method ‘difference‘: # GADF sin(φ_i - φ_j) sin_diff np.sin(np.subtract.outer(phi, -phi)) # 注意符号处理 return sin_diff # 批量处理 image_size 224 gaf_images [] for segment in all_segments[:100]: # 示例先处理100个 gaf_matrix gramian_angular_field(segment) # 缩放 resized_matrix cv2.resize(gaf_matrix, (image_size, image_size), interpolationcv2.INTER_LINEAR) # 伪彩色和归一化到[0,255] normalized ((resized_matrix 1) * 127.5).astype(np.uint8) colored cv2.applyColorMap(normalized, cv2.COLORMAP_VIRIDIS) gaf_images.append(colored) gaf_images np.array(gaf_images)4.4 步骤四数据集划分与保存将生成的图像数据集划分为训练集、验证集和测试集并保存为文件如TFRecord或直接NumPy数组标签方便模型加载。from sklearn.model_selection import train_test_split import pickle # 划分数据集先划分索引避免数据混乱 indices np.arange(len(gaf_images)) X_train_idx, X_temp_idx, y_train_idx, y_temp_idx train_test_split( indices, all_labels[:len(gaf_images)], test_size0.3, random_state42, stratifyall_labels[:len(gaf_images)] ) X_val_idx, X_test_idx, y_val_idx, y_test_idx train_test_split( X_temp_idx, y_temp_idx, test_size0.5, random_state42, stratifyy_temp_idx ) # 根据索引获取数据 X_train, y_train gaf_images[X_train_idx], all_labels[X_train_idx] X_val, y_val gaf_images[X_val_idx], all_labels[X_val_idx] X_test, y_test gaf_images[X_test_idx], all_labels[X_test_idx] print(f“训练集 {X_train.shape} 验证集 {X_val.shape} 测试集 {X_test.shape}“) # 保存数据集 save_dict { ‘X_train‘: X_train, ‘y_train‘: y_train, ‘X_val‘: X_val, ‘y_val‘: y_val, ‘X_test‘: X_test, ‘y_test‘: y_test, ‘label_names‘: {0: ‘正常‘, 1: ‘内圈故障‘, 2: ‘外圈故障‘, 3: ‘滚动体故障‘} } with open(‘seu_bearing_gaf_dataset.pkl‘, ‘wb‘) as f: pickle.dump(save_dict, f) print(“数据集已保存为 ‘seu_bearing_gaf_dataset.pkl‘“)5. 关键注意事项与避坑指南在实际操作中有几个细节如果不注意很容易导致模型效果不佳或结论错误。5.1 数据泄露问题这是时序数据划分中最常见的坑。绝对不能用随机打乱后再划分的方法错误做法将all_segments和all_labels用sklearn.model_selection.train_test_split直接随机划分。因为重叠切片的存在同一个原始样本的不同切片可能被分到了训练集和测试集导致模型通过“记忆邻居”就能做出正确预测严重高估泛化能力。正确做法按“样本源文件”划分。即将所有.mat文件列表打乱然后按比例分配给训练集、验证集和测试集。确保来自同一个原始数据文件的所有切片都只出现在同一个集合中。这样才能模拟现实场景模型用一批机器历史数据训练去预测另一批机器未来的状态。5.2 类别不平衡处理轴承故障数据中正常状态的数据往往远多于各种故障状态的数据因为故障是偶发事件。直接训练会导致模型偏向于预测“正常”类。对策在数据加载或训练时进行处理。过采样对少数类样本进行复制或使用SMOTE需谨慎对于图像数据简单的复制可能导致过拟合。欠采样随机丢弃一部分多数类样本可能损失有用信息。类别权重在损失函数中为少数类赋予更高的权重。这是最常用且有效的方法。在PyTorch或TensorFlow中可以方便地设置class_weight参数。# 以sklearn为例计算类别权重 from sklearn.utils.class_weight import compute_class_weight class_weights compute_class_weight(‘balanced‘, classesnp.unique(y_train), yy_train) # 在训练时将class_weights传递给损失函数5.3 GAF图像的可视化与检查生成GAF图像后一定要抽样可视化检查转换是否合理。import matplotlib.pyplot as plt fig, axes plt.subplots(2, 4, figsize(16, 8)) for i in range(4): # 展示4个类别每个类别2个样本 # 选取某个类别的样本索引 class_idx np.where(y_train i)[0] sample_idx class_idx[0] axes[0, i].imshow(X_train[sample_idx]) axes[0, i].set_title(f‘{label_names[i]} - Sample 1‘) axes[0, i].axis(‘off‘) sample_idx class_idx[1] axes[1, i].imshow(X_train[sample_idx]) axes[1, i].set_title(f‘{label_names[i]} - Sample 2‘) axes[1, i].axis(‘off‘) plt.tight_layout() plt.show()你需要观察不同类别的图像在纹理、颜色分布上是否有肉眼可辨的差异同一类别的不同样本是否具有相似的模式如果看起来都是杂乱无章的噪声可能需要检查GAF转换的参数如归一化是否正确或考虑原始信号是否信噪比过低。5.4 计算资源与效率优化GAF计算是计算密集型操作尤其是当数据量大、切片长度长时。向量化操作务必使用NumPy的广播和向量化函数如np.add.outer,np.subtract.outer来替代Python循环速度可提升成百上千倍。分批处理与存储不要试图一次性将所有数据转为GAF图像并加载进内存。应该设计一个生成器Generator或使用tf.data.Dataset/torch.utils.data.DataLoader的预处理管道在需要时动态生成批次数据或者预处理后存储为.tfrecord或.h5格式便于流式读取。并行处理利用multiprocessing库或joblib对多个数据文件或片段进行并行GAF转换。6. 进阶思考从SEU数据集到工业实际应用基于SEU数据集和GAF的方法在实验室环境下可以取得很高的准确率99%但这离真正的工业应用还有距离。在解读代码、复现结果之后我们更应该思考以下几个问题这也是项目价值能否延伸的关键6.1 模型学到了什么——可解释性分析当你的CNN模型在测试集上达到高精度后不要满足于此。我们需要知道模型是根据图像的哪些部分做出判断的这有助于验证模型是否真的学到了有物理意义的故障特征而不是一些无关的伪影。梯度加权类激活映射Grad-CAM这是最常用的可视化方法。它可以生成一个热力图叠加在原始GAF图像上显示哪些区域对模型决策的贡献最大。分析观察正常样本和故障样本的热力图差异。对于故障样本高亮区域是否集中在图像的某些特定对角线或区块这些区域可能对应着原始振动信号中周期性冲击发生的时间点。如果热力图总是集中在图像边缘或无关区域那模型的决策依据可能是可疑的。6.2 如何处理变工况与噪声——数据增强与域适应SEU数据集的工况是固定的。但现实中设备的转速、负载会变化环境噪声也更大。数据增强在生成GAF图像后或之前可以引入增强策略。时域增强对原始信号添加随机缩放、轻微抖动、添加高斯白噪声。图像域增强对GAF图像使用标准的图像增强如随机旋转小角度、水平/垂直翻转需谨慎GAF图像具有对称性翻转可能改变物理意义、亮度对比度微调。域适应Domain Adaptation如果目标工况与SEU数据集差异巨大可以考虑使用迁移学习或域自适应算法如DANN利用SEU的带标签数据和目标域的无标签数据让模型学会工况不变的特征。6.3 除了GAF还有哪些编码方式GAF不是唯一将时序信号转为图像的方法。了解其他方法有助于你在不同场景下做出选择。马尔可夫变迁场Markov Transition Field, MTF将时序数据离散化后计算状态间转移概率的矩阵。对刻画信号的状态变迁规律有效。递归图Recurrence Plots, RP展示信号在相空间中哪些时刻回到相似状态。对非线性、非平稳信号分析有优势。连续小波变换时频谱图CWT Scalogram这是更经典的方法能同时提供时间和频率信息。对于轴承故障冲击成分会在特定频带产生“条纹”非常直观。如何选择一个实用的建议是从时频谱图CWT开始。因为它最符合振动分析工程师的看图习惯特征物理意义明确。如果追求极致的端到端自动化再尝试GAF、MTF等编码方法并与CWT进行效果对比。6.4 部署考量轻量化与实时性实验室代码往往不考虑推理速度。但工业部署要求模型轻量化、推理快。模型压缩考虑使用MobileNet、EfficientNet等轻量级CNN backbone或者对训练好的模型进行剪枝、量化。推理流水线优化将GAF转换和模型推理集成到C或嵌入式环境中使用ONNX Runtime、TensorRT等推理引擎进行加速。边缘部署对于实时监测可以考虑在设备边缘如工控机、带算力的网关直接运行轻量级模型只将报警结果或特征上传到云端而非原始振动数据。解读“格拉姆角场东南大学轴承故障诊断代码”的最终目的绝不是为了在SEU数据集上刷出一个漂亮的数字。而是通过这个标准的“实验室样板”掌握“数据理解 - 特征工程GAF- 模型构建 - 评估分析”的完整方法论链条。当你吃透了数据集的每一个字节理解了GAF矩阵中每一个像素的来源与含义你才真正拥有了将这套方法迁移到其他设备、其他故障诊断场景的能力。数据是地基特征工程是骨架模型只是血肉。把地基打牢骨架搭正项目才能立得住走得远。