# Remaining Useful Life Estimation of Fan Slewing Bearings in Nonlinear Wiener Process with Random Covariate Effect

Authors: Mingjun Liu; Zengshou Dong; Hui Shi; Yujia Zheng; Min Zhang
Journal: Shock and Vibration
 (2022)

---

## Abstract

Since the degradation process of fan slewing bearings is easily influenced by the external environment, it is difficult to estimate its remaining useful life (RUL) accurately. A nonlinear Wiener degradation model considering the influence of random covariate is established for the prediction of the RUL of fan slewing bearings in this study. Firstly, considering the nonlinear and nonmonotonic properties of the operation process of fan slewing bearings, the degradation model of the nonlinear Wiener process of fan slewing bearings is established. Secondly, the combination of random covariate models and the nonlinear Wiener degradation process is researched. The stress effect which is used as a random covariate is introduced into the nonlinear Wiener degradation model in the form of the additive hazard model. Moreover, the closed expression for the RUL probability density function(PDF) is derived for the random variation of drift coefficients, the individual differences and the random variation of covariates. Thirdly, the maximum likelihood estimation algorithm is used to estimate the RUL parameters depending on the historical degradation data. Finally, the vibration data of fan slewing bearings monitored by sensors are used to verify the effectiveness of the proposed method. The results show that the proposed method can be used to improve the fitting degree of the model and the accuracy of RUL estimation.


---

## Body

## 1. Introduction

In recent years, with increasing concern to environmental and energy issues, wind energy as renewable energy has attracted more and more attention. Wind power generation is an important new power source in the world. Since 2005, the average annual growth rate of wind power capacity has reached 20%. It is estimated that wind power will account for 12% of the global power supply by 2030 [1]. In Europe, about 7.7% of total electricity consumption will be generated by offshore wind turbines, and the installed capacity is 66 GW [2]. The slewing bearings is one of the key components of the fan, which is used for the pitch and yaw system. For driving the pinion and gearbox, the gear is designed on the inner ring or outer ring (see Figure 1).Figure 1
Fan structure with (i) pitch bearing and (ii) yaw bearing.Unlike general-purpose industrial bearings, the large-scale fan slewing bearings is designed to be operated under harsh conditions with a high failure rate and maintenance cost [3]. During the operation of slewing bearings, the parts are subjected to alternate stress for a long time, which will cause wear and damage. At the same time, the fan slewing bearings are with low speed and under heavy load operation conditions. It is difficult to identify the faults formed and some faults can cause malignant accidents. Therefore, it is important to estimate the RUL of slewing bearings accurately [4, 5].The widely used RUL estimation methods of slewing bearings mainly include the physical model method and the data-driven method. The physical model method is mainly based on a large number of tests and statistical analysis of engineering data. It is usually based on stress and strain, critical plane model, material structure, etc. However, the physical model method is based on a large number of tests and statistical analysis of physical fault model data. Due to the large volume, various specifications as well as the complex load of slewing bearings, it is difficult to build the test platform and the cost is also very high. These factors hinder researchers from carrying out the experimental test.For data-driven method, attempts to derive models directly from collected degradation data or life data, which does not depend on physical or engineering principles completely. Monitoring signal changes with the development of wear. Through researching the change of monitoring signal, the degradation process can be analyzed.Zhang et al. [6] introduces the mixed distribution of Gaussian distribution as a distribution describing random effects in the degradation model of the diffusion process. Its drift coefficient is a linear combination of some known time-varying functions. Finally, the RUL PDF is derived. Zhang et al. [4] proposes that the features of temperature, torque, and vibration signal of the service sample and reference sample are extracted separately. PCA-based multiple sensitive features are used to establish performance decline indicators. Then, similarity between samples and reference samples can be calculated. The RUL of the sample is predicted based on the similarity. Feng et al. [7] establishes the relationship between bearing life and maximum load through the small sample bearing parameter estimation method. The RUL model of slewing bearings is established based on the Weibull distribution method. Then, the experimental analysis is performed on QNA-730-22 slewing bearings. Lu et al. [8] fuses the characteristics of vibration and other signals through principal component analysis, including the health state of operation, azimuth, peak value, wavelet energy entropy, and inherent mode function energy. Then, particle swarm optimization is used to optimize the degradation model of the least squares supporting vector machine. Aye and Heyns [9] obtains the optimal Gaussian process regression(GPR) by combining the existing simple mean and covariance functions to obtain the irregularity in the data of bearing degradation. GPR is improved to realize the estimation of the low error rate for the RUL of the low-speed bearing.In the above research process, random and dynamic characteristics are not considered for the degradation process of fan slewing bearings. The random process can be better used to describe its degradation state [10–12]. Therefore, the random process is used to describe the degradation process of fan slewing bearings in this study.Stochastic processes mainly include Wiener, Gamma, Markov, and other methods [13]. Since fan slewing bearing is influenced by internal or external factors, it has nonmonotonic characteristics [14]. The Wiener process is suitable to describe the nonmonotonic nonlinear degradation process. Therefore, the Wiener process is established to estimate the RUL of fan slewing bearings. Considering the effect of historical degradation data on the degradation model, Si et al. [15] proposes a Wiener degradation process with a recursive model, in which, the drift coefficient is updated by the recursive filter and other parameters are updated by the expectation maximization (EM) algorithm. At the same time, the distribution of the drift coefficient is also considered. Finally, the RUL distribution with high precision is obtained.Man and Zhou [16] establishes the Wiener degradation modeling with drift and uses the nonparametric baseline hazard model to build a joint modeling framework. Then, it is used to estimate the RUL of the system. However, the effect of multiple monitoring signals on the RUL is not considered. Paroissin [17] regards the Wiener process as a degradation model that starts at any time and regards the degradation process as a random delayed Wiener process. Assuming that the sample path is observed instantaneously under the same rules, statistical judgment is made according to the sampling scheme, and some progressive results are obtained. To simulate the degradation trajectory of industrial equipment, Huang et al. [18] establishes an adaptive skew-Wiener model. Making full use of prior knowledge and historical information, an online filtering algorithm for state estimation is proposed. The two-stage algorithm is used to estimate the unknown parameters. Finally, it is applied to motor bearings. Wang et al. [19] analyzes the impact of failure on the degradation process, establishes a degradation model based on the composite process and predicts the RUL distribution of the system without measuring noise online.In these models, parametric random variables and individual differences are considered only. The effect of the external environment acting on the system degradation is ignored. Therefore, these models cannot fully reflect the process of system degradation. However, the operating environment of wind power equipment is changeable and complex, including the continuous changes of strong wind, tropical high temperature, lightning, snow, and so on [20]. In this study, the degradation effect of the external environment acting on the slewing bearing through the blade propeller is considered (such as wind speed and change of wind direction). The variety of external environments leads to changes in the speed and rotation direction of the blade. Consequently, it produces the stress effect on the slewing bearing device which is connected to the blade. When the external environment changes, the stress effect changes randomly. And, the degenerate state of slewing bearing is influenced. Then, the stress effect, which is used as a random covariate [21], is introduced into the nonlinear Wiener degradation model.Scholars have also conducted some research on the degradation model based on covariates. Sun et al. [22] proposed that the effect of an external factor on degradation, which is used as covariates, is introduced into the nonlinear Wiener process in the form of proportional hazard models. The proportional hazards model requires that the relationship between the failure rate function and covariates is proportional. However, in general, the degradation process of the device is randomly varying, and the relationship between the failure rate function and covariates is not strictly proportional. Therefore, proportional hazards models are not suitable for engineering practice. Si et al. [13] analyzed Wiener process-based degradation models with covariates. The combination of covariates, which is in the form of a proportional hazards model and diffusion coefficient, is studied. In [21], the proportional hazards models of fixed covariate models and time-varying covariate models are analyzed. The proportional hazard models and the additive hazard models is compared. However, it is necessary to make further research to the combination of random covariate models and the nonlinear Wiener degradation process.Moreover, in many literature, when the effect of an external factor on degradation is considered, it is introduced into the degradation model in form of additive hazard models. Sun et al. [23] considers the impact of random shocks on the increments and rate of a degradation process, and the impact of random shocks is taken into a nonlinear Wiener process model in form of additive model. Chen et al. [24] proposes a nonlinear adaptive inverse Gaussian process along with the corresponding state space model considering measurement errors. The measurement errors is introduced into the degradation model in form of additive model. Therefore, the stress effect of the external environment acting on the fan slewing bearing through the blade propeller, which is used as the random covariates, is introduced into the nonlinear Wiener degradation model in the form of additive hazard models in this study. And, it provides a better estimation of the RUL of the slewing bearing.Based on the above analysis, a nonlinear Wiener degradation model considering the effect of random covariate is established for the estimation of the RUL of fan slewing bearing in this study. Firstly, the nonlinear Wiener process is established to estimate the RUL of fan slewing bearing. Next, the combination of random covariate models and the nonlinear Wiener degradation process is researched in this study. The degradation effect of the external environment acting on the slewing bearing through the blade propeller is considered. To improve the accuracy of model estimation, the stress effect, which is used as a random covariate, is introduced into the nonlinear Wiener degradation model in the form of additive hazard models. Moreover, a closed expression for the RUL probability density function (PDF) is derived for the random variation of drift coefficients, the individual differences, and the random variation of covariates. Thirdly, the maximum likelihood estimation algorithm is used to estimate the parameters of the PDF depending on the historical degradation data. Finally, the vibration data of the fan slewing bearings monitored by sensors are used to verify the effectiveness of the proposed method.


