高端旋转机械剩余使用寿命预测及其不确定性量化评估方法

Prediction of RUL and Uncertainty Quantification Evaluation Methods for High-end Rotating Machinery

摘要

针对高端旋转机械剩余使用寿命预测中的准确性和不确定性量化问题，提出了基于变分深度高斯过程（VDGP）的预测方法。通过构建深度高斯过程更新模型实现不确定性的递推量化，并采用诱导点和变分推断提高大数据的处理能力。C-MAPSS和风机行星齿轮数据集的实验表明，VDGP比高斯过程方法具有更高的预测准确度和更窄的置信区间，在C-MAPSS的 FD002数据集上，均方根误差、评分函数分别比现有最佳的对比方法减小0.21%和45.3%。

Abstract

To address the challenges of accuracy and uncertainty quantification in RUL prediction of high-end rotating machinery， a prediction method was proposed based on VDGP. The method achieved recursive uncertainty quantification by constructing a deep Gaussian process update model， and enhanced large-scale data processing capability through the use of inducing points and variational inference. Experiments on the C-MAPSS and wind turbine planetary gearbox datasets demonstrate that VDGP achieves higher prediction accuracy and narrower confidence intervals compared to the standard Gaussian process methods. On the C-MAPSS FD002 dataset， the root mean square error and scoring function are reduced by 0.21% and 45.3% respectively， relative to the best baseline method.

Graphical abstract

关键词

剩余使用寿命预测 / 变分深度高斯过程 / 不确定性量化 / 旋转机械 / 核函数

Key words

remaining useful life（RUL） prediction / variational deep Gaussian process（VDGP） / uncertainty quantification / rotating machinery / kernel function

0 引言

如今，旋转机械设备的故障预测与健康管理受到越来越多的关注^［1-2］。民用航空、风力发电等领域的高端旋转机械的可靠性和效率尤为重要。通过准确的剩余使用寿命（remaining useful life， RUL）预测实施预防性维护，可大幅降低安全事故的发生率，防止不必要的人员伤亡和财产损失^［3-4］。

随着信息物理系统的广泛应用，基于深度学习的RUL预测方法受到广泛关注。深度学习方法能有效提取状态监测数据中的设备健康特征并实现RUL预测。数据增强与关系建模方面，CHEN等^［5］利用带有分布外数据增强的收缩Mamba关系网络增强了多故障模式和复杂工况下的检测与定位能力；ZHU等^［6］提出基于多尺度长短期记忆网络（long short-term memory，LSTM）和多头自注意力机制的自适应分阶段预测模型，通过融合生产与操作数据动态更新权重，提高预测精度。在序列建模与注意力机制方面，CAO等^［7］构建了TCN-Transformer框架，并结合卡尔曼滤波降噪和非线性平滑，对潜在退化状态进行估计与趋势提取，以增强特征鲁棒性；PENG等^［8］基于健康状态与轨迹相似度，通过K-means聚类分类健康状态，并采用加权匹配策略提高了预测精度。在卷积网络与多模型融合方面，CHEN等^［9］采用多尺度卷积神经网络进行特征预训练，并将其输出输入到LSTM中，同时提取全局特征与局部特征；夏然等^［10］提出一种结合健康因子筛选和遗传算法优化的广义回归神经网络与BiLSTM-NAR的混合模型。在优化与迁移学习方面，WANG等^［11］通过仿真域迁移和贝叶斯优化减少对绝对RUL的依赖，提升预测性能；LUO等^［12］采用改进鲸鱼优化算法处理多故障特征；曹现刚等^［13］提出一种基于核主成分分析-长短期记忆网络的方法来提高旋转机械RUL预测的准确度。在自注意力与卷积改进方面，自注意力ConvLSTM模型^［14］用卷积算子替代LSTM内部连接以减少模型冗余并预测风机齿轮寿命；刘华开等^［15］将卷积先验、注意力机制和BiLSTM相结合，有效提取故障特征、抑制噪声并降低预测误差、提高预测精度。文献［16-17］ 在卷积神经网络中引入灵敏度权重分配、变分推断与多尺度结构，对模型进行增强，从而提高RUL预测性能。

尽管这些研究在特征提取、降噪和退化趋势建模等方面取得了显著进展，提高了预测准确度和模型的泛化能力，但在预测结果的不确定性量化方面仍存在不足。大型旋转机电装备结构复杂、工况多变，对预测精度要求高，还需要模型能提供可靠的不确定性评估。风机齿轮箱中的行星齿轮寿命预测存在较强不确定性：①人工特征提取依赖专家经验、具有主观性；②转速和负载变化引起啮合及故障特征频率波动，可能产生伪故障信号；③复杂传动结构导致冲击信号与其他部件的振动分量耦合。这些因素均加大了RUL预测的难度。

针对上述问题，本文构建了变分深度高斯过程（variational deep Gaussian process， VDGP）模型，用于实现高端旋转机械及其关键部件的RUL预测，并提供预测结果的不确定性量化，提高寿命预测的准确度。该方法通过构建不确定量分布函数的深度高斯过程（DGP）来更新模型，随着新监测信息的不断获取，持续更新并修正不确定量的分布函数，逐步精确量化预测结果的不确定性。为解决高斯过程（GP）在大规模监测数据下计算成本高的问题，在DGP中引入变分推断，并采用诱导点技术降低计算复杂度，提高RUL的预测精度。本文还对比了不同核函数和模型层数下的性能表现，为所提方法在工业场景下的应用提供参考。

1 高斯过程

高斯过程（GP）是一种非参数贝叶斯模型，用于对潜在函数f（ x ）的先验分布建模，在给定输入集合时诱导出输出的联合高斯分布。输入 x 与输出y的关系为

y=f（ x ）+ε（1）

式中：f（ x ）表示噪声去除后的真实映射；ε是均值为0、方差为{\sigma }_{\epsilon }^{\mathrm{2}}的高斯噪声。

