用Python和Matlab搞定东南大学齿轮箱数据集：从数据读取到故障分类实战-编程实验室

东南大学齿轮箱数据集实战：Python与Matlab双视角解析

引言：机械故障诊断的数据科学之旅

机械故障诊断领域正经历着从传统经验判断向数据驱动决策的转变。东南大学齿轮箱数据集作为国内知名的开放研究资源，包含了轴承和齿轮两大类共10种故障状态，为算法验证提供了宝贵素材。但面对8通道振动信号、多工况CSV文件时，许多研究者常陷入"数据在手，无从下手"的困境。

本文将采用工程化思维，以Python生态为主力工具链，Matlab为辅助参考，构建从原始振动信号到分类模型的完整分析流程。不同于简单的代码展示，我们会重点探讨：

多源异构振动信号的标准化处理方法
时频域特征工程的构建逻辑
分类模型在机械信号中的适配技巧
两种编程环境的协同工作模式

无论您是刚开始接触机械信号分析的在校学生，还是需要快速验证算法的工程人员，这套方法论都能帮助您建立系统化的分析框架。

1. 数据环境搭建与预处理

1.1 跨平台数据加载方案

东南大学数据集采用CSV格式存储，但包含特殊配置信息头。我们推荐以下跨语言处理方案：

Python方案（Pandas优先）

import pandas as pd def load_gear_data(filepath): # 跳过前两行配置信息 raw_data = pd.read_csv(filepath, skiprows=2, header=None) # 8通道信号标准化命名 columns = ['motor_vib', 'planet_x', 'planet_y', 'planet_z', 'motor_torque', 'reducer_x', 'reducer_y', 'reducer_z'] return raw_data.set_axis(columns, axis=1)

Matlab对照方案

function data = loadGearData(filename) opts = detectImportOptions(filename); opts.DataLines = [3 Inf]; data = readtable(filename, opts); data.Properties.VariableNames = {'motor_vib', 'planet_x',... 'planet_y', 'planet_z', 'motor_torque',... 'reducer_x', 'reducer_y', 'reducer_z'}; end

注意：原始数据采样频率为5120Hz，建议在加载阶段就记录该元数据，后续分析会频繁用到这个关键参数。

1.2 多文件批量处理框架

面对数十个CSV文件，需要建立系统化的批处理框架：

from pathlib import Path def build_dataset(root_dir): dataset = {} for csv_file in Path(root_dir).glob('*.csv'): # 从文件名解析故障类型和工况 parts = csv_file.stem.split('_') fault_type = parts[0] condition = parts[1] if len(parts)>1 else 'default' # 统一数据存储结构 if fault_type not in dataset: dataset[fault_type] = {} dataset[fault_type][condition] = load_gear_data(csv_file) return dataset

关键处理技巧：

使用文件命名规范自动提取元信息
构建嵌套字典保存原始信号
保留完整的信号通道关系

2. 振动信号可视化与探索

2.1 多通道信号对比分析

典型的多通道信号可视化方案：

import matplotlib.pyplot as plt def plot_multi_channel(data, start=0, duration=0.1): samples = int(duration * 5120) # 计算采样点数 subset = data.iloc[start:start+samples] fig, axes = plt.subplots(4, 2, figsize=(12, 8)) for idx, col in enumerate(data.columns): ax = axes[idx//2, idx%2] ax.plot(subset[col]) ax.set_title(f'Channel {col}') plt.tight_layout() return fig

信号观察要点：

电机振动（motor_vib）通常包含丰富的高频成分
行星齿轮XYZ三向振动呈现明显相关性
扭矩信号往往具有不同的量纲需要单独处理

2.2 时频域联合分析

结合快速傅里叶变换（FFT）和短时傅里叶变换（STFT）的综合观察：

from scipy import signal def time_freq_analysis(channel_data, fs=5120): # 时域特征 time_features = { 'peak': channel_data.max(), 'rms': np.sqrt(np.mean(channel_data**2)), 'kurtosis': channel_data.kurtosis() } # 频域特征 freqs, psd = signal.welch(channel_data, fs) dominant_freq = freqs[np.argmax(psd)] return {'time': time_features, 'main_freq': dominant_freq}

提示：齿轮故障常在啮合频率及其谐波处表现出异常，轴承故障则多在外圈/内圈特征频率附近出现峰值。

3. 特征工程构建策略

3.1 时域特征组合

推荐提取的时域特征集：

特征类别	具体指标	物理意义
幅值特征	峰值、峰峰值、RMS	振动能量水平
波形指标	偏度、峭度、波形因子	信号分布特性
脉冲指标	脉冲因子、裕度因子	冲击成分强度
相关性特征	通道间互相关系数	系统耦合程度

Python实现示例：

def extract_time_features(data_frame): features = {} for col in data_frame.columns: x = data_frame[col].values features.update({ f'{col}_peak': np.max(x), f'{col}_kurtosis': stats.kurtosis(x), f'{col}_crest': np.max(x)/np.sqrt(np.mean(x**2)) }) return features

