A time-varying copula-based prognostics method for bivariate accelerated degradation testing

Abstract

Accelerated degradation testing (ADT) has been widely used to accelerate failure/degradation processes and to quickly evaluate the reliability and lifetime of products. In particular, the application of copula function provides a convenient and efficient way to model the ADT data of products that have two or more s-dependent degradation measures. However, little effort has focused on the pointwise infimum and supremum of the multivariate joint-distribution function. For this paper, a novel prognostics method was developed for bivariate ADT data on the basis of Brownian motion and time-varying copula method, which can estimate the pointwise best-possible bounds on bivariate joint reliable life function with a given measure of association, such as Kendall’s τ or Spearman’s ρ. The proposed model is applied to the real ADT data of microwave assembly to illustrate its performance and effectiveness.

Keywords

Prognostics ADT s-dependency time-varying copula reliable life bounds

1 Introduction

The development of design and manufacturing technology and the improvement of materials have greatly enhanced the reliability and quality of products. For these products to be highly reliable and to have a long lifespan, it has become far more difficult to obtain sufficient amounts of time-to-failure data through testing products under the normal operating environments and sometimes even under harsher conditions [4, 27]. The traditional failure-time data-based reliability analysis method is obviously not suitable.

However, degradation (e.g., wear, crack growth, corrosion, fatigue, shocks, material aging) is a common phenomenon for most systems or components [6]. The degradation-threshold failure mechanism provides a natural link between the degradation process and product failure time. One can obtain failure time distribution by analyzing the first passage time of the degradation process to a fixed threshold. On the basis of this rationale, the degradation information obtained from historical data or degradation tests can be used to evaluate reliability and lifetime of such products. The accelerated degradation testing (ADT) approach, especially, has been presented for collecting degradation data within a short time, in which test specimens are exposed to harsher-than-normal conditions to accelerate degradation [22].

In the past decade, the ADT method has become an effective tool for reliable verification and lifetime evaluation of products [30]. ADT is also widely used in the field of intelligent monitoring and prognostics, which can provide significant prior information for failure diagnostic and remaining useful life (RUL) prediction. Based on the ADT, the degradation process for degradation sensitive features can be monitored, which can be applied for the RUL prediction. Liu et al. [17] employed ADT data to establish a degradation model, which was used as a prior for the field RUL prediction of target product. Wang et al. [16] integrated the ADT data and field data using Bayesian method to predict the product’s actual field reliability. Ye et al. [25] conducted ADT to demonstrate that miller platform voltage can be applied to online condition monitoring of power MOSFET gate oxide degradation.

In ADT data analysis, a degradation model is assumed to describe the degradation paths of the samples tested at different stress levels. Therefore, the essence of ADT modeling is to develop appropriate probability models that can describe the degradation process. Generally speaking, existing degradation models are often classified into two large classes: stochastic process models and general path models [30]. Stochastic process models have attracted more attention because of their properties of time-dependent structures, which can capture the degradation dynamics of a system. Typical models include Brownian motion (i.e., Wiener process), the Gamma process, and the Inverse Gaussian (IG) process. The Gamma and IG processes are suitable for modeling degradation processes, which are always positive and strictly increasing, whereas, Brownian motion does not have this restriction. Hence, the Brownian motion model can provide more flexibility to describe various degradation paths and has become the most popular stochastic degradation model. Doksum and Hóyland [13] assumed that degradation of cable insulation follows the Brownian motion and that the failure data can be analyzed by using IG distribution. Ye and Xie [30] conducted a comprehensive literature review on Brownian motion on the basis of the degradation model and its variants, such as covariates, random effects, and measurement errors. Moreover, Liao and Elsayed [7] developed an ADT model on the basis of Brownian motion to predict the field reliability of LED when considering stress variations. Ye et al. [31] introduced the new random effects Wiener process model such that the degradation rate changed along with the degradation volatility of a unit.

The aforementioned studies consider only the cases that involve a single degradation measure. However, modern engineering products usually have complex structures with different functions, which means that the products might have multiple degradation failure mechanisms. As a result, a product might have two or more degradation measures, and any of them could cause product failure. Therefore, it is necessary to develop a multivariate, or at least a bivariate, degradation model to estimate the reliability of products accurately. Some progress has been made on bivariate or multivariate degradation modeling [20 , 28]. These works considered either the s-independent assumption of the multiple degradation measures, or multivariate normal distribution, to address the degradation analysis of systems with multiple performance indicators. However, the s-independent assumptions cannot match the engineering practice, and multivariate normal distribution might not be suitable for all conditions [8].

Copulas provide an alternative way to estimate the joint distribution of multivariate data, which model separately the marginal distributions and the joint-dependence structure [18]. Moreover, copulas have no constraints on univariate marginal distributions. Because of these advantages, the copulas method has attracted increasing interest in ADT analysis. Sari et al. [11] used a copula function to describe the correlation between two performance characteristics of LED and combined it with a generalized linear regression model for bivariate degradation data. Similarly, Zhang et al. [12], Rodríguezpicón et al. [15], Peng et al. [23], Zhang et al. [24], Pan et al. [29] also used the copulas method to model the dependency between bivariate degradation measures. Sun et al. [5] proposed a new s-dependent nonlinear ADT model for products with multiple performance parameters by using the general Wiener process and copulas. However, previous studies are limited to using copulas with constant parameters, which are not appropriate in many scenarios, because the degradation process is time-dependent and the degradation dependency is also time-variant. There are two methods to capture the time-variant dependence between the degradation processes, including time functional form of copulas [10] and fixed copulas with time-variant parameters [1]. The second method is more popular for its convenient application. Patton et al. [1] proposed several time-varying copulas, whose evolution equation describes the time-varying of copula parameters and follows auto-regressive moving average ARMA(1,10). Wang and Pham [26] developed an s-dependent competing risk model for systems subjecting to multiple degradation processes and random shocks by using time-varying copulas.

