本发明属于可靠性工程领域，涉及一种基于逆高斯退化模型的剩余寿命贝叶斯预测方法。
背景技术：
当代社会对产品的可靠性要求越来越高，产品的质量监控与健康管理也朝着精细化方向发展，因此，如何提前准确预测出产品的剩余寿命成为了目前的一个研究热点。逆高斯退化模型具有优良的统计特性，适合对单调性能退化过程建模，已经广泛用于退化失效型产品的剩余寿命预测。为了有效融合多源退化数据提高剩余寿命的预测准确度和可信度，发明一种基于逆高斯退化模型的剩余寿命贝叶斯预测方法。
技术实现要素：
发明一种基于逆高斯退化模型的剩余寿命贝叶斯预测方法，该方法的具体技术方案为：
步骤一：建立逆高斯退化模型与剩余寿命预测模型
设产品的性能退化过程{Y(t),t≥0}服从逆高斯退化模型，则Y(t)～IG(μΛ(t),λΛ(t)2)，Y(t)的概率密度函数表示为
式中，μ为均值参数，λ为尺度参数，Λ(t)＝tΛ为时间函数。设产品性能指标的失效阈值为D，则产品的寿命ξ为Y(t)首次到达D的时间，ξ的累积分布函数为
式中，Φ(·)为标准Normal分布的累积分布函数。产品的剩余寿命L(t)是指从时刻t的性能退化量Y(t)到首次超过D的时间，可以表示为L(t)＝inf{x|Y(t+x)≥D,x>0}，则产品寿命ξ与剩余寿命L(t)之间的关系为ξ＝L(t)+t。
步骤二：设逆高斯退化模型的尺度参数与均值参数为随机参数，构建随机参数的共轭先验分布函数
为了有效融合先验退化数据有现场退化数据，提高剩余寿命预测值的准确度与可信度，设μ,λ作为随机参数。并且为了便于统计分析，采用μ,λ的共轭先验分布：设λ服从Gamma分布λ～Ga(a,b)，其概率密度函数表示为
设δ＝1/μ服从条件正态分布δ|λ～N(c,d/λ)，其概率密度函数为
其中，a,b,c,d是随机参数的超参数。
步骤三：根据Bayes理论推导随机参数的后验分布函数
设yi,j为第i个产品的第j次性能退化测量值，△yi,j为退化增量，△Λi,j为时间增量，根据逆Gaussian过程的统计特性建立似然函数为
建立完全对数似然函数为
设f(δi,λi)是随机参数δi,λi的联合先验概率密度函数，则有f(δi,λi)＝f(λi)·f(δi|λi)，通过Bayesian公式f(δi,λi|△yi)∝L(δi,λi)·f(δi,λi)推导出联合后验概率密度函数f(δi,λi|△yi)，得到随机参数δi,λi的后验分布为
其中，
步骤四：设计EM算法估计随机参数先验分布函数的超参数值
设计EM算法估计超参数Ω＝(a,b,c,d)，EM算法每一步迭代包含E步和M步。E步的任务是求取隐含数据项的期望值，设Ω(l)为第l次迭代后的估计值向量，则在第l+1次迭代中，隐含数据项λi,lnλi,λiδi,的期望值为
式中，ψ(·)为digamma函数。将完全似然函数式(1)中的各隐含数据项利用对应的期望值代替后，M步的任务是极大化式(1)，解得c(l+1)，d(l+1)，b(l+1)及a(l+1)的表达式如下
式中，ψ-1(·)为逆digamma函数。
EM算法的执行过程为：
初始化:设l＝0,Ω(0)＝(1,1,1,1)；
第l+1次迭代:
E步:计算E(λi|yi,Ω(l)),E(lnλi|yi,Ω(l)),E(λiδi|yi,Ω(l))及
M步:解得c(l+1)，d(l+1)，b(l+1)及a(l+1)，将Ω(l)更新为Ω(l+1)；
结束条件:max(Ω(l+1)-Ω(l))<10-3或l达到最大迭代数。
步骤五：估计随机参数的后验期望值
设产品在时间tj的现场性能退化数据为△xj，j＝1,2,…,k，首先将随机参数在时间tj的后验分布函数更新为
其中，
然后，根据Gamma函数与Normal函数的统计特性估计随机参数的后验期望值为
步骤六：预测产品的剩余寿命
利用随机参数的后验期望值更新累积失效分布函数
利用Fξ(t|△x)计算出E(ξ|△x)，产品在t时刻的剩余寿命预测值为L(t)＝E(ξ|△x)-t。
附图说明
图1基于逆高斯退化模型的剩余寿命贝叶斯预测方法。
具体实施方式
下面结合附图对本发明实现步骤进行进一步说明。
实施例：电连接器的主要失效模式有机械失效，电气失效，绝缘失效三种，机械失效主要由接插件应力松弛造成，性能退化量y为接插件应力值x相对于初始应力值x0的百分比变化y＝(x-x0)/x0×％，每个样品在0时刻的性能退化量为0，失效阈值为D＝30％。电连接器的6组先验性能退化数据及测量时间如表1所示，1组现场性能退化数据如表2所示。
表1电连接器的先验性能退化数据
表2电连接器的现场性能退化数据
步骤一：建立逆高斯退化模型与剩余寿命预测模型，产品性能指标的失效阈值为D＝30％，则产品的寿命ξ的累积分布函数为
产品的剩余寿命L(t)是指从时刻t的性能退化量Y(t)到首次超过D的时间L(t)＝inf{x|Y(t+x)≥0.3,x>0}。
步骤二：设逆高斯退化模型的尺度参数与均值参数为随机参数，构建随机参数的共轭先验分布函数。设μ,λ作为随机参数。并且为了便于统计分析，采用μ,λ的共轭先验分布：设λ服从Gamma分布λ～Ga(a,b)，其概率密度函数表示为
设δ＝1/μ服从条件正态分布δ|λ～N(c,d/λ)，其概率密度函数为
其中，a,b,c,d是随机参数的超参数。
步骤三：根据Bayes理论推导随机参数的后验分布函数。设yi,j为第i个产品的第j次性能退化测量值，△yi，j为退化增量，△Λi，j为时间增量，根据逆Gaussian过程的统计特性建立似然函数为
建立完全对数似然函数为
设f(δi,λi)是随机参数δi,λi的联合先验概率密度函数，则有f(δi,λi)＝f(λi)·f(δi|λi)，通过Bayesian公式f(δi,λi|△yi)∝L(δi,λi)·f(δi,λi)推导出联合后验概率密度函数f(δi,λi|△yi)，得到随机参数δi,λi的后验分布为
其中，
步骤四：设计EM算法估计随机参数先验分布函数的超参数值，EM算法的执行过程为：
初始化:设l＝0,Ω(0)＝(1,1,1,1)；
第l+1次迭代:
E步:计算E(λi|yi,Ω(l)),E(lnλi|yi,Ω(l)),E(λiδi|yi,Ω(l))及
M步:解得c(l+1)，d(l+1)，b(l+1)及a(l+1)，将Ω(l)更新为Ω(l+1)；
结束条件:max(Ω(l+1)-Ω(l))<10-3或l达到最大迭代数。
解得
步骤五：估计随机参数的后验期望值
利用产品在10个测量时刻的现场性能退化数据，结合先验参数估计值分别对超参数的后验估计值进行更新，更新结果如表3所示
表3超参数值更新结果
根据Gamma函数与Normal函数的统计特性估计出，产品在时刻1800h的随机参数后验期望值为E(λ|△x)＝0.647，E(μ|△x)＝0.601。
步骤六：预测产品的剩余寿命
利用随机参数的后验期望值更新累积失效分布函数
Fξ(t|△x)＝Φ(0.147(t0.509-49.917))-exp(2.153t0.509)Φ(-0.147((t0.509+49.917)))，
进而计算出E(ξ|△x)＝2169.906h，产品在1800h的剩余寿命预测值为L(t)＝369.906h，利用Bootstrap自助抽样法建立剩余寿命预测值的90％置信度区间为[158h,623h]。
如果不利用Bayes理论融合先验性能退化数据进行剩余寿命预测，只利用现场性能退化数据预测出
Fξ(t|△x)＝Φ(0.153(t0.482-61.004))-exp(2.869t0.482)Φ(-0.153((t0.482+61.004)))，
剩余寿命值为3259h，原大于实际剩余寿命560h；利用Bootstrap自助抽样法建立剩余寿命预测值的90％置信度区间为[626h,6755h]，区间长度远大于发明方法所建立的置信区间长度。可见，发明方法提高了剩余寿命预测值的准确度与置信度。