企业官网建设流程全解析

说实话,轴承和齿轮这类旋转机械的故障诊断,做了几年现场项目之后,最大的感受就是:数据好采,特征难提。传感器一贴,采集卡一开,几万几十万个点瞬间就来了,但真正到了“判断设备到底是正常还是早期磨损”这一步,很多人就卡住了。这几年国内设备健康管理、预测性维护的项目越来越多,十个需求里有八个都绕不开特征提取这个环节。今天我就以Matlab为工具,把这套在项目里反复打磨过的思路完整梳理一遍——36种覆盖时域、频域、时频域的多征兆域特征提取方法,以及它们在轴承故障诊断、齿轮磨损评估和寿命预测里的实际用法。内容偏实操,从原理到代码,再到排查过程中踩过的坑一次讲透,希望能给正在做同类项目的朋友省下几个月的摸索时间。

1. 多征兆域特征到底解决什么问题,36种特征是怎么分层的

在开始写代码之前,先把这个核心问题想清楚:为什么单靠一两个指标很难诊断故障。

举个例子,滚动轴承出现外圈点蚀时,振动信号的均方根值会升高,峭度也会变大。但如果设备转速不稳定,或者负载忽高忽低,均方根值照样波动。拿单一指标去做判断,误报率会高得让人崩溃。这是我在一个减速机监测项目里用血的教训换来的结论:必须从多个征兆域提取特征,把设备运行状态“投影”到高维空间里,才能做出稳定可靠的判断。

所谓“多征兆域”,我通常分成三个层次:

• 时域特征:直接对原始振动信号做统计计算,包括幅值、能量、分布形态等指标。这类特征计算最快,对冲击类故障敏感,适合在线粗筛。
• 频域特征:对信号做傅里叶变换或包络谱分析,提取频谱形状、谱峰、边频带等参数,专门捕捉轴承、齿轮的固有频率及其调制特征。
• 时频域特征:对非平稳信号做小波包分解、经验模态分解或熵类计算,刻画信号频率成分随时间变化的情况,这是早期微弱故障的最佳表征手段。

36种特征算下来,经典分布是:时域特征20种,频域特征10种,时频域特征6种。这个比例是经过工程验证的——时域特征计算成本低、数值稳定性好,用于批量处理海量历史数据非常友好;频域特征做故障机理解释;时频域特征做精细诊断和趋势预估。三层特征融合后在故障分类和寿命预测上的鲁棒性,远比单域特征高出几个量级。

注意一个容易被忽略的工程事实:特征的数量不是越多越好,覆盖的面越广越好。有些特征之间高度线性相关,比如均方根值和标准差在零均值信号下几乎等价。实际落地时先提取36种全量特征,再通过相关性分析和主成分分析做降维,才是比较稳妥的项目路线。

还有一类十分重要的特征需要单独提出来——用于寿命预测的退化特征。不同于故障诊断关注“有没有坏”“坏在哪”,寿命预测关注“还剩多久坏”。这需要特征具备单调的退化趋势,比如振动能量、谱峰能量、样本熵等,这些指标的数值会随着磨损程度稳定上升或下降。后面我会专门讲如何从这36种特征里筛选出这类“趋势候选者”。

2. 数据准备与信号预处理:特征提取之前必须先做的5步操作

特征提取的精度上限,在数据预处理阶段就已经框死了。预处理做得不到位,后面再花哨的算法也是白搭。以下这5步是我在每个项目里都会标准化执行的流程,建议直接抄。

2.1 从原始振动数据到纯净分析信号的完整流程

采集到的原始信号一般有三个问题:直流偏置、工频干扰、高频噪声。

第一步:去直流偏置。

振动加速度传感器输出的信号经常会叠加一个直流分量,这个分量如果不除掉,时域特征的均值、方根幅值、波形因子会全线失真。Matlab里一条线就能解决:

% x为原始振动信号 x_detrend = x - mean(x);

第二步:带通滤波去噪。

轴承故障的特征频率通常在1kHz到20kHz之间(具体取决于转速和几何尺寸),齿轮啮合频率则跟齿数相关。先用带通滤波器把关注频段框出来,滤掉低频干扰和高频噪声,信噪比会大幅提升。推荐用零相移滤波filtfilt,避免相位偏移影响后续特征频率的定位:

