Q-Learning-Incorporated Robust Relevance Vector Machine for Remaining Useful Life Prediction

Abstract

Accurate and reliable remaining useful life (RUL) prediction is crucial for improving equipment reliability and safety, realizing predictive maintenance. The relevance vector machine (RVM) method is commonly utilized for RUL prediction, profiting from its sparse property under a Bayesian framework. However, the RVM faces the issue of poor robustness, which is mainly manifested as poor prediction accuracy and difficulty in fitting when the predicted data fluctuate greatly. This is due to weights and random errors following Gaussian distributions, which are highly sensitive to outliers. Also, the traditional model training process heavily relies on an additional feature extraction process, which suffers from the problem of effective data loss as well as the risk of overfitting. Thus, a robust regression framework against outliers is developed by incorporating t-distribution into the RVM. And a Q-learning (QL) algorithm is embedded into the constructed robust RVM model to replace the feature extraction process. In addition, this paper firstly predicts the degradation trend of RUL to enhance the accuracy and interpretability of RUL prediction. Finally, a comparative experiment on the performance degradation of capacitors in the traction system is designed, and the root mean square errors for the QL-RRVM, QL-RVM, RRVM, and RVM models are obtained as 0.751, 8.599, 38.316, and 41.892, respectively. The experimental results confirm the superiority of the proposed method.

1. Introduction

Remaining Useful Life (RUL) prediction plays a crucial role in industrial systems for achieving efficient equipment maintenance and management [1]. Accurately predicting the RUL of equipment or components enables engineers to effectively avoid unexpected equipment failures in advance, and to strategically implement preventive maintenance, minimizing downtime and reducing production costs [2,3].

RUL prediction methods can be broadly classified into two main categories: model-based methods and data-driven methods. The model-based methods utilize physics models [4,5,6], mechanistic models [7,8] or stochastic process [9] models for RUL prediction. A physical equation was established and converted into a multivariate linear regression model to improve generalization capability of cold load prediction in Reference [10]. The storage life of a certain electronic component was predicted by analyzing its failure process based on a reaction-rate model in Reference [11]. A nonlinear dynamic system was established by the tensor product-based model transformation technique in Reference [12]. However, heavy reliance on prior knowledge and expert experience poses a significant obstacle to model-based approaches. It is extremely difficult to accurately obtain the physical models of failure mechanism with the increasing complexity of engineering equipment [13].

Data-driven approaches are booming in RUL prediction, as they do not require physical knowledge to infer the degradation process. Data-driven RUL prediction methods, such as neural networks [14,15,16] and the relevance vector machine (RVM) [17], aim to analyze historical and monitoring data to uncover the underlying relationships among data, enabling the prediction of degradation trends and subsequent analysis of RUL [18]. Neural network methods play a dominant role in the data-driven approach [19]. For instance, an improved Convolutional Neural Network (CNN) and Long Short-Term Memory (LSTM) were combined to enhance the prediction inaccuracy of lithium battery-life prediction in Reference [15]. A deep Bayesian learning method was proposed for fast and uncertainty-aware RUL prediction in Reference [16]. However, the prediction accuracy derived from neural networks is limited, due to the number of available data. In contrast, the RVM method, a small-sample prediction approach, demonstrates excellent performance when faced with insufficient data. The RVM method is under a Bayesian framework, where the weights for each input are controlled by a set of parameters [17,20]. RVM is widely applied in regression tasks, due to its sparsity and excellent prediction performance [21,22,23]. For example, Ref. [24] utilized an RVM method to describe a degradation process from initial fault to failure. Additionally, a degradation path-based RUL prediction framework is established using a dynamic multivariate RVM model in Reference [25].

While RVMs demonstrate considerable promise in achieving high accuracy and sparsity, their predictive prowess can falter in the presence of substantial data variability and outliers [26]. This limitation stems from the fact that the weights and random noise in RVMs adhere to Gaussian distributions, which are inherently susceptible to the influence of outliers and atypical data points. Furthermore, the application of RVM models for predictive or classification tasks typically necessitates an auxiliary feature-extraction procedure. For example, in Reference [27], features are derived from the battery-charging voltage curve through the application of random forest and Gaussian regression techniques. Similarly, in Reference [28], Principal Component Analysis (PCA) is employed to both reduce data dimensionality and extract pertinent features. However, these supplementary feature extraction techniques are not only contingent upon human expertise but also carry the risk of inadvertently leading to the loss or degradation of crucial information.

Motivated by these discussions, the RVM is extended into a robust one by integrating t-distribution into the RVM framework, for striving against outliers. And a Q-learning algorithm from the realm of reinforcement learning is introduced into the newly formulated RRVM model. In Q-learning, the state is the candidate feature, which transfers to the next state via random action. After all actions are performed, the optimal feature is selected corresponding to each time point. This integration culminates in the development of a QL-RRVM model that is capable of directly processing raw data and delivering precise predictions without the need for intermediate feature extraction steps. Specifically, the contributions of this study are summarized as follows.

• It is inferred that the weights and random errors of the proposed model obey t-distributions, with the assumption that hyperparameters are gamma distributions, which are more robust for outliers than classical RVM.
• The QL algorithm is introduced into the RRVM model for feature extraction to improve the model prediction accuracy.
• In the field of RUL prediction, the model transitions from static to dynamic, employing a forward-rolling prediction within a time-series context. It initially projects degradation trends essential for RUL estimation and subsequently refines these projections to forecast RUL, thereby improving the interpretability of the predictions.

