Estimation of the lifespan distribution of gold nanoparticles stabilized with lipoic acid by accelerated degradation tests and wiener process

Based on the bottom-up approach, particularly the colloidal method; gold nanoparticles with an average size of 2.5 nm were synthesized and stabilized using gold chloride as a precursor, sodium borohydride as a reducer, and lipoic acid as a stabilizer, following the methodology reported by Cornejo-Monroy et al [19].

Gold colloid stability implies that solid nanoparticles do not settle or aggregate at a significant speed [20]. When nanoparticles lose their stability by aggregation, particle size increases and creates agglomerates, losing their interesting properties.

One way to measure how these properties are affected according to their size is through their characterization by UV–vis spectroscopy [21], which is a simple and reliable method to monitor the stability of the gold colloids. As the nanoparticles become destabilized, the original characteristic peak will decrease in intensity due to the depletion of stable nanoparticles, and often the peak will be broadened to longer wavelengths due to the formation of aggregates or agglomerates [22]. The shape and peak position of UV–vis spectra are related to the morphology and size of the nanoparticle, as well as, the dispersion/aggregation of gold colloids [23]. Gold colloids present electronic transitions of bands in the visible range between 450 nm and 550 nm. Therefore, some visible wavelengths are absorbed, emitting a characteristic color that can be characterized and related to morphological changes in the nanomaterial [23]. In figure 1, J Martínez et al [24] show different sizes of gold nanoparticles as a function of the characteristic peak of the plasmon band, it can be observed that their characteristic peak is red-shift as the GNPs diameter increases. This is because the optical properties of gold nanoparticles change when the particles aggregate and the conduction electrons near the surface of each particle are delocalized and shared among neighboring nanoparticles [25]. When this occurs, the surface plasmon resonance shifts to lower energies, causing the characteristic absorption and scattering peak shift to longer wavelengths [24].

Zoom In Zoom Out Reset image size

In addition, it is known that charge repulsion effects between particles can be affected by the NaCl concentration of the solution [18]. This occurs because charges can be removed or neutralized by protonated or unprotonated ionizable groups or by the concentration of ions in solution. There is evidence that a high NaCl concentration can effectively mask the charge character of a carboxylate particle by having too many positively charged ions associated with surface charges thus causing aggregates in the material [25].

For this study, and as an example to apply an ADT analysis, to determine the failure rate and useful life of GNPs@LA, we considered that gold nanoparticles with UV–vis spectra peak greater than 525 nm failed and NaCl was used as a factor that accelerates the degradation.

One objective of the reliability analysis is to estimate the useful life of the product through the life distribution. To obtain the life distribution of a product using degradation data, the central step is to set up a model that describes the degradation process, called an accelerated degradation model. An accelerated degradation model is the combination of an accelerated model and a degradation model based on physics and statistical models.

An accelerated model shows the relationship between life and effort to establish the connection between degradation data and product life, it is essential to establish a suitable probability model to describe the behavior of collected degradation data, also known as degradation trajectories. Two types of degradation models are commonly used, which are the general trajectory models and the stochastic models.

General trajectory models are described as simple and easy to use, but they lack the ability to capture system dynamics. In contrast, stochastic models have great potential to capture random dynamics within degradation processes. The Wiener process, the Gamma process, and the Inverse Gaussian process are three common stochastic processes that have received many applications in degradation modeling [9, 26]. However, it should be noted that both the Gamma process and the Inverse Gaussian process are only suitable for modeling monotonous degradation trajectories. In comparison, the Wiener process applies to non-monotonous degradation processes that are frequently encountered in practice as it provides a satisfactory and flexible description of the degradation data [2]. The Winner process has been widely applied to degradation data analysis for example to light-emitting diodes [27], fatigue of metals [28], aluminum reduction cells [29], and micro-electromechanical systems [30], among others.

Since components of systems deteriorate over time and fail when the degradation level reaches a certain threshold. The Degradation information can be measured in a non-destructive way, after which an appropriate degradation model is chosen to describe the process through analysis of the data. Among the degradation models, the Wiener process with positive drift is well-established method due to its mathematical properties, it is expressed as follows,