## 2. Fan Slewing Bearing Degradation Model

### 2.1. Nonlinear Wiener Degradation Model

Assuming that the degradation value of the parameter at timet is Xt, the random Wiener degradation process can be expressed as follows:(1)Xt=X0+αΛt,θ+σBτt,γ,where X0 is the initial degenerate state. When X0=0, the system is in a healthy state. α is the drift coefficient. Λt,θ=∫0tλu,θdu, λt,θ, and τt,γ is the continuous nondecreasing function about time t. θ and γ are the parameter vectors. σ is the diffusion coefficient. Bτt,γ is the Brownian motion and it follows Gaussian distribution.The forms ofΛt,θ are t, tr, and expbt. Since the degradation of slewing bearings is a nonlinear process, and tris widely used because of its flexibility in describing linear degradation paths with r=1, nonlinear concave paths with r>1, and nonlinear convex paths with r<1. Therefore, it is selected as Λt,θ=tr. The nonlinear Wiener degradation process is as follows:(2)Xt=X0+αtr+σBτt,γ.Due to differences in the production materials and production processes of the system, the degradation rate in the degradation process of each system is different. Generally, the random parameterα represents different degradation rates, α∼Nμα,σα2. Suppose, r follows Gaussian distribution, N0,1.Property 1.
the properties of the Wiener process are given as follows:(1)
Increment,Xt1−X0, Xt2−Xt1,…,Xtn−Xtn−1 on different intervals is independent of each other(2)
Xti−Xti−1 follows Gaussian distribution, with expectation αΛti,θ−Λti−1,θ and variance σ2τti;γ−τti−1;γ
According to Property1, in the Wiener process, Xt is not a Gaussian process, and it is a degradation value.


### 2.2. Nonlinear Wiener Degradation Model Based on Random Covariate