The remainder of this paper is organized as follows. The RRVM model is given in Section 2. RUL prediction is performed in Section 3. QL and RRVM models are combined in Section 4. In Section 5, the effectiveness of the proposed approach is validified by a case study. Finally, a conclusion is summarized in Section 6.

2. Robust RVM Modelling

2.1. Modelling Process

For a nonlinear regression task, it is common to use the kernel function to represent the nonlinear relationship between the input vectors ${\left\{{\mathit{x}}_i\right\}}_{i=1}^N$ and corresponding targets ${\left\{t_i\right\}}_{i=1}^N$. Typically, predictions are based on a function y(x) defined over the input space, with “learning” being the process of inferring this function (potentially its parameters). A flexible and popular candidate set for y(x) is as follows:

$$y(\mathit{x})={\displaystyle \sum _{i=1}^N{w}_i{\phi }_i(\mathit{x})}=\mathit{\Phi }\mathit{W}$$

(1)

where $\mathit{\Phi }={[\phi ({\mathit{x}}_1),\phi ({\mathit{x}}_2),\cdots ,\phi ({\mathit{x}}_N)]}^{\mathrm{T}}\in {ℝ}^{N\times N+1}$ is a design matrix, in which $\phi ({\mathit{x}}_i)=[1,k({\mathit{x}}_i,{\mathit{x}}_1),k({\mathit{x}}_i,{\mathit{x}}_2),\cdots ,k({\mathit{x}}_i,{\mathit{x}}_N)]\in {ℝ}^{1\times (N+1)}$ and $k({\mathit{x}}_i,{\mathit{x}}_k)$ is a Gaussian kernel function between the vector ${\mathit{x}}_i$ and ${\mathit{x}}_k(k=1,\cdots ,N)$ in this study; $\mathit{W}={[w_1,w_2,\cdots ,w_{N+1}]}^{\mathrm{T}}\in {ℝ}^{(N+1)\times 1}$ denotes the weight matrix.

In the context of the RVM model, a time-series dataset ${\left\{{\mathit{x}}_i\right\}}_{i=1}\in {ℝ}^{1\times (N+L)}$ is provided, with ${\mathit{x}}_i=({\mathit{x}}_i,{\mathit{x}}_{i+1},\cdots ,{\mathit{x}}_{i+L-1})$$\in {ℝ}^{1\times L}$, where N denotes the number of target groups, and L represents the stride length of the sliding window. The target value can be represented as the sum of additive noise and the actual value

$$t_i=y({\mathit{x}}_i)+\epsilon$$

(2)

where $\epsilon$ is assumed to follow a Gaussian distribution with zero mean and variance ${\sigma }^2$. The precision $\beta$ is introduced as ${\beta }^{-1}={\sigma }^2$. Then, the probability density function (PDF) of the predicted value $t_i$ with respect to $w_i$ and ${\beta }^{-1}$ can be expressed as

$$p(t_i{\left|{\mathit{x}}_i,w_i,\beta \right.}^{-1})\sim \mathcal{N}(y({\mathit{x}}_i),{\beta }^{-1})={(2\pi {\beta }^{-1})}^{-{\displaystyle \frac{1}{2}}}\exp \left\{-{\displaystyle \frac{1}{2{\beta }^{-1}}}{(t_i-w_i\phi ({\mathit{x}}_i))}^2\right\}$$

(3)

N groups of data are considered in the following format:

$$p(\mathit{T}{\left|\mathit{X},\mathit{W},\beta \right.}^{-1})\sim {\displaystyle \prod _{i=1}^N\mathcal{N}(\mathit{\Phi }\mathit{W},{\beta }^{-1})}={(2\pi {\beta }^{-1})}^{-{\displaystyle \frac{N}{2}}}\exp \left\{-{\displaystyle \frac{1}{2{\beta }^{-1}}}{‖\mathit{T}-\mathit{\Phi }(\mathit{X})\mathit{W}‖}^2\right\}$$

(4)

where $\mathit{T}$, $\mathit{X}$, $\mathit{W}$ correspond to the matrix forms of $t_i$, $x_i$, $w_i$, respectively.

The parameter $\mathit{W}$ is defined as obeying a Gaussian distribution with zero mean and different variances to guarantee sparsity prediction in RVM. Its PDF is given as follows.

$$p(\mathit{W}|\mathit{\alpha })={\displaystyle \prod _{i=1}^{N+1}\mathcal{N}(w_i|0,{\mathit{\alpha }}_i^{-1})}$$

(5)

where $\mathit{\alpha }={({\alpha }_1,\dots ,{\alpha }_{N+1})}^{\mathrm{T}}\in {ℝ}^{N+1}$ and ${\alpha }_i,\hspace{0.33em}i\in [1,N+1]$ denote precision corresponding with $\mathit{W}$.

The training objective is to estimate the parameters $\mathit{W},\beta \hspace{0.33em}\mathrm{and}\hspace{0.33em}\mathit{\alpha }$ by maximizing the conditional PDF $p(\mathit{W},\mathit{\alpha },\beta |t_i)$ with the target $t_i$

$$p(\mathit{W},\mathit{\alpha },\beta |t_i)=p(\mathit{W}|t_i,\mathit{\alpha },\beta )p(\mathit{\alpha },\beta |t_i)$$

(6)

The posterior probability distribution of the parameter $\mathit{W}$ remains Gaussian. Then, there is