Besides, most of the previous work can only calculate the point estimator of system reliability for products with two or more degradation measures, and little effort has focused on the lower and upper bounds for the bivariate joint copula distribution. From this analysis, it is clear that the reliability evaluation of the system with two degradation measures has not been studied thoroughly. In this paper, a modified bivariate ADT model on the basis of Brownian motion and copula function is developed, which can estimate the joint-reliability bounds of a system in terms of a given measure of association, such as Kendall’s τ or Spearman’s ρ. The remainder of this paper is organized as follows. In the Section 2, the univariate ADT modeling method, on the basis of Brownian motion, is discussed in detail. The Section 3 introduces the bivariate time-varying copulas and measures of association. The Section 4 elaborates on the time-varying copula-based bivariate ADT modeling method, including the basic assumption, the bivariate dependent accelerated degradation model, and the parameter estimation method. The Section 5 section illustrates the validity of the proposed model and is followed by the Section 6.

2 Univariate accelerated degradation model

In engineering applications, two ADT methods have been widely used: constant-stress ADT (CSADT) and step-stress ADT (SSADT). CSADT divides the samples into several groups. Each group is tested under a certain accelerated stress level. In SSADT, all specimens are tested under the accelerated stress increasing step by step. For illustration purposes, this paper focuses on the bivariate CSADT model for a system with two performance parameters, which is the basis for studying more complex models.

Among the many stochastic process models, Brownian motion is the most widely used in degradation modeling and analysis.

Definition 1. The Brownian motion Y (t) is characterized by the following properties [3]:

Y (t) has independent increments in disjoint time intervals.

The increments are Gaussian random variables: ∀s, t > 0, Y (s + t) - Y (s) is normally distributed with the expected value 0 and variance σ²t, namely Y (s + t) - Y (s) ∼ N (0, σ²t); where σ is the diffusion coefficient.

The motion Y (t) is almost surely continuous in t.

When σ = 1, Y (t) is called standard Brownian motion, denoted by B (t).

There are several reasons to use the Brownian motion as a degradation model [30]. From the physical point of view, for many products, the degradation increment might be approximately and normally distributed because of the central limit theorem. External effects, such as shocks, are often independent, and the resulting degradation increments in disjoint time intervals are also independent. In this regard, Brownian motion, whose increment is normally distributed, is a good model for the degradation.

Therefore, Brownian motion is adopted to model the degradation process Y (t) under both the normal stress level S₀, which refers to the normal operational condition, and the accelerated stress levels S₁ < S₂ < … < S_K, where K is the total number of stress levels. $Y (t) = σ B (t) + μ \cdot t + y_{0}$ (1) where t is time; B (t) is a standard Brownian motion; y₀ is the initial degradation level at time zero; μ > 0 is the drift coefficient; and σ > 0 is the diffusion coefficient describing the variability attributable to unit-to-unit variation and the variation of operating and environmental conditions. Because of the degradation modes and the mechanism constant principle in ADT, the assumption is that the diffusion coefficient σ does not change with stress, i.e., σ₀ = σ₁ = … = σ_K = σ.

To extrapolate the lifetime or degradation rate of a product at the use conditions, some specific parameters of the degradation model should be assumed to be stress-related first, as described by a given stress-acceleration model. In Model (1), The drift coefficient μ, commonly known as the degradation rate, is regarded as a function of stress levels. Usually, such acceleration relationships are obtained from physical mechanism analysis or empirical experiences [14], e.g., the temperature-Arrhenius model or the voltage-Eyring model. After proper transformation, a stress-acceleration model for the degradation rate μ of the product under stress level S can be expressed as following general log-linear form [31]: $ln μ = ln d (S) = a + b \cdot φ (S)$ (2) where a and b are unknown parameters; d (S) is the degradation rate under stress level S; and φ (S) is a function of stress level depending on the stress type, e.g., φ (S_i) =1/S_i for absolute temperature and φ (S_i) = ln S_i for voltage.

According to the independent increment property of Brownian motion, the degradation increment Δy over the non-overlapped interval Δt is normally distributed with the mean μΔt and variance σ²Δt. Its probability density function (PDF) is $f (Δ y) = \frac{1}{\sqrt{2 π σ^{2} Δ t}} exp {- \frac{{[Δ y - μ \cdot Δ t]}^{2}}{2 σ^{2} Δ t}}$ (3)

Generally, a product has failed if its performance degradation exceeds a critical threshold for the first time. In other words, the first passage time (FPT) distribution is used as the lifetime distribution of the product. For a given critical threshold D, the lifetime T_D of the product is the instant when the degradation process Y (t) exceeds D for the first time, i.e.: $T_{D} = inf {t > 0 | Y (t) \geq D}$ (4)

It is well known that the first passage time of Brownian motion to a fixed threshold follows an inverse Gaussian distribution with the following PDF [13]:

$\begin{matrix} f (t; y_{0}, D) \\ = \frac{D - y_{0}}{σ \sqrt{2 π t^{3}}} exp {- \frac{{[(D - y_{0}) - μ \cdot t]}^{2}}{2 σ^{2} t}} \end{matrix}$ (5)

