CN106202792A - 一种融合电流和轴温的动量轮剩余寿命预测方法 - Google Patents

Abstract

本发明一种融合电流和轴温的动量轮剩余寿命预测方法，实现步骤为：首先分别建立在电流作为单退化量和轴承温度作为单退化量时的寿命分布函数，然后选取合适的二元Copula分布函数作为确定联合基于电流和基于轴承温度的联合分布函数，对Copula分布函数求取偏导数之后根据规范化的极大似然估计法对Copula概率密度函数中的未知参数进行估计，最后根据概率论中二维分布函数性质，获得融合电流和轴温的剩余寿命分布函数，本发明解决了现有技术中存在的动量轮寿命预测方法大多仅考虑单一退化量对剩余寿命预测不准确的问题。

Description

一种融合电流和轴温的动量轮剩余寿命预测方法

技术领域

本发明属于设备故障预测与健康管理技术领域，具体涉及一种融合电流和轴温的动量轮剩余寿命预测方法。

背景技术

传统剩余寿命预测方法大多是以产品某个单一物理量的历史寿命数据为基础，通过对其单变量的退化数据分析确定产品的剩余寿命分布轨迹。然而各种新型的行业，如机械工业、化学工业、电子工业、航天工业等众多领域涌现出了大量使用寿命长、可靠性高的产品。传统的寿命实验和加速寿命实验都很难在较短的时间内获得产品的失效数据，导致传统的剩余寿命预测不能得到可靠的预测结果。并且该类产品往往成本高昂、数量稀少、失效机理复杂，难以进行大样本加速寿命试验，因此目前没有公认的、成熟的动量轮失效模型。但是高可靠性、长寿命产品的失效与其性能退化存在着密切的联系，因此可以通过选取合适的退化参数，采用得到的性能退化数据进行可靠性建模。

动量轮是航空产品上必不可少的复杂设备，研究动量轮退化机理发现其失效原因与轴承温度、润滑剂系统、电流、转速等因素有关。现有动量轮的可靠性分析与寿命预测方法中，大多都是以动量轮单一的退化量作为关键性能参数，进行建模分析，估计动量轮的可靠性与剩余寿命。实际上，动量轮寿命受多个因素的共同影响。在外太空中，太阳的辐射强度周期性变化导致温度测量值会存在较大的波动。由于动量轮电路中的电流较小，当存在较大的测量噪声时，电流的真实值会被测量噪声所淹没，难以得出较为准确的电流值。因此，仅考虑单一退化量对剩余寿命的影响，不能得到较准确的退化模型。针对上述这些问题，我们提出将动量轮的轴温和电流作为关键性能参数，进而对动量轮可靠性分析与剩余寿命预测。

发明内容

本发明的目的是提供一种融合电流和轴温的动量轮剩余寿命预测方法，解决了现有技术中存在的动量轮寿命预测方法大多仅考虑单一退化量对剩余寿命预测不准确的问题。

本发明所采用的技术方案是，一种融合电流和轴温的动量轮剩余寿命预测方法，具体按照以下步骤实施：

步骤1、建立基于电流退化数据的失效模型，获得基于电流的动量轮剩余寿命预测函数；

步骤2、建立基于轴温退化数据的失效模型，得到基于轴温的动量轮剩余寿命预测函数；

步骤3、采用Copula函数融合电流和轴温的退化寿命。

本发明的特点还在于，

步骤1具体按照以下步骤实施：

步骤(1.1)、首先建立基于电流退化数据的失效模型y₁(t)＝α₁+β₁t，参数α₁、β₁分别代表电流的初始值及耗损率，y₁(t)表示电流单退化量下动量轮的性能退化函数,t为动量轮运行时间变量，其中参数α₁、β₁之间相互独立，和分别表示α₁的均值和方差，同理和分别表示β₁的均值和方差，因此y₁(t)的均值估计值和方差估计值表示为μ_α1，μ_β1，根据多个动量轮样本的退化量在不同时刻的观测值使用线性回归的方法估计得到，最终得到基于电流退化数据的失效模型其中，和分别表示y₁(t)均值和方差的估计值；

步骤(1.2)、根据步骤(1.1)得到的基于电流退化数据的失效模型求解基于电流的动量轮寿命失效分布函数F₁(t)：

$F\left(t\right)=P\left(y_1\right(t)>D_1)=\frac{1}{\sqrt{2{\mathrm{\πD}}_1}}{\∫}_{-\mathrm{\∞}}^{\stackrel{\‾}{x}}{e}^{-\frac{x^2}{2}}dx$

上式中，x为数学意义上的广义积分项，D₁为电流失效阈值，P(y₁(t)＞D₁)表示电流大于失效阈值的概率，表示时失效分布函数积分上限，中μ₁(t)，σ₁(t)分别为所述步骤(1.1)中的估计值