To estimate the RUL of a random degenerate system accurately, in practice, random variation of parameters, individual differences, and covariates effect are mainly considered [13, 25]. Parameter random variation means that the system degradation process is usually characterized by a random process. The individual difference means that the degradation paths of different equipment in the same system are different due to differences in working environment and load. The covariate effect means that the randomness of the external environment changes strongly, and it impacts on equipment degradation process through other components [20].The degradation process of the fan slewing bearing is influenced by the external environment (such as the change of wind speed and wind direction). The variety of external environments leads to changes in the speed and rotation direction of the blade. Consequently, it produces the stress effect on the slewing bearing device that is connected to the blade. It influences the degenerate state of slewing bearing. The stress effect is a random variation.Moreover, the proportional hazards model of covariates requires that the relationship between the failure rate function and covariates is proportional. However, in this study, the degradation process of slewing bearing is a random variation, and the relationship between the failure rate function and covariates is not strictly proportional. Generally, when the effect of an external factor on degradation is considered, it is introduced in an additive model into the degradation model. Therefore, considering the stress effect of the external environment acting on the slewing bearing through the blade propeller, as a random covariate, it is introduced into the nonlinear Wiener degradation process in the form of additive hazard models in this study. Considering the individual difference of slewing bearing samples, the random variation of drift coefficient, and the effect of a random covariate, the closed expression of RUL PDF of fan slewing bearing is derived.The stress effect of the external environment acting on the degradation process is set asM, and the total degradation state of equipment is Yt. The actual degradation process of slewing bearing at time t consists of two parts: self-degradation at time t and the degradation of stress effect at time t. The degradation model of slewing bearings is as follows:(3)Yt=Xt+M.In this study, the stress effect of the external environment acting on the slewing bearing through the blade propeller is considered. When the external environment changes, the stress effect changes randomly. Then, the degenerate state of slewing bearing is influenced, and the degradation of slewing bearing itself is determined by the inherent factors of the material. Therefore, the degeneration caused by stress and the degradation of slewing bearing itself is independent to each other. Since the Gamma distribution describes a monotonically increasing change process with time and the stress effect is a random variation, it may be nonmonotonous and not suitable for the Gamma distribution. Most of the external interference is assumed to be a Gaussian distribution [26]. According to the central limit theorem, the stress effect is supposed to be a Gaussian distribution, M∼Nμm,σm2. Since, the data-driven method is used to estimate RUL, based on the feature data, expectation μm and variance σm are estimated by using the maximum likelihood function, and then the stress effect is estimated.Combining with equations (2) and (3) in this study, the nonlinear Wiener degradation model of slewing bearing considering external environmental effects is established.(4)Yt=X0+αtr+M+σBτt,γ.The PDF is given as follows. Supposing that the degradation data obtained from the condition monitoring of the product at a discrete time pointtktk>0 is Ytk. If the service life of the product is T. Lk which is the RUL at time tk, can be expressed as Lk=T−tk. The first reaching failure threshold of the product is set as D. If Yt exceeds D, the product is judged to be invalid, and the definition of Lk is given.(5)Lk=inflk;Ylk+tk≥DYtk<D.Firstly, when the drift coefficientα is given, without considering the effect of external covariates, RUL PDF of Xt at time tk is analyzed. Once the degradation state Xt,t≥0 reaches a failure threshold D, the system is considered faulty. According to the definition of first passage time(FPT), the RUL can be defined as follows:(6)T=inft:Xt≥D|X0<D.Since it is difficult to obtain an accurate expression of the failure PDF for the nonlinear Wiener degradation model, the following equation is used to approximate the PDF according to literature [27].(7)flX1:k,θk≅12πlSll+λl,θkσexp−S2l2l.AmongSl=1/σD−Λl;θ, Λl;θ=∫tkl+tkλt;θdt, λl,θk=αrtr−1, Λl;θ=αl+tkr−tkr, and Dtk=D−Xtk, Dk is D-value between the failure threshold and the degenerate state value at time tk.Hence,Sl=1/σDtk−αl+tkr+αtkr.Letθk=μtk,σtk,r′, the PDF approximate expression of RUL Lk of Xt at tk is as follows:(8)flkX1:k,θk≅12πlkDtk−αlk+tkr+αlkrlk+tkr−1+αtkrlkσ•exp−Dtk−αlk+tkr+αtkr22lkσ2.When the degradation effect of the external environment acting on the slewing bearing is considered, according to the definition of FPT, the system is considered faulty, when the degradation stateYt,t≥0 reaches threshold D. According to Property 1, ΔYtk=Ytk−Ytk−1 are independent of each other and follow Gaussian distribution.Proof.
(9)Ytk=X0+αktkr+Mk+σkBkτtk,γk,Ytk−1=X0+αk−1tk−1r+Mk−1+σk−1Bk−1τtk−1,γk−1,ΔYtk=Ytk−Ytk−1=αktkr−αk+1tk+1r+Mk−Mk−1+σkBkτtk,γk−σk−1Bk−1τtk−1,γk−1,where(10)Xk−Xk−1=αktkr−αk−1tk−1r+σkBkτtk,γk−σk−1Bk−1τtk−1,γk−1.
Xk−Xk−1 follows Gaussian distribution.
Since,M∼Nμm,σm2, Mk−Mk−1 also follows Gaussian distribution.
Hence,ΔYtk=Ytk−Ytk−1 also follows Gaussian distribution.
Since the effect of the external environment is a random variation, it is difficult to obtain PDF directly. Therefore,Yt reaching threshold D can be equivalent to Xt reaching threshold De. RUL is defined as follows:(11)T=inft:Yt≥D|Y0<D=inft:Xt≥De|X0<De,where De is the time-varying random threshold reached by Xt, De=D−Mt, De∼ND−μm,σm2.
RULLk at tk is as follows:(12)Lk=inflk:Ytk+lk≥D=inflk:Xtk+lk≥De.
Next, the equivalent PDF is converted into a PDF in which historical degradation data is used to estimate the system RUL. It is assumed that the total degradation data obtained in the degradation process isY1:k=y1,y2,…,yk. The random covariates are supposed to be Gaussian distribution, M∼Nμm,σm2. μm and σm2 are included in equation (14). Based on the data, μm and σm2 are estimated. The degradation state caused by external stress is M1:k=m1,m2,…,mk, and the degradation state without considering the effect of external stress is X1:k=x1,x2,…,xk. Since, ytk=xtk+mtk, the RUL PDF can be calculated according to Y1:k=y1,y2,…,yk, and PDF is fLkY1:klk|Y1:k. Since α and De are independent random variables, α∼Nμα,σα2 and De∼ND−μm,σm2, the randomness of α and De should be considered in the calculation process. The full probability equation can be calculated according to the following equation:(13)fLk|Y1:klk|Y1:k=∫−∞+∞∫−∞+∞fLk|α,De,Y1:klk|α,De,Y1:kpDe|α,Y1:kpα|Y1:kdDedα=Eα|Y1:kEDe|α,Y1:kfLk|α,De,Y1:klk|α,De,Y1:k.
To calculatefLk|Y1:klk|Y1:k, Theorem 1 is given in [28].Theorem 1.
ifZ∼Nμz,σz2 and a,b,c,d,γ∈R, then:(14)EZa−bZ⋅exp−c−dZ22γ=γd2σz2+γ×a−bdσz2c+μzγd2σz2+γexp−c−dμz22d2σz2+γ.According to equation (13) and Theorem 1, fLk|Y1:klk|Y1:k can be obtained. Assuming α∼Nμα,σα2, De∼ND−μm,σm2, then:(15)fLk|Y1:klk|Y1:k=1q2πlk21lk+tkr−tkr2δα2+δm2+lkδ2D−μmlkδ2+1−δm2yk−δm2−1lk+tkr−tkr−lkrlk+tkr−1lk+tkr−tkrδα2D−yk−μm+μαδm2+lkδ2lk+tkr−tkr2δα2+δm2+lkδ2exp−D−yk−μm−lk+tkr−tkrμα22lk+tkr−tkr2δα2+δm2+lkδ2=1q2πlk21lk+tkr−tkr2σα2+qD−yk−μmlkσ2−lk+tkr−tkrlkσ2−qlkrlk+tkr−1lk+tkr−tkrσα2D−yk−μm+μαqlk+tkr−tkr2σα2+qexp−D−yk−μm−lk+tkr−tkrμα22lk+tkr−tkr2σα2+q,where q=σm2+lkσ2. Next, equation (15) is proved.Proof.
In equation (8), Dtk is replaced with De,k. According to equation (11), we can get, De,k=De−yk. Suppose De∼ND−μm,σm2, the result is as follows:(16)EDe|α,Y1:kflk|α,De,Y1:k=EDe|α,Y1:k−1σ2πlk3yk−De+αlk+tkr−αlkrlk+tkr−1−αtkr⋅exp−yk+αlk+tkr−αtkr−De22lkσ2=−12πlk2q3yk−D+μmlkσ2+αlkσ2lk+tkr−tkr−αqlkrlk+tkr−1exp−D−p−μm22q,where a=p−αlkrlk+tkr−1, c=p, b=d=1, γ=lkδ2, q=δm2+lkδ2, and p=yk+αlk+tkr−αtkr.
Furthermore, whenα∼Nμα,δα2, the result is as follows:(17)EαY1:kEDe|α,Y1:kflk|α,De,Y1:k=EαY1:k=−12πlk2q3yk−D+μmlkσ2+αlkσ2lk+tkr−tkr−αqlkrlk+tkr−1exp−D−p−μm22q=−1q2πlk21lk+tkr−tkr2σα2+qyk−D+μmlkσ2+lk+tkr−tkrlkσ2−qlkrlk+tkr−1lk+tkr−tkrσα2D−yk−μm+μαqlk+tkr−tkr2σα2+qexp−D−yk−μm−lk+tkr−tkrμα22lk+tkr−tkr2δα2+q,.where a=D−yk−μmlkσ2, b=lk+tkr−tkrlkσ2−qlkrlk+tkr−1, c=D−yk−μm, d=lk+tkr−tkr, and γ=q=σm2+lkσ2.


## 2.1. Nonlinear Wiener Degradation Model

Assuming that the degradation value of the parameter at timet is Xt, the random Wiener degradation process can be expressed as follows:(1)Xt=X0+αΛt,θ+σBτt,γ,where X0 is the initial degenerate state. When X0=0, the system is in a healthy state. α is the drift coefficient. Λt,θ=∫0tλu,θdu, λt,θ, and τt,γ is the continuous nondecreasing function about time t. θ and γ are the parameter vectors. σ is the diffusion coefficient. Bτt,γ is the Brownian motion and it follows Gaussian distribution.The forms ofΛt,θ are t, tr, and expbt. Since the degradation of slewing bearings is a nonlinear process, and tris widely used because of its flexibility in describing linear degradation paths with r=1, nonlinear concave paths with r>1, and nonlinear convex paths with r<1. Therefore, it is selected as Λt,θ=tr. The nonlinear Wiener degradation process is as follows:(2)Xt=X0+αtr+σBτt,γ.Due to differences in the production materials and production processes of the system, the degradation rate in the degradation process of each system is different. Generally, the random parameterα represents different degradation rates, α∼Nμα,σα2. Suppose, r follows Gaussian distribution, N0,1.Property 1.
the properties of the Wiener process are given as follows:(1)
Increment,Xt1−X0, Xt2−Xt1,…,Xtn−Xtn−1 on different intervals is independent of each other(2)
Xti−Xti−1 follows Gaussian distribution, with expectation αΛti,θ−Λti−1,θ and variance σ2τti;γ−τti−1;γ
According to Property1, in the Wiener process, Xt is not a Gaussian process, and it is a degradation value.


## 2.2. Nonlinear Wiener Degradation Model Based on Random Covariate

