Aeroengine performance degradation prediction method considering operating conditions

Abstract

It is significant to predict the performance degradation of complex electromechanical systems. Among the existing performance degradation prediction models, belief rule base (BRB) is a model that deal with quantitative data and qualitative information with uncertainty. However, when analyzing dynamic systems where observable indicators change frequently over time and working conditions, the traditional belief rule base (BRB) can not adapt to frequent changes in working conditions, such as the prediction of aeroengine performance degradation considering working condition. For the sake of settling this problem, this paper puts forward a new hidden belief rule base (HBRB) prediction method, in which the performance of aeroengines is regarded as hidden behavior, and operating conditions are used as observable indicators of the HBRB model to describe the hidden behavior to solve the problem of performance degradation prediction under different times and operating conditions. The performance degradation prediction case study of turbofan aeroengine simulation experiments proves the advantages of HBRB model, and the results testify the effectiveness and practicability of this method. Furthermore, it is compared with other advanced forecasting methods. The results testify this model can generate better predictions in aspects of accuracy and interpretability.

1. Introduction

With the rapid development of technology, some complex electromechanical systems with high technical and high integration have been widely used in engineering practice. These electromechanical systems will be affected by various mechanical stresses and thermal stresses during operation, and their performance will continue to decline. However, some complex electromechanical systems are affected by factors such as high temperature, high speed, load, and impact due to the uncertainty of the working environment, these factors will accelerate the rapid decline of the performance of the entire system [1]. When the performance of the electromechanical system degrades and suddenly fails due to these factors in an unknown state, it may cause huge economic losses and destructive consequences, and bring irreparable adverse effects to the living environment. Therefore, the secure and reliable operation of the complex system has great significance to protect the life and property safety of human, and the social economy. So it is vital to study the health status and performance degradation of complex electromechanical systems. By predicting the performance degradation of the evaluation results, experts can determine the best maintenance scheme.

In the field of forecasting, a large of researchers have put forward many methods ranging from simple linear regression to advance nonlinear forecasting [2]. Performance degradation prediction theory is becoming more and more mature, and the actual application effect of engineering is also improving. The existing performance degradation prediction modeling methods of complex electromechanical systems include the following three categories.

1) Physical-model-based methods, including Kalman predictor[3], particle predictor[4], and strong tracking predictor[5] etc. This kind of method uses some engineering principles to build predictive models. First, the established physical model must be accurate, select the state or parameters that are able to reflect the systematic behavior, and then use the data to forecast the system behavior. Only in this way can these methods be very effective. However, due to the widespread non-linearity and uncertainty in physical models, accurate mathematical models are often difficult to obtain, which increases the difficulty of predicting and modeling the performance degradation of complex electromechanical systems.

2) Qualitative-knowledge-based methods, such as expert system[6], Petri net[7] etc. Expert system (ES) can handle the uncertain information, and its reasoning process is interpretable and flexible. However, as we all know, the knowledge from experts in a specific field is usually inaccurate. Because of the dynamic complexity of expert systems and the limitations of expert knowledge, when meeting a new performance degradation prediction problem, the knowledge base is incomplete. In this case, methods based on qualitative knowledge have limited capabilities in predicting performance degradation and providing more robust problem solving.

3) Data-driven methods, this kind of forecasting method has become a popular health management and forecasting method [8,9]. These methods are mainly aimed at building a data model, adjusting the internal parameters of the model by inputting training data samples, and making predictions by testing data samples. Nowadays, many prediction methods have been proposed to realize system behavior prediction, including improved autoregressive moving-average-model (ARMA) methods [10], support-vector machine (SVM) [11-13], neural networks [14-15], etc. In general, data-driven methods have shown great advantages in some fields, but most of these methods belong to black box models and lack of interpretability. Therefore, it may be impossible to use this type of black box model in a prediction application system that requires systematic reasoning and provides an explanation of its results. In addition, the traditional data-driven methods need number of data for training are difficult to effectively deal with the performance degradation prediction modeling problem of complex electromechanical systems.

The prediction method based on logical rules is interpretable, and it can also enhance the system reliability. Based on this, Yang et al. [16] put forward a rule-based evidence reasoning approach, also called BRB system. It is compared with the IF-THEN rule base, the BRB method provides richer and truer solutions. Parameters and output of BRB system are comprehensible, and experts can intervene in its structural adjustment. At present, it has been applied to offshore systems, oil pipelines, multi-attribute decision analysis, etc [17-20].