Combined with Equation (2), the associated reliability function under normal operational condition S₀ can be expressed as

$\begin{matrix} R (t) & = & Φ [\frac{D - y_{0} - d (S_{0}) t}{σ \sqrt{t}}] \\ - exp (\frac{2 d (S_{0}) (D - y_{0})}{σ^{2}}) \\ Φ [- \frac{D - y_{0} + d (S_{0}) t}{σ \sqrt{t}}] \end{matrix}$ (6) where Φ (·) is the cumulative distribution function (CDF) of standard normal distribution.

Furthermore, the reliable life derived from Equation (6) can be used for maintenance decision making or for verifying the lifetime and reliability levels of the tested equipment.

3 Time-varying copula method

3.1 Bivariate time-varying copula

Copula function provides a convenient approach to identifying multivariate probability distribution by separating the learning of the marginal distributions from the learning of the multivariate dependence structure that links the marginal distributions to form a joint distribution [5].

On the basis of Sklar’s theorem [2], a bivariate copula function can be defined as follows.

Definition 2. A two-dimensional copula (or briefly, a copula) is a function C from I² to I = [0, 1] with the following properties [18]:

Boundary conditions: for every u in I, $C (u, 0) = C (0, u) = 0$ (7) and $C (u, 1) = C (1, u) = u$ (8)

Which is 2-increasing on I², i.e., for every u₁, u₂, v₁, v₂ in I such that u₁ ≤ u₂ and v₁ ≤ v₂,

$\begin{matrix} C (u_{2}, v_{2}) - C (u_{2}, v_{1}) \\ - C (u_{1}, v_{2}) + C (u_{1}, v_{1}) \geq 0 \end{matrix}$ (9)

Furthermore, if c (u, v) = ∂²C (u, v)/∂u∂v is the density of copula function C, the PDF corresponding to joint distribution function H (x, y) can be calculated as $h (x, y) = c (F (x), G (y)) \cdot f (x) \cdot g (y)$ (10) where f (x) and g (y) are the PDFs of margins F (x) and G (y), respectively.

Generally, the parameter of the copula function, θ, is a constant, i.e., the strength of dependence is unchanged throughout the degradation process. To capture the time-variance of the dependence, time-varying copula assumes θ to be a function of time, called dynamic evolution equation. In this study, the following ARMA(1, 10) process is employed to be the dynamic evolution equation for θ_t: $θ_{t} = Λ (ω + β θ_{t - 1} + α \frac{1}{10} \sum_{i = 1}^{10} | u_{t - i} - v_{t - i} |)$ (11) where Λ (·) is a modified logistic transformation function to ensure that θ_t is always within the domain. $\frac{1}{10} \sum_{i = 1}^{10} | u_{t - i} - v_{t - i} |$ is used to capture the variability in the dependence.

Table 1

Several time-varying Copulas

Time-Varying Copulas	C (u, v)	Transformation function Λ (x)
Time-Varying Normal	$\int_{- \infty}^{Φ^{- 1} (u)} \int_{- \infty}^{Φ^{- 1} (v)} \frac{1}{2 π \sqrt{1 - ρ_{t}^{2}}} exp {- \frac{s^{2} - 2 ρ_{t} sr + r^{2}}{2 (1 - ρ_{t}^{2})}} dsdr$	(1 - e^-x)/(1 + e^-x)
Time-Varying Rotated Gumbel	exp{ - ((- ln(1 - u)) ^{ρ _t} + (- ln(1 - v)) ^{ρ _t}) ^1/ρ_t }	1 + x²
Time-Varying SJC	$0.5 \times (C_{JC} (u, v; τ_{t}^{U}, τ_{t}^{L}) + C_{JC} (1 - u, 1 - v; τ_{t}^{U}, τ_{t}^{L}) + u + v - 1)$	(1 + e^-x) ^-1

Where, Φ is the distribution function of a standard normal random variable; $C_{J C} (u, v; τ_{t}^{U}, τ_{t}^{L}) = 1 - {(1 - {{[1 - {(1 - u)}^{κ t}]}^{- y t} + 1 {[1 - {(1 - v)}^{κ t}]}^{- y t} - 1}^{- 1 / y t})}^{1 / κ t}, τ_{t}^{U}$ and $τ_{t}^{L}$ are the upper and lower tail correlations of Time-Varying SJC copula respectively, and $κ_{t} = 1 / \log_{2} (2 - τ_{t}^{U}), γ_{t} = 1 / \log_{2} (2 - τ_{t}^{L})$ .

Table 1 lists three time-varying copulas developed by Patton [1], where u and v denote margins.

3.2 Measures of association

The most widely known measures of dependence or association between two random variables are Kendall’s τ and Spearman’s ρ. When these measures are applied to the bivariate degradation processes, the Kendall’s τ and Spearman’s ρ with a given copula function can be performed as [18] $τ = 4 \iint_{I^{2}} C (u, v) dC (u, v) - 1$ (12) $ρ = 12 \iint_{I^{2}} C (u, v) dudv - 3$ (13) where u and v are the CDFs of margins.

Although both Kendall’s τ and Spearman’s ρ measure the probability of concordance between random variables with a given copula, the values of τ and ρ are often quite different. The relationship between can be stated as $- 1 \leq 3 τ - 2 ρ \leq 1$ (14) $\frac{1 + ρ}{2} \geq {(\frac{1 + τ}{2})}^{2}, \frac{1 - ρ}{2} \geq {(\frac{1 - τ}{2})}^{2}$ (15)

