简介:这份资源是2022年第二届天府杯全国大学生数学建模竞赛A题的完整参赛论文,作者为参赛队#A485成员,面向备战数学建模竞赛的高校学生及对仪器故障智能诊断感兴趣的机器学习初学者。论文围绕基于图像识别与多模型建模的故障检测技术展开,系统覆盖信号去噪、时域特征提取、无监督与有监督分类建模以及迁移学习图像分类等核心环节。资源包内含1个PDF文件,大小约1.29MB,完整呈现从问题重述、模型假设到各问求解与评价指标分析的论文全貌。目前已有274人学习下载。读者可从中获取滑动平均、Savitzky-Golay滤波与小波变换的去噪对比思路,10类时域特征的提取方法,K-means、DBSCAN、Spectral聚类与SVM、Xgboost、MLP等八种有监督算法的建模流程,以及AlexNet、ResNet-18迁移学习在故障图像分类中的应用,为同类赛题提供可复用的方案参考与排错思路。
1. 天府杯A题仪器故障智能诊断:从赛题到可复现方案的拆解
仪器故障智能诊断这个方向,在数学建模竞赛里出现的频率越来越高。2022年天府杯数学建模A题把它作为核心命题,要求参赛队在有限时间内完成从数据预处理到诊断模型构建的全流程。这道题真正难的地方不在于模型有多复杂,而在于仪器信号本身信噪比低、故障样本稀缺、多类故障边界模糊——这三点恰好也是工业现场设备诊断的通用痛点。如果你正在做数学建模赛题、或者手头有设备振动/温度/电流信号需要做异常识别,这套思路可以直接迁移。我下面按“信号怎么读→特征怎么提→模型怎么选→结果怎么验”的顺序,把当时踩过的坑和最终跑通的路径讲清楚,代码基于Python,依赖numpy、scipy、sklearn,不需要GPU也能复现。
2. 仪器故障信号的数据规范化处理:从原始采样到可训练矩阵
2.1 为什么规范化处理是绕不过去的第一步
仪器故障诊断的数据来源通常是传感器采集的时序信号,采样率从几千赫兹到几十千赫兹不等。原始数据直接丢进模型,十有八九会翻车。原因有三个:不同传感器通道量纲差异大(振动加速度是m/s²,温度是℃,电流是A),不做归一化的话梯度下降会被大量纲特征主导;采样频率不一致导致时间轴对不齐,无法直接拼接;信号中混有工频干扰和随机噪声,不滤掉的话特征提取会提取出一堆噪声的统计量。
常见做法是分三步走:重采样对齐时间轴、带通滤波去工频和高频噪声、按通道做Z-score标准化。这三步做完,数据才算“干净”。很多参赛队在这一步偷懒,后面模型效果上不去,回头查半天以为是模型问题,其实是数据没洗干净。
2.2 用Python做重采样、滤波与标准化
下面这段代码处理的是典型的多通道仪器信号,假设原始数据是CSV格式,每列一个通道,最后一列是故障标签。
import numpy as np import pandas as pd from scipy import signal from sklearn.preprocessing import StandardScaler # 读取原始数据 df = pd.read_csv('instrument_data.csv') X_raw = df.iloc[:, :-1].values # 信号通道 y = df.iloc[:, -1].values # 故障标签 # 第一步:重采样到统一频率(假设目标为1000Hz) fs_target = 1000 fs_original = 5000 # 原始采样率,需根据实际数据修改 num_samples = int(X_raw.shape[0] * fs_target / fs_original) X_resampled = signal.resample(X_raw, num_samples, axis=0) # 第二步:带通滤波,保留5-200Hz频段(典型机械故障特征频段) b, a = signal.butter(4, [5, 200], btype='band', fs=fs_target) X_filtered = signal.filtfilt(b, a, X_resampled, axis=0) # 第三步:按通道做Z-score标准化 scaler = StandardScaler() X_scaled = scaler.fit_transform(X_filtered) print(f"处理后数据维度: {X_scaled.shape}") print(f"各通道均值: {X_scaled.mean(axis=0).round(4)}") print(f"各通道标准差: {X_scaled.std(axis=0).round(4)}")这段代码的逻辑链条是:signal.resample用傅里叶方法做重采样,比简单插值更保频谱;signal.butter设计4阶巴特沃斯带通滤波器,filtfilt做零相位滤波避免时延;StandardScaler按列减均值除标准差。参数方面,fs_target要根据你的故障特征频率来定——如果轴承故障特征频率在几百赫兹,采样率至少要到2kHz以上;[5, 200]这个频段是旋转机械的通用选择,但如果是电气类故障,需要往下调到0.5-50Hz。滤波阶数4阶是折中,阶数太高会振铃,太低则过渡带太宽。
注意:
filtfilt要求信号长度至少是滤波器阶数的3倍以上,数据太短会报错。如果单段信号不足,先做分段再滤波。
2.3 数据分段与样本增强的实操参数
时序信号做分类,不能整段丢进去,要切成固定长度的窗口。窗口长度怎么定?经验公式是至少覆盖3-5个故障特征周期。比如故障特征频率50Hz,周期20ms,窗口至少60-100ms,对应1000Hz采样率就是60-100个点。重叠率一般取50%,既能增加样本量又不至于引入太多冗余。
def segment_signal(X, y, window_size=128, overlap=0.5): step = int(window_size * (1 - overlap)) segments = [] labels = [] for i in range(0, X.shape[0] - window_size, step): seg = X[i:i + window_size, :] # 取窗口内多数标签作为该段标签 label = np.bincount(y[i:i + window_size].astype(int)).argmax() segments.append(seg) labels.append(label) return np.array(segments), np.array(labels) X_seg, y_seg = segment_signal(X_scaled, y, window_size=128, overlap=0.5) print(f"分段后样本数: {X_seg.shape[0]}, 每样本形状: {X_seg.shape[1:]}")window_size=128在1000Hz采样率下对应128ms,覆盖2-3个工频周期,对大多数旋转机械够用。overlap=0.5是常用值,如果样本严重不足可以提到0.75,但要注意验证集和训练集的分段不能来自同一段原始信号,否则数据泄漏会让准确率虚高。
3. 故障特征提取:时域、频域与时频域的选型逻辑
3.1 三种特征域的适用边界
时域特征(均值、方差、峭度、裕度)计算快,对冲击性故障敏感,但区分不了故障类型。频域特征(FFT谱峰、谱重心、边频带能量)能定位故障频率,但要求信号平稳。时频域特征(小波包能量、短时傅里叶)适合非平稳信号,但维度高、计算慢。
我的选型习惯是:先看信号平稳性。用ADF检验判断,p值小于0.05就用时域+频域组合,否则上小波包。天府杯A题的数据我印象里是变工况下的轴承信号,非平稳性明显,所以最终用的是小波包能量特征。
3.2 小波包能量特征提取的完整代码
import pywt def wavelet_packet_features(segment, wavelet='db4', level=3): """对单段信号做小波包分解,提取各频带能量""" features = [] for ch in range(segment.shape[1]): sig = segment[:, ch] wp = pywt.WaveletPacket(data=sig, wavelet=wavelet, mode='symmetric', maxlevel=level) # 获取第level层的所有节点 nodes = [node.path for node in wp.get_level(level, 'natural')] energies = [] for n in nodes: coeffs = wp[n].data energies.append(np.sum(coeffs ** 2)) energies = np.array(energies) # 归一化为能量占比 energies = energies / (energies.sum() + 1e-12) features.extend(energies) return np.array(features) # 对所有分段样本提取特征 feature_list = [wavelet_packet_features(seg) for seg in X_seg] X_features = np.array(feature_list) print(f"特征矩阵维度: {X_features.shape}")wavelet='db4'是Daubechies小波,对机械冲击信号匹配度好。level=3把信号分成8个频带,每个通道8维特征,如果原始有4个通道就是32维。能量占比归一化是为了消除量纲影响。如果特征维度还是太高,可以用PCA降到10-15维,保留95%方差即可。
提示:小波包分解前要确保信号长度是2的level次方的整数倍,否则
pywt会自动填充,填充方式选symmetric比zero更不容易引入边界突变。
3.3 特征筛选:别把所有特征都喂给模型
32维特征不算多,但里面肯定有冗余。我用的是互信息+递归特征消除的组合:先算每个特征与标签的互信息,筛掉低于阈值的,再用RFE交叉验证选最优子集。
from sklearn.feature_selection import mutual_info_classif, RFE from sklearn.svm import SVC # 互信息筛选 mi = mutual_info_classif(X_features, y_seg, random_state=42) mi_threshold = np.percentile(mi, 30) # 去掉最低30% selected_idx = np.where(mi > mi_threshold)[0] X_selected = X_features[:, selected_idx] print(f"互信息筛选后维度: {X_selected.shape[1]}") # RFE进一步筛选 estimator = SVC(kernel='linear', random_state=42) rfe = RFE(estimator, n_features_to_select=12, step=1) X_rfe = rfe.fit_transform(X_selected, y_seg) print(f"RFE筛选后维度: {X_rfe.shape[1]}")互信息阈值取30%分位数是经验值,如果特征本身就不多可以放宽到20%。RFE的n_features_to_select=12是根据样本量定的——样本数除以10大致是安全上限,样本500个左右选12-15个特征比较稳。
4. 诊断模型选型与训练:从SVM到一维CNN的取舍
4.1 小样本下为什么优先选SVM而不是深度学习
数学建模竞赛的数据量通常不大,天府杯A题给的训练样本我印象里在千级以内。这种量级下,一维CNN很容易过拟合,除非你做大量数据增强。SVM在特征工程到位的前提下,小样本表现更稳,训练也快,调参维度低。如果样本超过5000,再考虑上CNN。
SVM的关键参数是C和gamma。C控制惩罚力度,gamma控制核函数影响范围。网格搜索范围:C取[0.1, 1, 10, 100],gamma取[0.001, 0.01, 0.1, 1],用5折交叉验证选最优。
4.2 SVM训练与交叉验证的完整流程
from sklearn.model_selection import GridSearchCV, StratifiedKFold from sklearn.metrics import classification_report, confusion_matrix from sklearn.pipeline import Pipeline # 构建管道:标准化 + SVM pipe = Pipeline([ ('scaler', StandardScaler()), ('svm', SVC(kernel='rbf', probability=True, random_state=42)) ]) # 网格搜索参数 param_grid = { 'svm__C': [0.1, 1, 10, 100], 'svm__gamma': [0.001, 0.01, 0.1, 1] } cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=42) grid = GridSearchCV(pipe, param_grid, cv=cv, scoring='f1_macro', n_jobs=-1) grid.fit(X_rfe, y_seg) print(f"最优参数: {grid.best_params_}") print(f"最优F1: {grid.best_score_:.4f}") # 用最优模型做预测 y_pred = grid.predict(X_rfe) print(classification_report(y_seg, y_pred)) print(confusion_matrix(y_seg, y_pred))scoring='f1_macro'而不是accuracy,是因为故障类别通常不均衡,宏平均F1能反映小类表现。probability=True是为了后续画ROC曲线,但会稍微增加训练时间。如果类别极度不均衡,可以在SVC里加class_weight='balanced'。
4.3 一维CNN的适用场景与最小实现
如果样本量够(每类至少500段),可以试1D-CNN。结构不用太深,两层卷积+两层全连接足够。
import torch import torch.nn as nn class FaultCNN(nn.Module): def __init__(self, in_channels, num_classes): super().__init__() self.conv1 = nn.Conv1d(in_channels, 32, kernel_size=7, padding=3) self.conv2 = nn.Conv1d(32, 64, kernel_size=5, padding=2) self.pool = nn.MaxPool1d(2) self.fc1 = nn.Linear(64 * 32, 128) # 128点输入经过两次pool变32 self.fc2 = nn.Linear(128, num_classes) self.relu = nn.ReLU() self.dropout = nn.Dropout(0.3) def forward(self, x): x = self.pool(self.relu(self.conv1(x))) x = self.pool(self.relu(self.conv2(x))) x = x.view(x.size(0), -1) x = self.dropout(self.relu(self.fc1(x))) return self.fc2(x) # 输入形状: (batch, channels, length) model = FaultCNN(in_channels=X_seg.shape[2], num_classes=len(np.unique(y_seg))) print(model)卷积核7和5对应不同时间尺度,第一层抓短时冲击,第二层抓更宽的模式。Dropout(0.3)是防过拟合的关键,如果训练集准确率远高于验证集,先加大dropout再考虑减层。
5. 避坑与排查:仪器故障诊断里最容易翻车的五个点
5.1 验证集准确率99%但测试集崩了
现象:交叉验证F1到0.98,换一组数据直接掉到0.6。原因:分段时训练集和验证集来自同一段原始信号,窗口重叠导致信息泄漏。解决:按原始信号文件划分训练/验证/测试,同一文件的段只出现在一个集合里。用GroupKFold而不是StratifiedKFold。
5.2 所有样本被预测为多数类
现象:混淆矩阵里少数类全错。原因:类别不均衡+没设class_weight。解决:SVC加class_weight='balanced',CNN的损失函数用带权重的CrossEntropyLoss,权重取类别频率的倒数。
5.3 小波包特征维度爆炸且训练极慢
现象:4通道level=5产生128维特征,SVM训练要十几分钟。原因:分解层数太高,特征冗余严重。解决:level降到3,或者每层只取能量最大的前几个节点,再配合PCA降维。
5.4 滤波后信号幅值异常增大
现象:filtfilt之后信号幅值比原始大好几倍。原因:滤波器阶数太高导致数值不稳定,或者通带设置过窄引起谐振。解决:降低阶数到2-4阶,检查通带范围是否覆盖了信号主频,用signal.freqz看幅频响应。
5.5 标签对齐错误导致模型学反
现象:模型准确率始终在50%左右,像随机猜。原因:分段时标签取了窗口起始点的标签,但故障可能发生在窗口中间。解决:取窗口内多数标签,或者用滑窗内标签变化点做边界检测,确保每段标签一致。
6. 诊断结果的可信度验证:混淆矩阵之外还要看什么
模型跑出高准确率不等于方案可靠。我一般会补三个验证:一是看混淆矩阵的误判方向,如果A类全被判成B类,说明这两类特征空间重叠严重,需要回去补区分性特征;二是画t-SNE看特征聚类,类内散度大说明特征不稳;三是做对抗验证,故意加高斯噪声看F1下降幅度,下降超过15%说明模型鲁棒性不够。
from sklearn.manifold import TSNE import matplotlib.pyplot as plt # t-SNE可视化 tsne = TSNE(n_components=2, perplexity=30, random_state=42) X_tsne = tsne.fit_transform(X_rfe) plt.figure(figsize=(8, 6)) for label in np.unique(y_seg): mask = y_seg == label plt.scatter(X_tsne[mask, 0], X_tsne[mask, 1], label=f'Fault {label}', alpha=0.6) plt.legend() plt.title('t-SNE Feature Visualization') plt.show() # 噪声鲁棒性测试 noise_levels = [0.01, 0.05, 0.1, 0.2] for nl in noise_levels: X_noisy = X_rfe + np.random.normal(0, nl, X_rfe.shape) score = grid.best_estimator_.score(X_noisy, y_seg) print(f"Noise level {nl}: accuracy = {score:.4f}")perplexity=30适合样本量几百到几千的情况,太小会聚成团,太大则散开。噪声测试里,如果0.1噪声下准确率还能保持85%以上,这个模型拿到新数据上才比较放心。我自己的习惯是,任何诊断模型上线前必须过这三关,少一关都不踏实。希望帮到你。
本文还有配套的精品资源,点击获取