$$\begin{eqnarray}&&X\left(t\right)=\lambda {\rm{\Lambda }}\left(t\right)+\sigma B\left({\rm{\Lambda }}\left(t\right)\right),\end{eqnarray} \tag{ 1 }$$

where ${\rm{\Lambda }}\left(t\right)$ represents the transformed time scale and it is a monotonic continuous function that explains the non-linearization of the data, typical examples are: ${\rm{\Lambda }}\left(t\right)={{\rm{e}}}^{{bt}},$ ${\rm{\Lambda }}\left(t\right)={blnt},{\rm{\Lambda }}\left(t\right)={t}^{b},$ where $0\lt b\lt 1$ is a parameter to be estimated [31]. The parameters $\lambda$ and $\sigma$ stand for the drift and diffusion parameters, respectively. $B\left(t\right)$ corresponds to the standard Brownian motion which satisfies the following properties:

• i.  
  $B\left(0\right)=0.$
• ii.  
  $B\left(t\right)$ has normal distribution with mean 0 and variance $t.$
• iii.  
  ${\rm{B}}\left({\rm{t}}\right)$ has independent increments, for all ${{\rm{t}}}_{1}\lt {{\rm{t}}}_{2}\,\lt \ldots \,\lt \,{{\rm{t}}}_{{\rm{n}}}$

  $$\begin{eqnarray*}&&\begin{array}{l}B\left({t}_{n}\right)-B\left({t}_{n-1}\right),\,B\left({t}_{n-1}\right)-B\left({t}_{n-2}\right),\,\ldots ,\,B\left({t}_{2}\right)-\,B\left({t}_{1}\right),B\left({t}_{1}\right)\end{array}\end{eqnarray*}$$

  are independent.
• iv.  
  $\left\{B\left(t\right),t\geqslant 0\right\}$ has stationary increments. That is, the distribution of $B\left(t+s\right)-B\left(t\right)$ does not depend on $t.$

From the above properties, it can be deduced that random vector $\left(B\left({t}_{1}\right),B\left({t}_{2}\right),\ldots ,\,B\left({t}_{n}\right)\right){\sim }{\mathscr{N}}\left(0,{\rm{\Sigma }}\right),$ where the covariance matrix ${\rm{\Sigma }}={\left({\rm{\min }}\left\{{t}_{i},{t}_{j}\right\}\right)}_{i,j},i,j=1\,,\ldots ,\,n.$ The previous properties on $B\left(t\right)$ entail that the Wiener degradation process $X\left(t\right)$ has the same properties of the Brownian motion $B\left(t\right)$ but the first property. Additional is straightforward to see that $X\left(t\right){\sim }{\mathscr{N}}\left(\lambda {\rm{\Lambda }}\left(t\right),{\sigma }^{2}{\rm{\Lambda }}\left(t\right)\right).$

Due to imperfect instruments, random environments and among other factors, measurements errors are inescapably introduced. Thus, error of measurements $\epsilon ,$ with $\epsilon {\sim }{\mathscr{N}}\left(0,{\sigma }_{\epsilon }^{2}\right)$ is introduced, leading to an observe degradation process as follows

$$\begin{eqnarray}&&Y\left(t\right)=X\left(t\right)+{\epsilon }.\end{eqnarray} \tag{ 2 }$$

Since stress factors such as voltage, humidity, temperature, vibration, etc affect the performance of the degradation process, then an acceleration model can be used to integrate the covariate into the Wiener process. The most common way to incorporate the acceleration model into the Wiener process is to consider some model parameters as a covariate function which is typically called a link function $h\left(\cdot \right).$ The choice of the form of this function will depend on the way the acceleration factor influences the model parameters. Some accelerated models are the Arrhenius model, Inverse power model, Eyring model, and linear and quadratic model, which are summarized in the following table 1.

Table 1. Acceleration models and their link function.