步骤(1.3)、根据可靠度函数的性质可知，F₁(t)+R₁(t)＝1，则可获得基于电流的动量轮剩余寿命预测函数R₁(t)：

$R_1\left(t\right)=1-F_1\left(t\right)=1-\frac{1}{\sqrt{2\mathrm{\π}}{D}_1}{\∫}_{-\mathrm{\∞}}^{\stackrel{\‾}{x}}{e}^{-\frac{x^2}{2}}dx$

式中F₁(t)表示基于电流的动量轮寿命失效分布函数。

步骤2具体按照以下步骤实施：

步骤(2.1)、鉴于轴承温度的特点，使用有漂移的维纳过程来对基于轴温的动量轮寿命退化模型进行建模，利用线性函数描述其漂移特性：

y₂(t)＝a₀+a₁t+σ_wW(t)

上式中，a₀为轴承温度遥测数据的初值，a₁为轴承温度遥测值漂移参数，W(t)为标准的维纳过程，σ_w为扩散参数；

步骤(2.2)、以相同时间间隔在动量轮轴承温度遥测数据中取n₁+1个点：其中根据所述步骤(2.1)中的模型把相邻两个采样点数据代入，相减得到Δy_i＝a₁Δt_i+σ_wΔW(t_i)，Δt_i＝(t_i-t_i-1)，ΔW(t_i)＝[W(t_i)-W(t_i-1)]，Δy_i＝y_i-y_i-1，由维纳过程定义知，

步骤(2.3)、由维纳过程平稳独立增量可得步骤(2.2)中的 的联合概率密度，使用极大似然法估计未知参数，构造样本似然函数L(a₁,σ_w)：

$L(a_1,{\mathrm{\σ}}_w)=f({\mathrm{\Δy}}_1,{\mathrm{\Δy}}_2,...,{\mathrm{\Δy}}_{n_1})=f\left({\mathrm{\Δy}}_1\right)f\left({\mathrm{\Δy}}_2\right)...f\left({\mathrm{\Δy}}_{n_1}\right);$

步骤(2.4)、求解步骤(2.3)中构造的似然函数，得a₁的极大似然估计值和σ_w的极大似然估计值如下：

${\hat{a}}_1=\frac{1}{n_1}\underset{i=1}{\overset{n_1}{\Σ}}\frac{{\mathrm{\Δy}}_i}{{\mathrm{\Δt}}_i},$

${\widehat{\mathrm{\σ}}}_w={\[\frac{1}{n_1}\underset{i=1}{\overset{n_1}{\Σ}}\frac{{({\mathrm{\Δy}}_i-{\hat{a}}_1{\mathrm{\Δt}}_i)}^2}{{\mathrm{\Δt}}_i}\]}^{\frac{1}{2}};$

步骤(2.5)、设动量轮样本数为n₂，各动量轮轴承温度初值利用步骤(2.4)的方法得到各动量轮的漂移参数扩散参数根据估计得到其相应的极大似然估计值：

${\widehat{\mathrm{\μ}}}_{a_0}=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[a_0\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_0}\]}^2,$

${\widehat{\mathrm{\σ}}}_{a_0}^2=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[a_0\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_0}\]}^2,$

${\widehat{\mathrm{\μ}}}_{a_1}=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\hat{a}}_1\left(j\right),$

${\widehat{\mathrm{\σ}}}_{a_1}^2=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[{\hat{a}}_1\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_1}\]}^2;$

分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差。分别为的估计值。

步骤(2.6)、由于物理失效和性能失效都发生在性能首次超过失效阈值的时间，即一旦产品的劣化过程Y(t)达到某个特定的水平ζ时产品失效，ζ为失效阈值，因此，定义产品性能退化的样本轨迹模型Y(t)为超过ζ的函数，数学描述如下：

T＝inf{Y(t)≥ζ}＝{t|Y(t)≥ζ,Y(s)≤ζ,0≤s≤t}

式中，Y(t)表示产品的劣化过程退化变量函数模型，ζ为失效阈值，t表示时间，s为小于t时刻的某一时间，动量轮轴承温度性能退化过程Y(t)＝y₂(t)是维纳过程，y₂(t)为步骤(2.1)中温度特性描述函数，给定a₀，a₁，σ_w值时，动量轮的寿命分布概率密度函数，首次达到失效阈值的时间为逆高斯分布，数学描述为：

$f\left(t\right|a_0,a_1,\mathrm{\ζ})=\frac{\mathrm{\ζ}}{\sqrt{2{\mathrm{\π\σ}}_w^2{t}^3}}{e}^{\frac{{(\mathrm{\ζ}-a_0-a_{1}t)}^2}{2{\mathrm{\σ}}_w^{2}t}}$

式中a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数。ζ为失效阈值，t表示时间；

