assessment foR systems with two peRfoRmance chaRacteRistics based on gamma pRocesses with maRginal heteRogeneous Random effects

In recent years, important developments have been presented in the area of reliability inference of products and systems based on degradation models. Such models are important tools to obtain reliability information when few failure data is available [11], and consists in analyzing the gradual deterioration in performance of a performance characteristic (PC), also known as degradation process, in terms of the accumulated damage over time [27]. For some products, a failure would be defined at a specified critical level of degradation, which means that the product may not stop working completely as in the case of hard failures, but be defined when the cumulative degradation path crosses the critical level of degradation; such failures are known as soft failures. In general, if a failure can be defined in terms of a specified critical level of degradation it is possible to obtain a reliability assessment based on degradation process models [12]. Based on this, a modeling approach for degradation processes may consist in relating the degradation over time with a continuous stochastic process such that it is possible to describe the failure generating mechanisms. For certain PC, the desirable properties that a model must have to describe its degradation process are that the degradation process should always be positive and strictly increasing. In this paper, the gamma process is considered as a model to govern the degradation process of certain PC, given the characteristics that its increments are independent and non-negative having a gamma distribution that results in an always positive, strictly increasing stochastic process. As performance can only decrease over time, this is why it is considered to be suitable to model wear, crack growth, corrosion, consumption, fatigue, erosion, or any PC [16]. Some important applications of the gamma process in the reliability assessment of products can be found in Bagdonavicius and Nikulin [2], Park and Padgett [22], and Bagdonavicius and Nikulin [3]. RodRíguez-Picón LA. Reliability assessment for systems with two performance characteristics based on gamma processes with marginal heterogeneous random effects. eksploatacja i niezawodnosc – Maintenance and Reliability 2017; 19 (1): 8–18, http://dx.doi.org/10.17531/ ein.2017.1.2.


Introduction
In recent years, important developments have been presented in the area of reliability inference of products and systems based on degradation models.Such models are important tools to obtain reliability information when few failure data is available [11], and consists in analyzing the gradual deterioration in performance of a performance characteristic (PC), also known as degradation process, in terms of the accumulated damage over time [27].For some products, a failure would be defined at a specified critical level of degradation, which means that the product may not stop working completely as in the case of hard failures, but be defined when the cumulative degradation path crosses the critical level of degradation; such failures are known as soft failures.In general, if a failure can be defined in terms of a specified critical level of degradation it is possible to obtain a reliability assessment based on degradation process models [12].Based on this, a modeling approach for degradation processes may consist in relating the degradation over time with a continuous stochastic process such that it is possible to describe the failure generating mechanisms.
For certain PC, the desirable properties that a model must have to describe its degradation process are that the degradation process should always be positive and strictly increasing.In this paper, the gamma process is considered as a model to govern the degradation process of certain PC, given the characteristics that its increments are independent and non-negative having a gamma distribution that results in an always positive, strictly increasing stochastic process.As performance can only decrease over time, this is why it is considered to be suitable to model wear, crack growth, corrosion, consumption, fatigue, erosion, or any PC [16].Some important applications of the gamma process in the reliability assessment of products can be found in Bagdonavicius and Nikulin [2], Park and Padgett [22], and Bagdonavicius and Nikulin [3].