3.3 Selection among copula families

The classical Akaike information criterion (AIC) and the Bayesian information criterion (BIC) are used to compute the goodness-of-fit score for each time-varying copula under consideration and to decide which of the time-varying copulas best fits the given data.

The AIC is defined as $AIC = 2 q - 2 ln (L)$ (16)

The BIC is defined as $BIC = ln (n) \times q - 2 ln (L)$ (17)

Where q denotes the number of unknown parameters in the model; and ln(L) is the maximum value of the log-likelihood function and n represents the sample size.

4 Bivariate dependent ADT modeling method

4.1 Assumptions

To analyze the bivariate CSADT data, the following assumptions are made:

The test specimens are independent, and no catastrophic failures occur during the test.

The degradation observations of all specimens are collected at the same predetermined times.

A specimen is considered to have failed if one of the degradation measures reaches its corresponding failure threshold for the first time.

The degradation processes of all performance parameters in CSADT can be depicted by a Brownian motion with drift.

The two degradation measures are not independent, and the dependency between two performance parameters can be characterized by a time-varying copula.

4.2 System reliability model

If a product has two performance parameters and the degradation process under operating condition S₀ are random variables Y (t) = (Y₁ (t), Y₂ (t)), according to the aforementioned assumptions, the specimen is considered to have failed if any one of degradation measures reaches its corresponding failure threshold for the first time, which is known as D = (D₁, D₂). If the failure time of the kth performance parameter is T_k, and the lifetime of the system is T = min(T₁, T₂), then system reliability under normal condition S₀ can be expressed as

$\begin{matrix} R (t) & = & P (T > t) = P (T_{1} > t, T_{2} > t) \\ = & P (Y_{1} (t) < D_{1}, Y_{2} (t) < D_{2}) \end{matrix}$ (18)

If the two degradation measures are assumed to be independent, Equation (18) can be rewritten as

$\begin{matrix} R (t) & = & P (Y_{1} (t) < D_{1}) \times P (Y_{2} (t) < D_{2}) \\ = & R_{1} (t) \times R_{2} (t) \end{matrix}$ (19) where R_k (t), k = 1, 2 is the reliability function of kth performance parameter under normal operating condition S₀, as defined in Equation (6).

However, the various degradation measures are usually not independent of each other in engineering practice. In this case, Equation (19) will underestimate the reliability of the product. To model and measure the dependent structure between the two degradation measures, a copula function can be utilized in this paper.

If F_k (t) =1 - R_k (t) is the CDF of kth of the lifetime of performance parameter T_k under normal operating condition S₀, then, in accordance with

Sklar’s theorem and the definition of copula, the system joint CDF F (t) of T₁ and T₂ at time t can be calculated by $F (t) = P (T_{1} \leq t, T_{2} \leq t) = C (F_{1} (t), F_{2} (t); θ_{t})$ (20) where, θ_t is the parameter of a time-varying copula.

Therefore, the reliability of the product at time t under normal operating condition S₀ can be stated as [11 , 29]:

$\begin{matrix} R (t) & = & P (T_{1} > t, T_{2} > t) \\ = & 1 - P (T_{1} \leq t) - P (T_{2} \leq t) \\ + P (T_{1} \leq t, T_{2} \leq t) \\ = & 1 - F_{1} (t) - F_{2} (t) + C (F_{1} (t), F_{2} (t); θ_{t}) \\ = & R_{1} (t) + R_{2} (t) - 1 + C (F_{1} (t), F_{2} (t); θ_{t}) \end{matrix}$ (21)

4.3 Bounds on bivariate reliability functions

Fréchet and Hoeffding presented the fundamental best-possible bounds inequality for bivariate distribution functions with given margins [19]. For all u, v in I, the Fréchet-Hoeffding bounds for a bivariate joint copula distribution are defined as [18] $max (u + v - 1, 0) \leq C (u, v) \leq min (u, v)$ (22) where u and v are the CDFs of marginal random variables.

Furthermore, Nelsen et al. [19] improved the best-possible bounds on the bivariate copula distribution of continuous random variables with given margins and measures of association and derived the further narrowed pointwise infimum and supremum in terms of Kendall’s τ and Spearman’s ρ, as follows.

Theorem 2.When Kendall’sτis given, if ${\underline{K}}_{τ}$ and ${\bar{K}}_{τ}$ denote the pointwise infimum and supremum of a copula function withτin [−1,1], then for any (u, v) inI², ${\underline{K}}_{τ} - (u, v) = \max (0, u + v - 1, \frac{1}{2} [(u + v) - \sqrt{{(u - v)}^{2} + 1 - τ}])$ (23) ${\bar{K}}_{τ} (u, v) = min (u, v, \frac{1}{2} [(u + v - 1) + \sqrt{{(u + v - 1)}^{2} + 1 + τ}])$ (24)

Consequently, if X and Y are continuous random variables with joint distribution function H (x, y) and marginal distributions F (x) and G (y), respectively, then the best-possible bounds for H (x, y) are

${\bar{K}}_{τ} - (F (x), G (y)) \leq H (x, y) \leq {\bar{K}}_{τ} (F (x), G (y))$ (25)

Theorem 3.When Spearman’sρis given, if ${\underline{P}}_{ρ}$ and ${\bar{P}}_{ρ}$ denote the pointwise infimum and supremum of a copula function withρin [−1,1], then for any (u, v) inI²,

