[发明专利]滚动轴承集成期望最大化和粒子滤波的寿命预测方法

摘要

滚动轴承集成期望最大化和粒子滤波的寿命预测模型，首先采用峭度指标对轴承健康状态进行实时监测，确定寿命预测起始时刻；当满足预测起始条件后，采用有效值对轴承剩余寿命进行预测；在预测阶段，采用期望最大化方法对模型参数进行评估，同时采用粒子滤波方法对轴承状态进行评估，通过对模型参数和轴承状态的准确评估，提高剩余寿命预测精度，本发明能够实现对模型参数和滚动轴承状态的准确评估，并且在滚动轴承寿命预测中表现出了比传统指数模型更好的预测效果。

主权项

滚动轴承集成期望最大化和粒子滤波的寿命预测方法，其特征在于，包括以下步骤：1)实时监测并采集滚动轴承振动信号，计算滚动轴承峭度指标和有效值；2)计算滚动轴承健康时刻峭度指标的均值μ和标准差v，以确定健康状态下峭度指标的3v区间[μ‑3v,μ+3v]；3)判断滚动轴承连续l+1个时刻的峭度指标{mp+i}i＝0：l是否满足{|mp+i‑μ|＞3v}i＝0:l，如果满足该条件，则确定mp所对应的时刻tp为寿命预测的起始时刻；4)从寿命预测起始时刻开始，将滚动轴承振动信号有效值代入衰退模型：其中，xi为ti时刻状态值，是已知常数，θ′，β′和σ为三个未知参数，σB(ti)～N(0,σ2ti)服从布朗运动，对上式求对数得到以下变形形式：其中，θ＝ln(θ′)服从分布β＝β′‑σ2/2服从分布且规定以简化计算；对模型参数进行初始化，得初始参数μ0＝μθ,0，μ1＝μβ,0和5)从概率密度函数中进行随机采样，得到Ns个初始粒子粒子权值为其中Δt＝tj‑tj‑1为时间间隔；6)tk时刻得到观测序列S1:k＝{s1,s2,...,sk}，在参数已知条件下，观测序列的条件概率密度为：p(S1:k|θ,β)=(12πσ2Δt)kexp[-(s1-θ-βt1)22σ2t1-Σj=2k(sj-sj-1-βΔt)22σ2Δt];---(3)]]>根据贝叶斯理论得到参数θ和β的联合概率密度函数为：p(θ,β|S1:k)∝p(S1:k|θ,β)p(θ,β)∝exp[-(s1-θ-βt1)22σ2t1-Σj=2k(sj-sj-1-βΔt)22σ2Δt]exp[-(θ-μ0)22σ02]exp[-(β-μ1)22σ12]∝12πσθ,kσβ,k1-ρk2exp[-12(1-ρk2)((θ-μθ,k)σθ,k2-2ρk(θ-μθ,k)(β-μβ,k)σθ,kσβ,k+(β-μβ,k)2σβ,k2)]---(4)]]>由此得到参数更新公式如下：μθ,k=(s1σ02+μ0σ2t1)(σ2+σ12tk)-σ02t1(skσ12+μ1σ2-0.5σ4)(σ02+σ2t1)(σ12tk+σ2)-σ02σ12t1σθ,k2=σ02σ2t1(σ2+σ12tk)(σ02+σ2t1)(σ12tk+σ2)-σ02σ12t1μβ,k=(skσ12+μ1σ2-0.5σ4)(σ02+σ2t1)-σ12(s1σ02+μ0σ2t1)(σ02+σ2t1)(σ12tk+σ2)-σ02σ12t1σβ,k2=σ12σ2t1(σ02+σ2t1)(σ02+σ2t1)(σ12tk+σ2)-σ02σ12t1ρk=-σ0σ1t1(σ02+σ2t1)(σ12tk+σ2)---(5)]]>采用以上公式对参数μθ,k,μβ,k和进行更新；7)为待评估的模型参数，计算似然函数：ln[p(S1:k,θ,β|Θk)]=ln[p(S1:k|θ,β,Θk)]+ln[p(θ,β|Θk)]=-k2ln(Δt)-k+22ln(2π)-k2ln(σk2)-(s1-θ-βt1)22σk2t1-Σj=2k(sj-sj-1-βΔt)22σk2Δt-12ln(σ0,k2)-12ln(σ1,k2)-(θ-μ0,k)22σ0,k2-(β-μ1,k)22σ1,k2---(6)]]>似然函数的期望值为：其中为第i次评估的结果，计算满足条件的参数为第i+1次评估的结果，σ^k2(i+1)=1k(s12-2s1(μθ,k+μβ,kt1)+σθ,k2+σβ,k2+2t1(ρkσθ,kσβ,k+μθ,kμβ,k)+t12(μβ,k2+σβ,k2)t1+Σj=2k(sj-sj-1)2-(sj-sj-1)Δtμβ,k+(Δt)2(μβ,k2+σβ,k2)Δt)μ^0,k(i+1)=μθ,k,σ^0,k2(i+1)=σ^θ,k2,μ^1,k(i+1)=μβ,k,σ^1,k2(i+1)=σβ,k2;---(8)]]>8)将衰退模型改写为以下形式：sk＝sk‑1+βΔt+σW(Δt)， (9)其中W(Δt)＝B(tk)‑B(tk‑1)，由此得到重要密度函数为：p(sk|sk-1)=12πΔt(σβ,k2Δt+σk2)exp[-(sk-sk-1-μβ,kΔt)22Δt(σβ,k2Δt+σk2)],---(10)]]>从以上重要密度函数中进行重要性采样，得到粒子集9)采用tk时刻的观测值Sk对粒子权值进行更新，wki=wk-1i12πσk2tkexp[-(sk-ski)22σk2tk],wki=wkiΣi=1Nswki,---(11)]]>采用下式计算有效粒子数，N^eff=1Σi=1Ns(wki)2,---(12)]]>如果有效粒子数小于阈值NT，需要根据粒子权值大小进行重采样，得到新的粒子集使其满足粒子权值重置为10)采用粒子集对滚动轴承状态进行评估，s^k=Σi=1Ns(wkiski),---(13)]]>然后将状态评估结果代入下式对滚动轴承剩余寿命的概率密度函数进行预测，fLk|S1:k(lk|S1:k)=γ-s^k2πlk3(σβ,k2lk+σ2)exp[-(γ-s^k-μβ,klk)22lk(σβ,k2lk+σ2)],lk≥0,---(14)]]>其中lk为tk时刻的剩余寿命，γ为滚动轴承失效阈值。