Reliability assessment foR systems with two peRfoRmance chaRacteRistics based on gamma pRocesses with maRginal heteRogeneous Random effects ocena niezawodności systemów o dwóch paRametRach użytkowych opaRta na pRocesach gamma z bRzegowymi niejednoRodnymi efektami losowymi
In this paper, a reliability modeling approach for products with two performance characteristics related to two degradation processes is developed.The joint modeling of such processes is performed by using a copula function in order to consider the dependence structure between degradation processes.The proposed approach considers that different random effects affect the stochastic behavior of each performance characteristic.For such approach, different bivariate models with marginal gamma processes with heterogeneous random effects as marginal distributions are considered.As the random effects may differ between performance characteristics, different modifications of the structure of the parameters of the gamma process are proposed.Such that the random effects affect both the drift and diffusion, just the drift, and just the diffusion of the marginal gamma processes.The statistical inference of the joint bivariate models is performed via Bayesian approach.The obtained results show that a bivariate model with heterogeneous random effects has a slight better performance among the proposed models.Which implies that the bivariate heterogeneous random effects gamma process models may provide a better approach to model multivariate degradation data, and thus a better reliability assessment of the product under study.Keywords: degradation, gamma process, random effect, heterogeneous effect, copula function.

sciENcE aNd tEchNology
However, in most of the degradation processes, it can be found that the degradation of a product's characteristic is affected by different sources of variation.Which implies that the degradation in a product population has a large variation due to some unobservable effects.These effects are described by the variation in the degradation increments over time for every specific product's degradation path and the different behaviors of the degradation paths for every product, i.e., the degradation rate.The simple gamma process is unable to capture such variations.However, these variations have been well modeled by incorporating random effects into the gamma process.Lawless and Crowder [10] considered that the scale parameter of the gamma process have different realizations among products, and let the parameter to be random following a gamma distribution.Different applications of this model can be found in Tsai et al. [29], Wang et al. [30], Wang [32], and Pulcini [25].In the classical random effects gamma process model it can be noted that both the mean and variance of the process are affected by the random effect parameter, which results in a degradation process with random drift and random diffusion.Nevertheless, for some products it may be the case that the degradation paths are just characterized by a random drift or just a random diffusion.
Another important aspect to consider is that the functionality of a product may be related to multiple PC.In such cases, it is important to consider multivariate models to obtain robust reliability estimates.Some important applications of multivariate gamma processes in reliability analysis have been presented by Hao et al. [9], Pan and Balakrishnan [18], Pan and Balakrishnan [19], Pan et al. [20], Pan et al. [21], Park and Padgett [23], Sari et al. [26], Wang et al. [31], and Zhou et al. [36].In most of the multivariate degradation models, it is considered that each of the multiple PC are governed by a univariate stochastic gamma process and then the joint model is obtained via copula functions.By modeling with copula functions, the independence assumption is not assumed.Indeed it can be tested by considering the association parameter of the copula function [8], which makes the copula modeling approach quite attractive.In addition, most of the developed works do not consider random effects in the modeling [18-21, 26, 36].Random effects are an important aspect to consider in the degradation modeling of almost any product under study, given that most of the time there is a substantial subject-to-subject variability among the degradation processes of different individuals [33].Which, accounts to describe the individual variability that determines the heterogeneity among the degradation paths of different product units, also known as individual variability [35].
Considering that a product may consist in multiple PC and that random effects affect in different ways the multiple PC.It result important to develop multivariate models with heterogeneous random effects.Furthermore, each PC may experience different sources of variation, which means that the random effects may be characterized in different ways in the degradation paths.For example, the degradation paths of a PC 1 may be characterized by a large variation of the degradation rate and a low variation of the within degradation increments, and the degradation paths of a second PC may be characterized by both a large variation of the degradation rate and the within degradation increments.It is important to consider such scenarios when dealing with the degradation modeling.The importance relies on describing the heterogeneous behavior of the multiple PC degradation in terms of their best fitting stochastic gamma process with random effects to obtain robust reliability assessments.Hao et al. [9] and Wang et al. [31] developed bivariate gamma processes with random effects, however they only considered the classical gamma process with random effects as marginal distributions.The classical gamma process with random effects consider that the scale parameter is random, which means that the mean and variance of the gamma process are affected by the random effects parameter.This results in processes with random drift and random diffusion.In this paper, we model the degradation processes of two PC considering an Archimedean copula function and different gamma processes with heterogeneous random effects as marginal distributions.This, by proposing different modifications of the structure of the parameters of the marginal gamma processes, such that a random drift and a random diffusion gamma process models are obtained.A time-scale monotone transformation is considered to assure that degradation is a linear function of time [34].As the joint distributions are complex, the estimation of the parameters is performed via Gibbs sampling and Markov chain Monte Carlo (MCMC) implemented in OpenBUGS.The models are illustrated with the reliability assessment of a case study that consists of crack propagation data of two terminals of an electronic device.
The rest of the paper is organized as follows.Section 2 presents the simple gamma process and the bivariate modeling based on two PC.In Section 3, the different gamma processes with random effects for two PC are introduced.Section 4 presents the bivariate modeling based on the Frank copula function, and the heterogeneous random effects models are defined.Section 5 deals with inference method for the bivariate degradation model with random effects.Section 6 addresses the implementation of the proposed models in a fatigue-crack growth dataset.Finally, in Section 7 some concluding remarks are provided.