However, when predicting complex industrial systems, the traditional BRB is not suitable for directly predicting the implicit behavior of the system. Zhou [21] et al proposed the HBRB prediction model in 2015 to address the above problems.

At present, HBRB prediction method has been applied in the industrial field. Hu [22] et al. proposed a HBRB prediction model with multidimensional input for network security conditions. In the paper [23], it is found that there are many factors affecting the implicit behavior of gyroscope in navigation platform system, but the relationship between them is difficult to establish, so a delay implicit confidence rule base model considering the delay history information is proposed.

The performance of an aeroengine is an hidden behavior, and the observed data of aeroengine operation can be used to predict its performance. Therefore, this paper put forward a hidden belief rule base (HBRB) prediction model. The main contributions are as follows:

1. In view of the problem that traditional BRB can not adapt to frequent changes in aero engine performance when analyzing the observable index changes frequently with time and working conditions, a HBRB prediction model for aeroengine performance degradation considering hidden behavior is proposed.

2. The HBRB prediction model regards the performance of aeroengine at different times as the hidden behavior, and the working condition as the observable index of the model to describe its hidden behavior, which solves the problem of performance degradation prediction at different times and under different working conditions.

The remaining structure of this paper is as follows. In Section Ⅱ, introduce the traditional BRB model and performance degradation modeling. In Section Ⅲ, a performance degradation prediction model considering working conditions is put forward. In Section Ⅳ, conduct simulation experiment. In Section V, summarize the work of this paper.

2. Problem Formulation

Since the performance of aeroengine is a kind of implicit behavior, it is impossible to predict its performance directly, so this paper proposes an aeroengine performance degradation prediction method considering operating conditions. The model first regards the aero engine performance x(t) at time t as an hidden behavior, and uses operating conditions as an observable index O_t of performance x(t) to solve the problem of performance degradation prediction under different times and operating conditions. In the newly proposed HBRB prediction model, usually experts provide the initial parameters, which may be inaccurate, resulting in the decline of prediction results. So this paper adopts the projection covariance matrix adaption evolution strategy (P-CMA-ES) algorithm to optimize the HBRB parameters. Using two public simulation data sets published by NASA to construct case studies to testify the model effectiveness when observable indicators change frequently over time and working conditions. The HBRB model is interpretable and transparent in the model parameters and prediction results. In order to test the model advantages, compareing with several other mature forecasting methods, such as BP neural network, convolution neural network (CNN), fuzzy C-means clustering (FCM) and hidden Markov model (HMM).

This rest of this section will briefly introduce the traditional BRB model and the basic principle of performance degradation prediction modeling of complex electromechanical systems.

2.1 Brief description of BRB model in HBRB modeling

For the sake of establishing an hidden behavior prediction model, constructing a basic BRB model. The confidence rules are described as follows:

R_k = If x(t -1) is D_j, Then x(t) is

{(D₁, β_1,k), ⋯, (D_N, β_N,k), (D, β_D,k)}

with rule weight θ_k and attribute weight δ_k       (1)

where R_k(k = 1, ..., L) represents kth belief. x(t - 1) represents the system behavior and is the prerequisite attribute of BRB at t −1 time. θ_k is the kth rule weight. D_j is the referential value of x(t - 1) and D_j ∈ D, D = {D₁, D₂, ..., D_N} is a set of the prerequisite values of x(t - 1). δ_k is the prerequisite attribute. Because tn the BRB model only has a input, the premise attribute weight is δ_k = 1 . D_j(j = 1, ..., N) represents the jth consequent. β_j,k(j = 1, 2, ..., N, k = 1, 2, ..., N) represents the confidence relative to the jth evaluation result D_j. β_D,k(t) is the belief that is not assigned to any D_j, which can reflect the uncertainty of expert knowledge. In addition \(\begin{aligned}\beta_{D, k}(t)+\sum_{j=1}^{N} \beta_{j, k}(t)=1\end{aligned}\).

The belief rules in (1) can use the following four methods to establishe: 1) Extracte belief rules from expert experience, 2) Analyse existing historical data to extracte belief rules, 3) If available, use the previous belief rule base, 4) Generate rules randomly without any prior knowledge.

Therefore, the prediction model can be expressed as the following formula:

x(t + 1) =g(φ₁, x(t))       (2)

Equation (2) above is the system equation, where g represents the nonlinear mapping function, φ₁ represents the parameter vector, and φ₁ = [θ₁, θ₂, ..., θ_N, β_1,1, ..., β_N,N, β_D,1, ..., β_D,N]^T. For the sake of predicting the future behavior of complex electromechanical equipment, knowing the current information is necessary. That is, the current information should be observable. Therefore, the prediction model shown in (2) can't predict the hidden behavior. For the sake of settling this problem, a new HBRB prediction model is proposed.

2.2 HBRB prediction model

When analyzing the complex electromechanical system, because the system indexes can not be completely tested by sensors, such as the performance of aeroengine, so the traditional BRB model can not fully mine the information in the data. Therefore, a BRB prediction model considering hidden behavior is proposed. Its functions are as follows.

As described above, the hidden behavior x(t) can be reflected by observation O(t), assuming that the observation equation is:

O(t) = y(φ₂, x(t), η(t))       (3)

The above equation is the observation equation. Where y represents the nonlinear mapping function to be determined, which means the relationship between hidden behavior and observed indicators. φ₂ represents the parameter vector. η(t) represents the part of the test data contaminated by external factors. O(t) represents the observation information at time t , and O(t) = [o(1),o(2),...,o(t)]^T. In particular, we assume that η(t) obey the Gaussian distribution. Usually experts give initial parameter vectors, which may be inaccurate. In order to achieve accurate prediction, it is crucial to optimize these parameters.

3. Prediction model of brb considering hidden behavior

This section first introduces the objective function, and then introduces the reasoning process of HBRB model. Finally, using P-CMA-ES algorithmto optimize these parameters.

3.1 The construction process of likelihood function structure

When the observation information o(t) is a one-dimensional vector, it can be regarded as a scalar. Let O(t) = [o(1),o(2),...,o(t)]^T, the likelihood function is constructed as follows:

\(\begin{aligned}L(\Omega)=\prod_{t=2}^{T} f(o(t) \mid O(t-1))\end{aligned}\)       (4)

where Ω indicates the parameter vector, Ω = [φ^T₁,φ^T₂] is consist of rule weight θ_k, belief degree β_k and parameters used by (2).

(1) Calculation process of f(o(t)|O(t - 1))

It can be found from (4), which f(o(t)|O(t - 1)) requires to be calculated. The calculation process of f(o(t)|O(t - 1)) is given below, that is:

\(\begin{aligned}f(o(t) \mid O(t-1))=\sum_{x(t)=U\left(D_{1}\right)}^{U\left(D_{N}\right)} f(o(t) \mid x(t)) f(x(t) \mid O(t-1))\end{aligned}\)       (5)

where U(D_N) represents the utility value of D_j(j = 1, 2, ..., N).

(2) Calculation process of f(x(t)|O(t - 1))

For the sake of determine (5), f(x(t)|O(t - 1)) can be calculated as follows:

\(\begin{aligned}f(x(t) \mid O(t-1))=\sum_{x(t-1)=U\left(D_{1}\right)}^{U\left(D_{N}\right)} f(x(t-1) \mid O(t-1)) f(x(t) \mid x(t-1))\end{aligned}\)       (6)

Substituting (6) into (5) yields the following equation, that is:

\(\begin{aligned}f(o(t) \mid O(t-1))=\sum_{x(t)=U\left(D_{1}\right)}^{U\left(D_{N}\right)} f(o(t) \mid x(t)) \sum_{x(t-1)=U\left(D_{1}\right)}^{U\left(D_{N}\right)} f(x(t-1) \mid O(t-1)) f(x(t) \mid x(t-1))\end{aligned}\)       (7)

In formula (7), f(x(t - 1)|O(t - 1)) can be calculated by the following formula, that is:

\(\begin{aligned} f(x(t-1) \mid O(t-1)) & =f(x(t-1) \mid o(t-1), O(t-2)) \\ & =\frac{f(x(t-1), o(t-1) \mid O(t-2))}{f(o(t-1) \mid O(t-2))}\end{aligned}\)       (8)

Substituting (6) and (7) into (8) gives the following expression, that is:

\(\begin{aligned} f(x(t-1) \mid O(t-1))= & \frac{f(o(t-1) \mid x(t-1))}{\sum_{x(t-1)=U\left(D_{1}\right)}^{U\left(D_{N}\right)} f(o(t-1) \mid x(t-1))} \\ & \times \frac{\sum_{x(t-2)=U\left(D_{1}\right)}^{U\left(D_{N}\right)} f(x(t-2) \mid O(t-2)) f(x(t-1) \mid x(t-2))}{\sum_{x(t-2)=U\left(D_{1}\right)}^{U\left(D_{N}\right)} f(x(t-2) \mid O(t-2)) f(x(t-1) \mid x(t-2))}\end{aligned}\)       (9)

(3) Calculation process of f(x(t)|x(t - 1))

where f(x(t)|x(t - 1)) can be calculated by:

\(\begin{aligned}f\left(x(t)=U\left(D_{j}\right) \mid x(t-1)\right)=\frac{\beta_{j}(t)+\beta_{D}(t)}{\sum_{i=1}^{N}\left(\beta_{j}(t)+\beta_{D}(t)\right)}\end{aligned}\)       (10)

It should be noted that β_j(t) represents the lower bound likelihood value relative to D_j and (β_j(t) + β_D(t)) represents the upper bound likelihood value relative to D_j. these settings can reflect the uncertainty of the evaluation results.

(4) Reasoning process based on ER rules

In (6)-(9), you need to know f(x(t)|x(t - 1)), and its calculation process is as follows. First, x(t) is expressed as:

S(x(t)) = {(D_k, a_k(t)), k = 1, 2, ..., N}       (11)

where

\(\begin{aligned}\begin{array}{l}\alpha_{k}(t)=\frac{\mu\left(D_{k+1}\right)-x(t)}{\mu\left(D_{k+1}\right)-\mu\left(D_{k}\right)} \\ \alpha_{k+1}(t)=1-\alpha_{k}(t), U\left(D_{k}\right) \leq x(t) \leq U\left(D_{k+1}\right) \\ \alpha_{j}(t)=0, j=1,2, \ldots, N, j \neq k, k+1\end{array}\end{aligned}\)       (12)

where U(D_j) represents the utility value of D_j(j = 1, 2, ..., N).

In order to determine (10), a confidence level needs to be obtained β_j(t), which can be calculated by ER rule [24-26]. Where ω_k(t) represents the activation weight, it is expressed as:

\(\begin{aligned}\omega_{k}(t)=\frac{\theta_{k}\left(\alpha_{k}(t)\right)^{\delta}}{\sum_{l=1}^{N} \theta_{l}\left(\alpha_{l}(t)\right)^{\delta}}\end{aligned}\)       (13)

The model output can be expressed as:

S(x(t)) = {(D_j, β_j(t)), (D, β_D(t)), j = 1, 2, ..., N}       (14)

where β_j(t) is the belief of D_j at time t, and β_D(t) is the belief not set to any D_j, which reflects the uncertainty of the evaluation process.

\(\begin{aligned}\begin{array}{c}\hat{\beta}_{j}(t)=\frac{\mu \times\left[\prod_{k=1}^{N}\left(\omega_{k}(t) \beta_{j, k}+1-\omega_{k}(t) \sum_{i=1}^{N} \beta_{i, k}\right)-\prod_{k=1}^{N}\left(1-\omega_{k}(t) \sum_{i=1}^{N} \beta_{i, k}\right)\right]}{1-\mu \times\left[\prod_{k=1}^{N}\left(1-\omega_{k}(t)\right)\right]} \\ \beta_{D}(t)=1-\sum_{j=1}^{N} \beta_{j}(t)\end{array}\end{aligned}\)       (15)

\(\begin{aligned}\beta_{D}(t)=1-\sum_{j=1}^{N} \beta_{j}(t)\end{aligned}\)       (16)

where

\(\begin{aligned}\mu=\left[\sum_{j=1}^{N} \prod_{k=1}^{N}\left(\omega_{k}(t) \beta_{j, k}+1-\omega_{k}(t) \sum_{i=1}^{N} \beta_{i, k}\right)-(N-1) \prod_{k=1}^{N}\left(1-\omega_{k}(t) \sum_{i=1}^{N} \beta_{i, k}\right)\right]^{-1}\end{aligned}\)       (17)

When o(t) is multidimensional observation information, o(t) = [o₁(t), o₂(t), ..., o_n(t)]^T, assuming that o₁, o₂, ..., o_n is independent, the following formula can be obtained.