$P {0.3 em -}_{0.3 em ρ} (u, v) = max (0, u + v - 1, \frac{u + v}{2} - φ (u - v, 1 - ρ))$ (26) ${\bar{P}}_{ρ} (u, v) = min (u, v, \frac{u + v - 1}{2} + φ (u + v - 1, 1 + ρ))$ (27) where,

$\begin{matrix} φ (a, b) & = & \frac{1}{6} [{(9 b + 3 \sqrt{9 b^{2} - 3 a^{6}})}^{1 / 3} \\ + {(9 b - 3 \sqrt{9 b^{2} - 3 a^{6}})}^{1 / 3}] . \end{matrix}$

Therefore, the best-possible bounds for the joint copula distribution H (x, y) are

$P {0.3 em -}_{0.3 em ρ} (F (x), G (y)) \leq H (x, y) \leq {\bar{P}}_{ρ} (F (x), G (y))$ (28)

The proofs of Theorems 2 and 3 can be found in Nelsen et al. [19]. They also point out that the new infimum with a positive Kendall’s τ can improve the Fréchet-Hoeffding lower bound, whereas the new supremum ${\bar{K}}_{τ}$ , with a negative Kendall’s τ, is an improvement of upper bound. For the case of Spearman’s ρ, the Fréchet-Hoeffding lower bound will be improved when ρ ≥ -1/2, and the upper bound will be improved when ρ ≤ 1/2.

Therefore, if Kendall’s τ is known, the upper and lower bounds for system joint copula reliability (21) at time t under normal operating condition S₀ can be represented as

$\begin{matrix} R_{1} (t) + R_{2} (t) - 1 + K {0.3 em -}_{0.3 em τ} (F_{1} (t), F_{2} (t)) \leq R (t) \\ \leq R_{1} (t) + R_{2} (t) - 1 + {\bar{K}}_{τ} (F_{1} (t), F_{2} (t)) \end{matrix}$ (29)

Whereas, the pointwise infimum and supremum of system reliability (21) with Spearman’s ρ can be given by

$\begin{matrix} R_{1} (t) + R_{2} (t) - 1 + P {0.3 em -}_{0.3 em ρ} (F_{1} (t), F_{2} (t)) \leq R (t) \\ \leq R_{1} (t) + R_{2} (t) - 1 + {\bar{P}}_{ρ} (F_{1} (t), F_{2} (t)) \end{matrix}$ (30)

4.4 Parameter estimation

After obtaining the data from ADT, statistical inferences should be conducted to obtain estimates of the unknown parameters in the bivariate accelerated degradation model. Then, considering the given use conditions, system reliability assessment and lifetime evaluation can be performed with the estimated parameters.

As previously mentioned, the bivariate joint distribution (20) is specified by two margins with CDF F_k (t) and PDF f_k (t), k = 1, 2, and a copula function C with parameter set θ_t. If z _k = (a_k, b_k, σ_k) is the parameter set in the ADT model of kth performance parameter, then the joint PDF of copula-based bivariate distribution in Equation (20) can be expressed as

$\begin{matrix} f (t; z_{1}, z_{2}, θ_{t}) \\ = c (F_{1} (t; z_{1}), F_{2} (t; z_{2}); θ_{t}) \times \prod_{k = 1}^{2} f_{k} (t; z_{k}) \end{matrix}$ (31) where c(θ_t) is the density of the time-varying copula function C(θ_t) with dynamic evolution equation θ_t.

For a sample of size n, the full log-likelihood function can be expressed as

$\begin{matrix} ln L (z_{1}, z_{2}, θ_{t}) \\ = \sum_{i = 1}^{n} ln c (F_{1} (t; z_{1}), F_{2} (t; z_{2}); θ_{t}) ︸ dependence \begin{matrix} \end{matrix} structure \begin{matrix} \end{matrix} L_{C} \\ + \sum_{k = 1}^{2} \sum_{i = 1}^{n} ln f_{k} (t; z_{k}) ︸ m arg inals \begin{matrix} \end{matrix} \sum_{k = 1}^{2} L_{k} \end{matrix}$ (32)

The maximum likelihood estimates (MLEs) of the parameters (z₁, z₂, θ_t) can be obtained by maximizing Equation (32). However, as the dimension of parameters increases, it is difficult to obtain the optimal solution. Joe [9] proposed the inference functions for margins (IFM) method as a computationally attractive alternative to estimate the parameters in multivariate copula models. Technically, Equation (32) is first deconstructed into two parts. The first part, denoted by L_C, is the contribution from the dependence structure in the data, and the other part is the contributions from each margin, L_k, k = 1, 2. Then, the IFM method separately estimates the parameters of marginal distributions and the parameters of the copula in two stages. In the first stage, the estimates of the marginal distribution parameters are obtained from the log-likelihood L_k of each margin. In the second stage, the copula parameter estimate is computed by maximizing L_C on the basis of the marginal distribution parameters estimated in the first stage. In this paper, this method was used to estimate the parameters of the proposed bivariate ADT model.

Stage 1. Marginal distribution parameters estimation

Suppose that the sample size of a CSADT with K stress levels is n, and there are n_l specimens under stress level S_l. During the CSADT, all specimens are measured once every Δt time interval, and there are M_l inspections under S_l. Then, the observation of the kth performance parameter at time t_lij is y_k (t_lij), k = 1, 2; l = 1, …, K; i = 1, …, n_l; j = 1, …, M_l, where t_lij is the time of the jth measurement of the ith unit under the lth stress level. According to Equation (3), the likelihood function of the kth degradation measure is given by