Gamma process for two performance characteristics
Considering a non-negative-valued process {Z(t),t>0}, where Z(t) represents the measured degradation for an individual unit at time t, then the gamma process has the following properties: (a) Z(t+∆t)−Z(t)=∆Z(t) follows a gamma distribution Ga(v[τ(t+∆t)−τ(t) ],u), and (b) Z(t) has independent increments, Z(t 4 )−Z(t 3 ) and Z(t 2 )−Z(t 1 ) are independent ∀t 1 <t 2 <t 3 <t 4 .Now, let v(τ(t)) be a non-negative shape parameter with a time scale transformation in the form of τ(t,γ)= t γ , thus ∆τ(t)=τ(t+∆t)−τ(t), t≥0, v(0)≡0, and u>0 be a scale parameter.Then, Z(t),t>0 is governed by a gamma process with the parameters described above.Thus, the gamma process Ga(v(∆τ(t)),u) describes the degradation level of some characteristic at time t, and has a mean v(∆τ(t))/u and variance v(∆τ(t)) ⁄ u 2 , the probability density function (PDF) of ∆Z(t) is given by: where ( ) is the gamma function.
Considering that the degradation process of a certain product is governed by a gamma process, a failure of the product is said to have occurred the moment when the degradation path Z(t) crosses a critical level of degradation ω.Known also as the first-passage time.Thus, the first-passage time is defined as T ω = inf{t ω :Z(t)≥ω}.The cumulative distribution function (CDF) of t ω con be obtained as where Γ(v(τ(tω )),ωu) is the upper incomplete gamma function Now consider that a product has two PC and that are marginally governed by a gamma process with a time-scale transformation.

sciENcE aNd tEchNology
Also considering a degradation test (DT) with the next characteristics: N units are tested and M measurements for all the units are observed up to the termination time T, which results in degradation measurements Z ik (t j ) of the ith unit at the corresponding time t j , i=1,2,… ,N, j=1,2 ,…,M, and k=1,2 PC.Then, the degradation data can be presented as follows: According to the independent increment property of the gamma process, and , ,, for i=1,2,…,N,j=1,2,…,M, and k=1,2 PC.
Thus, it is possible to obtain independent random variables ∆Z ik (t j ) ~ Ga(v k (∆τ(t j ,γ k )),u k ), with the next PDF and CDF: where is the lower incomplete gamma function, i=1,2,…,N,j=1,2,…,M, and k=1,2 PC.If ω k is the critical degradation level of each PC for k=1,2.
Then the CDF of k t ω can be obtained from (2) as P(Z ik (t j )≥ω k ).