fs = 25600; % 采样率,根据实际采集情况设定 [b, a] = butter(4, [500, 10000]/(fs/2), 'bandpass'); x_filtered = filtfilt(b, a, x_detrend);

这里有个经验参数:带通下限不要低于500Hz,否则容易把转频和倍频的干扰引进来;上限设在10000Hz左右,能保留轴承外圈、内圈故障激发的高频共振成分。

第三步:去除趋势项。

转子升降速过程中信号会带有缓慢变化的趋势,这种趋势会污染频域特征里低频段的计算。用多项式拟合去掉趋势项:

x_poly = polyfit(1:length(x_filtered), x_filtered, 3); x_clean = x_filtered - polyval(x_poly, 1:length(x_filtered))';

第四步:滑窗分段。

故障诊断场景下,整段长信号没法直接提取一个特征集合——因为设备状态在变化。我会把连续信号切成固定长度的窗口,每个窗口独立计算特征。

window_len = 4096; % 窗口长度,对应0.16秒(采样率25600Hz) step = 1024; % 滑窗步长 num_windows = floor((length(x_clean) - window_len) / step) + 1; segments = zeros(num_windows, window_len); for i = 1:num_windows start_idx = (i - 1) * step + 1; segments(i, :) = x_clean(start_idx : start_idx + window_len - 1); end

窗口长度的选取有个原则:窗口内至少要包含20个以上旋转周期。低速设备(比如300rpm, 5Hz)意味着窗口至少1秒;高速设备(比如3000rpm, 50Hz)窗口可以短到0.2秒。我用的是4096个点,在25600Hz采样率下大约0.16秒,对3000rpm的电机来说包含约8个旋转周期,做故障诊断勉强够,做趋势监测建议窗口放大到8192点。

第五步:信号质量校验。

在进入特征提取前做一次“信号体检”:计算信噪比、检查是否存在削顶、断点。削顶信号会导致峭度严重失真。

if max(abs(x_clean)) >= 0.95 * max(abs(x)) warning('信号存在接近满量程的峰值,注意削顶风险'); end

2.2 关于转速基准和角度域重采样

如果你做的是变转速工况下的故障诊断,必须注意:转速波动会让等时间间隔采样的信号发生“频率模糊”,固定的频谱特征会漂移。此时需要做阶比分析,把时间域信号按转速重采样到角度域。在Matlab里可以用orderanalyze相关的工具箱,也可以自己用编码器脉冲信号做角度重采样。

不过我自己做项目时发现,绝大多数场景下设备要么是恒转速、要么是缓慢变速。对于恒转速或缓慢变速的设备,等时间间隔采样加上包络谱分析就够用了;只有快速升降速场景(比如变频驱动启停)才需要真正的阶比跟踪。别一开始就把方案做复杂,先把基础的36种特征在恒转速下跑通,再考虑变转速增强。

3. 36种特征逐域拆解:公式、代码和工程要点

这是全文的核心。我把36种特征按三个域逐一拆开,每个特征不仅给出Matlab实现,还会解释其物理意义、适用场景和常见的误用情况。这块内容当作“特征字典”来用比较合适,写到哪个特征时不记得了回来翻一下。

3.1 时域20种特征:最直接也最容易被低估

先补一个基础概念:对于一组长度为N的离散信号x(n),n = 1,2,...,N,记信号的均值、绝对均值、均方根值、方根幅值分别如下。

• 均值mean:信号的平均水平,反映直流分量,去均值后一般接近0。
• 绝对均值: mean(|x|):全波整流后的平均水平。
• 均方根值rms:信号能量水平,对稳态磨损敏感,公式为sqrt(sum(x.^2)/N)。
• 方根幅值: (mean(sqrt(abs(x))))^2:对小幅值冲击敏感,是一种轻量化的幅值特征。

然后是峰值类特征:

• 最大值max(x),最小值min(x)。
• 峰峰值max(x) - min(x):反映振动摆幅。
• 峰值max(abs(x)):冲击幅度的直接度量,对瞬时冲击类故障敏感,但对噪声也极为敏感。所以提取峰值前一定要先滤波。