To estimate the RUL of a random degenerate system accurately, in practice, random variation of parameters, individual differences, and covariates effect are mainly considered [13, 25]. Parameter random variation means that the system degradation process is usually characterized by a random process. The individual difference means that the degradation paths of different equipment in the same system are different due to differences in working environment and load. The covariate effect means that the randomness of the external environment changes strongly, and it impacts on equipment degradation process through other components [20].The degradation process of the fan slewing bearing is influenced by the external environment (such as the change of wind speed and wind direction). The variety of external environments leads to changes in the speed and rotation direction of the blade. Consequently, it produces the stress effect on the slewing bearing device that is connected to the blade. It influences the degenerate state of slewing bearing. The stress effect is a random variation.Moreover, the proportional hazards model of covariates requires that the relationship between the failure rate function and covariates is proportional. However, in this study, the degradation process of slewing bearing is a random variation, and the relationship between the failure rate function and covariates is not strictly proportional. Generally, when the effect of an external factor on degradation is considered, it is introduced in an additive model into the degradation model. Therefore, considering the stress effect of the external environment acting on the slewing bearing through the blade propeller, as a random covariate, it is introduced into the nonlinear Wiener degradation process in the form of additive hazard models in this study. Considering the individual difference of slewing bearing samples, the random variation of drift coefficient, and the effect of a random covariate, the closed expression of RUL PDF of fan slewing bearing is derived.The stress effect of the external environment acting on the degradation process is set asM, and the total degradation state of equipment is Yt. The actual degradation process of slewing bearing at time t consists of two parts: self-degradation at time t and the degradation of stress effect at time t. The degradation model of slewing bearings is as follows:(3)Yt=Xt+M.In this study, the stress effect of the external environment acting on the slewing bearing through the blade propeller is considered. When the external environment changes, the stress effect changes randomly. Then, the degenerate state of slewing bearing is influenced, and the degradation of slewing bearing itself is determined by the inherent factors of the material. Therefore, the degeneration caused by stress and the degradation of slewing bearing itself is independent to each other. Since the Gamma distribution describes a monotonically increasing change process with time and the stress effect is a random variation, it may be nonmonotonous and not suitable for the Gamma distribution. Most of the external interference is assumed to be a Gaussian distribution [26]. According to the central limit theorem, the stress effect is supposed to be a Gaussian distribution, M∼Nμm,σm2. Since, the data-driven method is used to estimate RUL, based on the feature data, expectation μm and variance σm are estimated by using the maximum likelihood function, and then the stress effect is estimated.Combining with equations (2) and (3) in this study, the nonlinear Wiener degradation model of slewing bearing considering external environmental effects is established.(4)Yt=X0+αtr+M+σBτt,γ.The PDF is given as follows. Supposing that the degradation data obtained from the condition monitoring of the product at a discrete time pointtktk>0 is Ytk. If the service life of the product is T. Lk which is the RUL at time tk, can be expressed as Lk=T−tk. The first reaching failure threshold of the product is set as D. If Yt exceeds D, the product is judged to be invalid, and the definition of Lk is given.(5)Lk=inflk;Ylk+tk≥DYtk<D.Firstly, when the drift coefficientα is given, without considering the effect of external covariates, RUL PDF of Xt at time tk is analyzed. Once the degradation state Xt,t≥0 reaches a failure threshold D, the system is considered faulty. According to the definition of first passage time(FPT), the RUL can be defined as follows:(6)T=inft:Xt≥D|X0<D.Since it is difficult to obtain an accurate expression of the failure PDF for the nonlinear Wiener degradation model, the following equation is used to approximate the PDF according to literature [27].(7)flX1:k,θk≅12πlSll+λl,θkσexp−S2l2l.AmongSl=1/σD−Λl;θ, Λl;θ=∫tkl+tkλt;θdt, λl,θk=αrtr−1, Λl;θ=αl+tkr−tkr, and Dtk=D−Xtk, Dk is D-value between the failure threshold and the degenerate state value at time tk.Hence,Sl=1/σDtk−αl+tkr+αtkr.Letθk=μtk,σtk,r′, the PDF approximate expression of RUL Lk of Xt at tk is as follows:(8)flkX1:k,θk≅12πlkDtk−αlk+tkr+αlkrlk+tkr−1+αtkrlkσ•exp−Dtk−αlk+tkr+αtkr22lkσ2.When the degradation effect of the external environment acting on the slewing bearing is considered, according to the definition of FPT, the system is considered faulty, when the degradation stateYt,t≥0 reaches threshold D. According to Property 1, ΔYtk=Ytk−Ytk−1 are independent of each other and follow Gaussian distribution.Proof.
(9)Ytk=X0+αktkr+Mk+σkBkτtk,γk,Ytk−1=X0+αk−1tk−1r+Mk−1+σk−1Bk−1τtk−1,γk−1,ΔYtk=Ytk−Ytk−1=αktkr−αk+1tk+1r+Mk−Mk−1+σkBkτtk,γk−σk−1Bk−1τtk−1,γk−1,where(10)Xk−Xk−1=αktkr−αk−1tk−1r+σkBkτtk,γk−σk−1Bk−1τtk−1,γk−1.
Xk−Xk−1 follows Gaussian distribution.
Since,M∼Nμm,σm2, Mk−Mk−1 also follows Gaussian distribution.
Hence,ΔYtk=Ytk−Ytk−1 also follows Gaussian distribution.
Since the effect of the external environment is a random variation, it is difficult to obtain PDF directly. Therefore,Yt reaching threshold D can be equivalent to Xt reaching threshold De. RUL is defined as follows:(11)T=inft:Yt≥D|Y0<D=inft:Xt≥De|X0<De,where De is the time-varying random threshold reached by Xt, De=D−Mt, De∼ND−μm,σm2.
RULLk at tk is as follows:(12)Lk=inflk:Ytk+lk≥D=inflk:Xtk+lk≥De.
Next, the equivalent PDF is converted into a PDF in which historical degradation data is used to estimate the system RUL. It is assumed that the total degradation data obtained in the degradation process isY1:k=y1,y2,…,yk. The random covariates are supposed to be Gaussian distribution, M∼Nμm,σm2. μm and σm2 are included in equation (14). Based on the data, μm and σm2 are estimated. The degradation state caused by external stress is M1:k=m1,m2,…,mk, and the degradation state without considering the effect of external stress is X1:k=x1,x2,…,xk. Since, ytk=xtk+mtk, the RUL PDF can be calculated according to Y1:k=y1,y2,…,yk, and PDF is fLkY1:klk|Y1:k. Since α and De are independent random variables, α∼Nμα,σα2 and De∼ND−μm,σm2, the randomness of α and De should be considered in the calculation process. The full probability equation can be calculated according to the following equation:(13)fLk|Y1:klk|Y1:k=∫−∞+∞∫−∞+∞fLk|α,De,Y1:klk|α,De,Y1:kpDe|α,Y1:kpα|Y1:kdDedα=Eα|Y1:kEDe|α,Y1:kfLk|α,De,Y1:klk|α,De,Y1:k.
To calculatefLk|Y1:klk|Y1:k, Theorem 1 is given in [28].Theorem 1.
ifZ∼Nμz,σz2 and a,b,c,d,γ∈R, then:(14)EZa−bZ⋅exp−c−dZ22γ=γd2σz2+γ×a−bdσz2c+μzγd2σz2+γexp−c−dμz22d2σz2+γ.According to equation (13) and Theorem 1, fLk|Y1:klk|Y1:k can be obtained. Assuming α∼Nμα,σα2, De∼ND−μm,σm2, then:(15)fLk|Y1:klk|Y1:k=1q2πlk21lk+tkr−tkr2δα2+δm2+lkδ2D−μmlkδ2+1−δm2yk−δm2−1lk+tkr−tkr−lkrlk+tkr−1lk+tkr−tkrδα2D−yk−μm+μαδm2+lkδ2lk+tkr−tkr2δα2+δm2+lkδ2exp−D−yk−μm−lk+tkr−tkrμα22lk+tkr−tkr2δα2+δm2+lkδ2=1q2πlk21lk+tkr−tkr2σα2+qD−yk−μmlkσ2−lk+tkr−tkrlkσ2−qlkrlk+tkr−1lk+tkr−tkrσα2D−yk−μm+μαqlk+tkr−tkr2σα2+qexp−D−yk−μm−lk+tkr−tkrμα22lk+tkr−tkr2σα2+q,where q=σm2+lkσ2. Next, equation (15) is proved.Proof.
In equation (8), Dtk is replaced with De,k. According to equation (11), we can get, De,k=De−yk. Suppose De∼ND−μm,σm2, the result is as follows:(16)EDe|α,Y1:kflk|α,De,Y1:k=EDe|α,Y1:k−1σ2πlk3yk−De+αlk+tkr−αlkrlk+tkr−1−αtkr⋅exp−yk+αlk+tkr−αtkr−De22lkσ2=−12πlk2q3yk−D+μmlkσ2+αlkσ2lk+tkr−tkr−αqlkrlk+tkr−1exp−D−p−μm22q,where a=p−αlkrlk+tkr−1, c=p, b=d=1, γ=lkδ2, q=δm2+lkδ2, and p=yk+αlk+tkr−αtkr.
Furthermore, whenα∼Nμα,δα2, the result is as follows:(17)EαY1:kEDe|α,Y1:kflk|α,De,Y1:k=EαY1:k=−12πlk2q3yk−D+μmlkσ2+αlkσ2lk+tkr−tkr−αqlkrlk+tkr−1exp−D−p−μm22q=−1q2πlk21lk+tkr−tkr2σα2+qyk−D+μmlkσ2+lk+tkr−tkrlkσ2−qlkrlk+tkr−1lk+tkr−tkrσα2D−yk−μm+μαqlk+tkr−tkr2σα2+qexp−D−yk−μm−lk+tkr−tkrμα22lk+tkr−tkr2δα2+q,.where a=D−yk−μmlkσ2, b=lk+tkr−tkrlkσ2−qlkrlk+tkr−1, c=D−yk−μm, d=lk+tkr−tkr, and γ=q=σm2+lkσ2.


