Bayes estimation of residual life by fusing multisource information

Abstract

Residual life estimation is essential for reliability engineering. Traditional methods may experience difficulties in estimating the residual life of products with high reliability, long life, and small sample. The Bayes model provides a feasible solution and can be a useful tool for fusing multisource information. In this study, a Bayes model is proposed to estimate the residual life of products by fusing expert knowledge, degradation data, and lifetime data. The linear Wiener process is used to model degradation data, whereas lifetime data are described via the inverse Gaussian distribution. Therefore, the joint maximum likelihood (ML) function can be obtained by combining lifetime and degradation data. Expert knowledge is used according to the maximum entropy method to determine the prior distributions of parameters, thereby making this work different from existing studies that use non-informative prior. The discussion and analysis of different types of expert knowledge also distinguish our research from others. Expert knowledge can be classified into three categories according to practical engineering. Methods for determining prior distribution by using the aforementioned three types of data are presented. The Markov chain Monte Carlo is applied to obtain samples of the parameters and to estimate the residual life of products due to the complexity of the joint ML function and the posterior distribution of parameters. Finally, a numerical example is presented. The effectiveness and practicability of the proposed method are validated by comparing it with residual life estimation that uses non-informative prior. Then, its accuracy and correctness are proven via simulation experiments.

Keywords

residual life estimation / Bayes model / linear Wiener process

Introduction

Residual life estimation is currently essential for reliability engineering. However, given the development in the manufacturing process, satellite tests exhibit the characte-ristics of high reliability and long life, and collected lifetime data are typically zero-failure or few-failure data. These situations result in considerable difficulty in residual life estimation. In practical engineering, various types of reliability information, including lifetime data, degradation data, and expert information, can be collected. However, traditional methods are no longer effective for fusing information. Fully utilizing the aforementioned reliability information has received sufficient attention and presents a solution for the residual life estimation of products with high reliability, long life, and small sample. The Bayes model can be a useful tool for solving this problem and is widely used to fuse multisource information.

The two most commonly used methods are as follows. The first method involves assuming that the degradation process conforms to a certain distribution, and the maximum likelihood (ML) function is established by fusing multisource data, followed by the estimation of parameters and residual life. ^{Pettit and Young (1999)} assumed that degradation is a Wiener process, and they fused failure lifetime and degradation data to estimate product lifetime. ^{Wang et al. (2012)} proposed an effective method that considers prior degradation information in evaluating reliability. Updating parameters using the Bayes method is one of the major advantages of their method. ^{Gebraeel et al. (2009)} replaced the degradation information of similar products with lifetime data and the estimated residual life by using lifetime data and field degradation information. ^{Padgett et al. (2004)} described degradation via the Gaussian process and achieved reliability evalu-ation by fusing failure lifetime and degradation data. Additional information in this field is available in ^{Yang and Yang (2011)} and ^{Peng et al. (2009)}.

Another method for fusing multisource information obtains prior distributions according to information sources and then fuses these distributions, which is the most important step. Zhang (2001) proposed a method based on the confidence level of prior information. Feng presented several methods for fusing prior information. These methods exhibit the following characteristics: A correlation information fusion method (^{Feng et al., 2003}), an adequacy measurement-based method (^{Feng et al., 2004}), a multisource information fusion method based on the maximum entropy Markov model (^{Feng et al., 2003}), an information fusion approach based on the Bayes fuzzy logistic operator (^{Feng and Zhou., 2008}), and the fusion of information from multiple sources based on ML-II theory in Bayesian analysis (^{Feng et al., 2006}). ^{Chai et al. (2005)} proposed a method for determining fusion weight based on the ML principle.

The conclusions that can be drawn after an extensive literature review of the aforementioned developments are as follows. (1) The literature on the information fusion of expert knowledge, lifetime data, and degradation data remains scarce. (2) Methods for determining prior distribution according to different types of expert information are limited. (3) The effect of fusing information from different types of sources requires further study. These conclusions have motivated our research.

This study aims to propose a Bayes residual life estimation method by fusing multisource information. First, prior distributions of parameters are determined using the maximum entropy method. Then, the joint ML function is established. Parameter estimation and residual life estimation can be achieved through the Gibbs sampling method. The effects of fusing different types of information sources are compared and discussed through a numerical study. In particular, expert knowledge is classified into three categories in practical engineering: (1) estimation of reliability, (2) estimation of lifetime, and (3) estimation of reliability and lifetime. In this study, we propose methods for determining the prior distribution of parameters under the three conditions.

The remainder of this paper is organized as follows. Section 2 presents the cumulative distribution function (CDF) and probability density function (PDF) of residual life based on the Wiener process. Section 3 introduces the process of obtaining the prior distributions of parameters according to different types of expert information. The ML function is established and the posterior distribution is obtained by fusing lifetime and degradation data in Section 4. Parameter estimation and residual life prediction are demonstrated in Section 5. A numerical example is provided in Section 6.