无量纲指标类,这类特征在变负载工况下特别好用,因为负载变化会对有量纲指标产生明显影响,无量纲指标相对稳健:

• 波形因子波形因子 = rms / mean(|x|)。
• 峰值因子峰值因子 = max(abs(x)) / rms。
• 脉冲因子脉冲因子 = max(abs(x)) / mean(|x|)。
• 裕度因子裕度因子 = max(abs(x)) / 方根幅值。
• 峭度峭度 = mean((x - mean(x)).^4) / (std(x).^4):对冲击特征极其敏感,是轴承早期故障最经典的指标之一。正常轴承振动近似高斯分布,峭度在3附近;出现局部损伤时信号中出现大量冲击成分,峭度会明显增大。
• 偏度偏度 = mean((x - mean(x)).^3) / (std(x).^3):反映幅值分布的非对称性,齿轮磨损导致的单向冲击会让偏度偏离0。

实测下来,峭度有两个值得注意的陷阱:一是在没有带通滤波的情况下,背景噪声会稀释冲击成分,峭度增长不明显;二是当冲击变成密集群时,振动信号趋于平稳的周期调制,峭度反而可能回落到一个相对较低的水平。所以现场诊断时要结合多个特征来判断,不能死盯峭度。

剩余有量纲幅值统计特征包括峰值指标、波形指标这类衍生参数,还有能量(即sum(x.^2))和一些高阶矩。我把20种时域特征完整罗列如下:

Matlab里计算时域特征写得非常紧凑:

function feat = time_domain_features(x) N = length(x); feat(1) = mean(x); feat(2) = mean(abs(x)); feat(3) = sqrt(mean(x.^2)); feat(4) = (mean(sqrt(abs(x))))^2; feat(5) = std(x); feat(6) = var(x); feat(7) = max(abs(x)); feat(8) = max(x); feat(9) = min(x); feat(10) = max(x) - min(x); feat(11) = feat(3) / feat(2); feat(12) = feat(7) / feat(3); feat(13) = feat(7) / feat(2); feat(14) = feat(7) / feat(4); feat(15) = kurtosis(x); feat(16) = skewness(x); feat(17) = feat(4) / feat(2); feat(18) = sum(x.^2); feat(19) = sum(x.^2) / N; feat(20) = mean((x - mean(x)).^6); end

3.2 频域10种特征:把故障“频率指纹”量化出来

频域特征的思路是,把时域波形通过FFT变换到频谱,再对频谱的形状做统计。这组特征特别擅长刻画与旋转频率、啮合频率相关的故障模式。

基础频域特征计算:

function feat = freq_domain_features(x, fs) N = length(x); X = fft(x); amp = abs(X(1:floor(N/2)+1)); amp(1) = 0; % 去直流 f = linspace(0, fs/2, floor(N/2)+1); % 谱能量相关 s = amp .^ 2; total_energy = sum(s); % 重心频率(频谱的“质量中心”) fc = sum(f .* s) / total_energy; % 均方频率 msf = sum((f.^2) .* s) / total_energy; % 频率方差 vf = sum(((f - fc).^2) .* s) / total_energy; feat(1) = total_energy; feat(2) = fc; feat(3) = msf; feat(4) = vf; feat(5) = sqrt(vf); % 频率标准差 feat(6) = max(amp) / sum(amp); % 主谱峰占比 feat(7) = sum(amp(f > 1000)) / sum(amp); % 高频能量占比 feat(8) = kurtosis(amp); % 频谱峭度 feat(9) = skewness(amp); % 频谱偏度 feat(10) = -sum(s ./ total_energy .* log(s ./ total_energy + eps)); % 谱熵 end

这段代码里包含了10种频域特征,工程上逐个拆解一下:

• 总谱能量:反映了振动总强度,跟时域的均方根值高度相关,作为全局退化指标。
• 重心频率:频谱能量的加权平均频率。低频故障(比如不平衡)会让重心频率下移,高频故障(比如轴承疲劳剥落)会让重心频率上移。
• 均方频率:频率平方的加权平均,强化高频成分的贡献,对早期高频冲击敏感。
• 频率方差:频谱能量在频率轴上展开的程度。频率方差变大说明信号频率成分发散,通常意味着故障产生了新的频率成分。
• 主谱峰占比:最大幅值与总幅值的比值。正常设备的谱峰集中在转频,这个值高;故障设备分散到多个故障特征频率后,这个值会降低。
• 高频能量占比:高频段能量占总能量的比例,直接指示冲击类故障。
• 频谱峭度、频谱偏度:量化频谱分布形态,配合谱熵使用。
• 谱熵:频谱能量分布的均匀程度,反映信号“杂乱程度”。早期故障引入新的频率成分时,谱熵会升高。