\(\begin{aligned}f(o(t) \mid x(t))=\prod_{j=1}^{n} f\left(o_{j}(t) \mid x(t)\right)\end{aligned}\)       (18)

When o(t) is calculated, the likelihood function shown in (4) is obtained, the calculation and implementation process of the likelihood function will be introduced below.

3.2 Objective function to be optimized of HBRB model

By calculating (4)-(18) can get the likelihood function shown in (4), and constructing the following nonlinear objective function:

\(\begin{aligned}\begin{array}{ll}\max & \{(\Omega)\} \\ \text { s.t. } & 0 \leq \theta_{k} \leq 1, k=1, \cdots, L \\ & 0 \leq \beta_{j, k} \leq 1, j=1,2, \cdots, N, ; k=1,2, \cdots, L \\ & 0 \leq \beta_{D, k} \leq 1, k=1,2, \ldots, N \\ & \beta_{D, k}+\sum_{j=1}^{N} \beta_{j, k}=1 \\ & c\left(\varphi_{2}\right) \leq 0 \\ & \operatorname{ceq}\left(\varphi_{2}\right)=0\end{array}\end{aligned}\)       (19)

where the definition of parameters θ_k, β_j,k, and β_D,k, has been given in (1), 2 ceq(φ₂) = 0 and c(φ₂) ≤ 0 represent the equality constraints and inequality constraints satisfied by φ₂ respectively, and the specific forms of ceq(φ₂) and c(φ₂) can be determined according to the equation (3).

3.3 Parameter optimization of HBRB prediction model based on P-CMA-ES

P-CMA-ES [27-29] algorithm is used to optimize the relevant parameters in this paper.

The specific steps are as follows:

(1) Set initial values:

The initial population is generated by taking a solution in the solution space as the central point, also known as the expectation. Give the initial average m⁰ = Ω⁰, initial covariance matrix C⁰, initial step size σ, overall size λ.

(2) Generate initial population:

The initial parameter vector Ω⁰ of the HBRB model is selected as the expectation. The resulting population is as follows:

Ω^g+1_k ~ m^g + σ^gℕ(0, C^g), for k = 1, ..., λ       (20)

where Ω^g+1_k is the ith solution of the generation (g + 1)th, m is the overall mean, σ is the step size, ℕ is normal distribution, C^g is the covariance matrix, its geometric meaning is the elliptic distribution of the population in the solution space.

(3) Projection operation:

The solution produced from the sampling operation can't satisfy the constraints. So a projection operation is put forward. The projection operation is obtained from the following formula.

Ω^g+1_k(1 + n_e x (j - 1) : n_e x j) = Ω^g+1_k(1 + n_e x (j - 1) : n_e x j) - A^T_e x (A_e x A^T_e)^-1 x Ω^g+1_k(1 + n_e x (j - 1) : n_e x j) x A_e       (21)

The above formula is to project the equality constraints in Ω^g+1_k one by one, and finally map them to the feasible domain hyperplane. Compared with other methods, projection operation uses analytical expression to solve equality constraints, which has high operation efficiency.

(4) Select operation:

Select some optimal solutions with better fitness value from the population. The larger the value, the better the solution. And by using the weighted average method to update the new center point of the population.

(5) Reorganization operation:

Update the expectations and guide the population to move towards a better solution. The update method is obtained by the following formula.

\(\begin{aligned}m^{g+1}=\sum_{i=1}^{\varepsilon} \gamma_{i} \Omega_{i: \lambda}^{g+1}\end{aligned}\)       (22)

where γ_i is the weight of the ith solution. λ is the number of solutions. Ω^g+1_i:λ is the ith solution of the gth generation.

(6) Update C operation:

The covariance matrix C is the equal probability density elliptic surface of the population distribution in the solution space. The initial C⁰ is the identity matrix I. In the process of evolution, the direction of the long axis of the covariance matrix always points to the optimal solution, the change of direction controls the evolutionary trend, and the length of the long axis controls the search range. The covariance matrix C is updated refer to (23).