CDF and PDF of residual life based on the linear Wiener process

The Wiener process can describe the degradation process of a product and exhibits good calculation and analysis properties. Therefore, the linear drift Wiener process model is selected for this study.

Let X\left(t\right) represent the degradation value at time t. If X\left(t\right) satisfies

(1)‚X\left(\mathrm{0}\right)=\mathrm{0};

(2)‚\left\{X\right(t),t>\mathrm{0}\} demonstrates independent and steady increments;

(3)‚X\left(t\right) is subject to a normal distribution with a mean value of \mu t and a variance of {\sigma }^{\mathrm{2}}t;

then \left\{X\right(t),t>\mathrm{0}\} is defined as a linear drift Wiener process with parameters \mu and \sigma. Evidently, X\left(t\right) exhibits the following form according to the preceding definition (^{Wang et al., 2015}):

(1)X(t)=\mu t+\sigma B\left(t\right),

where B\left(t\right) is the standard Brownian motion.

The CDF and PDF of lifetime, denoted by F\left(t\right) and f\left(t\right), respectively, can be presented as follows (^{Chen, 2016}):

(2)F\left(t\right)=\mathrm{1}-P(T\ge t)=\mathrm{\Phi }\left({\displaystyle \frac{\mu t-D}{\sigma \sqrt{t}}}\right)+\exp \left({\displaystyle \frac{\mathrm{2}\mu D}{{\sigma }^{\mathrm{2}}}}\right)\mathrm{\Phi }\left({\displaystyle \frac{-D-\mu t}{\sigma \sqrt{t}}}\right),

(3)f\left(t\right)={\displaystyle \frac{dF}{dt}}={\displaystyle \frac{D}{\sqrt{\mathrm{2}\mathrm{\pi }{\sigma }^{\mathrm{2}}{t}^{\mathrm{3}}}}}\exp [-{\displaystyle \frac{{(D-\mu t)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}t}}].

Then, the PDF of residual life at time \tau is given as follows:

(4)f\left(t\right)={\displaystyle \frac{dF}{dt}}={\displaystyle \frac{D-X_{\tau }}{\sqrt{\mathrm{2}\mathrm{\pi }{\sigma }^{\mathrm{2}}{t}^{\mathrm{3}}}}}\exp [-{\displaystyle \frac{{(D-\mu t)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}t}}],

where X_{\tau } represents the degradation data at \tau.

Therefore, the key step in lifetime and residual life estimation is to determine \mu and \sigma. The Bayes model can be used to fuse multisource information to improve the accuracy of parameter estimation. In the next section, the method for determining the prior distributions of parameters is presented.

Prior distributions of parameters

We assume that the lifetime of a certain product on a satellite platform conforms to the linear Wiener distribution. In practical engineering, the types of expert information can be divided into three categories: (1) point estimation of reliability at a specific time, (2) point estimation of lifetime, and (3) point estimation of reliability and lifetime. Methods for obtaining the prior distributions of parameters according to different types of expert information are proposed in this section.

Prior distribution based on expert estimation of reliability

Expert information about the parameters is as follows: The point estimation of reliability is R_{\mathrm{0}} at mission time t_{\mathrm{0}}. For this form of description, we can use information entropy as an objective function.

In accordance with the principle of maximum entropy:

(5)\mathrm{m}\mathrm{a}\mathrm{x}\mathrm{}H\left[\mathrm{\pi }\left(\mu \right)\right]=-{\displaystyle {\int }_{\mathrm{0}}^{+\infty }\mathrm{\pi }\left(\mu \right)\mathrm{l}\mathrm{n}\left[\mathrm{\pi }\left(\mu \right)\right]d\mu }.

Then, the extreme constraint problem can be expressed as follows:

(6)\{\begin{array}{l}\max \mathrm{}H\left[\mathrm{\pi }\left(\mu \right)\right]=-{\displaystyle {\int }_{\mathrm{0}}^{+\infty }\mathrm{\pi }\left(\mu \right)\ln \left[\mathrm{\pi }\left(\mu \right)\right]d\mu }\hfill \\ s.t.{\displaystyle {\int }_{\Theta }R\left(\mu ,{\sigma }^{\mathrm{2}},t_{\mathrm{0}}\right)\mathrm{\pi }\left(\mu \right)\mathrm{\pi }\left({\sigma }^{\mathrm{2}}\right)d\mu d{\sigma }^{\mathrm{2}}=R_{\mathrm{0}}}\hfill \end{array},

where \Theta represents the ranges of \mu and {\sigma }^{\mathrm{2}}. We use the Gamma distribution as the prior distributions of \mu and {\sigma }^{\mathrm{2}}:

(7)\mathrm{\pi }\left(x|a,b\right)={\displaystyle \frac{b^a}{\Gamma \left(a\right)}}{x}^{a-\mathrm{1}}\mathrm{e}\mathrm{x}\mathrm{p}\left(-bx\right).

Then, a and b can be calculated by solving this nonlinear programming problem with a nonlinear equality constraint, and the prior distributions of \mu and {\sigma }^{\mathrm{2}} can be determined. This programming problem can be difficult to solve, and a genetic algorithm can be a feasible approach to address it. We can also assume that {\sigma }^{\mathrm{2}} is uniformly distributed to simplify this problem.

Prior distribution based on expert estimation of lifetime

Expert information about the parameters is as follows: The point estimation of lifetime is L_{\mathrm{0}}.

In accordance with the principle of maximum entropy, the problem is expressed as follows (^{Liu, 2017}):

(8)\{\begin{array}{c}\hfill \max H\left[\mathrm{\pi }\left(\mu \right)\right]=-{\displaystyle {\int }_{\Theta }^{+\infty }\mathrm{\pi }\left(\mu \right)\ln \left[\mathrm{\pi }\left(\mu \right)\right]d\mu }\hfill \\ \hfill s.t.{\displaystyle {\int }_{{\Theta }_{\mathrm{1}}}\mu \mathrm{\pi }\left(\mu \right)d\mu ={\mu }_{\mathrm{0}}}\hfill \end{array},

where {\mu }_{\mathrm{0}} represents the estimation of \mu.

In accordance with Eq. (8), let \lambda ={\displaystyle \frac{D^{\mathrm{2}}}{{\sigma }^{\mathrm{2}}}} and \beta ={\displaystyle \frac{D}{\mu }}. Then,

f\left(t\right)={\displaystyle \frac{D}{\sqrt{\mathrm{2}\mathrm{\pi }{\sigma }^{\mathrm{2}}{t}^{\mathrm{3}}}}}\exp [-{\displaystyle \frac{{(D-\mu t)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}t}}]

=\sqrt{{\displaystyle \frac{D^{\mathrm{2}}}{\mathrm{2}\mathrm{\pi }{\sigma }^{\mathrm{2}}{t}^{\mathrm{3}}}}}\exp \left[{\displaystyle -\frac{D^{\mathrm{2}}}{{\sigma }^{\mathrm{2}}}·}{\displaystyle \frac{{(\mathrm{1}-\frac{t}{\beta })}^{\mathrm{2}}}{\mathrm{2}t}}\right]

(9)\begin{array}{l}=\sqrt{{\displaystyle \frac{\lambda }{\mathrm{2}\mathrm{\pi }{t}^{\mathrm{3}}}}}\exp [-{\displaystyle \frac{\lambda {(t-\beta )}^{\mathrm{2}}}{\mathrm{2}{\beta }^{\mathrm{2}}t}}],\end{array}

which is the standard form of the Wiener process. Let x={\displaystyle \frac{\mathrm{1}}{t}}. Then, dx=-{\displaystyle \frac{\mathrm{1}}{t^{\mathrm{2}}}}dt. Its expectation E\left(t\right) can be solved as follows:

\begin{array}{l}E\left(t\right)={\displaystyle {\int }_{\mathrm{0}}^{+\infty }tf\left(t\right)dt=}{\displaystyle {\int }_{\mathrm{0}}^{+\infty }t\sqrt{{\displaystyle \frac{\lambda }{\mathrm{2}\mathrm{\pi }{t}^{\mathrm{3}}}}}\exp [-{\displaystyle \frac{\lambda {(t-\beta )}^{\mathrm{2}}}{\mathrm{2}{\beta }^{\mathrm{2}}t}}]dt}\hfill \end{array}

={\displaystyle {\int }_{\mathrm{0}}^{+\infty }\sqrt{{\displaystyle \frac{\lambda }{\mathrm{2}\mathrm{\pi }{\displaystyle \frac{\mathrm{1}}{x}}}}}\exp [-{\displaystyle \frac{\lambda {({\displaystyle \frac{\mathrm{1}}{x}}-\beta )}^{\mathrm{2}}}{\mathrm{2}{\beta }^{\mathrm{2}}{\displaystyle \frac{\mathrm{1}}{x}}}}]·(-{\displaystyle \frac{\mathrm{1}}{x^{\mathrm{2}}}})dx}

{\displaystyle =}\sqrt{{\displaystyle \frac{\lambda }{\mathrm{2}\mathrm{\pi }}}}{\displaystyle {\int }_{\mathrm{0}}^{+\infty }(-x^{-\frac{\mathrm{3}}{\mathrm{2}}})\exp [-{\displaystyle \frac{\lambda {(\mathrm{1}-\beta x)}^{\mathrm{2}}}{\mathrm{2}{\beta }^{\mathrm{2}}x}}]dx}

=\sqrt{{\displaystyle \frac{\lambda }{\mathrm{2}\mathrm{\pi }}}}{\displaystyle {\int }_{\mathrm{0}}^{+\infty }[-\beta {\left({\beta }^{\mathrm{2}}x\right)}^{-\frac{\mathrm{3}}{\mathrm{2}}}]}

(10)·\exp [-{\displaystyle \frac{\lambda }{\mathrm{2}{\beta }^{\mathrm{2}}}}·{\displaystyle \frac{{(\mathrm{1}-{\beta }^{\mathrm{2}}x)}^{\mathrm{2}}}{{\beta }^{\mathrm{2}}x}}]d\left({\beta }^{\mathrm{2}}x\right)=\beta .

Therefore, we obtain {\mu }_{\mathrm{0}}={\displaystyle \frac{D}{L_{\mathrm{0}}}}. Assume that the distribution of \mu is the Gamma distribution as follows:

(11)\mathrm{\pi }\left(\mu \right)={\displaystyle \frac{b^a}{\Gamma \left(a\right)}}{\mu }^{a-\mathrm{1}}{e}^{-b\mu }.

Then, the expectation of \mu can be calculated as follows:

E\left(\mu \right)={\displaystyle {\int }_{\mathrm{0}}^{+\infty }\mu \mathrm{\pi }\left(\mu \right)d\mu }={\displaystyle {\int }_{\mathrm{0}}^{+\infty }{\displaystyle \frac{b^a}{\Gamma \left(a\right)}}{\mu }^a{e}^{-b\mu }d\mu }

={\displaystyle {\int }_{\mathrm{0}}^{+\infty }{\displaystyle \frac{\mathrm{1}}{b\Gamma \left(a\right)}}{\left(b\mu \right)}^a{e}^{-b\mu }d\left(b\mu \right)=}{\displaystyle {\int }_{\mathrm{0}}^{+\infty }{\displaystyle \frac{\mathrm{1}}{b\Gamma \left(a\right)}}·t^a{e}^{-t}dt}

(12)={\displaystyle \frac{\mathrm{1}}{b\Gamma \left(a\right)}}·\Gamma \left(a+\mathrm{1}\right)={\displaystyle \frac{a}{b}}.

The obtained constraint condition is {\displaystyle \frac{a}{b}}={\mu }_{\mathrm{0}}={\displaystyle \frac{D}{L_{\mathrm{0}}}}. Hence, we have

(13)b=a{\displaystyle \frac{L_{\mathrm{0}}}{D}}.

We can transform the constrained nonlinear programming problem into an unconstrained nonlinear programming problem by using the preceding equation. Finally, MATLAB is a powerful tool for solving this problem. The obtained accuracy and efficiency are satisfactory. Given the lack of information, we also assume that {\sigma }^{\mathrm{2}} is uniformly distributed.

Prior distribution based on expert estimation of lifetime and reliability

To simplify the calculation, we assume that {\sigma }^{\mathrm{2}} conforms to the uniform distribution presented in Sections 3.1 and 3.2. However, after fusing expert data in the form of reliability at a specific time, the Gamma distribution can be used as a prior distribution of {\sigma }^{\mathrm{2}}. The expert data of both types of point estimation can be used to obtain a more accurate distribution of the residual life of products by replacing the original uniform distribution. The specific method is described as follows.

Expert information about the parameters is as follows: the point estimation of reliability at mission time t_{\mathrm{0}} is R_{\mathrm{0}} and that of lifetime is L_{\mathrm{0}}.

The Gamma distribution is selected as the prior distribution of {\sigma }^{\mathrm{2}}. We derive {\mu }_{\mathrm{0}}={\displaystyle \frac{D}{L_{\mathrm{0}}}} in Section 2.2. The point estimation of {\sigma }^{\mathrm{2}}, which is denoted by {\sigma }_{\mathrm{0}}^{\mathrm{2}}, can be obtained by solving R\left({\mu }_{\mathrm{0}},{\sigma }^{\mathrm{2}},t\right)=R_{\mathrm{0}}.

Similarly, the programming problem is expressed as follows according to the principle of maximum entropy:

(14)\{\begin{array}{c}\hfill \max H\left[\mathrm{\pi }\left({\sigma }^{\mathrm{2}}\right)\right]=-{\displaystyle {\int }_{\mathrm{0}}^{+\infty }\mathrm{\pi }\left({\sigma }^{\mathrm{2}}\right)\ln \left[\mathrm{\pi }\left({\sigma }^{\mathrm{2}}\right)\right]d{\sigma }^{\mathrm{2}}}\hfill \\ \hfill s.t.{\displaystyle \underset{\mathrm{0}}{\overset{\infty }{\int }}{\sigma }^{\mathrm{2}}\mathrm{\pi }\left({\sigma }^{\mathrm{2}}\right)d{\sigma }^{\mathrm{2}}={\sigma }_{\mathrm{0}}^{\mathrm{2}}}\hfill \end{array}.

After obtaining the values of the parameters of the Gamma distribution and substituting them into Eq. (14), the prior distribution of {\sigma }^{\mathrm{2}} is derived. The method for determining the prior distribution of \mu is presented in Section 3.2.

Given that \mu and {\sigma }^{\mathrm{2}} are mutually independent, the joint distribution function of \mu and {\sigma }^{\mathrm{2}}, which can also be regarded as the joint prior distribution, is expressed as follows:

(15)\mathrm{\pi }(\mu ,{\sigma }^{\mathrm{2}})=\mathrm{\pi }\left(\mu \right)\mathrm{\pi }\left({\sigma }^{\mathrm{2}}\right).

When the prior distributions of \mu and {\sigma }^{\mathrm{2}} that we obtained under the three conditions are substituted, we can derive the joint prior distribution \mathrm{\pi }(\mu ,{\sigma }^{\mathrm{2}}).

Posterior distributions of the parameters

We assume that we have collected m lifetime data, including k\left(k\le m\right) failure lifetime data T_{\mathrm{1}},T_{\mathrm{2}},...,T_k and m-k censored lifetime data T_{C_{\mathrm{1}}},T_{C_{\mathrm{2}}},...,T_{C_{m-k}}, where D is the failure threshold of degradation. Given the drift parameter \mu and the diffusion coefficient \sigma, the ML function based on lifetime data can be given as

L\left(\mu ,{\sigma }^{\mathrm{2}}|T\right)={\displaystyle \prod _{i=\mathrm{1}}^k{\displaystyle \frac{D}{\sqrt{\mathrm{2}\mathrm{\pi }{\sigma }^{\mathrm{2}}{T}_i^{\mathrm{3}}}}}\mathrm{e}\mathrm{x}\mathrm{p}\left[-{\displaystyle \frac{{\left(D-\mu T_i\right)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}{T}_i}}\right]}

(16)\times {\displaystyle \prod _{i=\mathrm{1}}^{m-k}\left[\Phi \left({\displaystyle \frac{D-\mu T_{C_i}}{\sigma \sqrt{T_{C_i}}}}\right)+\mathrm{e}\mathrm{x}\mathrm{p}\left({\displaystyle \frac{\mathrm{2}\mu D}{{\sigma }^{\mathrm{2}}}}\right)\Phi \left({\displaystyle \frac{-D-\mu T_{C_i}}{\sigma \sqrt{T_{C_i}}}}\right)\right]}.

Assume that the performance degradation data of n products were collected. At the initial time t_{i\mathrm{0}}, the performance degradation value for sample i is X_{i\mathrm{0}}=\mathrm{0}. The degradation values at times t_{i\mathrm{1}},t_{i\mathrm{2}},\cdots ,t_{imi} are X_{i\mathrm{1}},X_{i\mathrm{2}},\cdots ,X_{imi}. The change in the performance degradation of product i between times t_{i(j-\mathrm{1})} and tij is \Delta x_{ij}=Xij-X_{i(j-\mathrm{1})}. Hence,

(17)\Delta x_{ij}\sim N\left(\mu \Delta t_{ij},{\sigma }^{\mathrm{2}}\Delta t_{ij}\right),

where \Delta t_{ij}=t_{ij}-t_{i(j-\mathrm{1})},i=\mathrm{1},\mathrm{2},...,n,j=\mathrm{1},\mathrm{2},...,m_i.

The ML function of the parameters established based on the performance degradation data are

(18)L\left(\mu ,{\sigma }^{\mathrm{2}}|X\right)={\displaystyle \prod _{i=\mathrm{1}}^n{\displaystyle \prod _{j=\mathrm{1}}^{mi}{\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{2}{\sigma }^{\mathrm{2}}\Delta t_{ij}}}}\mathrm{e}\mathrm{x}\mathrm{p}\left[-{\displaystyle \frac{{\left(\Delta x_{ij}-\mu \Delta t_{ij}\right)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}\Delta t_{ij}}}\right]}}.

When lifetime and performance degradation data are combined, the ML function is obtained as follows:

L\left(\mu ,{\sigma }^{\mathrm{2}}|X,T\right)={\displaystyle \prod _{i=\mathrm{1}}^n{\displaystyle \prod _{j=\mathrm{1}}^{m_i}{\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{2}{\sigma }^{\mathrm{2}}\Delta t_{ij}}}}\mathrm{e}\mathrm{x}\mathrm{p}\left[-{\displaystyle \frac{{\left(x_{ij}-\mu \Delta t_{ij}\right)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}\Delta t_{ij}}}\right]}}