$$p(\mathit{W}|t_i,\mathit{\alpha },\beta )=\mathcal{N}(\mathit{W}|\mathit{\Lambda },\mathit{\Sigma })$$

(7)

where $\mathit{\Lambda }=\beta \mathit{\Sigma }{\mathit{\Phi }}^{\mathrm{T}}\mathit{T}$ is the mean; $\hspace{0.33em}\mathit{\Sigma }={(\mathit{A}+\beta {\mathit{\Phi }}^{\mathrm{T}}\mathit{\Phi })}^{-1}$ denotes the variance, in which $\mathit{A}=diag({\alpha }_1,{\alpha }_2,\dots ,{\alpha }_{N+1})\in {ℝ}^{N+1}$.

The poor robustness of the traditional RVM model stems from the assumption that both the weights and random errors follow Gaussian distributions, which are sensitive to outliers and anomalous data points. Thus, the parameters $\mathit{\alpha }$ and $\beta$ are assumed to follow a gamma distribution in this study. Then, the predicted values follow the t-distribution, which is more robust than the Gaussian distribution. A gamma distribution is introduced for the parameters $\mathit{\alpha }$ and $\beta$. And their corresponding PDFs are given below.

$$p(\mathit{\alpha })={\displaystyle \prod _{i=0}^l\Gamma ({\alpha }_i|a,b)},p(\beta )=\Gamma (\beta |c,d)$$

(8)

where a, c are shape parameters; b and d are scale parameters; $\Gamma (x|y,z)=\Gamma {(y)}^{-1}{z}^y{x}^{y-1}{e}^{-zx}$ in which $\Gamma (x)={\displaystyle {\int }_0^{\infty }{t}^{x-1}{e}^{-t}dt}.$

Then, the posterior PDFs of weights $\mathit{W}$ and random errors $\epsilon$ are

$$p(\mathit{W})={\displaystyle \int p(\mathit{W}|\mathit{\alpha })p(\mathit{\alpha })d\mathit{\alpha }}={\displaystyle \frac{\Gamma \left(\frac{\nu }{2}+\frac{1}{2}\right)}{\Gamma \left(\frac{\nu }{2}\right)}}{\left({\displaystyle \frac{\lambda }{\pi \nu }}\right)}^{{\displaystyle \frac{1}{2}}}{\left(1+{\displaystyle \frac{\lambda {\mathit{W}}^2}{\nu }}\right)}^{-\left({\displaystyle \frac{\nu }{2}}+{\displaystyle \frac{1}{2}}\right)}=t(\mathit{W})$$

(9)

$$p(\epsilon )={\displaystyle \int p(\epsilon |\beta )p(\beta )d\beta }={\displaystyle \frac{\mathsf{\Gamma }\left(\frac{v^{\prime }}{2}+\frac{1}{2}\right)}{\mathsf{\Gamma }\left(\frac{v^{\prime }}{2}\right)}}{\left({\displaystyle \frac{{\lambda }^{\prime }}{\pi v^{\prime }}}\right)}^{{\displaystyle \frac{1}{2}}}{\left(1+{\displaystyle \frac{{\lambda }^{\prime }{\epsilon }^2}{{\nu }^{\prime }}}\right)}^{-\left({\displaystyle \frac{v^{\prime }}{2}}+{\displaystyle \frac{1}{2}}\right)}=t(\epsilon )$$

(10)

where the parameters $\lambda =a/b$ (${\lambda }^{\prime }=c/d$) and $\nu =2a$ ($v^{\prime }=2c$) present the precision and degrees of freedom of the weights $\mathit{W}$ (error $\epsilon$) obeying the t-distribution.

The Equations (9) and (10) suggest that both the weights $\mathit{W}$ and the random errors $\epsilon$ adhere to t-distribution. It is noteworthy that the t-distribution is distinguished by the presence of extended tails, a characteristic that sets it apart from the Gaussian distribution. Thus, the established regression model (2) is robust by assuming that the regression model coefficients $\mathit{W}$ and the random errors $\epsilon$ follow a long-tailed t-distribution.

The subsequent section deduces the forecasts of the model, which is predicated on the completion of parameter estimation (the parameter estimation process is detailed in Section 3), assuming optimal estimation of parameters $\mathit{\alpha }\mathrm{and}\beta$ has been achieved. Initially, the model is trained assuming that under the condition at time i, the most fitting $\mathit{\alpha }\hspace{0.33em}\mathrm{and}\beta$ for the current data, obtained by training the model with data from time i − 1, is denoted as ${\mathit{\alpha }}_{i-1}^{best}\hspace{0.33em}\mathrm{and}{\beta }_{i-1}^{best}$. Then, we calculate the predictive distribution of data $X_{i-1}$ by combining Equations (10) and (11) as follows.

$$p(t_i|t_{i-1})={\displaystyle \int p(t_i|\mathit{W},\mathit{\alpha },\beta )p(\mathit{W},\mathit{\alpha },\beta |t_{i-1})}d\mathit{W}d\mathit{\alpha }d\beta$$

(11)

$$p(t_i|t_{i-1},{\mathit{\alpha }}_{i-1}^{best},{\beta }_{i-1}^{best})={\displaystyle \int p(t_i|\mathit{W},{\beta }_{i-1}^{best})p(\mathit{W}|t_{i-1},{\mathit{\alpha }}_{i-1}^{best},{\beta }_{i-1}^{best})}d\mathit{W}$$