\(\begin{aligned} \boldsymbol{C}^{g+1}= & \left(1-a_{1}-a_{\varepsilon}\right) \boldsymbol{C}^{g}+a_{1} q^{g+1}\left(q^{g+1}\right)^{\mathrm{T}} \\ & +a_{\varepsilon} \sum_{i=1}^{\varepsilon} \gamma_{i}\left(\frac{\left(\Omega_{1: \lambda}^{g+1}-\operatorname{mean}^{g}\right)}{\eta^{g}}\right)\left(\frac{\left(\Omega_{1: \lambda}^{g+1}-\text { mean }^{g}\right)}{\eta^{g}}\right)^{\mathrm{T}}\end{aligned}\)       (23)

where α₁ and α_ε denote the learning factor. q denotes the evolution path. The rule is updated as below:

\(\begin{aligned}q^{g+1}=\left(1-a_{q}\right) q^{g}+\sqrt{a_{q}\left(2-a_{q}\right)\left(\sum_{i=1}^{\varepsilon} \gamma_{i}^{2}\right)^{-1}} \frac{\text { mean }^{g+1}-\text { mean }^{g}}{\eta^{g}}\end{aligned}\)       (24)

where α_q ≤ 1 denotes the evolution path parameter. The step size η is updated refer to (25).

\(\begin{aligned}\eta^{g+1}=\eta^{g} \exp \left(\frac{a_{\eta}}{d_{\eta}}\left(\frac{\left\|q_{\eta}^{g+1}\right\|}{E\|\mathbb{N}(0, \boldsymbol{I})\|}-1\right)\right)\end{aligned}\)       (25)

where d_η is the damping coefficient, E||ℕ(0, I)|| is the expectation of Euclidean paradigm ℕ(0, I). q_η is the conjugate evolution path. a_η is the q_η parameter. q_η update reference following formula.

\(\begin{aligned}q_{\eta}^{g+1}=\left(1-a_{\eta}\right) q_{\eta}^{g}+\sqrt{a_{\eta}\left(2-a_{\eta}\right)\left(\sum_{i=1}^{\varepsilon} \gamma_{i}^{2}\right)^{-1}} C^{(g)-\frac{1}{2}} \frac{\text { mean }^{g+1}-\text { mean }^{g}}{\eta^{g}}\end{aligned}\)       (26)

Repeat steps (1)-(6) until the accuracy requirements are met, and then output the best parameters Ω_best of HBRB model.

3.4 Hidden behavior prediction algorithm steps by HBRB model

To sum up, the prediction process of HBRB model is as follows.

Step 1: Set the initial parameters x⁰ and Ω⁰, and the parameters in Ω⁰ meet the constraints shown in (19).

Step 2: When the observation information O(t) = [o(1),o(2),...,o(t)]^T arrives, when o(t) is a one-dimensional vector, use (5)-(17) to calculate the likelihood function as shown in (4). When o(t) is a multi-dimensional vector, use (5)-(18) to establish the likelihood function.

Step 3: Calculate the optimal parameters through P-CMA-ES algorithm. The steps are as follows in Fig. 1.

Fig. 1. Optimization algorithm specific steps

Step 4: The HBRB prediction model is established to predict the hidden behavior.

4. Case Studies

In this section, turbofan engine simulation test will be constructed to examine the effectiveness of the put forward HBRB model.

Turbofan aeroengine is one of the important components of aircraft, which has great significance to the stable operation of aircraft. The structure diagram of turbofan aeroengine is shown in Fig. 2. Due to the existence of complex system characteristics, such as strong nonlinearity, non-stationarity and uncertainty of turbofan aeroengines, condition monitoring based on multi-source heterogeneous data has become a very challenging task.

Fig. 2. The structure diagram of turbofan aeroengine

4.1 Data set description

This section uses performance parameter measurements from turbofan engine general test data to evaluate the performance degradation prediction status, which is NASA's public data. Including four sets of data, and each set of data contains the training and testing set. The training set contains complete monitoring data of several turbofan engines, and each turbofan engine has recorded monitoring data of all flight cycles from normal operation to failure. The monitoring data of each flight cycle includes 24 monitoring parameters, three of which are operating conditions. The test set contains incomplete data for several turbofan engines, and the termination state of each engine in the test set is not a performance failure state. The four groups of data in the public test data set have their own characteristics: there is only 1 working condition in group 1 and group 3, and there are 6 working conditions in group 2 and group 4. There is only 1 failure mode for group 1 and group 2, and there are 2 failure modes for group 3 and group 4.

4.2 Data preparation