注意一个关键工程点:频域特征计算前,时域信号必须经过加窗处理,否则FFT的频谱泄漏会让计算出来的重心频率、均方频率完全失真。Matlab的fft本身没有内建窗函数,需要手动加汉宁窗或者汉明窗:

win = hann(N, 'periodic'); x_win = x .* win; X = fft(x_win);

3.3 时频域6种特征:捕捉故障的“出生和成长”

纯时域和纯频域特征处理平稳信号很有效,但故障早期信号往往是非平稳的——冲击事件随时序随机出现,频率成分随时变化。时频域特征,尤其是熵类指标,在这个环节价值巨大。

我用实战中最常用、效果最好的6种,给出计算公式和Matlab代码。

特征F31 小波包能量熵

小波包分解比小波分解多了一层对高频段的再分解,能把信号精细地划分到多个频带。每个频带的能量分布会直接反映冲击故障落到了哪个频带。

function e = wp_energy_entropy(x) wpt = wpdec(x, 3, 'db4'); % 3层小波包分解 E = wprcoef(wpt, 1:8); % 重建8个子频带 energy = sum(E.^2, 2); total = sum(energy); p = energy / total; e = -sum(p .* log(p + eps)); % 小波包能量熵 end

正常轴承信号的能量分布比较均匀,能量熵较高;故障发生后,冲击能量集中在特定频带,能量熵显著降低。这个特征对早期故障识别非常灵敏,而且对噪声有一定的鲁棒性。

特征F32 EMD能量熵

经验模态分解能把信号分解成若干个本征模态函数(IMF),每个IMF对应不同的振荡模式。做完分解,按IMF能量分布计算熵值:

function e = emd_energy_entropy(x) imfs = emd(x, 'Display', 0); % Matlab自带emd energy = sum(imfs.^2, 2); total = sum(energy); p = energy / total; e = -sum(p .* log(p + eps)); end

这里有个需要提醒的点:EMD有个叫“端点效应”的数学问题,在数据两端容易产生飞翼。实测中把信号按窗口切完后直接做EMD会得到部分失真结果,建议每个窗口在两端各扩展50个点,做完EMD后再裁掉。

特征F33 排列熵

排列熵是检测时间序列复杂度和动力学突变的重要指标,计算的是相邻数据点的排序模式的重复程度。它计算速度快,抗噪能力强,特别适合在线监测场景。

function pe = permutation_entropy(x, m, tau) N = length(x); perm_patterns = zeros(N - (m-1)*tau, 1); for i = 1:(N - (m-1)*tau) pattern = x(i : tau : i + (m-1)*tau); [~, idx] = sort(pattern); perm_patterns(i) = poly2num(idx); % 具体编码方式可自定义 end [counts, ~] = histcounts(perm_patterns, 'BinMethod', 'integers'); p = counts / sum(counts); p = p(p > 0); pe = -sum(p .* log(p)); end

排序模式编码的细节在Matlab中可以直接用sort函数的索引转数字编码实现,关键是参数选择:嵌入维度m通常取4到6,m太小则区分度低,m太大则计算量大且需要长数据;时间延迟τ取1即可。排列熵在齿轮磨损趋势上的表现相当好——磨损加剧时,振动信号的规律性增强,排列熵呈下降趋势。

特征F34 样本熵

样本熵衡量时间序列中新信息出现的速率,样本熵越高,系统越随机。轴承进入故障期后,振动信号中的冲击会交替出现,样本熵会出现先升后降的过程。计算要写循环,我用的是已有的开源自定义函数SampEn,核心代码段如下:

function se = sample_entropy(x, m, r) N = length(x); % 构造m维向量 template = zeros(N-m, m); for i = 1:(N-m) template(i, :) = x(i : i+m-1); end % 统计匹配对数,细节略 end