Acceleration model	Link function
Arrhenius Model	$h\left({S}_{k}\right)={\rm{\exp }}\left(-\frac{\beta }{{S}_{k}}\right)$
Inverse Power	$h\left({S}_{k}\right)={\rm{\exp }}\left(-\beta {\rm{ln}}\left({S}_{k}\right)\right)$
Linear Model	$h\left({S}_{k}\right)={\rm{\exp }}\left(-\beta \left({S}_{k}\right)\right)$
Quadratic Model	$h\left({S}_{k}\right)={\rm{\exp }}\left(-{\beta }_{1}\left({S}_{k}\right)+{\beta }_{2}{{(S}_{k})}^{2}\right)$

Therefore, the acceleration model of the drift parameter and diffusion parameter depends on the stress factor employing the link function as follows

$$\begin{eqnarray}&&{\lambda }_{k}=\eta h\left({S}_{k}\right),\,{\sigma }_{k}^{2}=\kappa h\left({S}_{k}\right),\end{eqnarray} \tag{ 3 }$$

where ${S}_{k}$ denotes the stress level and $\eta$ represents a variability parameter and $k$ is a constant factor associated with the diffusion. From now on $h\left({S}_{k}\right)$ will be represented as ${h}_{k}.$

It is common to find differences between the degradation trajectories from unit to unit of the population. This type of difference is the result of non-observable random effects. In order to express this, some model parameter will be specifics for each unity, obtaining a process with certain parametric distribution [28]. To model this, the drift parameter ${\lambda }_{k}$ will be specific for each unit and follows a normal distribution, on the other hand the parameter ${\sigma }_{k}$ will be taken as a constant, Peng and Tseng [32]. Si et al [33, 34] and Tsai et al [35]. Thus, it is assumed that the variability parameter $\eta$ is a random variable with normal distribution ${\mathscr{N}}\left({\mu }_{\eta },{\sigma }_{\eta }^{2}\right).$ Putting all this together into (2), it is obtained

$$\begin{eqnarray}&&{Y}_{k}\left(t\right)=Y\left(t,| ,{S}_{k}\right)=\eta {h}_{k}{\rm{\Lambda }}\left(t\right)+\sqrt{\kappa {h}_{k}}B\left({\rm{\Lambda }}\left(t\right)\right)+\epsilon ,\end{eqnarray} \tag{ 4 }$$

where (4) models the ADT with random effects, error measurements, and covariates. Some studies that configure more than one variant in the Wiener process, for instance, Li Sun et al [36] describe a methodology for the model and estimation of parameters through a constant stress ADT (CSADT) applying the non-linear Wiener process, with covariates, random effects, and measurement errors. A CSADT is a test plan consisting of three to four levels of tests with different proportions of units in each one, where mainly at the low effort level more samples run than at a high level and this type of plan can provide accurate estimates for an ADT. Following the notation in [36], the increasing applied stress level is ${{\rm{S}}}_{1}\,,\ldots {{\rm{S}}}_{{\rm{k}}}\,,\ldots {{\rm{S}}}_{{\rm{K}}}\,,$ where $K$ denotes the maximum stress level. Also, there are ${N}_{k}$ units tested under each constant stress ${S}_{k}$ and each unit is measured ${M}_{{ki}}$ times at the $k$ stress level with $i=1,2,\,\ldots ,\,{N}_{k}.$ The transformed time will be expressed as ${\omega }_{{kij}}={\rm{\Lambda }}\left({t}_{{kij}}\right)$ with $j=\mathrm{1,2},\,\ldots ,\,{M}_{{ki}}.$

Therefore, the degradation observed under the Wiener process with its four variants is shown as follows

$$\begin{eqnarray}&&{\,y}_{{kij}}\left({{\rm{\omega }}}_{{kij}},| ,{S}_{k}\right)={\eta }_{i}{{h}_{k}\omega }_{{kij}}+\sqrt{\kappa {h}_{k}}B\left({\omega }_{{kij}}\right)+{\varepsilon }_{{kij}},\end{eqnarray} \tag{ 5 }$$

with ${\eta }_{i}\sim N\left({\mu }_{\eta },{\sigma }_{\eta }^{2}\right)\,,$ ${\varepsilon }_{{kij}}\sim N\left(0,{\sigma }_{\varepsilon }^{2}\right)$ and ${h}_{k}$ the link function of the accelerated model.

Model (5) comprises the following set of unknown parameters,