This paper uses the "train_FD001" and "train_FD002" files provided by the dataset to build the performance degradation curve. In "train_FD001", 21 sensor parameters are continuously monitored under single working condition, and the turbofan engine is only affected by time factor. In "train_FD002", 21 sensor parameters are continuously monitored under variable working conditions, and the turbofan engine is affected by time factors and working conditions, variable operation refers to the mode of aeroengine switching according to the requirements or randomly under multiple operating conditions, under this condition, the uncertainty information of aeroengine will greatly increase. When the working modes are constantly switched, the aeroengine needs to adjust its flight parameters at any time. In many cases, it is overloaded, and the working conditions are extremely complicated, which will accelerate the performance degradation of the aeroengine. The time series length of each turbofan engine unit is different. "Column FD001" and "column FD002" include 100 monitoring units and 260 monitoring units respectively. "FD001 train" and "FD002 train" include 20631 and 53759 engine cycles, respectively. All turbofan engine units are reduced from a slight wear condition to a failure threshold, each turbofan engine records monitoring data for all flight cycles from normal operation to failure. The monitoring data of each flight cycle includes 24 monitoring parameters, three of which are operating conditions, namely Height, Mach number, Throttle Rotation transformer Angle, which are used as the input of the observation equation in this paper. The test set contains incomplete data of multiple turbofan engines, and the termination state of each engine in the test set is not a performance failure state. One aeroengine unit data in "train_FD001" and "train_FD002" provided by each selected data set is used for performance degradation prediction modeling, and the results of health status assessment are used as the input of HBRB model, that is, the performance of aeroengine is regarded as hidden behavior and input into HBRB model.

4.3 Prediction of aeroengine performance degradation with hidden behavior under single working condition

4.3.1 The HBRB prediction model

Single condition operation means that the aeroengine works under an operating setting. Under this operating condition, the main reason that affects the degradation speed of the aeroengine is the time factor. And some external factors lead to some uncertain information in the process of aeroengine operation and data acquisition.

An aeroengine unit data in "train_FD001" is used for performance degradation prediction modeling, and the input of this model is shown in Fig. 3.

Fig. 3. Input of performance degradation prediction model

For the sake of constructing the aeroengine prediction model, four reference values are selected for the performance degradation factor, which are Small (S), Medium (M), Large (L) and Very Large (VL). The above reference points need to be quantified, the results are given in Table 1.

Table 1. Reference value of performance degradation factor

In order to predict the performance of aeroengine, it is crucial to build HBRB model. The kth rule is described as:

R_k : If x(t) is D_k, Then x(t + 1) is {(D₁, β_1,k), ⋯, (D_N, β_N,k), (D, β_D,k)}

with rule weight θ_k and attribute weight δ_i,k       (27)

Among them, x(t) represents the performance degradation value of the current aeroengine, x(t + 1) represents predicted performance degradation value of the aeroengine at the next time, and the range is [0,1]. The larger the value, the more serious the performance degradation. R_k(k = 1, 2, ⋯, L) stands for rule k and L stands for the number of rules. D_k ∈ D(D = {D₁, ⋯, D_N}) represents the prediction result of rule k, and also represents the division of reference level of aeroengine health status. For example, it can be divided into four levels: excellent, good, medium and poor. D represents the set of levels. β_j,k(j = 1, 2, ⋯, N) is the belief degree of the jth result in rule k, and N is the number of results. β_D,k represents residual confidence and represents the ignorance of expert knowledge. θ_k and δ_i are rule weight and premise attribute weight respectively.

The observation information in this paper is three sensor parameters representing different working conditions, namely Height, Mach number, Throttle Rotation transformer Angle. As shown in Fig. 4.

Fig. 4. Data collection physical map

The initial parameters are shown in Table 2 and Table 3. The initial belief given by expert experience is given in Table 2. It is assumed that the initial values of θ_k are all 1. The initial values of the observed parameters are given in Table 3.

Table 2. Initial rule weight and belief of HBRB prediction model

Table 3. Initial parameters of observation equation

When the above parameters in Table 1-3 are given, can gain initial parameter vector Ω⁰. Then, the initial HBRB parameters are optimized.

Because the number of aeroengine cycles in train_FD001 is 195, that is, the total data sample in this paper is 195, this paper takes the 45 as testing samples. Based on the initial prediction model and training data, the parameter Ω⁰ is trained. Table 4 shows the rule weight and belief after training. Table 5 shows the observation equation parameter after training.