样本熵算起来比较慢,当处理上万条窗口时尤其明显。工程上建议对样本熵的窗口做降采样,或者改用快速近似熵。

特征F35 模糊熵

模糊熵是样本熵的改进版,用指数函数模糊化相似度的度量,对噪声和参数选择的敏感度更低。模糊熵在轴承早期故障识别中与峭度配合使用效果显著。Matlab代码可以直接搜索开源的FuzzyEntropy函数,参数通常设为m=2, r=0.2*std(x),实测效果稳定。

特征F36 瞬时能量均值

计算信号经带通滤波后的包络(用希尔伯特变换取包络),再对包络序列求均值:

envelope = abs(hilbert(x_band)); instant_energy = envelope .^ 2; feat36 = mean(instant_energy);

这个特征与轴承故障的直接相关性最高——故障脉冲对轴承结构产生激励,包络能量随故障程度加深稳步上升,是寿命预测里最常用的健康因子之一。

3.4 特征体系的构建原则:不要遗漏,也不要重复

36种特征全部提取完后,数据结构应该是一张二维表:每一行是一个时间窗口,每一列是一种特征。实际做项目时,我通常会按下面的逻辑组织特征工程管线:

1. 先分窗,计算36种全量特征。
2. 做标准化(去均值和单位方差),因为峭度、偏度这类无量纲特征与其他有量纲特征在数值尺度上差异很大,直接灌给模型会被中等幅值的大数值特征支配。
3. 计算特征相关矩阵,剔除高度相关的特征(相关系数绝对值>0.95),保留物理意义更明确的那一个。
4. 利用主成分分析或t-SNE做降维和可视化,确认故障样本在特征空间里的可分性。

这里放一个实操细节:标准化选用zscore就行,但对于寿命预测任务,注意rolling窗口的标准化方式——如果对整个历史数据一次性做zscore,相当于“偷看未来”,对未来样本不公平。正确的做法是滑动窗口在线zscore,即用当前时刻之前的数据计算均值和方差来标准化当前样本。

4. 从特征向量到故障诊断与寿命预测:完整的落地流程

拿到36种特征后,后面的路径就清晰了。我在这里给出一套经过数次现场测试的完整流程,同时注明参数选择的原因。

4.1 滑动窗口特征矩阵构建

以某齿轮箱振动数据为例(实测采样率25600Hz,持续5分钟,电机转速约1450rpm)。实际操作中我选择窗口长度8192点(0.32秒),步长2048点,得到约375个窗口。每个窗口提取36维特征,特征矩阵X大小是375行乘36列。

选择8192点窗口有两大原因:一是在1450rpm(转频24Hz)下,一个窗口包含约8个旋转周期,足够捕捉到轴承外圈故障的特征频率及其边频带;二是8192是2的13次方,FFT计算效率高,而且频率分辨率达到3.125Hz,足以分辨轴承故障特征频率与相邻谱峰的间距。

4.2 二次特征:包络谱与边频带系数

除了上述36种基础特征,在故障诊断环节我通常还会额外计算一组“二次特征”——包络谱峰值因子和边频带能量占比。这类特征严格来说可以并入频域特征族,但它们对齿轮和轴承故障的机理表征能力太强,单独列出来讲更有助于理解多征兆域融合的含义。

% 包络谱计算 envelope = abs(hilbert(x_band)); env_fft = abs(fft(envelope)); f_axis = linspace(0, fs/2, length(env_fft) / 2);

在包络谱上,如果存在轴承外圈故障,会在外圈故障特征频率BPFO及其倍频处出现明显峰值。BPFO的计算公式是:

BPFO = (n / 2) * f_r * (1 - (d / D) * cos(alpha))

其中n是滚动体个数,fr是转频,d是滚动体直径,D是节圆直径,alpha是接触角。这些参数在轴承型号手册里可以查到。齿轮故障则更关注啮合频率f_m = z * f_r(z为齿数)及两侧的边频带。以啮合频率为中心,两侧间隔为故障轴转频的边频带幅值比值,就是“边频带能量占比”特征。

诊断时,把各故障特征频率处的幅值、边频带能量、基础频域特征组合起来,输入分类器,效果会远超只用单一频谱峰值。