Gamma processes with random effects
Considering that the sampled product under study may experience different sources of variation during its operation, it results appropriate to incorporate product-to-product variability in the modeling of the degradation processes.In this case, it is assumed that γ k and v k are fixed parameters that are common to all products.The random effects are introduced by letting that u k be a random parameter that follows a gamma distribution representing the heterogeneity among different products.Considering that, a product has two performance characteristics whose degradation have been observed during a DT with degradation measurements Z ik (t j ) of the ith unit at the corresponding time t j , i=1,2,…,N,j=1,2,…,M, and k=1,2 PC.Then, the gamma process with random effects (RE) for k=1,2 can be written as / 2 with PDF defined as follows: Thus, the PDF of the degradation increments ∆Z ik (t j ) of the ith unit at the corresponding time t j , i=1,2,…,N,j=1,2,…,M,k=1,2 is given by: According to Lawless and Crowder [10], for any fixed t, the ran- ) follows an F -distribution with 2v k ∆τ(t j ,γ k ) and 2δ k degrees of freedom, thus, the CDF of ∆Z ik (t j ) is defined as: Taking into account (8), then the CDF of the lifetime when any of the degradation paths of the k =1,2 PC reach the respective critical con be obtained as: Given that the means (v k (∆τ(t j ,γ k ))) ⁄ u k and variances (v k (∆τ(t j ,γ k ))) ⁄(u k 2 ) of the respective k=1,2 gamma processes are affected by the random effects parameters u k .It is expected that, the degradation rates of the degradation paths tend to have a larger dispersion in both processes.In addition, it is also expected that the variances of the degradation observations within each unit of both processes tend to be large.Figure 1 provides pseudo-random paths of the gamma degradation process under the RE model.However, it may be the case that for some products, both PC present only large dispersions of the degradation rates or that only present large variances of the degradation observations within each unit.
A modification of the parameters' structure of the gamma processes is proposed as Z ik (t j )~Ga(η k (∆τ (t j ,γ k )) u k ,u k ) with u k following a gamma distribution.The degradation means and variances of this model can be obtained as η k (∆τ (t j ,γ k )) and (η k (∆τ (t j ,γ k ))) ⁄u k , respectively.It can be noted that only the variances are affected by the random parameter u k .Thus, the gamma process with random diffusion (10) sciENcE aNd tEchNology The CDF of the lifetime when the degradation paths of any of the PC for k=1,2 reach a critical level ω k is defined as } and it is obtained in the same way as: As random effects are involved in the scale parameter u k of the RV model, it is expected that the variance of the degradation observations within each unit to be significant.However, a low level of variation in the degradation rates may be observed.Thus, the RV model is suitable for the degradation modeling of products for which overall degradation rate is low and a large unit-specific degradation variation exists.Figure 2 provides pseudo-random paths of the gamma degradation process under the RV model.
A second modification of the structure of the parameters of the gamma process is proposed as the degradation means and variances of this model are defined as ζ k (∆τ(t j ,γ k ))u k and ζ k ∆(τ(t j ,γ k )), respectively.Again, if u k follows a gamma distribution, only the means are affected by the random effects parameters.Thus, the random drift gamma process (RD) for k=1,2 PC results in with PDF described as: The CDF of the lifetime when the degradation paths of any of the degradation paths for k=1,2 reach a critical level ω k is defined as } and it is obtained in the same way as: Because random effects are considered in the parameters u k , the means of the RD gamma process varies to a certain level.It leads to a larger dispersion of the degradation rates, which is reflected in larger variations of the first-time passage distributions.Thus, this model is appropriate for the modeling of the degradation of products in which significant variation of the degradation rate within the products' samples is observed.In Figure 3, a set of pseudo-random paths of the gamma degradation process under RD model are illustrated.