$\begin{matrix} L_{k} = \prod_{l = 1}^{K} \prod_{i = 1}^{n_{l}} \prod_{j = 1}^{M_{l} - 1} \frac{1}{\sqrt{2 π σ_{k}^{2} Δ t}} \\ \times exp {- \frac{{[(y_{k} (t_{li (j + 1)}) - y_{k} (t_{lij})) - μ_{kl} \cdot Δ t]}^{2}}{2 σ_{k}^{2} Δ t}}, \\ k = 1, 2 \end{matrix}$ (33)

Considering the acceleration model given by Equation (2), μ_kl can be expressed as μ_kl = d_k (S_l) = exp(a_k + b_kφ (S_l)); where a_k and b_k are unknown parameters of the kth degradation measure. Therefore, the model parameter set of the degradation measure is z _k = (a_k, b_k, σ_k) and k = 1, 2. If Δy_klij = y_k (t_li(j+1)) - y_k (t_lij), then Equation (33) can be expressed as

$\begin{matrix} ln L = - \frac{1}{2} \sum_{l = 1}^{K} \sum_{i = 1}^{n_{l}} \sum_{j = 1}^{M_{l} - 1} {ln (2 π σ_{k}^{2} Δ t) \\ + \frac{{[Δ y_{klij} - exp (a_{k} + b_{k} φ (S_{l})) \cdot Δ t]}^{2}}{σ_{k}^{2} Δ t}}, \\ k = 1, 2 \end{matrix}$ (34)

Then, the first partial derivatives for each parameter are calculated and set to be 0. The parameter set z _k = (a_k, b_k, σ_k) of each marginal distribution can be estimated using MLE method.

Stage 2. Time-varying copula parameter estimation

By replacing the marginal parameters z _k, k = 1, 2, with their first-stage estimate ${\hat{Z}}_{k}$ , the time-varying copula parameters θ_t = (ω, β, α) can be estimated by maximizing the following log-likelihood function:

${\hat{θ}}_{t} = arg max \sum_{k = 1}^{2} \sum_{i = 1}^{n} ln c (F_{1} (t; {\hat{z}}_{1}), F_{2} (t; {\hat{z}}_{2}); θ_{t})$ (35)

The time-varying copulas contain more parameters than constant copulas and the log-likelihood function may take a non-convex form. Therefore, we adopt the Matlab copula toolbox developed by Andrew Patton [1] to estimate the time-varying copula parameters θ_t.

In addition, the AIC and BIC methods are used to quantitatively select the best-fitting model from candidate copula functions.

Fig.1

The accelerated degradation data.

Fig.3

Marginal reliability functions and CDF for power and noise of microwave assembly.

Fig.2

Scatter plots of the two degradation processes of a specimen.

5 Case study: Microwave assembly

Here, a case study demonstrates the effectiveness of the proposed method in real-world industrial applications.

Table 2
Estimators of univariate ADT model parameters

Parameters $\hat{a}$ $\hat{b}$ μ (S₀) $\hat{σ}$

Power 6.8739 −4639.6119 1.2929×10^–4 0.0876

Noise 18.2298 −9230.4836 1.7485×10^–6 0.0824

Parameters	$\hat{a}$	$\hat{b}$	μ (S₀)	$\hat{σ}$
Power	6.8739	−4639.6119	1.2929×10^–4	0.0876
Noise	18.2298	−9230.4836	1.7485×10^–6	0.0824

The microwave assembly is used to transmit/receive signals with high reliability and a long lifespan. Fault-mechanism analysis and historical information show that temperature is the major variable influencing the performance degradation of microwave assembly. To evaluate the lifetime and reliability of microwave assembly, a temperature CSADT with three accelerated stress levels was conducted using a thermal test chamber. And, two key performance parameters, power and noise, were measured. The normal operational condition is S₀ = 25°C, and the accelerated stress levels are S₁ = 70°C; S₂ = 80°C; and S₃ = 85°C. The sample size under each stress level is 3. The power and noise of microwave assembly were measured every hour by a computerized measuring system, and the monitoring time in accelerated stress levels S₁, S₂, and S₃ were 174, 155, and 114 hours, respectively. The accelerated degradation processes of power and noise are illustrated in Fig. 1.

Scatter plots of the two degradation processes of a specimen are depicted in Fig. 2, which visualize the dependent relationship between power and noise. The degradation data are normalized into the range [0, 1] prior to scatter plotting. Figure 2 shows that the two degradation measures are obvious correlated.

The proposed bivariate ADT modeling method can address different dependence situations through the copula function. Consequently, the proposed method is implemented to model and analyze the CSADT data of microwave assembly directly.

First, the Brownian motion model is used to establish the univariate ADT models for power and noise, respectively. The model parameters were estimated by using the first-stage procedure of the IFM method, and the results are presented in Table 2. By considering that the relative failure threshold D - y₀ of power and noise are 12 and 9, respectively, the reliability and CDF functions under normal operational condition S₀= 25°C for two marginal degradation processes are obtained, as shown in Fig. 3.

Table 3

Parameter estimation and goodness-of-fit for candidate copulas

Copulas	Parameter(s) Estimation	AIC	BIC	Rank
Time-varying Normal	[4.9999, −0.4624, 0.5337]	−31223	−31223	5
Time-varying Rotated Gumbel	[2.0732, 0.1035, −4.9999]	−40899	−40899	1
Time-varying SJC	Upper tail:[1.9097, −24.9999, 2.1328]	−38468	−38468	3
	Lower tail:[−10.9233, −9.6070, −18.8507]
