去年接手某个轴承退化趋势预测任务时我被一个现象卡了很久振动特征在健康阶段波动极小各种深度学习模型的预测都贴着实测走可一旦设备进入早期退化段模型的预测曲线就开始“犹豫”要么滞后几个采样周期要么把退化趋势平滑成一条近似直线等到真正发生严重磨损时RUL还停留在“基本健康”的状态。数据量并不少标签也不至于完全缺失但趋势就是追不上去。后来我把黎曼几何和物理信息神经网络的思路结合起来做成了Riemannian-PINN在同一批数据上用Pytorch重跑退化中后期的趋势才算真正稳了下来。这篇文章会拆开讲三件事为什么传统数据驱动模型在退化趋势预测上会翻车为什么要把特征放到SPD流形上用黎曼距离度量以及一个完整的、能跑的Riemannian-PINN训练管线长什么样。适合正在做设备健康管理、故障预测与健康管理、或者对PINN变体感兴趣的朋友参考。这里面的实现细节和踩坑经验我会尽量写透方便你直接迁移到自己的数据集上。1. 退化趋势预测难在哪三个看似平常却致命的假设1.1 训练集与测试集的分布漂移滚动轴承退化预测本质上是一个时间序列回归问题但和股票预测、气象预报不一样的是轴承的退化过程几乎不存在“重复样本”。同一型号轴承在相同工况下由于材料微观缺陷、装配误差、润滑状态不同退化轨迹的起始时间、加速点、寿命长度都差异很大。更麻烦的是训练集里往往是某一台设备从健康到失效的全周期数据测试集却是另一台设备从健康到失效的数据两者的特征分布天然不一致。传统LSTM、GRU这类模型默认训练数据和测试数据同分布。它们会把训练集里“第300天、RMS2.1”这种具体数值记住一旦测试设备的早期退化特征和训练设备不同模型就会拿健康阶段的模式去解释预测结果自然偏向保守或滞后。分布漂移是所有数据驱动方法在做健康管理时绕不过去的墙单纯堆模型容量并不能解决它。1.2 标签稀疏且定义模糊退化趋势预测的标签通常有两种形式一种是剩余使用寿命RUL百分比从1递减到0一种是健康指标HI从0递增到1。问题在于标签在哪个时间点开始“算退化”并没有客观标准。有人用工况阈值切分有人用3sigma准则检测突变点但不同方法切出来的退化起点可能相差几十个周期。更棘手的是大部分退化数据集只有末端附近有相对可信的标签前中期几乎都是模糊的“健康”或“轻微退化”。数据驱动模型在这种标签分布下很容易学到一种偷懒的解把所有时刻都预测成健康值因为这样在大部分标签点上误差都很小。这也是很多模型验证集指标不错、实际部署却完全失灵的原因之一。1.3 高维特征在欧氏空间里的距离失真做退化预测通常不会直接用原始振动波形而是提取时域、频域、时频域统计量比如RMS、峭度、峰值因子、谱峭度、包络谱特征等。这些特征在高维空间里构成一个可能存在非线性流形结构的点云。大多数模型默认用欧氏距离度量样本间差异但真实退化轨迹往往不是直线变化的欧氏距离会低估或高估某些阶段的真实退化程度。举个直观的例子两个时间窗的振动特征一个来自健康阶段一个来自早期退化阶段它们的欧氏距离可能很小而两个都来自退化末期的特征欧氏距离反而很大。如果模型完全依赖欧氏距离判断状态就容易把早期退化当成健康抖动把晚期状态差异当成不同故障模式。这就是为什么需要在特征表示层面引入更贴合数据几何结构的度量方式。2. 为什么选择SPD流形把退化特征放进正确的几何空间2.1 协方差矩阵为什么是SPD矩阵在滚动轴承振动监测中我们通常采集多个通道的振动信号常见的有加速度计的x、y、z三向或者同一测点的多组传感器。对每个滑动时间窗计算各通道之间的协方差矩阵会得到一个对称正定SPD矩阵。协方差矩阵描述的是通道间联合波动特性天然携带了振动能量分布和相位关系的信息比单一统计量更完整。SPD矩阵构成的集合不是普通的平坦空间而是一个黎曼流形。意思是在这个空间里两点之间最短路径不是直线而是测地线直接拿矩阵元素做减法和欧氏距离会得到和流形真实结构不一致的结果。比如协方差矩阵A和B它们的欧氏均值可能不再是正定矩阵这在物理上根本解释不通。2.2 欧氏距离在SPD流形上的失效案例我用一个实际算例说明这个问题。假设某轴承健康窗口的协方差矩阵近似为对角线矩阵数值都集中在1附近退化窗口的协方差矩阵则出现明显的能量集中某个方向的特征值增大到10。此时用欧氏距离计算两者差异会受矩阵元素绝对值的影响可能把“健康与早期退化”的距离算得比“早期退化与严重退化”还大因为前者在多个元素上都有小幅变化后者只在一个元素上大幅变化。Log-Euclidean距离是SPD流形上一种常用的黎曼度量做法是先对矩阵取矩阵对数再计算欧氏距离。这个操作的物理含义是它把矩阵间的差异从“绝对数值差异”变成“对数域的相对差异”避免了大数值元素主导距离计算的问题。实测中同一对SPD矩阵欧氏距离和Log-Euclidean距离给出的大小排序经常完全不同这直接影响到后续聚类、分类和趋势预测的效果。2.3 工程上构造SPD特征的三个细节第一协方差矩阵要保证数值正定。实际计算中由于窗口长度有限或数值精度问题协方差矩阵可能出现极小的负特征值解决办法是加上一个很小的单位阵系数比如1e-5或1e-6。第二时间窗长度的选择要兼顾稳定性和敏感性。窗口太短协方差估计噪声大窗口太长退化细节被平滑掉。经验值是取设备旋转周期的10到20倍具体需要结合采样率实验。第三矩阵对数可以作为一种输入特征变换。由于SPD流形不是线性的直接把协方差矩阵的元素拼成向量喂给神经网络网络很难学到流形结构。先把矩阵映射到对数域再取上三角元素展开成向量相当于把流形上的点投影到其切空间这样既保留了黎曼距离的结构又让特征向量进入网络时是“友好”的欧氏表示。这段代码我在后面第4章会详细给出。3. 把退化动力学写进损失函数PINN约束设计3.1 物理信息神经网络为什么能补数据不足物理信息神经网络PINN的核心思想很朴素如果某个系统遵循已知的动力学规律那就在损失函数里加一项“违背规律”的惩罚项强迫网络的输出不仅拟合有限标签还要满足微分方程约束。对于滚动轴承退化问题观测数据难免稀疏、带噪声但退化趋势总体上是有规律可循的比如早期缓慢、中后期加速不会完全随机。有了物理约束网络就相当于被限制在一个合理的解空间里不再可能把所有时刻都预测成健康值。因为退化曲线一旦保持水平物理残差会立刻放大模型必须让输出随时间呈现符合退化方程的变化。这在标签稀疏的前中期尤其重要物理约束起到了伪标签的作用。3.2 退化控制方程怎么选受限增长模型轴承退化曲线大致呈S形初期退化非常缓慢中期逐渐加速后期接近寿命终点时变化速率再次放缓受限于最大退化程度。这种形态用受限增长模型描述最自然dy/dt r * y * (1 - y / K)其中y是归一化退化指标范围在0到1r是退化速率系数K是退化上限通常设为1。这个方程的直观含义是退化指标越小增长越慢越接近上限增长也被抑制中间阶段增长最快。相比纯指数模型受限增长模型不会出现退化指标无限膨胀的离谱预测更适合做寿命末期校准。实际使用时r需要根据设备历史数据大致估计。有一个快速标定的办法取退化中段两个标签点用差分近似dy/dt再除以对应时刻的y*(1-y)就能得到初始r值。代码里我会把它设成一个可训练的小参数或者直接作为常数传入。3.3 自动微分求物理残差PINN实现的关键一步是用自动微分计算网络输出对时间的导数。Pytorch的torch.autograd.grad可以方便地完成这一点。假设网络输入包含连续时间t输出退化指标y那么对t求梯度就能得到dy/dt然后和方程右侧的ry(1-y/K)计算MSE就是物理残差损失。这里有个细节t在输入网络前要做归一化处理让数值落在0到1之间否则梯度会被时间量纲拉偏。同时求导时要设置create_graphTrue因为物理损失本身还要参与反向传播需要保留计算图。这个点新手很容易漏漏了之后梯度链断裂loss能算出来但模型参数完全不更新。3.4 黎曼几何正则让预测速率对齐特征流形变化只加物理约束还不够因为物理方程里的r是全局常数无法刻画不同设备、不同退化阶段之间的个体差异。这时候就需要几何信息来“校准”退化速率。我的做法是对相邻时间窗的SPD特征矩阵计算它们的Log-Euclidean距离作为该时段特征流形上的实际变化幅度同时计算网络预测的退化指标在该时段的增量让两者成比例对齐。这个约束的直观解释是特征空间里走得越快退化指标就应该涨得越快。它把“物理规律”和“特征几何”两套信息绑在同一个约束里既避免了纯物理方程过于僵化又弥补了纯数据模型缺乏机制约束的问题。实际数据里只要SPD特征构造合理这个几何正则能明显改善早期退化阶段的趋势捕捉能力。4. Pytorch实现Riemannian-PINN可运行的训练管线4.1 数据管线滑动窗口、SPD特征提取与标签生成先把数据结构设计清楚。假设原始数据是多通道振动信号形状为(num_samples, n_channels)采样率固定。先用滑动窗口切段每个窗口计算一个SPD协方差矩阵再通过矩阵对数映射成特征向量。标签按照设备全寿命周期归一化把健康起点标为0失效终点标为1中间线性插值。下面是SPD特征提取的核心代码import torch import torch.nn as nn import numpy as np def sym_logm(X, eps1e-6): # X: (batch, n, n) SPD矩阵 evals, evecs torch.linalg.eigh(X) evals torch.clamp(evals, mineps) log_evals torch.log(evals) # 利用 X V diag(λ) V^TlogX V diag(logλ) V^T return (evecs * log_evals.unsqueeze(-2)) evecs.transpose(-2, -1) def extract_spd_features(windows): # windows: (num_windows, n_channels, window_len) cov_mats [] for X in windows: Xc X - X.mean(dim-1, keepdimTrue) C (Xc Xc.transpose(-1, -2)) / (X.shape[-1] - 1) C C 1e-5 * torch.eye(X.shape[-2]) cov_mats.append(C) cov_mats torch.stack(cov_mats) # (N, n, n) log_mats sym_logm(cov_mats) # (N, n, n) n log_mats.shape[-1] idx torch.triu_indices(n, n) feat_vec log_mats[:, idx[0], idx[1]] # 取上三角展平 return cov_mats, feat_vecsym_logm里用特征分解实现矩阵对数比直接调用SciPy稳定得多而且支持反向传播。取上三角是为了避免重复输入冗余的下三角元素同时保留了SPD矩阵对数域的完整信息。如果你的数据是多通道振动这一步基本可以直接复用。4.2 模型定义时间与Log-Euclidean特征的双输入结构网络结构采用简单的多层感知机输入拼接连续时间t和Log-Euclidean特征向量z输出归一化退化指标y。为什么不让网络直接吃cov_mats因为SPD流形的非线性几何结构会让普通全连接层学得很吃力而经过矩阵对数投影后特征向量已经落在线性切空间里网络只需要拟合一个相对平滑的映射。class RiemPINN(nn.Module): def __init__(self, feat_dim, hidden128, num_layers4): super().__init__() dims [1 feat_dim] [hidden] * num_layers [1] layers [] for i in range(len(dims) - 1): layers.append(nn.Linear(dims[i], dims[i 1])) if i len(dims) - 2: layers.append(nn.Tanh()) else: layers.append(nn.Sigmoid()) # 退化指标限制在[0,1] self.net nn.Sequential(*layers) def forward(self, t, z): # t: (batch, 1)z: (batch, feat_dim) x torch.cat([t, z], dim-1) return self.net(x)激活函数选择Tanh而不是ReLU是因为退化曲线是连续平滑的Tanh的平滑性更适合做自动微分求导。最后一层用Sigmoid把输出限制在0到1和归一化标签对齐。隐藏层宽度128、层数4是经验值数据量小的场景可以缩到64和3数据量大的话适度加宽。4.3 组合损失函数数据拟合、物理残差、黎曼正则这是整个模型的核心。损失包含三部分有监督数据点的MSE、退化方程残差、SPD流形变化速率对齐。物理残差依赖自动微分黎曼正则依赖Log-Euclidean距离。三个损失各自解决一个问题最后的组合权重需要细心调度。def log_euclidean_dist(A, B): # A, B: (batch, n, n) SPD矩阵 diff sym_logm(A) - sym_logm(B) return torch.norm(diff, dim(-2, -1)) # (batch,) def pinn_loss(model, batch_data, batch_phys, batch_geom, r0.02, K1.0): # batch_data: (t_data, z_data, y_data) t_data, z_data, y_data batch_data y_pred model(t_data, z_data) loss_data nn.functional.mse_loss(y_pred, y_data) # batch_phys: (t_phys, z_phys) t_phys, z_phys batch_phys t_phys.requires_grad_(True) y_phys model(t_phys, z_phys) dy_dt torch.autograd.grad( y_phys.sum(), t_phys, grad_outputstorch.ones_like(y_phys), create_graphTrue )[0] f_phys r * y_phys * (1 - y_phys / K) loss_phys nn.functional.mse_loss(dy_dt, f_phys) # batch_geom: (A_spd, B_spd, z_a, z_b, t_a, t_b, scale) A_spd, B_spd, z_a, z_b, t_a, t_b, scale batch_geom d_riemann log_euclidean_dist(A_spd, B_spd) y_a model(t_a, z_a) y_b model(t_b, z_b) delta_y (y_b - y_a).squeeze(-1) target_delta scale * d_riemann loss_geom nn.functional.mse_loss(delta_y, target_delta) return loss_data, loss_phys, loss_geomscale是退化指标增量与SPD特征距离之间的比例常数。实践中可以用带标签的样本做一个最小二乘预标定先算历史相邻窗口的特征距离再算对应的标签增量拟合一个线性系数。代码注释里建议先跑一遍带标签段标定好scale再进入正式训练。关于物理损失我特别提醒一点r的取值会影响退化曲线的整体陡峭度。如果r设得过大物理残差会强行把曲线拉成陡峭的S形即使数据标签不支持如果r过小退化趋势又会被过度平滑。稳妥的做法是把r定在0.01到0.05之间观察物理残差的数量级再做微调。4.4 训练主循环与损失权重调度训练主循环要注意两个问题一是类别不平衡三个损失的量纲不同需要分别用权重配平二是早期训练阶段网络还没把数据拟合好物理和几何约束加入太早会干扰收敛。我的做法是数据损失权重固定为1物理损失权重随训练进程从0.01线性增加到0.5几何损失权重在第50个epoch之后才启用。def train_pinn(model, dl_data, dl_phys, dl_geom, epochs300): opt torch.optim.Adam(model.parameters(), lr1e-3) sched torch.optim.lr_scheduler.CosineAnnealingLR(opt, T_maxepochs) for epoch in range(epochs): w_phys 0.01 0.49 * min(epoch / 200, 1.0) w_geom 0.1 if epoch 50 else 0.0 for batch_data, batch_phys, batch_geom in zip(dl_data, dl_phys, dl_geom): loss_data, loss_phys, loss_geom pinn_loss( model, batch_data, batch_phys, batch_geom ) loss loss_data w_phys * loss_phys w_geom * loss_geom opt.zero_grad() loss.backward() nn.utils.clip_grad_norm_(model.parameters(), max_norm1.0) opt.step() sched.step() if epoch % 50 0: print(fepoch {epoch}: fdata{loss_data.item():.4f}, fphys{loss_phys.item():.4f}, fgeom{loss_geom.item():.4f})梯度裁剪非常重要。物理残差对时间求导之后梯度数值放大了不少如果不裁剪训练很容易在初期就发散。Clip到1.0是我在这个任务上试过比较稳的阈值如果你的数据尺度不同可以从0.5到2.0之间调。学习率用余弦退火可以让后期的损失更低避免在平台期来回震荡。评估时不要只看整体RMSE要分阶段看健康段、早期退化段、中后期退化段分别计算误差。退化趋势预测最容易出问题的是中后期预测偏保守如果分段误差显示中后期MAE明显偏大大概率是物理权重不够或者r参数偏小对应调整即可。5. 训练过程踩过的坑权重失衡、梯度爆炸与平台期5.1 物理残差权重过大导致拟合完全失效我第一次跑完整管线时把物理损失权重设成和数据损失相同结果训练到100个epoch时预测曲线虽然形态很“标准”的S形但数值整体偏低末期退化指标只有0.6左右和真实标签0.9完全对不上。原因是物理方程是一个全局约束它要求所有时刻的y都满足dy/dt ry(1-y/K)如果这个约束权重过大模型就会牺牲局部的数据拟合来换取整体曲线形态的“正确”。解决方式是物理权重从小到大退火。前20个epoch让数据损失主导网络先大致拟合标签然后逐渐增加物理约束让曲线形态从“噪声轨迹”慢慢被拉成平滑的S形。这条经验几乎可以迁移到所有PINN类任务物理约束是指导性的不是主导性的它应该帮助网络在数据稀疏处保持合理外推而不是取代数据拟合。5.2 矩阵对数特征值过小时数值不稳定SPD矩阵的特征值如果非常接近0sym_logm虽然做了clamp但梯度在log(λ)的导数位置会变得异常大导致反向传播不稳定。我遇到过训练十几个epoch后loss突然变成NaN的情况排查到最终原因是某个时间窗的振动信号接近静默协方差矩阵的所有元素都趋向0特征值数量级掉到1e-8log之后数值直接发散。两个补救手段一是把clamp的下限从1e-6提高到1e-4虽然损失一点精度但换来稳定二是在特征提取前对原始振动信号做标准化减去均值后除以标准差让协方差矩阵的元素量级更稳定。如果你用的是多传感器采集的原始波形这条尤其重要不同测点的量纲差异会让SPD特征矩阵奇异性明显增加。5.3 退化起点不一致导致几何正则失控黎曼几何正则的设计前提是“特征流形距离越大退化指标增量越大”。但不同设备进入退化的起点不一致有的设备在第100天就开始退化有的第300天才开始。如果几何约束直接套在所有时间窗对上健康阶段特征的自然漂移会被几何正则误判成退化增量预测曲线会出现早期虚高的虚假退化。我的处理方法是给几何正则加一个门控只有当SPD特征距离超过某个健康参考阈值时才激活几何对齐损失。这个阈值可以用全寿命数据里前10%窗口的Log-Euclidean距离分布来标定取95分位数作为“健康漂移上限”。低于阈值的样本对几何损失直接置0。加了门控之后早期退化段的虚高预测明显减少。5.4 训练后期进入平台期不是过拟合是约束打架训练到250个epoch后三个损失各自看起来都平稳了但预测曲线始终和真实标签隔着一层薄膜末期趋势有就是不够陡。后来我把三个损失的梯度方向单独拿出来分析发现数据损失在推着曲线上调物理损失却因为r参数固定而限制曲线增速几何损失又试图保持特征距离对齐三个梯度方向在后期相互抵消形成了一个隐性死锁。解法是分阶段解耦前150个epoch把r当成可训练参数让它在反向传播中学习150个epoch后冻结r只更新网络权重。实测下来末期曲线的陡峭度明显改善验证集分段误差比之前降低了大概三成。这种“先学参数、再固定约束”的思路在处理多约束联合训练时非常实用。5.5 关于评价指标的提醒RMSE不是一切最后说一个评估层面的经验。退化趋势预测的最终目标是给维护决策提供依据所以真正重要的不是所有点的平均误差而是退化加速点的定位准不准、末期RUL偏不偏。我曾经用RMSE作为唯一指标调参调到后期RMSE很漂亮但退化加速点比实际晚了将近50个周期这类模型在真实场景里根本不敢用。建议在验证阶段同时看三个指标整体RMSE、退化加速点定位误差、末期30%阶段的平均绝对误差。如果后两个指标不达标即使整体RMSE再好看也要回头调物理权重或几何门控。这套Riemannian-PINN框架的价值不在于它一定比LSTM的RMSE低多少而在于它用物理规律和特征几何把预测约束在“像真实退化”的轨迹上减少那些看似合理实则危险的误判。做工业健康管理的朋友应该都懂能让模型在关键时刻不犯错比让它平时更精确重要得多。