## 3. Parameter Estimation

Now, the parametersΘ=μα,σα2,μm,σm2,σ2,r in equation (14) are estimated. The parameter vector is Θ=μα,σα2,μm,σm2,σ2,r. Assuming that Ω=λ0,λ1,…,λk is a drift coefficient vector up to tk, to reflect the update process of Θ, Θk=μα,k,σα,k2,μm,k,σm,k2,σk2,rk is used to represent the parameter vector based on the degradation observation value Y1:k, and the parameter estimation vector is expressed as Θ^k=μ^α,k,σ^α,k2,μ^m,k,σ^m,k2,σ^k2,r^k. To estimate Θk, the maximum likelihood function is used to calculate.Assuming that the degradation state data of the operating system are measured in the order of time0<t0<t0<t1<...<tk and the measurement interval is the same. The corresponding degradation observation value is Y1:k=y1,y2,...,yk. According to equations (2) and (3), the degradation observation value can be expressed as: Yt=αΛt;θ+σBt+Mt. Assuming that, α follows Gaussian distribution, α∼Nμα,σα2, Mti follows Gaussian process, Mti∼μm,σm2. According to the property of the Wiener process, ΔYi follows Gaussian distribution, ΔYi∼NμαiΛti+μm1i,σi2Pi+σmi21i1i′+σαi2ΛtiΛti′, where 1i=1,1,…,1′ is a ni dimensional vector of 1 element.The maximum likelihood estimation method is used to determine the estimated value of parameters in this study.Let(18)Λt=tr,Qi=ti,1,ti,2,…,ti,ni′,Ψi=Λti,1,Λti,2,…,Λti,ni′,Yi=yiti,1,yiti,2,…,yiti,ni′,Y=Y1,Y2,…,Ym.The total degradation state information includes the measured values ofm products at different times. Product i is measured ni times, the cumulative degradation value corresponding at times ti,j is yiti,j, j=1,2,…,ni, and i=1,2,…,m.According to the property of the Wiener process, it is obtained that,(19)ΔYi∼Nηi,Σi,whereηi=μαiΨi+μmi1i, Σi=σi2Pi+σmi21i1i′+σαi2ΨiΨi′, Pi=ti,1ti,1⋯ti,1ti,1ti,2⋯ti,2⋮⋮⋮⋮ti,1ti,2...ti,ni.According to equation (19), the log-likelihood function of Yi can be obtained as follows:(20)lΘi,rYi=−ni2ln2π−12lnΣi−12Yi−ηi′∑i−1Yi−ηi.Furthermore, the log-likelihood function ofY is as follows:(21)lΘi,rY=−ln2π∑i=1mni2−12∑i=1mlnΣi−12∑i=1mYi−ηi′∑i−1Yi−ηi.where Θ=Θ1,Θ2,…,Θm′.r is fixed and lΘX is maximized the maximum likelihood estimation value Θ^ir of Θi is obtained.The results are as follows:(22)μ^αr=∑i=1mΨi′∑i−1Yi+Yi′∑i−1Ψi−Ψi′∑i−1μm1i−μm′1i′∑i−1Ψi∑i=1m2Ψi′∑i−1Ψi,(23)σ^α2r=∑i=1mYi−η^i′Σ/σα2−1Yi−η^i∑i=1mni,(24)μ^mr=∑i=1m1i′∑i−1Yi−1i′∑i−1μαΨi+Yi′∑i−11i−Ψi′μα′∑i−11i∑i=1m21i′∑i−11i,(25)σ^m2r=∑i=1mYi−η^i′Σ/σm2−1Yi−η^i∑i=1mni,(26)σ^2r=∑i=1mYi−η^i′Σ/σ2−1Yi−η^i∑i=1mni.LetΘ^ir,i=1,2,…,m be used in equations (21) and (27) to estimate r.(27)lrXi,Θ^i,i=1,2,…,mr∝−12∑i=1mlnΣi−12∑i=1mYi−ηi′∑i−1Yi−ηi.Using one-dimensional search method, the estimated valuer^ of r can be obtained, then r^ is used in equations (22)–(26) to determine Θ^i.


## 4. Case Analysis