Bivariate modeling with copula functions
Supposing that the degradation process of k=1,2 PC are observed for i=1,2,…,N units and j=1,2,…,M observations are recorded at specified times t_j for each degradation process.Now consider that, stochastic gamma processes govern the degradation process of each ) .In this paper, it is considered that the depend- ence structure between degradation processes k=1,2 is described by a copula function in the form H(x,y)=C{F(x),G(y)}. Where, F(x),G(y) are two marginal CDF and x,y are two random variables.If H(x,y) is a joint distribution with marginal distributions F(x) and G(y), then there exists a unique copula C for x,y [15].Thus, a bivariate copula is a CDF defined in [0,1] 2 with uniform marginal distributions [0,1].In this case, an Archimedean Frank copula is considered to model two degradation measurements Z ik (t j ) with k=1,2 PC, i=1,2,…,N,j=1,2,…,M as , where F(∆Z ik (t j )) is the CDF of the respective marginal stochastic process for k=1,2 and θ is the association parameter of the copula function.It should be noted that the distributions described in (7 -13) refer to the case where the degradation levels before t j are unknown.Instead, when it is assumed that the degradation levels measured before t j are known, the conditional CDF are considered as F(∆Z ik (t j )) | Z i1 (t 1 ), …,Z i1 (t( j-1 )) for k=1,2, i=1,2,…,N, and j=1,2,…,M.Such that the bivariate copula results in With the marginal CDF defined as: , ,Z t du ik j k 1 (14) where the conditional PDF of the random parameter is given by: The bivariate model is denoted as H(∆Z i1 (t j ), ∆Z i2 (t j ))=C(U ij1 , V ij2 ), where U ij1 =F(∆Z i1 (t j ) | Z i2 (t 1 ),…,Z i1 (t( j-1 )) and V ij2 =F(∆Z i2 (t j ) | Z i2 (t 1 ),…,Z i2 (t( j-1 )) are the CDF of the marginal gamma process.If a simple gamma process is considered then the CDF is represented by (5).If random effects models are considered in the marginal processes and the degradation levels measured before t j are known, then the respective CDF can be obtained from (14).If the degradation levels measured before t j are unknown then the CDF can be obtained from ( 10) and ( 12) for the RV and RD models, respectively.Parametric copula functions impose strong assumptions or restrictions on the correlation structure among CDF.This denotes the need to determine how a specific copula function is selected, and how the chosen copula function can be validated by data.Such a topic has been widely studied in the literature, for example, by comparing the sample empirical copula and theoretical copula.Practitioners may be referred to Durrleman et al. [5] and Melchiori [14].In this paper, the Frank copula function is considered as C(U ij1 ,V ij2 | θ), to illustrate the joint modeling of the proposed models with random effects for the two degradation processes.In the literature, it has been found that the Frank copula works well when dealing with systems with two PC and when dealing with crack propagation data as can be seen in the works of [9,31,36], as the case study presented in this paper.The Frank copula function is presented in (15).Then, the log-likelihood function of the set of parameters ϑ η ζ γ δ ϕ θ , , , , , , for k=1,2, and i=1,2,…, is described in (16). where The product is considered to have failed if any of the PC reach the critical degradation level ω k , k=1,2.Thus, the reliability function can be described as Different combinations of the marginal processes for the bivariate modeling are proposed, these combinations are presented in Table 1.The bivariate modeling implies that the dependence structure between the degradation processes of the two characteristics can be tested by using the association parameter (θ), which justifies the use of the copula function.The models B 1 ,B 2 , and B 2 consider bivariate gamma models with homogeneous random effects for the RE, RV, and RD models, respectively.However, another important aspect to consider is that the degradation of each characteristic may be different in terms of the random effects that affect the behavior of the degradation paths.This means that, e.g., the degradation paths of the characteristic 1 may exhibit a large variation of the degradation rate, and a low level of variation of the within degradation increments in the paths, while the characteristic 2 may exhibit a large variation in the within degradation increments in the paths and a low variation of the degradation rates.For such cases, it may be important to consider heterogeneous marginal models with random effects in the bivariate modeling.In this paper, two heterogeneous models are considered and denoted as B 4 and B 5 in Table 1, where in the model B 4 it is considered that the degradation of the characteristic 1 is governed by a gamma process with RD, while the degradation of the characteristic 2 is governed by a gamma process with RV.In the model B 5 , it is considered that U ij1 =F RE (∆Z i1 (t j ) ) and V ij2 =F RV (∆Z i2 (t j )).Further bivariate combinations of gamma processes with heterogeneous random effects can be considered depending on the observed behavior of the degradation paths of the degradation processes under study.

Parameters estimation
As the bivariate functions consider different combinations of gamma processes with random effects in terms of the scale parameter of the gamma process, it is easy to note that the joint bivariate distribution function for any of the models B i ,i=1,2,3,4,5 may result in a non-standard complex form.

Table 1. Bivariate models and bivariate heterogeneous Models
sciENcE aNd tEchNology However, the MCMC can be utilized to estimate the parameters ϑ i of interest from B i , ,i=1,2,3,4,5.The MCMC procedure consists in generating samples from the joint posterior distribution.In this case, Gibbs sampling algorithm is utilized to obtain such samples from the joint distribution.Generally, the algorithm consists in dividing the parameter vector into d subvectors, ϑ=(ϑ 1 ,…,ϑ d ), such that each iteration of the algorithm cycles through the subvectors of δ, drawing each subset conditional on the value of all vectors.This process can be seen as generating a realization of a Markov chain that is built from a set of base transition probabilities.When the base transition probabilities are applied in sequence the algorithm can be described as simulating a homogeneous Markov Chain ϑ (1) ,ϑ (2) ,ϑ (3) ,…,.It is assumed that the prior distributions for all the parameters of interest are non-informative and there is prior independence among the parameters of interest.The implementation of the Gibbs sampling algorithm for the estimation of the parameters is implemented using the OpenBUGS package software [13].Zeros trick is used in OpenBUGS given that the loglikelihood is not a standard distribution [17].Important information about this algorithm can be found in Gelfand and Smith [6], Casella and George [4], Smith and Roberts [28].
Considering that degradation measurements Z ik (t j ) have been observed of the ith unit at the corresponding time t j , i=1,2,…,N,j=1,2,… ,M for k=1,2 PC.Thus, ∆Z ik (t j )=Z ik (t j )−Z ik (t j−1 )) is the degradation increment at the interval time t j −t j−1 .Then, the degradation increments ∆Z ik (t j ) are independent random variables that follow a gamma distribution as Ga(v k (∆τ (t j ,γ k )), u k ) for the RE model, Ga(η k (∆τ (t j ,γ k )), u k ) for the RV model, and Ga(ζ k (∆τ (t j ,γ k )), u k ) for the RD model, in all cases with u k following a gamma distribution.It should be noted that, under a gamma process with random effects, the degradation increments ∆Z ik (t j ) are independent random variables only conditional to the value of (v k ,η k ,ζ k ) for any of the described random effects models.Different combinations of the bivariate models with different random effects as marginal distributions are considered as described in Section 4. Given that random effects are generated by in all the models, it is considered that u k,i follows the same prior distributions with hyper-parameters (δ k ,φ k ), which accounts for pooling the information of the random effects among the different degradation trajectories for k=1,2.The non-informative prior distributions for these hyper-parameters are considered as δ .In addition, a non-informative gamma prior distribution is considered for the shape parameters of the gamma processes ) , in order to avoid neg- ative values of the time-scale transformation.
The bivariate joint Frank copula function specified in (15) is implemented in OpenBUGS with the marginal CDF described in Table 1 as B i by considering the proposed prior distributions for every parameter of interest.