记得在调试时先把包络谱画出来人工看一眼——有一次程序误把工频谐波当成了边频带,自动特征提取算法完全无误,但画图后发现峰值是电网干扰的谐波成分。算法判断前的可视化验证永远不能省略。

4.3 故障模式识别:基于特征的分类模型搭建

特征全部提取好之后,故障分类就是一个标准的机器学习问题。在工程中我用下来,几种模型按优先级排列:

• 随机森林:对特征重要性有内置评估,能容忍特征间的非线性关系,小样本情况下表现稳定。在第一个项目里我用它实现了96%以上的准确率。
• 支持向量机(RBF核):样本量小于1000时经典且强,但需要调C和gamma,建议配合网格搜索。
• 轻量梯度提升树(LightGBM):数据量大、特征维度高时的首选,速度快,精度高。
• 深度神经网络:仅当样本量超过10万时才推荐,否则容易过拟合且解释性差。

分类模型本身不是重点,重点是特征的判别力是否足够。我会先用主成分分析把36种特征降到2维或3维,直接可视化看故障样本的聚类情况。如果散点图上正常、外圈故障、内圈故障、齿轮磨损各自抱团,然后再上分类器;如果类别重叠严重,说明特征集不够,需要回到包络谱或时频域层面补充特征。

[coeff, score, latent] = pca(feat_matrix); figure; gscatter(score(:,1), score(:,2), labels);

这个可视化简直是个照妖镜。有一次模拟混合故障数据,时域特征在主成分散点图上几乎完全混叠,后来加了边频带能量占比和样本熵之后,类别才分离开。多征兆域特征的价值,在这一刻体现得最直观。

4.4 寿命预测:健康因子构造与剩余寿命拟合

寿命预测和故障诊断在思路上有个关键差异:诊断关心的是“现在处于哪个状态”,预测关心的是“这个状态在什么时间到达极限”。

我把36种特征全部按时间顺序绘制成曲线后发现:均值、均方根、峰值、总谱能量、包络能量、小波包能量熵这6个指标呈现出不错的单调上升趋势,适合构造健康因子。实际操作中,一种可靠的做法是:

第一步:构造健康因子。把上述趋势特征做加权融合,权重由各特征与“真实退化程度”的相关性决定。如果样本有从安装到失效的完整生命周期数据,可以按时间点标记健康度从1到0,然后拟合回归得到权重。

第二步:提取退化起始点。用滑窗的方法计算特征数列的斜率均值,超过阈值时标记为退化开始。这个点定得准不准,直接影响寿命预测误差。

第三步:趋势外推。我建议用基于指数退化模型的非线性拟合,因为机械磨损的典型曲线是指数型的——早期平缓,后期急剧上升。

% 假设hi为健康因子时间序列,t为时间(单位:小时) ft = fittype('a * exp(b * t) + c', 'independent', 't'); model = fit(t, hi, ft, 'StartPoint', [0.1, 0.001, 0.9]); future_t = t(end):0.1:(t(end) + 2000); predicted_hi = model(future_t); failure_time = future_t(find(predicted_hi >= 0.99, 1)); rul = failure_time - t(end);

这里有个很容易踩的坑:直接对全局数据拟合指数模型,会把早期平稳段的微小波动也纳入拟合,导致外推偏差。我的做法是只取“退化起始点”到当前时刻的数据段进行拟合,起始点之前的数据当成正常段忽略。这样可以大幅提升预测稳定性。

此外,寿命预测一定要给出不确定性区间。简单的方式是用贝叶斯线性回归或高斯过程回归代替普通的最小二乘拟合,这样不仅能得到期望值,还能得到置信区间,现场做决策时更有参考意义。

5. 常见问题与排查技巧实录

这套流程跑了多个项目之后,我把高频问题整理成一个速查表,按频率排序,每一条都融合了实际的排查经验。

5.1 最容易被忽略的坑:工况混叠

实际现场数据往往包含多个工况:负载变化、转速波动、温度变化。同一个轴承在空载和满载下的特征幅值差异,可能比“正常与故障”的差异还大。训练集如果只包含空载数据,测试时一旦加载到满载,模型分分钟就会失效。