To verify the effectiveness of the methods, the degradation data of fan slewing bearings obtained from the actual measurement are used for analysis. Original data comes from the acceleration data of fan slewing bearings provided by Yang et al. [29]. According to Lin et al. [30], the material of the yaw and pitch bearing is 42CrMo. The sampling frequency of the sensor is 25.6 kHz. The sampling time and the sampling interval are 1 min and 20 min, respectively. By converting the number of sampling points into the detection time, the curve of degradation data can be obtained.During the test, the a1 data, a2 data, and a3 data are vibration data collected by three sensors. The 3 groups of data have a small initial amplitude, gentle growth in the middle period, and a sharp increase in the later period. These are consistent with the general trend of signal change during device degradation. Therefore, the 3 groups of data are selected for the experimental validation. The 3 groups of data are fitted (see Figure2). The characteristics of the 3 groups of data are extracted separately, then the characteristics are used for the parameter estimation. a1 characteristics value is used to test for RUL estimation.Figure 2
The 3 groups of data are collected. (a) The a1 original signal. (b) The a2 original signal. (c) The a3 original signal.
(a)(b)(c)The key to describing the RUL of slewing bearing is the extraction of the feature of the degraded state. In practice, it is difficult to obtain direct data to characterize the operating state. Usually, indirect data are used for analysis. The vibration data monitored by the sensor is used to characterize the running state in this study. At present, the commonly used characteristics [31] of signals include maximum, peak value, variance, root mean square (RMS), and kurtosis (see Table 1). The sensitivity of each signal characteristic value to the change of fan slewing bearings operating state is different. By calculating the correlation coefficient between each characteristic value and the original signal, the characteristic value characterizing the original information is selected. According to equation (28), the correlation coefficient is calculated, and the correlation coefficient of the different signal characteristic values and the original signal is shown (see Table 1).(28)Corr=∑i=1nXi−X¯xi−x¯∑i=1nXi−X¯2•∑i=1nxi−x¯2.Table 1
Time domain characteristic indexes and correlation coefficient.
FeatureFeature expressionCorrMaximumXmax=maxxi0.87Peak valueXp−p=maxxi−minxi0.82Varianceσ2=1/n∑i=1nxi−X¯20.86RMSXRMSΔt=1/n∑i=1nxi21.00KurtosisK=1/n∑i=1nxi−X¯4/1/n∑i=1nxi−X¯220.58These formulations are used to obtain the characteristic value of the original signal at pointii=1,…,n. n is the total number of data. xi represents the original data signal value at point i, X¯ is the mean of an characteristic value in the time domain, Xi is a characteristic value of the original signal at point i, and Δt is the sampling time.Considering that the vibration signal obtained has the effect of positive and negative values, the RMS is one of the important indicators to characterize the degree of the signal, which can better reflect the change of slewing bearing performance state, and it is easy to calculate. Therefore, the RMS is only selected for feature extraction in this study.The characteristic values extracted from the 3 groups of degradation data are shown (see Figure3). The curve contains the degradation data of 12000 min from the beginning of the run-in stage to the completion of the fatigue life test.Figure 3
Feature extraction. (a) a1 characteristic value. (b) a2 characteristic value. (c) a3 characteristic value.
(a)(b)(c)The overall change trend of the whole curve after removing individual singular points reflects the corresponding relationship between the wear of the monitoring points of the experimental slewing bearings and the vibration energy. The whole curve reaches a stable phase at the beginning of the run-in stages. After 9000 min, the gear wear and the amplitude increase significantly until failure occurs.(1)
When the monitoring timet is 0–1000 min, the slewing bearings are in the run-in stages.(2)
When the monitoring timet is 1000–9000 min, the slewing bearings are in the normal wear stage.(3)
When the monitoring timet is 9000–12000 min, the characteristic value increases dramatically. It can be seen from the results after the fatigue failure of the final slewing bearings that the wear intensifies until the fault occurs at 12000 min.According to the degradation characteristics, when the monitoring timet is 0.3×103 min, the equipment starts to enter the degradation state. Therefore, the RUL estimation begins at 0.3×103 min. In [8], many destructive experiments are used to obtain the failure threshold. For the same class of fan slewing bearing, the average value of the failure threshold which is calculated based on the failure threshold of the samples is 2.3×10−3mm⋅s−2.


## 5. Remaining Useful Life Estimation

### 5.1. PDF ofYt Degradation Incremental Simulation

To verify the rationality of equation (4), Monte Carlo simulation is used to fit the PDF of Yt model. The a1 characteristic value is used as test data to verify. Then, ytk−ytk−1 is calculated, Δytk=ytk−ytk−1. Firstly, according to Δytk, the occurrence frequency of Δytk in each interval is used to replace the probability, and the histogram is obtained. To facilitate the analysis, the histogram is represented by b1. Secondly, the curve is fitted from the probability histogram, which is represented by b2. Thirdly, the expectation and variance of Δytk are calculated as 0.0047 and 0.0021, respectively. The PDF curve obeying the Gaussian distribution of expectation and variance is obtained, and the PDF curve is represented by b3 (see Figure 4).Figure 4
PDF of Monte Carlo simulation results.As can be seen from Figure4, through comparative analysis, b2 and b3 curves are close to each other. Therefore, the fitted curves follow the Gaussian distribution, and equation (4) proposed in this study can fit the degradation process of slewing bearing.


### 5.2. Experimental Verification

To facilitate analysis, the method proposed in this study is recorded as the M1 model and the traditional nonlinear Wiener degradation model is recorded as the M2 model [32]. The random Wiener degradation process of the M2 model can be expressed as follows:(29)Xt=X0+μtr+σBt.The parameter that is needed to be estimated isΘ′=μ,σ2,r. In [14, 23], the nonlinear Wiener model with proportional hazard models is recorded as the M3 model. The RUL estimation of the fan slewing bearing begins at the monitoring time t. In the early stage, the degradation data is limited. To ensure the accuracy of the fitting effect, the a1, a2, and a3 data are used for training the initial parameter estimation before the monitoring time t, Θ=μα,σα2,μm,σm2,σ2,r, and the accuracy of the model is improved. The remaining data of a1 data are selected for testing the RUL estimation after the monitoring time t. Therefore, the generalization abilities of the algorithms is actually assessed. The estimated parameter values at different monitoring points are shown in Tables 2 and 3.Table 2
Parameter estimation of the M1model and M2model.
Measurement time (×103 min)M1 ModelM2 Modelμασα2μmσm2σ2rμσ2r0.33.65910.53313.68260.45020.43160.89632.87230.99890.551021.90720.20331.88210.17520.25300.58531.48160.95610.680841.50770.18271.50750.17520.24610.20030.87640.82490.5099Table 3
Parameter estimation of the M3 model.
Measurement time (×103 min)M3 Modelμασα2μmσm2σ2r0.33.60300.44133.57340.48780.50891.115621.70300.09731.70340.10180.10190.170641.35300.08431.35340.08580.08790.1106Meanwhile, to compare the performance of the M1 model, M2 model, and M3 model, the Akaike information criterion (AIC) and Log-LF are introduced as the criteria to evaluate the fitting degree of modeling methods. The calculation equation is as follows:(30)AIC=2γ−lnLθY,where γ is the number of unknown parameters in the model and LθY is the likelihood function value.The Log-LF can skew the observations when overfitting occurs. In contrast to the Log-LF model, the AIC tries to find the most fitting model based on the model complexity and the correlation simultaneously. This would avoid overfitting issues, which are caused mainly by the high nonlinear and complex data [33]. Therefore, in this study, the Log-LF and AIC are calculated at 300 min, 2000 min, and 4000 min, respectively. And the Log-LF and AIC are used simultaneously to validate the fitting accuracy of three models. The AIC and Log-LF comparison of the three models is shown in Table 4.Table 4
AIC comparison between M1 model, M2 model, and M3 model.
Measurement time (×103 min)Modellog-LFAIC0.3M1−109.8223231.6446M2−131.4121268.8242M3−117.8287246.87382M1−537.43691086.8738M2−599.21911204.4382M3−551.43591114.87184M1−881.81271775.6254M2−985.00851976.0170M3−912.28391836.5678In Tables2 and 3, it can be concluded that the parameters, Θ=μα,σα2,μm,σm2,σ2,r, are constantly updated with the increase of the monitoring time and degradation data. At each time, the AIC value of the M1 model is less than that of the M2 model and M3 model, and the Log-LF value of the M1 model is greater than that of the M2 model and M3 model. At the same time, the AIC value of the M3 model is less than that of the M2 model, and the Log-LF value of the M3 model is greater than that of the M2 model. Therefore, the fitting effect of the M1 model is better than that of the M2 model and M3 model and the fitting effect of the M2 model is the worst.To verify the effectiveness of the M1 model in this study, RUL estimation results at different observation times are compared between the M1 model, M2 model, and M3 model (see Figure5).Figure 5
Comparison of RUL between the M1 model, M2 model, and M3 model at different times. (a)tk = 0.3 RUL estimation. (b) tk = 2 RUL estimation. (c) tk = 4 RUL estimation. (d) tk = 6 RUL estimation. (e) tk = 7 RUL estimation. (f) tk = 8 RUL estimation. (g) tk = 9 RUL estimation. (h) tk = 10 RUL estimation. (i) tk = 11 RUL estimation.
(a)(b)(c)(d)(e)(f)(g)(h)(i)In Figure5, the RUL estimation of the M1 model is compared with the RUL estimation of the M2 model, the M3 model and the real RUL. At the same, compared with the M2 models without considering the effects of the covariates, the RUL of the M1 model and M3 model considering the effects of the covariates has a higher estimation accuracy. The RUL estimation of the M1 model introducing the additive hazard model is more accurate than the RUL estimation of the M3 model introducing the proportional hazards model. The RUL estimation of the M2 model has the largest error. The estimated value of the M1 model that is proposed in this study is closer to the real value and has higher accuracy. In the early stage, less data is used, and the error is large. As the system running time increases, more monitoring data will be obtained, the estimated results will be more converged. The curve of RUL PDF becomes narrower and higher, and the variance becomes smaller. Hence, the uncertainty of estimation becomes smaller, and RUL estimation results are more accurate. It indicates that the RUL estimation of the M1 model is closer to the real value, and it has the best estimation effect.In addition, the accuracy of the M1 model is verified. The RUL estimation results of the M1 model at different times are displayed in Figure6. At the same time, it is compared with the RUL estimation results of the M2 model in Figure 7 and the RUL estimation results of the M3 model in Figure 8. From the results, it can be concluded that the RUL estimation accuracy of slewing bearings is significantly improved by using the M1 model method.Figure 6
RUL estimation of the M1 model in a1 sample.Figure 7
RUL estimation of the M2 model in a1 sample.Figure 8
RUL estimation of the M3 model in a1 sample.The accuracy of the M1 method is illustrated further by the mean square error (MSE) of RUL estimation (see Tables5, 6, and Figure 9).(31)MSEk=∫0∞lk−l˜k2fLklkY1:kdlk.Table 5
Error analysis of RUL estimation value of the M1model and M2 model.
Measurement time (×104 min)M1 Model (103 mm·s−2)Relative error (%)MSE (10−3)M2 Model (103 mm·s−2)Relative error (%)Real values (103 mm·s−2)MSE (10−3)0.374022.093.557011.766.426.63411.564.24581033.245.23357.673.5156820.264.77211.513.144868.1873.89221.232.684655.3883.16210.712.314242.8692.55150.202.013330.98101.79110.041.452820.30110.83170.030.633710.14Table 6
Error analysis of RUL estimation value of the M1 model and M3 model.
Measurement time (×104 min)M1 Model (103 mm·s−2)Relative error (%)MSE (10−3)M3 Model (103 mm·s−2)Relative error (%)Real values (103 mm·s−2)MSE (10−3)0.374022.095.435411.739.3126.63411.565.13491023.7245.23357.673.8752817.0664.77211.513.574165.9073.89221.233.043953.8483.16210.712.573642.0492.55150.202.222630.61101.79110.041.681620.10110.83170.030.732710.07Figure 9
RUL MSE comparison.Equation (32) is a discrete form of equation (31). In the experimental verification, MSE is calculated using equation (32).(32)MSEk=1n∑k=1nlk−l˜k2,where lk represents the estimated RUL at tk, l˜k represents the real RUL at time tk, and n is the total number of predictions at time tk.According to equation (32), RUL MSE can be calculated at each time. In Figure 9 and Tables 5 and 6, it can be obtained from that when the observation time is 0.3 × 104 min, the relative error and MSE are relatively large. With increasing observation data, the relative error and MSE gradually decrease. It indicates that RUL estimation is more accurate. At time 11 × 104 min, since, the slewing bearings is on the verge of failure, the relative error increases. But MSE still tends to decrease. At the same time, compared to the estimation results of the M2 model and the M3 model, the M1 method has higher accuracy and can be effectively used to fit the degenerative process of fan slewing bearing.To further illustrate the effectiveness of the M1 model in predicting the RUL of fan slewing bearing, a prediction model of fan slewing bearing based on convolutional neural networks (CNN) is established and compared with the prediction method proposed in this paper. Firstly, in the fatigue test of slewing bearing, the degenerate state characteristic value before the monitoring timet is selected as the input vector of network training. Secondly, the three-layer CNN neural network is established and used to train data. Finally, the prediction curve of the degradation state at the current time is obtained by using the trained network (see Figure 10).Figure 10
The degenerate states of the CNN model at different times. (a)tk = 9. (b) tk = 10. (c) tk = 11.
(a)(b)(c)The predicted value is compared with the actual value through a simulation example (see Figure10). Compared with different monitoring points at 9, 10, and 11, with the increase of training data, the error of RUL of the CNN model decreases gradually. In Figure 10(c), at the monitoring time 11×103 min, the average failure time of the slewing bearing predicted by the CNN model is 12.28×103 min. Compared with the CNN model, the average failure time predicted by the M1 model is 11.78×103 min, and the average actual failure time of the slewing bearing is 12×103 min. It can be seen that at the same monitoring point, the prediction accuracy of RUL based on the M1 model proposed in this paper is higher than that using the CNN model.