$$\begin{eqnarray}&&{\boldsymbol{\Theta }}=\left\{{\mu }_{\eta },{\sigma }_{\eta }^{2},\kappa ,\beta ,{\sigma }_{\varepsilon }^{2}\right\}\end{eqnarray} \tag{ 6 }$$

that will be estimated in the next section.

In the previous section a model for CSADT was formulated in (5), to estimate the unknown parameters set $\tilde{{\boldsymbol{\Theta }}}$ let's consider the following vectors ${{\boldsymbol{T}}}_{{ki}}={\left({t}_{{ki}1},\,\ldots ,\,{t}_{{ki}{M}_{{ki}}}\right)}^{T},$ ${t}_{{kij}}={\omega }_{{kij}},$ ${{\boldsymbol{y}}}_{{ki}}={\left({y}_{{ki}1},\,\ldots ,\,{y}_{{ki}{M}_{{ki}}}\right)}^{T},$ ${{\boldsymbol{y}}}_{k}={\left({{\boldsymbol{y}}}_{k1},\,\ldots ,\,{{\boldsymbol{y}}}_{k{N}_{k}}\right)}^{T}$ and ${\boldsymbol{Y}}=({{\boldsymbol{y}}}_{1},\,\ldots ,\,{{\boldsymbol{y}}}_{K\,}$) for $k=\mathrm{1,2},\,\ldots \,K,$ $i=\mathrm{1,2},\,\ldots ,\,{N}_{k}.$ The properties ii and iii of the Wiener processes implies that the Brownian movement ${{\boldsymbol{y}}}_{{ki}}$ follows a multivariate normal distribution with a mean,

$$\begin{eqnarray}&&{\mu }_{{ki}}={\mu }_{\eta }{h}_{k}{{\boldsymbol{T}}}_{{ki}}\end{eqnarray} \tag{ 7 }$$

and covariance

$$\begin{eqnarray*}&&{{\boldsymbol{\Omega }}}_{{ki}}={\sigma }_{\eta }^{2}{\tilde{{\boldsymbol{\Omega }}}}_{{ki},}\,{\tilde{{\boldsymbol{\Omega }}}}_{{\boldsymbol{ki}}}={\tilde{{\boldsymbol{{\rm H}}}}}_{{ki}}+{h}_{k}^{2}{{\boldsymbol{T}}}_{{ki}}{{\boldsymbol{T}}}_{{ki}}^{{\bf{T}}},\end{eqnarray*}$$

$$\begin{eqnarray}&&{\tilde{{\boldsymbol{{\rm H}}}}}_{{ki}}=\tilde{\kappa }{h}_{k}\left[\begin{array}{cc}\begin{array}{cc}\begin{array}{c}{\omega }_{{ki}1}\\ {\omega }_{{ki}1}\end{array} & \begin{array}{c}{\omega }_{{ki}1}\\ {\omega }_{{ki}2}\end{array}\end{array} & \begin{array}{cc}\begin{array}{c}\ldots \\ \ldots \end{array} & \begin{array}{c}{\omega }_{{ki}1}\\ {\omega }_{{ki}2}\end{array}\end{array}\\ \begin{array}{c}\begin{array}{cc}\vdots & \vdots \end{array}\\ \begin{array}{cc}{\omega }_{{ki}1} & {\omega }_{{ki}2}\end{array}\end{array} & \begin{array}{cc}\begin{array}{c}\ddots \\ \ldots \end{array} & \begin{array}{c}\vdots \\ {\omega }_{{ki}{M}_{{ki}}}\end{array}\end{array}\end{array}\right]+{\tilde{\sigma }}_{\varepsilon }^{2}{{\boldsymbol{I}}}_{{ki}}\end{eqnarray} \tag{ 8 }$$

where ${{\boldsymbol{I}}}_{{ki}}$ is the identity matrix of order ${M}_{{ki}},$ where $\tilde{k}=\frac{k}{{\sigma }_{\eta }^{2}}$ and ${\tilde{\sigma }}_{\epsilon }^{2}=\frac{{\sigma }_{\epsilon }^{2}}{{\sigma }_{\eta }^{2}}$ are a reparameterization to facilitate the inference. The log-likelihood function of unknown parameters $\mathop{{\boldsymbol{\Theta }}}\limits^{\sim }=\left\{{\mu }_{\eta },{\sigma }_{\eta }^{2},\tilde{\kappa },\beta ,{\tilde{\sigma }}_{\varepsilon }^{2}\right\}$ is expressed as