Case study
In the following, a case study is presented, which consisted in the fatigue-crack propagation of two cracks in two terminals of an electronic device.Each device has two terminals, whose function is to transfer a signal to a receptor.Some cracks may be present in both terminals.The propagation of the cracks to a certain critical length can lead to failure of the device given the inability of transferring the signal to the receptor.A DT based on vibration was carried out in order to study the propagation of cracks in the terminals of ten devices.As every device has two terminals, two sets of fatigue-crack growth data were obtained.The length increments of the cracks for both terminals were measured every 0.1 hundred thousand cycles until 0.9 hundred thousand cycles.The measurements for every crack were performed at the same measurement times, and by considering equally distanced inspection times.The increments of the degradation were measure by considering a vision system with special software applications to measure crack propagations.The obtained data are presented in Table 2, the units are in millimeters.In Figure 4, the crack degradation paths for the two terminals of every device are illustrated.
The propagation of the cracks can be seen as a degradation process and therefore a stochastic process.As the device counts with two terminals, the crack propagation of every terminal can be seen as degradation processes.Considering that the cracks are from different positions, it is important to assess the dependence in the cracks propagation of the terminals for every device.It is considered that the device has failed if the length of any of the two terminal' cracks cross the critical limit of 0.663 mm.Such critical level of degradation was obtained by considering the total width of the terminal defined by the customer of the product.In this case, if a crack length exceeds the total width of the terminal, a failure of the system can be obtained given that the inability of the product to transfer a signal the receptor.As can be seen from the cumulative degradation paths in Figure 4, there are several differences among the degradation paths of a single terminal and between terminals.In the degradation paths of terminal 1, a large variation of the degradation rate can be observed which may indicates that the random effects affect the degradation rate, and the opposite can be seen in the degradation paths of terminal 2, where the variation is rather small.From both terminals, it can be noted that a large variation of the within increments is observed, which may indicates that the random effects affect the diffusion of the degradation process of both terminals.From these behaviors, it can be noted that the random effects have different influences on the degradation processes as described above.The different bivariate models presented in Table 1 are considered, these models were estimated via Gibbs sampling by using OpenBUGS.