步骤(2.8)、设动量轮的失效阈值ζ近似服从u_ζ＝75℃，σ_ζ＝3℃的正态分布，即ζ～N(75,3²)，则在考虑随机失效阈值的情况下，动量轮的寿命分布概率密度函数为：

$f\left(t\right)={\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数。分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差；

步骤(2.9)、对步骤(2.8)得到的动量轮的寿命分布概率密度函数积分即可得到寿命分布函数：

$F_2\left(t\right)={\∫}_{-\mathrm{\∞}}^{t}f\left(\mathrm{\τ}\right)d\mathrm{\τ}={\∫}_{-\mathrm{\∞}}^t{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}d\mathrm{\τ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；τ是对t积分的时间变量。分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差；

步骤(2.10)、根据可靠度函数的性质可知，F₂(t)+R₂(t)＝1，则得基于轴温的动量轮剩余寿命预测函数R₂(t)为：

$R_2\left(t\right)=1-F_2\left(t\right)={\∫}_{-\mathrm{\∞}}^t{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}d\mathrm{\τ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；τ是对t积分的时间变量。分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差。

步骤3具体按照以下步骤实施：

步骤(3.1)、选取二元Copula函数形式为Clayton：

$C_{\mathrm{\γ}}(u,v)={(u^{-\mathrm{\γ}}+v^{-\mathrm{\γ}}-1)}^{-\frac{1}{\mathrm{\γ}}},$

其中u＝F₁(t)，v＝F₂(t)，分别表示基于电流和轴温的单变量退化数据对应的动量轮寿命分布函数，γ为未知参数，由所述步骤(1.2)和所述步骤(2.9)得F₁(t)，F₂(t)；

步骤(3.2)、通过对步骤(3.1)中Copula函数求偏导数进而求得Copula分布函数的密度函数：

$c(u,v;\mathrm{\α})=\frac{{\∂}^{2}C(u,v;\mathrm{\α})}{\∂u\∂v}$

步骤(3.3)、使用规范化的极大似然法求解步骤(3.2)的Copula分布函数的密度函数，得到未知参数γ的估计值

步骤(3.4)、根据联合概率分布计算公式，得到使用Copula分布函数的动量轮剩余寿命概率函数R(t)如下：

R(t)＝P(y₁(t)≤D₁,y₂(t)≤D₂)

＝1-P(y₁(t)＞D₁)-P(y₂(t)＞D₂)+P(y₁(t)＞D₁,y₂(t)＞D₂)

＝1-F₁(t)-F₂(t)+C_γ(u,v)

＝R₁(t)+R₂(t)+C_γ(u,v)-1

其中，u＝F₁(t)，v＝F₂(t)分别表示基于电流和轴温的单变量退化数据对应的动量轮寿命分布函数，D₁，D₂分别表示动量轮失效的电流和轴承温度阈值，R₁(t)＝1-F₁(t)，R₂(t)＝1-F₂(t)分别为相对应的单退化量动量轮剩余寿命概率分布函数，Copula分布函数表示为C_γ(u,v)。

本发明的有益效果是，一种融合电流和轴温的动量轮剩余寿命预测方法，根据动量轮中电流数据的近似线性特点，确定线性退化轨迹函数，建立退化轨迹与剩余寿命模型之间的函数关系，得到基于电流的动量轮的寿命预测曲线和基于轴温的动量轮寿命预测曲线。鉴于两个退化量相关性未知，本方法采用二元Copula函数联合两个边缘分布得到动量轮的剩余寿命分布，融合两种不同的退化量，从数据相关的角度，得到了温度和电流两个边缘分布的相关性。

附图说明

图1是本发明一种融合电流和轴温的动量轮剩余寿命预测方法中基于电流的动量轮寿命估计流程图；

图2是本发明一种融合电流和轴温的动量轮剩余寿命预测方法中基于轴温的动量轮寿命估计流程图；

图3是本发明一种融合电流和轴温的动量轮剩余寿命预测方法中采用Copula函数融合电流和轴温的动量轮寿命估计流程图。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

本发明一种融合电流和轴温的动量轮剩余寿命预测方法，具体按照以下步骤实施：

步骤1、建立基于电流退化数据的失效模型，获得基于电流的动量轮剩余寿命预测函数，如图1所示，具体按照以下步骤实施：

定义基于电流单退化量下动量轮的性能退化函数为y₁(t)，由于电流的退化数据呈现近似线性，所以采用线性退化模型中的单调回归模型来描述退化过程。在估计电流模型相关参数的过程中，已经存在误差及不确定性，因此不再选择测量误差。

步骤(1.1)、首先建立基于电流退化数据的失效模型y₁(t)＝α₁+β₁t，参数α₁、β₁分别代表电流的初始值及耗损率，y₁(t)表示电流单退化量下动量轮的性能退化函数，t为动量轮运行时间变量，其中参数α₁、β₁之间相互独立，和分别表示α₁的均值和方差，同理和分别表示β₁的均值和方差，因此y₁(t)的均值估计值和方差估计值表示为μ_α1，μ_β1，根据多个动量轮样本的退化量在不同时刻的观测值使用线性回归的方法估计得到，最终得到基于电流退化数据的失效模型其中，和分别表示y₁(t)的均值和方差估计值；

步骤(1.2)、根据步骤(1.1)得到的基于电流退化数据的失效模型求解基于电流的动量轮寿命失效分布函数F₁(t)：

$F\left(t\right)=P\left(y_1\right(t)>D_1)=\frac{1}{\sqrt{2{\mathrm{\πD}}_1}}{\∫}_{-\mathrm{\∞}}^{\stackrel{\‾}{x}}{e}^{-\frac{x^2}{2}}dx$

上式中，x为数学意义上的广义积分项，D₁为电流失效阈值，P(y₁(t)＞D₁)表示电流大于失效阈值的概率，表示时失效分布函数积分上限，中μ₁(t)，σ₁(t)分别为所述步骤(1.1)中的估计值

步骤(1.3)、根据可靠度函数的性质可知，F₁(t)+R₁(t)＝1，则可获得基于电流的动量轮剩余寿命预测函数R₁(t)：

$R_1\left(t\right)=1-F_1\left(t\right)=1-\frac{1}{\sqrt{2\mathrm{\π}}{D}_1}{\∫}_{-\mathrm{\∞}}^{\stackrel{\‾}{x}}{e}^{-\frac{x^2}{2}}dx$

式中F₁(t)表示基于电流的动量轮寿命失效分布函数；

步骤2、建立基于轴温退化数据的失效模型，得到基于轴温的动量轮剩余寿命预测函数，如图2所示，具体按照以下步骤实施：

步骤(2.1)、动量轮轴承温度随时间有缓慢的波动上升趋势，轴承温度随时间上升反映了动量轮性能的退化。然而，由于太阳入射角、功耗、噪声干扰等某些随机不确定因素的影响，遥测数据在上升的同时也存在一定的波动。鉴于轴承温度的特点，使用有漂移的维纳过程来对基于轴温的动量轮寿命退化模型进行建模，利用线性函数描述其漂移特性：

y₂(t)＝a₀+a₁t+σ_wW(t)

上式中，a₀为轴承温度遥测数据的初值，a₁为轴承温度遥测值漂移参数，W(t)为标准的维纳过程，σ_w为扩散参数；

步骤(2.2)、以相同时间间隔在动量轮轴承温度遥测数据中取n₁+1个点：其中根据所述步骤(2.1)中的模型把相邻两个采样点数据代入，相减得到Δy_i＝a₁Δt_i+σ_wΔW(t_i)，Δt_i＝(t_i-t_i-1)，ΔW(t_i)＝[W(t_i)-W(t_i-1)]，Δy_i＝y_i-y_i-1，由维纳过程定义知，

步骤(2.3)、由维纳过程平稳独立增量可得步骤(2.2)中的 的联合概率密度，使用极大似然法(Maximum Likelihood,ML)估计未知参数，构造样本似然函数L(a₁,σ_w)：

$L(a_1,{\mathrm{\σ}}_w)=f({\mathrm{\Δy}}_1,{\mathrm{\Δy}}_2,...,{\mathrm{\Δy}}_{n_1})=f\left({\mathrm{\Δy}}_1\right)f\left({\mathrm{\Δy}}_2\right)...f\left({\mathrm{\Δy}}_{n_1}\right);$

步骤(2.4)、求解步骤(2.3)中构造的似然函数，得a₁的极大似然估计值和σ_w的极大似然估计值如下：

${\hat{a}}_1=\frac{1}{n_1}\underset{i=1}{\overset{n_1}{\Σ}}\frac{{\mathrm{\Δy}}_i}{{\mathrm{\Δt}}_i},$

${\widehat{\mathrm{\σ}}}_w={\[\frac{1}{n_1}\underset{i=1}{\overset{n_1}{\Σ}}\frac{{({\mathrm{\Δy}}_i-{\hat{a}}_1{\mathrm{\Δt}}_i)}^2}{{\mathrm{\Δt}}_i}\]}^{\frac{1}{2}};$

步骤(2.5)、设动量轮样本数为n₂，各动量轮轴承温度初值利用步骤(2.4)的方法得到各动量轮的漂移参数扩散参数根据估计得到其相应的极大似然估计值：

${\widehat{\mathrm{\μ}}}_{a_0}=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[a_0\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_0}\]}^2,$

${\widehat{\mathrm{\σ}}}_{a_0}^2=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[a_0\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_0}\]}^2,$