{\displaystyle \prod _{i=\mathrm{1}}^k{\displaystyle \frac{D}{\sqrt{\mathrm{2}\mathrm{\pi }{\sigma }^{\mathrm{2}}{T}_i^{\mathrm{3}}}}}\mathrm{e}\mathrm{x}\mathrm{p}\left[-{\displaystyle \frac{{(D-\mu T_i)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}{T}_i}}\right]}{\displaystyle \prod _{i=\mathrm{1}}^{m-k}[F\left({\displaystyle \frac{D-\mu T_{C_i}}{\sigma \sqrt{T_{C_i}}}}\right)}

(19){\displaystyle +\mathrm{e}\mathrm{x}\mathrm{p}\left({\displaystyle \frac{\mathrm{2}\mu D}{{\sigma }^{\mathrm{2}}}}\right)F\left({\displaystyle \frac{-D-\mu T_{C_i}}{\sigma \sqrt{T_{C_i}}}}\right)]}.

Then, the posterior distributions of the parameters are shown according to Bayesian theory as follows:

(20)\mathrm{\pi }(\mu ,{\sigma }^{\mathrm{2}}|X,T)={\displaystyle \frac{\mathrm{\pi }(\mu ,{\sigma }^{\mathrm{2}})L\left(\mu ,{\sigma }^{\mathrm{2}}|X,T\right)}{{\displaystyle \iint \mathrm{\pi }(\mu ,{\sigma }^{\mathrm{2}})L\left(\mu ,{\sigma }^{\mathrm{2}}|X,T\right)d\mu d{\sigma }^{\mathrm{2}}}}}.