解决办法是在建训练集时有意识地覆盖不同工况,或者在特征层面进行工况归一化——把特征除以对应工况下的历史基线值。这是我在某风电齿轮箱项目里踩过最深的一个坑,之后凡是现场数据都先做工况标记。

5.2 那个让所有人都困惑的“样本不平衡”问题

故障样本在自然工况中极少出现。训练分类器时,健康样本占了95%以上,故障样本只有几十条。这种情况下,随机森林等模型会直接倾向于把样本全部判为“健康”,准确率“虚高”到95%以上——但模型其实什么都没学到。

处理办法:

• 先用无监督异常检测(如一分类支持向量机、孤立森林)划分“正常”与“异常”,再做有监督分类区分故障类型。
• 对故障样本做合成少数类过采样技术(SMOTE),生成合成样本来平衡类别。
• 评估指标不要用准确率,改用F1-score、召回率、精确率。

5.3 现场部署的轻量化考虑

特征提取算法的计算耗时大量集中在EMD和小波包分解上。在一套现场在线监测系统里,如果每个窗口都做EMD加小波包,普通工控机会非常吃力。

我的优化做法是:在线粗筛阶段只用20种时域特征加10种频域特征(FFT很快);当异常评分超过阈值时,才触发时频域特征做精细分析。这样既保住了实时性,又保留了诊断深度。实测下来,30种特征在单核处理器上处理一帧8192点信号,耗时从每秒上百帧压缩到几个毫秒,完全够用。

6. 案例复盘:从时域特征开始直到寿命预测的完整实战

最后用一个具体的模拟项目把整个流程串起来。项目背景是某旋转设备试验台,轴承数据包含正常、外圈故障、内圈故障、滚动体故障四种状态,另有3组从早期磨损到失效的完整生命周期数据用于寿命预测。

我按以下步骤执行:

第一步:预处理。原始采样率48kHz,带通滤波设为1000Hz到16000Hz,去直流、去趋势。窗口长度设为12000点(0.25秒),步长3000点,共得到约800个窗口。

第二步:特征提取。写了一个主循环,遍历所有窗口,每个窗口提取36种基础特征,另加包络谱的故障特征频率幅值和边频带能量占比,合计42维特征。

feature_matrix = zeros(num_windows, 42); for i = 1:num_windows seg = segments(i, :); idx_tf = 1:20; feature_matrix(i, idx_tf) = time_domain_features(seg); feature_matrix(i, idx_tf(end)+1 : idx_tf(end)+10) = freq_domain_features(seg, fs); feature_matrix(i, 31:36) = tfr_domain_features(seg); % 自定义封装时频域6特征 feature_matrix(i, 37:42) = envelope_additional_features(seg, fs); end

第三步:特征筛选。做zscore标准化,计算相关矩阵,剔除相关系数大于0.95的特征后剩29维;再做主成分分析,取累计贡献率95%以上的主成分,实际保留了7个主成分,信息保留度很好。

第四步:训练模型。划分训练集与测试集,用随机森林进行五折交叉验证,准确率98.5%,F1 scores平均0.98。比只用时域特征(准确率91%)提升明显,增幅主要来自频域和高频包络特征对轴承故障的识别力。

第五步:寿命预测。对3组生命周期数据,分别计算健康因子序列(由能量类特征加权融合),用指数退化模型拟合退化段,预测结果与真实寿命的误差在10%以内,满足工程指导需求。

这套流程的可复用性很高,在没有公开数据集的情况下也可以利用仿真信号先把整个pipeline跑通——在Matlab里用bearing相关的工具箱生成仿真故障信号,或者自己用sin加随机冲击脉冲构造模拟数据。重点是先把流程架构跑通,再替换真实数据时只需要调整参数,不改变整体框架。

最后单独说一点我个人的体会:特征提取这套东西,代码写得再漂亮也只是工程的一半。另一半是理解数据从哪里来、设备是什么结构、故障机理是怎样的。有一次程序判断某轴承为内圈故障,但现场拆解后发现是外圈故障——因为当时负载方向导致外圈故障被遮挡,解调谱上内圈故障特征频率反而更明显。这种案例提醒我,特征数据永远要和设备原理、工况条件联合解读。如果引导读者养成“先画图验证、再建模型决策”的习惯,这篇分享的最大价值就算兑现了。