${\widehat{\mathrm{\μ}}}_{a_1}=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\hat{a}}_1\left(j\right),$

${\widehat{\mathrm{\σ}}}_{a_1}^2=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[{\hat{a}}_1\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_1}\]}^2;$

分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差。分别为的估计值。

步骤(2.6)、由于物理失效和性能失效都发生在性能首次超过失效阈值的时间，即一旦产品的劣化过程Y(t)达到某个特定的水平ζ时产品失效，ζ为失效阈值，因此，定义产品性能退化的样本轨迹模型Y(t)为超过ζ的函数，数学描述如下：

T＝inf{Y(t)≥ζ}＝{t|Y(t)≥ζ,Y(s)≤ζ,0≤s≤t}

式中，Y(t)表示产品的劣化过程退化变量函数模型，ζ为失效阈值，t表示时间，s为小于t时刻的某一时间，动量轮轴承温度性能退化过程Y(t)＝y₂(t)是维纳过程，y₂(t)为步骤(2.1)中温度特性描述函数，给定a₀，a₁，σ_w值时，动量轮的寿命分布概率密度函数，首次达到失效阈值的时间为逆高斯分布，数学描述为：

$f\left(t\right|a_0,a_1,\mathrm{\ζ})=\frac{\mathrm{\ζ}}{\sqrt{2{\mathrm{\π\σ}}_w^2{t}^3}}{e}^{\frac{{(\mathrm{\ζ}-a_0-a_{1}t)}^2}{2{\mathrm{\σ}}_w^{2}t}}$