Parameter estimation and residual life prediction

Estimating the parameters of a distribution is a difficult task because of the complexity of the ML function and the posterior distribution of parameters \mu and \sigma. Therefore, the Markov chain Monte Carlo sampling method is selected to address this problem.

The conditional density function of \mu is

\mathrm{\pi }\left(\mu |{\sigma }^{\mathrm{2}},T_h,X_h\right)\propto \mathrm{\pi }\left(\mu ,{\sigma }^{\mathrm{2}}\right)

(21)\begin{array}{l}\mathrm{e}\mathrm{x}\mathrm{p}\left[-{\displaystyle \sum _{i=\mathrm{1}}^n{\displaystyle \sum _{j=\mathrm{1}}^{m_i}{\displaystyle \frac{{\left(\Delta x_{ij}-\mu \Delta t_{ij}\right)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}\Delta t_{ij}}}-{\displaystyle \sum _{i=\mathrm{1}}^k\frac{{\left(D-\mu T_i\right)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}{T}_i}}}}\right]\hfill \\ {\displaystyle \prod _{i=\mathrm{1}}^{m-k}\left[\Phi \left({\displaystyle \frac{D-\mu T_{C_i}}{\sigma \sqrt{T_{C_i}}}}\right)+\mathrm{e}\mathrm{x}\mathrm{p}\left({\displaystyle \frac{\mathrm{2}\mu D}{{\sigma }^{\mathrm{2}}}}\right)\Phi \left({\displaystyle \frac{-D-\mu T_{C_i}}{\sigma \sqrt{T_{C_i}}}}\right)\right]}.\hfill \end{array}