## 5.1. PDF ofYt Degradation Incremental Simulation

To verify the rationality of equation (4), Monte Carlo simulation is used to fit the PDF of Yt model. The a1 characteristic value is used as test data to verify. Then, ytk−ytk−1 is calculated, Δytk=ytk−ytk−1. Firstly, according to Δytk, the occurrence frequency of Δytk in each interval is used to replace the probability, and the histogram is obtained. To facilitate the analysis, the histogram is represented by b1. Secondly, the curve is fitted from the probability histogram, which is represented by b2. Thirdly, the expectation and variance of Δytk are calculated as 0.0047 and 0.0021, respectively. The PDF curve obeying the Gaussian distribution of expectation and variance is obtained, and the PDF curve is represented by b3 (see Figure 4).Figure 4
PDF of Monte Carlo simulation results.As can be seen from Figure4, through comparative analysis, b2 and b3 curves are close to each other. Therefore, the fitted curves follow the Gaussian distribution, and equation (4) proposed in this study can fit the degradation process of slewing bearing.


## 5.2. Experimental Verification

To facilitate analysis, the method proposed in this study is recorded as the M1 model and the traditional nonlinear Wiener degradation model is recorded as the M2 model [32]. The random Wiener degradation process of the M2 model can be expressed as follows:(29)Xt=X0+μtr+σBt.The parameter that is needed to be estimated isΘ′=μ,σ2,r. In [14, 23], the nonlinear Wiener model with proportional hazard models is recorded as the M3 model. The RUL estimation of the fan slewing bearing begins at the monitoring time t. In the early stage, the degradation data is limited. To ensure the accuracy of the fitting effect, the a1, a2, and a3 data are used for training the initial parameter estimation before the monitoring time t, Θ=μα,σα2,μm,σm2,σ2,r, and the accuracy of the model is improved. The remaining data of a1 data are selected for testing the RUL estimation after the monitoring time t. Therefore, the generalization abilities of the algorithms is actually assessed. The estimated parameter values at different monitoring points are shown in Tables 2 and 3.Table 2
Parameter estimation of the M1model and M2model.
Measurement time (×103 min)M1 ModelM2 Modelμασα2μmσm2σ2rμσ2r0.33.65910.53313.68260.45020.43160.89632.87230.99890.551021.90720.20331.88210.17520.25300.58531.48160.95610.680841.50770.18271.50750.17520.24610.20030.87640.82490.5099Table 3
Parameter estimation of the M3 model.
Measurement time (×103 min)M3 Modelμασα2μmσm2σ2r0.33.60300.44133.57340.48780.50891.115621.70300.09731.70340.10180.10190.170641.35300.08431.35340.08580.08790.1106Meanwhile, to compare the performance of the M1 model, M2 model, and M3 model, the Akaike information criterion (AIC) and Log-LF are introduced as the criteria to evaluate the fitting degree of modeling methods. The calculation equation is as follows:(30)AIC=2γ−lnLθY,where γ is the number of unknown parameters in the model and LθY is the likelihood function value.The Log-LF can skew the observations when overfitting occurs. In contrast to the Log-LF model, the AIC tries to find the most fitting model based on the model complexity and the correlation simultaneously. This would avoid overfitting issues, which are caused mainly by the high nonlinear and complex data [33]. Therefore, in this study, the Log-LF and AIC are calculated at 300 min, 2000 min, and 4000 min, respectively. And the Log-LF and AIC are used simultaneously to validate the fitting accuracy of three models. The AIC and Log-LF comparison of the three models is shown in Table 4.Table 4
AIC comparison between M1 model, M2 model, and M3 model.
Measurement time (×103 min)Modellog-LFAIC0.3M1−109.8223231.6446M2−131.4121268.8242M3−117.8287246.87382M1−537.43691086.8738M2−599.21911204.4382M3−551.43591114.87184M1−881.81271775.6254M2−985.00851976.0170M3−912.28391836.5678In Tables2 and 3, it can be concluded that the parameters, Θ=μα,σα2,μm,σm2,σ2,r, are constantly updated with the increase of the monitoring time and degradation data. At each time, the AIC value of the M1 model is less than that of the M2 model and M3 model, and the Log-LF value of the M1 model is greater than that of the M2 model and M3 model. At the same time, the AIC value of the M3 model is less than that of the M2 model, and the Log-LF value of the M3 model is greater than that of the M2 model. Therefore, the fitting effect of the M1 model is better than that of the M2 model and M3 model and the fitting effect of the M2 model is the worst.To verify the effectiveness of the M1 model in this study, RUL estimation results at different observation times are compared between the M1 model, M2 model, and M3 model (see Figure5).Figure 5
Comparison of RUL between the M1 model, M2 model, and M3 model at different times. (a)tk = 0.3 RUL estimation. (b) tk = 2 RUL estimation. (c) tk = 4 RUL estimation. (d) tk = 6 RUL estimation. (e) tk = 7 RUL estimation. (f) tk = 8 RUL estimation. (g) tk = 9 RUL estimation. (h) tk = 10 RUL estimation. (i) tk = 11 RUL estimation.
(a)(b)(c)(d)(e)(f)(g)(h)(i)In Figure5, the RUL estimation of the M1 model is compared with the RUL estimation of the M2 model, the M3 model and the real RUL. At the same, compared with the M2 models without considering the effects of the covariates, the RUL of the M1 model and M3 model considering the effects of the covariates has a higher estimation accuracy. The RUL estimation of the M1 model introducing the additive hazard model is more accurate than the RUL estimation of the M3 model introducing the proportional hazards model. The RUL estimation of the M2 model has the largest error. The estimated value of the M1 model that is proposed in this study is closer to the real value and has higher accuracy. In the early stage, less data is used, and the error is large. As the system running time increases, more monitoring data will be obtained, the estimated results will be more converged. The curve of RUL PDF becomes narrower and higher, and the variance becomes smaller. Hence, the uncertainty of estimation becomes smaller, and RUL estimation results are more accurate. It indicates that the RUL estimation of the M1 model is closer to the real value, and it has the best estimation effect.In addition, the accuracy of the M1 model is verified. The RUL estimation results of the M1 model at different times are displayed in Figure6. At the same time, it is compared with the RUL estimation results of the M2 model in Figure 7 and the RUL estimation results of the M3 model in Figure 8. From the results, it can be concluded that the RUL estimation accuracy of slewing bearings is significantly improved by using the M1 model method.Figure 6
RUL estimation of the M1 model in a1 sample.Figure 7
RUL estimation of the M2 model in a1 sample.Figure 8
RUL estimation of the M3 model in a1 sample.The accuracy of the M1 method is illustrated further by the mean square error (MSE) of RUL estimation (see Tables5, 6, and Figure 9).(31)MSEk=∫0∞lk−l˜k2fLklkY1:kdlk.Table 5
Error analysis of RUL estimation value of the M1model and M2 model.
Measurement time (×104 min)M1 Model (103 mm·s−2)Relative error (%)MSE (10−3)M2 Model (103 mm·s−2)Relative error (%)Real values (103 mm·s−2)MSE (10−3)0.374022.093.557011.766.426.63411.564.24581033.245.23357.673.5156820.264.77211.513.144868.1873.89221.232.684655.3883.16210.712.314242.8692.55150.202.013330.98101.79110.041.452820.30110.83170.030.633710.14Table 6
Error analysis of RUL estimation value of the M1 model and M3 model.
Measurement time (×104 min)M1 Model (103 mm·s−2)Relative error (%)MSE (10−3)M3 Model (103 mm·s−2)Relative error (%)Real values (103 mm·s−2)MSE (10−3)0.374022.095.435411.739.3126.63411.565.13491023.7245.23357.673.8752817.0664.77211.513.574165.9073.89221.233.043953.8483.16210.712.573642.0492.55150.202.222630.61101.79110.041.681620.10110.83170.030.732710.07Figure 9
RUL MSE comparison.Equation (32) is a discrete form of equation (31). In the experimental verification, MSE is calculated using equation (32).(32)MSEk=1n∑k=1nlk−l˜k2,where lk represents the estimated RUL at tk, l˜k represents the real RUL at time tk, and n is the total number of predictions at time tk.According to equation (32), RUL MSE can be calculated at each time. In Figure 9 and Tables 5 and 6, it can be obtained from that when the observation time is 0.3 × 104 min, the relative error and MSE are relatively large. With increasing observation data, the relative error and MSE gradually decrease. It indicates that RUL estimation is more accurate. At time 11 × 104 min, since, the slewing bearings is on the verge of failure, the relative error increases. But MSE still tends to decrease. At the same time, compared to the estimation results of the M2 model and the M3 model, the M1 method has higher accuracy and can be effectively used to fit the degenerative process of fan slewing bearing.To further illustrate the effectiveness of the M1 model in predicting the RUL of fan slewing bearing, a prediction model of fan slewing bearing based on convolutional neural networks (CNN) is established and compared with the prediction method proposed in this paper. Firstly, in the fatigue test of slewing bearing, the degenerate state characteristic value before the monitoring timet is selected as the input vector of network training. Secondly, the three-layer CNN neural network is established and used to train data. Finally, the prediction curve of the degradation state at the current time is obtained by using the trained network (see Figure 10).Figure 10
The degenerate states of the CNN model at different times. (a)tk = 9. (b) tk = 10. (c) tk = 11.
(a)(b)(c)The predicted value is compared with the actual value through a simulation example (see Figure10). Compared with different monitoring points at 9, 10, and 11, with the increase of training data, the error of RUL of the CNN model decreases gradually. In Figure 10(c), at the monitoring time 11×103 min, the average failure time of the slewing bearing predicted by the CNN model is 12.28×103 min. Compared with the CNN model, the average failure time predicted by the M1 model is 11.78×103 min, and the average actual failure time of the slewing bearing is 12×103 min. It can be seen that at the same monitoring point, the prediction accuracy of RUL based on the M1 model proposed in this paper is higher than that using the CNN model.