3.2 频域特征提取

关键频域特征包括：

频谱重心：反映能量分布位置

def spectral_centroid(psd, freqs): return np.sum(freqs*psd) / np.sum(psd)

频带能量比：划分5个典型频段计算能量占比

% Matlab实现 function ratio = band_energy_ratio(psd, f_edges) total = sum(psd); ratio = []; for i = 1:length(f_edges)-1 mask = (f >= f_edges(i)) & (f < f_edges(i+1)); ratio(end+1) = sum(psd(mask)) / total; end end

谐波成分检测：针对齿轮特征频率的谐波分析

4. 故障分类模型构建

4.1 特征选择与数据集构建

典型特征矩阵构建流程：

from sklearn.pipeline import make_pipeline from sklearn.preprocessing import StandardScaler # 构建特征矩阵 def build_feature_matrix(dataset): X, y = [], [] for fault_type, conditions in dataset.items(): for condition, data in conditions.items(): # 合并时域频域特征 features = extract_time_features(data) features.update(extract_freq_features(data)) X.append(list(features.values())) y.append(fault_type) return np.array(X), np.array(y) # 标准化管道 preprocessor = make_pipeline( StandardScaler(), PCA(n_components=0.95) )

4.2 模型训练与评估

推荐模型对比表：

模型类型	优点	注意事项
SVM	小样本效果好，泛化能力强	需要仔细调参
随机森林	自动特征选择，抗噪声	可能过拟合
1D CNN	自动特征学习，端到端	需要大量数据
XGBoost	处理不平衡数据优秀	超参数较多

Python实现示例：

from sklearn.ensemble import RandomForestClassifier from sklearn.model_selection import cross_val_score # 随机森林实现 model = RandomForestClassifier( n_estimators=200, max_depth=10, class_weight='balanced' ) cv_scores = cross_val_score( model, X_preprocessed, y, cv=5, scoring='accuracy' ) print(f"Mean CV Accuracy: {cv_scores.mean():.2%}")

4.3 模型解释与故障分析

SHAP值分析示例：

import shap # 训练最终模型 model.fit(X_train, y_train) # 解释模型预测 explainer = shap.TreeExplainer(model) shap_values = explainer.shap_values(X_test) # 可视化特定样本分析 shap.force_plot( explainer.expected_value[0], shap_values[0][0,:], X_test[0,:], feature_names=feature_names )

关键分析维度：

各通道振动信号的贡献度排序
时域与频域特征的相对重要性
不同故障类型的判别特征差异

5. 工程实践中的进阶技巧

5.1 数据增强策略

机械信号特有的数据增强方法：

时域窗口滑动：通过重叠采样增加样本量

def sliding_window(data, window_size, step): n_samples = len(data) return [data[i:i+window_size] for i in range(0, n_samples-window_size, step)]

噪声注入：添加符合实际工况的高斯噪声

def add_industrial_noise(signal, snr_db=20): rms_signal = np.sqrt(np.mean(signal**2)) noise = np.random.normal(0, rms_signal/(10**(snr_db/20)), len(signal)) return signal + noise

5.2 跨语言协作模式

Python与Matlab混合编程方案：

数据交换格式：
- 使用HDF5(.h5)实现大型数组交换
- 通过MAT文件接口scipy.io.savemat/loadmat

实时通信：

% Matlab作为计算引擎 eng = py.matlab.engine.start_matlab() result = eng.sqrt(4.0)

性能分工：
- Matlab处理信号处理专用算法
- Python负责机器学习流程

5.3 实际部署考量

边缘设备部署优化策略：

特征压缩：使用PCA保留95%方差

from sklearn.decomposition import PCA pca = PCA(n_components=0.95) X_compressed = pca.fit_transform(X)

模型轻量化：
- 决策树剪枝
- 神经网络量化
实时性保障：
- 固定长度滑动窗口
- 增量特征计算

常见问题解决方案

Q1 数据量不足怎么办？

采用迁移学习，复用ImageNet预训练模型
使用生成对抗网络(GAN)进行数据增强

Q2 不同工况数据分布不一致？

添加工况标签作为额外特征
采用领域自适应(Domain Adaptation)技术

Q3 如何选择采样长度？

齿轮分析至少包含5个完整啮合周期
轴承故障建议0.5-1秒时长

Q4 特征重要性如何评估？

排列重要性(Permutation Importance)
基于SHAP值的全局解释

延伸方向与资源推荐

进阶学习路径：

时频分析：小波变换、Hilbert-Huang变换
深度学习：1D CNN、LSTM、Transformer
异常检测：One-Class SVM、Autoencoder

参考工具库：

Python: PyBearings, SciPy Signal
Matlab: Predictive Maintenance Toolbox
开源项目: SKF @ptitude Analyst

公开数据集扩展：

NASA轴承数据集
凯斯西储大学数据集
Paderborn大学轴承数据集

用Python和Matlab搞定东南大学齿轮箱数据集：从数据读取到故障分类实战