\mathrm{G}\mathrm{P}的定义为f(\mathit{x})~GP(m(\mathit{x}),k(\mathit{x},{\mathit{x}}^{\text{'}}))，其中，m（ x ）为均值函数，k（ x， x ′）为核函数， x 、 x ′为输入变量。使用较为广泛的RBF内核

\begin{array}{c}k\left(\mathit{x},{\mathit{x}}^{\text{'}}\right)={\sigma }_{\mathrm{f}}^{\mathrm{2}}\mathrm{e}\mathrm{x}\mathrm{p}(-\frac{\mathit{x}-({\mathit{x}}^{\text{'}}{)}^{\mathrm{2}}}{\mathrm{2}{l}^{\mathrm{2}}})\end{array}（2）

式中：{\sigma }_{\mathrm{f}}^{\mathrm{2}}为RBF核函数的信号方差；l为长度尺度参数。

图1展示了GP回归的基本工作流程。GP回归的目标是利用已有的观测数据来学习输入与输出之间的潜在函数关系，并据此对新的输入样本进行预测。

给定训练数据集D_t =｛X_t，y_t ｝，下标t表示第t个训练样本，在此基础上，根据GP的定义，训练集的观测输出y_t 与待预测输入 X_*对应的函数值 f ^* 共同服从联合多元高斯分布，即

\begin{array}{c}\left[\begin{array}{c}{y}_t\\ {\mathit{f}}^{\mathbf{*}}\end{array}\right]~N(\mathrm{0},\left[\begin{array}{cc}{\mathit{\Sigma }}_t& k({\mathit{X}}_t,{\mathit{X}}_{\mathrm{*}})\\ k({\mathit{X}}_{\mathrm{*}},{\mathit{X}}_t)& k({\mathit{X}}_{\mathrm{*}},{\mathit{X}}_{\mathrm{*}})\end{array}\right])\end{array}（3）

其中， Σ_t 为由训练输入计算得到的协方差矩阵，{\mathit{\Sigma }}_t=k({\mathit{X}}_t,{\mathit{X}}_t)+{\sigma }_{\epsilon }^{\mathrm{2}}\mathit{I}；k（ X_*， X_*）为新点之间的协方差矩阵；k（ X_*， X_t ）衡量新输入点与训练输入点之间的相关性（相似性），反映它们在函数空间中的关系。由联合多元高斯分布的性质可知，条件分布p({\mathit{f}}^{\mathbf{*}}|{\mathit{X}}_t,{\mathit{y}}_t,{\mathit{X}}_{\mathrm{*}})也是一个多元高斯分布，该条件分布的均值为k({\mathit{X}}_{\mathrm{*}},{\mathit{X}}_t){\mathit{\Sigma }}_t^{-\mathrm{1}}{\mathit{y}}_t，预测点{\mathit{X}}_{\mathrm{*}}的不确定性协方差矩阵为k({\mathit{X}}_{\mathrm{*}},{\mathit{X}}_{\mathrm{*}})-k({\mathit{X}}_{\mathrm{*}},{\mathit{X}}_t)·

{\mathit{\Sigma }}_t^{-\mathrm{1}}k({\mathit{X}}_t,{\mathit{X}}_{\mathrm{*}})。

因此，GP的后验仍然是一个GP，则后验分布均值函数和协方差函数分别为

\begin{array}{c}{m}_t(\mathit{x})=k(\mathit{x},{\mathit{X}}_t){\mathit{\Sigma }}_t^{-\mathrm{1}}{\mathit{y}}_t\end{array}（4）

\begin{array}{c}{k}_t(\mathit{x},{\mathit{x}}^{\text{'}})=k(\mathit{x},{\mathit{x}}^{\text{'}})-k(\mathit{x},{\mathit{X}}_t){\mathit{\Sigma }}_t^{-\mathrm{1}}k({\mathit{X}}_t,{\mathit{x}}^{\text{'}})\end{array}（5）

其中，m_t （ x ）为观测数据条件下 x 的预测均值；k_t （ x， x '）为 x 的预测方差即不确定性大小^［18］。

2 基于VDGP的PUL预测及不确定性量化评估

构建DGP模型，通过堆叠多个GP层来捕捉更复杂的非线性特征。在单层GP的基础上，将GP推广为深度结构，DGP将每层GP的输出作为下一层的输入，原理如图2所示。第l（l=1，2，…，L）层的输入和输出的关系为{\mathit{X}}^{(l)}=\mathrm{G}{\mathrm{P}}_l(\mathrm{0},k_l({\mathit{X}}^{(l-\mathrm{1})},{\mathit{X}}^{(l-\mathrm{1})}))+{\mathit{\epsilon }}^{(l)}，其中，{\mathit{X}}^{\left(\mathrm{0}\right)}为原始输入数据， ε^（l）是第l层的高斯噪声。

DGP模型的边缘似然为

p(\mathit{Y}|{\mathit{X}}^{\left(\mathrm{0}\right)})={\int }_{\mathit{\xi }}[p(\mathit{Y}|{\mathit{F}}^{\left(L\right)})p({\mathit{F}}^{\left(L\right)}|{\mathit{X}}^{\left(L-\mathrm{1}\right)})\cdot

{\displaystyle \prod _{l-\mathrm{1}}^{L-\mathrm{1}}}p({\mathit{X}}^{\left(l\right)}|{\mathit{F}}^{\left(l\right)})p({\mathit{F}}^{\left(l\right)}|{\mathit{X}}^{\left(l-\mathrm{1}\right)})]\mathrm{d}\mathit{\xi }（6）

训练DGP模型等价于最大化输出 Y 的边缘似然。为简化符号表示，设 F^（l）为对应 X^（l）（l=1，2，…，L-1）、 Y 的无噪声变量即 X^（l）= F^（l）+ ε^（l）（l=1，2，…，L-1）和 Y = F^（L）+ ε^（L）。由于DGP每层的输出仅依赖其上一层，因此具有条件独立性性质。基于条件独立性属性，隐藏层 Y 的条件概率密度p({\mathit{X}}^{(l)},{\mathit{F}}^{(l)}|{\mathit{X}}^{(l)},{\mathit{X}}^{(l-\mathrm{1})},\cdots ,{\mathit{X}}^{(\mathrm{0})},{\mathit{F}}^{(l)},{\mathit{F}}^{(l-\mathrm{1})},\cdots ,

F^{(\mathrm{1})})可简化为p({\mathit{X}}^{\left(l\right)},{\mathit{F}}^{\left(l\right)}|{\mathit{X}}^{\left(l-\mathrm{1}\right)})，并进一步分解为p({\mathit{X}}^{\left(l\right)}|{\mathit{F}}^{\left(l\right)})p({\mathit{F}}^{\left(l\right)}|{\mathit{X}}^{\left(l-\mathrm{1}\right)})。令\xi =\{{\mathit{F}}^{\left(\mathrm{1}\right)},{\mathit{F}}^{\left(\mathrm{2}\right)},\cdots ,

{\mathit{F}}^{\left(L\right)},{\mathit{X}}^{\left(\mathrm{1}\right)},{\mathit{X}}^{\left(\mathrm{2}\right)},\cdots ,{\mathit{X}}^{\left(L-\mathrm{1}\right)}\}表示需要边缘化的中间变量集合。因此，模型似然可表示为向量积分，如式（6）所示。式（6）的积分在解析上没有封闭解。积分的主要困难在于：协方差矩阵依赖随机隐变量，且 X^（l）以非线性形式出现在条件分布p({\mathit{F}}^{\left(l\right)}|{\mathit{X}}^{\left(l-\mathrm{1}\right)})的协方差矩阵的逆中，使得解析求解不可行。

2.1 变分推断方法

传统的高斯模型计算成本高。对于旋转机械长时间运行监测的大规模数据集，矩阵求逆会导致计算复杂度迅速增加。为降低计算复杂度并使后续推断可行，引入诱导点近似来逼近原始GP。选取较少的诱导点能显著降低计算量，从而提高大规模数据集下的计算效率。在该近似框架下，每层 GP 可由诱导点 Z^（l）及其诱导变量 u^（l）表示：

{\mathit{f}}^{\left(l\right)}({\mathit{X}}^{\left(l-\mathrm{1}\right)})\approx \mathit{K}({\mathit{X}}^{\left(l-\mathrm{1}\right)},{\mathit{z}}^{\left(l\right)}){\mathit{K}}^{-\mathrm{1}}({\mathit{z}}^{\left(l\right)},{\mathit{z}}^{\left(l\right)}){\mathit{u}}^{\left(l\right)}（7）

式中：\mathit{K}({\mathit{X}}^{\left(l-\mathrm{1}\right)},{\mathit{z}}^{\left(l\right)})为输出{\mathit{X}}^{\left(l-\mathrm{1}\right)}与诱导点集合 Z^（l）的协方差矩阵；\mathit{K}({\mathit{z}}^{\left(l\right)},{\mathit{z}}^{\left(l\right)})为诱导点之间的协方差矩阵。

基于诱导点近似的GP采用变分推断来近似真实后验分布，并通过最大化证据下界（evidence lower bound， ELBO）来优化变分参数和模型超参数。记 f_t 是与 y_t 对应的无噪声函数值，q({\mathit{f}}_t,\mathit{u})=p({\mathit{f}}_t|\mathit{u})\varphi (\mathit{u})为变分分布，其中，\varphi (\mathit{u})为对诱导变量\mathit{u}的变分后验分布。为便于推导，符号中省略了输入 X_t 与诱导点 Z 的显式表示（隐含在协方差矩阵的计算中），则对数边际似然为

\mathrm{l}\mathrm{n}p({\mathit{y}}_t)=\mathrm{l}\mathrm{n}{\int }_{\mathit{u}}^{}{\int }_{{\mathit{f}}_t}^{}q({\mathit{f}}_t,\mathit{u})\frac{p({\mathit{y}}_t|{\mathit{f}}_t)}{q({\mathit{f}}_t,\mathit{u})}p({\mathit{f}}_t|\mathit{u})p(\mathit{u})\mathrm{d}{\mathit{f}}_t\mathrm{d}\mathit{u}=

\mathrm{l}\mathrm{n}{\int }_{\mathit{u}}^{}{\int }_{{\mathit{f}}_t}^{}q({\mathit{f}}_t,\mathit{u})\frac{p({\mathit{y}}_t|{\mathit{f}}_t)p(\mathit{u})}{\varphi (\mathit{u})}\mathrm{d}{\mathit{f}}_t\mathrm{d}\mathit{u}（8）

通过应用詹森不等式得到变分下界：

\mathrm{l}\mathrm{n}p({\mathit{y}}_t)\ge F_v(\mathit{Z},\varphi )={\int }_{\mathit{u}}^{}{\int }_{{\mathit{f}}_t}^{}q({\mathit{f}}_t,\mathit{u})\frac{p({\mathit{y}}_t|{\mathit{f}}_t)p(\mathit{u})}{\varphi (\mathit{u})}\mathrm{d}{\mathit{f}}_t\mathrm{d}\mathit{u}=

{\int }_{\mathit{u}}^{}\varphi (\mathit{u})\mathrm{l}\mathrm{n}\frac{\mathrm{N}({\mathit{y}}_t|{\mathit{a}}_t,{\sigma }_{\epsilon }^{\mathrm{2}}{\mathit{I}}_N)p(\mathit{u})}{\varphi (\mathit{u})}\mathrm{d}\mathit{u}（9）

{\mathit{a}}_t={\mathit{K}}_{NM}{\mathit{K}}_{NM}^{-\mathrm{1}}\mathit{u}\tilde{\mathit{K}}={\mathit{K}}_{NN}-{\mathit{K}}_{NM}{\mathit{K}}_{MM}^{-\mathrm{1}}{\mathit{K}}_{MN}

式中：N（*）为高斯分布的概率密度函数； K_NN 为p（ f_d | X ）的协方差矩阵； K_NM 为 f_d 和 u_d 之间的协方差矩阵，下标d表示多输出（多维）高斯过程中的第d个输出维度； K_MM 为诱导点集合Z上的协方差矩阵；\varphi为变分分布对应的变分参数集合。

此处使用GP诱导变量对数似然积分定理，即在稀疏高斯过程（SGP）框架下，给定一个GP映射 Y ～GP（ X ），该GP的诱导变量为 U，诱导点为 Z，对于特定的输出维度d，有

\int p({\mathit{f}}_d|{\mathit{u}}_d,\mathit{X})\mathrm{l}\mathrm{n}p({\mathit{y}}_d|{\mathit{f}}_d)\mathrm{d}{\mathit{f}}_d=\mathrm{l}\mathrm{n}N({\mathit{y}}_d|{\mathit{a}}_d,{\sigma }_{\epsilon }^{\mathrm{2}}{\mathit{I}}_N)-\frac{\mathrm{t}\mathrm{r}\tilde{\mathit{K}}}{\mathrm{2}{\sigma }_{\epsilon }^{\mathrm{2}}}（10）

{\mathit{a}}_d={\mathit{K}}_{NM}{\mathit{K}}_{MM}^{-\mathrm{1}}{\mathit{u}}_d

通过对变分分布\varphi \left(\mathit{u}\right)进行最优设置，可使变分下界F_v达到最大值。结果表明，最优分布实际上是高斯分布。最大化后的下界为

F_v(\mathit{Z})=\mathrm{l}\mathrm{n}N({\mathit{y}}_t|\mathrm{0},{\sigma }_{\epsilon }^{\mathrm{2}}{\mathit{I}}_N+{\mathit{Q}}_{NN})-\frac{\mathrm{t}\mathrm{r}\tilde{\mathit{K}}}{\mathrm{2}{\sigma }_{\epsilon }^{\mathrm{2}}}（11）

{\mathit{Q}}_{NN}={\mathit{K}}_{NM}{\mathit{K}}_{MM}^{-\mathrm{1}}{\mathit{K}}_{MN}

基于前文推导得到的证据下界公式（式（11）），可通过最大化该下界来训练SGP模型。训练完成后，通过{\int }_{\mathit{u}}^{}p({\mathit{f}}^{\mathrm{*}}|\mathit{u})q(\mathit{u})\mathrm{d}\mathit{u}获得{\mathit{f}}^{\mathrm{*}}的预测分布，其中，诱导变量 u 可视为对 f_t 的充分近似，从而在保证预测精度的同时有效降低计算复杂度。

采用诱导点和变分推理获得可处理的下界。图3所示为具有诱导变量和2个隐藏层的有监督VDGP结构，其中，横向并列的节点属于同一层，纵向相邻的层通过条件依赖相连，下标i表示数据集中的第i个样本。

引入诱导变量{\mathit{U}}^{\left(l\right)}\in {\mathbf{R}}^{M_l\times D_l}和诱导点{\mathit{Z}}^{\left(l\right)}\in {\mathbf{R}}^{M_l\times D_{l-\mathrm{1}}}来扩展概率空间。给定输入{\mathit{X}}^{\left(\mathrm{0}\right)}和诱导点{\mathit{Z}}^{\left(l\right)}后，所有变量的联合概率密度为

\begin{array}{c}P=p(\mathit{Y},\{{\mathit{F}}^{(l)},{\mathit{U}}^{(l)}{\}}_{l=\mathrm{1}}^L,\{{\mathit{X}}^{(l)}{\}}_{l=\mathrm{1}}^{L-\mathrm{1}}|{\mathit{X}}^{(\mathrm{0})},\{{\mathit{Z}}^{(l)}{\}}_{l=\mathrm{1}}^L)=\\ {\displaystyle \prod _{l=\mathrm{1}}^{L-\mathrm{1}}}p({\mathit{X}}^{(l)}|{\mathit{F}}^{(l)})p({\mathit{F}}^{(l)},{\mathit{U}}^{(l)}|{\mathit{X}}^{(l-\mathrm{1})},{\mathit{Z}}^{(l)})\cdot \end{array}

p(\mathit{Y}|{\mathit{F}}^{(L)})p({\mathit{F}}^{(L)},{\mathit{U}}^{(L)}|{\mathit{X}}^{(L-\mathrm{1})},{\mathit{Z}}^{(L)})（12）

其中， U^（l）为第l层高斯过程的诱导变量，各个诱导点对应的函数取值服从与该层潜在函数 F^（l）相同的高斯过程先验分布； Z^（l）为诱导变量对应的输入位置，被视为可优化的变分参数。为简化符号，在不引起歧义时省略 Z 的显式书写（ Z 仍参与核矩阵的计算），则联合概率密度为p（ U^（l））与p（ F^（l）| U^（l）， X^（l-1））的乘积。由于核函数使协方差矩阵的元素成为 X^（l-1）的非线性函数，故协方差矩阵及其逆出现在对全部隐变量进行边缘化的积分内部，从而令解析边缘化不可得。为获得可计算的近似，使用变分推断的方法解决计算问题。其中，变分分布q(\{{\mathit{F}}^{(l)},{\mathit{U}}^{(l)}{\}}_{l=\mathrm{1}}^L,\{{\mathit{X}}^{(l)}{\}}_{l=\mathrm{1}}^{L-\mathrm{1}})用于近似真实的后验分布p(\{{\mathit{F}}^{(l)},{\mathit{U}}^{(l)}{\}}_{l=\mathrm{1}}^L,\{{\mathit{X}}^{(l)}{\}}_{l=\mathrm{1}}^{L-\mathrm{1}}|\mathit{Y})，即有

\begin{array}{l}Q=q(\{{\mathit{F}}^{(l)},{\mathit{U}}^{(l)}{\}}_{l=\mathrm{1}}^L,\{{\mathit{X}}^{(l)}{\}}_{l=\mathrm{1}}^{L-\mathrm{1}})=\\ {\displaystyle \prod _{l=\mathrm{2}}^L}p({\mathit{F}}^{(l)}|{\mathit{U}}^{(l)},{\mathit{X}}^{(l-\mathrm{1})})q({\mathit{U}}^{(l)})q({\mathit{X}}^{(l-\mathrm{1})})\cdot \end{array}

p({\mathit{F}}^{(\mathrm{1})}|{\mathit{U}}^{(\mathrm{1})},{\mathit{X}}^{(\mathrm{0})})q({\mathit{U}}^{(l)})（13）

q({\mathit{X}}^{(l)})={\displaystyle \prod _{d=\mathrm{1}}^{D_l}}q({\mathit{x}}_d^{(l)})q({\mathit{U}}^{(l)})={\displaystyle \prod _{d=\mathrm{1}}^{D_l}}q({\mathit{u}}_d^{(l)})

通过詹森不等式可得lnp(\mathit{Y}|{\mathit{X}}^{\left(\mathrm{0}\right)})的变分下界F_v。令\mathit{\zeta }=({\mathit{F}}^{(\mathrm{1})},{\mathit{F}}^{(\mathrm{2})},\cdots ,{\mathit{F}}^{(L)},{\mathit{U}}^{(\mathrm{1})},{\mathit{U}}^{(\mathrm{2})},\cdots ,{\mathit{U}}^{(L)},{\mathit{X}}^{(\mathrm{1})},

{\mathit{X}}^{(\mathrm{2})},\cdots ,{\mathit{X}}^{(L-\mathrm{1})})，得到

\begin{array}{c}\mathrm{l}\mathrm{n}p(\mathit{Y}|{\mathit{X}}^{(\mathrm{0})})=\mathrm{l}\mathrm{n}{\int }_{\mathit{\zeta }}^{}Q\frac{P}{Q}\mathrm{d}\mathit{\zeta }\ge F_v={\int }_{\mathit{\zeta }}^{}Q\mathrm{l}\mathrm{n}\frac{P}{Q}\mathrm{d}\mathit{\zeta }=\\ \begin{array}{c}{\int }_{\mathit{\zeta }}^{}Q\mathrm{l}\mathrm{n}\frac{p(\mathit{Y}|{\mathit{F}}^{(L)})p({\mathit{U}}^{(L)})}{q({\mathit{U}}^{(L)})}\mathrm{d}\mathit{\zeta }\end{array}+{\displaystyle \sum _{l=\mathrm{1}}^{L-\mathrm{1}}}{\int }_{\mathit{\zeta }}^{}Q\mathrm{l}\mathrm{n}\frac{\mathrm{1}}{q({\mathit{X}}^{(l)})}\mathrm{d}\mathit{\zeta }+\\ {\displaystyle \sum _{l=\mathrm{1}}^{L-\mathrm{1}}}{\int }_{{\mathit{\theta }}_{\mathrm{2}}}^{}Q\mathrm{l}\mathrm{n}\frac{p({\mathit{X}}^{(l)}|{\mathit{F}}^{(l)})p({\mathit{U}}^{(l)})}{q({\mathit{U}}^{(l)})}\mathrm{d}\mathit{\zeta }=\end{array}F_{\mathit{Y}}+{\displaystyle \sum _{l=\mathrm{1}}^{L-\mathrm{1}}}H(q({\mathit{X}}^{(l)}))+{\displaystyle \sum _{l=\mathrm{1}}^{L-\mathrm{1}}}F({\mathit{X}}^{(l)})（14）

式中：H(q({\mathit{X}}^{\left(l\right)}))为分布q({\mathit{X}}^{\left(l\right)})的熵。

根据变分下界的维度分解性质，给定GP映射 Y ～GP（ X ），若其诱导变量为 U 、诱导点为 Z，则变分期望E_{q\left(\mathit{F},\mathit{U}\right)}\{\mathrm{l}\mathrm{n}p(\mathit{Y}|\mathit{F})+\mathrm{l}\mathrm{n}\frac{p(\mathit{U})}{q(\mathit{U})}\}可按输出维度d进行线性分解：

\sum _d\int q({\mathit{f}}_d,{\mathit{u}}_d)\mathrm{l}\mathrm{n}\frac{p({\mathit{u}}_d)p({\mathit{y}}_d|{\mathit{f}}_d)}{q({\mathit{u}}_d)}\mathrm{d}{\mathit{f}}_d\mathrm{d}{\mathit{u}}_d

q(\mathit{F},\mathit{U})=p(\mathit{F}|\mathit{U},\mathit{X})q(\mathit{U})

\begin{array}{c}{F}_{\mathit{Y}}=E_{q({\mathit{F}}^{(L)},{\mathit{U}}^{(L)},{\mathit{X}}^{(L-\mathrm{1})})}\{\mathrm{l}\mathrm{n}p(\mathit{Y}|{\mathit{F}}^{(L)})+\mathrm{l}\mathrm{n}\frac{p({\mathit{U}}^{(L)})}{q({\mathit{U}}^{(L)})}\}\end{array}（15）

\begin{array}{c}{F}_{{\mathit{X}}^{(l)}}=E_{q({\mathit{X}}^{(l)},{\mathit{F}}^{(l)},{\mathit{U}}^{(l)},{\mathit{X}}^{(l-\mathrm{1})})}\{\mathrm{l}\mathrm{n}p({\mathit{X}}^{(l)}|{\mathit{F}}^{(l)})+\mathrm{l}\mathrm{n}\frac{p({\mathit{U}}^{(l)})}{q({\mathit{U}}^{(l)})}\}\end{array}（16）

\begin{array}{c}{F}_{{\mathit{X}}^{(\mathrm{1})}}=E_{q({\mathit{X}}^{(\mathrm{1})},{\mathit{F}}^{(\mathrm{1})},{\mathit{U}}^{\left((\mathrm{1})\right)})}\{\mathrm{l}\mathrm{n}p({\mathit{X}}^{(\mathrm{1})}|{\mathit{F}}^{(\mathrm{1})})+\mathrm{l}\mathrm{n}\frac{p({\mathit{U}}^{(\mathrm{1})})}{q({\mathit{U}}^{(\mathrm{1})})}\}\end{array}（17）

其中，E_{q\left(·\right)}\{\mathit{f}\}表示在变分分布q(·)下对函数 f 取期望，q(·)表示由联合变分分布Q对联合变量中除括号所指变量以外的其余潜在变量积分所得的边缘变分分布。因此，需找到q({\mathit{U}}^{(l)})的最优形式，并证明式（15）～式（17）是可解的。

2.2 诱导变量优化

在式（15）中应用GP诱导变量最优变分分布定理，给定第l层的诱导点 Z^（l）、核函数超参数（长度尺度σ_l 等）及核函数k_l(·,·)，若q({\mathit{u}}_d^{\left(L\right)})\propto p({\mathit{u}}_d^{\left(L\right)})·

\mathrm{e}\mathrm{x}\mathrm{p}(E_{q({\mathit{X}}^{\left(L-\mathrm{1}\right)})}\{\mathrm{l}\mathrm{n}N({\mathit{y}}_d|{\mathit{a}}_d^{\left(L\right)},{\sigma }_L^{\mathrm{2}}{\mathit{I}}_N)\})，则变分目标函数F_Y 达到最大。令{\mathit{\epsilon }}_L=E_{q({\mathit{X}}^{\left(L-\mathrm{1}\right)})}\{\mathrm{l}\mathrm{n}N({\mathit{y}}_d|{\mathit{a}}_d^{\left(L\right)},{\sigma }_L^{\mathrm{2}}{\mathit{I}}_N)\}， 则{\mathit{\epsilon }}_L的最大值为

{\displaystyle \sum _{d=\mathrm{1}}^{D_L}}ln\int p({\mathit{u}}_d^{\left(L\right)})\mathrm{e}\mathrm{x}\mathrm{p}({\mathit{\epsilon }}_L)\mathrm{d}{\mathit{u}}_d^{\left(L\right)}-\frac{D_L}{\mathrm{2}{\sigma }_L^{\mathrm{2}}}\mathrm{t}\mathrm{r}(E_{q({\mathit{X}}^{\left(L-\mathrm{1}\right)})}\{{\tilde{\mathit{K}}}_L\})

因此，通过最大化F_Y，可求得q({\mathit{u}}_d^{\left(L\right)})的条件最优分布。由于{\epsilon }_L是{\mathit{u}}_d^{\left(L\right)}的二次函数，它与先验q({\mathit{u}}_d^{\left(L\right)})相乘后仍为高斯形式，因此最优分布是高斯分布。此时，出现在最大值表达式中的{\mathit{u}}_d^{\left(L\right)}的积分变成标准的高斯积分，从而可求出封闭形式的解析解。

在式（16）中，应用中间层最优诱导变量分布推导，具体而言，对于任意一层 l，在固定其他层的近似后验分布和模型参数不变的条件下，对q({\mathit{X}}^{(l)})进行更新。F_{{\mathit{X}}^{\left(l\right)}}对\forall l\in \left\{\mathrm{2,3},\cdots ,L-\mathrm{1}\right\}最大化时，q({\mathit{u}}_d^{\left(L\right)})与p({\mathit{u}}_d^{\left(L\right)})\mathrm{e}\mathrm{x}\mathrm{p}(E_{q({\mathit{X}}^{\left(l\right)},{\mathit{X}}^{\left(l-\mathrm{1}\right)})}·\{\mathrm{l}\mathrm{n}N({\mathit{x}}_d^{\left(l\right)}|{\mathit{a}}_d^{\left(l\right)},{\sigma }_l^{\mathrm{2}}{\mathit{I}}_N)\})成正比。令{\epsilon }_l=E_{q\left({\mathit{X}}^{\left(l\right)},{\mathit{X}}^{\left(l-\mathrm{1}\right)}\right)}·\{\mathrm{l}\mathrm{n}N({\mathit{x}}_d^{\left(l\right)}|{\mathit{a}}_d^{\left(l\right)},{\sigma }_l^{\mathrm{2}}{\mathit{I}}_N)\}，则相应的变分目标函数F_{{\mathit{X}}^{\left(l\right)}}的最大值为{\displaystyle \sum _{d=\mathrm{1}}^{D_l}}\int p({\mathit{u}}_d^{\left(l\right)})\mathrm{e}\mathrm{x}\mathrm{p}({\mathit{\epsilon }}_{\mathit{l}})\mathrm{d}{\mathit{u}}_d^{\left(l\right)}-\frac{D_l}{\mathrm{2}{\sigma }_l^{\mathrm{2}}}\mathrm{t}\mathrm{r}(E_{q({\mathit{X}}^{\left(l-\mathrm{1}\right)})}\{{\tilde{\mathit{K}}}_l\})。

在式（17）中，应用第一层最优诱导变量分布推导，即给定q(X^{\left(\mathrm{1}\right)})对\foralll以及所有其他参数，F_{{\mathit{X}}^{\left(\mathrm{1}\right)}}最大化时，q({\mathit{u}}_d^{\left(\mathrm{1}\right)})与p({\mathit{u}}_d^{\left(\mathrm{1}\right)})\mathrm{e}\mathrm{x}\mathrm{p}(E_{q({\mathit{X}}^{\left(\mathrm{1}\right)})}·

\{\mathrm{l}\mathrm{n}N({\mathit{x}}_d^{\left(\mathrm{1}\right)}|{\mathit{a}}_d^{\left(\mathrm{1}\right)},{\sigma }_{\mathrm{1}}^{\mathrm{2}}{\mathit{I}}_N)\})成正比。令{\mathit{\epsilon }}_{\mathrm{1}}=E_{q({\mathit{X}}^{\left(\mathrm{1}\right)})}\{\mathrm{l}\mathrm{n}N({\mathit{x}}_d^{\left(\mathrm{1}\right)}|{\mathit{a}}_d^{\left(\mathrm{1}\right)},{\sigma }_{\mathrm{1}}^{\mathrm{2}}{\mathit{I}}_N)\}，则与中间层相关的下界项的最大值为{\displaystyle \sum _{d=\mathrm{1}}^{D_{\mathrm{1}}}}\int p({\mathit{u}}_d^{\left(\mathrm{1}\right)})\mathrm{e}\mathrm{x}\mathrm{p}({\mathit{\epsilon }}_{\mathrm{1}})\mathrm{d}{\mathit{u}}_d^{\left(\mathrm{1}\right)}-\frac{D_{\mathrm{1}}}{\mathrm{2}{\sigma }_{\mathrm{1}}^{\mathrm{2}}}\mathrm{t}\mathrm{r}(\{{\tilde{\mathit{K}}}_{\mathrm{1}}\})。

在式（15）～式（17）中利用诱导变量的条件最优分布消去q({\mathit{U}}^{\left(l\right)})，将F_v中与诱导变量相关的项化简为封闭形式，从而使下界中的各项都能解析计算。因此，将F_v作为证据下界目标函数，通过基于梯度的算法优化模型^［19］最大化边缘似然。优化过程包含两组参数，一组是模型参数集Θ₁，包括{\sigma }_l和k_l(\mathit{x},{\mathit{x}}^{\text{'}})的参数；另一组是变分参数集Θ₂，包括{\mathit{Z}}^{\left(l\right)}和q({\mathit{X}}^{\left(l\right)})的参数。

2.3 量化不确定性

为逼近难以直接求解的真实后验分布p(f|X,y)，变分推断通过构建可优化的变分分布q(f)，将后验近似问题转化为最大化证据下界：

L(q)=E_{q\left(f\right)}[\mathrm{l}\mathrm{n}p(y|f)]-KL(q(f)||p(f))（18）

KL(q(f)||p(f))=\int q(f)\mathrm{l}\mathrm{n}\frac{q(f)}{p(f)}\mathrm{d}f

其中，L(q)为最大化证据下界（evidence lower bound，ELBO）。给定变分分布q（f）下的ELBO，模型能更好地拟合观测数y，这一优化过程等价于最小化KL（q（f）|| p（f）），使变分分布q（f）逼近真实的后验分布。训练完成后，可直接从近似后验q（f）中计算预测的均值和方差来预测模型并量化其不确定性。

每个高斯层不仅输出预测的均值和方差，还将不确定性通过层与层之间传递，使得最终输出层能反映整个网络中逐层传播的不确定性。VDGP的第一层会根据输入 x 计算出对应的预测分布f_{\mathrm{1}}(\mathit{x})~N({\mu }_{\mathrm{1}}(\mathit{x}),{\sigma }_{\mathrm{1}}^{\mathrm{2}}(\mathit{x}))。第二层将该分布作为输入并进行推理，输出新的预测分布f_{\mathrm{2}}(\mathit{x})~N({\mu }_{\mathrm{2}}(\mathit{x}),{\sigma }_{\mathrm{2}}^{\mathrm{2}}(\mathit{x}))。通过多层递归计算，VDGP能捕捉并传播每层的不确定性，使最终输出的预测结果不仅具有更高的准确性，也能全面量化整体的不确定性。

不确定性量化使用后验分布的均值和方差来构建预测值的置信区间。对于给定的输入 x_*，在VDGP中，95%置信区间的上界B_U和下界B_L分别为

B_{\mathrm{U}}=\mu ({\mathit{x}}_{\mathit{*}})+z_{\alpha /\mathrm{2}}\sigma ({\mathit{x}}_{\mathit{*}})（19）

B_{\mathrm{L}}=\mu ({\mathit{x}}_{\mathit{*}})-z_{\alpha /\mathrm{2}}\sigma ({\mathit{x}}_{\mathit{*}})（20）

式中：μ（ x_*）为后验分布的均值；σ（ x_*）为后验分布的标准差；z_α/2为对应标准正态分布的临界值，95%置信水平下，z_α/2=1.96。

2.4 RUL预测的评价指标与预测流程

为评估RUL预测模型在测试数据上的性能，采用如下评分函数来衡量模型的质量^［20］：

（21）

h_i=E_i^{\text{'}}-E_i

式中：n为测试集中的样本总数；E_i^{\text{'}}为测试样本i的测试RUL；E_i 为测试样本i的真实RUL。

均方根误差E_{\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{E}}=\sqrt[]{\frac{\mathrm{1}}{n}{\displaystyle \sum _{i=\mathrm{1}}^n}{h}_i^{\mathrm{2}}}常用于评估RUL预测。较小的E_SCORE和E_RMSE均表明预测结果与真实寿命更接近，模型性能更优。

采用平均绝对误差E_MAE衡量预测准确度。

提出的VDGP方法框架如图4所示。首先，预处理原始数据，并将其分为训练集和测试集。训练过程中，引入诱导变量 U^（l）和诱导点 Z^（l）以减小计算复杂度并进行变分推断。通过变分推断建立并计算目标函数F_v，然后利用梯度下降算法不断更新模型参数以最大化F_v，直到达到收敛条件。模型训练完成后，使用测试集进行RUL预测。

\mathrm{l}\mathrm{n}p({\mathit{y}}_t)\ge F_v(\mathit{Z},\varphi )=

{\int }_u^{}\varphi (\mathit{u})\mathrm{l}\mathrm{n}\frac{N({\mathit{y}}_t|a_t,{\sigma }_{\epsilon }^{\mathrm{2}}{\mathit{I}}_N)p(\mathit{u})}{\varphi (\mathit{u})}\mathrm{d}\mathit{u}-\frac{\mathrm{1}}{\mathrm{2}{\sigma }_{\epsilon }^{\mathrm{2}}}\mathrm{t}\mathrm{r}({\tilde{\mathit{K}}}_L)

{\theta }_{t+\mathrm{1}}={\theta }_t-\eta \nabla F_v({\theta }_t)

3 试验验证

3.1 C-MPASS数据集验证

NASA的C-MAPSS数据集^［21］是分析涡扇发动机退化的仿真数据集。该数据集（含690万个数据点）分为4个子数据集（FD001～FD004），每个子数据集包含多个不同引擎的多变量时间序列，并已划分为训练集和测试集。

C-MAPSS数据集包含来自21个传感器的多元时间数据。为消除不同传感器输出的数值大小和单位差异，在预处理阶段将每个传感器的原始读数映射至［-1，1］内。第i个样本的第j个传感器的读数x_i，j 归一化后为

x_{i,j}^{\text{'}}=\frac{\mathrm{2}(x_{i,j}-x_j^{\mathrm{m}\mathrm{i}\mathrm{n}})}{x_j^{\mathrm{m}\mathrm{a}\mathrm{x}}-x_j^{\mathrm{m}\mathrm{i}\mathrm{n}}}-\mathrm{1}

式中：x_j^{\mathrm{m}\mathrm{a}\mathrm{x}}、x_j^{\mathrm{m}\mathrm{i}\mathrm{n}}分别为训练数据集中第j个传感器读数的最大值和最小值。

RUL标签的建模直接影响预测的准确度。发动机在开始时运行正常，然后逐渐恶化，因此采用分段线性RUL代替实际的RUL^［22-23］，如图5所示。假设RUL标签RUL_early在开始时保持不变（如图5，［0，50 min］的RUL是保持不变的），然后RUL线性退化，直到C-MAPSS发动机失效。分段RUL方案中，将RUL_early设为125^［22-23］。

3.1.1 RUL预测核函数分析

核函数和高斯层数对预测的准确度具有重大影响。通过选择合适的核函数和层数，VDGP 模型能有效地捕获复杂的非线性退化特征，提高 RUL 预测的精度与模型的灵活性和通用性^［24］。本文首先构建两层VDGP模型，并在子数据集 FD001上通过应用多种常用核函数筛选表现优异的核函数。模型两层采用相同的核函数，隐藏单元数设为10，每层的诱导变量数设为100，每个核函数的初始长度为10，初始隐噪声方差为0.01。本文采用梯度下降法优化VDGP模型。梯度下降法在非凸目标函数和高维参数空间中容易陷入局部最优解，且对学习率敏感，导致收敛慢或参数更新振荡。为解决这些问题，采用以下方法进行改进：①引入Adam优化算法，Adam首先利用一阶与二阶矩的估计对每个参数动态调整学习率，再通过自适应学习率机制提高稳定性和收敛速度；②合理初始化变分的分布参数，并在ELBO目标函数中加入正则化项来增强优化的稳定性；③后期采用学习率衰减策略，避免过大步长跳过全局最优解。

诱导点位置 Z^（l）、诱导变量均值 m^（l）、协方差矩阵 S^（l）等变分参数按文献［25］的默认方法随机初始化。采用表1中的核函数进行实验，所得结果如表2、表3所示。后期（退化阶段）是离故障更近的阶段，关系设备的故障和维护，因此后期的性能更受关注。除在每个测试引擎的最后一个周期评估指标外，表2、表3中还展示了在RUL标签不大于15的评估指标。

由表2可知，RBF和Matern5/2的总体表现最优，其次是RQ核函数。为进一步分析不同核函数组合的使用效果，在两层VDGP中按所有可能的组合方式配置，并比较其性能，结果如表3所示。表现最优的组合核函数：①第一层用RBF核，第二层用Matern5/2核；②第一层用Matern5/2核，第二层用RBF核。不同层使用核函数的顺序会显著影响预测结果。

在两层VDGP方法中，第一层用于对输入特征进行初步变换，第二层对第一层的输出进行建模。第一层的数学表达式为

H~GP(\mathrm{0},k_{\mathrm{1}}(X,X^{\text{'}}))

其中，H为隐藏层的输出，核函数k₁描述输入特征X之间的相似性。第二层的数学表达式为

Y~GP(\mathrm{0},k_{\mathrm{2}}(H,H^{\text{'}}))

其中，k₂是第二层的核函数，用于刻画隐藏层输出H与预测输出和Y之间的函数映射关系。 VDGP中，隐藏层输出H的分布由第一层核函数k₁和变分分布q（u）共同决定，因此有

q(H)=\int p(H|X,u)q(u)\mathrm{d}u

k₁决定了H的平滑性和尺度。当数据包含快速变化或高频退化特征时，选择较低平滑度（更“非平滑”）的核函数能增强模型对这些复杂行为的表达能力。对于变化较平缓的退化过程，高平滑度核函数通常表现更稳定。第二层核函数k₂决定了隐藏表示H向最终输出Y的映射方式，即控制H与Y之间的相关性结构与平滑性。当两层核函数具有一致的特性时，输入到输出的映射将更连续，隐藏表示的变化更平稳，从而提高模型性能、增强梯度计算的稳定性，并使优化目标ELBO更易收敛。核函数特性不同时，第一层生成的隐藏层输出H较为平滑，第二层的高灵活性可能引入剧烈的梯度波动，导致优化器难以找到全局最优解，模型难以正确捕捉输入X和输出Y之间的关系。

3.1.2 VDGP层数分析

对于深度学习及其扩展模型，更深的结构在某些情况下可提高预测精度，但会带来过拟合和计算复杂度增加的问题。层数大于3时^［26］，VDGP模型的性能提升并不显著，计算复杂度却线性增加，因此，需通过实验进一步确定VDGP模型的最优层数。为确定最优层数，基于核函数RBF、Matern5/2和RQ，在3层VDGP模型上进行实验，并将实验结果（表4）与2层VDGP模型的结果（表2）进行对比，从而确定最优层数。

2层VDGP和3层VDGP的结果如图6所示，在真实E_RUL≤15的测试子集上评估指标并绘制于图7，3层的VDGP核函数在名称后面加入后缀-3。

评价指标E_RMSE、E_SCORE、E_MAE越小，效果越好。由图6可看出，2层的VDGP效果均优于3层的VDGP。由图7可看出最优的指标依然是2层的VDGP。因此，2层的VDGP效果最佳。

3.1.3 RUL预测方法对比分析

基于核函数的选取实验，确定表5所示的4个VDGP模型。根据4个子数据集不同的复杂程度，将子数据集FD001～FD004对应VDGP的隐藏单元数分别设置为10、120、20、140。

VDGP的实验结果和其他方法的实验结果^［27-41］如表6所示，VDGP在4个子数据集上的预测性能均显著优于其他模型。复杂多工况数据集FD002和FD004上，VDGP充分展现了动态工况下建模的能力。VDGP采用多层高斯过程构成的分层架构，各层可使用不同核函数以递归方式提取和表达数据的复杂特征，将原始输入逐步映射到更高维的特征空间，从而增强对多变工况的建模能力。

进一步分析发现，所有模型在单一工况数据集FD001、FD003上的预测性能普遍优于多工况数据集FD002和FD004。FD002和FD004包含6种不同工况，导致退化特征与RUL的关系更为复杂，增加了模型学习的难度。为验证这一猜测，选取低压压缩机出口温度（s_2）和高压压缩机出口温度（s_3）与RUL的关系进行分析。本文仅展示FD001和FD002中的3个发动机数据，如图8、图9所示。FD003的数据趋势与FD001相似，FD004的数据趋势与FD002相似，因此不再展示。

对于数据集FD001，s_2、s_3在发动机退化初期（RUL在［175，275］）变化不明显；在退化中期（RUL在［75，175］），变化幅度显著增大；到了退化后期（RUL在［0，75］），变量值随RUL减少而非线性加速增大。与FD001相比，FD002工况更多且变化复杂，呈现出幅度更大、形态更不规则的波动。这表明工况变化对发动机的退化特征有显著影响，导致各模型对FD001、FD003的预测效果优于FD002和FD004。

为进一步得出最优的VDGP模型，对四种 VDGP 结构进行对比实验，并在E_RUL≤15的后期退化阶段对模型性能进行了重点比较。由表7知，VDGP-M-M总体表现最优；所有类型的VDGP在后期退化阶段的预测准确度更高，再次说明VDGP模型具有更大的实际应用价值。

另外，将VDGP和GP进行对比实验来验证VDGP方法的预测准确度与置信区间。GP在大规模数据集上模型训练容易超出计算机内存。为避免超出内存限制，从数据集FD001~FD004里分别取前100个发动机训练GP。FD001和FD003训练集恰好各有100个发动机。基于表7所示的实验结果，选择最优的VDGP-R-R和VDGP-M-M，并与GP模型使用相同的数据集进行训练。实验结果如表8所示。

分析表8所示结果可知，对于“最后周期”和“E_RUL≤15”，VDGP-R-R和VDGP-M-M的预测精度均超过GP，这表明VDGP方法的整体预测性能优于GP。因此，使用完整训练集的模型预测结果优于使用部分训练集的模型。为进一步比较置信区间效果，选取表8中表现最优的VDGP-M-M与最优的GP（RBF）方法，并在相同的FD001数据集上使用相同的训练与测试进行对比。由图10可以看出，VDGP的置信区间比GP的置信区间更窄。这表明VDGP在预测中的不确定性更小，预测结果的可信度更高。

3.1.4 VDGP的RUL预测性能分析

模型VDGP-M-M在FD001～FD004上的预测结果如图11～图14所示。为进一步展示VDGP模型在整个退化过程中对单个发动机RUL预测的性能，从所有子数据集中随机选取4个测试发动机，并给出其从初始运行到失效的完整RUL预测轨迹，如图15～图18所示。

3.2 风机行星齿轮数据集验证

数据集来源于笔者开展的风机齿轮箱行星齿轮故障加速寿命实验。

3.2.1 试验数据获取

1）试验设备。如图19所示，试验系统包括行星齿轮箱、电机、制动器及其控制器、数据采集系统（加速度传感器、电涡流传感器等）。加速度传感器分别安装在测点1和测点3，其中，测点1位于行星齿轮箱的第二级，传感器安装方向与水平方向成45°；测点2、测点3分别安装在第三级齿轮箱的垂直和水平位置。电涡流传感器安装在测点4，通过测速齿轮盘获取输出轴的实时转速。

为提高试验效率并降低试验成本，在不改变磨损失效机制的前提下，采用了2种加速疲劳寿命试验的方法：①利用制动器施加一定的载荷，加速老化进程；②将第三级轮系中的某个行星齿轮替换为未进行表面热处理的齿轮，加快磨损。

设置采样频率为20 480 Hz，单次采样时长为10 s，每隔10 min进行一次采样。行星齿轮箱高速轴的转频为600 Hz，磁粉制动器的转矩设定为8.33 N·m。

齿轮在早期阶段健康，但性能逐渐下降，直到发生故障。整个生命周期内，4个传感器对齿轮进行数据采样，记录了从健康到失效的完整退化轨迹。

在齿轮数据集上，将所提方法与代表性基线模型进行对比评估，并采用与C-MAPSS发动机实验相同的数据预处理流程。该数据集包含4个传感器的多元时间数据，试验的齿轮数量为1，以较高频率连续采样，形成高密度时序记录。因此，间隔抽取一半的数据用于训练，将剩下的数据用于模型的性能测试。归一化方面沿用C-MAPSS发动机实验部分所示的Min–Max方法，将所有传感器通道映射至区间［-1，1］，并使用训练集统计量对测试集进行同一变换。

3.2.2 RUL预测方法的对比分析

表5中的4种VDGP变体模型在C-MAPSS发动机数据集上的应用效果已被验证，但这些模型在新数据集上的泛化能力有待验证。因此，选择LSTM^［28］、MLP^［42］、CNN^［27］、LSTM with attention^［29］、Bi-level LSTM Scheme^［35］、CBLSTM^［40］、LSTM-Z^［41］与4种VDGP变体模型进行对比分析，它们的性能结果见表9。

由表9可得以下结论：

1）VDGP在齿轮数据集上的表现优于其他模型，泛化能力和鲁棒性更强，后期退化阶段的预测准确度更高。这说明VDGP具有更大的应用价值。

2）4种VDGP模型中，VDGP-R-R和VDGP-M-M的性能最优。在整个测试集的RUL预测中，VDGP-M-M表现最优，VDGP-R-R在后期退化阶段表现最优，但VDGP-R-M和VDGP-M-R组合的表现较差。因此不同层使用不同核函数的VDGP模型泛化能力较差，这与3.1.1节的不同核函数组合理论的分析一致。

基于表9，将VDGP-R-R和VDGP-M-M与GP进行比较。实验结果如表10和图20～图23所示，可得如下结论：

1）VDGP的预测准确度显著优于GP。与GP相比，VDGP的预测置信区间更小。基于齿轮数据集上的对比预测结果，VDGP的置信区间始终较小。这证明了VDGP的不确定性量化效果好，且在处理复杂数据集时有更好的鲁棒性和稳定性。

2）VDGP-M-M在整体上的RUL预测效果优于VDGP-R-R，但VDGP-R-R在后期的预测效果更好。因此，可根据实际需求选择应用这两种模型，以实现更为准确的RUL预测。

4 结论

1）提出了基于变分深度高斯过程（VDGP）的预测方法，即通过对原始数据的预处理与多层次建模精准预测剩余使用寿命（RUL），有效刻画了预测结果的概率分布。采用引入诱导点和变分推断的方式求解，提高了模型处理大规模状态监测数据的性能。VDGP不仅能提高预测准确度，还能提供置信区间量化的预测不确定性，为结果的可信性评估提供重要支持。

2）为验证模型性能，进行了多层次及不同核函数组合的对比研究，确定了核函数和模型层数等关键超参数。其中，2层VDGP-M-M变体在数据集FD001上的 均方根误差、评分函数分别为16.80和578.33，较其他核函数平均减小3.69% 和2.39%。与高斯过程（GP）方法相比，VDGP提高了预测准确度，置信区间也更小，展现了更高的稳定性和可靠性。

3）C-MAPSS和风机行星齿轮数据集的测试验证了VDGP方法的有效性和通用性。C-MAPSS数据集的对比结果表明，VDGP 模型整体性能优异，预测准确度较高，在不同规模的训练数据上展现出良好的稳健性。风机齿轮数据集的实验验证了模型的通用性。实验结果显示，VDGP在预测准确度、稳定性和鲁棒性上均优于RCNN+ABi-LSTM方法。

参考文献

[1] GAO Pengjie， WANG Junliang， SHI Ziqi， et al. Long-term Temporal Attention Neural Network with Adaptive Stage Division for Remaining Useful Life Prediction of Rolling Bearings［J］. Reliability Engineering & System Safety， 2024， 251： 110218.

[2] LIN Ruiguan， YU Yaowei， WANG Huawei， et al. Remaining Useful Life Prediction in Prognostics Using Multi-scale Sequence and Long Short-term Memory Network［J］. Journal of Computational Science， 2022， 57： 101508.

[3] CHEN Chuang， LU Ningyun， JIANG Bin， et al. Prediction Interval Estimation of Aeroengine Remaining Useful Life Based on Bidirectional Long Short-term Memory Network［J］. IEEE Transactions on Instrumentation and Measurement， 2021， 70： 3527213.

[4] 刘小峰， 亢莹莹， 柏林. 轴承自驱式独立退化轨迹构建与剩余寿命灰色预测［J］. 中国机械工程， 2024， 35（9）： 1613-1621.

[5] LIU Xiaofeng， KANG Yingying， BO Lin. Self-driven Independent Degradation Trajectory Construction and Remaining Life Gray Prediction for Bearings［J］. China Mechanical Engineering， 2024， 35（9）： 1613-1621.

[6] CHEN Zuoyi， HUANG Hong-Zhong， DENG Zhongwei， et al. Shrinkage Mamba Relation Network with Out-of-distribution Data Augmentation for Rotating Machinery Fault Detection and Localization under Zero-faulty Data［J］. Mechanical Systems and Signal Processing， 2025， 224： 112145.

[7] ZHU Ting， CHEN Zhen， ZHOU Di， et al. Adaptive Staged Remaining Useful Life Prediction of Roller in a Hot Strip Mill Based on Multi-scale LSTM with Multi-head Attention［J］. Reliability Engineering & System Safety， 2024， 248： 110161.

[8] CAO Wei， MENG Zong， LI Jimeng， et al. A Remaining Useful Life Prediction Method for Rolling Bearing Based on TCN-transformer［J］. IEEE Transactions on Instrumentation and Measurement， 2024， 74： 3501309.

[9] PENG Peng， LI Yonghua， GUO Zhongyi. High-performance Remaining Useful Life Prediction for Aeroengine Based on Combining Health States and Trajectory Similarity［J］. Engineering Applications of Artificial Intelligence， 2025， 141： 109799.

[10] CHEN Zesheng， TU Xiaotong， HU Yue， et al. Real-time Bearing Remaining Useful Life Estimation Based on the Frozen Convolutional and Activated Memory Neural Network［J］. IEEE Access， 2019， 7： 96583-96593.

[11] 夏然， 苏春. 基于健康因子和混合Bi-LSTM-NAR模型的锂离子电池剩余寿命预测［J］. 中国机械工程， 2024， 35（5）： 851-859.

[12] XIA Ran， SU Chun. Remaining Useful Life Prediction for Lithium-ion Batteries Based on Health Indicators and Hybrid Bi-LSTM-NAR Model［J］. China Mechanical Engineering， 2024， 35（5）： 851-859.

[13] WANG Zhe， YANG Lechang， FANG Xiaolei， et al. Image-based Remaining Useful Life Prediction through Adaptation from Simulation to Experimental Domain［J］. Reliability Engineering & System Safety， 2025， 255： 110668.

[14] LUO Peien， YIN Zhonggang， ZHANG Yanqing， et al. A Novel Whale Optimization Algorithm Based on Music Theory Knowledge for RUL Prediction of Motor Bearing［J］. IEEE Transactions on Instrumentation and Measurement， 2023， 72： 3534211.

[15] 曹现刚， 叶煜， 赵友军， 等. 基于KPCA-LSTM的旋转机械剩余使用寿命预测［J］. 振动与冲击， 2023， 42（24）： 81-91.

[16] CAO Xiangang， YE Yu， ZHAO Youjun， et al. Remaining Useful Life Prediction of Rotating Machinery Based on KPCA-LSTM［J］. Journal of Vibration and Shock， 2023， 42（24）： 81-91.

[17] LI Biao， TANG Baoping， DENG Lei， et al. Self-attention ConvLSTM and Its Application in RUL Prediction of Rolling Bearings［J］. IEEE Transactions on Instrumentation and Measurement， 2021， 70： 3518811.

[18] 刘华开， 丁康， 何国林， 等. 联合故障机理和卷积神经网络的齿轮剩余使用寿命预测方法研究［J］. 机械工程学报， 2024， 60（12）： 116-125.

[19] LIU Huakai， DING Kang， HE Guolin， et al. Research on the Prediction Method of Remaining Useful Life of Gears by Combining Fault Mechanism and Convolutional Neural Network［J］. Journal of Mechanical Engineering， 2024， 60（12）： 116-125.

[20] LIU Xiaofei， LEI Yaguo， LI Naipeng， et al. RUL Prediction of Machinery Using Convolutional-vector Fusion Network through Multi-feature Dynamic Weighting［J］. Mechanical Systems and Signal Processing， 2023， 185： 109788.

[21] WANG Biao， LEI Yaguo， LI Naipeng， et al. Multiscale Convolutional Attention Network for Predicting Remaining Useful Life of Machinery［J］. IEEE Transactions on Industrial Electronics， 2020， 68（8）： 7496-7504.

[22] DUVENAUD D. Automatic model construction with Gaussian processes［D］. Cambridge， UK： University of Cambridge， 2014.

[23] LI Yanhui， YANG Pengfei， BAI Lu， et al. BRDF Modeling and Optimization of a Target Surface Based on the Gradient Descent Algorithm［J］. Applied Optics， 2023， 62（36）： 9486-9492.

[24] HEIMES F O. Recurrent Neural Networks for Remaining Useful Life Estimation［C］∥2008 International Conference on Prognostics and Health Management. Denver： IEEE， 2008： 1-6..

[25] SAXENA A， GOEBEL K， SIMON D， et al. Damage Propagation Modeling for Aircraft Engine Run-to-failure Simulation［C］∥2008 International Conference on Prognostics and Health Management. Denver： IEEE， 2008： 1-9.

[26] LI Han， ZHAO Wei， ZHANG Yuxi， et al. Remaining Useful Life Prediction Using Multi-scale Deep Convolutional Neural Network［J］. Applied Soft Computing， 2020， 89： 106113.

[27] LIU Hui， LIU Zhenyu， JIA Weiqiang， et al. A Novel Deep Learning-based Encoder-decoder Model for Remaining Useful Life Prediction［C］∥2019 International Joint Conference on Neural Networks （IJCNN）. Budapest： IEEE， 2019： 1-8.

[28] THANG B， HERNANDEZ-LOBATO D， HERNANDEZ-LOBATO J， et al. Deep Gaussian Processes for Regression using Approximate Expectation Propagation ［C］∥Proceedings of the 33rd International Conference on Machine Learning （ICML）. New York， 2016： 1472-1481.

[29] DAMIANOU A. Deep Gaussian Processes and Variational Propagation of Uncertainty［D］. Sheffield： University of Sheffield， 2015.

[30] XIE Yucen， ZOU Jianxiao， PENG Chao， et al. A Novel PEM Fuel Cell Remaining Useful Life Prediction Method Based on Singular Spectrum Analysis and Deep Gaussian Processes［J］. International Journal of Hydrogen Energy， 2020， 45（55）： 30942-30956.

[31] SATEESH BABU G， ZHAO Peilin， LI Xiaoli. Deep Convolutional Neural Network Based Regression Approach for Estimation of Remaining Useful Life［C］∥Database Systems for Advanced Applications. Cham： Springer， 2016： 214-228.

[32] ZHENG Shuai， RISTOVSKI K， FARAHAT A， et al. Long Short-term Memory Network for Remaining Useful Life Estimation［C］∥2017 IEEE International Conference on Prognostics and Health Management （ICPHM）. Dallas： IEEE， 2017： 88-95.

[33] CHEN Zhenghua， WU Min， ZHAO Rui， et al. Machine Remaining Useful Life Prediction via an Attention-based Deep Learning Approach［J］. IEEE Transactions on Industrial Electronics， 2020， 68（3）： 2521-2531.

[34] BENKER M， FURTNER L， SEMM T， et al. Utilizing Uncertainty Information in Remaining Useful Life Estimation via Bayesian Neural Networks and Hamiltonian Monte Carlo［J］. Journal of Manufacturing Systems， 2021， 61： 799-807.

[35] GAO Hui， LI Yibin， ZHAO Ying， et al. Dual Channel Feature Attention-based Approach for RUL Prediction Considering the Spatiotemporal Difference of Multisensor Data［J］. IEEE Sensors Journal， 2023， 23（8）： 8514-8525.

[36] ElDALI M， KUMAR K D. Fault Diagnosis and Prognosis of Aerospace Systems Using Growing Recurrent Neural Networks and LSTM［C］∥2021 IEEE Aerospace Conference （50100. SkyBig： IEEE， 2021： 1-20.

[37] DUAN Yuhang， LI Honghui， HE Mengqi， et al. A BiGRU Autoencoder Remaining Useful Life Prediction Scheme with Attention Mechanism and Skip Connection［J］. IEEE Sensors Journal， 2021， 21（9）： 10905-10914.

[38] XIA Tangbin， SHU Junqing， XU Yuhui， et al. Multiscale Similarity Ensemble Framework for Remaining Useful Life Prediction［J］. Measurement， 2022， 188： 110565.

[39] SONG Tao， LIU Chao， WU Rui， et al. A Hierarchical Scheme for Remaining Useful Life Prediction with Long Short-term Memory Networks［J］. Neurocomputing， 2022， 487： 22-33.

[40] LIN Lin， WU Jinlei， FU Song， et al. Channel Attention & Temporal Attention Based Temporal Convolutional Network： a Dual Attention Framework for Remaining Useful Life Prediction of the Aircraft Engines［J］. Advanced Engineering Informatics， 2024， 60： 102372.

[41] LI Yasong， ZHOU Zheng， SUN Chuang， et al. Life-cycle Modeling Driven by Coupling Competition Degradation for Remaining Useful Life Prediction［J］. Reliability Engineering & System Safety， 2023， 238： 109480.

[42] 闫啸家， 梁伟阁， 张钢， 等. 基于RCNN-ABiLSTM的机械设备剩余寿命预测方法［J］. 系统工程与电子技术， 2023， 45（3）： 931-940.

[43] YAN Xiaojia， LIANG Weige， ZHANG Gang， et al. Prediction Method for Mechanical Equipment Based on RCNN-ABiLSTM［J］. Systems Engineering and Electronics， 2023， 45（3）： 931-940.

[44] WANG Huan， LIU Zhiliang， PENG Dandan， et al. Understanding and Learning Discriminant Features Based on Multiattention 1DCNN for Wheelset Bearing Fault Diagnosis［J］. IEEE Transactions on Industrial Informatics， 2020， 16（9）： 5735-5745.

[45] YOU Dazhang， CHEN Linbo， LIU Fei， et al. Intelligent Fault Diagnosis of Bearing Based on Convolutional Neural Network and Bidirectional Long Short-term Memory［J］. Shock and Vibration， 2021， 2021（1）： 7346352.

[46] 庄雨璇， 李奇， 杨冰如， 等. 基于LSTM的轴承故障诊断端到端方法［J］. 噪声与振动控制， 2019， 39（6）： 187-193.

[47] ZHUANG Yuxuan， LI Qi， YANG Bingru， et al. An End-to-end Approach for Bearing Fault Diagnosis Based on LSTM［J］. Noise and Vibration Control， 2019， 39（6）： 187-193.

基金资助

国家自然科学基金(62303065)

北京信息科技大学勤信人才项目(QXTCPC202120)