Table 4. Rule weight and confidence of HBRB prediction model after training

Table 5. The parameters of observation equation after training

After training, the last 45 groups of data are used to verify the HBRB prediction model. In order to demonstrate the stability of the model, ten prediction experiments are conducted, each of which uses the same data. The results are shown in Fig. 5, it is a comparison between the trained HBRB prediction model and the test data. The average value is 0.0057 and the variance is 6.1032 × 10^-6. It can be seen from the results in Fig. 5 that HBRB model has a good prediction effect on aeroengine.

Fig. 5. Prediction results of HBRB model

As can be seen from Fig. 5, the values of the ten experiments all fitted, indicating that the HBRB prediction model can predict the aeroengine performance stably and accurately under a single working condition.

4.3.2 Comparison test verification under single working condition

For the sake of examining the model superiority, other models are used for comparison. The comparison models include initial BRB, BP neural network, CNN, FCM and HMM. The mean square error (MSE) can indicate the performance of the algorithm. Fig. 6 shows the specific results.

Fig. 6. Test results of other models under single working condition

According to the analysis of the simulation Fig. 6, the overall volatility of the 10 tests of BP neural network and CNN is very small, but it can't simulate the test data well. The overall trend of HMM gradually decreases, deviating from the actual value of aero-engine. The whole FCM curve fluctuates greatly, which can't reflect the actual situation well. The initial BRB curve does not follow the test data well. The specific mean square error is shown in Table 6 below.

Table 6. Test results of other models under single working condition

4.4 Hidden behavior prediction of aeroengine performance degradation under variable operating conditions

In the same way as the performance degradation prediction under single working condition, an aeroengine unit data in "train_FD002" is used for performance degradation prediction modeling, and the results of health state assessment are used as the input of HBRB model. The performance degradation factor x(t) and observation information vector are input into the HBRB model, and the algorithm flow steps 1-4 are used for modeling and optimization. There are 271 aeroengine cycles in "train_FD002", so 50 test samples are set, and in order to demonstrate the stability of the model, ten prediction experiments are conducted, each of which uses the same data. It can be seen from the results in Fig. 7 that HBRB model has good predictability. According to statistics, the average value of 10 experiments is 0.0014 and the variance is 1.0440 × 10^-7.

Fig. 7. Prediction results of HBRB model

As can be seen from Fig. 7, the values of each experiment fit well with the original performance values. The trend of HBRB prediction results is very similar to the original performance trend, which indicates that the HBRB prediction model can predict the aeroengine performance well under varying working conditions.

For the sake of examining the advantage of this method, which is compared with initial BRB, BP neural network, CNN, FCM and HMM model. Fig. 8 shows the specific experimental results.

Fig. 8. Test results of other models under working-condition division

According to the analysis of simulation Fig. 8, the overall trend of BP is similar to the actual value, indicating that the performance of aviation act tends to deteriorate in winter, but it is far from the actual value. The overall trend of Hidden Markov Model (HMM) gradually decreases, which deviates from the actual situation of aeroengine. The overall fluctuation of FCM curve is large, which can not well reflect the actual situation. The overall change trend of 10 groups of tests of CNN and Initial BRB is not obvious, and it can not greatly reflect the performance degradation state of aeroengine. The specific mean square error is shown in Table 7 below.

Table 7. Test results under working conditions of other models

It can be seen from the specific MSE values in Table 7 that HBRB model solves the problem of unstable performance curve prediction modeling caused by frequent switching of working conditions, which reflects the advantages of this method.

5. Conclusion

In this paper, an aeroengine performance degradation HBRB prediction model considering hidden behavior is proposed. In this model, the performance of aeroengine is regarded as hidden behavior, and the working state is regarded as the observable index of HBRB model to describe its hidden behavior, which can solve the problem of performance degradation prediction in different time and working states. Taking the simulation experiment of turbofan engine as an example, two examples of performance degradation prediction of turbofan engine system under single and variable operating conditions are established. For the sake of gaining the optimal parameters of HBRB model, using P-CMA-ES to optimize the model parameters. Finally, by comparing with other models, the superiority of HBRB model is proved, and its results can produce better prediction accuracy, which verifies the effectiveness and practicability of this method. This method has a wide engineering application prospect.

The model is also applicable to other systems, such as welding robot systems. In the future work, how to establish a prediction model that is closer to the actual engineering is still a problem that we need to solve.