(12)

Upon incorporating the PDF of the respective distribution into the aforementioned equation and performing simplifications, it is determined that the resulting form aligns with the properties of the Gaussian distribution, i.e., $p(t_i|t_{i-1},{\mathit{\alpha }}_{i-1}^{best},{\beta }_{i-1}^{best})=\mathcal{N}(t_i|y_i,{\beta }_i)$, where $y_i={\mathit{W}}^{\mathrm{T}}\mathit{\Phi }({\mathit{x}}_{i-1})$ and ${\beta }_i={\beta }_{i-1}^{best}+\mathit{\Phi }{({\mathit{x}}_{i-1})}^{\mathrm{T}}\mathit{\Sigma }\mathit{\Phi }({\mathit{x}}_{i-1})$. In this study, the mean value $y_i$ is considered as the predicted value for the data ${\mathit{x}}_{i-1}$. Similarly, a prediction for data at moment i + 1 has $p(t_{i+1}|t_i,{\mathit{\alpha }}_i^{best},{\beta }_i^{best})=\mathcal{N}(t_{i+1}|y_{i+1},{\beta }_{i+1})$, where $y_{i+1}={\mathit{W}}^{\mathrm{T}}\mathit{\Phi }({\mathit{x}}_i)$ and ${\beta }_{i+1}={\beta }_i^{best}+\mathit{\Phi }{({\mathit{x}}_i)}^{\mathrm{T}}\mathit{\Sigma }\mathit{\Phi }({\mathit{x}}_i)$.

Employing this method, the time-series data is sequentially forecasted, and the actual data at time i + L + 1 are utilized as the training target for $X_i$. This approach aims to progressively reduce the error, ultimately achieving accurate rolling forecasts.

2.2. Parameter Estimation

In this section, parameters $\mathit{\alpha }$ and $\beta$ in the regression model are estimated. The optimal parameters are sought by maximizing the posterior probability, $\max (p(\mathit{\alpha },\beta |t_i))$.

The PDF of the predicted values with respect to the parameters $\mathit{\alpha },\beta$ are calculated according to the convolution formula, as follows [29].

$$p(t_i|\mathit{\alpha },\beta )={\displaystyle {\int }_{-\infty }^{+\infty }p(t_i|\mathit{W},\beta )p(\mathit{W}|\mathit{\alpha })d\mathit{W}}={({\displaystyle \frac{\beta }{2\pi }})}^{{\displaystyle \frac{1}{2}}}{({\displaystyle \frac{1}{\beta }}+\mathit{\Phi }{\mathit{A}}^{-1}{\mathit{\Phi }}^{\mathrm{T}})}^{-{\displaystyle \frac{1}{2}}}e\mathit{x}p(-{\displaystyle \frac{1}{2}}({\displaystyle \frac{1}{\beta }}+\mathit{\Phi }{\mathit{A}}^{-1}{\mathit{\Phi }}^{\mathrm{T}})t_i^2)$$

(13)

Taking the logarithm of the above relevant equation, the following log-likelihood function is given,

$$LF=lnp(t_i|ln\mathit{\alpha },ln\beta )+lnp(ln\mathit{\alpha })+lnp(ln\beta )$$

(14)

As $p(ln\mathit{\alpha })=\mathit{\alpha }p(\mathit{\alpha }),\hspace{0.33em}\mathrm{and}\hspace{0.33em}p(t_i|ln\mathit{\alpha },ln\beta )=p(t_i|\mathit{\alpha },\beta )$, the likelihood function (14), ignoring the the terms unrelated to $\mathit{\alpha },\beta$, is converted as follows.

$$LF=-{\displaystyle \frac{1}{2}}[ln|{\displaystyle \frac{1}{\beta }}\mathit{I}+{\mathit{\Phi }}^{\mathrm{T}}{\mathit{A}}^{-1}\mathit{\Phi }|+t_i{({\displaystyle \frac{1}{\beta }}\mathit{I}+{\mathit{\Phi }}^{\mathrm{T}}{\mathit{A}}^{-1}\mathit{\Phi })}^{-1}{t}_i]+{\displaystyle \sum _{i=1}^{N+1}(a-ln{\alpha }_i-b{\alpha }_i)+cln\beta -d\beta }$$

(15)

Then, take the partial derivatives of $\mathit{\alpha },\beta$ respectively and let them be 0.

$${\displaystyle \frac{\partial LF}{\partial ln\mathit{\alpha }}}={\displaystyle \frac{1}{2}}[1-\mathit{\alpha }({\mathit{\Lambda }}^2+\mathit{\Sigma })]+a-b\mathit{\alpha }=0$$

(16)

$${\displaystyle \frac{\partial LF}{\partial ln\beta }}={\displaystyle \frac{1}{2}}[{\displaystyle \frac{N}{\beta }}-{(t_i-\mathit{\Phi }\mathit{\Lambda })}^2-tr(\mathit{\Sigma }{\mathit{\Phi }}^{\mathrm{T}}\mathit{\Phi })]+c-d\beta =0$$

(17)

An implicit expression for the parameters $\mathit{\alpha }\hspace{0.33em}\mathrm{and}\hspace{0.33em}\beta$ can be obtained.

$${\alpha }_i={\displaystyle \frac{{\gamma }_i+2a}{{\mathit{\Lambda }}_i^2+2b}},\hspace{0.33em}\beta ={\displaystyle \frac{N-\mathit{\Sigma }\mathit{\Lambda }+2c}{{(t_i-\mathit{\Phi }\mathit{\Lambda })}^2+2d}},i=1,\dots ,N+1.$$