## 6. Conclusions

For the nonlinear and nonmonotonic nature of the degradation process of the fan slewing bearing, and to improve the accuracy of RUL estimation of the fan slewing bearing, the nonlinear Wiener process is selected to describe the degradation process of the fan slewing bearing. Then, at the same time, the effect of the external environment acting on the degradation of fan slewing bearing is analyzed. Therefore, the stress effect of the external environment (such as wind speed and change of wind direction) acting on the slewing bearing through the blade propeller, which is used as a random covariate, is introduced into the nonlinear Wiener degradation process in the form of additive hazard model. Meanwhile, an approximate expression is derived for RUL PDF based on the first-reach time by considering the random variation in the drift coefficients, the individual variance, and the random variation of the covariate. The parameters of the degenerate models,Θ=μα,δα2,μm,δm2,δ2,r, are estimated by maximum likelihood estimation using vibration data. In this study, the a1 sample data are tested for RUL. It is shown that as the observation time and observation data increases, the parameters, Θ=μα,δα2,μm,δm2,δ2,r, are constantly updated, and the error of RUL estimation results gradually decreases. Compared with the nonlinear Wiener model that does not consider the effect of the external environment, and the nonlinear Wiener model with proportional hazard models, the RUL estimation results of the proposed method have high accuracy and small relative errors and MSE. Therefore, the proposed method improves the reliability of the model estimation, and it can effectively be used to model the degenerative process of fan slewing bearing.



---

*Source: 5441760-2022-11-15.xml*