The conditional density function of {\sigma }^{\mathrm{2}} is

\mathrm{\pi }\left({\sigma }^{\mathrm{2}}|\mu ,T,X\right)\propto \mathrm{\pi }\left(\mu ,{\sigma }^{\mathrm{2}}\right).{\sigma }^{-\sum _{i=\mathrm{1}}^n{m}_i+m}

\times \mathrm{e}\mathrm{x}\mathrm{p}\left[-{\displaystyle \sum _{i=\mathrm{1}}^n{\displaystyle \sum _{j=\mathrm{1}}^{m_i}{\displaystyle \frac{{\left(\Delta x_{ij}-\mu \Delta t_{ij}\right)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}\Delta t_{ij}}}-{\displaystyle \sum _{i=\mathrm{1}}^k\frac{{\left(D-\mu T_i\right)}^{\mathrm{2}}}{\mathrm{2}{\sigma }^{\mathrm{2}}{T}_i}}}}\right]

(22){\displaystyle \times \prod _{i=\mathrm{1}}^{m-k}\left[\Phi \left({\displaystyle \frac{D-\mu T_{C_i}}{\sigma \sqrt{T_{C_i}}}}\right)+exp\left({\displaystyle \frac{\mathrm{2}\mu D}{{\sigma }^{\mathrm{2}}}}\right)\Phi \left({\displaystyle \frac{-D-\mu T_{C_i}}{\sigma \sqrt{T_{C_i}}}}\right)\right]}.

Given sample size T. The steps of the Gibbs sampling method are as follows.

(1) Let t=1.

