当前位置：首页 > news >正文

【数据-模型融合驱动】基于指数退化模型和自适应加权损失函数的机械健康趋势预测（Python）

news 2026/4/12 21:48:17

机械设备在用久了之后，总会慢慢老化、磨损，最终出故障。如果能提前知道它啥时候会坏，就能提前安排维修，避免突然停机造成损失。但预测这事儿不好办，因为每台设备的工况不一样，磨损速度也不一样。传统的办法要么靠老师傅的经验，要么设定一个固定的检修周期，但都不够准。现在有了数据驱动的方法，可以用传感器监测设备的振动、温度这些信号，从中提取能反映健康状况的特征，比如振动的剧烈程度。然后，用一个数学模型（比如指数函数）来描述设备从开始出问题到完全失效这个过程。为了让预测更准，这个模型会重点学习最近一段时间的数据，因为越靠近现在，趋势越能代表未来的走向。这样，就能算出设备大概还能撑多久，让维修计划更有依据。

这个算法其实就是给机械设备建立一个健康衰减曲线，利用历史数据拟合一个指数模型，重点学习最近的变化趋势，然后往前推算出设备还能正常工作多久。它把数据、数学和实际经验结合起来，既考虑了设备退化的普遍规律，又能捕捉到当下的变化速度，从而给出一个动态更新的剩余寿命预测，帮助更科学的安排维修，减少意外停机。

算法步骤

收集信号，提取特征：先装上传感器，采集设备运行时的振动数据。然后把数据切成一段一段的，算每一段的振动剧烈程度（比如均方根值和最大值），得到一条能反映设备健康状态的时间序列。

找准“开始变坏”的点：从这条序列里找出设备从正常转入退化状态的时刻，也就是故障起始点，后边的数据才用来建模。

划分历史数据和测试数据：把起始点之后的数据分成两部分，前面一段用来训练模型，后面一段留着验证预测准不准。

搭个指数模型：假设退化趋势是指数上升的，先根据历史数据猜一下模型的初始参数，比如初始退化水平和上升速度。

训练模型，让近期数据更重要：用历史数据去拟合这个指数模型，但拟合的时候给最近的数据更高的权重，这样模型会更关注当前的变化趋势，而不是平均地看整个历史。训练过程就是反复调整参数，让模型的预测值和实际值之间的误差（加权的）越来越小。

监控训练，防止过头：一直训练直到误差不再明显下降，或者达到预设的训练次数，避免模型学过头了。

外推预测，算出剩余寿命：用训练好的模型向未来推演，看看按照这个趋势，设备的振动水平什么时候会达到设定的故障阈值（比如经验值或标准值），那个时间点就是预测的故障时间。

画图展示，对比实际：把原始数据、历史数据、模型拟合曲线、预测曲线、阈值线都画在一张图上，直观看出预测得怎么样。

评估准不准：拿预测的故障时间和实际发生的故障时间比较，算算误差有多大，就知道这模型靠不靠谱了。

# 导入必要的模块

import osimport timeimport scipy.ioimport scipy.statsfrom matplotlib import pyplot as pltimport numpy as npimport pandas as pdimport torchimport torch.nn as nnimport torch.optim as optimprint("PyTorch版本:", torch.__version__)# 载入数据文件PHM_path = 'PHM'PHM_bearing_files = [os.path.join(PHM_path, file) for file in os.listdir(PHM_path)]# 将mat文件转换为数组的函数def mat_to_arr(file):"""将MATLAB的.mat文件转换为Python数组参数:file: 文件路径返回:h: 原始振动信号FPT: 故障起始时间kurtosis: 峭度指标rms: 均方根值ma: 最大值"""# 加载mat文件并重塑为1维数组h = scipy.io.loadmat(file)['h'].reshape(-1)# 将信号重塑为多段，每段2560个点h2 = h.reshape(-1, 2560)# 计算每段的峭度（反映信号尖峰程度的统计量）kurtosis = np.array([scipy.stats.kurtosis(i) for i in h2])# 计算每段的均方根值（反映信号能量大小的指标）rms = np.array([np.mean(i**2)**0.5 for i in h2])# 对RMS进行平滑处理，使用[0.3,0.4,0.3]的卷积核rms = np.convolve(rms, [0.3, 0.4, 0.3], mode='same')# 计算每段信号的绝对最大值ma = np.array([np.max(np.abs(i)) for i in h2])# 计算故障起始时间（Fault Progression Time）FPT = int(len(h2)) * 1700 / 2560print("故障起始时间FPT:", FPT)return h, FPT, kurtosis, rms, ma