Estimation of Proposed Models
Considering the dataset described in Based on the Frank copula function (15) and the log-likelihood function in (16) with ϑ i ,i=1,2,3,4,5, the estimated parameters were obtained by using OpenBUGS.It is assumed that there are prior independence among the parameters of interest.Estimation of the parameters was performed by using the zeros trick in the developed algorithm to specify the bivariate log-likelihood function.Two sets of initial values are considered in the algorithm to assess the convergence of the parameters of interest with the Brooks-Gelman-Rubin (BGR) statistic [7].A total of 50,000 iterations were considered for burn-in and 100,000 were considered for estimation purposes.As two sets of initial values were determined for every parameter, the BGR statistic was calculated for the parameters of interest.In general, it was found that convergence is achieved in every parameter in accordance with the BGR graphs obtained from OpenBUGS.The estimations obtained for the parameters of every model B i , along with the standard deviation, Monte Carlo error and some percentiles are presented in Table 3.
Information criteria was used to select the best fitting random effects bivariate model.The Akaike Information Criterion (AIC) is used for such purpose, which is an appealing tool for model selection based on information.The model with the lowest value of AIC is considered as the best fitting model [1].This criterion is defined as AIC=−2×l(ϑ 1 )+2R; where, l(ϑ 1 ) is the evaluated log-likelihood function from (16) for any of the models i=1,2,3,4,5, for the Frank copula and R is the number of parameters.In addition, the Kendall coefficient was estimated from the copula parameter (θ), this coefficient provides a good alternative to measure the level of dependence between the marginal distributions from the copula.The Kendall coefficient can be obtained from the Frank copula as  4.
As can be noted from Table 4, the best fitting model is the heterogeneous model B 5 , given that it has the lowest value of AIC.This model considers a RE gamma process model for the terminal 1 and a RV gamma process model for the terminal 2. The result makes sense, given that as can be seen from the cumulative degradation paths in Figure 4, the paths of terminal 1 have large variations of the degradation rate and the within degradation increments.Such characteristics are described when the RE model is considered given that both the mean and variance of the gamma process are affected by the random effects parameter.On the other hand, the degradation paths of terminal 2 have a lower level of variation of the degradation rates, and a large variation of the within degradation increments, such characteristics are described when the RV model is considered given that only the variance of the gamma process is affected by the random effects parameter, as described in Section 3.However, it must be noted that the differences among the AIC values for the proposed models is slight.Furthermore, the B 3 model, which considers that both terminals 1 and 2 are governed by a gamma process with RD, is the model with the poorest performance.The level of dependence described by the Kendall coefficient from the bivariate RE-RV model is 0.0411, which may indicates a low level of dependence, however, it must be considered when dealing with the reliability assessment of the product as a whole.In addition, from the estimates in Table 3 it can be noted that in the model B 5 the credible interval of θ include the value 0. Hence, the degradation processes seem to be independent, because for θ→0 the marginal CDF are independent, which explains why the obtained Kendall coefficients are small.In Figure 5, the contour and density plots are presented for the best fitting heterogeneous model via Frank copula.variance for terminal 2 under RV in B 5 is 0.007221.This means that, as the terminal 1 has a larger variance, some of the degradation paths will have large first-passage times, and for terminal 2 the first-time passage times may occur earlier with more probability.In addition, the reliability functions for the dependent and independent scenarios constantly overlap, which can be explained by the fact the considered level of dependence in terms of the Kendall coefficient is quite small 0.0411.