式中a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数。ζ为失效阈值，t表示时间；

步骤(2.8)、设动量轮的失效阈值ζ近似服从u_ζ＝75℃，σ_ζ＝3℃的正态分布，即ζ～N(75,3²)，则在考虑随机失效阈值的情况下，动量轮的寿命分布概率密度函数为：

$f\left(t\right|\mathrm{\ζ})={\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{2{\mathrm{\π\σ}}_w^2{t}^2}}{e}^{-\frac{{(\mathrm{\ζ}-a_0-a_{1}t)}^2}{2{\mathrm{\σ}}_w^{2}t}}\·\frac{1}{\sqrt{2{\mathrm{\π\σ}}_{a_0}^2}}{e}^{-\frac{{(a_0-{\mathrm{\μ}}_{a_0}t)}^2}{2{\mathrm{\σ}}_{a_0}^{2}t}}\·\frac{1}{\sqrt{2{\mathrm{\π\σ}}_{a_1}^2}}{e}^{-\frac{{(a_1-{\mathrm{\μ}}_{a_1}t)}^2}{2{\mathrm{\σ}}_{a_1}^{2}t}}{\mathrm{da}}_0{\mathrm{da}}_1$

对上式积分计算，化简整理得：

$f\left(t\right)={\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差；

步骤(2.9)、实际上由动量轮的失效机理分析可知，一般认为当动量轮轴温超过了一定的水平则轴承就不能正常工作。本方法采用轴温失效阈值ζ＝75℃。但受动量轮设计、生产误差、工作环境、外界干扰、功耗不恒定等众多因素影响，使用随机失效阈值更加符合动量轮真实的情况。参考动量轮设计要求和工程经验，动量轮失效阈值通常在75℃左右波动，实际上其失效值95％以上都是在70～80℃之间。基于上述事实，设动量轮的失效阈值近似服从u_ζ＝75℃，σ_ζ＝3℃的正态分布，即ζ～N(75,3²)。由此可得在考虑随机失效阈值的情况下，动量轮的寿命分布概率密度函数为

$\begin{array}{c}f\left(t\right)={\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}N\left(\mathrm{\ζ}\right|{\mathrm{\μ}}_{\mathrm{\ζ}},{\mathrm{\σ}}_{\mathrm{\ζ}}^2)f\left(t\right|\mathrm{\ζ})d\mathrm{\ζ}\\ ={\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}\end{array}$

对上面得到的动量轮的寿命分布概率密度函数积分即可得到寿命分布函数：

$F_2\left(t\right)={\∫}_{-\mathrm{\∞}}^{t}f\left(\mathrm{\τ}\right)d\mathrm{\τ}={\∫}_{-\mathrm{\∞}}^t{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}d\mathrm{\τ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；τ是对t积分的时间变量。分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差；

步骤(2.10)、根据可靠度函数的性质可知，F₂(t)+R₂(t)＝1，则得基于轴温的动量轮剩余寿命预测函数R₂(t)为：

$R_2\left(t\right)=1-F_2\left(t\right)={\∫}_{-\mathrm{\∞}}^t{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}d\mathrm{\τ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；τ是对t积分的时间变量。分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差；

步骤3、由概率论与数理统计理论可知，如果设备的各退化量之间相互独立，只要令设备在单退化量基础上计算得到的边缘分布函数相乘，即可得到其剩余寿命联合分布函数。若设备的各个单退化量相关性已知，即相关系数及协方差已知，可以通过数理统计中的协方差矩阵得到设备剩余寿命的联合分布。然而，由于目前动量轮失效机理模型尚不明确，很难从中提取某一个退化量作为表征动量轮寿命的决定性因素。

卫星动量轮中电流与轴温之间存在相关性，但是相关系数未知。因此，采用Copula函数融合电流和轴温的退化寿命，如图3所示，具体按照以下步骤实施：

步骤(3.1)、选取二元Copula函数形式为Clayton：

$C_{\mathrm{\γ}}(u,v)={(u^{-\mathrm{\γ}}+v^{-\mathrm{\γ}}-1)}^{-\frac{1}{\mathrm{\γ}}},$

其中u＝F₁(t)，v＝F₂(t)，分别表示基于电流和轴温的单变量退化数据对应的动量轮寿命分布函数，γ为未知参数，由所述步骤(1.2)和所述步骤(2.9)得F₁(t)，F₂(t)；

步骤(3.2)、通过对步骤(3.1)中Copula函数求偏导数进而求得Copula分布函数的密度函数：

$c(u,v;\mathrm{\α})=\frac{{\∂}^{2}C(u,v;\mathrm{\α})}{\∂u\∂v}$

步骤(3.3)、使用规范化的极大似然法(Canonical Maximum Likelihood method,CML)求解步骤(3.2)的Copula分布函数的密度函数，得到未知参数γ的估计值

步骤(3.4)、根据联合概率分布计算公式，得到使用Copula分布函数的动量轮剩余寿命概率函数R(t)如下：

R(t)＝P(y₁(t)≤D₁,y₂(t)≤D₂)

＝1-P(y₁(t)＞D₁)-P(y₂(t)＞D₂)+P(y₁(t)＞D₁,y₂(t)＞D₂)

＝1-F₁(t)-F₂(t)+C_γ(u,v)

＝R₁(t)+R₂(t)+C_γ(u,v)-1

其中，u＝F₁(t)，v＝F₂(t)分别表示基于电流和轴温的单变量退化数据对应的动量轮寿命分布函数，D₁，D₂分别表示动量轮失效的电流和轴承温度阈值，R₁(t)＝1-F₁(t)，R₂(t)＝1-F₂(t)分别为相对应的单退化量动量轮剩余寿命概率分布函数，Copula分布函数表示为C_γ(u,v)。

本发明一种融合电流和轴温的动量轮剩余寿命预测方法，针对动量轮的长寿命、价格高昂、样本数量少、可用测量数据少等特点，并且考虑到动量轮测量存在噪声、设计误差，工艺偏差，运行环境不恒定等不确定性因素，采用融合基于电流的单变量退化模型和基于轴承温度的单变量退化模型的动量轮寿命预测方法，克服了基于单变量退化模型的寿命预测数据利用不充分的缺陷，在一定程度上提高了模型的可信度，为工程上提供了一种新型寿命预测方法。

Claims (4)

1.一种融合电流和轴温的动量轮剩余寿命预测方法，其特征在于，具体按照以下步骤实施：

步骤1、建立基于电流退化数据的失效模型，获得基于电流的动量轮剩余寿命预测函数；

步骤2、建立基于轴温退化数据的失效模型，得到基于轴温的动量轮剩余寿命预测函数；

步骤3、采用Copula函数融合电流和轴温的退化寿命。

2.根据权利要求1所述的一种融合电流和轴温的动量轮剩余寿命预测方法，其特征在于，所述步骤1具体按照以下步骤实施：

步骤(1.1)、首先建立基于电流退化数据的失效模型y₁(t)＝α₁+β₁t，参数α₁、β₁分别代表电流的初始值及耗损率，y₁(t)表示电流单退化量下动量轮的性能退化函数，t为动量轮运行时间变量，其中参数α₁、β₁之间相互独立，和分别表示α₁的均值和方差，同理和分别表示β₁的均值和方差，因此y₁(t)的均值估计值和方差估计值表示为μ_α1，μ_β1，根据多个动量轮样本的退化量在不同时刻的观测值使用线性回归的方法估计得到，最终得到基于电流退化数据的失效模型其中，和分别表示y₁(t)的均值和方差的估计值；

步骤(1.2)、根据所述步骤(1.1)得到的基于电流退化数据的失效模型求解基于电流的动量轮寿命失效分布函数F₁(t)：

$F_1\left(t\right)=P\left(y_1\right(t)>D_1)=\frac{1}{\sqrt{2{\mathrm{\πD}}_1}}{\∫}_{-\mathrm{\∞}}^{\stackrel{\‾}{x}}{e}^{-\frac{x^2}{2}}dx$

上式中，x为数学意义上的广义积分项，D₁为电流失效阈值，P(y₁(t)＞D₁)表示电流大于失效阈值的概率，表示时失效分布函数积分上限，中μ₁(t)，σ₁(t)分别为所述步骤(1.1)中的估计值

步骤(1.3)、根据可靠度函数的性质可知，F₁(t)+R₁(t)＝1，则可获得基于电流的动量轮剩余寿命预测函数R₁(t)：

$R_1\left(t\right)=1-F_1\left(t\right)=1-\frac{1}{\sqrt{2{\mathrm{\πD}}_1}}{\∫}_{-\mathrm{\∞}}^{\stackrel{\‾}{x}}{e}^{-\frac{x^2}{2}}dx$

式中F₁(t)表示基于电流的动量轮寿命失效分布函数。

3.根据权利要求1所述的一种融合电流和轴温的动量轮剩余寿命预测方法，其特征在于，所述步骤2具体按照以下步骤实施：

步骤(2.1)、鉴于轴承温度的特点，使用有漂移的维纳过程来对基于轴温的动量轮寿命退化模型进行建模，利用线性函数描述其漂移特性：

y₂(t)＝a₀+a₁t+σ_wW(t)

上式中，a₀为轴承温度遥测数据的初值，a₁为轴承温度遥测值漂移参数，W(t)为标准的维纳过程，σ_w为扩散参数；

步骤(2.2)、以相同时间间隔在动量轮轴承温度遥测数据中取n₁+1个点：(t₀,y₀)，(t₁,y₁)，…，其中根据所述步骤(2.1)中的模型把相邻两个采样点数据代入，相减得到Δy_i＝a₁Δt_i+σ_wΔW(t_i)，Δt_i＝(t_i-t_i-1)，ΔW(t_i)＝[W(t_i)-W(t_i-1)]，Δy_i＝y_i-y_i-1，由维纳过程定义知，i＝1,2,…n₁；

步骤(2.3)、由维纳过程平稳独立增量可得所述步骤(2.2)中的Δy₁，Δy₂，…，的联合概率密度，使用极大似然法估计未知参数，构造样本似然函数L(a₁,σ_w)：

$L(a_1,{\mathrm{\σ}}_w)=f({\mathrm{\Δy}}_1,{\mathrm{\Δy}}_2,...,{\mathrm{\Δy}}_{n_1})=f\left({\mathrm{\Δy}}_1\right)f\left({\mathrm{\Δy}}_2\right)...f\left({\mathrm{\Δy}}_{n_1}\right);$

步骤(2.4)、求解所述步骤(2.3)中构造的似然函数，得a₁的极大似然估计值和σ_w的极大似然估计值如下：

${\hat{a}}_1=\frac{1}{n_1}\underset{i=1}{\overset{n_1}{\Σ}}\frac{{\mathrm{\Δy}}_i}{{\mathrm{\Δt}}_i},$

${\widehat{\mathrm{\σ}}}_w={\[\frac{1}{n_1}\underset{i=1}{\overset{n_1}{\Σ}}\frac{{({\mathrm{\Δy}}_i-{\hat{a}}_1{\mathrm{\Δt}}_i)}^2}{{\mathrm{\Δt}}_i}\]}^{\frac{1}{2}};$

步骤(2.5)、设动量轮样本数为n₂，各动量轮轴承温度初值利用所述步骤(2.4)的方法得到各动量轮的漂移参数扩散参数根据 估计得到其相应的极大似然估计值：

${\widehat{\mathrm{\μ}}}_{a_0}=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[a_0\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_0}\]}^2,$

${\widehat{\mathrm{\σ}}}_{a_0}^2=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[a_0\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_0}\]}^2,$