$$\begin{eqnarray}&&l\left(\mathop{{\boldsymbol{\Theta }}}\limits^{\sim },| ,{\boldsymbol{Y}}\right)\,=\,-\frac{1}{2}{\rm{ln}}\left(2\pi \right)\displaystyle \sum _{k=1}^{K}\displaystyle \sum _{i=1}^{{N}_{k}}{M}_{{ki}}\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray*}&&-\frac{1}{2}{\rm{ln}}\left({\sigma }_{\eta }^{2}\right)\displaystyle \sum _{k=1}^{K}\displaystyle \sum _{i=1}^{{N}_{k}}{M}_{{ki}}-\frac{1}{2}\displaystyle \sum _{k=1}^{K}\displaystyle \sum _{i=1}^{{N}_{k}}{ln}\left|{\tilde{{\boldsymbol{\Omega }}}}_{{ki}}\right|\end{eqnarray*}$$

$$\begin{eqnarray*}&&-\frac{1}{2{\sigma }_{\eta }^{2}}\displaystyle \sum _{k=1}^{K}\displaystyle \sum _{i=1}^{{N}_{k}}{\left({{\boldsymbol{y}}}_{{ki}}-{\mu }_{\eta }{h}_{k}{{\boldsymbol{T}}}_{{ki}}\right)}^{T}{\tilde{{\boldsymbol{\Omega }}}}_{{ki}}^{-{\bf{1}}}\left({{\boldsymbol{y}}}_{{ki}}-{\mu }_{\eta }{h}_{k}{{\boldsymbol{T}}}_{{ki}}\right).\end{eqnarray*}$$

To estimate the parameters ${\mu }_{\eta }$ and ${\sigma }_{\eta }^{2},$ the values of $\left(\tilde{\kappa },\beta ,{\tilde{\sigma }}_{\varepsilon }^{2}\right)$ are fixed to specific values, then the derivatives $\frac{\partial l\left(\tilde{{\boldsymbol{\Theta }}},| ,{\boldsymbol{Y}}\right)\,}{\partial {\mu }_{\eta }}\,{and}\,\frac{\partial l\left(\tilde{{\boldsymbol{\Theta }}},| ,{\boldsymbol{Y}}\right)\,}{\partial {\sigma }_{\eta }^{2}}$ of (9) concerning the parameters ${\mu }_{\eta }$ and ${\sigma }_{\eta }^{2}$ are computed and equated to zero, obtaining the following close -form expression to the maximum likelihood estimator MLE ${\hat{\mu }}_{\eta }$ and ${\hat{\sigma }}_{\eta }^{2}$ of ${\mu }_{\eta }$ and ${\sigma }_{\eta }^{2}$

$$\begin{eqnarray}&&{\hat{\mu }}_{\eta }=\frac{{\sum }_{k=1}^{K}{\sum }_{i=1}^{{N}_{k}}{h}_{k}{{\boldsymbol{T}}}_{{ki}}^{{\boldsymbol{T}}}{\tilde{{\boldsymbol{\Omega }}}}_{{ki}}^{-{\bf{1}}}{{\boldsymbol{y}}}_{{ki}}}{{\sum }_{k=1}^{K}{\sum }_{i=1}^{{N}_{k}}{h}_{k}^{2}{{{\boldsymbol{T}}}_{{ki}}^{{\boldsymbol{T}}}{\tilde{{\boldsymbol{\Omega }}}}_{{ki}}^{-{\bf{1}}}{\boldsymbol{T}}}_{{ki}}}\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}&&{\hat{\sigma }}_{\eta }^{2}=\frac{1}{{\sum }_{k=1}^{K}{\sum }_{i=1}^{{N}_{k}}{M}_{{ki}}}\,\displaystyle \sum _{k=1}^{K}\displaystyle \sum _{i=1}^{{N}_{k}}{\left({{\boldsymbol{y}}}_{{ki}}-{\hat{\mu }}_{\eta }{h}_{k}{{\boldsymbol{T}}}_{{ki}}\right)}^{T}{\tilde{{\boldsymbol{\Omega }}}}_{{ki}}^{-{\bf{1}}}\left({{\boldsymbol{y}}}_{{ki}}-{\hat{\mu }}_{\eta }{h}_{k}{{\boldsymbol{T}}}_{{ki}}\right).\end{eqnarray} \tag{ 11 }$$

