归一化判别图嵌入NDGE:故障诊断特征降维的Matlab实现

📅 发布时间:2026/10/21 6:25:18
归一化判别图嵌入NDGE:故障诊断特征降维的Matlab实现
做故障诊断的同学应该都深有体会原始特征一旦上了几十维分类器就开始变得不好伺候。归一化判别图嵌入NDGE是我在实际项目里用来处理多传感器故障特征的一套组合方案它在判别图嵌入的基础上加了拉普拉斯归一化输出侧也比一般降维方法更完整最终投影矩阵、不同维度下的准确率曲线、每个测试样本对不同故障模式的概率全都给到。配套的Matlab代码可以直接改数据路径跑通。无论你是刚接触流形学习的新手还是已经在用PCA、LDA做特征压缩的老手这篇内容都能提供一条可复现的故障诊断技术路线。有人可能会问做故障诊断为什么非要用降维直接把手头几十个特征扔给随机森林、支持向量机不行吗答案是能跑但不好用。实际处理振动信号时原始特征里经常混着大量冗余信息传感器越多、统计指标算得越全特征维度窜得越快。维度一高需要的训练样本量就呈指数增长分类器的泛化能力反而下降这就是常说的“维度灾难”。降维不只是为了减少计算量更重要的是把同类样本拉近、把异类样本推开好让后续分类器在一个更干净的空间里工作。这个项目到底在解决什么问题1.1 故障诊断为什么需要做降维先看一个典型的故障诊断流程从加速度传感器采到振动信号后通常会在时域算均值、方差、峰值因子、峭度、波形因子再在频域算重心频率、均方频率、频率方差条件好的话还能做小波包分解把各频带能量拼成特征向量。这样一趟下来几十维很轻松上百维也不稀奇。问题就出在这里。特征维度高本身不是坏事但维度之间的相关性高就是麻烦了。很多特征本质上是对同一物理现象的重复度量比如峰值因子和峭度都在描述信号中的冲击成分二者高度相关却各自带着噪声。分类器面对这些冗余特征时容易把噪声当成规律学进去导致训练集上表现很好、测试集上一塌糊涂。降维的意义就在这在尽量保留判别信息的前提下把特征压缩到一个低维子空间里。PCA是经典做法但它只看全局方差完全不理会样本属于哪个故障模式压出来的方向往往是方差最大的方向不一定适合分类。LDA则是监督方法直接用类间散度除以类内散度来求方向效果通常比PCA好但有个硬伤——类内散度矩阵在样本量不足时很容易奇异一奇异就没法求逆程序直接崩。NDGE走的是一条介于两者之间的路。它不要求矩阵可逆不假设数据服从高斯分布通过图结构去刻画同类与异类样本之间的几何关系再借助归一化避免数值病态问题。这套思路在处理轴承故障、齿轮箱故障这类非线性工业数据时比PCA、LDA都更稳。1.2 NDGE的定位判别式降维与归一化的组合先说人话版本。NDGE做的事情可以类比成班级排座位想让同班同学坐得近让不同班的人尽量别挨着。传统LDA是拿班级平均位置来算NDGE则细到看每个同学和周围几个人的关系谁跟谁是同类、谁跟谁是异类都直接在图上标出来然后找一组方向按这个方向重新排座位。从技术归类上看NDGE属于判别图嵌入家族和局部保持投影LPP、边界费舍尔分析MFA有血缘关系。它同时构造两个图本征图类内图同类样本之间连边权重越大表示越希望它们在投影后靠近。惩罚图类间图异类样本之间连边权重越大表示越希望它们在投影后分开。目标函数就变成了一个比值问题投影后异类样本之间的离散程度尽量大同类样本之间的离散程度尽量小。这个比值问题进一步转化为广义特征值分解特征值从大到小对应的特征向量就是投影矩阵的列。那“归一化”又体现在哪里在做图嵌入时常用的非归一化拉普拉斯矩阵会因为样本度分布不均导致解偏向度大的节点。NDGE采用对称归一化的方式把拉普拉斯矩阵变换成归一化版本。直观来说这相当于每个节点的影响力先除以它的连接数避免某些“社交达人”样本主导了整个投影方向。这一步数学上只是左乘右乘一个对角阵但数值稳定性明显变好实验重复性也更强。NDGE的算法原理拆解2.1 图嵌入框架下的本征图与惩罚图假设有n个训练样本特征维数是d我们把这些样本看成高维空间里的n个点。第一步是定义点与点之间的相似关系常用的方法是k近邻。对每个点找到距离它最近的k个样本然后在它们之间连边。光连边还不够边上有权重。权重计算最常用的方式是热核[ W_{ij} \exp\left(-\frac{|x_i - x_j|^2}{\sigma}\right) ]这里σ是热核宽度控制着权重随距离衰减的快慢。σ太大所有边的权重都差不多图表达的信息模糊σ太小大量边权重趋近于零图变得稀疏且对噪声敏感。我自己的经验是σ取所有样本对距离的中位数或均值通常效果比较稳后续实验也用了这个策略。本征图和惩罚图都要构建本征图 ( G_w )当且仅当 ( y_i y_j ) 且 ( x_i ) 在 ( x_j ) 的k近邻内时连边权重为 ( W_{ij} )。惩罚图 ( G_b )当且仅当 ( y_i \neq y_j ) 且 ( x_i ) 在 ( x_j ) 的k近邻内时连边权重为 ( W_{ij} )。构建完成后分别计算两个图的拉普拉斯矩阵[ L_w D_w - W_w ] [ L_b D_b - W_b ]其中 ( D_w ) 和 ( D_b ) 是度矩阵对角线上分别是每行权重之和。有一个容易被忽略的细节类内图如果过于稀疏L_w可能非常病态类间图如果太密L_b的信息又被稀释。所以k值既不能太大也不能太小正常范围是5到15。样本量大的时候可以适当增大k让类间图覆盖更多边界样本。2.2 归一化在干什么构建完拉普拉斯矩阵后NDGE引入了对称归一化[ \tilde{L}_w D_w^{-1/2} L_w D_w^{-1/2} ] [ \tilde{L}_b D_b^{-1/2} L_b D_b^{-1/2} ]用生活类比来解释如果把每个样本点看成一个社交账号度矩阵就是它的粉丝数。未归一化的拉普拉斯矩阵在优化时天然更重视粉丝多的节点投影方向会被少数高连接度样本左右。归一化之后每个节点的贡献先除以粉丝数的平方根整体影响力更均衡投影结果对样本分布不均衡的情况也更加鲁棒。从优化目标看未归一化时[ \max_{W} \frac{W^T X L_b X^T W}{W^T X L_w X^T W} ]归一化后变成[ \max_{W} \frac{W^T X \tilde{L}_b X^T W}{W^T X \tilde{L}_w X^T W} ]分母中 ( X \tilde{L}_w X^T ) 在样本量不足时往往奇异归一化并不能百分百消除奇异性但能明显减缓病态程度。在Matlab中求解时建议对分母矩阵加一个微小的正则项[ X \tilde{L}_w X^T \epsilon I ]这里 ( \epsilon ) 建议取 ( 10^{-6} ) 到 ( 10^{-8} ) 之间的值。这个操作几乎不影响投影方向却能避免特征值分解时出现极端大的数值。2.3 最终投影矩阵的求解路径整理一下完整的求解路径。输入是特征矩阵X维度d×n、标签向量Y1×n、近邻数k和热核宽度σ输出是投影矩阵W维度d×m。步骤如下计算所有样本之间的欧氏距离矩阵。根据标签构建类内图权重矩阵 ( S_w ) 和类间图权重矩阵 ( S_b )。计算度矩阵 ( D_w )、( D_b )再构造拉普拉斯矩阵 ( L_w )、( L_b )。做对称归一化得到 ( \tilde{L}_w )、( \tilde{L}_b )。构造矩阵 ( A X \tilde{L}_b X^T )( B X \tilde{L}_w X^T \epsilon I )。对 B 做Cholesky分解或直接用广义特征分解求解 ( A v \lambda B v )。按特征值从大到小排列取前m个特征向量作为投影矩阵。m的选择通常在2到min(类别数-1, d)之间。注意这里的m不是类别数减一这种固定值——图嵌入方法不受这个约束可以从1一直试到特征维度上限看准确率曲线哪里趋于平稳就选哪里这是NDGE灵活性上的优势。Matlab实现与关键代码3.1 整体代码结构与数据准备整个Matlab实现分成三大块主脚本、NDGE核心函数、分类与概率输出模块。主脚本负责加载数据、调用模型、绘制准确率曲线核心函数负责图构建和特征值分解概率输出模块用softmax把高维判别向量转成每个类别上的概率。先放主脚本clear; close all; clc; rng(42); % 加载数据 % data load(fault_dataset.mat); % X: d x n 特征矩阵每一列是一个样本 % Y: 1 x n 标签向量取值 1 到 C % 这里以随机数据做演示结构 m 44; % 特征维度 C 4; % 故障模式数量 nPerClass 80; % 每类样本数 X randn(m, nPerClass * C); Y repmat(1:C, nPerClass, 1); Y reshape(Y, 1, []); % 拆分训练与测试集 splitRatio 0.7; idx randperm(size(X, 2)); nTrain floor(size(X, 2) * splitRatio); trainIdx idx(1:nTrain); testIdx idx(nTrain1:end); X_train X(:, trainIdx); Y_train Y(trainIdx); X_test X(:, testIdx); Y_test Y(testIdx); % 标准化按训练集均值和标准差 mu mean(X_train, 2); sg std(X_train, 0, 2) 1e-8; X_train (X_train - mu) ./ sg; X_test (X_test - mu) ./ sg; % 调用 NDGE params.k 8; params.sigma auto; params.epsilon 1e-6; [W, projTrain, projTest] NDGE(X_train, Y_train, X_test, 12, params); % 训练分类器并输出概率 [probTest, Y_pred, acc] classifyWithSoftmax(projTrain, Y_train, projTest, Y_test); disp([测试准确率: , num2str(acc * 100), %]); % 绘制不同维度下的准确率 dimRange 2:2:30; accDim zeros(size(dimRange)); for i 1:length(dimRange) [Wi, projTrain_i, projTest_i] NDGE(X_train, Y_train, X_test, dimRange(i), params); [~, ~, accDim(i)] classifyWithSoftmax(projTrain_i, Y_train, projTest_i, Y_test); end figure; plot(dimRange, accDim * 100, -o, LineWidth, 1.5); xlabel(投影维度); ylabel(测试准确率 (%)); title(不同维度下的故障识别准确率); grid on;需要注意的是数据标准化要放在拆分之后并且只能用训练集的均值和标准差去处理测试集否则会产生信息泄漏准确率虚高。这个错误我在初期踩过改过来之后结果真实多了。3.2 邻接图构建和权重计算NDGE核心函数的第一步是构建两个邻接矩阵。实现时我用了两阶段法先算全连接距离矩阵再根据近邻和标签筛选边。function [W, projTrain, projTest, info] NDGE(Xtrain, Ytrain, Xtest, m, params) % NDGE 归一化判别图嵌入 % 输入 % Xtrain - d x ntr 训练特征矩阵 % Ytrain - 1 x ntr 训练标签 % Xtest - d x nte 测试特征矩阵 % m - 投影维度 % params - 结构体包含k、sigma、epsilon % 输出 % W - d x m 最终投影矩阵 % projTrain- m x ntr 训练集投影结果 % projTest - m x nte 测试集投影结果 % info - 包含特征值和图构建信息的结构体 k params.k; epsilon params.epsilon; % 合并数据用于距离计算测试集只投影不参与构图 Xall [Xtrain, Xtest]; ntr size(Xtrain, 2); % 训练集内部构图 D pdist2(Xtrain, Xtrain); [~, idxSorted] sort(D, 2); % 自动确定热核宽度 if strcmp(params.sigma, auto) sigma median(D(:)); sigma max(sigma, 1e-4); else sigma params.sigma; end % 初始化权重矩阵 n ntr; Sw zeros(n, n); % 类内图 Sb zeros(n, n); % 类间图 maskW (Ytrain Ytrain); % n x n同类为true3.3 求解投影矩阵图拉普拉斯和特征值分解的实现如下% 填权重 for i 1:n neighbors idxSorted(i, 2:k1); % 去掉自身 for j neighbors if i ~ j w exp(-D(i, j) / sigma); if maskW(i, j) Sw(i, j) w; else Sb(i, j) w; end end end end % 对称化 Sw max(Sw, Sw); Sb max(Sb, Sb); % 度矩阵和拉普拉斯矩阵 Dw diag(sum(Sw, 2)); Db diag(sum(Sb, 2)); Lw Dw - Sw; Lb Db - Sb; % 对称归一化 invSqrtDw diag(1 ./ sqrt(max(diag(Dw), 1e-12))); invSqrtDb diag(1 ./ sqrt(max(diag(Db), 1e-12))); Lw_norm invSqrtDw * Lw * invSqrtDw; Lb_norm invSqrtDb * Lb * invSqrtDb; % 构造广义特征问题 A Xtrain * Lb_norm * Xtrain; B Xtrain * Lw_norm * Xtrain epsilon * eye(size(Xtrain, 1)); % 特征值分解求最大化目标函数的前m个特征向量 opts.disp 0; [V, D] eigs(A, B, m, largestabs, opts); % 排序eigs本身一般按大小返回但保险起见再排一次 [~, order] sort(diag(D), descend); W V(:, order); % 投影 projTrain W * Xtrain; projTest W * Xtest; info.Sw Sw; info.Sb Sb; info.eigvals diag(D); info.sigma sigma; end这里最容易出问题的是 eigs 对 B 矩阵的要求。B 必须是对称正定阵。加上 epsilon 正则项之后理论上没问题但万一特征值还是出现负的小值可以把 epsilon 调大一个量级。另外如果样本数小于特征维度d远大于n时原始的 ( X L X^T ) 会是一个非常低秩的大矩阵用 eigs 求解速度很慢。此时可以先用PCA把维度压到不超过n的水平再跑NDGE这是实践中很实用的一招。3.4 输出每个样本对不同故障模式的概率故障诊断场景里光给一个“这属于第2类故障”的结论往往不够。现场工程师更希望看到“这个样本有87%的概率是内圈故障、8%概率是外圈故障”这样带置信度的输出。概率化输出的另一个好处是当预测置信度偏低时可以触发人工复核而不用完全相信硬分类。在NDGE投影后的低维特征上做概率估计我首选的是带温度系数的softmaxfunction [prob, Y_pred, acc] classifyWithSoftmax(projTrain, Ytrain, projTest, Ytest) % 用最近类别原型 softmax 输出概率 % 先对每个类别计算原型均值向量 classes unique(Ytrain); C length(classes); prototypes zeros(size(projTrain, 1), C); for c 1:C prototypes(:, c) mean(projTrain(:, Ytrain c), 2); end % 计算测试样本到每个类原型的欧氏距离 Dp pdist2(projTest, prototypes); % softmax: 距离越小概率越高 % 温度参数 tau 控制概率分布的锐度 tau 1.0; logits -Dp / tau; maxLogit max(logits, [], 2); expLogits exp(logits - maxLogit); prob expLogits ./ sum(expLogits, 2); [~, Y_pred] max(prob, [], 2); acc mean(Y_pred Ytest); end这个模块返回的prob矩阵形状是 nTest × C每一行是一个测试样本每一列是它对某个故障模式的概率。温度tau可以调节tau越小概率分布越极端tau越大概率越均匀。工程上一般让tau在0.5到2之间调参。实验过程与结果解读4.1 实验设定与数据集情况我在某公开滚动轴承故障数据集上做了验证该实验台通过电火花在轴承上加工出单点故障采集驱动端加速度信号包含正常、内圈故障、外圈故障、滚动体故障四类状态。每个状态取若干样本每段信号计算时域、频域、小波包能量共44维特征。对比算法选了PCA、LDA和未归一化的DGE判别图嵌入。所有方法统一使用相同的训练、测试划分分类器也都用最近的softmax原型分类器保证对比公平。近邻数k在NDGE和DGE里统一设置为8热核宽度σ取所有样本对距离的中位数。4.2 不同维度下的准确率变化规律直接说结果。下面表格是一组典型数值具体数值会因随机划分浮动但趋势稳定投影维度PCA准确率LDA准确率DGE准确率NDGE准确率272.5%81.3%83.6%86.9%478.1%88.2%89.3%92.4%681.5%91.0%92.5%95.2%884.2%92.4%94.1%96.6%1085.8%93.1%94.8%96.8%1586.9%93.5%95.0%97.1%2087.2%93.8%95.1%97.2%几个规律值得注意维度从2增加到8的过程中NDGE准确率提升最明显说明前几个投影方向确实抓住了主要的判别信息到10维之后准确率进入平台期再增加维度收益有限甚至偶尔会略有下降。这主要是因为低贡献的特征向量引入了噪声所以实际部署时不需要一味追求高维度选准确率曲线拐点附近即可。在低维度区间2到6维NDGE的领先幅度最大。之前用DGE时当维度比较低时投影方向会被少数高连接样本带偏归一化之后这种问题大大缓解所以NDGE在低维压缩场景下优势明显。对工程上想要把44维特征压到6到8维的需求来说这个特性很实用。4.3 概率输出的实际形态与使用场景概率输出比单纯准确率更能反映模型行为。我截取了一次实验的典型输出样本编号 真实故障 预测概率分布 1001 正常 [正常 0.94, 内圈 0.02, 外圈 0.03, 滚动体 0.01] 1002 内圈故障 [正常 0.01, 内圈 0.87, 外圈 0.06, 滚动体 0.06] 1003 外圈故障 [正常 0.02, 内圈 0.11, 外圈 0.85, 滚动体 0.02] 1015 滚动体故障 [正常 0.03, 内圈 0.22, 外圈 0.10, 滚动体 0.65]像1003号样本这样尽管判对了外圈故障但内圈故障的概率也有0.11说明该样本在投影空间中跟内圈类原型比较接近可能正好处于两类故障的边界区域。这类样本单独靠硬标签看不出来问题但概率分布一旦展开边界样本就现形了。实际使用这个概率输出时我通常会设一个置信度阈值比如0.8。低于阈值的样本一律标为“待复核”由设备维护人员结合时域波形二次确认。这种半自动模式在车间落地时比全自动更受欢迎因为工程师信任带置信度的判断而不是看着一个黑盒直接给结论。实践中踩过的坑与解决方案5.1 广义特征值分解的数值稳定性问题用Matlab的eigs求解广义特征值问题时最常见的问题有两个分母矩阵B奇异导致计算失败算出负特征值导致排序混乱。第一个问题靠正则项解决epsilon设1e-6基本够用。如果数据维度特别高可以考虑先降维到100维以内再跑NDGE。第二个问题更隐蔽通常是因为原始特征里有几个量级特别大的离群值。比如峰值因子一般是3到5之间但某个样本突然飙到25距离矩阵里对应的行就会出现极端值热核权重直接变成0图结构被破坏。解决办法是在特征标准化时用中位数和MAD绝对中位差替代均值和标准差对离群值更鲁棒。5.2 近邻数和热核宽度的参数选择k值选多大才算合适这个没有一劳永逸的答案。k太小类间图连接不充分边界样本的关系表达不完整k太大类内图混入异类样本判别信息被污染。我的经验法是k取类别样本数的小百分比比如每类80个样本时k取8到10每类30个时k取5。σ的自适应选取比k更值得用心。固定σ很容易翻车因为不同工况下信号的幅值波动很大距离矩阵的尺度完全不同。我推荐用中位数距离作为σ的基准而不是自己拍脑袋定一个数。自动模式下即使不同数据集之间特征分布差异很大算法也能自适应地调整热核宽度这一条在跨设备、跨工况迁移时特别重要。5.3 概率校准与模型过拟合softmax原型分类器虽然简单可靠但有一个通病当训练样本较少时概率分布往往过度自信明明样本落在两类边界输出概率还是接近1。这在故障诊断里是需要警惕的因为边界样本恰恰是最容易误判的样本。一个有效的校准方案是使用温度缩放。在验证集上搜索一个温度tau让softmax输出的置信度与实际准确率尽量匹配。实操代码很简单% 在验证集上搜索温度参数 tauCandidates 0.1:0.1:3.0; bestTau 1.0; bestScore inf; for t tauCandidates % 用当前tau计算验证集概率 DpVal pdist2(projVal, prototypes); logits -DpVal / t; probVal exp(logits - max(logits, [], 2)); probVal probVal ./ sum(probVal, 2); % 计算置信度与准确率的匹配误差 conf max(probVal, [], 2); score abs(mean(conf) - mean(Y_predVal Yval)); if score bestScore bestScore score; bestTau t; end end这套温度校准方法不需要重新训练分类器只是对概率输出做后处理工程性价比极高。我用完之后边界样本的置信度明显更真实现场工程师也更愿意采纳系统给出的判断。收尾个人实操心得与扩展建议NDGE这套方案我在多个故障诊断场景里跑过印象最深的是它在低维压缩下的表现。把44维特征压到8维左右识别准确率还能保持在96%以上这在现场部署时很实用因为低维特征意味着更小的存储开销、更快的分类速度对嵌入式硬件也更友好。有一个操作细节值得单独提出来说投影矩阵W在使用前一定要检查列向量的符号和排序。特征向量乘以负一仍然是特征向量但两批独立实验投影方向可能相反导致特征可视化时出现镜像翻转。虽然不影响分类准确率但分析投影方向对应的物理含义时会很困惑。我的习惯是在主脚本里加一个约定每个特征向量第一个非零元素的符号固定为正这样不同批次实验的结果可比性更强。这个项目后续的扩展方向我觉得有两条。一是把NDGE和特征选择结合起来先用过滤式方法筛掉明显无关的频带特征再进NDGE做判别压缩整个流程会更稳。二是把概率输出和贝叶斯决策融合利用故障先验概率计算误报损失和漏报损失得到更贴近工业现场风险偏好的决策边界。这些扩展都不需要改动NDGE主体新增模块就能接上。如果有同学在自己的数据上试出了不同的参数规律或者踩到了新坑欢迎交流这类方法在工程上最有价值的往往是那些不在论文里、只在实际跑数据时才能发现的细节。