${\widehat{\mathrm{\μ}}}_{a_1}=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\hat{a}}_1\left(j\right),$

${\widehat{\mathrm{\σ}}}_{a_1}^2=\frac{1}{n_2}\underset{j=1}{\overset{n_2}{\Σ}}{\[{\hat{a}}_1\left(j\right)-{\widehat{\mathrm{\μ}}}_{a_1}\]}^2;$

分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差。分别为的估计值。

步骤(2.6)、由于物理失效和性能失效都发生在性能首次超过失效阈值的时间，即一旦产品的劣化过程Y(t)达到某个特定的水平ζ时产品失效，ζ为失效阈值，因此，定义产品性能退化的样本轨迹模型Y(t)为超过ζ的函数，数学描述如下：

T＝inf{Y(t)≥ζ}＝{t|Y(t)≥ζ,Y(s)≤ζ,0≤s≤t}

式中，Y(t)表示产品的劣化过程退化变量函数模型，ζ为失效阈值，t表示时间，s为小于t时刻的某一时间，动量轮轴承温度性能退化过程Y(t)＝y₂(t)是维纳过程，y₂(t)为步骤(2.1)中温度特性描述函数，给定a₀，a₁，σ_w值时，动量轮的寿命分布概率密度函数，首次达到失效阈值的时间为逆高斯分布，数学描述为：