(2) The initial values {\mu }^{\left(\mathrm{1}\right)} and {\left({\sigma }^{\mathrm{2}}\right)}^{\left(\mathrm{1}\right)} of parameters \mu and {\sigma }^{\mathrm{2}} are randomly generated.

(3) Let t=t+1. {\mu }^{\left(t\right)} and {\left({\sigma }^{\mathrm{2}}\right)}^{\left(t\right)} are sampled according to the conditional PDF.

(4) Repeat Step 3 until all the values with size T are obtained.

After the data set is obtained via Gibbs sampling, we estimate the parameters by calculating the mean value of samples as follows:

\widehat{\mu }=E\left[\mu \right]=\underset{T\to \infty }{\mathrm{l}\mathrm{i}\mathrm{m}}{\displaystyle \frac{\mathrm{1}}{T}}{\displaystyle \sum _{i=\mathrm{1}}^t{\mu }^{\left(i\right)}}

(23)\widehat{{\sigma }^{\mathrm{2}}}=E\left[{\sigma }^{\mathrm{2}}\right]=\underset{T\to \infty }{\mathrm{l}\mathrm{i}\mathrm{m}}{\displaystyle \frac{\mathrm{1}}{T}}{\displaystyle \sum _{i=\mathrm{1}}^t{\left({\sigma }^{\mathrm{2}}\right)}^{\left(i\right)}}.

Finally, the PDF and CDF of residual life can be obtained after estimating the parameters. Then, we can predict the residual life of the product.

Numerical example

In this section, a certain product of a satellite platform is selected as the object of analysis to achieve residual life estimation based on the Bayesian model by fusing multisource information. The feasibility and practicability of the method are validated. Then, the accuracy of the proposed method is proven through simulation experiments by comparing the proposed method with methods that use non-informative prior and a degradation model.

The reliability information collected for this product is presented as follows.

(1) Expert knowledge: The point estimate of product lifetime is L_{\mathrm{0}}=\mathrm{40} months.

(2) The lifetime data are shown in Fig. 1, all of which are zero-failure lifetime data.

(3) The degradation data are the temperatures of this product with a failure threshold of 42°C, which are roughly subject to the linear Wiener process.

Data processing

The median filtering method is used to remove outliers from the performance degradation data (^{Liu, 2017}), and the degradation data are converted to once a day. The curve of the temperature is plotted in Fig. 2(a). To facilitate data analysis, temperature data are converted into once a half-month, as shown in Fig. 2(b).

Prior distribution of the parameters with information

The maximum entropy method and the optimization tool in MATLAB are used according to the point estimate of lifetime by experts to obtain a=\mathrm{1.0045},b=\mathrm{2.009}.

When combining lifetime and degradation data, we use Gibbs sampling to obtain the sample trajectory of \mu (Fig. 3) and the posterior distribution of \mu (Fig. 4).

Then, we can get the point estimates of the parameters from the prior distribution with information.

\mu =\mathrm{0.5441},\sigma =\mathrm{0.5939}

The initial performance degradation value is X_{t\mathrm{0}}=\mathrm{19.34}°C. When \tau =\mathrm{200} days, X_{\tau }=\mathrm{22.43}°C. Therefore, the PDFs of lifetime and residual lifetime are shown in Figs. 5(a) and 5(b), respectively.

Residual life estimation based on non-informative prior

When expert knowledge is unavailable, we assume that the prior distributions of \mu and {\sigma }^{\mathrm{2}} are uniform distributions. After fusing lifetime and degradation data, the sample trajectory of \mu that is obtained via sampling is shown in Fig. 6 and the posterior distribution of \mu is shown in Fig. 7.

Then, we can obtain the point estimates of the parameters.

{\mu }_{\mathrm{0}}=\mathrm{0.6031},\hspace{0.17em}\hspace{0.17em}{\sigma }_{\mathrm{0}}=\mathrm{0}\mathrm{.5745}

The PDFs of lifetime and residual life can be determined and are presented and compared in Section 6.4.

Comparative analysis

Table 1 provides three cases: prior distribution of parameters with information (Case 1), prior distribution of parameters without information (Case 2), and modeling of performance degradation data (Case 3). The lifetime and residual life estimations of the three cases are compared.

Posterior distributions of the three cases are depicted in Fig. 8.

A conclusion can be drawn that the degradation model can solve the problem of estimating the residual life of products with high reliability, long life, and small sample. Nevertheless, limitations also exist. Determining the threshold and goodness of fit of the degradation model influences the accuracy of estimation. The comparison of the posterior distributions of the parameters shows that the distribution of a parameter becomes more concentrated by fusing multisource information, which can improve the precision of prediction. The difference between the estimation of \mu obtained in Cases 1 and 2 has drawn our attention because expert knowledge does not entirely correspond to the degradation data. Therefore, the accuracy of expert information is high in estimating residual life. If the attitude of experts is too optimistic, then the prediction of residual life may be longer.

Simulation experiment