(18)

where ${\gamma }_i=1-{\alpha }_i{\Sigma }_{ii},i=1,\dots ,N+1.$ and ${\Sigma }_{ii}$ is the corresponding element in the i-th row and i-th column of matrix $\mathit{\Sigma }$.

Specifically, the step-by-step parameter estimation process according to Equation (18) is as follows: firstly, substitute the optimal parameters obtained in Section 2.2 into the model and provide initial values ${\alpha }_0$ and ${\beta }_0$, then calculate the values of A, $\mathit{\Sigma }$ and $\mathit{\Lambda }$ based on ${\alpha }_0$, ${\beta }_0$ using $\mathit{A}=diag({\alpha }_1,{\alpha }_2,\dots ,{\alpha }_M)$, $\mathit{\Sigma }={(\mathit{A}+\beta {\mathit{\Phi }}^{\mathrm{T}}\mathit{\Phi })}^{-1}$ and $\mathit{\Lambda }=\beta \mathit{\Sigma }{\mathit{\Phi }}^{\mathrm{T}}{t}_n$. Next, substitute the values of A, $\mathit{\Sigma }$, $\mathit{\Lambda }$, ${\alpha }_0$ and ${\beta }_0$ into Equation (18) to calculate the values for the next time step, denoted as ${\alpha }_1$ and ${\beta }_1$.Thirdly, use the newly obtained ${\alpha }_1$ and ${\beta }_1$ to recalculate the values of A, $\mathit{\Sigma }$ and $\mathit{\Lambda }$. After obtaining the three new values for A, $\mathit{\Sigma }$ and $\mathit{\Lambda }$, substitute them back into Equation (18) to obtain the values of ${\alpha }_2$ and ${\beta }_2$ for the third time step. By iteratively repeating the above steps until traversing the data, a set of values for $\mathit{\alpha }\hspace{0.33em}\mathrm{and}\hspace{0.33em}\beta$ can be obtained, achieving the objective of parameter estimation for two parameters in this section.

3. RUL Prediction

Assuming that historical and current observations are $\{{\mathit{x}}_1,{\mathit{x}}_2,\cdots ,{\mathit{x}}_k\}$, the lifetime L of moment $tim{e}_l$ can be defined according to First Hitting Time (FHT) as

$$L=\inf \{tim{e}_l:t(tim{e}_k+tim{e}_l)\equiv t_{k+l}\in {\rm B}|X_{1:k}\}$$

(19)

where $tim{e}_l$ denotes the length of remaining time from the moment $tim{e}_k$ to the moment of failure, and l refers to the prediction step. $X_{1:k}$ denotes the measurement data from moment $tim{e}_1$ to moment $tim{e}_k$. $t(tim{e}_k+tim{e}_l)$ denotes the predicted mean value at moment $tim{e}_k+tim{e}_l$. $\mathsf{{\rm B}}$ is a known bounded set with elements of fault thresholds, the exact values of which are determined based on expert experience. The physical meaning of Equation (19) is that the moment when the predicted degradation trajectory exceeds the given threshold for the first time is defined as the failure moment, and the time length from the current moment to the failure moment is the remaining lifetime at the current moment.

Assuming that failure occurs when the degradation trajectory reaches the threshold $H\in \mathsf{{\rm B}}$, the distribution of lifetimes at each moment before the variable reaches H is predicted next. The Cumulative Distribution Function (CDF) of the lifetime L, given the observed data $\{{\mathit{x}}_1,{\mathit{x}}_2,\cdots ,{\mathit{x}}_k\}$, is calculated by the following equation.

$$p\{L\le tim{e}_l|X_{1:k}\}=1-p\{z\le {\displaystyle \frac{H-t_{k+l}}{{\beta }_{k+l}}}\}=p\{z\ge {\displaystyle \frac{H-t_{k+l}}{{\beta }_{k+l}}}\}={\displaystyle {\int }_{-\infty }^{g_{k+l}}z}dz\triangleq \rho (g_{k+l})$$

(20)

where ${\beta }_{k+l}$ denotes the predicted variance at moment $tim{e}_k+tim{e}_l$; z is a random variable obeying a standard normal distribution; and $g_{k+l}=(t_{k+l}-H)/{\beta }_{k+l}$ and $\rho (\cdot )$ is the CDF of the random variable z.

From a practical application point of view, for lifetimes there is always $L\ge 0$. Therefore, the truncation condition CDF is introduced based on $L\ge 0$. The truncation CDF for the lifetime L is given by the following equation.

$$p\{L\le tim{e}_l|X_{1:k},L\ge 0\}={\displaystyle \frac{p\{0\le L\le tim{e}_l|X_{1:k}\}}{p\{L\ge 0|X_{1:k}\}}}={\displaystyle \frac{\rho (g_{k+l})-\rho (g_k)}{1-\rho (g_k)}}$$

(21)

Differentiating on $tim{e}_l$ yields the truncation condition PDF for the remaining lifetime L as

$$p_{L|X_{1:k},L\ge 0}={\displaystyle \frac{\mathcal{N}(g_{k+l})\Delta g_{k+l}}{1-\rho (g_k)}}$$

(22)

