## 2 Methodology
### Parameter Initialization and Updating Algorithm
\right.. \tag{9}\]  
Let \(\partial\ln l\Big{(}\Theta|\mathbf{y_{c}}\Big{)}/\partial\mu_{\lambda}=\partial\ln l \Big{(}\Theta|\mathbf{y_{c}}\Big{)}/\partial\sigma_{\lambda}^{2}=0\), then the maximum likelihood estimations (MLEs) of \(\mu_{\lambda}\) and \(\sigma_{\lambda}^{2}\) based on the monitoring data up to \(t_{c}\), denoted as \(\mu_{\lambda,0}\) and \(\sigma_{\lambda,0}^{2}\), are derived as  
\[\left\{ \begin{array}{c}\mu_{\lambda,0}=\frac{\mathbf{t_{c}}^{\prime}\big{(}\mathbf{t _{c}}\mathbf{t_{c}}^{\prime}+\widetilde{\Omega}\big{)}^{-1}\mathbf{y_{c}}}{\mathbf{t_{c}} \big{(}\mathbf{t_{c}}\mathbf{t_{c}}^{\prime}+\widetilde{\Omega}\big{)}^{-1}\mathbf{t_{c}}} \\ \sigma_{\lambda,0}^{2}=\frac{\big{(}\mathbf{y_{c}}-\mu_{\lambda,0}\mathbf{t_{c}}\big{)} ^{\prime}\big{(}\mathbf{t_{c}}\mathbf{t_{c}}^{\prime}+\widetilde{\Omega}\big{)}^{-1} \big{(}\mathbf{y_{c}}-\mu_{\lambda,0}\mathbf{t_{c}}\big{)}}{c}\end{array} \right.. \tag{10}\]  
Substitute Equation (10) into Equation (8), we obtained a simplified version of the log-likelihood function with model parameters \(\theta\), \(m\) and \(\sigma_{m}^{2}\). In this way, the MLEs of \(\theta\), \(m\) and \(\sigma_{m}^{2}\), denoted as \(\theta_{0}\), \(m_{0}\) and \(\sigma_{m,0}^{2}\), can be obtained using some numerical optimization algorithms. Finally, we could obtain \(\mu_{\lambda,0}\) and \(\sigma_{\lambda,0}^{2}\) by substituting \(\theta_{0}\), \(m_{0}\) and \(\sigma_{m,0}^{2}\) into Equation (10). After offline training of model parameters based on the extracted RMS values up to \(t_{c}\), we could further carry out online RUL prediction when a new monitoring data is available. Based on the proposed model, the RUL of the unit at \(t_{i}\), denoted as \(r_{i}\), is defined as  
\[R_{i}=\inf\{r_{i}:x(t_{i}+r_{i})\geq d\big{|}\mathbf{x_{i}}\}=\inf\{r_{i}:y(t_{i}+r_ {i})\geq d+\varepsilon\big{|}\mathbf{y_{i}}\}, \tag{11}\]  
where \(d\) represents the failure threshold, \(\mathbf{x_{i}}=\{x_{1},\cdots,x_{i}\}\), \(\mathbf{y_{i}}=\{y_{1},\cdots,y_{i}\}\) and \(i=c+1,c+2,\cdots,n\). According to the definition, the conditional probability density function (PDF) of \(r_{i}\) can be derived as  
\[f_{kj|\lambda,y_{i}}\big{(}r_{i}|\lambda,\mathbf{y_{i}}\big{)}=\frac{d_{i}+ \varepsilon-\lambda(p_{i}-r_{i}q_{i})}{\sqrt{2\pi m\lambda r_{i}^{3}}}\cdot \exp\!\left[-\frac{\left(d_{i}+\varepsilon-\lambda p_{i}\right)^{2}}{2m\lambda r _{i}}\right], \tag{12}\]  
where \(d_{i}=d-y_{i}\), \(p_{i}=\int_{t_{i}}^{t_{i}+r_{i}}\varphi(\tau;\theta)d\tau\) and \(q_{i}=\varphi(t_{i}+r_{i}\theta)\). Equation (12) states that the conditional PDF of \(r_{i}\) uses the latest monitoring data instead of the entire historical data. This is conflict with the degradation evolution pattern introduced in Equation (11). To incorporate the history of observations while maintain the nice property of Wiener process model, we considered an updating procedure for parameter \(\lambda\) by using a linear random walk model [32], i.e.,  
\[\lambda_{i}=\lambda_{i-1}+\eta, \tag{13}\]  
where \(\eta\sim N(0,\kappa)\). This model is chosen because of its simplicity and practical rationality. According to the property of random walk, the hidden drift at time point \(t_{i}\), \(\lambda_{i}\) is a time-dependent variable dependent on \(\lambda_{i-1}\) and adjusted by a noise term. Accordingly, the state space model for observation \(y_{i}\) can be constructed as  
\[y_{i}=y_{i-1}+\lambda_{i-1}\int_{t_{i-1}}^{t_{i}}\varphi(\tau;\theta)d\tau+ \sqrt{m\lambda_{i-1}}\Delta B_{i}+\Delta\varepsilon_{i}, \tag{14}\]  
where \(t_{0}=0\), \(y_{0}=0\), \(\Delta B_{i}=B(t_{i})-B(t_{i-1})\sim N(0,t_{i}-t_{i-1})\) and \(\Delta\varepsilon_{i}=\varepsilon_{i}-\varepsilon_{i-1}\sim N\big{(}0,2\sigma _{m}^{2}\big{)}\). Since the existence of hidden states in both state transition function Equation (13) and observation function Equation (14), we could not utilize maximum likelihood estimation method to update the model parameters and system states directly.