Concluding remarks
This paper introduces a bivariate degradation modeling approach for two PC based on marginal gamma processes with heterogeneous random effects.As the random effects may affect in different ways the degradation of each PC, different modifications of the structure of the parameters of the marginal gamma processes were performed.Such modifications are considered in the bivariate modeling, which allows to describe the heterogeneous behaviors of the degradation processes due to random effects.Different combinations of marginal gamma processes with random effects joined by a Frank copula were considered, and it was found that a bivariate model with marginal RE and RV models for the two degradation processes is slightly the best fitting model according to information criteria.The selection of the best fitting bivariate model can be explained by the behaviors of the degradation paths of each PC.The terminal 1 has a large variation of the degradation rate and the within degradation increments, such characteristics are described by the gamma process with RE.While the terminal 2, has only a large variation of the within degradation increments, which is a characteristic of the gamma process with RV.Under the selected model, it is expected that the estimations of the reliability and the remaining useful life of the product under study to be more accurate, as noted in the comparison of reliability functions in Figure 6.The comparison of the reliability functions in Figure 7 shows that terminal 2 should decide the reliability of the product.The proposed models can be extended to consider more than two marginal degradation processes.Future research can be conducted by considering different distributions for the random effects scale parameter of the gamma process, and by letting the shape parameter to be a random effects parameter.
neous RD-RV model.This denotes that there may be some differences in the reliability assessment if homogeneous or heterogeneous random effects are considered in the modeling.Although not so obvious, the same differences can be noted in the reliability functions in Figure 6b for the RV, RV (with RD), and RV (with RE) models.The device is considered to have failed if any of the two terminal's cracks exceeds the critical degradation level ω k =0.663.As can be noted from Figure 7 the reliability of terminal 1 is higher than the reliability of terminal 2. Which can be explained by considering that the degradation mean for terminal 1 under the RE model in B 5 is 0.3715, and the degradation mean for terminal 2 under the RV model in B 5 is 0.3855.As can be noted the degradation mean for terminal 2 is greater than the degradation mean for terminal 2, which means that the first-passage times of terminal 2 will occur earlier than the ones from terminal 1.In addition, the variance of terminal 1 under RE in B 5 results in 0.007379, and the

Fig. 1 .
Fig. 1.Pseudo-random paths of the gamma degradation process under the RE model

Fig. 7 .
Fig. 7. Comparison of reliability functions for terminal 1 and 2 with RE andRV models.

Table 4 .
Estimation of information criteria and dependence structure for bivariate models

Table 3 .
Parameters estimates for the considered bivariate models