where $\Delta g_{k+l}$ denotes the derivative of $g_{k+l}$ with respect to $tim{e}_{k+l}$. $\mathcal{N}(\cdot )$ denotes the PDF of a random variable obeying a standard normal distribution. Further, the desired value of RUL can be obtained as follows.

$$E_{L|{\mathit{X}}_{1:k},L\ge 0}={\displaystyle {\int }_{-\infty }^{+\infty }{L}_{p_{L|X_{1:K,L\ge 0}}}}dL$$

(23)

4. RRVM with Embedded Q-Learning

4.1. Q-Learning Algorithm

Within the RL paradigm, the QL technique aligns with the structure of Markov Decision Processes (MDPs) [30]. It accomplishes this by mastering a Q-function designed to gauge the prospective returns of selecting a particular action within a defined state. Employing temporal difference learning to refine these estimates, the method unearths the optimal policy devoid of any pre-existing knowledge of the environment’s underlying dynamics. In the throes of decision-making, the agent discerns its environs and carries out actions from a spectrum of potential choices.

At time step t, with the current environmental state $S_t$, the Agent selects and performs action $A{C}_t$, causing the environment to transition from $S_t$ to $S_{t+1}$, and simultaneously provides the reward $RE(S_t,A{C}_t)$ to the Agent. The Agent continuously repeats this process until the training and learning process is concluded. In the QL algorithm [31], the action-value function $Q(S_t,A{C}_t)$ represents the maximum cumulative reward that the Agent will receive after choosing and taking action $A{C}_t$ in state $S_t$. This value is determined by the immediate reward received by the Agent after choosing and executing the action, plus the value obtained by subsequently following the optimal policy. The value of $Q(S_t,A{C}_t)$ can be expressed by the formula

$$Q(S_t,A{C}_t)=RE(S_t,A{C}_t)+\gamma {\max }_{ac}Q(S_{t+1},A{C}_{t+1})$$

(24)

where AC represents any action from the action set, and the constant parameter $\gamma (0\le \gamma \le 1)$ is known as the discount factor. During the Agent’s training and learning process, it always selects the action corresponding to the state with the maximum Q-value, and then iteratively trains, based on this policy. Through multiple training sessions, the Q-table, which stores the Q-values, is continuously updated. To ensure that QL converges at the appropriate time, an appropriate learning rate, often denoted as $\chi$, is incorporated into the formula. With the learning rate, the update rule for $Q(S_t,A{C}_t)$ can be expressed as

$$Q(S_t,A{C}_t)=(1-\chi )Q(S_t,A{C}_t)+\chi (RE(S_t,A{C}_t)+\gamma \max Q(S_{t+1},A{C}_{t+1}))$$

(25)

4.2. QL-Based Feature Extraction Process

In response to the issues of simplicity in the feature extraction process, susceptibility to loss of effective information, and difficulty in extracting features with high variability in traditional machine learning methods, a feature extraction process is designed based on the QL algorithm under intelligent decision-making. This algorithm is integrated into the front of the RRVM model, providing the model with accurate feature information, which serves as an essential foundation for subsequent precise prediction processes.

The training dataset is initialized as an empty set, and a QL action list with n actions is established to determine the processing of original data. For example, action $A{C}_{3,i}=[0.6\hspace{0.33em}\hspace{0.33em}0.4]$ signifies that at the i-th sampling point, the third action is executed to combine the original data’s median and mean in a 6:4 ratio, thereby obtaining the feature data for that particular sampling point. Using the predictive error of the RRVM model in its current state as the reward function, the specific calculation process involves first using the data formed by randomly selecting actions at corresponding points to train and validate the RRVM. Then, the training error and prediction error are combined in the same proportion as the corresponding action, to obtain the final reward function. The specific framework of the QL-RRVM model is illustrated in Figure 1.

The specific processing steps represented by steps ➀ to ⑬ are as follows:

①②: The agent receives the original signal and selects an action from the action set according to a random selection strategy.

③: The agent processes the original signal at the current moment according to the randomly selected action, using this as the current state $S_i$.

④⑤⑥: Using the RRVM model for training, first calculate matrices A, $\mathit{\Lambda }$ and $\mathit{\Sigma }$ based on the current state $S_i$ and the Gaussian radial basis kernel function. Then, determine the optimal parameters ${\mathit{\alpha }}_i^{best}\mathrm{a}\mathrm{n}\mathrm{d}{\beta }_i^{best}$ for the current moment i. Finally, calculate the predicted value $y_{i+1}={\mathit{W}}^{\mathrm{T}}\mathit{\Phi }({\mathit{x}}_i)$ based on $p(t_{i+1}|t_i,{\mathit{\alpha }}_i^{best},{\beta }_i^{best})\sim \mathcal{N}(t_{i+1}|y_{i+1},{\beta }_{i+1})$.

⑦: Combine the true value ${\tilde{y}}_{i+1}$ corresponding to the current moment’s state with the predicted value $y_{i+1}$ at the same moment to calculate the absolute value of the prediction error. In the specific calculation, both the training set and test set errors are considered: $E{r}_i=A{C}_{j,i}(|{\tilde{y}}_{i+1}^{train}-y_{i+1}^{train}|+|{\tilde{y}}_{i+1}^{test}-y_{i+1}^{test}|)$, where $A{C}_{j,i}$ is the j-th action randomly selected at the i-th moment, with $1\le j\le 6$; $y_{i+1}^{train}$ is the predicted value for the i-th point in the training set, and ${\tilde{y}}_{i+1}^{train}$ is its corresponding true value. The same applies to the test set.

