1. 項目背景與核心挑戰去年參加天府杯數學建模競賽的經歷現在回想起來依然覺得收獲頗豐。我們團隊當時選的是A題關于儀器故障智能診斷。這個題目乍一看感覺像是傳統工業領域的問題但組委會給的數據集和問題描述一下子就把我們拉到了數據科學和智能算法的前沿。題目要求我們基于給定的傳感器時序數據構建一個能夠自動、精準識別多種潛在故障模式的智能診斷系統。這不僅僅是套用一個現成的分類模型那么簡單它涉及到信號處理、特征工程、模型選擇、結果可解釋性等一系列環環相扣的挑戰。當時我們面臨的核心痛點非常明確第一數據是典型的多維時間序列包含了振動、溫度、壓力等多種傳感器在不同時間點的讀數噪聲大、維度高直接喂給模型效果肯定不好。第二故障模式并非獨立發生早期故障信號極其微弱容易被噪聲淹沒如何從海量數據中提取出對故障敏感的特征是診斷準確率提升的關鍵。第三賽題不僅要求診斷出故障類型還希望我們對故障的嚴重程度或發展階段做出評估這要求模型具備一定的回歸或排序能力。最終我們團隊通過一套融合了信號處理、傳統機器學習與深度學習的混合策略成功解決了這些問題并拿到了一等獎。這篇文章我就把我們的解題思路、技術細節以及用Python實現的核心代碼毫無保留地分享出來。無論你是正在備戰數學建模競賽還是對工業預測性維護、時序數據分析感興趣相信都能從中獲得直接的啟發和可復用的代碼。2. 解題總覽從問題定義到技術路線圖面對“儀器故障智能診斷”這樣一個開放性問題第一步也是最關鍵的一步就是明確我們要解決的具體是什么問題。組委會提供的數據通常是一個包含多個csv文件的數據包每個文件可能對應一臺設備、一段時間內的運行數據或者不同故障模式下的樣本。列通常包括時間戳、若干傳感器通道如acc_x,acc_y,temp,pressure等以及一個label或fault_type列在訓練集中。我們的目標可以拆解為三個層次故障檢測判斷設備在某個時間窗口內是否發生了故障二分類正常 vs 異常。故障識別如果發生故障具體是哪種類型多分類如軸承內圈故障、外圈故障、齒輪磨損等。故障程度評估量化故障的嚴重性回歸或有序分類如輕微、中等、嚴重。技術路線上我們沒有押寶單一模型而是設計了一個分階段的流水線Pipeline這樣既能保證基礎模型的穩健性又能利用深度模型挖掘深層特征。整體流程如下原始時序數據 - 數據預處理與清洗 - 時域/頻域/時頻域特征提取 - 特征選擇 - (路徑A)傳統機器學習模型 - (路徑B)深度學習模型 - 模型融合與決策 - 結果輸出與可視化這個雙路徑設計是我們的核心策略。路徑A傳統機器學習依賴精心設計的特征工程模型如XGBoost、LightGBM解釋性強訓練快。路徑B深度學習如1D-CNN、LSTM能自動學習特征對原始數據中的復雜模式捕捉能力更強。兩者優勢互補通過加權投票或堆疊Stacking方式融合最終診斷的魯棒性和準確性顯著提升。3. 數據預處理為模型提供“干凈”的燃料原始工業傳感器數據幾乎不可能是完美無缺的。直接建模等于讓模型在噪音中學習事倍功半。我們的預處理步驟主要解決以下四個問題3.1 缺失值與異常值處理傳感器可能短暫失靈產生NaN或明顯超出物理量程的異常值。import pandas as pd import numpy as np def handle_missing_and_outliers(df, sensor_columns): 處理缺失值和基于標準差/分位數的異常值。 df: 包含傳感器數據的DataFrame sensor_columns: 傳感器列名的列表 df_filled df.copy() # 1. 缺失值處理對于時間序列用前后時刻的均值填充更合理 for col in sensor_columns: df_filled[col] df_filled[col].interpolate(methodlinear) # 線性插值 # 如果開頭或結尾還有NaN用最近的有效值填充 df_filled[col] df_filled[col].fillna(methodbfill).fillna(methodffill) # 2. 異常值處理使用基于IQR四分位距的方法 for col in sensor_columns: Q1 df_filled[col].quantile(0.25) Q3 df_filled[col].quantile(0.75) IQR Q3 - Q1 lower_bound Q1 - 1.5 * IQR upper_bound Q3 1.5 * IQR # 將超出邊界的值替換為邊界值或視為缺失值再填充 df_filled[col] np.where((df_filled[col] lower_bound) | (df_filled[col] upper_bound), np.nan, df_filled[col]) # 再次填充因異常值替換產生的NaN df_filled[col] df_filled[col].interpolate(methodlinear).fillna(methodbfill).fillna(methodffill) return df_filled注意對于高頻振動信號簡單的插值可能會引入虛假頻率成分。在要求極高的場景下需要考慮更專業的信號處理方法如基于模型預測的插值。但在數學建模的有限時間內IQR插值是穩健且高效的選擇。3.2 數據標準化與平滑不同傳感器量綱和量級差異巨大例如加速度單位是g溫度是攝氏度。必須進行標準化防止量級大的特征主導模型。我們通常使用StandardScaler減去均值除以標準差因為它能保留數據的分布形狀對后續的PCA等線性變換友好。同時為了抑制高頻噪聲可以對信號進行滑動平均濾波。from sklearn.preprocessing import StandardScaler def normalize_and_smooth(df, sensor_columns, window_size5): 標準化并應用簡單的移動平均平滑。 window_size: 滑動窗口大小需為奇數。 df_processed df.copy() scaler StandardScaler() # 先平滑再標準化順序有時有影響可根據實驗調整 for col in sensor_columns: # 滑動平均平滑 df_processed[col] df_processed[col].rolling(windowwindow_size, centerTrue, min_periods1).mean() # 標準化 df_processed[sensor_columns] scaler.fit_transform(df_processed[sensor_columns]) # 保存scaler用于后續的測試數據轉換 return df_processed, scaler3.3 樣本切片與標簽對齊原始數據是長序列但我們需要將其切割成固定長度的時間窗口作為模型的一個個“樣本”。這里的關鍵是標簽對齊一個時間窗口對應一個故障標簽。通常我們假設在一個短時間窗口內故障類型是穩定的。采用滑動窗口方法進行切片并可以設置重疊overlap以增加樣本量。def create_samples(data_sequence, labels, window_size, step_size): 將長序列切割成重疊的時間窗口樣本。 data_sequence: 形狀為 (n_timesteps, n_features) 的傳感器數據數組 labels: 形狀為 (n_timesteps,) 的標簽數組每個時間點一個標簽 window_size: 窗口長度時間步數 step_size: 滑動步長 X, y [], [] n_samples len(data_sequence) for start in range(0, n_samples - window_size 1, step_size): end start window_size X.append(data_sequence[start:end]) # 取窗口內最主要的標簽作為該樣本的標簽對于分類問題 window_labels labels[start:end] from scipy.stats import mode label, _ mode(window_labels, keepdimsFalse) y.append(label) return np.array(X), np.array(y)實操心得window_size和step_size是超參數。window_size要足夠長以包含故障特征周期可通過分析故障頻率初步估算但太長會混入不同狀態的信息。step_size小于window_size會產生重疊樣本能有效增加數據量防止切割時恰好切掉關鍵特征但也會引入樣本相關性。我們通常設置重疊率為50%。4. 特征工程從原始信號中“榨取”信息這是傳統機器學習路徑路徑A的成敗關鍵。好的特征應該對故障敏感同時對工況變化如轉速、負載相對魯棒。我們從三個域進行特征提取4.1 時域特征直接從時間序列的幅值統計信息中提取計算簡單物理意義明確。import numpy as np from scipy import stats def extract_time_domain_features(signal): 提取單個傳感器通道在一個時間窗口內的時域特征。 features {} features[mean] np.mean(signal) features[std] np.std(signal) features[rms] np.sqrt(np.mean(signal**2)) # 均方根值反映能量 features[peak] np.max(np.abs(signal)) # 峰值 features[skewness] stats.skew(signal) # 偏度衡量分布不對稱性 features[kurtosis] stats.kurtosis(signal) # 峰度衡量分布尖銳程度 features[crest_factor] features[peak] / features[rms] if features[rms] ! 0 else 0 # 峰值因子 features[clearance_factor] features[peak] / (np.mean(np.sqrt(np.abs(signal)))**2) if np.mean(np.sqrt(np.abs(signal))) ! 0 else 0 # 裕度因子 # 還可以增加波形因子、脈沖因子等 return features4.2 頻域特征故障常常在振動信號的頻譜中表現出特定的頻率成分如軸承的故障特征頻率。通過快速傅里葉變換FFT將信號轉換到頻域。from scipy.fft import fft, fftfreq def extract_freq_domain_features(signal, sampling_rate): 提取頻域特征。 signal: 時間窗口信號 sampling_rate: 采樣頻率 (Hz) n len(signal) yf fft(signal) # 取單邊頻譜 yf_abs 2.0/n * np.abs(yf[:n//2]) xf fftfreq(n, 1/sampling_rate)[:n//2] features {} features[dominant_freq] xf[np.argmax(yf_abs)] # 主頻 features[dominant_amp] np.max(yf_abs) # 主頻幅值 # 計算頻譜重心、均方頻率、頻率方差等 features[spectral_centroid] np.sum(xf * yf_abs) / np.sum(yf_abs) if np.sum(yf_abs) ! 0 else 0 features[spectral_rms] np.sqrt(np.sum((xf**2) * yf_abs) / np.sum(yf_abs)) if np.sum(yf_abs) ! 0 else 0 # 可以計算特定頻帶如故障特征頻率附近的能量占比 return features4.3 時頻域特征對于非平穩信號即統計特性隨時間變化的信號單純的頻域分析會丟失時間信息。短時傅里葉變換STFT或小波變換能提供聯合時頻信息。我們常用小波包變換WPT因為它能對高頻部分進行更精細的分解適合提取故障引起的瞬態沖擊特征。import pywt # 需要安裝PyWavelets def extract_wavelet_features(signal, waveletdb4, level3): 進行小波包分解并計算各節點子頻帶的能量作為特征。 wp pywt.WaveletPacket(datasignal, waveletwavelet, modesymmetric, maxlevellevel) # 獲取第level層所有節點的名稱如 aaa, aad, ada, ... nodes [node.path for node in wp.get_level(level, natural)] energy_features [] for node_name in nodes: node_coeffs wp[node_name].data node_energy np.sum(node_coeffs**2) energy_features.append(node_energy) # 通常將能量歸一化構成能量分布向量 total_energy np.sum(energy_features) energy_features_norm [e/total_energy for e in energy_features] if total_energy ! 0 else energy_features return energy_features_norm將所有傳感器通道、所有域的特征拼接起來會得到一個高維特征向量。接下來必須進行特征選擇去除冗余和無關特征防止“維數災難”。我們使用了基于樹模型如XGBoost的特征重要性排序結合遞歸特征消除RFE來選擇Top-N個最重要的特征。5. 模型構建雙路徑融合策略5.1 路徑A基于特征工程的機器學習模型我們選擇了LightGBM作為主力模型。它訓練速度快對類別不平衡有一定處理能力并且能輸出特征重要性與我們的特征工程流程完美契合。import lightgbm as lgb from sklearn.model_selection import train_test_split, StratifiedKFold from sklearn.metrics import accuracy_score, classification_report, confusion_matrix def train_lightgbm(X_features, y, paramsNone): X_features: 特征工程后得到的特征矩陣 (n_samples, n_features) y: 標簽 if params is None: params { objective: multiclass, # 多分類 num_class: len(np.unique(y)), metric: multi_logloss, boosting_type: gbdt, num_leaves: 31, learning_rate: 0.05, feature_fraction: 0.9, bagging_fraction: 0.8, bagging_freq: 5, verbose: -1, seed: 42 } # 劃分訓練集和驗證集 X_train, X_val, y_train, y_val train_test_split(X_features, y, test_size0.2, stratifyy, random_state42) # 創建Dataset train_data lgb.Dataset(X_train, labely_train) val_data lgb.Dataset(X_val, labely_val, referencetrain_data) # 訓練使用早停法防止過擬合 model lgb.train(params, train_data, valid_sets[val_data], num_boost_round1000, callbacks[lgb.early_stopping(stopping_rounds50), lgb.log_evaluation(period100)]) # 驗證集評估 y_pred model.predict(X_val, num_iterationmodel.best_iteration) y_pred_class np.argmax(y_pred, axis1) print(fValidation Accuracy: {accuracy_score(y_val, y_pred_class):.4f}) print(classification_report(y_val, y_pred_class)) # 可視化特征重要性 lgb.plot_importance(model, max_num_features20, figsize(10,6)) return model5.2 路徑B基于原始信號的深度學習模型我們設計了一個結合1D-CNN和LSTM的混合網絡。CNN擅長提取局部空間特征如振動信號中的沖擊波形LSTM擅長捕捉時間依賴關系。模型直接輸入標準化后的原始時序窗口數據(window_size, n_sensors)。import tensorflow as tf from tensorflow.keras import layers, models, callbacks def build_hybrid_cnn_lstm(input_shape, num_classes): 構建1D-CNN LSTM混合模型。 input_shape: (window_size, n_sensors) model models.Sequential([ # 第一部分1D-CNN 提取局部特征 layers.Input(shapeinput_shape), layers.Conv1D(filters64, kernel_size3, activationrelu, paddingsame), layers.BatchNormalization(), layers.MaxPooling1D(pool_size2), layers.Conv1D(filters128, kernel_size3, activationrelu, paddingsame), layers.BatchNormalization(), layers.MaxPooling1D(pool_size2), layers.Dropout(0.3), # 第二部分LSTM 捕捉時序依賴 # 將CNN輸出的序列輸入到LSTM。return_sequencesTrue表示輸出每個時間步的狀態。 layers.LSTM(units64, return_sequencesTrue), layers.Dropout(0.3), layers.LSTM(units32), layers.Dropout(0.3), # 第三部分全連接層分類 layers.Dense(units64, activationrelu), layers.Dense(unitsnum_classes, activationsoftmax) ]) model.compile(optimizertf.keras.optimizers.Adam(learning_rate0.001), losssparse_categorical_crossentropy, metrics[accuracy]) model.summary() return model # 訓練深度學習模型 def train_deep_model(model, X_train_seq, y_train, X_val_seq, y_val, epochs50): X_train_seq: 原始序列樣本形狀 (n_samples, window_size, n_sensors) early_stopping callbacks.EarlyStopping(monitorval_loss, patience10, restore_best_weightsTrue) reduce_lr callbacks.ReduceLROnPlateau(monitorval_loss, factor0.5, patience5, min_lr1e-6) history model.fit(X_train_seq, y_train, validation_data(X_val_seq, y_val), epochsepochs, batch_size32, callbacks[early_stopping, reduce_lr], verbose1) return model, history踩坑實錄直接訓練這個混合網絡很容易過擬合尤其是在數據量有限的情況下。我們采用了強力的正則化策略除了網絡結構中的Dropout和BatchNorm還在數據上做了隨機縮放、添加高斯噪聲、時間軸輕微扭曲等數據增強顯著提升了模型的泛化能力。另外LSTM層對輸入數據的標準化非常敏感務必確保輸入數據已標準化。5.3 模型融合112的策略我們采用了加權投票法進行融合。兩個模型在驗證集上的準確率作為其權重的基礎。def weighted_ensemble_predict(model_lgb, model_dl, X_feat, X_seq, weightsNone): 加權投票融合。 model_lgb: LightGBM模型輸入特征工程后的數據X_feat model_dl: 深度學習模型輸入原始序列數據X_seq weights: 兩個模型的權重列表如 [0.4, 0.6]。默認為None則根據驗證集準確率自動計算。 proba_lgb model_lgb.predict(X_feat, num_iterationmodel_lgb.best_iteration) # 已經是概率形式 proba_dl model_dl.predict(X_seq) if weights is None: # 這里假設我們已經有了兩個模型在某個驗證集上的準確率 acc_lgb, acc_dl # 例如acc_lgb 0.92, acc_dl 0.94 acc_lgb, acc_dl 0.92, 0.94 total_acc acc_lgb acc_dl weights [acc_lgb/total_acc, acc_dl/total_acc] # 加權平均概率 weighted_proba weights[0] * proba_lgb weights[1] * proba_dl final_pred np.argmax(weighted_proba, axis1) return final_pred, weighted_proba融合后我們在測試集上的準確率比單一的最佳模型通常是深度學習模型提升了約1-2個百分點更重要的是對于某些單一模型容易混淆的故障類別融合模型的判斷更加穩定。6. 故障嚴重程度評估與結果可視化對于故障程度評估我們將其建模為一個**有序分類Ordinal Regression**問題而不是簡單的多分類或回歸。因為“輕微”、“中等”、“嚴重”之間存在明確的順序關系。我們使用了“序數邏輯回歸”的思想將其轉化為多個二分類問題例如模型1區分“無/輕微” vs “中等/嚴重”模型2區分“無/輕微/中等” vs “嚴重”或者直接使用支持有序分類的損失函數如CORAL損失函數在神經網絡中的實現。結果可視化對于診斷系統的可解釋性至關重要。我們主要做了以下幾類圖混淆矩陣熱力圖清晰展示模型在各類別上的混淆情況。特征重要性條形圖從LightGBM模型獲取告訴我們哪些傳感器、哪些特征對診斷貢獻最大這對于后續的傳感器優化布置有指導意義。t-SNE/PCA降維圖將高維特征或深度學習模型最后一層隱藏層的輸出降到2維或3維進行可視化觀察不同故障類別的樣本在特征空間是否能夠被良好區分。關鍵傳感器信號對比圖將正常狀態和不同故障狀態下的關鍵傳感器如振動最大的那個原始信號或頻譜圖畫在一起直觀展示故障特征。import matplotlib.pyplot as plt import seaborn as sns from sklearn.manifold import TSNE def visualize_tsne(features, labels, titlet-SNE Visualization of Features): 使用t-SNE對高維特征進行降維可視化。 tsne TSNE(n_components2, random_state42, perplexity30) features_2d tsne.fit_transform(features) plt.figure(figsize(10,8)) scatter plt.scatter(features_2d[:,0], features_2d[:,1], clabels, cmaptab20, alpha0.7, s10) plt.colorbar(scatter) plt.title(title) plt.xlabel(t-SNE Component 1) plt.ylabel(t-SNE Component 2) plt.tight_layout() plt.show()7. 參賽總結與可復現性建議回顧整個項目拿到一等獎的關鍵在于系統性的問題拆解和務實的技術選型。我們沒有追求最花哨的模型而是確保數據預處理、特征工程、基礎模型訓練每個環節都扎實可靠最后用融合策略提升天花板。有幾個特別重要的點想分享關于數據數學建模競賽給的數據往往“不完美”可能存在標簽噪聲、傳感器漂移等問題。我們花了近三分之一的時間在數據探索和清洗上這是后續所有工作的基石。可視化每一類故障的典型信號波形和頻譜能建立直觀認識甚至能發現數據本身可能存在的問題。關于特征時域、頻域、時頻域特征各有千秋。對于周期性明顯的故障如軸承頻域特征非常有效對于瞬態沖擊故障如齒輪斷齒小波包能量特征可能更好。不要盲目堆砌特征一定要結合特征重要性分析進行篩選。關于模型LightGBM這類樹模型對特征工程的質量要求高但訓練快、調參相對簡單、解釋性強非常適合作為基線模型和提供特征重要性。深度學習模型潛力大但依賴大量數據和高超的調參技巧防止過擬合。雙路徑并行的策略讓我們在有限時間內既能有一個穩健的保底方案又能沖擊更高的性能。關于代碼在競賽中代碼的可復現性和模塊化至關重要。我們將整個流程封裝成多個函數和類數據加載、預處理、特征提取、模型訓練、評估可視化使得調整參數、更換模型、交叉驗證變得非常方便。最終提交的論文中清晰的流程圖和核心代碼片段也是加分項。如果你想在自己的項目或未來的競賽中復現這套方法我的建議是從理解數據開始畫出數據分布聽聽“數據的聲音”。先搭建一個簡單的基線系統比如只用時域特征LightGBM確保整個Pipeline能跑通。迭代優化在此基礎上逐步加入頻域特征、嘗試深度學習模型、調整融合策略。每次只改變一個變量評估其效果。重視驗證策略使用分層K折交叉驗證來更穩健地評估模型性能避免因為數據劃分的偶然性導致過擬合。這個項目讓我深刻體會到解決一個復雜的工程問題往往不是靠一個“銀彈”算法而是靠對問題的深刻理解、扎實的基礎工作以及將多種工具巧妙組合的系統性思維。希望這份詳細的總結和代碼能為你打開一扇門助你在智能診斷或相關的數據科學道路上走得更遠。