$f\left(t\right|a_0,a_1,\mathrm{\ζ})=\frac{\mathrm{\ζ}}{\sqrt{2{\mathrm{\π\σ}}_w^2{t}^3}}{e}^{-\frac{{(\mathrm{\ζ}-a_0-a_{1}t)}^2}{2{\mathrm{\σ}}_w^{2}t}}$

式中a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；ζ为失效阈值，t表示时间；

步骤(2.8)、设动量轮的失效阈值ζ近似服从u_ζ＝75℃，σ_ζ＝3℃的正态分布，即ζ～N(75,3²)，则在考虑随机失效阈值的情况下，动量轮的寿命分布概率密度函数为：

$f\left(t\right)={\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_1}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差；

步骤(2.9)、对所述步骤(2.8)得到的动量轮的寿命分布概率密度函数积分即可得到寿命分布函数：

$F_2\left(t\right)={\∫}_{-\mathrm{\∞}}^{t}f\left(\mathrm{\τ}\right)d\mathrm{\τ}={\∫}_{-\mathrm{\∞}}^t{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_1}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}d\mathrm{\τ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；τ是对t积分的时间变量。分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差；

步骤(2.10)、根据可靠度函数的性质可知，F₂(t)+R₂(t)＝1，则得基于轴温的动量轮剩余寿命预测函数R₂(t)为：

$R_2\left(t\right)=1-F_2\left(t\right)=1-{\∫}_{-\mathrm{\∞}}^t{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{\mathrm{\ζ}}{\sqrt{4{\mathrm{\π}}^2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_0}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)t^2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}{e}^{-\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{a_0}-{\mathrm{\μ}}_{a_1}t)}^2}{2({\mathrm{\σ}}_w^{2}t+{\mathrm{\σ}}_{a_1}^2+{\mathrm{\σ}}_{a_1}^2{t}^2)}+\frac{{(\mathrm{\ζ}-{\mathrm{\μ}}_{\mathrm{\ζ}})}^2}{2{\mathrm{\σ}}_{\mathrm{\ζ}}^2}}d\mathrm{\ζ}d\mathrm{\τ}$