Gaussian	0.9994	203268	203275	8
Student’s t	[0.9855, 22.4676]	−10	−35603	7
Frank	35.7014	−39081	−39074	2
Gumbel	8.3701	−36317	−36309	4
Clayton	10.7212	−23038	−23030	6

Fig.4

Time-vary Normal copula parameter estimation.

Fig.5

Time-vary Rotated Gumbel copula parameter estimation.

Fig.6

Time-vary SJC copula parameter estimation.

Then, three time-varying copulas are employed to describe the dependence between two marginal degradation processes, including the time-varying Normal copula, time-varying Rotated Gumbel copula, and time-varying SJC copula. The copula toolbox developed by Patton [1] is used to estimate the time-varying copula parameters θ_t = (ω, β, α). We also compared the time-varying copulas with some common used constant copulas. The parameter estimation results for the constant copulas and time-varying copulas are shown in Table 3. The parameter estimation curves for the time-varying copulas are illustrated in Figs. 4–6, where the solid lines indicate the dynamic evolution processes of time-varying copulas, and the dashed lines display the estimation of corresponding constant copulas.

Fig.7

Joint CDF of microwave assembly with different copulas.

To determine the most-suitably fitting copula function, AIC and BIC are used as the criteria. From Table 3, one can see that both the AIC and BIC values of time-varying Rotated Gumbel Copula are the smallest among constant copulas and time-varying copulas. Hence, it might be the best choice among the candidates.

Plugging the marginal CDF into the fitted copulas, we can obtain the CDFs of the joint distribution. Figure 7 displays the joint CDF curves of microwave assembly with the best five copulas (time-varying Rotated Gumbel, Frank, time-varying SJC, Gumbel, and time-varying Normal copula), and they are also compared with s-independent situation in the figure.

Next, the pointwise infimum and supremum of the joint CDF with Kendall’s τ and Spearman’s ρ of the microwave assembly can be obtained based on time-varying rotated Gumbel copula by using Equations 28), as shown in Fig. 8. As above-mentioned, Spearman’s ρ lower bound is a modified limit for the Fréchet-Hoeffding lower limit. Therefore, the upper and lower bounds, in terms of Spearman’s ρ, are adopted.

Fig.8

Upper and lower bounds of the joint CDF of the microwave assembly.

Finally, the system joint reliability and its upper and lower bounds of microwave assembly are calculated by Equations (21) and (30). Figure 9 shows the marginal reliability and system reliability functions under both s-dependent and s-independent assumptions. One can see that the reliability function, by assuming s-independent degradation features, tends to be much lower than that established on the basis of the s-dependent assumption. Clearly, dependence among the multiple features cannot be ignored. This conclusion could provide more scientific decision-making support for the life evaluation, maintenance scheduling, and warranty policy of the products with multiple performance parameters.

6 Conclusions

Fig.9

Comparison of system reliability with bounds under s-dependent assumption.

This paper investigates a time-varying copula-based reliability bounds modeling method for bivariate CSADT. To establish the bivariate dependent accelerated degradation model, the Brownian motion, whose first passage time to a fixed threshold follows an inverse Gaussian distribution, was used to characterize the degradation path of each performance parameter, and the time-varying copula method was adopted to model and measure the dependence structure between two degradation measures. On the basis of the best-possible bounds on the bivariate copula distribution of continuous random variables with given margins and measures of association, the pointwise infimum and supremum of the joint reliability function was derived for the systems with two performance parameters. To reduce the computational difficulty, the IFM method was used to estimate the unknown parameters of the proposed model. Then, the AIC and BIC were used as the criteria to determine the best-fitting copula family from several candidates. Finally, a case study on CSADT data of microwave assembly was implemented to demonstrate the usability and validity of the proposed method. The results reveal that this method can serve as an efficient tool to evaluate system reliability bounds under an s-dependent assumption.

Footnotes

Acknowledgments

This work was supported in part by the National Natural Science Foundation of China (Grant Nos. 61603018, 51775020 and 61104182) and the Fundamental Research Funds for the Central Universities (No. YWF-16-JCTD-A-02-06).

References

Patton

Modelling asymmetric exchange rate dependence, International Economic Review47 (2006), 527–556.

Sklar

, Fonctions de répartition à n dimensions et leursmarges, Publications de l’Institut de Statistique de L’Université de Paris8 (1959), 229–231.

Revuz

, Yor

, Continuous martingales and Brownianmotion, Springer Science & Business Media, New York, USA2013.

Zio

, Some challenges and opportunities in reliabilityengineering, IEEE Transactions on Reliability65 (2016), 1769–1782.

Sun

, Liu

, Li

and Liao

, Stochastic modeling and analysisof multiple nonlinear accelerated degradation processes throughinformation fusion, Sensors16 (2016), 1242.

Yang

, Environmental-stress-screening using degradationmeasurements, IEEE Transactions on Reliability51 (2002), 288–293.

Liao

and Elsayed

E.A.

, Reliability inference for field conditionsfrom accelerated degradation testing, Naval Research Logistics53 (2006), 576–587.

Kat

H.M.

, The dangers of using correlation to measure dependence, The Journal of Alternative Investments6 (2003), 54–58.

Joe

, Asymptotic efficiency of the two-stage estimation method forcopula-based models, Journal of Multivariate Analysis94 (2005), 401–419.

10.

Rodriguez

J.C.