To test the accuracy of the proposed method, we assume that µ = 1.5, s = 0.4, and D = 90°C. Then, the expected lifetime is 60 months by using Eq. (3). Therefore, the expert information is as follows: the estimation of lifetime is 60 months. The degradation data between 1 and 20 months are generated via simulation and shown in Fig. 9.

Failure lifetime data (T_fail) and zero-failure lifetime data (T_censor) are provided in Table 2.

Therefore, when current time \tau =\mathrm{10} months, the true value of residual life is 49.1141 months by integrating Eq. (4). Posterior distributions of three cases are depicted in Fig. 10.

Additional details are provided in Table 3.

The accuracy of lifetime and residual life estimations is improved by fusing expert knowledge. Moreover, the results are more precise if more information is fused and used. We should emphasize that the advantage of this method will be more evident if lifetime and degradation data are not easily tractable.

Conclusions

A Bayes method for residual life estimation by fusing multisource information is proposed in this study. The determination of the prior distributions of parameters by using different types of expert information is discussed. The proposed method is proven to be valid and accurate by comparing with methods that use non-informative prior and the degradation model. Moreover, the residual life of a product is approximately estimated. The accuracy of expert information is important. Inaccurate expert information may lead to considerable errors. Nevertheless, the proposed method is limited. Guaranteeing or testing the accuracy of expert information should be studied further. Moreover, additional efforts should be exerted to improve this method when prior information is changed or the degradation process is different.

References

[1] Chai J, Shi Y M, Wei J Q, Li X C (2005). The fusion method for prior dictribution in multi-sources of prior information. Science Techno-logy and Engineering, 5(20): 1479–1481 (in Chinese)

[2] Chen H (2016). Research on Residual Life Prediction for Typical Satellite Platform Subsystem Based on Multi-source Information Fusion. Dissertation for the Master’s Degree. Changsha: National University of Defense Technology (in Chinese)

[3] Feng J, Dong C, Liu Q, Zhou J L (2004). Fusion of information from multiple sources based on adequacy measure in Bayesian analysis. Mini-Micro Systems, 25(7): 1356–1358

[4] Feng J, Liu Q, Zhou J L, Dong C (2003). Correlation information fusion method and application in reliability analysis. Journal of Systems Engineering and Electronics, 25(6): 682–684

[5] Feng J, Liu Q, Zhou J L, Dong C (2003). Multi-source information fusion based on MEMM in Bayes analysis. Quality and Reliability, (6): 31–34

[6] Feng J, Zhou J L (2008). Small-sample reliability information fusion approach based on Bayes-fuzzy logistic operator. Journal of Aerospace Power, 23(9): 1633–1636

[7] Feng J, Zhou J L, Sun Q (2006). Fusion of information of multiple sources based on ML-II theory in Bayesian analysis. Mathematics in Practice and Theory, 36(6)

[8] Gebraeel N, Elwany A, Pan J (2009). Residual life predictions in the absence of prior degradation knowledge. IEEE Transactions on Reliability, 58(1): 106–117

[9] Liu S, Chen H, Bo G, Jia X, Qi J (2017). Residual life estimation by fusing few failure lifetime and degradation data from real-time updating. IEEE International Conference on Software Quality, 177–184

[10] Liu S Q (2017). Residual Life of Satellite Platform System Fusing Multiple-souce Information. Dissertation for the Master’s Degree.Changsha: National University of Defense Technology (in Chinese)

[11] Padgett W J, Tomlinson M A (2004). Inference from accelerated degradation and failure data based on Gaussian process models. Lifetime Data Analysis, 10(2): 191–206

[12] Peng B H, Zhou J L, Jin G (2009). Reliability assessment of metallized film capacitor using multiple reliability information sources. High Power Laser and Particle Beams, 21(8): 1271–1275

[13] Pettit L I, Young K D S (1999). Bayesian analysis for inverse gaussian lifetime data with measures of degradation. Journal of Statistical Computation and Simulation, 63(3): 217–234

[14] Wang X L, Guo B, Cheng Z J (2012). Reliability assessment of products with wiener process degradation by fusing multiple information. Tien Tzu Hsueh Pao, 40(5): 977–982

[15] Wang X L, Jiang P, Xing Y Y, Guo B (2015). Residual Life Estimation for Nonlinear-Deterioriate Products. Beijing: National Defense Industry Press (in Chinese)

[16] Yang K, Yang G B (2011). Degradation reliability assessment using severe critical values. International Journal of Reliability Quality and Safety Engineering, 5(1): 85–95

[17] Zhang J H (2001). Accuracy detection method using bayesian multi-sensor data fusing technique. Journal of National Vniversity of Defense Technology, 23(3): 93–97 (in Chinese)

RIGHTS & PERMISSIONS

The Author(s) 2018. Published by Higher Education Press. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0)