⑧: The current Q-Value is saved and transmitted to the agent, which then iteratively performs steps ① to ⑦ until all actions have been traversed.

⑨⑩⑪: The Q-Values obtained from executing all actions are summarized in a Q-Table, and the original data is processed accordingly, based on the data in this table, to obtain the final training dataset, referred to as Data.

⑫⑬: Using the obtained dataset, the model is trained and validated. Finally, unknown labeled signals are input into the model to achieve prediction.

5. Case Study

5.1. Platform Description and Data Preprocessing

In this section, an experimental setup jointly developed by the China CRRZ Zhuzhou Institute and Central South University [32,33] is employed to substantiate the QL-RRVM method. Intermediate DC-link support capacitor-performance degradation data are extracted to predict the RUL of the capacitor. Figure 2 illustrates the primary components of the testbed, encompassing four main sections: (a) data acquisition and real-time monitoring; (b) signal processing; (c) traction control system; and (d) dSAPCE. It is important to note that the immediate DC-link is equipped with a pair of DC-link capacitors, namely the upper and lower DC-link capacitors. These capacitors are instrumental in sustaining voltage stability for the high-speed train’s traction systems. Given that two voltage sensors are installed within the DC-link, the upper and lower terminal voltages are read, to accurately mirror the capacitors’ health status.

Initially, the experimental platform generates eight distinct sets of original voltage signals under uniform operational conditions, differentiated by their sampling batches. These are sequentially labeled as original dataset 1 through original dataset 8. Each dataset encompasses both the upper and lower terminal voltages. Figure 3 provides a visual representation, encapsulating the entire lifespan of the system, highlighting the stages of extended normal operation, progressive degradation, and eventual shutdown. As this study focuses on the data pertaining to the degradation process, the original data are segmented, to clearly exhibit degradation patterns, which are illustrated in Figure 4.

5.2. QL Algorithm

In this section, the process of feature data extraction is implemented using QL, corresponding to the steps from ① to ⑪ in Figure 1. During this process, the parameters of the QL algorithm are set as follows: the learning rate alpha is set to 0.1, the discount factor gamma = 0.9, and the exploration rate epsilon = 0.1, with the number of training rounds defined as 100. The actions to be taken are defined in a matrix named ‘actions’, which consists of six sets of vectors, each representing a data processing action. For instance, the vector [0.6 0.4] signifies that the current data point is processed as 0.6 times the median plus 0.4 times the mean. The Q-table is initialized as a zero-vector matrix, and the exploration rate is designed to be greater than 0.01 and gradually decreasing. After the process, multiple sets of original data (as shown in Figure 5 (1)) are processed, resulting in the new data as depicted in Figure 5 (2).

5.3. Voltage Degradation-Trend Prediction

In this subsection, a prediction of the voltage degradation trend is implemented. Initially, the RRVM model is trained using the final operational data obtained in Section 2. The hyperparameters are set with the following initial values: scalar hyperparameter to 0.2, inverse noise variance to 0.5, kernel function bandwidth to 3, and the maximum number of iterations to 50. The prediction lag L is set to 10, and the prediction step is set to 1, meaning that for each training instance, every 10 data points are used to predict the corresponding subsequent data points. The final training results of the QL-RRVM model are depicted in Figure 6, where the upper and lower curves in Figure 6a correspond to the degradation trends of the upper and lower voltages, respectively. Figure 6b,c represent the corresponding training errors. From this figure, it can be observed that the QL-RRVM model exhibits a satisfactory prediction performance for this dataset.

After the training is completed, the original data from Figure 5a are used as the prediction set to verify the predictive capability of the model. Additionally, a comparative experiment is designed based on the QL-RRVM and QL-RVM models on the original dataset 6, to show the effectiveness of the robustness on outliers. Figure 7a–c correspond to the overall prediction comparison, the lower voltage-prediction error comparison, and the upper voltage-prediction error comparison, respectively. By comparing the prediction effects of the QL-RRVM model with the QL-RVM model, it can be observed that the predictive accuracy of the improved model presented in this paper is higher, and the prediction trend is smoother, indicating better robustness.

To validate the superiority of incorporating QL in the proposed QL-RRVM model, the QL algorithm component is removed, and the following comparative experiments are conducted. Figure 8 illustrates the degradation trend prediction experiment using the RRVM model and the RVM model. It can be observed that significant prediction errors occur at both the beginning and the end of the dataset. This is attributed to the unprocessed dataset lacking representativeness, which makes it difficult to make accurate predictions for other data.

In this study, four evaluation criteria—Mean Absolute Error (MAE), Root Mean Squared Error (RMSE), Mean Absolute Percentage Error (MAPE) and Pearson Correlation Coefficient—are employed to further assess the prediction performance. The criteria are given below.

$$MAE={\displaystyle \frac{1}{M}}{\displaystyle \sum _{i=1}^M|{\tilde{y}}_i-y_i|}$$

(26)

$$RMSE=\sqrt{{\displaystyle \frac{1}{N}}{\displaystyle \sum _{i=1}^N{({\tilde{y}}_i-y_i)}^2}}$$

(27)

$$MAPE={\displaystyle \frac{1}{N}}{\displaystyle \sum _{i=1}^N\left|\frac{{\tilde{y}}_i-y_i}{y_i}\right|}\times 100$$

(28)