Substituting (10) and (11) in (9) and simplifying the log-likelihood function, it is obtained

$$\begin{eqnarray*}&&l\left(\tilde{{\boldsymbol{\Theta }}},| ,{\boldsymbol{Y}}\right)=-\frac{1}{2}\left({\rm{ln}}\left(2\pi \right)+1\right)\displaystyle \sum _{k=1}^{K}\displaystyle \sum _{i=1}^{{N}_{k}}{M}_{{ki}}\end{eqnarray*}$$

$$\begin{eqnarray}&&-\frac{1}{2}{ln}\left({\hat{\sigma }}_{\eta }^{2}\right)\displaystyle \sum _{k=1}^{K}\displaystyle \sum _{i=1}^{{N}_{k}}{M}_{{ki}}-\frac{1}{2}\displaystyle \sum _{k=1}^{K}\displaystyle \sum _{i=1}^{{N}_{k}}{ln}\left|{\tilde{{\boldsymbol{{\rm H}}}}}_{{ki}}+{h}_{k}^{2}{{\boldsymbol{T}}}_{{ki}}{{\boldsymbol{T}}}_{{ki}}^{{\bf{T}}}\right|.\,\end{eqnarray} \tag{ 12 }$$

Note that the matrix ${\tilde{{\boldsymbol{{\rm H}}}}}_{{ki}}$ depends on the parameters $\tilde{\kappa },\beta$ and ${\tilde{\sigma }}_{\varepsilon }^{2}$ and ${h}_{k}^{2}$ depends on $\beta .$ Therefore, the maximum likelihood estimates of $\tilde{\kappa },\beta ,{\tilde{\sigma }}_{\varepsilon }^{2}$ can be obtained by maximizing the log-likelihood function (12) by employing the L-BFGS-B quasi-Newton optimization method that can be found in the R-project packages [37]. The value of $\hat{\kappa }$ can be obtained by the following equations:

$$\begin{eqnarray}&&\hat{\kappa }=\,\hat{\tilde{\kappa }}\cdot {\hat{\sigma }}_{\eta }^{2},\end{eqnarray} \tag{ 13 }$$

$$\begin{eqnarray}&&{\hat{\sigma }}_{\varepsilon }^{2}={\hat{\tilde{\sigma }}}_{\varepsilon }^{2}\cdot {\hat{\sigma }}_{\eta }^{2},\end{eqnarray} \tag{ 14 }$$

One objective of the reliability analysis is to estimate the useful life of the product through the life distribution. To obtain the life distribution of a product using degradation data Li et al [31] incorporate the measurement errors into the deduction of the expression for the CDF and PDF of the failure time

$$\begin{eqnarray}&&T={\inf }\left\{t:\,Y\left(t\right)\geqslant w\right\},\end{eqnarray} \tag{ 15 }$$

where $T$ corresponds to the first time that the degradation process $Y$ hits a failure threshold $w.$ They deduced the life PDF expression for each stress level ${S}_{k},$ which can be found in equation (12) at [31] Note that in this formulation they used two different time scales ${\rm{\Lambda }}\left(t\right)\,{\rm{and}}\tau \left(t\right),$ when $\tau \left(t\right)={\rm{\Lambda }}(t),$ this PDF is reduced to,