, Measuring financial contagion: A copula approach, Journal of Empirical Finance14 (2007), 401–423.

11.

Sari

J.K.

, Newby

M.J.

, Brombacher

A.C.

and Tang

L.C.

, Bivariateconstant stress degradation model: LED lighting system reliabilityestimation with two-stage modelling, Quality and ReliabilityEngineering International25 (2009), 1067–1084.

12.

Zhang

, Ma

and Zhao

, A stress-strength time-varyingcorrelation interference model for structural re-liability analysisusing copulas, IEEE Transactions on Reliability66 (2017), 351–365.

13.

Doksum

K.A.

and Hóyland

, Models for variable-stressaccelerated life testing experiments based on Wiener processes andthe Inverse Gaussian distribution, Technometrics34 (1992), 74–82.

14.

Escobar

L.A.

and Meeker

W.Q.

, A review of accelerated test models, Statistical Science21 (2006), 552–577.

15.

Rodríguezpicón

L.A.

and Floresochoa

V.H.

, L.C.Mén-dezgonzález and M.A. Rodríguezmedina, Bivariatedegradation modelling with marginal heterogeneousstochastic processes, Journal ofStatistical Computation & Simulation (2017), 1–20.

16.

Wang

, Pan

, Li

and Jiang

, A Bayesian reliabilityevaluation method with integrated accelerated degradation testingand field information, Reliability Engineering & SystemSafety112 (2013), 38–47.

17.

Liu

, Li

X.-Y.

, Jiang

T.-M.

and Sun

F.-Q.

, Utilizing accelerateddegradation and field data for life prediction of highly reliableproducts, Quality and Reliability Engineering International32 (2016), 2281–2297.

18.

Nelsen

R.B.

An introduction to copulasSpringer Science & Business MediaNew York, USA, 2007.

19.

Nelsen

R.B.

, Quesada-Molina

J.J.

, Rodriíguez-Lallena

J.A.

and Úbeda-Flores

, Bounds on bivariate distribution functions withgiven margins and measures of association, Communications inStatistics - Theory and Methods30 (2001), 1055–1062.

20.

, Lu

and Kolarik

W.J.

, Multivariate performance reliabilityprediction in real-time, Reliability Engineering & SystemSafety72 (2001), 39–45.

21.

Huang

and Askin

R.G.

, Reliability analysis of electronic deviceswith multiple competing failure modes involving performance agingdegradation, Quality and Reliability Engineering International19 (2003), 241–254.

22.

Meeker

W.Q.

, Escobar

L.A.

and Lu

C.J.

, Accelerated degradationtests: Modeling and analysis, Technometrics40 (1998), 89–99.

23.

Peng

W.W.

, Li

Y.F.

, Yang

Y.J.

, Zhu

S.P.

and Huang

H.Z.

, Bivariateanalysis of incomplete degradation observations based on InverseGaussian processes and copulas, IEEE Transactions on Reliability65 (2016), 624–639.

24.

Zhang

X.P.

, Shang

J.Z.

, Chen

, Zhang

C.H.

and WangStatistical

Y.S.

, Statistical inference of accelerated life testing with dependentcompeting failures based on copula theory, IEEE Transactions on Reliability63 (2014), 764–780.

25.

X.R.

, Chen

, Wang

Y.X.

, Zhai

G.F.

and Vachtsevanos

G.J.

, Onlinecondition monitoring of power MOSFET gate oxide degradation based onmiller platform voltage, IEEE Transactions on Power Electronics32 (2017), 4776–4784.

26.

Wang

Y.P.

and Pham

, Modeling the dependent competing risks withmultiple degradation processes and random shock using time-varyingcopulas, IEEE Transactions on Reliability61 (2012), 13–22.

27.

Song

, Liu

, Yang

and Peng

, Data-driven hybrid remaininguseful life estimation approach for spacecraft lithium-ion battery, Microelectronics Reliability75 (2017), 142–153.

28.

Pan

Z.Q.

and Balakrishnan

, Reliability modeling of degradation ofproducts with multiple performance characteristics based on gammaprocesses, Reliability Engineering & System Safety96 (2011), 949–957.

29.

Pan

Z.Q.

, Balakrishnan

and Sun

Q.A.

, Bivariate constant-stressaccelerated degradation model and inference, Communications in Statistics - Simulation and Computation40 (2011), 259–269.

30.

Z.S.

and Xie

, Stochastic modelling and analysis of degradationfor highly reliable products, Applied Stochastic Models in Business and Industry31 (2015), 16–32.

31.

Z.S.

, Chen

and Shen

, A new class of Wiener process modelsfor degradation analysis, Reliability Engineering & System Safety139 (2015), 58–67.

A time-varying copula-based prognostics method for bivariate accelerated degradation testing

Abstract

Keywords

1 Introduction

2 Univariate accelerated degradation model

3.1 Bivariate time-varying copula

4.1 Assumptions

4.2 System reliability model

Table 2 Estimators of univariate ADT model parameters Parameters a ˆ b ˆ μ (S0) σ ˆ Power 6.8739 −4639.6119 1.2929×10–4 0.0876 Noise 18.2298 −9230.4836 1.7485×10–6 0.0824

Footnotes

Acknowledgments

References

Table 2
Estimators of univariate ADT model parameters

Parameters $\hat{a}$ $\hat{b}$ μ (S₀) $\hat{σ}$

Power 6.8739 −4639.6119 1.2929×10^–4 0.0876

Noise 18.2298 −9230.4836 1.7485×10^–6 0.0824