$$R={\displaystyle \frac{{\displaystyle \sum _{i=1}^N(y_i-{\overline{y}}_{*})({\tilde{y}}_i-\overline{y})}}{\sqrt{{\displaystyle \sum _{i=1}^N{(y_i-{\overline{y}}_{*})}^2{({\tilde{y}}_i-\overline{y})}^2}}}}$$

(29)

where N, ${\tilde{y}}_i$, $y_i$, ${\overline{y}}_{*}$ and $\overline{y}$ represent the number of samples, the true values, the predicted values, the true value means and the predicted value means, respectively.

Utilizing four evaluation metrics, trend predictions are conducted among multiple models, which are shown in

Table 1. Evaluation of the results of four models for predicting degradation trends on different datasets.

DATA	MODEL	MAE (V)	RMSE (V)	MAPE (%)	R
DATA SET 1 (Up and Down)	QL-RRVM	0.0109	0.0393	0.0008	0.99
		0.0576	0.751	0.0039	0.99
	QL-RVM	5.161	6.518	0.319	0.983
		6.749	8.599	0.501	0.969
	RRVM	15.790	41.714	0.968	0.818
		14.270	38.316	1.039	0.846
	RVM	24.317	50.758	1.483	0.729
		19.432	41.892	1.399	0.815
DATA SET 6 (Up and Down)	QL-RRVM	0.203	3.441	0.0137	0.997
		0.057	0.751	0.003	0.998
	QL-RVM	5.652	7.861	0.359	0.986
		6.967	9.041	0.504	0.982
	RRVM	1.417	10.835	0.095	0.974
		4.829	24.082	0.328	0.872
	RVM	8.051	12.581	0.513	0.964
		11.064	29.578	0.775	0.806

. By comparing the data in the table, it is evident that the QL-RRVM model provides the best prediction for the degradation trend, followed by the RRVM model. This further demonstrates the effectiveness of the introduction of the QL algorithm and the robustness improvements made in this paper.

5.4. RUL Prediction

Finally, the predicted degradation trends obtained in Section 5.4 are utilized to forecast the RUL of the components. In this paper, the RUL of a component is linked to its corresponding degradation trend, with the final degradation moment defined as the failure state, i.e., RUL = 0. Figure 9 presents the training results of the RUL prediction for the components using the QL-RRVM model based on the degradation trends derived from Figure 7. It is observed that the model presented in this paper accurately predicts the trend of RUL changes.

Similarly, multiple comparative experiments are designed to verify the superiority of the model proposed in this paper. Figure 10 compares the RUL prediction performance of the RRVM and RVM models when the QL algorithm is incorporated. Figure 10a,b represent the prediction results of the two models under the lower- and upper-voltage data, respectively, while Figure 10c–f correspond to the respective prediction errors. From Figure 10, it can be seen that the RRVM model proposed in this paper provides more accurate RUL predictions compared to the RVM model, with less fluctuation.

Figure 11 presents the RUL predictions for the same dataset by four different models. Figure 11a–d correspond to the prediction results of the four models, respectively, with the bar charts in the figures representing the absolute values of the prediction errors. From Figure 11, one can visually discern the differences in the prediction results among the four models, where the QL-RRVM model yields the best results, and the RVM model the worst. Similarly, comparative experiments are conducted on multiple datasets, and the results are accurately assessed using three evaluation metrics, as shown in

Table 2. Evaluation of RUL prediction results for four models on different datasets.

DATA	MODEL	MAE (RUL)	RMSE (RUL)	R
DATA SET 1	QL-RRVM	63.023	77.391	0.996
	QL-RVM	86.729	108.418	0.983
	RRVM	131.399	157.154	0.984
	RVM	110.987	147.024	0.907
DATA SET 3	QL-RRVM	53.386	75.552	0.967
	QL-RVM	78.646	133.706	0.890
	RRVM	149.002	172.939	0.819
	RVM	197.572	235.4965	0.654
DATA SET 6	QL-RRVM	44.581	55.7997	0.983
	QL-RVM	96.39	119.922	0.969
	RRVM	86.472	112.719	0.958
	RVM	110.987	147.024	0.907

. The comparison of corresponding metrics between the QL-RRVM and QL-RVM models demonstrates the effectiveness of the robustness improvements made in this paper. The comparison of corresponding metrics between the QL-RRVM and RRVM models proves the efficacy of incorporating the QL algorithm in this study.

6. Conclusions

This paper proposes a robust RVM model based on QL for RUL prediction. The robustness issue of the traditional RVM model is addressed by introducing a t-distribution; the QL algorithm is incorporated into the model, which processes the original signals by selecting different actions to obtain the features data. The article achieves RUL prediction based on the rolling prediction of degradation trends, enhancing the prediction interpretability.

The effectiveness and superiority of the proposed method are verified under identical data and operational conditions using four models: QL-RRVM, QL-RVM, RRVM, and RVM. As the corresponding RMSE values obtained are 75.552, 133.706, 172.939, and 235.4965, respectively, it indicates that the QL-RRVM model exhibits an approximate 67% improvement in prediction performance. Furthermore, it is evident that the present study significantly enhances feature extraction and robustness.

The proposed QL-RRVM model involves multiple nested loops and complex parameter optimization, which leads to substantial time consumption during training and prediction. In future research, the QL algorithm is planned to be optimized to address this issue. Additionally, within the RRVM module, the large number of kernel functions corresponding to each variable affects computational efficiency. Therefore, future focuses include implementing pruning on the kernel functions to enhance the model’s computational efficiency.