$$\begin{eqnarray}&&{f}_{{S}_{k}}\left(t\right)=\frac{{\sigma }_{\varepsilon }^{2}{h}_{k}{\mu }_{\eta }+\left(\kappa {h}_{k}+{{h}_{k}}^{2}{\sigma }_{\eta }^{2}\right)w}{\sqrt{2\pi {\left({h}_{k}\omega \left(\kappa +{h}_{k}{\sigma }_{\eta }^{2}\omega \right)+{\sigma }_{\varepsilon }^{2}\right)}^{3}}}\cdot \frac{d\omega }{{dt}}\cdot {\rm{\exp }}\left\{\frac{-{\left(w-{h}_{k}{\mu }_{\eta }\omega \right)}^{2}}{2\left({h}_{k}\omega \left(\kappa +{h}_{k}{\sigma }_{\eta }^{2}\omega \right)+{\sigma }_{\varepsilon }^{2}\right)}\right\}\end{eqnarray} \tag{ 16 }$$

which is the one used in this work. Note ${f}_{{S}_{k}}\left(t\right)$ depends on the parameters set $\tilde{{\boldsymbol{\Theta }}}.$

To obtain degradation data of GNPs@LA under a CSADT, several samples were synthesized following the methodology reported by Cornejo et al [19] and using NaCl as an acceleration factor. The stress test levels were based on an exploratory study, leaving three different levels of effort for three different populations. UV–vis absorption spectra of colloids were carried out every third-day generating degradation signals for low, medium, and high levels. In this study, the material degradation was determined considering the absorbance between 450 and 550 nm, and the maximum characteristic peak of gold. Gold colloids with a characteristic peak greater than 525 were considered as a failure.

Since this work proposes to estimate the life distribution of GNPs@LA applying a CSADT. Accordingly, the following can be defined:

• The percentage levels of NaCl are indexed as $k=1,\,\ldots ,\,K$
• The population of samples is indexed as $i=1,\,\ldots ,\,{N}_{k}$
• Measurement times are indexed as $j=1,\,\ldots ,\,{M}_{{ki}}$ Once the degradation signals have been defined as degradation data, we proceed to obtain the configuration of the CSADT that describes the degradation trajectories. Thus, this work proposes the non-linear Wiener process with drift parameter, random effects, measurement errors, and different link functions in the covariability (table 1, section 2.3) is used. More specifically the model assumes:
• A drift parameter with transformed times ${{\rm{\omega }}}_{{\rm{kij}}}={\rm{\Lambda }}\left({{\rm{t}}}_{{\rm{kij}}}\right)={\left({{\rm{t}}}_{{\rm{kij}}}\right)}^{{\rm{b}}}$ that explain the non-linearity of the data.
• Random effects in the drift parameter. Different units have different drift parameters while all diffusion parameters have the same value under a certain effort. So, the parameter ${\rm{\eta }}$ in the drift term at (5) is assumed to has a normal distribution $\eta \sim {\rm{N}}\left({\mu }_{\eta },{\sigma }_{\eta }^{2}\right).$
• Measurement errors with ${\varepsilon }_{{\rm{kij}}}\sim {\rm{N}}\left(0,{\sigma }_{\varepsilon }^{2}\right).$
• Covariability, which is expressed via the link function ${{\rm{h}}}_{{\rm{k}}},$ where the accelerated model is immersed. It is highlighted three accelerated models that were evaluated in this work corresponding, inverse power, linear model, and the quadratic model, see (table 1, section 2.3).

Now, proposing the CSADT and the Wiener process with these characteristics, the degradation process can be established as formula (5). Thus, the following parameter set $\tilde{{\rm{\Theta }}}=\left\{{\mu }_{\eta },{\sigma }_{\eta }^{2},\tilde{{\kappa }},\beta ,{\tilde{{\sigma }}}_{\varepsilon }^{2}\right\}$ will be estimated to get the life distribution. We remark that $\beta =({\beta }_{1},{\beta }_{2})$ in the quadratic model case and $\beta ={\beta }_{1}$ in the remaining cases.

All the above were programmed in the statistical software R, where the MLE was applied to obtain the parameters set $\tilde{{\boldsymbol{\Theta }}}.$ Figure 2 shows the proposed methodology.

Zoom In Zoom Out Reset image size

As can be seen in figure 2, different degradation models have been proposed. To select the best model will be used the AIC criterion [38] and for the evaluation and validation, the estimated parameters of the model will be employed and the Bootstrapping distribution [39] will be calculated with the construction of confidence intervals.