u_ζ和σ_ζ分别表示失效阈值的均值和标准差，ζ为失效阈值，t表示时间，a₀，a₁，σ_w分别为步骤(2.1)定义的轴承温度遥测数据的初值，轴承温度遥测值漂移参数，扩散参数；τ是对t积分的时间变量。分别为轴承温度遥测数据的初值a₀和轴承温度遥测值漂移参数a₁的均值和方差。

4.根据权利要求1所述的一种融合电流和轴温的动量轮剩余寿命预测方法，其特征在于，所述步骤3具体按照以下步骤实施：

步骤(3.1)、选取二元Copula函数形式为Clayton：

$C_{\mathrm{\γ}}(u,v)={(u^{-\mathrm{\γ}}+v^{-\mathrm{\γ}}-1)}^{-\frac{1}{\mathrm{\γ}}},$

其中u＝F₁(t)，v＝F₂(t)，分别表示基于电流和轴温的单变量退化数据对应的动量轮寿命分布函数，γ为未知参数，由所述步骤(1.2)和所述步骤(2.9)得F₁(t)，F₂(t)；

步骤(3.2)、通过对所述步骤(3.1)中Copula函数求偏导数进而求得Copula分布函数的密度函数：

$c(u,v;\mathrm{\α})=\frac{{\∂}^{2}C(u,v;\mathrm{\α})}{\∂u\∂v}$

步骤(3.3)、使用规范化的极大似然法求解所述步骤(3.2)的Copula分布函数的密度函数，得到未知参数γ的估计值

步骤(3.4)、根据联合概率分布计算公式，得到使用Copula分布函数的动量轮剩余寿命概率函数R(t)如下：

R(t)＝P(y₁(t)≤D₁,y₂(t)≤D₂)

＝1-P(y₁(t)＞D₁)-P(y₂(t)＞D₂)+P(y₁(t)＞D₁，y₂(t)＞D₂)

＝1-F₁(t)-F₂(t)+C_γ(u,v)

＝R₁(t)+R₂(t)+C_γ(u,v)-1

其中，u＝F₁(t)，v＝F₂(t)，分别表示基于电流和轴温的单变量退化数据对应的动量轮寿命分布函数，D₁，D₂分别表示动量轮失效的电流和轴承温度阈值，R₁(t)＝1-F₁(t)，R₂(t)＝1-F₂(t)分别为相对应单退化量的动量轮剩余寿命概率分布函数，Copula分布函数表示为C_γ(u,v)。